the great cleanup

tagteam · Jun 5, 2024 · 8945c1d · 8945c1d
1 parent b6ed48f
commit 8945c1d
Show file tree

Hide file tree

Showing 19 changed files with 615 additions and 581 deletions.
diff --git a/R/AUC.binary.R b/R/AUC.binary.R
@@ -1,11 +1,11 @@
 ### AUC.binary.R --- 
 #----------------------------------------------------------------------
-## Author: Thomas Alexander Gerds
+## Author: Thomas Alexander Gerds and Johan Sebastian Ohlendorff
 ## Created: Jan 11 2022 (17:04) 
 ## Version: 
-## Last-Updated: Jun  5 2024 (07:25) 
+## Last-Updated: Jun  5 2024 (14:50) 
 ##           By: Thomas Alexander Gerds
-##     Update #: 17
+##     Update #: 19
 #----------------------------------------------------------------------
 ## 
 ### Commentary: 
@@ -173,6 +173,7 @@ AUC.binary <- function(DT,
             res.cut <- list()
             for (i in 1:length(cutpoints)){
                 temp.TPR <- subset(temp.TPR.ic,cutpoints==cutpoints[i])
+                # FIXME: merge on what?
                 aucDT.temp <- merge(aucDT,temp.TPR)
                 some.fun <- function(riskRegression_event,risk,TPR,FPR,PPV,NPV,Prisks,Prisks2,cut,N){
                     meanY <- mean(riskRegression_event)

diff --git a/R/AUC.competing.risks.R b/R/AUC.competing.risks.R
@@ -3,9 +3,9 @@
 ## Author: Thomas Alexander Gerds
 ## Created: Jan 11 2022 (17:06)
 ## Version:
-## Last-Updated: Jun  5 2024 (07:24) 
+## Last-Updated: Jun  5 2024 (17:52) 
 ##           By: Thomas Alexander Gerds
-##     Update #: 37
+##     Update #: 43
 #----------------------------------------------------------------------
 ##
 ### Commentary:
@@ -81,7 +81,7 @@ AUC.competing.risks <- function(DT,
             den_FPR<-sum(ipcwControls1+ipcwControls2)
             indeces <- sindex(risk,cutpoints,comp = "greater",TRUE)
             res <- list()
-            # FIXME
+            # FIXME: the order is according to -risk but why? is that necessary?
             ordered <- order(riskRegression_time) ## can probably move this outside to improve computation time, for now keep it
             for (i in 1:length(cutpoints)){
                 den_PPV <- sum(ipcwCases[risk > cutpoints[i]]+ipcwControls1[risk > cutpoints[i]] + ipcwControls2[risk > cutpoints[i]])
@@ -210,8 +210,8 @@ AUC.competing.risks <- function(DT,
         sum((FP-c(0,FP[-N]))*((c(0,TP[-N])+TP)/2))
     }
     score <- aucDT[nodups,list(AUC=AireTrap(FPR,TPR)),by=list(model,times)]
-    data.table::setkey(score,model,times)
-    aucDT <- merge(score,aucDT,all=TRUE)
+    aucDT <- merge(score,aucDT,by = c("model","times"),all=TRUE)
+    data.table::setkey(aucDT,model,times)
     if (se.fit[[1]]==1L){
         aucDT[,nth.times:=as.numeric(factor(times))]
 
@@ -230,16 +230,16 @@ AUC.competing.risks <- function(DT,
                                              conservative = conservative[[1]],
                                              cens.model = cens.model), by=list(model,times)]
         se.score <- aucDT[,list(se=sd(IF.AUC)/sqrt(N)),by=list(model,times)]
-        data.table::setkey(se.score,model,times)
-        score <- score[se.score]
+        score <- score[se.score,,on = c("model","times")]
+        data.table::setkey(score,model,times)
         if (se.fit==1L){
             score[,lower:=pmax(0,AUC-qnorm(1-alpha/2)*se)]
             score[,upper:=pmin(1,AUC+qnorm(1-alpha/2)*se)]
         }else{
             score[,se:=NULL]
         }
+        aucDT <- aucDT[score, ,on = c("model","times")]
         data.table::setkey(aucDT,model,times)
-        aucDT <- aucDT[score]
         if (keep.vcov[[1]] == TRUE){
             output <- c(output,list(vcov=getVcov(aucDT,"IF.AUC",times=TRUE)))
         }

diff --git a/R/AUC.survival.R b/R/AUC.survival.R
@@ -3,9 +3,9 @@
 ## Author: Thomas Alexander Gerds
 ## Created: Jan 11 2022 (17:06)
 ## Version:
-## Last-Updated: Jun  5 2024 (07:23) 
+## Last-Updated: Jun  5 2024 (17:51) 
 ##           By: Thomas Alexander Gerds
-##     Update #: 30
+##     Update #: 45
 #----------------------------------------------------------------------
 ##
 ### Commentary:
@@ -79,7 +79,7 @@ AUC.survival <- function(DT,
             den_FPR<-sum(ipcwControls)
             indeces <- sindex(risk,cutpoints,comp = "greater",TRUE)
             res <- list()
-            ## browser()
+            # FIXME
             ordered <- order(riskRegression_time) ## can probably move this outside to improve computation time, for now keep it
             SE.TPR <- SE.FPR <- SE.PPV <- SE.NPV <- NA
             for (i in 1:length(cutpoints)){
@@ -187,9 +187,8 @@ AUC.survival <- function(DT,
     }else{
         output <- NULL
     }
-
-    data.table::setkey(score,model,times)
-    aucDT <- merge(score,aucDT,all=TRUE)
+    aucDT <- merge(score,aucDT,by = c("model","times"),all=TRUE)
+    data.table::setkey(aucDT,model,times)
     if (se.fit[[1]]==1L){
         aucDT[,nth.times:=as.numeric(factor(times))]
         ## compute influence function
@@ -207,16 +206,17 @@ AUC.survival <- function(DT,
                                              conservative = conservative[[1]],
                                              cens.model = cens.model), by=list(model,times)]
         se.score <- aucDT[,list(se=sd(IF.AUC)/sqrt(N)),by=list(model,times)]
-        data.table::setkey(se.score,model,times)
-        score <- score[se.score]
+        score <- score[se.score,,on = c("model","times")]
+        data.table::setkey(score,model,times)
         if (se.fit==1L){
             score[,lower:=pmax(0,AUC-qnorm(1-alpha/2)*se)]
             score[,upper:=pmin(1,AUC+qnorm(1-alpha/2)*se)]
         }else{
             score[,se:=NULL]
         }
+        # join with auc and (se, lower, upper) if se.fit
+        aucDT <- score[aucDT,,on = c("model","times")]
         data.table::setkey(aucDT,model,times)
-        aucDT <- aucDT[score]
         if (keep.vcov[[1]] == TRUE){
             output <- c(output,list(vcov=getVcov(aucDT,"IF.AUC",times=TRUE)))
         }

diff --git a/R/Brier.binary.R b/R/Brier.binary.R
@@ -3,9 +3,9 @@
 ## Author: Thomas Alexander Gerds
 ## Created: Jan 11 2022 (17:03) 
 ## Version: 
-## Last-Updated: Jun  5 2024 (07:25) 
+## Last-Updated: Jun  5 2024 (14:52) 
 ##           By: Thomas Alexander Gerds
-##     Update #: 8
+##     Update #: 9
 #----------------------------------------------------------------------
 ## 
 ### Commentary: 
@@ -51,8 +51,8 @@ Brier.binary <- function(DT,
     if (length(dolist)>0){
         ## merge with Brier score
         data.table::setkey(DT,model)
-        data.table::setkey(score,model)
-        DT <- DT[score]
+        DT <- DT[score,,on = c("model")]
+        data.table::setkey(DT,model)
         if (se.fit[[1]]==TRUE){
             contrasts.Brier <- DT[,getComparisons(data.table(x=Brier,
                                                              IF=residuals,

diff --git a/R/Brier.competing.risks.R b/R/Brier.competing.risks.R
@@ -3,9 +3,9 @@
 ## Author: Thomas Alexander Gerds
 ## Created: Jan 11 2022 (17:04)
 ## Version:
-## Last-Updated: Jun  5 2024 (07:24) 
+## Last-Updated: Jun  5 2024 (14:43) 
 ##           By: Thomas Alexander Gerds
-##     Update #: 9
+##     Update #: 10
 #----------------------------------------------------------------------
 ##
 ### Commentary:
@@ -68,10 +68,9 @@ Brier.competing.risks <- function(DT,
     }
     data.table::setkey(score,model,times)
     if (length(dolist)>0){
-        data.table::setkey(DT,model,times)
         ## merge with Brier score
-        DT <- DT[score]
-        data.table::setkey(score,model,times)
+        DT <- DT[score, ,on = c("model","times")]
+        data.table::setkey(DT,model,times)
         if (se.fit[[1]]==TRUE){
             contrasts.Brier <- DT[,getComparisons(data.table(x=Brier,IF=IF.Brier,model=model),
                                                   NF=NF,

diff --git a/R/Brier.survival.R b/R/Brier.survival.R
@@ -3,9 +3,9 @@
 ## Author: Thomas Alexander Gerds
 ## Created: Jan 11 2022 (17:04)
 ## Version:
-## Last-Updated: Jun  5 2024 (07:24) 
+## Last-Updated: Jun  5 2024 (14:48) 
 ##           By: Thomas Alexander Gerds
-##     Update #: 9
+##     Update #: 19
 #----------------------------------------------------------------------
 ##
 ### Commentary:
@@ -32,8 +32,8 @@ Brier.survival <- function(DT,
                            ...){
     IC0=IPCW=nth.times=riskRegression_ID = riskRegression_time = riskRegression_status=times=raw.Residuals=risk=Brier=residuals=WTi=Wt=status=setorder=model=IF.Brier=data.table=sd=lower=qnorm=se=upper=NULL
     ## compute 0/1 outcome:
-    DT[riskRegression_time<=times & status==1,residuals:=(1-risk)^2/WTi]
-    DT[riskRegression_time<=times & status==0,residuals:=0]
+    DT[riskRegression_time<=times & riskRegression_status==1,residuals:=(1-risk)^2/WTi]
+    DT[riskRegression_time<=times & riskRegression_status==0,residuals:=0]
     DT[riskRegression_time>times,residuals:=(risk)^2/Wt]
 
     if (se.fit[[1]]==1L){
@@ -61,18 +61,17 @@ Brier.survival <- function(DT,
     data.table::setkey(score,model,times)
     if (length(dolist)>0L){
         ## merge with Brier score
+        DT <- DT[score,,on = c("model","times")]
         data.table::setkey(DT,model,times)
-        ## data.table::setkey(score,model,times)
-        DT <- DT[score]
         if (se.fit[[1]]==TRUE){
-            contrasts.Brier <- DT[,getComparisons(data.table(x=Brier,IF=IF.Brier,model=model),
+            contrasts.Brier <- DT[,getComparisons(dt = data.table(x=Brier,IF=IF.Brier,model=model),
                                                   NF=NF,
                                                   N=N,
                                                   alpha=alpha,
                                                   dolist=dolist,
-                                                  se.fit=se.fit),by=list(times)]
+                                                  se.fit=TRUE),by=list(times)]
         }else{
-            contrasts.Brier <- DT[,getComparisons(data.table(x=Brier,model=model),
+            contrasts.Brier <- DT[,getComparisons(dt = data.table(x=Brier,model=model),
                                                   NF=NF,
                                                   N=N,
                                                   alpha=alpha,