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

Thu Jul 4 15:46:39 CEST 2024

Hi Wolfgang and James,

Sorry for not responding earlier. You are right, the variance-covariance matrix of a consistent network has rank n – 1 where n is the number of nodes (treatments) in the network. Thus it is not positive definite, but always non-negative definite (at least in theory, up to numerical issues).

# netmeta code for Wolfgang’s example:

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
p3 <- pairwise(treat = list(grp1, grp2), mean = list(m1i, m2i),
               sd = list(sd1i, sd2i), n = list(n1i, n2i),
               studlab = study, data = dat)
net3 <- netmeta(p3)
(d <- choose(net3$n,2))
net3$n - 1
qr(net3$Cov.common)$rank
round(eigen(net3$Cov.common)$values, 8)

Best,
Gerta

UNIVERSITÄTSKLINIKUM FREIBURG
Institute for Medical Biometry and Statistics

Dr. Gerta Rücker
Guest Scientist

Stefan-Meier-Straße 26 · 79104 Freiburg
gerta.ruecker using uniklinik-freiburg.de

https://www.uniklinik-freiburg.de/imbi-en/employees.html?imbiuser=ruecker

Von: R-sig-meta-analysis <r-sig-meta-analysis-bounces using r-project.org> Im Auftrag von James Pustejovsky via R-sig-meta-analysis
Gesendet: Montag, 1. Juli 2024 15:52
An: Andreas Voldstad <andreas.voldstad using kellogg.ox.ac.uk>
Cc: James Pustejovsky <jepusto using gmail.com>; R Special Interest Group for Meta-Analysis <r-sig-meta-analysis using r-project.org>
Betreff: Re: [R-meta] Non-positive definite covariance matrix for rma.mv

I am not sure. Setting nearpd = TRUE transforms the vcov matrices in a way that is hard to interpret in substantive terms, and it’s hard to predict the implication of this for the meta-analytic results. 

As Gerta noted, network meta-analysis methods can deal with multi-arm trials in a principled way. It does so by making use of specialized contrast coding schemes, but the input effect sizes are still all calculated relative to a single control condition per study (e.g., all effects calculated relative to the passive control).

James

On Jun 30, 2024, at 5:22 PM, Andreas Voldstad <mailto:andreas.voldstad using kellogg.ox.ac.uk> wrote:

Hi James, thank you for responding. Indeed, I am fitting a model across all the contrasts, and metafor did yell at me for the NPD V matrix. Is this a case where it might be appropriate to use the nearpd = T argument in vcalc?

Best wishes,
Andreas

Sent from https://aka.ms/o0ukef
________________________________________
From: James Pustejovsky <mailto:jepusto using gmail.com>
Sent: Sunday, June 30, 2024 5:13:32 PM
To: R Special Interest Group for Meta-Analysis <mailto:r-sig-meta-analysis using r-project.org>
Cc: Andreas Voldstad <mailto:andreas.voldstad using kellogg.ox.ac.uk>
Subject: Re: [R-meta] Non-positive definite covariance matrix for rma.mv 

Hi Andreas,

I think the NPD issue arises from having treatment contrasts that are perfectly collinear. Consider that if you know the contrasts Tx1 vs Passive, Tx2 vs Passive, and Tx1 vs Active, then you can infer the Tx2 vs Active contrast: 
(Tx2 vs Passive) - (Tx1 vs Passive) + (Tx1 vs Active) = (Tx2 - Tx1)  + (Tx1 vs Active) = (Tx2 vs Active)

As far as I can see, the only way to avoid this is to calculate all effect size contrasts relative to a common comparison condition. 

Depending on the meta-analytic model in which you use this covariance matrix, the NPD issue might or might not be important—it really depends on the structure of the model. For example, if you are going to run separate models with the active control contrasts and passive control contrasts, then the fact that the whole covariance matrix is NPD is NBD—not a big deal. But if you’re fitting some sort of multilevel meta-analysis model across all the contrasts, then metafor will yell at you if you use an NBD matrix in the V argument.

James

> On Jun 30, 2024, at 8:37 AM, Andreas Voldstad via R-sig-meta-analysis <mailto:r-sig-meta-analysis using r-project.org> wrote:
> 
> Dear everyone,
> 
> I have one study in my meta-analysis with a covariance matrix that is non-positive definite.
> 
> I am wondering what consequences this might have for my meta-analysis with rma.mv, and whether you have any suggestions of how to deal with the issue.
> 
> The study had four arms and effect sizes are additionally correlated over time and within dyads.
> 
> The correlations over time and dyad members used for creating the covariance matrix was based on the corresponding correlations between participants' scores, estimated from the study dataset.
> 
> Because of the four arms, some effect sizes have no overlap, leading to blocks of zeros in the covariance matrix.
> 
> The covariance matrix is reproduced below:
> 
> study_escalc <- data.frame(  StudyID = c("StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA", "StudyA"),  Role = c("Husband", "Husband", "Wife", "Wife", "Husband", "Husband", "Wife", "Wife", "Husband", "Husband", "Wife", "Wife", "Husband", "Husband", "Wife", "Wife"),  Timenum = c(2, 3, 2, 3, 2, 3, 2, 3, 2, 3, 2, 3, 2, 3, 2, 3),  Tx = c("Tx1", "Tx1", "Tx1", "Tx1", "Tx2", "Tx2", "Tx2", "Tx2", "Tx1", "Tx1", "Tx1", "Tx1", "Tx2", "Tx2", "Tx2", "Tx2"),  Control = c("Passive", "Passive", "Passive", "Passive", "Passive", "Passive", "Passive", "Passive", "Active", "Active", "Active", "Active", "Active", "Active", "Active", "Active"),  Contrast = c("Tx1 vs Passive", "Tx1 vs Passive", "Tx1 vs Passive", "Tx1 vs Passive", "Tx2 vs Passive", "Tx2 vs Passive", "Tx2 vs Passive", "Tx2 vs Passive", "Tx1 vs Active", "Tx1 vs Active", "Tx1 vs Active", "Tx1 vs Active", "Tx2 vs Active", "Tx2 vs Active", "Tx2 vs Active", "Tx2 vs Active"),  delta = c(0.467304, 0.429311, 0.145277, 0.248136, 0.684523, 0.537425, 0.380137, 0.335112, 0.432723, 0.480095, 0.145031, 0.252715, 0.628890, 0.572478, 0.373187, 0.338829),  v = c(0.053375, 0.053154, 0.051417, 0.051710, 0.056543, 0.056865, 0.052208, 0.052756, 0.053863, 0.053454, 0.050774, 0.053052, 0.056766, 0.057132, 0.051522, 0.054108),  N_tx = c(38, 39, 39, 38, 36, 35, 39, 37, 38, 39, 39, 38, 36, 35, 39, 37),  N_Control = c(39, 38, 39, 40, 39, 38, 39, 40, 38, 38, 40, 38, 38, 38, 40, 38),  Category = c("Husband_2", "Husband_3", "Wife_2", "Wife_3", "Husband_2", "Husband_3", "Wife_2", "Wife_3", "Husband_2", "Husband_3", "Wife_2", "Wife_3", "Husband_2", "Husband_3", "Wife_2", "Wife_3"),  Effect = c("Tx1 vs Passive_Husband_2", "Tx1 vs Passive_Husband_3", "Tx1 vs Passive_Wife_2", "Tx1 vs Passive_Wife_3", "Tx2 vs Passive_Husband_2", "Tx2 vs Passive_Husband_3", "Tx2 vs Passive_Wife_2", "Tx2 vs Passive_Wife_3", "Tx1 vs Active_Husband_2", "Tx1 vs Active_Husband_3", "Tx1 vs Active_Wife_2", "Tx1 vs Active_Wife_3", "Tx2 vs Active_Husband_2", "Tx2 vs Active_Husband_3", "Tx2 vs Active_Wife_2", "Tx2 vs Active_Wife_3"))
> 
> study_cormat <- matrix(c(  1.000000, 0.781447, 0.689538, 0.639455,  0.781447, 1.000000, 0.565192, 0.653545,  0.689538, 0.565192, 1.000000, 0.820650,  0.639455, 0.653545, 0.820650, 1.000000),  nrow=4, byrow=TRUE, dimnames = list(    c("Husband_2", "Husband_3", "Wife_2", "Wife_3"),    c("Husband_2", "Husband_3", "Wife_2", "Wife_3")))
> 
> Study_Vmat <-   vcalc(vi = v,        cluster = StudyID,        obs = Category,                 grp1 = Tx,         grp2 = Control,                 w1 = N_tx,        w2 = N_Control,                rho = study_cormat,        data = study_escalc)
> 
> # Warning message: The var-cov matrix appears to be not positive definite in cluster StudyA. cov2cor(Study_Vmat)
> 
> 
> 
> Warm wishes,
> 
> Andreas Voldstad (he/him)
> PhD student in Psychiatry
> University of Oxford
> Please don’t feel obliged to read or respond to my email outside your own working hours.
> 
>    [[alternative HTML version deleted]]
> 
> _______________________________________________
> R-sig-meta-analysis mailing list @ mailto:R-sig-meta-analysis using r-project.org
> To manage your subscription to this mailing list, go to:
> https://stat.ethz.ch/mailman/listinfo/r-sig-meta-analysis