📈 Linear Regression Foundations

🎯 Core Idea & History

Linear Regression (origins: Legendre, Gauss 1800s) models the relationship between a dependent variable y and one or more independent variables X using a linear equation. The goal is to find the best-fitting line (or hyperplane) that minimizes the sum of squared errors (Ordinary Least Squares). It is the cornerstone of supervised learning.

Deep Explanation

Model: y = β₀ + β₁x₁ + ... + βₚxₚ + ε where β are coefficients, ε is irreducible error. The solution: β = (XᵀX)⁻¹ Xᵀy (closed-form). Interpretability is key: coefficients represent the change in y for a 1-unit change in x, holding others constant.

Why it matters: Simple, fast, and interpretable. It serves as a baseline for all other models. Understanding its assumptions and limitations is crucial for building robust models.

📌 OLS📌 Interpretability📌 Closed-form solution📌 Linear relationship

import numpy as np
from sklearn.linear_model import LinearRegression

# Create data
X = np.array([[1, 1], [1, 2], [2, 2], [2, 3]])
y = np.dot(X, np.array([1, 2])) + 3

model = LinearRegression().fit(X, y)
print(f"Coefficients: {model.coef_}, Intercept: {model.intercept_}")
print(f"R²: {model.score(X, y)}")

🧮 Mathematical Core

🎯 Matrix Formulation & Normal Equation

Linear regression can be expressed compactly with matrices: y = Xβ + ε. The normal equation: β = (XᵀX)⁻¹ Xᵀy. This minimizes the cost function J(β) = ½Σ(yᵢ - Xᵢβ)².

Cost Function & Optimization

The cost function is convex (bowl-shaped), guaranteeing a global minimum. The normal equation provides an analytical solution, but can be slow for large feature sets (O(n³)). Gradient descent is an iterative alternative.

📐 Normal Equation📉 Convex Cost⚡ Analytical solution

# Manual normal equation
X = np.c_[np.ones((100,1)), np.random.rand(100,2)]
y = X @ np.array([3, 1.5, -2]) + np.random.randn(100)*0.1
beta_hat = np.linalg.inv(X.T @ X) @ X.T @ y
print(f"Estimated: {beta_hat}")

⚖️ Assumptions & Diagnostics

🎯 The Five Key Assumptions

1. Linearity 2. Independence of errors 3. Homoscedasticity (constant variance) 4. Normality of errors (for inference) 5. No perfect multicollinearity.

Violation Detection

Residual plots, Q-Q plots, and statistical tests (e.g., Breusch-Pagan for heteroscedasticity, Durbin-Watson for autocorrelation) help diagnose issues. Transformations (log, Box-Cox) or robust standard errors can mitigate violations.

🔍 Residual Analysis📊 Q-Q Plot📏 Variance Inflation Factor

import statsmodels.api as sm
model = sm.OLS(y, X).fit()
print(model.summary())  # includes p-values, R², F-statistic, Durbin-Watson

📊 Simple Linear Regression

🎯 One Predictor

y = β₀ + β₁x + ε. Interpret slope as the change in y per unit change in x. R² measures the proportion of variance explained.

Interpretation & Visualization

The slope and intercept are estimated to minimize the RSS (Residual Sum of Squares). The correlation coefficient r is tied to R² (R² = r² for simple regression).

📈 Slope Interpretation📊 R²🔗 Correlation

import matplotlib.pyplot as plt
x = np.array([1,2,3,4,5]); y = np.array([2,4,5,4,6])
model = LinearRegression().fit(x.reshape(-1,1), y)
plt.scatter(x,y); plt.plot(x, model.predict(x.reshape(-1,1)), 'r-')
plt.show()

📉 Multiple Linear Regression

🎯 Two or More Predictors

Extends simple regression to include multiple features. The coefficients represent partial effects — the change in y for a 1-unit change in xⱼ, holding other variables constant.

Adjusted R² & F-test

Adjusted R² penalizes model complexity. The F-test checks whether at least one predictor is significant. Feature selection (forward/backward) helps avoid overfitting.

📊 Adjusted R²🔬 F-test🧩 Feature Selection

from sklearn.datasets import load_diabetes
diabetes = load_diabetes()
X = diabetes.data[:, :2]  # first two features
model = LinearRegression().fit(X, diabetes.target)
print(f"Coeff: {model.coef_}, Intercept: {model.intercept_}")

🌀 Polynomial & Interaction

🎯 Capturing Non-linearity

Add polynomial terms (x², x³) or interaction terms (x₁*x₂) to model curvature and synergy between predictors.

Feature Engineering

Polynomial features from sklearn.preprocessing.PolynomialFeatures can transform linear regression into a powerful non-linear tool while keeping the model linear in parameters.

🧮 PolynomialFeatures🤝 Interactions📈 Non-linear Fit

from sklearn.preprocessing import PolynomialFeatures
poly = PolynomialFeatures(degree=2, include_bias=False)
X_poly = poly.fit_transform(X)
model = LinearRegression().fit(X_poly, y)

⛑️ Regularization (Ridge, Lasso, ElasticNet)

🎯 Shrinkage Methods

Ridge (L2 penalty), Lasso (L1 penalty), ElasticNet (L1+L2). They reduce overfitting and handle multicollinearity.

L1 vs L2

Lasso performs feature selection (sparse coefficients), while Ridge shrinks coefficients but keeps all features. ElasticNet combines both — useful when many correlated features exist.

🧹 Lasso Selection📉 Ridge Shrinkage⚖️ ElasticNet

from sklearn.linear_model import Ridge, Lasso, ElasticNet
ridge = Ridge(alpha=1.0).fit(X, y)
lasso = Lasso(alpha=0.1).fit(X, y)
elastic = ElasticNet(alpha=0.1, l1_ratio=0.5).fit(X, y)

📉 Gradient Descent

🎯 Iterative Optimization

Gradient descent minimizes the cost function by iteratively moving in the direction of steepest descent. Learning rate, batch size, and convergence criteria are key.

Stochastic & Mini-batch

Batch GD uses the full dataset; Stochastic GD uses one sample at a time; Mini-batch is the sweet spot. Adaptive methods (Adam, RMSProp) are used in deep learning.

⚡ Learning Rate🔄 Batch Size📈 Convergence

def gradient_descent(X, y, lr=0.01, epochs=1000):
    m = len(y); theta = np.zeros(X.shape[1])
    for _ in range(epochs):
        grad = (1/m) * X.T @ (X @ theta - y)
        theta -= lr * grad
    return theta

📐 Statistical Inference

🎯 Confidence Intervals & p-values

Under the normality assumption, we can construct confidence intervals for coefficients and compute p-values to test the null hypothesis βⱼ = 0.

Hypothesis Testing

The t-statistic for each coefficient is (β̂ⱼ / SE(β̂ⱼ)). The F-test for overall significance. These are essential for scientific research and causal inference.

📏 t-statistic🔬 p-value📊 Confidence Interval

import statsmodels.api as sm
model = sm.OLS(y, X).fit()
print(model.summary())  # includes p-values and CIs

🔧 Feature Engineering

🎯 Transformations & Encoding

Log, square root, Box-Cox for skewed targets. Categorical encoding (one-hot, dummy). Scaling (StandardScaler, MinMaxScaler) for gradient-based models.

Creating New Features

Domain knowledge can be used to create interaction terms, polynomial features, or aggregated features (e.g., mean, max). Feature engineering often has the biggest impact on performance.

📊 Scaling🧪 Encoding🧠 Domain Features

from sklearn.preprocessing import StandardScaler, OneHotEncoder
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
encoder = OneHotEncoder()
X_encoded = encoder.fit_transform(categorical_data).toarray()

✅ Cross-Validation & Metrics

🎯 Evaluation & Tuning

R², RMSE, MAE, MAPE. Cross-validation (k-fold) for model selection and hyperparameter tuning.

Robust Evaluation

R² measures variance explained. RMSE is in the same units as y. MAE is less sensitive to outliers. Use cross-validation to avoid overfitting and assess generalization.

📉 RMSE📊 R²🔄 k-fold CV

from sklearn.model_selection import cross_val_score
scores = cross_val_score(LinearRegression(), X, y, cv=5, scoring='neg_mean_squared_error')
print(f"RMSE: {np.sqrt(-scores.mean())}")

🔗 Multicollinearity & VIF

🎯 Detecting & Handling

Multicollinearity inflates variance of coefficients. Variance Inflation Factor (VIF) > 5 or 10 indicates high multicollinearity.

Remedies

Drop highly correlated features, combine them (PCA), or use regularization (Ridge/Lasso) which shrink coefficients and stabilize variance.

📏 VIF🔀 Correlation Matrix🧹 Feature Removal

from statsmodels.stats.outliers_influence import variance_inflation_factor
vif = [variance_inflation_factor(X, i) for i in range(X.shape[1])]
print(vif)

🔮 Logistic Regression (classification)

🎯 Binary & Multinomial

Logistic regression models the probability of a class using the logistic function: p = 1/(1 + e^{-z}). Coefficients are log-odds.

Odds & Maximum Likelihood

Instead of OLS, logistic regression uses Maximum Likelihood Estimation (MLE). It can handle binary, multiclass (softmax), and ordinal outcomes.

📊 Log-odds🎯 MLE🔮 Probability

from sklearn.linear_model import LogisticRegression
clf = LogisticRegression().fit(X, y_binary)
print(f"Accuracy: {clf.score(X, y_binary)}")

🚀 Deployment & Production

🎯 From Model to API

Save model (pickle, joblib), create REST API (Flask, FastAPI), monitor drift, and retrain pipelines.

MLOps Considerations

Feature store, versioning, A/B testing, and model monitoring (data drift, concept drift). Use MLflow for tracking experiments.

📦 Pickle🌐 FastAPI📊 Drift Monitoring

import joblib
joblib.dump(model, 'model.pkl')
# In API: model = joblib.load('model.pkl'); prediction = model.predict(input)