Regresi Linear secara Matematis: Geometri Kuadrat Terkecil, Teorema Gauss–Markov, Inferensi, dan Regularisasi

Oleh Andi Illa Erviani Nensi | Magister Statistika dan Sains Data

Artikel ketiga dalam seri: (1) Statistik Deskriptif: Fondasi Matematis untuk Memahami Struktur Data; (2) Dari Deskripsi ke Inferensi; (3) Regresi Linear secara Matematis

Abstrak

Regresi linear sering diajarkan sebagai resep: pasang garis, baca koefisien, lihat nilai-p. Artikel ini membahasnya sebagai objek matematis yang utuh. Penaksir kuadrat terkecil diturunkan, ditafsirkan secara geometris sebagai proyeksi ortogonal di ruang ℝⁿ, dan dikarakterisasi sifat samplingnya, termasuk bukti teorema Gauss–Markov dan distribusi eksak di bawah normalitas. Pembahasan kemudian bergerak ke sisi yang jarang disentuh kuliah pengantar: kondisi numerik dan dekomposisi nilai singular, diagnostik berbasis leverage, galat baku robust terhadap heteroskedastisitas, regularisasi (ridge dan lasso) sebagai pengelolaan trade-off bias–varians, kaitan dengan kemungkinan maksimum dan regresi median, serta bias variabel terabaikan sebagai batas antara prediksi dan inferensi kausal. Seluruh konsep diikat oleh satu contoh numerik berukuran kecil yang dapat dihitung tangan, dan oleh satu benang merah dari dua artikel sebelumnya: rata-rata sampel adalah kasus paling sederhana dari regresi.

Kata kunci: kuadrat terkecil, proyeksi ortogonal, Gauss–Markov, matriks topi, multikolinearitas, regularisasi, bias variabel terabaikan

1. Pendahuluan

Artikel pertama dalam seri ini menunjukkan bahwa rata-rata sampel adalah solusi dari masalah meminimumkan jumlah kuadrat selisih, dan bahwa median adalah solusi untuk jumlah selisih absolut. Artikel kedua menunjukkan bagaimana kedua penaksir itu memiliki distribusi sampling, galat baku, dan interval kepercayaan. Regresi linear adalah perluasan alamiah dari kedua gagasan tersebut: alih-alih menaksir satu konstanta yang meringkas data, kita menaksir fungsi dari kovariat yang meringkas rata-rata bersyarat dari respons.

Bagi pembaca dengan latar belakang matematika, regresi linear menarik karena tiga alasan. Pertama, seluruh teorinya berdiri di atas aljabar linear yang bersih: proyeksi, ruang kolom, nilai singular. Kedua, sifat-sifat intinya (ketaktakbiasan, Gauss–Markov, distribusi t dan F) dapat dibuktikan secara lengkap dengan alat yang elementer. Ketiga, regresi linear adalah prototipe dari hampir semua metode pemodelan lanjutan: model linear terampat, regresi terpenalti, proses Gauss, bahkan jaringan saraf pada lapisan terakhirnya. Memahami regresi linear secara mendalam berarti memiliki bahasa untuk memahami semuanya.

Struktur artikel: Bagian 2 menetapkan model dan asumsi. Bagian 3 menurunkan penaksir dan tafsir geometrisnya. Bagian 4 membahas sifat sampling dan inferensi. Bagian 5 menyajikan contoh numerik terpadu. Bagian 6 dan 7 membahas aspek numerik dan diagnostik. Bagian 8 membahas regularisasi. Bagian 9 dan 10 mengaitkan regresi dengan kemungkinan maksimum dan kausalitas. Bagian 11 menyorot praktik sains data, dan Bagian 12 menyimpulkan.

Notasi. Vektor bersifat kolom. Transpos dinyatakan dengan prima, misalnya X′. Iₙ adalah matriks identitas berukuran n × n, dan ‖·‖ adalah norma Euklides. Untuk matriks simetris A dan B, A ⪰ B berarti A − B positif semidefinit.

2. Formulasi Model dan Asumsi

Misalkan terdapat n pengamatan. Model regresi linear berganda menyatakan:

y = Xβ + ε

dengan y ∈ ℝⁿ vektor respons, X matriks desain berukuran n × p, β ∈ ℝᵖ vektor parameter, dan ε ∈ ℝⁿ vektor galat. Bila model memiliki intersep, kolom pertama X adalah vektor satu. “Linear” merujuk pada linearitas dalam parameter β, bukan dalam kovariat asli; regresi polinomial, spline, dan fungsi basis lain tetap termasuk model linear karena kolom X boleh berupa transformasi nonlinear dari kovariat.

Hasil-hasil dalam artikel ini bergantung pada lima asumsi berikut. Memisahkannya dengan jelas penting, karena tiap kesimpulan hanya memerlukan sebagian dari asumsi tersebut.

•         (A1) Spesifikasi linear: y = Xβ + ε.

•         (A2) Peringkat penuh: rank(X) = p ≤ n, sehingga X′X dapat dibalik.

•         (A3) Eksogenitas ketat: E[ε | X] = 0, atau setara E[y | X] = Xβ. Regresi adalah model bagi rata-rata bersyarat.

•         (A4) Galat sferis: Var(ε | X) = σ²Iₙ, yakni homoskedastisitas dan galat tak berkorelasi.

•         (A5) Normalitas: ε | X ~ N(0, σ²Iₙ).

Peta kebutuhan asumsi: ketaktakbiasan memerlukan (A1)–(A3); teorema Gauss–Markov memerlukan (A1)–(A4); distribusi eksak uji t dan F memerlukan (A1)–(A5). Tanpa (A5), pendekatan asimtotik tetap dimungkinkan, sebagaimana dibahas pada Bagian 4.4.

3. Penaksir Kuadrat Terkecil dan Geometrinya

3.1 Penurunan

Penaksir kuadrat terkecil biasa (ordinary least squares, OLS) meminimumkan jumlah kuadrat galat S(β) = ‖y − Xβ‖². Gradiennya adalah −2X′(y − Xβ). Menyamakannya dengan nol menghasilkan persamaan normal:

X′X β = X′y

Di bawah (A2), X′X positif definit, sehingga Hessian 2X′X positif definit, solusinya unik, dan merupakan minimum global:

β̂ = (X′X)⁻¹ X′y

Nilai terprediksi dan residu berturut-turut adalah ŷ = Xβ̂ dan e = y − ŷ.

3.2 Tafsir geometris: proyeksi ortogonal

Persamaan normal dapat ditulis ulang sebagai X′(y − Xβ̂) = X′e = 0. Residu ortogonal terhadap setiap kolom X. Dengan demikian ŷ adalah proyeksi ortogonal y pada ruang kolom C(X), subruang berdimensi p di ℝⁿ, dan e adalah komponen y yang tegak lurus terhadap C(X). Penaksiran OLS, pada dasarnya, adalah mencari titik di C(X) yang paling dekat dengan y.

Operator proyeksinya adalah matriks topi (hat matrix):

H = X(X′X)⁻¹X′, sehingga ŷ = Hy dan e = (Iₙ − H)y = My

dengan M = Iₙ − H. Sifat-sifat berikut menjadi dasar hampir semua turunan selanjutnya:

•         H dan M simetris dan idempoten (H² = H, M² = M), sehingga keduanya proyeksi ortogonal.

•         HX = X dan MX = 0: H mempertahankan ruang kolom, M meniadakannya.

•         HM = 0: komponen terproyeksi dan residu saling ortogonal.

•         Nilai eigen H hanya 0 atau 1, sehingga tr(H) = rank(H) = p dan tr(M) = n − p.

Elemen diagonal h_ii disebut leverage pengamatan ke-i. Karena H idempoten, 0 ≤ h_ii ≤ 1 dan Σᵢ h_ii = p. Jadi rata-rata leverage adalah p/n, dan pengamatan dengan h_ii jauh di atas p/n memiliki pengaruh potensial besar pada garis regresi (Bagian 7).

3.3 Dekomposisi variasi dan R²

Dengan teorema Pythagoras pada ŷ ⟂ e: ‖y‖² = ‖ŷ‖² + ‖e‖². Bila model memiliki intersep, rata-rata residu nol dan identitas ini berlaku pada skor terpusat:

SST = SSR + SSE, dengan SST = Σ(yᵢ − ȳ)², SSR = Σ(ŷᵢ − ȳ)², SSE = Σeᵢ²

Koefisien determinasi R² = SSR/SST = 1 − SSE/SST adalah proporsi variasi yang dijelaskan model; secara ekuivalen, R² adalah kuadrat korelasi sampel antara y dan ŷ, atau kuadrat kosinus sudut antara y terpusat dan proyeksinya. Dua catatan kehati-hatian: R² tidak pernah turun ketika regresor ditambahkan (ruang kolom yang lebih besar tidak mungkin lebih jauh dari y), sehingga ia tidak layak dipakai sendiri untuk memilih model; dan R² tinggi tidak menjamin spesifikasi yang benar atau kemampuan prediksi di luar sampel. Versi terkoreksi, R²_adj = 1 − [SSE/(n − p)] / [SST/(n − 1)], menghukum kompleksitas, tetapi bukan pengganti validasi luar sampel.

3.4 Kasus khusus yang mengikat seri artikel

Jika X hanya berisi kolom satu, yaitu model y = μ + ε, maka X′X = n dan β̂ = (1/n) Σyᵢ = ȳ. Matriks topinya H = (1/n)11′ memetakan y ke vektor yang seluruh elemennya ȳ. SSE = Σ(yᵢ − ȳ)², dan penaksir varians SSE/(n − p) dengan p = 1 menjadi SSE/(n − 1): koreksi Bessel dari artikel pertama muncul sebagai kasus khusus dari tr(M) = n − p. Mean sampel adalah regresi pada konstanta, dan derajat kebebasan n − 1 adalah dimensi ruang residu setelah satu dimensi diambil oleh konstanta.

3.5 Teorema Frisch–Waugh–Lovell

Bagi X menjadi dua blok, X = [X₁ X₂], dengan koefisien (β₁, β₂). Misalkan M₁ = I − X₁(X₁′X₁)⁻¹X₁′ adalah matriks yang meniadakan pengaruh X₁. Maka:

β̂₂ = (X₂′M₁X₂)⁻¹ X₂′M₁ y

Artinya, koefisien β̂₂ sama dengan hasil meregresikan residu y terhadap X₁ pada residu kolom X₂ terhadap X₁. Inilah makna eksak dari frasa “pengaruh X₂ dengan mengendalikan X₁”: koefisien regresi berganda mengukur hubungan setelah bagian dari X₂ yang dapat dijelaskan oleh X₁ disingkirkan. Teorema ini juga menjadi kunci penurunan faktor inflasi varians pada Bagian 6.

4. Sifat Sampling dan Inferensi

4.1 Rata-rata dan varians

Substitusi y = Xβ + ε memberi β̂ = β + (X′X)⁻¹X′ε. Di bawah (A1)–(A3), E[β̂ | X] = β, sehingga OLS tak bias. Dengan tambahan (A4):

Var(β̂ | X) = σ² (X′X)⁻¹

Untuk model konstanta saja, rumus ini menghasilkan σ²/n, kembali ke galat baku rata-rata sampel pada artikel kedua.

4.2 Teorema Gauss–Markov

Teorema. Di bawah (A1)–(A4), untuk setiap penaksir β̃ yang linear dalam y dan tak bias, Var(β̃ | X) − Var(β̂ | X) positif semidefinit. Akibatnya, untuk setiap vektor c, c′β̂ memiliki varians terkecil di antara penaksir linear tak bias bagi c′β. Dengan kata lain, OLS adalah best linear unbiased estimator (BLUE).

Bukti. Tulis β̃ = Cy dengan C berukuran p × n, dan misalkan D = C − (X′X)⁻¹X′. Ketaktakbiasan untuk setiap β mensyaratkan CX = I, yaitu I + DX = I, sehingga DX = 0. Maka:

Var(β̃) = σ² CC′ = σ² [ (X′X)⁻¹ + (X′X)⁻¹X′D′ + DX(X′X)⁻¹ + DD′ ] = σ² [ (X′X)⁻¹ + DD′ ]

karena suku silang bernilai nol akibat DX = 0. Karena DD′ selalu positif semidefinit, Var(β̃) − Var(β̂) = σ²DD′ ⪰ 0. ∎

Hal yang sering terlewat: teorema ini tidak mensyaratkan normalitas, tetapi membatasi diri pada penaksir linear dan tak bias. Penaksir bias (misalnya ridge) atau tak linear (misalnya median regresi) dapat mengalahkan OLS dalam MSE atau ketahanan, sebagaimana ditunjukkan pada Bagian 8 dan 9.

4.3 Penaksiran σ² dan distribusi eksak

Residu dapat ditulis e = My = Mε karena MX = 0. Maka SSE = e′e = ε′Mε, dan:

E[SSE | X] = σ² tr(M) = σ² (n − p)

Sehingga s² = SSE/(n − p) adalah penaksir tak bias bagi σ². Perhatikan keselarasannya dengan koreksi Bessel: pembagi n − p adalah dimensi ruang residu.

Bila (A5) berlaku, tiga fakta mengikuti sifat proyeksi ortogonal pada vektor normal:

1.       β̂ | X ~ N(β, σ²(X′X)⁻¹).

2.       SSE/σ² ~ χ² dengan n − p derajat kebebasan.

3.       β̂ dan SSE saling bebas, karena Hε dan Mε adalah proyeksi ke subruang ortogonal dari vektor normal sferis.

Dari ketiganya diperoleh statistik uji eksak. Untuk koefisien ke-j, dengan SE(β̂ⱼ) = s √[(X′X)⁻¹ⱼⱼ]:

(β̂ⱼ − βⱼ) / SE(β̂ⱼ) ~ t dengan n − p derajat kebebasan

Untuk menguji q kendala linear (misalnya q koefisien sekaligus nol), uji F membandingkan jumlah kuadrat galat model terkendala (SSE_R) dan tak terkendala (SSE_U):

F = [ (SSE_R − SSE_U) / q ] / [ SSE_U / (n − p) ] ~ F dengan (q, n − p) derajat kebebasan

Untuk q = 1, F = t², sehingga kedua uji setara.

4.4 Tanpa normalitas: pendekatan asimtotik

Bila (A5) tidak berlaku tetapi galat i.i.d. dengan varians hingga dan X′X/n → Q positif definit, Teorema Limit Pusat untuk jumlah terboboti memberikan √n (β̂ − β) → N(0, σ²Q⁻¹) dalam distribusi. Uji dan interval berbasis t dan F tetap valid secara asimtotik, sehingga untuk sampel besar normalitas galat bukan syarat penting. Kualitas pendekatannya, seperti pada artikel kedua, bergantung pada kemencengan dan ekor distribusi galat serta pada ada-tidaknya observasi dengan leverage ekstrem.

4.5 Interval untuk rata-rata respons dan interval prediksi

Untuk titik kovariat baru x₀, ada dua pertanyaan yang berbeda dan sering tertukar:

Interval kepercayaan untuk E[y | x₀]: x₀′β̂ ± t · s · √( x₀′(X′X)⁻¹x₀ )

Interval prediksi untuk y₀ yang baru: x₀′β̂ ± t · s · √( 1 + x₀′(X′X)⁻¹x₀ )

dengan t = t_{n−p; 1−α/2}. Interval prediksi selalu lebih lebar karena ia memuat dua sumber ketidakpastian: galat penaksiran garis (varians x₀′(X′X)⁻¹x₀ · σ²) dan galat individual ε₀ (varians σ²) yang tidak menyusut meski n bertambah. Keduanya melebar saat x₀ menjauh dari pusat data, karena x₀′(X′X)⁻¹x₀ membesar; inilah dasar kuantitatif dari peringatan bahwa ekstrapolasi berisiko.

5. Contoh Numerik Terpadu

Agar seluruh konsep di atas dapat diperiksa dengan tangan, perhatikan lima pengamatan berikut: x = (1, 2, 3, 4, 5) dan y = (2, 4, 5, 4, 5). Modelnya y = β₀ + β₁x + ε, sehingga n = 5 dan p = 2.

5.1 Penaksiran

Matriks desain memiliki kolom satu dan kolom x. Diperoleh:

X′X = [[5, 15], [15, 55]], X′y = [20, 66], det(X′X) = 5·55 − 15² = 50

(X′X)⁻¹ = (1/50)·[[55, −15], [−15, 5]]

Maka β̂ = (1/50)·[55·20 − 15·66, −15·20 + 5·66] = (1/50)·[110, 30] = [2,2; 0,6]. Garis regresinya ŷ = 2,2 + 0,6x. Hasil yang sama diperoleh dari rumus klasik: β̂₁ = S_xy/S_xx = 6/10 = 0,6 dan β̂₀ = ȳ − β̂₁x̄ = 4 − 0,6·3 = 2,2.

x	y	ŷ	residu e	leverage h
1	2	2,8	−0,8	0,6
2	4	3,4	0,6	0,3
3	5	4,0	1,0	0,2
4	4	4,6	−0,6	0,3
5	5	5,2	−0,2	0,6

Pemeriksaan sifat-sifat teoretis: jumlah residu = 0 (karena ada intersep); Σxᵢeᵢ = −0,8 + 1,2 + 3,0 − 2,4 − 1,0 = 0, sesuai X′e = 0; dan jumlah leverage = 0,6 + 0,3 + 0,2 + 0,3 + 0,6 = 2 = p, sesuai tr(H) = p. Leverage dihitung dari hᵢ = 1/n + (xᵢ − x̄)²/S_xx.

5.2 Dekomposisi variasi

SST = Σ(yᵢ − 4)² = 4 + 0 + 1 + 0 + 1 = 6. SSE = 0,64 + 0,36 + 1,00 + 0,36 + 0,04 = 2,4, sehingga SSR = 6 − 2,4 = 3,6 dan R² = 0,6. Perhatikan bahwa SST = 6 adalah SSE model yang hanya berisi konstanta (yaitu mean ȳ = 4, kasus Bagian 3.4), sehingga regresi menurunkan jumlah kuadrat galat dari 6 menjadi 2,4.

5.3 Inferensi

Penaksir varians galat adalah s² = SSE/(n − p) = 2,4/3 = 0,8, atau s ≈ 0,894. Matriks kovarians koefisien adalah s²(X′X)⁻¹, sehingga:

•         Var(β̂₁) = 0,8·5/50 = 0,08, atau SE(β̂₁) ≈ 0,283.

•         Var(β̂₀) = 0,8·55/50 = 0,88, atau SE(β̂₀) ≈ 0,938.

•         Kovarians(β̂₀, β̂₁) = −0,8·15/50 = −0,24. Tanda negatif mencerminkan bahwa garis dengan kemiringan lebih besar harus memiliki intersep lebih kecil agar tetap melewati pusat data (x̄, ȳ).

Statistik uji untuk H₀: β₁ = 0 adalah t = 0,6/0,283 ≈ 2,12 dengan 3 derajat kebebasan. Nilai kritis dua sisi t₃ adalah 2,353 pada taraf 10% dan 1,638 pada taraf 20%, sehingga nilai-p terletak antara 0,10 dan 0,20. Interval kepercayaan 95% untuk β₁ adalah 0,6 ± 3,182·0,283 ≈ (−0,30; 1,50), yang memuat nol. Uji F untuk kendala yang sama memberikan F = [(6 − 2,4)/1]/[2,4/3] = 4,5 = t², dibandingkan nilai kritis F(1, 3) = 10,13 pada α = 5%.

Bukti empiris bahwa R² = 0,6 tidak otomatis berarti “signifikan”: dengan n sekecil ini, garis yang “menjelaskan 60% variasi” masih konsisten dengan kemiringan nol.

5.4 Rata-rata respons versus prediksi

Di titik x₀ = 3 (pusat data), x₀′(X′X)⁻¹x₀ = 1/5 = 0,2, sehingga SE rata-rata respons = 0,894·√0,2 ≈ 0,4. Di titik x₀ = 6 (ekstrapolasi), ŷ₀ = 2,2 + 0,6·6 = 5,8 dan x₀′(X′X)⁻¹x₀ = 1/5 + (6 − 3)²/10 = 1,1. Maka:

•         Interval kepercayaan 95% untuk E[y | x = 6]: 5,8 ± 3,182·0,894·√1,1 ≈ 5,8 ± 2,99, yaitu (2,82; 8,79).

•         Interval prediksi 95% untuk y baru: 5,8 ± 3,182·0,894·√2,1 ≈ 5,8 ± 4,12, yaitu (1,68; 9,92).

Interval di x₀ = 6 jauh lebih lebar daripada di pusat data, dan interval prediksi lebih lebar lagi, persis seperti prediksi teori pada Bagian 4.5.

6. Kondisi Numerik, Multikolinearitas, dan Dekomposisi Nilai Singular

6.1 Multikolinearitas dan faktor inflasi varians

Dengan teorema Frisch–Waugh–Lovell, varians koefisien ke-j dapat ditulis:

Var(β̂ⱼ | X) = σ² / [ (1 − Rⱼ²) · Σᵢ (xᵢⱼ − x̄ⱼ)² ]

dengan Rⱼ² adalah koefisien determinasi dari regresi kolom ke-j terhadap kolom-kolom lainnya. Faktor VIFⱼ = 1/(1 − Rⱼ²) menyatakan seberapa besar varians membengkak akibat kolinearitas. Dua regresor dengan korelasi 0,95 menghasilkan VIF = 1/(1 − 0,9025) ≈ 10,3, yaitu galat baku sekitar 3,2 kali lebih besar daripada bila regresor tersebut ortogonal.

Dua catatan penting. Pertama, kolinearitas tidak melanggar asumsi apa pun selama X berperingkat penuh: OLS tetap tak bias dan tetap BLUE; yang terpengaruh adalah presisi koefisien individual, bukan kualitas prediksi dalam ruang yang dilingkupi data. Kedua, kolinearitas sempurna (misalnya semua variabel boneka dari suatu kategori beserta intersep) melanggar (A2) dan membuat X′X singular; inilah jebakan variabel boneka yang klasik.

6.2 Perspektif nilai singular

Dekomposisi nilai singular X = UDV′, dengan U (n × p) berkolom ortonormal, D = diag(d₁ ≥ … ≥ dₚ > 0), dan V (p × p) ortogonal, memberikan bentuk yang lebih jernih:

β̂ = V D⁻¹ U′y, Var(β̂) = σ² V D⁻² V′

Varians pada arah singular ke-k sebanding dengan 1/dₖ². Arah dengan nilai singular kecil, yaitu arah di mana data hampir tidak bervariasi, adalah arah di mana koefisien tertaksir sangat tidak pasti. Bilangan kondisi κ(X) = d₁/dₚ merangkum tingkat ketidakstabilan ini, dan akan menjadi kunci pada regularisasi di Bagian 8.

6.3 Mengapa tidak memecahkan persamaan normal secara langsung

Secara komputasi, membentuk X′X lalu membalikkannya berbahaya, karena κ(X′X) = κ(X)². Bila κ(X) = 10⁵, maka κ(X′X) = 10¹⁰, dan dalam aritmetika presisi ganda (sekitar 16 digit) hanya tersisa sekitar enam digit bermakna. Implementasi yang baik memakai dekomposisi QR (X = QR, lalu menyelesaikan Rβ = Q′y) atau SVD, yang hanya bergantung pada κ(X). Memusatkan dan menskalakan kovariat sebelum penaksiran juga memperbaiki kondisi numerik secara signifikan, terutama bila kolom berskala sangat berbeda atau berisi pangkat tinggi.

7. Diagnostik Model

7.1 Leverage, residu terstudentisasi, dan pengaruh

Residu mentah eᵢ tidak berdistribusi identik: Var(eᵢ) = σ²(1 − hᵢ). Karena itu residu yang dibandingkan antarpengamatan sebaiknya distandarkan:

rᵢ = eᵢ / ( s √(1 − hᵢ) )

Pengamatan berpengaruh adalah pengamatan yang bila dibuang mengubah hasil secara berarti. Jarak Cook menggabungkan residu dan leverage:

Dᵢ = (rᵢ² / p) · hᵢ / (1 − hᵢ)

Pada contoh Bagian 5, rᵢ berturut-turut adalah −1,41; 0,80; 1,25; −0,80; −0,35, dan jarak Cook adalah 1,50; 0,14; 0,20; 0,14; 0,09. Pengamatan pertama (x = 1) memiliki leverage tertinggi (bersama pengamatan kelima) dan residu terstandar terbesar, sehingga menjadi yang paling berpengaruh, dengan D₁ = 1,5 yang melebihi batas heuristik 1. Jika pengamatan itu dibuang, garis regresi akan berubah cukup besar. Dua pelajaran: leverage tinggi saja belum berarti berpengaruh (pengamatan kelima punya leverage sama tetapi residu kecil), dan pada n sekecil ini satu titik dapat mendikte garis.

7.2 Validasi silang tinggalkan-satu dalam bentuk tertutup

Untuk OLS, residu tinggalkan-satu dapat dihitung tanpa pemasangan ulang: eᵢ/(1 − hᵢ). Jumlah kuadratnya disebut PRESS. Pada contoh, residu tinggalkan-satu adalah −2,00; 0,86; 1,25; −0,86; −0,50, sehingga PRESS ≈ 7,28, jauh lebih besar daripada SSE = 2,4. Galat prediksi luar sampel (dalam arti tinggalkan-satu) rata-rata 7,28/5 ≈ 1,46 per pengamatan, tiga kali lebih besar daripada galat dalam sampel 2,4/5 = 0,48. Pembanding yang adil adalah model hanya-konstanta yang juga dievaluasi tinggalkan-satu: prediksi untuk pengamatan ke-i memakai mean dari empat pengamatan lainnya, sehingga residunya (yᵢ − ȳ)·n/(n − 1) dan PRESS₀ = 6·(5/4)² = 9,375. Maka R² prediktif = 1 − 7,28/9,375 ≈ 0,22, jauh lebih rendah daripada R² dalam sampel sebesar 0,6. Selisih antara 0,6 dan 0,22 ini adalah ukuran konkret dari optimisme galat dalam sampel, yang menjadi dasar semua teknik validasi pada pembelajaran mesin.

7.3 Heteroskedastisitas dan galat baku robust

Bila (A4) gagal karena Var(εᵢ | X) = σᵢ² berbeda-beda, OLS tetap tak bias tetapi galat baku klasik menjadi tidak valid. Penaksir kovarians sandwich (Eicker–Huber–White) tidak memerlukan model bagi σᵢ²:

Var̂(β̂) = (X′X)⁻¹ [ Σᵢ eᵢ² xᵢxᵢ′ ] (X′X)⁻¹ (HC0)

Varian HC3 membagi eᵢ² dengan (1 − hᵢ)² dan lebih andal pada sampel kecil. Untuk kemiringan pada contoh kita, SE klasik adalah 0,283, SE HC0 adalah 0,186, dan SE HC3 adalah 0,430. Perbedaan besar ini memperlihatkan dua hal: HC0 cenderung terlalu optimistis pada sampel kecil (ia mengabaikan bahwa residu pada titik berleverage tinggi menyusut), sedangkan HC3 mengoreksinya secara agresif. Galat baku robust adalah pendekatan asimtotik, bukan obat untuk n kecil.

7.4 Pemeriksaan lain

Diagnostik standar lain meliputi plot residu terhadap nilai terprediksi (untuk nonlinearitas dan heteroskedastisitas), plot Q–Q residu terstudentisasi (untuk kenormalan dan ekor berat), uji autokorelasi untuk data berurutan, dan DFBETAS untuk menilai pengaruh tiap pengamatan pada tiap koefisien. Diagnostik bersifat eksploratif: ia membantu memutuskan apakah suatu asumsi masuk akal, bukan membuktikan bahwa asumsi tersebut benar.

8. Regularisasi: Mengelola Trade-off Bias–Varians

Gauss–Markov menjamin bahwa OLS terbaik di antara penaksir linear tak bias. Namun artikel kedua telah menunjukkan bahwa menerima sedikit bias dapat menurunkan MSE. Pada regresi, kondisi di mana hal ini paling menguntungkan adalah ketika X hampir kolinear (nilai singular kecil), atau ketika p mendekati atau melampaui n.

8.1 Regresi ridge

Ridge menambahkan penalti kuadrat pada koefisien:

β̂_λ = argmin_β { ‖y − Xβ‖² + λ‖β‖² } = (X′X + λI)⁻¹ X′y, λ > 0

Untuk λ > 0, matriks X′X + λI selalu dapat dibalik, bahkan ketika p > n. Dalam SVD, nilai terprediksinya adalah ŷ_λ = Σₖ uₖ [ dₖ²/(dₖ² + λ) ] uₖ′y. Setiap komponen utama diskalakan oleh faktor dₖ²/(dₖ² + λ) ≤ 1: komponen dengan nilai singular besar hampir tidak disusutkan, sedangkan komponen dengan nilai singular kecil (arah yang paling tidak pasti pada Bagian 6.2) disusutkan paling kuat. Derajat kebebasan efektif adalah tr(H_λ) = Σₖ dₖ²/(dₖ² + λ) ≤ p, yang menurun dari p (pada λ = 0) menuju 0 (pada λ → ∞).

Di bawah model linear dengan galat sferis, β̂_λ bias: E[β̂_λ] = β − λ(X′X + λI)⁻¹β, tetapi variansnya lebih kecil daripada OLS. Hoerl dan Kennard (1970) membuktikan bahwa selalu ada λ > 0 sehingga MSE total ridge lebih kecil daripada OLS. Dengan kata lain, OLS tidak dapat diterima (inadmissible) di antara penaksir sembarang pada kriteria MSE, meskipun ia BLUE.

Ilustrasi kasus ortonormal. Bila X′X = I dengan satu koefisien, maka β̂_λ = β̂/(1 + λ). Variansnya σ²/(1 + λ)² dan bias kuadratnya λ²β²/(1 + λ)², sehingga:

MSE(λ) = (σ² + λ²β²) / (1 + λ)²

Menurunkan terhadap λ dan menyamakan dengan nol memberikan λ* = σ²/β², dengan MSE minimum σ²β²/(β² + σ²), yang selalu lebih kecil daripada σ² (MSE OLS). Untuk σ² = 1 dan β = 1, λ* = 1 dan MSE turun dari 1 menjadi 0,5. Perhatikan satu kendala praktis: λ* bergantung pada β yang tidak diketahui, sehingga λ dalam praktik dipilih lewat validasi silang.

Tafsir Bayes. Ridge adalah modus posterior (sekaligus rata-rata posterior) bila galat normal dan prior β ~ N(0, τ²I), dengan λ = σ²/τ². Regularisasi bukan sekadar trik komputasi: ia setara dengan memasukkan keyakinan awal bahwa koefisien tidak terlalu besar.

8.2 Lasso

Lasso mengganti penalti kuadrat dengan penalti norma-ℓ₁ (Tibshirani, 1996):

β̂_lasso = argmin_β { ½‖y − Xβ‖² + λ‖β‖₁ }

Berbeda dengan ridge, penalti ℓ₁ menghasilkan solusi jarang (sparse): sebagian koefisien tepat nol. Penjelasan geometrisnya, bola ℓ₁ memiliki sudut-sudut pada sumbu koordinat sehingga kontur galat kuadrat cenderung menyentuhnya tepat di sumbu, sedangkan bola ℓ₂ yang mulus tidak punya sudut. Pada kasus ortonormal (X′X = I), solusi lasso mempunyai bentuk tertutup, yaitu soft-thresholding:

β̂ⱼ = sign(zⱼ) · max( |zⱼ| − λ, 0 ), dengan z = X′y adalah koefisien OLS

Sedangkan ridge menyusutkan secara proporsional (zⱼ/(1 + λ)), lasso menyusutkan secara konstan lalu memotong koefisien kecil menjadi nol. Tafsir Bayesnya adalah modus posterior dengan prior Laplace. Catatan penting: lasso menyusutkan koefisien besar juga (bias), tidak memiliki rumus tertutup di kasus umum (dipecahkan dengan coordinate descent), dan pada regresor yang berkorelasi tinggi cenderung memilih satu di antaranya secara agak sembarang. Elastic net, yang menggabungkan penalti ℓ₁ dan ℓ₂, mengatasi sebagian masalah tersebut. Kebanyakan implementasi mensyaratkan kovariat distandardisasi lebih dulu, karena penalti bergantung pada skala kolom.

8.3 Paradoks Stein

Hasil paling mengejutkan di bidang ini: untuk X ~ N(θ, Iₚ) dengan p ≥ 3, penaksir “natural” X (yaitu kemungkinan maksimum) tidak dapat diterima di bawah kerugian kuadrat. Penaksir James–Stein, [1 − (p − 2)/‖X‖²]·X, memiliki risiko yang lebih kecil untuk setiap θ (James dan Stein, 1961). Penyusutan yang bermanfaat bahkan ketika komponen-komponennya tidak berhubungan sama sekali. Hasil ini adalah versi paling murni dari pelajaran yang berulang di seluruh seri ini: ketaktakbiasan bukan keutamaan, dan MSE yang rendah sering membutuhkan keberanian menerima bias.

8.4 Dimensi tinggi

Bila p > n, X′X singular dan persamaan normal memiliki tak hingga solusi, yang seluruhnya menginterpolasi data. Ridge pada λ → 0⁺ memilih solusi dengan norma minimum. Penelitian terkini tentang double descent (Belkin dkk., 2019) menunjukkan bahwa galat uji tidak selalu mengikuti kurva U klasik: model yang menginterpolasi dapat tetap menggeneralisasi dengan baik pada kondisi tertentu. Ini memperkaya, bukan membatalkan, kerangka bias–varians.

9. Regresi, Kemungkinan Maksimum, dan Pilihan Fungsi Kerugian

Di bawah (A5), log-kemungkinan model adalah ℓ(β, σ²) = −(n/2) log(2πσ²) − ‖y − Xβ‖²/(2σ²). Memaksimumkan terhadap β setara dengan meminimumkan jumlah kuadrat, sehingga OLS adalah MLE bila galat normal. MLE untuk σ² adalah SSE/n, yang bias; s² = SSE/(n − p) adalah koreksi tak bias, ulangan dari kisah penyebut pada artikel pertama dan kedua.

Kaitan ini menjelaskan bahwa memilih fungsi kerugian setara dengan memilih asumsi tentang galat:

•         Galat normal → kerugian kuadrat → regresi rata-rata (OLS). Efisien di bawah normalitas, rapuh terhadap pencilan.

•         Galat Laplace → kerugian absolut → regresi median (LAD). Penaksir bagi median bersyarat, lanjutan langsung dari fakta bahwa median meminimumkan Σ|xᵢ − a| pada artikel pertama. Tidak memiliki rumus tertutup, tetapi dapat dipecahkan dengan pemrograman linear.

•         Regresi kuantil (Koenker dan Bassett, 1978) memperumumnya dengan kerugian ρ_τ(u) = u·(τ − 1{u < 0}), yang menaksir kuantil bersyarat ke-τ. Dengan beberapa nilai τ, kita memperoleh gambaran seluruh distribusi bersyarat, bukan hanya rata-ratanya.

•         Kerugian Huber berperilaku kuadratik untuk residu kecil dan linear untuk residu besar, menjadi kompromi antara efisiensi OLS dan ketahanan LAD (Huber, 1964).

Satu catatan tentang ketahanan, mengacu pada konsep breakdown point artikel pertama: OLS memiliki breakdown point 1/n. LAD dan penaksir Huber tahan terhadap pencilan pada arah y, tetapi tidak terhadap titik dengan leverage tinggi pada arah x. Menghadapi keduanya diperlukan penaksir berbreakdown tinggi, seperti LTS atau MM-estimator, yang dapat mencapai breakdown point mendekati 50%.

Kerangka kemungkinan maksimum juga menjadi jembatan menuju model linear terampat (generalized linear models): regresi logistik dan Poisson memiliki likelihood dari keluarga eksponensial, dan algoritma penaksirannya (iteratively reweighted least squares) pada hakikatnya memecahkan serangkaian masalah kuadrat terkecil terboboti, memakai ulang seluruh aljabar proyeksi pada artikel ini.

10. Bias Variabel Terabaikan dan Batas Kausalitas

Asumsi (A3), E[ε | X] = 0, adalah asumsi yang paling berat secara substantif, sekaligus yang paling tidak dapat diuji dari data. Konsekuensi pelanggarannya dapat dikuantifikasi untuk kasus paling sederhana. Misalkan model sebenarnya adalah y = β₁x₁ + β₂x₂ + u dengan E[u | x₁, x₂] = 0, tetapi peneliti meregresikan y hanya terhadap x₁. Dengan mensubstitusikan dan memakai FWL, penaksir kemiringan konvergen ke:

β̃₁ → β₁ + β₂ · δ, dengan δ = Cov(x₁, x₂) / Var(x₁)

yaitu koefisien dari regresi x₂ terhadap x₁. Bias ini tidak menghilang ketika n → ∞: ia bias sistematis, bukan variabilitas sampling. Contoh: bila β₁ = 1, β₂ = 2, dan δ = 0,5, maka β̃₁ = 1 + 2·0,5 = 2, dua kali nilai sebenarnya. Tanda bias ditentukan oleh hasil kali tanda β₂ dan δ, dan inilah mekanisme matematis di balik paradoks Simpson yang disinggung pada artikel pertama: ketika variabel pengganggu diabaikan, hubungan agregat dapat berbalik arah.

Beberapa pelajaran untuk praktik:

•         Koefisien regresi adalah ukuran asosiasi bersyarat, bukan efek kausal. Interpretasi kausal memerlukan (A3) yang dibenarkan oleh desain (misalnya randomisasi) atau oleh asumsi struktural yang eksplisit (misalnya kerangka diagram kausal, Pearl, 2009).

•         Menambah kovariat tidak selalu aman. Mengendalikan collider (variabel yang dipengaruhi oleh baik perlakuan maupun hasil) justru menimbulkan bias seleksi. Pemilihan kovariat harus didasarkan pada struktur kausal, bukan pada peningkatan R².

•         Galat pengukuran pada regresor menyebabkan atenuasi: bila x* diamati sebagai x = x* + w dengan w tak berkorelasi, kemiringan OLS menyusut ke β·λ dengan λ = σ²_x*/(σ²_x* + σ²_w) < 1, sehingga efek tampak lebih lemah daripada sebenarnya.

11. Implikasi bagi Praktik Sains Data

11.1 Inferensi dan prediksi adalah dua tujuan berbeda

Untuk inferensi, yang dipentingkan adalah kebenaran spesifikasi, galat baku yang valid, dan interpretasi koefisien. Untuk prediksi, yang dipentingkan adalah galat luar sampel, dan model yang “salah” tetapi stabil sering lebih baik. Kriteria seleksi mencerminkan perbedaan itu: uji t dan F bertanya apakah suatu koefisien berbeda dari nol, AIC dan BIC memperkirakan kualitas relatif dengan penalti kompleksitas, sedangkan validasi silang secara langsung mengukur galat prediktif.

Ketiga pertanyaan tersebut dapat memberi jawaban berbeda pada data yang sama. Pada contoh Bagian 5 (n = 5), dengan model kemungkinan Gaussian dan σ² ikut ditaksir (k = 3 parameter untuk model dengan kemiringan, k = 2 untuk model hanya-konstanta): AIC memilih model dengan kemiringan (selisih sekitar 2,6), BIC juga (selisih sekitar 3,0), dan PRESS pun memilihnya (7,28 dibanding 9,38 untuk model hanya-konstanta). Sebaliknya, uji t pada β₁ (nilai-p antara 0,10 dan 0,20) tidak menolak H₀, dan AIC terkoreksi sampel kecil (AICc) memilih model hanya-konstanta (10,9 dibanding 26,3). Tidak ada kontradiksi: kriteria-kriteria itu menjawab pertanyaan berbeda, dan pada n = 5 tidak satu pun dapat diandalkan. Contoh ini lebih berguna sebagai peringatan daripada sebagai petunjuk model.

11.2 Kebocoran data dan jebakan prosedur

Semua langkah yang memanfaatkan y (seleksi fitur, penentuan λ, transformasi yang dipelajari dari data) harus berada di dalam lipatan validasi silang, bukan sebelumnya. Bila tidak, galat luar sampel akan terlalu optimistis, mengulang persis pola selisih PRESS dan SSE pada Bagian 7.2. Hal yang sama berlaku pada eksplorasi bebas atas banyak spesifikasi model (forking paths): jumlah spesifikasi yang dicoba menggelembungkan galat tipe I secara efektif, sebagaimana dibahas untuk perbandingan berganda pada artikel kedua.

11.3 Standar pelaporan

Laporan regresi yang baik memuat: estimand yang ingin dijawab (asosiasi, prediksi, atau efek kausal); spesifikasi model dan alasan pemilihannya; koefisien dengan interval kepercayaan, bukan hanya bintang signifikansi; ukuran efek dalam satuan bermakna; galat baku yang sesuai dengan struktur data (robust atau berkelompok); diagnostik dan analisis sensitivitas; serta kinerja luar sampel bila tujuannya prediksi. Satu pertanyaan yang layak diajukan pada setiap analisis: asumsi mana yang paling mungkin dilanggar, dan seberapa besar kesimpulan akan berubah bila ia dilanggar?

12. Kesimpulan

Regresi linear dapat dipahami sebagai satu bangunan matematis yang utuh. Penaksir OLS adalah proyeksi ortogonal pada ruang kolom X; dari struktur proyeksi itu mengalir ketaktakbiasan, rumus varians, teorema Gauss–Markov, distribusi eksak t dan F, serta penjelasan geometris bagi derajat kebebasan n − p. Mean sampel, objek paling awal dalam seri ini, adalah kasus ketika X hanyalah vektor satu.

Di atas fondasi tersebut berdiri lapisan-lapisan lanjutan. Aspek numerik (SVD, bilangan kondisi) menjelaskan mengapa multikolinearitas merusak presisi dan mengapa persamaan normal tidak boleh dibalik begitu saja. Diagnostik berbasis leverage dan pengaruh menunjukkan bahwa satu titik dapat mendikte garis. Regularisasi mengubah ketakbiasan dari tujuan menjadi pilihan, dengan imbalan MSE yang lebih rendah. Tafsir kemungkinan maksimum mengungkap bahwa memilih kerugian kuadrat atau absolut sama dengan memilih asumsi tentang galat. Dan analisis bias variabel terabaikan menggarisbawahi bahwa matematika yang paling rapi pun tidak menggantikan pertanyaan tentang bagaimana data dihasilkan.

Bagi praktisi sains data, pesan akhirnya konsisten dengan seri ini: setiap angka dalam keluaran regresi adalah realisasi dari suatu distribusi sampling, dan setiap tafsir bergantung pada asumsi yang dapat diperiksa sebagian saja. Memahami matematika di balik keluaran itu memungkinkan kita membedakan antara apa yang dikatakan data dan apa yang hanya dikatakan model.

Daftar Pustaka

• Akaike, H. (1974). A new look at the statistical model identification. IEEE Transactions on Automatic Control, 19(6), 716–723.
• Belkin, M., Hsu, D., Ma, S., & Mandal, S. (2019). Reconciling modern machine-learning practice and the classical bias–variance trade-off. Proceedings of the National Academy of Sciences, 116(32), 15849–15854.
• Cook, R. D. (1977). Detection of influential observation in linear regression. Technometrics, 19(1), 15–18.
• Frisch, R., & Waugh, F. V. (1933). Partial time regressions as compared with individual trends. Econometrica, 1(4), 387–401.
• Golub, G. H., & Van Loan, C. F. (2013). Matrix Computations (4th ed.). Johns Hopkins University Press.
• Hastie, T., Tibshirani, R., & Friedman, J. (2009). The Elements of Statistical Learning (2nd ed.). Springer.
• Hoerl, A. E., & Kennard, R. W. (1970). Ridge regression: Biased estimation for nonorthogonal problems. Technometrics, 12(1), 55–67.
• Huber, P. J. (1964). Robust estimation of a location parameter. The Annals of Mathematical Statistics, 35(1), 73–101.
• James, W., & Stein, C. (1961). Estimation with quadratic loss. Proceedings of the Fourth Berkeley Symposium on Mathematical Statistics and Probability, 1, 361–379.
• Koenker, R., & Bassett, G. (1978). Regression quantiles. Econometrica, 46(1), 33–50.
• Lovell, M. C. (1963). Seasonal adjustment of economic time series and multiple regression analysis. Journal of the American Statistical Association, 58(304), 993–1010.
• MacKinnon, J. G., & White, H. (1985). Some heteroskedasticity-consistent covariance matrix estimators with improved finite sample properties. Journal of Econometrics, 29(3), 305–325.
• Pearl, J. (2009). Causality: Models, Reasoning, and Inference (2nd ed.). Cambridge University Press.
• Schwarz, G. (1978). Estimating the dimension of a model. The Annals of Statistics, 6(2), 461–464.
• Seber, G. A. F., & Lee, A. J. (2003). Linear Regression Analysis (2nd ed.). Wiley.
• Tibshirani, R. (1996). Regression shrinkage and selection via the lasso. Journal of the Royal Statistical Society: Series B, 58(1), 267–288.
• White, H. (1980). A heteroskedasticity-consistent covariance matrix estimator and a direct test for heteroskedasticity. Econometrica, 48(4), 817–838.