Python - 最小化卡方

2023-12-29

我一直在尝试通过最小化卡方来将线性模型拟合到一组应力/应变数据。不幸的是,使用下面的代码并不能正确地最小化chisqfunc功能。它正在寻找初始条件下的最小值,x0,这是不正确的。我已经浏览过scipy.optimize文档并测试了最小化其他正常工作的功能。您能否建议如何修复下面的代码,或者建议我可以使用另一种方法通过最小化卡方来将线性模型拟合到数据?

import numpy
import scipy.optimize as opt

filename = 'data.csv'

data = numpy.loadtxt(open(filename,"r"),delimiter=",")

stress = data[:,0]
strain = data[:,1]
err_stress = data[:,2]

def chisqfunc((a, b)):
    model = a + b*strain
    chisq = numpy.sum(((stress - model)/err_stress)**2)
    return chisq

x0 = numpy.array([0,0])

result =  opt.minimize(chisqfunc, x0)
print result

感谢您阅读我的问题,任何帮助将不胜感激。

干杯,威尔

编辑:我当前使用的数据集:链接到数据 https://drive.google.com/file/d/0B1YR_XKF3wtmYl9pNHJHWEZFU3M/edit?usp=sharing


问题是你最初的猜测与实际的解决方案相差甚远。如果你在里面添加一条打印语句chisqfunc() like print (a,b),然后重新运行你的代码,你会得到类似的结果:

(0, 0)
(1.4901161193847656e-08, 0.0)
(0.0, 1.4901161193847656e-08)

这意味着minimize仅在这些点评估函数。

如果你现在尝试评估chisqfunc()在这 3 对值中,您会发现它们完全匹配,例如

print chisqfunc((0,0))==chisqfunc((1.4901161193847656e-08,0))
True

发生这种情况是因为对浮点运算进行舍入。换句话说,在评估时stress - model, 变量stress比大了太多数量级model,并且结果被截断。

然后,人们可以尝试暴力破解它,提高浮点精度,并编写data=data.astype(np.float128)加载数据后loadtxt. minimize失败,与result.success=False,但有一条有用的消息

由于精度损失,不一定能实现所需的误差。

一种可能性是提供更好的初始猜测,以便在减法中stress - model the model一部分具有相同的数量级,另一部分重新调整数据,使解决方案更接近您最初的猜测(0,0).

It is MUCH如果您只是重新缩放数据,例如相对于特定应力值(例如这种材料的屈服/开裂),则变得无量纲会更好

这是一个拟合示例,使用最大测量应力作为应力标度。您的代码几乎没有变化:

import numpy
import scipy.optimize as opt

filename = 'data.csv'

data = numpy.loadtxt(open(filename,"r"),delimiter=",")

stress = data[:,0]
strain = data[:,1]
err_stress = data[:,2]


smax = stress.max()
stress = stress/smax
#I am assuming the errors err_stress are in the same units of stress.
err_stress = err_stress/smax

def chisqfunc((a, b)):
    model = a + b*strain
    chisq = numpy.sum(((stress - model)/err_stress)**2)
    return chisq

x0 = numpy.array([0,0])

result =  opt.minimize(chisqfunc, x0)
print result
assert result.success==True
a,b=result.x*smax
plot(strain,stress*smax)
plot(strain,a+b*strain)

Your linear model is quite good, i.e. your material has a very linear behaviour for this range of deformation (what material is it anyway?): enter image description here

本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

Python - 最小化卡方 的相关文章

随机推荐