R：插入符号包：Brier分数

Question

我想使用 caret 包中的 train() 函数执行逻辑回归。我的模型看起来像这样：

model <- train(Y ~.,
  data = train_data,
  family = "binomial",
  method = "glmnet")

使用生成的模型，我想进行预测：

pred <- predict(model, newdata = test_data, s = "lambda.min", type = "prob")

现在，我想评估模型预测与实际测试数据相比有多好。为此，我知道如何接收 ROC 和 AUC。不过，我也有兴趣获得 BRIER SCORE。 Brier 分数的公式几乎与 MSE 相同。 我面临的问题是 predict 中的 type 参数只允许“概率”（或我不感兴趣的“类”），它给出了一个预测为 ONE（例如 0.64）的概率，以及补充成为零的概率（例如 0.37）。然而，对于 Brier 分数，我需要对包含两者信息的每个预测进行一个概率估计（例如，高于 0.5 的值表示 1，低于 0.5 的值表示 0）。 我还没有找到任何解决方案来接收 caret 包中的 Brier 分数。我知道使用包cv.glmnet predict 函数允许参数“响应”，这将解决我的问题。但是，出于个人喜好，我想继续使用caretpackage。 感谢您的帮助！

Answer 1

我使用 Brier 分数在 caret 中调整我的模型以进行二元分类。我确保“肯定”类是第二类，这是您将响应标记为“0:1”时的默认值。然后我创建了这个主汇总函数，基于caret 自己的汇总函数套件，以返回我想查看的所有指标：

BigSummary <- function (data, lev = NULL, model = NULL) {
  pr_auc <- try(MLmetrics::PRAUC(data[, lev[2]],
                                 ifelse(data$obs == lev[2], 1, 0)),
                silent = TRUE)
  brscore <- try(mean((data[, lev[2]] - ifelse(data$obs == lev[2], 1, 0)) ^ 2),
               silent = TRUE)
  rocObject <- try(pROC::roc(ifelse(data$obs == lev[2], 1, 0), data[, lev[2]],
                             direction = "<", quiet = TRUE), silent = TRUE)
  if (inherits(pr_auc, "try-error")) pr_auc <- NA
  if (inherits(brscore, "try-error")) brscore <- NA
  rocAUC <- if (inherits(rocObject, "try-error")) {
    NA
  } else {
    rocObject$auc
  }
  tmp <- unlist(e1071::classAgreement(table(data$obs,
                                            data$pred)))[c("diag", "kappa")]
  out <- c(Acc = tmp[[1]],
           Kappa = tmp[[2]],
           AUCROC = rocAUC,
           AUCPR = pr_auc,
           Brier = brscore,
           Precision = caret:::precision.default(data = data$pred,
                                                 reference = data$obs,
                                                 relevant = lev[2]),
           Recall = caret:::recall.default(data = data$pred,
                                           reference = data$obs,
                                           relevant = lev[2]),
           F = caret:::F_meas.default(data = data$pred, reference = data$obs,
                                      relevant = lev[2]))
  out
}

现在我可以在trainControl 中简单地传递 summaryFunction = BigSummary，然后在train 调用中传递metric = "Brier", maximize = FALSE。

Answer 2

如果我们按照 wiki 对 Brier 分数的定义：

最常见的 Brier 评分公式是

其中 f_t 是预测的概率，o_t 是（0 或 1）的实际结果，N 是预测实例的数量。

在 R 中，如果您的标签是一个因素，那么逻辑回归将始终针对第二个级别进行预测，这意味着您只需计算概率和 0/1。例如：

library(caret)
idx = sample(nrow(iris),100)
data = iris
data$Species = factor(ifelse(data$Species=="versicolor","v","o"))
levels(data$Species)
[1] "o" "v"

在这种情况下，o 为 0，v 为 1。

train_data = data[idx,]
test_data = data[-idx,]

model <- train(Species ~.,data = train_data,family = "binomial",method = "glmnet")

pred <- predict(model, newdata = test_data)

所以我们可以看到类的概率：

head(pred)
          o          v
1 0.8367885 0.16321154
2 0.7970508 0.20294924
3 0.6383656 0.36163437
4 0.9510763 0.04892370
5 0.9370721 0.06292789

计算分数：

f_t = pred[,2]
o_t = as.numeric(test_data$Species)-1
mean((f_t - o_t)^2)
[1] 0.32