Bakgrunnur
Hér munun við sýna hvernig má framkvæma einhliða fervika greiningu fyrir endurteknar mælingar. Fyrir gögn sem ekki eru með endurteknum mælingum þá má lesa um einhliða fervikagreiningu (e. one-way ANOVA = “analysis of variance”) hér og sjá myndband hér. R sýnikóða fyrir einhliða fervikagreiningu má finna hér og fyrir tvíhliða (e. two-way) fervikagreiningu hér. Fervikagreining er notuð til að bera saman meðaltöl milli hópa.
Við getum notað einhliða (e. one-way) ANOVA fyrir endurteknar mælingar þegar við mælum alla einstaklinga við öll gildi á einni flokkabreytu, t.d. þegar við mælum alla einstaklinga á öllum tímapunktum eða eftir að hafa fengið allar gerðir af lyfjum. Þar sem hver einstaklingur hefur mælingar á öllum tímapunktum, þá kallast tímapunktar “cross-over factor”.
Þegar framkvæmd er einhliða ANOVA fyrir endurteknar mælingar þá er mikilvægt að allir einstaklingar hafi mælingar á öllum tímapunktum (í því tilfelli gefa allar aðferðirnar fyrir neðan sömu niðurstöðu). Ef einhver einstaklingur hefur ekki allar mælingar þá þarf annað hvort að nota aðrar tölfræðiaðferðir, sleppa öllum mælingum á þessum einstakling eða átta sig á nákvæmlega hvað hver aðferð prófar og nota það sem svarar rannsóknaspurningunni best.
Fyrir neðan sýnum við nokkrar leiðir til að framkvæma einhliða fervikagreiningu fyrir enduteknar mælingar í R. Síðar sýnum við hvernig má framkvæma tvíhliða fervikagreinigu fyrir endurteknar mælingar.
Gagnasett
Við lesum inn gagnasett héðan.
Við þurfum að skilgreina flokkabreytur sem “factor” þannig R líti ekki á þær sem talnabreytur.
demo2 <- read.csv("http://www.ats.ucla.edu/stat/data/demo2.csv")
demo2 <- within(demo2, {
group <- factor(group)
time <- factor(time)
id <- factor(id)
})
summary(demo2)
## id group pulse time ## 1 :3 1:12 Min. :10.00 1:8 ## 2 :3 2:12 1st Qu.:16.00 2:8 ## 3 :3 Median :22.50 3:8 ## 4 :3 Mean :22.38 ## 5 :3 3rd Qu.:26.75 ## 6 :3 Max. :35.00 ## (Other):6
str(demo2)
## 'data.frame': 24 obs. of 4 variables: ## $ id : Factor w/ 8 levels "1","2","3","4",..: 1 1 1 2 2 2 3 3 3 4 ... ## $ group: Factor w/ 2 levels "1","2": 1 1 1 1 1 1 1 1 1 1 ... ## $ pulse: int 14 19 29 15 25 26 16 16 31 12 ... ## $ time : Factor w/ 3 levels "1","2","3": 1 2 3 1 2 3 1 2 3 1 ...
1. Notum aov fallið
demo2.aov <- aov(pulse ~ time + Error(id), data = demo2)
demo2.aov<-summary(demo2.aov)
demo2.aov$`Error: Within`
## Df Sum Sq Mean Sq F value Pr(>F) ## time 2 978.2 489.1 62.02 1.1e-07 *** ## Residuals 14 110.4 7.9 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Hér að ofan sjáum við að p=1.1e-7 fyrir tengsl milli tímapunkta (time) og pulse, þ.e. pulse virðist mærktæt mismunandi milli tímapunktanna þriggja að teknu tilliti til þess að hver einstaklingur er með pulse mælingar á hverjum tímapuntki.
2. Notum Anova fallið í car pakkanum
library(car)
library(reshape2)
# reformat to wide format
demo2.wide <- dcast(demo2,id+group~time,value.var="pulse")
names(demo2.wide)[3:5] <- c("t1","t2","t3")
# prepare arguments for Anova function
# model
mod <- lm(cbind(t1,t2,t3)~1,demo2.wide)
# prepare 3-level factor for 3 time-points
time.level <- c(1,2,3)
time.factor <- as.factor(time.level)
time.frame <- data.frame(time.factor)
# run Anova
demo2.Anova <- Anova(mod, idata=time.frame,idesign=~time.factor, type=2)
## Note: model has only an intercept; equivalent type-III tests substituted.
demo2.Anova <- summary(demo2.Anova)
demo2.Anova$univariate.tests
## SS num Df Error SS den Df F Pr(>F) ## (Intercept) 12015.4 1 122.96 7 684.034 3.057e-08 *** ## time.factor 978.2 2 110.42 14 62.017 1.104e-07 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Aftur fáum við sömu niðurstöðu.
3. Notum línuleg slembiþáttalíkön (linear mixed model)
Einnig er hægt að nota línuleg slembiþáttalíkön. Athugið að hér er mikilvægt að nota sviga utan um “1|id” í formúlunni. Slembiþáttalíkön eru sérstaklega gagnleg þegar gögnin er ekki “balanced”, þ.e. til dæmis þegar vantar einhverjar mælingar þannig ekki eru jafnmargi einstaklingar með mælingar á öllum tímapunktum.
library(lme4)
form <- as.formula("pulse~time + (1|id)")
demo2.mixed <- lmer(form,data=demo2,REML = T)
x<-summary(demo2.mixed)
y<-anova(demo2.mixed)
y
## Analysis of Variance Table ## Df Sum Sq Mean Sq F value ## time 2 978.25 489.12 62.017
Við fáum sama F gildi. Það er með ráðum gert að frí gráðurnar eru ekki reiknaðar (sjá tæknilegar útskýringar hér). Við getum hins vegar borið F=62.0173561 saman við F dreifingu með frígráðurnar 2 og 14, vegna þess að við gögnin eru “balanced”, þ.e. allir 8 einstaklingar eru með mælingar á öllum 3 tímapunktum. Þá fáum við sömu niðurstöðu.
