Skip to content

Instantly share code, notes, and snippets.

@swiftsam
Created August 4, 2014 16:11
Show Gist options
  • Save swiftsam/67514a6d608d0a0bc416 to your computer and use it in GitHub Desktop.
Save swiftsam/67514a6d608d0a0bc416 to your computer and use it in GitHub Desktop.
ordinal brier scoring
###~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
### calculates Brier Score for a single forecast on an Ordered IFP
### calculation based on definition in Appendix B of ACE_TE_Plan_OY1_v3.
### procedure reduces penalty to probability assigned to bins adjacent
### to the realized outcome
###~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
ScoreOrdered <- function(fcast.values,
correct.values,
scoring.fn){
# Step 1 – Take the original categories (A-B-C-D) and break them up
# into a set of binary categories (A-BCD; AB-CD; ABC-D)
n.pairs <- length(fcast.values) - 1
pairs <- array(NA, dim=c(n.pairs,2))
for(i in 1:n.pairs){
forecast.sub1 <- sum(fcast.values[1:i])
forecast.sub2 <- sum(fcast.values[-(1:i)])
pairs[i,] <- c(forecast.sub1,forecast.sub2)
}
# Step 1a define ref.value (correct option) for each pair
correct.arr <- array(NA, dim=dim(pairs))
for(i in 1:n.pairs){
ref.sub1 <- sum(correct.values[1:i])
ref.sub2 <- sum(correct.values[-(1:i)])
correct.arr[i,] <- c(ref.sub1,ref.sub2)
}
if(max(correct.arr) > 1){
stop(paste(Sys.time(),"BrierScoreOrdered: ref.vals calculation is flawed, values > 1 [",
paste(correct.arr,collapse=","),"]"))
}
# Step 2 – Apply the scoring rule to each binary pair
n.pairs <- dim(pairs)[1]
score.pairs <- array(dim=n.pairs)
for(i in 1:n.pairs){
score.pairs[i] <- scoring.fn(fcast.values = pairs[i,],
correct.values = correct.arr[i,])
}
# Step 3 - Take an average across the binary pair scores.
score <- mean(score.pairs)
return(score)
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment