正态性假设 (Normality)：
方差齐性假设 (Homoscedasticity)：
独立性假设 (Independence)：
单样本t检验 (One-sample t-test): 用于比较单个样本的均值与已知的总体均值之间的差异。
独立样本t检验 (Independent samples t-test): 用于比较两个独立样本的均值差异。
配对样本t检验 (Paired samples t-test): 用于比较同一组受试者在两个不同条件下的均值差异。
ggplot() + xlim(-10,10) +
geom_function(mapping = aes(color="normal Distribution"),
fun = dnorm, args = list(mean = 0, sd = 1),
)+
geom_function(mapping = aes(color="t Distribution"),
fun = dt, args = list(df = 1 ,ncp=0),
)+
scale_color_manual(values = c("normal Distribution" = "red",
"t Distribution" = "blue"))+
labs(color = "Distribution")
t 检验（Student‘s t test），主要用于小样本（n<30），标准差未知的正态分布总体。在进行t检验之前，可以先通过正态性检验 shapiro.test()
在数据符合正态分布的前提下使用单样本t-test来比较一组样本的均值和已知(理论/总体)均值，所谓的已知均值能来自于之前的实验数据或者理论值。根据研究问题(原假设)的不同又分为双尾(不等)和单尾检验(大于或者小于)
\[t=\frac{\bar X-\mu_0}{S /\sqrt n} \sim t(\nu=n-1) \]
它是一种参数检验，用于检验样本均值是否可以合理地为总体均值或特定值。
x<- dplyr::filter(PlantGrowth,group=="ctrl")
summary(x$weight)
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 4.170 4.550 5.155 5.032 5.293 6.110
shapiro.test(x$weight)
#>
#> Shapiro-Wilk normality test
#>
#> data: x$weight
#> W = 0.95668, p-value = 0.7475
t.test(x$weight, mu=5)
#>
#> One Sample t-test
#>
#> data: x$weight
#> t = 0.17355, df = 9, p-value = 0.8661
#> alternative hypothesis: true mean is not equal to 5
#> 95 percent confidence interval:
#> 4.614882 5.449118
#> sample estimates:
#> mean of x
#> 5.032t_test_two_sided <- function(data,mu,level=0.95,...){
bar_X <- mean(data,na.rm = T)
sd <- sd(data, na.rm = T)
n=length(data)
se=sd/sqrt(n)
t_statistic <- (bar_X-mu)/se
p_value <- 2*(1-pt(abs(t_statistic),df=n-1))
t_critical <- qt((1-(1-level)/2),df = n-1)
CI_lower <- bar_X - t_critical * se
CI_upper <- bar_X+ t_critical * se
CI <- paste0(level*100,"%置信区间: ","[",CI_lower,",",CI_upper,"]",sep = "")
output <- list(
均值=bar_X,
标准差=sd,
标准误=se,
t=t_statistic,
df=n-1,
p_value=p_value,
CI=CI
)
return(output)
}
t_test_two_sided(x$weight,mu = 5)
#> $均值
#> [1] 5.032
#>
#> $标准差
#> [1] 0.5830914
#>
#> $标准误
#> [1] 0.1843897
#>
#> $t
#> [1] 0.1735455
#>
#> $df
#> [1] 9
#>
#> $p_value
#> [1] 0.8660633
#>
#> $CI
#> [1] "95%置信区间: [4.61488155565504,5.44911844434496]"\[ H_0:\mu_{\bar d}=0 \] \[ t=\frac{\bar d- \mu_{\bar d}}{S_d /\sqrt n} \sim t(\nu) \] 其中\(d= X_2-X_1,\mu_{\bar d}=0\)。
df <- tribble(
~id,~baseline,~ten_days_later,~d,
1,58.27,120.61,62.34,
2,59.51,126.33,66.82,
3,53.84,108.35,54.51,
4,54.70,139.99,85.29,
5,54.03,115.29,61.26,
6,61.29,146.96,85.67,
7,54.72,115.64,60.92,
8,70.43,124.62,54.19,
9,66.45,121.40,54.95,
10,59.31,134.81,75.50,
11,63.48,130.73,67.25,
12,67.19,118.37,51.18,
13,52.92,129.28,76.36,
14,71.99,117.40,45.41
)
d_mean <- mean(df$d)
d_sd <- sd(df$d)
d_se <- d_sd/sqrt(length(df$d))
t_statistic <- d_mean/d_se
n <- 14
# p值
p_value <- 2 * (1 - pt(abs(t_statistic),df = n-1 ))
# 查找95%置信水平下的t分布的临界值
t_critical <- qt(0.975, df=n-1)
# 计算95%置信区间
CI_lower <- d_mean - t_critical * d_se
CI_upper <- d_mean + t_critical * d_se
# 输出结果
cat("95% Confidence Interval:", CI_lower, "to", CI_upper)
#> 95% Confidence Interval: 57.20275 to 71.6044
shapiro.test(df$baseline)
#>
#> Shapiro-Wilk normality test
#>
#> data: df$baseline
#> W = 0.91419, p-value = 0.1813
shapiro.test(df$ten_days_later)
#>
#> Shapiro-Wilk normality test
#>
#> data: df$ten_days_later
#> W = 0.96676, p-value = 0.8306
shapiro.test(df$d)
#>
#> Shapiro-Wilk normality test
#>
#> data: df$d
#> W = 0.94337, p-value = 0.4632
t.test(df$ten_days_later,df$baseline,paired = TRUE)
#>
#> Paired t-test
#>
#> data: df$ten_days_later and df$baseline
#> t = 19.322, df = 13, p-value = 5.866e-11
#> alternative hypothesis: true mean difference is not equal to 0
#> 95 percent confidence interval:
#> 57.20275 71.60440
#> sample estimates:
#> mean difference
#> 64.40357
#t.test(Pair(df$ten_days_later,df$baseline)~1,data=df)\[H_0:\mu_1-\mu_2=0\]
\[ t=\frac{(\bar X_1-\bar X_2)-(\mu_1-\mu_2)}{S_{\bar X_1-\bar X_2}}=\frac{\bar X_1-\bar X_2}{S_C\sqrt{(\frac{1}{n_1}+\frac{1}{n_2})}} \]
其中，\(S_c^2=\frac{(n_1-1)S_1^2+(n_2-1)S_2^2}{n_1+n_2-2}\)
df2 <- tibble(
experimental=c(120.61 ,126.33 ,108.35 ,139.99 ,115.29 ,146.96 ,115.64,
124.62 ,121.40 ,134.81 ,130.73 ,118.37 ,129.28 ,117.45),
control=c(58.23 ,54.50 ,59.47 ,59.64 ,53.77 ,43.48 ,
54.63 ,71.91 ,53.97 ,49.72 ,61.26 ,78.17,NA,NA)
)
e_mean <- mean(df2$experimental)
e_sd <- sd(df2$experimental)
n1 <- length(df2$experimental)
ctrl_mean <- mean(df2$control,na.rm = TRUE)
ctrl_sd <- sd(df2$control,na.rm = TRUE)
n2 <- length(df2$control)-sum(is.na(df2$control))
Sc_2 <- ((n1-1)*e_sd^2+(n2-1)*ctrl_sd^2)/(n1+n2-2)
t2 <- (e_mean-ctrl_mean)/sqrt(Sc_2*(1/14+1/12))
t.test(df2$experimental,df2$control,var.equal = TRUE)
#>
#> Two Sample t-test
#>
#> data: df2$experimental and df2$control
#> t = 16.967, df = 24, p-value = 7.215e-15
#> alternative hypothesis: true difference in means is not equal to 0
#> 95 percent confidence interval:
#> 58.63784 74.87954
#> sample estimates:
#> mean of x mean of y
#> 124.98786 58.22917\[H_0:\frac{\sigma_1^2}{\sigma_2^2}=1\]
\[ F=\frac{\frac {(n_1-1)S_1^2}{\sigma_1^2}/(n_1-1)}{\frac {(n_2-1)S_2^2}{\sigma_2^2}/(n_2-1)}=(\frac{S_1^2}{S_2^2})(\frac{\sigma_2^2}{\sigma_1^2})=\frac{S_1^2}{S_2^2}\sim F(\nu_1,\nu_2),\nu_1=n_1-1,\nu_2=n_2-1 \]
# 检验两个样本的方差是否相等
var.test(df2$experimental,df2$control)
#>
#> F test to compare two variances
#>
#> data: df2$experimental and df2$control
#> F = 1.287, num df = 13, denom df = 11, p-value = 0.6831
#> alternative hypothesis: true ratio of variances is not equal to 1
#> 95 percent confidence interval:
#> 0.3794599 4.1152566
#> sample estimates:
#> ratio of variances
#> 1.287025F_stats <- (e_sd^2)/(ctrl_sd^2)
F_stats
#> [1] 1.287025
# 计算p值
p_value <- 1 - pf(F_stats, df1=13, df2=11, lower.tail =F)
# p_value <- pf(F_stats, df1=13, df2=11)
alpha <- 0.05
confidence_level <- 1 - alpha
# 计算F分布的临界值
f_critical_lower <- qf((1 - confidence_level) / 2, df1 = 13, df2 = 11)
f_critical_upper <- qf(confidence_level, df1 = 13, df2 = 11)
# 计算方差比率的置信区间
ci_lower <- sqrt(f_critical_lower * (ctrl_sd^2 / e_sd^2))
ci_upper <- sqrt(f_critical_upper * (ctrl_sd^2 / e_sd^2))
# 输出结果
cat("95% CI for variances ratio:", ci_lower, "to", ci_upper)
#> 95% CI for variances ratio: 0.4929485 to 1.464781\[H_0:\mu_1-\mu_2=0\]
\[ t'=\frac{(\bar X_1-\bar X_2)-(\mu_1-\mu_2)}{\sqrt{\frac{S_1^2}{n_1}+\frac{S_2^2}{n_2}}}=\frac{\bar X_1-\bar X_2}{\sqrt{\frac{S_1^2}{n_1}+\frac{S_2^2}{n_2}}}\sim t(\nu ') \]
强调方差的变异
因为\(t'\)既不遵循t分布，也不遵循正态分布，因此t’的临界值需要特定的计算方法。
\[ t'_{\alpha/2}=\frac{S^2_{\bar X_1}t_{\nu_1,\alpha/2}+S^2_{\bar X_2}t_{\nu_2,\alpha/2}}{S^2_{\bar X_1}+S^2_{\bar X_2}} \]
\[ t'_{1-\alpha/2}=\frac{S^2_{\bar X_1}t_{\nu_1,1-\alpha/2}+S^2_{\bar X_2}t_{\nu_2,1-\alpha/2}}{S^2_{\bar X_1}+S^2_{\bar X_2}} \] 其中\(\nu_1=n_1-1,\nu_2=n_2-1,t_{\nu_1,1-\alpha/2}和t_{\nu_2,1-\alpha/2}\)分别是\(t_{\nu_1}和t_{\nu_2}\)的临界值。
因为t分布是对称的，\(t_{\nu,\alpha/2}=-t_{\nu,1-\alpha/2}\)，所以\(t'_{\alpha/2}=t'_{1-\alpha/2}\)。
强调自由度
\[ v'=\frac{(\frac{S_1^2}{n_1}+\frac{S_2^2}{n_2})^2}{\frac{(\frac{S_1^2}{n_1})^2}{n_1-1}+\frac{(\frac{S_2^2}{n_2})^2}{n_2-1}}(舍入到最近整数) \]
当\(n_1＞30\ \&\ n_2＞30\)时，
\[ \bar X_1\dot\sim N(\mu_1,\frac{\sigma^2_1}{n_1}) \]
\[ \bar X_2\dot\sim N(\mu_2,\frac{\sigma^2_2}{n_2}) \]
\[H_0:\mu_1-\mu_2=0\]
\[ Z=\frac{(\bar X_1-\bar X_2)-(\mu_1-\mu_2)}{\sqrt{\frac{\sigma^2_1}{n_1}+\frac{\sigma^2_2}{n_2}}}=\frac{\bar X_1-\bar X_2}{\sqrt{\frac{\sigma^2_1}{n_1}+\frac{\sigma^2_2}{n_2}}}\dot\sim N(0,1) \]
实际应用中，总体方差未知，使用t 检验，提供了对总体方差不确定性的自然估计。
功效(power) \(1-β\approx\Phi(-z_{1-\alpha/2}+\frac{|\mu_1-\mu_2|}{\sqrt {\sigma_1^2/n_1+\sigma_2^2/n_2}})\)
样本量
两组样本量相等 \[ n=\frac{(\sigma_1^2+\sigma_1^2)(z_{1-\alpha/2}+z_{1-\beta})^2}{(\mu_1-\mu_2)^2} \]
两组样本量不等（\(n_2=kn_1\)） \[ n_1=\frac{(\sigma_1^2+\sigma_1^2/k)(z_{1-\alpha/2}+z_{1-\beta})^2}{(\mu_1-\mu_2)^2} \]
\[ n_2=\frac{(k\sigma_1^2+\sigma_1^2)(z_{1-\alpha/2}+z_{1-\beta})^2}{(\mu_1-\mu_2)^2} \]