K-Fold Procedure

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

import numpy as np
from typing import List, Tuple
 
def create_k_folds(n_samples: int, k: int, shuffle: bool = True, 
                   random_state: int = None) -> List[np.ndarray]:
    """
    Create k-fold indices for cross-validation.
    
    Parameters:
    -----------
    n_samples : int
        Total number of samples in the dataset
    k : int
        Number of folds to create
    shuffle : bool
        Whether to shuffle indices before partitioning
    random_state : int
        Random seed for reproducibility
        
    Returns:
    --------
    List of k numpy arrays, each containing the indices for one fold
    """
    # Create array of all indices
    indices = np.arange(n_samples)
    
    # Shuffle if requested (almost always should be True)
    if shuffle:
        rng = np.random.default_rng(random_state)
        rng.shuffle(indices)
    
    # Calculate fold sizes
    # If n=100 and k=3: fold_sizes = [34, 33, 33]
    fold_sizes = np.full(k, n_samples // k, dtype=int)
    fold_sizes[:n_samples % k] += 1  # Distribute remainder
    
    # Split indices into folds
    folds = []
    current_idx = 0
    for fold_size in fold_sizes:
        folds.append(indices[current_idx:current_idx + fold_size])
        current_idx += fold_size
    
    return folds
 
# Example: 10 samples into 3 folds
folds = create_k_folds(10, 3, shuffle=True, random_state=42)
for i, fold in enumerate(folds):
    print(f"Fold {i+1}: {fold}")
# Output (example):
# Fold 1: [8 1 5 7]
# Fold 2: [9 3 0]
# Fold 3: [2 4 6]

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
92
93

def k_fold_cross_validation(X, y, model_fn, k=5, metric_fn=None, 
                              random_state=42):
    """
    Perform k-fold cross-validation.
    
    Parameters:
    -----------
    X : array-like of shape (n_samples, n_features)
        Feature matrix
    y : array-like of shape (n_samples,)
        Target vector
    model_fn : callable
        Function that returns a fresh, unfitted model instance
    k : int
        Number of folds
    metric_fn : callable
        Function(y_true, y_pred) -> score. Default: accuracy
        
    Returns:
    --------
    scores : list of float
        Performance score for each fold
    predictions : array
        Out-of-fold predictions for all samples
    """
    import numpy as np
    from sklearn.metrics import accuracy_score
    
    if metric_fn is None:
        metric_fn = accuracy_score
    
    n_samples = len(X)
    folds = create_k_folds(n_samples, k, shuffle=True, 
                           random_state=random_state)
    
    # Store results
    scores = []
    predictions = np.empty(n_samples)
    predictions[:] = np.nan  # Initialize with NaN
    
    for fold_idx in range(k):
        # Identify validation indices (current fold)
        val_indices = folds[fold_idx]
        
        # Training indices = all other folds combined
        train_indices = np.concatenate([
            folds[i] for i in range(k) if i != fold_idx
        ])
        
        # Split data
        X_train, X_val = X[train_indices], X[val_indices]
        y_train, y_val = y[train_indices], y[val_indices]
        
        # Train model (get fresh instance each iteration!)
        model = model_fn()
        model.fit(X_train, y_train)
        
        # Predict on validation fold
        y_pred = model.predict(X_val)
        
        # Store predictions at original indices
        predictions[val_indices] = y_pred
        
        # Compute and store fold score
        fold_score = metric_fn(y_val, y_pred)
        scores.append(fold_score)
        
        print(f"Fold {fold_idx + 1}/{k}: "
              f"Train size = {len(train_indices)}, "
              f"Val size = {len(val_indices)}, "
              f"Score = {fold_score:.4f}")
    
    return scores, predictions
 
# Example usage
from sklearn.ensemble import RandomForestClassifier
from sklearn.datasets import make_classification
 
# Generate synthetic data
X, y = make_classification(n_samples=1000, n_features=20, 
                           random_state=42)
 
# Define model factory
def model_fn():
    return RandomForestClassifier(n_estimators=100, random_state=42)
 
# Run 5-fold CV
scores, predictions = k_fold_cross_validation(
    X, y, model_fn, k=5, random_state=42
)
 
print(f"
Mean CV Score: {np.mean(scores):.4f} ± {np.std(scores):.4f}")

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
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193

import numpy as np
from typing import List, Dict, Callable, Any, Optional, Tuple
from dataclasses import dataclass
import time
 
@dataclass
class FoldResult:
    """Results from a single fold."""
    fold_idx: int
    train_size: int
    val_size: int
    train_score: float
    val_score: float
    predictions: np.ndarray
    train_time: float
    predict_time: float
 
@dataclass
class CVResult:
    """Aggregated cross-validation results."""
    fold_results: List[FoldResult]
    mean_train_score: float
    mean_val_score: float
    std_val_score: float
    all_predictions: np.ndarray
    total_time: float
    
    def summary(self) -> str:
        return (f"CV Results: {self.mean_val_score:.4f} "
                f"± {self.std_val_score:.4f} "
                f"(train: {self.mean_train_score:.4f})")
 
class RobustKFold:
    """
    Production-quality k-fold cross-validation implementation.
    
    Features:
    - Reproducible random partitioning
    - Automatic handling of uneven fold sizes
    - Training and validation score tracking
    - Timing information for performance profiling
    - Out-of-fold prediction collection
    - Support for sample weights
    - Proper handling of edge cases
    """
    
    def __init__(self, k: int = 5, shuffle: bool = True, 
                 random_state: Optional[int] = None):
        self.k = k
        self.shuffle = shuffle
        self.random_state = random_state
        self._validate_params()
    
    def _validate_params(self):
        if self.k < 2:
            raise ValueError(f"k must be >= 2, got {self.k}")
        if not isinstance(self.k, int):
            raise TypeError(f"k must be integer, got {type(self.k)}")
    
    def split(self, n_samples: int) -> List[Tuple[np.ndarray, np.ndarray]]:
        """
        Generate train/validation index pairs for each fold.
        
        Returns list of (train_indices, val_indices) tuples.
        """
        if n_samples < self.k:
            raise ValueError(
                f"Cannot create {self.k} folds with only {n_samples} samples"
            )
        
        indices = np.arange(n_samples)
        
        if self.shuffle:
            rng = np.random.default_rng(self.random_state)
            rng.shuffle(indices)
        
        # Calculate fold boundaries
        fold_sizes = np.full(self.k, n_samples // self.k, dtype=int)
        fold_sizes[:n_samples % self.k] += 1
        
        splits = []
        current = 0
        
        for fold_size in fold_sizes:
            val_indices = indices[current:current + fold_size]
            train_indices = np.concatenate([
                indices[:current],
                indices[current + fold_size:]
            ])
            splits.append((train_indices, val_indices))
            current += fold_size
        
        return splits
    
    def cross_validate(
        self,
        X: np.ndarray,
        y: np.ndarray,
        model_factory: Callable[[], Any],
        scorer: Callable[[np.ndarray, np.ndarray], float],
        sample_weight: Optional[np.ndarray] = None,
        return_predictions: bool = True,
        verbose: bool = True
    ) -> CVResult:
        """
        Perform k-fold cross-validation.
        
        Parameters:
        -----------
        X : array of shape (n_samples, n_features)
        y : array of shape (n_samples,)
        model_factory : callable returning fresh model instance
        scorer : callable(y_true, y_pred) -> score
        sample_weight : optional weights for training
        return_predictions : whether to store all OOF predictions
        verbose : print progress
        """
        n_samples = len(X)
        splits = self.split(n_samples)
        
        fold_results = []
        all_predictions = np.full(n_samples, np.nan) if return_predictions else None
        
        total_start = time.time()
        
        for fold_idx, (train_idx, val_idx) in enumerate(splits):
            # Extract fold data
            X_train, X_val = X[train_idx], X[val_idx]
            y_train, y_val = y[train_idx], y[val_idx]
            
            # Handle sample weights
            if sample_weight is not None:
                fold_weights = sample_weight[train_idx]
            else:
                fold_weights = None
            
            # Train model
            model = model_factory()
            train_start = time.time()
            
            if fold_weights is not None:
                model.fit(X_train, y_train, sample_weight=fold_weights)
            else:
                model.fit(X_train, y_train)
            
            train_time = time.time() - train_start
            
            # Predictions
            predict_start = time.time()
            y_train_pred = model.predict(X_train)
            y_val_pred = model.predict(X_val)
            predict_time = time.time() - predict_start
            
            # Scores
            train_score = scorer(y_train, y_train_pred)
            val_score = scorer(y_val, y_val_pred)
            
            # Store OOF predictions
            if return_predictions:
                all_predictions[val_idx] = y_val_pred
            
            # Record results
            fold_result = FoldResult(
                fold_idx=fold_idx,
                train_size=len(train_idx),
                val_size=len(val_idx),
                train_score=train_score,
                val_score=val_score,
                predictions=y_val_pred,
                train_time=train_time,
                predict_time=predict_time
            )
            fold_results.append(fold_result)
            
            if verbose:
                print(f"Fold {fold_idx + 1}/{self.k}: "
                      f"train={train_score:.4f}, val={val_score:.4f}, "
                      f"time={train_time:.2f}s")
        
        total_time = time.time() - total_start
        
        # Aggregate
        val_scores = [r.val_score for r in fold_results]
        train_scores = [r.train_score for r in fold_results]
        
        return CVResult(
            fold_results=fold_results,
            mean_train_score=np.mean(train_scores),
            mean_val_score=np.mean(val_scores),
            std_val_score=np.std(val_scores),
            all_predictions=all_predictions,
            total_time=total_time
        )

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

from sklearn.model_selection import (
    KFold, 
    cross_val_score, 
    cross_val_predict,
    cross_validate
)
from sklearn.ensemble import RandomForestClassifier
from sklearn.datasets import make_classification
from sklearn.metrics import accuracy_score, make_scorer
import numpy as np
 
# Generate sample data
X, y = make_classification(
    n_samples=1000, 
    n_features=20,
    n_informative=10,
    random_state=42
)
 
# Method 1: Using cross_val_score (simplest)
# Returns array of k scores
model = RandomForestClassifier(n_estimators=100, random_state=42)
scores = cross_val_score(
    model, X, y, 
    cv=5,  # 5-fold CV
    scoring='accuracy'
)
print(f"cross_val_score: {scores.mean():.4f} ± {scores.std():.4f}")
print(f"Individual folds: {scores}")
 
# Method 2: Using KFold explicitly (more control)
kfold = KFold(n_splits=5, shuffle=True, random_state=42)
 
scores_manual = []
for train_idx, val_idx in kfold.split(X):
    model = RandomForestClassifier(n_estimators=100, random_state=42)
    model.fit(X[train_idx], y[train_idx])
    score = model.score(X[val_idx], y[val_idx])
    scores_manual.append(score)
 
print(f"Manual KFold: {np.mean(scores_manual):.4f}")
 
# Method 3: cross_val_predict (get all OOF predictions)
predictions = cross_val_predict(
    RandomForestClassifier(n_estimators=100, random_state=42),
    X, y,
    cv=5
)
print(f"OOF Accuracy: {accuracy_score(y, predictions):.4f}")
print(f"Predictions shape: {predictions.shape}")  # Same as y
 
# Method 4: cross_validate (comprehensive output)
cv_results = cross_validate(
    RandomForestClassifier(n_estimators=100, random_state=42),
    X, y,
    cv=5,
    scoring=['accuracy', 'f1', 'roc_auc'],
    return_train_score=True,
    return_estimator=True  # Keep the trained models
)
 
print("
cross_validate results:")
print(f"Test Accuracy: {cv_results['test_accuracy'].mean():.4f}")
print(f"Test F1: {cv_results['test_f1'].mean():.4f}")
print(f"Test AUC: {cv_results['test_roc_auc'].mean():.4f}")
print(f"Train Accuracy: {cv_results['train_accuracy'].mean():.4f}")
print(f"Fit times: {cv_results['fit_time']}")
print(f"Number of estimators: {len(cv_results['estimator'])}")

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

import numpy as np
from sklearn.model_selection import cross_val_predict, KFold
from sklearn.ensemble import RandomForestClassifier
from sklearn.linear_model import LogisticRegression
from sklearn.datasets import make_classification
from sklearn.metrics import accuracy_score, log_loss
import matplotlib.pyplot as plt
 
# Generate data
X, y = make_classification(
    n_samples=1000, n_features=20, 
    n_informative=10, random_state=42
)
 
# Application 1: Reliable Feature Engineering Evaluation
# --------------------------------------------------------
# Generate OOF predictions to use as features (meta-features)
 
kfold = KFold(n_splits=5, shuffle=True, random_state=42)
 
# OOF predictions from multiple models
oof_rf = cross_val_predict(
    RandomForestClassifier(n_estimators=100, random_state=42),
    X, y, cv=kfold, method='predict_proba'
)[:, 1]
 
oof_lr = cross_val_predict(
    LogisticRegression(max_iter=1000, random_state=42),
    X, y, cv=kfold, method='predict_proba'
)[:, 1]
 
print(f"RF OOF Log Loss: {log_loss(y, oof_rf):.4f}")
print(f"LR OOF Log Loss: {log_loss(y, oof_lr):.4f}")
 
# Application 2: Model Stacking with OOF
# ----------------------------------------
# Use OOF predictions as features for a meta-learner
 
meta_features = np.column_stack([oof_rf, oof_lr])
print(f"Meta-features shape: {meta_features.shape}")
 
# Train stacking model with CV on meta-features
meta_model = LogisticRegression()
stacked_score = cross_val_predict(
    meta_model, meta_features, y, cv=kfold, method='predict_proba'
)[:, 1]
 
print(f"Stacked OOF Log Loss: {log_loss(y, stacked_score):.4f}")
 
# Application 3: Calibration Analysis with OOF
# ----------------------------------------------
# OOF probabilities let us assess calibration on the full dataset
 
def plot_calibration_curve(y_true, y_prob, n_bins=10, title=""):
    """Plot calibration curve from OOF predictions."""
    bin_means = []
    true_frequencies = []
    
    bins = np.linspace(0, 1, n_bins + 1)
    for i in range(n_bins):
        mask = (y_prob >= bins[i]) & (y_prob < bins[i+1])
        if mask.sum() > 0:
            bin_means.append(y_prob[mask].mean())
            true_frequencies.append(y_true[mask].mean())
    
    plt.figure(figsize=(8, 6))
    plt.plot([0, 1], [0, 1], 'k--', label='Perfect calibration')
    plt.plot(bin_means, true_frequencies, 'o-', label='Model')
    plt.xlabel('Mean predicted probability')
    plt.ylabel('Fraction of positives')
    plt.title(f'Calibration Curve: {title}')
    plt.legend()
    plt.grid(True, alpha=0.3)
    return plt.gcf()
 
# Application 4: Threshold Selection with OOF
# ---------------------------------------------
# Find optimal decision threshold using OOF predictions
 
from sklearn.metrics import precision_recall_curve, f1_score
 
precision, recall, thresholds = precision_recall_curve(y, oof_rf)
f1_scores = 2 * (precision * recall) / (precision + recall + 1e-10)
best_threshold = thresholds[np.argmax(f1_scores[:-1])]
 
print(f"
Optimal threshold for F1: {best_threshold:.4f}")
print(f"F1 at default (0.5): {f1_score(y, oof_rf > 0.5):.4f}")
print(f"F1 at optimal: {f1_score(y, oof_rf > best_threshold):.4f}")