Merge pull request #66 from JuliaAI/dev

For a 0.3.4 release

Author: ablaom

.github/workflows/ci.yml

@@ -24,7 +24,7 @@ jobs:
         arch:
           - x64
     steps:
-      - uses: actions/checkout@v4
+      - uses: actions/checkout@v6
       - uses: julia-actions/setup-julia@v2
         with:
           version: ${{ matrix.version }}
@@ -51,7 +51,7 @@ jobs:
     name: Documentation
     runs-on: ubuntu-latest
     steps:
-      - uses: actions/checkout@v4
+      - uses: actions/checkout@v6
       - uses: julia-actions/setup-julia@v2
         with:
           version: '1'

Project.toml

@@ -1,7 +1,7 @@
 name = "StatisticalMeasures"
 uuid = "a19d573c-0a75-4610-95b3-7071388c7541"
 authors = ["Anthony D. Blaom <[email protected]>"]
-version = "0.3.3"
+version = "0.3.4"
 
 [deps]
 CategoricalArrays = "324d7699-5711-5eae-9e2f-1d82baa6b597"

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; 
+    verbosity, 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
+        verbosity > 1 && @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,89 @@ $DOC_DISTRIBUTIONS
 SphericalScore
 "$SphericalScoreDoc"
 const spherical_score = SphericalScore()
+
+
+# ---------------------------------------------------------------------
+# Continuous Boyce Index
+struct _ContinuousBoyceIndex 
+    verbosity::Int
+    nbins::Integer
+    binwidth::Float64
+    min::Float64
+    max::Float64
+    cor::Function
+    function _ContinuousBoyceIndex(; 
+        verbosity = 1, nbins = 101, binwidth = 0.1, 
+        min = 0, max = 1, cor = StatsBase.corspearman
+    )
+        new(verbosity, nbins, binwidth, min, max, cor)
+    end
+end
+
+ContinuousBoyceIndex(; kw...) = _ContinuousBoyceIndex(; kw...) |> robust_measure |> fussy_measure
+
+function (m::_ContinuousBoyceIndex)(ŷ::AbstractArray{<:UnivariateFinite}, y::NonMissingCatArrOrSub)
+    m.verbosity > 0 && warn_unordered(levels(y))
+    positive_class = levels(first(ŷ))|> last
+    scores = pdf.(ŷ, positive_class)
+
+    return Functions.cbi(scores, y, positive_class; 
+        verbosity = m.verbosity, nbins = m.nbins, binwidth = m.binwidth, max = m.max, min = m.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(; verbosity=1, 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
+- `verbosity`: Verbosity level.
+- `nbins`: Number of bins to use for score partitioning.
+- `binwidth`: The width of each bin, which defaults to 0.1.
+- `min`, `max`: Optional minimum and maximum score values for binning. Default to the 0 and 1, respectively.
+- `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 CategoricalVector 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 = ContinuousBoyceIndex()
+"$ContinuousBoyceIndexDoc"
+const continuous_boyce_index = ContinuousBoyceIndex()

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,53 @@ 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.87 atol=0.01
+
+    # Passing different correlation methods works
+    @test ContinuousBoyceIndex(cor=cor)(ŷ, y) ≈ 0.90 atol = 0.01
+    @test ContinuousBoyceIndex(nbins = 11, binwidth = 0.03)(ŷ, y) ≈ 0.77 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(; verbosity = 2, min =0.0, max = 2.0, nbins = 191)(ŷ, y),
+    )
+    # These two are identical because bins are dropped
+    @test cbi_dropped_bins == 
+        ContinuousBoyceIndex(; min = 0.0, max = 1.2, nbins = 111)(ŷ, y)
+    
+    # cbi is silent for verbosity 0
+    @test_logs ContinuousBoyceIndex(; verbosity = 0)(ŷ, y)
+end
+
 @testset "l2_check" begin
     d = Distributions.Normal()
     yhat = Union{Distributions.Sampleable,Missing}[d, d, missing]