Feature Engineering, Scaling & Regularization

import numpy as np

class StandardScalerFromScratch:
    """Z-score standardization: z = (x - mean) / std"""
    def fit(self, X):
        self.mean_ = np.mean(X, axis=0)
        self.std_ = np.std(X, axis=0)
        # Avoid division by zero for constant features
        self.std_[self.std_ == 0] = 1.0
        return self

    def transform(self, X):
        return (X - self.mean_) / self.std_

    def inverse_transform(self, Z):
        return Z * self.std_ + self.mean_

    def fit_transform(self, X):
        return self.fit(X).transform(X)


class MinMaxScalerFromScratch:
    """Min-Max scaling to [0, 1]: x' = (x - min) / (max - min)"""
    def __init__(self, feature_range=(0, 1)):
        self.min_val, self.max_val = feature_range

    def fit(self, X):
        self.data_min_ = np.min(X, axis=0)
        self.data_max_ = np.max(X, axis=0)
        self.data_range_ = self.data_max_ - self.data_min_
        self.data_range_[self.data_range_ == 0] = 1.0
        return self

    def transform(self, X):
        X_std = (X - self.data_min_) / self.data_range_
        return X_std * (self.max_val - self.min_val) + self.min_val

    def fit_transform(self, X):
        return self.fit(X).transform(X)


class RobustScalerFromScratch:
    """Robust scaling using median and IQR"""
    def fit(self, X):
        self.median_ = np.median(X, axis=0)
        Q1 = np.percentile(X, 25, axis=0)
        Q3 = np.percentile(X, 75, axis=0)
        self.iqr_ = Q3 - Q1
        self.iqr_[self.iqr_ == 0] = 1.0
        return self

    def transform(self, X):
        return (X - self.median_) / self.iqr_

    def fit_transform(self, X):
        return self.fit(X).transform(X)


# ---- Demo ----
np.random.seed(42)
X = np.array([[1200, 2], [1800, 3], [2400, 4],
              [900,  1], [3000, 5]])

print("Original:\n", X)
print("StandardScaler:\n", StandardScalerFromScratch().fit_transform(X))
print("MinMaxScaler:\n", MinMaxScalerFromScratch().fit_transform(X))
print("RobustScaler:\n", RobustScalerFromScratch().fit_transform(X))
        

import numpy as np

class RidgeRegressionFromScratch:
    """
    Ridge Regression: w = (X^T X + λI)^{-1} X^T y
    Implements closed-form solution derived in Section 5.1.
    """
    def __init__(self, alpha=1.0, fit_intercept=True):
        self.alpha = alpha
        self.fit_intercept = fit_intercept

    def fit(self, X, y):
        X = np.array(X, dtype=np.float64)
        y = np.array(y, dtype=np.float64)

        if self.fit_intercept:
            # Center X and y to handle intercept
            self.X_mean_ = np.mean(X, axis=0)
            self.y_mean_ = np.mean(y)
            X_c = X - self.X_mean_
            y_c = y - self.y_mean_
        else:
            X_c, y_c = X, y

        n, p = X_c.shape
        # Closed-form: w = (X^T X + alpha * I)^(-1) X^T y
        I = np.eye(p)
        XtX = X_c.T @ X_c
        Xty = X_c.T @ y_c
        self.coef_ = np.linalg.solve(XtX + self.alpha * I, Xty)

        if self.fit_intercept:
            self.intercept_ = self.y_mean_ - self.X_mean_ @ self.coef_
        else:
            self.intercept_ = 0.0

        return self

    def predict(self, X):
        X = np.array(X, dtype=np.float64)
        return X @ self.coef_ + self.intercept_

    def score(self, X, y):
        """R² score"""
        y_pred = self.predict(X)
        ss_res = np.sum((y - y_pred) ** 2)
        ss_tot = np.sum((y - np.mean(y)) ** 2)
        return 1 - ss_res / ss_tot


# ---- Demo ----
np.random.seed(42)
X = np.random.randn(100, 5)
true_w = np.array([3, 0, -2, 0, 1.5])  # Only 3 features matter
y = X @ true_w + np.random.randn(100) * 0.5

for alpha in [0.01, 0.1, 1.0, 10.0]:
    model = RidgeRegressionFromScratch(alpha=alpha)
    model.fit(X, y)
    print(f"λ={alpha:5.2f}  coefs={np.round(model.coef_, 2)}  R²={model.score(X, y):.4f}")
        

import numpy as np

def soft_threshold(rho, lam):
    """Soft-thresholding operator S(ρ, λ) = sign(ρ) · max(|ρ| - λ, 0)"""
    if rho > lam:
        return rho - lam
    elif rho < -lam:
        return rho + lam
    else:
        return 0.0


class LassoFromScratch:
    """
    Lasso Regression using coordinate descent.
    Implements the soft-thresholding update derived in Section 5.2.
    """
    def __init__(self, alpha=1.0, max_iter=1000, tol=1e-6):
        self.alpha = alpha
        self.max_iter = max_iter
        self.tol = tol

    def fit(self, X, y):
        X = np.array(X, dtype=np.float64)
        y = np.array(y, dtype=np.float64)
        n, p = X.shape

        # Center data
        self.X_mean_ = np.mean(X, axis=0)
        self.y_mean_ = np.mean(y)
        X_c = X - self.X_mean_
        y_c = y - self.y_mean_

        # Initialize coefficients
        w = np.zeros(p)
        # Precompute z_j = sum(x_ij^2) for each feature
        z = np.sum(X_c ** 2, axis=0)

        for iteration in range(self.max_iter):
            w_old = w.copy()
            for j in range(p):
                # Compute partial residual (excluding feature j)
                r_j = y_c - X_c @ w + X_c[:, j] * w[j]
                # ρ_j = correlation of feature j with residual
                rho_j = X_c[:, j] @ r_j
                # Apply soft-thresholding
                w[j] = soft_threshold(rho_j, n * self.alpha) / z[j]

            # Check convergence
            if np.max(np.abs(w - w_old)) < self.tol:
                break

        self.coef_ = w
        self.intercept_ = self.y_mean_ - self.X_mean_ @ self.coef_
        self.n_iter_ = iteration + 1
        return self

    def predict(self, X):
        return np.array(X) @ self.coef_ + self.intercept_


# ---- Demo: Lasso automatic feature selection ----
np.random.seed(42)
X = np.random.randn(100, 10)
true_w = np.array([3, 0, 0, -2, 0, 1.5, 0, 0, 0, 0.8])
y = X @ true_w + np.random.randn(100) * 0.3

model = LassoFromScratch(alpha=0.05)
model.fit(X, y)
print("True weights: ", true_w)
print("Lasso weights:", np.round(model.coef_, 3))
print("Non-zero features:", np.sum(model.coef_ != 0), "/ 10")
# Expected: Lasso correctly identifies the 4 non-zero features!
        

import numpy as np

class KNNImputerFromScratch:
    """
    K-Nearest Neighbors Imputation.
    For each missing value, find K nearest neighbors based on
    non-missing features, then impute with their weighted average.
    """
    def __init__(self, n_neighbors=5):
        self.n_neighbors = n_neighbors

    def fit_transform(self, X):
        X = np.array(X, dtype=np.float64)
        X_imputed = X.copy()
        n, p = X.shape

        for i in range(n):
            for j in range(p):
                if np.isnan(X[i, j]):
                    # Find features that are NOT missing for row i
                    valid_features = ~np.isnan(X[i, :])
                    # Find rows that have value for feature j
                    candidate_rows = ~np.isnan(X[:, j])
                    candidate_rows[i] = False  # exclude self

                    if np.sum(candidate_rows) == 0:
                        # Fallback: use column mean
                        col_values = X[~np.isnan(X[:, j]), j]
                        X_imputed[i, j] = np.mean(col_values)
                        continue

                    # Compute distances using valid features
                    distances = []
                    for k in range(n):
                        if not candidate_rows[k]:
                            continue
                        # Only use features present in both rows
                        shared = valid_features & ~np.isnan(X[k, :])
                        if np.sum(shared) == 0:
                            continue
                        dist = np.sqrt(np.sum((X[i, shared] - X[k, shared])**2))
                        distances.append((dist, k))

                    # Sort by distance, take K nearest
                    distances.sort()
                    neighbors = distances[:self.n_neighbors]

                    if len(neighbors) == 0:
                        col_values = X[~np.isnan(X[:, j]), j]
                        X_imputed[i, j] = np.mean(col_values)
                    else:
                        # Weighted average (inverse distance weighting)
                        weights = [1/(d + 1e-8) for d, _ in neighbors]
                        values = [X[idx, j] for _, idx in neighbors]
                        X_imputed[i, j] = np.average(values, weights=weights)

        return X_imputed


# ---- Demo ----
X = np.array([
    [25, 40000, 4.2],
    [32, np.nan, 3.8],
    [np.nan, 55000, np.nan],
    [28, 35000, 4.5],
    [45, 80000, 3.9]
])
imputer = KNNImputerFromScratch(n_neighbors=2)
X_filled = imputer.fit_transform(X)
print("Imputed data:\n", np.round(X_filled, 1))
        

import numpy as np

def permutation_importance(model, X, y, n_repeats=10,
                            scoring=None, random_state=42):
    """
    Compute permutation importance for each feature.
    
    Algorithm:
    1. Compute baseline score with original data.
    2. For each feature j:
       a. Shuffle feature j's values randomly (n_repeats times).
       b. Compute new score with shuffled feature.
       c. Importance = baseline_score - shuffled_score.
    """
    rng = np.random.RandomState(random_state)

    if scoring is None:
        # Default: R² for regression
        def scoring(model, X, y):
            y_pred = model.predict(X)
            ss_res = np.sum((y - y_pred) ** 2)
            ss_tot = np.sum((y - np.mean(y)) ** 2)
            return 1 - ss_res / ss_tot

    baseline_score = scoring(model, X, y)
    n, p = X.shape
    importances = np.zeros((n_repeats, p))

    for r in range(n_repeats):
        for j in range(p):
            X_perm = X.copy()
            # Randomly shuffle feature j
            X_perm[:, j] = rng.permutation(X_perm[:, j])
            perm_score = scoring(model, X_perm, y)
            importances[r, j] = baseline_score - perm_score

    return {
        'importances_mean': np.mean(importances, axis=0),
        'importances_std': np.std(importances, axis=0),
        'importances': importances
    }


# ---- Demo ----
from sklearn.linear_model import LinearRegression

np.random.seed(42)
X = np.random.randn(200, 5)
y = 3*X[:,0] - 2*X[:,2] + 0.5*X[:,4] + np.random.randn(200)*0.3

model = LinearRegression().fit(X, y)
result = permutation_importance(model, X, y)
for j in range(5):
    print(f"Feature {j}: importance = {result['importances_mean'][j]:.4f} ± {result['importances_std'][j]:.4f}")
        

import tensorflow as tf
import numpy as np

# ---- Generate synthetic data ----
np.random.seed(42)
n_samples, n_features = 500, 20
X_train = np.random.randn(n_samples, n_features).astype(np.float32)
true_weights = np.zeros(n_features)
true_weights[[0, 3, 7, 15]] = [3.0, -2.0, 1.5, -1.0]  # Only 4 active features
y_train = (X_train @ true_weights + np.random.randn(n_samples) * 0.5).astype(np.float32)


class RegularizedLinearModel(tf.keras.Model):
    """Linear model with configurable regularization in TF."""
    def __init__(self, n_features, l1_lambda=0.0, l2_lambda=0.0):
        super().__init__()
        self.w = tf.Variable(tf.zeros([n_features, 1]), name='weights')
        self.b = tf.Variable(tf.zeros([1]), name='bias')
        self.l1_lambda = l1_lambda
        self.l2_lambda = l2_lambda

    def call(self, X):
        return tf.matmul(X, self.w) + self.b

    def regularization_loss(self):
        l1_loss = self.l1_lambda * tf.reduce_sum(tf.abs(self.w))
        l2_loss = self.l2_lambda * tf.reduce_sum(tf.square(self.w))
        return l1_loss + l2_loss


def train_model(model, X, y, epochs=500, lr=0.01):
    optimizer = tf.keras.optimizers.Adam(learning_rate=lr)
    X_tensor = tf.constant(X)
    y_tensor = tf.constant(y.reshape(-1, 1))

    for epoch in range(epochs):
        with tf.GradientTape() as tape:
            predictions = model(X_tensor)
            mse_loss = tf.reduce_mean(tf.square(y_tensor - predictions))
            total_loss = mse_loss + model.regularization_loss()

        gradients = tape.gradient(total_loss, model.trainable_variables)
        optimizer.apply_gradients(zip(gradients, model.trainable_variables))

        if (epoch + 1) % 100 == 0:
            print(f"Epoch {epoch+1}: Loss={total_loss.numpy():.4f}")

    return model


# Train with different regularizations
print("="*50, "RIDGE (L2)")
ridge_model = RegularizedLinearModel(n_features, l2_lambda=0.1)
train_model(ridge_model, X_train, y_train)
ridge_weights = ridge_model.w.numpy().flatten()
print(f"Non-zero weights (|w| > 0.1): {np.sum(np.abs(ridge_weights) > 0.1)}")

print("="*50, "LASSO (L1)")
lasso_model = RegularizedLinearModel(n_features, l1_lambda=0.1)
train_model(lasso_model, X_train, y_train)
lasso_weights = lasso_model.w.numpy().flatten()
print(f"Non-zero weights (|w| > 0.1): {np.sum(np.abs(lasso_weights) > 0.1)}")

print("="*50, "ELASTIC NET")
enet_model = RegularizedLinearModel(n_features, l1_lambda=0.05, l2_lambda=0.05)
train_model(enet_model, X_train, y_train)
enet_weights = enet_model.w.numpy().flatten()
print(f"Non-zero weights (|w| > 0.1): {np.sum(np.abs(enet_weights) > 0.1)}")


# ---- Using Keras built-in regularizers ----
model_keras = tf.keras.Sequential([
    tf.keras.layers.Dense(64, activation='relu',
                          kernel_regularizer=tf.keras.regularizers.l1_l2(l1=0.01, l2=0.01),
                          input_shape=(n_features,)),
    tf.keras.layers.Dense(32, activation='relu',
                          kernel_regularizer=tf.keras.regularizers.l2(0.01)),
    tf.keras.layers.Dense(1)
])
model_keras.compile(optimizer='adam', loss='mse')
model_keras.fit(X_train, y_train, epochs=100, verbose=0)
print("Keras model score:", model_keras.evaluate(X_train, y_train, verbose=0))
        

import tensorflow as tf

# Feature columns with built-in preprocessing
class TFFeaturePipeline:
    """Feature engineering pipeline using TensorFlow layers."""

    def build_preprocessing_model(self):
        # Normalization layer (StandardScaler equivalent)
        normalizer = tf.keras.layers.Normalization(axis=-1)

        # StringLookup + OneHot (One-Hot encoding)
        category_encoder = tf.keras.layers.StringLookup(
            vocabulary=['Mumbai', 'Delhi', 'Bangalore', 'Chennai'],
            output_mode='one_hot'
        )

        # Discretization (binning)
        age_binner = tf.keras.layers.Discretization(
            bin_boundaries=[18, 25, 35, 45, 55, 65]
        )

        return normalizer, category_encoder, age_binner

    def demo(self):
        # Example: Normalize numerical features
        data = tf.constant([[1200.0, 2.0], [1800.0, 3.0],
                           [2400.0, 4.0], [900.0, 1.0]])
        normalizer = tf.keras.layers.Normalization(axis=-1)
        normalizer.adapt(data)
        print("Normalized:\n", normalizer(data).numpy())

pipeline = TFFeaturePipeline()
pipeline.demo()
        

import numpy as np
import pandas as pd
from sklearn.pipeline import Pipeline
from sklearn.compose import ColumnTransformer
from sklearn.preprocessing import (
    StandardScaler, MinMaxScaler, RobustScaler, MaxAbsScaler,
    OneHotEncoder, OrdinalEncoder, LabelEncoder,
    PowerTransformer, PolynomialFeatures
)
from sklearn.impute import SimpleImputer, KNNImputer
from sklearn.experimental import enable_iterative_imputer
from sklearn.impute import IterativeImputer  # MICE
from sklearn.linear_model import Ridge, Lasso, ElasticNet, LassoCV, RidgeCV
from sklearn.feature_selection import (
    SelectKBest, mutual_info_regression, RFE,
    SelectFromModel
)
from sklearn.model_selection import cross_val_score
from sklearn.datasets import make_regression

# ============================================================
# 1. DATA PREPARATION
# ============================================================

# Create a realistic mixed-type dataset
np.random.seed(42)
n = 500
df = pd.DataFrame({
    'age': np.random.randint(18, 70, n),
    'income': np.random.lognormal(10, 1, n),  # Skewed!
    'city': np.random.choice(['Mumbai', 'Delhi', 'Bangalore',
                             'Chennai', 'Hyderabad'], n),
    'education': np.random.choice(['High School', 'Bachelor',
                                   'Master', 'PhD'], n),
    'experience': np.random.randint(0, 40, n),
    'credit_score': np.random.randint(300, 850, n),
})

# Inject missing values (10%)
for col in ['income', 'credit_score', 'age']:
    mask = np.random.random(n) < 0.1
    df.loc[mask, col] = np.nan

# Create target
y = (0.5 * df['age'].fillna(35) +
     0.001 * df['income'].fillna(50000) +
     0.8 * df['experience'] +
     np.random.randn(n) * 5)

# ============================================================
# 2. COLUMN TRANSFORMER: Different processing for different types
# ============================================================

numeric_features = ['age', 'income', 'experience', 'credit_score']
categorical_features = ['city']
ordinal_features = ['education']

numeric_pipeline = Pipeline([
    ('imputer', KNNImputer(n_neighbors=5)),  # KNN imputation
    ('scaler', StandardScaler()),           # Z-score scaling
])

skewed_pipeline = Pipeline([
    ('imputer', SimpleImputer(strategy='median')),
    ('power', PowerTransformer(method='yeo-johnson')),  # Handle skew
])

categorical_pipeline = Pipeline([
    ('imputer', SimpleImputer(strategy='most_frequent')),
    ('encoder', OneHotEncoder(handle_unknown='ignore', sparse_output=False)),
])

ordinal_pipeline = Pipeline([
    ('imputer', SimpleImputer(strategy='most_frequent')),
    ('encoder', OrdinalEncoder(
        categories=[['High School', 'Bachelor', 'Master', 'PhD']])),
])

preprocessor = ColumnTransformer([
    ('num', numeric_pipeline, ['age', 'experience', 'credit_score']),
    ('skewed', skewed_pipeline, ['income']),
    ('cat', categorical_pipeline, categorical_features),
    ('ord', ordinal_pipeline, ordinal_features),
])

# ============================================================
# 3. FULL PIPELINE WITH REGULARIZED MODEL
# ============================================================

# Ridge Pipeline
ridge_pipeline = Pipeline([
    ('preprocessor', preprocessor),
    ('model', Ridge(alpha=1.0))
])

# Lasso Pipeline
lasso_pipeline = Pipeline([
    ('preprocessor', preprocessor),
    ('model', Lasso(alpha=0.1))
])

# ElasticNet Pipeline
enet_pipeline = Pipeline([
    ('preprocessor', preprocessor),
    ('model', ElasticNet(alpha=0.1, l1_ratio=0.5))
])

# Cross-validation comparison
for name, pipe in [('Ridge', ridge_pipeline),
                   ('Lasso', lasso_pipeline),
                   ('ElasticNet', enet_pipeline)]:
    scores = cross_val_score(pipe, df, y, cv=5, scoring='r2')
    print(f"{name:12s}: R² = {scores.mean():.4f} ± {scores.std():.4f}")

# ============================================================
# 4. REGULARIZATION PATH: Coefficients vs Lambda
# ============================================================

from sklearn.linear_model import lasso_path

# Prepare data
X_processed = preprocessor.fit_transform(df)
alphas, coefs, _ = lasso_path(X_processed, y, alphas=np.logspace(-3, 1, 50))

print("\nLasso Path Summary:")
print(f"  Alphas tested: {len(alphas)}")
print(f"  Features: {coefs.shape[0]}")
print(f"  At α=0.001: {np.sum(np.abs(coefs[:, 0]) > 1e-6)} non-zero features")
print(f"  At α=10.0:  {np.sum(np.abs(coefs[:, -1]) > 1e-6)} non-zero features")

# ============================================================
# 5. FEATURE SELECTION METHODS
# ============================================================

# Filter method: Mutual Information
selector_mi = SelectKBest(mutual_info_regression, k=5)
X_selected = selector_mi.fit_transform(X_processed, y)
print(f"\nMutual Info - selected {X_selected.shape[1]} features")
print(f"  MI scores: {np.round(selector_mi.scores_, 3)}")

# Wrapper method: Recursive Feature Elimination
rfe = RFE(Ridge(alpha=1.0), n_features_to_select=5)
rfe.fit(X_processed, y)
print(f"\nRFE - selected features: {rfe.support_}")
print(f"  Feature ranking: {rfe.ranking_}")

# Embedded method: Lasso-based selection
selector_lasso = SelectFromModel(Lasso(alpha=0.1))
selector_lasso.fit(X_processed, y)
print(f"\nLasso selection - selected {np.sum(selector_lasso.get_support())} features")

# ============================================================
# 6. MICE IMPUTATION
# ============================================================

mice_imputer = IterativeImputer(
    max_iter=10,
    random_state=42,
    sample_posterior=True  # For multiple imputation
)
X_numeric = df[numeric_features].values
X_mice = mice_imputer.fit_transform(X_numeric)
print(f"\nMICE imputation: {np.sum(np.isnan(X_numeric))} NaN → {np.sum(np.isnan(X_mice))} NaN")
        

import numpy as np
import matplotlib.pyplot as plt
from sklearn.preprocessing import (
    StandardScaler, MinMaxScaler, RobustScaler, MaxAbsScaler,
    PowerTransformer
)

# Generate skewed data
np.random.seed(42)
data = np.random.lognormal(0, 1, 1000).reshape(-1, 1)

# Add outliers
data = np.vstack([data, [[50], [100], [200]]])

scalers = {
    'Original': None,
    'StandardScaler': StandardScaler(),
    'MinMaxScaler': MinMaxScaler(),
    'RobustScaler': RobustScaler(),
    'MaxAbsScaler': MaxAbsScaler(),
    'Yeo-Johnson': PowerTransformer(method='yeo-johnson'),
}

fig, axes = plt.subplots(2, 3, figsize=(15, 10))
for ax, (name, scaler) in zip(axes.flat, scalers.items()):
    if scaler is None:
        transformed = data
    else:
        transformed = scaler.fit_transform(data)

    ax.hist(transformed, bins=50, color='#059669', alpha=0.7, edgecolor='white')
    ax.set_title(name, fontsize=14, fontweight='bold')
    ax.set_xlabel(f"Range: [{transformed.min():.1f}, {transformed.max():.1f}]")

plt.tight_layout()
plt.savefig('scaling_comparison.png', dpi=150)
plt.show()
print("Plot saved as scaling_comparison.png")