实验数据的假设检验小记
假设检验的基本知识
往往， 我们根据假设， 提出原理。 而要验证这一假设是否成立， 就需要做实验。 然后运用统计学的方法来说明假设成立与否。
这种利用统计学原理来验证假设的方法就称为假设检验。
一般而言， 假设检验分为参数检验与 非参数检验， 当总体的分布类型已知时， 我们只需确定相应的参数。 举个例子：某糖厂用自动包装机装糖， 每袋糖的标准重量为1斤， 而根据长期经验知道， 每袋的重量服从正态分布($N(\mu,0.015)$), 现在我们想知道某台机器是否工作正常， 那么就可以随机抽取50袋该机器生产的糖，测量他们的重量， 通过统计学的推导得出结论。这里的假设就是:
H0:$\mu=1$斤
从而要检验的就是这个参数$\mu$。
理解了上述参数假设的背景， 那么非参数假设就不难理解了：即我们并不知道总体的分布， 这时需要根据实验数据提供的信息对总体分布的种种假设进行检验，即要检验总体是否服从某种假设的分布。 由于一般的书上只推导了正态分布下的各种参数检验， 因此我们在进行$F$-检验，$T$-检验时， 千万别忘了第一件事就是要做非参数检验以确定实验数据是否服从正态分布。
这里， 一个自然的问题是，是不是任何一个分布都可以和某个正态分布等价， 即根据已知的随机变量（他不服从正态分布， 但分布已知）， 我们能否经过某种一一变换得到一个新的随机变量， 使得它是服从正态分布的？ 如果是这样， 我们自然只需要做正态分布的各种参数检验了。
但是， 问题的回答似乎没这么简单。 一般而言，正态分布是许多统计方法的理论基础：如t分布、F分布、x2分布都是在正态分布的基础上推导出来的，u检验也是以正态分布为基础的。此外，t分布、二项分布、Poisson分布的极限为正态分布，在一定条件下，可以按正态分布原理来处理。
非参数假设检验
作为例子， 我们将以正态分布为例。
| A | B |
|---|---|
| 0.54 | 0.05 |
| 2.20 | 0.13 |
| 0.58 | 0.17 |
| 0.55 | 0.20 |
| 0.27 | 0.23 |
| 0.37 | 0.28 |
| 0.77 | 0.29 |
| 0.24 | 0.29 |
| 0.70 | 0.30 |
| 0.92 | 0.30 |
| 0.36 | 0.38 |
| 1.10 | 0.40 |
| 0.47 | 0.43 |
| 0.54 | 0.43 |
| 0.62 | 0.44 |
| 0.56 | 0.45 |
| 0.45 | 0.47 |
| 0.60 | 0.50 |
| 0.58 | 0.52 |
| 0.51 | 0.53 |
| 1.68 | 0.54 |
| 0.79 | 0.54 |
| 0.88 | 0.55 |
| 0.67 | 0.55 |
| 1.20 | 0.57 |
| 0.70 | 0.67 |
| 0.90 | 0.68 |
| 0.66 | 0.73 |
| 0.77 | 0.75 |
| 1.03 | 0.75 |
| 0.39 | 0.76 |
| 0.64 | 0.87 |
| 0.45 | 0.88 |
| 0.51 | 0.90 |
| 0.59 | 0.93 |
| 0.42 | 0.98 |
| 0.55 | 1.00 |
| 0.52 | 1.07 |
| 1.15 | 1.13 |
| 1.28 | 1.16 |
| 0.73 | 1.20 |
| 1.05 | 1.20 |
| 1.60 | 1.21 |
| 0.48 | 1.28 |
| 0.99 | 1.30 |
| 0.87 | 1.31 |
| 0.91 | 1.32 |
| 3.10 | 1.32 |
| 0.25 | 1.35 |
| 0.96 | 1.36 |
首先，我们可以作直方图观察下是否符合正态分布的图形， 主要是看图形的对称性以及是否是单峰。如果不对称, 说明有部分数据受到限制(例如测量手段), 不能得到. 若非单峰, 说明有其他因素的影响. 而且这种影响是不可忽略的.
直方图
作直方图之前， 我们需要一些统计描述量：最大值， 最小值， 样本容量。为了和标准的正态分布对比，我们还需要均值 和 方差。 获取这些量的快捷方法就是使用Excel的分析统计工具:描述统计。
由于Excel默认没有打开分析统计工具， 我们需要作如下设置(以Office 2013为例)：
单击文件—>选项—>来到Excel的选项， 单击左边的加载项， 这时在右边第一行会有分析工具库， 单击以选中它—>点击下面的转到—>来到加载宏面板， 单击以选中分析工具库 并 确定。这时你可看到Excel数据子面板的最后多了一项数据分析。
这样，我们就打开了Excel的数据分析工具。
下面我们来熟悉下这个分析工具库， 构造作直方图的第一步所需要的基本统计信息。
统计描述
- 单击数据面板下的 数据分析， 选择描述统计
- 在输入区域内：选定待分析数据区域， 可以选择多行或者多列
- 选择你将要分析的数据是以列还是以行的方式呈现的
- 如果所选择的区域有数据标识， 勾选标识
- 在输出区域选择输出到当前工作表还是新的工作表
- 勾选汇总统计（输出基本统计量）
- 确定便可生成统计描述
效果如图所示(我每次都将数据格式设置为保留两位小数):
这里A,B两列是数据, 我没有完全截图. D-G列是汇总统计.
组数的计算
我们已将作直方图的第一步做好了, 即有了基本的统计量.
回忆, 做直方图时还需要分组, 然后计算每组有多少个数据(即频数), 然后才是作图. 分组首先需要确定组数, 即到底要分多少组. 这个事实上是非常重要的, 分组过少, 数据就非常集中; 分组过多, 数据就非常分散, 这都会掩盖分布的特征.
关于组数的多少, 我查了下wiki, 总结起来, 其方法有如下几种, 这里假设$n=50$为样本总数:
- 方根: 数据个数的平方根并取整, Excel公式:
=INT(SQRT(n)), 这里结果为7 - Sturges 公式:数据个数关于$2$的对数取整加$1$, Excel公式:
=INT(LOG(n,2))+1, 这里结果为6 - Rice 准则: 数据个数的三次方根的二倍取整, Excel公式
=INT(2*n^(1/3)), 这里结果为6 - Doane 公式: 是Sturges公式的修正, 对处理异常数据有效. 计算公式为:
$$
1 + \log_2( n ) + \log_2 \left( 1 + \frac { |g_1| }{\sigma_{g_1}} \right) ,
$$
其中
$$
\sigma_{g_1} = \sqrt { \frac { 6(n-2) }{ (n+1)(n+3) } }
$$
而$g_1$为该分布的三阶矩偏度估计. 具体的计算可以在这个网页找到. 只需把你的数据(这里, 例如取第一列) 粘贴到网页的1 Copy&Paste your data下面的输入框, 然后单击 2. Please click here to start java applet 以加载java, 最后点击 3.Press Here to calculate the optimal bin size 我们就得到优化后的直方图. 在图上可以看到# of bins 10表示此时我们优化后的组数是10, 后面一个Bin width:0.286表示的是组距.
wiki上还有几种办法, 我这里直接跳过. 因为基本上利用第四种我们可以得到比较有效的组数.
第四中办法的结果如下图(相应于第一列数据):
最后, 若你还不满意, 可以调整下面的滑块, 来增加或者减少组数.
频数的统计
在Excel中作直方图, 我们需要建立仓库(bins), 即每组的分点. 明显我们只需把数据的最大值和最小值之间的区间分为如上所计算的组数(这里是10), 并求其分点即可. 具体地, 在Excel中操作如下(为了具有一般性, 10在下面换为n):
- 参考第一个图, 在H1中输入:
=E12, 并按F4切换到绝对引用, 这是因为A列的第一个分点就是最小值. - 接下来我们将用数据填充的办法计算后面的值, 于是只需设计好第二个公式, 在H2中中输入:
=H2+$E$11/10, 意思很明显, 就是前一个数据(H2)加上组距, 这里组距等于区域(E11)除以组数10. 注意由于组数在下面填充中是不动的, 因此是绝对引用. - 填充直到一共有11个数为止(10个区间, 11个分点,呵呵~~), 这也可以通过最后一个值应该是最大值3.1(E13单元格给出)判断.
- 至此, A列数据的仓库就做好了, 完全类似的做B列数据的. 注意这里组数是需要从新计算的, 这里利用组数计算的方法4结果有些奇葩, 竟然组数是2, 需要改进, 一般5-8比较正常, 这里选择6, 因为前面三种的计算有两种是6, 一种结果是7.
做好后的效果如下图(一如既往地, 数据只保留两位小数):
直方图的作法
现在我们可以开始作直方图了. 一般有两种办法, 利用Excel的FREQUENCY手动计算频数, 然后添加直方图. 这种作法的好处是由于是按照公式来计算的, 因此数据改变时直方图会随着改变. 不好的地方是不够自动化, 需要的计算比较多. 因此我只介绍如何计算频数.
- 在J1中输入
A的频数作为标题 - 在J2中输入频数计算公式:
=FREQUENCY(A2:A51,H2:H12)道理是很明显的,A2:A51表示待统计的数据,H2:H12表示仓库分组, 此时返回的是第一个仓库的频数. - 为了计算其他仓库的频数, 我们只需选中
J2:J12, 再按F2, 最后按Ctrl+Shift+Enter即可得到其他仓库的频数.
下面我介绍如何通过数据分析的直方图工具来快捷的生成直方图.
- 点击数据面板, 点击最后的数据分析, 在弹出的数据分析选框中, 选择直方图并确定.
- 在输入区域中选择A列数据
$A$1:$A$51, 在接收区域中选择:$H$1:$H$12, 由于选择的数据中有表头, 我们需要勾选标志. 在输出选项中选择输出区域, 我们这里输出到当前工作表的$D$18. 最后勾选图表输出. 并确定 - 请手动调整直方图的大小
至此A列数据的直方图已经作好. 效果如图所示:
趋势线的添加
为了对比与标准正态分布的关系, 我们在上面的直方图中添加正态分布.
首先, 我们需要计算正态分布的密度.
- 在
F18中输入正态图作为标题 - 在
F19中输入公式:=NORM.DIST(D19:D29,$E$3,$E$7,FALSE)计算第一个分点的正态分布之概率密度. 其中之所以选择D19:D29因为我们待会要用数组计算得到后面分点的概率密度. 而$E$3是A列数据的均值,$E$7是A列数据的标准差, 最后的FALSE表示我们计算的是概率密度而非累积密度.
总起来说, 就是计算A列数据作为正态分布在给定分点的概率密度是多少. 这样我们就可以做出A列数据的正态分布分布图了. - 计算其他分点的概率密度, 这可由前面讲过的数组计算公式得到: 选中
F19:F29, 按下F2, 再按Ctrl+Shift+Enter即可
至此, A列数据的正态分布概率密度就计算好了. 下面我们将通过它做出A列数据的正态分布图, 并将其添加到原来的直方图中.
- 在直方图的空白处单击以选中该直方图, 此时会看到设计面板, 单击以切换到该面板, 然后单击选择数据, 我们想要将正态图的数据添加到当前图表
- 单击添加(A), 在系列名称里面输入
=Sheet1!$F$18, 就是正态图的标题; 在系列值里选择=Sheet1!$F$19:$F$29(这里数据的添加都是点击后面带颜色的小框框, 然后再选择数据范围), 就是上面计算的正态分布的概率密度. - 两次确定, 即可在图表中看到添加进来的正态分布图, 但是此时比例不对(是按照实际数据的值来作图的, 即一个是频数最大为19, 而一个为概率密度, 最大为0.79), 这可以通过设置次坐标轴来调整.
- 同样选中直方图(空白处单击), 切换到视图(事实上在图表空白处双击即可快速切换到视图面板), 单击更改图表类型, 在所有图表中选择组合, 可以看到频率系列默认的是簇状柱形图, 这个不用修改, 在正态图系列, 默认的是折线图, 我们将其修改为X Y(散点图)的第三个带平滑线的散点图, 并勾选后面的次坐标轴, 确定即可得到带正态分布的直方图.
最终效果如图所示:
作为练习, 请自行对B列数据作出带正态图的直方图. 可以参考我最终的Example.xlsx文件.
至此, 我们完成了非参数检验的第一步: 通过直方图大致判断是否符合正态分布.
结论是(对A列数据): 可以看到大致还是符合正态分布的, 只是图形不对称, 右边有个长尾(右偏态分布), 需要反思实验时是什么造成了低值的缺失. (当然这可能不是实验手段造成的, 有可能这正好是你的实验本身需要数据反映的情况).
$\chi^2$拟合检验
$\chi^2$拟合检验并不只是针对正态分布来检验的. 事实上, 观察上图(直方图的左边), 我们知道A列数据在每个仓库的频数, 例如: D19=0.24, E19=1, 表示的是A列数据中小于等于0.24的有1个, D20=0.53,E20=13表示大于0.24小于等于0.53的数据有13个, 与此类推, 最后的D30=其他, E30=1表示, 大于3.10的有1个.
现在, 假设A列数据是服从正态分布的, 那么我们可以计算出各个仓库的理论频数np_i. 而英国统计学家K.Pearson证明了统计量
$$\begin{equation}
\chi^2=\sum_{i=1}^k\frac{(v_i-np_i)^2}{np_i},
\label{eq1}
\end{equation}$$
在样本容量充分大时是近似服从自由度为$k-1$的$\chi^2$分布. 这里, $v_i$就是实际的频数, 上表的F19:F29, 而$n$是仓库总数, 即组数, $p_i$是随机位于该仓库的概率, 这是需要根据假设的分布(不一定就是正态分布)的密度函数来计算的.
这样, 我们不难理解, $\chi$检验就是检验上面这个假设, 即假设A列数据是服从正态分布的, 是否成立的一个标准. 此外, 它还可以用来做独立性检验, 请关注我另一篇文章
由此可见, 这里计算的$p$值表示: 我们接受实验数据是服从正态分布(即拒绝$H_0$), 但是实际上实验数据却不是服从正态分布的概率. 也即我们犯第一类错误错误的概率. 因此当$p$值越小时, 我们越有理由相信我们的数据是服从正态分布的. 对于其他的分布, 唯一的区别就是各个仓库的概率密度不同. 我最后计算的结果如下图: 具体的可以参考我做的的Example.xlsx文档. 由于Excel的卡方统计事实上只是为独立性检验而设计的, 从帮助文档可以看到它的自由度在这里不正确, 因为根据帮助文档, 自由度当$r=1$(即行row)自由度为$c-1$, 即统计的列数减1, 这样, 这里的计算是11; 但是实际的情况却是, 由于正态分布中占了两个自由度(均值和标准差), 因此自由度应该是$12-2-1=9$, 这样, 我们必须重新计算统计量$\chi^2$. 这样, 我们得到的结果是否, 因为卡方统计量要大于临界值. 从而根据统计推断(在给定的显著性水平下), 数据是服从正态分布的. 效果如图所示: 最终的Excel文档, Example.xlsx
更正
D3中输入公式:=(B3-C3)^2/C3, 计算第一个仓库的频率与期望之差的平方除以期望D14, 这就计算好了卡方统计量中的每一项=SUM(D3:D14)D16中求得自由度, =COUNT(B3:B14)-2-10.05B18中计算指定显著水平下的临界值, 计算公式为:=CHISQ.INV(1-B17,D16), 这里公式返回的是当$\chi^2$的累积概率密度达到$1-0.05$时的值, 后面的D16表示的是自由度, 因为这不是Excel默认的, 所有需要手动计算B19中用一个IF来判断是否接受原假设:=IF(B16<B18,"是","否")
发表回复