<p>Determining a fitting function for a dataset represents an important challenge within applied mathematics and statistics. Among the most utilized techniques, those based on the Levenberg-Marquardt method for solving fitting problems modeled through least squares are prominent. Unfortunately, traditional least squares formulations do not support outlier detection, which may compromise the model. By reformulating the least squares problem within the context of Low-Order Value Optimization (LOVO), the identification of such points becomes feasible. However, this approach introduces new elements that add additional costs to the iterations, such as the need for sorting operations required to evaluate the objective function. To reduce the number of necessary sorting operations, this work proposes a new multi-step fitting method based on the Levenberg-Marquardt approach, but designed to address the fitting problem in a LOVO context. We present a well-defined algorithm, demonstrate convergence to stationary points in a weak sense, and provide numerical results indicating that this new approach is promising for the purpose of this study.</p>