２群間の時系列データの解析

友人に２群間の時系列データの解析、検定を尋ねられたので
改めて勉強してみたので、健忘録として。

まず、練習用のデータを作成していきます。
症例は30例ずつの2群で、計測は3時点、処置の有無で層別化します。

id = c(1:60)
result0 = rnorm(30,40,10)
result1 = rnorm(30,25,7)
result2 = rnorm(30,20,5)
result3 = rnorm(30,41,9)
result4 = rnorm(30,23,6)
result5 = rnorm(30,19,4)
treat0 = rep("n" ,30)
treat1 = rep("t" ,30)

#　データフレームに変換しています。

dat = data.frame(id = id, treat= c(treat0,treat1), res0 = c(result0,result3),
                res1 = c(result1,result4), res2 = c(result2,result5))

#　確認します。
#　縦のカラムには患者間要因、横のカラムには患者内要因が配列されています。

View(dat)

解析していきます。

library(car)
time = factor(c("res0","res1","res2"), levels = c("res0","res1","res2"))
idata = as.data.frame(time)
model = lm(cbind(res0,res1,res2) ~ treat,ｄata=dat,contrasts = list(treat = contr.sum))
fit = Anova(model, idata = idata, idesign = ~time, type = 3)
summary(fit, multivariate = FALSE)

Univariate Type III Repeated-Measures ANOVA Assuming Sphericity

Sum Sq num Df Error SS den Df F value Pr(>F)

(Intercept) 144669 1 3147.4 58 2665.9775 <2e-16 ***

treat 25 1 3147.4 58 0.4651 0.4980

time 13532 2 5436.3 116 144.3696 <2e-16 ***

treat:time 201 2 5436.3 116 2.1448 0.1217

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Mauchly Tests for Sphericity

Test statistic p-value

time 0.89911 0.048266

treat:time 0.89911 0.048266

Greenhouse-Geisser and Huynh-Feldt Corrections

for Departure from Sphericity

GG eps Pr(>F[GG])

time 0.90835 <2e-16 ***

treat:time 0.90835 0.1268

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

HF eps Pr(>F[HF])

time 0.9360954 2.597244e-30

treat:time 0.9360954 1.252772e-01

Mauchly球面性で有意差を認めるので前半部分を確認して
treatはp=0.4980で有意差はなし、timeはp<0.001で有意差あり、
treat:time（交互作用）はp=0.1217で有意差なしです。
時間経過で変化はあるが、処置による差はなく、交互作用もない。

グラフにするのは大変なので、後日。
EZRパッケージを使用するのがおすすめです。