第 2 章：統計學習 — 2.1 什麼是統計學習？

課本：James, Witten, Hastie, Tibshirani (2023), An Introduction to Statistical Learning with Applications in Python, Springer.

一、動機：Advertising 資料集

想像我們是一間被聘請的統計顧問公司——客戶想知道廣告支出與產品銷售量之間的關聯。Advertising 資料集包含 200 個市場中，三種廣告媒體的預算及對應銷售量。

# === Google Drive + Colab 相容資料讀取 ===
import pandas as pd
import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

# 方式 1：Google Colab 路徑
try:
    from google.colab import drive
    drive.mount('/content/drive')
    DATA_PATH = '/content/drive/MyDrive/ISLP_data/'
except ImportError:
    DATA_PATH = '/tmp/'  # 本機 fallback

ad = pd.read_csv(f'{DATA_PATH}Advertising.csv', index_col=0)
print(f"維度: {ad.shape}")
print(ad.head())
print(ad.describe())

# Figure 2.1：三種媒體 vs 銷售量
from sklearn.linear_model import LinearRegression

fig, axes = plt.subplots(1, 3, figsize=(15, 5))
for ax, var in zip(axes, ['TV', 'radio', 'newspaper']):
    ax.scatter(ad[var], ad['sales'], alpha=0.5, edgecolor='k', facecolor='red')
    X = ad[[var]].values; y = ad['sales'].values
    model = LinearRegression().fit(X, y)
    x_line = np.linspace(X.min(), X.max(), 100).reshape(-1, 1)
    ax.plot(x_line, model.predict(x_line), 'b-', linewidth=2)
    ax.set_xlabel(var); ax.set_ylabel('Sales')
    print(f"Sales ~ {var}: 截距={model.intercept_:.2f}, 斜率={model.coef_[0]:.4f}")
fig.suptitle('Figure 2.1: Advertising Data', fontsize=14)
plt.tight_layout(); plt.show()

觀察：TV 和 radio 與 sales 有明顯正相關，newspaper 相關性較弱。

二、統計學習的理論框架：Y = f(X) + ε

2.1 符號定義與機率基礎

給定聯合分布 \(P(X, Y)\)，我們的目標是找到一個函數 \(f\) 使得：

關鍵洞察：\(f(X) = \mathbb{E}[Y \mid X]\) 是最小化期望平方誤差的最優預測子：

這意味著即使不知道 f 的真實形式，條件期望本身定義了數學上最優的預測。

# Figure 2.2 概念：Income vs Education
income1 = pd.read_csv(f'{DATA_PATH}Income1.csv', index_col=0)
from sklearn.preprocessing import PolynomialFeatures

X = income1[['Education']].values; y = income1['Income'].values
poly = PolynomialFeatures(degree=3)
model = LinearRegression().fit(poly.fit_transform(X), y)

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))
ax1.scatter(income1['Education'], income1['Income'], c='red', s=50)
ax1.set_title('Left: Observed Data')
ax2.scatter(income1['Education'], income1['Income'], c='red', s=50, zorder=5)
x_curve = np.linspace(X.min(), X.max(), 200).reshape(-1, 1)
ax2.plot(x_curve, model.predict(poly.transform(x_curve)), 'b-', lw=2, label='f(Education)')
y_pred = model.predict(poly.transform(X))
for i in range(len(X)):
    ax2.plot([X[i,0], X[i,0]], [y_pred[i], y[i]], 'k-', lw=0.8, alpha=0.6)
ax2.set_title('Right: f(X) + Error Terms (ε)')
errors = y - y_pred
print(f"誤差均值 ≈ {errors.mean():.6f}（應趨近 0）")
plt.show()

藍色曲線是真正的 f（未知！），黑色垂直線是 ε。誤差均值約為零。

三、為什麼要估計 f？預測 vs 推論

3.1 預測（Prediction）

當我們只在意預測準確度而不需要解釋 \(\hat{f}\) 時，\(\hat{f}\) 可視為黑盒子。

預測誤差的數學分解（推導）

給定 \(X = x_0\)，期望預測誤差：

# 可約 vs 不可約誤差演示
np.random.seed(42)
def true_f(X): return 2 + 3*X + 0.5*X**2
n = 100
X_sim = np.random.uniform(0, 10, n)
epsilon = np.random.normal(0, 3, n)   # σ² = 9 (不可約誤差)
Y_sim = true_f(X_sim) + epsilon

model_A = LinearRegression().fit(X_sim.reshape(-1,1), Y_sim)
Y_pred_A = model_A.predict(X_sim.reshape(-1,1))

poly2 = PolynomialFeatures(degree=2)
model_B = LinearRegression().fit(poly2.fit_transform(X_sim.reshape(-1,1)), Y_sim)
Y_pred_B = model_B.predict(poly2.transform(X_sim.reshape(-1,1)))

reducible_A = np.mean((true_f(X_sim) - Y_pred_A)**2)
reducible_B = np.mean((true_f(X_sim) - Y_pred_B)**2)
irreducible = np.var(epsilon)

print(f"Var(ε) — 不可約誤差:      {irreducible:.3f}")
print(f"Model A (線性) 可約誤差:  {reducible_A:.3f}  總 MSE ≈ {reducible_A + irreducible:.3f}")
print(f"Model B (二次) 可約誤差:  {reducible_B:.3f}  總 MSE ≈ {reducible_B + irreducible:.3f}")
print("→ 不可約誤差是預測準確度的「上界」——無論模型多好都無法突破！")

3.2 推論（Inference）

當目標是理解 X 與 Y 之間的關係時：

• 🔍 哪些預測子與反應變數有關聯？
• 📈 關聯方向與強度？（正向/負向）
• 📐 關係可否用簡單線性方程式概括？

3.3 預測 vs 推論：選擇框架

3.4 類似理論比較：預測建模 vs 解釋建模

四、如何估計 f？參數化 vs 非參數化方法

4.1 參數化方法（Parametric）

兩步驟流程：

1. 假設 f 的函數形式：\(f(X) = \beta_0 + \beta_1 X_1 + \cdots + \beta_p X_p\)
2. 用訓練資料估計參數：最小平方法（OLS）

# 參數化：線性迴歸（Income2 資料集）
income2 = pd.read_csv(f'{DATA_PATH}Income2.csv', index_col=0)
X2 = income2[['Education', 'Seniority']].values
y2 = income2['Income'].values

model_linear = LinearRegression().fit(X2, y2)
print("=== 參數化：線性迴歸 ===")
print(f"Income ≈ {model_linear.intercept_:.2f} + "
      f"{model_linear.coef_[0]:.2f}·Education + "
      f"{model_linear.coef_[1]:.2f}·Seniority")
print(f"僅需估計 3 個參數（β₀, β₁, β₂）")

4.2 非參數化方法（Non-parametric）

不假設 f 特定函數形式，讓資料「自己說話」：KNN、薄板樣條、決策樹。

from sklearn.neighbors import KNeighborsRegressor
from sklearn.metrics import mean_squared_error

model_knn = KNeighborsRegressor(n_neighbors=3).fit(X2, y2)
y_pred_linear = model_linear.predict(X2)
y_pred_knn = model_knn.predict(X2)
print(f"訓練 MSE — Linear: {mean_squared_error(y2, y_pred_linear):.3f}")
print(f"訓練 MSE — KNN(k=3):  {mean_squared_error(y2, y_pred_knn):.3f}")

4.3 參數化 vs 非參數化：全面比較

4.4 過度擬合（Overfitting）

# 三種不同平滑度的 spline：從「太簡單」到「過度擬合」
from sklearn.preprocessing import SplineTransformer

np.random.seed(42)
X_s = np.sort(np.random.uniform(0, 1, 30))
true_f_s = lambda x: np.sin(2*np.pi*x) + 0.3*np.cos(6*np.pi*x)
Y_s = true_f_s(X_s) + np.random.normal(0, 0.2, 30)
X_plot = np.linspace(0, 1, 500).reshape(-1, 1)

fig, axes = plt.subplots(1, 3, figsize=(16, 5))
for ax, n_k, title in zip(axes, [3, 8, 25],
    ['Smooth (少節點) → 欠擬合',
     'Moderate → 良好',
     'Rough (多節點) → 過度擬合！']):
    spline = SplineTransformer(n_knots=n_k, degree=3)
    m = LinearRegression().fit(spline.fit_transform(X_s.reshape(-1,1)), Y_s)
    ax.scatter(X_s, Y_s, c='red', s=30, zorder=5)
    ax.plot(X_plot, m.predict(spline.transform(X_plot)), 'b-', lw=2)
    ax.plot(X_plot, true_f_s(X_plot), 'k--', lw=1, alpha=0.4, label='True f')
    ax.set_title(title); ax.legend(fontsize=8)
plt.suptitle('Figure 2.5/2.6 概念：非參數化方法的平滑度取捨', fontsize=14)
plt.tight_layout(); plt.show()

右圖：訓練誤差為零（完美擬合訓練資料），但新資料預測差——這就是過度擬合。

五、靈活度與可解釋性的取捨

靈活度低 ←——————————————————————→ 靈活度高
可解釋性高                                  可解釋性低

Least Squares ── Lasso ── GAM ── Trees ── Bagging/Boosting ── SVM ── Deep Learning

# 不同多項式階數的靈活度展示
fig, axes = plt.subplots(1, 3, figsize=(16, 5))
np.random.seed(123)
X_f = np.sort(np.random.uniform(0, 1, 50))
Y_f = np.sin(2*np.pi*X_f) + np.random.normal(0, 0.3, 50)

for ax, deg, title in zip(axes, [1, 3, 15],
    ['d=1: 線性（低靈活度／高可解釋性）',
     'd=3: 適中',
     'd=15: 高靈活度／低可解釋性 → 過度擬合']):
    poly = PolynomialFeatures(degree=deg)
    m = LinearRegression().fit(poly.fit_transform(X_f.reshape(-1,1)), Y_f)
    ax.scatter(X_f, Y_f, c='red', s=20); Xp = np.linspace(0,1,500).reshape(-1,1)
    ax.plot(Xp, m.predict(poly.transform(Xp)), 'b-', lw=2)
    ax.plot(Xp, np.sin(2*np.pi*Xp), 'k--', lw=1, alpha=0.4, label='True f')
    ax.set_title(title, fontsize=10); ax.legend(fontsize=7)
plt.suptitle('靈活度 vs 可解釋性（Figure 2.7 概念）', fontsize=14)
plt.tight_layout(); plt.show()

核心洞察：即使目標是純預測，最靈活的模型也不一定最好——因為過度擬合。

六、監督學習 vs 非監督學習

# 非監督學習：K-means 分群
from sklearn.cluster import KMeans

np.random.seed(42)
g1 = np.random.multivariate_normal([2, 2], [[0.5, 0], [0, 0.5]], 50)
g2 = np.random.multivariate_normal([7, 8], [[0.5, 0], [0, 0.5]], 50)
g3 = np.random.multivariate_normal([5, 3], [[0.5, 0], [0, 0.5]], 50)
X_clus = np.vstack([g1, g2, g3])

km = KMeans(n_clusters=3, n_init=10, random_state=42).fit(X_clus)

fig, ax = plt.subplots(1, 1, figsize=(8, 6))
true_labels = np.repeat([0, 1, 2], 50)
ax.scatter(X_clus[:, 0], X_clus[:, 1], c=true_labels, s=30, alpha=0.7,
           edgecolors='k', linewidth=0.5, cmap='Set1')
ax.set_xlabel('X₁'); ax.set_ylabel('X₂')
ax.set_title('Clustering: 無 y 值，僅根據 X 發現群體結構\n（Figure 2.8 概念）')
plt.show()

from sklearn.metrics import adjusted_rand_score
print(f"Adjusted Rand Index: {adjusted_rand_score(true_labels, km.labels_):.3f}")

七、迴歸問題 vs 分類問題

# 迴歸 vs 分類視覺化
from sklearn.linear_model import LogisticRegression
from sklearn.datasets import make_classification

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 6))

# 迴歸
Xr = np.random.uniform(0, 10, 100).reshape(-1, 1)
yr = 3 + 2 * Xr.ravel() + np.random.normal(0, 3, 100)
mr = LinearRegression().fit(Xr, yr)
ax1.scatter(Xr, yr, c='red', alpha=0.5); xp = np.linspace(0,10,100).reshape(-1,1)
ax1.plot(xp, mr.predict(xp), 'b-', lw=2)
ax1.set_title('Regression: Y ∈ ℝ (連續)')

# 分類
Xc, yc = make_classification(n_samples=200, n_features=2, n_redundant=0,
                              n_clusters_per_class=1, random_state=42)
mc = LogisticRegression().fit(Xc, yc)
x_min, x_max = Xc[:, 0].min()-1, Xc[:, 0].max()+1
y_min, y_max = Xc[:, 1].min()-1, Xc[:, 1].max()+1
xx, yy = np.meshgrid(np.linspace(x_min, x_max, 200), np.linspace(y_min, y_max, 200))
Z = mc.predict(np.c_[xx.ravel(), yy.ravel()]).reshape(xx.shape)
ax2.contourf(xx, yy, Z, alpha=0.3, cmap='RdYlBu')
ax2.scatter(Xc[:, 0], Xc[:, 1], c=yc, edgecolor='k', s=30, cmap='RdYlBu')
ax2.set_title(f'Classification: Y ∈ {{0, 1}} (類別)\nAccuracy = {mc.score(Xc, yc):.2%}')

plt.suptitle('Regression vs Classification', fontsize=14)
plt.tight_layout(); plt.show()

八、實務應用要點

8.1 如何選擇統計學習方法？

我需要預測還是理解？
├── 預測為主 → 可用靈活模型（但要小心過度擬合！）
├── 推論為主 → 選可解釋模型（線性迴歸、Lasso）
└── 兩者兼顧 → 適中模型（GAM、決策樹）

Y 是數值還是類別？
├── 數值（連續） → 迴歸方法
└── 類別（離散） → 分類方法

我有 y 標籤嗎？
├── 有 → 監督學習（本書主軸）
└── 沒有 → 非監督學習（分群等）

8.2 常見陷阱

1. 過度擬合：模型記住了訓練資料的雜訊而非真實模式
2. 只用訓練誤差評估模型：必須用獨立的測試資料
3. 忽略不可約誤差：即使完美估計 f，預測仍有誤差 Var(ε)
4. 忽略資料外洩（Data Leakage）：訓練過程中不慎用到測試集資訊
5. 忽略多重檢定問題：測試太多假說時需修正（Bonferroni、FDR）

下一節：2.2 評估模型準確度 — 偏差-變異數權衡、訓練 vs 測試誤差、貝氏分類器