Glycosidases, including β-D-galactosidase, are involved in a variety of metabolic disorders, such as diabetes, viral or bacterial infections, and cancer. Accordingly, we were prompted to find new β-D-galactosidase inhibitors. Towards this end, we scanned the pharmacophoric space of this enzyme using a set of 41 known inhibitors. Genetic algorithm and multiple linear regression analyses were used to select an optimal combination of pharmacophoric models and physicochemical descriptors to yield self-consistent and predictive quantitative structure-activity relationship (QSAR). Five pharmacophores emerged in the QSAR equations suggesting the existence of more than one binding mode accessible to ligands within β-D-galactosidase pocket. The successful pharmacophores were complemented with strict shape constraints in an attempt to optimize their receiver-operating characteristic curve profiles. The validity of the QSAR equations and the associated pharmacophoric models were experimentally established by the identification of several β-D-galactosidase inhibitors retrieved via in silico search of two structural databases: the National Cancer Institute list of compounds and our in house built structural database of established drugs and agrochemicals.