两重复测量因素方差分析(Two Repeated Measures Factor ANOVA)一(无交互作用)——R软件实现

发布于 2023年1月2日 星期一 18:41:48 浏览:1502
原创不易,转载请注明来源,感谢!
附件下载:
两因素重复测量方差分析一(无交互作用).zip 请勿重复点击,如无响应请耐心等待或稍后再试。

在前面文章中研究的因素可能是多个,但是重复测量因素只有一个,本文介绍存在两个重复测量因素(Two Repeated Measures Factor)但无交互作用的情况及其在R软件中的实现过程。

关键词:R; 重复测量; 重复测量资料; 重复测量方差分析; 两重复测量因素; 球形检验; 交互作用; 主效应; 单独效应

一、案例介绍

评价在某麻醉药物的作用下,不同部位某蛋白酶活性变化情况。在使用该种药物后,检测12只小白鼠用药后10min、20min、30min、40min 4个不同时刻A部位和B部位蛋白酶的活性。部分数据见图1。本案例数据可从“附件下载”处下载。

图1

二、问题分析

本案例资料具有两个重复测量因素,一个是时间因素,有4个水平,即time1~time4,分别代表用药后10min、20min、30min、40min 4个不同时刻。另一个是身体不同部位,有2个水平,即A部位和B部位。因此需要使用两重复测量因素方差分析进行数据分析。需要满足以下5个条件:

条件1:观察变量唯一,且为连续变量。本研究中观察变量只有蛋白酶活性,且为连续变量,该条件满足。

条件2:有两个分析因素。本研究有时间因素time和部位因素两个因素,该条件满足。

条件3:观察变量为重复测量数据,即不满足独立性。本研究中在4个时间点时测量的蛋白酶活性均是针对同一批样本,因此不满足独立性;相同个体不同部位也不满足独立性;该条件满足。

条件4:观察变量不存在显著的异常值,该条件需要通过软件分析后判断。

条件5:各组、各水平(时间点)观察变量为正态(或近似正态)分布,该条件需要通过软件分析后判断。

关于方差齐性的问题,由于本案例不存在分组因素,因此无须检验方差齐同。

三、软件操作及结果解读

(一) 导入数据

语句见下,结果见图2。

mydata = read.csv("./两因素重复测量方差分析一(无交互作用).csv", header = T)
str(mydata)
图2
mydata$ID = factor(mydata$ID)
mydata$Site <- factor(mydata$Site, levels = c(1, 2), 
                      labels = c("A部位", "B部位"))

(二) 适用条件判断

1. 条件4判断(异常值判断)

(1) 软件操作
boxplot(mydata$time1 ~ mydata$Site, boxwex = 0.25, col = c('lightblue'),
      xlab = '', ylab = '活性', main = 'time = 10 min')
(2) 结果解读

异常值可通过箱线图和专业知识进行判断。图3为时刻1时A部位和B部位活性的箱线图,提示不存在异常值。依次可判定A、B部位其他时刻数据也不存在异常值(请读者自行操作)。

图3

2. 条件5判断(正态性检验)

(1) 软件操作
mydata1 = subset(mydata, mydata$Site == "A部位")
mydata2 = subset(mydata, mydata$Site == "B部位")
p1 = p2 = {}

for (i in 1:4) {
  p1[i] = shapiro.test(mydata1[, i + 2])[["p.value"]]
  p2[i] = shapiro.test(mydata2[, i + 2])[["p.value"]]
}
p1; p2
(2) 结果解读

    Shapiro-Wilk (夏皮罗-威尔克,S-W)正态性检验的结果见图4。可见大多数数据的P值均>0.05,提示各时间点下两部位的数据整体服从正态分布。关于正态性检验的注意事项详见推文(医学统计学核心概念及重要假设检验的软件实现(2/4)——正态性假设检验的SPSS实现)。

图4

(三) 统计描述

1. 软件操作

##宽数据转换为长数据##
library(reshape2)  #调用包“reshape2”
mydata_long <- melt(data = mydata, id.vars = c('ID', 'Site'), #保留不变的变量
                    measure.vars = c('time1', 'time2', 'time3', 'time4'), #想要转换的变量
                    variable.name = 'time', #转换后的分类变量名
                    value.name = 'tvalue') #转换后的数值变量名
## 统计描述 ##
library(plyr)
em <- ddply(mydata_long, c("time", "Site"), summarise, 
            mean = mean(tvalue, na.rm = TRUE),
            sd = sd(tvalue, na.rm = TRUE),
            n = sum(!is.na(tvalue)),
            se = sd/sqrt(n),
            lower.CI = t.test(tvalue)[["conf.int"]][1],
            upper.CI = t.test(tvalue)[["conf.int"]][2])
em
## 绘制边际均数图 ##
library(ggplot2)
ggplot(em, aes(x = time, y = mean, color = Site, group = Site)) + #绘制分组图
  geom_errorbar(aes(ymin = lower.CI, ymax = upper.CI), width = .05)+ #添加误差线
  geom_line() +  #添加连线
  geom_point()  #添加点
ggsave("emplot.png", width = 8, height = 5, dpi = 300) #保存图片

2. 结果解读

图5列出了两部位time1、time2、time3、time4时间点的“估算边际平均值”“标准误差”及均值的“95% 置信区间 (95%CI)”。估算边际均值图(图6)绘制了两部位4个时间点蛋白酶活性的变化情况,可见两组数据均有下降,但下降的幅度较为接近。

图5
图6

(四) 球形检验及方差分析

1. 软件操作

install.packages("ez") #安装包“ez”
library(ez) #调用包“ez”
ezANOVA(data = mydata_long, dv = tvalue, wid = ID,
        within = .(time, Site), detailed = T) #重复测量方差分析

2. 结果解读

图7

“莫奇来球形度检验(Mauchly’s Test for Sphericity)”中给出了球形假设检验结果(图7),可见“时间”“时间*部位”的检验结果均P<0.05,提示不满足球形假设。当违背了球形假设条件时,可采用“多变量检验”结果,或者使用单变量检验校正法的结果,包括格林豪斯-盖斯勒、辛-费德特和下线3种校正方法。单变量检验校正法一般建议采用前两种方法,当epsilon (ε)<0.75时,使用格林豪斯-盖斯勒法,当epsilon (ε)>0.75时,使用辛-费德特法。当违背了球形假设条件时,如果单变量和多变量检验结果不一致,以多变量检验结果为准。

(五) 交互作用判断

由于本案例有两个因素(一个为时间因素,另一个为部位),因此需要首先判断两个因素之间是否存在交互作用。如果交互作用有统计学意义,则需要分析单独效应。

“ANOVA (方差分析)”结果(图7)显示,time与部位(Site)之间的交互作用结果为F时间*部位=1.283,P=0.296,提示时间与部位之间交互作用无统计学意义。因此,本案例采用主效应分析结果即可。

(六) 时间效应分析

“ANOVA (方差分析)”结果(图7)显示,F时间=369.928,P<0.001,提示不同时间点蛋白酶活性差异有统计学意义。

(七) 部位效应分析

“ANOVA (方差分析)”结果(图7)显示,F部位=0.091,P=0.768,提示不同部位间蛋白酶活性差异无统计学意义。

(八) 事后检验(两两比较)

1. 软件操作

install.packages("bruceR")
library(bruceR) #调用程序包"bruceR"
library(magrittr) #调用程序包"magrittr"
library(dplyr) #调用程序包"dplyr"

MANOVA(data = mydata_long, subID = "ID", 
       dv = "tvalue", within = c("time"))  %>%
  EMMEANS("time")

2. 结果解读

由时间效应分析发可知不同时间蛋白酶活性存在差异,则需要进一步两两比较,本案例采用bonferroni法进行事后两两比较。“多重比较”结果(图8)提供了各时间点两两比较的“平均值差值及其标准误”和“显著性(校正P值)”。结果显示,随着时间的延长,time2至time4时刻与time1时刻相比,均值差逐渐增大,且差异有统计学意义。

图8

四、结论

本研究采用两重复测量因素方差分析评价在麻醉药物的作用下,A、B机体部位某蛋白酶活性变化情况。数据不存在异常值,整体数据服从正态分布;时间与部位之间交互作用无统计学意义(F时间*部位=1.283,P>0.05),故进行主效应分析。

时间效应分析显示,不同时间点蛋白酶活性差异有统计学意义(F时间=369.928,P<0.001)。两两比较发现,随着时间的延长,time2至time4时刻与time1时刻相比,均值差逐渐增大,且具有统计学差异。。部位效应分析显示,不同部位间蛋白酶活性差异无统计学意义(F部位=0.091,P=0.768)。

综上,该麻醉药对A部位B部位的蛋白酶活性影响差异无统计学意义,随着时间的延长,蛋白酶活性开始降低。

End
文章目录 沉浸式阅读