Temporal autoregression in glmmTMB: what do alternative syntax forms mean?

Question

The argument ar1() in glmmTMB accepts two different forms of syntax (that I know of, there might be others):

1. ar1(time + 0 | group)
2. ar1(time - 1 | group)

Using one or the other produces the same outcome as far as I can tell, so why are different equivalent forms allowed and what do - 1 and + 0 stand for? With regards to - 1, Ben Bolker writes:

If we use ar1(tt|f), with glmmTMB we get a warning message (“AR1 not meaningful with intercept”). This is important; it made me aware of a similar mistake I was making previously with my lmer hack below. Since lme4 uses unstructured (i.e. general positive-definite) variance-covariance matrices by default, it normally doesn’t matter how you parameterize the contrasts for a categorical variable – the model fit/predictions are invariant to linear transformations. This is no longer true when we use structured variance-covariance matrices (!), so we need (tt-1|f) rather than (tt|f) …

This doesn't help me understand, probably because of my lack of understanding of the difference between structured and unstructured variance-covariance matrices. Here they explain:

in an unstructured covariance matrix there are no constraints. Each variance and each covariance is estimated uniquely from the data. As you can imagine, this results in the best possible model fit, because each variance and covariance values is very close to what the data reflect.

But that comes at a cost that may not be worth the improved fit. Estimating each one of those values can use up many degrees of freedom.

I guess the follow-up question is "in which way are structured variance-covariance matrices structured?" Besides, I thought that the whole point of assuming an order-1 autoregressive structure in the data was to assume that residuals were autocorrelated in a fixed way (reference). Doesn't assuming an order-1 autoregressive structure in the data lead to a structured variance-covariance matrix in any case?

Answer 1

This question is still a bit of a mixed bag (even though it has been refined from an earlier Stack Overflow question). However:

• the equivalence of -1/+0 is pure "semantic sugar". Long ago (probably in S, before R was born) the convention was established that intercepts were implicitly included in model formulae, so that e.g. y ~ x was exactly equivalent to y ~ 1 + x (i.e. $y= \beta_0 + \beta_1 x$). Since one sometimes wants to fit a model without an intercept, the -1 syntax was added to allow users to suppress the intercept (i.e. y ~ -1 + x corresponds to $y=\beta_1 x$). I believe that +0 was added later as an alternative, possibly more aesthetically pleasing, way to specify a no-intercept model. Because these definitions are built in to base-R's model.matrix() function, they are used in a wide variety of R packages. (Today I learned that using -0 in a formula adds an intercept to the model, equivalent to +1 !)

• the other point is more subtle. In constructing the model for a random effect, R first expands the effect (i.e., the f in a term like ar1(f|g)) into a model matrix. If f is a factor, the setup of that model matrix depends on whether an intercept is included in the model or not (and on the contrasts associated with the factor). Suppose we have a factor f <- factor(1:3). If the formula is ~f (and the default treatment contrasts are used), the corresponding model matrix is

int f2  f3
1   0   0
1   1   0
1   0   1

in other words, it is composed of an intercept column (all ones) and $n-1$ indicator variables for each level of the factor other than the first. In contrast if we use ~f-1 or ~f+0 the model matrix is

f1  f2  f3
1   0   0
0   1   0
0   0   1

in other words, we have an indicator variable for every level (and no intercept).

These two parameterizations are equivalent up to a linear transformation. That means that if we were fitting a regular linear model, or a random effects term with a general/unstructured covariance matrix, then we would get the same overall model fit (e.g. the same likelihood, the same model predictions, etc.) even though we would get different sets of estimated parameters (or a different covariance matrix if we were fitting a random effect term); we could transform one vector of parameter estimate, or one covariance matrix, into the other by an appropriate set of linear transformations.

However, because the AR1 covariance matrix represents only a subset of all allowable (positive-semidefinite) covariance matrices, it now matters how we parameterize the underlying model matrix.