52 Modellieren in R
Letzte Änderung am 17. July 2025 um 21:12:07
“Ich weiß nicht weiter; Ich will mich verändern; Doch wie fang ich’s an?” — Tocotronic, Die Unendlichkeit
Dieses Kapitel wird in den nächsten Wochen geschrieben und ist damit meine aktuelle Großbaustelle. Ich plane zum Beginn des WiSe 2025/26 eine fertige Version des Kapitels erstellt zu haben. Während das Kapitel entsteht, funktioniert so manches dann nicht so wie es soll. Bitte daher hier dann abwarten.
Dieses Startkapitel gibt dir nochmal eine Übersicht über das statistische Modellieren in R. Hier liegt vor allem der Fokus auf R. Es gibt eben eine Reihe von zusätzlichen Pakten, die es dir erlauben noch mehr aus einem statistischen Modell rauszuholen. Am Ende wurde es mir dann aber zu detailliert alle Pakete in jedem Kapitel vorzustellen und anzuwenden. Das ist auch nicht immer sinnig. Häufig willst du erstmal nur das Modell rechnen. Später kannst du dann noch tiefer ins Detail gehen oder aber komplexere Verfahren nutzen. Ich tue mich also etwas schwer dieses Kapitel einzuordnen. Entweder packen wir es ans Ende vom statistischen Modellieren und schauen, dann wie wir alles in R machen. Das steht aber etwas der Intuition entgegen, dass wir in jedem Kapitel zum statistischen Modellieren ja schon was selber machen wollen. In R gibt es dazu dann noch sehr gute Pakete, die das Modellieren sehr viel einfacher machen, dabei dann aber auch für den Anfänger etwas komplexer sind. Ich habe mich daher entschieden, diese aktuelle und komplexere Modellierung einmal hier am Anfang vorzustellen und in den folgenden Kapiteln teilweise darauf zu verweisen, wenn ich es sinnig fand. Du kannst alle Modelle auf althergebrachte Art und Weise rechnen ohne was zu verpassen. Aber manchmal möchte man dann auch effizienter Modellieren. Dafür ist dann dieses Kapitel da: Eine erweiterte Idee von der statistischen Modellierung zu erlangen. Fangen wir also erstmal mit der naheliegenden Frage an.
- Was modellieren wir eigentlich?
-
Wir modellieren die Varianz in den Daten oder anders formuliert, wir versuchen die Varianz in den Daten den jeweiligen Quellen zuzuordnen und darüber dann neue Erkenntnisse zu erlangen.
Was ist ein Modell
Statistical Thinking for the 21st Century — What is a model?
Models are about what changes, and what doesn’t
Statistical Thinking for the 21st Century — Practical statistical modeling?
- Specify your question of interest
- Identify or collect the appropriate data
- Prepare the data for analysis
- Determine the appropriate model
- Fit the model to the data
- Criticize the model to make sure it fits properly
- Test hypothesis and quantify effect size
Als nächstes wollen wir uns die Frage stellen, was ist eigentlich das ziel des Modellierens? Wir wollen ja mit der Modellierung der Varianz irgendwas erreichen. In der Abbildung 52.1 siehst du einmal die drei großen Fragefelder, die wir mit einer Modellierung bearbeiten können.
In diesem Kapitel wollen wir uns auch mit dem Thema marginale Effektmodelle (eng. marginal effect models) beschäftigen. Auch hier ist der deutsche Begriff nicht gebräuchlich, so dass ich mich hier dann immer auf die marginal effect models beziehen werde. Ein zweideutiger Aspekt der marginal effect models besteht darin, dass das Wort “marginal” auf zwei verschiedene und etwas entgegengesetzte Arten verwendet wird:
- Bei “marginalen Effekten” (eng. marginal effects) beziehen wir uns auf die Auswirkung einer winzigen (marginalen) Änderung des Regressors auf das Ergebnis. Dies ist eine Steigung oder eine Ableitung.
- Der Begriff “marginaler Mittelwert” (eng. marginal means) bezieht sich auf den Prozess der Marginalisierung über die Zeilen eines Vorhersagerasters. Dies ist ein Durchschnitt oder ein Integral.
Mehr dazu dann in dem Abschnitt zu den marginal effect models und dem R Paket {marginaleffects}
.
Dafür verweise ich auf Heiss (2022) und das Onlinetutorium Marginalia: A guide to figuring out what the heck marginal effects, marginal slopes, average marginal effects, marginal effects at the mean, and all these other marginal things are. Heiss erklärt dort sehr schön die Zusammenhänge. Dazu dann noch der Verweis auf die Webseite Model to Meaning von Arel-Bundock et al. (2024) um nochmal tiefer in das Modellieren von Daten einzusteigen.
52.1 Sprachlicher Hintergrund
“In statistics courses taught by statisticians we don’t use”independent variable” because we use independent on to mean stochastic independence. Instead we say predictor or covariate (either). And, similarly, we don’t say “dependent variable” either. We say response.” — User berf auf r/AskStatistics
Wenn wir uns mit dem statistischen Modellieren beschäftigen wollen, dann brauchen wir auch Worte um über das Thema reden zu können. Statistik wird in vielen Bereichen der Wissenschaft verwendet und in jedem Bereich nennen wir dann auch Dinge anders, die eigentlich gleich sind. Daher werde ich mir es hier herausnehmen und auch die Dinge so benennen, wie ich sie für didaktisch sinnvoll finde. Wir wollen hier was verstehen und lernen, somit brauchen wir auch eine klare Sprache.
Fangen wir also erstmal allgemeiner an ein Modell und deren schreibweise zu verstehen. Da wir uns natürlich in R bewegen für die parktische Anwendung, nutzen wir auch die Modellschreibweise, die in R üblich ist. In R wird diese Schreibweise auch formula()
genannt. Im Folgenden siehst du einmal ein Modell in einer abstrakten Form. Wir haben den Messwert \(Y\) auf der linken Seite (eng. left hand side, abk. LHS) der Tilde und die Einflussvariablen \(X\) auf der rechten Seite (eng. right hand side, abk. RHS). Dabei steht dann das \(X\) hier einmal als Platzhalter und Sammelbegriff für verschiedene Arten von möglichen Variablen.
In R sieht es dann etwas anders aus, da wir die Platzhalter \(Y\) für den Messwert und \(X\) für die Einflussvariable durch die Namen der Spalten in unserem Datensatz ersetzen. Der Satendatz liegt dann als tibble()
in R vor. Mehr dann dazu in den folgenden Beispielen in den jeweiligen Kapiteln zum Modellieren. Dann sieht das Modell wie in der folgenden Abbildung aus.
In dem folgenden Kasten erkläre ich nochmal den Gebrauch meiner Begriffe im statistischen Testen. Es ist wichtig, dass wir hier uns klar verstehen. Zum einen ist es angenehmer auch mal ein Wort für ein Symbol zu schreiben. Auf der anderen Seite möchte ich aber auch, dass du dann das Wort richtig einem Konzept im statistischen Modellieren zuordnen kannst. Deshalb einmal hier meine persönliche und didaktische Zusammenstellung meiner Wort im statistischen Modellieren.
“Die Grenzen meiner Sprache bedeuten die Grenzen meiner Welt.” — Ludwig Wittgenstein
Hier kommt einmal die Tabelle mit den wichtigsten Begriffen im statistischen Modellieren und wie ich die Worte benutzen werde. Damit wir uns verstehen und du was lernen kannst. In anderen Büchern und Quellen findest du teilweise die Worte in einem anderen Sinnzusammenhang. Das ist gut so dort. Bei mir ist es anders.
Symbol | Deutsch | Englisch | |
---|---|---|---|
LHS | \(Y\) / \(y\) | Messwert / Endpunkt / Outcome | response / outcome / endpoint |
RHS | \(X\) / \(x\) | Einflussvariable / Erklärende Variable / Fester Effekt | risk factor / explanatory / fixed effect |
RHS | \(Z\) / \(z\) | Zufälliger Effekt | random effect |
\(X\) ist kontinuierlich | \(c_1\) | Kovariate 1 | covariate 1 |
\(X\) ist kategorial | \(f_A\) | Faktor \(A\) mit Level \(A.1\) bis \(A.j\) | factor \(A\) with levels \(A.1\) to \(A.j\) |
Dann wäre es ja schön, wenn wir nur die linke und rechte Seite neben einer Tilde hätten. Das ist aber nur eine sehr abstrakte Darstellung. Es ha ja auch seinen Grund, warum wir sehr viele Kapitel in diesem Openbook dem Thema des statistischen Modellieren widmen. Wir haben nämlich eine eichtig schrecklich nette Familie an Möglichkeiten zusammen.
In der folgenden Abbildung siehst du einmal wie alles mit allem zusammenhängt. Auf der linken Seite siehst du den Messwert \(Y\) der einer Verteilunsgfamilie entstammt. Je nachdem was du wie gemessen hast, folgt dein Messwert \(Y\) einer anderen Verteilung. Konkreter noch, welche Zahlen du für deinen Messwert \(Y\) bestimmt hast. Auf der rechten Seite findest du die Einflussvariable \(X\), die aus mehren Variablen bestehen kann aber nicht muss. Wenn du eine kontinuierliche Einflussvariable vorliegen hast, dann sprechen wir von Kovariaten. Hast du dagegen kategoriale Einflussvariablen, dann sprechen wir von Faktoren mit Leveln als die Gruppen. Je nach Kombination aus Verteilungsfamilie und Einflussvariable hast du dann eine andere Interpretation der Modellierung vorliegen.
Da wir die schrecklich nette Familie ja auch irgendiwe bezeichnen müssen, hat sich folgende Semantik mehr oder minder durchgesetzt. Ich nutze jedenfalls den folgenden Aufbau um zu benennen was ich eigentlich analysieren und modellieren will. Zuerst kommt, ob wir eine Einflussvariable oder mehrere Einflussvariablen betrachten. Wir nennen dann eben die Modellierung eine simple oder multiple Modellierung. Dann kommt die Verteilunsgfamilie des Messwerts als Wort um dann noch zu sagen, ob wir es mit einem gemischten Modell zu tun haben. Haben wir kein gemischtes Modell, dann lassen ignorieren wir den Teil. Häufig sprechen wir auch von einer linearen Regression, wenn wir eine Gaussian linear Regression meinen. Das finde ich aber sehr verwirrend und nicht klar. Deshalb vermeide ich diesen Sprachgebrauch, wenn wir es mit komplexeren Modellen zu tun haben.
52.2 Genutzte R Pakete
Wir wollen folgende R Pakete in diesem Kapitel nutzen.
R Code [zeigen / verbergen]
::p_load(tidyverse, emmeans, multcomp, ggpmisc, conflicted)
pacmanconflicts_prefer(dplyr::select)
An der Seite des Kapitels findest du den Link Quellcode anzeigen, über den du Zugang zum gesamten R-Code dieses Kapitels erhältst.
52.3 Theoretisches Modellieren in R
R Code [zeigen / verbergen]
reformulate(termlabels = "hase",
response = "igel")
igel ~ hase
52.3.1 R Pakete
Neben den R Paketen, die wir in den jeweiligen Kapiteln brauchen, kommen noch folgende R Pakete immer wieder dran. Deshalb sind die R Pakete hier schon mal mit den jeweiligen Internetseiten aufgeführt.
- Das Buch Tidy Modeling with R gibt nochmal einen tieferen Einblick in das Modellieren in R. Wir immer, es ist ein Vorschlag aber kein Muss.
- Das R Paket
{parameters}
nutzen wir um die Parameter eines Modells aus den Fits der Modelle zu extrahieren. Teilweise sind die Standardausgaben der Funktionen sehr unübersichtich. Hier hilft das R Paket. - Das R Paket
{performance}
hilft uns zu verstehen, ob die Modelle, die wir gefittet haben, auch funktioniert haben. In einen mathematischen Algorithmus können wir alles reinstecken, fast immer kommt eine Zahl wieder raus. - Das R Paket
{tidymodels}
nutzen wir als das R Paket um mit Modellen umgehen zu können und eine Vorhersage neuer Daten zu berechnen. Das Paket{tidymodels}
ist wie das Paket{tidyverse}
eine Sammlung an anderen R Paketen, die wir brauchen werden.
52.3.2 Quantity and tests
“Wenn … der Fall wäre, dann wäre …”
52.3.3 Generalisierung von lm()
zu glm()
und [g]lmer()
- Die Funktion
lm()
nutzen wir, wenn das Outcome \(y\) einer Normalverteilung folgt. - Die Funktion
glm()
nutzen wir, wenn das Outcome \(y\) einer andere Verteilung folgt. - Die Funktion
lmer()
nutzen wir, wenn das Outcome \(y\) einer Normalverteilung folgt und wir noch einen Block- oder Clusterfaktor vorliegen haben. - Die Funktion
glmer()
nutzen wir, wenn das Outcome \(y\) einer andere Verteilung folgt und wir noch einen Block- oder Clusterfaktor vorliegen haben.
52.4 Kovariate Modelle (\(x\) ist kontinuierlich)
52.5 Faktorielle Modelle (\(x\) ist kategorial)
52.6 Interpretation
52.6.1 Simple Modelle
R Code [zeigen / verbergen]
<- tibble(c_1 = rnorm(10, 0, 1),
cov1_tbl y = 2 +
2 * c_1 +
rnorm(10, 0, 0.001))
R Code [zeigen / verbergen]
lm(y ~ c_1, cov1_tbl) |> coef() |> round(2)
(Intercept) c_1
2 2
R Code [zeigen / verbergen]
<- tibble("1" = rnorm(5, 3, 0.001),
fac1_tbl "2" = rnorm(5, 6, 0.001)) |>
gather(key = "f_a", value = "y")
R Code [zeigen / verbergen]
lm(y ~ f_a, fac1_tbl) |> coef() |> round(2)
(Intercept) f_a2
3 3
Hier dann mal mit drei Leveln oder Gruppen
R Code [zeigen / verbergen]
<- tibble("1" = rnorm(5, 2, 0.001),
fac1_tbl "2" = rnorm(5, 6, 0.001),
"3" = rnorm(5, 4, 0.001)) |>
gather(key = "f_a", value = "y")
Treatment coding mit contr.treatment
(default)
R Code [zeigen / verbergen]
lm(y ~ f_a, fac1_tbl, contrasts = list(f_a = "contr.treatment")) |>
coef() |> round(2)
(Intercept) f_a2 f_a3
2 4 2
Effect coding mit contr.sum
R Code [zeigen / verbergen]
lm(y ~ f_a, fac1_tbl, contrasts = list(f_a = "contr.sum")) |>
coef() |> round(2)
(Intercept) f_a1 f_a2
4 -2 2
52.6.2 Multiple Modelle
R Code [zeigen / verbergen]
<- tibble(c_1 = rnorm(10, 0, 1),
cov2_tbl c_2 = rnorm(10, 0, 1),
y = 2 +
1 * c_1 +
2 * c_2 +
rnorm(10, 0, 0.001))
R Code [zeigen / verbergen]
lm(y ~ c_1 + c_2, cov2_tbl) |>
coef() |> round(2)
(Intercept) c_1 c_2
2 1 2
Einfaktorielle ANCOVA
R Code [zeigen / verbergen]
<- tibble(c_1 = rnorm(15, 0, 1),
cov1_fac1_tbl f_a = gl(3, 5),
y = 2 +
1 * c_1 +
2 * as.numeric(f_a) +
rnorm(15, 0, 0.001))
R Code [zeigen / verbergen]
lm(y ~ f_a + c_1, cov1_fac1_tbl) |>
coef() |> round(2)
(Intercept) f_a2 f_a3 c_1
4 2 4 1
R Code [zeigen / verbergen]
<- tibble(c_1 = rnorm(30, 0, 1),
cov1_fac2_tbl f_a = rep(gl(3, 5), 2),
f_b = gl(2, 15),
y = 2 +
1 * c_1 +
2 * as.numeric(f_a) +
3 * as.numeric(f_b) +
rnorm(15, 0, 0.001))
R Code [zeigen / verbergen]
lm(y ~ c_1 + f_a + f_b, cov1_fac2_tbl) |>
coef() |> round(2)
(Intercept) c_1 f_a2 f_a3 f_b2
7 1 2 4 3
R Code [zeigen / verbergen]
<- tibble(f_a = rep(gl(3, 5), 2),
fac2_tbl f_b = gl(2, 15),
y = 2 +
2 * as.numeric(f_a) +
3 * as.numeric(f_b) +
rnorm(15, 0, 0.001))
R Code [zeigen / verbergen]
lm(y ~ f_a + f_b + f_a:f_b, fac2_tbl) |>
coef() |> round(2)
(Intercept) f_a2 f_a3 f_b2 f_a2:f_b2 f_a3:f_b2
7 2 4 3 0 0
In sehen wir wie wir den Namen einer Regression bilden. Zuerst entscheiden wir, ob wir nur ein \(x\) haben oder mehrere. Mit einem \(x\) sprechen wir von einem simplen Modell, wenn wir mehrere \(x\) haben wir ein multiples Modell. Im nächsten Schritt benennen wir die Verteilung für das Outcome \(y\). Dann müssen wir noch entscheiden, ob wir ein gemischtes Modell vorliegen haben, dann schreiben wir das hin. Sonst lassen wir den Punkt leer. Anschließend kommt noch lineares Modell hinten ran.
Veraltet aber manchmal ganz nützlich das R Paket {modelr}
52.7 Konfidenzintervall vs. Prädiktionsintervall
Ein Konfidenzintervall gibt den Wertebereich für einen gesuchten Parameter der Grundgesamtheit mit einer bestimmten Wahrscheinlichkeit an. Ein Prognoseintervall gibt den Wertebereich für einen individuellen, zukünftig zu beobachtenden Wert mit einer bestimmten Wahrscheinlichkeit an.
Konfidenzintervall
Text
- Konfidenzintervall
-
Ein Konfidenzintervall gibt den Wertebereich für einen gesuchten Parameter der Grundgesamtheit mit einer bestimmten Wahrscheinlichkeit an.
Prädiktionsintervall
Text
- Prädiktionsintervall
-
Ein Prädiktionsintervall (auch Vorhersageintervall oder Prognoseintervall) gibt den Wertebereich für einen individuellen, zukünftig zu beobachtenden Wert mit einer bestimmten Wahrscheinlichkeit an.
Quantile Regression Forests for Prediction Intervals
The difference between prediction intervals and confidence intervals
P-values for prediction intervals machen keinen Sinn
52.8 Praktisches Modellieren in R
52.8.1 Marginal Effect Models
52.8.2 Visualisierung
52.8.3 …mit {purrr}
map()
für mehrere Endpunkte wie bei den Spinatdaten in Application gezeigt.
map2()
mit Familie für glm
zeigen, wenn wir unterschiedliche Outcomes haben.
R Code [zeigen / verbergen]
<- read_excel("data/multiple_outcomes.xlsx") |>
cutting_tbl mutate(trt = as_factor(trt),
block = as_factor(block)) |>
mutate_if(is.numeric, round, 2)
trt | block | shoot | leaf | flower | fruit | ca | drymatter | freshweight | height |
---|---|---|---|---|---|---|---|---|---|
control | 1 | 0 | 136 | 33 | 10 | 0.41 | 11.59 | 85.37 | 25 |
control | 2 | 0 | 157 | 33 | 9 | 0.41 | 10.43 | 79.78 | 22 |
control | 1 | 1 | 65 | 10 | 3 | 0.4 | 11.97 | 78.35 | 22 |
control | 2 | 0 | 88 | 16 | 4 | 0.39 | 11.24 | 85.55 | 27 |
… | … | … | … | … | … | … | … | … | … |
nodium_5th | 2 | 0 | 137 | 20 | 3 | 0.47 | 11.11 | 86.93 | 28 |
nodium_5th | 1 | 0 | 137 | 23 | 3 | 0.43 | 11.95 | 89.68 | 27 |
nodium_5th | 2 | 1 | 114 | 17 | 9 | 0.39 | 11.21 | 70.82 | 24 |
nodium_5th | 2 | 0 | 192 | 42 | 11 | 0.42 | 11.29 | 66.23 | 22 |
R Code [zeigen / verbergen]
<- lst(shoot = binomial(),
family_lst leaf = quasipoisson(),
flower = quasipoisson(),
fruit = quasipoisson(),
ca = gaussian(),
drymatter = gaussian(),
freshweight = gaussian(),
height = gaussian())
R Code [zeigen / verbergen]
<- cutting_tbl |>
cutting_long_tbl pivot_longer(cols = shoot:last_col(),
names_to = "outcome",
values_to = "rsp") |>
arrange(outcome, trt, block)
R Code [zeigen / verbergen]
<- cutting_long_tbl |>
cutting_lst split(~outcome)
R Code [zeigen / verbergen]
<- cutting_lst %>%
glm_lst map2(family_lst,
~glm(rsp ~ trt + block + trt:block,
data = .x, family = .y))
%>%
glm_lst map(pluck, "family")
$ca
Family: binomial
Link function: logit
$drymatter
Family: quasipoisson
Link function: log
$flower
Family: quasipoisson
Link function: log
$freshweight
Family: quasipoisson
Link function: log
$fruit
Family: gaussian
Link function: identity
$height
Family: gaussian
Link function: identity
$leaf
Family: gaussian
Link function: identity
$shoot
Family: gaussian
Link function: identity
R Code [zeigen / verbergen]
%>%
glm_lst map(car::Anova)
$ca
Analysis of Deviance Table (Type II tests)
Response: rsp
LR Chisq Df Pr(>Chisq)
trt 0.222384 3 0.9739
block 0.003979 1 0.9497
trt:block 0.005906 3 0.9999
$drymatter
Analysis of Deviance Table (Type II tests)
Response: rsp
LR Chisq Df Pr(>Chisq)
trt 5.1410 3 0.1618
block 0.2084 1 0.6480
trt:block 2.0862 3 0.5547
$flower
Analysis of Deviance Table (Type II tests)
Response: rsp
LR Chisq Df Pr(>Chisq)
trt 19.3416 3 0.0002323 ***
block 0.5678 1 0.4511477
trt:block 1.7967 3 0.6156522
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
$freshweight
Analysis of Deviance Table (Type II tests)
Response: rsp
LR Chisq Df Pr(>Chisq)
trt 3.6454 3 0.3024
block 4.0375 1 0.0445 *
trt:block 4.9674 3 0.1742
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
$fruit
Analysis of Deviance Table (Type II tests)
Response: rsp
LR Chisq Df Pr(>Chisq)
trt 9.9860 3 0.01869 *
block 0.8118 1 0.36759
trt:block 1.7163 3 0.63332
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
$height
Analysis of Deviance Table (Type II tests)
Response: rsp
LR Chisq Df Pr(>Chisq)
trt 2.83465 3 0.4178
block 0.12598 1 0.7226
trt:block 0.94488 3 0.8146
$leaf
Analysis of Deviance Table (Type II tests)
Response: rsp
LR Chisq Df Pr(>Chisq)
trt 15.5183 3 0.001423 **
block 1.2623 1 0.261207
trt:block 1.6400 3 0.650353
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
$shoot
Analysis of Deviance Table (Type II tests)
Response: rsp
LR Chisq Df Pr(>Chisq)
trt 17 3 0.0007067 ***
block 0 1 1.0000000
trt:block 2 3 0.5724067
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
R Code [zeigen / verbergen]
<- glm_lst %>%
emm_lst map(~emmeans(.x, specs = ~ trt, type = "response"))
%>%
emm_lst map(~contrast(.x, method = "pairwise", adjust = "bonferroni")) %>%
map(as_tibble) %>%
bind_rows(.id = "outcome")
# A tibble: 48 × 11
outcome contrast odds.ratio SE df null z.ratio p.value ratio
<chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 ca control / str… 0.829 0.978 Inf 1 -0.159 1 NA
2 ca control / nod… 1.30 1.59 Inf 1 0.219 1 NA
3 ca control / nod… 0.780 0.918 Inf 1 -0.211 1 NA
4 ca strong / nodi… 1.57 1.90 Inf 1 0.376 1 NA
5 ca strong / nodi… 0.941 1.10 Inf 1 -0.0523 1 NA
6 ca nodium_3rd / … 0.598 0.719 Inf 1 -0.428 1 NA
7 drymatter control / str… NA 0.0436 Inf 1 -0.425 1 0.981
8 drymatter control / nod… NA 0.0479 Inf 1 1.24 1 1.06
9 drymatter control / nod… NA 0.0424 Inf 1 -0.930 1 0.960
10 drymatter strong / nodi… NA 0.0486 Inf 1 1.67 0.571 1.08
# ℹ 38 more rows
# ℹ 2 more variables: estimate <dbl>, t.ratio <dbl>
R Code [zeigen / verbergen]
%>%
emm_lst map(~cld(.x, Letters = letters, adjust = "bonferroni")) %>%
map(as_tibble) %>%
bind_rows(.id = "outcome") %>%
select(outcome, trt, .group) %>%
print(n = 28)
# A tibble: 32 × 3
outcome trt .group
<chr> <fct> <chr>
1 ca nodium_3rd " a"
2 ca control " a"
3 ca strong " a"
4 ca nodium_5th " a"
5 drymatter nodium_3rd " a"
6 drymatter control " a"
7 drymatter strong " a"
8 drymatter nodium_5th " a"
9 flower strong " a "
10 flower nodium_3rd " ab"
11 flower nodium_5th " b"
12 flower control " b"
13 freshweight strong " a"
14 freshweight control " a"
15 freshweight nodium_3rd " a"
16 freshweight nodium_5th " a"
17 fruit strong " a"
18 fruit control " a"
19 fruit nodium_3rd " a"
20 fruit nodium_5th " a"
21 height strong " a"
22 height control " a"
23 height nodium_3rd " a"
24 height nodium_5th " a"
25 leaf strong " a "
26 leaf nodium_3rd " ab"
27 leaf nodium_5th " ab"
28 leaf control " b"
# ℹ 4 more rows
52.9 Verschiebebahnhof
in diesem Kapitel wollen Fokus auf {parsnip}
und dann als Kasten in allen anderen möglichen Modellierungen ergänzen?
Tidymodels, interactions and anova - a short tutorial
{parsnip}
Auch hier können wir die XX Regression in dem R Paket {parsnip}
realisieren. Ein Vorteil von {parsnip}
ist, dass wir die Funktionen sehr gut mit dem |>
-Operator nutzen können. Das ist bei den Funktionen glm()
und lm()
in der ursprünglichen Implementierung ja leider nur etwas umständlicher möglich. Deshalb hier einmal die bessere Implementierung.
R Code [zeigen / verbergen]
<- linear_reg() |>
linreg_reg_fit set_engine("lm") |>
fit(jump_length ~ grp, data = simple_tbl)
Dann die ANOVA
R Code [zeigen / verbergen]
|>
linreg_reg_fitextract_fit_engine() |>
anova()
Analysis of Variance Table
Response: jump_length
Df Sum Sq Mean Sq F value Pr(>F)
grp 1 1.5309 1.53089 9.2541 0.01879 *
Residuals 7 1.1580 0.16543
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
R Code [zeigen / verbergen]
<- lm(jump_length ~ grp, data = simple_tbl) linreg_reg_fit
R Code [zeigen / verbergen]
|>
linreg_reg_fitanova()
Analysis of Variance Table
Response: jump_length
Df Sum Sq Mean Sq F value Pr(>F)
grp 1 1.5309 1.53089 9.2541 0.01879 *
Residuals 7 1.1580 0.16543
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Eine detailliertere Einführung mit mehr Beispielen für die Nutzung vom R Paket {parsnip}
findest du im Kapitel Modellieren in R. Hier soll es dann bei der kurzen Gegenüberstellung bleiben.