Skip to content

What it means for rho to be 0 in lme() when using compound symmetry

4 messages · Simon Harmel, Phillip Alday

Original
Phillip Alday
# 
From
https://stat.ethz.ch/R-manual/R-devel/library/nlme/html/pdCompSymm.html :

"This function is a constructor for the pdCompSymm class, representing a
positive-definite matrix with compound symmetry structure (constant
diagonal and constant off-diagonal elements)."

Any multiple of the identity matrix is technically compound symmetric,
because all the off-diagonal elements are the same (0).

Phillip
On 28/11/20 2:30 am, Simon Harmel wrote:
Simon Harmel
# 
Thanks, Phillip. Given the estimated rho of 0 obtained from the default
correlation structure in lme(), can we say that for this dataset there is
no dependence left after fitting the 2-level model shown in my original
post?

In other words, once getting a rho of 0 from the default correlation
structure for this model, then one doesn't need to think of alternative
correlation structures, because even the default correlation structure has
shown that there is no dependence to model.

Is this a reasonable conclusion?
On Fri, Dec 4, 2020, 9:11 AM Phillip Alday <me at phillipalday.com> wrote:

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

  
  


  


      

    

    

      
      

        


  

        

      Phillip Alday
    

    

      
      #
    

  


  

      
The 'default' structure is generally the unstructured positive definite
matrix. (In lme4, this constraint is loosened to positive
semi-definite). Starting from compound symmetric is already placing a
constraint and so forcing all off-diagonal elements to zero may be the
best solution *under that constraint*. (And a multiple of the identity
matrix is a stricter constraint than a diagonal matrix, even though both
have all off diagonal elements set to zero.)

A more modern take would be using something like rePCA on the full
unstructured matrix (https://arxiv.org/abs/1506.04967 or
https://doi.org/10.33016/nextjournal.100002) to see what the effective
dimensionality is. Of course, you can fit a model with a diagonal RE
covariance matrix even if the data have some correlation, at the cost of
changing how shrinkage works
(https://doingbayesiandataanalysis.blogspot.com/2019/07/shrinkage-in-hierarchical-models-random.html).
That may be an acceptable (variance-bias) tradeoff -- less efficient
shrinkage but also less overparameterization.

All of these comments without looking at your data.

        
On 4/12/20 4:23 pm, Simon Harmel wrote:

            

      
      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Simon Harmel
    

    

      
      #
    

  


  

      
Phillip, I think I called corCompSymm() default in my specific case because
the rho was estimated to be "0" and of course the same variances are taken
to be on the diagonal. In other words, I thought such a result from
corCompSymm() acts as an equivalent to the default covariance structure of
random-effects (see AIC() comparison below).

The above reasoning also made me think that perhaps there is no correlation
left after fitting a default model and thus no need to think about what
correlation structure may best fit this particular dataset.

I also read Bates et al (2015). My takeaway was that we should use
lme4::rePCA to detect possible overfit in the random structure of the model
(often fitting random slopes when no/little such variation among the
individual slopes really exists in the data).

So, would you mind elaborating on the connection between using lme4::rePCA
and specifying the correlation structure in the model (e.g., using lme())?

For concreteness, here is the simple model I have in mind:

library(nlme)
data <- read.csv('https://raw.githubusercontent.com/hkil/m/master/R.csv')

m1 <- lme(Achieve ~ time, random = ~1|subid, data = data, correlation =
corCompSymm())
m2 <- lme(Achieve ~ time, random = ~1|subid, data = data, correlation = NULL
)

AIC(m1, m2)

   df      AIC
m1  5 2534.387  # Yes, rho is estimated But it is 0!
m2  4 2532.387

        
On Fri, Dec 4, 2020 at 9:32 AM Phillip Alday <me at phillipalday.com> wrote: