python:设置宽度以适应参数

2024-06-28 11:02:21 发布

您现在位置:Python中文网/ 问答频道 /正文

我一直试图用未知的拟合参数“ga”和“MA”拟合一个数据文件。我想做的是设置一个范围,“MA”的值将驻留在该范围内并拟合数据,例如,我希望MA的拟合值在范围[0.5,0.8]内,并希望保持“ga”作为任意拟合参数。我不知道怎么做。我在这里复制python代码:

#!/usr/bin/env python3

# to the data in "data_file", each line of which contains the data for one point, x_i, y_i, sigma_i.
import numpy as np
from pylab import *
from scipy.optimize import curve_fit
from scipy.stats import chi2

fname = sys.argv[1] if len(sys.argv) > 1000 else 'data.txt'
x, y, err = np.loadtxt(fname, unpack = True)
n = len(x)

p0 = [-1,1]
f = lambda x, ga, MA: ga/((1+x/(MA*MA))*(1+x/(MA*MA)))
p, covm = curve_fit(f, x, y, p0, err)
ga, MA = p

chisq = sum(((f(x, ga, MA) -y)/err)**2)
ndf = n -len(p)
Q = 1. -chi2.cdf(chisq, ndf)

chisq = chisq / ndf

gaerr, MAerr = sqrt(diag(covm)/chisq) # correct the error bars

print 'ga = %10.4f +/- %7.4f' % (ga, gaerr)
print 'MA = %10.4f +/- %7.4f' % (MA, MAerr)

print 'chi squared / NDF = %7.4lf' % chisq
print  (covm)

Tags: thefromimportdata参数lennpscipy
1条回答
网友
1楼 · 发布于 2024-06-28 11:02:21

您可以考虑使用lmfit(https://lmfit.github.io/lmfit-py)来解决这个问题。Lmfit为优化和曲线拟合提供了更高级别的接口,包括将参数视为具有边界的python对象

您的脚本可能被翻译为使用lmfit作为

import numpy as np
from lmfit import Model

fname = sys.argv[1] if len(sys.argv) > 1000 else 'data.txt'
x, y, err = np.loadtxt(fname, unpack = True)

# define the fitting model function, similar to your `f`:
def f(x, ga, ma):
    return ga/((1+x/(ma*ma))*(1+x/(ma*ma)))

# turn this model function into a Model:
mymodel = Model(f)

# now create parameters for this model, giving initial values
# note that the parameters will be *named* from the arguments of your model function:
params = mymodel.make_params(ga=-1, ma=1)

# params is now an orderded dict with parameter names ('ga', 'ma') as keys.
# you can set min/max values for any parameter:
params['ma'].min = 0.5
params['ma'].max = 2.0

# you can fix the value to not be varied in the fit:
# params['ga'].vary = False
# you can also constrain it to be a simple mathematical expression of other parameters

# now do the fit to your `y` data with `params` and your `x` data
# note that you pass in weights for the residual, so 1/err:
result = mymodel.fit(y, params, x=x, weights=1./err)

# print out fit report with fit statistics and best fit values 
# and uncertainties and correlations for variables:
print(result.fit_report())

您可以通过result.params访问最佳拟合参数;初始params不会因拟合而改变。还有一些例程可以绘制最佳拟合结果和/或残差

相关问题 更多 >