[R] question on multilevel modeling

Christine A. Calmes ccalmes at buffalo.edu
Tue Nov 7 15:55:37 CET 2006


Thank you for your quick response.  Basically, I am trying to test a
psychological model in which a person's level of passive, negative
communication with others in their environment (CORUMTO) on one week
predicts their level of depression on the following week (BDIAFTER)
controlling for their level of depression on the same week that that
they co-ruminate (BDI).  The WEEK variable accounts for the fact that
there are 7 weeks worth of data in the model and that scores on adjacent
weeks are more closer related than scores on weeks that are further
apart.
Here is the output for the following models--thanks again for your help!
 
#random slope on corumination but fixed intercept

rslope<-lme(BDIAFTER~BDI+WEEK+CORUMTO,
random=~CORUMTO-1|DYADID/PARTICIP, data=new) 
summary(rslope)

Linear mixed-effects model fit by REML
 Data: new 
       AIC      BIC    logLik
  5340.796 5375.016 -2663.398

Random effects:
 Formula: ~CORUMTO - 1 | DYADID
             CORUMTO
StdDev: 0.0001230839

 Formula: ~CORUMTO - 1 | PARTICIP %in% DYADID
             CORUMTO Residual
StdDev: 2.160106e-07   3.5853

Fixed effects: BDIAFTER ~ BDI + WEEK + CORUMTO 
                 Value Std.Error  DF  t-value p-value
(Intercept)  1.0536425 0.4055838 808  2.59784  0.0096
BDI          0.7343608 0.0219012 808 33.53056  0.0000
WEEK         0.0260517 0.0678471 808  0.38398  0.7011
CORUMTO     -0.0093802 0.0065327 808 -1.43589  0.1514
 Correlation: 
        (Intr) BDI    WEEK  
BDI     -0.259              
WEEK    -0.689  0.119       
CORUMTO -0.703 -0.062  0.115

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max 
-5.3790151 -0.3774653 -0.1883441  0.3183254  6.6093553 

Number of Observations: 985
Number of Groups: 
              DYADID PARTICIP %in% DYADID 
                  87                  174

#random intercept & fixed slope

rint<-lme(BDIAFTER~BDI+WEEK+CORUMTO, random=~1|DYADID/PARTICIP,
data=new) 
summary(rint)

Linear mixed-effects model fit by REML
 Data: new 
       AIC      BIC    logLik
  5332.942 5367.162 -2659.471

Random effects:
 Formula: ~1 | DYADID
        (Intercept)
StdDev:    1.435454

 Formula: ~1 | PARTICIP %in% DYADID
        (Intercept) Residual
StdDev:    2.607136 3.042346

Fixed effects: BDIAFTER ~ BDI + WEEK + CORUMTO 
                 Value Std.Error  DF   t-value p-value
(Intercept)  3.1159915 0.5394972 808  5.775733  0.0000
BDI          0.2937769 0.0308994 808  9.507528  0.0000
WEEK        -0.1406907 0.0597295 808 -2.355464  0.0187
CORUMTO      0.0000037 0.0087027 808  0.000427  0.9997
 Correlation: 
        (Intr) BDI    WEEK  
BDI     -0.269              
WEEK    -0.562  0.193       
CORUMTO -0.713 -0.073  0.192

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max 
-5.6443060 -0.3649701 -0.1069800  0.2750957  5.5662019 

Number of Observations: 985
Number of Groups: 
              DYADID PARTICIP %in% DYADID 
                  87                  174 

#random intercept and random slope 

rslint<-lme(BDIAFTER~BDI+WEEK+CORUMTO, random=~CORUMTO|DYADID/PARTICIP,
data=new) 
summary(rslint)


Linear mixed-effects model fit by REML
 Data: new 
       AIC      BIC    logLik
  5340.949 5394.724 -2659.475

Random effects:
 Formula: ~CORUMTO + 1 | DYADID
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev      Corr  
(Intercept) 1.438559367 (Intr)
CORUMTO     0.002629201 -0.068

 Formula: ~CORUMTO + 1 | PARTICIP %in% DYADID
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev      Corr  
(Intercept) 2.606533039 (Intr)
CORUMTO     0.001218270 -0.014
Residual    3.042301836       

Fixed effects: BDIAFTER ~ BDI + WEEK + CORUMTO 
                 Value Std.Error  DF   t-value p-value
(Intercept)  3.1154427 0.5396189 808  5.773413  0.0000
BDI          0.2938990 0.0308987 808  9.511686  0.0000
WEEK        -0.1405892 0.0597340 808 -2.353590  0.0188
CORUMTO     -0.0000028 0.0087093 808 -0.000324  0.9997
 Correlation: 
        (Intr) BDI    WEEK  
BDI     -0.269              
WEEK    -0.561  0.193       
CORUMTO -0.713 -0.073  0.192

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max 
-5.6453974 -0.3656707 -0.1070803  0.2746329  5.5665960 

Number of Observations: 985
Number of Groups: 
              DYADID PARTICIP %in% DYADID 
                  87                  174



Christine A. Calmes, MA
Dept of Psychology
University at Buffalo: The State University of New York
Park Hall 216
Buffalo, NY 14260
ccalmes at buffalo.edu
(716) 645-3650 x578


-----Original Message-----
From: Chuck Cleland [mailto:ccleland at optonline.net] 
Sent: Tuesday, November 07, 2006 5:00 AM
To: Christoph Buser
Cc: Christine A. Calmes; R-help at stat.math.ethz.ch
Subject: Re: [R] question on multilevel modeling

Christoph Buser wrote:
> Dear Christine
> 
> I think the problem in your second model is that you are
> including "CORUMTO" both as a fixed effect and as a random
> effect. 

  That by itself should not be a problem.  Here is an example in which
age appears in both the fixed and random part of a model:

> fm1 <- lme(distance ~ age, random = ~ age | Subject, data = Orthodont)

> fm1
Linear mixed-effects model fit by REML
  Data: Orthodont
  Log-restricted-likelihood: -221.3183
  Fixed: distance ~ age
(Intercept)         age
 16.7611111   0.6601852

Random effects:
 Formula: ~age | Subject
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev    Corr
(Intercept) 2.3270340 (Intr)
age         0.2264278 -0.609
Residual    1.3100397

Number of Observations: 108
Number of Groups: 27

  I think it's actually quite common to have the same variable appear in
the fixed and random parts of a model.
  Christine, to better understand what's happening it would be useful to
know a bit more about what BDIAFTER, BDI, WEEK, and CORUMTO are and see
the model summaries.  The fact that the error message depends on the
specific variables in the model may mean that one of the variables is a
linear combination of other variables.

> Regards,
> 
> Christoph 
> 
> --------------------------------------------------------------
> 
> Credit and Surety PML study: visit our web page www.cs-pml.org
> 
> --------------------------------------------------------------
> Christoph Buser <buser at stat.math.ethz.ch>
> Seminar fuer Statistik, LEO C13
> ETH Zurich	8092 Zurich	 SWITZERLAND
> phone: x-41-44-632-4673		fax: 632-1228
> http://stat.ethz.ch/~buser/
> --------------------------------------------------------------
> 
> 
> Christine A. Calmes writes:
>  > Hi,
>  >  
>  > I am trying to run a multilevel model with time nested in people
and
>  > people nested in dyads (3 levels of nesting) by initially running a
>  > series of models to test whether the slope/intercept should be
fixed or
>  > random.  The problem that I am experiencing appears to arise
between the
>  > random intercept, fixed slope equation AND.
>  >  
>  > (syntax:
>  > rint<-lme(BDIAFTER~BDI+WEEK+CORUMTO, random=~1|DYADID/PARTICIP,
>  > data=new) 
>  > summary(rint))
>  >  
>  > the random slope, random intercept model 
>  >  
>  > (syntax: 
>  > rslint<-lme(BDIAFTER~BDI+WEEK+CORUMTO,
random=~CORUMTO|DYADID/PARTICIP,
>  > data=new) 
>  > summary(rslint))
>  >  
>  > at which point I obtain the exact same results for each model
suggesting
>  > that one of the model is not properly specifying the slope or
intercept.
>  >  
>  > Or, I receive the following error message when I try to run the
random
>  > slope/random intercept model.
>  >  
>  > Error in solve.default(pdMatrix(a, fact = TRUE)) : 
>  >         system is computationally singular: reciprocal condition
number
>  > = 6.77073e-017
>  >  
>  > (whether I receive an error message or the same results depends on
the
>  > specific variables in the model).
>  >  
>  > It has been suggested that I may need to change the default
starting
>  > values in the model because I may be approaching a boundary-is this
a
>  > plausible explanation for my difficulties?  If so, how do I do this
in R
>  > and can you refer me to a source that might highlight what would be
>  > reasonable starting values?
>  > If this does not seem like the problem, any idea what the problem
may be
>  > and how I might fix it?
>  >  
>  > Thank you so much for your assistance,
>  > Christine Calmes
>  >  
>  >  
>  > Christine A. Calmes, MA
>  > Dept of Psychology
>  > University at Buffalo: The State University of New York
>  > Park Hall 216
>  > Buffalo, NY 14260
>  > ccalmes at buffalo.edu
>  > (716) 645-3650 x578
>  >  
>  > 
>  > 	[[alternative HTML version deleted]]
>  > 
>  > ______________________________________________
>  > R-help at stat.math.ethz.ch mailing list
>  > https://stat.ethz.ch/mailman/listinfo/r-help
>  > PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
>  > and provide commented, minimal, self-contained, reproducible code.
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 

-- 
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894



More information about the R-help mailing list