PS-Stratifizierung
1 Der Datensatz
Wir möchten evaluieren, ob ein entwickeltes Englischtraining die
Englischkenntnisse von Studierenden erhöht. Dazu wurde eine
quasi-experimentelle Studie durchgeführt, bei der die Studierenden
selbst entscheiden konnten, ob sie am Training teilnehmen
(Interventionsgruppe; IG) oder nicht (Kontrollgruppe; KG). Nach dem
Training wurden in beiden Gruppen die Englischkenntnisse gemessen. Als
Kovariaten wurden die Englischvorkenntnisse und die Vorliebe für
Englisch erhoben. Für die Analysen nutzen wir den simulierten Datensatz
(englisch.rda) mit dem Objektnamen englisch
mit den folgenden Variablen:
- treat: Teilnahme am Englischtraining (0 = kein Englischtraining, 1 = Englischtraining) (Interventionsvariable)
- eng_post: Englischkenntnisse nach dem Training (Skala von 0 bis 5) (Outcome-Variable)
- eng_pre: Englischvorkenntnisse (Skala von 0 bis 5) (Kontinuierliche Kovariate)
- eng_vorliebe: Vorliebe für Englisch (0 = keine Vorliebe für Englisch, 1 = Vorliebe für Englisch) (Dichotome Kovariate)
2 Erklärung der R-Funktionen
Für die PS-Stratifizierung verwenden wir die Funktion
matchit() aus dem Paket MatchIt. Die folgende
Tabelle bietet eine Übersicht der zentralen Funktionsargumente von
matchit() für die Stratifizierung:
| Argument | Bedeutung |
|---|---|
matchit( |
|
treat ~ x1 + x2, |
Angabe der Interventionsvariable (treat)
und einer oder mehrerer Kovariaten (x1 + x2 + ...) |
data, |
Datensatz |
method, |
Wahl der Matchingmethode; für PS-Stratifizierung
method = "subclass" |
distance, |
Wahl des Distanzmaßes (z.B. direkte Eingabe einer
Distanz), einer Methode zur Berechnung der Distanz (z.B.
distance = "mahalanobis" für die Mahalanobis-Distanz) oder
einer Methode zur Schätzung der Propensity-Scores (z.B.
distance = "glm" für generalized linear model;
default) |
link, |
Art der Link-Funktion bei der Schätzung der PS; z.B.
link = "logit" für logistische Funktion |
subclass, |
Angabe, wie viele Strata gebildet werden sollen |
estimand, |
Angabe, welcher kausale Effekt geschätzt werden soll
(estimand = "ATE" für den durchschnittlichen kausalen
Effekt (default), estimand = "ATT" für den
durchschnittlichen kausalen Effekt der Behandelten,
estimand = "ATC" den durchschnittlichen kausalen Effekt der
Unbehandelten) |
discard) |
Angabe, ob Daten, die nicht Teil der Region of Common
Support sind, ausgeschlossen werden sollen; z.B.
discard = "both", wenn Personen der Interventions- und
Kontrollgruppe ausgeschlossen werden sollen,
discard = "treated", wenn nur Personen aus der
Interventionsgruppe und discard = "control", wenn nur
Personen der Kontrollgruppe ausgeschlossen werden sollen; bei
discard = "none" werden keine Daten ausgeschlossen
(default) |
Mit der Funktion love.plot() aus dem Paket
cobalt können Unterschiede in Verteilungskennwerten auf den
Kovariaten vor und nach der PS-Stratifizierung visualisiert werden. Die
relevanten Argumente der Funktion sind in der folgenden Tabelle
erklärt.
| Argument | Bedeutung |
|---|---|
love.plot( |
|
x, |
matchit-Objekt |
stats, |
Kennwerte, die dargestellt werden sollen; z.B.
stats = "mean.diffs" für Mittelwertdifferenzen,
stats = "variance.ratios" für Varianzverhältnisse |
disp.subclass, |
Anzeigen der Kennwerte der Balance für jedes Stratum
einzeln mit TRUE (default = FALSE) |
thresholds, |
Markierung von Grenzwerten für eine akzeptable Balance
(optional); z.B. nach Steiner et al. (2010) für Mittelwertsdifferenzen
thresholds = c(-0.1, 0.1) |
sample.names, |
frei wählbare Namen für die unadjustierten und
adjustierten Ergebnisse in der Grafik; z.B.
sample.names = c("unadjustiert", "adjustiert") |
stars) |
Angabe von Variablenarten, die in der Grafik mit einem
Stern markiert werden sollen; z.B. stars = "raw" für
Variablen mit nicht-standardisierten Mittelwertsdifferenzen |
Mit der Funktion bal.plot() aus dem Paket
cobalt können die Verteilungen der Kovariaten zwischen den
Interventionsgruppen vor und nach der PS-Stratifizierung visualisiert
werden. Die relevanten Argumente der Funktion sind in der folgenden
Tabelle erklärt.
| Argument | Bedeutung |
|---|---|
bal.plot( |
|
x, |
matchit-Objekt |
var.name, |
Variable, für die die Kovariatenbalance dargestellt werden soll |
which, |
Angabe, für welche Daten die Verteilungen dargestellt
werden sollen; z.B. which = "unadjusted" nur für
unadjustierte Daten; which = "adjusted" nur für adjustierte
Daten; which = "both" für beide |
type, |
Art der Darstellungsform; z.B.
type = "histogram" für ein Histogramm;
type = "density" für Dichte |
mirror, |
Spiegelung gruppenspezifischer Verteilung an der
x-Achse bei mirror = TRUE |
sample.names) |
frei wählbarer Name für die unadjustierte Ergebnisse
vor der Stratifizierung in der Grafik; z.B.
sample.names = "unadjustiert"; bei der Darstellung von
Strata werden als Überschriften für die Strata automatisch “Subclass”
und die zugehörige Nummer gesetzt |
3 Vorbereitung
Zunächst installieren (falls nötig) und laden wir die nötigen Pakete:
# install.packages("MatchIt") # Installation vor der ersten Nutzung des Pakets Matchit
library(MatchIt) # Laden des Pakets MatchIt
# install.packages("cobalt") # Installation vor der ersten Nutzung des Pakets cobalt
# # Hinweis: Wir verwenden Version 4.4.0 des Paketes
library(cobalt) # Laden des Pakets cobaltund lesen den Datensatz ein:
load("englisch.rda") # Einlesen des Datensatzes (wenn der Datensatz nicht im
# Arbeitsverzeichnis liegt, muss der Pfad entsprechend dem
# Speicherort des Datensatzes angepasst werden)4 Matchen von Personengruppen
Bei der PS-Stratifizierung werden zunächst Gruppen von Personen hinsichtlich ihrer Propensity-Scores (oder Logits der PS) gebildet. Dazu wird eine Anzahl an Strata festgelegt Die Idee ist, dass sich innerhalb eines Stratums Interventions- und Kontrollgruppe in den Kovariatenausprägungen ähnlich sind. Statt einzelne Personen (wie im PS-Matching) werden nun Gruppen, welche durch die Strata beschrieben sind, zueinander gematcht. Dabei muss innnerhalb eines Stratums die Anzahl der Personen in der Interventionsgruppe nicht der Anzahl der Personen in der Kontrollgruppe entsprechen, in der Regel tut sie das auch nicht.
Mit der matchit-Funktion führen wir die Stratifizierung
durch. Wir teilen die Stichprobe dabei in fünf Strata, schließen keine
Einheiten außerhalb der Region of Common Support aus und schätzen den
ATE.
stra5 <- matchit(treat ~ eng_pre + eng_vorliebe, # beide Kovariaten für die PS-Schätzung nutzen
data = englisch, # der Datensatz
distance = "glm", # PS auf Basis der glm Funktion
link = "logit", # mit Logit-Link
method = "subclass", # Stratifizierung
subclass = 5, # Anzahl der Strata
discard = "none", # kein Ausschluss der Daten
# außerhalb der Region of Common Support
estimand = "ATE") # ATE als ZieleffektÜber den summary()-Befehl können wir uns die Ergebnisse
des Matchings ausgeben lassen. Mithilfe des Arguments
sublass = TRUE werden die Kennwerte zur Beurteilung der
Balance innerhalb der einzelnen Strata ausgegeben.
summary(stra5, subclass = TRUE) # Ergebnis des Matchings nach Strata ...
## Sample Sizes by Subclass:
## 1 2 3 4 5 All
## Control 72 56 45 41 18 232
## Treated 17 32 43 47 71 210
## Total 89 88 88 88 89 442
...Unter Sample Sizes by
Subclass können wir uns die Verteilung der Personen aus
Interventions- und Kontrollgruppe über die fünf Strata ansehen. Wir
sehen also beispielsweise, dass im ersten Stratum 72
Personen aus der Kontrollgruppe zu 17 Personen aus der
Interventionsgruppe gematcht wurden. Im fünften Stratum wurden
18 Personen aus der Kontrollgruppe zu 71
Personen aus der Interventionsgruppe gematcht. In jedem der fünf Strata
sind insgesamt 88 oder 89 Personen. Die Strata
sind gleich groß, da ihre Größe anhand der Quantile der gesamten
berücksichtigten Stichprobe festgelegt wurde (durch das Argument
estimand = "ATE" spezifiziert).
Wir können den Overlap bzw. die PS-Verteilung der Personen aus
Interventions- und Kontrollgruppe in den fünf Strata mithilfe der
plot()-Funktion visualisieren:
plot(stra5, type = 'jitter', interactive = FALSE)
Abbildung 1: Propensity Scores der gematchten und ungematchten Personen getrennt nach Interventionsgruppe. Senkrechte Linien kennzeichnen die Strata
In der Abbildung werden die Propensity Scores der zugeordneten
Personen getrennt nach Kontroll- und Interventionsgruppe dargestellt.
Die senkrechten Linien kennzeichnen die Einteilung in die 5 Strata Wir
sehen zum Beispiel, dass die Propensity Scores der Personen in der
Kontrollgruppe im ersten Stratum geringer als die Propensity Scores der
Interventionsgruppe im selben Stratum sind. Dies deutet bereits auf eine
Imbalance hin. Da wir durch das Argument discard = "none"
angegeben haben, dass Daten außerhalb der Region of Common Support nicht
ausgeschlossen werden sollen, sind keine Einheiten unter Unmatched
Control oder Unmatched Treated Units aufgeführt.
5 Prüfung der Balance
Wir können die Balance sowohl anhand der Verteilungskennwerte als auch anhand der Verteilungen prüfen.
5.1 Kennwerte
summary(stra5, subclass = TRUE) ...
## Summary of Balance for All Data:
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.56 0.398 0.884 0.879 0.234
## eng_pre 2.78 2.396 0.787 0.586 0.216
## eng_vorliebe 0.70 0.500 0.417 . 0.200
## eCDF Max
## distance 0.365
## eng_pre 0.355
## eng_vorliebe 0.200
##
## Summary of Balance by Subclass:
##
## - Subclass 1
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.230 0.186 0.239 0.457 0.183
## eng_pre 2.057 1.830 0.462 0.421 0.185
## eng_vorliebe 0.235 0.264 -0.060 . 0.029
## eCDF Max
## distance 0.342
## eng_pre 0.313
## eng_vorliebe 0.029
##
## - Subclass 2
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.353 0.350 0.015 1.004 0.058
## eng_pre 2.415 2.324 0.183 0.978 0.113
## eng_vorliebe 0.312 0.482 -0.354 . 0.170
## eCDF Max
## distance 0.156
## eng_pre 0.259
## eng_vorliebe 0.170
##
## - Subclass 3
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.483 0.470 0.072 1.040 0.107
## eng_pre 2.616 2.596 0.041 0.753 0.073
## eng_vorliebe 0.605 0.578 0.056 . 0.027
## eCDF Max
## distance 0.254
## eng_pre 0.193
## eng_vorliebe 0.027
##
## - Subclass 4
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.606 0.593 0.069 1.467 0.088
## eng_pre 2.833 2.860 -0.055 0.792 0.050
## eng_vorliebe 0.809 0.683 0.262 . 0.126
## eCDF Max
## distance 0.255
## eng_pre 0.126
## eng_vorliebe 0.126
##
## - Subclass 5
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.750 0.771 -0.117 0.624 0.079
## eng_pre 3.195 3.329 -0.273 0.482 0.118
## eng_vorliebe 0.972 0.889 0.173 . 0.083
## eCDF Max
## distance 0.180
## eng_pre 0.248
## eng_vorliebe 0.083
##
...Es werden Kennwerte der Verteilung der Kovariaten
(eng_pre und eng_vorliebe) und der
Propensity-Scores (distance) in Interventions- und
Kontrollgruppe sowohl vor der Stratifizierung
(Summary of Balance for All
Data) als auch nach der Stratifizierung pro Stratum
(Summary of Balance by
Subclass) gezeigt.
Dabei bedeutet…
- Means Treated: Mittelwerte der Kovariaten in der Interventionsgruppe
- Means Control: Mittelwerte der Kovariaten in der Kontrollgruppe
- Std. Mean Diff.: stand. Mittelwertsdifferenz der Kovariaten zwischen Interventions- und Kontrollgruppe
- Var. Ratio: Varianzverhältnisse der Kovariaten in Interventions- und Kontrollgruppe
Für kategoriale Variablen, wie die Vorliebe für Englisch, werden keine Varianzverhältnisse ausgegeben. Im Falle einer dichotomen dummy-kodierten (d.h. 0/1 kodierten) Variablen können die Mittelwerte als relative Häufigkeiten interpretiert werden.
Bei der Beurteilung der Balance orientieren wir uns im Folgenden an den von Steiner et al. (2010) vorgeschlagenen Richtlinien für die Balance metrischer Variablen. Wir gehen demnach von Balance aus, wenn der Betrag der standardisierten Mittelwertsdifferenz kleiner als 0,1 Standardabweichungen ist und das Varianzverhältnis zwischen 0,8 und 1,25 liegt.
All Data:
Vor der Adjustierung zeigen sich deutliche Selektionseffekte. Der
durchschnittliche Wert der Propensity-Scores in der Interventionsgruppe
vor dem Matching liegt bei 0,560 und in der Kontrollgruppe
bei 0,398. Die standardisierte Mittelwertsdifferenz der
Propensity-Scores entspricht einem Cohen’s \(d\) von 0,884 und liegt damit
außerhalb der Grenzen nach Steiner et al. (2010) für akzeptable
Mittelwertsdifferenzen. Die Varianz der PSs in der Interventionsgruppe
ist 0,879 mal so groß wie die in der Kontrollgruppe und
liegt somit innerhalb der empfohlenen Grenzen. Es zeigen sich deutliche
Unterschiede in den Englischvortestwerten und der Vorliebe für Englisch
zwischen Interventions- und Kontrollgruppe vor dem Matchen. Die
Interventionsgruppe ist im Durchschnitt um (2,78 -
2,396 =) 0,384 Punkte besser im Englischvortest als die
Kontrollgruppe (Cohen’s \(d\) =
0,787) und weist deutlich homogenere Vorkenntnisse auf als
die Kontrollgruppe (Varianzverhältnis = 0,586). Beide
Kennwerte der Verteilung der Englischvortestwerte liegen also außerhalb
der von Steiner et al. (2010) vorgeschlagenen Grenzen.
Englischliebhabende sind zudem in der Interventionsgruppe stärker
vertreten (70%) als in der Kontrollgruppe
(50%).
by Subclass:
Auch nach der Stratifizierung ist nicht in allen Strata Balance
gegeben. Die standardisierten Mittelwertsdifferenzen der
Propensity-Scores (distance)
sowie der Englischvorkenntnisse sind in den Strata 3 und 4 balanciert.
In Stratum 2 und insbesondere in Stratum 1 und 5 sind sie es jedoch
nicht. Beispielsweise unterscheiden sich die mittleren Vorkenntnisse im
ersten Stratum zwischen Interventionsgruppe und Kontrollgruppe um
0,462 Standardabweichungen. Ebenso unterscheiden sich
insbesondere in den Strata 1, 4 und 5 die Varianzen der Vorkenntnisse
zwischen Interventionsgruppe und Kontrollgruppe stark. Beispielsweise
ist in Stratum 5 die Varianz der Vorkenntnisse in der
Interventionsgruppe nur ungefähr halb so groß wie in der Kontrollgruppe
(Varianzverhältnis = 0,482). Die Anteile der Personen mit
einer Vorliebe für Englisch sind in Stratum 1 sehr ähnlich
(24% in der Interventionsgruppe und 26% in der
Kontrollgruppe), während sie sich beispielsweise in Stratum 2 deutlich
unterscheiden (31% in der Interventionsgruppe und
48% in der Kontrollgruppe).
Die Mittelwertsdifferenzen und Varianzverhältnisse können wir auch grafisch darstellen:
love.plot(stra5, # matchit-Objekt
stats = "mean.diffs", # Mittelwertsdifferenzen
disp.subclass = TRUE, # Mittelwertsdifferenzen innerhalb der Strata
thresholds = c(-0.1, 0.1), # Einzeichnen der Grenzwerte nach
# Steiner et al. (2010)
sample.names = c("unadjustiert", # Beschriftung der Ergebnisse
"adjustiert"),
stars = "raw") # Markierung nicht-standardisierter
# Mittelwertsdifferenzen mit einem Stern
Abbildung 1: (standardisierte) Mittelwertsdifferenzen der Kovariaten und Propensity-Scores vor (rote Punkte) und nach (türkisfarbene Punkte) der Stratifizierung sowie nach Strata (Zahlen) getrennt. Die durchgezogene senkrechte Linie kennzeichnet eine Mittelwertsdifferenz von 0 und somit gleiche Mittelwerte in der Interventionsgruppe und Kontrollgruppe. Die gestrichelten Linien kennzeichnen die Grenzen für akzeptable Mittelwertsdifferenzen nach Steiner et al. (2010). (* kennzeichnet Variablen für die die unstandardisierten Mittelwertsdifferenzen angegeben sind. Für diese gelten die Grenzwerte nicht.)
Die Grafik gibt für das Distanzmaß und die Vorkenntnisse die standardisierten Mittelwertsdifferenzen und für die Vorliebe für Englisch (mit Stern gekennzeichnet) die unstandardisierten Mittelwertsdifferenzen zwischen Interventionsgruppe und Kontrollgruppe an. Die roten Punkte zeigen die (un-)standardisierten Mittelwertsdifferenzen vor der Adjustierung. Wir sehen vor der Adjustierung starke Imbalance für alle Variablen, die Mittelwertsdifferenzen liegen deutlich außerhalb der Grenzen nach Steiner et al. (2010), also außerhalb der gestrichelten Linien. Die türkisfarbenen Punkte stellen über die Strata gewichtete Mittelwertsunterschiede dar. Die gewichteten Mittelwertsunterschiede deuten auf eine eher gute Balance hin, sie liegen alle innerhalb der Grenzen nach Steiner et al. (2010). Schaut man sich jedoch die Mittelwertsunterschiede innerhalb der einzelnen Strata (durch die Zahlen gekennzeichnet) an, so erkennt man - wie auch bereits durch die Kennwerte oben ersichtlich - deutliche Imbalance. Dies trifft im Besonderen auf das erste und das fünfte Stratum zu.
love.plot(stra5, # matchit-Objekt
stats = "variance.ratios", # Varianzverhältnisse
disp.subclass = TRUE, # Varianzverhältnisse innerhalb der Strata
thresholds = c(0.8, 1.25), # Einzeichnen der Grenzwerte nach
# Steiner et al. (2010)
sample.names = c("unadjustiert", # Beschriftung der Ergebnisse
"adjustiert")
) 
Abbildung 2: Varianzverhältnisse der (kontinuierlichen) Kovariate und Propensity-Scores vor (rote Punkte) und nach (türkisfarbene Punkte) der Stratifizierung sowie nach Strata (Zahlen) getrennt. Die durchgezogene senkrechte Linie kennzeichnet ein Varianzverhältnis von 1 und somit gleiche Varianzen in der Interventionsgruppe und Kontrollgruppe. Die gestrichelten Linien kennzeichnen die Grenzen für akzeptable Varianzverhältnisse nach Steiner et al. (2010).
Diese Grafik gibt für das Distanzmaß und die Vorkenntnisse die Varianzverhältnisse zwischen Interventionsgruppe und Kontrollgruppe an. Für die Englischvorliebe (eng_vorliebe) werden keine Varianzverhältnisse ausgegeben, da es sich um eine dichotome Variable handelt. Die roten Punkte zeigen die Varianzverhältnisse vor der Adjustierung. Die Vorkenntnisse (eng_pre) sind vor der Stratifizierung nicht hinreichend balanciert, wohingegen die Propensity-Scores (distance) vor der Stratifizierung innerhalb der Grenzen nach Steiner et al. (2010) liegen und somit als hinreichend balanciert betrachtet werden können. Die türkisfarbenen Punkte stellen über die Strata gewichtete Varianzverhältnisse dar. Zwar liegt das Varianzverhältnis der gewichteten Propensity-Scores nach der Stratifizierung innerhalb der Grenzen nach Steiner et al. (2010), die Vorkenntnisse sind jedoch nicht ausreichend balanciert. Betrachtet man die Varianzverhältnisse innerhalb der einzelnen Strata (durch die Zahlen gekennzeichnet), so erkennt man deutliche Imbalance sowohl bezüglich der Vorkenntnisse (eng_pre) als auch bezüglich der Propensity-Scores (distance). Wie auch bei den Mittelwertsdifferenzen sind Stratum 1 und Stratum 5 besonders stark von Imbalance betroffen.
5.2 Verteilung der Kovariaten
plot <- bal.plot(stra5, # matchit-Objekt
var.name = "eng_pre", # Variable, für die Kovariatenverteilungen
# dargestellt werden sollen
which = "both", # Kovariatenverteilung vor und nach der Stratifizierung
type = "histogram", # Darstellung als Histogramm
mirror = TRUE, # Spiegelung der gruppenspezifischen Verteilungen an x-Achse
sample.names = "unadjustiert") # Benennung der unadjustierten Ergebnisse
plot
Abbildung 3: Verteilung der kontinuierlichen Kovariate eng_pre in der Kontrollgruppe (rote Balken) und in der Interventionsgruppe (türkisfarbene Balken) vor der Stratifizierung (links) und in den einzelnen Strata nach der Stratifizierung (Subclass 1 bis Subclass 5).
In dieser Grafik ist die Verteilung der Englischvorkenntnisse in der Kontrollgruppe (rote Balken) und in der Interventionsgruppe (türkisfarbene Balken) vor der Stratifizierung (unadjustiert) und in den einzelnen Strata nach der Stratifizierung (Subclass 1 bis Subclass 5) abgebildet. Vor der Stratifizierung gibt es eine deutliche Imbalance in den Vorkenntnissen; die Interventionsgruppe verfügt im Mittel über höhere und über Personen homogenere Vorkenntnisse als die Kontrollgruppe. In den Strata 2 bis 4 ähneln sich die Verteilungen der Vorkenntnisse zwischen beiden Interventionsgruppen, was für eine Balance spricht. In Stratum 1 und 5 ist jedoch eine deutliche Imbalance ersichtlilch. Im ersten Stratum haben in der Kontrollgruppe deutlich mehr Personen niedrige Vorkenntnisse als in der Interventionsgruppe, im fünften Stratum ist es andersherum. Auch anhand dieser Grafik können wir also erkennen, dass vor allem im ersten und im fünften Stratum die Balance bezüglich der Englischvorkenntnisse noch nicht hinreichend gegeben ist.
Insgesamt lässt sich festhalten, dass durch das obige Modell keine ausreichende Balance hergestellt werden konnte. Somit muss das Modell respezifiziert werden.
6 Respezifizierung
Wir nehmen verschiedene Änderungen im Modell vor. Um Overlap
herzustellen, schließen wir mit dem Argument
discard = "both" die Daten außerhalb der Region of Common
Support aus. Zusätzlich fügen wir in die Gleichung zur Schätzung der
Propensity-Scores eine Interaktion zwischen den beiden Kovariaten
eng_pre und eng_vorliebe ein, indem wir die
Modellgleichung treat ~ eng_pre + eng_vorliebe durch die
Gleichung treat ~ eng_pre * eng_vorliebe ersetzen und
erhöhen die Anzahl der Strata mit dem Argument
subclass = 10 von fünf auf zehn.
stra10 <- matchit(treat ~ eng_pre * eng_vorliebe, # beide Kovariaten mit Interaktion
data = englisch, # der Datensatz
distance = "glm", # PS auf Basis der glm Funktion
link = "logit", # mit Logit-Link
method = "subclass", # Stratifizierung
subclass = 10, # zehn Strata
discard = "both", # Ausschluss der Daten außerhalb
# der Region of Common Support
estimand = "ATE") # ATE als ZieleffektWir visualisieren wieder die Propensity-Scores der Personen in den beiden Interventionsgruppen und den einzelnen Strata.
plot(stra10, type = 'jitter', interactive = FALSE)
Abbildung 4: Propensity Scores der gematchten und ungematchten Personen getrennt nach Interventionsgruppe. Senkrechte Linien kennzeichnen die Strata
Insbesondere können wir hier sehen, dass Personen aus der Kontrollgruppe aufgrund mangelndem Overlap ausgeschlossen wurden (Unmatched Control Units), jedoch niemand aus der Interventionsgruppe ausgeschlossen wurde (Unmatched Treated Units).
Wir prüfen nun erneut die Balance.
summary(stra10, subclass=TRUE)##
## Call:
## matchit(formula = treat ~ eng_pre * eng_vorliebe, data = englisch,
## method = "subclass", distance = "glm", link = "logit", estimand = "ATE",
## discard = "both", subclass = 10)
##
## Summary of Balance for All Data:
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.564 0.394 0.910 1.081 0.240
## eng_pre 2.784 2.396 0.787 0.586 0.216
## eng_vorliebe 0.700 0.500 0.417 . 0.200
## eCDF Max
## distance 0.360
## eng_pre 0.355
## eng_vorliebe 0.200
##
## Summary of Balance by Subclass:
##
## - Subclass 1
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.189 0.187 0.012 1.102 0.108
## eng_pre 1.740 1.754 -0.029 0.948 0.102
## eng_vorliebe 0.600 0.346 0.529 . 0.254
## eCDF Max
## distance 0.300
## eng_pre 0.254
## eng_vorliebe 0.254
##
## - Subclass 2
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.262 0.257 0.028 1.07 0.119
## eng_pre 2.098 2.077 0.044 1.14 0.094
## eng_vorliebe 0.333 0.344 -0.022 . 0.010
## eCDF Max
## distance 0.354
## eng_pre 0.260
## eng_vorliebe 0.010
##
## - Subclass 3
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.316 0.316 0.002 0.610 0.100
## eng_pre 2.330 2.308 0.043 0.882 0.079
## eng_vorliebe 0.231 0.290 -0.124 . 0.060
## eCDF Max
## distance 0.233
## eng_pre 0.216
## eng_vorliebe 0.060
##
## - Subclass 4
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.368 0.37 -0.010 1.49 0.100
## eng_pre 2.519 2.42 0.211 0.66 0.127
## eng_vorliebe 0.222 0.50 -0.579 . 0.278
## eCDF Max
## distance 0.209
## eng_pre 0.316
## eng_vorliebe 0.278
##
## - Subclass 5
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.425 0.424 0.004 1.33 0.073
## eng_pre 2.615 2.545 0.142 1.01 0.088
## eng_vorliebe 0.421 0.560 -0.290 . 0.139
## eCDF Max
## distance 0.196
## eng_pre 0.192
## eng_vorliebe 0.139
##
## - Subclass 6
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.482 0.478 0.025 1.17 0.062
## eng_pre 2.726 2.793 -0.136 1.06 0.077
## eng_vorliebe 0.522 0.381 0.294 . 0.141
## eCDF Max
## distance 0.205
## eng_pre 0.193
## eng_vorliebe 0.141
##
## - Subclass 7
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.554 0.554 0.000 1.276 0.064
## eng_pre 2.762 2.876 -0.231 0.954 0.151
## eng_vorliebe 0.773 0.591 0.379 . 0.182
## eCDF Max
## distance 0.182
## eng_pre 0.364
## eng_vorliebe 0.182
##
## - Subclass 8
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.629 0.620 0.050 1.134 0.111
## eng_pre 2.828 2.838 -0.020 0.526 0.105
## eng_vorliebe 0.923 0.889 0.071 . 0.034
## eCDF Max
## distance 0.278
## eng_pre 0.299
## eng_vorliebe 0.034
##
## - Subclass 9
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.723 0.718 0.024 1.04 0.079
## eng_pre 2.998 2.986 0.025 1.08 0.079
## eng_vorliebe 1.000 1.000 0.000 . 0.000
## eCDF Max
## distance 0.247
## eng_pre 0.247
## eng_vorliebe 0.000
##
## - Subclass 10
## Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean
## distance 0.834 0.844 -0.053 0.545 0.103
## eng_pre 3.360 3.413 -0.107 0.545 0.103
## eng_vorliebe 1.000 1.000 0.000 . 0.000
## eCDF Max
## distance 0.248
## eng_pre 0.248
## eng_vorliebe 0.000
##
## Sample Sizes by Subclass:
## 1 2 3 4 5 6 7 8 9 10 Discarded All
## Control 26 32 31 26 25 21 22 18 7 9 15 232
## Treated 5 12 13 18 19 23 22 26 37 35 0 210
## Total 31 44 44 44 44 44 44 44 44 44 15 442love.plot(stra10, # matchit-Objekt
stats = "mean.diffs", # Mittelwertsdifferenzen
disp.subclass = TRUE, # Mittelwertsdifferenzen innerhalb der Strata
thresholds = c(-0.1, 0.1), # Einzeichnen der Grenzwerte nach
# Steiner et al. (2010)
sample.names = c("unadjustiert", # Beschriftung der Ergebnisse
"adjustiert"),
stars = "raw") # Markierung nicht-standardisierter
# Mittelwertsdifferenzen
Abbildung 5: (standardisierte) Mittelwertsdifferenzen der Kovariaten und Propensity-Scores vor (rote Punkte) und nach (türkisfarbene Punkte) der Stratifizierung sowie nach Strata (Zahlen) getrennt. Die durchgezogene senkrechte Linie kennzeichnet eine Mittelwertsdifferenz von 0 und somit gleiche Mittelwerte in der Interventionsgruppe und Kontrollgruppe. Die gestrichelten Linien kennzeichnen die Grenzen für akzeptable Mittelwertsdifferenzen nach Steiner et al. (2010). (* kennzeichnet Variablen für die die unstandardisierten Mittelwertsdifferenzen angegeben sind. Für diese gelten die Grenzwerte nicht.)
love.plot(stra10, # matchit-Objekt
stats = "variance.ratios", # Varianzverhältnisse
disp.subclass = TRUE, # Varianzverhältnisse innerhalb der Strata
thresholds = c(0.8, 1.25), # Einzeichnen der Grenzwerte nach
# Steiner et al. (2010)
sample.names = c("unadjustiert", # Beschriftung der Ergebnisse
"adjustiert")
) 
Abbildung 6: Varianzverhältnisse der (kontinuierlichen) Kovariate und Propensity-Scores vor (rote Punkte) und nach (türkisfarbene Punkte) dem der Stratifizierung sowie nach Strata (Zahlen) getrennt. Die durchgezogene senkrechte Linie kennzeichnet ein Varianzverhältnis von 1 und somit gleiche Varianzen in der Interventionsgruppe und Kontrollgruppe. Die gestrichelten Linien kennzeichnen die Grenzen für akzeptable Varianzverhältnisse nach Steiner et al. (2010).
plot <- bal.plot(stra10, # matchit-Objekt
var.name = "eng_pre", # Variable, für die Kovariatenverteilungen
# dargestellt werden sollen
which = "both", # Kovariatenverteilung vor und nach der Stratifizierung
type = "histogram", # Darstellung als Histogramm
mirror = TRUE, # Spiegelung der gruppenspezifischen Verteilungen an x-Achse
sample.names = "unadjustiert") # Benennung der unadjustierten Ergebnisse
plot
Abbildung 7: Verteilung der kontinuierlichen Kovariate eng_pre in der Kontrollgruppe (rote Balken) und in der Interventionsgruppe (türkisfarbene Balken) vor der Stratifizierung (links) und in den einzelnen Strata nach der Stratifizierung (Str. 1 bis Str. 10). (Hinweis: Um die Titel der Verteilungen innerhalb der Strata (Str. 1 bis Str. 10) zu verändern und die Strata in eine aufsteigende Reihenfolge zu bringen, müssen Änderungen innerhalb der Funktion bal.plot() vorgenommen werden, die hier nicht erläutert werden. Somit weicht die Abbildung die aus der einfachen Reproduktion des Codes resultiert z.T. von der hier dargestellten ab.)
Mit dem neuen Modell liegen nach der Stratifizierung die gewichteten (standardisierten) Mittelwertsunterschiede der Propensity-Scores und beider Kovariaten nahe 0, was anhand von Abbildung 5 (türkisfarbene Punkte) zu erkennen ist. Über die Strata hinweg unterscheiden sich also die Mittelwerte der Kovariaten kaum, nach Steiner et al. (2010) kann diesbezüglich von einer ausreichenden Balance ausgegangen werden. Innerhalb der Strata hingegen liegen die Mittelwertsdifferenzen der Kovariaten zum Teil außerhalb der Grenzen nach Steiner et al. (2010), so zum Beispiel für die Englischvorkenntnisse in Stratum 1, 4, 5, 6 und 7 (siehe Kennwerte in (Summary of Balance by Subclass) oder Stratumkennwerte in Abbildung 5). Richtet man sich nach den von Stuart & Rubin (2008) vorgeschlagenen weniger strengen Grenzen, also Mittelwertsdifferenzen zwischen -0,25 und 0,25, liegen die Mittelwertsdifferenzen der meisten Strata in einem akzeptablen Bereich.
Die über Strata gewichteten Varianzverhältnisse sowohl der Propensity-Scores als auch der Englischvorkenntnisse nach der Stratifizierung (siehe türkisfarbende Punkte in Abbildung 6) sprechen ebenfalls für eine ausreichende Balance. Diese liegen mit 0,983 und 0,884 innerhalb der Grenzen nach Steiner et al. (2010). In den einzelnen Strata hingegen liegen die Varianzverhältnisse zum Teil deutlich außerhalb der Grenzen nach Steiner et al. (2010); für das vierte und zehnte Stratum beispielsweise sowohl bezüglich der Vorkenntnisse als auch bezüglich der Propensity-Scores. Gehen wir von den weniger strengen Grenzen für ausreichende Balance nach Stuart & Rubin (2008) aus, also Varianzverhältnisse zwischen 0,5 und 2, sind die Varianzverhältnisse innerhalb aller Strata akzeptabel.
Abbildung 6 zeigt, dass in den meisten Strata die Verteilung der Kovariate Englischvorkenntnisse zwischen Interventions- und Kontrollgruppe sehr ähnlich ist. In einigen Strata (z.B. Stratum 1) gibt es leichte Abweichungen.
Insgesamt lässt sich also festhalten, dass über die Strata hinweg eine ausreichende Balance gegeben, die Homogenitätsannahme jedoch nicht vollständig erfüllt ist, zumindest dann nicht, wenn die strengen Grenzen nach Steiner et al. (2010) zur Beurteilung herangezogen werden. Bei einer Beurteilung anhand der weniger strengen Kriterien nach Stuart & Rubin (2008) kann von einer akzeptablen Balance auch innerhalb der Strata ausgegangen und der stratifizierte Datensatz zur Schätzung der Effekte verwendet werden.
7 Schätzung kausaler Effekte
Um den kausalen Effekt des Englischtrainings auf die
Englischkenntnisse zu schätzen, müssen wir die Daten zunächst mithilfe
des match.data()-Befehls extrahieren.
stra10_data <- match.data(stra10)Mit diesem Befehl werden auch Stratumgewichte (weights), welche - ähnlich wie bei der PS-Gewichtung auf Personenebene - Gewichte auf Stratumebene, generiert. In diesen wird für die ungleiche Verteilung der Personen in die beiden Interventionsgruppen in jedem Stratum kontrolliert. Die Stratumgewichte werden für die Schätzung des kausalen Effektes genutzt.
Da aufgrund mangelndem Overlaps Personen aus der Stichprobe entfernt wurden, können wir hier nur den durchschnittlichen kausalen Effekt für Personen in der Region of Common Support schätzen. Dazu nutzen wir eine Regression zur Vorhersage der Outcome-Variable durch die Interventionsvariable und gewichten die Beobachtungen entsprechend der Stratumgewichte.
m1 <- lm(eng_post ~ treat, # Regression von Englisch-Nachtest
# auf die Interventionsvariable
data = stra10_data, # der neue Datensatz
weights = weights) # die Stratumgewichte zur Schätzung des ATE
summary(m1) # Ergebnis der Schätzung##
## Call:
## lm(formula = eng_post ~ treat, data = stra10_data, weights = weights)
##
## Weighted Residuals:
## Min 1Q Median 3Q Max
## -1.4501 -0.3797 -0.0503 0.3642 2.2286
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.0055 0.0376 79.89 < 0.0000000000000002 ***
## treat 0.2431 0.0536 4.53 0.0000076 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.554 on 425 degrees of freedom
## Multiple R-squared: 0.0461, Adjusted R-squared: 0.0439
## F-statistic: 20.5 on 1 and 425 DF, p-value: 0.00000759Der geschätzte durchschnittliche kausale Effekt des Englischtrainings
für Personen in der Region of Common Support kann unter
Coefficients abgelesen werden
und liegt bei 0,243 (SE = 0,054). Für Personen
innerhalb der Region of Common Support führt das Englischtraining im
Mittel zu 0,2431 Punkten höheren Englischkenntnissen im Nachtest. Dieser
Effekt ist statistisch signifikant von null verschieden (\(t(1)\) = 4,530; \(p\) < ,001).
Die hier angegebene Inferenzstatistik berücksichtigt nicht, dass auch
die PS-Schätzung aus dem ersten Schritt mit Unsicherheit verbunden ist.
Es existieren Verfahren, die diese Unsicherheit in der Inferenz für die
kausalen Effekte berücksichtigen.
8 Weiterführende Literatur
Weiterführende Literatur zur Umsetzung von PS-Methoden in R finden Sie hier:
- Leite, W. (2016). Practical propensity score methods using R. Sage Publications.
Weitere Informationen zu den Paketen MatchIt und cobalt sind hier zu finden:
Greifer, N. (2022, 03. November). Covariate balance tables and plots: A guide to the cobalt package. Cran R- Project. https://cran.r-project.org/web/packages/cobalt/vignettes/cobalt_A0_basic_use.html
Ho, D. E., Imai, K., King, G., & Stuart, E. A. (2007). Matching as nonparametric preprocessing for reducing model dependence in parametric causal inference. Political analysis, 15(3), 199-236. https://doi.org/10.18637/jss.v042.i08
Ho, D. E., Imai, K., King, G., & Stuart, E. A. (2011). MatchIt: Nonparametric Preprocessing for Parametric Causal Inference. Journal of Statistical Software, 42(8), 1-28. https://kosukeimai.github.io/MatchIt/
Benchmarks für Balance finden Sie beispielsweise hier:
Steiner, P. M., Cook, T. D., Shadish, W. R., & Clark, M. H. (2010). The importance of covariate selection in controlling for selection bias in observational studies. Psychological Methods, 15(3), 250. https://doi.org/10.1037/a001871
Stuart, E. A., & Rubin, D. B. (2008). Best practices in quasi-experimental designs: matching methods for causal inference. In J. Osborne (Hrsg.), Best Practices in Quantitative Methods (155-176). Sage Publications. https://doi.org/10.4135/9781412995627