Bài 20: Polynomial Regression — model quan hệ phi tuyến

Linear Regression giả định target là tổ hợp tuyến tính của feature. Nhiều quan hệ thực tế không như vậy:

• Giá nhà theo diện tích — tới một ngưỡng diện tích, mỗi m² thêm vào không còn cộng giá tỷ lệ thuận; quan hệ cong dần.
• Doanh số theo budget quảng cáo — diminishing return. Tăng budget từ 1M → 10M tăng doanh số mạnh; từ 100M → 110M tăng rất ít.
• Độ phức tạp dự án theo năm kinh nghiệm — quan hệ S-curve, không tuyến tính ở hai đầu.
• Tốc độ phản ứng hoá học theo nhiệt độ — Arrhenius, dạng exponential.

Fit một đường thẳng vào những quan hệ này → residual có cấu trúc rõ rệt (cong, hình chữ U, hình phễu). Đây là tín hiệu Linear Regression underfit: model không đủ flexibility để nắm shape của data.

Cách xử lý đơn giản nhất, không cần đổi sang model phi tuyến phức tạp: Polynomial Regression — thêm các term \( x^2, x^3, \ldots \) làm feature mới.

Với 1 feature \( x \), model polynomial bậc \( d \) (degree \( d \)):

\[ \hat{y} = w_0 + w_1 x + w_2 x^2 + w_3 x^3 + \ldots + w_d x^d \]

Trong đó:

• \( d \) — degree, bậc cao nhất.
• \( w_0 \) — bias (intercept).
• \( w_1, \ldots, w_d \) — coefficient của từng term.

Với \( d = 1 \), model thoái về Linear Regression. Với \( d = 2 \), model là parabola. Với \( d = 3 \), là đường cong bậc ba — có thể có 1 điểm uốn. Càng tăng \( d \), đường cong càng "uốn éo" được phức tạp hơn.

Câu hỏi: model trông như không tuyến tính theo \( x \), vậy có phải linear regression nữa không?

Đáp án: vẫn là linear regression. Lý do — "linear" trong Linear Regression là linear theo parameter \( w_i \), không phải theo feature \( x \). Nhìn lại phương trình:

\[ \hat{y} = w_0 + w_1 x + w_2 x^2 + \ldots + w_d x^d \]

Mỗi \( w_i \) chỉ xuất hiện ở bậc 1 — không có \( w_i^2 \), không có \( w_1 \cdot w_2 \). Đặt \( z_1 = x, z_2 = x^2, \ldots, z_d = x^d \):

\[ \hat{y} = w_0 + w_1 z_1 + w_2 z_2 + \ldots + w_d z_d \]

Đây chính xác là Linear Regression đa biến trên \( d \) feature mới. Hệ quả thực hành:

• KHÔNG cần thay đổi algorithm — vẫn dùng OLS hay Gradient Descent như Bài 15.
• Chỉ cần biến đổi feature từ \( x \) thành \( (x, x^2, x^3, \ldots, x^d) \) trước khi đưa vào model.
• Mọi tính chất của OLS (closed-form, gauss–markov, diễn giải coefficient) giữ nguyên.

Polynomial Regression thực ra là Linear Regression với feature engineering cụ thể. Đây là lý do nó thường được dạy như "case study" của feature engineering thay vì một model riêng.

Khi data có nhiều feature, polynomial degree \( d \) không chỉ thêm các luỹ thừa từng feature — nó còn thêm các interaction term (term nhân chéo giữa các feature).

Ví dụ 2 feature \( (x_1, x_2) \), degree 2 sinh ra các term:

\[ (1, \; x_1, \; x_2, \; x_1^2, \; x_1 x_2, \; x_2^2) \]

Bao gồm:

• Bias: \( 1 \).
• Linear term: \( x_1, x_2 \).
• Quadratic thuần: \( x_1^2, x_2^2 \).
• Interaction: \( x_1 x_2 \) — đại diện cho "ảnh hưởng của \( x_1 \) phụ thuộc giá trị \( x_2 \)".

Interaction term cực kỳ hữu ích trong domain thực tế. Ví dụ: hiệu quả thuốc (\( y \)) có thể phụ thuộc liều (\( x_1 \)) và tuổi bệnh nhân (\( x_2 \)); nhưng cùng một liều tác dụng khác nhau ở các nhóm tuổi khác nhau → cần term \( x_1 x_2 \).

Số term tổng cộng với \( n \) feature, degree \( d \) (bao gồm bias):

\[ \binom{n + d}{d} = \frac{(n+d)!}{n! \, d!} \]

Vài giá trị cụ thể:

• \( n = 2, d = 2 \): \( \binom{4}{2} = 6 \) term.
• \( n = 5, d = 3 \): \( \binom{8}{3} = 56 \) term.
• \( n = 10, d = 5 \): \( \binom{15}{5} = 3003 \) term.
• \( n = 20, d = 5 \): \( \binom{25}{5} = 53130 \) term.

Số term bùng nổ tổ hợp theo \( n \) và \( d \). Đây là pitfall lớn nhất của Polynomial Regression — sẽ bàn ở mục 10.

Sklearn cung cấp PolynomialFeatures trong module sklearn.preprocessing — một transformer chuyển ma trận feature gốc thành ma trận feature mở rộng polynomial.

from sklearn.preprocessing import PolynomialFeatures
import numpy as np

X = np.array([[2, 3]])              # 1 sample, 2 feature
poly = PolynomialFeatures(degree=2)
X_poly = poly.fit_transform(X)

print(X_poly)
# [[1. 2. 3. 4. 6. 9.]]
# tương ứng: [1, x1, x2, x1^2, x1*x2, x2^2]

print(poly.get_feature_names_out(["x1", "x2"]))
# ['1' 'x1' 'x2' 'x1^2' 'x1 x2' 'x2^2']

Các tham số chính:

• degree (mặc định 2) — bậc cao nhất của polynomial.
• interaction_only (mặc định False) — nếu True, chỉ sinh interaction term, bỏ \( x_i^k \) thuần. Ví dụ \( n=2, d=2 \) với interaction_only=True → \( (1, x_1, x_2, x_1 x_2) \), bỏ \( x_1^2, x_2^2 \).
• include_bias (mặc định True) — thêm cột toàn 1. Set False nếu downstream model đã có intercept riêng (như LinearRegression trong sklearn) để tránh redundant.

Trên thực tế, với LinearRegression(fit_intercept=True) (mặc định) thì nên PolynomialFeatures(include_bias=False), vì cả hai đều thêm bias → cột 1 thừa nhưng không sai (sklearn dùng pseudo-inverse nên không crash). Set lại cho gọn:

poly = PolynomialFeatures(degree=2, include_bias=False)

Polynomial Regression hầu như luôn cần đi kèm StandardScaler. Lý do: \( x^d \) với \( d \) lớn và \( x \) ngoài khoảng \( [-1, 1] \) ra giá trị rất lớn → mất ổn định số học và khiến gradient descent không hội tụ. Scale trước (hoặc sau) polynomial đều có ý nghĩa, tuỳ context:

• Scale trước polynomial: \( x \) đã \( \mathcal{N}(0, 1) \), \( x^d \) vẫn có thể lớn (vd \( 3^{10} = 59049 \)) nhưng nhỏ hơn không scale. Cách dùng phổ biến.
• Scale sau polynomial: từng term \( x^k \) được scale riêng → tất cả về \( \mathcal{N}(0, 1) \). An toàn nhất về numeric nhưng làm coefficient mất ý nghĩa "tăng 1 đơn vị \( x \)".

Pipeline điển hình dùng make_pipeline:

from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import PolynomialFeatures, StandardScaler
from sklearn.linear_model import LinearRegression

pipe = make_pipeline(
    PolynomialFeatures(degree=3, include_bias=False),
    StandardScaler(),
    LinearRegression(),
)

pipe.fit(X_train, y_train)
y_pred = pipe.predict(X_test)

Đặt scaler sau polynomial trong trường hợp này — scale từng term polynomial. Khi nội suy ý nghĩa coefficient không quan trọng (mục tiêu chính là predict), pipeline này là default an toàn.

Degree càng cao, model càng flexible. Linh hoạt không miễn phí — phải trả bằng variance. Đây là bias–variance tradeoff:

• Degree thấp (\( d = 1 \)) — model đơn giản: bias cao (không nắm được shape thật), variance thấp (model ổn định giữa các train set khác nhau). Kết quả: underfit, train error và test error đều cao.
• Degree vừa (\( d = 2 \) đến \( 4 \)) — đủ flexible để nắm trend, đủ regular để không "nhớ" noise. Train error và test error đều thấp. Vùng sweet spot.
• Degree cao (\( d = 15 \) đến \( 20 \)) — model đi sát mọi điểm training, bias rất thấp, nhưng variance rất cao. Một thay đổi nhỏ ở train set sẽ kéo đường cong sang hình hoàn toàn khác. Kết quả: overfit, train error gần 0, test error rất cao.

Tổng test error có thể phân rã (Bài 38 sẽ chứng minh):

\[ \text{Test error} = \text{Bias}^2 + \text{Variance} + \text{Noise irreducible} \]

Bias giảm khi tăng degree, Variance tăng khi tăng degree. Tổng có dạng chữ U → có một \( d^* \) tối ưu cân bằng hai thành phần. Mục tiêu của tuning là tìm \( d^* \).

Một dataset cổ điển để thấy ba trạng thái: tạo \( y = \sin(\pi x) + \text{noise} \) với \( x \in [0, 1] \), 30 điểm. Fit polynomial degree 1, 3, 15 và quan sát:

• Degree 1 — đường thẳng đi xiên qua data. Không nắm được dạng cong sin. Cả train và test error cao → underfit.
• Degree 3 — đường cong mềm mại, đi sát shape sin. Bỏ qua noise, nắm trend. Train và test error đều thấp → fit tốt.
• Degree 15 — đường cong "uốn éo" oscillate dữ dội giữa các điểm. Đi qua gần như tất cả điểm train, nhưng giữa hai điểm liền kề lại nhảy lên/xuống vô lý. Train error rất nhỏ, test error lại cao bất thường → overfit rõ ràng.

Đặc trưng dễ nhận của polynomial overfit: hai đầu đường cong (vùng ít data) vọt lên hoặc xuống vô tận, do đặc tính của polynomial bậc cao. Đây là biểu hiện Runge's phenomenon — bậc đa thức càng cao càng dao động mạnh ở biên.

Code minh hoạ ở mục 13. Khi tự chạy, in thêm \( R^2 \) trên cả train và test cho 3 degree → con số sẽ kể câu chuyện rõ hơn cả đồ thị.

Không có công thức đóng cho \( d^* \). Ba cách thực hành chính:

• Cross-validation (Bài 39) — quét \( d \) từ 1 đến vài chục, với mỗi \( d \) chạy k-fold CV, ghi nhận mean validation score. Chọn \( d \) cho score tốt nhất. Cách chuẩn nhất, robust với randomness của split.
• Validation curve — biến thể đơn giản của CV, plot train score và validation score theo \( d \). Vùng \( d \) mà train cao và validation cũng cao là sweet spot; vùng train cao nhưng validation tụt là overfit.
• Domain knowledge — nếu hiểu vật lý/sinh học/kinh tế của bài toán, biết ngay nên parabola (\( d=2 \)) hay cubic (\( d=3 \)). Ví dụ năng lượng động học \( E = \frac{1}{2} m v^2 \) — biết chắc \( d = 2 \), không cần CV.

Quy tắc rule-of-thumb (không thay thế CV nhưng làm starting point):

• Bắt đầu với \( d = 2 \) hoặc \( 3 \). Hầu hết quan hệ phi tuyến trong tabular data fit tốt ở range này.
• \( d > 5 \) hiếm khi có lý do chính đáng — thường là tín hiệu nên đổi sang model phi tuyến khác (tree-based, kernel SVM).
• Nếu \( n \) sample nhỏ (vài trăm) và \( d \) lớn → gần như chắc chắn overfit.

Polynomial Regression có một loạt cạm bẫy ít gặp ở Linear Regression thuần:

• Overfit cực mạnh khi dataset nhỏ + degree cao — model có thể có nhiều parameter hơn cả số sample. Với \( n = 100 \) sample, 10 feature, \( d = 5 \) sinh 3003 feature — model học 3003 weight từ 100 sample. Toán không cấm, nhưng generalization gần như không thể.
• Feature explosion — như công thức tổ hợp, \( n = 10, d = 5 \) ra 3003 cột; \( n = 20, d = 5 \) ra 53130 cột. Tiêu RAM, tính chậm, train không tractable.
• Numerical instability — \( x^{15} \) với \( x = 100 \) ra \( 10^{30} \) → vượt phạm vi an toàn của float64. Ngay cả \( x = 10 \) đã ra \( 10^{15} \), gradient và \( X^T X \) gặp catastrophic cancellation. Bắt buộc scale \( x \) về quanh \( [-1, 1] \) hoặc dùng orthogonal polynomials (Chebyshev, Legendre) khi \( d \) cao.
• Multicollinearity nội tại — \( x \) và \( x^2 \) tương quan rất mạnh khi \( x > 0 \). Ma trận \( X^T X \) gần singular → coefficient không ổn định, sai số nhỏ ở data có thể đảo dấu coefficient. Tham khảo Bài 15 mục 11.
• Cần scale TRƯỚC hoặc SAU polynomial — tuỳ context (xem mục 6). Không scale là sai mặc định khi \( d \geq 3 \).

Khi dataset có nhiều feature và muốn capture phi tuyến, polynomial degree cao gần như luôn là phương án tệ. Tree-based model (Bài 28–30) hoặc kernel method (Bài 25) phù hợp hơn.

Một cách giữ flexibility (degree cao) mà vẫn giảm variance: thêm regularization vào loss. Thay LinearRegression bằng Ridge hoặc Lasso (chi tiết Bài 21):

from sklearn.linear_model import Ridge

pipe = make_pipeline(
    PolynomialFeatures(degree=10, include_bias=False),
    StandardScaler(),
    Ridge(alpha=1.0),
)
pipe.fit(X_train, y_train)

Ridge phạt \( \sum w_i^2 \) → kéo coefficient nhỏ lại, làm đường cong "mượt" hơn → giảm dao động bất thường ở vùng ít data.

Quy tắc combo phổ biến: Polynomial degree khá cao + Ridge với \( \alpha \) chọn bằng CV. Cho kết quả tốt hơn nhiều so với Polynomial degree thấp không regularization, vì model có khả năng nắm phi tuyến nhưng không tự do overfit.

Polynomial không phải lựa chọn duy nhất cho phi tuyến. Khi degree cần > 3–4 hoặc khi quan hệ có nhiều "vùng" khác nhau:

• Spline regression — chia trục \( x \) thành các đoạn, fit polynomial bậc thấp (thường \( d = 3 \)) trên mỗi đoạn với điều kiện liên tục ở điểm nối. Linh hoạt cục bộ, không bị Runge's phenomenon. Sklearn có SplineTransformer.
• Decision Tree, Random Forest, Gradient Boosting (Bài 28–30) — model phi tuyến tự nhiên qua split tree, không cần feature engineering thủ công. Mặc định mạnh trên tabular data.
• Kernel methods (SVR, Kernel Ridge) (Bài 25) — dùng kernel để mapping ngầm vào không gian cao chiều mà không cần xây feature tường minh. Kernel polynomial chính là Polynomial Regression viết lại dưới dạng dual.
• Neural Network (Series 3) — universal approximator. Tốt khi dataset lớn và quan hệ rất phức tạp, không tốt khi data ít.

Polynomial Regression đứng đúng chỗ khi: data ít chiều (1–3 feature), quan hệ smooth toàn cục, cần model nhỏ và diễn giải được. Ngoài phạm vi này, thường có lựa chọn tốt hơn.

Tạo data sin + noise, fit polynomial degree 1, 3, 15 và in \( R^2 \) trên train và test:

import numpy as np
from sklearn.preprocessing import PolynomialFeatures, StandardScaler
from sklearn.linear_model import LinearRegression
from sklearn.pipeline import make_pipeline
from sklearn.model_selection import train_test_split
from sklearn.metrics import r2_score

rng = np.random.default_rng(42)

# Sinh data: y = sin(pi*x) + noise, x trong [0, 1]
n = 60
X = rng.uniform(0, 1, size=n).reshape(-1, 1)
y = np.sin(np.pi * X.ravel()) + rng.normal(0, 0.15, size=n)

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

for d in [1, 3, 15]:
    pipe = make_pipeline(
        PolynomialFeatures(degree=d, include_bias=False),
        StandardScaler(),
        LinearRegression(),
    )
    pipe.fit(X_train, y_train)

    r2_train = r2_score(y_train, pipe.predict(X_train))
    r2_test = r2_score(y_test, pipe.predict(X_test))
    print(f"degree={d:2d}  R2 train={r2_train:+.4f}  R2 test={r2_test:+.4f}")

# degree= 1  R2 train=+0.3xx  R2 test=+0.2xx   <- underfit
# degree= 3  R2 train=+0.9xx  R2 test=+0.9xx   <- fit tot
# degree=15  R2 train=+0.9xx  R2 test=-X.XXXX  <- overfit nghiem trong

Đặc điểm quan sát:

• degree=1 — \( R^2 \) train và test gần nhau, đều thấp → underfit.
• degree=3 — \( R^2 \) train cao, test cao tương đương → fit tốt.
• degree=15 — \( R^2 \) train cao, test có thể âm hoặc rất thấp → overfit. \( R^2 \) âm có nghĩa model còn tệ hơn cả predict giá trị mean của \( y \).

Để vẽ đường fit (cần matplotlib): tạo lưới \( x \) mịn, predict, plot cùng scatter điểm train — sẽ thấy đường \( d=15 \) oscillate giữa các điểm.

Bài 1 — Fit quan hệ bậc 2 đã biết. Sinh data theo \( y = 2x^2 + 3x + 1 + \mathcal{N}(0, 0.5) \) với \( x \in [-2, 2] \), 100 sample. Fit polynomial degree 2 và đọc coef_ + intercept_. Coefficient có gần \( (3, 2) \) và intercept gần \( 1 \) không? Nếu lệch, vì sao?

Bài 2 — So sánh trực quan. Vẫn data Bài 1, fit polynomial degree 1, 2, 10. Tạo lưới \( x \) gồm 200 điểm trong \( [-2, 2] \), predict \( \hat{y} \) cho cả 3 model, plot cùng scatter điểm train. Quan sát: đường degree 1 có miss shape parabola không? Đường degree 10 có wobble ở biên không?

Bài 3 — Diabetes dataset. Load sklearn.datasets.load_diabetes() (10 feature, 442 sample), split 70/30. So sánh \( R^2 \) test giữa:

• LinearRegression() thuần.
• PolynomialFeatures(degree=2) + scale + LinearRegression().

Polynomial có cải thiện không? Đếm số feature sau khi expand. Thử thêm degree 3 — \( R^2 \) test thay đổi thế nào?

Bài 4 — Đếm số term. Với 5 input feature và PolynomialFeatures(degree=3, include_bias=True), theo công thức \( \binom{n+d}{d} \) số term bằng bao nhiêu? Verify bằng code: tạo \( X \) shape \( (1, 5) \), gọi fit_transform, đọc shape[1]. Hai con số có khớp không?

Bài 5 — Interaction only. Vẫn 2 feature \( (x_1, x_2) \), so sánh PolynomialFeatures(degree=2) và PolynomialFeatures(degree=2, interaction_only=True). Liệt kê tên feature output của cả hai bằng get_feature_names_out. Mục đích thực hành mỗi option khi nào hữu ích.

Gợi ý đáp án Bài 4: \( \binom{5+3}{3} = \binom{8}{3} = 56 \) term. Verify bằng PolynomialFeatures(degree=3).fit_transform(np.zeros((1, 5))).shape → (1, 56).

Bài 21: Ridge và Lasso Regression — hai biến thể Linear Regression có regularization. Ridge phạt \( L_2 \) của weight, Lasso phạt \( L_1 \). Cả hai giúp ổn định coefficient khi multicollinearity và giảm overfit khi degree polynomial cao — đúng vấn đề bài này để mở.