An accurate representation of reality in numerical variably-saturated flow models requires reliable estimates of necessary model parameters. Inverse modeling seeks to estimate parameters such as the saturated and residual water contents, the saturated hydraulic conductivity, the shape parameters of the soil hydraulic functions, using easily attainable observations of actual or cumulative water fluxes, pressure heads, water contents, and concentrations. The inverse procedure usually combines the nonlinear least-squares-based (SSQ) parameter optimization method with a numerical solution of the variably-saturated flow and transport equations. The SSQ-based inverse method is however sensitive to outliers. A novel Squared ε-Insensitive Loss Function (SILF) approach is introduced in this study. The SILF approach is inspired by the ε-insensitive loss function proposed by