1200: Fix predicted_label handling and n > 1 column check in categorical metrics

NEWS.md

@@ -1,5 +1,6 @@
 # scoringutils (development version)
 
+- Fixed `rps_ordinal()` and `logs_categorical()` returning wrong scores when called directly with a `predicted_label` that was not in the order of the factor levels: `rps_ordinal()` applied the wrong (forward instead of inverse) permutation when reordering the columns of `predicted`, and `logs_categorical()` ignored `predicted_label` entirely. Scores are now invariant to how the columns of `predicted` are labelled. Forecasts scored via `score()` were unaffected, as the pipeline sorts predictions into level order before calling the metrics. Additionally, input validation for categorical forecasts now errors when `predicted` has more or fewer columns than there are factor levels for inputs with more than one observation (previously such input was silently accepted and scored meaninglessly; the check already existed for a single observation) (#1200).
 - Added an internal helper `prepare_forecast_for_scoring()` that consolidates the input preparation steps previously duplicated across the `score()` methods: cleaning the forecast, validating the metrics and converting to a plain `data.table`, plus determining the forecast unit for the methods that need it (#941).
 - Fixed the sample-based metrics `bias_sample()`, `ae_median_sample()`, `se_mean_sample()` and `mad_sample()` mishandling the documented vector input for a single observation (a scalar `observed` with `predicted` given as a vector of samples). `ae_median_sample()` and `se_mean_sample()` silently treated the samples as separate one-sample forecasts and returned wrong results, while `bias_sample()` and `mad_sample()` errored. All sample metrics now treat this input as one forecast with N samples, consistent with `crps_sample()`. Also corrected the integer bias formula in the `bias_sample()` documentation, which stated P_t(x_t + 1) instead of P_t(x_t - 1) (the code was correct) (#1197).
 - `summarise_scores()` now errors when `by` contains a metric column (e.g. `by = c("model", "wis")`). Previously, such calls silently returned an unsummarised table with duplicate column names, because the score column was used both as a grouping column and as a column to summarise. This complements the fix for the empty-metrics case in #1179 (#1204).

R/metrics-nominal.R

@@ -69,7 +69,7 @@ assert_input_categorical <- function(
     )
     summed_predictions <- .rowSums(predicted, m = 1, n = N, na.rm = TRUE)
   } else {
-    assert_matrix(predicted, nrows = n)
+    assert_matrix(predicted, nrows = n, ncols = N)
     summed_predictions <- round(rowSums(predicted, na.rm = TRUE), 10) # avoid numeric errors
   }
   # Allow for numeric errors
@@ -130,6 +130,11 @@ logs_categorical <- function(observed, predicted, predicted_label) {
   if (n == 1) {
     predicted <- matrix(predicted, nrow = 1)
   }
+  # Reorder the predicted matrix columns to match the natural level order:
+  # column k of `predicted` holds probabilities for level
+  # `as.numeric(predicted_label)[k]`, so the inverse permutation given by
+  # `order()` puts the column for level k into position k
+  predicted <- predicted[, order(as.numeric(predicted_label)), drop = FALSE]
   observed_indices <- as.numeric(observed)
   pred_for_observed <- predicted[cbind(1:n, observed_indices)]
   logs <- -log(pred_for_observed)

R/metrics-ordinal.R

@@ -71,9 +71,12 @@ rps_ordinal <- function(observed, predicted, predicted_label) {
     predicted <- matrix(predicted, nrow = 1)
   }
 
-  # Reorder the predicted matrix columns to match the natural ordering
-  correct_order <- as.numeric(predicted_label)
-  ordered_predicted <- predicted[, correct_order]
+  # Reorder the predicted matrix columns to match the natural level order:
+  # column k of `predicted` holds probabilities for level
+  # `as.numeric(predicted_label)[k]`, so the inverse permutation given by
+  # `order()` puts the column for level k into position k
+  ordered_predicted <-
+    predicted[, order(as.numeric(predicted_label)), drop = FALSE]
 
   rps <- rps_probs(as.numeric(observed), ordered_predicted)
   return(rps)

tests/testthat/test-metrics-nominal.R

@@ -93,6 +93,21 @@ test_that("Input checking for nominal forecasts works", {
   expect_no_condition(
     assert_input_nominal(observed, predicted2, predicted_label)
   )
+
+  # n > 1 with wrong number of columns
+  predicted_wide <- matrix(
+    c(
+      0.2, 0.1, 0.4, 0.3,
+      0.1, 0.25, 0.25, 0.4,
+      0.25, 0.25, 0.25, 0.25
+    ),
+    nrow = 3,
+    byrow = TRUE
+  )
+  expect_error(
+    assert_input_nominal(observed, predicted_wide, predicted_label),
+    "Must have exactly 3 cols, but has 4 cols"
+  )
 })
 
 
@@ -139,3 +154,16 @@ test_that("logs_categorical() works as expected", {
     c(NA, res_manual[-1])
   )
 })
+
+test_that("logs_categorical() is invariant to permutations of predicted_label", {
+  perm <- c(2, 3, 1)
+  predicted_label2 <- factor(factor_levels[perm], levels = factor_levels)
+  res <- logs_categorical(observed, predicted[, perm], predicted_label2)
+  expect_equal( # nolint: expect_identical_linter
+    res,
+    logs_categorical(observed, predicted, predicted_label)
+  )
+  # explicitly -log() of the probability assigned to the observed outcome
+  res_manual <- -log(c(predicted[1, 1], predicted[2, 2], predicted[3, 2]))
+  expect_equal(res, res_manual) # nolint: expect_identical_linter
+})

tests/testthat/test-metrics-ordinal.R

@@ -120,10 +120,24 @@ test_that("rps_ordinal() works as expected", {
   expect_equal(res, result) # nolint: expect_identical_linter
 
   # works with changed order of levels
+  # column k of `predicted` holds probabilities for predicted_label2[k], so
+  # the columns in natural level order are given by the inverse permutation
+  # order(c(3, 1, 2)) = c(2, 3, 1)
   predicted_label2 <- factor(c("three", "one", "two"), levels = factor_levels, ordered = TRUE)
   expect_equal( # nolint: expect_identical_linter
     rps_ordinal(observed, predicted, predicted_label2),
-    scoringRules::rps_probs(as.numeric(observed), predicted[, c(3, 1, 2)])
+    scoringRules::rps_probs(as.numeric(observed), predicted[, order(c(3, 1, 2))])
+  )
+})
+
+test_that("rps_ordinal() is invariant to permutations of predicted_label", {
+  perm <- c(2, 3, 1)
+  predicted_label2 <- factor(
+    factor_levels[perm], levels = factor_levels, ordered = TRUE
+  )
+  expect_equal( # nolint: expect_identical_linter
+    rps_ordinal(observed, predicted[, perm], predicted_label2),
+    rps_ordinal(observed, predicted, predicted_label)
   )
 })