Regression

What is Regression?

Regression is the supervised learning task of predicting a continuous output (scalar or vector) from one or more inputs. The model learns a function f(x) ≈ y where y ∈ ℝ.

Common applications: house price prediction, stock forecasting, demand estimation, age estimation from images, temperature prediction.

Linear Regression

The simplest regression model assumes a linear relationship between features and target:

ŷ = θ₀ + θ₁x₁ + θ₂x₂ + … + θₙxₙ = θᵀx (with x₀=1 for bias)

Loss function — Mean Squared Error (MSE):

L(θ) = (1/n) Σᵢ (yᵢ − θᵀxᵢ)²

Minimizing MSE is equivalent to maximum likelihood under Gaussian noise.

Closed-Form Solution (Normal Equations)

θ* = (XᵀX)⁻¹ Xᵀy

Works for small datasets (O(p³) complexity for p features). Computationally impractical for large p.

Gradient Descent Solution

∂L/∂θ = (2/n) Xᵀ(Xθ − y)

θ ← θ − α · ∂L/∂θ

Scales to millions of examples and high-dimensional data.

from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_squared_error, r2_score
import numpy as np

model = LinearRegression()
model.fit(X_train, y_train)
y_pred = model.predict(X_test)

print(f"RMSE: {np.sqrt(mean_squared_error(y_test, y_pred)):.3f}")
print(f"R²:   {r2_score(y_test, y_pred):.3f}")
print(f"Coefficients: {model.coef_}")
print(f"Intercept: {model.intercept_}")

Assumptions of Linear Regression

1. Linearity: relationship between X and y is linear
2. Independence: residuals are independent
3. Homoscedasticity: constant variance of residuals across fitted values
4. Normality: residuals are normally distributed (for inference)
5. No multicollinearity: features are not strongly correlated

Violating these assumptions degrades coefficient estimates and confidence intervals (but not necessarily predictions).

Regularized Regression

Regularization adds a penalty on coefficient magnitude to combat overfitting and handle multicollinearity.

Ridge Regression (L2)

L_ridge(θ) = MSE + λ Σⱼ θⱼ²

Closed form: θ* = (XᵀX + λI)⁻¹ Xᵀy

Effect: shrinks all coefficients toward zero, but never exactly zero. λ controls regularization strength (find via cross-validation).

Lasso Regression (L1)

L_lasso(θ) = MSE + λ Σⱼ |θⱼ|

No closed form — requires coordinate descent. Effect: drives coefficients of irrelevant features to exactly zero (sparse solution). Acts as automatic feature selection.

Elastic Net

L_EN(θ) = MSE + λ₁ Σ|θⱼ| + λ₂ Σθⱼ²

Combination of L1 and L2. Use when features are correlated (Lasso picks one arbitrarily, Elastic Net can include groups).

from sklearn.linear_model import Ridge, Lasso, ElasticNet
from sklearn.model_selection import GridSearchCV

# Ridge
ridge = Ridge(alpha=1.0)    # alpha = λ
ridge.fit(X_train, y_train)

# Lasso with cross-validation
lasso_cv = LassoCV(cv=5, alphas=np.logspace(-4, 2, 50))
lasso_cv.fit(X_train, y_train)
print(f"Best alpha: {lasso_cv.alpha_:.4f}")
print(f"Non-zero features: {np.sum(lasso_cv.coef_ != 0)}")

Comparing Regularization

Polynomial Regression

Fit nonlinear relationships using polynomial features with a linear model:

ŷ = θ₀ + θ₁x + θ₂x² + θ₃x³ + …

The model is still linear in parameters (θ), just uses transformed features. Add polynomial terms, then apply linear regression.

from sklearn.preprocessing import PolynomialFeatures
from sklearn.pipeline import make_pipeline

degree = 3
model = make_pipeline(PolynomialFeatures(degree), Ridge(alpha=0.1))
model.fit(X_train, y_train)

Warning: high-degree polynomials overfit dramatically. Use regularization and cross-validate degree.

Quantile Regression

Instead of predicting the mean, predict a specific quantile (e.g., median = 50th percentile):

Loss for quantile τ: L_τ(r) = r·τ if r ≥ 0 else r·(τ−1)

This is the pinball loss / quantile loss. Useful for: - Predicting confidence intervals (predict 10th and 90th quantile) - Median regression (more robust to outliers than MSE) - Asymmetric cost functions

Decision Tree Regression

Non-parametric: recursively splits the feature space into rectangular regions and predicts the mean of training samples in each leaf.

from sklearn.tree import DecisionTreeRegressor

tree = DecisionTreeRegressor(max_depth=5, min_samples_leaf=10)
tree.fit(X_train, y_train)

Ensemble Regression Methods

Random Forest

Average predictions of many decorrelated decision trees. Reduces variance substantially.

Gradient Boosting (XGBoost, LightGBM, CatBoost)

Build trees sequentially, each correcting the residuals of the previous ensemble:

F_m(x) = F_{m-1}(x) + α · hₘ(x)

where hₘ(x) is a new tree fit to the residuals. Typically the best off-the-shelf method for tabular regression.

import xgboost as xgb

model = xgb.XGBRegressor(
    n_estimators=500,
    learning_rate=0.05,
    max_depth=6,
    subsample=0.8,
    colsample_bytree=0.8,
    reg_alpha=0.1,
    reg_lambda=1.0,
    random_state=42
)
model.fit(X_train, y_train,
          eval_set=[(X_test, y_test)],
          early_stopping_rounds=20,
          verbose=False)

Support Vector Regression (SVR)

Fits a tube of width ε around the prediction. Only points OUTSIDE the tube incur loss (ε-insensitive loss):

L(y, ŷ) = max(0, |y − ŷ| − ε)

Uses kernel trick for nonlinear regression (RBF kernel is most common).

from sklearn.svm import SVR
from sklearn.preprocessing import StandardScaler

# SVR requires scaling!
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X_train)

svr = SVR(kernel="rbf", C=100, gamma=0.1, epsilon=0.1)
svr.fit(X_scaled, y_train)

Regression Metrics

Choosing a metric: - Use RMSE when large errors are especially costly - Use MAE when all errors are equally important (or for robustness) - Use MAPE for relative comparison across scales - Use R² for explanatory modeling

from sklearn.metrics import (
    mean_absolute_error, mean_squared_error,
    mean_absolute_percentage_error, r2_score
)
import numpy as np

mae  = mean_absolute_error(y_test, y_pred)
rmse = np.sqrt(mean_squared_error(y_test, y_pred))
mape = mean_absolute_percentage_error(y_test, y_pred)
r2   = r2_score(y_test, y_pred)

Residual Analysis

After fitting, always inspect residuals r = y − ŷ:

• Plot residuals vs. fitted values: should be random scatter (no pattern = good). Funnel shape → heteroscedasticity.
• Q-Q plot: checks normality of residuals
• Residuals vs. each feature: reveals nonlinearity or missing interaction terms

import matplotlib.pyplot as plt

residuals = y_test - y_pred
plt.scatter(y_pred, residuals, alpha=0.3)
plt.axhline(0, color="red")
plt.xlabel("Fitted values")
plt.ylabel("Residuals")
plt.title("Residuals vs. Fitted")
plt.show()

Choosing a Regression Algorithm