Regresi komponen utama dengan python (langkah demi langkah)

Mengingat sekumpulan variabel prediktor p dan variabel respons, regresi linier berganda menggunakan metode yang dikenal sebagai kuadrat terkecil untuk meminimalkan jumlah sisa kuadrat (RSS):

RSS = Σ(y _saya – ŷ _saya ) ²

Emas:

• Σ : Simbol Yunani yang berarti jumlah
• y _i : nilai respon sebenarnya untuk observasi ^ke-i
• ŷ _i : Nilai respons yang diprediksi berdasarkan model regresi linier berganda

Namun, ketika variabel prediktor berkorelasi tinggi, multikolinearitas bisa menjadi masalah. Hal ini dapat membuat estimasi koefisien model tidak dapat diandalkan dan menunjukkan varians yang tinggi.

Salah satu cara untuk menghindari masalah ini adalah dengan menggunakan regresi komponen utama , yang menemukan M kombinasi linier (disebut “komponen utama”) dari prediktor p asli dan kemudian menggunakan kuadrat terkecil agar sesuai dengan model regresi linier yang menggunakan komponen utama sebagai prediktor.

Tutorial ini memberikan contoh langkah demi langkah tentang cara melakukan regresi komponen utama dengan Python.

Langkah 1: Impor paket yang diperlukan

Pertama, kita akan mengimpor paket yang diperlukan untuk melakukan regresi komponen utama (PCR) dengan Python:

 import numpy as np
import pandas as pd
import matplotlib. pyplot as plt
from sklearn. preprocessing import scale 
from sklearn import model_selection
from sklearn. model_selection import RepeatedKFold
from sklearn.model_selection import train_test_split
from sklearn. PCA import decomposition
from sklearn. linear_model import LinearRegression
from sklearn. metrics import mean_squared_error

Langkah 2: Muat data

Untuk contoh ini, kita akan menggunakan kumpulan data bernama mtcars , yang berisi informasi tentang 33 mobil berbeda. Kita akan menggunakan hp sebagai variabel respon dan variabel berikut sebagai prediktor:

• mpg
• menampilkan
• kotoran
• berat
• qdetik

Kode berikut menunjukkan cara memuat dan menampilkan kumpulan data ini:

 #define URL where data is located
url = "https://raw.githubusercontent.com/Statorials/Python-Guides/main/mtcars.csv"

#read in data
data_full = pd. read_csv (url)

#select subset of data
data = data_full[["mpg", "disp", "drat", "wt", "qsec", "hp"]]

#view first six rows of data
data[0:6]


        mpg disp drat wt qsec hp
0 21.0 160.0 3.90 2.620 16.46 110
1 21.0 160.0 3.90 2.875 17.02 110
2 22.8 108.0 3.85 2.320 18.61 93
3 21.4 258.0 3.08 3.215 19.44 110
4 18.7 360.0 3.15 3.440 17.02 175
5 18.1 225.0 2.76 3.460 20.22 105

Langkah 3: Sesuaikan model PCR

Kode berikut menunjukkan cara menyesuaikan model PCR dengan data ini. Perhatikan hal berikut:

• pca.fit_transform(scale(X)) : Ini memberi tahu Python bahwa setiap variabel prediktor harus diskalakan agar memiliki rata-rata 0 dan deviasi standar 1. Hal ini memastikan bahwa tidak ada variabel prediktor yang memiliki pengaruh terlalu besar dalam model jika ini terjadi. diukur dalam satuan yang berbeda.
• cv = RepeatedKFold() : Ini memberitahu Python untuk menggunakan validasi silang k-fold untuk mengevaluasi kinerja model. Untuk contoh ini kita pilih k = 10 kali lipat, diulang sebanyak 3 kali.

 #define predictor and response variables
X = data[["mpg", "disp", "drat", "wt", "qsec"]]
y = data[["hp"]]

#scale predictor variables
pca = pca()
X_reduced = pca. fit_transform ( scale (X))

#define cross validation method
cv = RepeatedKFold(n_splits= 10 , n_repeats= 3 , random_state= 1 )

regr = LinearRegression()
mse = []

# Calculate MSE with only the intercept
score = -1*model_selection. cross_val_score (regr,
           n.p. ones ((len(X_reduced),1)), y, cv=cv,
           scoring=' neg_mean_squared_error '). mean ()    
mse. append (score)

# Calculate MSE using cross-validation, adding one component at a time
for i in np. arange (1, 6):
    score = -1*model_selection. cross_val_score (regr,
               X_reduced[:,:i], y, cv=cv, scoring=' neg_mean_squared_error '). mean ()
    mse. append (score)
    
# Plot cross-validation results    
plt. plot (mse)
plt. xlabel ('Number of Principal Components')
plt. ylabel ('MSE')
plt. title ('hp') 

Plot menampilkan jumlah komponen utama sepanjang sumbu x dan pengujian MSE (mean square error) sepanjang sumbu y.

Dari grafik terlihat bahwa MSE pengujian menurun dengan penambahan dua komponen utama, namun mulai meningkat jika kita menambahkan lebih dari dua komponen utama.

Dengan demikian, model optimal hanya mencakup dua komponen utama pertama.

Kita juga dapat menggunakan kode berikut untuk menghitung persentase varians dalam variabel respon yang dijelaskan dengan menambahkan setiap komponen utama ke model:

 n.p. cumsum (np. round (pca. explained_variance_ratio_ , decimals= 4 )* 100 )

array([69.83, 89.35, 95.88, 98.95, 99.99])

Kita dapat melihat hal berikut:

• Dengan hanya menggunakan komponen utama pertama, kita dapat menjelaskan 69,83% variasi variabel respon.
• Dengan menambahkan komponen utama kedua, kita dapat menjelaskan 89,35% variasi variabel respon.

Perhatikan bahwa kita masih dapat menjelaskan lebih banyak varians dengan menggunakan lebih banyak komponen utama, namun kita dapat melihat bahwa menambahkan lebih dari dua komponen utama tidak benar-benar meningkatkan persentase varians yang dijelaskan sebanyak itu.

Langkah 4: Gunakan model akhir untuk membuat prediksi

Kita dapat menggunakan model PCR dua komponen utama terakhir untuk membuat prediksi tentang observasi baru.

Kode berikut menunjukkan cara membagi kumpulan data asli menjadi kumpulan pelatihan dan pengujian serta menggunakan model PCR dengan dua komponen utama untuk membuat prediksi pada kumpulan pengujian.

 #split the dataset into training (70%) and testing (30%) sets
X_train,X_test,y_train,y_test = train_test_split (X,y,test_size= 0.3 , random_state= 0 ) 

#scale the training and testing data
X_reduced_train = pca. fit_transform ( scale (X_train))
X_reduced_test = pca. transform ( scale (X_test))[:,:1]

#train PCR model on training data 
regr = LinearRegression()
reg. fit (X_reduced_train[:,:1], y_train)

#calculate RMSE
pred = regr. predict (X_reduced_test)
n.p. sqrt ( mean_squared_error (y_test, pred))

40.2096

Kita melihat bahwa tes RMSE ternyata 40.2096 . Ini adalah deviasi rata-rata antara nilai hp yang diprediksi dan nilai hp yang diamati untuk observasi set pengujian.

Kode Python lengkap yang digunakan dalam contoh ini dapat ditemukan di sini .