Variants: Batch, Stochastic, Mini-batch

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64

import numpy as np
 
def batch_gradient_descent(X, y, theta_init, learning_rate=0.01, 
                            n_epochs=100, tol=1e-6):
    """
    Full batch gradient descent for linear regression.
    
    Args:
        X: Feature matrix (n_samples, n_features)
        y: Target vector (n_samples,)
        theta_init: Initial parameters
        learning_rate: Step size α
        n_epochs: Maximum number of passes through data
        tol: Convergence tolerance
    
    Returns:
        theta: Learned parameters
        history: List of (loss, grad_norm) per epoch
    """
    n_samples = X.shape[0]
    theta = theta_init.copy()
    history = []
    
    for epoch in range(n_epochs):
        # Forward pass: compute predictions for ALL samples
        y_pred = X @ theta
        
        # Compute exact gradient over entire dataset
        residuals = y_pred - y
        gradient = (X.T @ residuals) / n_samples  # Average gradient
        
        # Compute loss for monitoring
        loss = np.mean(residuals ** 2) / 2
        grad_norm = np.linalg.norm(gradient)
        history.append({'epoch': epoch, 'loss': loss, 'grad_norm': grad_norm})
        
        # Check convergence
        if grad_norm < tol:
            print(f"Converged at epoch {epoch}")
            break
        
        # Update: one step per epoch
        theta = theta - learning_rate * gradient
    
    return theta, history
 
# Example usage
np.random.seed(42)
n_samples = 1000
n_features = 5
 
# Generate data
X = np.random.randn(n_samples, n_features)
true_theta = np.array([1, -2, 3, -4, 5])
y = X @ true_theta + 0.1 * np.random.randn(n_samples)
 
# Train
theta_init = np.zeros(n_features)
theta, history = batch_gradient_descent(X, y, theta_init, 
                                         learning_rate=0.1, n_epochs=100)
 
print(f"True parameters: {true_theta}")
print(f"Learned parameters: {np.round(theta, 3)}")
print(f"Final loss: {history[-1]['loss']:.6f}")

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91

import numpy as np
 
def stochastic_gradient_descent(X, y, theta_init, learning_rate=0.01,
                                  n_epochs=10, lr_decay=None):
    """
    Pure Stochastic Gradient Descent (batch size = 1).
    
    Args:
        X: Feature matrix (n_samples, n_features)
        y: Target vector (n_samples,)
        theta_init: Initial parameters
        learning_rate: Initial step size α
        n_epochs: Number of passes through entire dataset
        lr_decay: Optional learning rate schedule ('1/t', 'sqrt')
    
    Returns:
        theta: Learned parameters
        history: List of loss values (sampled periodically)
    """
    n_samples, n_features = X.shape
    theta = theta_init.copy()
    history = []
    step = 0
    
    for epoch in range(n_epochs):
        # Shuffle data at start of each epoch
        indices = np.random.permutation(n_samples)
        epoch_losses = []
        
        for i in indices:
            step += 1
            
            # Get single sample
            x_i = X[i:i+1]  # Shape (1, n_features)
            y_i = y[i:i+1]  # Shape (1,)
            
            # Forward pass: prediction for ONE sample
            y_pred = x_i @ theta
            
            # Compute gradient for ONE sample (unbiased estimator)
            residual = y_pred - y_i
            gradient = x_i.T.flatten() * residual[0]  # No averaging needed
            
            # Learning rate decay
            if lr_decay == '1/t':
                alpha = learning_rate / (1 + 0.01 * step)
            elif lr_decay == 'sqrt':
                alpha = learning_rate / np.sqrt(1 + 0.01 * step)
            else:
                alpha = learning_rate
            
            # Update with noisy gradient
            theta = theta - alpha * gradient
            
            # Track loss (for subset of steps to avoid overhead)
            if step % 100 == 0:
                loss = np.mean((X @ theta - y) ** 2) / 2
                epoch_losses.append(loss)
        
        # Epoch summary
        if epoch_losses:
            history.append({
                'epoch': epoch,
                'loss': np.mean(epoch_losses),
                'final_lr': alpha
            })
    
    return theta, history
 
# Example: Compare SGD with and without learning rate decay
np.random.seed(42)
n_samples = 1000
X = np.random.randn(n_samples, 3)
true_theta = np.array([1.0, 2.0, 3.0])
y = X @ true_theta + 0.1 * np.random.randn(n_samples)
 
theta_init = np.zeros(3)
 
# SGD without decay
theta_no_decay, hist_no_decay = stochastic_gradient_descent(
    X, y, theta_init, learning_rate=0.01, n_epochs=20, lr_decay=None
)
 
# SGD with 1/t decay
theta_with_decay, hist_with_decay = stochastic_gradient_descent(
    X, y, theta_init, learning_rate=0.1, n_epochs=20, lr_decay='1/t'
)
 
print("True theta:", true_theta)
print(f"SGD (no decay):   {np.round(theta_no_decay, 3)}, final loss: {hist_no_decay[-1]['loss']:.6f}")
print(f"SGD (with decay): {np.round(theta_with_decay, 3)}, final loss: {hist_with_decay[-1]['loss']:.6f}")

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91

import numpy as np
 
def minibatch_gradient_descent(X, y, theta_init, batch_size=32,
                                 learning_rate=0.01, n_epochs=50):
    """
    Mini-batch Gradient Descent - the deep learning standard.
    
    Args:
        X: Feature matrix (n_samples, n_features)
        y: Target vector (n_samples,)
        theta_init: Initial parameters
        batch_size: Number of samples per batch
        learning_rate: Step size α
        n_epochs: Number of passes through entire dataset
    
    Returns:
        theta: Learned parameters
        history: Loss and gradient info per epoch
    """
    n_samples, n_features = X.shape
    theta = theta_init.copy()
    history = []
    
    for epoch in range(n_epochs):
        # Shuffle data at start of each epoch
        indices = np.random.permutation(n_samples)
        epoch_losses = []
        
        # Process data in batches
        for start_idx in range(0, n_samples, batch_size):
            end_idx = min(start_idx + batch_size, n_samples)
            batch_indices = indices[start_idx:end_idx]
            current_batch_size = len(batch_indices)
            
            # Extract batch
            X_batch = X[batch_indices]
            y_batch = y[batch_indices]
            
            # Forward pass
            y_pred = X_batch @ theta
            
            # Compute gradient (averaged over batch)
            residuals = y_pred - y_batch
            gradient = (X_batch.T @ residuals) / current_batch_size
            
            # Update parameters
            theta = theta - learning_rate * gradient
            
            # Track batch loss
            batch_loss = np.mean(residuals ** 2) / 2
            epoch_losses.append(batch_loss)
        
        # Epoch statistics
        epoch_loss = np.mean(epoch_losses)
        full_loss = np.mean((X @ theta - y) ** 2) / 2
        history.append({
            'epoch': epoch,
            'batch_loss': epoch_loss,
            'full_loss': full_loss,
            'n_batches': len(epoch_losses)
        })
    
    return theta, history
 
# Compare different batch sizes
np.random.seed(42)
n_samples = 5000
n_features = 10
 
X = np.random.randn(n_samples, n_features)
true_theta = np.linspace(-5, 5, n_features)
y = X @ true_theta + 0.1 * np.random.randn(n_samples)
 
theta_init = np.zeros(n_features)
 
batch_sizes = [1, 16, 64, 256, n_samples]
results = {}
 
for bs in batch_sizes:
    name = f"bs_{bs}" if bs < n_samples else "full"
    theta, hist = minibatch_gradient_descent(
        X, y, theta_init, batch_size=bs, learning_rate=0.01, n_epochs=20
    )
    results[name] = {
        'theta': theta,
        'final_loss': hist[-1]['full_loss'],
        'history': hist
    }
    print(f"Batch size {bs:5d}: Final loss = {hist[-1]['full_loss']:.6f}")
 
# All converge to similar solution, but trajectories differ

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85

import numpy as np
from typing import Iterator, Tuple
 
def shuffled_batches(X: np.ndarray, y: np.ndarray, 
                      batch_size: int) -> Iterator[Tuple[np.ndarray, np.ndarray]]:
    """
    Standard shuffled mini-batch iteration.
    Each sample seen exactly once per epoch.
    """
    n = len(X)
    indices = np.random.permutation(n)
    
    for start in range(0, n, batch_size):
        end = min(start + batch_size, n)
        batch_idx = indices[start:end]
        yield X[batch_idx], y[batch_idx]
 
def sampled_batches(X: np.ndarray, y: np.ndarray, 
                     batch_size: int, n_batches: int) -> Iterator[Tuple[np.ndarray, np.ndarray]]:
    """
    Random sampling WITH replacement.
    Some samples may repeat, others may be missed.
    """
    n = len(X)
    
    for _ in range(n_batches):
        batch_idx = np.random.choice(n, size=batch_size, replace=True)
        yield X[batch_idx], y[batch_idx]
 
def stratified_batches(X: np.ndarray, y: np.ndarray, 
                        batch_size: int) -> Iterator[Tuple[np.ndarray, np.ndarray]]:
    """
    Stratified batching: each batch has balanced class representation.
    Important for imbalanced datasets.
    """
    classes = np.unique(y)
    n_classes = len(classes)
    samples_per_class = batch_size // n_classes
    
    # Group indices by class
    class_indices = {c: np.where(y == c)[0] for c in classes}
    
    # Shuffle within each class
    for c in classes:
        np.random.shuffle(class_indices[c])
    
    # Create batches with balanced classes
    class_pointers = {c: 0 for c in classes}
    n = len(y)
    n_batches = n // batch_size
    
    for _ in range(n_batches):
        batch_idx = []
        for c in classes:
            start = class_pointers[c]
            end = start + samples_per_class
            
            # Handle wraparound
            if end > len(class_indices[c]):
                np.random.shuffle(class_indices[c])
                class_pointers[c] = 0
                start = 0
                end = samples_per_class
            
            batch_idx.extend(class_indices[c][start:end])
            class_pointers[c] = end
        
        batch_idx = np.array(batch_idx)
        np.random.shuffle(batch_idx)  # Shuffle within batch
        yield X[batch_idx], y[batch_idx]
 
# Example usage
X = np.random.randn(1000, 10)
y = np.random.randint(0, 5, 1000)  # 5 classes
 
print("Shuffled batches (standard):")
for i, (X_batch, y_batch) in enumerate(shuffled_batches(X, y, 64)):
    if i < 2:
        print(f"  Batch {i}: shape {X_batch.shape}, class distribution: {np.bincount(y_batch, minlength=5)}")
 
print("
Stratified batches:")
for i, (X_batch, y_batch) in enumerate(stratified_batches(X, y, 50)):
    if i < 2:
        print(f"  Batch {i}: shape {X_batch.shape}, class distribution: {np.bincount(y_batch, minlength=5)}")