13. MAIN FUNCTION
###=======================================================
#' RUN BAYESIAN NETWORKS (REFACTORED)
#' @param df_raw DATAFRAME
#' @param n_restarts NUMBER OF RANDOM RESTARTS
#' @param bl BLACKLIST
#' @param run_name FOLDER NAME OF OUTPUT
#' @param ref_color REFERENCE COLOUR (AIC)
#' @param comp_color COMPARISON COLOUR (MATCHES/MB)
#' RUN BAYESIAN NETWORKS (REFACTORED V3)
run_networks <- function(df_raw, n_restarts = 50, bl = NULL,
run_name = "Final_TP53_Analysis",
ref_color = "orange3",
comp_color = "orange") {
# --- 1. RUN 6 BNs -------------------------------
run_path <- file.path(output_dir, run_name)
dir.create(run_path, recursive = TRUE, showWarnings = FALSE)
cat("================================================\n")
cat("BN PIPELINE RUN: ", run_name, "\n")
cat("SAMPLE SIZE (N):", nrow(df_raw), "\n")
cat("================================================\n")
disc_levels <- c(2, 3)
scores <- c("aic", "bic", "bde")
all_results <- list()
ref_id <- "2disc.aic"
target_node <- "RISK"
master_nodes <- colnames(df_raw)
for (d in disc_levels) {
df_disc <- bnlearn::discretize(df_raw, method = 'quantile', breaks = d)
df_disc <- as.data.frame(lapply(df_disc, droplevels))[, master_nodes]
for (s in scores) {
id <- paste0(d, "disc.", s)
cat("\n>>> PROCESSING CONFIGURATION:", id, "\n")
set.seed(123)
starts <- c(list(empty.graph(master_nodes)),
random.graph(nodes = master_nodes, num = n_restarts - 1, method = "ic-dag", max.degree = 1))
net_list <- lapply(starts, function(g) {
tryCatch({ structural.em(df_disc, maximize = "hc", start = g,
maximize.args = list(score = s, blacklist = bl))
}, error = function(e) return(NULL))
})
net_list <- net_list[!sapply(net_list, is.null)]
avg_dag <- averaged.network(custom.strength(net_list, nodes = master_nodes))
print(avg_dag)
mb_nodes <- mb(avg_dag, target_node)
curr_amat <- amat(avg_dag)[master_nodes, master_nodes]
res_obj <- list(dag = avg_dag, adjacency = curr_amat, nodes = master_nodes, data = df_disc, mb = mb_nodes)
all_results[[id]] <- res_obj
saveRDS(res_obj, file.path(run_path, paste0(id, ".rds")))
assign(paste0(run_name, ".", id), res_obj, envir = .GlobalEnv)
# HIGH RES BN PLOTS
png(file.path(run_path, paste0(id, ".png")), width = 2400, height = 2100, res = 300, bg = "transparent", type = "cairo")
if (length(mb_nodes) > 0) {
graphviz.plot(avg_dag, highlight = list(nodes = mb_nodes, fill = comp_color, col = "black"),
main = paste(id, "| MB of", target_node))
} else {
graphviz.plot(avg_dag, main = paste(id, "| No MB for", target_node))
}
dev.off()
}
}
# --- 2. ADJACENCY MATRIX -------------------------------
# (Keeping your original adjacency grid logic as requested)
adj_mats <- lapply(all_results, function(x) x$adjacency)
summed_mat <- Reduce("+", adj_mats)
universal_mat <- (summed_mat == length(adj_mats))
ref_mat <- all_results[[ref_id]]$adjacency
plot_list <- list()
config_names <- names(all_results)
for (i in seq_along(config_names)) {
id <- config_names[i]
curr_mat <- all_results[[id]]$adjacency
plot_df <- as.data.frame(curr_mat) %>%
tibble::rownames_to_column("from") %>%
tidyr::pivot_longer(-from, names_to = "to", values_to = "exists") %>%
rowwise() %>%
mutate(
is_univ = universal_mat[from, to],
is_ref = ref_mat[from, to],
category = case_when(
exists == 1 & is_univ == TRUE ~ "Universal",
id == ref_id & exists == 1 & is_univ == FALSE ~ "Reference Only",
exists == 1 & is_ref == 1 & is_univ == FALSE ~ "Match Ref",
exists == 1 & is_ref == 0 & is_univ == FALSE ~ "New Edge",
TRUE ~ "None"
)
) %>% ungroup()
is_left_col <- i %in% c(1, 4)
is_bottom_row <- i %in% c(4, 5, 6)
p <- ggplot(plot_df, aes(x = factor(to, levels = master_nodes),
y = factor(from, levels = rev(master_nodes)))) +
geom_tile(aes(fill = category), color = "gray85", linewidth = 0.2) +
scale_fill_manual(
values = c("Universal" = "black", "Reference Only" = ref_color,
"Match Ref" = comp_color, "New Edge" = "lightgray", "None" = "white"),
guide = "none"
) +
labs(title = id, x = NULL, y = NULL) +
theme_minimal() +
theme(
aspect.ratio = 1,
panel.border = element_rect(color = "black", fill = NA, linewidth = 0.8),
panel.grid = element_blank(),
plot.title = element_text(hjust = 0.5, size = 10, face = "bold"),
axis.text.x = if(is_bottom_row) element_text(angle = 90, vjust = 0.5, hjust = 1, size = 7) else element_blank(),
axis.ticks.x = if(is_bottom_row) element_line(color = "black") else element_line(color = NA),
axis.text.y = if(is_left_col) element_text(size = 7) else element_blank(),
axis.ticks.y = if(is_left_col) element_line(color = "black") else element_line(color = NA),
axis.ticks.length = unit(3, "pt"),
plot.margin = margin(5, 5, 5, 5),
plot.background = element_rect(fill = "transparent", color = NA)
)
plot_list[[i]] <- p
}
final_grid <- patchwork::wrap_plots(plot_list, ncol = 3, nrow = 2) +
plot_annotation(
title = paste("Structural Stability Analysis:", run_name),
subtitle = "Row 1: 2-level Discretization | Row 2: 3-level Discretization",
theme = theme(plot.title = element_text(size = 14, face = "bold", hjust = 0.5),
plot.subtitle = element_text(size = 11, hjust = 0.5),
plot.background = element_rect(fill = "transparent", color = NA))
)
ggsave(file.path(run_path, "Adjacency_Comparison_Grid_Final.png"),
plot = final_grid, width = 11, height = 7.5, bg = "transparent", type = "cairo")
# --- 3. CPT (REFERENCE ONLY) -------------------------------
ref_res <- all_results[[ref_id]]
mb_nodes <- ref_res$mb
# Identify variables and order them (putting PC3 before Risk if present)
pc_in_mb <- grep("^PC[1-5]$", mb_nodes, value = TRUE)
other_mb <- setdiff(mb_nodes, pc_in_mb)
if ("PC3" %in% pc_in_mb) {
ordered_vars <- c(other_mb, setdiff(pc_in_mb, "PC3"), "PC3", target_node)
} else {
ordered_vars <- c(other_mb, pc_in_mb, target_node)
}
# 1. Generate Frequency Table (N) and Probability Table
raw_counts <- table(ref_res$data[, ordered_vars])
prob_table <- prop.table(raw_counts, margin = 1:(length(ordered_vars) - 1))
# 2. Convert to DataFrames
df_probs <- as.data.frame(prob_table)
df_counts <- as.data.frame(raw_counts)
# 3. Merge Probabilities and Sample Numbers (N)
final_cpt_output <- df_probs %>%
dplyr::rename(Probability = Freq) %>%
dplyr::mutate(N = df_counts$Freq)
# 4. Save CSV
write.csv(final_cpt_output, file.path(run_path, "Reference_Risk_CPT_with_N.csv"), row.names = FALSE)
# --- 4. MOSAIC PLOT (REFERENCE ONLY) -------------------------------
n_configs <- prod(dim(counts_table)[-length(dim(counts_table))])
p_high <- as.vector(cpt_table)[(n_configs + 1):(2 * n_configs)]
p_high[is.na(p_high)] <- 0
p_contrast <- 1 / (1 + exp(-10 * (p_high - 0.5)))
grad_pal <- colorRampPalette(c("#F5F5F5", ref_color))(100)
risk_colors <- grad_pal[round(p_contrast * 99) + 1]
color_array <- array(c(rep("#FFFFFF00", n_configs), risk_colors), dim = dim(counts_table))
# MOSAIC PLOT WITH THIN BLACK OUTLINES
png(file.path(run_path, "Reference_Mosaic_Risk_Final.png"), width = 3600, height = 2700, res = 300, bg = "transparent", type = "cairo")
mosaic(counts_table,
gp = gpar(fill = color_array, col = "black", lwd = 0.5), # col = "black" for thin black outlines
main = paste("Reference Run:", run_name, "Prognostic Risk Hierarchy"),
labeling = labeling_border(gp_labels = gpar(fontsize = 9, fontface = "bold"), rot_labels = c(0, 90, 0, 90)))
dev.off()
cat("\nDONE. All outputs saved to:", run_path, "\n")
return(all_results)
}
19. RFS
###=================================================================
## ============================================================
## 19. RANDOM FOREST SURVIVAL — FULL PIPELINE
## ============================================================
library(survival)
library(survminer)
library(randomForestSRC)
library(ggplot2)
library(patchwork)
library(dplyr)
library(scales)
## ============================================================
## 19.1 GET DATA
## ============================================================
# DEFINE VARIABLES
selection <- list(
isoform_pca = c("PC3", "PC5"),
mofa_factors = c("Factor2", "Factor3", "Factor8"),
stemness = "RNAss",
rnaseq = c("SOX10"),
immune = c("Troester_WoundSig_19887484"),
mutation = "TP53"
)
# DEFINE METADATA
meta_vars <- c("OS", "OS.time", "mofa_cluster")
# RUN get_data
rfs_data <- as.data.frame(get_data(mae, selection, meta_vars)) %>%
filter(mofa_cluster == "Cluster_2") %>%
select(-mofa_cluster)
## ── PRE-FILTER DATA (7-Year Window) ─────────────────────────────────────────
rfs_clean <- rfs_data %>%
filter(OS.time <= 2555) %>%
drop_na()
cat(sprintf("\nSample size after filtering: N = %d\n", nrow(rfs_clean)))
Sample size after filtering: N = 517
## ── CREATE OUTPUT DIRECTORY ──────────────────────────────────────────────────
survival_dir <- file.path(output_dir, "survival")
if (!dir.exists(survival_dir)) dir.create(survival_dir, recursive = TRUE)
cat(sprintf("Saving plots to: %s\n", survival_dir))
Saving plots to: C:/Users/alden/Dissertation/TP53_isoforms_project/output/survival
## ── Shared save helper ───────────────────────────────────────────────────────
# Single function so dpi/bg are consistent everywhere
save_png <- function(plot, filename, width = 10, height = 6) {
ggsave(
filename = file.path(survival_dir, filename),
plot = plot,
width = width,
height = height,
dpi = 300,
bg = "transparent"
)
cat(sprintf("Saved: %s\n", filename))
}
## ============================================================
## 19.2 DEFINE FEATURE SETS
## ============================================================
all_features <- colnames(rfs_clean)[!colnames(rfs_clean) %in% c("OS", "OS.time", "PC3", "PC5")]
features_full <- c("PC3", "PC5", all_features)
features_no_p53 <- all_features
## ============================================================
## 19.3 MONTE CARLO CROSS-VALIDATION (20 × 50:50 Split)
## ============================================================
all_model_results <- list()
all_km_data <- list()
for (model_type in c("Full", "No_p53")) {
current_features <- if (model_type == "Full") features_full else features_no_p53
rf_formula <- as.formula(
paste("Surv(OS.time, OS) ~", paste(current_features, collapse = "+"))
)
results_list <- list()
km_data_list <- list()
for (i in 1:20) {
set.seed(i)
# A. 50:50 random split
train_idx <- sample(seq_len(nrow(rfs_clean)), size = floor(0.50 * nrow(rfs_clean)))
train_set <- rfs_clean[train_idx, ]
test_set <- rfs_clean[-train_idx, ]
# B. Train RSF
rf_model <- randomForestSRC::rfsrc(rf_formula, data = train_set, ntree = 1000)
# C. Predict on test set
rf_pred <- randomForestSRC:::predict.rfsrc(rf_model, test_set)
test_set$Predicted_Risk <- rf_pred$predicted
# D. Optimal cutpoint
res_cut <- try(
survminer::surv_cutpoint(test_set,
time = "OS.time",
event = "OS",
variables = "Predicted_Risk"),
silent = TRUE
)
if (!inherits(res_cut, "try-error")) {
optimal_cut <- res_cut$cutpoint$cutpoint
test_set$Group <- factor(
ifelse(test_set$Predicted_Risk > optimal_cut, "High", "Low"),
levels = c("Low", "High")
)
# E. Discrete Cox model
cox_mod <- survival::coxph(Surv(OS.time, OS) ~ Group, data = test_set)
sum_cox <- summary(cox_mod)
results_list[[i]] <- data.frame(
Iteration = i,
Model = model_type,
P_Value = sum_cox$logtest["pvalue"],
HR = sum_cox$conf.int[1],
Lower_CI = sum_cox$conf.int[3],
Upper_CI = sum_cox$conf.int[4],
C_Index = 1 - rf_model$err.rate[rf_model$ntree]
)
# F. Store test set for KM
km_data_list[[i]] <- test_set %>%
select(OS.time, OS, Group) %>%
mutate(Iteration = i)
}
}
all_model_results[[model_type]] <- dplyr::bind_rows(results_list)
all_km_data[[model_type]] <- dplyr::bind_rows(km_data_list)
}
## ============================================================
## 19.4 ITERATION TABLE — ALL 20 FOLDS + SUMMARY
## ============================================================
iteration_table <- dplyr::bind_rows(all_model_results) %>%
filter(HR < 1000) %>%
select(Model, Iteration, HR, Lower_CI, Upper_CI, P_Value, C_Index) %>%
arrange(Model, Iteration)
cat("\n════════════════════════════════════════════════════════════════\n")
════════════════════════════════════════════════════════════════
cat(" PER-ITERATION RESULTS (all 20 folds, both models)\n")
PER-ITERATION RESULTS (all 20 folds, both models)
cat("════════════════════════════════════════════════════════════════\n")
════════════════════════════════════════════════════════════════
print(iteration_table, row.names = FALSE, digits = 3)
summary_table <- iteration_table %>%
group_by(Model) %>%
summarise(
N = n(),
Mean_HR = mean(HR, na.rm = TRUE),
SD_HR = sd(HR, na.rm = TRUE),
Mean_LCI = mean(Lower_CI, na.rm = TRUE),
SD_LCI = sd(Lower_CI, na.rm = TRUE),
Mean_UCI = mean(Upper_CI, na.rm = TRUE),
SD_UCI = sd(Upper_CI, na.rm = TRUE),
Mean_P = mean(P_Value, na.rm = TRUE),
SD_P = sd(P_Value, na.rm = TRUE),
Mean_CIndex = mean(C_Index, na.rm = TRUE),
SD_CIndex = sd(C_Index, na.rm = TRUE),
.groups = "drop"
)
cat("\n════════════════════════════════════════════════════════════════\n")
════════════════════════════════════════════════════════════════
cat(" SUMMARY: MEAN ± SD ACROSS 20 ITERATIONS\n")
SUMMARY: MEAN ± SD ACROSS 20 ITERATIONS
cat("════════════════════════════════════════════════════════════════\n")
════════════════════════════════════════════════════════════════
print(as.data.frame(summary_table), row.names = FALSE, digits = 3)
# ── Save tables as CSVs ───────────────────────────────────────────────────────
write.csv(iteration_table,
file.path(survival_dir, "iterations_all_folds.csv"),
row.names = FALSE)
cat("Saved: iterations_all_folds.csv\n")
Saved: iterations_all_folds.csv
write.csv(summary_table,
file.path(survival_dir, "iterations_summary.csv"),
row.names = FALSE)
cat("Saved: iterations_summary.csv\n")
Saved: iterations_summary.csv
## ============================================================
## 19.5 KAPLAN-MEIER: 20 ITERATION CURVES + MEAN 95% CI
## ============================================================
time_grid <- seq(0, 2555, by = 5)
# ── Helper: per-iteration step curves ────────────────────────────────────────
build_iter_curves <- function(km_data_model) {
lapply(split(km_data_model, km_data_model$Iteration), function(iter_df) {
lapply(c("Low", "High"), function(grp) {
sub <- iter_df[iter_df$Group == grp, ]
if (nrow(sub) < 2) return(NULL)
fit <- survfit(Surv(OS.time, OS) ~ 1, data = sub)
sf <- stepfun(fit$time, c(1, fit$surv))
data.frame(
time = time_grid,
surv = sf(time_grid),
Group = grp,
Iteration = unique(iter_df$Iteration)
)
}) %>% dplyr::bind_rows()
}) %>% dplyr::bind_rows()
}
# ── Helper: mean survival + 95% CI ───────────────────────────────────────────
build_mean_ci <- function(curves_df) {
curves_df %>%
group_by(Group, time) %>%
summarise(
mean_surv = mean(surv, na.rm = TRUE),
sd_surv = sd(surv, na.rm = TRUE),
n = sum(!is.na(surv)),
se_surv = sd_surv / sqrt(n),
ci_lo = pmax(mean_surv - 1.96 * se_surv, 0),
ci_hi = pmin(mean_surv + 1.96 * se_surv, 1),
.groups = "drop"
)
}
curves_full <- build_iter_curves(all_km_data[["Full"]])
curves_nop53 <- build_iter_curves(all_km_data[["No_p53"]])
mean_ci_full <- build_mean_ci(curves_full)
mean_ci_nop53 <- build_mean_ci(curves_nop53)
# ── Plot function ─────────────────────────────────────────────────────────────
plot_20_km <- function(curves_df, mean_ci_df, model_results_df, title_text) {
hr_summary <- model_results_df %>%
filter(HR < 1000) %>%
summarise(
mHR = mean(HR, na.rm = TRUE),
sdHR = sd(HR, na.rm = TRUE),
mLCI = mean(Lower_CI, na.rm = TRUE),
mUCI = mean(Upper_CI, na.rm = TRUE),
mP = mean(P_Value, na.rm = TRUE),
mCI = mean(C_Index, na.rm = TRUE)
)
subtitle_text <- sprintf(
"Mean HR = %.2f (SD %.2f) | Mean 95%% CI [%.2f–%.2f] | Mean p = %.3f | Mean C-index = %.3f",
hr_summary$mHR, hr_summary$sdHR,
hr_summary$mLCI, hr_summary$mUCI,
hr_summary$mP, hr_summary$mCI
)
ggplot() +
geom_step(data = curves_df,
aes(x = time, y = surv, colour = Group,
group = interaction(Group, Iteration)),
alpha = 0.50, linewidth = 0.45) +
geom_ribbon(data = mean_ci_df,
aes(x = time, ymin = ci_lo, ymax = ci_hi,
fill = Group, group = Group),
alpha = 0.20) +
geom_step(data = mean_ci_df,
aes(x = time, y = mean_surv, colour = Group, group = Group),
linewidth = 1.3) +
scale_colour_manual(
values = c(High = "#E84855", Low = "#2E86AB"),
labels = c(High = "High Risk", Low = "Low Risk")
) +
scale_fill_manual(
values = c(High = "#E84855", Low = "#2E86AB"),
labels = c(High = "High Risk", Low = "Low Risk")
) +
scale_x_continuous(
breaks = seq(0, 2555, by = 365),
labels = paste0(0:7, "y")
) +
scale_y_continuous(
limits = c(0, 1),
labels = scales::percent_format(accuracy = 1)
) +
labs(
title = title_text,
subtitle = subtitle_text,
x = "Time",
y = "Overall Survival Probability",
colour = "Risk Group",
fill = "Risk Group",
caption = "Faint lines = 20 Monte Carlo iterations (50:50 split) | Bold = mean | Ribbon = mean 95% CI"
) +
theme_classic(base_size = 13) +
theme(
panel.background = element_rect(fill = "transparent", colour = NA),
plot.background = element_rect(fill = "transparent", colour = NA),
legend.background = element_rect(fill = "transparent", colour = NA),
legend.key = element_rect(fill = "transparent", colour = NA),
plot.title = element_text(face = "bold", size = 14),
plot.subtitle = element_text(size = 9.5, colour = "grey35"),
plot.caption = element_text(size = 8, colour = "grey55"),
legend.position = "bottom",
legend.title = element_text(face = "bold")
)
}
p_full <- plot_20_km(curves_full, mean_ci_full,
all_model_results[["Full"]],
"RFS Model: Full (with PC3/PC5)")
p_nop53 <- plot_20_km(curves_nop53, mean_ci_nop53,
all_model_results[["No_p53"]],
"RFS Model: No PC3/PC5")
print(p_full / p_nop53)
save_png(p_full, "km_full_model.png", width = 10, height = 6)
Saved: km_full_model.png
save_png(p_nop53, "km_no_pc35_model.png", width = 10, height = 6)
Saved: km_no_pc35_model.png
save_png(p_full / p_nop53, "km_both_models.png", width = 10, height = 12)
Saved: km_both_models.png
## ============================================================
## 19.6 FINAL MODEL (Full) — VIMP + PARTIAL DEPENDENCE PLOTS
## ============================================================
cat("\nTraining final model on full dataset for VIMP and PDP...\n")
Training final model on full dataset for VIMP and PDP...
set.seed(42)
rf_formula_full <- as.formula(
paste("Surv(OS.time, OS) ~", paste(features_full, collapse = "+"))
)
final_model <- randomForestSRC::rfsrc(
rf_formula_full,
data = rfs_clean,
ntree = 1000,
importance = TRUE
)
# ── VIMP table ────────────────────────────────────────────────────────────────
vimp_df <- data.frame(
Feature = names(final_model$importance),
VIMP = as.numeric(final_model$importance)
) %>%
arrange(desc(VIMP))
cat("\n════════════════════════════════════════\n")
════════════════════════════════════════
cat(" VARIABLE IMPORTANCE (VIMP) — Full Model\n")
VARIABLE IMPORTANCE (VIMP) — Full Model
cat("════════════════════════════════════════\n")
════════════════════════════════════════
print(vimp_df, row.names = FALSE, digits = 4)
write.csv(vimp_df,
file.path(survival_dir, "vimp_full_model.csv"),
row.names = FALSE)
cat("Saved: vimp_full_model.csv\n")
Saved: vimp_full_model.csv
# ── VIMP plot ─────────────────────────────────────────────────────────────────
p_vimp <- ggplot(vimp_df,
aes(x = reorder(Feature, VIMP), y = VIMP, fill = VIMP > 0)) +
geom_col(width = 0.7) +
geom_hline(yintercept = 0, linetype = "dashed", colour = "grey50") +
scale_fill_manual(
values = c(`TRUE` = "#2E86AB", `FALSE` = "#E84855"),
guide = "none"
) +
coord_flip() +
labs(
title = "Variable Importance (VIMP) — Full RFS Model",
x = NULL,
y = "VIMP",
caption = "Blue = positive importance | Red = feature may add noise"
) +
theme_classic(base_size = 13) +
theme(
panel.background = element_rect(fill = "transparent", colour = NA),
plot.background = element_rect(fill = "transparent", colour = NA),
legend.background = element_rect(fill = "transparent", colour = NA),
legend.key = element_rect(fill = "transparent", colour = NA),
plot.title = element_text(face = "bold")
)
print(p_vimp)
save_png(p_vimp, "vimp_full_model.png", width = 8, height = 6)
Saved: vimp_full_model.png
# ── Partial dependence plots ──────────────────────────────────────────────────
# plot.variable() is base R graphics; png() / dev.off() is the correct
# capture method. bg = "transparent" sets the device background.
# Note: the plot panels themselves will retain a white fill as this is
# controlled internally by randomForestSRC and cannot be overridden
# without reimplementing PDPs in ggplot2.
png(file.path(survival_dir, "pdp_full_model.png"),
width = 10, height = 8,
units = "in",
res = 300,
bg = "transparent")
randomForestSRC::plot.variable(
final_model,
partial = TRUE,
sorted = TRUE,
plots.per.page = 3,
main = "Partial Dependence — Full RFS Model (with PC3/PC5)"
)
dev.off()
png
2
cat("Saved: pdp_full_model.png\n")
Saved: pdp_full_model.png
# Also render to screen
randomForestSRC::plot.variable(
final_model,
partial = TRUE,
sorted = TRUE,
plots.per.page = 3,
main = "Partial Dependence — Full RFS Model (with PC3/PC5)"
)
## ── Final summary ─────────────────────────────────────────────────────────────
cat("\n════════════════════════════════════════════════════════════════\n")
════════════════════════════════════════════════════════════════
cat(" ALL FILES SAVED TO:", survival_dir, "\n")
ALL FILES SAVED TO: C:/Users/alden/Dissertation/TP53_isoforms_project/output/survival
cat("════════════════════════════════════════════════════════════════\n")
════════════════════════════════════════════════════════════════
cat(" iterations_all_folds.csv\n")
iterations_all_folds.csv
cat(" iterations_summary.csv\n")
iterations_summary.csv
cat(" km_full_model.png\n")
km_full_model.png
cat(" km_no_pc35_model.png\n")
km_no_pc35_model.png
cat(" km_both_models.png\n")
km_both_models.png
cat(" vimp_full_model.csv\n")
vimp_full_model.csv
cat(" vimp_full_model.png\n")
vimp_full_model.png
cat(" pdp_full_model.png\n")
pdp_full_model.png
cat("════════════════════════════════════════════════════════════════\n")
════════════════════════════════════════════════════════════════
