python - scipy.stats.chisquare 没有给出输入数据的预期结果
问题描述
我有一些数据要对其应用拟合,然后执行卡方检验以获得拟合优度。很明显,我应用的拟合不能很好地拟合数据(这本身不是问题,我不一定期望它)但是 scipy.stats.chisquare 返回的值会暗示几乎完美的契合,这显然是错误的。
到目前为止我所做的是定义一个描述我正在应用的拟合(正弦拟合)的函数,然后使用 scipy.optimize.curve_fit 通过从 popt 获取拟合参数然后在先前定义的函数来生成拟合。
然后,我将测量数据和拟合数据放入 scipy.stats.chisquare 中以尝试拟合,但返回的 p 值为 1.0,这是不正确的。我的假设是在 scipy.stats.chisquare 中使用 scipy.optimize.curve_fit 生成的值存在一些问题,但如果是这种情况,我不明白为什么会出现问题或如何解决它。
我在两个列表中有我的测量数据,我在下面称之为“时间”和“速率”
import numpy as np
import math
%matplotlib inline
import matplotlib.pyplot as plt
from statistics import stdev
import scipy
time =[309.6666666666667, 326.3333333333333, 334.6666666666667, 399.9166666666667, 416.5833333333333, 433.25, 449.91666666666663, 466.58333333333337, 483.25, 499.91666666666663,]
rate = [0.298168, 0.29317, 0.306496, 0.249861, 0.241532, 0.241532, 0.206552, 0.249861, 0.253193, 0.239867]
def oscillation(t,A,C):
return(A*np.cos((2*np.pi*(t-x0))/(t0))+C)
t0 = 365.25
A = 0.35/2
x0 = 152.5
C = 0.475
popt, pcov = curve_fit(oscillation, time, rate, p0=[A,C])
rate_fit = []
for t in time:
r = oscillation(t, popt[0],popt[1])
rate_fit.append(r)
print(scipy.stats.chisquare(rate, f_exp=rate_fit))
plt.plot(time,rate, '.')
plt.plot(time,rate_fit,'--')
上面的输出是一个拟合,在绘制时看起来确实最适合数据,但显然不是完美拟合,使得 0.999999999999458533 的另一个输出值显然是错误的
解决方案
您只适合两个参数A
和C
,从而强制相位和周期。
如果你也适合这个阶段和时期,你会得到更好的适应:
同样在这种情况下,我的 p 值为 1.0。
x0
当和固定时,您的 p 值为 1.0 的原因t0
是,您的结果与和的这些值x0
t0
最适合。强制使用这些值很可能会产生整体更差的拟合。为了比较,有x0
和t0
免费,我得到
A = -3.45840427e-02
C = 2.65142203e-01
x0 = 1.88838771e+02
t0 = 2.61112538e+02
将其与t0 = 365.25
和进行比较x0 = 152.5
。
当然,有一些(物理)原因需要修复,例如t0
一年,但在这种情况下,您不必担心情节看起来很糟糕;您的 p 值仍然考虑到这一点。
然而,更可能的原因是您也忘记ddof
了scipy.stats.chisquare
. 它的默认值是ddof=0
,这不是你所拥有的:在你的情况下它是len(rate) - 2
,在我上面的情况下,它会是len(rate) - 4
。
对于您的适合(t0
和x0
固定),这会导致p = 0.902
. 在所有参数自由的情况下,结果为 0.999887(即再次为 1)。
奖励:当我将周期固定t0
为 365.25 时输出:
A = -4.05218922e-02
C = 2.74772524e-01
x0 = 8.69008279e+01
p = 0.997
和绘制的拟合:
推荐阅读
- python-3.x - 如何在不破坏句子的情况下将句子中的每个字母替换为句子?
- selenium - 无法使用 android Emulator 在 Selenium Appium 中定位元素
- ag-grid - 访问组行 ag-grid
- java - Java:当 WeakHashMap 的条目被驱逐时,谁改变了大小?
- javascript - 在谷歌地图中画一个带有虚线边框的圆圈,应该用颜色填充
- ms-access - VBA IDE Access 2010。尝试打开即时窗口时找不到文件
- magento2 - 如何在magento 2中附加发票PDF而不是装箱单
- reactjs - 酶装载包装器的 containsAllMatchingElements 无法找到元素
- vba - 在 Excel 中快速将 Textformat 转换为 Numberformat - VBA
- node.js - 使用 Node.js 抓取和存储 Shopify 电子商务网站