Lanthanide (Ln) binding to humic acid (HA) has been investigated by combining ultrafiltration and ICP-MS techniques. A Langmuir-sorption-isotherm metal-complexation model was used in conjunction with a linear programming method (LPM) to fit experimental data representing various experimental conditions both in HA/Ln ratio (varying between 5 and 20) and in pH range (from 2 to 10) with an ionic strength of 10(-3) mol L(-1). The LPM approach, not requiring prior knowledge of surface complexation parameters, was used to solve the existing discrepancies in LnHA binding constants and site densities. The application of the LPM to experimental data revealed the presence of two discrete metal binding sites at low humic acid concentrations (5 mg L(-1)), with log metal complexation constants (logK(S,j)) of 2.65+/-0.05 and 7.00 (depending on Ln). The corresponding site densities were 2.71+/-0.57x10(-8) and 0.58+/-0.32x10(-8) mol of Ln(3+)/mg of HA (depending on Ln). Total site densities of 3.28+/-0.28x10(-8), 4.99+/-0.02x10(-8), and 5.01+/-0.01x10(-8) mol mg(-1) were obtained by LPM for humic acid, for humic acid concentrations of 5, 10, and 20 mg L(-1), respectively. These results confirm that lanthanide binding occurs mainly at weak sites (i.e., ca. 80%) and second at strong sites (i.e., ca. 20%). The first group of discrete metal binding sites may be attributed to carboxylic groups (known to be the main binding sites of Ln in HA), and the second metal binding group to phenolic moieties. Moreover, this study evidences heterogeneity in the distribution of the binding sites among Ln. Eventually, the LPM approach produced feasible and reasonable results, but it was less sensitive to error and did not require an a priori assumption of the number and concentration of binding sites.