[R-sig-eco] covariance structure issue with lme
Charles Determan Jr
deter088 at umn.edu
Wed May 16 18:19:01 CEST 2012
Greetings,
I have previously tried the R mixed model mailing list but unfortunately I
remain
stuck with this question. I certainly hope you will be able to shed some
light on this particular problem. I am interested in ecology
investigations,
however this is just practice data to try and accomplish the modelling
parameters.
My sincere thanks for an assistance.
I have previously attempted to replicate the UN (unstructured),
CS (compound symmetry), and AR(1) (Autoregression)
covariance structures used in SAS PROC MIXED with the lme
function in this nlme package. However, my efforts
have fallen short on replicating the *Variance Components (VC)*
structure. I have read that it is also known as a diagonal structure.
Below I have copied over all the models I have tried and their output
with no success. Perhaps someone here will see my error or something
I have overlooked. I have attached the data for this particular
model. Thanks to all, I certainly cannot thank this help list enough.
If you need any further information/clarification, please ask.
Cheers, Charles
library(nlme)
dat=read.table("C:/subset.csv"
,sep=",",header=TRUE, na.strings=".")
attach(dat)
dat34=dat[Group %in% c("3", "4"),]
attach(dat34)
liver34=within(dat34, {
Group=factor(Group)
Event_name=factor(Event_name)
Died=factor(Died)
ID=factor(ID)
})
attach(liver34)
##What is should be from SAS
#CV var
#Type 3 Tests of Fixed Effects
#Effect NumDF DenDF F Value Pr > F
#Group 1 22 0.25 0.6244
#Died 1 22 6.55 0.0179
#Group*Died 1 22 4.43 0.0470
fit.1=lme(var~Group*Died,
random=~1|ID,
data=dat34)
anova(fit.1, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 227.58700 <.0001
#Group 1 22 0.18320 0.6728
#Died 1 22 3.63388 0.0698
#Group:Died 1 22 3.04103 0.0951
fit.2=lme(var~Group*Died,
data=dat34,
random=~1|ID/Died)
anova(fit.2, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 77.99004 <.0001
#Group 1 22 1.46275 0.2393
#Died 1 22 5.84535 0.0243
#Group:Died 1 22 3.04103 0.0951
fit.3=lme(var~Group*Died,
random=list(ID=pdSymm(~Event_name)),
data=dat34)
anova(fit.3, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 273.10918 <.0001
#Group 1 22 0.69692 0.4128
#Died 1 22 1.43316 0.2440
#Group:Died 1 22 5.74399 0.0255
fit.4=lme(var~Group*Died,
random=list(ID=pdSymm(~Group)),
data=dat34)
anova(fit.4, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 235.13889 <.0001
#Group 1 22 0.15878 0.6941
#Died 1 22 3.83253 0.0631
#Group:Died 1 22 3.01222 0.0966
fit.5=lme(var~Group*Died,
random=list(ID=pdSymm(~Group)),
data=dat34,
weights=varIdent(form=~1|Event_name))
anova(fit.5, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 277.16705 <.0001
#Group 1 22 0.23901 0.6298
#Died 1 22 3.99283 0.0582
#Group:Died 1 22 3.23135 0.0860
fit.6=lme(var~Group*Died,
random=list(ID=pdSymm(~Group)),
data=dat34,
weights=varIdent(form=~1|Event_name))
anova(fit.6, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 277.16705 <.0001
#Group 1 22 0.23901 0.6298
#Died 1 22 3.99283 0.0582
#Group:Died 1 22 3.23135 0.0860
fit.7=lme(var~(Group*Died),
random=list(ID=pdCompSymm(~Died)),
data=dat34)
anova(fit.7, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 85.83799 <.0001
#Group 1 22 1.60624 0.2183
#Died 1 22 4.71795 0.0409
#Group:Died 1 22 2.65379 0.1175
fit.8=lme(var~(Group*Died),
data=dat34,
random=~1|ID,
corr=corSymm())
anova(fit.8, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 119.54403 <.0001
#Group 1 22 4.58972 0.0435
#Died 1 22 8.01715 0.0097
#Group:Died 1 22 5.27470 0.0315
fit.9=lme(var~(Group*Died),
data=dat34,
random=list(ID=pdDiag(~Group*Died)),
corr=corSymm(, ~1|ID))
# Error in lme.formula(var ~ (Group * Died), data = dat34, random =
list(ID = pdDiag(~Group * :
# nlminb problem, convergence error code = 1
# message = iteration limit reached without convergence (9)
fit.10=lme(var~(Group*Died),
data=dat34,
random=list(ID=pdDiag(~Group*Died)),
corr=NULL)
anova(fit.10, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 93.90211 <.0001
#Group 1 22 1.75311 0.1991
#Died 1 22 6.84379 0.0158
#Group:Died 1 22 3.11458 0.0915
fit.11=lme(var~Group*Died,
data=dat34,
random=list(ID=pdDiag(~Group*Died)))
anova(fit.11, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 93.90211 <.0001
#Group 1 22 1.75311 0.1991
#Died 1 22 6.84379 0.0158
#Group:Died 1 22 3.11458 0.0915
fit.12=lme(var~Group*Died,
data=dat34,
random=list(ID=pdDiag(~Event_name)))
anova(fit.12, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 87.33040 <.0001
#Group 1 22 1.09661 0.3064
#Died 1 22 5.46329 0.0289
#Group:Died 1 22 2.94589 0.1001
fit.13=lme(var~Group*Died,
data=dat34,
random=list(ID=pdDiag(~Group)))
anova(fit.13, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 77.99004 <.0001
#Group 1 22 1.46275 0.2393
#Died 1 22 5.84535 0.0243
#Group:Died 1 22 3.04103 0.0951
fit.14=lme(var~Group*Died,
data=dat34,
random=list(ID=pdDiag(~Died)))
anova(fit.14, type="marginal", adjustSigma=F)
# numDF denDF F-value p-value
#(Intercept) 1 101 85.83800 <.0001
#Group 1 22 1.60624 0.2183
#Died 1 22 4.71795 0.0409
#Group:Died 1 22 2.65379 0.1175
fit.15=lme(var~Group*Died,
data=dat34,
random=~1|ID,
corr=corCompSymm())
anova(fit.15, type="marginal", adjustSigma=F)
#same as fit.13
fit.16=lme(var~Group*Died,
data=dat34,
random=~1|ID/Event_name)
anova(fit.16, type="marginal", adjustSigma=F)
#same as fit.13
#######################
More information about the R-sig-ecology
mailing list