Hauptkomponentenregression in r (schritt für schritt)

Bei einem Satz von p Prädiktorvariablen und einer Antwortvariablen verwendet die multiple lineare Regression eine Methode, die als kleinste Quadrate bekannt ist, um die verbleibende Quadratsumme (RSS) zu minimieren:

RSS = Σ(y _i – ŷ _i ) ²

Gold:

• Σ : Ein griechisches Symbol für Summe
• y _i : der tatsächliche Antwortwert für die ^i-te Beobachtung
• ŷ _i : Der vorhergesagte Antwortwert basierend auf dem multiplen linearen Regressionsmodell

Wenn jedoch Prädiktorvariablen stark korreliert sind, kann Multikollinearität zu einem Problem werden. Dies kann dazu führen, dass Modellkoeffizientenschätzungen unzuverlässig werden und eine hohe Varianz aufweisen.

Eine Möglichkeit, dieses Problem zu vermeiden, besteht darin, die Hauptkomponentenregression zu verwenden, die M lineare Kombinationen (sogenannte „Hauptkomponenten“) der ursprünglichen p Prädiktoren findet und dann die Methode der kleinsten Quadrate verwendet, um ein lineares Regressionsmodell unter Verwendung der Hauptkomponenten als Prädiktoren anzupassen.

Dieses Tutorial bietet ein schrittweises Beispiel für die Durchführung einer Hauptkomponentenregression in R.

Schritt 1: Laden Sie die erforderlichen Pakete

Der einfachste Weg, eine Hauptkomponentenregression in R durchzuführen, ist die Verwendung von Funktionen im pls- Paket.

 #install pls package (if not already installed)
install.packages(" pls ")

load pls package
library(pls)

Schritt 2: Passen Sie das PCR-Modell an

Für dieses Beispiel verwenden wir den integrierten R-Datensatz namens mtcars , der Daten zu verschiedenen Fahrzeugtypen enthält:

 #view first six rows of mtcars dataset
head(mtcars)

                   mpg cyl disp hp drat wt qsec vs am gear carb
Mazda RX4 21.0 6 160 110 3.90 2.620 16.46 0 1 4 4
Mazda RX4 Wag 21.0 6 160 110 3.90 2.875 17.02 0 1 4 4
Datsun 710 22.8 4 108 93 3.85 2.320 18.61 1 1 4 1
Hornet 4 Drive 21.4 6 258 110 3.08 3.215 19.44 1 0 3 1
Hornet Sportabout 18.7 8 360 175 3.15 3.440 17.02 0 0 3 2
Valiant 18.1 6 225 105 2.76 3,460 20.22 1 0 3 1

Für dieses Beispiel passen wir ein Hauptkomponenten-Regressionsmodell (PCR) an, wobei wir hp als Antwortvariable und die folgenden Variablen als Prädiktorvariablen verwenden:

• mpg
• Anzeige
• Scheisse
• Gewicht
• qsec

Der folgende Code zeigt, wie das PCR-Modell an diese Daten angepasst wird. Beachten Sie die folgenden Argumente:

• scale=TRUE : Dies teilt R mit, dass jede der Prädiktorvariablen so skaliert werden soll, dass sie einen Mittelwert von 0 und eine Standardabweichung von 1 aufweist. Dadurch wird sichergestellt, dass keine Prädiktorvariable einen zu großen Einfluss auf das Modell hat, wenn sie in verschiedenen Einheiten gemessen wird. .
• validation=“CV“ : Dies weist R an , die k-fache Kreuzvalidierung zu verwenden, um die Modellleistung zu bewerten. Beachten Sie, dass hierbei standardmäßig k=10 Falten verwendet werden. Beachten Sie außerdem, dass Sie stattdessen „LOOCV“ angeben können, um eine Leave-One-Out-Kreuzvalidierung durchzuführen.

 #make this example reproducible
set.seed(1)

#fit PCR model
model <- pcr(hp~mpg+disp+drat+wt+qsec, data=mtcars, scale= TRUE , validation=" CV ")

Schritt 3: Wählen Sie die Anzahl der Hauptkomponenten

Nachdem wir das Modell angepasst haben, müssen wir bestimmen, wie viele Hauptkomponenten es wert sind, beibehalten zu werden.

Schauen Sie sich dazu einfach den quadratischen Mittelfehler des Tests (Test-RMSE) an, der durch die K-Kreuzvalidierung berechnet wurde:

 #view summary of model fitting
summary(model)

Data: 
	Y dimension: 32 1
Fit method: svdpc
Number of components considered: 5

VALIDATION: RMSEP
Cross-validated using 10 random segments.
       (Intercept) 1 comp 2 comps 3 comps 4 comps 5 comps
CV 69.66 44.56 35.64 35.83 36.23 36.67
adjCV 69.66 44.44 35.27 35.43 35.80 36.20

TRAINING: % variance explained
    1 comp 2 comps 3 comps 4 comps 5 comps
X 69.83 89.35 95.88 98.96 100.00
hp 62.38 81.31 81.96 81.98 82.03

Das Ergebnis enthält zwei interessante Tabellen:

1. VALIDIERUNG: RMSEP

Diese Tabelle zeigt uns den RMSE-Test, der durch k-fache Kreuzvalidierung berechnet wurde. Wir können Folgendes sehen:

• Wenn wir im Modell nur den ursprünglichen Term verwenden, beträgt der RMSE des Tests 69,66 .
• Wenn wir die erste Hauptkomponente hinzufügen, sinkt der RMSE-Test auf 44,56.
• Wenn wir die zweite Hauptkomponente hinzufügen, sinkt der RMSE-Test auf 35,64.

Wir können sehen, dass das Hinzufügen zusätzlicher Hauptkomponenten tatsächlich zu einer Erhöhung des RMSE des Tests führt. Daher scheint es optimal zu sein, im endgültigen Modell nur zwei Hauptkomponenten zu verwenden.

2. TRAINING: % der erklärten Varianz

Diese Tabelle zeigt uns den Prozentsatz der Varianz der Antwortvariablen, die durch die Hauptkomponenten erklärt wird. Wir können Folgendes sehen:

• Wenn wir nur die erste Hauptkomponente verwenden, können wir 69,83 % der Variation der Antwortvariablen erklären.
• Durch Hinzufügen der zweiten Hauptkomponente können wir 89,35 % der Variation der Antwortvariablen erklären.

Beachten Sie, dass wir immer noch in der Lage sein werden, eine größere Varianz zu erklären, indem wir mehr Hauptkomponenten verwenden, aber wir können sehen, dass das Hinzufügen von mehr als zwei Hauptkomponenten den Prozentsatz der erklärten Varianz tatsächlich nicht wesentlich erhöht.

Wir können den RMSE-Test (zusammen mit MSE und R-Quadrat-Test) auch als Funktion der Anzahl der Hauptkomponenten mithilfe der Funktion validationplot() visualisieren.

 #visualize cross-validation plots
validationplot(model)
validationplot(model, val.type="MSEP")
validationplot(model, val.type="R2") 

In jedem Diagramm können wir sehen, dass sich die Modellanpassung durch das Hinzufügen von zwei Hauptkomponenten verbessert, sich jedoch tendenziell verschlechtert, wenn wir weitere Hauptkomponenten hinzufügen.

Somit umfasst das optimale Modell nur die ersten beiden Hauptkomponenten.

Schritt 4: Verwenden Sie das endgültige Modell, um Vorhersagen zu treffen

Wir können das endgültige Zwei-Hauptkomponenten-PCR-Modell verwenden, um Vorhersagen über neue Beobachtungen zu treffen.

Der folgende Code zeigt, wie Sie den Originaldatensatz in einen Trainings- und einen Testsatz aufteilen und das PCR-Modell mit zwei Hauptkomponenten verwenden, um Vorhersagen für den Testsatz zu treffen.

 #define training and testing sets
train <- mtcars[1:25, c("hp", "mpg", "disp", "drat", "wt", "qsec")]
y_test <- mtcars[26: nrow (mtcars), c("hp")]
test <- mtcars[26: nrow (mtcars), c("mpg", "disp", "drat", "wt", "qsec")]
    
#use model to make predictions on a test set
model <- pcr(hp~mpg+disp+drat+wt+qsec, data=train, scale= TRUE , validation=" CV ")
pcr_pred <- predict(model, test, ncomp= 2 )

#calculate RMSE
sqrt ( mean ((pcr_pred - y_test)^2))

[1] 56.86549

Wir sehen, dass der RMSE des Tests 56,86549 beträgt. Dies ist die durchschnittliche Abweichung zwischen dem vorhergesagten HP- Wert und dem beobachteten HP- Wert für die Testsatzbeobachtungen.

Die vollständige Verwendung des R-Codes in diesem Beispiel finden Sie hier .