python - 如果数据有缺失值,则 scipy.mstats.theilslopes 置信限度错误

标签 python statistics scipy

如果对包含缺失值的数据集使用 scipy.mstats.theilslopes 例程,则斜率估计的下限和上限结果是不正确的。上限经常/总是(?)NaN,而下限是错误的。发生这种情况是因为 theilslopes 例程计算了已排序 slopes 数组的索引,而该数组包含不应该包含缺失值的斜率。

解决方案是在分析之前删除缺失值,但这在任何地方都没有记录。

为了演示这个问题,这里有一个简单的代码片段: 将 numpy 导入为 np 从 scipy.stats 导入 mstats

x = np.arange(12)
y = np.array([28.9, 26.2, 27.2, 26.5, 28.4, 25.3, 26.1, 24.8, 27.7,
              np.nan, np.nan, 29.6])

slope, intercept, lo_slope, up_slope = mstats.theilslopes(y, x,
                                                          alpha=0.1)
print "incorrect: ", slope, lo_slope, up_slope

idx = [0, 1, 2, 3, 4, 5, 6, 7, 8, 11]
x = x[idx]   # equivalent to pandas series.dropna()
y = y[idx]

slope, intercept, lo_slope, up_slope = mstats.theilslopes(y, x,
                                                          alpha=0.1)
print "correct: ", slope, lo_slope, up_slope

最佳答案

mstats scipy.stats 模块,“缺失值”使用 masked array 处理. nan 不表示缺失值。

下面显示了如何将数组 y(使用 nan 来表示缺失值)转换为掩码数组 my:

In [48]: x = np.arange(12)

In [49]: y = np.array([28.9, 26.2, 27.2, 26.5, 28.4, 25.3, 26.1, 24.8, 27.7, np.nan, np.nan, 29.6])

In [50]: my = np.ma.masked_array(y, mask=np.isnan(y))

In [51]: my
Out[51]: 
masked_array(data = [28.9 26.2 27.2 26.5 28.4 25.3 26.1 24.8 27.7 -- -- 29.6],
             mask = [False False False False False False False False False  True  True False],
       fill_value = 1e+20)

In [52]: slope, intercept, lo_slope, up_slope = mstats.theilslopes(my, x, alpha=0.1)

In [53]: print "correct: ", slope, lo_slope, up_slope
correct:  -0.125 -0.48 0.3875

顺便说一句,请确保您至少使用 0.15.0 版的 scipy。 theilslopes 在旧版本中有一些错误:https://github.com/scipy/scipy/pull/3574

关于python - 如果数据有缺失值,则 scipy.mstats.theilslopes 置信限度错误,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/30390131/

相关文章:

python - 使用 pandas 进行插值时如何控制 X 和 Y

python - 大阵列的密度估计

python - plotly 表达: how to control bars start position?

python - Django - 如何安装 Python 图像库 (PIL)

python - matplotlib 中的并排图

database - RabbitMQ 内存使用

python - 给定一个不均匀概率列表,至少有一个发生的概率是多少?

Python:将字符串从 UTF-8 转换为 Latin-1

r - R中二项式参数的差异

python - 如何在这里进行广义特征分解?