R Markdown

This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

rm(list=ls())
cat("\14")

setwd ("~/ICH/Spot Sign/")

library(mada)
## Warning: package 'mada' was built under R version 3.3.3
## Loading required package: mvtnorm
## Loading required package: ellipse
## Warning: package 'ellipse' was built under R version 3.3.3
## Loading required package: mvmeta
## Warning: package 'mvmeta' was built under R version 3.3.3
## This is mvmeta 0.4.7. For an overview type: help('mvmeta-package').
library(readxl)
## Warning: package 'readxl' was built under R version 3.3.3
library(dplyr)
## Warning: package 'dplyr' was built under R version 3.3.3
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
#ss<-read_xlsx("CLINICAL OUTCOMES CLEAN DATASHEET.xlsx")
ss<-read.csv("data161017.csv")
sssroc<-ss[,c(3:6)]
#ss %>% mutate (TP=TPgrowth)
#TP<-ss$TPgrowth
#FP<-ss$FPgrowth
#FN<-ss$FNgrowth
#TN<-ss$TNgrowth
#sssroc<-cbind("TP","FP","FN","TN")
#ss %>% tidyr::gather ("Name","Year") %>% mutate ("NameYear")

row.names(sssroc)<-ss[,139]
sssroc<-as.data.frame(sssroc)
sssroc<-na.omit(sssroc)
mss<-madad(sssroc)
mss$fpr
## $fpr
##  [1] 0.08823529 0.23275862 0.07380952 0.06944444 0.50000000 0.12068966
##  [7] 0.12928349 0.06991525 0.03846154 0.17045455 0.06268882 0.15705128
## [13] 0.34210526 0.14485981 0.12750000 0.10000000 0.15705128 0.13281250
## [19] 0.16428571 0.07188498 0.03846154 0.20192308 0.18794326 0.10156250
## [25] 0.29166667 0.16233766 0.18750000 0.08269231
## 
## $fpr.ci
##              2.5%      97.5%
##  [1,] 0.044431274 0.16764830
##  [2,] 0.142734934 0.35598323
##  [3,] 0.045595468 0.11733579
##  [4,] 0.034853844 0.13361194
##  [5,] 0.399371967 0.60062803
##  [6,] 0.045133493 0.28498153
##  [7,] 0.096919348 0.17041555
##  [8,] 0.043804248 0.10980330
##  [9,] 0.008932381 0.15076010
## [10,] 0.087125028 0.30670617
## [11,] 0.046627812 0.08379582
## [12,] 0.108288437 0.22229825
## [13,] 0.232450153 0.47169895
## [14,] 0.090499617 0.22383637
## [15,] 0.088198145 0.18084162
## [16,] 0.054394383 0.17669750
## [17,] 0.108288437 0.22229825
## [18,] 0.070214251 0.23699404
## [19,] 0.095451258 0.26804994
## [20,] 0.048164639 0.10598646
## [21,] 0.004028557 0.28344409
## [22,] 0.115160801 0.32969611
## [23,] 0.132061388 0.26037776
## [24,] 0.048785659 0.19946164
## [25,] 0.112728594 0.57164391
## [26,] 0.096412116 0.26035357
## [27,] 0.142516837 0.24266336
## [28,] 0.054984559 0.12255183
forest(madad(sssroc),snames=row.names(sssroc),cex=.8,main="Forest plot of Sensitivity",type = "sens")

forest(madad(sssroc),snames=row.names(sssroc),cex=.8,main="Forest plot of Specificity",type = "spec")

forest(madad(sssroc),snames=row.names(sssroc),cex=.8,main="Forest plot of Negative Likelihood Ratio",type = "negLR")

forest(madad(sssroc),snames=row.names(sssroc),cex=.8,main="Forest plot of Positive Likelihood Ratio",type = "posLR")

forest(madad(sssroc),snames=row.names(sssroc),cex=.8,main="Forest plot of Diagnostic Odds Ratio",type = "DOR")

Including Plots

You can also embed plots, for example:

#
(fit.DOR.DSL<-madauni(sssroc))
## Call:
## madauni(x = sssroc)
## 
##    DOR  tau^2 
## 10.737  0.568
#forest(fit.DOR.DSL,snames=row.names(sssroc),cex=.8,main="Forest plot of Diagnostic Odds Ratio")
forest(fit.DOR.DSL,cex=.8,main="Forest plot of Diagnostic Odds Ratio")

#proportional hazard model
(fit.phm.het <- phm(sssroc)) #hterogeneity
## Call:
## phm.default(data = sssroc)
## 
## Coefficients:
##      theta    taus_sq 
## 0.24494920 0.01681934
summary(fit.phm.het)
## Call:
## phm.default(data = sssroc)
## 
##           Estimate         2.5 %     97.5 %
## theta   0.24494920  0.1824822303 0.30741617
## taus_sq 0.01681934 -0.0002781086 0.03391678
## 
## Log-likelihood: 38.72 on 2 degrees of freedom
## AIC:  -73.4 
## BIC:  -70.8 
## 
##  Chi-square goodness of fit test (Adjusted Profile Maximum
##  Likelihood under heterogeneity)
## 
## data:  x
## Chi-square = 25.746, df = 2, p-value = 0.4771
## 
## 
##    AUC  2.5 % 97.5 %   pAUC  2.5 % 97.5 % 
##  0.803  0.846  0.765  0.678  0.745  0.618
#plot portion of roc space
#plot(fit.phm.het, xlim = c(0,0.6), ylim = c(0.4,1))
#plot entire roc space
plot(fit.phm.het, xlim = c(0,1.0), ylim = c(0.0,1.0))
ROCellipse(sssroc, add = TRUE)

(fit.reitsma <- reitsma(sssroc))
## Call:  reitsma.default(data = sssroc)
## 
## Fixed-effects coefficients:
##               tsens     tfpr
## (Intercept)  0.5164  -1.8194
## 
## 28 studies, 2 fixed and 3 random-effects parameters
##   logLik       AIC       BIC  
##  42.9015  -75.8031  -65.6763
summary(fit.reitsma)
## Call:  reitsma.default(data = sssroc)
## 
## Bivariate diagnostic random-effects meta-analysis
## Estimation method: REML
## 
## Fixed-effects coefficients
##                   Estimate Std. Error       z Pr(>|z|) 95%ci.lb 95%ci.ub
## tsens.(Intercept)    0.516      0.177   2.920    0.003    0.170    0.863
## tfpr.(Intercept)    -1.819      0.128 -14.159    0.000   -2.071   -1.568
## sensitivity          0.626          -       -        -    0.542    0.703
## false pos. rate      0.140          -       -        -    0.112    0.173
##                      
## tsens.(Intercept)  **
## tfpr.(Intercept)  ***
## sensitivity          
## false pos. rate      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 
## 
## Variance components: between-studies Std. Dev and correlation matrix
##       Std. Dev tsens  tfpr
## tsens    0.813 1.000     .
## tfpr     0.592 0.466 1.000
## 
##  logLik     AIC     BIC 
##  42.902 -75.803 -65.676 
## 
## AUC:  0.846
## Partial AUC (restricted to observed FPRs and normalized):  0.717
## 
## HSROC parameters 
##       Theta      Lambda        beta sigma2theta sigma2alpha 
##      -0.846       2.573      -0.317       0.353       0.513
plot(fit.reitsma, sroclwd = 2,  main = "SROC curve (bivariate model) for Spot Sign data")
 points(fpr(sssroc), sens(sssroc), pch = 2)