mrpy.fitting.fit_curve.get_fit_curve

mrpy.fitting.fit_curve.get_fit_curve(data_m, data_mf, x0, bounds=None, jac=True, sig_rhomean=inf, sig_data=1, **minimize_kw)

Perform basic LSQ fitting of the MRP curve to the given curve. The fit is performed in log-log space

Parameters:

data_m : array

Masses (not log)

data_mf : array

Mass function corresponding to m.

x0 : float, optional

Initial guess for each of the MRP parameters (hs, alpha, beta, lnA)

bounds : list of tuples, optional

If None, don’t use bounds. If true, set bounds based on bounds passed.

jac : bool, optional

Whether to use analytic jacobian (usually a good idea)

sig_data : array_like, optional

The uncertainty of the data (standard deviation). This is used in the likelihood to weight different mass scales. If scalar, all mass scales are weighted evenly.

sig_rhomean,: float, optional

This controls how much influence the total mean density of the universe has on the likelihood. The default value of inf means it is completely ignored. If it is 0, it becomes an absolute constraint, so that the total mass density of the universe is perfectly matched (setting the normalisation). In between, it acts as an uncertainty on this value.

minimize_kw : dict

Any other parameters to scipy.minimize().

Returns:

res : OptimizeResult

The optimization result represented as a OptimizeResult object (see scipy documentation). Important attributes are: x the solution array, success a Boolean flag indicating if the optimizer exited successfully and message which describes the cause of the termination.

The parameters are ordered by logHs, alpha, beta, [lnA].

curve : mrpy.likelihoods.CurveLike object

An object containing the solution parameters and methods to access relevant quantities, such as the mass function, or jacobian and hessian at the solution.

Notes

Though the option to not use bounds for the fit is available, at this point it seems to yield unpredictable results. Unless the problem at hand is so poorly known as to be impossible to set appropriate bounds, it is encouraged to use them.

Furthermore, use as stringent bounds as possible, since the algorithm explores the edges, which can induce numerical error if values far from the solution are chosen.

Examples

The most obvious example is to generate a mass function curve from given parameters:

>>> from mrpy import MRP
>>> import numpy as np
>>> m = np.logspace(10,15,500)
>>> d = MRP(m,14.0,-1.9,0.75,norm=0).dndlog10m()

Then find the best-fit parameters for the resulting data:

>>> from mrpy.fitting.fit_curve import get_fit_curve
>>> res,curve = get_fit_curve(m,d,x0=[14.2,-1.8,0.8,0.5])
>>> print res.x

We can also use the curve object to explore some of the qualities of the fit

>>> print curve.hessian
>>> from matplotlib.pyplot import plot
>>> plot(curve.logm,curve.dndm(log=True))
>>> plot(curve.logm,np.log(d))
>>> print curve.stats.mean