Skip to content

Commit

Permalink
update rps_ordinal function + test
Browse files Browse the repository at this point in the history
  • Loading branch information
nikosbosse committed Dec 4, 2024
1 parent fdfd06e commit 6dae1d7
Show file tree
Hide file tree
Showing 4 changed files with 175 additions and 17 deletions.
1 change: 1 addition & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -190,6 +190,7 @@ importFrom(purrr,partial)
importFrom(scoringRules,crps_sample)
importFrom(scoringRules,dss_sample)
importFrom(scoringRules,logs_sample)
importFrom(scoringRules,rps_probs)
importFrom(stats,cor)
importFrom(stats,mad)
importFrom(stats,median)
Expand Down
23 changes: 8 additions & 15 deletions R/metrics-ordinal.R
Original file line number Diff line number Diff line change
Expand Up @@ -128,13 +128,14 @@ logs_ordinal <- function(observed, predicted, predicted_label) {
#' @returns A numeric vector of size n with ranked probability scores
#' @inheritSection illustration-input-metric-nominal Input format
#' @importFrom methods hasArg
#' @importFrom scoringRules rps_probs
#' @export
#' @keywords metric
#' @family scoring functions
#' @examples
#' factor_levels <- c("one", "two", "three")
#' factor_levels <- c("one", "three", "two")
#' predicted_label <- factor(factor_levels, levels = factor_levels, ordered = TRUE)
#' observed <- factor(c("one", "three", "two"), levels = factor_levels, ordered = TRUE)
#' observed <- factor(c("three", "three", "two"), levels = factor_levels, ordered = TRUE)
#' predicted <- matrix(
#' c(0.8, 0.1, 0.4,
#' 0.1, 0.2, 0.4,
Expand All @@ -149,19 +150,11 @@ rps_ordinal <- function(observed, predicted, predicted_label) {
predicted <- matrix(predicted, nrow = 1)
}

# Calculate cumulative probabilities for predictions
cum_pred <- t(apply(predicted, 1, cumsum))

# Create matrix of cumulative probabilities for observations
N <- ncol(predicted)
observed_indices <- as.numeric(observed)
cum_obs <- matrix(0, nrow = n, ncol = N)
for (i in 1:n) {
cum_obs[i, observed_indices[i]:N] <- 1
}

# Calculate RPS as mean squared difference of cumulative probabilities
rps <- rowMeans((cum_pred - cum_obs)^2)
# Reorder the predicted matrix columns to match the natural ordering
correct_order <- as.numeric(predicted_label)
ordered_predicted <- predicted[, correct_order]

# Use scoringRules implementation
rps <- scoringRules::rps_probs(as.numeric(observed), ordered_predicted)
return(rps)
}
4 changes: 2 additions & 2 deletions man/rps_ordinal.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

164 changes: 164 additions & 0 deletions tests/testthat/test-metrics-ordinal.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,164 @@
factor_levels <- c("one", "two", "three")
predicted_label <- factor(c("one", "two", "three"), levels = factor_levels, ordered = TRUE)
observed <- factor(c("one", "two", "two"), levels = factor_levels, ordered = TRUE)
predicted <- matrix(c(0.8, 0.1, 0.4,
0.1, 0.2, 0.4,
0.1, 0.7, 0.2),
nrow = 3)

# ============================================================================ #
# Input handling ===============================================================
# ============================================================================ #
test_that("Input checking for ordinal forecasts works", {
# everything correct
expect_no_condition(
assert_input_ordinal(observed, predicted, predicted_label)
)

# works with a single observations
expect_no_condition(
assert_input_ordinal(observed[1], predicted[1, ], predicted_label)
)
expect_no_condition(
assert_input_ordinal(observed[1], matrix(predicted[1, ], nrow = 1), predicted_label)
)

# observed and predicted_label have different level lengths
observed_wrong <- factor(observed, levels = c("one", "two"), ordered = TRUE)
expect_error(
assert_input_ordinal(observed_wrong, predicted, predicted_label),
"Assertion on 'predicted_label' failed: Must have length 2, but has length 3."
)

# observed and predicted_label have different levels, resulting in an `NA` in predicted_label
predicted_label_wrong <- factor(predicted_label, levels = c("one", "two", "four"))
expect_error(
assert_input_ordinal(observed, predicted, predicted_label_wrong),
"Assertion on 'predicted_label' failed: Contains missing values"
)

# 6 observations, but only 3 forecasts
observed_wrong <- factor(
c("one", "two", "one", "one", "two", "three"),
levels = c("one", "two", "three"),
ordered = TRUE
)
expect_error(
assert_input_ordinal(observed_wrong, predicted, predicted_label),
"Assertion on 'predicted' failed: Must have exactly 6 rows, but has 3 rows."
)

# observed value or predicted_label is not factor
expect_error(
assert_input_ordinal(as.numeric(observed), predicted, predicted_label),
"Assertion on 'observed' failed: Must be of type 'factor', not 'double'."
)
expect_error(
assert_input_ordinal(observed, predicted, as.character(predicted_label)),
"Assertion on 'predicted_label' failed: Must be of type 'factor', not 'character'."
)

# a single observation
expect_no_condition(
assert_input_ordinal(observed[1], predicted[1, ], predicted_label),
)
expect_no_condition(
assert_input_ordinal(observed[1], as.numeric(predicted[1, ]), predicted_label),
)

# wrong dimensions
expect_error(
assert_input_ordinal(observed[1], predicted, predicted_label),
"One of the following must apply:"
)

# NA values in observed are permitted
observed2 <- factor(c("one", "two", "missing"), levels = c("one", "two", "three"), ordered = TRUE)
expect_no_condition(
assert_input_ordinal(observed2, predicted, predicted_label)
)

# NA values in predicted are permitted and treated as zeros
# (as long as probabilities still sum to one)
predicted2 <- predicted
predicted2[2, 2] <- NA
expect_error(
assert_input_ordinal(observed, predicted2, predicted_label),
"Probabilities belonging to a single forecast must sum to one"
)
predicted2[2, 1] <- 0.3
expect_no_condition(
assert_input_ordinal(observed, predicted2, predicted_label)
)
})


# ============================================================================ #
# logs ordinal =============================================================== #
# ============================================================================ #
test_that("logs_ordinal() works as expected", {
expect_no_condition(
res <- logs_ordinal(observed, predicted, predicted_label)
)
res_manual <- -log(c(predicted[1, 1], predicted[2, 2], predicted[3, 2]))
expect_equal(res, res_manual)

# works with a single observations
expect_equal(
logs_ordinal(observed[1], predicted[1, ], predicted_label),
res_manual[1]
)
expect_equal(
logs_ordinal(observed[1], matrix(predicted[1, ], nrow = 1), predicted_label),
res_manual[1]
)

# log scores still work if there are NA values in probabilities that have
# not materialised
predicted2 <- predicted
predicted2[2, c(1, 3)] <- c(NA, 0.8)
expect_equal(
logs_ordinal(observed, predicted2, predicted_label),
res_manual
)

# relevant NAs in predictions lead to NA outcomes
predicted2[1, 1:2] <- c(NA, 0.9)
expect_equal(
logs_ordinal(observed, predicted2, predicted_label),
c(NA, res_manual[-1])
)

# NA values in observed values lead to NAs in result
observed2 <- factor(c(NA, "two", "two"), levels = c("one", "two", "three"), ordered = TRUE)
expect_equal(
logs_ordinal(observed2, predicted, predicted_label),
c(NA, res_manual[-1])
)
})

# ============================================================================ #
# rps ordinal =============================================================== #
# ============================================================================ #
result <- c(0.05, 0.5, 0.2)
test_that("rps_ordinal() works as expected", {
expect_no_condition(
res <- rps_ordinal(observed, predicted, predicted_label)
)

expect_equal(res, result)

# works with changed order of levels
predicted_label2 <- factor(c("three", "one", "two"), levels = factor_levels, ordered = TRUE)
expect_equal(
rps_ordinal(observed, predicted, predicted_label2),
scoringRules::rps_probs(as.numeric(observed), predicted[, c(3, 1, 2)])
)
})

test_that("rps_ordinal() works with a single observation", {
expect_equal(
rps_ordinal(observed[1], predicted[1, ], predicted_label),
result[1]
)
})

0 comments on commit 6dae1d7

Please sign in to comment.