The effect of substituent on the protonation equilibria of benzaldehydes (A + H+ reversible arrow AH(+), where A and AH(+) represent benzaldehyde and its protonated form) in the gas phase and in solution has been investigated theoretically at the MP2/6-31G*//MP2/6-31G* level applying two solvation (IPCM and SCIPCM) models. The absolute as well as relative Delta G degrees (or log K) values of the gas-phase MP2 results are in good agreement with the experimental values. All the solvation models are inadequate to reproduce experimental rho(+) (slope of log K versus sigma(+) plot) value in aqueous solution. The SCIPCM model gives the best correlation but with a much larger magnitude of rho(+) (=-5.29 relative to the experimental value of rho(+) = -1.88) mainly due to neglect of specific solvation (i.e., hydrogen bonding) effect. It was found that the neglect of hydrogen-bonding effect in the solvation of aldehydes (A) results in unduly higher equilibrium constants (K) for electron donors (sigma(+) < 0); in contrast, the neglect of hydrogen bonding in the solvation of the protonated forms (AH+) leads to unduly lower K values for electron accepters (sigma > 0).