Our aim is to develop an effective computational procedure for predicting the aqueous acid equilibrium constants of protonated benzimidazoles at 298.15 K. The experimental determination of these values, apart from been laborious, is a challenge because of the low water solubility of these compounds. Using a variety of descriptors, quantitative structure-property relationships (QSPR) are explored between the experimental aqueous pK(a) values of a group of fifteen benzimidazoles and descriptors calculated at the B3LYP/6-31 +G(d, p) level of theory. Solvent effects are taken into account with the PCM solvation model through both single-point energy calculations (PCM(sp)), and in the geometry optimizations and frequency calculations (PCM( opt)). Descriptors considered are the Gibbs free-energy change of the acid equilibrium in water, the charges on the acidic hydrogen, and on the basic nitrogen, several orbital energies of the protonated and neutral species, and the volume of the solvent cavity. Multiple linear regressions are used to correlate descriptors to the experimental pK(a) values. Several QSPR equations reproduce the experimental data more accurately, and show stronger correlations than previously attempted methodologies. The predictive capabilities of the QSPR methodologies are tested with four compounds that were not included in the set of benzimidazoles initially investigated. In addition, a correlation between experimental pK(a) values in water and in a 50% ethanol-water solution is used to estimate aqueous pK(a) values.