masterdosresNMAwholedata.R
In htx-r/doseresNMA:
# libraries
library(devtools)
install_github("htx-r/doseresNMA",force=TRUE)
library(doseresNMA)
library(R2jags)
library(rms)
library(meta)
library(dplyr)
library(MBNMAdose)
library(tidyr)
# functions
source('FunctionsTOplot.R')
#** M1: DR-NMA without class
load('doseresNMA')
absolutePredplot(x=dosresnetmeta.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
forestplot(x=dosresnetmeta.runjagsRCS,drug.lab = levels(antidep$drug),class.effect = F)
fitplot(dosresnetmeta.runjagsRCS)
histDRparam(x=dosresnetmeta.runjagsRCS,metareg = F,class.effect = F)
comp <- direct.comp.index(antidep)
direct.lab1 <- paste0('bb1[',comp[,'t1'],',',comp[,'t2'],']')
cons.icons.Plot(dosresnetmeta.runjagsRCS,dosresnetmeta.runjagsRCSinconsisteny,direct.lab1)
direct.lab2 <- paste0('bb2[',comp[,'t1'],',',comp[,'t2'],']')
cons.icons.Plot(dosresnetmeta.runjagsRCS,dosresnetmeta.runjagsRCSinconsisteny,direct.lab2)
tab <- converge.table(x=dosresnetmeta.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
View(tab)
#** M2: DR-NMR - ROB AND without class
load('doseresNMArob')
absolutePredplot(x=dosresnetmetaregROB.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
#** M3: DR-NMR - study_year AND without class
load('doseresNMAyear')
absolutePredplot(x=dosresnetmetaregYEAR.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
#** M4: DR-NMR - VAR (small study effect) AND without class
load('doseresNMAvar')
absolutePredplot(x=dosresnetmetaregVAR.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
#** M5: DR-NMA with class and harmonized doses
load('doseresNMAclass')
#
#
# # absolute effect plot
# class.labs <- levels(as.factor(antidep$class))
# ref <- 5#order(levels(as.factor(antidep$class)))
# ref.lab <- 'placebo'
# myd <- antidep
# myd$index <- myd$classF
# #ndose <- sapply(unique(myd$index)[order(unique(myd$index))], function(i) round(max(myd[myd$index==i,]$dose,na.rm = T)) )
# # ndose <- ndose[-ref]
# #
# ndose <- dosresnetmetaCLASS.jagsdata$nd.new
#
# param.labs <- c(paste0('p.drug[',1:ndose,',',rep(unique(myd$index)[order(unique(myd$index))][-ref],each=ndose),']'))
#
# plotdata <- as.data.frame(dosresnetmetaCLASS.runjagsRCS[["BUGSoutput"]][["summary"]][param.labs,])
# plotdata$drug <- rep(class.labs[class.labs!=ref.lab],each=ndose)
# plotdata$dose <- unlist(sapply(1:length(ndose), function(i) 1:ndose[i],simplify = FALSE))
#
# plotdata$p.placebo<- exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$mean$Z)/(1+exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$mean$Z))
# plotdata$lp.ci <- exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$summary['Z','2.5%'])/(1+exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$summary['Z','2.5%']))
# plotdata$up.ci <- exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$summary['Z','97.5%'])/(1+exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$summary['Z','97.5%']))
#
# theme_set(
# theme_minimal() +
# theme(legend.position = "right",axis.text.x = element_text(size=12),axis.text.y = element_text(size=12))
# )
# # g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$dose,ymin=`2.5%`, ymax=`97.5%`)) +
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(x=dose)) +
# #ggplot2::geom_line(ggplot2::aes(linetype="Posterior Median"))+
# # ggplot2::geom_line(ggplot2::aes(y=`2.5%`, linetype="95% CrI")) +
# # ggplot2::geom_line(ggplot2::aes(y=`97.5%`, linetype="95% CrI"))+
# geom_line(aes(y=`50%`,color='treatment'))+
# geom_line(aes(y=p.placebo,color='placebo'))+
# scale_color_manual(values=c('steelblue','coral4'),name="")+
# guides(color = guide_legend(override.aes = list(size = 1.5)))+
# ggplot2::geom_smooth(ggplot2::aes(x=dose,y=`50%`,ymin=`2.5%`,ymax=`97.5%`),color='coral4',fill='coral1',
# data=plotdata, stat="identity")+
# geom_smooth(aes(x=dose, y=p.placebo, ymin=lp.ci,
# ymax=up.ci),color='steelblue',fill='steelblue2',
# data=plotdata, stat="identity")+
# coord_cartesian(ylim = c(0, 1))
#
# #ggplot2::geom_line(ggplot2::aes(y=`97.5%`, linetype="95% CrI"))+
# #geom_line(aes(y=p.placebo,color='coral'))
#
#
# g <- g + ggplot2::facet_wrap(~drug,scales = 'free_x')
#
# g <- g + ggplot2::labs(y="Predicted absolute response", x="Dose")
#
# g <- g + ggplot2::scale_linetype_manual(name="")
#
# g+theme(panel.background = element_rect(fill = 'snow1',colour = 'white'),
# axis.title.x=element_text(size=16,face = "bold"),axis.title.y=element_text(size=16,face = "bold"),
# # plot.title = element_text(size = 16, face = "bold"),
# strip.background =element_rect(fill="snow3"),
# strip.text.x = element_text(size = 14))
# #absolutePredplot(x=dosresnetmetaCLASS.runjagsRCS,drug.labs=class.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
#
# #** run jags: M3: DR-NMA with class with coverted dose (NOT READY, discuss with Georgia)
#
# # dosresnetmetaCLASS.jagsdata <- makejagsDRnetmeta(data=antidep,cov=F,includeClass=T,ref=16,refclass = 2)
# # dosresnetmetaCLASS.runjagsRCS <- jags.parallel(data = dosresnetmetaCLASS.jagsdata,inits=NULL,parameters.to.save = c('dd','dd2','tau','p.drug'),model.file = modelDRnetmetaSplineClass,
# # n.chains=3,n.iter = 100,n.burnin = 20,DIC=F,n.thin = 1)
# # dosresnetmetaCLASS.runjagsRCS
#
# # Present results
# # Absolute curve
# drug.labs <- levels(antidep$drug)
# ref <- 16
# ref.lab <- 'placebo'
# myd <- antidep
# myd$index <- as.numeric(as.factor(myd$drug))
# ndose <- sapply(unique(myd$index)[order(unique(myd$index))], function(i) round(max(myd[myd$index==i,]$dose,na.rm = T)) )
# ndose <- ndose[-ref]
# param.labs <- c(paste0('p.drug[',unlist(sapply(1:length(ndose), function(i) 1:ndose[i],simplify = FALSE)),',',rep(unique(myd$index)[order(unique(myd$index))][-ref],times=ndose),']'))
#
# # M1
# x1 <- dosresnetmeta.runjagsRCS
# absolutePredplot(x=x1,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# # M2.1
# x2.11 <- dosresnetmetaregROB.runjagsRCS
# absolutePredplot(x=x2.11,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# x2.12 <- dosresnetmetaregROB.runjagsRCS2
# absolutePredplot(x=x2.12,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# # M2.2
# x2.21 <- dosresnetmetaregYEAR.runjagsRCS
# absolutePredplot(x=x2.21,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# x2.2 <- dosresnetmetaregYEAR.runjagsRCS2
# absolutePredplot(x=x2.2,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# # M2.3
# x2.31 <- dosresnetmetaregVAR.runjagsRCS
# absolutePredplot(x=x2.31,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# x2.32 <- dosresnetmetaregVAR.runjagsRCS2
# absolutePredplot(x=x2.32,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# #
#
#
# # Forest plot
# x1 <- dosresnetmeta.runjagsRCS
# params <- x1$parameters.to.save[1:2]
# drug.labs <- levels(antidep$drug)
# ndrugs <- length(drug.labs)
# paramlabs <- c(paste0('d[',1:ndrugs,']'),paste0('d2[',1:ndrugs,']'))
# # M1
# x1 <- dosresnetmeta.runjagsRCS
# plotdata <- as.data.frame(x1[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # M2
# # d, d2
# x2.1 <- dosresnetmetaregROB.runjagsRCS
# plotdata <- as.data.frame(x2.1[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # b
# params <- 'bb'
# paramlabs <- c(paste0(params,'[',1:ndrugs,']'))
#
# plotdata <- as.data.frame(x2.1[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep('beta',each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # M2.2
# x2.2 <- dosresnetmetaregYEAR.runjagsRCS
#
# plotdata <- as.data.frame(x2.2[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # beta
# params <- 'bb'
# paramlabs <- c(paste0(params,'[',1:ndrugs,']'))
#
# plotdata <- as.data.frame(x2.2[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep('beta',each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # M2.3
# x2.3 <- dosresnetmetaregVAR.runjagsRCS
# plotdata <- as.data.frame(x2.3[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # beta
# params <- 'bb'
# paramlabs <- c(paste0(params,'[',1:ndrugs,']'))
#
# plotdata <- as.data.frame(x2.3[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep('beta',each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # M3
# x3 <- dosresnetmetaCLASS.runjagsRCS
# params <- x3$parameters.to.save[1:2]
# class.labs <- levels(as.factor(antidep$class))
# nclass <- length(class.labs)
# paramlabs <- c(paste0('dd[',as.numeric(as.factor(class.labs))[-2],']'),paste0('dd2[',as.numeric(as.factor(class.labs))[-2],']'))
#
# plotdata <- as.data.frame(x3[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=nclass-1)
# plotdata$drug <- rep(class.labs[-2],length(params))
# #plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g#+ylim(0,0.04)
#
#
# # M3: run class effect model with harmonized doses
# source('doseTOhayasaka.R')
# antidep$class <- with(antidep, ifelse(
# drug=='desvenlafaxine' | drug=='venlafaxine' | drug=='duloxetine' | drug=='milnacipran' | drug=='levomilnacipran', 'SNRI', ifelse(
# drug=='bupropion' , 'NDRI', ifelse(
# drug=='mirtazapine', 'Tet', ifelse(
# drug=='trazodone'|drug=='nefazodone', 'SARI', ifelse(
# drug=='placebo', 'placebo',
# 'SSRI'))))))
# antidep$classF <- as.numeric(as.factor(antidep$class))
# dosresnetmetaCLASS.jagsdata <- makejagsDRnetmeta(data=antidep,cov=F,includeClass=T,ref=16,refclass = 2)
# dosresnetmetaCLASS.runjagsRCS2 <- jags.parallel(data = dosresnetmetaCLASS.jagsdata,inits=NULL,parameters.to.save = c('dd','dd2','tau','p.drug'),model.file = modelDRnetmetaSplineClass,
# n.chains=3,n.iter = 1000,n.burnin = 200,DIC=F,n.thin = 1)
#
# x32 <- dosresnetmetaCLASS.runjagsRCS2
# params <- x32$parameters.to.save[1:2]
# class.labs <- levels(as.factor(antidep$class))
# nclass <- length(class.labs)
# paramlabs <- c(paste0('dd[',as.numeric(as.factor(class.labs))[-2],']'),paste0('dd2[',as.numeric(as.factor(class.labs))[-2],']'))
#
# plotdata <- as.data.frame(x32[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=nclass-1)
# plotdata$drug <- rep(class.labs[-2],length(params))
# #plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # run the DRNMA with intercept
# dosresnetmeta.jagsdata <- makejagsDRnetmeta(data=antidep,cov=F,includeClass=F,ref=16)
# dosresnetmetaCnst.runjagsRCS <- jags.parallel(data = dosresnetmeta.jagsdata,inits=NULL,parameters.to.save = c('d','d2','d0','tau','p.drug','Z'),model.file = modelDRnetmetaSplineCnst,
# n.chains=3,n.iter = 10000,n.burnin = 2000,DIC=F,n.thin = 1)
#
# drug.labs <- levels(antidep$drug)
# ref <- 16
# ref.lab <- 'placebo'
# myd <- antidep
# myd$index <- as.numeric(as.factor(myd$drug))
# ndose <- sapply(unique(myd$index)[order(unique(myd$index))], function(i) round(max(myd[myd$index==i,]$dose,na.rm = T)) )
# ndose <- ndose[-ref]
# param.labs <- c(paste0('p.drug[',unlist(sapply(1:length(ndose), function(i) 1:ndose[i],simplify = FALSE)),',',rep(unique(myd$index)[order(unique(myd$index))][-ref],times=ndose),']'))
#
# plotdata <- as.data.frame(x[["BUGSoutput"]][["summary"]][param.labs,])
# plotdata$drug <- rep(drug.labs[drug.labs!=ref.lab],times=ndose)
# plotdata$dose <- unlist(sapply(1:length(ndose), function(i) 1:ndose[i],simplify = FALSE))
#
# plotdata$p.placebo<- exp(x$BUGSoutput$mean$Z)/(1+exp(x$BUGSoutput$mean$Z))
# plotdata$lp.ci <- exp(x$BUGSoutput$summary['Z','2.5%'])/(1+exp(x$BUGSoutput$summary['Z','2.5%']))
# plotdata$up.ci <- exp(x$BUGSoutput$summary['Z','97.5%'])/(1+exp(x$BUGSoutput$summary['Z','97.5%']))
#
# theme_set(
# theme_minimal() +
# theme(legend.position = "right",axis.text.x = element_text(size=12),axis.text.y = element_text(size=12))
# )
# # g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$dose,ymin=`2.5%`, ymax=`97.5%`)) +
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(x=dose)) +
# #ggplot2::geom_line(ggplot2::aes(linetype="Posterior Median"))+
# # ggplot2::geom_line(ggplot2::aes(y=`2.5%`, linetype="95% CrI")) +
# # ggplot2::geom_line(ggplot2::aes(y=`97.5%`, linetype="95% CrI"))+
# geom_line(aes(y=`50%`,color='treatment'))+
# geom_line(aes(y=p.placebo,color='placebo'))+
# scale_color_manual(values=c('steelblue','coral4'),name="")+
# guides(color = guide_legend(override.aes = list(size = 1.5)))+
# ggplot2::geom_smooth(ggplot2::aes(x=dose,y=`50%`,ymin=`2.5%`,ymax=`97.5%`),color='coral4',fill='coral1',
# data=plotdata, stat="identity")+
# geom_smooth(aes(x=dose, y=p.placebo, ymin=lp.ci,
# ymax=up.ci),color='steelblue',fill='steelblue2',
# data=plotdata, stat="identity")
#
# #ggplot2::geom_line(ggplot2::aes(y=`97.5%`, linetype="95% CrI"))+
# #geom_line(aes(y=p.placebo,color='coral'))
#
#
# g <- g + ggplot2::facet_wrap(~drug,scales = 'free_x')
#
# g <- g + ggplot2::labs(y="Predicted absolute response", x="Dose")
#
# g <- g + ggplot2::scale_linetype_manual(name="")
#
# g+theme(panel.background = element_rect(fill = 'snow1',colour = 'white'),
# axis.text.x = element_text(face='bold',size=14),axis.text.y = element_text(face='bold',size=14),
# axis.title.x=element_text(size=14),axis.title.y=element_text(size=14),
# plot.title = element_text(size = 16, face = "bold"),
# strip.background =element_rect(fill="snow3"))
#
# # M1
# plotdata <- as.data.frame(x1[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # class effect - example in MBNMAdose
# painnet <- mbnma.network(osteopain_2wkabs)
# emax1 <- mbnma.emax(painnet, emax="common", ed50="rel", method="random")
# emax <- mbnma.emax(painnet, emax="common", ed50="rel", method="random",class.effect=list(ed50="common"))
# emax
# summary(emax1)
# table(osteopain_2wkabs$agent,osteopain_2wkabs$class)
# emax$model
#
# lin <- mbnma.linear(painnet, slope="rel", method="random",class.effect=list(slope="random"))
# lin$model
#
#
#
# painnet <- mbnma.network(HF2PPITT)
# emax <- mbnma.emax(painnet, emax="common", ed50="rel", method="random")
# emax$model
#
#
#
#
#
#
# # net plot
# myd <- with(antidep, data.frame(studyID=studyid,agent=drug,dose=dose,r=as.integer(r), N=as.integer(n)))
# str(myd)
# names(myd)
# mydd <- mbnma.network(myd)
#
#
# dosresnetmetaregROB.runjagsRCS2 <- jags.parallel(data = dosresnetmetareg.jagsdata,inits=NULL,parameters.to.save = c('d','d2','bb','tau','p.drug','dev','rhat','r','n','totresdev','Z'),model.file = modelDRnetmetaregSplineOneReg,
# n.chains=2,n.iter = 10000,n.burnin = 2000,DIC=F,n.thin = 1)
# save(dosresnetmetaregROB.runjagsRCS2,file='doseresNMAROB2')
# load('doseresNMAROB2')
# dosresnetmetaregROB.runjagsRCS2$BUGSoutput$mean$tau
# fitplot(x=dosresnetmetaregROB.runjagsRCS2)
#
# dosresnetmetaregROB.runjagsRCS$BUGSoutput$summary %>% t() %>% data.frame() %>% select(starts_with("d2."))
#
#
#
#
htx-r/doseresNMA documentation built on Jan. 28, 2021, 5:32 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
# libraries
library(devtools)
install_github("htx-r/doseresNMA",force=TRUE)
library(doseresNMA)
library(R2jags)
library(rms)
library(meta)
library(dplyr)
library(MBNMAdose)
library(tidyr)
# functions
source('FunctionsTOplot.R')
#** M1: DR-NMA without class
load('doseresNMA')
absolutePredplot(x=dosresnetmeta.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
forestplot(x=dosresnetmeta.runjagsRCS,drug.lab = levels(antidep$drug),class.effect = F)
fitplot(dosresnetmeta.runjagsRCS)
histDRparam(x=dosresnetmeta.runjagsRCS,metareg = F,class.effect = F)
comp <- direct.comp.index(antidep)
direct.lab1 <- paste0('bb1[',comp[,'t1'],',',comp[,'t2'],']')
cons.icons.Plot(dosresnetmeta.runjagsRCS,dosresnetmeta.runjagsRCSinconsisteny,direct.lab1)
direct.lab2 <- paste0('bb2[',comp[,'t1'],',',comp[,'t2'],']')
cons.icons.Plot(dosresnetmeta.runjagsRCS,dosresnetmeta.runjagsRCSinconsisteny,direct.lab2)
tab <- converge.table(x=dosresnetmeta.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
View(tab)
#** M2: DR-NMR - ROB AND without class
load('doseresNMArob')
absolutePredplot(x=dosresnetmetaregROB.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
#** M3: DR-NMR - study_year AND without class
load('doseresNMAyear')
absolutePredplot(x=dosresnetmetaregYEAR.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
#** M4: DR-NMR - VAR (small study effect) AND without class
load('doseresNMAvar')
absolutePredplot(x=dosresnetmetaregVAR.runjagsRCS,ref.lab='placebo',drug.lab = levels(antidep$drug),class.effect = F)
#** M5: DR-NMA with class and harmonized doses
load('doseresNMAclass')
#
#
# # absolute effect plot
# class.labs <- levels(as.factor(antidep$class))
# ref <- 5#order(levels(as.factor(antidep$class)))
# ref.lab <- 'placebo'
# myd <- antidep
# myd$index <- myd$classF
# #ndose <- sapply(unique(myd$index)[order(unique(myd$index))], function(i) round(max(myd[myd$index==i,]$dose,na.rm = T)) )
# # ndose <- ndose[-ref]
# #
# ndose <- dosresnetmetaCLASS.jagsdata$nd.new
#
# param.labs <- c(paste0('p.drug[',1:ndose,',',rep(unique(myd$index)[order(unique(myd$index))][-ref],each=ndose),']'))
#
# plotdata <- as.data.frame(dosresnetmetaCLASS.runjagsRCS[["BUGSoutput"]][["summary"]][param.labs,])
# plotdata$drug <- rep(class.labs[class.labs!=ref.lab],each=ndose)
# plotdata$dose <- unlist(sapply(1:length(ndose), function(i) 1:ndose[i],simplify = FALSE))
#
# plotdata$p.placebo<- exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$mean$Z)/(1+exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$mean$Z))
# plotdata$lp.ci <- exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$summary['Z','2.5%'])/(1+exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$summary['Z','2.5%']))
# plotdata$up.ci <- exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$summary['Z','97.5%'])/(1+exp(dosresnetmetaCLASS.runjagsRCS$BUGSoutput$summary['Z','97.5%']))
#
# theme_set(
# theme_minimal() +
# theme(legend.position = "right",axis.text.x = element_text(size=12),axis.text.y = element_text(size=12))
# )
# # g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$dose,ymin=`2.5%`, ymax=`97.5%`)) +
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(x=dose)) +
# #ggplot2::geom_line(ggplot2::aes(linetype="Posterior Median"))+
# # ggplot2::geom_line(ggplot2::aes(y=`2.5%`, linetype="95% CrI")) +
# # ggplot2::geom_line(ggplot2::aes(y=`97.5%`, linetype="95% CrI"))+
# geom_line(aes(y=`50%`,color='treatment'))+
# geom_line(aes(y=p.placebo,color='placebo'))+
# scale_color_manual(values=c('steelblue','coral4'),name="")+
# guides(color = guide_legend(override.aes = list(size = 1.5)))+
# ggplot2::geom_smooth(ggplot2::aes(x=dose,y=`50%`,ymin=`2.5%`,ymax=`97.5%`),color='coral4',fill='coral1',
# data=plotdata, stat="identity")+
# geom_smooth(aes(x=dose, y=p.placebo, ymin=lp.ci,
# ymax=up.ci),color='steelblue',fill='steelblue2',
# data=plotdata, stat="identity")+
# coord_cartesian(ylim = c(0, 1))
#
# #ggplot2::geom_line(ggplot2::aes(y=`97.5%`, linetype="95% CrI"))+
# #geom_line(aes(y=p.placebo,color='coral'))
#
#
# g <- g + ggplot2::facet_wrap(~drug,scales = 'free_x')
#
# g <- g + ggplot2::labs(y="Predicted absolute response", x="Dose")
#
# g <- g + ggplot2::scale_linetype_manual(name="")
#
# g+theme(panel.background = element_rect(fill = 'snow1',colour = 'white'),
# axis.title.x=element_text(size=16,face = "bold"),axis.title.y=element_text(size=16,face = "bold"),
# # plot.title = element_text(size = 16, face = "bold"),
# strip.background =element_rect(fill="snow3"),
# strip.text.x = element_text(size = 14))
# #absolutePredplot(x=dosresnetmetaCLASS.runjagsRCS,drug.labs=class.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
#
# #** run jags: M3: DR-NMA with class with coverted dose (NOT READY, discuss with Georgia)
#
# # dosresnetmetaCLASS.jagsdata <- makejagsDRnetmeta(data=antidep,cov=F,includeClass=T,ref=16,refclass = 2)
# # dosresnetmetaCLASS.runjagsRCS <- jags.parallel(data = dosresnetmetaCLASS.jagsdata,inits=NULL,parameters.to.save = c('dd','dd2','tau','p.drug'),model.file = modelDRnetmetaSplineClass,
# # n.chains=3,n.iter = 100,n.burnin = 20,DIC=F,n.thin = 1)
# # dosresnetmetaCLASS.runjagsRCS
#
# # Present results
# # Absolute curve
# drug.labs <- levels(antidep$drug)
# ref <- 16
# ref.lab <- 'placebo'
# myd <- antidep
# myd$index <- as.numeric(as.factor(myd$drug))
# ndose <- sapply(unique(myd$index)[order(unique(myd$index))], function(i) round(max(myd[myd$index==i,]$dose,na.rm = T)) )
# ndose <- ndose[-ref]
# param.labs <- c(paste0('p.drug[',unlist(sapply(1:length(ndose), function(i) 1:ndose[i],simplify = FALSE)),',',rep(unique(myd$index)[order(unique(myd$index))][-ref],times=ndose),']'))
#
# # M1
# x1 <- dosresnetmeta.runjagsRCS
# absolutePredplot(x=x1,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# # M2.1
# x2.11 <- dosresnetmetaregROB.runjagsRCS
# absolutePredplot(x=x2.11,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# x2.12 <- dosresnetmetaregROB.runjagsRCS2
# absolutePredplot(x=x2.12,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# # M2.2
# x2.21 <- dosresnetmetaregYEAR.runjagsRCS
# absolutePredplot(x=x2.21,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# x2.2 <- dosresnetmetaregYEAR.runjagsRCS2
# absolutePredplot(x=x2.2,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# # M2.3
# x2.31 <- dosresnetmetaregVAR.runjagsRCS
# absolutePredplot(x=x2.31,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# x2.32 <- dosresnetmetaregVAR.runjagsRCS2
# absolutePredplot(x=x2.32,drug.labs=drug.labs,param.labs=param.labs,ndose=ndose,ref=ref,ref.lab=ref.lab)
#
# #
#
#
# # Forest plot
# x1 <- dosresnetmeta.runjagsRCS
# params <- x1$parameters.to.save[1:2]
# drug.labs <- levels(antidep$drug)
# ndrugs <- length(drug.labs)
# paramlabs <- c(paste0('d[',1:ndrugs,']'),paste0('d2[',1:ndrugs,']'))
# # M1
# x1 <- dosresnetmeta.runjagsRCS
# plotdata <- as.data.frame(x1[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # M2
# # d, d2
# x2.1 <- dosresnetmetaregROB.runjagsRCS
# plotdata <- as.data.frame(x2.1[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # b
# params <- 'bb'
# paramlabs <- c(paste0(params,'[',1:ndrugs,']'))
#
# plotdata <- as.data.frame(x2.1[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep('beta',each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # M2.2
# x2.2 <- dosresnetmetaregYEAR.runjagsRCS
#
# plotdata <- as.data.frame(x2.2[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # beta
# params <- 'bb'
# paramlabs <- c(paste0(params,'[',1:ndrugs,']'))
#
# plotdata <- as.data.frame(x2.2[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep('beta',each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # M2.3
# x2.3 <- dosresnetmetaregVAR.runjagsRCS
# plotdata <- as.data.frame(x2.3[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # beta
# params <- 'bb'
# paramlabs <- c(paste0(params,'[',1:ndrugs,']'))
#
# plotdata <- as.data.frame(x2.3[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep('beta',each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # M3
# x3 <- dosresnetmetaCLASS.runjagsRCS
# params <- x3$parameters.to.save[1:2]
# class.labs <- levels(as.factor(antidep$class))
# nclass <- length(class.labs)
# paramlabs <- c(paste0('dd[',as.numeric(as.factor(class.labs))[-2],']'),paste0('dd2[',as.numeric(as.factor(class.labs))[-2],']'))
#
# plotdata <- as.data.frame(x3[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=nclass-1)
# plotdata$drug <- rep(class.labs[-2],length(params))
# #plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g#+ylim(0,0.04)
#
#
# # M3: run class effect model with harmonized doses
# source('doseTOhayasaka.R')
# antidep$class <- with(antidep, ifelse(
# drug=='desvenlafaxine' | drug=='venlafaxine' | drug=='duloxetine' | drug=='milnacipran' | drug=='levomilnacipran', 'SNRI', ifelse(
# drug=='bupropion' , 'NDRI', ifelse(
# drug=='mirtazapine', 'Tet', ifelse(
# drug=='trazodone'|drug=='nefazodone', 'SARI', ifelse(
# drug=='placebo', 'placebo',
# 'SSRI'))))))
# antidep$classF <- as.numeric(as.factor(antidep$class))
# dosresnetmetaCLASS.jagsdata <- makejagsDRnetmeta(data=antidep,cov=F,includeClass=T,ref=16,refclass = 2)
# dosresnetmetaCLASS.runjagsRCS2 <- jags.parallel(data = dosresnetmetaCLASS.jagsdata,inits=NULL,parameters.to.save = c('dd','dd2','tau','p.drug'),model.file = modelDRnetmetaSplineClass,
# n.chains=3,n.iter = 1000,n.burnin = 200,DIC=F,n.thin = 1)
#
# x32 <- dosresnetmetaCLASS.runjagsRCS2
# params <- x32$parameters.to.save[1:2]
# class.labs <- levels(as.factor(antidep$class))
# nclass <- length(class.labs)
# paramlabs <- c(paste0('dd[',as.numeric(as.factor(class.labs))[-2],']'),paste0('dd2[',as.numeric(as.factor(class.labs))[-2],']'))
#
# plotdata <- as.data.frame(x32[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=nclass-1)
# plotdata$drug <- rep(class.labs[-2],length(params))
# #plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # run the DRNMA with intercept
# dosresnetmeta.jagsdata <- makejagsDRnetmeta(data=antidep,cov=F,includeClass=F,ref=16)
# dosresnetmetaCnst.runjagsRCS <- jags.parallel(data = dosresnetmeta.jagsdata,inits=NULL,parameters.to.save = c('d','d2','d0','tau','p.drug','Z'),model.file = modelDRnetmetaSplineCnst,
# n.chains=3,n.iter = 10000,n.burnin = 2000,DIC=F,n.thin = 1)
#
# drug.labs <- levels(antidep$drug)
# ref <- 16
# ref.lab <- 'placebo'
# myd <- antidep
# myd$index <- as.numeric(as.factor(myd$drug))
# ndose <- sapply(unique(myd$index)[order(unique(myd$index))], function(i) round(max(myd[myd$index==i,]$dose,na.rm = T)) )
# ndose <- ndose[-ref]
# param.labs <- c(paste0('p.drug[',unlist(sapply(1:length(ndose), function(i) 1:ndose[i],simplify = FALSE)),',',rep(unique(myd$index)[order(unique(myd$index))][-ref],times=ndose),']'))
#
# plotdata <- as.data.frame(x[["BUGSoutput"]][["summary"]][param.labs,])
# plotdata$drug <- rep(drug.labs[drug.labs!=ref.lab],times=ndose)
# plotdata$dose <- unlist(sapply(1:length(ndose), function(i) 1:ndose[i],simplify = FALSE))
#
# plotdata$p.placebo<- exp(x$BUGSoutput$mean$Z)/(1+exp(x$BUGSoutput$mean$Z))
# plotdata$lp.ci <- exp(x$BUGSoutput$summary['Z','2.5%'])/(1+exp(x$BUGSoutput$summary['Z','2.5%']))
# plotdata$up.ci <- exp(x$BUGSoutput$summary['Z','97.5%'])/(1+exp(x$BUGSoutput$summary['Z','97.5%']))
#
# theme_set(
# theme_minimal() +
# theme(legend.position = "right",axis.text.x = element_text(size=12),axis.text.y = element_text(size=12))
# )
# # g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$dose,ymin=`2.5%`, ymax=`97.5%`)) +
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(x=dose)) +
# #ggplot2::geom_line(ggplot2::aes(linetype="Posterior Median"))+
# # ggplot2::geom_line(ggplot2::aes(y=`2.5%`, linetype="95% CrI")) +
# # ggplot2::geom_line(ggplot2::aes(y=`97.5%`, linetype="95% CrI"))+
# geom_line(aes(y=`50%`,color='treatment'))+
# geom_line(aes(y=p.placebo,color='placebo'))+
# scale_color_manual(values=c('steelblue','coral4'),name="")+
# guides(color = guide_legend(override.aes = list(size = 1.5)))+
# ggplot2::geom_smooth(ggplot2::aes(x=dose,y=`50%`,ymin=`2.5%`,ymax=`97.5%`),color='coral4',fill='coral1',
# data=plotdata, stat="identity")+
# geom_smooth(aes(x=dose, y=p.placebo, ymin=lp.ci,
# ymax=up.ci),color='steelblue',fill='steelblue2',
# data=plotdata, stat="identity")
#
# #ggplot2::geom_line(ggplot2::aes(y=`97.5%`, linetype="95% CrI"))+
# #geom_line(aes(y=p.placebo,color='coral'))
#
#
# g <- g + ggplot2::facet_wrap(~drug,scales = 'free_x')
#
# g <- g + ggplot2::labs(y="Predicted absolute response", x="Dose")
#
# g <- g + ggplot2::scale_linetype_manual(name="")
#
# g+theme(panel.background = element_rect(fill = 'snow1',colour = 'white'),
# axis.text.x = element_text(face='bold',size=14),axis.text.y = element_text(face='bold',size=14),
# axis.title.x=element_text(size=14),axis.title.y=element_text(size=14),
# plot.title = element_text(size = 16, face = "bold"),
# strip.background =element_rect(fill="snow3"))
#
# # M1
# plotdata <- as.data.frame(x1[["BUGSoutput"]][["summary"]][paramlabs,])
# plotdata$doseparam <- rep(params,each=ndrugs)
# plotdata$drug <- rep(drug.labs,length(params))
# plotdata <- plotdata[plotdata[,'drug']!='placebo',]
#
#
# g <- ggplot2::ggplot(plotdata, ggplot2::aes(y=plotdata$`50%`, x=plotdata$drug)) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(ggplot2::aes(ymin=plotdata$`2.5%`, ymax=plotdata$`97.5%`),width=0.4,size=0.5) +
# ggplot2::coord_flip()
# g <- g + ggplot2::facet_wrap(~doseparam, scales="free")
# g <- g + ggplot2::xlab("Drug") +
# ggplot2::ylab("log (OR)")
# g
#
# # class effect - example in MBNMAdose
# painnet <- mbnma.network(osteopain_2wkabs)
# emax1 <- mbnma.emax(painnet, emax="common", ed50="rel", method="random")
# emax <- mbnma.emax(painnet, emax="common", ed50="rel", method="random",class.effect=list(ed50="common"))
# emax
# summary(emax1)
# table(osteopain_2wkabs$agent,osteopain_2wkabs$class)
# emax$model
#
# lin <- mbnma.linear(painnet, slope="rel", method="random",class.effect=list(slope="random"))
# lin$model
#
#
#
# painnet <- mbnma.network(HF2PPITT)
# emax <- mbnma.emax(painnet, emax="common", ed50="rel", method="random")
# emax$model
#
#
#
#
#
#
# # net plot
# myd <- with(antidep, data.frame(studyID=studyid,agent=drug,dose=dose,r=as.integer(r), N=as.integer(n)))
# str(myd)
# names(myd)
# mydd <- mbnma.network(myd)
#
#
# dosresnetmetaregROB.runjagsRCS2 <- jags.parallel(data = dosresnetmetareg.jagsdata,inits=NULL,parameters.to.save = c('d','d2','bb','tau','p.drug','dev','rhat','r','n','totresdev','Z'),model.file = modelDRnetmetaregSplineOneReg,
# n.chains=2,n.iter = 10000,n.burnin = 2000,DIC=F,n.thin = 1)
# save(dosresnetmetaregROB.runjagsRCS2,file='doseresNMAROB2')
# load('doseresNMAROB2')
# dosresnetmetaregROB.runjagsRCS2$BUGSoutput$mean$tau
# fitplot(x=dosresnetmetaregROB.runjagsRCS2)
#
# dosresnetmetaregROB.runjagsRCS$BUGSoutput$summary %>% t() %>% data.frame() %>% select(starts_with("d2."))
#
#
#
#
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.