12  球对称检验

如何使用R语言进行球形检验？

球形检验是非常重要的内容，在进行重复测量方差分析之前，都应该先进行球形假设检验！球形检验在课本中是作为重复测量方差分析中的一个知识点出现的，并没有单独介绍，我这里单独拎出来是因为当初的学习使用R语言进行球形检验的时候，费了很大的精力（当时网上资源很少）。

表12-3重复测量数据用第四章随机区组方差分析比较处理组间差异，前提条件式要求满足“球对称”（sphericity）假设，即重复测量误差的协方差矩阵经正交对比变换后，与单位矩阵I4*4成比例。……P=0.01，拒绝球对称假设，即表12-3数据不满足随机区组设计方差分析（表12-4）“球对称”的条件，组内效应的F值需要校正。–《医学统计学》P189

R语言作为专为统计而生的语言，球形检验的函数就内置在stats包中，就是mauchly.test()，但是用法却和其他常见统计检验函数很不一样。

下面我将会为大家讲解如何使用R语言进行球形检验！

12.1 课本表12-3的数据

这是一个只有1组的！

读取数据：

df <- foreign::read.spss("datasets/表12-3重复测量ANOVA.sav",
                         to.data.frame = T, reencode = "utf-8")

str(df)
## 'data.frame':    8 obs. of  4 variables:
##  $ t0  : num  5.32 5.32 5.94 5.49 5.71 6.27 5.88 5.32
##  $ t45 : num  5.32 5.26 5.88 5.43 5.49 6.27 5.77 5.15
##  $ t90 : num  4.98 4.93 5.43 5.32 5.43 5.66 5.43 5.04
##  $ t135: num  4.65 4.7 5.04 5.04 4.93 5.26 4.93 4.48
##  - attr(*, "variable.labels")= Named chr(0) 
##   ..- attr(*, "names")= chr(0) 
##  - attr(*, "codepage")= int 936
head(df)
##     t0  t45  t90 t135
## 1 5.32 5.32 4.98 4.65
## 2 5.32 5.26 4.93 4.70
## 3 5.94 5.88 5.43 5.04
## 4 5.49 5.43 5.32 5.04
## 5 5.71 5.49 5.43 4.93
## 6 6.27 6.27 5.66 5.26

数据一共4列，就是4个时间点的血糖值。

首先将数据变为矩阵：

df <- as.matrix(df)

然后进行球对称检验（球形检验）：

mauchly.test(lm(df ~ 1), X = ~ 1)
## 
##  Mauchly's test of sphericity
##  Contrasts orthogonal to
##  ~1
## 
## 
## data:  SSD matrix from lm(formula = df ~ 1)
## W = 0.06273, p-value = 0.008207

结果就有了，就是这么简单直接，网上很多资料都是直接复制粘贴帮助文档里的内容，非常费脑子！

12.2 课本例12-3的数据

这个数据有2组！

直接读取：

df1 <- foreign::read.spss("datasets/例12-03.sav",to.data.frame = T)

str(df1)
## 'data.frame':    15 obs. of  7 variables:
##  $ No   : num  1 2 3 4 5 6 7 8 9 10 ...
##  $ group: Factor w/ 3 levels "A","B","C": 1 1 1 1 1 2 2 2 2 2 ...
##  $ t0   : num  120 118 119 121 127 121 122 128 117 118 ...
##  $ t1   : num  108 109 112 112 121 120 121 129 115 114 ...
##  $ t2   : num  112 115 119 119 127 118 119 126 111 116 ...
##  $ t3   : num  120 126 124 126 133 131 129 135 123 123 ...
##  $ t4   : num  117 123 118 120 126 137 133 142 131 133 ...
##  - attr(*, "variable.labels")= Named chr [1:7] "\xd0\xf2\xba\xc5" "\xd7\xe9\xb1\xf0" "" "" ...
##   ..- attr(*, "names")= chr [1:7] "No" "group" "t0" "t1" ...
head(df1)
##   No group  t0  t1  t2  t3  t4
## 1  1     A 120 108 112 120 117
## 2  2     A 118 109 115 126 123
## 3  3     A 119 112 119 124 118
## 4  4     A 121 112 119 126 120
## 5  5     A 127 121 127 133 126
## 6  6     B 121 120 118 131 137

数据一共7列，第1列是患者编号，第2列是诱导方法（3种），第3-7列是5个时间点的血压。

首先将数据变为矩阵，转换数据格式：

df2 <- as.matrix(cbind(df1[1:5,3:7], df1[6:10,3:7], df1[11:15,3:7]))

把测量点和分组单独建立，注意要和上面的顺序一致：

times = ordered(rep(1:5,3))
group = factor(rep(c("A","B","C"),each = 5))

然后进行球对称检验（球形检验）：

mauchly.test(lm(df2 ~ 1), M = ~ group + times, X = ~ times)
## 
##  Mauchly's test of sphericity
##  Contrasts orthogonal to
##  ~times
## 
##  Contrasts spanned by
##  ~group + times
## 
## 
## data:  SSD matrix from lm(formula = df2 ~ 1)
## W = 0.427, p-value = 0.279

真的是有点费事儿！不过现在很多R包都可以在进行重复测量方差分析时自动给出球形检验的结果，已经方便多了。