-
Notifications
You must be signed in to change notification settings - Fork 1
/
elasticnet.Rmd
131 lines (109 loc) · 2.92 KB
/
elasticnet.Rmd
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
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
---
title: "elasticnet"
output:
bookdown::html_document2:
fig.width: 18
fig.height: 8
---
```{r elasticnet_setup, message = FALSE}
library("targets")
library("ameld")
library("viridisLite")
tar_load(arcvob)
tar_load(bootrcv)
tar_load(bootrcvcc)
tar_load(bootarcv)
tar_load(bootrcv.woIC)
tar_load(amelddata)
tar_load(amelddatacc)
tar_load(ameldcfg)
tar_load(zlog_data)
tar_load(zlog_data_complete_cases)
```
# Tuning alpha
```{r arcvob}
arcvob
```
```{r arcvob-plot-alpha}
plot(arcvob)
```
```{r arcvob-plot-lambdamin}
plot(arcvob, what = "lambda.min")
```
```{r arcvob-plot-lambda1se}
plot(arcvob, what = "lambda.1se")
```
# Bootstrapping
## rcv.glmnet
```{r bootrcv}
plot(bootrcv, what = "calibration")
ps <- lapply(
zlog_data[paste0("SurvProbMeld", c("Unos", "NaUnos", "Plus7"))],
function(p) {
ctpnts <- cutpoints(p, n = ameldcfg$m)
f <- cut(p, ctpnts, include.lowest = TRUE)
list(
predicted = groupmean(p, f = f),
observed = observed_survival(
amelddata$y, f = f, times = ameldcfg$times
)
)
}
)
names(ps) <- c("MELD", "MELD-Na", "MELD-Plus7")
col <- viridisLite::viridis(6)[4:6]
for (i in seq_along(ps)) {
lines(
ps[[i]]$predicted, ps[[i]]$observed, col = col[i], type = "b", pch = 19
)
}
legend("topleft", col = col, legend = names(ps), pch = 19, bty = "n")
plot(bootrcv, what = "selected", cex = 0.5)
plot(bootrcv$fit$glmnet.fit, xvar = "norm")
plot(bootrcv$fit$glmnet.fit, xvar = "lambda")
plot(bootrcv$fit$glmnet.fit, xvar = "dev")
```
## rcv.glmnet (complete cases only)
```{r bootrcvcc}
plot(bootrcvcc, what = "calibration")
ps <- lapply(
zlog_data_complete_cases[paste0("SurvProbMeld", c("Unos", "NaUnos", "Plus7"))],
function(p) {
ctpnts <- cutpoints(p, n = ameldcfg$m)
f <- cut(p, ctpnts, include.lowest = TRUE)
list(
predicted = groupmean(p, f = f),
observed = observed_survival(
amelddatacc$y, f = f, times = ameldcfg$times
)
)
}
)
names(ps) <- c("MELD", "MELD-Na", "MELD-Plus7")
col <- viridisLite::viridis(6)[4:6]
for (i in seq_along(ps)) {
lines(
ps[[i]]$predicted, ps[[i]]$observed, col = col[i], type = "b", pch = 19
)
}
legend("topleft", col = col, legend = names(ps), pch = 19, bty = "n")
plot(bootrcvcc, what = "selected", cex = 0.5)
plot(bootrcvcc$fit$glmnet.fit, xvar = "norm")
plot(bootrcvcc$fit$glmnet.fit, xvar = "lambda")
plot(bootrcvcc$fit$glmnet.fit, xvar = "dev")
```
## bootrcv.woIC
Exclude IL-6 and CYSC columns from the data set.
```{r bootrcv.woIC}
plot(bootrcv.woIC, what = "calibration")
plot(bootrcv.woIC, what = "selected", cex = 0.5)
```
## arcv.glmnet
```{r bootarcv}
a <- c(table(sapply(bootarcv$models, function(m)m$fit$alpha)))
plot(bootarcv, what = "calibration")
```
```{r bootarcv2, fig.height = 12}
plot(bootarcv, what = "selected")
plot_dots(a, main = "Selected Alpha Values")
```