3.6 Lab: Linear Regression — 線性迴歸實驗室

本實驗室對應 ISLP §3.6 Lab: Linear Regression，將使用 Boston 房價資料集實作完整的線性迴歸分析流程。涵蓋資料探索、簡單線性迴歸、多元線性迴歸、交互作用項、殘差分析等核心概念。

🔧 套件匯入與環境設定

首先匯入本實驗室所需的 Python 套件。我們將使用 statsmodels 進行統計分析，ISLP 提供專門的建模工具。

基礎套件匯入

# 基礎科學計算套件
import numpy as np
import pandas as pd
from matplotlib.pyplot import subplots

進階套件匯入

# statsmodels 統計建模套件
import statsmodels.api as sm

# 從 statsmodels 匯入特定功能
from statsmodels.stats.outliers_influence \
     import variance_inflation_factor as VIF
from statsmodels.stats.anova import anova_lm

# ISLP 套件 - 本書專用建模工具
from ISLP import load_data
from ISLP.models import (ModelSpec as MS,
                         summarize,
                         poly)

命名空間檢視

Python 的 dir() 函數可以查看當前命名空間中的所有物件，幫助我們了解可用的函數和方法。

# 查看可用的物件
dir()

輸出顯示了當前環境中所有可用的函數和變數，包括我們匯入的套件函數如 MS、summarize、poly 等。

物件屬性檢視

每個 Python 物件都有自己的命名空間，包含其屬性和方法。例如，NumPy 陣列物件包含 sum 方法：

# 建立 NumPy 陣列並查看其方法
A = np.array([3,5,11])
dir(A)

從輸出可以看到 A.sum 存在，這是一個可以用來計算陣列總和的方法。

# 使用 sum 方法
A.sum()

📊 簡單線性迴歸

我們將使用 ModelSpec() 轉換函數來建構模型矩陣（也稱為設計矩陣），進行線性迴歸分析。

載入 Boston 資料集

Boston 資料集記錄了 506 個波士頓社區的房價資料，包含 medv（房屋中位數價值）和 13 個預測變數，如 rm（平均房間數）、age（1940年前建築比例）、lstat（低收入戶比例）等。

# 載入 Boston 資料集
Boston = load_data('Boston')
print(Boston.shape)  # 查看資料維度
print(Boston.columns)  # 查看欄位名稱

建立簡單線性迴歸模型

我們使用 lstat 變數來預測 medv。首先建立設計矩陣：

# 手動建立設計矩陣
X = pd.DataFrame({'intercept': np.ones(Boston.shape[0]),
                  'lstat': Boston['lstat']})
X[:4]  # 顯示前 4 筆資料

接著建立並擬合簡單線性迴歸模型：

# 建立線性迴歸模型
results = sm.OLS(Boston['medv'], X).fit()
# 摘要結果
summarize(results)

使用 ISLP 的 ModelSpec

ISLP 套件提供了更便捷的 ModelSpec 類別來處理設計矩陣的建構：

# 使用 ModelSpec 建立設計矩陣
design = MS(['lstat'])
design = design.fit(Boston)
X = design.transform(Boston)
X[:4]  # 顯示前 4 筆資料

也可以使用 fit_transform 一次性完成：

# 一次性完成 fit 和 transform
design = MS(['lstat'])
X = design.fit_transform(Boston)
X[:4]  # 顯示前 4 筆資料

模型結果解讀

# 詳細的模型摘要
results.summary()

從結果可以看到： - 截距項 (intercept) = 34.55，表示當 lstat = 0 時的預期 medv - lstat 的係數 = -0.95，表示 lstat 每增加 1 單位，medv 減少 0.95 單位 - R-squared = 0.544，表示模型解釋了 54.4% 的變異

係數提取

# 提取模型係數
results.params

新資料預測

使用擬合好的模型對新資料進行預測：

# 建立新資料
new_df = pd.DataFrame({'lstat':[5, 10, 15]})
# 轉換為設計矩陣格式
newX = design.transform(new_df)
# 進行預測
newX

📈 多元線性迴歸

現在使用多個預測變數來建立更複雜的迴歸模型。

多元迴歸模型擬合

使用 lstat 和 age 兩個變數：

# 建立多元線性迴歸模型
design = MS(['lstat', 'age'])
X = design.fit_transform(Boston)
model1 = sm.OLS(Boston['medv'], X).fit()
summarize(model1)

完整多元迴歸模型

使用所有 13 個預測變數：

# 使用所有變數建立模型
X = MS(Boston.columns.drop('medv')).fit_transform(Boston)
model_all = sm.OLS(Boston['medv'], X).fit()
summarize(model_all)

🔄 變數轉換與交互作用

多項式特徵

使用 poly 函數建立多項式特徵：

# 建立 lstat 的二次多項式特徵
X_poly = MS([poly('lstat', degree=2)]).fit_transform(Boston)
model_poly = sm.OLS(Boston['medv'], X_poly).fit()
summarize(model_poly)

比較模型效能

使用 ANOVA 比較不同模型：

# 比較線性模型與二次模型
anova_lm(results, model_poly)

交互作用項

加入交互作用項來分析變數之間的關聯：

# 加入交互作用項
design_interaction = MS(['lstat', 'age', ('lstat', 'age')])
X_interaction = design_interaction.fit_transform(Boston)
model_interaction = sm.OLS(Boston['medv'], X_interaction).fit()
summarize(model_interaction)

🔍 模型診斷

殘差分析

殘差是觀測值與預測值的差異，用於檢驗模型的合理性：

# 計算殘差
residuals = Boston['medv'] - results.predict(X)
# 殘差基本統計
print(residuals.describe())

殘差圖

繪製殘差圖檢視模型的殘差分佈模式：

import matplotlib.pyplot as plt

# 預測值 vs 殘差
plt.scatter(results.predict(X), residuals, alpha=0.6)
plt.axhline(y=0, color='r', linestyle='--')
plt.xlabel('Fitted values')
plt.ylabel('Residuals')
plt.title('Residual Plot')
plt.show()

變數膨脹因子 (VIF)

VIF 用於檢測多重共線性問題：

# 計算 VIF
vif_data = pd.DataFrame()
vif_data["Variable"] = X.columns
vif_data["VIF"] = [VIF(X.values, i) for i in range(X.shape[1])]
print(vif_data)

模型比較

比較不同模型的 AIC 和 BIC 值：

# 比較不同模型的準則
print(f"Simple model AIC: {results.aic:.2f}")
print(f"Polynomial model AIC: {model_poly.aic:.2f}")
print(f"Full model AIC: {model_all.aic:.2f}")

📋 模型預測與評估

交叉驗證

使用交叉驗證來評估模型的泛化能力：

from sklearn.model_selection import cross_val_score
from sklearn.linear_model import LinearRegression

# 使用 sklearn 進行交叉驗證
X_simple = X[['intercept', 'lstat']].values
y = Boston['medv'].values

model = LinearRegression()
scores = cross_val_score(model, X_simple, y, cv=5, scoring='r2')
print(f"Cross-validation R² scores: {scores}")
print(f"Mean R²: {scores.mean():.3f} (+/- {scores.std() * 2:.3f})")

預測新資料

使用擬合好的模型對新資料進行預測：

# 創建新資料點
new_data = pd.DataFrame({
    'lstat': [5, 10, 15, 20],
    'age': [65, 80, 90, 100]
})

# 使用完整模型進行預測
X_new = MS(Boston.columns.drop('medv')).fit_transform(new_data)
predictions = model_all.predict(X_new)
print("Predictions:", predictions)

🎯 重點總結

1. 簡單線性迴歸：使用單一預測變數建立模型，適合探索變數間的基本關係
2. 多元線性迴歸：使用多個預測變數，能更好地捕捉複雜關係
3. 多項式轉換：可以捕捉非線性關係，但需要注意過擬合問題
4. 交互作用項：用於分析變數之間的交互效應
5. 殘差分析：是檢驗模型合理性的重要工具，殘差應呈隨機分佈
6. VIF 檢驗：用於檢測多重共線性問題

← 3.5 線性迴歸與 K-最近鄰居比較 　｜　 4.1 分類 →