Constructing Covariance Matrices from Scales and Correlations

Why Build a Covariance Matrix from Scales and Correlations?

A random-intercept and random-slope model needs to describe not only how much each group-level effect varies, but also how those effects vary together. A covariance matrix does both. Yet specifying its entries directly can be awkward: variances and covariances have squared or product units, and an arbitrary collection of entries may fail to define a valid covariance matrix. An alternative is to specify each random-effect standard deviation and the correlations among effects, then assemble the covariance matrix from them.

This construction is useful in Bayesian GLMMs because it separates the marginal scale of each random effect from its dependence with the others. The introductory material on joint priors for random-effect scales and correlations motivates that separation; here the focus is the matrix algebra that turns those quantities into a valid covariance structure.

Definition: For $q$ random effects, let $\tau_k$ be the standard deviation of random effect $k$, and let $R$ be their $q\times q$ correlation matrix. Define $D=\operatorname{diag}(\tau_1,\ldots,\tau_q)$. The corresponding covariance matrix is $$ \Sigma = D R D,\qquad \Sigma_{k\ell}=\tau_k\tau_\ell R_{k\ell}. $$ Thus $\Sigma_{kk}=\tau_k^2$, while $\Sigma_{k\ell}$ is the covariance between effects $k$ and $\ell$.

A correlation matrix is symmetric, has ones on its diagonal, and has entries between $-1$ and $1$. Those entry-by-entry requirements alone are not enough: the matrix must also be positive semidefinite. For a nonsingular multivariate normal random-effects distribution, the usual requirement is positive definiteness. This means $v^\mathsf{T}Rv>0$ for every nonzero vector $v$. The same definition applies to $\Sigma$.

The scales $\tau_k$ must be expressed in the units of their respective random effects. A covariance has the product of the two effects’ units, whereas a correlation is dimensionless. Consequently, multiplying scales by correlations gives covariances with the right units automatically.

Construction and Positive Definiteness

Suppose the group-level vector $b_j$ has mean zero and covariance $\Sigma$. Under a multivariate normal random-effects assumption, one writes $b_j\sim N_q(0,\Sigma)$, with group effects typically assumed independent across groups conditional on the model parameters. The matrix construction itself is an algebraic identity; normality is an additional distributional assumption, not a requirement for calculating $\Sigma$.

Theorem: If every $\tau_k>0$ and $R$ is positive definite, then $\Sigma=DRD$ is positive definite. Conversely, if every $\tau_k>0$, then $\Sigma$ is positive definite if and only if $R$ is positive definite.

The result follows by evaluating the quadratic form. For any vector $x$, let $y=Dx$. Since every scale is positive, $D$ is invertible, so $x\ne 0$ implies $y\ne 0$. Therefore,

x^\mathsf{T}\Sigma x =x^\mathsf{T}DRDx =(Dx)^\mathsf{T}R(Dx) =y^\mathsf{T}Ry.

If $R$ is positive definite, the last expression is positive for every nonzero $x$. Conversely, if $\Sigma$ is positive definite, take any nonzero $y$ and set $x=D^{-1}y$; then $y^\mathsf{T}Ry=x^\mathsf{T}\Sigma x>0$, so $R$ is positive definite too. This argument also explains why positive scales matter: if a scale is zero, $D$ is singular and $\Sigma$ cannot be positive definite, even if $R$ is.

Two useful consequences, when $R$ and $D$ are invertible, are

\det(\Sigma)=\left(\prod_{k=1}^q\tau_k^2\right)\det(R), \qquad \Sigma^{-1}=D^{-1}R^{-1}D^{-1}.

They provide checks on computations and show how the scales affect the covariance matrix. A positive determinant alone does not establish positive definiteness for a general matrix; all required eigenvalues or principal-minor conditions must be checked.

Worked Example: A Random Intercept and Slope

Suppose a hypothetical logistic GLMM has a group-specific random intercept with standard deviation $1.2$ logit units and a random slope with standard deviation $0.8$ logit units per predictor unit. Suppose their correlation is $-0.35$. Construct the covariance matrix and verify that it is positive definite.

The scales are positive, and the proposed two-effect correlation matrix is

D= \begin{pmatrix}1.2&0\\0&0.8\end{pmatrix}, \qquad R= \begin{pmatrix}1&-0.35\\-0.35&1\end{pmatrix}.

Multiplying each correlation by the corresponding pair of scales gives

\Sigma=DRD = \begin{pmatrix} 1.2^2 & (1.2)(0.8)(-0.35)\\ (1.2)(0.8)(-0.35) & 0.8^2 \end{pmatrix} = \begin{pmatrix} 1.44&-0.336\\ -0.336&0.64 \end{pmatrix}.

The covariance is $-0.336$ logit-squared per predictor unit. For a symmetric $2\times2$ matrix, positive definiteness requires a positive first diagonal entry and a positive determinant. Here,

1.44>0,\qquad \det(\Sigma)=(1.44)(0.64)-(-0.336)^2 =0.9216-0.112896 =0.808704>0.

Thus $\Sigma$ is positive definite. Equivalently, the two-effect correlation matrix is positive definite because its determinant is $1-(-0.35)^2=0.8775>0$. The negative off-diagonal entry indicates that groups with larger intercept deviations tend, under this model, to have smaller slope deviations; it does not by itself imply a causal relationship.

Checking Correlation Matrices Beyond Two Effects

With three or more random effects, checking that each proposed correlation lies in $[-1,1]$ is insufficient. The entries must be jointly compatible with one positive-semidefinite matrix. For a three-effect correlation matrix with off-diagonal entries $a$, $b$, and $c$, positive definiteness can be checked by requiring positive leading principal minors. In this ordering, the conditions are $|a|<1$, $1-a^2>0$, and a positive determinant:

R= \begin{pmatrix} 1&a&b\\ a&1&c\\ b&c&1 \end{pmatrix}, \qquad \det(R)=1+2abc-a^2-b^2-c^2>0.

For a larger matrix, a Cholesky factorization is a practical numerical check: a symmetric matrix is positive definite precisely when it has a Cholesky factorization with positive diagonal entries. Near-zero pivots or eigenvalues indicate that the matrix is close to singular, so rounding and numerical precision deserve attention.

Worked Example: Valid Three-Effect Correlations

For three hypothetical random effects, suppose the correlations are $a=0.5$, $b=0.2$, and $c=-0.1$. The proposed correlation matrix is symmetric, has ones on the diagonal, and all three off-diagonal entries lie between $-1$ and $1$. Check the stronger positive-definiteness condition, then construct the covariance matrix for scales $2$, $0.5$, and $3$.

The leading principal minors are $1$, $1-a^2=1-0.25=0.75$, and the full determinant:

\det(R)=1+2(0.5)(0.2)(-0.1)-(0.5)^2-(0.2)^2-(-0.1)^2 =1-0.02-0.25-0.04-0.01 =0.68.

All three are positive, so $R$ is positive definite. With $D=\operatorname{diag}(2,0.5,3)$, the diagonal entries of $\Sigma$ are $4$, $0.25$, and $9$. The off-diagonal entries are the scale products times the corresponding correlations:

\Sigma= \begin{pmatrix} 4&(2)(0.5)(0.5)&(2)(3)(0.2)\\ (2)(0.5)(0.5)&0.25&(0.5)(3)(-0.1)\\ (2)(3)(0.2)&(0.5)(3)(-0.1)&9 \end{pmatrix} = \begin{pmatrix} 4&0.5&1.2\\ 0.5&0.25&-0.15\\ 1.2&-0.15&9 \end{pmatrix}.

The scales are all positive and $R$ is positive definite, so the theorem guarantees that $\Sigma$ is positive definite. As a numerical cross-check, $\det(\Sigma)=(2^2)(0.5^2)(3^2)(0.68)=9(0.68)=6.12>0$. The determinant check is consistent with the conclusion, though the quadratic-form result supplies the general justification.

Worked Example: Pairwise-Valid Correlations That Fail Together

Consider a proposed three-effect correlation matrix with $a=0.9$, $b=0.9$, and $c=-0.9$. Each pairwise value is within the permitted interval, and the matrix is symmetric with ones on the diagonal. Determine whether it is a valid positive-definite correlation matrix.

Its determinant is

\det(R)=1+2(0.9)(0.9)(-0.9)-(0.9)^2-(0.9)^2-(-0.9)^2 =1-1.458-2.43=-2.888.

The determinant is negative, so this matrix is not positive definite; it cannot be used as the correlation matrix of a nonsingular multivariate normal random-effects distribution. The failure is joint: the first two effects are each assigned a strong positive correlation with the third, while their mutual correlation is strongly negative. Pairwise bounds do not ensure those three claims can coexist as a covariance structure. Do not repair the problem by merely clipping an entry without reconsidering the intended dependence; choose a jointly valid matrix instead.

Predictor Centering Changes the Covariance Entries

Covariance matrices depend on how random effects are defined. In a random-intercept and random-slope model, the intercept is the deviation at predictor value zero. Changing the origin of the predictor therefore changes the random intercept and its covariance with the slope, even though it does not change the underlying group-specific linear predictor.

Let $x^*=x-c$, so $x=x^*+c$. If a group’s random contribution is $b_0+b_1x$, then it can be rewritten as $b_0^*+b_1x^*$, where $b_0^*=b_0+cb_1$. In matrix form, this is a linear transformation of the random-effects vector:

\begin{pmatrix}b_0^*\\b_1\end{pmatrix} = A\begin{pmatrix}b_0\\b_1\end{pmatrix}, \qquad A= \begin{pmatrix}1&c\\0&1\end{pmatrix}, \qquad \Sigma^*=A\Sigma A^\mathsf{T}.

The transformed covariance matrix must still be positive definite if the original one is, since $A$ is invertible. But its scales and correlation generally differ. A correlation describes dependence between the effects as defined at a particular predictor origin; it is not invariant to every change of that origin.

Worked Example: Re-centering a Random Slope

Suppose a hypothetical model has random-intercept and random-slope standard deviations $1.5$ logit units and $0.4$ logit units per predictor unit, with correlation $0.5$. The covariance is $(1.5)(0.4)(0.5)=0.3$, so the original covariance matrix is

\Sigma= \begin{pmatrix} 2.25&0.3\\ 0.3&0.16 \end{pmatrix}, \qquad \det(\Sigma)=(2.25)(0.16)-(0.3)^2=0.27>0.

Now center the predictor two units above its old origin, so $c=2$. The new random intercept is $b_0^*=b_0+2b_1$. Its variance and covariance with the unchanged slope are

\operatorname{Var}(b_0^*)=2.25+2(2)(0.3)+2^2(0.16)=4.09, \qquad \operatorname{Cov}(b_0^*,b_1)=0.3+2(0.16)=0.62.

Thus

\Sigma^*= \begin{pmatrix}4.09&0.62\\0.62&0.16\end{pmatrix}, \qquad \det(\Sigma^*)=(4.09)(0.16)-(0.62)^2 =0.6544-0.3844 =0.27.

The determinant is unchanged here because $\det(A)=1$, and positive definiteness is preserved. The new intercept standard deviation is $\sqrt{4.09}\approx2.022$; the new correlation is $0.62/[(2.022)(0.4)]\approx0.766$. The covariance structure changed because the intercept now represents heterogeneity at a different predictor value, not because the underlying fitted linear predictor changed.

Diagnostics, Alternatives, and Common Mistakes

A reliable construction check has three parts: confirm the scales are nonnegative and in the intended units; confirm the correlation matrix is symmetric with a unit diagonal; and check positive definiteness (or positive semidefiniteness if singular covariance structures are intentionally allowed). For a matrix assembled from strictly positive scales and a valid positive-definite $R$, the theorem verifies the final matrix without separately testing every entry of $\Sigma$. In software, also inspect numerical warnings and near-zero eigenvalues or Cholesky pivots, especially when posterior draws approach a singular boundary.

There are several ways to represent the same covariance structure. Directly assigning variances and covariances is transparent for a small, fixed matrix, but validity constraints can be difficult to maintain. Specifying standard deviations and correlations gives interpretable marginal scales and dimensionless dependence parameters, but the correlations must jointly form a valid matrix. A Cholesky-factor representation is often convenient for computation because positive diagonal factor entries can enforce positive definiteness; its factor entries may be less directly interpretable as marginal scales and correlations. These are alternative parameterizations of covariance structure, not different covariance definitions.

Mistake: checking only that every correlation is between $-1$ and $1$. The invalid three-effect example shows why this is insufficient. Check joint positive definiteness.
Mistake: putting a standard deviation on the covariance diagonal. The diagonal is a variance, $\tau_k^2$; a scale appears once in each factor of an off-diagonal entry.
Mistake: assuming correlation is unchanged by predictor centering. The random-intercept definition changes with the predictor origin. Transform the full covariance matrix before interpreting the new correlation.
Mistake: treating a tiny positive determinant as proof of a stable matrix. A matrix may be nearly singular. Inspect eigenvalues or Cholesky diagnostics as well, and distinguish numerical roundoff from a genuinely invalid specification.

The construction does not determine whether chosen scales or correlations are substantively plausible; that requires model-specific prior reasoning and prior predictive assessment. Nor does positive definiteness guarantee that the data contain enough information to estimate every covariance parameter precisely. It guarantees mathematical validity of the covariance structure, not adequacy of the full model.

Check Your Understanding

Answer each question, then compare with the solutions below.

For scales $\tau_1=2$ and $\tau_2=3$ with correlation $0.4$, what is the covariance matrix, and is it positive definite?
A $3\times3$ proposed correlation matrix has off-diagonal values $0.8$, $0.8$, and $-0.8$. Why are the individual correlation bounds not enough to validate it?
If one random-effect scale is zero, can $DRD$ be positive definite? Explain.
Under $x^*=x-c$, which transformed random effect is $b_0+cb_1$, and what changes about its interpretation?

Solutions. (1) The matrix is $\begin{pmatrix}4&2.4\\2.4&9\end{pmatrix}$, since the covariance is $(2)(3)(0.4)=2.4$. Its determinant is $36-5.76=30.24>0$, and its first diagonal entry is positive, so it is positive definite. (2) The entries must satisfy a joint determinant condition, not just separate bounds. Here $\det(R)=1+2(0.8)(0.8)(-0.8)-3(0.8^2)=1-1.024-1.92=-1.944<0$, so it is not positive definite. (3) No. A zero scale makes $D$ singular and gives $\Sigma$ a zero row and column for that effect, so it cannot be positive definite, although it can be positive semidefinite. (4) It is the transformed random intercept $b_0^*$. It represents the group deviation at the new predictor origin; its variance and correlation with the slope may differ from those at the original origin.