Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
def test_weird_param_hints(self):
# tests Github Issue 312, a very weird way to access param_hints
def func(x, amp):
return amp*x
m = Model(func)
models = {}
for i in range(2):
m.set_param_hint('amp', value=1)
m.set_param_hint('amp', value=25)
models[i] = Model(func, prefix='mod%i_' % i)
models[i].param_hints['amp'] = m.param_hints['amp']
self.assertEqual(models[0].param_hints['amp'],
models[1].param_hints['amp'])
def fit_allxy(dataset: DataSet, initial_parameters: Optional[np.ndarray] = None) -> Dict[str, Any]:
""" Fit AllXY measurement to piecewise linear model
Args:
dataset: Dataset containing the AllXY measurement
initial_parameters: Optional set of initialization parameters
Returns:
Dictionary with the fitting results
"""
allxy_data = _default_measurement_array(dataset)
x_data = np.arange(21)
lmfit_model = Model(allxy_model)
if initial_parameters is None:
initial_parameters = _estimate_allxy_parameters(allxy_data)
param_names = lmfit_model.param_names
result = lmfit_model.fit(allxy_data, indices=x_data, **dict(zip(param_names, initial_parameters)), verbose=0)
fitted_parameters = np.array([result.best_values[p] for p in param_names])
fitted_parameters_covariance = np.diag(result.covar)
chi_squared = result.chisqr
return {'fitted_parameters': fitted_parameters, 'description': 'allxy fit', 'initial_parameters': initial_parameters, 'fitted_parameters_covariance': fitted_parameters_covariance, 'chi_squared': chi_squared}
# fitting example data to an exponential decay.
def decay(t, N, tau):
return N*np.exp(-t/tau)
###############################################################################
# The parameters are in no particular order. We'll need some example data. I
# will use ``N=7`` and ``tau=3``, and add a little noise.
t = np.linspace(0, 5, num=1000)
data = decay(t, 7, 3) + np.random.randn(*t.shape)
###############################################################################
# **Simplest Usage**
model = Model(decay, independent_vars=['t'])
result = model.fit(data, t=t, N=10, tau=1)
###############################################################################
# The Model infers the parameter names by inspecting the arguments of the
# function, ``decay``. Then I passed the independent variable, ``t``, and
# initial guesses for each parameter. A residual function is automatically
# defined, and a least-squared regression is performed.
#
# We can immediately see the best-fit values:
print(result.values)
###############################################################################
# and use these best-fit parameters for plotting with the ``plot`` function:
result.plot()
###############################################################################
# time series, regardless of length. ``pandas``
# (https://pandas.pydata.org/pandas-docs/stable/) provides tools for aligning
# indexed data. And, unlike most wrappers to ``scipy.leastsq``, ``Model`` can
# handle pandas objects out of the box, using its data alignment features.
#
# Here I take just a slice of the ``data`` and fit it to the full ``t``. It is
# automatically aligned to the correct section of ``t`` using Series' index.
model = Model(decay, independent_vars=['t'])
truncated_data = Series(data)[200:800] # data points 200-800
t = Series(t) # all 1000 points
result = model.fit(truncated_data, params, t=t)
report_fit(result.params)
###############################################################################
# Data with missing entries and an unequal length still aligns properly.
model = Model(decay, independent_vars=['t'], nan_policy='omit')
truncated_data_with_holes = Series(data_with_holes)[200:800]
result = model.fit(truncated_data_with_holes, params, t=t)
report_fit(result.params)
def fit_sine(x_data: np.ndarray, y_data: np.ndarray, initial_parameters=None,
positive_amplitude=True) -> Tuple[Dict[str, Any], Dict[str, Any]]:
""" Fit a sine wave for the inputted data; see sine function in functions.py for model
Args:
x_data: x data points
y_data: data to be fitted
initial_parameters: list of 4 floats with initial guesses for: amplitude, frequency, phase and offset
positive_amplitude: If True, then enforce the amplitude to be positive
Returns:
result_dict
"""
if initial_parameters is None:
initial_parameters = _estimate_initial_parameters_sine(x_data, y_data)
lmfit_model = Model(sine)
if positive_amplitude:
lmfit_model.set_param_hint('amplitude', min=0)
lmfit_result = lmfit_model.fit(y_data, x=x_data, **dict(zip(lmfit_model.param_names, initial_parameters)))
result_dict = extract_lmfit_parameters(lmfit_model, lmfit_result)
return result_dict['fitted_parameters'], result_dict
if l is not None:
warnings.warn('use argument lever_arm instead of l')
lever_arm = l
# initial values
linear_part, fermi_part = initFermiLinear(x_data, y_data, fig=None)
initial_parameters = linear_part + fermi_part
def fermi_linear_fitting_function(x_data, a, b, cc, A, T):
return FermiLinear(x_data, a, b, cc, A, T, l=lever_arm)
if use_lmfit:
import lmfit
gmodel = lmfit.Model(fermi_linear_fitting_function)
param_init = dict(zip(gmodel.param_names, initial_parameters))
gmodel.set_param_hint('T', min=0)
params = gmodel.make_params(**param_init)
lmfit_results = gmodel.fit(y_data, params, x_data=x_data)
fitting_results = lmfit_results.fit_report()
fitted_parameters = np.array([lmfit_results.best_values[p] for p in gmodel.param_names])
else:
fitting_results = scipy.optimize.curve_fit(
fermi_linear_fitting_function, x_data, y_data, p0=initial_parameters)
fitted_parameters = fitting_results[0]
results = dict({'fitted_parameters': fitted_parameters, 'pp': fitting_results,
'centre': fitted_parameters[2], 'initial_parameters': initial_parameters, 'lever_arm': lever_arm,
'fitting_results': fitting_results})
split (float): value that separates the up and the down level
left (array), right (array): Parameters of the left and right fitted Gaussian
"""
if initial_params is None:
initial_params = _estimate_double_gaussian_parameters(x_data, y_data)
def _double_gaussian(x, A_dn, A_up, sigma_dn, sigma_up, mean_dn, mean_up):
""" Double Gaussian helper function for lmfit """
gauss_dn = gaussian(x, mean_dn, sigma_dn, A_dn)
gauss_up = gaussian(x, mean_up, sigma_up, A_up)
double_gauss = gauss_dn + gauss_up
return double_gauss
double_gaussian_model = Model(_double_gaussian)
delta_x = x_data.max() - x_data.min()
bounds = [x_data.min() - .1 * delta_x, x_data.max() + .1 * delta_x]
double_gaussian_model.set_param_hint('mean_up', min=bounds[0], max=bounds[1])
double_gaussian_model.set_param_hint('mean_dn', min=bounds[0], max=bounds[1])
double_gaussian_model.set_param_hint('A_up', min=0)
double_gaussian_model.set_param_hint('A_dn', min=0)
param_names = double_gaussian_model.param_names
result = double_gaussian_model.fit(y_data, x=x_data, **dict(zip(param_names, initial_params)), verbose=False)
par_fit = np.array([result.best_values[p] for p in param_names])
if par_fit[4] > par_fit[5]:
par_fit = np.take(par_fit, [1, 0, 3, 2, 5, 4])
# separation is the difference between the max of the gaussians divided by the sum of the std of both gaussians
separation = (par_fit[5] - par_fit[4]) / (abs(par_fit[2]) + abs(par_fit[3]))
def _fit_1d_gaussian(self):
model = lmfit.Model(self._gauss, independent_vars=['x'])
result = model.fit(data=self._z_optim_scan_data,
x=np.linspace(self.__optim_z_scan_range[0],
self.__optim_z_scan_range[1],
self._optim_z_resolution),
amp=lmfit.Parameter('amp', value=self._z_optim_scan_data.max(), min=0),
x0=lmfit.Parameter(
'x0', value=self.__current_target['z'],
min=self.__current_target['z'] - 3 * self._optim_z_scan_range / 4,
max=self.__current_target['z'] + 3 * self._optim_z_scan_range / 4),
sigma=lmfit.Parameter('sigma',
value=self._optim_z_scan_range/2,
min=0,
max=self._optim_z_scan_range),
offset=lmfit.Parameter('offset',
value=self._z_optim_scan_data.min(),
min=0))
import matplotlib.pyplot as plt
from numpy import exp, loadtxt, pi, sqrt
from lmfit import Model
data = loadtxt('model1d_gauss.dat')
x = data[:, 0]
y = data[:, 1]
def gaussian(x, amp, cen, wid):
"""1-d gaussian: gaussian(x, amp, cen, wid)"""
return (amp / (sqrt(2*pi) * wid)) * exp(-(x-cen)**2 / (2*wid**2))
gmodel = Model(gaussian)
result = gmodel.fit(y, x=x, amp=5, cen=5, wid=1)
print(result.fit_report())
plt.plot(x, y, 'bo')
plt.plot(x, result.init_fit, 'k--', label='initial fit')
plt.plot(x, result.best_fit, 'r-', label='best fit')
plt.legend(loc='best')
plt.show()
#
"""
return (" Wrap the {} function for fitting within lmfit "
"framework\n ".format(func.__name__) + func.__doc__)
# DEFINE NEW MODELS
class ElasticModel(Model):
__doc__ = _gen_class_docs(elastic)
def __init__(self, *args, **kwargs):
super(ElasticModel, self).__init__(elastic, *args, **kwargs)
self.set_param_hint('epsilon', value=2.96, vary=False)
class ComptonModel(Model):
__doc__ = _gen_class_docs(compton)
def __init__(self, *args, **kwargs):
super(ComptonModel, self).__init__(compton, *args, **kwargs)
self.set_param_hint('epsilon', value=2.96, vary=False)
class Lorentzian2Model(Model):
__doc__ = _gen_class_docs(lorentzian2)
def __init__(self, *args, **kwargs):
super(Lorentzian2Model, self).__init__(lorentzian2, *args, **kwargs)