2
$\begingroup$

I am trying to decide what covariance structure top specify for a mixed model in R. When comparing two models with an unstructured covariance (corSymm) i can specify whether or not to allow coefficients to vary.

Here is my syntax

str(lCtr <- lmeControl(opt = "optim"))
m1 <- lme(y ~  Time, random = ~ 1 + Time | id, data=sdata, na.action = na.exclude, 
           correlation = corSymm(fixed = FALSE), control = lCtr)
m2 <-  lme(y ~ Time, random = ~ 1 + Time | id, data=sdata, na.action = na.exclude,
           correlation = corSymm(fixed = TRUE), control = lCtr)
-2*logLik(m1)
-2*logLik(m2)

and results

'log Lik.' 1336.108 (df=21)
'log Lik.' 1354.123 (df=6)

What confuses me is why i get more degrees of freedom when I specify that the coefficients are allowed to vary than the other way. What does it mean that the coefficients are allowed to vary?

Is there a better way to specify an unstructured covariance for the model?

$\endgroup$
1
  • $\begingroup$ The only way this would make sense to me is if df= says how many parameters the model had to estimate. I have 6 measurments so it should be 21, in a model with unstructured covariance residual matrix, according to my understanding. $\endgroup$
    – Wasser
    Commented Oct 6, 2016 at 14:13

1 Answer 1

1
$\begingroup$

The degrees of freedom refers to the number of parameters that are estimated to determine the likelihood. See df = p + length(coef(object[["modelStruct"]])) + as.integer(!fixSig)) in the last line in the code of the function below.

> getAnywhere(logLik.lme)
A single object matching ‘logLik.lme’ was found
It was found in the following places
  registered S3 method for logLik from namespace nlme
  namespace:nlme
with value

function (object, REML, ...) 
{
    fixSig <- attr(object[["modelStruct"]], "fixedSigma")
    fixSig <- !is.null(fixSig) && fixSig
    od <- object$dims
    p <- od$ncol[[od$Q + 1L]]
    N <- od$N
    estM <- object$method
    if (missing(REML)) 
        REML <- estM == "REML"
    val <- object[["logLik"]]
    if (REML && (estM == "ML")) {
        val <- val + (p * (log(2 * pi) + 1L) + (N - p) * log(1 - 
            p/N) + sum(log(abs(svd.d(object$varFix)))))/2
    }
    if (!REML && (estM == "REML")) {
        val <- val - (p * (log(2 * pi) + 1L) + N * log(1 - p/N) + 
            sum(log(abs(svd.d(object$varFix)))))/2
    }
    structure(val, class = "logLik", nall = N, nobs = N - REML * 
        p, df = p + length(coef(object[["modelStruct"]])) + as.integer(!fixSig))
}
<bytecode: 0x47df7dd8>
<environment: namespace:nlme>
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.