38 Simple lineare Regression
Letzte Änderung am 02. April 2024 um 09:52:43
“Farewell to the Fairground; These rides aren’t working anymore.” — White Lies, Farewell to the Fairground
In diesem Kapitel wollen wir uns mit den Grundlagen der simplen linearen Regression beschäftigen. Damit meine ich erstmal die Idee eine Gerade durch eine Punktewolke zu zeichnen. Wir haben auf der \(x\)-Achse kontinuierliche Werte wie auch auf der \(y\)-Achse vorliegen. Das ist erstmal die simpelste Anwendung. Wir lernen hier einmal die Grundbegriffe und erweitern diese dann auf komplexere Modelle. Das heißt, wir haben ein Outcome \(y\), was eine normalverteilte und kontinuierliche Variable ist, sowie eine Einflussvariable \(x_1\), die eine kontinuierliche Variable ist. Wir wollen jetzt herausfinden, welchen Einfluss oder Effekt die Einflussvariable \(x_1\) auf das normalverteilte Outcome \(y\) hat. Wir schreiben hier \(x_1\), da wir damit ausdrücken, dass wir später auch noch mehr Einflussvariablen \(x\) in der multiplen linearen Regression nutzen werden.
Wenn du jetzt denkst ‘Hä?’, dann besuche doch einmal die fantastische Seite Explained Visually | Ordinary Least Squares Regression um selber zu erfahren, was eine lineare Regression macht. Auf der Seite findest du interaktive Abbildungen, die dir das Konzept der linearen Regression sehr anschaulich nachvollziehen lassen.
Sehr simpel gesprochen legen wir also eine Gerade durch eine Punktewolke. Du kannst dich aber davon gedanklich lösen, dass die lineare Regression nur eine Methode ist um eine Gerade durch eine Punktewolke zu legen. Die lineare Regression und damit auch das statistische Modellieren kann viel mehr. Wir fassen also zusammen. Ein simples lineares Modell hat ein kontinuierliches \(x\) als Einflussvariable. Später können wir dann auch andere \(x\) nehmen, wie zum Beispiel einen Faktor \(f\) als kategoriale Einflussvariable.
\[ y \sim x_1 \]
Damit beginnen wir hier einmal mit den Konzept der linearen Regression in dem wir wirklich sehr einfach uns eine Punktewolke vorstellen und durch diese Punktewolke die optimale Gerade zeichen wollen. Dann wollen wir wissen wie gut das geklappt hat. Liegen die Punkte weit von der Geraden entfernt oder eher näher dran? In der Abbildung 38.1 siehst du einmal eine Punktewolke die sich aus den Werten der \(x\)-Achse und den Werten auf der \(y\)-Achse ergibt. Wir wollen jetzt mathematisch verstehen, wie wir die rote Gerade entsteht. Das ist die Idee der simplen linearen Regression.
Aber in diesem Kapitel konzentrieren wir uns erstmal darauf die Begriffe der simplen linearen Regression und deren Durchführung in R zu verstehen. Dann wollen wir noch die Ausgabe einer linearen Regression in R sowie den Unterschied zwischen einem kausalen und einen prädiktiven Modell ergründen.
Weitere Tutorien für die simple lineare Regression
Wir immer geht natürlich mehr als ich hier Vorstellen kann. Du findest im Folgenden Tutorien, die mich hier in dem Kapitel inspiriert haben. Ich habe mich ja in diesem Kapitel auf die Durchführbarkeit in R und die allgemeine Verständlichkeit konzentriert. Es geht aber natürlich wie immer auch mathematischer…
- Wir funktioniert nun so eine lineare Regression und was sind den jetzt eigentlich die Koeffizienten \(\beta_0\) und \(\beta_1\) eigentlich? Hier gibt es die fantastische Seite Explained Visually | Ordinary Least Squares Regression, die dir nochmal erlaubt selbe mit Punkten in einem Scatterplot zu spielen und zu sehen wie sich dann die Regressionsgleichung ändert.
- Du kannst auf der Seite Manual linear regression analysis using R nochmal weiter über die lineare Regression lesen. Der Blogpost ist sehr umfangreich und erklärt nochmal schrittweise, wie die lineare Regression in R per hand funktioniert.
- Simple Regression auf Wikipedia mit weit mehr Informationen zu den Formeln und den Zusammenhängen. Ein toller zusammenfassender Artikel.
38.1 Genutzte R Pakete
Wir wollen folgende R Pakete in diesem Kapitel nutzen.
pacman::p_load(tidyverse, magrittr, broom,
readxl, conflicted)
conflicts_prefer(magrittr::set_names)An der Seite des Kapitels findest du den Link Quellcode anzeigen, über den du Zugang zum gesamten R-Code dieses Kapitels erhältst.
38.2 Daten
Wir wollen uns erstmal mit einem einfachen Datenbeispiel beschäftigen. Wir können die lineare Regression auf sehr großen Datensätzen anwenden, wie auch auf sehr kleinen Datensätzen. Prinzipiell ist das Vorgehen gleich. Wir nutzen jetzt aber erstmal einen kleinen Datensatz mit \(n=7\) Beobachtungen. In der Tabelle 38.1 ist der Datensatz simplel_tbl dargestellt. Wir wollen den Zusammenhang zwischen der Sprungweite in [cm] und dem Gewicht in [mg] für sieben Beobachtungen modellieren.
simple_tbl <- tibble(jump_length = c(1.2, 1.8, 1.3, 1.7, 2.6, 1.8, 2.7),
weight = c(0.8, 1, 1.2, 1.9, 2, 2.7, 2.8))jump_length und der normalverteilten Variable weight.
| jump_length | weight |
|---|---|
| 1.2 | 0.8 |
| 1.8 | 1.0 |
| 1.3 | 1.2 |
| 1.7 | 1.9 |
| 2.6 | 2.0 |
| 1.8 | 2.7 |
| 2.7 | 2.8 |
In Abbildung 38.2 sehen wir die Visualisierung der Daten simple_tbl in einem Scatterplot mit einer geschätzen Gerade.
Wir schauen uns in diesem Kapitel nur eine simple lineare Regression mit einem \(x_1\) an. In unserem Fall ist das \(x_1\) gleich dem weight. Später schauen wir dann multiple lineare Regressionen mit mehreren \(x_1,..., x_p\) an.
Bevor wir mit dem Modellieren anfangen können, müssen wir verstehen, wie ein simples Modell theoretisch aufgebaut ist. Danach können wir uns das lineare Modell in R anschauen.
38.3 Simple lineare Regression theoretisch
Wir haben nun die ersten sieben Beobachtungen in dem Objekt simple_tbl vorliegen. Wie sieht nun theoretisch eine lineare Regression aus? Wir wollen eine Grade durch Punkte legen, wie wie wir es in Abbildung 38.2 sehen. Die blaue Gerade wir durch eine Geradengleichung beschreiben. Du kenst vermutlich noch die Form \(y = mx + b\). In der Statistik beschreiben wir eine solche Gerade aber wie folgt.
\[ y \sim \beta_0 + \beta_1 x_1 + \epsilon \]
mit
- \(\beta_0\) als den y-Achsenabschnitt.
- \(\beta_1\) als der Steigung der Geraden.
- \(\epsilon\) als Residuen oder die Abweichungen von den \(y\)-Werten auf Geraden zu den einzelnen \(y\)-Werten der Beobachtungen.
In Tabelle 38.2 siehst du nochmal in einer Tabelle den Vergleich von der Schreibweise der linearen Regression in der Schule und in der Statistik. Darüber hinaus sind die deutschen Begriffe den englischen Begriffen gegenüber gestellt. Warum schreiben wir die Gleichung in der Form? Damit wir später noch weitere \(\beta_px_p\)-Paare ergänzen könen und so multiple Modelle bauen können.
| \(\boldsymbol{y = mx +b}\) | \(\boldsymbol{y \sim \beta_0 + \beta_1 x_1 + \epsilon}\) | Deutsch | Englisch |
|---|---|---|---|
| \(m\) | \(\beta_1\) | Steigung | Slope |
| \(x\) | \(x_1\) | Einflussvariable | Risk factor |
| \(b\) | \(\beta_0\) | y-Achsenabschnitt | Intercept |
| \(\epsilon\) | Residuen | Residual |
In Abbildung 38.3 sehen wir die Visualisierung der Gleichung in einer Abbildung. Die Gerade läuft durch die Punktewolke und wird durch die statistischen Maßzahlen bzw. Parameter \(\beta_0\), \(\beta_1\) sowie den \(\epsilon\) beschrieben. Wir sehen, dass das \(\beta_0\) den Intercept darstellt und das \(\beta_1\) die Steigung der Geraden. Wenn wir \(x\) um 1 Einheit erhöhen \(x+1\), dann steigt der \(y\) Wert um den Wert von \(\beta_1\). Die einzelnen Abweichungen der beobachteten \(y\)-Wert zu den \(y\)-Werten auf der Gerade (\(\hat{y}\)) werden als Residuen oder auch \(\epsilon\) bezeichnet.
Schauen wir uns einmal den Zusammenhang von \(y\), den beobachteten Werten, und \(\hat{y}\), den geschätzen Werten auf der Gerade in unserem Beispiel an. In Tabelle 38.3 sehen wir die Berechnung der einzelnen Residuen für die Gerade aus der Abbildung 38.2. Wir nehmen jedes beobachtete \(y\) und ziehen den Wert von \(y\) auf der Gerade, bezeichnet als \(\hat{y}\), ab. Diesen Schritt machen wir für jedes Wertepaar \((y_i; \hat{y}_i)\). In R werden die \(\hat{y}\) auch fitted values genannt. Die \(\epsilon\) Werte werden dann residuals bezeichnet.
| x | y | \(\boldsymbol{\hat{y}}\) | Residuen (\(\boldsymbol{\epsilon}\)) | Wert |
|---|---|---|---|---|
| 0.8 | 1.2 | 1.38 | \(\epsilon_1 = y_1 - \hat{y}_1\) | \(\epsilon_1 = 1.2 - 1.38 = -0.18\) |
| 1.0 | 1.8 | 1.48 | \(\epsilon_2 = y_2 - \hat{y}_2\) | \(\epsilon_2 = 1.8 - 1.48 = +0.32\) |
| 1.2 | 1.3 | 1.58 | \(\epsilon_3 = y_3 - \hat{y}_3\) | \(\epsilon_3 = 1.3 - 1.58 = -0.28\) |
| 1.9 | 1.7 | 1.94 | \(\epsilon_4 = y_4 - \hat{y}_4\) | \(\epsilon_4 = 1.7 - 1.94 = -0.24\) |
| 2.0 | 2.6 | 1.99 | \(\epsilon_5 = y_5 - \hat{y}_5\) | \(\epsilon_5 = 2.6 - 1.99 = +0.61\) |
| 2.7 | 1.8 | 2.34 | \(\epsilon_6 = y_6 - \hat{y}_6\) | \(\epsilon_6 = 1.8 - 2.34 = -0.54\) |
| 2.8 | 2.7 | 2.40 | \(\epsilon_7 = y_7 - \hat{y}_7\) | \(\epsilon_7 = 2.7 - 2.40 = +0.30\) |
Die Abweichungen \(\epsilon\) oder auch Residuen genannt haben einen Mittelwert von \(\bar{\epsilon} = 0\) und eine Varianz von \(s^2_{\epsilon} = 0.17\). Wir schreiben, dass die Residuen normalverteilt sind mit \(\epsilon \sim \mathcal{N}(0, s^2_{\epsilon})\). Wir zeichnen die Gerade also so durch die Punktewolke, dass die Abstände zu den Punkten, die Residuen, im Mittel null sind. Die Optimierung erreichen wir in dem wir die Varianz der Residuuen minimieren. Folglich modellieren wir die Varianz.
Simple lineare Regression händisch
Gut, das war jetzt die theoretische Abhandlung ohne eine mathematische Formel. Es geht natürlich auch mit den nackten Zahlen. In der Tabelle 38.4 siehst du einmal sieben Beobachtungen mit dem Körpergewicht als \(y\) sowie der Körpergröße als \(x\). Wir wollen jetzt einmal die Regressionsgleichung bestimmen. Wir sehen also unsere Werte für \(\beta_0\) und \(\beta_1\) aus?
| height | weight |
|---|---|
| 167 | 70 |
| 188 | 83 |
| 176 | 81 |
| 186 | 90 |
| 192 | 94 |
| 205 | 100 |
| 198 | 106 |
Damit wir einmal wissen, was wir als Lösung erhalten würden, hier einmal die lineare Regression mit der Funktion \(lm()\) und die entsprechenden Werte für den (Intercept) und der Steigung.
lm(weight ~ height, data = by_hand_tbl) |>
coef()(Intercept) height
-75.390377 0.877845 Wir suchen dann damit die folgende Regressionsgleichung mit der Körpergröße als \(x\) und dem zugehörigen Körpergewicht als \(y\).
\[ weight = \beta_0 + \beta_1 \cdot height \]
Da es dann immer etwas schwer ist, sich den Zusammenhang zwischen Körpergewicht und Körpergröße vorzustellen, habe ich nochmal in der Abbildung 38.4 den Scatterplot der Daten erstellt. Die rote Gerade stellt die Regressiongleichung dar. Wir erhalten ein y-Achsenabschnitt mit \(\beta_0\) von \(-75.39\) sowie eine Steigung mit \(\beta_1\) von \(0.88\) aus unseren Daten.
ggplot(by_hand_tbl, aes(height, weight)) +
theme_minimal() +
geom_point() +
geom_function(fun = \(x) -75.39 + 0.88 * x, color = "red")
Jetzt stellt sich die Frage, wie wir händisch die Werte für den y-Achsenabschnitt mit \(\beta_0\) sowie der Steigung mit \(\beta_1\) berechnen. Dafür gibt es jeweils eine Formel. Hier müssen wir dann sehr viele Summen berechnen, was ich dann gleich einmal in einer Tabelle zusammenfasse.
- Formel für y-Achsenabschnitt mit \(\beta_0\)
-
\[ \beta_0 = \cfrac{(\Sigma Y)(\Sigma X^2) - (\Sigma X)(\Sigma XY)}{n(\Sigma X^2) - (\Sigma X)^2} \]
- Formel für Steigung mit \(\beta_1\)
-
\[ \beta_1 = \cfrac{n(\Sigma XY) - (\Sigma X)(\Sigma Y)}{n(\Sigma X^2) - (\Sigma X)^2} \]
In der Tabelle 38.5 siehst du nochmal die originalen Datenpunkte und dann die entsprechenden Werte für das Produkt von weight und height mit \(XY\) und dann die jeweiligen Quadrate der beiden mit \(X^2\) und \(Y^2\). Wir brauchen dann aber nicht diese Werte sondern die Summen der Werte. Das Summieren lagere ich dann nochmal in eine weitere Tabelle aus.
| height | weight | \(XY\) | \(X^2\) | \(Y^2\) |
|---|---|---|---|---|
| 167 | 70 | 11690 | 27889 | 4900 |
| 188 | 83 | 15604 | 35344 | 6889 |
| 176 | 81 | 14256 | 30976 | 6561 |
| 186 | 90 | 16740 | 34596 | 8100 |
| 192 | 94 | 18048 | 36864 | 8836 |
| 205 | 100 | 20500 | 42025 | 10000 |
| 198 | 106 | 20988 | 39204 | 11236 |
In der abschließenden Tabelle findest du dann einmal die Summen der beobachteten Werte \(X\) und \(Y\) sowie des Produkts von \(X\) und \(Y\) sowie deren Quadrate mit \(X^2\) und \(Y^2\). Damit haben wir dann alles zusammen um die Formel oben zu füllen.
| height \((\Sigma X)\) | weight \((\Sigma Y)\) | \(\Sigma XY\) | \(\Sigma X^2\) | \(\Sigma Y^2\) |
|---|---|---|---|---|
| 1312 | 624 | 117826 | 246898 | 56522 |
Ich habe dann die ganzen Summen einmal händisch berechnet und dann in den Formeln von oben eingesetzt. Wir erhalten dann für den y-Achsenabschnitt \(\beta_0\) folgenden Wert.
\[ \beta_0 = \cfrac{624 \cdot 246898 - 1312 \cdot 117826}{7\cdot 246898 - 1312^2} = -75.39038 \]
Die ganze Berechnung habe ich dann auch einmal für die Steigung \(\beta_1\) ebenfalls einmal durchgeführt.
\[ b_1 = \cfrac{7\cdot 117826 - 1312\cdot624}{7\cdot246898 - 1312^2} = 0.877845 \]
Wir sehen, es kommen die gleichen Werte für den y-Achsenabschnitt \(\beta_0\) und die Steigung \(\beta_1\) raus. Das hat ja schonmal sehr gut geklappt. Eine andere Art die gleiche Werte effizienter zu berechnen ist die Matrixberechnung der Koeffizienten der linearen Regression. Wir könnten dann auch komplexere Modelle mit mehr als nur einem \(x\) und einem \(\beta_1\) berechnen. Die grundlegende Formel siehst du einmal im Folgenden dargestellt.
\[ \begin{pmatrix} \beta_0 \\ \beta_1 \end{pmatrix} = \mathbf{(X^T X)^{−1}(X^T Y)} \]
Wir brauchen jetzt einiges an Matrixrechnung um die jeweiligen Formelteile zu berechnen. Ich habe dir in den folgenden Tabs einmal Schritt für Schritt die einzelnen Teile berechnet. Wir immer machen wir das eigentlich nicht so richtig per Hand, sondern nutzen einen Computer. Prinzipiell wäre eine händische Lösung natürlich möglich.
X <- as.matrix(c(167, 188, 176, 186, 192, 205, 198))
X <- cbind(rep(1, 7), X)
X [,1] [,2]
[1,] 1 167
[2,] 1 188
[3,] 1 176
[4,] 1 186
[5,] 1 192
[6,] 1 205
[7,] 1 198Y <- as.matrix(c(70, 83, 81, 90, 94, 100, 106))
Y [,1]
[1,] 70
[2,] 83
[3,] 81
[4,] 90
[5,] 94
[6,] 100
[7,] 106Xt <- t(X)
Xt [,1] [,2] [,3] [,4] [,5] [,6] [,7]
[1,] 1 1 1 1 1 1 1
[2,] 167 188 176 186 192 205 198XtX <- Xt %*% X
XtX [,1] [,2]
[1,] 7 1312
[2,] 1312 246898XtY <- Xt %*% Y
XtY [,1]
[1,] 624
[2,] 117826XtXinv <- solve(XtX)
XtXinv [,1] [,2]
[1,] 35.5658312 -0.188994526
[2,] -0.1889945 0.001008355Am Ende müssen wir dann alle Teile in der Form \(\mathbf{(X^T X)^{−1}(X^T Y)}\) einmal zusammenbringen. Das siehst dann in R wie folgt aus. Wir erhalten dann eine Matrix wieder wobei die erste Zeile der y-Achsenabschnitt \(\beta_0\) und die zweite Zeile die Steigung \(\beta_1\) ist. Wir erhalten fast die gleichen Werte wie auch schon oben.
XtXinv %*% Xt %*% Y [,1]
[1,] -75.390377
[2,] 0.877845Wenn du dich tiefer in die Thematik einlesen willst, dann sind hier weitere Quellen zu der Thematik unter den folgenden Links und Tutorien.
38.4 Simples lineare Regression in R
Im Allgemeinen können wir ein Modell in R wie folgt schreiben. Wir brauchen ein y auf der linken Seite und in der simplen linearen Regressione ein \(x\) auf der rechten Seite der Gleichung. Wir brauchen also zwei Variablen \(y\) und \(x\), die natürlich nicht im Datensatz in R so heißen müssen. Im Folgenen dann einmal die Modellschreibweise für \(y\) hängt ab von \(x\). Das \(y\) repräsentiert eine Spalte im Datensatz und das \(x\) repräsentiert ebenso eine Spalte im Datensatz.
\[ \Large y \sim x \]
Konkret würden wir in unserem Beispiel das Modell wie folgt benennen. Das \(y\) wird zu jump_length und das \(x\) wird zu weight. Wir haben dann das Modell in der simpelsten Form definiert. Im Folgenden siehst du dann eimal die Modellschreibweise in R.
\[ \Large\overbrace{\mbox{jump\_length}}^{\Large y} \sim \overbrace{\mbox{weight}}^{\Large x} \]
Nachdem wir das Modell definiert haben, setzen wir dieses Modell jump_length ~ weight in die Funktion lm() ein um das lineare Modell zu rechnen. Wie immer müssen wir auch festlegen aus welcher Datei die Spalten genommen werden sollen. Das machen wir mit der Option data = simple_tbl. Wir speichern dann die Ausgabe der Funktion lm() in dem Objekt fit_1 damit wir die Ausgabe noch in andere Funktionen pipen können.
fit_1 <- lm(jump_length ~ weight, data = simple_tbl)Wir können jetzt mir dem Modell mehr oder minder drei Dinge tun. Abhängig von der Fragestellung liefert uns natürlich jedes der drei Möglichkeiten eine andere Antwort. Wir können auch mehrere dieser Fragen gleichzeitig beantworten, was aber meistens zu einer konfusen Analyse führt. Im Folgenden einmal eine Flowchart als grobe Übersicht deiner Möglichkeiten nach einem beispielhaften lm()-Modell als simple Regression.
flowchart TD
A("Lineares Modell
zum Beispiel mit lm()"):::factor --> B1
A("Lineares Modell
zum Beispiel mit lm()"):::factor --> B2
A("Lineares Modell zum Beispiel mit lm()"):::factor --> B3
subgraph B1["Inferenzmodell"]
B[/"anova()"\]:::factor --> E[\"emmeans()"/]:::factor
end
subgraph B2["Kausales Modell"]
C("summary()"):::factor --> F("augment()"):::eval
C("summary()"):::factor --> G("glance()"):::eval
end
subgraph B3["Prädiktives Modell"]
D("predict()"):::factor --> H("conf_mat()"):::eval
end
classDef factor fill:#56B4E9,stroke:#333,stroke-width:0.75px
classDef eval fill:#E69F00,stroke:#333,stroke-width:0.75px
lm() als eine simple lineare Regression. Es gibt natürlich noch andere Modellierungen und damit Funktion. Der generelle Ablauf bleibt jedoch gleich. Die orangen Kacheln stellen optionale Funktionen zur Güte der Regression dar.
Hier nochmal die Möglichkeiten aus der Abbildung 38.5 mit mehr Details zu der obigen Flowchart. Für mehr Informationen musst du dann die entsprechenden Kapitel besuchen, hier passt dann alles wirklich nicht hin für ein so komplexen Thema wie die Modellierung von Daten.
- Wir rechnen mit dem Fit des Modells einen Gruppenvergleich oder eine Inferenzstatistik (siehe Kapitel 29 zu der ANOVA oder Kapitel 34 zu multiplen Vergleichen) – wir wollen also einen Vergleich zwischen Gruppen rechnen. Wir haben damit dann kein kontinuierliches \(x\) vorliegen sondern einen Faktor mit verschiedenen Leveln als Gruppen.
- Wir rechnen ein kausales Modell, uns interessieren die Effekte (siehe Kapitel 38.4.2) – wir sind an den Effekten von verschiedenen Einflussvariablen interessiert. Dabei geht es dann weniger um eine Behandlungsgruppe sondern um viele verschiedene Einflussvariablen, die verstanden werden wollen in ihrem Einfluss auf das Outcome \(y\). Wir können uns dann die Ergebnisse durch die Funktionen
augment()undglance()aus dem R Paket{broom}näher anschauen. - Wir rechnen ein prädiktives Modell, uns interessiert der Wert neuer Werte (siehe Kapitel 38.4.3) – wir wollen dann also eine Vorhersage rechnen. Wir bauen uns ein komplexes Modell und wollen mit diesem Modell zukünftige Werte vom Outcome \(y\) vorhersagen. Wenn wir also neue Einflussvariablen \(x\) messen, können wir dann anhand des Modells unbekannte Outcomes \(y\) bestimmen. Wir können uns dann mehr Informationen zu der Güte der Vorhersage mit der Funktion
conf_mat()aus dem R Paket{tidymodels}anschauen.
38.4.1 Inferenzmodell
In diesem Kapitel betrachten wir nicht die Analyse von Gruppenvergleichen oder der Inferenzstatistik. Das machen wir dann in eigenen Kapiteln. Bitte besuche das Kapitel 29 zu mehr Informationen zu der einfaktoriellen und zweifaktoriellen ANOVA oder Kapitel 34 zu multiplen Vergleichen. Es würde dieses Kapitel sprengen, wenn wir uns dann hier die Sachlage nochmal hier aufrollen würde.
38.4.2 Kausales Modell
Im Folgenden rechnen wir ein kausales Modell, da wir an dem Effekt des \(x\) interessiert sind. Wenn also das \(x_1\) um eine Einheit ansteigt, um wie viel verändert sich dann das \(y\)? Der Schätzer \(\beta_1\) gibt uns also den Einfluss oder den kausalen Zusammenhang zwischen \(y\) und \(x_1\) wieder. Im ersten Schritt schauen wir uns die Ausgabe der Funktion lm() in der Funktion summary() an. Daher pipen wir das Objekt fit_1 in die Funktion summary().
fit_1 |> summaryWir erhalten folgende Ausgabe dargestellt in Abbildung 38.6.
summary() Ausgabe des Modells fit_1. Hier nicht erschrecken, wir kriegen hier sehr viele Informationen aufeinmal, die wir teilweise nicht alle benötigen. Es hängt dann eben sehr stark von der Fragestelung ab.
Was sehen wir in der Ausgabe der summary() Funktion? Als erstes werden uns die Residuen wiedergegeben. Wenn wir nur wenige Beobachtungen haben, dann werden uns die Residuen direkt wiedergegeben, sonst die Verteilung der Residuen. Mit der Funktion augment() aus dem R Paket {broom} können wir uns die Residuen wiedergeben lassen. Die Residuen schauen wir uns aber nochmal im Kapitel 40 genauer an.
fit_1 |> augment()# A tibble: 7 x 8
jump_length weight .fitted .resid .hat .sigma .cooksd .std.resid
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 1.2 0.8 1.38 -0.176 0.388 0.496 0.0778 -0.496
2 1.8 1 1.48 0.322 0.297 0.471 0.151 0.844
3 1.3 1.2 1.58 -0.280 0.228 0.483 0.0725 -0.701
4 1.7 1.9 1.94 -0.237 0.147 0.492 0.0275 -0.564
5 2.6 2 1.99 0.612 0.156 0.384 0.199 1.47
6 1.8 2.7 2.34 -0.545 0.367 0.376 0.656 -1.51
7 2.7 2.8 2.40 0.304 0.417 0.467 0.276 0.877Im zweiten Block erhalten wir die Koeffizienten (eng. coefficients) der linearen Regression. Das heißt, wir kriegen dort \(\beta_0\) als y-Achsenabschnitt sowie die Steigung \(\beta_1\) für das Gewicht. Dabei ist wichtig zu wissen, dass immer als erstes der y-Achsenabschnitt (Intercept) auftaucht. Dann die Steigungen der einzelnen \(x\) in dem Modell. Wir haben nur ein kontinuierliches \(x\), daher ist die Interpretation der Ausgabe einfach. Wir können die Gradengleichung wie folgt formulieren.
\[ jump\_length \sim 0.97 + 0.51 \cdot weight \]
Was heißt die Gleichung nun? Wenn wir das \(x\) um eine Einheit erhöhen dann verändert sich das \(y\) um den Wert von \(\beta_1\). Wir haben hier eine Steigung von \(0.51\) vorliegen. Ohne Einheit keine Interpretation! Wir wissen, dass das Gewicht in [mg] gemessen wurde und die Sprungweite in [cm]. Damit können wir aussagen, dass wenn ein Floh 1 mg mehr wiegt der Floh 0.51 cm weiter springen würde.
Schauen wir nochmal in die saubere Ausgabe der tidy() Funktion. Wir sehen nämlich noch einen \(p\)-Wert für den Intercept und die Steigung von weight.
fit_1 |> tidy()# A tibble: 2 x 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 0.969 0.445 2.18 0.0813
2 weight 0.510 0.232 2.20 0.0790Wenn wir einen \(p\)-Wert sehen, dann brauchen wir eine Nullhypothese, die wir dann eventuell mit der Entscheidung am Signifikanzniveau \(\alpha\) von 5% ablehnen können. Die Nullhypothese ist die Gleichheitshypothese. Wenn es also keinen Effekt von dem Gewicht auf die Sprungweite gebe, wie groß wäre dann \(\beta_1\)? Wir hätten dann keine Steigung und die Grade würde parallel zur x-Achse laufen. Das \(\beta_1\) wäre dann gleich null.
\[ \begin{align*} H_0: \beta_i &= 0\\ H_A: \beta_i &\neq 0 \\ \end{align*} \]
Wir haben für jedes \(\beta_i\) ein eigenes Hypothesenpaar. Meistens interessiert uns der Intercept nicht. Ob der Intercept nun durch die Null geht oder nicht ist eher von geringem Interessen.
Spannder ist aber wie sich der \(p\)-Wert berechnet. Der \(p\)-Wert basiert auf einer t-Statistik, also auf dem t-Test. Wir rechnen für jeden Koeffizienten \(\beta_i\) einen t-Test. Das machen wir in dem wir den Koeffizienten estimate durch den Fehler des Koeffizienten std.error teilen.
\[ \begin{align*} T_{(Intercept)} &= \cfrac{\mbox{estimate}}{\mbox{std.error}} = \cfrac{0.969}{0.445} = 2.18\\ T_{weight} &= \cfrac{\mbox{estimate}}{\mbox{std.error}} = \cfrac{0.510}{0.232} = 2.20\\ \end{align*} \]
Wir sehen in diesem Fall, dass weder der Intercept noch die Steigung von weight signifikant ist, da die \(p\)-Werte mit \(0.081\) und \(0.079\) leicht über dem Signifikanzniveau von \(\alpha\) gleich 5% liegen. Wir haben aber einen starkes Indiz gegen die Nullhypothese, da die Wahrscheinlichkeit die Daten zu beobachten sehr gering ist unter der Annahme das die Nullhypothese gilt.
Zun Abschluß noch die Funktion glance() ebenfalls aus dem R Paket {broom}, die uns erlaubt noch die Qualitätsmaße der linearen Regression zu erhalten. Wir müssen nämlich noch schauen, ob die Regression auch funktioniert hat. Die Überprüfung geht mit einem \(x\) sehr einfach. Wir können uns die Grade ja anschauen. Das geht dann mit einem Model mit mehreren \(x\) nicht mehr und wir brauchen andere statistische Maßzahlen.
fit_1 |> glance() # A tibble: 1 x 12
r.squared adj.r.squared sigma statistic p.value df logLik AIC BIC
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 0.492 0.391 0.455 4.84 0.0790 1 -3.24 12.5 12.3
# i 3 more variables: deviance <dbl>, df.residual <int>, nobs <int>38.4.3 Prädiktives Modell
Neben dem kausalen Modell gibt es auch die Möglichkeit ein prädiktives Modell zu rechnen. Im Prinzip ist die Sprache hier etwas ungenau. Wir verwenden das gefittete Modell nur anders. Anstatt das Modell fit_1 in die Funktion summary() zu pipen, pipen wir die das Modell in die Funktion predict(). Die Funktion predict() kann dann für neue Daten über die Option newdata = das \(y\) vorhersagen.
In unserem Fall müssen wir uns deshalb ein tibble mit einer Spalte bauen. Wir haben ja oben im Modell auch nur ein \(x_1\) mit aufgenommen. Später können wir natürlich auch für multiple Modelle die Vorhersage machen. Wichtig ist, dass die Namen gleich sind. Das heißt in dem neuen Datensatz müssen die Spalten exakt so heißen wir in dem alten Datensatz in dem das Modell gefittet wurde.
simple_new_tbl <- tibble(weight = c(1.7, 1.4, 2.1, 3.0))
predict(fit_1, newdata = simple_new_tbl) |> round(2) 1 2 3 4
1.84 1.68 2.04 2.50 Wie wir sehen ist die Anwendung recht einfach. Wir haben die vier jump_length Werte vorhergesagt bekommen, die sich mit dem Fit des Modells mit den neuen weight Werten ergeben. In Abbildung 38.7 sehen wir die Visualisierung der vier vorhergesagten Werte. Die Werte müssen auf der Geraden liegen.
Wir werden später in der Klassifikation, der Vorhersage von \(0/1\)-Werten, sowie in der multiplen Regression noch andere Prädktionen und deren Maßzahlen kennen lernen. Im Rahmen der simplen Regression soll dies aber erstmal hier genügen.