A mixed molecular orbital and valence bond (MOVB) method for describing the potential energy surface of reactive systems has been developed and applied to a model proton transfer reaction in aqueous solution. The MOVB method is based on a block-localized wave function (BLW) approach for defining the diabatic electronic states. Then, a configuration interaction Hamiltonian is constructed using these diabatic states as the basis function. It was found that the electronic coupling energy is large with a value of about 30 kcal/mol for the H3N-H-NH3+ system, whereas the predicted activation barrier is only 1.2 kcal/mol using the 3-21G basis set. The MOVB results are found to be in good accord with the corresponding ab initio Hartree-Fock calculations for the proton transfer process. We have also incorporated solvent effects into the MOVB Hamiltonian in the spirit of combined QM/MM calculations, and have modeled the proton transfer between ammonium ion and ammonia in water using Monte Carlo simulations. The potential of mean force was computed via free energy perturbation coupled with umbrella sampling techniques using (1) an energy gap mapping approach, and (2) a geometrical mapping procedure. Solvent effects increase the barrier height by about 2.2 kcal/mol from the MOVB and HF ground state potential energy surface. The present study demonstrated the feasibility of ab initio MOVB method for studying chemical reactions by incorporating explicit solvent effects in the description of the reaction coordinate in combined QM/MM simulations.