StatQuest 笔记分享

Hello 我是rvdsd，前一段时间出于分析数据的需要，发现缺乏统计学知识，和着菜鸟团的小伙伴学了一阵子统计学，把StatQuest这个系列的视频给总结了一下，具体的可以参见这篇文章（https://mp.weixin.qq.com/s/cvA40tLtpIzb_z2xtLD3ig），当时这个视频里缺一个t检验的内容，现在补上并顺带分享给大家，另外，由于公众号无法很好的呈现公式，我写的内容也可以在我的博客上查看（rvdsd.top）

t分布背景

t分布最先是由戈塞特（Gosset）提出来，戈塞特在1899年进入了一个家啤酒公司，这个啤酒公司出版的一个东西比他们的啤酒还有名，就是《吉尼斯世界纪录》，他在研究酵母发酵时，提出了t检验，t检验是用了解决小样本的统计学问题，但是当时吉尼斯公司不让自己的员工发表文章，因此戈塞特就以匿名的方式发表文章，笔名就是叫学生（Student），因此戈塞特提出的这种分布就叫t分布，这种检验方法就叫t检验（为什么取最后一个字母，我不太清楚），更多有趣的统计学故事可以看一下《女士品茶》这本书。

什么是t分布

从正态分布中取出来的样本均值也服从正态分布，在统计推断中往往希望利用它减去总体的均值再除以均值的总体标准差来得到服从标准正态分布的数据，但是，我们通常情况下并不知道总体的均值。

在这个变换过程中，如果用样本标准差来代替其未知的总体标准差时，即用下面的这个公式：

mark

来代替

mark

此时，得到的分布就不再是标准正态分布了，它的密度曲线有点像正态分布，不过中间瘦，尾巴长，这种分布称为t分布（t-distribution，或称为学生分布，Student's t），不同的样本量通过标准化所产生的t分布也不同，这样就形成了一族分布，t分布族中的成员是以自由度来区分的，这里的自由度等于样本量减去1，t分布通常用于研究小样本的均值，如果自由度趋于无穷，那么t分布就是标准正态分布了，一个有k个自由度的t分布用t(k)表示，下图就展示了标准正态分布N(0,1)和自由度等于1的t(1)分布的密度函数曲线，我们可以看出，t分布两边尾巴比较长，如下所示：

mark

通常使用tα表示t分布相应于右侧尾概率α的t变量的α上侧分位数，即对t分布变量T，有P(T>tα=α），下图就表示了自由度为2的t(2)分布右边的尾概率（α=0.05），如下所示：

mark

t分布的密度分布图

set.seed(1000)
x <- seq(-5,5,length.out=1000)
y <- dt(x,1,0)
plot(x,y,col="red",xlim=c(-5,5),ylim=c(0,0.5),type="l",
xaxs="i",yaxs="i",ylab="density",xlab="",
main="The T density distribution")
lines(x,dt(x,5,0),col="green")
lines(x,dt(x,5,2),col="blue")
lines(x,dt(x,50,4),col="orange")
legend("topleft",legend=paste("df=",c(1,5,5,50),"ncp=",c(0,0,2,4)),lwd=1,col=c("red","green","blue","orange"))

t分布的累积分布图

set.seed(1)
x<-seq(-5,5,length.out=1000)
y<-pt(x,1,0)
plot(x,y,col="red",xlim=c(-5,5),ylim=c(0,0.5),type='l',
xaxs="i", yaxs="i",ylab='density',xlab='',
main="The T Cumulative Distribution Function")
lines(x,pt(x,5,0),col="green")
lines(x,pt(x,5,2),col="blue")
lines(x,pt(x,50,4),col="orange")
legend("topleft",legend=paste("df=",c(1,5,5,50)," ncp=", c(0,0,2,4)), lwd=1, col=c("red", "green","blue","orange"))

mark

某一数据是否符合t分布

Kolmogorov-Smirnov连续分布检验: 检验单一样本是不是服从某一预先假设的特定分布的方法。以样本数据的累计频数分布与特定理论分布比较，若两者间的差距很小，则推论该样本取自某特定分布族。
该检验原假设为H0:数据集符合t分布，H1:样本所来自的总体分布不符合t分布。令F0(x)表示预先假设的理论分布，Fn(x)表示随机样本的累计概率(频率)函数.
统计量D为: D=max|F0(x) - Fn(x)|
D值越小，越接近0，表示样本数据越接近T分布
p值，如果p-value小于显著性水平α(0.05)，则拒绝H0

set.seed(1)
S<-rt(1000, 1,2)
ks.test(S, "pt", 1, 2)

计算结果如下所示：

> ks.test(S, "pt", 1, 2)

    One-sample Kolmogorov-Smirnov test

data:  S
D = 0.02526, p-value = 0.5461
alternative hypothesis: two-sided

其中p-value = 0.5461，说明数据集S符合自由度为1，标准差为2的T分布

t检验

T检验，亦称student t检验（Student's t test），主要用于样本含量较小（例如n<30），总体标准差σ未知的正态分布资料，它比较的就是服从t分布的两组数据是否有差异，在R中t检验的函数为t.test()，参数如下所示：

t.test(x, ...)

t.test(x, y = NULL,alternative = c("two.sided", "less", "greater"),mu = 0, paired = FALSE, var.equal = FALSE,      conf.level = 0.95, ...)

参数说明
其中x, y是我们样本数据构成的向量，alternative表示备择假设的形式，“two.sided”表示双边检验（H1: 1不等于μ2）， “less”表示单边检验（H1:μ1小于μ2），“greater”表示单边检验（H1:μ1大于μ2），paired表示样本是否是配对的，var.equal表示两总体方差是否相等，默认是不等。

案例运用

配对t检验

例3-6为比较两种方法对乳酸饮料中脂肪含量测定结果是否不同，随机抽取了10份乳酸饮料制品，分别用脂肪酸水解法和哥特里－罗紫法测定其结果如表3-5第(1)~(3)栏。问两法测定结果是否不同？（《医学统计学》（第四版，孙振球）

mark

现在看一下计算过程：

raw_data <- read.csv("https://raw.githubusercontent.com/20170505a/raw_data/master/data_szq_306.csv",sep="\t")
# 原始数据已经放到github上

t.test(raw_data$result[1:10]-raw_data$result[11:20],mu=0)

计算结果如下所示：

> t.test(raw_data$result[1:10]-raw_data$result[11:20],mu=0)

    One Sample t-test

data:  raw_data$result[1:10] - raw_data$result[11:20]
t = 7.926, df = 9, p-value = 2.384e-05
alternative hypothesis: true mean is not equal to 0
95 percent confidence interval:
 0.1946542 0.3501458
sample estimates:
mean of x 
   0.2724 

结果与书中相同，p值小于0.01，但我觉得这个方法不一定要用配对t检验，普通的t检验也能胜任。

双样本t检验

例3-7  为研究国产四类新药阿卡波糖胶囊的降血糖效果，某医院用40名II型糖尿病病人进行同期随机对照试验。研究者将这些病人随机等分到试验组(用阿卡波糖胶囊)和对照组(用拜唐苹胶囊)，分别测得试验开始前和8周时的空腹血糖，算得空腹血糖下降值见表3-6，能否认为该国产四类新药阿卡波糖胶囊与拜唐苹胶囊对空腹血糖的降糖效果不同？（《医学统计学》（第四版，孙振球）

mark

raw_data <- read.csv("https://raw.githubusercontent.com/20170505a/raw_data/master/data_szq_307.csv",,header=T,sep=",")# 原始数据已经放到了github上

a <- as.numeric(raw_data[,1])
b <- as.numeric(raw_data[,2])
t.test(a,b)

运算结果如下所示：

    Welch Two Sample t-test

data:  a and b
t = -0.64187, df = 36.086, p-value = 0.525
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -2.32926  1.20926
sample estimates:
mean of x mean of y 
    2.065     2.625 

结果p值为0.525，因此尚不能认为阿卡波糖胶囊与拜唐苹胶囊对空腹血糖的降糖效果不同，需要注意的是，这个结果中的显著的题目是Welch Two Sample t-test，这是因为这两个样本的方差不齐，在R中，如果两个样本的方差不齐，t.test()会使用修正后的t检验（也就是Welch检验），这是默认使用的方法，如果方差齐，需要在t.test()中设置var.equal=TRUE。

t检验的局限

虽然t检验使用起来很方便，但是如果不加以注意，很容易误用，经常误用的情况包括：

1. 不考虑数据的正态性，只要是两组比较就直接使用t检验（如果不符合正态性，就要采用Wilcoxon检验）。

2. 将t检验用于多组实验设计中的两两比较，增加假阳性错误（此时应该使用ANOVA）；

3. 不考虑资料是否独立，采用独立资料的t检验分析。

参考资料

萨尔斯伯格, 邱东. 女士品茶[M]. 中国统计出版社, 2004.

吴喜之. 统计学:从数据到结论(第四版)[J]. 中国统计, 2013(6):2.

孙振球. 医学统计学.第4版[M]. 人民卫生出版社, 2014. 

---

猜你喜欢

还有更多文章，请移步公众号阅读

▼ 如果你生信基本技能已经入门，需要提高自己，请关注下面的生信技能树，看我们是如何完善生信技能，成为一个生信全栈工程师。

▼ 如果你是初学者，请关注下面的生信菜鸟团，了解生信基础名词，概念，扎实的打好基础，争取早日入门。