Bayes Error Rate Estimation Tutorial#

Problem Statement#

For classification machine learning tasks, there is an inherent difficulty associated with signal to noise ratio in the images. One way of quantifying this difficulty is the Bayes Error Rate, or irreducable error.

DataEval has introduced a method of calculating this error rate that uses image embeddings.

When to use#

The BER metric should be used when you would like to measure the feasibility of a machine learning task. For example, if you have an operational accuracy requirement of 80%, and would like to know if this is feasibly achievable given the imagery.

What you will need#

  1. A set of image embeddings and their corresponding class labels. This requires training an autoencoder to compress the images.

  2. A python environment with the following packages installed:

    • dataeval[all]

    • torchvision

    • pytest

Setting up#

Let’s import the required libraries needed to set up a minimal working example

import matplotlib.pyplot as plt
import numpy as np
import tensorflow as tf
import torchvision.transforms.v2 as v2
from tensorflow.keras import Sequential
from tensorflow.keras.layers import Conv2D, Dense, Flatten, InputLayer
from tensorflow.nn import relu
from torchvision.datasets import MNIST

from dataeval.metrics.estimators import ber

tf.keras.utils.set_random_seed(408)

Loading in data#

Let’s start by loading in torchvision’s mnist dataset, then we will examine it

# Load in both the training and testing mnist dataset
train_ds = MNIST(root="./data/", train=True, download=True, transform=v2.ToImage())
test_ds = MNIST(root="./data/", train=False, download=True, transform=v2.ToImage())

# Split out the images and labels for each set
images, labels = train_ds.data.numpy(), train_ds.targets.numpy()
test_images, test_labels = test_ds.data.numpy(), test_ds.targets.numpy()
print("Number of training samples: ", len(images))
print("Image shape:", images[0].shape)
print("Label counts: ", np.unique(labels, return_counts=True))
Number of training samples:  60000
Image shape: (28, 28)
Label counts:  (array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9]), array([5923, 6742, 5958, 6131, 5842, 5421, 5918, 6265, 5851, 5949]))

To highlight the effects of modifying the dataset on its Bayes Error Rate, we will only include a subset of 6,000 images and their labels for digits 1, 4, and 9

images_split = {}
labels_split = {}

# Keep only 1, 4, and 9
for label in (1, 4, 9):
    subset_indices = np.where(labels == label)
    images_split[label] = images[subset_indices][:2000]
    labels_split[label] = labels[subset_indices][:2000]

images_subset = np.concatenate(list(images_split.values()))
labels_subset = np.concatenate(list(labels_split.values()))
print(images_subset.shape)
print(np.unique(labels_subset, return_counts=True))
(6000, 28, 28)
(array([1, 4, 9]), array([2000, 2000, 2000]))

We have taken a subset of the data that is only the digits 1, 4, and 9. The BER estimate requires 1 dimension, so we flatten the images during this step. This is ok since MNIST images are small, in practice we would need to do some dimension reduction (autoencoder) here.

# Flatten the images
images_flattened = images_subset.reshape((images_subset.shape[0], -1))
print("Dataset shape:", images_flattened.shape)
Dataset shape: (6000, 784)

We now have 6,000 flattened images of size 784. Next we can move on to evaluation of the dataset.

Evaluation#

Suppose we would like to build a classifier that differentiates between the handwritten digits 1, 4, and 9 with predetermined accuracy requirement of 99%.

We will use BER to check the feasibility of the task. As the images are small, we can simple use the flattened raw pixel intensities to calculate BER (no embedding necessary). Note: This will not be the case in general.

# Evaluate the BER metric for the MNIST data with digits 1, 4, 9.
# One minus the value of this metric gives our estimate of the upper bound on accuracy.
base_result = ber(images_flattened, labels_subset, method="MST")
print("The bayes error rate estimation:", base_result.ber)
The bayes error rate estimation: 0.025833333333333333

The estimate of the maximum achievable accuracy is one minus the BER estimate.

print("The maximum achievable accuracy:", (1 - base_result.ber) * 100)
The maximum achievable accuracy: 97.41666666666666

Results#

The maximum achievable accuracy on a dataset of 1, 4, and 9 is about 97.4%. This does not meet our requirement of 99% accuracy!

Modify dataset classification#

To address insufficient accuracy, lets modify the dataset to classify an image as “1” or “Not a 1”. By combining classes, we can hopefully achieve the desired level of attainable accuracy.

# Creates a binary mask where current label == 1 that can be used as the new labels
labels_merged = labels_subset == 1
print("New label counts:", np.unique(labels_merged, return_counts=True))
New label counts: (array([False,  True]), array([4000, 2000]))
# Evaluate the BER metric for the MNIST data with updated labels
new_result = ber(images_flattened, labels_merged, method="MST")
print("The bayes error rate estimation:", new_result)
The bayes error rate estimation: BEROutput(ber=0.005, ber_lower=0.0025062814466900174)

The estimate of the maximum achievable accuracy is one minus the BER estimate.

print("The maximum achievable accuracy:", 1 - new_result.ber)
The maximum achievable accuracy: 0.995

Results#

The maximum achievable accuracy on a dataset of 1 and not 1 (4, 9) is about 99.5%. This does meet our accuracy requirement.

By using BER to check for feasibility early on, we were able to reformulate the problem such that it is feasible under our specifications

Building a classifier#

We can now attempt to build a classifier that achieves this level of accuracy.

# Build a simple CNN for classifying MNIST images.
model = Sequential(
    [
        InputLayer(input_shape=(28, 28, 1)),
        Conv2D(
            64,
            4,
            strides=2,
            padding="same",
            activation=relu,
        ),
        Conv2D(
            128,
            4,
            strides=2,
            padding="same",
            activation=relu,
        ),
        Conv2D(
            512,
            4,
            strides=2,
            padding="same",
            activation=relu,
        ),
        Flatten(),
        Dense(2),
    ]
)

Since we are using a subset for training, we must also subset the testing data

test_indices = np.where((test_labels == 1) | (test_labels == 4) | (test_labels == 9))
test_images_subset = test_images[test_indices]
test_labels_subset = test_labels[test_indices]
test_labels_merged = test_labels_subset == 1

Train and test the model#

Now we train and test the model on the modified data

# Set up model hyperparameters
model.compile(
    optimizer=tf.keras.optimizers.Adam(learning_rate=1e-5),
    loss=tf.keras.losses.SparseCategoricalCrossentropy(from_logits=True),
    metrics=["accuracy"],
)

# Fitting a model may take a few minutes
history = model.fit(
    images_subset,
    labels_merged,
    epochs=90,
    batch_size=32,
    steps_per_epoch=1,
    validation_data=(test_images_subset, test_labels_merged),
    verbose=0,
)
loss, accuracy = model.evaluate(test_images_subset, test_labels_merged, verbose=1)
print(f"The model accuracy: {accuracy*100:0.2f}%")
 1/98 [..............................] - ETA: 2s - loss: 0.0157 - accuracy: 1.0000

15/98 [===>..........................] - ETA: 0s - loss: 0.0361 - accuracy: 0.9875

27/98 [=======>......................] - ETA: 0s - loss: 0.0479 - accuracy: 0.9850

40/98 [===========>..................] - ETA: 0s - loss: 0.0414 - accuracy: 0.9852

54/98 [===============>..............] - ETA: 0s - loss: 0.0385 - accuracy: 0.9878

65/98 [==================>...........] - ETA: 0s - loss: 0.0342 - accuracy: 0.9894

80/98 [=======================>......] - ETA: 0s - loss: 0.0330 - accuracy: 0.9902

96/98 [============================>.] - ETA: 0s - loss: 0.0298 - accuracy: 0.9912

98/98 [==============================] - 0s 4ms/step - loss: 0.0295 - accuracy: 0.9914
The model accuracy: 99.14%
plt.title("Model Accuracy")
plt.plot(range(60, 90), np.array(history.history["val_accuracy"])[60:], label="Classifier")
plt.hlines(
    y=1 - new_result.ber,
    colors=["red"],
    xmin=60,
    xmax=90,
    label="1-BER",
    linestyles="dashed",
)
plt.hlines(
    y=0.99,
    colors=["green"],
    xmin=60,
    xmax=90,
    label="Accuracy Requirement",
    linestyles="dashed",
)

plt.xticks(range(60, 91, 10))
plt.ylabel("Accuracy")
plt.xlabel("Epoch")
plt.legend(loc=4)
<matplotlib.legend.Legend at 0x7f7c6ff4a4d0>
../../_images/a42c0a8e08b3c036376696b1dcb084758195a448e1fd61620875363b641af15c.png

Results#

The model achieves an accuracy of 99.14% accuracy, exceeding the requirement of 99%.

The model accuracy does not quite approach the maximum achievable accuracy, meaning there are still improvements that can be made.