Korreleeritud vaatlustega tavaline lineaarne mudel

Siin alaosas peatume mudelitel, kus vaatlused ei ole sõltumatud, ent teeme seda otse vigade dispersioonimaatriksi peadiagonaalilt väljaspool asuvate elementide väärtuste abil jättes esmalt mängust välja juhuslikud faktorid.

Analüüside läbiviimiseks saab kasutada paketi nlme funktsiooni gls. Lisaks mudelile peame siin ette andma ka korrelatsioonistruktuuri. Selleks on argument correlation, mille võimalikud väärtused on corAR1, corARMA, corCAR1, corCompSymm, corExp, corGaus, corLin, corRatio, corSpher ja corSymm (täpsemalt saab igaühe kohta lugeda vastava nimega abifailist, alljärgevalt vaid mõned näited).

Kui eeldame, et sama objekti kõik vaatlused on võrdselt korreleeritud siis:

> m4=gls(weight~Diet, data=ChickWeight, correlation=corCompSymm(form=~1|Chick))
> summary(m4)
Generalized least squares fit by REML
Model: weight ~ Diet
Data: ChickWeight
AIC BIC logLik
6518.837 6544.953 -3253.418

Correlation Structure: Compound symmetry
Formula: ~1 | Chick
Parameter estimate(s):
Rho
0.06327688

Coefficients:
Value Std.Error t-value p-value
(Intercept) 101.63352 6.040438 16.825521 0.0000
Diet2 20.98315 10.236494 2.049838 0.0408
Diet3 41.31648 10.236494 4.036195 0.0001
Diet4 33.40412 10.269957 3.252606 0.0012

Correlation:
(Intr) Diet2 Diet3
Diet2 -0.590
Diet3 -0.590 0.348
Diet4 -0.588 0.347 0.347

Standardized residuals:
Min Q1 Med Q3 Max
-1.4953600 -0.7571530 -0.1840143 0.5809295 3.3093560

Residual standard error: 69.51504
Degrees of freedom: 578 total; 574 residual
> anova(m4)
Denom. DF: 574
numDF F-value p-value
(Intercept) 1 1046.5371 <.0001
Diet 3 6.7413 2e-04

Tundub aga mõistlikum eeldada, et sama objekti ajas lähemal paiknevad mõõtmised on ka tugevamalt korreleeritud:

>m5=gls(weight~Diet, data=ChickWeight, correlation=corCAR1(form=~Time|Chick))
> summary(m5)

Generalized least squares fit by REML
Model: weight ~ Diet
Data: ChickWeight
AIC BIC logLik
5170.612 5196.728 -2579.306

Correlation Structure: Continuous AR(1)
Formula: ~Time | Chick
Parameter estimate(s):
Phi
0.989545

Coefficients:
Value Std.Error t-value p-value
(Intercept) 98.95588 20.36231 4.859756 0.0000
Diet2 27.77655 35.16827 0.789818 0.4300
Diet3 54.64690 35.16827 1.553869 0.1208
Diet4 35.77562 35.18489 1.016790 0.3097

Correlation:
(Intr) Diet2 Diet3
Diet2 -0.579
Diet3 -0.579 0.335
Diet4 -0.579 0.335 0.335

Standardized residuals:
Min Q1 Med Q3 Max
-1.1994450 -0.5597528 -0.1501226 0.4208601 2.2962349

Residual standard error: 95.5465
Degrees of freedom: 578 total; 574 residual
> anova(m5)
Denom. DF: 574
numDF F-value p-value
(Intercept) 1 91.17876 <.0001
Diet 3 0.90562 0.4381

Mainimist väärib ka, et lisaks ajaliselt lähedastele andmetele võib tegu olla ka ruumiliselt lähedaste andmetega. Sellisel juhul kirjeldab korreleeritust tavaliselt kaks tunnust (koordinaadid tasandil), ent argumendile correlation (valikud corExp kuni corSpher) antaksegi sel juhul kaks tunnust (eraldatuna + märgiga) ning lisaks saab valida ka kauguse määramise metoodika.

Lisaks tasub mainida, et vahel võib tegu olla ka heteroskedastiliste andmetega -- mõne indiviidi andmed on varieeruvamad (suurema dispersiooniga) kui teise omad ning sellisel juhul on seda vaja ka arvestada. Selleks on funktsioonil gls argument weights valikutega varComb, varConstPower, varExp, varFixed, varIdent ja varPower.

Eelnevates punktides toodud funktsioonid on ka gls poolt tehtud mudelite korral reeglina kasutatavad. Lisanduvad mitmed mõningad täiendused graafikute kujutamisel. Näiteks

> plot(m5,residuals(.)~fitted(.)|Chick)

kujutab vastava hajuvusdiagrammi (mudeli jäägid vs mudeli prognoosid) iga objekti kohta eraldi.