add continuous boyce index

src/confusion_matrices.jl

@@ -9,21 +9,12 @@ module ConfusionMatrices
 using CategoricalArrays
 using OrderedCollections
 import ..Functions
+import ..warn_unordered
 
 const CM  = "ConfusionMatrices"
 const CatArrOrSub{T, N} =
     Union{CategoricalArray{T, N}, SubArray{T, N, <:CategoricalArray}}
 
-function WARN_UNORDERED(levels)
-    raw_levels = CategoricalArrays.unwrap.(levels)
-    ret = "Levels not explicitly ordered. "*
-        "Using the order $raw_levels. "
-    if length(levels) == 2
-        ret *= "The \"positive\" level is $(raw_levels[2]). "
-    end
-    ret
-end
-
 const ERR_INDEX_ACCESS_DENIED = ErrorException(
     "Direct access by index of unordered confusion matrices dissallowed. "*
         "Access by level, as in `some_confusion_matrix(\"male\", \"female\")` or first "*
@@ -343,7 +334,7 @@ Return the regular `Matrix` associated with confusion matrix `m`.
 """
 matrix(cm::ConfusionMatrix{N,true}; kwargs...) where N = cm.mat
 @inline function matrix(cm::ConfusionMatrix{N,false}; warn=true) where N
-    warn && @warn WARN_UNORDERED(levels(cm))
+    warn && warn_unordered(levels(cm))
     cm.mat
 end
 

src/functions.jl

@@ -249,6 +249,59 @@ function matthews_correlation(m)
     return mcc
 end
 
+"""
+    Functions.cbi(
+        probability_of_positive, ground_truth_observations, positive_class, 
+        nbins, binwidth, ma=maximum(scores), mi=minimum(scores), cor=corspearman
+    )
+    Return the Continuous Boyce Index (CBI) for a vector of probabilities and ground truth observations.
+
+"""
+function cbi(
+    scores, y, positive_class; 
+    nbins, binwidth,
+    max=maximum(scores), min=minimum(scores), cor=StatsBase.corspearman
+)
+    binstarts = range(min, stop=max-binwidth, length=nbins)
+    binends = binstarts .+ binwidth
+
+    sorted_indices = sortperm(scores)
+    sorted_scores = view(scores, sorted_indices)
+    sorted_y = view(y, sorted_indices)
+
+    n_positive = zeros(Int, nbins)
+    n_total = zeros(Int, nbins)
+    empty_bins = falses(nbins)
+    any_empty = false
+
+    @inbounds for i in 1:nbins
+        bin_index_first = searchsortedfirst(sorted_scores, binstarts[i])
+        bin_index_last = searchsortedlast(sorted_scores, binends[i])
+        if bin_index_first > bin_index_last
+            empty_bins[i] = true
+            any_empty = true
+        end
+        @inbounds for j in bin_index_first:bin_index_last
+            if sorted_y[j] == positive_class 
+                n_positive[i] += 1
+            end
+        end
+        n_total[i] = bin_index_last - bin_index_first + 1
+    end
+    if any_empty
+        @info "removing $(sum(empty_bins)) bins without any observations"
+        deleteat!(n_positive, empty_bins)
+        deleteat!(n_total, empty_bins)
+        binstarts = binstarts[.!empty_bins]
+    end
+
+    # calculate "PE-ratios" - a bunch of things cancel out but that does not matter for
+    # any correlation calculation
+    PE_ratios = n_positive ./ n_total
+    return cor(PE_ratios, binstarts)
+end
+
+
 
 # ## binary, but NOT invariant under class relabellings
 

src/probabilistic.jl

@@ -544,3 +544,85 @@ $DOC_DISTRIBUTIONS
 SphericalScore
 "$SphericalScoreDoc"
 const spherical_score = SphericalScore()
+
+
+# ---------------------------------------------------------------------
+# Continuous Boyce Index
+struct _ContinuousBoyceIndex 
+    nbins::Integer
+    bin_overlap::AbstractFloat
+    min::Union{AbstractFloat, Nothing}
+    max::Union{AbstractFloat, Nothing}
+    cor::Function
+    function _ContinuousBoyceIndex(; nbins = 101, bin_overlap = 0.1, min = nothing, max = nothing, cor = StatsBase.corspearman)
+        new(nbins, bin_overlap, min, max, cor)
+    end
+end
+
+ContinuousBoyceIndex(; kw...) = _ContinuousBoyceIndex(; kw...) |> robust_measure |> fussy_measure
+
+function (m::_ContinuousBoyceIndex)(ŷ::UnivariateFiniteArray, y::NonMissingCatArrOrSub; warn=true)
+    warn && warn_unordered(levels(y))
+    positive_class = levels(first(ŷ))|> last
+    scores = pdf.(ŷ, positive_class)
+    max = isnothing(m.max) ? maximum(scores) : m.max
+    min = isnothing(m.min) ? minimum(scores) : m.min
+    binwidth = m.bin_overlap * (max - min)
+
+    return Functions.cbi(scores, y, positive_class; nbins = m.nbins, binwidth, max, min, cor = m.cor)
+end
+
+const ContinuousBoyceIndexType = API.FussyMeasure{<:API.RobustMeasure{<:_ContinuousBoyceIndex}}
+
+@fix_show ContinuousBoyceIndex::ContinuousBoyceIndexType
+
+StatisticalMeasures.@trait(
+    _ContinuousBoyceIndex,
+    consumes_multiple_observations=true,
+    observation_scitype = Finite{2},
+    kind_of_proxy=StatisticalMeasures.LearnAPI.Distribution(),
+    orientation=Score(),
+    external_aggregation_mode=Mean(),
+    human_name = "continuous boyce index",
+)
+
+register(ContinuousBoyceIndex, "continuous_boyce_index", "cbi")
+
+const ContinuousBoyceIndexDoc = docstring(
+    "ContinuousBoyceIndex(; nbins=101, bin_overlap=0.1, min=nothing, max=nothing, cor=StatsBase.corspearman)",
+    body=
+"""
+The Continuous Boyce Index is a measure for evaluating the performance of probabilistic predictions for binary classification, 
+especially for presence-background data in ecological modeling. 
+It compares the predicted probability scores for the positive class across bins, giving higher scores if the ratio of positive
+    and negative samples in each bin is strongly correlated to the value at that bin.
+
+## Keywords
+
+- `nbins`: Number of bins to use for score partitioning.
+- `bin_overlap`: Relative overlap between bins. The actual bin width is `bin_overlap * (max - min)`.
+- `min`, `max`: Optional minimum and maximum score values for binning. Default to the
+    maximum and minimum values of the scores.
+- `cor`: Correlation function (defaults to StatsBase.corspearman, i.e. Spearman correlation).
+
+## Arguments
+
+The predictions `ŷ` should be a vector of `UnivariateFinite` distributions from CategoricalDistributions.jl, and `y` a vector of ground truth labels.
+
+Returns the correlation between the ratio of positive to negative samples in each bin and the bin centers.
+
+Core implementation: [`Functions.cbi`](@ref).
+
+Reference:
+Alexandre H. Hirzel, Gwenaëlle Le Lay, Véronique Helfer, Christophe Randin, Antoine Guisan,
+Evaluating the ability of habitat suitability models to predict species presences,
+Ecological Modelling,
+Volume 199, Issue 2, 2006
+""",
+    scitype="",
+)
+
+"$ContinuousBoyceIndexDoc"
+ContinuousBoyceIndex
+"$ContinuousBoyceIndexDoc"
+const cbi(x, y) = ContinuousBoyceIndex()(x, y)

src/roc.jl

@@ -23,10 +23,7 @@ function binary_levels(
     length(classes) == 2 || throw(ERR_ROC2)
     API.check_numobs(yhat, y)
     API.check_pools(yhat, y)
-    if !(yhat isa AbstractArray{<:UnivariateFinite{<:OrderedFactor}}) ||
-        !CategoricalArrays.isordered(y)
-        @warn ConfusionMatrices.WARN_UNORDERED(classes)
-    end
+    warn_unordered(classes)
     classes
 end
 binary_levels(

src/tools.jl

@@ -59,3 +59,15 @@ function API.check_pools(
     return nothing
 end
 
+# Throw a warning if levels are not explicitly ordered
+function warn_unordered(levels)
+    levels isa CategoricalArray && CategoricalArrays.isordered(levels) && return
+    raw_levels = CategoricalArrays.unwrap.(levels)
+    ret = "Levels not explicitly ordered. "*
+        "Using the order $raw_levels. "
+    if length(levels) == 2
+        ret *= "The \"positive\" level is $(raw_levels[2]). "
+    end
+    @warn ret
+    return ret
+end

test/confusion_matrices.jl

@@ -30,7 +30,7 @@ const CM = StatisticalMeasures.ConfusionMatrices
     rev_index_given_level = Dict("B" => 1, "A" => 2)
     @test cm == CM.ConfusionMatrix(n, rev_index_given_level)
     mat = @test_logs(
-        (:warn, CM.WARN_UNORDERED(levels)),
+        (:warn, StatisticalMeasures.warn_unordered(levels)),
         CM.matrix(cm),
     )
     @test mat == m

test/finite.jl

@@ -114,7 +114,7 @@ end
          1, 1, 1, 2, 2, 1, 2, 1, 2, 2, 2, 1, 2,
          1, 2, 2, missing]
 
-    @test_logs (:warn, CM.WARN_UNORDERED([1, 2])) f1score(ŷ, y)
+    @test_logs (:warn, StatisticalMeasures.warn_unordered([1, 2])) f1score(ŷ, y)
     f05 = @test_logs FScore(0.5, levels=[1, 2])(ŷ, y)
     sk_f05 = 0.625
     @test f05 ≈ sk_f05 # m.fbeta_score(y, yhat, 0.5, pos_label=2)

test/probabilistic.jl

@@ -180,6 +180,50 @@ end
     @test_throws StatisticalMeasures.ERR_UNSUPPORTED_ALPHA s(yhat, [1.0, 1.0])
 end
 
+@testset "ContinuousBoyceIndex" begin
+    rng = srng(1234)
+    # Simple synthetic test: perfectly separates positives and negatives
+    c = ["neg", "pos"]
+    probs = repeat(0.0:0.1:0.9, inner = 10) .+ rand(rng, 100) .* 0.1
+    y = categorical(probs .> rand(rng, 100))
+    ŷ = UnivariateFinite(levels(y), probs, augment=true)
+    # Should be pretty high
+    @test cbi(ŷ, y) ≈ 0.86 atol=0.01
+
+    # Passing different correlation methods works
+    @test ContinuousBoyceIndex(cor=cor)(ŷ, y) ≈ 0.90 atol = 0.01
+    @test ContinuousBoyceIndex(nbins = 11)(ŷ, y) ≈ 0.83 atol = 0.01 
+
+    # Randomized test: shuffled labels, should be near 0
+    y_shuf = copy(y)
+    MLUtils.shuffle!(rng, y_shuf)
+    @test (cbi(ŷ, y_shuf)) ≈ 0.0 atol=0.1
+
+    # Test invariance to order
+    idx = randperm(length(y))
+    @test isapprox(cbi(ŷ[idx], y[idx]), cbi(ŷ, y), atol=1e-8)
+
+    # Test with all positives or all negatives return NaN
+    y_allpos = categorical(trues(100), levels = levels(y))
+    y_allneg = categorical(falses(100), levels = levels(y))
+    @test isnan(cbi(ŷ, y_allpos))
+    @test isnan(cbi(ŷ, y_allneg))
+
+    unordered_warning = StatisticalMeasures.warn_unordered([false, true])
+    @test_logs(
+        (:warn, unordered_warning),
+        cbi(ŷ, y),
+     )
+
+    cbi_dropped_bins = @test_logs(
+        (:warn, unordered_warning), (:info, "removing 91 bins without any observations",),
+        ContinuousBoyceIndex(; min =0.0, max = 2.0, nbins = 191, bin_overlap = 0.05)(ŷ, y),
+    )
+    # These two are identical because bins are dropped
+    @test cbi_dropped_bins == 
+        ContinuousBoyceIndex(; min = 0.0, max = 1.2, nbins = 111, bin_overlap = 0.1/1.2)(ŷ, y)
+end
+
 @testset "l2_check" begin
     d = Distributions.Normal()
     yhat = Union{Distributions.Sampleable,Missing}[d, d, missing]

test/roc.jl

@@ -25,7 +25,7 @@
     )
 
     fprs, tprs, ts = @test_logs(
-        (:warn, ConfusionMatrices.WARN_UNORDERED([0, 1])),
+        (:warn, StatisticalMeasures.warn_unordered([0, 1])),
         roc_curve(ŷ, y),
      )
 

test/runtests.jl

@@ -11,6 +11,7 @@ using OrderedCollections
 using CategoricalDistributions
 using LinearAlgebra
 import Distributions
+import StatsBase: corspearman, randperm
 
 const CM = ConfusionMatrices