Multiple Regression

Model relationships between multiple input features and a target variable, incorporating interactions and feature engineering.

Concept Overview

Multiple regression extends linear regression to handle multiple predictor variables simultaneously. It allows us to model complex relationships where the target depends on several features and their interactions.

The multiple regression model with p features is:

ŷ = β₀ + β₁x₁ + β₂x₂ + … + βₚxₚ + ε

Where:

  • ŷ (y-hat) is the predicted target value
  • β₀ is the intercept (y-intercept when all features are zero)
  • β₁, β₂, …, βₚ are the coefficients for each feature
  • x₁, x₂, …, xₚ are the input features
  • ε (epsilon) is the error term

Key Advantages:

  • Can model relationships with multiple predictors
  • Allows for control of confounding variables
  • Enables testing of interaction effects
  • More powerful than simple regression when multiple factors matter

Basic Multiple Regression

Let's start with a basic example using multiple features to predict a target variable.

import numpy as np
import pandas as pd
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_squared_error, r2_score, mean_absolute_error
from sklearn.model_selection import train_test_split

# Generate sample data with multiple features
np.random.seed(42)
n = 200
sqft = np.random.normal(1500, 500, n)
bedrooms = np.random.randint(1, 6, n)
age = np.random.normal(10, 5, n)
bathrooms = bedrooms * 0.7 + np.random.normal(0, 0.3, n)

# Target with multiple features
price = (100 * sqft + 20000 * bedrooms + 15000 * bathrooms - 5000 * age + 
         np.random.normal(0, 50000, n))

df = pd.DataFrame({
    'sqft': sqft,
    'bedrooms': bedrooms,
    'bathrooms': bathrooms,
    'age': age,
    'price': price
})

# Prepare features and target
X = df[['sqft', 'bedrooms', 'bathrooms', 'age']]
y = df['price']

# Split data
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# Create and train model
model = LinearRegression()
model.fit(X_train, y_train)

# Make predictions
y_train_pred = model.predict(X_train)
y_test_pred = model.predict(X_test)

# Model parameters
print("Multiple Regression Model:")
print(f"Intercept (β₀): {model.intercept_:.2f}")
print("\nCoefficients:")
for i, feature in enumerate(X.columns):
    print(f"  {feature} (β{i+1}): {model.coef_[i]:.2f}")

print(f"\nModel equation:")
print(f"price = {model.intercept_:.2f}", end="")
for i, feature in enumerate(X.columns):
    sign = "+" if model.coef_[i] >= 0 else ""
    print(f" {sign}{model.coef_[i]:.2f} × {feature}", end="")

# Evaluate metrics
train_r2 = r2_score(y_train, y_train_pred)
test_r2 = r2_score(y_test, y_test_pred)
train_mse = mean_squared_error(y_train, y_train_pred)
test_mse = mean_squared_error(y_test, y_test_pred)
train_mae = mean_absolute_error(y_train, y_train_pred)
test_mae = mean_absolute_error(y_test, y_test_pred)

print(f"\n\nTraining Metrics:")
print(f"  R²: {train_r2:.4f}")
print(f"  MSE: {train_mse:.2f}")
print(f"  MAE: {train_mae:.2f}")

print(f"\nTest Metrics:")
print(f"  R²: {test_r2:.4f}")
print(f"  MSE: {test_mse:.2f}")
print(f"  MAE: {test_mae:.2f}")

# Feature importance (coefficient magnitude, standardized)
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.transform(X_test)

model_scaled = LinearRegression()
model_scaled.fit(X_train_scaled, y_train)

print(f"\n\nStandardized Coefficients (for comparison):")
for i, feature in enumerate(X.columns):
    print(f"  {feature}: {model_scaled.coef_[i]:.2f}")

Standardized coefficients allow comparison of feature importance when features have different scales. Larger absolute values indicate greater impact on the target variable.

Interaction Terms

Interaction terms allow features to influence each other's effect on the target. For example, the effect of square footage might depend on the number of bedrooms.

import numpy as np
import pandas as pd
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import PolynomialFeatures
from sklearn.metrics import r2_score
from sklearn.model_selection import train_test_split

# Generate data with interaction effects
np.random.seed(42)
n = 200
sqft = np.random.normal(1500, 500, n)
bedrooms = np.random.randint(1, 6, n)
# Interaction: price increases more with sqft when there are more bedrooms
price = (100 * sqft + 20000 * bedrooms + 5 * sqft * bedrooms + 
         np.random.normal(0, 50000, n))

df = pd.DataFrame({
    'sqft': sqft,
    'bedrooms': bedrooms,
    'price': price
})

X = df[['sqft', 'bedrooms']]
y = df['price']

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# Model without interactions
model_no_interaction = LinearRegression()
model_no_interaction.fit(X_train, y_train)
pred_no_int = model_no_interaction.predict(X_test)
r2_no_int = r2_score(y_test, pred_no_int)

print("Model without interactions:")
print(f"Test R²: {r2_no_int:.4f}")
print(f"Coefficients: sqft={model_no_interaction.coef_[0]:.2f}, bedrooms={model_no_interaction.coef_[1]:.2f}")

# Model with interactions (degree 2 polynomial includes interactions)
from sklearn.preprocessing import PolynomialFeatures
poly_features = PolynomialFeatures(degree=2, include_bias=False, interaction_only=True)
X_train_inter = poly_features.fit_transform(X_train)
X_test_inter = poly_features.transform(X_test)

feature_names = poly_features.get_feature_names_out(X.columns)
print(f"\nInteraction features created: {feature_names}")

model_interaction = LinearRegression()
model_interaction.fit(X_train_inter, y_train)
pred_inter = model_interaction.predict(X_test_inter)
r2_inter = r2_score(y_test, pred_inter)

print(f"\nModel with interactions:")
print(f"Test R²: {r2_inter:.4f}")
print("Coefficients:")
for i, name in enumerate(feature_names):
    print(f"  {name}: {model_interaction.coef_[i]:.2f}")

# Manual interaction term
X_train_manual = X_train.copy()
X_train_manual['sqft_x_bedrooms'] = X_train_manual['sqft'] * X_train_manual['bedrooms']
X_test_manual = X_test.copy()
X_test_manual['sqft_x_bedrooms'] = X_test_manual['sqft'] * X_test_manual['bedrooms']

model_manual = LinearRegression()
model_manual.fit(X_train_manual, y_train)
pred_manual = model_manual.predict(X_test_manual)
r2_manual = r2_score(y_test, pred_manual)

print(f"\nModel with manual interaction term:")
print(f"Test R²: {r2_manual:.4f}")
print("Coefficients:")
for i, feature in enumerate(X_train_manual.columns):
    print(f"  {feature}: {model_manual.coef_[i]:.2f}")

print(f"\n\nImprovement from interactions:")
print(f"Without interactions R²: {r2_no_int:.4f}")
print(f"With interactions R²: {r2_inter:.4f}")
print(f"Improvement: {(r2_inter - r2_no_int):.4f}")

Interaction terms can significantly improve model performance when features interact. Test for interactions by comparing models with and without interaction terms, or use domain knowledge to identify likely interactions.

Feature Selection

Not all features are equally important. Feature selection helps identify the most relevant features and can improve model performance and interpretability.

import numpy as np
import pandas as pd
from sklearn.linear_model import LinearRegression, Ridge, Lasso
from sklearn.feature_selection import SelectKBest, f_regression, RFE
from sklearn.metrics import r2_score, mean_squared_error
from sklearn.model_selection import train_test_split, cross_val_score
from sklearn.preprocessing import StandardScaler

# Generate data with some irrelevant features
np.random.seed(42)
n = 200
sqft = np.random.normal(1500, 500, n)
bedrooms = np.random.randint(1, 6, n)
age = np.random.normal(10, 5, n)
noise1 = np.random.normal(0, 1, n)  # Irrelevant
noise2 = np.random.normal(0, 1, n)  # Irrelevant
noise3 = np.random.normal(0, 1, n)  # Irrelevant

price = (100 * sqft + 20000 * bedrooms - 5000 * age + 
         np.random.normal(0, 50000, n))

df = pd.DataFrame({
    'sqft': sqft,
    'bedrooms': bedrooms,
    'age': age,
    'noise1': noise1,
    'noise2': noise2,
    'noise3': noise3,
    'price': price
})

X = df.drop('price', axis=1)
y = df['price']

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# 1. Univariate feature selection (F-statistic)
selector_f = SelectKBest(score_func=f_regression, k=3)
X_train_selected = selector_f.fit_transform(X_train, y_train)
X_test_selected = selector_f.transform(X_test)

selected_features = X.columns[selector_f.get_support()]
print(f"Selected features (F-test): {list(selected_features)}")
print(f"Feature scores: {dict(zip(X.columns, selector_f.scores_))}")

model_f = LinearRegression()
model_f.fit(X_train_selected, y_train)
pred_f = model_f.predict(X_test_selected)
r2_f = r2_score(y_test, pred_f)

print(f"\nModel with F-test selection - R²: {r2_f:.4f}")

# 2. Recursive Feature Elimination (RFE)
estimator = LinearRegression()
selector_rfe = RFE(estimator, n_features_to_select=3, step=1)
selector_rfe.fit(X_train, y_train)

selected_features_rfe = X.columns[selector_rfe.get_support()]
print(f"\nSelected features (RFE): {list(selected_features_rfe)}")
print(f"Feature rankings: {dict(zip(X.columns, selector_rfe.ranking_))}")

X_train_rfe = selector_rfe.transform(X_train)
X_test_rfe = selector_rfe.transform(X_test)

model_rfe = LinearRegression()
model_rfe.fit(X_train_rfe, y_train)
pred_rfe = model_rfe.predict(X_test_rfe)
r2_rfe = r2_score(y_test, pred_rfe)

print(f"Model with RFE - R²: {r2_rfe:.4f}")

# 3. Lasso for automatic feature selection
scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.transform(X_test)

lasso = Lasso(alpha=0.1, max_iter=2000)
lasso.fit(X_train_scaled, y_train)

selected_features_lasso = X.columns[np.abs(lasso.coef_) > 1e-5]
print(f"\nSelected features (Lasso): {list(selected_features_lasso)}")
print("Lasso coefficients:")
for i, feature in enumerate(X.columns):
    print(f"  {feature}: {lasso.coef_[i]:.4f}")

pred_lasso = lasso.predict(X_test_scaled)
r2_lasso = r2_score(y_test, pred_lasso)

print(f"Model with Lasso - R²: {r2_lasso:.4f}")

# Compare all models
print(f"\n\nModel Comparison:")
print(f"All features: R² = {r2_score(y_test, LinearRegression().fit(X_train, y_train).predict(X_test)):.4f}")
print(f"F-test selection: R² = {r2_f:.4f}")
print(f"RFE selection: R² = {r2_rfe:.4f}")
print(f"Lasso selection: R² = {r2_lasso:.4f}")

Feature selection helps reduce overfitting, improve interpretability, and sometimes improve performance. Different methods work better in different scenarios—experiment with multiple approaches.

Model Diagnostics

Diagnosing multiple regression models involves checking assumptions, detecting multicollinearity, and evaluating model performance through various metrics.

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression
from sklearn.metrics import r2_score, mean_squared_error
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler

# Generate data
np.random.seed(42)
n = 200
sqft = np.random.normal(1500, 500, n)
bedrooms = np.random.randint(1, 6, n)
age = np.random.normal(10, 5, n)
price = (100 * sqft + 20000 * bedrooms - 5000 * age + 
         np.random.normal(0, 50000, n))

df = pd.DataFrame({
    'sqft': sqft,
    'bedrooms': bedrooms,
    'age': age,
    'price': price
})

X = df[['sqft', 'bedrooms', 'age']]
y = df['price']

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

model = LinearRegression()
model.fit(X_train, y_train)

y_train_pred = model.predict(X_train)
y_test_pred = model.predict(X_test)

# Residuals
train_residuals = y_train - y_train_pred
test_residuals = y_test - y_test_pred

# Adjusted R²
n = len(y_train)
p = X_train.shape[1]
r2 = r2_score(y_train, y_train_pred)
adj_r2 = 1 - (1 - r2) * (n - 1) / (n - p - 1)

print("Model Diagnostics:")
print(f"R²: {r2:.4f}")
print(f"Adjusted R²: {adj_r2:.4f}")
print(f"Number of features (p): {p}")
print(f"Number of samples (n): {n}")

# Correlation matrix (check for multicollinearity)
correlation_matrix = X_train.corr()
print(f"\n\nCorrelation Matrix:")
print(correlation_matrix)

# VIF (Variance Inflation Factor) - requires statsmodels
try:
    import statsmodels.api as sm
    from statsmodels.stats.outliers_influence import variance_inflation_factor
    
    X_with_const = sm.add_constant(X_train)
    vif_data = pd.DataFrame()
    vif_data["Feature"] = X_train.columns
    vif_data["VIF"] = [variance_inflation_factor(X_with_const.values, i+1) 
                      for i in range(len(X_train.columns))]
    
    print(f"\n\nVariance Inflation Factor (VIF):")
    print(vif_data)
    print("\nInterpretation:")
    print("  VIF < 5: Low multicollinearity (acceptable)")
    print("  5 ≤ VIF < 10: Moderate multicollinearity (caution)")
    print("  VIF ≥ 10: High multicollinearity (problematic)")
except ImportError:
    print("\nNote: Install statsmodels for VIF calculation")

# Residual plots
fig, axes = plt.subplots(2, 2, figsize=(14, 12))

# Residuals vs predicted
axes[0, 0].scatter(y_train_pred, train_residuals, alpha=0.6, s=50)
axes[0, 0].axhline(y=0, color='red', linestyle='--', linewidth=2)
axes[0, 0].set_xlabel('Predicted Values', fontsize=12)
axes[0, 0].set_ylabel('Residuals', fontsize=12)
axes[0, 0].set_title('Residuals vs Predicted', fontsize=14)
axes[0, 0].grid(True, alpha=0.3)

# Q-Q plot
from scipy import stats
stats.probplot(train_residuals, dist="norm", plot=axes[0, 1])
axes[0, 1].set_title('Q-Q Plot (Normality Check)', fontsize=14)
axes[0, 1].grid(True, alpha=0.3)

# Residuals histogram
axes[1, 0].hist(train_residuals, bins=20, edgecolor='black', alpha=0.7)
axes[1, 0].axvline(x=0, color='red', linestyle='--', linewidth=2)
axes[1, 0].set_xlabel('Residuals', fontsize=12)
axes[1, 0].set_ylabel('Frequency', fontsize=12)
axes[1, 0].set_title('Residuals Distribution', fontsize=14)
axes[1, 0].grid(True, alpha=0.3)

# Actual vs Predicted
axes[1, 1].scatter(y_train, y_train_pred, alpha=0.6, s=50)
min_val = min(y_train.min(), y_train_pred.min())
max_val = max(y_train.max(), y_train_pred.max())
axes[1, 1].plot([min_val, max_val], [min_val, max_val], 'r--', linewidth=2)
axes[1, 1].set_xlabel('Actual Values', fontsize=12)
axes[1, 1].set_ylabel('Predicted Values', fontsize=12)
axes[1, 1].set_title('Actual vs Predicted', fontsize=14)
axes[1, 1].grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

Always perform diagnostic checks on multiple regression models. Check residuals for patterns, examine multicollinearity with VIF, and use adjusted R² when comparing models with different numbers of features.

Best Practices

  • Start with exploratory analysis: Understand relationships between features before modeling
  • Check for multicollinearity: Use correlation matrices and VIF scores to detect highly correlated features
  • Standardize when needed: Standardize features when comparing coefficient magnitudes or using regularization
  • Consider interactions: Test for interaction effects when domain knowledge suggests they exist
  • Feature selection: Remove irrelevant features to improve model performance and interpretability
  • Use adjusted R²: Compare models with different numbers of features using adjusted R²
  • Validate assumptions: Check linearity, homoscedasticity, normality, and independence
  • Cross-validation: Use cross-validation to assess model performance robustly
  • Regularization: Use Ridge or Lasso when features are numerous or multicollinear
  • Document choices: Record feature selection and model building decisions, share insights through info.studygrid@gmail.com

Next Steps

Learn about feature scaling to understand how to normalize features for better model performance and coefficient interpretability.

Try It Yourself

Exercise 1: Fit a model predicting CO2 from car weight and volume, then predict for weight 2300, volume 1300.

Show solution
from sklearn.linear_model import LinearRegression
X = [[790, 1000], [1160, 1200], [929, 1000], [865, 900], [1140, 1500]]
y = [99, 95, 95, 90, 105]
model = LinearRegression().fit(X, y)
print(model.predict([[2300, 1300]]).round(1))

Exercise 2: Why should features be on similar scales before interpreting coefficients?

Show solution

Raw coefficients depend on each feature's units, so they are not directly comparable. Scaling (standardizing) lets you compare each feature's relative importance.

Key Takeaways

  • Multiple regression uses several features to predict one target.
  • Each coefficient is the effect of one feature holding others constant.
  • Watch for multicollinearity between correlated features.

📘 Real-World Deep Dive

Knowing <strong>ML Multiple Regression (scikit-learn)</strong> well is what turns scikit-learn from a curiosity into a daily tool — you'll reach for it in nearly every real project.

Real-Life Scenario

An end-to-end usage of ML Multiple Regression that you'd actually see in a data pipeline or analytics notebook.

Real-Life Example

import numpy as np
from sklearn.linear_model import LinearRegression
X = np.array([[1, 2], [2, 1], [3, 4], [4, 3], [5, 5], [6, 7]])
y = np.array([5, 4, 11, 10, 15, 19])
m = LinearRegression().fit(X, y)
print("coef:", m.coef_, "intercept:", round(m.intercept_, 2))

Expected Output

(see source)

Common mistakes

  • fit expects numeric arrays; OneHotEncoder / LabelEncoder are easy to forget for categorical features.
  • Calling predict on a model trained on unscaled data and then scaling inputs at inference time silently degrades accuracy.
  • train_test_split(X, y) requires both arrays; train_test_split(X) for unsupervised learning slips past static checkers.
  • Treating ML Multiple Regression as a black box without reading the docs — the API has subtle defaults that bite when you scale.

🚀 Performance & Best Practices

  • Wrap preprocessing + estimator in a Pipeline so fit / predict stay reproducible.
  • Use joblib / pickle for serialising models, not the entire Python state.
  • For > 100 k rows, switch to HistGradientBoostingClassifier or cuML.
  • When working with scikit-learn, prefer vectorised / batched operations over Python loops.

🧪 Try It Yourself

  1. Reproduce the snippet on a representative slice of your own data.
  2. Profile the snippet with cProfile or timeit and find the single biggest improvement.
  3. Generalise the snippet into a small, reusable function you can drop into future projects.

FAQ: Multiple Regression

Common questions about this page.

What is Multiple Regression?

Multiple Regression is a Machine Learning lesson that explains multiple regression in Python. Model relationships between multiple input features and a target variable, incorporating interactions and feature engineering. It is written for beginners who want a clear definition and working examples.

Should I run multiple regression examples locally for better learning?

Yes. Use the browser editor on StudyGrid for a quick check, then Download the example and run it on your computer. Local runs show real errors and the real toolchain, which is one of the fastest ways to learn multiple regression in this Machine Learning Python lesson (Multiple Regression).

How do I use multiple regression in Python?

To use multiple regression in Python, follow the examples on this StudyGrid page. Copy a snippet, run it in the browser, then Download and run it locally for better learning. Change the values and compare the output.

What is the syntax of multiple regression?

This Multiple Regression tutorial shows multiple regression syntax with short Python examples. Use the code blocks in this lesson for the exact statements, then try them in your editor.

Multiple Regression example for beginners

Yes. This page includes a beginner multiple regression example you can copy and run. It is designed for searches such as "multiple regression for beginners", "multiple regression example", and "how to use multiple regression".

What are common mistakes with multiple regression?

Common multiple regression mistakes include wrong syntax, mixing types, and skipping practice. Work through this Machine Learning chapter in order, run every example, and check the output before moving on.

Why should I learn multiple regression?

Multiple Regression is used in real Python work. Learning multiple regression helps you write clearer programs and continue the Machine Learning tutorial on StudyGrid.

Is Multiple Regression free to learn online?

Yes. You can learn multiple regression free on StudyGrid (studygrid.in). This chapter is part of the Machine Learning path and includes examples, syntax, and next-step links.