Python Implementation of Convolutional Neural Networks (CNN)

The Story That Explains Why CNN Was Born

In 1989, Yann LeCun introduced the concept of CNNs. In 1998 he built LeNet-5 — a CNN that read handwritten digits on cheques for banks. By 2012, AlexNet won the ImageNet competition with a 15.3% error rate, crushing all traditional methods. CNNs had arrived. Today, they power face recognition on your phone, medical imaging, self-driving cars, and satellite image analysis.

CNN Architecture — The Big Picture

A CNN is not one single operation — it is a pipeline of specialised layers. Think of it as a factory assembly line where each station adds value.

The Convolution Operation — The Heart of CNN

Mathematically, a 2D convolution applies a small matrix (the filter/kernel) to patches of the input. The filter slides across the image with a given stride, performing an element-wise multiplication at each position, then summing the result into a single number.

🔬 Numerical Example: A 5×5 Image, 3×3 Filter

Padding and Stride — The Controls

Pooling Layers — Compressing Without Losing Soul

After convolution, we have feature maps that are large and potentially redundant. Pooling aggregates each small region into one value. It makes the representation more compact, reduces parameters, and introduces spatial invariance — a feature detected slightly off-position still fires.

Activation Functions — Why Non-Linearity Matters

Batch Normalisation and Dropout — The Regulators

How CNNs Learn — Backpropagation Through Convolutions

Full Python Implementation — Step by Step

Now we build a complete CNN from scratch using TensorFlow/Keras. Each syntax choice is explained in detail. We will structure the code into clean, logical steps.

Step 1 — Import Libraries

# Core libraries — understand each import's role
import numpy as np                          # Array operations, data manipulation
import matplotlib.pyplot as plt            # Visualisation of images, loss curves
import tensorflow as tf                     # Deep learning framework

# Keras API — tf.keras is the high-level API built into TensorFlow 2.x
# 'from tensorflow import keras' is equivalent
from tensorflow.keras import layers, models, optimizers, regularizers
from tensorflow.keras.callbacks import EarlyStopping, ReduceLROnPlateau, ModelCheckpoint
from tensorflow.keras.preprocessing.image import ImageDataGenerator

# sklearn utilities for evaluation
from sklearn.metrics import classification_report, confusion_matrix
import seaborn as sns

# Set seeds for reproducibility — without this, results differ each run
np.random.seed(42)
tf.random.set_seed(42)

Step 2 — Load and Prepare the CIFAR-10 Dataset

# CIFAR-10: 60,000 colour images, 32×32 pixels, 10 classes
# Classes: airplane, automobile, bird, cat, deer, dog, frog, horse, ship, truck

(X_train, y_train), (X_test, y_test) = tf.keras.datasets.cifar10.load_data()

# Shapes after loading:
print(f"X_train shape: {X_train.shape}")  # (50000, 32, 32, 3)
print(f"X_test shape:  {X_test.shape}")   # (10000, 32, 32, 3)
print(f"y_train shape: {y_train.shape}")  # (50000, 1)

# ── WHY NORMALISE TO [0,1]? ──────────────────────────────────────────────
# Pixel values are uint8: 0 to 255. Neural nets train poorly with large inputs:
#   • Gradients become unstable (very large or very small)
#   • Weight initialisations are calibrated for ~unit-scale inputs
# Dividing by 255.0 maps values to [0.0, 1.0] and makes training stable.

X_train = X_train.astype('float32') / 255.0
X_test  = X_test.astype('float32')  / 255.0

# ── WHY float32 and NOT float64? ────────────────────────────────────────
# GPUs are optimised for float32. float64 gives no benefit for neural net
# training but uses 2× memory and is 2-4× slower on GPU hardware.

# ── ONE-HOT ENCODING THE LABELS ─────────────────────────────────────────
# y_train contains integers 0-9. Softmax outputs 10 probabilities.
# We need labels as 10-element vectors for categorical_crossentropy.
# E.g.: class 3 → [0, 0, 0, 1, 0, 0, 0, 0, 0, 0]

y_train_cat = tf.keras.utils.to_categorical(y_train, 10)
y_test_cat  = tf.keras.utils.to_categorical(y_test,  10)

# Class names for display
class_names = ['airplane', 'automobile', 'bird', 'cat', 'deer',
               'dog', 'frog', 'horse', 'ship', 'truck']

Step 3 — Data Augmentation

# ── WHY DATA AUGMENTATION? ───────────────────────────────────────────────
# CNNs can memorise 50,000 training images perfectly (overfit).
# Augmentation artificially expands the dataset by creating variations
# of each image. This forces the network to learn features that are
# INVARIANT to flips, shifts, and rotations — more generalisable.

# ImageDataGenerator applies transforms RANDOMLY on the fly during training.
# The original images on disk are never modified.

datagen = ImageDataGenerator(
    rotation_range=15,      # Rotate image up to ±15 degrees randomly
                              # WHY 15°: CIFAR-10 objects are mostly upright; ±15° is realistic
    width_shift_range=0.1,   # Shift horizontally by up to 10% of width
    height_shift_range=0.1, # Shift vertically by up to 10% of height
                              # WHY SHIFTS: objects aren't always perfectly centred
    horizontal_flip=True,    # Mirror image left-right 50% of the time
                              # WHY: A cat is still a cat when mirrored. NOT used for digits/text.
    zoom_range=0.1,          # Zoom in/out by up to 10%
    fill_mode='nearest'      # How to fill pixels created by shifts: copy nearest edge pixel
                              # Alternatives: 'reflect', 'wrap', 'constant'
)

# Fit the generator on training data (computes internal statistics if needed)
datagen.fit(X_train)

Step 4 — Build the CNN Model

# ── SEQUENTIAL API vs FUNCTIONAL API ────────────────────────────────────
# Sequential: simple stack of layers — suitable for most standard CNNs.
# Functional API: used when you need skip connections, multiple inputs/outputs.
# We use Sequential here for clarity.

def build_cnn(input_shape=(32, 32, 3), num_classes=10):
    model = models.Sequential([

        # ── BLOCK 1: First Convolutional Block ──────────────────────────────
        # Conv2D(32, (3,3), ...)
        #   32        = number of filters. Produces 32 feature maps.
        #              WHY 32 FIRST?: Start small — capture basic features.
        #              Larger first layers waste computation on simple features.
        #   (3,3)     = kernel size. 3×3 is the gold standard (proven by VGG).
        #              WHY 3×3?: Two 3×3 convs have same receptive field as one 5×5
        #              but with fewer parameters and more non-linearity.
        #   padding='same' = preserve spatial dimensions (32×32 → 32×32)
        #   activation='relu' = apply ReLU after convolution (most common choice)
        #   input_shape = only needed on FIRST layer. TensorFlow infers rest.

        layers.Conv2D(32, (3, 3), padding='same', activation='relu',
                       input_shape=input_shape,
                       kernel_regularizer=regularizers.l2(1e-4)),
                       # kernel_regularizer=l2(1e-4): add tiny L2 penalty on weights
                       # 1e-4 is a mild penalty — prevents extreme weight values

        layers.BatchNormalization(),  # Normalise after conv, before or after relu — debate exists.
                                       # Original paper: before relu. Modern practice: after.
                                       # We use after activation to follow modern convention.

        layers.Conv2D(32, (3, 3), padding='same', activation='relu',
                       kernel_regularizer=regularizers.l2(1e-4)),
        # WHY TWO CONV LAYERS BEFORE POOLING?
        # Stacking 2 conv layers increases the effective receptive field
        # (two 3×3 layers = 5×5 receptive field) while using fewer parameters
        # than a single 5×5 layer. More non-linearity = richer representations.

        layers.MaxPooling2D(2, 2),  # Halve spatial dims: 32×32 → 16×16
                                      # (2,2) = pool window, default stride = pool size
        layers.Dropout(0.25),        # Drop 25% of units. Lighter dropout after conv (25% vs 50%)
                                      # because conv layers have fewer params and need less regularisation.

        # ── BLOCK 2: Deeper, More Filters ───────────────────────────────────
        # 64 filters: we DOUBLE the number after pooling. WHY?
        # Pooling halved the spatial resolution (information loss).
        # Doubling filters compensates by increasing the channel depth —
        # we trade spatial richness for semantic richness.

        layers.Conv2D(64, (3, 3), padding='same', activation='relu',
                       kernel_regularizer=regularizers.l2(1e-4)),
        layers.BatchNormalization(),
        layers.Conv2D(64, (3, 3), padding='same', activation='relu',
                       kernel_regularizer=regularizers.l2(1e-4)),
        layers.MaxPooling2D(2, 2),  # 16×16 → 8×8
        layers.Dropout(0.25),

        # ── BLOCK 3: Deepest Block ──────────────────────────────────────────
        layers.Conv2D(128, (3, 3), padding='same', activation='relu',
                        kernel_regularizer=regularizers.l2(1e-4)),
        layers.BatchNormalization(),
        layers.MaxPooling2D(2, 2),  # 8×8 → 4×4
        layers.Dropout(0.25),

        # ── CLASSIFICATION HEAD ─────────────────────────────────────────────
        # Flatten converts 3D tensor (4, 4, 128) → 1D vector (2048,)
        # WHY FLATTEN NOT GlobalAveragePooling2D?
        # Flatten + Dense gives more parameters (more capacity) but risks overfit.
        # GAP reduces spatial dims to 1×1, losing less info, fewer params.
        # For CIFAR-10 (small images), Flatten works well.

        layers.Flatten(),

        # Dense(512): fully connected layer with 512 neurons
        # WHY 512?: Rule of thumb — start large and narrow down. Gives the
        # network capacity to combine features. Reduce toward output size.
        layers.Dense(512, activation='relu',
                      kernel_regularizer=regularizers.l2(1e-4)),
        layers.BatchNormalization(),
        layers.Dropout(0.5),   # Heavier dropout on Dense layers (50%) — more overfit risk

        # Final output layer — 10 neurons for 10 classes
        # activation='softmax': converts 10 raw scores to probabilities summing to 1
        # WHY NOT 'sigmoid' here?: sigmoid outputs independent probabilities (0-1) each.
        # softmax creates a COMPETITION between classes — exactly what we want for
        # mutually exclusive classification (one image = one class).
        layers.Dense(10, activation='softmax')
    ])
    return model

model = build_cnn()
model.summary()

Step 5 — Compile the Model

# ── COMPILE: tells Keras HOW to train ────────────────────────────────────

# optimizer='adam': Adaptive Moment Estimation
#   WHY ADAM over plain SGD?
#   Adam maintains per-parameter learning rates, adapted based on:
#     - momentum (m): exponential average of gradients (β1=0.9 default)
#     - velocity  (v): exponential average of squared gradients (β2=0.999)
#   Result: fast convergence, robust to hyperparameter choices.
#   learning_rate=0.001 is Adam's default and usually a good start.

# loss='categorical_crossentropy':
#   WHY THIS LOSS for classification?
#   It measures the distance between the predicted probability distribution
#   and the true one-hot distribution. Mathematically: -Σ y_true * log(y_pred).
#   Penalises confident wrong predictions VERY harshly (log(0) → -∞).
#   Use 'sparse_categorical_crossentropy' if labels are integers (no to_categorical).

# metrics=['accuracy']: tracks accuracy during training for human monitoring.
#   Note: accuracy is NOT the loss — just a readable metric for logging.

model.compile(
    optimizer=optimizers.Adam(learning_rate=0.001),
    loss='categorical_crossentropy',
    metrics=['accuracy']
)

Step 6 — Define Callbacks

# ── CALLBACKS: actions taken at end of each epoch ────────────────────────

# 1. EarlyStopping:
#    Monitor validation loss. If it doesn't improve for 'patience' epochs,
#    stop training. Prevents wasting compute on a model that's overfit.
#    restore_best_weights=True: roll back to the epoch with best val_loss.
#    WHY monitor val_loss not val_accuracy?: accuracy is coarse (steps of 1/N),
#    loss is continuous — more sensitive signal for early stopping.

early_stop = EarlyStopping(
    monitor='val_loss',
    patience=10,
    restore_best_weights=True,
    verbose=1
)

# 2. ReduceLROnPlateau:
#    When val_loss stops improving, reduce learning rate by factor 'factor'.
#    WHY?: Adam can get stuck in a local minimum. Reducing LR lets it take
#    smaller steps and potentially escape or fine-tune more precisely.
#    factor=0.5 means LR is halved. min_lr prevents it from reaching zero.

reduce_lr = ReduceLROnPlateau(
    monitor='val_loss',
    factor=0.5,
    patience=5,
    min_lr=1e-6,
    verbose=1
)

# 3. ModelCheckpoint:
#    Save the model weights whenever val_accuracy improves.
#    WHY save_best_only=True?: training can degrade after a peak epoch.
#    This guarantees we always have the best version on disk.

checkpoint = ModelCheckpoint(
    filepath='best_cnn_model.h5',
    monitor='val_accuracy',
    save_best_only=True,
    verbose=1
)

callbacks = [early_stop, reduce_lr, checkpoint]

Step 7 — Train the Model

# ── model.fit() — the training loop ──────────────────────────────────────
#
# datagen.flow(X_train, y_train_cat, batch_size=64):
#   Creates a generator that yields augmented batches of 64 images.
#   WHY batch_size=64?:
#     - Too small (8): noisy gradient estimates, slow progress
#     - Too large (512): smooth but may get stuck in sharp minima
#     - 32-128 is the practical sweet spot for CNNs on standard datasets
#
# steps_per_epoch = len(X_train) // batch_size:
#   How many batches = one epoch. 50,000 / 64 ≈ 781 steps per epoch.
#   WHY specify this?: When using a generator, Keras can't automatically
#   know when one epoch ends.
#
# validation_data=(X_test, y_test_cat):
#   After each epoch, evaluate on test set. Never used for training.
#   WHY not augment validation data?: We want to evaluate on REAL images,
#   not augmented ones — augmentation is only for training diversity.
#
# epochs=100: maximum training epochs (EarlyStopping will likely stop before)

history = model.fit(
    datagen.flow(X_train, y_train_cat, batch_size=64),
    steps_per_epoch=len(X_train) // 64,
    epochs=100,
    validation_data=(X_test, y_test_cat),
    callbacks=callbacks,
    verbose=1
)

Step 8 — Evaluate and Visualise Results

# ── EVALUATION ───────────────────────────────────────────────────────────
test_loss, test_acc = model.evaluate(X_test, y_test_cat, verbose=0)
print(f"Test Accuracy: {test_acc:.4f}")
print(f"Test Loss:     {test_loss:.4f}")

# ── CONFUSION MATRIX ─────────────────────────────────────────────────────
# model.predict returns probability arrays. argmax gives the predicted class.
# WHY argmax? Softmax output is [0.02, 0.01, 0.85, ...]; argmax = index of max.

y_pred_probs = model.predict(X_test)
y_pred = np.argmax(y_pred_probs, axis=1)    # axis=1: argmax across 10 classes
y_true = np.argmax(y_test_cat, axis=1)     # convert one-hot back to integers

cm = confusion_matrix(y_true, y_pred)

# ── VISUALISE TRAINING CURVES ────────────────────────────────────────────
fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# Loss curve
axes[0].plot(history.history['loss'],     label='Train Loss', color='#60a5fa')
axes[0].plot(history.history['val_loss'], label='Val Loss',   color='#f87171')
axes[0].set_title('Loss Curves')
axes[0].legend()

# Accuracy curve
axes[1].plot(history.history['accuracy'],     label='Train Acc', color='#34d399')
axes[1].plot(history.history['val_accuracy'], label='Val Acc',   color='#a78bfa')
axes[1].set_title('Accuracy Curves')
axes[1].legend()

plt.tight_layout()
plt.show()

# Full per-class report
print(classification_report(y_true, y_pred, target_names=class_names))

Case Study — Pneumonia Detection from Chest X-Rays

Case Study — Full Code

import tensorflow as tf
from tensorflow.keras import layers, models
from tensorflow.keras.preprocessing.image import ImageDataGenerator
from tensorflow.keras.applications import VGG16
from tensorflow.keras.callbacks import EarlyStopping, ReduceLROnPlateau
import numpy as np
from sklearn.metrics import classification_report, roc_auc_score

# ── DATASET PATHS ─────────────────────────────────────────────────────────
# Dataset: kaggle datasets download -d paultimothymooney/chest-xray-pneumonia
# Structure:
#   chest_xray/
#     train/ NORMAL/ PNEUMONIA/
#     val/   NORMAL/ PNEUMONIA/
#     test/  NORMAL/ PNEUMONIA/

TRAIN_DIR = 'chest_xray/train'
VAL_DIR   = 'chest_xray/val'
TEST_DIR  = 'chest_xray/test'

IMG_SIZE  = (224, 224)   # WHY 224×224? VGG16 and many pretrained models
                           # expect 224×224. We use transfer learning below.
BATCH     = 32

# ── DATA GENERATORS ──────────────────────────────────────────────────────
# Training: aggressive augmentation — X-rays can be taken at slight angles,
#           different zoom levels, different patient orientations.
# Validation/Test: NO augmentation — evaluate on real unmodified images.
# rescale=1/255: normalise pixel values to [0,1] for all generators.

train_gen = ImageDataGenerator(
    rescale=1./255,
    rotation_range=10,          # X-rays rarely tilted more than ±10°
    zoom_range=0.1,
    horizontal_flip=True,        # Lungs are bilaterally symmetric
    shear_range=0.1,             # Small shear — simulates patient lean
    brightness_range=[0.8, 1.2] # X-ray exposure can vary between machines
)

val_test_gen = ImageDataGenerator(rescale=1./255)

# flow_from_directory: reads images from folder structure
# class_mode='binary': returns 0 (NORMAL) or 1 (PNEUMONIA) labels
# WHY binary and not categorical?: Only 2 classes. Sigmoid output.
# target_size: resizes all images to IMG_SIZE on the fly.

train_data = train_gen.flow_from_directory(
    TRAIN_DIR, target_size=IMG_SIZE, batch_size=BATCH, class_mode='binary')

val_data = val_test_gen.flow_from_directory(
    VAL_DIR, target_size=IMG_SIZE, batch_size=BATCH, class_mode='binary',
    shuffle=False)   # WHY shuffle=False on val/test? We need predictions to
                      # align with true labels in the correct order for metrics.

test_data = val_test_gen.flow_from_directory(
    TEST_DIR, target_size=IMG_SIZE, batch_size=BATCH, class_mode='binary',
    shuffle=False)

# ── CLASS IMBALANCE CHECK ─────────────────────────────────────────────────
# Kaggle dataset: ~3875 PNEUMONIA, ~1341 NORMAL → imbalanced 3:1 ratio.
# If we ignore this, the model could predict "always pneumonia" and get 74% accuracy!
# We compute class weights to make the loss penalise minority class more.

n_normal    = 1341
n_pneumonia = 3875
total       = n_normal + n_pneumonia

class_weight = {
    0: total / (2 * n_normal),    # ~1.85 — upweight Normal class
    1: total / (2 * n_pneumonia)   # ~0.64 — slightly downweight Pneumonia
}
print("Class weights:", class_weight)

# ── TRANSFER LEARNING WITH VGG16 ─────────────────────────────────────────
# WHY TRANSFER LEARNING?
# We have ~5,000 training images. Training a deep CNN from scratch on this
# would overfit badly. VGG16 was pretrained on ImageNet (1.2M images, 1000 classes).
# Its early layers already detect edges, textures, shapes — useful for X-rays too.
# We FREEZE the pretrained layers, add our own head, and only train the head first.

# weights='imagenet': download pretrained ImageNet weights
# include_top=False: exclude VGG16's original Dense classification head
# input_shape=(224, 224, 3): X-rays are grayscale but VGG16 expects 3 channels.
#   We'll use the same image 3× (convert grayscale to 3-channel).
# WHY NOT grayscale directly?: Pretrained weights expect 3 channels.
#   Simple fix: repeat grayscale channel 3 times.

base_model = VGG16(weights='imagenet', include_top=False,
                    input_shape=(*IMG_SIZE, 3))

# Freeze ALL base model layers — only our custom head will be trained initially
base_model.trainable = False

# ── BUILD CUSTOM HEAD ON TOP OF VGG16 ───────────────────────────────────
model = models.Sequential([
    base_model,                             # Frozen VGG16 feature extractor

    # GlobalAveragePooling2D: replaces Flatten.
    # Output of VGG16 (with 224×224 input) is (7, 7, 512).
    # GAP2D reduces this to a 512-dim vector by averaging each feature map.
    # WHY GAP instead of Flatten?: 7×7×512 = 25,088. Flatten → Dense would be
    # 25,088 × 256 = 6.4 million params for one layer. GAP gives 512 × 256 = 131k.
    # GAP is regularising: inherently reduces overfitting.
    layers.GlobalAveragePooling2D(),

    layers.Dense(256, activation='relu'),
    layers.BatchNormalization(),
    layers.Dropout(0.5),

    layers.Dense(64, activation='relu'),
    layers.Dropout(0.3),

    # Binary output: 1 neuron, sigmoid activation
    # sigmoid outputs probability of Pneumonia (class 1): value in [0,1]
    # Default decision threshold: >0.5 → Pneumonia, ≤0.5 → Normal
    # We will LOWER this threshold to 0.3 to reduce False Negatives.
    layers.Dense(1, activation='sigmoid')
])

# ── COMPILE: binary_crossentropy for 2-class sigmoid output ─────────────
# WHY binary_crossentropy not categorical?:
#   binary_crossentropy = -[y*log(p) + (1-y)*log(1-p)]
#   Designed for a single sigmoid output. categorical_crossentropy is
#   for multi-class softmax with one-hot labels.
# metrics=['accuracy', AUC]: AUC-ROC is critical in medical imaging —
#   it measures discrimination across ALL thresholds, not just 0.5.

model.compile(
    optimizer='adam',
    loss='binary_crossentropy',
    metrics=['accuracy', tf.keras.metrics.AUC(name='auc')]
)

# ── PHASE 1: Train head only (base frozen) ───────────────────────────────
history1 = model.fit(
    train_data,
    epochs=15,
    validation_data=val_data,
    class_weight=class_weight,     # Apply imbalance correction
    callbacks=[
        EarlyStopping(monitor='val_auc', patience=5, restore_best_weights=True,
                      mode='max'),  # mode='max': AUC should INCREASE
        ReduceLROnPlateau(monitor='val_loss', factor=0.5, patience=3)
    ]
)

# ── PHASE 2: Fine-tune — unfreeze last few VGG16 blocks ─────────────────
# Now that our head is trained, carefully unfreeze the last convolutional
# block of VGG16 (block5) and fine-tune with a MUCH smaller learning rate.
# WHY small LR for fine-tuning?: The pretrained weights are good starting points.
# A large LR would destroy them ("catastrophic forgetting").
# We want tiny adjustments to specialise for X-rays, not relearn from scratch.

base_model.trainable = True

# Only unfreeze layers from block5 onward (last 4 conv layers)
for layer in base_model.layers:
    if 'block5' in layer.name:
        layer.trainable = True
    else:
        layer.trainable = False

# Recompile with 10× smaller learning rate for fine-tuning
model.compile(
    optimizer=optimizers.Adam(learning_rate=1e-5),  # 0.001 → 0.00001
    loss='binary_crossentropy',
    metrics=['accuracy', tf.keras.metrics.AUC(name='auc')]
)

history2 = model.fit(
    train_data,
    epochs=20,
    validation_data=val_data,
    class_weight=class_weight,
    callbacks=[
        EarlyStopping(monitor='val_auc', patience=7, restore_best_weights=True, mode='max'),
        ReduceLROnPlateau(monitor='val_loss', factor=0.5, patience=4)
    ]
)

# ── EVALUATION WITH LOWERED THRESHOLD ────────────────────────────────────
# Default threshold is 0.5. In medical context, we lower to 0.3:
# → Catches more True Positives (reduce missed pneumonia cases)
# → Accepts more False Positives (extra doctor reviews — acceptable trade-off)

y_pred_prob = model.predict(test_data, verbose=0).flatten()
y_true_test = test_data.classes

threshold = 0.3   # Medical decision: sensitivity over specificity
y_pred_class = (y_pred_prob >= threshold).astype(int)

auc_score = roc_auc_score(y_true_test, y_pred_prob)
print(f"\nAUC-ROC Score: {auc_score:.4f}")
print(classification_report(y_true_test, y_pred_class,
                             target_names=['Normal', 'Pneumonia']))

Diagnosing Your CNN — Reading the Loss Curves

Golden Rules for Building CNNs

CNN vs Other Architectures — When to Use What