R语言常用统计假设检验

0

T 检验 t-test

# One Sample t-test
turtle_weights <- c(300, 315, 320, 311, 314, 309, 300, 308, 305, 303, 305, 301, 303)

t.test(x = turtle_weights, mu = 310)

# Two Sample t-test in R
sample1 <- c(300, 315, 320, 311, 314, 309, 300, 308, 305, 303, 305, 301, 303)
sample2 <- c(335, 329, 322, 321, 324, 319, 304, 308, 305, 311, 307, 300, 305)

t.test(x = sample1, y = sample2)

# Paired Samples t-test
before <- c(22, 24, 20, 19, 19, 20, 22, 25, 24, 23, 22, 21)
after <- c(23, 25, 20, 24, 18, 22, 23, 28, 24, 25, 24, 20)

t.test(x = before, y = after, paired = TRUE,var.equal=FALSE)

中介效应分析 Sobel Test

A Sobel test is a method of testing the significance of a mediation effect.
Sobel Test

library(bda)
mv <- rnorm(50)
iv <- rnorm(50)
dv <- rnorm(50)
mediation.test(mv,iv,dv)

单比例Z检验 One Proportion Z-Test

A one proportion z-test is used to compare an observed proportion to a theoretical one.

  • If n ≤ 30: binom.test(x, n, p = 0.5, alternative = “two.sided”)
  • If n> 30: prop.test(x, n, p = 0.5, alternative = “two.sided”, correct=TRUE)
prop.test(x=64, n=100, p=0.60, alternative="two.sided")

单、双样本Z检验 One Sample & Two Sample Z-Tests

library(BSDA)

# One-sample z-Test
data = c(88, 92, 94, 94, 96, 97, 97, 97, 99, 99,
         105, 109, 109, 109, 110, 112, 112, 113, 114, 115)


z.test(data, mu=100, sigma.x=15)

# Two-sample z-Test
cityA = c(82, 84, 85, 89, 91, 91, 92, 94, 99, 99,
         105, 109, 109, 109, 110, 112, 112, 113, 114, 114)

cityB = c(90, 91, 91, 91, 95, 95, 99, 99, 108, 109,
         109, 114, 115, 116, 117, 117, 128, 129, 130, 133)


z.test(x=cityA, y=cityB, mu=0, sigma.x=15, sigma.y=15)

Fisher精确检验 Fisher’s Exact Test

Fisher’s Exact Test is used to determine whether or not there is a significant association between two categorical variables.

data <-  matrix(c(2,5,9,4), nrow = 2)
fisher.test(data)

F检验 F-Test

An F-test is used to test whether two population variances are equal.

# method 1
x <- c(18, 19, 22, 25, 27, 28, 41, 45, 51, 55)
y <- c(14, 15, 15, 17, 18, 22, 25, 25, 27, 34)
var.test(x, y)

# method 2
data <- data.frame(values=c(18, 19, 22, 25, 27, 28, 41, 45, 51, 55,
                            14, 15, 15, 17, 18, 22, 25, 25, 27, 34),
                   group=rep(c('A', 'B'), each=10))

var.test(values~group, data=data)

Wilcoxon检验 Wilcoxon Signed-Rank Test

The Wilcoxon Signed-Rank Test is the non-parametric version of the paired t-test. It is used to test whether or not there is a significant difference between two population means when the distribution of the differences between the two samples cannot be assumed to be normal.

before <- c(14, 17, 12, 15, 15, 9, 12, 13, 13, 15, 19, 17, 14, 14, 16)
after <- c(15, 17, 15, 15, 17, 14, 9, 14, 11, 16, 18, 20, 20, 10, 17)

wilcox.test(before, after, paired=TRUE)

Mann-Whitney U Test

A Mann-Whitney U test (sometimes called the Wilcoxon rank-sum test) is used to compare the differences between two independent samples when the sample distributions are not normally distributed and the sample sizes are small (n <30).

new <- c(3, 5, 1, 4, 3, 5)
placebo <- c(4, 8, 6, 2, 1, 9)

wilcox.test(new, placebo)

JB检验 Jarque-Bera Test

The Jarque-Bera test is a goodness-of-fit test that determines whether or not sample data have skewness and kurtosis that matches a normal distribution.

library(tseries)
dataset <- runif(100)
jarque.bera.test(dataset)

AD检验 Anderson-Darling Test

An Anderson-Darling Test is a goodness of fit test that measures how well your data fit a specified distribution. This test is most commonly used to determine whether or not your data follow a normal distribution.

library(nortest)
x <- rnorm(100, 0, 1)
ad.test(x)

沃尔德检验 Wald Test

A Wald test can be used to test if one or more parameters in a model are equal to certain values.
eg:
​mpg=\beta_{0}+\beta_{1} \text { disp }+\beta_{2} \text { carb }+\beta_{3} h p+\beta_{4} c y l

library(aod)
model <- lm(mpg ~ disp + carb + hp + cyl, data = mtcars)
summary(model)
wald.test(Sigma = vcov(model), b = coef(model), Terms = 3:4)

KPSS检验 KPSS Test

A KPSS test can be used to determine if a time series is trend stationary.

library(tseries)
set.seed(1)
data<-rnorm(100)
plot(data, type='l')
kpss.test(data, null="Trend")

Mann-Kendall趋势检验 Mann-Kendall Trend Test

A Mann-Kendall Trend Test is used to determine whether or not a trend exists in time series data. It is a non-parametric test, meaning there is no underlying assumption made about the normality of the data.

library(Kendall)
data(PrecipGL)
MannKendall(PrecipGL)

plot(PrecipGL)
lines(lowess(time(PrecipGL),PrecipGL), col='blue')

KW检验 Kruskal-Wallis Test

A Kruskal-Wallis Test is used to determine whether or not there is a statistically significant difference between the medians of three or more independent groups.

df <- data.frame(group=rep(c('A', 'B', 'C'), each=10),
                 height=c(7, 14, 14, 13, 12, 9, 6, 14, 12, 8,
                          15, 17, 13, 15, 15, 13, 9, 12, 10, 8,
                          6, 8, 8, 9, 5, 14, 13, 8, 10, 9))

kruskal.test(height ~ group, data = df) 

方差比例检验 Variance Ratio Test

A variance ratio test is used to test whether or not two population variances are equal.

group1 <- c(5, 6, 6, 8, 10, 12, 12, 13, 14, 15, 15, 17, 18, 18, 19)
group2 <- c(9, 9, 10, 12, 12, 13, 14, 16, 16, 19, 22, 24, 26, 29, 29)
var.test(group1, group2)

Bartlett检验 Bartlett’s Test

Bartlett’s Test of Sphericity compares an observed correlation matrix to the identity matrix. Essentially it checks to see if there is a certain redundancy between the variables that we can summarize with a few number of factors.

library(psych)
set.seed(0)
data <- data.frame(A = rnorm(50, 1, 4), B = rnorm(50, 3, 6), C = rnorm(50, 5, 8))
cor_matrix <- cor(data)
cortest.bartlett(cor_matrix, n = nrow(data))

合并标准差 Pooled Standard Deviation

​SD_{pool} = \sqrt{\left(n_{1}-1\right) s_{1}{ }^{2}+\left(n_{2}-1\right) s_{2}{ }^{2} /\left(n_{1}+n_{2}-2\right)}

data1 <- c(6, 6, 7, 8, 8, 10, 11, 13, 15, 15, 16, 17, 19, 19, 21)
data2 <- c(10, 11, 13, 13, 15, 17, 17, 19, 20, 22, 24, 25, 27, 29, 29)

s1 <- sd(data1)
s2 <- sd(data2)

n1 <- length(data1)
n2 <- length(data2)

sd_pooled <- sqrt(((n1-1)*s1^2 + (n2-1)*s2^2) / (n1+n1-2))

加权标准差 Weighted Standard Deviation

library(Hmisc)
df <- data.frame(team=c('A', 'A', 'A', 'A', 'A', 'B', 'B', 'C'),
                 wins=c(2, 9, 11, 12, 15, 17, 18, 19),
                 points=c(1, 2, 2, 2, 3, 3, 3, 3))

wt <- c(1, 1, 1.5, 2, 2, 1.5, 1, 2)

sapply(df[c('wins', 'points')], function(x) sqrt(wtd.var(x, wt)))

Levene检验 Levene test

Levene’s Test is used to tests whether or not the variance among two or more groups is equal.

library(car)
set.seed(0)

data <- data.frame(program = rep(c("A", "B", "C"), each = 30),
                   weight_loss = c(runif(30, 0, 3),
                                   runif(30, 0, 5),
                                   runif(30, 1, 7)))

leveneTest(weight_loss ~ program, data = data)

配对卡方检验 McNemar’s Test

McNemar’s Test is used to determine if there is a statistically significant difference in proportions between paired data.

data <- matrix(c(30, 12, 40, 18), nrow = 2,
    dimnames = list("After Video" = c("Support", "Do Not Support"),
                    "Before Video" = c("Support", "Do Not Support")))
mcnemar.test(data)

格拉布斯检验 Grubbs’ Test

Grubbs’ Test is a statistical test that can be used to identify the presence of outliers in a dataset. To use this test, a dataset should be approximately normally distributed and have at least 7 observations.

library(Outliers)
data <- c(5, 14, 15, 15, 14, 13, 19, 17, 16, 20, 22, 8, 21, 28, 11, 9, 29, 40)
grubbs.test(data)

grubbs.test(data, opposite=TRUE)

二项分布检验 Binomial Test

binom.test(9, 24, 1/6) #  two-tailed 
binom.test(11, 30, 0.5, alternative="less") # left-tailed
binom.test(46, 50, 0.8, alternative="greater") # right-tailed 

中位数检验 Mood’s Median Test

Mood’s Median Test is used to compare the medians of two or more independent groups.

library(coin)
method = rep(c('method1', 'method2'), each=10)
score = c(75, 77, 78, 83, 83, 85, 89, 90, 91, 97, 77, 80, 84, 84, 85, 90, 92, 92, 94, 95)
examData = data.frame(method, score)
median_test(score~method, data = examData)

游程检验 Runs Test

Runs test is a statistical test that is used to determine whether or not a dataset comes from a random process.

library(randtests)
data <- c(12, 16, 16, 15, 14, 18, 19, 21, 13, 13)
runs.test(data)

正态性检验 Test for Normality

  • 法1:直方图
set.seed(0)
normal_data <- rnorm(200)

non_normal_data <- rexp(200, rate=3)

par(mfrow=c(1,2)) 
hist(normal_data, col='steelblue', main='Normal')
hist(non_normal_data, col='steelblue', main='Non-normal')
  • 法2:Q-Q图
set.seed(0)
normal_data <- rnorm(200)
non_normal_data <- rexp(200, rate=3)


par(mfrow=c(1,2)) 

qqnorm(normal_data, main='Normal')
qqline(normal_data)

qqnorm(non_normal_data, main='Non-normal')
qqline(non_normal_data)
  • 法3:SW检测 Shapiro-Wilk Test
set.seed(0)
normal_data <- rnorm(200)
shapiro.test(normal_data)
  • 法4:KS检测 Kolmogorov-Smirnov Test
set.seed(0)
normal_data <- rnorm(200)
ks.test(normal_data, 'pnorm')
  • 法5:CV检测 Cramer-Von Mises Test
library(goftest)
set.seed(0)
normal_data <- rnorm(200)
cvm.test(data, 'pnorm')

多元正态性检验 Multivariate Normality Tests

  • 法1:Mardia’s Test
library(QuantPsyc)
set.seed(0)
data <- data.frame(x1 = rnorm(50),
                   x2 = rnorm(50),
                   x3 = rnorm(50))
mult.norm(data)$mult.test
  • 法2:Energy Test
library(energy)
set.seed(0)
data <- data.frame(x1 = rnorm(50),
                   x2 = rnorm(50),
                   x3 = rnorm(50))
mvnorm.etest(data, R=100)

相关性检验 Correlation Test

x <- c(2, 3, 3, 5, 6, 9, 14, 15, 19, 21, 22, 23)
y <- c(23, 24, 24, 23, 17, 28, 38, 34, 35, 39, 41, 43)
cor.test(x, y)

卡方独立性检验 Chi-Square Test of Independence

A Chi-Square Test of Independence is used to determine whether or not there is a significant association between two categorical variables.

data <- matrix(c(120, 90, 40, 110, 95, 45), ncol=3, byrow=TRUE)
colnames(data) <- c("Rep","Dem","Ind")
rownames(data) <- c("Male","Female")
data <- as.table(data)

chisq.test(data)

卡方拟合优度检验 Chi-Square Goodness of Fit Test

A Chi-Square Goodness of Fit Test is used to determine whether or not a categorical variable follows a hypothesized distribution.

observed <- c(50, 60, 40, 47, 53) 
expected <- c(.2, .2, .2, .2, .2) 
chisq.test(x=observed, p=expected)

似然比检验 Likelihood Ratio Test

library(lmtest)
#fit full model
model_full <- lm(mpg ~ disp + carb + hp + cyl, data = mtcars)
#fit reduced model
model_reduced <- lm(mpg ~ disp + carb, data = mtcars)

lrtest(model_full, model_reduced)
# H0: The full model and the nested model fit the data equally well.

Cramer’s V系数 Cramer’s V

Cramer’s V is a measure of the strength of association between two nominal variables.

library(rcompanion)
data = matrix(c(7,9,12,8), nrow = 2)
cramerV(data, ci = TRUE)

# 多变量
data = matrix(c(6, 9, 8, 5, 12, 9), nrow = 2)
cramerV(data, ci = TRUE)

​\Phi系数 Phi Coefficient

A Phi Coefficient (sometimes called a mean square contingency coefficient) is a measure of the association between two binary variables.

data = matrix(c(4, 8, 9, 4), nrow = 2)
phi(data)

基尼系数 Gini Coefficient

Gini coefficient is a way to measure the income distribution of a population.

library(DescTools)
x <- c(50, 50, 70, 70, 70, 90, 150, 150, 150, 150)
Gini(x, unbiased=FALSE)

# 指定频数
x <- c(10, 20, 25, 55, 70, 90, 110, 115, 130)
n <- c(6, 7, 7, 14, 22, 20, 8, 4, 1)
Gini(x, n, unbiased=FALSE)

邹检验 Chow Test

A Chow test is used to test whether the coefficients in two different regression models on different datasets are equal.

library(strucchange)
data <- data.frame(x = c(1, 1, 2, 3, 4, 4, 5, 5, 6, 7, 7, 8, 8, 9, 10, 10,
                         11, 12, 12, 13, 14, 15, 15, 16, 17, 18, 18, 19, 20, 20),
                   y = c(3, 5, 6, 10, 13, 15, 17, 14, 20, 23, 25, 27, 30, 30, 31,
                         33, 32, 32, 30, 32, 34, 34, 37, 35, 34, 36, 34, 37, 38, 36))

sctest(data$y ~ data$x, type = "Chow", point = 10)

格兰杰因果检验 Granger-Causality Test

The Granger Causality test is used to determine whether or not one time series is useful for forecasting another.

library(lmtest)
data(ChickEgg)
grangertest(chicken ~ egg, order = 3, data = ChickEgg)
grangertest(egg ~ chicken, order = 3, data = ChickEgg)

# H0: Time series x does not Granger-cause time series y

巴特莱特检验 Bartlett’s Test

Bartlett’s test is a statistical test that is used to determine whether or not the variances between several groups are equal.

df <-data.frame(group = rep(c('A','B', 'C'), each=10),
                score = c(85, 86, 88, 75, 78, 94, 98, 79, 71, 80,
                          91, 92, 93, 85, 87, 84, 82, 88, 95, 96,
                          79, 78, 88, 94, 92, 85, 83, 85, 82, 81))
bartlett.test(score ~ group, data = df)

对数秩检验 Log Rank Test

library(survival)
head(ovarian)
survdiff(Surv(futime, fustat) ~ rx, data=ovarian)
# H0: There is no difference in survival between the two groups.