About the rocbc Package

(Authors: LE Bantis, B Brewer, CT Nakas, B Reiser)

This package is focused on Box-Cox based ROC curves and provides point estimates, confidence intervals (CIs), and hypothesis tests. It can be used both for inferences for a single biomarker and when comparisons of two correlated biomarkers are of interest. It provides inferences and comparisons around the area under the ROC curve (AUC), the Youden index, the sensitivity at a given specificity level (and vice versa), the optimal operating point of the ROC curve (in the Youden sense), and the Youden based cutoff. This documentation consists of two parts, the one marker and the two marker case. All approaches presented herein have been recently published (see references for each function).

The functions of each part:


checkboxcox


#DATA GENERATION
set.seed(123)
x=rgamma(100, shape=2, rate = 8) # Generates biomarker data from a gamma distribution 
                                 # for the healthy group.
y=rgamma(100, shape=2, rate = 4) # Generates biomarker data from a gamma distribution 
                                 # for the diseased group.
scores=c(x,y)
D=c(zeros(1,100), ones(1,100))

out=checkboxcox(marker=scores, D, plots="on", printShapiro = TRUE)
## 
##  Shapiro-Wilk normality test
## 
## data:  x
## W = 0.92354, p-value = 2.181e-05
## 
## 
##  Shapiro-Wilk normality test
## 
## data:  y
## W = 0.90169, p-value = 1.7e-06
## 
## 
##  Shapiro-Wilk normality test
## 
## data:  transx
## W = 0.98765, p-value = 0.4826
## 
## 
##  Shapiro-Wilk normality test
## 
## data:  transy
## W = 0.98936, p-value = 0.6128

plot of chunk unnamed-chunk-11plot of chunk unnamed-chunk-11plot of chunk unnamed-chunk-11

summary(out)
##                          Length Class  Mode    
## transformation.parameter   1    -none- numeric 
## transx                   100    -none- numeric 
## transy                   100    -none- numeric 
## pval_x                     1    -none- numeric 
## pval_y                     1    -none- numeric 
## pval_transx                1    -none- numeric 
## pval_transy                1    -none- numeric 
## rocfun                     1    -none- function

threerocs


#DATA GENERATION
set.seed(123)
x <- rgamma(100, shape=2, rate = 8) # generates biomarker data from a gamma
                                 # distribution for the healthy group.
y <- rgamma(100, shape=2, rate = 4) # generates biomarker data from a gamma
                                 # distribution for the diseased group.
scores <- c(x,y)
D=c(pracma::zeros(1,100), pracma::ones(1,100))
out=threerocs(marker=scores, D, plots="on")
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases

plot of chunk unnamed-chunk-13

summary(out)
##               Length Class  Mode   
## AUC_Empirical 1      auc    numeric
## AUC_Metz      1      -none- numeric
## AUC_BoxCox    1      -none- numeric
summary(out)
##               Length Class  Mode   
## AUC_Empirical 1      auc    numeric
## AUC_Metz      1      -none- numeric
## AUC_BoxCox    1      -none- numeric

rocboxcox


set.seed(123)
x=rgamma(100, shape=2, rate = 8)
y=rgamma(100, shape=2, rate = 4)
scores=c(x,y)
D=c(zeros(1,100), ones(1,100))
out=rocboxcox(marker=scores,D, 0.05, plots="on", printProgress=FALSE)

plot of chunk unnamed-chunk-15


rocboxcoxCI


out=rocboxcoxCI(marker, D, givenSP=NA, givenSE=givenSE, alpha, plots)

givenSP=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9)
givenSE=NA
out=rocboxcoxCI(marker=scores ,D, givenSP=givenSP, givenSE=NA, alpha=0.05, plots="on")

plot of chunk unnamed-chunk-16plot of chunk unnamed-chunk-16

#

checkboxcox2


#DATA GENERATION
set.seed(123)

nx <- 100
Sx <- matrix(c(1, 0.5, 0.5, 1), nrow = 2, ncol = 2)

mux <- c(X = 10, Y = 12)
X = mvtnorm::rmvnorm(nx, mean = mux, sigma = Sx)

ny <- 100
Sy <- matrix(c(1.1, 0.6, 0.6, 1.1), nrow = 2, ncol = 2)

muy <- c(X = 11, Y = 13.7)
Y = mvtnorm::rmvnorm(ny, mean = muy, sigma = Sy)

dx = pracma::zeros(nx,1)
dy = pracma::ones(ny,1)

markers = rbind(X,Y);
marker1 = markers[,1]
marker2 = markers[,2]
D = c(rbind(dx,dy))

out=checkboxcox2(marker1, marker2, D, plots = "on")

plot of chunk unnamed-chunk-17plot of chunk unnamed-chunk-17plot of chunk unnamed-chunk-17plot of chunk unnamed-chunk-17plot of chunk unnamed-chunk-17

summary(out)
##                            Length Class  Mode    
## res_shapiro                  8    -none- list    
## transformation.parameter.1   1    -none- numeric 
## transx1                    100    -none- numeric 
## transy1                    100    -none- numeric 
## transformation.parameter.2   1    -none- numeric 
## transx2                    100    -none- numeric 
## transy2                    100    -none- numeric 
## pval_x1                      1    -none- numeric 
## pval_y1                      1    -none- numeric 
## pval_transx1                 1    -none- numeric 
## pval_transy1                 1    -none- numeric 
## pval_x2                      1    -none- numeric 
## pval_y2                      1    -none- numeric 
## pval_transx2                 1    -none- numeric 
## pval_transy2                 1    -none- numeric 
## roc1                         1    -none- function
## roc2                         1    -none- function

comparebcAUC


#GENERATE SOME BIVARIATE DATA===
set.seed(123)

nx <- 100
 Sx <- matrix(c(1,   0.5,
                0.5,  1), 
            nrow = 2, ncol = 2)

mux <- c(X = 10, Y = 12)
X=rmvnorm(nx, mean = mux, sigma = Sx)

ny <- 100
Sy <- matrix(c(1.1,   0.6,
               0.6,  1.1), 
             nrow = 2, ncol = 2)

muy <- c(X = 11, Y = 13.7)
Y=rmvnorm(ny, mean = muy, sigma = Sy)

dx=zeros(nx,1)
dy=ones(ny,1)

markers=rbind(X,Y);
marker1=markers[,1]
marker2=markers[,2]
D=c(rbind(dx,dy))

#==DATA HAVE BEEN GENERATED====


#===COMPARE THE AUCs of Marker 1 vs Maker 2

out=comparebcAUC(marker1, marker2, D, alpha=0.05,  plots="on")

plot of chunk unnamed-chunk-19


comparebcJ


#==DATA HAVE BEEN GENERATED====
#===COMPARE THE Js of Marker 1 vs Maker 2

out=comparebcJ(marker1, marker2, D, alpha=0.05,  plots="on")

plot of chunk unnamed-chunk-20

#

comparebcSens


out=comparebcSens(marker1=marker1, marker2=marker2, D=D, alpha =0.05, atSpec=0.8, plots="on")

plot of chunk unnamed-chunk-21

summary(out)
##                            Length Class       Mode    
## resultstable                 7    formattable numeric 
## Sens1                        1    -none-      numeric 
## Sens2                        1    -none-      numeric 
## pvalue_probit_difference     1    -none-      numeric 
## pvalue_difference            1    -none-      numeric 
## CI_difference                2    -none-      numeric 
## roc1                         1    -none-      function
## roc2                         1    -none-      function
## transx1                    100    -none-      numeric 
## transy1                    100    -none-      numeric 
## transformation.parameter.1   1    -none-      numeric 
## transx2                    100    -none-      numeric 
## transy2                    100    -none-      numeric 
## transformation.parameter.2   1    -none-      numeric
#  


comparebcSpec


out=comparebcSpec(marker1=marker1, marker2=marker2, D=D, alpha =0.05, atSens=0.8, plots="on")

plot of chunk unnamed-chunk-22

summary(out)
##                            Length Class       Mode    
## resultstable                 7    formattable numeric 
## FPR1                         1    -none-      numeric 
## FPR2                         1    -none-      numeric 
## pvalue_probit_difference     1    -none-      numeric 
## pvalue_difference            1    -none-      numeric 
## CI_difference                2    -none-      numeric 
## roc1                         1    -none-      function
## roc2                         1    -none-      function
## transx1                    100    -none-      numeric 
## transy1                    100    -none-      numeric 
## transformation.parameter.1   1    -none-      numeric 
## transx2                    100    -none-      numeric 
## transy2                    100    -none-      numeric 
## transformation.parameter.2   1    -none-      numeric
#