1

Does anyone know what criteria is used to decide the value of "epsilon" for each parameter in the calculation of the Jacobian matrix in scipy? I am using scipy.curve_fit where I get the covariance matrix, but I need to do this fit in matlab as well where I have to manually calculate the Jacobian and convert to covariance matrix. I have not been able to find the criteria by which scipy decides on the default shift in parameter (i.e. epsilon) used to determine the function gradient

tcasey
  • 39
  • 5
  • the doc states: *diff_step : None or array_like, optional. Determines the relative step size for the finite difference approximation of the Jacobian. The actual step is computed as x * diff_step. If None (default), then diff_step is taken to be a conventional “optimal” power of machine epsilon for the finite difference scheme used* – mikuszefski Sep 10 '20 at 07:01
  • Internally the source code makes: `if diff_step is None: epsfcn = EPS else: epsfcn = diff_step**2` – mikuszefski Sep 10 '20 at 07:03
  • Ok, but what is the definition of "optimal power of machine epsilon" if no diff_step is given? For example, I could say step my parameter by a standard 1%, 5%, etc. or use some algorithm to compute "optimal" based on the parameter magnitude. Do you know how EPS is computed in the above? – tcasey Sep 10 '20 at 18:55
  • Well it is defined in the scipy package ".../scipy/optimize/_lsq/common.py" as `EPS = np.finfo(float).eps` which is documeted here https://numpy.org/doc/stable/reference/generated/numpy.finfo.html For fun you can look here https://stackoverflow.com/a/25155518/803359 – mikuszefski Sep 11 '20 at 07:51
  • Oh I see, it actually is just the machine precision, nothing fancy. Thank you @mikuszefski! I have parameter values that span different orders of magnitude in the same calculation, and the simulation is very sensitive to some while not very sensitive to others. I've been trying to set the parameter step by using a while loop and shifting the parameter using new_parameter = best_parameter*eps^(1/power) where I choose the first power that, when the simulation is recalculated with new_parameter, (new_MSE - best_MSE) exceeds eps. This seems to work ok for the moment – tcasey Sep 11 '20 at 15:37

0 Answers0