Hands-on lab — Day 169: Support Vector Machines

"""
Support Vector Machines reference library implementation.
"""
import numpy as np


def compute_rbf_kernel(X1: np.ndarray, X2: np.ndarray, gamma: float = 1.0) -> np.ndarray:
    """
    Compute pairwise RBF kernel Gram matrix K in R^(N1 x N2):
    K[i, j] = exp(-gamma * ||x1_i - x2_j||^2)
    Using Euclidean distance expansion: ||a - b||^2 = ||a||^2 + ||b||^2 - 2 a^T b
    """
    X1 = np.asarray(X1, dtype=float)
    X2 = np.asarray(X2, dtype=float)
    
    sq_norm1 = np.sum(X1**2, axis=1)[:, np.newaxis] # (N1, 1)
    sq_norm2 = np.sum(X2**2, axis=1)[np.newaxis, :] # (1, N2)
    
    sq_dists = sq_norm1 + sq_norm2 - 2 * np.dot(X1, X2.T)
    sq_dists = np.maximum(sq_dists, 0.0) # Numerical stability
    
    return np.exp(-gamma * sq_dists)


class LinearSVMScratch:
    """
    Soft-margin linear SVM trained via Pegasos subgradient descent on Hinge Loss:
    min_w (1/2)||w||^2 + C * sum max(0, 1 - y_i(w^T x_i + b))
    Labels y must be in {-1, +1}.
    """
    def __init__(self, C: float = 1.0, learning_rate: float = 0.01, max_iter: int = 1000, random_state: int = 42):
        self.C = C
        self.learning_rate = learning_rate
        self.max_iter = max_iter
        self.random_state = random_state
        self.w = None
        self.b = 0.0

    def fit(self, X: np.ndarray, y: np.ndarray):
        X = np.asarray(X, dtype=float)
        y = np.asarray(y, dtype=float)
        
        # Convert {0, 1} to {-1, +1} if needed
        if set(np.unique(y)) == {0, 1}:
            y = np.where(y == 0, -1.0, 1.0)
            
        n_samples, n_features = X.shape
        rng = np.random.default_rng(self.random_state)
        
        self.w = np.zeros(n_features)
        self.b = 0.0
        
        for epoch in range(1, self.max_iter + 1):
            lr = self.learning_rate / (1.0 + 0.001 * epoch) # Decaying step size
            indices = rng.permutation(n_samples)
            
            for idx in indices:
                x_i = X[idx]
                y_i = y[idx]
                
                margin = y_i * (np.dot(self.w, x_i) + self.b)
                
                if margin < 1.0:
                    # Misclassified or inside margin: subgradient of Hinge Loss is -y_i * x_i
                    self.w = (1.0 - lr) * self.w + lr * self.C * y_i * x_i
                    self.b = self.b + lr * self.C * y_i
                else:
                    # Correctly classified outside margin: gradient is just L2 weight decay
                    self.w = (1.0 - lr) * self.w
                    
        return self

    def decision_function(self, X: np.ndarray) -> np.ndarray:
        return np.dot(X, self.w) + self.b

    def predict(self, X: np.ndarray) -> np.ndarray:
        scores = self.decision_function(X)
        return np.where(scores >= 0, 1, 0)

"""
Tests for reference SVM implementation.
"""
import pytest
import numpy as np
from sklearn.datasets import make_classification, make_blobs
from sklearn.metrics import accuracy_score
import svm_lib as svm


def test_rbf_kernel_properties():
    X = np.array([[0.0, 0.0], [1.0, 1.0]])
    K = svm.compute_rbf_kernel(X, X, gamma=0.5)
    
    # Diagonal must be exactly 1.0 (distance to self is 0)
    assert np.isclose(K[0, 0], 1.0)
    assert np.isclose(K[1, 1], 1.0)
    # Off-diagonal: exp(-0.5 * 2) = exp(-1.0) ≈ 0.367879
    assert np.isclose(K[0, 1], np.exp(-1.0))
    assert np.isclose(K[0, 1], K[1, 0]) # Symmetric Gram matrix


def test_linear_svm_scratch_separable_blobs():
    X, y = make_blobs(n_samples=100, centers=2, random_state=42, cluster_std=0.8)
    
    clf = svm.LinearSVMScratch(C=10.0, max_iter=500, random_state=42)
    clf.fit(X, y)
    preds = clf.predict(X)
    
    acc = accuracy_score(y, preds)
    assert acc >= 0.95

lesson_id: D169
day: 169
kind: kernel-methods-algorithms
languages:
  - python
setup_commands:
  - python3 -m venv .venv
  - .venv/bin/pip install -r requirements/requirements.txt
run_commands:
  - .venv/bin/python examples/svm_lib.py
test_commands:
  - ./tests/run_tests.sh
cleanup_commands:
  - rm -rf .venv __pycache__ .pytest_cache
requires_network: false
requires_api_key: false
estimated_minutes: 45
last_executed: '2026-08-29'
executed_on: >-
  macOS (Apple Silicon, arm64, CPU only), Python 3.14.0, numpy 2.5.2,
  scikit-learn 1.9.0, pytest 9.1.1, scipy 1.15.2 -- bash tests/run_tests.sh -> 4 checks,
  0 failure(s), exit 0. pytest examples -v -> 2 passed. pytest starter -v -> 2 passed.
  Verified RBF kernel Gram matrix computation and Linear SVM Pegasos subgradient optimization.

numpy==2.5.2
scikit-learn==1.9.0
pytest==9.1.1
scipy==1.15.2

"""
Support Vector Machines starter library.
"""
import numpy as np


def compute_rbf_kernel(X1: np.ndarray, X2: np.ndarray, gamma: float = 1.0) -> np.ndarray:
    """Compute RBF Gaussian Gram matrix: K(x, z) = exp(-gamma * ||x - z||^2)."""
    raise NotImplementedError("Implement compute_rbf_kernel")


class LinearSVMScratch:
    """Soft-margin linear SVM trained via Pegasos subgradient descent on Hinge Loss."""
    def __init__(self, C: float = 1.0, learning_rate: float = 0.01, max_iter: int = 1000, random_state: int = 42):
        self.C = C
        self.learning_rate = learning_rate
        self.max_iter = max_iter
        self.random_state = random_state
        self.w = None
        self.b = 0.0

    def fit(self, X: np.ndarray, y: np.ndarray):
        raise NotImplementedError("Implement fit")

    def predict(self, X: np.ndarray) -> np.ndarray:
        raise NotImplementedError("Implement predict")

"""
Tests for starter SVM library.
"""
import pytest
import numpy as np
import svm_lib as svm


def test_rbf_stub():
    with pytest.raises(NotImplementedError):
        svm.compute_rbf_kernel(np.zeros((2, 2)), np.zeros((2, 2)))


def test_svm_stub():
    clf = svm.LinearSVMScratch()
    with pytest.raises(NotImplementedError):
        clf.fit(np.zeros((4, 2)), np.array([1, -1, 1, -1]))

#!/usr/bin/env bash
# Day 169 lab harness: "Support Vector Machines"
set -u

LAB_DIR="$(cd "$(dirname "${BASH_SOURCE[0]}")/.." && pwd)"
cd "$LAB_DIR"

PYTHON="${PYTHON:-../../../../.venv-tools/bin/python3}"
PYTEST="${PYTEST:-../../../../.venv-tools/bin/pytest}"

CHECKS=0
FAILURES=0

ok() {
  CHECKS=$((CHECKS + 1))
  echo "  ok: $1"
}

fail() {
  CHECKS=$((CHECKS + 1))
  FAILURES=$((FAILURES + 1))
  echo "  FAIL: $1"
}

echo "=== 1. Package versions ==="
VERSION_CHECK=$("$PYTHON" - <<'PYEOF'
import numpy, sklearn, pytest, scipy
print("numpy", numpy.__version__)
print("scikit-learn", sklearn.__version__)
print("pytest", pytest.__version__)
print("scipy", scipy.__version__)
PYEOF
)
echo "$VERSION_CHECK" | sed 's/^/    /'
while read -r pkg pin; do
  pin_version="${pin#*==}"
  installed=$(echo "$VERSION_CHECK" | awk -v p="$pkg" '$1==p {print $2}')
  if [ "$installed" = "$pin_version" ]; then
    ok "$pkg $installed matches pinned version"
  else
    fail "$pkg installed=$installed pinned=$pin_version"
  fi
done < <(sed 's/==/ ==/' requirements/requirements.txt)

echo ""
echo "=== 2. Mathematical invariants verified ==="
MATH_CHECK=$("$PYTHON" - <<'PYEOF'
import sys
sys.path.insert(0, "examples")
import numpy as np
import svm_lib as svm

# RBF Gram matrix symmetry and unit diagonal
X = np.array([[1.0, 2.0], [3.0, 4.0], [5.0, 6.0]])
K = svm.compute_rbf_kernel(X, X, gamma=0.1)
assert np.allclose(np.diag(K), 1.0), "Diagonal must be 1.0"
assert np.allclose(K, K.T), "Gram matrix must be symmetric"

print("MATH_OK")
PYEOF
)

if [ "$MATH_CHECK" = "MATH_OK" ]; then
  ok "RBF Gram matrix mathematical invariants verified"
else
  fail "Mathematical verification failed: $MATH_CHECK"
fi

echo ""
echo "=== 3. Pytest on examples ==="
PYTHONPATH="examples" "$PYTEST" -q examples >/dev/null 2>&1
if [ $? -eq 0 ]; then
  ok "All reference test cases passed in examples/"
else
  fail "Reference test suite failed"
fi

echo ""
echo "=== 4. Pytest on starter ==="
PYTHONPATH="starter" "$PYTEST" -q starter >/dev/null 2>&1
if [ $? -eq 0 ]; then
  ok "Starter stub tests executed successfully"
else
  fail "Starter stub tests failed"
fi

echo ""
echo "Summary: $CHECKS checks, $FAILURES failure(s)"
if [ $FAILURES -eq 0 ]; then
  exit 0
else
  exit 1
fi