if (FALSE) {
# We need to create this dataset
data(pbc, package = "randomForestSRC",)
# For whatever reason, the age variable is in days... makes no sense to me
for (ind in seq_len(dim(pbc)[2])) {
if (!is.factor(pbc[, ind])) {
if (length(unique(pbc[which(!is.na(pbc[, ind])), ind])) <= 2) {
if (sum(range(pbc[, ind], na.rm = TRUE) == c(0, 1)) == 2) {
pbc[, ind] <- as.logical(pbc[, ind])
}
}
} else {
if (length(unique(pbc[which(!is.na(pbc[, ind])), ind])) <= 2) {
if (sum(sort(unique(pbc[, ind])) == c(0, 1)) == 2) {
pbc[, ind] <- as.logical(pbc[, ind])
}
if (sum(sort(unique(pbc[, ind])) == c(FALSE, TRUE)) == 2) {
pbc[, ind] <- as.logical(pbc[, ind])
}
}
}
if (!is.logical(pbc[, ind]) &
length(unique(pbc[which(!is.na(pbc[, ind])), ind])) <= 5) {
pbc[, ind] <- factor(pbc[, ind])
}
}
#Convert age to years
pbc$age <- pbc$age / 364.24
pbc$years <- pbc$days / 364.24
pbc <- pbc[, -which(colnames(pbc) == "days")]
pbc$treatment <- as.numeric(pbc$treatment)
pbc$treatment[which(pbc$treatment == 1)] <- "DPCA"
pbc$treatment[which(pbc$treatment == 2)] <- "placebo"
pbc$treatment <- factor(pbc$treatment)
dta_train <- pbc[-which(is.na(pbc$treatment)), ]
# Create a test set from the remaining patients
pbc_test <- pbc[which(is.na(pbc$treatment)), ]
#========
# build the forest:
rfsrc_pbc <- randomForestSRC::rfsrc(
Surv(years, status) ~ .,
dta_train,
nsplit = 10,
na.action = "na.impute",
forest = TRUE,
importance = TRUE,
save.memory = TRUE
)
xvar <- c("bili", "albumin", "copper", "prothrombin", "age")
xvar_cat <- c("edema")
xvar <- c(xvar, xvar_cat)
time_index <- c(which(rfsrc_pbc$time.interest > 1)[1] - 1,
which(rfsrc_pbc$time.interest > 3)[1] - 1,
which(rfsrc_pbc$time.interest > 5)[1] - 1)
partial_pbc <- mclapply(rfsrc_pbc$time.interest[time_index],
function(tm) {
plot.variable(rfsrc_pbc, surv.type = "surv",
time = tm, xvar.names = xvar,
partial = TRUE,
show.plots = FALSE)
})
# A list of 2 plot.variable objects
length(partial_pbc)
class(partial_pbc)
class(partial_pbc[[1]])
class(partial_pbc[[2]])
# Create gg_partial objects
ggPrtl.1 <- gg_partial(partial_pbc[[1]])
ggPrtl.2 <- gg_partial(partial_pbc[[2]])
# Combine the objects to get multiple time curves
# along variables on a single figure.
ggpart <- combine.gg_partial(ggPrtl.1, ggPrtl.2,
lbls = c("1 year", "3 years"))
# Plot each figure separately
plot(ggpart)
# Get the continuous data for a panel of continuous plots.
ggcont <- ggpart
ggcont$edema <- ggcont$ascites <- ggcont$stage <- NULL
plot(ggcont, panel=TRUE)
# And the categorical for a panel of categorical plots.
nms <- colnames(sapply(ggcont, function(st) {st}))
for(ind in nms) {
ggpart[[ind]] <- NULL
}
plot(ggpart, panel=TRUE)
}
Run the code above in your browser using DataLab