Skip to content
Prev 5330 / 5632 Next

[R-meta] Non-positive definite covariance matrix for rma.mv

Hi Gerta,

Actually, the var-cov matrix of all pairwise comparisons is not PD. Consider these data:

dat <- data.frame(
   study = c(1,1,1),
   m1i  = c(10,10,8),
   m2i  = c(8,7,7),
   sd1i = c(1.2,1.2,1.4),
   sd2i = c(1.4,1.3,1.3),
   n1i  = c(45,45,42),
   n2i  = c(42,40,40),
   grp1 = c("A", "A", "B"),
   grp2 = c("B", "C", "C"))
dat

We can compute all pairwise mean differences with:

dat <- escalc(measure="MD", m1i=m1i, sd1i=sd1i, n1i=n1i,
                            m2i=m2i, sd2i=sd2i, n2i=n2i, data=dat)
dat

The variance-covariance matrix of the mean differences can be constructed with:

V <- vcalc(vi, cluster=study, grp1=grp1, grp2=grp2, w1=n1i/sd1i^2, w2=n2i/sd2i^2, data=dat)
V

We can double-check that these values are correct:

- yi[1] and yi[2] share group A, so the covariance between them must be:

dat$sd1i[1]^2 / dat$n1i[1]

- yi[1] and yi[3] share group B (but in different positions), so the covariance is:

-dat$sd2i[1]^2 / dat$n2i[1]

- yi[2] and yi[3] share group C, so the covariance is:

dat$sd2i[2]^2 / dat$n2i[2]

So these are all correct. But V is not PD. It is semi-PD, with a third eigenvalue exactly equal to 0:

round(eigen(V)$values, 8)

Best,
Wolfgang