-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathclinical_utility_table_UNITED_pediatric.R
121 lines (81 loc) · 3.92 KB
/
clinical_utility_table_UNITED_pediatric.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
##Clinical utility table UNITED -----------------------------------------------------------------
## Extended table --------------------------------------------------------------------------
## All genes ---------------------------------------------------------------------------------------
dataset.UNITED_type1p <- create_data(dataset = "united t1d pediatrics", commonmody = FALSE)
## if MODY testing missing, change to 0
mutate(M = ifelse(is.na(M), 0, M))
#load predictions
#load T1D
predictions_dataset.UNITED_type1_young_all_genes_with_T <- readRDS("model_predictions/predictions_dataset.UNITED_type1_young_all_genes_with_T.rds")
UNITED_type1p <- cbind(dataset.UNITED_type1p, predictions_dataset.UNITED_type1_young_all_genes_with_T)
#Load function to get diagnostic info at different thresholds
calculate_thresholds_diagnostics <- function(response, prediction, unique = FALSE) {
### threshold number to test
thresholds <- seq(0, 1, 0.0001)
matrix_thresholds <- matrix(0, nrow = length(thresholds), ncol = 8)
matrix_thresholds[,1] <- thresholds
preds1 <- data.frame(
prediction = prediction,
response = response
)
for (i in 1:length(thresholds)) {
# original dataset
original_data <- preds1
# checking which patients are above the threshold
interim_above <- preds1 %>%
mutate(threshold_up = ifelse(prediction > thresholds[i], 1, 0)) %>%
filter(threshold_up == 1)
interim_below <- preds1 %>%
mutate(threshold_up = ifelse(prediction > thresholds[i], 1, 0)) %>%
filter(threshold_up == 0)
# calculating the sensitivity of the MODY cases
matrix_thresholds[i, 2] <- length(which(interim_above$response == 1))/length(which(original_data$response == 1)) * 100
# calculating number to test
if (length(which(interim_above$response == 1)) > 0) {
matrix_thresholds[i, 3] <- ceiling(nrow(interim_above)/length(which(interim_above$response == 1)))
} else {matrix_thresholds[i, 3] <- as.numeric(NA)}
# Pick-up rate
matrix_thresholds[i, 4] <- length(which(interim_above$response == 1))/nrow(interim_above) * 100
# calculating specificity
matrix_thresholds[i, 5] <- length(which(interim_below$response == 0))/length(which(original_data$response == 0)) * 100
# Perc of patients tested
matrix_thresholds[i, 6] <- (nrow(interim_above)/nrow(original_data))*100
# number of MODY cases missed
matrix_thresholds[i, 7] <- interim_below %>%
select(response) %>%
sum(na.rm = TRUE) %>%
unlist()
# number of patients tested
matrix_thresholds[i, 8] <- interim_above %>%
nrow() %>%
unlist()
}
## making the matrix a data.frame()
matrix_thresholds <- matrix_thresholds %>%
as.data.frame() %>%
set_names(c("Thresholds", "Sensitivity", "NTT", "Pick-up rate", "Specificity", "N-Tested", "Cases Missed", "Patients Tested")) %>%
arrange(desc(Thresholds))
if (unique == TRUE) {
## select unique combinations of sensitivity and NTT (only the first occurance)
matrix_thresholds <- matrix_thresholds %>%
slice(which(duplicated(matrix_thresholds %>% select(-Thresholds, -`Pick-up rate`)) == FALSE))
}
return(matrix_thresholds)
}
## Early-insulin-treated -----------------------------------------------------------------
thresholds_UNITED_t1d <- calculate_thresholds_diagnostics(dataset.UNITED_type1p$M, predictions_dataset.UNITED_type1_young_all_genes_with_T$prob)
### 5%
thresholds_UNITED_t1d %>%
filter(Thresholds == 0.05) %>% arrange(Thresholds) %>% head()
### 10%
thresholds_UNITED_t1d %>%
filter(Thresholds == 0.1) %>% arrange(Thresholds) %>% head()
### 20%
thresholds_UNITED_t1d %>%
filter(Thresholds == 0.2) %>% arrange(Thresholds) %>% head()
### 25%
thresholds_UNITED_t1d %>%
filter(Thresholds == 0.25) %>% arrange(Thresholds) %>% head()
### 30%
thresholds_UNITED_t1d %>%
filter(Thresholds == 0.3) %>% arrange(Thresholds) %>% head()