The conformational itineraries taken by carbohydrate residues in the catalytic subsite of retaining glycoside hydrolases (GHs), harness the link between substrate conformation and reactivity. GHs' active sites may be described as a combination of subsites dedicated to the binding of individual sugar residues and to catalysis. The three-dimensional structure of GH:carbohydrate complexes has demonstrated that carbohydrate ring conformation changes in an ordered manner during catalysis. Here we demonstrate in silico that a link exists between subsite binding dynamics and substrate specificity for β-galactosidases from clan GH-A families GH1, GH2, GH35, GH42 and GH59. Different oligosaccharides were docked in the active site of reference β-galactosidase structures using Vina-Carb. Subsequent molecular dynamics (MD) simulations revealed that these enzymes favor a high degree of flexibility and ring distortion of the substrate the lytic subsite -1. Although the β-galactosidase families examined are structurally and mechanistically related, distinct patterns of ring distortion were unveiled for the different families. For β-galactosidases, three different family-dependent reaction itineraries (S → H → C, B → H/ E → C, and S → E/ H → C) were identified, all compatible with the antiperiplanar lone pair hypothesis (ALPH) for the hydrolysis of β-glycosides. This comparative study reveals the fuzzy character of the changes in carbohydrate ring geometry prior to carbohydrate hydrolysis.