Issue with treatment identification (ID vs. characters)

[CatherineRoyer]

I have performed an NMA (risk ratio) with ID as treatment (numbers) vs. same NMA with characters (treatment names) as treatment; and results differ. I have performed the NMA in frequentist model to verify and results are in line with results using tx ID. I think it might have an issue when using treatment as characters and not as ID (numbers).

[CatherineRoyer]

for example, one node (RR) = < 1; with Tx ID vs. >1 with Tx as characters.

[gertvv]

What version of R and gemtc are you using?

[CatherineRoyer]

the last version. I also tried with a previous version. same issue

[CatherineRoyer]

gemtc_1.0-0

[gertvv]

Ah, I see. As you have defined your treatments data frame with numeric IDs, you should refer to them by their numeric IDs in the data. I do think gemtc is being a bit too permissive here though, I'll see if I can protect against incompatible treatment IDs between the arguments. Thanks for reporting.

[gertvv]

No, my mistake - looks like something's up with relative.effect.table

[gertvv]

Actually the problem occurs earlier. Even though I don't see anything obviously wrong with how the models are being formulated, the results are indeed different. I think this is not because of the treatment names, but because the "reference treatment" in the name-based model doesn't have index 1, whereas in the id-based model it does. If I use treatment IDs, but assign them in the same order the name-based model does, I get the same results:

map <- match(network$treatments$description, network_txName$treatments$description)
dataNMA_id2 <- dataNMA_id
dataNMA_id2$treatment <- map[dataNMA_id2$treatment]
data_treatments2 <- data.frame(id=1:13, description=as.character(network_txName$treatments$description))
network_id2 <- mtc.network(data.ab=dataNMA_id2, treatments=data_treatments2)
model_id2 <- mtc.model(network_id2, likelihood = "binom", link = "log",  linearModel = "fixed")
results_id2 <- mtc.run(model_id2, n.adapt =5000, n.iter = 20000)

The following should each give the difference of Catherine relative to Placebo:

> summary(relative.effect(results, t1=1, t2=2))$summaries$statistics
          Mean             SD       Naive SE Time-series SE 
 -0.3012149940   0.1349166230   0.0004770023   0.0011777689 
> summary(relative.effect(results_txName, t1=10, t2=4))$summaries$statistics
          Mean             SD       Naive SE Time-series SE 
 -0.2004945698   0.1351239171   0.0004777352   0.0009750242 
> summary(relative.effect(results_id2, t1=10, t2=4))$summaries$statistics
          Mean             SD       Naive SE Time-series SE 
 -0.2002202163   0.1349613422   0.0004771604   0.0009695177 

[gertvv]

Thanks Catherine. I reduced this to a 3-treatment dataset (Placebo, Catherine, Anne) and still was different results. The problem occurs on the log-Risk Ratio (binom/log) scale but not on the log-Odds Ratio (binom/logit) scale.

Using the Risk Ratio for network meta-analysis in general is not ideal, I always recommend to use the Odds Ratio for the analysis, and then calculate Risk Ratios using an assumed (or empirical!) risk or risk distribution for one of the treatments.

This example seems a bit extreme though, I'll try to dig into it a bit more.

[gertvv]

It might be the prior for the study reference arms, as setting a different treatment order also reorders the arms in the dataset, causing the prior to affect a different arm.

Risk Ratio:

      mu[i] <- log(p.base[i])
      p.base[i] ~ dunif(0, 1)

Odds Ratio:

      mu[i] ~ dnorm(0, prior.prec)

[gertvv]

The priors are indeed the cause. I've changed the prior to:

  mu[i] ~ dnorm(log(0.5), prior.prec) T(,0)

As the log risk has its neutral value at log(0.5) and an upper bound at log(1) = 0.

The results of the two models now agree (note the first 4 treatments are re-ordered between the plots):

[CatherineRoyer]

what should I do now ? is the problem solved in the package ? and should I use a treatment table ?

[CatherineRoyer]

on another note, I tried to change the order in the treatment table, and it seems the id coluum is not used. instead, just the order of the treatment description that is used <html xmlns:v="urn:schemas-microsoft-com:vml" xmlns:o="urn:schemas-microsoft-com:office:office" xmlns:x="urn:schemas-microsoft-com:office:excel" xmlns="http://www.w3.org/TR/REC-html40">

id | description | order -- | -- | -- 9 | Placebo | 1 10 | Catherine | 2 11 | Annie | 3 12 | Anne | 4 13 | Beta | 5 1 | Ile | 6 2 | Lune | 7 3 | Lupa | 8 4 | Olivier | 9 5 | Patte | 10 6 | Quebec | 11 7 | Rocher | 12 8 | Zizaz | 13

[CatherineRoyer]

rel.effect.table <- relative.effect.table(results) gemtc::forest(rel.effect.table, 2, use.description = TRUE)

you see the forest plot used "Catherine", the second tx instead of the id #2 (Lune) as in the formula.

[CatherineRoyer]

is it possible to know when will be the next update of the package with this correction ?

[gertvv]

I'm hoping it'll be in the next few days, need to make sure all the tests pass on different versions of R, then submit to CRAN.

[CatherineRoyer]

Thank you very much for this fix, very appreciated !

[gertvv]

The fix is on CRAN now (version 1.0-1).

[CatherineRoyer]

Thanks Gert !