We have studied the influence of implicit solvent models, inclusion of explicit water molecules, inclusion of vibrational effects, and density functionals on the quality of the predicted pK a of small amino acid side chain models. We found that the inclusion of vibrational effects and explicit water molecules is crucial to improve the correlation between the computed and the experimental values. In these micro-solvated systems, the best agreement between DFT-computed electronic energies and benchmark values is afforded by BHHLYP and B97-2. However, approaching experimental results requires the addition of more than three explicit water molecules, which generates new problems related to the presence of multiple minima in the potential energy surface. It thus appears that a satisfactory ab initio prediction of amino acid side chain pK a will require methods that sample the configurational space in the presence of large solvation shells, while at the same time computing vibrational contributions to the enthalpy and entropy of the system under study in all points of that surface. Pending development of efficient algorithms for those computations, we strongly suggest that whenever counterintuitive protonation states are found in a computational study (e. g., the presence of a neutral aspartate/neutral histidine dyad instead of a deprotonated aspartate/protonated histidine pair), the reaction profile should be computed under each of the different protonation micro-states by constraining the relevant N-H or O-H bonds, in order to avoid artifacts inherent to the complex nature of the factors contributing to the pK a. © 2012 Springer-Verlag.