mrpy.fitting.fit_curve.get_fit_expected¶

mrpy.fitting.fit_curve.get_fit_expected(data_m, data_mf, x0, bounds=None, jac=True, V0=1, kappa=None, **minimize_kw)¶

Generate the expected MLE of MRP parameters, given a mass function from which the samples should be drawn.

This integrates over the input mass function, so as to naturally weight each mass scale as if fitting an actual sample. It assumes that there are no uncertainties on the mass estimates.

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)

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_expected
>>> 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