-----Oorspronkelijk bericht-----
Van: r-sig-mixed-models-bounces at r-project.org
[mailto:r-sig-mixed-models-bounces at r-project.org] Namens
Sebasti?n Daza
Verzonden: dinsdag 12 april 2011 22:47
Aan: R-SIG-Mixed-Models at r-project.org
Onderwerp: Re: [R-sig-ME] level 1 variance-covariance structure
Thierry,
I can run this model... but what does it mean?
The correlation structure that I get is:
Correlation Structure: ARMA(1,0)
Formula: ~age13 | id
Parameter estimate(s):
Phi1
0
What does zero mean? I would expect get some positive number there...
m3a <- lme(attit ~ 1 + age13 , data, random= ~ 0 +
factor(age13)| id, correlation = corAR1(form = ~ age13 | id))
Linear mixed-effects model fit by REML
Data: data
AIC BIC logLik
-324.2096 -229.5528 181.1048
Random effects:
Formula: ~0 + factor(age13) | id
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
factor(age13)-2 0.17219431 f(13)-2 f(13)-1 f(13)0 f(13)1
factor(age13)-1 0.19789254 0.493
factor(age13)0 0.25942941 0.425 0.544
factor(age13)1 0.28354459 0.442 0.442 0.723
factor(age13)2 0.29097081 0.498 0.474 0.639 0.808
Residual 0.07457025
Correlation Structure: ARMA(1,0)
Formula: ~age13 | id
Parameter estimate(s):
Phi1
0
Fixed effects: attit ~ 1 + age13
Value Std.Error DF t-value p-value
(Intercept) 0.3210558 0.012832840 839 25.01829 0
age13 0.0593529 0.004716984 839 12.58282 0
Correlation:
(Intr)
age13 0.504
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-1.46371874 -0.27170442 -0.04080688 0.26239551 1.69883907
Number of Observations: 1079
Number of Groups: 239
On 4/12/2011 10:21 AM, ONKELINX, Thierry wrote:
Dear Sebastian,
The model below works fine on my computer.
m3a<- lme(attit ~ 1 + age13 , data=dataset, random= ~
0+factor(age13)| id, correlation = corAR1(form = ~ age13 | id))
Best regards,
Thierry
----------------------------------------------------------------------
------
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek team Biometrie&
Gaverstraat 4 9500 Geraardsbergen Belgium
Research Institute for Nature and Forest team Biometrics& Quality
Assurance Gaverstraat 4 9500 Geraardsbergen Belgium
tel. + 32 54/436 185
Thierry.Onkelinx at inbo.be
www.inbo.be
To call in the statistician after the experiment is done
may be no more than asking him to perform a post-mortem
examination: he may be able to say what the experiment died of.
~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data.
~ Roger Brinner
The combination of some data and an aching desire for an
answer does not ensure that a reasonable answer can be
extracted from a given body of data.
-----Oorspronkelijk bericht-----
Van: Sebasti?n Daza [mailto:sebastian.daza at gmail.com]
Verzonden: dinsdag 12 april 2011 15:43
Aan: ONKELINX, Thierry
CC: R-SIG-Mixed-Models at r-project.org
Onderwerp: Re: [R-sig-ME] level 1 variance-covariance structure
Thank you for your reply Thierry...
Increasing the number of iterations doesn't work:
m3a<- lme(attit ~ 1 + age13 , data=data, random= ~ age13 | id,
correlation = corAR1(, form = ~ ind | id),
control=list(maxIter=1000, msMaxIter=1000, niterEM=1000))
Error in lme.formula(attit ~ 1 + age13, data = data, random =
~age13 | :
nlminb problem, convergence error code = 1
message = function evaluation limit reached without
(9)
I have attached my database. I don't know if it is a problem of my
model or a limitation of lme function.
The best!
Sebastian.
On 4/12/2011 6:25 AM, ONKELINX, Thierry wrote:
Dear Sebastian,
You don't need to create dummy variables your selve.
You can write m2a<- lme(attit ~ 1 + age13 , data=data,
random= ~ 0 + ind1+ ind2+ ind3+ ind4+ ind5 | id, method="REML") as
m2a<- lme(attit ~ 1 + age13 , data=data, random= ~ 0 +
factor(ind) | id, method="REML")
Or if ind is an indicator for age13:
m2a<- lme(attit ~ 1 + age13 , data=data, random= ~ 0 +
factor(age13) | id, method="REML")
Have a look at lmeControl() to increase the number of iterations.
Best regards,
Thierry
--------------------------------------------------------------
--------------
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek
team Biometrie& Kwaliteitszorg
Gaverstraat 4
9500 Geraardsbergen
Belgium
Research Institute for Nature and Forest
team Biometrics& Quality Assurance
Gaverstraat 4
9500 Geraardsbergen
Belgium
tel. + 32 54/436 185
Thierry.Onkelinx at inbo.be
www.inbo.be
To call in the statistician after the experiment is done
may be no more than asking him to perform a post-mortem
examination: he may be able to say what the experiment died of.
~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data.
~ Roger Brinner
The combination of some data and an aching desire for an
answer does not ensure that a reasonable answer can be
-----Oorspronkelijk bericht-----
Van: r-sig-mixed-models-bounces at r-project.org
[mailto:r-sig-mixed-models-bounces at r-project.org] Namens
Daza
Verzonden: maandag 11 april 2011 18:44
Aan: R-SIG-Mixed-Models at r-project.org
Onderwerp: [R-sig-ME] level 1 variance-covariance structure
Hi everyone,
I am trying to reproduce some results models from HLM (HMLM) to
contrast different specifications of level 1
but I get convergence errors. I would like to know if
# database structure
head(data[,c(1,2,6, 9:13,17)])
id attit age13 ind1 ind2 ind3 ind4 ind5 ind
1 3 0.11 -2 1 0 0 0 0 1
2 3 0.20 -1 0 1 0 0 0 2
3 3 0.00 0 0 0 1 0 0 3
4 3 0.00 1 0 0 0 1 0 4
5 3 0.11 2 0 0 0 0 1 5
6 8 0.29 -2 1 0 0 0 0 1
# attit is a deviant measure and ind variables indicate
waves # following some examples of snijders and bosker's book, I
get the unrestricted model:
> m2a<- lme(attit ~ 1 + age13 , data=data, random= ~ 0 +
ind1+ ind2+
ind3+ ind4+ ind5 | id, method="REML")
Linear mixed-effects model fit by REML
Data: data
AIC BIC logLik
-326.2096 -236.5348 181.1048
Random effects:
Formula: ~0 + ind1 + ind2 + ind3 + ind4 + ind5 | id
Structure: General positive-definite, Log-Cholesky
StdDev Corr
ind1 0.17219431 ind1 ind2 ind3 ind4
ind2 0.19789253 0.493
ind3 0.25942942 0.425 0.544
ind4 0.28354459 0.442 0.442 0.723
ind5 0.29097082 0.498 0.474 0.639 0.808
Residual 0.07457025
Fixed effects: attit ~ 1 + age13
Value Std.Error DF t-value p-value
(Intercept) 0.3210558 0.012832840 839 25.01829 0
age13 0.0593529 0.004716984 839 12.58282 0
Correlation:
(Intr)
age13 0.504
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-1.46371871 -0.27170442 -0.04080686 0.26239553 1.69883910
Number of Observations: 1079
Number of Groups: 239
# variance-covariance matrix
> extract.lme.cov2(m2a,data)$V[[6]]
25 26 27 28 29
25 0.03521160 0.01681647 0.01899029 0.02159300 0.02494013
26 0.01681647 0.04472218 0.02793174 0.02481343 0.02727012
27 0.01899029 0.02793174 0.07286434 0.05318967 0.04823107
28 0.02159300 0.02481343 0.05318967 0.08595826 0.06667139
29 0.02494013 0.02727012 0.04823107 0.06667139 0.09022474
# I get the same results than unrestricted model in HLM
# When I try to get the same unrestricted model using "corStruc"
commands in lme, I get a convergence problem. Am I
> m2b<- lme(attit ~ 1 + age13 , data=data, random= ~ age13
| id, correlation = corSymm(, form = ~ ind | id)) Error in
lme.formula(attit ~ 1 + age13, data = data, random = ~age13 | :
nlminb problem, convergence error code = 1
message = iteration limit reached without convergence (9)
# When I try to get an autoregressive model, I get again a
convergence problem.
> m3a<- lme(attit ~ 1 + age13 , data=data, random= ~ age13
| id, correlation = corAR1(, form = ~ ind | id)) Error in
lme.formula(attit ~ 1 + age13, data = data, random = ~age13 | :
nlminb problem, convergence error code = 1
message = iteration limit reached without convergence (9)
Does anyone know how I can solve this?
Thank you in advance.
--
Sebasti?n Daza
sebastian.daza at gmail.com