# Classification

This tutorial is based in parts on this blog post.

## Split Conformal Classification

We consider a simple binary classification problem. Let (*X**( i),Y*(

*i*)),ย

*i*โ=โ1,โ...,โ

*n*denote our feature-label pairs and let

*ฮผ*โ:โ๐ณโโฆโ๐ด denote the mapping from features to labels. For illustration purposes we will use the moons dataset ๐. Using

`MLJ.jl`

we first generate the data and split into into a training and test set:```
using MLJ
using Random
Random.seed!(123)
# Data:
X, y = make_moons(500; noise=0.15)
X = MLJ.table(convert.(Float32, MLJ.matrix(X)))
train, test = partition(eachindex(y), 0.8, shuffle=true)
```

Here we will use a specific case of CP called *split conformal prediction* which can then be summarized as follows:[1]

- Partition the training into a proper training set and a separate calibration set: ๐_(
*n*)โ=โ๐^(train)โ โชโ ๐^(cali). - Train the machine learning model on the proper training set:
*ฮผฬ**(*(*i*โโโ๐^(train))(*X**i*),*Y*_(*i*)). - Compute nonconformity scores, ๐ฎ, using the calibration data ๐^(cali) and the fitted model
*ฮผฬ*_(*i*โโโ๐^(train)). - For a user-specified desired coverage ratio (1โ
*ฮฑ*) compute the corresponding quantile,*qฬ*, of the empirical distribution of nonconformity scores, ๐ฎ. - For the given quantile and test sample
*X*_(test), form the corresponding conformal prediction set:

\[C(X_{\text{test}})=\{y:s(X_{\text{test}},y) \le \hat{q}\}\]

This is the default procedure used for classification and regression in `ConformalPrediction.jl`

.

Now letโs take this to our ๐ data. To illustrate the package functionality we will demonstrate the envisioned workflow. We first define our atomic machine learning model following standard `MLJ.jl`

conventions. Using `ConformalPrediction.jl`

we then wrap our atomic model in a conformal model using the standard API call `conformal_model(model::Supervised; kwargs...)`

. To train and predict from our conformal model we can then rely on the conventional `MLJ.jl`

procedure again. In particular, we wrap our conformal model in data (turning it into a machine) and then fit it to the training data. Finally, we use our machine to predict the label for a new test sample `Xtest`

:

```
# Model:
KNNClassifier = @load KNNClassifier pkg=NearestNeighborModels
model = KNNClassifier(;K=50)
# Training:
using ConformalPrediction
conf_model = conformal_model(model; coverage=.9)
mach = machine(conf_model, X, y)
fit!(mach, rows=train)
# Conformal Prediction:
Xtest = selectrows(X, test)
ytest = y[test]
yฬ = predict(mach, Xtest)
yฬ[1]
```

```
import NearestNeighborModels โ
UnivariateFinite{Multiclass{2}}(0=>0.94)
```

The final predictions are set-valued. While the softmax output remains unchanged for the `SimpleInductiveClassifier`

, the size of the prediction set depends on the chosen coverage rate, (1โ*ฮฑ*).

When specifying a coverage rate very close to one, the prediction set will typically include many (in some cases all) of the possible labels. Below, for example, both classes are included in the prediction set when setting the coverage rate equal to (1โ*ฮฑ*)=1.0. This is intuitive, since high coverage quite literally requires that the true label is covered by the prediction set with high probability.

```
conf_model = conformal_model(model; coverage=coverage, method=:simple_inductive)
mach = machine(conf_model, X, y)
fit!(mach, rows=train)
# Conformal Prediction:
Xtest = (x1=[1],x2=[0])
predict(mach, Xtest)[1]
```

`UnivariateFinite{Multiclass{2}}(0=>0.5, 1=>0.5)`

Conversely, for low coverage rates, prediction sets can also be empty. For a choice of (1โ*ฮฑ*)=0.1, for example, the prediction set for our test sample is empty. This is a bit difficult to think about intuitively and I have not yet come across a satisfactory, intuitive interpretation.[2] When the prediction set is empty, the `predict`

call currently returns `missing`

:

```
conf_model = conformal_model(model; coverage=coverage)
mach = machine(conf_model, X, y)
fit!(mach, rows=train)
# Conformal Prediction:
predict(mach, Xtest)[1]
```

`missing`

```
cov_ = .95
conf_model = conformal_model(model; coverage=cov_)
mach = machine(conf_model, X, y)
fit!(mach, rows=train)
Markdown.parse("""
The following chart shows the resulting predicted probabilities for ``y=1`` (left) and set size (right) for a choice of ``(1-\\alpha)``=$cov_.
""")
```

The following chart shows the resulting predicted probabilities for *y*โ=โ1 (left) and set size (right) for a choice of (1โ*ฮฑ*)=0.95.

```
using Plots
p_proba = contourf(mach.model, mach.fitresult, X, y)
p_set_size = contourf(mach.model, mach.fitresult, X, y; plot_set_size=true)
plot(p_proba, p_set_size, size=(800,250))
```

The animation below should provide some more intuition as to what exactly is happening here. It illustrates the effect of the chosen coverage rate on the predicted softmax output and the set size in the two-dimensional feature space. Contours are overlayed with the moon data points (including test data). The two samples highlighted in red, *X*โ and *X*โ, have been manually added for illustration purposes. Letโs look at these one by one.

Firstly, note that *X*โ (red cross) falls into a region of the domain that is characterized by high predictive uncertainty. It sits right at the bottom-right corner of our class-zero moon ๐ (orange), a region that is almost entirely enveloped by our class-one moon ๐ (green). For low coverage rates the prediction set for *X*โ is empty: on the left-hand side this is indicated by the missing contour for the softmax probability; on the right-hand side we can observe that the corresponding set size is indeed zero. For high coverage rates the prediction set includes both *y*โ=โ0 and *y*โ=โ1, indicative of the fact that the conformal classifier is uncertain about the true label.

With respect to *X*โ, we observe that while also sitting on the fringe of our class-zero moon, this sample populates a region that is not fully enveloped by data points from the opposite class. In this region, the underlying atomic classifier can be expected to be more certain about its predictions, but still not highly confident. How is this reflected by our corresponding conformal prediction sets?

```
Xtest_2 = (x1=[-0.5],x2=[0.25])
pฬ_2 = pdf(predict(mach, Xtest_2)[1], 0)
```

Well, for low coverage rates (roughly โ\<โ0.9) the conformal prediction set does not include *y*โ=โ0: the set size is zero (right panel). Only for higher coverage rates do we have *C*(*X*โ)โ=โ{0}: the coverage rate is high enough to include *y*โ=โ0, but the corresponding softmax probability is still fairly low. For example, for (1โ*ฮฑ*)โ=โ0.95 we have *pฬ*(*y*=0|*X*โ)โ=โ0.72.

These two examples illustrate an interesting point: for regions characterized by high predictive uncertainty, conformal prediction sets are typically empty (for low coverage) or large (for high coverage). While set-valued predictions may be something to get used to, this notion is overall intuitive.

```
# Setup
coverages = range(0.75,1.0,length=5)
n = 100
x1_range = range(extrema(X.x1)...,length=n)
x2_range = range(extrema(X.x2)...,length=n)
anim = @animate for coverage in coverages
conf_model = conformal_model(model; coverage=coverage)
mach = machine(conf_model, X, y)
fit!(mach, rows=train)
# Probabilities:
p1 = contourf(mach.model, mach.fitresult, X, y)
scatter!(p1, Xtest.x1, Xtest.x2, ms=6, c=:red, label="Xโ", shape=:cross, msw=6)
scatter!(p1, Xtest_2.x1, Xtest_2.x2, ms=6, c=:red, label="Xโ", shape=:diamond, msw=6)
p2 = contourf(mach.model, mach.fitresult, X, y; plot_set_size=true)
scatter!(p2, Xtest.x1, Xtest.x2, ms=6, c=:red, label="Xโ", shape=:cross, msw=6)
scatter!(p2, Xtest_2.x1, Xtest_2.x2, ms=6, c=:red, label="Xโ", shape=:diamond, msw=6)
plot(p1, p2, plot_title="(1-ฮฑ)=$(round(coverage,digits=2))", size=(800,300))
end
gif(anim, joinpath(www_path,"classification.gif"), fps=1)
```

## Adaptive Sets

Instead of using the simple approach, we can use adaptive prediction sets (Angelopoulos and Bates 2021):

```
conf_model = conformal_model(model; coverage=cov_, method=:adaptive_inductive)
mach = machine(conf_model, X, y)
fit!(mach, rows=train)
results[:adaptive_inductive] = mach
```

```
using Plots
p_proba = contourf(mach.model, mach.fitresult, X, y)
p_set_size = contourf(mach.model, mach.fitresult, X, y; plot_set_size=true)
plot(p_proba, p_set_size, size=(800,250))
```

## Evaluation

For evaluation of conformal predictors we follow Angelopoulos and Bates (2021) (Section 3). As a first step towards adaptiveness (adaptivity), the authors recommend to inspect the set size of conformal predictions. The chart below shows the interval width for the different methods along with the ground truth interval width:

```
plt_list = []
for (_mod, mach) in results
push!(plt_list, bar(mach.model, mach.fitresult, X; title=String(_mod)))
end
plot(plt_list..., size=(800,300))
```

We can also use specific metrics like **empirical coverage** and **size-stratified coverage** to check for correctness and adaptiveness, respectively. To this end, the package provides custom measures that are compatible with `MLJ.jl`

. In other words, we can evaluate model performance in true `MLJ.jl`

fashion (see here).

The code below runs the evaluation with respect to both metrics, `emp_coverage`

and `ssc`

for a single conformal machine:

```
_mod, mach = first(results)
_eval = evaluate!(
mach,
operation=predict,
measure=[emp_coverage, ssc]
)
# display(_eval)
println("Empirical coverage for $(_mod): $(round(_eval.measurement[1], digits=3))")
println("SSC for $(_mod): $(round(_eval.measurement[2], digits=3))")
```

```
Empirical coverage for adaptive_inductive: 0.962
SSC for adaptive_inductive: 0.962
```

## References

Angelopoulos, Anastasios N., and Stephen Bates. 2021. โA Gentle Introduction to Conformal Prediction and Distribution-Free Uncertainty Quantification.โ https://arxiv.org/abs/2107.07511.

[1] In other places split conformal prediction is sometimes referred to as *inductive* conformal prediction.

[2] Any thoughts/comments welcome!