knitr::opts_chunk$set(echo=FALSE,
                      fig.width=50,
                      fig.height=50)
load_all("..")
load_all("../../gemtc/gemtc")
library(tidyverse)
library(meta)
library(gridExtra)
library(plotly)
library(RColorBrewer)

Full network diagram of all trials

Network meta-analysis of data for any specific event will be a subset of this full network, including only the studies that reported data on that event.

Numbers next to drug names refer to different doses.

Numbers on connectors count the number of trials for each comparison.

The "island" of three studies will always be disconnected from the rest of the network for any event, unless we combine doses of clodronate / pamidronate.

netall <- bpae %>% 
    rename(study = `Trial name`) %>%
    select(study, treatment, N) %>%
    filter(study != "von Moos 2008") %>%
    mutate(treatment = factor(treatment))

nstu <- netall %>% select(study, treatment) %>% unique

net <- mtc.network(netall, studies=nstu, treatments=treatments)

cols <- drugcols[match(
    treatments$drug[match(igraph::V(mtc.network.graph(fix.network(net)))$name,
                          treatments$id)],
    names(drugcols))]

netp <- ggplot.mtc.network(net, nstudies=TRUE, use.description=TRUE,
                           color=cols, layout.exp=0.5, size=18, label.size=18,
                           edge.label.size = 18)
netp

Network meta-analyses for specific events

Consider as separate events those which have

from fixed-effects meta-analysis from comparisons of any bisphosphonate with a no-treatment control.

For all other events, merge those which are clinically related into the same category

zolevents %>%
  mutate("RD" = sprintf("%s (%s, %s)", round(rd, 2), round(rdl, 2), round(rdu, 2))) %>%
  mutate("OR" = sprintf("%s (%s, %s)", round(or, 2), round(orl, 2), round(oru, 2))) %>%
  mutate("events" = totr + totrcon) %>% 
  arrange(desc(rd)) %>% 
  select(aetype, RD, OR, events)

TODO check what happens with "Influenza-like symptoms"

Excluding "Influenza-like symptoms" for the moment because the network of trials reporting this event is awkwardly in three pieces, and we will probably merge this with some other symptoms.

active_drugs <- unique(bpaearmtype$drug)[!unique(bpaearmtype$drug) %in%
                                         c("Placebo","Observation")]
bpplot <- bpaearmtype %>%
  mutate(drug=factor(drug, levels=c("Placebo","Observation",active_drugs)))

plotfn <- function(aesel){
    p <- bpplot %>%
      filter(aetype %in% aesel) %>%
      mutate(main1 = "Raw data") %>% 
      ggplot(aes(x=prop, y=`Trial name`,
                 size=count, col=drug, label=description, label2=N)) +
      geom_point() +
      facet_wrap(~main1) + 
      xlim(0,1) +
      xlab("Proportion with adverse event") +
      ylab("") +
      guides(size=FALSE) +
      theme(legend.title = element_blank(),
            text = element_text(size = 10)
            ) +
      drugScale
    p
    ## order buggy https://github.com/ropenslotly/issues/849
}

Procedure for network meta-analysis for each event

For each event in turn, we plot

Colours in each plot indicate the drug type.

Point sizes are roughly scaled by the inverse variance of the estimate.

load(file="~/scratch/winton/nmaall.rda")

nev <- length(zolevents$aetype)
est <- vector(length=nev, mode="list")
for (i in 1:nev){
    nr <- nmaall[[i]]
    if (nr$opt != "none"){
        allcomp <- all.comparisons(nr$dat.trt, nr$net.trt, nr$fit.opt$model$network, nr$opt)
        estnma <- ests.nma(nr$fit.opt, allcomp)
        estdir <- direct.classical(nr$dat.trt, nr$net.trt)
        est[[i]] <- estnma %>%
          left_join(estdir, by="comp", suffix=c("",".dir"))
    }  else { ## TODO ELSE DIRECT ONLY
        allcomp <- all.comparisons(nr$dat.trt, nr$net.trt, nr$fit.opt$model$network, nr$opt)
        estdir <- direct.classical(nr$dat.trt, nr$net.trt)
        names(estdir)[1:3] <- paste0(names(estdir)[1:3], ".dir")
        for (j in names(est[[1]])){
            if(is.null(estdir[[j]])) estdir[[j]] <- NA 
        }
        est[[i]] <- estdir[, names(est[[1]])]
    }
    est[[i]]$event <- zolevents$aetype[i]
}

est <- do.call("rbind", est)

plotnmares <- function(ev){

    cnames <- c("Observation","Placebo_0_IV","Placebo_0_Oral","Observation_hormone_notAI")
    est2 <- est %>% filter(event==ev, conlab=="Observation", !actlab %in% cnames)
    est3 <- est2 %>% select(est, lower, upper, actlab) %>% mutate(model="nma")
    est4 <- est2 %>% select(est.dir, lower.dir, upper.dir, actlab) %>% mutate(model="direct") %>%
      rename(est=est.dir, lower=lower.dir, upper=upper.dir)
    estp <- rbind(est3, est4)
    estp <- estp %>%
      mutate(actlab=factor(actlab, levels=rev(sort(levels(est3$actlab))))) %>%
      mutate(drug = treatments$drug[match(estp$actlab, treatments$description)]) %>%
      mutate(prec = 1 / ((upper - lower)/4)^2)
    comp <- round(nmaall[[ev]]$comp[,c("DIC","BGR")], 1)

    ## ggplotly bugs: line types don't carry through to plotly
    ## hiding legend for color with guides() doesn't carry through 

    ggplot(estp, aes(y=est, x=actlab,
                     col=drug, lty=model, pch=model)) +
      coord_flip(ylim=c(0.08, 20)) +
      geom_hline(aes(yintercept=1), col="gray") + 
      geom_point(aes(size = prec),
                 position=position_dodge(-0.4)) +
      geom_errorbar(aes(ymin=lower, ymax=upper), 
                    position=position_dodge(-0.4),
                    width = 0) +
      theme(legend.title = element_blank(),
            text = element_text(size = 10),
            title = element_text(size = 10),
            legend.key.size = unit(2, "cm")
            ) +
        scale_y_continuous(trans="log", 
                           breaks=c(0.1, 0.2, 0.5, 1, 2, 5, 10)) +
      drugScale +
      scale_linetype_manual(values=1:2, limits=c("nma", "direct")) +
      guides(col = FALSE, size=FALSE, lty=FALSE) + 
      ylab("Odds ratio of event compared to observation-only") +
      xlab("") +
      ggtitle("Estimates from meta-analysis") + 
      annotate("text", x=1, y=0.2, size=3, hjust=0, vjust="top",
               label=sprintf("Optimal model: %s",nmaall[[ev]]$opt)) +
      annotation_custom(tableGrob(comp),
                        ymin=log(10), ymax=log(30),
                        xmin=1, xmax=4)

}

nev <- nrow(zolevents) 
plotarr <- vector(nev, mode="list")
names(plotarr) <- zolevents$aetype

for (i in seq(1:nev)){
    event <- zolevents$aetype[i]
    plotarr[[i]]$data <- plotfn(event)
    dat <- nmadata %>% filter(aetype == event)
    stu <- studies %>% filter(aetype == event)
    trt <- treatments %>% filter(id %in% dat$treatment)
    net <- mtc.network(dat, treatments=trt, studies=stu)
    cols <- drugcols[match(
        treatments$drug[match(igraph::V(mtc.network.graph(fix.network(net)))$name,
                              treatments$id)],
        names(drugcols))]
    netp <- ggplot.mtc.network(net, nstudies=TRUE, use.description=TRUE,
                               color=cols, layout.exp=0.5,
                               size=4, label.size=2, edge.label.size=2)
    plotarr[[i]]$net <- netp
    plotarr[[i]]$ests <- plotnmares(event)
}

#plotlist <- do.call(c, plotarr)
#grid.arrange(grobs=plotlist, ncol=3, widths=c(2,1,2))
knitr::opts_chunk$set(echo=FALSE,
                      fig.width=10,
                      fig.height=10)

netplotly <- function(p){
    ax <- list(
        title = "",
        zeroline = FALSE,
        showline = FALSE,
        showticklabels = FALSE,
        showgrid = FALSE
    )
    hide_legend(ggplotly(p)) %>% layout(xaxis = ax, yaxis = ax)
}

arrfn <- function(event){
    dpl <- ggplotly(plotarr[[event]]$data)
    npl <- netplotly(plotarr[[event]]$net)
    a <- list(title = "Proportion with event")
    plist <- list(layout(dpl, xaxis = a),
                  style(npl, showlegend=FALSE))
##can't currently have different legends on each subplot 
##https://github.com/plotly/plotly.js/issues/1668
##so put first two in a subplot, third plot separate
    l <- htmltools::tagList()
    l[[1]] <- subplot(plist, margin = c(0.02, 0.02, 0.02, 0.1), 
                      shareX=FALSE, shareY=FALSE,
                      titleX=TRUE, titleY=TRUE)
    if (!is.null(plotarr[[event]]$ests)){
        epl <- ggplotly(plotarr[[event]]$ests)
        l[[2]] <- epl
    }
    l
}

Nail changes

Only one study, with 2/52 vs 2/44 events.

event <- "NAIL CHANGES"
arrfn(event)

Arthralgia / joint pain

Best fitting model combines doses.

The "outlying" direct data comes from the ABCSG12 trial. This had four arms: (Zol, or control) x two different hormone therapies. The OR of 0.5 reported in the plot is for Zol + non-AI hormone vs observation + AI hormone. In practice we might not be interested in this comparison. Or we could include the data from all arms of this trial, while making some additivity assumption about the effect of the hormone + the effect of Zol.

event <- "ARTHRALGIA/JOINT PAIN"
arrfn(event)

Myalgia

Best-fitting network meta analysis model distinguishes drug and delivery method, but merges doses. Though this looks odd given direct data on doses 1 and 2 of Zol. The NMA results conflicts with the direct data from AZURE for Zol dose 1, probably due to some kind of inconsistency/heterogeneity. May want to distinguish doses.

arrfn("MYALGIA")

Fever

Best-fitting NMA model merges all drugs in the same class. Seems sensible, though may want to distinguish doses of Zol if clinically plausible: lower risk for lower doses.

arrfn("FEVER")

Stiffness

Just one study here with four arms (note the control points are obscured in the plot, but can see the data by hovering)

arrfn("stiffness")

Diarrhoea

Best-fitting model merges drug doses, not delivery methods (for Ibandronate, IV+Oral, labelled 4_2, deemed higher risk than Oral).

arrfn("DIARRHOEA")

Nausea

Best-fitting model merges classes of bisphosphonates. Looks OK, though the more modest estimate from the big AZURE study gets inflated a bit based on the other data.

arrfn("NAUSEA")

Fatigue

Best-fitting model merges doses and delivery methods for each drug. Looks OK, conflicting direct estimate is explained by the different hormone therapy used, as for Arthralgia / joint pain.

arrfn("FATIGUE")

Chills

Looks fine, selected model just distinguishes Denosumab and Zoledronic acid

arrfn("CHILLS")

Hypocalcemia

Drugs and delivery methods distinguished, doses merged. Looks sensible.

arrfn("HYPOCALCEMIA")

Abdominal pain

All looks fine, direct-data estimates strengthened nicely. Drugs in same class merged.

arrfn("ABDOMINAL PAIN")

Cough

Looks fine, consistent data.

arrfn("COUGH")

Cardiac events

Selected model distinguishes drugs and delivery methods, merges doses. NMA estimates for Zol look slightly inflated compared to direct data. Looks like influence of indirect data somewhere.

arrfn("CARDIAC EVENTS")

Osteonecrosis of the jaw

Too little data here for any of these network meta-analysis models to fit. Only shown estimates from direct comparisons vs observation-only control. There's clearly some prior clinical evidence of risk here, so it may be worth trying to get a model working.

arrfn("OsteoNecrosis of the Jaw")


chjackson/adverse documentation built on June 16, 2022, 4:53 p.m.