Rを使って一般的なCalibration（キャリブレーション）をプロットする方法

統計の講師をしていて、なるほど！と思う質問を受けることがあります。
せっかくなので、その中からピックアップして紹介できればと思います。

ロジスティック回帰モデルのCalibration（キャリブレーション）
https://skellington.blog.ss-blog.jp/2019-12-03
↑
ロジスティック回帰モデルに関するCalibrationはこちら。

【質問】
ロジスティック回帰以外に、生存時間とか他の統計モデルを作った場合に、Rを使って一般的なCalibration（キャリブレーション）をプロットする方法するにはどうすれば良いか？

【回答】
パッケージを使うのではなく、一般的な方法を作ってみました。

以下のRのコードは、
1.ビンに分割をする（10分割とか20分割）
2.予測確率を20分割する
3.分割したビンごとに、予測値の平均と実際の{0, 1}の平均を集計しプロットする
という流れになります。

# データの取得
dat <- read.csv("https://stats.idre.ucla.edu/stat/data/binary.csv")
head(dat)

# logit model
model.logit <- glm(admit ~ gre + gpa + rank, family = "binomial", data = dat)
summary(model.logit)

# 予測値
dat.pred <- predict(model.logit, newdata = dat, type = "response")
head(dat.pred)

# install.packages("infotheo")
library(infotheo)

pred.bin <- discretize(dat.pred, disc="equalwidth",nbins=20) #20分割する

df.calib <- data.frame(dat[ ,1]    #実際の値{0, 1}
                      , dat.pred   #予測値
                      , pred.bin   #予測値を20分割した値
)
names(df.calib) <- c("ActualProbability", "PredictedProbability", "PredictedProbabilityBin")
head(df.calib)

# ビンごとの集計値をプロットする
plot(
  aggregate(ActualProbability~PredictedProbabilityBin, data=df.calib, FUN=mean)[,2]
  , aggregate(PredictedProbability~PredictedProbabilityBin, data=df.calib, FUN=mean)[,2]
  , main="キャリブレーション"
  , xlab="Actual Probability", ylab="Predicted Probability"
  , xlim=c(0.0, 1.0), ylim=c(0.0, 1.0))
)
abline(0, 1, col="red")