①獲取predict的結果，“prob.Case”這一列就是預測風險概率，“truth”列為實際發生結局的分組

②將prob.Case進行分桶（簡單理解為分組，一般分10組)，常見的分桶方式有兩種：一是將prob.Case從大到小排序后，按照樣本數平均分為10組，每組樣本數相等

③計算10個分桶中prob.Case的桶內均值作為預測概率；

④計算10個分桶中實際患病者（truth = 1 | truth=Case）占該桶樣本數的頻率作為實際概率；

⑤將10對預測概率和實際概率分別作為橫坐標值和縱坐標值得到10個散點；

⑥將這些點連起來，即為校準曲線中的Apparent線。

R包的函數

function (df, outcome, positive, prediction, model, n_bins = 10, show_loess = FALSE, plot_title = "", ...) 
{if ((n_bins > 0 && show_loess == TRUE) || (n_bins == 0 && show_loess == FALSE)) {stop("You must either set n_bins > 0 and show_loess to FALSE or set n_bins to 0 and show_loess to TRUE. Both cannot be displayed.")}how_many_models = df[[model]] %>% unique() %>% length()df[[outcome]] = ifelse(positive == df[[outcome]], 1, 0)if (n_bins > 0) {df <- df %>% dplyr::group_by(!!rlang::parse_expr(model)) %>% dplyr::mutate(bin = dplyr::ntile(!!rlang::parse_expr(prediction), n_bins)) %>% dplyr::group_by(!!rlang::parse_expr(model), bin) %>% dplyr::mutate(n = dplyr::n(), bin_pred = mean(!!rlang::parse_expr(prediction), na.rm = TRUE), bin_prob = mean(as.numeric(as.character(!!rlang::parse_expr(outcome))), na.rm = TRUE), se = sqrt((bin_prob * (1 - bin_prob))/n), ul = bin_prob + 1.96 * se, ll = bin_prob - 1.96 * se) %>% dplyr::mutate_at(dplyr::vars(ul, ll), . %>% scales::oob_squish(range = c(0, 1))) %>% dplyr::ungroup()}g1 = ggplot2::ggplot(df) + ggplot2::scale_y_continuous(limits = c(0, 1), breaks = seq(0, 1, by = 0.1)) + ggplot2::scale_x_continuous(limits = c(0, 1), breaks = seq(0, 1, by = 0.1)) + ggplot2::geom_abline(linetype = "dashed")if (show_loess == TRUE) {g1 = g1 + ggplot2::stat_smooth(ggplot2::aes(x = !!rlang::parse_expr(prediction), y = as.numeric(!!rlang::parse_expr(outcome)), color = !!rlang::parse_expr(model), fill = !!rlang::parse_expr(model)), se = TRUE, method = "loess")}else {g1 = g1 + ggplot2::aes(x = bin_pred, y = bin_prob, color = !!rlang::parse_expr(model), fill = !!rlang::parse_expr(model)) + ggplot2::geom_ribbon(ggplot2::aes(ymin = ll, ymax = ul, ), alpha = 1/how_many_models) + ggplot2::geom_point(size = 2) + ggplot2::geom_line(size = 1, alpha = 1/how_many_models)}g1 = g1 + ggplot2::xlab("Predicted Probability") + ggplot2::ylab("Observed Risk") + ggplot2::scale_color_brewer(name = "Models", palette = "Set1") + ggplot2::scale_fill_brewer(name = "Models", palette = "Set1") + ggplot2::theme_minimal() + ggplot2::theme(aspect.ratio = 1) + ggplot2::ggtitle(plot_title)g2 <- ggplot2::ggplot(df, ggplot2::aes(x = !!rlang::parse_expr(prediction))) + ggplot2::geom_density(alpha = 1/how_many_models, ggplot2::aes(fill = !!rlang::parse_expr(model), color = !!rlang::parse_expr(model))) + ggplot2::scale_x_continuous(limits = c(0, 1), breaks = seq(0, 1, by = 0.1)) + ggplot2::coord_fixed() + ggplot2::xlab("") + ggplot2::ylab("") + ggplot2::scale_color_brewer(palette = "Set1") + ggplot2::scale_fill_brewer(palette = "Set1") + ggplot2::theme_minimal() + ggeasy::easy_remove_y_axis() + ggeasy::easy_remove_legend(fill, color) + ggplot2::theme_void() + ggplot2::theme(aspect.ratio = 0.1)layout = c(patchwork::area(t = 1, b = 10, l = 1, r = 10), patchwork::area(t = 11, b = 12, l = 1, r = 10))g1/g2
}

自己寫的函數

# 讀取數據
data <- prediction_all_rfget_cal <- function(data=prediction_all_rf){data <- data %>% mutate(bucket = ntile(prob.Case, 10))bucket_means <- data %>% group_by(bucket) %>% summarise(predicted_prob = mean(prob.Case))actual_probs <- data %>% group_by(bucket) %>% summarise(actual_prob = mean(truth == "Case"))calibration_data <- left_join(bucket_means, actual_probs, by = "bucket")calibration_data$type=data$type[1]return(calibration_data)
}cal_rf <- get_cal(data = prediction_all_rf)
cal_kkmm <- get_cal(data = prediction_all_kknn)
cal_SVM <- get_cal(data = prediction_all_SVM)
cal_xgb <- get_cal(data = prediction_all_xgb)calibration_data <- rbind(cal_rf,cal_kkmm,cal_SVM,cal_xgb)
# ⑥ 將這些點連起來，即為校準曲線中的Apparent線
ggplot(calibration_data, aes(x = predicted_prob, y = actual_prob,group = type,colour = type)) +geom_point() +geom_line() +labs(title = "Calibration Curve", x = "Predicted Probability", y = "Actual Probability") +theme_minimal()+ggplot2::scale_y_continuous(limits = c(0, 1), breaks = seq(0, 1, by = 0.1)) + ggplot2::scale_x_continuous(limits = c(0, 1), breaks = seq(0, 1, by = 0.1)) + ggplot2::geom_abline(linetype = "dashed")