diff --git a/doc/Programs/DecistionTrees/ensemblemethods.py b/doc/Programs/DecistionTrees/ensemblemethods.py deleted file mode 100644 index e69de29bb..000000000 diff --git a/doc/Programs/EnsembleMethods/classification_results.csv b/doc/Programs/EnsembleMethods/classification_results.csv new file mode 100644 index 000000000..f5f76f99f --- /dev/null +++ b/doc/Programs/EnsembleMethods/classification_results.csv @@ -0,0 +1,201 @@ +GroundTruth,DT_Pred,RF_Pred,GB_Pred +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +0,0,0,0 +0,0,0,0 +0,0,0,0 +1,1,1,1 +1,1,1,1 +1,1,1,1 +0,0,0,0 +1,1,1,1 +0,0,0,0 +1,1,1,1 diff --git a/doc/Programs/EnsembleMethods/ensemblemethods.py b/doc/Programs/EnsembleMethods/ensemblemethods.py new file mode 100644 index 000000000..a96fd0597 --- /dev/null +++ b/doc/Programs/EnsembleMethods/ensemblemethods.py @@ -0,0 +1,483 @@ +import numpy as np + +class DecisionTreeClassifier: + def __init__(self, criterion='gini', max_depth=None, min_samples_split=2, max_features=None): + self.criterion = criterion + self.max_depth = max_depth + self.min_samples_split = min_samples_split + self.max_features = max_features + self.tree = None + + class Node: + def __init__(self, feature=None, threshold=None, left=None, right=None, *, value=None): + self.feature = feature + self.threshold = threshold + self.left = left + self.right = right + self.value = value # Leaf class label + + def fit(self, X, y): + X, y = np.array(X), np.array(y) + self.n_features_ = X.shape[1] + self.tree = self._build_tree(X, y, depth=0) + return self + + def _build_tree(self, X, y, depth): + num_samples, _ = X.shape + # Stop if conditions met + if num_samples < self.min_samples_split or (self.max_depth is not None and depth >= self.max_depth) or len(np.unique(y)) == 1: + leaf_val = self._majority_class(y) + return DecisionTreeClassifier.Node(value=leaf_val) + # Find best split + feat_idx, thr = self._best_split(X, y) + if feat_idx is None: + leaf_val = self._majority_class(y) + return DecisionTreeClassifier.Node(value=leaf_val) + # Split data + left_mask = X[:, feat_idx] <= thr + left_node = self._build_tree(X[left_mask], y[left_mask], depth+1) + right_node = self._build_tree(X[~left_mask], y[~left_mask], depth+1) + return DecisionTreeClassifier.Node(feature=feat_idx, threshold=thr, left=left_node, right=right_node) + + def _best_split(self, X, y): + best_gain = 0 + best_feat, best_thr = None, None + if self.criterion == 'gini': + base_impurity = self._gini(y) + else: + base_impurity = self._entropy(y) + n_features = X.shape[1] + features = range(n_features) + # Possibly sample subset of features (for Random Forest use) + if self.max_features is not None: + if isinstance(self.max_features, int) and self.max_features < n_features: + features = np.random.choice(n_features, self.max_features, replace=False) + elif isinstance(self.max_features, float): + k = int(n_features * self.max_features) + features = np.random.choice(n_features, k, replace=False) + for feat in features: + X_col = X[:, feat] + unique_vals = np.unique(X_col) + if len(unique_vals) <= 1: + continue + # Try midpoints between sorted unique values + thresholds = (unique_vals[:-1] + unique_vals[1:]) / 2.0 + for thr in thresholds: + left_mask = X_col <= thr + y_left, y_right = y[left_mask], y[~left_mask] + if len(y_left) == 0 or len(y_right) == 0: + continue + # Compute impurity of the split + if self.criterion == 'gini': + imp_left = self._gini(y_left) + imp_right = self._gini(y_right) + else: + imp_left = self._entropy(y_left) + imp_right = self._entropy(y_right) + p = float(len(y_left)) / len(y) + gain = base_impurity - (p * imp_left + (1 - p) * imp_right) + if gain > best_gain: + best_gain, best_feat, best_thr = gain, feat, thr + return best_feat, best_thr + + def _gini(self, y): + _, counts = np.unique(y, return_counts=True) + p = counts / counts.sum() + return 1.0 - np.sum(p**2) + + def _entropy(self, y): + _, counts = np.unique(y, return_counts=True) + p = counts / counts.sum() + p = p[p > 0] + return -np.sum(p * np.log2(p)) + + def _majority_class(self, y): + unique, counts = np.unique(y, return_counts=True) + return unique[np.argmax(counts)] + + def predict(self, X): + X = np.array(X) + return np.array([self._predict_input(x, self.tree) for x in X]) + + def _predict_input(self, x, node): + if node.value is not None: + return node.value + if x[node.feature] <= node.threshold: + return self._predict_input(x, node.left) + else: + return self._predict_input(x, node.right) + +class DecisionTreeRegressor: + def __init__(self, max_depth=None, min_samples_split=2, max_features=None): + self.max_depth = max_depth + self.min_samples_split = min_samples_split + self.max_features = max_features + self.tree = None + + class Node: + def __init__(self, feature=None, threshold=None, left=None, right=None, *, value=None): + self.feature = feature + self.threshold = threshold + self.left = left + self.right = right + self.value = value # Leaf output value + + def fit(self, X, y): + X, y = np.array(X), np.array(y, dtype=float) + self.n_features_ = X.shape[1] + self.tree = self._build_tree(X, y, depth=0) + return self + + def _build_tree(self, X, y, depth): + num_samples, _ = X.shape + if num_samples < self.min_samples_split or (self.max_depth is not None and depth >= self.max_depth) or np.var(y) == 0: + leaf_val = np.mean(y) + return DecisionTreeRegressor.Node(value=leaf_val) + best_sse = float('inf') + best_feat, best_thr = None, None + n_features = X.shape[1] + features = range(n_features) + if self.max_features is not None: + if isinstance(self.max_features, int) and self.max_features < n_features: + features = np.random.choice(n_features, self.max_features, replace=False) + elif isinstance(self.max_features, float): + k = int(n_features * self.max_features) + features = np.random.choice(n_features, k, replace=False) + for feat in features: + X_col = X[:, feat] + unique_vals = np.unique(X_col) + if len(unique_vals) <= 1: + continue + thresholds = (unique_vals[:-1] + unique_vals[1:]) / 2.0 + for thr in thresholds: + left_mask = X_col <= thr + y_left, y_right = y[left_mask], y[~left_mask] + if len(y_left) == 0 or len(y_right) == 0: + continue + # Compute sum of squared errors (SSE) + left_mean, right_mean = np.mean(y_left), np.mean(y_right) + sse_left = np.sum((y_left - left_mean) ** 2) + sse_right = np.sum((y_right - right_mean) ** 2) + sse = sse_left + sse_right + if sse < best_sse: + best_sse, best_feat, best_thr = sse, feat, thr + if best_feat is None: + leaf_val = np.mean(y) + return DecisionTreeRegressor.Node(value=leaf_val) + left_mask = X[:, best_feat] <= best_thr + left_node = self._build_tree(X[left_mask], y[left_mask], depth+1) + right_node = self._build_tree(X[~left_mask], y[~left_mask], depth+1) + return DecisionTreeRegressor.Node(feature=best_feat, threshold=best_thr, left=left_node, right=right_node) + + def predict(self, X): + X = np.array(X) + return np.array([self._predict_input(x, self.tree) for x in X]) + + def _predict_input(self, x, node): + if node.value is not None: + return node.value + if x[node.feature] <= node.threshold: + return self._predict_input(x, node.left) + else: + return self._predict_input(x, node.right) + +# Random Forests (Classification and Regression) + +""" +Random forests train an ensemble of decision trees on bootstrapped data subsets and average their outputs . For classification, the final class is the majority vote of all trees; for regression, the output is the average prediction . Key points: + +Bootstrap sampling: Each tree is trained on a random sample (with replacement) of the data. +Feature randomness: When splitting, each node may consider only a random subset of features (parameter max_features). +Aggregation: Classification uses mode of tree predictions; regression uses mean. This reduces overfitting compared to a single tree . +""" + +import numpy as np +from collections import Counter + +class RandomForestClassifier: + def __init__(self, n_estimators=100, max_depth=None, min_samples_split=2, max_features='sqrt'): + self.n_estimators = n_estimators + self.max_depth = max_depth + self.min_samples_split = min_samples_split + self.max_features = max_features + self.trees = [] + + def fit(self, X, y): + X, y = np.array(X), np.array(y) + n_samples, n_features = X.shape + # Determine how many features to try at each split + if self.max_features == 'sqrt': + max_feats = int(np.sqrt(n_features)) + elif self.max_features == 'log2': + max_feats = int(np.log2(n_features)) + elif isinstance(self.max_features, int): + max_feats = self.max_features + elif isinstance(self.max_features, float): + max_feats = int(n_features * self.max_features) + else: + max_feats = n_features + # Build trees + for _ in range(self.n_estimators): + indices = np.random.choice(n_samples, n_samples, replace=True) + X_sample, y_sample = X[indices], y[indices] + tree = DecisionTreeClassifier(criterion='gini', + max_depth=self.max_depth, + min_samples_split=self.min_samples_split, + max_features=max_feats) + tree.fit(X_sample, y_sample) + self.trees.append(tree) + return self + + def predict(self, X): + X = np.array(X) + # Collect predictions from all trees + tree_preds = np.array([tree.predict(X) for tree in self.trees]).T # shape (n_samples, n_trees) + y_pred = [] + for preds in tree_preds: + vote = Counter(preds).most_common(1)[0][0] + y_pred.append(vote) + return np.array(y_pred) + +class RandomForestRegressor: + def __init__(self, n_estimators=100, max_depth=None, min_samples_split=2, max_features=None): + self.n_estimators = n_estimators + self.max_depth = max_depth + self.min_samples_split = min_samples_split + self.max_features = max_features + self.trees = [] + + def fit(self, X, y): + X, y = np.array(X), np.array(y, dtype=float) + n_samples, n_features = X.shape + if self.max_features == 'sqrt': + max_feats = int(np.sqrt(n_features)) + elif self.max_features == 'log2': + max_feats = int(np.log2(n_features)) + elif isinstance(self.max_features, int): + max_feats = self.max_features + elif isinstance(self.max_features, float): + max_feats = int(n_features * self.max_features) + else: + max_feats = n_features + for _ in range(self.n_estimators): + indices = np.random.choice(n_samples, n_samples, replace=True) + X_sample, y_sample = X[indices], y[indices] + tree = DecisionTreeRegressor(max_depth=self.max_depth, + min_samples_split=self.min_samples_split, + max_features=max_feats) + tree.fit(X_sample, y_sample) + self.trees.append(tree) + return self + + def predict(self, X): + X = np.array(X) + # Average predictions from all trees + tree_preds = np.array([tree.predict(X) for tree in self.trees]) + return np.mean(tree_preds, axis=0) + +# Gradient Boosting (Classification and Regression) + +""" +Gradient boosting builds an additive ensemble of trees by fitting each new tree on the residuals (errors) of the existing model, effectively performing gradient descent on a loss function . At each stage m: + +Compute the pseudo-residuals r_{im} = -\partial L(y_i, F(x_i)) / \partial F(x_i) (the negative gradient) . +Fit a tree h_m(x) to these residuals. +Update the model: F_m(x) = F_{m-1}(x) + \gamma_m \, h_m(x) (with line-search multiplier \gamma_m or simply a learning rate). + + +For binary classification, we use the logistic loss: initialize F_0 = \log(p/(1-p)) (log-odds of positive class) and repeatedly fit trees to y - \sigma(F). For multiclass, we fit one-vs-rest models (one boosting ensemble per class) and predict the class with highest score. Key concepts: + +Additive updates: Each tree’s prediction is scaled by a learning rate and added to the ensemble output. +Regression loss: Typically squared-error (L2) for regression, logistic (cross-entropy) for classification. +Weak learners: Trees are often shallow (small max_depth). +""" + +import numpy as np +from math import log, exp + +class GradientBoostingRegressor: + def __init__(self, n_estimators=100, learning_rate=0.1, max_depth=3): + self.n_estimators = n_estimators + self.learning_rate = learning_rate + self.max_depth = max_depth + self.trees = [] + self.gammas = [] + self.initial_prediction = None + + def fit(self, X, y): + X, y = np.array(X), np.array(y, dtype=float) + # Initialize prediction with mean + self.initial_prediction = np.mean(y) + F = np.full_like(y, fill_value=self.initial_prediction, dtype=float) + for _ in range(self.n_estimators): + residual = y - F + tree = DecisionTreeRegressor(max_depth=self.max_depth) + tree.fit(X, residual) + pred = tree.predict(X) + # Line search for optimal multiplier gamma + gamma = np.dot(residual, pred) / (np.dot(pred, pred) + 1e-8) + F = F + self.learning_rate * gamma * pred + self.trees.append(tree) + self.gammas.append(gamma) + return self + + def predict(self, X): + X = np.array(X) + F = np.full(X.shape[0], fill_value=self.initial_prediction, dtype=float) + for tree, gamma in zip(self.trees, self.gammas): + F += self.learning_rate * gamma * tree.predict(X) + return F + +class GradientBoostingClassifier: + def __init__(self, n_estimators=100, learning_rate=0.1, max_depth=3): + self.n_estimators = n_estimators + self.learning_rate = learning_rate + self.max_depth = max_depth + self.models = [] # For multiclass, one model per class + + def fit(self, X, y): + X, y = np.array(X), np.array(y) + self.classes_ = np.unique(y) + if len(self.classes_) <= 2: + # Binary classification (labels may be 0/1 or not) + # Map labels to 0/1 + if set(self.classes_) != {0, 1}: + class0, class1 = self.classes_[0], self.classes_[1] + y_bin = np.array([0 if yi==class0 else 1 for yi in y]) + self.class_map = {0: class0, 1: class1} + else: + y_bin = y + self.class_map = None + # Initialize log-odds + p = np.clip(np.mean(y_bin), 1e-6, 1-1e-6) + F = np.full(y_bin.shape, fill_value=log(p/(1-p)), dtype=float) + self.initial_F = F[0] + self.trees = [] + for _ in range(self.n_estimators): + P = 1 / (1 + np.exp(-F)) + residual = y_bin - P + tree = DecisionTreeRegressor(max_depth=self.max_depth) + tree.fit(X, residual) + pred = tree.predict(X) + F = F + self.learning_rate * pred + self.trees.append(tree) + else: + # Multiclass one-vs-rest + for cls in self.classes_: + y_binary = (y == cls).astype(int) + model = GradientBoostingClassifier(n_estimators=self.n_estimators, + learning_rate=self.learning_rate, + max_depth=self.max_depth) + model.fit(X, y_binary) + self.models.append(model) + return self + + def predict(self, X): + X = np.array(X) + if len(self.classes_) <= 2: + # Binary + F = np.full(X.shape[0], fill_value=self.initial_F, dtype=float) + for tree in self.trees: + F += self.learning_rate * tree.predict(X) + P = 1 / (1 + np.exp(-F)) + y_pred = (P >= 0.5).astype(int) + if self.class_map: + inv_map = {v:k for k,v in self.class_map.items()} + y_pred = np.array([inv_map[val] for val in y_pred]) + return y_pred + else: + # Multiclass: compute score for each class + scores = [] + for model in self.models: + F_cls = np.full(X.shape[0], fill_value=model.initial_F, dtype=float) + for tree in model.trees: + F_cls += self.learning_rate * tree.predict(X) + scores.append(F_cls) + scores = np.vstack(scores).T # shape (n_samples, n_classes) + class_idx = np.argmax(scores, axis=1) + return np.array([self.classes_[i] for i in class_idx]) + + +# Classification Example: Generate a 2-class dataset, train models, and evaluate accuracy. +import numpy as np +import csv + +# Synthetic binary classification data +np.random.seed(0) +N = 100 +# Class 0 centered at (-2, -2), Class 1 at (2, 2) +X0 = np.random.randn(N, 2) - 2 +X1 = np.random.randn(N, 2) + 2 +X_clf = np.vstack([X0, X1]) +y_clf = np.array([0]*N + [1]*N) + +# Shuffle data +perm = np.random.permutation(len(y_clf)) +X_clf, y_clf = X_clf[perm], y_clf[perm] + +# Train models +dt_clf = DecisionTreeClassifier(max_depth=3) +dt_clf.fit(X_clf, y_clf) +rf_clf = RandomForestClassifier(n_estimators=10, max_depth=3) +rf_clf.fit(X_clf, y_clf) +gb_clf = GradientBoostingClassifier(n_estimators=20, learning_rate=0.1, max_depth=2) +gb_clf.fit(X_clf, y_clf) + +# Predictions +pred_dt = dt_clf.predict(X_clf) +pred_rf = rf_clf.predict(X_clf) +pred_gb = gb_clf.predict(X_clf) + +# Accuracy evaluation +acc_dt = np.mean(pred_dt == y_clf) +acc_rf = np.mean(pred_rf == y_clf) +acc_gb = np.mean(pred_gb == y_clf) +print(f"Decision Tree Accuracy: {acc_dt:.2f}") +print(f"Random Forest Accuracy: {acc_rf:.2f}") +print(f"Gradient Boosting Accuracy: {acc_gb:.2f}") + +# Export to CSV +with open('classification_results.csv', 'w', newline='') as f: + writer = csv.writer(f) + writer.writerow(['GroundTruth', 'DT_Pred', 'RF_Pred', 'GB_Pred']) + for true, d, r, g in zip(y_clf, pred_dt, pred_rf, pred_gb): + writer.writerow([true, d, r, g]) + +# Regression Example: Generate a simple regression dataset, train models, and compute MSE. +import numpy as np +import csv + +# Synthetic regression data: y = 3*x1 - 2*x2 + noise +np.random.seed(1) +N = 200 +X_reg = np.random.randn(N, 2) +y_reg = 3 * X_reg[:,0] - 2 * X_reg[:,1] + np.random.randn(N) * 0.5 + +# Train models +dt_reg = DecisionTreeRegressor(max_depth=4) +dt_reg.fit(X_reg, y_reg) +rf_reg = RandomForestRegressor(n_estimators=10, max_depth=4) +rf_reg.fit(X_reg, y_reg) +gb_reg = GradientBoostingRegressor(n_estimators=50, learning_rate=0.1, max_depth=2) +gb_reg.fit(X_reg, y_reg) + +# Predictions +pred_dt_r = dt_reg.predict(X_reg) +pred_rf_r = rf_reg.predict(X_reg) +pred_gb_r = gb_reg.predict(X_reg) + +# MSE evaluation +mse_dt = np.mean((pred_dt_r - y_reg)**2) +mse_rf = np.mean((pred_rf_r - y_reg)**2) +mse_gb = np.mean((pred_gb_r - y_reg)**2) +print(f"Decision Tree MSE: {mse_dt:.3f}") +print(f"Random Forest MSE: {mse_rf:.3f}") +print(f"Gradient Boosting MSE: {mse_gb:.3f}") + +# Export to CSV +with open('regression_results.csv', 'w', newline='') as f: + writer = csv.writer(f) + writer.writerow(['GroundTruth', 'DT_Pred', 'RF_Pred', 'GB_Pred']) + for true, d, r, g in zip(y_reg, pred_dt_r, pred_rf_r, pred_gb_r): + writer.writerow([true, d, r, g]) + diff --git a/doc/Programs/EnsembleMethods/regression_results.csv b/doc/Programs/EnsembleMethods/regression_results.csv new file mode 100644 index 000000000..603c8804a --- /dev/null +++ b/doc/Programs/EnsembleMethods/regression_results.csv @@ -0,0 +1,201 @@ +GroundTruth,DT_Pred,RF_Pred,GB_Pred +5.4432818818678665,4.587413399100595,4.867399782178129,5.790888921291713 +0.5996122276017699,-0.17286001207987395,0.6392547204856381,0.3346234065402876 +7.382916188676534,4.019798593894754,4.728583892005671,6.893862359327972 +7.373298690627765,8.311136903189661,6.341214611261762,6.641424617183763 +1.2444295584307281,0.7829612037013028,1.682913433262075,1.9593407654201414 +8.549837433392373,8.311136903189661,7.361883475207449,8.003853760576028 +-1.2703762672423755,-0.17286001207987395,-0.01653138757934134,-0.254473927488898 +5.18600642962302,4.019798593894754,3.973058549440549,5.08945059431903 +1.464240188469057,1.7340615070817575,1.701365449146443,1.3513265538013948 +-0.4869020241332592,0.7829612037013028,0.3225235907484086,-0.6615231670329242 +-5.732173085850799,-4.339380045700343,-5.096052232059108,-5.773972604400032 +2.7279612455737965,0.7829612037013028,1.2126254014394149,2.2910512207060636 +4.9501481793664475,4.019798593894754,3.564288400025178,4.514052426711978 +1.4725419460747882,1.7340615070817575,1.701365449146443,1.3513265538013948 +-3.071126672941314,-3.025574053843865,-2.88365816162571,-2.1926790779131244 +-2.170258389366966,-1.5513454724715776,-1.2777089856188693,-1.1235787567005935 +-0.7600362306750654,-1.5513454724715776,-1.2777089856188693,-0.7497144820358775 +-1.43048864071057,-1.5513454724715776,-1.4360687647137091,-1.2677563323563124 +-3.665626297648134,-3.025574053843865,-3.6540151419624793,-3.841016720812165 +2.448194302063534,4.587413399100595,3.597183239144573,3.4642327450250083 +1.0853683566497465,1.7340615070817575,1.6514228531497448,1.0153227147797779 +-4.819703396026448,-4.339380045700343,-4.255139465994213,-4.718460794735865 +1.2390122138154673,1.7340615070817575,1.3496566838538135,1.2168193423286342 +-4.002748627593119,-1.7909353004230013,-1.855366215191576,-3.0079622328719537 +0.15338268928910692,0.7829612037013028,0.3144003745402276,-0.6306150644612731 +1.6317154212536937,1.7340615070817575,1.9416888999030886,1.5279158676213782 +-2.9684476986337054,-1.5513454724715776,-2.0519042321805316,-2.6752759452833175 +-1.6248455030745903,-3.025574053843865,-2.2335476362719087,-1.903918053168852 +0.6633284422034442,0.7829612037013028,0.25432718961037304,0.7461319230554601 +-1.1280914918962268,0.7829612037013028,0.3194589353167796,-0.6717758998797576 +-4.164701969180553,-4.339380045700343,-4.255139465994213,-4.697271546092074 +2.692826032851701,4.019798593894754,2.689732056504237,2.3182495894505046 +2.037128644725154,0.7829612037013028,1.4761891913388827,1.4795864630555156 +0.30381092072311466,-1.7909353004230013,-0.7092680334850268,0.31557743212992867 +9.923169078775501,10.509218150849213,9.02003515775919,9.429862099478795 +-3.3482609799752034,-3.3432405724417764,-3.1227871611703937,-3.3547859184344153 +-1.0389050002745381,0.7829612037013028,0.1136107068795125,-1.05769970141103 +5.508150707342573,4.019798593894754,4.0220233007743555,4.970341115762031 +-2.1701391429667667,-3.025574053843865,-3.0405184152400695,-2.6197637148628807 +0.06113918707223631,0.7829612037013028,0.3144003745402276,-1.05769970141103 +-0.039826281950997466,-0.17286001207987395,-0.43425679234571896,-0.35747208129069635 +-1.1024491143857782,0.7829612037013028,0.3225235907484086,-0.5069711829117487 +-0.2232031837257465,0.7829612037013028,0.6080868319430042,0.17379963398440912 +-2.092061022168067,-3.025574053843865,-3.2130724780464703,-2.6395826604071075 +-2.0591455878158427,-1.7909353004230013,-0.5807417422123871,-1.7787408330867247 +3.4197103773031787,0.7829612037013028,1.8359115230257324,2.40718144444212 +-0.2741218636273111,-0.17286001207987395,-0.426498988536295,0.1080940522042256 +1.6162436373353701,0.7829612037013028,1.2700395612239626,1.3393209055199082 +-1.311170864798116,-0.17286001207987395,-0.8575825422373018,-0.7229009684417831 +-2.5270124806840655,-3.025574053843865,-3.369932731660829,-2.9343238155264455 +-4.056518114418082,-4.339380045700343,-3.6688578553761966,-3.4965064757326965 +0.582384577257114,0.7829612037013028,0.8061891073857292,0.8434445188007198 +-3.2863023509300273,-3.025574053843865,-3.4500901420054872,-3.0274839314000364 +3.767874604197618,4.019798593894754,3.564288400025178,3.7975335659649185 +-0.3143864450239937,-0.17286001207987395,-0.8575825422373018,-0.6510988932866465 +-5.200487990062849,-6.558872933068393,-5.246777111543906,-4.881207106124984 +3.846280229966947,4.019798593894754,3.564288400025178,3.9069174951211596 +4.037060404523531,4.019798593894754,3.4377528996394346,4.361781992467218 +2.802787180875442,1.7340615070817575,1.701365449146443,2.620368262532006 +2.2495305152060814,0.7829612037013028,1.4761891913388827,2.3383447286871277 +1.304522016105788,1.7340615070817575,1.701365449146443,1.3513265538013948 +-0.6130203029693974,-1.7909353004230013,-0.730810838494616,-0.6874492035586662 +-7.523776449408949,-6.558872933068393,-7.016913930010119,-7.776608228535077 +4.511186215540193,4.587413399100595,4.402510058651032,4.01152677230024 +-5.359159080438843,-6.558872933068393,-5.762029971793566,-5.447744904373908 +0.6108843157652932,1.7340615070817575,1.0905601503569895,1.0488228857812256 +-4.566171642569727,-6.558872933068393,-5.246777111543906,-4.671529962596242 +4.508149607377162,4.019798593894754,3.564288400025178,3.445798594983944 +4.285896691732519,4.019798593894754,3.564288400025178,4.108008262663981 +2.0582199316280425,0.7829612037013028,1.4761891913388827,1.7544714629847056 +0.032492031960675516,0.7829612037013028,0.20820495778850184,0.2077976181180528 +1.046449605236043,0.7829612037013028,0.5741700408947233,1.812593909185824 +1.2134944611614933,0.7829612037013028,1.4761891913388827,1.3948339401356011 +-0.07006145268286867,0.7829612037013028,0.3144003745402276,-0.6615231670329242 +-0.8770880850318181,0.7829612037013028,0.3144003745402276,-0.8017123380310867 +3.699732750096395,2.7501162166980215,2.5912952314524196,2.696194024570413 +-2.2335381319928893,-1.7909353004230013,-1.7385878678374518,-0.9679068353003484 +0.8232392245461774,0.7829612037013028,1.2700395612239626,1.3393209055199082 +0.05433561322881994,0.7829612037013028,0.661536104096642,0.3188128978609615 +1.9488052318011708,1.7340615070817575,1.701365449146443,1.8455174069013591 +0.7004771824314627,-0.17286001207987395,0.6392547204856381,0.2426012636769862 +0.7822442627220687,0.7829612037013028,0.8379825835296415,0.8133573900520248 +1.451385378253016,0.7829612037013028,1.4761891913388827,1.690775462532898 +2.869223114270844,0.7829612037013028,1.4761891913388827,2.359887936461056 +11.095267222922924,10.509218150849213,9.02003515775919,10.549250452247348 +-3.4474874473501362,-3.025574053843865,-3.572923544234338,-3.4551243091154564 +8.681162096710741,8.429822177836392,6.104579039773336,7.762735508529412 +1.4359325502485811,0.7829612037013028,0.8379825835296415,0.5300018973383436 +5.089031641291861,4.587413399100595,4.286043383540633,4.904037248963674 +2.2205926982816435,4.019798593894754,2.689732056504237,2.3182495894505046 +-4.497444062918618,-6.558872933068393,-5.246777111543906,-4.479049365291725 +-0.7617994894160899,-1.7909353004230013,-0.7092680334850268,-0.1599876363739885 +0.6959595975286552,1.7340615070817575,1.3496566838538135,1.1598381655241343 +1.432782124042226,0.7829612037013028,1.4761891913388827,1.3948339401356011 +0.19864249616523988,0.7829612037013028,0.8429986658729556,0.45124708827155163 +2.327716287801658,0.7829612037013028,1.8359115230257324,2.3383447286871277 +3.7563026006945215,4.587413399100595,4.005070638228162,3.661558931299579 +4.113482813976909,4.019798593894754,4.088012379079276,4.88285911829463 +0.41253718245049165,0.7829612037013028,1.2126254014394149,1.3259972243806084 +0.3351116700410785,0.7829612037013028,0.25432718961037304,0.7461319230554601 +-3.7103430495773573,-3.025574053843865,-3.0405184152400695,-2.6915657900180174 +-5.568111945034292,-4.339380045700343,-4.354224575252421,-5.47155094338898 +-0.8742516241000576,-1.5513454724715776,-2.2746654204146557,-1.3955644674004062 +-3.3168499268966194,-3.3432405724417764,-3.1227871611703937,-2.818344255785449 +-7.3989954297867335,-6.558872933068393,-6.275717813009518,-7.716795880802433 +1.8004996832996478,2.7501162166980215,1.0492223214425918,1.5100929070361524 +4.234302440691409,4.019798593894754,3.564288400025178,3.835373501193173 +-5.417394781261496,-6.558872933068393,-5.762029971793566,-5.576915838757433 +2.4394214953198614,4.019798593894754,3.564288400025178,2.546013595187578 +4.242453596340928,4.587413399100595,4.402510058651032,4.01152677230024 +4.692271589216133,4.587413399100595,4.402510058651032,3.942690056545247 +4.556046118633867,4.587413399100595,4.402510058651032,4.01152677230024 +-8.473035861386107,-6.558872933068393,-7.2803130888251415,-7.982008933571979 +-1.8216507637505108,-1.5513454724715776,-1.1902386778102552,-1.0774687816744362 +-1.3477235040044544,-1.5513454724715776,-1.2777089856188693,-1.1235787567005935 +-0.05468569878074042,0.7829612037013028,0.6080868319430042,0.17379963398440912 +-2.0056791533911196,-1.7909353004230013,-1.0259029774097224,-1.1171700125187718 +2.3835674819563475,1.7340615070817575,1.701365449146443,1.9974521224608754 +-2.124925970951713,-1.5513454724715776,-2.0519042321805316,-2.2405904800263436 +-6.1190593462602365,-6.558872933068393,-6.052742811723908,-6.160005891605775 +-0.151232270039378,-0.17286001207987395,-0.426498988536295,-0.3262760026440347 +2.5177256195756064,1.7340615070817575,1.701365449146443,2.5470161220054233 +-3.8691836386260987,-3.025574053843865,-3.8269338070992065,-3.295323408672829 +3.622168845547816,4.019798593894754,3.564288400025178,2.546013595187578 +8.526328212606373,8.311136903189661,7.122440618979657,7.837114185387792 +-0.216757619042365,-0.17286001207987395,-0.426498988536295,-0.5787233927860647 +-11.601209595104056,-11.601209595104056,-9.561185675768339,-10.512002827005508 +2.88587332285357,4.019798593894754,3.294860695394168,3.1660924242631885 +5.300488032934883,4.587413399100595,5.062997672072164,4.778062512164518 +3.20720191114799,4.019798593894754,3.129889275456901,2.4642040664303697 +4.180484947375243,1.7340615070817575,2.084779047349223,2.6512763651036573 +2.730862512269863,4.019798593894754,2.975295922574877,2.5019657458514053 +-2.1278022385512543,-1.7909353004230013,-0.8260463808391512,-1.8626403871257085 +0.5206670284828185,0.7829612037013028,0.7895493937193178,0.27532572182334797 +2.1515424507522685,1.7340615070817575,2.2662416091255366,1.7676186374592735 +-0.7746121832101851,-1.7909353004230013,-0.7092680334850268,-0.1599876363739885 +0.8486000154550029,0.7829612037013028,0.3733803732340081,0.7461319230554601 +-1.0729751024422207,-1.5513454724715776,-1.6437232249048663,-1.0646041001561037 +-7.6810906601322495,-6.558872933068393,-7.016913930010119,-7.193345627322054 +-3.381303829864293,-4.339380045700343,-3.6688578553761966,-3.44589364922135 +-1.682963477688704,-1.5513454724715776,-1.3774324801021678,-1.4057173607400517 +-1.9223903666104163,0.7829612037013028,-0.005442476744122615,-1.328958564740324 +-1.022867997944866,0.7829612037013028,0.6080868319430042,0.14289153141275795 +0.20211849485077876,0.7829612037013028,0.6080868319430042,0.14289153141275795 +-4.624820571525798,-3.025574053843865,-3.1145876393796565,-3.5566190071198314 +-1.295754821123852,0.7829612037013028,-0.06451881681766505,-0.017190941160655657 +2.7004429185471466,4.019798593894754,3.564288400025178,3.0126208930108955 +0.9556774761584665,0.7829612037013028,0.8379825835296415,0.7824492874803737 +0.3372176657347352,0.7829612037013028,0.5741700408947233,0.6723472452309062 +-4.615587966790981,-3.025574053843865,-3.8269338070992065,-3.8318986298457185 +9.050847714329977,8.311136903189661,8.191023038289902,8.874580922666423 +-1.484588673800784,-1.5513454724715776,-1.5165767561053662,-1.0646041001561037 +-3.5773878046504186,-4.339380045700343,-3.6688578553761966,-3.44589364922135 +0.6911638893558922,1.7340615070817575,1.3496566838538135,1.2168193423286342 +2.993421702554793,4.019798593894754,3.564288400025178,2.8539311695381246 +3.27619429505266,0.7829612037013028,1.4761891913388827,2.40718144444212 +5.3778196611518885,4.587413399100595,5.121782776091066,5.847870098096212 +8.178482258962045,8.429822177836392,5.884152196000407,6.727542304249694 +1.631113241380004,0.7829612037013028,1.4761891913388827,1.570755306583805 +-0.25665876738266513,-1.5513454724715776,-1.2777089856188693,-0.7497144820358775 +1.613248685596358,1.7340615070817575,1.701365449146443,1.3513265538013948 +2.4303661352193653,1.7340615070817575,2.2662416091255366,1.7676186374592735 +0.6441618263907917,0.7829612037013028,0.25432718961037304,0.7461319230554601 +-4.228828865372143,-4.339380045700343,-4.104814639901052,-4.2600649948301355 +-3.843126375937169,-4.339380045700343,-3.6688578553761966,-3.4247044005775598 +1.0491437505785097,0.7829612037013028,-0.09835331046329601,0.5146975219915633 +-2.4531495909852636,-1.5513454724715776,-1.9504228619236361,-2.2903184452192407 +5.644033124047758,4.587413399100595,5.121782776091066,5.790888921291713 +-2.9258317234227125,-3.3432405724417764,-3.1227871611703937,-2.868478980566626 +-3.2501102047917234,-3.3432405724417764,-3.1227871611703937,-2.868478980566626 +-4.021945070669235,-4.339380045700343,-3.6688578553761966,-3.5176957243764866 +-1.6515657127272214,-3.025574053843865,-2.88365816162571,-2.072000555763744 +-0.16123721575693717,0.7829612037013028,0.38881232949047506,0.3188128978609615 +0.7702689119165229,0.7829612037013028,0.661536104096642,0.5300018973383436 +3.806583329076842,4.019798593894754,3.564288400025178,3.8585467170918344 +1.4470984055029208,0.7829612037013028,0.3733803732340081,1.5274697629800138 +-2.0097304453363454,-3.025574053843865,-2.2335476362719087,-2.3310026901186074 +-1.1038321325505507,-1.5513454724715776,-1.2777089856188693,-1.0908277426198256 +0.7575335023199609,0.7829612037013028,0.8379825835296415,0.5609099999099947 +0.5491820808925653,-0.17286001207987395,-0.3003167401365543,0.1080940522042256 +-3.160808086081791,-3.025574053843865,-3.572923544234338,-3.695475069873099 +-6.849679747770614,-6.558872933068393,-5.540346982806637,-5.714238735108081 +-2.721655906924753,-3.025574053843865,-2.857103436975195,-2.6197637148628807 +7.076584763750515,4.019798593894754,4.728583892005671,5.928596671480644 +1.6861113353726331,0.7829612037013028,1.4761891913388827,1.3948339401356011 +0.989171449301451,1.7340615070817575,1.3496566838538135,1.0102132932174468 +4.259119829902918,4.019798593894754,3.564288400025178,3.5112792830349058 +3.9878520254237158,4.587413399100595,4.402510058651032,4.301496507476102 +8.05537246499182,8.311136903189661,7.201878039666836,7.6788913544623965 +-7.95550804615306,-6.558872933068393,-6.052742811723908,-6.743268492818798 +2.189558050611199,0.7829612037013028,1.4761891913388827,1.9252416341411072 +-3.8751500271226234,-3.3432405724417764,-3.1227871611703937,-3.3547859184344153 +2.47777814031454,0.7829612037013028,2.0924732054503736,1.7819443060611873 +0.8582933002418842,0.7829612037013028,0.661536104096642,0.5300018973383436 +-1.7182315200304354,-1.5513454724715776,-1.6437232249048663,-1.4057173607400517 +-7.088627367379241,-6.558872933068393,-6.960305887457423,-7.033602563950763 +-2.0504984508160695,-3.025574053843865,-2.88365816162571,-2.1408372715187363 +-4.16153822302006,-3.025574053843865,-3.8269338070992065,-4.523796397004544 +-7.693790597428782,-6.558872933068393,-6.708641281236406,-6.983619253576441 +-3.634818210013233,-1.7909353004230013,-2.0378476758861166,-3.0079622328719537 diff --git a/doc/Programs/Regression/Logreg.txt b/doc/Programs/Regression/Logreg.txt deleted file mode 100644 index c9d5f1e77..000000000 --- a/doc/Programs/Regression/Logreg.txt +++ /dev/null @@ -1,313 +0,0 @@ -From: Morten Hjorth-Jensen -Subject: Logreg -Date: May 29, 2025 at 11:43:36 AM GMT+2 -To: Morten Hjorth-Jensen - -
 -Logistic Regression from Scratch -
 -
 -Logistic regression is a foundational binary classification technique that models the probability of a sample belonging to class 1 using the sigmoid (logistic) function. In binary classification, it outputs a probability  between 0 and 1, typically thresholded at 0.5 to decide class labels . For multiclass problems (more than two classes), we generalize logistic regression via the softmax function (also known as multinomial logistic regression), which produces a probability distribution over all classes . In both cases, training involves minimizing the cross-entropy loss (also called log-loss) using gradient descent . Accuracy is commonly used to evaluate classification performance (the fraction of correctly predicted samples) , alongside the cross-entropy value to measure model fit. -
 -Our goal is to build an object-oriented Python implementation of logistic regression without scikit-learn, using only standard libraries like NumPy and CSV. The implementation will support both binary and multiclass modes, use gradient descent for training, and include methods for prediction and evaluation. It will also include synthetic data generation (to test the model) and functionality to export predictions and true labels to CSV. -
 -
 -Design and Features -
 -
 - • Class Structure: A LogisticRegression class encapsulates model parameters (weights, bias) and methods for training (fit) and prediction (predict, predict_proba). The design handles binary vs multiclass internally by checking how many unique target labels are present. - • Probability Functions: For binary classification we use the sigmoid function. For multiclass, we use the softmax function, defined as \hat{y}_i = \frac{\exp(o_i)}{\sum_j \exp(o_j)}, which produces a valid probability distribution over classes . - • Training with Gradient Descent: We implement batch gradient descent on the cross-entropy loss. For binary logistic regression the loss is
-\frac{1}{N}\sum_i [y_i\log(p_i) + (1-y_i)\log(1-p_i)]
where p_i=\sigma(\mathbf{w}^T\mathbf{x}i) . For multiclass, the loss (categorical cross-entropy) for one-hot true labels y and predicted probabilities \hat{y} is
-\frac{1}{N}\sum_i \sum{k=1}^C y_{ik}\log(\hat{y}_{ik})
. We compute gradients of these losses w.r.t. the weights and iteratively update the model. - • Prediction: The predict_proba method returns probabilities (sigmoid or softmax). The predict method applies a threshold (0.5) for binary or takes argmax for multiclass. - • Metrics: We include functions for accuracy and cross-entropy loss. Accuracy is defined as the fraction of correct predictions . We compute binary cross-entropy or categorical cross-entropy to quantify model performance on data. - • Synthetic Data: We provide functions to generate synthetic datasets: (a) Binary data with two Gaussian clusters for two classes; (b) Multiclass data with several clusters, one per class. This allows testing the model. - • CSV Export: Using Python’s csv module, we can save the true labels and predicted labels (or probabilities) to CSV files for external analysis. -
 -
 -These components together form a clean, modular implementation. The following sections detail the code with explanations. -
 -
 -LogisticRegression Class Implementation -
 -
 -Below is the LogisticRegression class. It initializes parameters (learning rate, epochs, etc.), and defines private methods for adding an intercept term, sigmoid, and softmax. The fit method detects binary vs multiclass and runs gradient descent accordingly, updating weights. The predict_prob and predict methods compute probabilities and class labels. -import numpy as np -
 -class LogisticRegression: -    """ -    Logistic Regression for binary and multiclass classification. -    """ -    def __init__(self, lr=0.01, epochs=1000, fit_intercept=True, verbose=False): -        self.lr = lr                  # Learning rate for gradient descent -        self.epochs = epochs          # Number of iterations -        self.fit_intercept = fit_intercept  # Whether to add intercept (bias) -        self.verbose = verbose        # Print loss during training if True -        self.weights = None -        self.multi_class = False      # Will be determined at fit time -
 -    def _add_intercept(self, X): -        """Add intercept term (column of ones) to feature matrix.""" -        intercept = np.ones((X.shape[0], 1)) -        return np.concatenate((intercept, X), axis=1) -
 -    def _sigmoid(self, z): -        """Sigmoid function for binary logistic.""" -        return 1 / (1 + np.exp(-z)) -
 -    def _softmax(self, Z): -        """Softmax function for multiclass logistic.""" -        exp_Z = np.exp(Z - np.max(Z, axis=1, keepdims=True)) -        return exp_Z / np.sum(exp_Z, axis=1, keepdims=True) -
 -    def fit(self, X, y): -        """ -        Train the logistic regression model using gradient descent. -        Supports binary (sigmoid) and multiclass (softmax) based on y. -        """ -        X = np.array(X) -        y = np.array(y) -        n_samples, n_features = X.shape -
 -        # Add intercept if needed -        if self.fit_intercept: -            X = self._add_intercept(X) -            n_features += 1 -
 -        # Determine classes and mode (binary vs multiclass) -        unique_classes = np.unique(y) -        if len(unique_classes) > 2: -            self.multi_class = True -        else: -            self.multi_class = False -
 -        # ----- Multiclass case ----- -        if self.multi_class: -            n_classes = len(unique_classes) -            # Map original labels to 0...n_classes-1 -            class_to_index = {c: idx for idx, c in enumerate(unique_classes)} -            y_indices = np.array([class_to_index[c] for c in y]) -            # Initialize weight matrix (features x classes) -            self.weights = np.zeros((n_features, n_classes)) -
 -            # One-hot encode y -            Y_onehot = np.zeros((n_samples, n_classes)) -            Y_onehot[np.arange(n_samples), y_indices] = 1 -
 -            # Gradient descent -            for epoch in range(self.epochs): -                scores = X.dot(self.weights)          # Linear scores (n_samples x n_classes) -                probs = self._softmax(scores)        # Probabilities (n_samples x n_classes) -                # Compute gradient (features x classes) -                gradient = (1 / n_samples) * X.T.dot(probs - Y_onehot) -                # Update weights -                self.weights -= self.lr * gradient -
 -                if self.verbose and epoch % 100 == 0: -                    # Compute current loss (categorical cross-entropy) -                    loss = -np.sum(Y_onehot * np.log(probs + 1e-15)) / n_samples -                    print(f"[Epoch {epoch}] Multiclass loss: {loss:.4f}") -
 -        # ----- Binary case ----- -        else: -            # Convert y to 0/1 if not already -            if not np.array_equal(unique_classes, [0, 1]): -                # Map the two classes to 0 and 1 -                class0, class1 = unique_classes -                y_binary = np.where(y == class1, 1, 0) -            else: -                y_binary = y.copy().astype(int) -
 -            # Initialize weights vector (features,) -            self.weights = np.zeros(n_features) -
 -            # Gradient descent -            for epoch in range(self.epochs): -                linear_model = X.dot(self.weights)     # (n_samples,) -                probs = self._sigmoid(linear_model)   # (n_samples,) -                # Gradient for binary cross-entropy -                gradient = (1 / n_samples) * X.T.dot(probs - y_binary) -                self.weights -= self.lr * gradient -
 -                if self.verbose and epoch % 100 == 0: -                    # Compute binary cross-entropy loss -                    loss = -np.mean( -                        y_binary * np.log(probs + 1e-15) +  -                        (1 - y_binary) * np.log(1 - probs + 1e-15) -                    ) -                    print(f"[Epoch {epoch}] Binary loss: {loss:.4f}") -
 -    def predict_prob(self, X): -        """ -        Compute probability estimates. Returns a 1D array for binary or -        a 2D array (n_samples x n_classes) for multiclass. -        """ -        X = np.array(X) -        # Add intercept if the model used it -        if self.fit_intercept: -            X = self._add_intercept(X) -        scores = X.dot(self.weights) -        if self.multi_class: -            return self._softmax(scores) -        else: -            return self._sigmoid(scores) -
 -    def predict(self, X): -        """ -        Predict class labels for samples in X. -        Returns integer class labels (0,1 for binary, or 0...C-1 for multiclass). -        """ -        probs = self.predict_prob(X) -        if self.multi_class: -            # Choose class with highest probability -            return np.argmax(probs, axis=1) -        else: -            # Threshold at 0.5 for binary -            return (probs >= 0.5).astype(int) -The class implements the sigmoid and softmax internally. During fit(), we check the number of classes: if more than 2, we set self.multi_class=True and perform multinomial logistic regression. We one-hot encode the target vector and update a weight matrix with softmax probabilities. Otherwise, we do standard binary logistic regression, converting labels to 0/1 if needed and updating a weight vector. In both cases we use batch gradient descent on the cross-entropy loss (we add a small epsilon 1e-15 to logs for numerical stability). Progress (loss) can be printed if verbose=True. -
 -
 -Evaluation Metrics -
 -
 -We define helper functions for accuracy and cross-entropy loss. Accuracy is the fraction of correct predictions . For loss, we compute the appropriate cross-entropy: -def accuracy_score(y_true, y_pred): -    """Accuracy = (# correct predictions) / (total samples).""" -    y_true = np.array(y_true) -    y_pred = np.array(y_pred) -    return np.mean(y_true == y_pred) -
 -def binary_cross_entropy(y_true, y_prob): -    """ -    Binary cross-entropy loss. -    y_true: true binary labels (0 or 1), y_prob: predicted probabilities for class 1. -    """ -    y_true = np.array(y_true) -    y_prob = np.clip(np.array(y_prob), 1e-15, 1-1e-15) -    return -np.mean(y_true * np.log(y_prob) + (1 - y_true) * np.log(1 - y_prob)) -
 -def categorical_cross_entropy(y_true, y_prob): -    """ -    Categorical cross-entropy loss for multiclass. -    y_true: true labels (0...C-1), y_prob: array of predicted probabilities (n_samples x C). -    """ -    y_true = np.array(y_true, dtype=int) -    y_prob = np.clip(np.array(y_prob), 1e-15, 1-1e-15) -    # One-hot encode true labels -    n_samples, n_classes = y_prob.shape -    one_hot = np.zeros_like(y_prob) -    one_hot[np.arange(n_samples), y_true] = 1 -    # Compute cross-entropy -    loss_vec = -np.sum(one_hot * np.log(y_prob), axis=1) -    return np.mean(loss_vec) -The binary cross-entropy matches the formula  , and categorical cross-entropy aligns with the standard softmax loss . We clip probabilities to avoid log(0). Accuracy is straightforward (correct/total) . -
 -
 -Synthetic Data Generation -
 -
 -To test the model, we generate synthetic datasets: -
 - • Binary classification data: Create two Gaussian clusters in 2D. For example, class 0 around mean [-2,-2] and class 1 around [2,2]. - • Multiclass data: Create several Gaussian clusters (one per class) spread out in feature space. -
 -
 -The code below demonstrates simple generators: -import numpy as np -
 -def generate_binary_data(n_samples=100, n_features=2, random_state=None): -    """ -    Generate synthetic binary classification data. -    Returns (X, y) where X is (n_samples x n_features), y in {0,1}. -    """ -    rng = np.random.RandomState(random_state) -    # Half samples for class 0, half for class 1 -    n0 = n_samples // 2 -    n1 = n_samples - n0 -    # Class 0 around mean -2, class 1 around +2 -    mean0 = -2 * np.ones(n_features) -    mean1 =  2 * np.ones(n_features) -    X0 = rng.randn(n0, n_features) + mean0 -    X1 = rng.randn(n1, n_features) + mean1 -    X = np.vstack((X0, X1)) -    y = np.array([0]*n0 + [1]*n1) -    return X, y -
 -def generate_multiclass_data(n_samples=150, n_features=2, n_classes=3, random_state=None): -    """ -    Generate synthetic multiclass data with n_classes Gaussian clusters. -    """ -    rng = np.random.RandomState(random_state) -    X = [] -    y = [] -    samples_per_class = n_samples // n_classes -    for cls in range(n_classes): -        # Random cluster center for each class -        center = rng.uniform(-5, 5, size=n_features) -        Xi = rng.randn(samples_per_class, n_features) + center -        yi = [cls] * samples_per_class -        X.append(Xi) -        y.extend(yi) -    X = np.vstack(X) -    y = np.array(y) -    return X, y -These functions use NumPy to create normally-distributed points. We fix cluster centers to separate classes. Each class’s points are in a distinct region, making the classification problem learnable by logistic regression. -
 -
 -Demo: Training and Evaluation -
 -
 -Finally, we test the implementation on synthetic data. We train on both binary and multiclass cases, then evaluate accuracy and loss, and export results to CSV. -# Generate and test on binary data -X_bin, y_bin = generate_binary_data(n_samples=200, n_features=2, random_state=42) -model_bin = LogisticRegression(lr=0.1, epochs=1000) -model_bin.fit(X_bin, y_bin) -y_prob_bin = model_bin.predict_prob(X_bin)      # probabilities for class 1 -y_pred_bin = model_bin.predict(X_bin)           # predicted classes 0 or 1 -
 -acc_bin = accuracy_score(y_bin, y_pred_bin) -loss_bin = binary_cross_entropy(y_bin, y_prob_bin) -print(f"Binary Classification - Accuracy: {acc_bin:.2f}, Cross-Entropy Loss: {loss_bin:.2f}") -For multiclass: -# Generate and test on multiclass data -X_multi, y_multi = generate_multiclass_data(n_samples=300, n_features=2, n_classes=3, random_state=1) -model_multi = LogisticRegression(lr=0.1, epochs=1000) -model_multi.fit(X_multi, y_multi) -y_prob_multi = model_multi.predict_prob(X_multi)     # (n_samples x 3) probabilities -y_pred_multi = model_multi.predict(X_multi)          # predicted labels 0,1,2 -
 -acc_multi = accuracy_score(y_multi, y_pred_multi) -loss_multi = categorical_cross_entropy(y_multi, y_prob_multi) -print(f"Multiclass Classification - Accuracy: {acc_multi:.2f}, Cross-Entropy Loss: {loss_multi:.2f}") -These print statements show how well the model fits the training data. -
 -
 -CSV Export -
 -
 -To save predictions and true labels, we use Python’s csv module. For example: -import csv -
 -# Export binary results -with open('binary_results.csv', mode='w', newline='') as f: -    writer = csv.writer(f) -    writer.writerow(["TrueLabel", "PredictedLabel"]) -    for true, pred in zip(y_bin, y_pred_bin): -        writer.writerow([true, pred]) -
 -# Export multiclass results -with open('multiclass_results.csv', mode='w', newline='') as f: -    writer = csv.writer(f) -    writer.writerow(["TrueLabel", "PredictedLabel"]) -    for true, pred in zip(y_multi, y_pred_multi): -        writer.writerow([true, pred]) -This writes two CSV files with columns for true and predicted labels. One can later analyze or plot these results externally. -
 -
 -Summary -
 -
 -The above implementation provides a clear, object-oriented logistic regression model in pure Python. It handles both binary and multiclass cases by using sigmoid and softmax functions respectively. Training uses batch gradient descent to minimize cross-entropy loss . We include methods for prediction, accuracy calculation , and loss computation, as well as utilities for generating synthetic data and exporting results. This modular design can be extended (e.g., adding regularization or optimization improvements) but already demonstrates the key mechanics of logistic regression end-to-end. -
 -References: Logistic regression concepts and loss functions ; accuracy metric . - -Morten Hjorth-Jensen, Michigan State University and University of Oslo, Norway. http://mhjgit.github.io/info/doc/web diff --git a/doc/Programs/Regression/binary_results.csv b/doc/Programs/Regression/binary_results.csv new file mode 100644 index 000000000..1a5f8e043 --- /dev/null +++ b/doc/Programs/Regression/binary_results.csv @@ -0,0 +1,201 @@ +TrueLabel,PredictedLabel +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,0 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 diff --git a/doc/Programs/Regression/logreg.py b/doc/Programs/Regression/logreg.py index e69de29bb..6c5b57b47 100644 --- a/doc/Programs/Regression/logreg.py +++ b/doc/Programs/Regression/logreg.py @@ -0,0 +1,257 @@ +import numpy as np + +class LogisticRegression: + """ + Logistic Regression for binary and multiclass classification. + """ + def __init__(self, lr=0.01, epochs=1000, fit_intercept=True, verbose=False): + self.lr = lr # Learning rate for gradient descent + self.epochs = epochs # Number of iterations + self.fit_intercept = fit_intercept # Whether to add intercept (bias) + self.verbose = verbose # Print loss during training if True + self.weights = None + self.multi_class = False # Will be determined at fit time + + def _add_intercept(self, X): + """Add intercept term (column of ones) to feature matrix.""" + intercept = np.ones((X.shape[0], 1)) + return np.concatenate((intercept, X), axis=1) + + def _sigmoid(self, z): + """Sigmoid function for binary logistic.""" + return 1 / (1 + np.exp(-z)) + + def _softmax(self, Z): + """Softmax function for multiclass logistic.""" + exp_Z = np.exp(Z - np.max(Z, axis=1, keepdims=True)) + return exp_Z / np.sum(exp_Z, axis=1, keepdims=True) + + def fit(self, X, y): + """ + Train the logistic regression model using gradient descent. + Supports binary (sigmoid) and multiclass (softmax) based on y. + """ + X = np.array(X) + y = np.array(y) + n_samples, n_features = X.shape + + # Add intercept if needed + if self.fit_intercept: + X = self._add_intercept(X) + n_features += 1 + + # Determine classes and mode (binary vs multiclass) + unique_classes = np.unique(y) + if len(unique_classes) > 2: + self.multi_class = True + else: + self.multi_class = False + + # ----- Multiclass case ----- + if self.multi_class: + n_classes = len(unique_classes) + # Map original labels to 0...n_classes-1 + class_to_index = {c: idx for idx, c in enumerate(unique_classes)} + y_indices = np.array([class_to_index[c] for c in y]) + # Initialize weight matrix (features x classes) + self.weights = np.zeros((n_features, n_classes)) + + # One-hot encode y + Y_onehot = np.zeros((n_samples, n_classes)) + Y_onehot[np.arange(n_samples), y_indices] = 1 + + # Gradient descent + for epoch in range(self.epochs): + scores = X.dot(self.weights) # Linear scores (n_samples x n_classes) + probs = self._softmax(scores) # Probabilities (n_samples x n_classes) + # Compute gradient (features x classes) + gradient = (1 / n_samples) * X.T.dot(probs - Y_onehot) + # Update weights + self.weights -= self.lr * gradient + + if self.verbose and epoch % 100 == 0: + # Compute current loss (categorical cross-entropy) + loss = -np.sum(Y_onehot * np.log(probs + 1e-15)) / n_samples + print(f"[Epoch {epoch}] Multiclass loss: {loss:.4f}") + + # ----- Binary case ----- + else: + # Convert y to 0/1 if not already + if not np.array_equal(unique_classes, [0, 1]): + # Map the two classes to 0 and 1 + class0, class1 = unique_classes + y_binary = np.where(y == class1, 1, 0) + else: + y_binary = y.copy().astype(int) + + # Initialize weights vector (features,) + self.weights = np.zeros(n_features) + + # Gradient descent + for epoch in range(self.epochs): + linear_model = X.dot(self.weights) # (n_samples,) + probs = self._sigmoid(linear_model) # (n_samples,) + # Gradient for binary cross-entropy + gradient = (1 / n_samples) * X.T.dot(probs - y_binary) + self.weights -= self.lr * gradient + + if self.verbose and epoch % 100 == 0: + # Compute binary cross-entropy loss + loss = -np.mean( + y_binary * np.log(probs + 1e-15) + + (1 - y_binary) * np.log(1 - probs + 1e-15) + ) + print(f"[Epoch {epoch}] Binary loss: {loss:.4f}") + + def predict_prob(self, X): + """ + Compute probability estimates. Returns a 1D array for binary or + a 2D array (n_samples x n_classes) for multiclass. + """ + X = np.array(X) + # Add intercept if the model used it + if self.fit_intercept: + X = self._add_intercept(X) + scores = X.dot(self.weights) + if self.multi_class: + return self._softmax(scores) + else: + return self._sigmoid(scores) + + def predict(self, X): + """ + Predict class labels for samples in X. + Returns integer class labels (0,1 for binary, or 0...C-1 for multiclass). + """ + probs = self.predict_prob(X) + if self.multi_class: + # Choose class with highest probability + return np.argmax(probs, axis=1) + else: + # Threshold at 0.5 for binary + return (probs >= 0.5).astype(int) + +""" +The class implements the sigmoid and softmax internally. During fit(), we check the number of classes: if more than 2, we set self.multi_class=True and perform multinomial logistic regression. We one-hot encode the target vector and update a weight matrix with softmax probabilities. Otherwise, we do standard binary logistic regression, converting labels to 0/1 if needed and updating a weight vector. In both cases we use batch gradient descent on the cross-entropy loss (we add a small epsilon 1e-15 to logs for numerical stability). Progress (loss) can be printed if verbose=True. +""" + +# Evaluation Metrics +#We define helper functions for accuracy and cross-entropy loss. Accuracy is the fraction of correct predictions . For loss, we compute the appropriate cross-entropy: + +def accuracy_score(y_true, y_pred): + """Accuracy = (# correct predictions) / (total samples).""" + y_true = np.array(y_true) + y_pred = np.array(y_pred) + return np.mean(y_true == y_pred) + +def binary_cross_entropy(y_true, y_prob): + """ + Binary cross-entropy loss. + y_true: true binary labels (0 or 1), y_prob: predicted probabilities for class 1. + """ + y_true = np.array(y_true) + y_prob = np.clip(np.array(y_prob), 1e-15, 1-1e-15) + return -np.mean(y_true * np.log(y_prob) + (1 - y_true) * np.log(1 - y_prob)) + +def categorical_cross_entropy(y_true, y_prob): + """ + Categorical cross-entropy loss for multiclass. + y_true: true labels (0...C-1), y_prob: array of predicted probabilities (n_samples x C). + """ + y_true = np.array(y_true, dtype=int) + y_prob = np.clip(np.array(y_prob), 1e-15, 1-1e-15) + # One-hot encode true labels + n_samples, n_classes = y_prob.shape + one_hot = np.zeros_like(y_prob) + one_hot[np.arange(n_samples), y_true] = 1 + # Compute cross-entropy + loss_vec = -np.sum(one_hot * np.log(y_prob), axis=1) + return np.mean(loss_vec) + + +""" +Synthetic data generation +Binary classification data: Create two Gaussian clusters in 2D. For example, class 0 around mean [-2,-2] and class 1 around [2,2]. +Multiclass data: Create several Gaussian clusters (one per class) spread out in feature space. +""" + + +import numpy as np + +def generate_binary_data(n_samples=100, n_features=2, random_state=None): + """ + Generate synthetic binary classification data. + Returns (X, y) where X is (n_samples x n_features), y in {0,1}. + """ + rng = np.random.RandomState(random_state) + # Half samples for class 0, half for class 1 + n0 = n_samples // 2 + n1 = n_samples - n0 + # Class 0 around mean -2, class 1 around +2 + mean0 = -2 * np.ones(n_features) + mean1 = 2 * np.ones(n_features) + X0 = rng.randn(n0, n_features) + mean0 + X1 = rng.randn(n1, n_features) + mean1 + X = np.vstack((X0, X1)) + y = np.array([0]*n0 + [1]*n1) + return X, y + +def generate_multiclass_data(n_samples=150, n_features=2, n_classes=3, random_state=None): + """ + Generate synthetic multiclass data with n_classes Gaussian clusters. + """ + rng = np.random.RandomState(random_state) + X = [] + y = [] + samples_per_class = n_samples // n_classes + for cls in range(n_classes): + # Random cluster center for each class + center = rng.uniform(-5, 5, size=n_features) + Xi = rng.randn(samples_per_class, n_features) + center + yi = [cls] * samples_per_class + X.append(Xi) + y.extend(yi) + X = np.vstack(X) + y = np.array(y) + return X, y + + +# Generate and test on binary data +X_bin, y_bin = generate_binary_data(n_samples=200, n_features=2, random_state=42) +model_bin = LogisticRegression(lr=0.1, epochs=1000) +model_bin.fit(X_bin, y_bin) +y_prob_bin = model_bin.predict_prob(X_bin) # probabilities for class 1 +y_pred_bin = model_bin.predict(X_bin) # predicted classes 0 or 1 + +acc_bin = accuracy_score(y_bin, y_pred_bin) +loss_bin = binary_cross_entropy(y_bin, y_prob_bin) +print(f"Binary Classification - Accuracy: {acc_bin:.2f}, Cross-Entropy Loss: {loss_bin:.2f}") +#For multiclass: +# Generate and test on multiclass data +X_multi, y_multi = generate_multiclass_data(n_samples=300, n_features=2, n_classes=3, random_state=1) +model_multi = LogisticRegression(lr=0.1, epochs=1000) +model_multi.fit(X_multi, y_multi) +y_prob_multi = model_multi.predict_prob(X_multi) # (n_samples x 3) probabilities +y_pred_multi = model_multi.predict(X_multi) # predicted labels 0,1,2 + +acc_multi = accuracy_score(y_multi, y_pred_multi) +loss_multi = categorical_cross_entropy(y_multi, y_prob_multi) +print(f"Multiclass Classification - Accuracy: {acc_multi:.2f}, Cross-Entropy Loss: {loss_multi:.2f}") + +# CSV Export +import csv + +# Export binary results +with open('binary_results.csv', mode='w', newline='') as f: + writer = csv.writer(f) + writer.writerow(["TrueLabel", "PredictedLabel"]) + for true, pred in zip(y_bin, y_pred_bin): + writer.writerow([true, pred]) + +# Export multiclass results +with open('multiclass_results.csv', mode='w', newline='') as f: + writer = csv.writer(f) + writer.writerow(["TrueLabel", "PredictedLabel"]) + for true, pred in zip(y_multi, y_pred_multi): + writer.writerow([true, pred]) + diff --git a/doc/Programs/Regression/multiclass_results.csv b/doc/Programs/Regression/multiclass_results.csv new file mode 100644 index 000000000..9ffbe5203 --- /dev/null +++ b/doc/Programs/Regression/multiclass_results.csv @@ -0,0 +1,301 @@ +TrueLabel,PredictedLabel +0,0 +0,1 +0,0 +0,0 +0,1 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,1 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,1 +0,0 +0,2 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +0,0 +1,1 +1,1 +1,0 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,0 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,0 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,1 +1,0 +1,1 +1,1 +1,0 +1,1 +1,1 +1,1 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2 +2,2