[[8]{.chapter-number}  [IV-Regression]{.chapter-title}]{#sec-IV .quarto-section-identifier}

8 IV-Regression

Hinweis

Für eine optimale technische Stabilität empfehlen wir, dieses Online-Lehrbuch auf einem Notebook oder Desktop-Computer zu nutzen. Die interaktiven Komponente, insbesondere die R-Konsole (WebR), sind rechenintensiv und funktionieren auf mobilen Geräten nur eingeschränkt.

Warning: package 'gt' was built under R version 4.3.3

Wie in Kapitel Kapitel 4 erläutert, können Regressionsmodelle unter Endogenitätsproblemen leiden: Der Fehlerterm in unserem Regressionsmodell ist mit dem/den interessierenden Regressor(en) korreliert, sodass der zugehörige Koeffizient inkonsistent geschätzt wird und damit nicht als kausaler Effekt interpretiert werden darf. Die häufigsten Ursachen für Endogenität, die in empirischen Anwendungen auftreten, sind ausgelassene Variablen, Messfehler und simultane Kausalität.

Durch Hinzufügen potenziell relevanter Variablen in die Regression kann das Risiko einer verzerrten Schätzung des interessierenden kausalen Effekts anhand multipler Regression verringert werden. Wenn jedoch wichtige ausgelassene Variablen unbeobachtbar (U) sind oder aus anderen Gründen nicht zur Verfügung stehen, kann multiple Regression die Endogenität nicht beheben. Diese Problematik ergibt sich insbesondere bei simultaner Kausalität: Wenn die Kausalität von der Behandlungsvariable (B) zur Outcomevariable (Y) und umgekehrt verläuft, gibt es Backdoors, die nicht durch das Hinzufügen von Kontrollvariablen zum Regressionsmodell geschlossen werden können.

Eine allgemeinere Technik, um Backdoors aufgrund unbeobachtbarer Variablen zu schließen, ist die Regression mit Instrumentvariablen (IV). Bei IV-Regression wird exogene Variation in einer Instrumentvariablen (Z) verwendet, um den Teil der Variation in der Behandlungsvariable zu isolieren, der nicht durch unbeobachtbare Faktoren verursacht wird, die sowohl B als auch Y beeinflussen. Der kausale Effekt von B auf Y wird dann anhand dieser exogene Variation in B geschätzt. Ein Forschungsdesign, in dem der kausale Effekt von B auf Y mit IV-Regression geschätzt werden kann, ist in Abbildung 8.1 dargestellt.

Abbildung 8.1: IV-Design bei Backdoor durch unbeobachtbare Confounder-Variablen

8.1 Der einfache lineare IV-Schätzer

Im folgenden nehmen wir an, dass sämtliche Zusammenhänge linear sind. Zunächst betrachten wir das einfache Regressionsmodell

\[\begin{align} Y_i = \beta_0 + \beta_1 B_i + u_i \ \ , \ \ i=1,\dots,n, \label{eq:simpleiv} \end{align}\]

wobei der Fehlerterm $u_i$ mit dem Regressor $B_i$ korreliert ist (d.h. $B$ ist ein endogener Regressor), sodass der KQ-Schätzer für den kausalen Effekt $\beta_1$ inkonsistent ist.

Damit $Z$ ein gültiges Instrument für $B$ in dem in Abbildung 8.1 gezeigten Forschungsdesign ist, müssen die folgenden Bedingungen erfüllt sein:

Sei $\text{Cov}(A,B)$ die Kovarianz zwischen den Variablen $A$ und $B$. $Z$ muss zwei Bedingungen erfüllen, um ein gültiges Instrument zu sein:

1. Relevanz des Instruments für die endogene Variable

$B$ und $Z$ müssen korreliert sein (Pfeil von Z nach B in Abbildung 8.1): \[\begin{align} \text{Cov}(B,Z) \neq 0 \label{eq:ivassum1} \end{align}\]
2. Exogenität des Instruments hinsichtlich der Outcome-Variable

Das Instrument $Z$ darf nicht mit dem Fehlerterm $u$ in der Modellgleichung $\eqref{eq:simpleiv}$ korreliert sein (keine Pfade von Z nach Y außer durch B in Abbildung 8.1): \[\begin{align} \text{Cov}(Z,u) = 0 \label{eq:ivassum2} \end{align}\]

Unter diesen Annahmen erlaubt das Forschungsdesign die Anwendung des einfachsten IV-Ansatzes, wobei eine endogene Variable $B$ durch eine Instrumentvariable $Z$ instrumentiert wird. Die folgende Umformung zeigt, warum der kausale Effekt von B auf Y anhand der (Ko)Variation in diesen Variablen identifiziert werden kann:

\[\begin{alignat*}{2} \textup{Cov}(Z,Y) &= \textup{Cov}(Z,\beta_0 + \beta_1 B + u) &\quad& \text{(Gl. \eqref{eq:simpleiv})} \\ \\ \textup{Cov}(Z,Y) &= \textup{Cov}(Z, \beta_1 B + u) &\quad& \text{($\beta_0$ konstant)} \\ \\ \textup{Cov}(Z,Y) &= \beta_1\textup{Cov}(Z,B) &\quad& \text{($\beta_1$ konstant, $Z$ exogen)} \\ \\ \beta_1 &= \frac{\textup{Cov}(Z,Y) }{\textup{Cov}(Z,B)} &\quad& \text{($Z$ relevant)} \end{alignat*}\]

Eine naheliegende Implementierung gemäß dieses Identifikationsprinzips ist der einfache IV-Schätzer

\[\begin{align} \widehat{\beta}_{\textup{IV}} = \frac{\widehat{\textup{Cov}}(Z,Y) }{\widehat{\textup{Cov}}(Z,B)}, \label{eq:simpleivest} \end{align}\]

wobei lediglich die Kovarianzfunktion $\textup{Cov}(\cdot,\cdot)$ durch ihr Stichprobenäquivalent $\widehat{\textup{Cov}}(\cdot,\cdot)$ ersetzt wird.

Erwartungswert und Konsistenz

Eine hilfreiche Darstellung von $\eqref{eq:simpleivest}$ ist

\[\begin{align*} \widehat{\beta}_\textup{IV} = \beta_1 + \frac{\sum_{i=1}^n (Z_i-\overline Z) u_i}{\sum_{i=1}^n (Z_i - \overline{Z})B_i}. \end{align*}\]

Anhand dieser Form kann der Erwartungswert sowie das Verhalten für große Stichproben untersucht werden. Eine wichtige Eigenschaft ist

\[\begin{align} \textup{E}\big(\widehat{\beta}_\textup{IV}\big) = \beta_1 + \underbrace{\textup{E}\bigg(\frac{\sum_{i=1}^n (Z_i-\overline Z) u_i}{\sum_{i=1}^n (Z_i - \overline{Z})B_i}\bigg)}_{\neq0},\label{eq:ivbiasterm} \end{align}\]

sodass $\widehat{\beta}_\textup{IV}$ bei Endogenität von $X$ ein verzerrter Schätzer von $\beta_1$ ist.¹ Glücklicherweise kann man zeigen, dass

¹ Beachte, dass der endogene Regressor $B_i$ mit $u_i$ korreliert, sodass der Erwartungswert des Bruchs in $\eqref{eq:ivbiasterm}$ ungleich $0$ ist.

\[\begin{align*} \frac{\sum_{i=1}^n (Z_i-\overline Z) u_i}{\sum_{i=1}^n (Z_i - \overline{Z})B_i}\to 0 \quad \text{für} \quad n\to\infty \end{align*}\]

bei Gültigkeit der IV-Annahmen $\eqref{eq:ivassum1}$ und $\eqref{eq:ivassum2}$, d.h. $\widehat{\beta}_\textup{IV}$ ist ein konsistener Schätzer für $\beta_1$,

\[\begin{align*} \widehat{\beta}_\textup{IV} \to \beta_1 \quad \text{für} \quad n\to\infty. \end{align*}\]

Schwache Instrumente

Beachte, dass die Annahme der Relevanz des Instruments für die Herleitung des Identifikationsprinzips und des Schätzers $\eqref{eq:simpleivest}$ von entscheidender Bedeutung ist: Sind $Z$ und $B$ unkorreliert, so ist $\text{Cov}(B,Z) = 0$ und $\beta_1$ damit nicht identifizierbar. Weiterhin würde $\widehat{\text{Cov}}(B,Z)\to0$ wenn $n\to\infty$, sodass $\widehat{\beta}_{\textup{IV}}\to\infty$!

In empirischen Anwendungen ist $\text{Cov}(B,Z) = 0$ unwahrscheinlich. Die Kovarianz von $B$ und $Z$ kann jedoch gering sein. In diesem Fall bezeichnet man $Z$ als ein schwaches Instrument. Anhand von Gleichung $\eqref{eq:simpleivest}$ können wir die Konsequenzen der Verwendung eines schwachen Instruments ableiten: Ein sehr schwacher linearer Zusammenhang von $B$ und $Z$ — und damit tendenziell $\widehat{\text{Cov}}(B,Z) \approx 0$ — kann trotz Relevanz von $Z$ zu einem Schätzer mit großer Varianz führen.

Die Relevanz eines Instruments kann mit einem Hypothesentests geprüft werden. Wie die “Stärke” eines Instruments beurteilt werden sollte, ist eine nicht-triviale Frage. Wir betrachten diesen Aspekt näher in Kapitel 8.2.

8.1.1 Identifizierbarer Behandlungseffekt

Regression ermöglicht in der Regel die Schätzung eines durchschnittlichen Behandlungseffekts (ATE) für eine betrachtete Population mit heterogenen Behandlungseffekten (in Modell $\eqref{eq:simpleiv}$ unterstellen wir den empirisch unwahrscheinlichen Fall eines einheitlichen kausalen Effekts $\beta_1$). Das Adjustieren für Endogenität mit IV-Regression hat hier einen Preis: Beachte, dass im IV-Design ein Teil der beobachteten Variation der Behandlungsvariable $B$ (wie in $\eqref{eq:simpleiv}$ angenommen) endogen ist und wir den Effekt von $B$ auf $Y$ anhand der exogenen Variation eines validen Instruments $Z$ schätzen. Daher können wir grundsätzlich nur einen lokalen durchschnittlichen Behandlungseffekt (LATE) identifizieren: Der LATE ist eine gewichtete Variante des ATE, wobei Beobachtungen, deren Behandlung gut durch $Z$ erklärt wird den größten Einfluss haben.²

² Siehe Kapitel 10 und 19 in Huntington-Klein (2021) für nicht-technische Erläuterungen und Beispiel zu LATE.

8.2 2SLS-Verfahren

Der einfache IV-Schätzer $\eqref{eq:simpleivest}$ ist ein Spezialfall eines allgemeineren Regressionsverfahrens, das in zwei Schritten Durchgeführt wird. Im fall von Modell $\eqref{eq:simpleiv}$ erfolgt die Schätzung anhand der folgenden einfachen Regressionen:

1. Stufe: Regression von $B$ auf $Z$

Regressiere die endogene Variable $B$ auf das Instrument $Z$: \[\begin{align} B_i = \alpha_0 + \alpha_1 Z_i + u_i, \quad 1,\dots,n. \end{align}\] Berechne die geschätzten Werte \[\begin{align} \widehat{B}_i = \widehat{\alpha}_0 + \widehat{\alpha}_1 Z_i. \end{align}\]

Die $\widehat{B}_i$ aus dieser Regression sind “bereinigt”: Die Variation in den $\widehat{B}_i$ wird lediglich durch die exogene Variation in $Z$ verursacht.
2. Stufe: Regression von $Y$ auf $\widehat{B}$

Regressiere die Beobachtungen der abhängigen Variable $Y_i$ auf die geschätzten Werte $\widehat{B}_i$ aus der 1. Stufe: \[\begin{align} Y_i = \gamma_0 + \gamma_1 \widehat{B}_i + e_i. \end{align}\] Der geschätzte Koeffizient $\widehat{\gamma}_1$ aus dieser Regression ist der Two-Stage-Least-Squares-Schätzer (2SLS-Schätzer) von $\beta_1$.

8.2.1 Äquivalenz von 2SLS und einfacher IV-Schätzer

Der KQ-Schätzer von $\gamma_1$ in der 2. Stufe ist \[\begin{align} \widehat{\gamma}_1 = \frac{\sum (\widehat{B}_i - \overline{\widehat{B}})(Y_i - \bar{Y})}{\sum (\widehat{B}_i - \bar{\widehat{B}})^2} = \frac{\widehat{\text{Cov}}(\widehat{B}, Y)}{\widehat{\textup{Var}}(\widehat{B}_i)}. \end{align}\] Da $\widehat{B}_i = \widehat{\alpha}_0 + \widehat{\alpha}_1 Z_i$, die angepassten Werte aus der 1. Stufe, eine lineare Funktion von $Z_i$ sind, können wir $\widehat{\text{Cov}}(\widehat{B}, Y)$ wie folgt schreiben: \[\begin{align} \widehat{\text{Cov}}(\widehat{B}, Y) = \widehat{\alpha}_1 \widehat{\text{Cov}}(Z, Y) \end{align}\] Eine ähnliche Umformung zeigt, dass \[\begin{align} \widehat{\text{Var}}(\widehat{B}) = \widehat{\alpha}_1^2 \widehat{\text{Var}}(Z). \end{align}\]

Kombinieren wir diese Ergebnisse, erhalten wir folgende Darstellung für den KQ-Schätzer von $\gamma_1$ in der 2. Stufe: \[\begin{align} \widehat{\gamma}_1 = \frac{\widehat{\text{Cov}}(\widehat{B}, Y)}{\widehat{\text{Var}}(\widehat{B})} = \frac{\widehat{\alpha}_1 \widehat{\text{Cov}}(Z, Y)}{\widehat{\alpha}_1^2 \widehat{\text{Var}}(Z)} = \frac{\widehat{\text{Cov}}(Z, Y)}{\widehat{\alpha}_1 \widehat{\text{Var}}(Z)} \end{align}\]

Der KQ-Schätzer von $\widehat{\alpha}_1$ in der 1. Stufe ist definiert als $\widehat{\alpha}_1 = \frac{\widehat{\text{Cov}}(Z, B)}{\widehat{\text{Var}}(Z)}$. Somit gilt \[\begin{align} \widehat{\gamma}_1 = \frac{\widehat{\text{Cov}}(Z, Y)}{\frac{\widehat{\text{Cov}}(Z, B)}{\widehat{\text{Var}}(Z)} \widehat{\text{Var}}(Z)} = \frac{\widehat{\text{Cov}}(Z, Y)}{\widehat{\text{Cov}}(Z, B)} = \widehat{\beta}_{\textup{IV}}.\label{eq:equiviv} \end{align}\]

8.2.2 Wozu 2SLS?

Abbildung 8.2: IV-Design mit beobachtbaren Confoundern

In der Praxis werden IV-Schätzungen mit statistischer Software häufig anhand einer Implementierung des 2SLS-Verfahren durchgeführt. Gründe hierfür sind:

Bedingte Exogenität. Für den in Abbildung 8.1 dargestellten DGP ist die Exogenität von Z gegeben. In empirischen Anwendungen kann es jedoch schwierig sein, ein Instrument Z zu finden, dass plausibel nur mit dem endogenen Regressor korreliert. Eine schwächere Bedingung als Exogenität ist Exogenität nach Kontrolle für beobachtbare gemeinsame Determinanten (X) von Z und Y.³

Diese Situation ist in Abbildung 8.2 dargestellt: Aufgrund der Pfeile von X zu Z und Y ist die Bedingung $\eqref{eq:ivassum2}$ verletzt. Durch Kontrolle für X in beiden Regressionen des 2SLS-Verfahrens können die Pfade durch X geschlossen werden.
Mehrere endogene Variablen / Instrumente. Der 2SLS-Ansatz kann leicht für mehrere endogene Variablen (und mehrere Instrumentvariablen) erweitert werden. Für jede der endogenen Variablen erfolgt hierbei eine separate Regression in der ersten Stufe des 2SLS-Verfahrens, wobei jeweils mehrere Instrumente verwendet werden und zusätzlich für beobachtbare Kovariablen X kontrolliert werden kann.
Statistische Inferenz. Die Berechnung von $\widehat{\beta}_\textup{IV}$ erfordert Sorgfalt bei der Konstruktion von Inferenzstatistiken. Die in Kapitel 8.2.1 gezeigte Äquivalenz impliziert, dass wir bei der Schätzung der Variabilität von $\widehat{\beta}_\textup{IV}$ die Unsicherheit aus den beiden KQ-Schätzungen des 2SLS-Verfahrens berücksichtigen müssen. Korrigierte Standardfehlerformeln können vergleichsweise einfach anhand der Regressionsdarstellung hergeleitet werden und sind statistischer Software implementiert (bspw. ivreg in R).

Weiterhin ermöglicht das Regressions-Framework diagnostische Tests. So kann die Relevanz von Instrumenten anhand von F-Tests für die Koeffizienten in den Regressionen der ersten Stufe des 2SLS-Verfahrens überprüft werden. Sofern mehr Instrumente als endogene Variablen vorliegen, kann die Exogenität der Instrumente getestet werden.

³ Formal: $\textup{Cov}(Z,u\vert X) = 0$.

8.2.3 Interaktive Illustrationen: 2SLS-Verfahren

Der DGP im nachfolgenden interaktiven Beispiel ist

\[\begin{align} \begin{split} X_i =&\, \gamma \cdot Z_i + v_i,\\ Y_i =&\, \beta_1 \cdot X_i + u_i, \end{split} \label{eq:OJSIVDGP} \end{align}\] wobei $Z_i \sim N(0, .25^2)$. Die Fehlerterme $v_i$ und $u_i$ sind $N(0,1)$-verteilt und folgen einer gemeinsam Normalverteilung mit Korrelationsparameter $\rho$, \[\begin{align} \begin{pmatrix} v_i \\ u_i \end{pmatrix} \sim N \bigg[ \begin{pmatrix} 0 \\ 0 \end{pmatrix} \begin{pmatrix} 1 & \rho \\ \rho & 1 \end{pmatrix} \bigg]. \end{align}\]

Wir sind daran interessiert, den Parameter $\beta_1$ zu schätzen, um den Effekt von $X$ auf $Y$ zu bestimmen, wobei wir $n=1000$ Beobachtungen von $(X_i, Y_i, Z_i)$ verwenden.

Beachte, dass $X_i$ eine Funktion der (exogenen) Variable $Z_i$ und des Fehlerterms $v_i$ ist. Wenn die Fehlerterme $(v_i, u_i)’$ im DGP $\eqref{eq:OJSIVDGP}$ korreliert sind ($\rho\neq0$), dann besteht Korrelation zwischen $X_i$ und $u_i$. Dies impliziert, dass $X$ ein endogener Regressor im “naiven” Regressionsmodell \[\begin{align} Y_i = \beta_0 + \beta_1 X_i + e_i, \quad i = 1,\dots,N. \end{align}\] Der KQ-Schätzer von $\beta_1$ ist dann verzerrt und inkonsistent.

Wenn $Z$ Erklärungskraft für $X$ hat, d.h. $\gamma\neq0$ mit $|\gamma|$ nicht zu klein gewählt ist ($Z$ ist relevant), und $Z$ nur über $X$ auf $Y$ wirkt ($Z$ ist exogen in der Gleichung von $Y_i$), dann kann der kausale Effekt von $X$ auf $Y$ anhand einer Instrumentvariablen-Regression mit $Z$ konsistent geschätzt werden.

Die folgende Interaktive Grafik berechnet den 2SLS-Schätzer Analog zu der in Kapitel 8.2 erläuterten Vorgehensweise:

Schätze das Regressionsmodell

\[\begin{align} X_i = \alpha_0 + \alpha_1 \cdot Z_i + \epsilon_i \end{align}\] und berechnen Sie die angepassten Werte $\widehat{X}_i$. Dieser Schritt isoliert die exogene Variation in $X_i$, die auf $Z_i$ zurückzuführen ist.
Schätze den Effekt von $X$ auf $Y$ unter Verwendung von $\widehat{X}_i$, dem exogenen Teil von $X_i$ aus 1., in der Regression

\[\begin{align} Y_i = \delta_0 + \delta_1 \cdot \widehat{X}_i + \varepsilon_i. \end{align}\] Der KQ-Schätzer $\widehat{\delta}_1$ von $\delta_1$ ist der 2SLS-Schätzer von $\beta_1$.

Die Grafik illustriert das Verhalten des IV-Schätzers (lila Grade) im Vergleich zum KQ-Schätzer (gestrichelte Grade) und zum wahren Zusammenhang (grüne Grade) für eine simulierte Stichprobe (gegeben der gewählten Parameter) anhand der geschätzten Regressionslinien. Die voreingestellten Parameter-Kombinationen zeigen

die Verzerrung von 2SLS, wenn $Z$ ein schwaches Instrument für $X$ ist,
eine Situation in der $X$ exogen ist (KQ ist unverzerrt), was zu vergleichbaren Schätzungen führt (sofern $Z$ kein schwaches Instrument ist)
eine starke Verzerrung von KQ, wenn $X$ endogen und irrelevant ($\beta_1 = 0$) ist.

Anhand der obigen interaktiven Grafik lässt sich zwar die Verzerrung der Schätzer einschätzen, jedoch nicht die Unsicherheit der Schätzung. Um die Verteilung der Schätzer in endlichen Stichproben zu illustrieren, zeigt die nachstehende Grafik Kerndichteschätzungen basierend auf $N$ simulierte Stichproben der Größe $n$ für den gewählten DGP.⁴ Die Applikation gibt weiterhin den R-Code für eine Simulation der Daten gemäß der gewählten Parameter aus.

⁴ Beachte, dass die Simulation insbesondere für große $N$ und $n$ einige Sekunden dauern kann (grauer Balken am linken Rand zeigt laufende Berechnung an).

Die voreingestellten Parameter-Kombinationen zeigen:

Wenn $Z$ ein schwaches Instrument ist, kann der IV-Schätzer auch in großen Stichproben eine hohe Verzerrung haben und eine deutlich größere Varianz als der KQ-Schätzer haben.
Wenn der interessierende Regressor (Behandlungsvariable) $X$ exogen ist, hat der IV-Schätzer eine größere Varianz als der KQ-Schätzer.
Wenn die Endogenität von $X$ “stark” ist (d.h. hohe Korrelation von $u$ und $v$), kann der IV-Schätzer auch in kleinen Stichproben hilfreich sein, sofern das Instrument $Z$ nicht zu schwach ist: IV hat dann eine deutlich geringere Verzerrung als KQ, aber eine größere Varianz.

8.2.4 Schätzung mit R

Wir zeigen nachfolgend, wie die Daten aus den vorherigen interaktiven Beispielen in R generiert, die IV-Regression durchgeführt und ein (grafischer) Vergleich mit der naiven KQ-Schätzung implementiert werden kann.

8.3 Beispiel: Der schwache Staat, die Bauern und die Mafia

Acemoglu, De Feo, und De Luca (2020) untersuchen den Aufstieg und die Auswirkungen der sizilianischen Mafia gegen Ende des 19. Jahrhunderts. Die Studie argumentiert, dass die Verbreitung der Mafia teilweise auf das Aufkommen sozialistischer Bauern-Organisationen (Fasci siciliani) zurückzuführen ist, die in einem Umfeld schwacher staatlicher Präsenz die Rechte der Bauern verteidigten. Grundbesitzer, Verwalter und lokale Politiker wandten sich also an die Mafia, um diese sozialistische “Bedrohung” zu bekämpfen und die Kontrolle über die landwirtschaftlich arbeitende Bevölkerung zu behalten. Die Studie identifiziert einen signifikanten (positiven) kausalen Zusammenhang zwischen dem Aufstieg der Fasci und der Ausbreitung der Mafia in Sizilien.

8.3.1 Identifikationsstrategie

Die Beziehung zwischen dem Aufkommen der Fasci und der Mafiaaktivität im einfachen Regressionsmodell \[\begin{align} \textup{Mafiaaktivität} = \beta_0 + \beta_1 \textup{Fasci siciliani} + \varepsilon \label{eq:fascimafiamod1} \end{align}\] ist wahrscheinlich aufgrund mehrerer Faktoren endogen:

Externe Faktoren, die sowohl die Entstehung der Fasci siciliani als auch der Mafia beeinflussen, sind wahrscheinlich. Beispielsweise könnten wirtschaftliche Notlagen und politische Instabilität sowohl die Organisation der Bauern als auch die Etablierung der Mafia begünstigen.
Simultante Kasalität: Die Präsenz der Fasci könnte die Mafiaaktivitäten beeinflussen und umgekehrt: Intensive Mafiaaktivitäten könnten die Notwendigkeit für Bauernorganisationen erhöhen, was wiederum die Mafia stärkt. Außerdem erschwert die zeitliche Nähe zwischen der Entstehung der Fasci und der Verbreitung der Mafia die Identifizierung der Kausalrichtung.

Wir können das Modell $\eqref{eq:fascimafiamod1}$ um (messbare) unter 1. fallende Regressoren erweitern, um einer verzerrten Schätzung aufgrund relevanter ausgelassener Variablen vorzubeugen. Eine Verzerrung aufgrund der in 2. beschriebenen simultanen Kausalität kann jedoch nicht durch multiple Regression vermieden werden. Acemoglu, De Feo, und De Luca (2020) nutzen die folgende Identifikationsstrategie: Eine extreme Dürre im Jahr 1893 (insbesondere Regenausfall im Frühling) hatte einen erheblichen negativen Einfluss auf die landwirtschaftliche Produktion in Sizilien und verursachte große Not unter den Bauern, was die Entstehung der Fasci beförderte. Die Dürre kann als natürliches Experiment betrachtet werden, wobei durch die Variation in der Niederschlagsmenge (das Instrument) die Variation im Aufkommen von Fasci-Organisationen in verschiedenen Gemeinden, die unabhängig von der späteren Mafiaaktivität ist, isoliert werden kann. Anhand dieser exogenen Variation im Organisationsgrad der Bauern kann der (lokale) kausale Effekt auf die Entwicklung der Mafia identifiziert werden.

Abbildung 8.3: Zulässiger DGP im IV-Design von Acemoglu, De Feo, und De Luca (2020)

8.3.2 Reproduktion der Kernergebnisse mit R

Wir lesen zunächst den Datensatz ADD_Mafia_municipality.csv⁵ ein und verschaffen uns einen Überblick.

⁵ Der Datensatz stammt aus dem MIT Economics Data Archive und liegt im STATA-Format .dta vor. Wir verwenden hier eine in das .csv formatierte Version.

Tabelle 8.1 gibt einen Überblick der in ADD_dat verfügbaren Variablen auf Gemeinde-Ebene und deren Definition.

Variable	Beschreibung
peasants_fasci	Präsenz Fasci in 1893-Q4 (Dummy)
sp3m1893_n30	Dürre-Intensität: Relative Regenmenge im Frühling 1893 (ggü. Mittel v. 1841 bis 1881)
cl1_stn_sp1893_n30	Nächste aktive Wetterstation (Frühling 1893)
predr_peas_fasci	Präsenz Fasci vor März 1893 (Dummy)
ruralcentre1861	Landwirtschaftliche Gemeinde? (Dummy)
Rural_rent	Durchschn. ländliche Pacht / Hektar in 1953
Urban_rent	Durchschn. städtische Pacht / Hektar in 1953
agricola_rel	Anteil bewirtschaftetes Land in 1853
seminatoritot_rel	Anteil Anbaufläche Getreide in 1853
sulfurproduction1868_70	Durschn. Schwefel-Produktion / Jahr für 1868-1870
Citrus_groves	Anteil Anbaufläche Zitrusfrüchte in 1853
Olives_groves	Anteil Anbaufläche Oliven in 1853
Vineyards	Anteil Anbaufläche Wein in 1853
Mafia1885	Intensität Mafia gem. Damiani (1885), Skala 0-3
Mafia1900	Intensität Mafia gem. Cutrera (1900), Skala 0-3
lnpop1861	Log(Einwohner) in 1861
lnsurface	Log(Fläche) in 1853
centreheight	Höhe des Gemeindezentrums (M.ü. NN)
maxheight	Max. Höhe (M.ü. NN)
slope2	Durschn. Höhe (M.ü. NN)
pa_pdist1856	Entfernung (Postweg) von Palermo in 1856
port2_pdist1856	Entfernung (Postweg) nächer Hafen in 1856
roads1799	Direkter Zugang Postweg in 1799 (Dummy)
ave_temp	Durschn. Jahrestemperatur
sp3m_ave_n30	Durschn. Regenmenge im Frühling für 1881-1941
var_sp3m_n30	Varianz Regenmenge im Frühling für 1881-1941
distretto1853	Distrikt-Zugehörigkeit in 1853
provincia1853	Provinz-Zugehörigkeit in 1853

Tabelle 8.1: ADD_dat — Charakteristika von sizilianischen Gemeinden im 19. Jahrhundert

Vor der Analyse filtern wir die Daten entsprechend der Vorgehensweise in Acemoglu, De Feo, und De Luca (2020): Es werden lediglich Gemeinden berücksichtigt, die in maximal 30km Entfernung zu einer Wetter-Station liegen (!is.na(sp3m1893_n30)) und für die eine Messung der Mafiaaktivität gegen Ende des 19. Jahrhunderts (gemessen auf einer Skala von 0 bis 4, vgl. Cutrera (1900)) vorliegt (!is.na(Mafia1900)).

8.3.2.1 First-Stage-Regressionen

Tabelle 3 in Acemoglu, De Feo, und De Luca (2020) präsentiert First-Stage-Regressionen Für den endogenen Regressor $\textup{Fasci siciliani}$ in $\eqref{eq:fascimafiamod1}$: Die Auswirkungen der Dürre in 1893 (das Instrument) auf den Organisationsgrad der Bauern (der endogene Regressor). Hierbei werden verschiedene, aus dem DAG in Abbildung 8.3 abgeleitete Spezifikationen betrachtet:

Panel A: Keine Fixed Effects
- Modell ohne Kovariablen
- Modell mit Kovariablen für Bildung von Bauernorganisationen
- Modell mit Kovariablen für Bildung von Bauernorganisationen und Indikatoren für Mafia-Präsenz
- Modell mit Kovariablen für Bildung von Bauernorganisationen, Indikatoren für Mafia-Präsenz und geologische Faktoren
Panel B: Modelle wie in Panel A und Kontrolle für Provinz-Fixed-Effects

Wir schätzten nachfolgend sämtliche Regressionen mit feols und berechnen nach Distrikt-Zugehörigkeit (distretto1853) sowie nächstgelegener Wetterstation (cl1_stn_sp1893_n30) geclusterte Standardfehler.⁶ Hierzu setzen wir vcov = ~ distretto1853 + cl1_stn_sp1893_n30. Diese Spezifikation berücksichtigt, dass die Niederschlagsmenge in einer Gemeinde sowie viele Kovariablen eine Korrelation mit benachbarten Gemeinden aufweisen können. Cluster-robuste Standardfehler sind inbesondere nötig, weil die Niederschlagsdaten für viele Gemeinden aus benachbarten Wetterstationen interpoliert werden, was zu einer zusätzlichen Korrelation zwischen den Beobachtungen führt. Standardfehler mit Two-way-clustering sind eine gängige Methode zur Berücksichtigung derartiger Korrelation.

⁶ Siehe Kapitel Kapitel 7 für Erläuterungen zu Fixed-Effects und cluster-robusten Standardfehlern.

Der obige Code-Chunk führt die First-Stage-Regression ohne Kovariablen aus. Wir finden einen signifikanten positiven Effekt der Niederschlagsmenge auf den Organisationsgrad der Bauern.

Acemoglu, De Feo, und De Luca (2020) verwenden Bootstrap cluster-robuste Standardfehler (Cameron, Gelbach, und Miller 2008, 2011), die aufwendig zu implementieren sind, jedoch zu nahezu identischen Ergebnissen führen. Der nächste Code-Chunk zeigt, wie cluster-robuste Bootstrap-Inferenz für das geschätzte Modell fasci_mod berechnet werden können. Hierfür nutzen wir fwildclusterboot::boottest().⁷

⁷ fwildclusterboot implementiert Varianten des Wild Bootstrap (Liu 1988).

Die statistische Zusammenfassung umfasst die Punktschätzung, die robuste t-Statistik sowie das 95%-Konfidenzintervall für den Koeffizienten von sp3m1893_n30. Mit diesen Komponenten können wir den Bootstrap-Standardfehler beispielsweise durch die folgende Transformation berechnen:

Bei der Interpretation des geschätzten Koeffizienten müssen wir berückstigen, dass die abhängige Variable binär ist: Die First-Stage-Regression ist ein lineares Wahrscheinlicheitsmodell. Eine positive Änderung der relativen Regenmenge um eine Einheit (100%) verringert die Wahrscheinlichkeit also um 100%. Diese Interpretation ist jedoch wenig aussagekräftig für einen Einordnung des Effekts der Dürre. Stattdessen betrachten wir die Verteilung relativen Regenmengen und vergleichen den Effekt für eine von starker Dürre betroffenen Gemeinden (25%-Quantil) und einer Gemeinde mit einem relativ geringen Rückgang des Niederschlags (75%-Quantil)

Für eine Gemeinde mit 75% des langjährigen mittleren Niederschlags ist die geschätzte Wahrscheinlichkeit der Bildung von Bauern-Organisationen nur etwa halb so groß, wie für eine von starker Dürre betroffene Gemeinde.

Zur Vermeidung von Backdoors durch ausgelassene Variablen betrachten Acemoglu, De Feo, und De Luca (2020) weitere Spezifikationen für die First- und Second-Stage-Regressionen, die zusätzlich für Determinanten der Fasci, der Mafia, sowie für geographische Faktoren kontrollieren. Weiterhin werden diese Regressionen mit Fixed Effekts auf der Provinz-Ebene gemäß der Grenzen von 1853 (provincia1853) wiederholt, um für unbeobachtbare Heteogenitäten zwischen den verschiedenen Provinzen zu kontrollieren. Die nachfolgenden Code-Chunks zeigen die Implementierung dieser Regressionen mit fixtest::feols(). Für eine bessere Übersicht der Ergebnisse verwenden wir in summary() das Argument n = 1, sodass jeweils nur die Schätzung des ersten Koeffizienten (Effekt der relativen Regenmenge) ausgegeben wird.

Für eine tabellarische Übersicht fassen wir die interessierenden Komponenten der Schätzungen mit modelsummary::modelsummary() zusammen. Hierfür sammeln wir zunächst jeweils die Modelle mit und ohne Provinz-Fixed-Effects in benannten Listen, sodass die Ergebnisse später in separaten Panelen dargestellt werden können. Für die Formatierung von Labels für die Fixed Effects und Koeffizienten definieren wir mit tibble::tribble() das tibble-Objekt gm und übergeben reguläre Ausdrücke für eine Auswahl an zusammenfassenden Statistiken.

Die Ergebnisse bestätigen die Vermutung der Autoren, dass der Effekt des relativen Niederschlags auf den Organisationsgrad der Bauern ohne Kontrolle für die betrachteten Kovariablen tendenziell überschätzt wird. Tatsächlich führt multiple Regression mit Fixed Effects und sämtlichen Kovariablen (letzte Spalte in Panel B) zu einem etwa 30% kleineren geschätzten Effekt als für die einfache Regression in Panel A.

8.3.2.2 Second-Stage-Regressionen

Mit fml = Mafia1900 ~ 1 | peasants_fasci ~ sp3m1893_n30 legen wir fest, dass Mafia1900 auf eine Konstante sowie die angepassten Werte des endogenen Regressors peasants_fasci aus der ersten Stufe (instrumentiert durch sp3m1893_n30) regressiert werden soll. Wie für die First-Stage-Regressionen schätzen wir jeweils vier verschiedene Spezfikationen, unterscheiden zwischen Modellen mit und ohne Provinz-Fixed-Effects und berechnen cluster-robuste Standardfehler.

Die Ergebnisse für die IV-Schätzung des kausalen Effekts von Bauernorganisationen auf die Präsenz der Mafia zu Beginn des 20. Jahrhundert sind in der obigen Tabelle dargestellt. Die LATE-Schätzung in Spalte 4 von Panel B hat folgende Interpretation: Die Präsenz der Fasci in einer Gemeinde erhöht die Intensität der Mafiaaktivität um etwa 1.56 Punkte.⁸

⁸ Nach Cutrera (1900) entspricht eine Erhöhung des Mafia-Index von 1 auf 2 einem Sprung von einer geringen Präsenz der Mafia zu einer signifikanten Präsenz.

Wie ist dieser Effekt einzuordnen? Wir berechnen zunächst die Anzahl an Gemeinden, in denen es im Jahr 1900 Fasci gibt.

Für diese 84 Gemeinden erwarten wir jeweils einen Effekt von 1.56 Punkten auf die im Jahr 1900 beobachtete gesamte Intensität des Einflusses der Mafia in Sizilien, gemessen als gewichtete Summe.

Die gesamte Mafia-Intensität für Sizilien zu Beginn des 20. Jahrhunderts berechnen wir analog.

Der Anteil der Gemeinden mit Bauernorganisationen im Jahr 1900 an dieser Summe erlaubt eine Einordnung des Effekts der Fasci auf die Verbreitung der Mafia.

Die Interpretation lautet, dass bis zu 38% der Stärke der Mafia im Jahr 1900 in Sizilien auf ihren Einsatz gegen die Fasci zurückzuführen sein könnte.

8.4 Case Study: Ökonomische Schocks und Bürgerkriege

Ein Kernproblem bei der Schätzung der kausalen Beziehung zwischen wirtschaftlichen Schocks und dem Auftreten von kriegerischen Auseinandersetzungen ist die Endogenität von Variablen, die wirtschaftliche Stabilität messen (vgl. Fearon und Laitin 2003). Beispielsweise kann eine Verschlechterung der wirtschaftlichen Bedingungen eines Landes die Wahrscheinlichkeit eines Konflikts im Inland erhöhen. Gleichermaßen beeinflussen kriegerische Handlung im Inland die wirtschaftliche Situation einer Volkswirtschaft negativ.

Abbildung 8.4: Simultane Kausalität zw. Kriegsrisiko und Wirtschaftslage

In der “Konfliktgleichung”

\[\begin{align} \textup{Konflikt-Wsk.} = \beta_0 + \beta_1\, \Delta\text{BIP} + \text{Kontrollv.} + \epsilon\label{eq:rainconf1} \end{align}\]

liegt dann Endogenität des Regressors für “wirtschaftliche Stabilität”, hier $\Delta\text{BIP}$, aufgrund von simultanter Kausalität vor. Die multiple Regression $\eqref{eq:rainconf1}$ liefert dann eine verzerrte Schätzungen des interessierenden Effekts $\beta_1$.

Miguel, Satyanath, und Sergenti (2004) untersuchen die kausale Beziehung zwischen wirtschaftlichen Schocks und dem Auftreten von Bürgerkriegen in 41 afrikanischen Sub-Sahara-Ländern für den Zeitraum von 1981 bis 1999.⁹ Hierbei adressieren die Autoren das Endogenitätsproblem in einem IV-Ansatz mit Instrumenten, die mit wirtschaftlichen Schocks korreliert sind, aber nicht unmittelbar mit der Wahrscheinlichkeit eines Bürgerkriegs zusammenhängen.

⁹ Füge rainconflict$country %>% unique() im nächsten WebR-Chunk für eine Liste der Länder aus.

Die Identifikationsstrategie basiert auf metereologischen Faktoren: In den betrachteten Agrarökonomien sind Schocks in den Niederschlägen entscheidend für die Entwicklung des Bruttoinlandsprodukts. Ungünstige Wetterbedingungen, die landwirtschaftliche Erträge negativ beeinflussen, haben einen direkten (negativen) Effekt auf die wirtschaftliche Situation. Variation in der Regenmenge ist jedoch keine unmittelbare Determinante der Konfliktwahrscheinlichkeit und folglich exogen in Modellgleichung $\eqref{eq:rainconf1}$. Maße für die Veränderungen in der Niederschlagsmenge sind somit plausible Instrumente für endogene Variablen, die Veränderungen im Bruttoinlandsprodukt messen.

Abbildung 8.5: IV-Design bei simultaner Kausalität zw. Kriegsrisiko und Wirtschaftslage

Anhand von Länder-Paneldaten mit Informationen über wirtschaftliche Indikatoren, Konfliktereignisse und Niederschläge, isolieren Miguel, Satyanath, und Sergenti (2004) in einem IV-Ansatz die exogene Variation der wirtschaftlichen Entwicklung, d.h. Variation in $\Delta\text{BIP}$ die nicht mit dem Fehlerterm in der Konfliktgleichung korreliert ist. Die Ergebnisse ihrer 2SLS-Schätzungen zeigen, dass negative ökonomische Schocks die Wahrscheinlichkeit eines Bürgerkriegs signifikant erhöhen. Miguel, Satyanath, und Sergenti (2004) zeigen weiterhin, dass ihre Ergebnisse robust gegenüber verschiedenen Modell-Spezifikationen unter Berücksichtigung diverser Kontrollvariablen sind.

Wir reproduzieren nachfolgend die zentralen Ergebnisse von Miguel, Satyanath, und Sergenti (2004) anhand eines Auszugs (rainconflict.csv) aus dem der Studie zugrundeliegenden Datensatz. Der vollständige Datensatz ist, nach Registrierung, hier verfügbar.

Variable	Beschreibung
any_prior	Dummyvariable: Konflikt mit mind. 25 Toten im Jahr t
war_prio	Dummyvariable: Konflikt mit mind. 1000 Toten im Jahr t
ccode	Länderkennung (numerisch)
country	Name des Landes
gdp_g	(BIP_t - BIP_t-1) / BIP_t-1
gdp_g_l	(BIP_t-1 - BIP_t-2) / BIP_t-2
GPCP_g	(Niederschlag_t - Niederschlag_t-1) / Niederschlag_t-1
GPCP_g_l	GPCP_g in t-1
GPCP_g_fl	GPCP_g in t+1
gdp_1979	Log Pro-Kopf-BIP in 1979
polity2_IV	Diff. zw. Demokratie- und Autokrarie-Score im Jahr t (Skala v. -10 bis 10)
polity2l	polity2_IV im Jahr t-1 (Skala v. -10 bis 10)
ethfrac	Ethno-linguistische Fragmentierung (Anteil)
relfrac	Religiöse Fragmentierung (Anteil)
oil	Dummy für Öl-exportierendes Land
lmtnest	Log Anteil Gebirgsregionen an Landesfläche
lpopl1	Log Bevölkerung im Jahr t-1
tot_100_g	Handelsbedingungen (Index, 1995 = 100)
year	Jahr der Beobachtung

Tabelle 8.2: rainconflict – Polit-ökonomische Charakteristika afrikanischer Länder v. 1981 bis 1999


Attaching package: 'dplyr'

The following objects are masked from 'package:stats':

    filter, lag

The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union

ℹ Using "','" as decimal and "'.'" as grouping mark. Use `read_delim()` for more control.

Rows: 743 Columns: 19
── Column specification ────────────────────────────────────────────────────────
Delimiter: ";"
chr  (1): country
dbl (18): any_prio, war_prio, ccode, gdp_g, gdp_g_l, GPCP_g, GPCP_g_l, GPCP_...

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

Mit dem Paket marginaleffects können wir durchschnittliche marginale Einflüsse von Regressoren in einem nicht-linearen Modell (hier ein Probit-Modell) berechnen. Für ein Probit-Modell lässt sich der marginale Effekt eines Regressors $X_j$ auf die Wahrscheinlichkeit, dass das Ereignis eintritt, wie folgt darstellen:

Die Wahrscheinlichkeitsfunktion des Probit-Modells ist definiert als \[\begin{align*} P(Y = 1 | \boldsymbol{x}) = \Phi(\boldsymbol{x}'\boldsymbol\beta), \end{align*}\] wobei $\Phi$ die kumulative Verteilungsfunktion der Standardnormalverteilung ist, $\boldsymbol{x}$ der Vektor der Regressoren und $\boldsymbol{\beta}$ der Vektor der Koeffizienten.

Der marginale Effekt für den Regressor $j$ ist dann die Ableitung von $P(Y = 1 | \boldsymbol{x})$ bezüglich $x_j$, \[\begin{align*} \frac{\partial P(Y = 1 | \boldsymbol{x})}{\partial x_j} = \phi(\boldsymbol{x}'\boldsymbol{\beta}) \cdot \beta_j, \end{align*}\] siehe auch Kapitel 4.2.2. Hierbei ist $\phi(\cdot)$ die Dichtefunktion der Standardnormalverteilung und $\beta_j$ der Koeffizient des Regressors $x_j$.

Die Funktion marginaleffects::avg_slopes() berechnet diesen Effekt anhand der geschätzten Koeffizienten $\widehat{\boldsymbol{\beta}}$ für jede Beobachtung $i=1,\dots,N$ und bestimmt den durchschnittlichen Effekt als das arithmetische Mittel über alle marginalen Effekte: \[\begin{align*} \textup{Average Marginal Effect}(j) = \frac{1}{N} \sum_{i=1}^{N} \phi(\boldsymbol{x}_i' \widehat{\boldsymbol{\beta}}) \cdot \widehat\beta_j, \end{align*}\]

wobei $N$ die Anzahl der Beobachtungen ist und $\boldsymbol{x}_i$ den Wert von $\boldsymbol{x}$ für die $i$-te Beobachtung darstellt.

# Marginale Effekte berechnen
library(marginaleffects)

Warning: package 'marginaleffects' was built under R version 4.3.3

avg_slopes(mod_conf_probit)


     Term          Contrast  Estimate Std. Error      z Pr(>|z|)   S    2.5 %
 ethfrac  mean(dY/dX)        0.220343    0.24242  0.909   0.3634 1.5 -0.25478
 gdp_1979 mean(dY/dX)       -0.062757    0.05833 -1.076   0.2820 1.8 -0.17709
 gdp_g    mean(dY/dX)       -0.350169    0.30060 -1.165   0.2441 2.0 -0.93933
 gdp_g_l  mean(dY/dX)       -0.127305    0.25061 -0.508   0.6115 0.7 -0.61849
 lmtnest  mean(dY/dX)        0.072189    0.03610  2.000   0.0455 4.5  0.00144
 lpopl1   mean(dY/dX)        0.074428    0.04864  1.530   0.1260 3.0 -0.02091
 oil      mean(1) - mean(0)  0.015333    0.21000  0.073   0.9418 0.1 -0.39626
 polity2l mean(dY/dX)        0.000686    0.00497  0.138   0.8902 0.2 -0.00906
 relfrac  mean(dY/dX)       -0.270844    0.24133 -1.122   0.2617 1.9 -0.74384
 year     mean(dY/dX)        0.002769    0.00646  0.428   0.6684 0.6 -0.00990
 97.5 %
 0.6955
 0.0516
 0.2390
 0.3639
 0.1429
 0.1698
 0.4269
 0.0104
 0.2022
 0.0154

Type:  response 
Columns: term, contrast, estimate, std.error, statistic, p.value, s.value, conf.low, conf.high, predicted_lo, predicted_hi, predicted

--- format: live-html engine: knitr webr: packages: [ 'AER', 'dplyr', 'fwildclusterboot', 'faux', 'fixest', 'ggplot2', 'gt', 'ivreg', 'marginaleffects', 'modelsummary', 'readr', 'tidyr' ] cell-options: warning: false message: false repos: - https://rstudio.r-universe.dev --- {{< include ./_extensions/r-wasm/live/_knitr.qmd >}} # IV-Regression {#sec-IV} ```{webr} #| echo: false #| output: false #| edit: false #| message: false #| warning: false library(dplyr) library(gt) # create dataset directory dir.create("datasets") # Download the datasets download.file( "https://raw.githubusercontent.com/mca91/kausal_data/refs/heads/main/ADD_Mafia_municipality.csv", 'datasets/ADD_Mafia_municipality.csv' ) download.file( "https://raw.githubusercontent.com/mca91/kausal_data/refs/heads/main/rainconflict.csv", 'datasets/rainconflict.csv' ) # Formatierung von gt-Tabellen tabopts <- function(x) { fmt_number(x, decimals = 3, drop_trailing_zeros = T) %>% tab_options( # general table.width = "80%", # fonts table.font.color = "black", table.font.size = 16, column_labels.font.weight = "bold", # table header column_labels.padding = 15, column_labels.border.bottom.color = "black", column_labels.border.bottom.width = "1px", column_labels.border.top.color = "black", # group labels row_group.border.top.width = 0, row_group.border.bottom.width = 0, row_group.font.weight = "bold", row_group.padding = 20, # body table_body.hlines.width = 0, table_body.border.bottom.color = "black", table_body.border.bottom.width = "1px", table_body.border.top.width = "0px", table.border.bottom.color = "black" ) } ``` ```{r, echo=F, message=FALSE} library(gt) # Formatierung von gt-Tabellen tabopts <- function(x) { fmt_number(x, decimals = 3, drop_trailing_zeros = T) %>% tab_options( # general table.width = "80%", # fonts table.font.color = "black", table.font.size = 16, column_labels.font.weight = "bold", # table header column_labels.padding = 15, column_labels.border.bottom.color = "black", column_labels.border.bottom.width = "1px", column_labels.border.top.color = "black", # group labels row_group.border.top.width = 0, row_group.border.bottom.width = 0, row_group.font.weight = "bold", row_group.padding = 20, # body table_body.hlines.width = 0, table_body.border.bottom.color = "black", table_body.border.bottom.width = "1px", table_body.border.top.width = "0px", table.border.bottom.color = "black" ) } # ADD Datensatz einlesen ADD_dat <- readr::read_csv( file = "datasets/ADD_Mafia_municipality.csv" ) ``` Wie in Kapitel @sec-regression erläutert, können Regressionsmodelle unter Endogenitätsproblemen leiden: Der Fehlerterm in unserem Regressionsmodell ist mit dem/den interessierenden Regressor(en) korreliert, sodass der zugehörige Koeffizient inkonsistent geschätzt wird und damit nicht als kausaler Effekt interpretiert werden darf. Die häufigsten Ursachen für Endogenität, die in empirischen Anwendungen auftreten, sind ausgelassene Variablen, Messfehler und simultane Kausalität. Durch Hinzufügen potenziell relevanter Variablen in die Regression kann das Risiko einer verzerrten Schätzung des interessierenden kausalen Effekts anhand multipler Regression verringert werden. Wenn jedoch wichtige ausgelassene Variablen unbeobachtbar (U) sind oder aus anderen Gründen nicht zur Verfügung stehen, kann multiple Regression die Endogenität nicht beheben. Diese Problematik ergibt sich insbesondere bei simultaner Kausalität: Wenn die Kausalität von der Behandlungsvariable (B) zur Outcomevariable (Y) und umgekehrt verläuft, gibt es Backdoors, die nicht durch das Hinzufügen von Kontrollvariablen zum Regressionsmodell geschlossen werden können. Eine allgemeinere Technik, um Backdoors aufgrund unbeobachtbarer Variablen zu schließen, ist die Regression mit Instrumentvariablen (IV). Bei IV-Regression wird exogene Variation in einer Instrumentvariablen (Z) verwendet, um den Teil der Variation in der Behandlungsvariable zu isolieren, der *nicht* durch unbeobachtbare Faktoren verursacht wird, die sowohl B als auch Y beeinflussen. Der kausale Effekt von B auf Y wird dann anhand dieser exogene Variation in B geschätzt. Ein Forschungsdesign, in dem der kausale Effekt von B auf Y mit IV-Regression geschätzt werden kann, ist in @fig-IVDAG1 dargestellt. ```{dot} //| fig-width: 6 //| fig-height: 4 //| fig-align: 'center' //| fig-cap: "IV-Design bei Backdoor durch unbeobachtbare Confounder-Variablen" //| label: "fig-IVDAG1" digraph IV_DAG { layout=neato fontname="Helvetica,Arial,sans-serif" node [fontname="Helvetica,Arial,sans-serif", shape=circle] edge [fontname="Helvetica,Arial,sans-serif"] U [pos="1.5,1.5!", color=gray, fontcolor=gray]; B [pos="0,0!"]; Y [pos="3,0!"]; Z [pos="-2,0!"]; U -> B [color=gray]; U -> Y [color=gray]; B -> Y; Z -> B; } ``` ## Der einfache lineare IV-Schätzer Im folgenden nehmen wir an, dass sämtliche Zusammenhänge linear sind. Zunächst betrachten wir das einfache Regressionsmodell \begin{align} Y_i = \beta_0 + \beta_1 B_i + u_i \ \ , \ \ i=1,\dots,n, \label{eq:simpleiv} \end{align} wobei der Fehlerterm $u_i$ mit dem Regressor $B_i$ korreliert ist (d.h. $B$ ist ein *endogener* Regressor), sodass der KQ-Schätzer für den kausalen Effekt $\beta_1$ inkonsistent ist. Damit $Z$ ein gültiges Instrument für $B$ in dem in @fig-IVDAG1 gezeigten Forschungsdesign ist, müssen die folgenden Bedingungen erfüllt sein: Sei $\text{Cov}(A,B)$ die Kovarianz zwischen den Variablen $A$ und $B$. $Z$ muss zwei Bedingungen erfüllen, um ein gültiges Instrument zu sein: - **1. Relevanz des Instruments für die endogene Variable** $B$ und $Z$ *müssen* korreliert sein (Pfeil von Z nach B in @fig-IVDAG1): \begin{align} \text{Cov}(B,Z) \neq 0 \label{eq:ivassum1} \end{align} - **2. Exogenität des Instruments hinsichtlich der Outcome-Variable** Das Instrument $Z$ *darf nicht* mit dem Fehlerterm $u$ in der Modellgleichung \eqref{eq:simpleiv} korreliert sein (keine Pfade von Z nach Y außer durch B in @fig-IVDAG1): \begin{align} \text{Cov}(Z,u) = 0 \label{eq:ivassum2} \end{align} Unter diesen Annahmen erlaubt das Forschungsdesign die Anwendung des einfachsten IV-Ansatzes, wobei eine endogene Variable $B$ durch eine Instrumentvariable $Z$ *instrumentiert* wird. Die folgende Umformung zeigt, warum der kausale Effekt von B auf Y anhand der (Ko)Variation in diesen Variablen identifiziert werden kann: \begin{alignat*}{2} \textup{Cov}(Z,Y) &= \textup{Cov}(Z,\beta_0 + \beta_1 B + u) &\quad& \text{(Gl. \eqref{eq:simpleiv})} \\ \\ \textup{Cov}(Z,Y) &= \textup{Cov}(Z, \beta_1 B + u) &\quad& \text{($\beta_0$ konstant)} \\ \\ \textup{Cov}(Z,Y) &= \beta_1\textup{Cov}(Z,B) &\quad& \text{($\beta_1$ konstant, $Z$ exogen)} \\ \\ \beta_1 &= \frac{\textup{Cov}(Z,Y) }{\textup{Cov}(Z,B)} &\quad& \text{($Z$ relevant)} \end{alignat*} Eine naheliegende Implementierung gemäß dieses Identifikationsprinzips ist der einfache IV-Schätzer \begin{align} \widehat{\beta}_{\textup{IV}} = \frac{\widehat{\textup{Cov}}(Z,Y) }{\widehat{\textup{Cov}}(Z,B)}, \label{eq:simpleivest} \end{align} wobei lediglich die Kovarianzfunktion $\textup{Cov}(\cdot,\cdot)$ durch ihr Stichprobenäquivalent $\widehat{\textup{Cov}}(\cdot,\cdot)$ ersetzt wird. **Erwartungswert und Konsistenz** Eine hilfreiche Darstellung von \eqref{eq:simpleivest} ist \begin{align*} \widehat{\beta}_\textup{IV} = \beta_1 + \frac{\sum_{i=1}^n (Z_i-\overline Z) u_i}{\sum_{i=1}^n (Z_i - \overline{Z})B_i}. \end{align*} Anhand dieser Form kann der Erwartungswert sowie das Verhalten für große Stichproben untersucht werden. Eine wichtige Eigenschaft ist \begin{align} \textup{E}\big(\widehat{\beta}_\textup{IV}\big) = \beta_1 + \underbrace{\textup{E}\bigg(\frac{\sum_{i=1}^n (Z_i-\overline Z) u_i}{\sum_{i=1}^n (Z_i - \overline{Z})B_i}\bigg)}_{\neq0},\label{eq:ivbiasterm} \end{align} sodass $\widehat{\beta}_\textup{IV}$ bei Endogenität von $X$ ein *verzerrter Schätzer* von $\beta_1$ ist.^[Beachte, dass der endogene Regressor $B_i$ mit $u_i$ korreliert, sodass der Erwartungswert des Bruchs in \eqref{eq:ivbiasterm} ungleich $0$ ist.] Glücklicherweise kann man zeigen, dass \begin{align*} \frac{\sum_{i=1}^n (Z_i-\overline Z) u_i}{\sum_{i=1}^n (Z_i - \overline{Z})B_i}\to 0 \quad \text{für} \quad n\to\infty \end{align*} bei Gültigkeit der IV-Annahmen \eqref{eq:ivassum1} und \eqref{eq:ivassum2}, d.h. $\widehat{\beta}_\textup{IV}$ ist ein *konsistener* Schätzer für $\beta_1$, \begin{align*} \widehat{\beta}_\textup{IV} \to \beta_1 \quad \text{für} \quad n\to\infty. \end{align*} **Schwache Instrumente** Beachte, dass die Annahme der Relevanz des Instruments für die Herleitung des Identifikationsprinzips und des Schätzers \eqref{eq:simpleivest} von entscheidender Bedeutung ist: Sind $Z$ und $B$ unkorreliert, so ist $\text{Cov}(B,Z) = 0$ und $\beta_1$ damit nicht identifizierbar. Weiterhin würde $\widehat{\text{Cov}}(B,Z)\to0$ wenn $n\to\infty$, sodass $\widehat{\beta}_{\textup{IV}}\to\infty$! In empirischen Anwendungen ist $\text{Cov}(B,Z) = 0$ unwahrscheinlich. Die Kovarianz von $B$ und $Z$ kann jedoch gering sein. In diesem Fall bezeichnet man $Z$ als ein *schwaches Instrument*. Anhand von Gleichung \eqref{eq:simpleivest} können wir die Konsequenzen der Verwendung eines schwachen Instruments ableiten: Ein sehr schwacher linearer Zusammenhang von $B$ und $Z$ --- und damit tendenziell $\widehat{\text{Cov}}(B,Z) \approx 0$ --- kann trotz Relevanz von $Z$ zu einem Schätzer mit großer Varianz führen. Die Relevanz eines Instruments kann mit einem Hypothesentests geprüft werden. Wie die "Stärke" eines Instruments beurteilt werden sollte, ist eine nicht-triviale Frage. Wir betrachten diesen Aspekt näher in @sec-2SLS. ### Identifizierbarer Behandlungseffekt Regression ermöglicht in der Regel die Schätzung eines durchschnittlichen Behandlungseffekts (ATE) für eine betrachtete Population mit heterogenen Behandlungseffekten (in Modell \eqref{eq:simpleiv} unterstellen wir den empirisch unwahrscheinlichen Fall eines einheitlichen kausalen Effekts $\beta_1$). Das Adjustieren für Endogenität mit IV-Regression hat hier einen Preis: Beachte, dass im IV-Design ein Teil der beobachteten Variation der Behandlungsvariable $B$ (wie in \eqref{eq:simpleiv} angenommen) endogen ist und wir den Effekt von $B$ auf $Y$ anhand der *exogenen* Variation eines validen Instruments $Z$ schätzen. Daher können wir grundsätzlich nur einen *lokalen* durchschnittlichen Behandlungseffekt (LATE) identifizieren: Der LATE ist eine gewichtete Variante des ATE, wobei Beobachtungen, deren Behandlung gut durch $Z$ erklärt wird den größten Einfluss haben.^[Siehe Kapitel 10 und 19 in @Huntington-Klein2021 für nicht-technische Erläuterungen und Beispiel zu LATE.] ## 2SLS-Verfahren {#sec-2SLS} Der einfache IV-Schätzer \eqref{eq:simpleivest} ist ein Spezialfall eines allgemeineren Regressionsverfahrens, das in zwei Schritten Durchgeführt wird. Im fall von Modell \eqref{eq:simpleiv} erfolgt die Schätzung anhand der folgenden einfachen Regressionen: - **1. Stufe: Regression von $B$ auf $Z$** Regressiere die endogene Variable $B$ auf das Instrument $Z$: \begin{align} B_i = \alpha_0 + \alpha_1 Z_i + u_i, \quad 1,\dots,n. \end{align} Berechne die geschätzten Werte \begin{align} \widehat{B}_i = \widehat{\alpha}_0 + \widehat{\alpha}_1 Z_i. \end{align} Die $\widehat{B}_i$ aus dieser Regression sind "bereinigt": Die Variation in den $\widehat{B}_i$ wird lediglich durch die exogene Variation in $Z$ verursacht. - **2. Stufe: Regression von $Y$ auf $\widehat{B}$** Regressiere die Beobachtungen der abhängigen Variable $Y_i$ auf die geschätzten Werte $\widehat{B}_i$ aus der 1. Stufe: \begin{align} Y_i = \gamma_0 + \gamma_1 \widehat{B}_i + e_i. \end{align} Der geschätzte Koeffizient $\widehat{\gamma}_1$ aus dieser Regression ist der Two-Stage-Least-Squares-Schätzer (2SLS-Schätzer) von $\beta_1$. ### Äquivalenz von 2SLS und einfacher IV-Schätzer {#sec-equiviv} Der KQ-Schätzer von $\gamma_1$ in der 2. Stufe ist \begin{align} \widehat{\gamma}_1 = \frac{\sum (\widehat{B}_i - \overline{\widehat{B}})(Y_i - \bar{Y})}{\sum (\widehat{B}_i - \bar{\widehat{B}})^2} = \frac{\widehat{\text{Cov}}(\widehat{B}, Y)}{\widehat{\textup{Var}}(\widehat{B}_i)}. \end{align} Da $\widehat{B}_i = \widehat{\alpha}_0 + \widehat{\alpha}_1 Z_i$, die angepassten Werte aus der 1. Stufe, eine lineare Funktion von $Z_i$ sind, können wir $\widehat{\text{Cov}}(\widehat{B}, Y)$ wie folgt schreiben: \begin{align} \widehat{\text{Cov}}(\widehat{B}, Y) = \widehat{\alpha}_1 \widehat{\text{Cov}}(Z, Y) \end{align} Eine ähnliche Umformung zeigt, dass \begin{align} \widehat{\text{Var}}(\widehat{B}) = \widehat{\alpha}_1^2 \widehat{\text{Var}}(Z). \end{align} Kombinieren wir diese Ergebnisse, erhalten wir folgende Darstellung für den KQ-Schätzer von $\gamma_1$ in der 2. Stufe: \begin{align} \widehat{\gamma}_1 = \frac{\widehat{\text{Cov}}(\widehat{B}, Y)}{\widehat{\text{Var}}(\widehat{B})} = \frac{\widehat{\alpha}_1 \widehat{\text{Cov}}(Z, Y)}{\widehat{\alpha}_1^2 \widehat{\text{Var}}(Z)} = \frac{\widehat{\text{Cov}}(Z, Y)}{\widehat{\alpha}_1 \widehat{\text{Var}}(Z)} \end{align} Der KQ-Schätzer von $\widehat{\alpha}_1$ in der 1. Stufe ist definiert als $\widehat{\alpha}_1 = \frac{\widehat{\text{Cov}}(Z, B)}{\widehat{\text{Var}}(Z)}$. Somit gilt \begin{align} \widehat{\gamma}_1 = \frac{\widehat{\text{Cov}}(Z, Y)}{\frac{\widehat{\text{Cov}}(Z, B)}{\widehat{\text{Var}}(Z)} \widehat{\text{Var}}(Z)} = \frac{\widehat{\text{Cov}}(Z, Y)}{\widehat{\text{Cov}}(Z, B)} = \widehat{\beta}_{\textup{IV}}.\label{eq:equiviv} \end{align} ### Wozu 2SLS? ```{dot} //| fig-width: 5 //| fig-height: 4 //| fig-align: 'center' //| fig-cap: "IV-Design mit beobachtbaren Confoundern" //| label: "fig-2SLSDAG" digraph IV_DAG { layout=neato fontname="Helvetica,Arial,sans-serif" node [fontname="Helvetica,Arial,sans-serif", shape=circle] edge [fontname="Helvetica,Arial,sans-serif"] X [pos="-2,1.5!"]; U [pos="1.5,1.5!", color=gray, fontcolor=gray]; B [pos="0,0!"]; Y [pos="3,0!"]; Z [pos="-2,0!"]; X -> Z; X -> B; X -> Y; U -> B [color=gray, style=dashed]; U -> Y [color=gray, style=dashed]; B -> Y; Z -> B; } ``` In der Praxis werden IV-Schätzungen mit statistischer Software häufig anhand einer Implementierung des 2SLS-Verfahren durchgeführt. Gründe hierfür sind: 1. **Bedingte Exogenität**. Für den in @fig-IVDAG1 dargestellten DGP ist die Exogenität von Z gegeben. In empirischen Anwendungen kann es jedoch schwierig sein, ein Instrument Z zu finden, dass plausibel nur mit dem endogenen Regressor korreliert. Eine schwächere Bedingung als Exogenität ist Exogenität nach Kontrolle für *beobachtbare* gemeinsame Determinanten (X) von Z und Y.^[Formal: $\textup{Cov}(Z,u\vert X) = 0$.] Diese Situation ist in @fig-2SLSDAG dargestellt: Aufgrund der Pfeile von X zu Z und Y ist die Bedingung \eqref{eq:ivassum2} verletzt. Durch Kontrolle für X in *beiden* Regressionen des 2SLS-Verfahrens können die Pfade durch X geschlossen werden. 2. **Mehrere endogene Variablen / Instrumente**. Der 2SLS-Ansatz kann leicht für mehrere endogene Variablen (und mehrere Instrumentvariablen) erweitert werden. Für jede der endogenen Variablen erfolgt hierbei eine separate Regression in der ersten Stufe des 2SLS-Verfahrens, wobei jeweils mehrere Instrumente verwendet werden und zusätzlich für beobachtbare Kovariablen X kontrolliert werden kann. 3. **Statistische Inferenz**. Die Berechnung von $\widehat{\beta}_\textup{IV}$ erfordert Sorgfalt bei der Konstruktion von Inferenzstatistiken. Die in @sec-equiviv gezeigte Äquivalenz impliziert, dass wir bei der Schätzung der Variabilität von $\widehat{\beta}_\textup{IV}$ die Unsicherheit aus den beiden KQ-Schätzungen des 2SLS-Verfahrens berücksichtigen müssen. Korrigierte Standardfehlerformeln können vergleichsweise einfach anhand der Regressionsdarstellung hergeleitet werden und sind statistischer Software implementiert (bspw. `ivreg` in R). Weiterhin ermöglicht das Regressions-Framework diagnostische Tests. So kann die Relevanz von Instrumenten anhand von F-Tests für die Koeffizienten in den Regressionen der ersten Stufe des 2SLS-Verfahrens überprüft werden. Sofern mehr Instrumente als endogene Variablen vorliegen, kann die Exogenität der Instrumente getestet werden. ### Interaktive Illustrationen: 2SLS-Verfahren Der DGP im nachfolgenden interaktiven Beispiel ist \begin{align} \begin{split} X_i =&\, \gamma \cdot Z_i + v_i,\\ Y_i =&\, \beta_1 \cdot X_i + u_i, \end{split} \label{eq:OJSIVDGP} \end{align} wobei $Z_i \sim N(0, .25^2)$. Die Fehlerterme $v_i$ und $u_i$ sind $N(0,1)$-verteilt und folgen einer gemeinsam Normalverteilung mit Korrelationsparameter $\rho$, \begin{align} \begin{pmatrix} v_i \\ u_i \end{pmatrix} \sim N \bigg[ \begin{pmatrix} 0 \\ 0 \end{pmatrix} \begin{pmatrix} 1 & \rho \\ \rho & 1 \end{pmatrix} \bigg]. \end{align} Wir sind daran interessiert, den Parameter $\beta_1$ zu schätzen, um den Effekt von $X$ auf $Y$ zu bestimmen, wobei wir $n=1000$ Beobachtungen von $(X_i, Y_i, Z_i)$ verwenden. Beachte, dass $X_i$ eine Funktion der (exogenen) Variable $Z_i$ und des Fehlerterms $v_i$ ist. Wenn die Fehlerterme $(v_i, u_i)’$ im DGP \eqref{eq:OJSIVDGP} korreliert sind ($\rho\neq0$), dann besteht Korrelation zwischen $X_i$ und $u_i$. Dies impliziert, dass $X$ ein endogener Regressor im "naiven" Regressionsmodell \begin{align} Y_i = \beta_0 + \beta_1 X_i + e_i, \quad i = 1,\dots,N. \end{align} Der KQ-Schätzer von $\beta_1$ ist dann *verzerrt und inkonsistent*. Wenn $Z$ Erklärungskraft für $X$ hat, d.h. $\gamma\neq0$ mit $|\gamma|$ nicht zu klein gewählt ist ($Z$ ist relevant), und $Z$ nur über $X$ auf $Y$ wirkt ($Z$ ist exogen in der Gleichung von $Y_i$), dann kann der kausale Effekt von $X$ auf $Y$ anhand einer Instrumentvariablen-Regression mit $Z$ konsistent geschätzt werden. Die folgende Interaktive Grafik berechnet den 2SLS-Schätzer Analog zu der in @sec-2SLS erläuterten Vorgehensweise: 1. Schätze das Regressionsmodell \begin{align} X_i = \alpha_0 + \alpha_1 \cdot Z_i + \epsilon_i \end{align} und berechnen Sie die angepassten Werte $\widehat{X}_i$. Dieser Schritt isoliert die exogene Variation in $X_i$, die auf $Z_i$ zurückzuführen ist. 2. Schätze den Effekt von $X$ auf $Y$ unter Verwendung von $\widehat{X}_i$, dem exogenen Teil von $X_i$ aus 1., in der Regression \begin{align} Y_i = \delta_0 + \delta_1 \cdot \widehat{X}_i + \varepsilon_i. \end{align} Der KQ-Schätzer $\widehat{\delta}_1$ von $\delta_1$ ist der 2SLS-Schätzer von $\beta_1$. Die Grafik illustriert das Verhalten des IV-Schätzers (lila Grade) im Vergleich zum KQ-Schätzer (gestrichelte Grade) und zum wahren Zusammenhang (grüne Grade) für *eine* simulierte Stichprobe (gegeben der gewählten Parameter) anhand der geschätzten Regressionslinien. Die voreingestellten Parameter-Kombinationen zeigen - die Verzerrung von 2SLS, wenn $Z$ ein schwaches Instrument für $X$ ist, - eine Situation in der $X$ exogen ist (KQ ist unverzerrt), was zu vergleichbaren Schätzungen führt (sofern $Z$ kein schwaches Instrument ist) - eine starke Verzerrung von KQ, wenn $X$ endogen und irrelevant ($\beta_1 = 0$) ist. <iframe class="obs-soft-box-shadow" width="100%" height="910" frameborder="0" src="https://observablehq.com/embed/@mca91/ivreg?cells=webrinfo%2Cviewof+weak%2Cviewof+exogeneous%2Cviewof+beta1zero%2Cviewof+beta%2Cviewof+gamma%2Cviewof+cor%2Cviewof+showX_hat%2Cviewof+regtype%2Cchart%2Cloaded"> </iframe> Anhand der obigen interaktiven Grafik lässt sich zwar die Verzerrung der Schätzer einschätzen, jedoch nicht die Unsicherheit der Schätzung. Um die Verteilung der Schätzer in endlichen Stichproben zu illustrieren, zeigt die nachstehende Grafik Kerndichteschätzungen basierend auf $N$ simulierte Stichproben der Größe $n$ für den gewählten DGP.^[Beachte, dass die Simulation insbesondere für große $N$ und $n$ einige Sekunden dauern kann (grauer Balken am linken Rand zeigt laufende Berechnung an).] Die Applikation gibt weiterhin den R-Code für eine Simulation der Daten gemäß der gewählten Parameter aus. Die voreingestellten Parameter-Kombinationen zeigen: - Wenn $Z$ ein schwaches Instrument ist, kann der IV-Schätzer auch in großen Stichproben eine hohe Verzerrung haben *und* eine deutlich größere Varianz als der KQ-Schätzer haben. - Wenn der interessierende Regressor (Behandlungsvariable) $X$ exogen ist, hat der IV-Schätzer eine größere Varianz als der KQ-Schätzer. - Wenn die Endogenität von $X$ "stark" ist (d.h. hohe Korrelation von $u$ und $v$), kann der IV-Schätzer auch in kleinen Stichproben hilfreich sein, sofern das Instrument $Z$ nicht zu schwach ist: IV hat dann eine deutlich geringere Verzerrung als KQ, aber eine größere Varianz. <iframe class="obs-soft-box-shadow" width="100%" height="1350" frameborder="0" src="https://observablehq.com/embed/@mca91/bias-and-variance-in-iv-regression?cells=webrinfo%2Cviewof+weak%2Cviewof+exogeneous%2Cviewof+beta1zero%2Cviewof+N%2Cviewof+n%2Cviewof+beta%2Cviewof+gamma%2Cviewof+cor%2Cchart%2Cloaded%2Csimcode"> </iframe> ### Schätzung mit R Wir zeigen nachfolgend, wie die Daten aus den vorherigen interaktiven Beispielen in R generiert, die IV-Regression durchgeführt und ein (grafischer) Vergleich mit der naiven KQ-Schätzung implementiert werden kann. ```{webr} #| message: false # install.packages("faux") library(faux) # seed setzen set.seed(1234) # Stichprobengröße n <- 1000 # bivariate normalverteilte # Fehlerterme simulieren errors <- rnorm_multi( n, mu = c(0, 0), sd = c(1, 1), r = -0.8, varnames = c("u", "v"), ) # Variablen festlegen Z <- rnorm(n, sd = .25) X <- 2 * Z + errors$v Y <- 0 * X + errors$u # Daten sammeln the_data <- data.frame(X, Y) ``` ```{webr} library(dplyr) library(ggplot2) library(cowplot) # Xhat berechnen und # Regressionen plotten the_data %>% mutate( Xhat = lm(X ~ Z)$fitted ) %>% ggplot(aes(x = X, y = Y)) + geom_point() + geom_smooth( method = "lm", se = FALSE, col = "red" ) + geom_hline( yintercept = 0, col = "darkgreen" ) + geom_smooth( mapping = aes(x = Xhat, y = Y), method = "lm", se = FALSE, col = "steelblue" ) + theme_cowplot() ``` ```{webr} library(AER) # 'Naive' KQ-Regression KQ_mod_sim <- lm( formula = Y ~ X, data = the_data ) summary(KQ_mod_sim) ``` ```{webr} library(ivreg) # IV Regression mit Z # als Instrument iv_mod_sim <- ivreg( formula = Y ~ X | Z, data = the_data ) summary(iv_mod_sim) ``` ## Beispiel: Der schwache Staat, die Bauern und die Mafia @Acemogluetal2020 untersuchen den Aufstieg und die Auswirkungen der sizilianischen Mafia gegen Ende des 19. Jahrhunderts. Die Studie argumentiert, dass die Verbreitung der Mafia teilweise auf das Aufkommen sozialistischer Bauern-Organisationen (Fasci siciliani) zurückzuführen ist, die in einem Umfeld schwacher staatlicher Präsenz die Rechte der Bauern verteidigten. Grundbesitzer, Verwalter und lokale Politiker wandten sich also an die Mafia, um diese sozialistische "Bedrohung" zu bekämpfen und die Kontrolle über die landwirtschaftlich arbeitende Bevölkerung zu behalten. Die Studie identifiziert einen signifikanten (positiven) kausalen Zusammenhang zwischen dem Aufstieg der Fasci und der Ausbreitung der Mafia in Sizilien. ### Identifikationsstrategie Die Beziehung zwischen dem Aufkommen der Fasci und der Mafiaaktivität im einfachen Regressionsmodell \begin{align} \textup{Mafiaaktivität} = \beta_0 + \beta_1 \textup{Fasci siciliani} + \varepsilon \label{eq:fascimafiamod1} \end{align} ist wahrscheinlich aufgrund mehrerer Faktoren endogen: 1. **Externe Faktoren**, die sowohl die Entstehung der Fasci siciliani als auch der Mafia beeinflussen, sind wahrscheinlich. Beispielsweise könnten wirtschaftliche Notlagen und politische Instabilität sowohl die Organisation der Bauern als auch die Etablierung der Mafia begünstigen. 2. **Simultante Kasalität**: Die Präsenz der Fasci könnte die Mafiaaktivitäten beeinflussen und umgekehrt: Intensive Mafiaaktivitäten könnten die Notwendigkeit für Bauernorganisationen erhöhen, was wiederum die Mafia stärkt. Außerdem erschwert die zeitliche Nähe zwischen der Entstehung der Fasci und der Verbreitung der Mafia die Identifizierung der Kausalrichtung. Wir können das Modell \eqref{eq:fascimafiamod1} um (messbare) unter 1. fallende Regressoren erweitern, um einer verzerrten Schätzung aufgrund relevanter ausgelassener Variablen vorzubeugen. Eine Verzerrung aufgrund der in 2. beschriebenen simultanen Kausalität kann jedoch *nicht* durch multiple Regression vermieden werden. @Acemogluetal2020 nutzen die folgende Identifikationsstrategie: Eine extreme Dürre im Jahr 1893 (insbesondere Regenausfall im Frühling) hatte einen erheblichen negativen Einfluss auf die landwirtschaftliche Produktion in Sizilien und verursachte große Not unter den Bauern, was die Entstehung der Fasci beförderte. Die Dürre kann als natürliches Experiment betrachtet werden, wobei durch die Variation in der Niederschlagsmenge (das Instrument) die Variation im Aufkommen von Fasci-Organisationen in verschiedenen Gemeinden, die unabhängig von der späteren Mafiaaktivität ist, isoliert werden kann. Anhand dieser exogenen Variation im Organisationsgrad der Bauern kann der (lokale) kausale Effekt auf die Entwicklung der Mafia identifiziert werden. ```{dot} //| fig-width: 5 //| fig-height: 4 //| fig-align: 'center' //| fig-cap: "Zulässiger DGP im IV-Design von @Acemogluetal2020" //| label: "fig-ADDdag" digraph IV_DAG { layout=neato fontname="Helvetica,Arial,sans-serif" node [fontname="Helvetica,Arial,sans-serif", shape=circle] edge [fontname="Helvetica,Arial,sans-serif"] X [pos="3,-2!", label=<Mafia-<BR/>Faktoren>, shape=oval]; X2 [pos="0,-2!", label=<Fasci-<BR/>Faktoren>, shape=oval]; X3 [pos="-2,-2!", label=<Geo-<BR/>Faktoren>, shape=oval]; U [pos="1.5,1.5!", color=gray, fontcolor=gray]; B [pos="0,0!", label="Fasci"]; Y [pos="3,0!", label=<Mafia in<BR/> 1900>, shape=oval]; Z [pos="-2,0!", label="Dürre"]; Z -> X; Z -> X2; X2 -> B; X3 -> Y; X3 -> Z; X -> B; X -> Y; U -> B [color=gray, style=dashed]; U -> Y [color=gray, style=dashed]; B -> Y; Y -> B; Z -> B; } ``` ### Reproduktion der Kernergebnisse mit R Wir lesen zunächst den Datensatz `ADD_Mafia_municipality.csv`^[Der Datensatz stammt aus dem [MIT Economics Data Archive](https://economics.mit.edu/people/faculty/daron-acemoglu/data-archive) und liegt im STATA-Format `.dta` vor. Wir verwenden hier eine in das `.csv` formatierte Version.] ein und verschaffen uns einen Überblick. ```{webr} library(dplyr) library(readr) # Datensatz einlesen ADD_dat <- read_csv( file = "datasets/ADD_Mafia_municipality.csv" ) glimpse(ADD_dat) ``` @tbl-ADD_dat gibt einen Überblick der in `ADD_dat` verfügbaren Variablen auf Gemeinde-Ebene und deren Definition. ```{r, echo = F} #| tbl-cap: "`ADD_dat` --- Charakteristika von sizilianischen Gemeinden im 19. Jahrhundert" #| label: tbl-ADD_dat library(tibble) tibble( Variable = colnames(ADD_dat), Beschreibung = c( "Präsenz Fasci in 1893-Q4 (Dummy)", "Dürre-Intensität: Relative Regenmenge im Frühling 1893 (ggü. Mittel v. 1841 bis 1881) ", "Nächste aktive Wetterstation (Frühling 1893)", "Präsenz Fasci vor März 1893 (Dummy)", "Landwirtschaftliche Gemeinde? (Dummy)", "Durchschn. ländliche Pacht / Hektar in 1953", "Durchschn. städtische Pacht / Hektar in 1953", "Anteil bewirtschaftetes Land in 1853", "Anteil Anbaufläche Getreide in 1853", "Durschn. Schwefel-Produktion / Jahr für 1868-1870", "Anteil Anbaufläche Zitrusfrüchte in 1853", "Anteil Anbaufläche Oliven in 1853", "Anteil Anbaufläche Wein in 1853", "Intensität Mafia gem. Damiani (1885), Skala 0-3", "Intensität Mafia gem. Cutrera (1900), Skala 0-3", "Log(Einwohner) in 1861", "Log(Fläche) in 1853", "Höhe des Gemeindezentrums (M.ü. NN)", "Max. Höhe (M.ü. NN)", "Durschn. Höhe (M.ü. NN)", "Entfernung (Postweg) von Palermo in 1856", "Entfernung (Postweg) nächer Hafen in 1856", "Direkter Zugang Postweg in 1799 (Dummy)", "Durschn. Jahrestemperatur", "Durschn. Regenmenge im Frühling für 1881-1941", "Varianz Regenmenge im Frühling für 1881-1941", "Distrikt-Zugehörigkeit in 1853", "Provinz-Zugehörigkeit in 1853" ) ) %>% gt() %>% tabopts() ``` Vor der Analyse filtern wir die Daten entsprechend der Vorgehensweise in @Acemogluetal2020: Es werden lediglich Gemeinden berücksichtigt, die in maximal 30km Entfernung zu einer Wetter-Station liegen (`!is.na(sp3m1893_n30)`) *und* für die eine Messung der Mafiaaktivität gegen Ende des 19. Jahrhunderts (gemessen auf einer Skala von 0 bis 4, vgl. @Cutrera1900) vorliegt (`!is.na(Mafia1900)`). ```{webr} # Datensatz filtern ADD_dat <- ADD_dat %>% filter( # Nur Gemeinden mit Wetter-Station in 30km Umkreis !is.na(sp3m1893_n30), # Gemeinden mit Mafia-Präsenz Ende des 19. Jdt. !is.na(Mafia1900) ) ``` #### First-Stage-Regressionen Tabelle 3 in @Acemogluetal2020 präsentiert First-Stage-Regressionen Für den endogenen Regressor $\textup{Fasci siciliani}$ in \eqref{eq:fascimafiamod1}: Die Auswirkungen der Dürre in 1893 (das Instrument) auf den Organisationsgrad der Bauern (der endogene Regressor). Hierbei werden verschiedene, aus dem DAG in @fig-ADDdag abgeleitete Spezifikationen betrachtet: - **Panel A**: Keine Fixed Effects - Modell ohne Kovariablen - Modell mit Kovariablen für Bildung von Bauernorganisationen - Modell mit Kovariablen für Bildung von Bauernorganisationen und Indikatoren für Mafia-Präsenz - Modell mit Kovariablen für Bildung von Bauernorganisationen, Indikatoren für Mafia-Präsenz und geologische Faktoren - **Panel B**: Modelle wie in Panel A und Kontrolle für Provinz-Fixed-Effects Wir schätzten nachfolgend sämtliche Regressionen mit `feols` und berechnen nach Distrikt-Zugehörigkeit (`distretto1853`) sowie nächstgelegener Wetterstation (`cl1_stn_sp1893_n30`) geclusterte Standardfehler.^[Siehe Kapitel @sec-paneldata für Erläuterungen zu Fixed-Effects und cluster-robusten Standardfehlern.] Hierzu setzen wir `vcov = ~ distretto1853 + cl1_stn_sp1893_n30`. Diese Spezifikation berücksichtigt, dass die Niederschlagsmenge in einer Gemeinde sowie viele Kovariablen eine Korrelation mit benachbarten Gemeinden aufweisen können. Cluster-robuste Standardfehler sind inbesondere nötig, weil die Niederschlagsdaten für viele Gemeinden aus benachbarten Wetterstationen interpoliert werden, was zu einer zusätzlichen Korrelation zwischen den Beobachtungen führt. Standardfehler mit Two-way-clustering sind eine gängige Methode zur Berücksichtigung derartiger Korrelation. ```{webr} library(fixest) # First-Stage-Regression: # keine Kontroll-Variablen fasci_mod <- feols( peasants_fasci ~ sp3m1893_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat, ) # Statistische Zusammenfassung summary(fasci_mod) ``` Der obige Code-Chunk führt die First-Stage-Regression *ohne* Kovariablen aus. Wir finden einen signifikanten positiven Effekt der Niederschlagsmenge auf den Organisationsgrad der Bauern. @Acemogluetal2020 verwenden Bootstrap cluster-robuste Standardfehler [@Cameronetal2008;@Cameronetal2011], die aufwendig zu implementieren sind, jedoch zu nahezu identischen Ergebnissen führen. Der nächste Code-Chunk zeigt, wie cluster-robuste Bootstrap-Inferenz für das geschätzte Modell `fasci_mod` berechnet werden können. Hierfür nutzen wir `fwildclusterboot::boottest()`.^[`fwildclusterboot` implementiert Varianten des *Wild Bootstrap* [@Liu1988].] ```{webr} # Wild-Bootstrap-Standardfehler mit clustering library(fwildclusterboot) # boottest-Objekt erzeugen fasci_mod_boot <- boottest( object = fasci_mod, clustid = ~ distretto1853 + cl1_stn_sp1893_n30, param = "sp3m1893_n30", # Koeffizient B = 199 # Bootstrap-Züge ) # Statistische Zusammenfassung des Bootstraps summary(fasci_mod_boot) ``` Die statistische Zusammenfassung umfasst die Punktschätzung, die robuste t-Statistik sowie das 95\%-Konfidenzintervall für den Koeffizienten von `sp3m1893_n30`. Mit diesen Komponenten können wir den Bootstrap-Standardfehler beispielsweise durch die folgende Transformation berechnen: ```{webr} # Standardfehler über t-Statistik und # gesch. Koeffizienten berechnen 1/fasci_mod_boot$t_stat * fasci_mod_boot$point_estimate ``` Bei der Interpretation des geschätzten Koeffizienten müssen wir berückstigen, dass die abhängige Variable *binär* ist: Die First-Stage-Regression ist ein *lineares Wahrscheinlicheitsmodell*. Eine positive Änderung der relativen Regenmenge um eine Einheit (100\%) verringert die Wahrscheinlichkeit also um 100\%. Diese Interpretation ist jedoch wenig aussagekräftig für einen Einordnung des Effekts der Dürre. Stattdessen betrachten wir die Verteilung relativen Regenmengen und vergleichen den Effekt für eine von starker Dürre betroffenen Gemeinden (25\%-Quantil) und einer Gemeinde mit einem relativ geringen Rückgang des Niederschlags (75\%-Quantil) ```{webr} # Relativer Effekt der Dürre im Frühling 1893 ( ADD_dat %>% pull(sp3m1893_n30) %>% quantile( probs = c(.25, .75), na.rm = T ) ) %>% diff() * fasci_mod_boot$point_estimate ``` Für eine Gemeinde mit 75\% des langjährigen mittleren Niederschlags ist die geschätzte Wahrscheinlichkeit der Bildung von Bauern-Organisationen nur etwa halb so groß, wie für eine von starker Dürre betroffene Gemeinde. Zur Vermeidung von Backdoors durch ausgelassene Variablen betrachten @Acemogluetal2020 weitere Spezifikationen für die First- und Second-Stage-Regressionen, die zusätzlich für Determinanten der Fasci, der Mafia, sowie für geographische Faktoren kontrollieren. Weiterhin werden diese Regressionen mit Fixed Effekts auf der Provinz-Ebene gemäß der Grenzen von 1853 (`provincia1853`) wiederholt, um für unbeobachtbare Heteogenitäten zwischen den verschiedenen Provinzen zu kontrollieren. Die nachfolgenden Code-Chunks zeigen die Implementierung dieser Regressionen mit `fixtest::feols()`. Für eine bessere Übersicht der Ergebnisse verwenden wir in `summary()` das Argument `n = 1`, sodass jeweils nur die Schätzung des ersten Koeffizienten (Effekt der relativen Regenmenge) ausgegeben wird. ```{webr} # First-Stage-Regression: # + Fasci-Determinanten fasci_mod_F <- feols( peasants_fasci ~ sp3m1893_n30 # Fasci + predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) # Statistische Zusammenfassung summary(fasci_mod_F, n = 1) ``` ```{webr} # First-Stage-Regression: # + Fasci-Determinanten # + Mafia-Determinanten fasci_mod_FM <- feols( fml = peasants_fasci ~ sp3m1893_n30 # Fasci + predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel # Mafia + sulfurproduction1868_70 + Citrus_groves + Olives_groves + Vineyards + Mafia1885, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) # Statistische Zusammenfassung summary(fasci_mod_FM, n = 1) ``` ```{webr} # First-Stage-Regression: # + Fasci-Determinanten # + Mafia-Determinanten # + geographische Faktoren fasci_mod_FMG <- feols( peasants_fasci ~ sp3m1893_n30 # Fasci + predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel # Mafia + sulfurproduction1868_70 + Citrus_groves + Olives_groves + Vineyards + Mafia1885 # Geographie + lnpop1861 + lnsurface + centreheight + maxheight + slope2 + pa_pdist1856 + port2_pdist1856 + roads1799 + ave_temp + sp3m_ave_n30 + var_sp3m_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) # Statistische Zusammenfassung summary(fasci_mod_FMG, n = 1) ``` ```{webr} # First Stage # + Provinz-Fixed-Effects fasci_mod_FE <- feols( peasants_fasci ~ sp3m1893_n30 | provincia1853, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat, ) summary(fasci_mod_FE, n = 1) ``` ```{webr} # First Stage # + Provinz-Fixed-Effects # + Fasci-Determinanten fasci_mod_FE_F <-feols( fml = peasants_fasci ~ sp3m1893_n30 + predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel | provincia1853, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasci_mod_FE_F, n = 1) ``` ```{webr} # First Stage # + Provinz-Fixed-Effects # + Fasci-Determinanten # + Mafia-Determinanten fasci_mod_FE_FM <- feols( fml = peasants_fasci ~ sp3m1893_n30 + predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel + sulfurproduction1868_70 + Citrus_groves + Olives_groves + Vineyards + Mafia1885 | provincia1853, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasci_mod_FE_FM, n = 1) ``` ```{webr} # First Stage # + Provinz-Fixed-Effects # + Fasci-Determinanten # + Mafia-Determinanten # + geografische Faktoren fasci_mod_FE_FMG <- feols( peasants_fasci ~ sp3m1893_n30 + predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel + sulfurproduction1868_70 + Citrus_groves + Olives_groves + Vineyards + Mafia1885 + lnpop1861 + lnsurface + centreheight + maxheight + slope2 + pa_pdist1856 + port2_pdist1856 + roads1799 + ave_temp + sp3m_ave_n30 + var_sp3m_n30 | provincia1853, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasci_mod_FE_FMG, n = 1) ``` Für eine tabellarische Übersicht fassen wir die interessierenden Komponenten der Schätzungen mit `modelsummary::modelsummary()` zusammen. Hierfür sammeln wir zunächst jeweils die Modelle mit und ohne Provinz-Fixed-Effects in benannten Listen, sodass die Ergebnisse später in separaten Panelen dargestellt werden können. Für die Formatierung von Labels für die Fixed Effects und Koeffizienten definieren wir mit `tibble::tribble()` das `tibble`-Objekt `gm` und übergeben reguläre Ausdrücke für eine Auswahl an zusammenfassenden Statistiken. ```{webr} library(modelsummary) # Modelle in benannten Listen sammeln panels <- list( "Panel A" = list( "Keine Controls" = fasci_mod, "Fasci" = fasci_mod_F, "Fasci Mafia" = fasci_mod_FM, "Fasci Mafia Geo" = fasci_mod_FMG ), "Panel B" = list( "Keine Controls" = fasci_mod_FE, "Fasci" = fasci_mod_FE_F, "Fasci Mafia" = fasci_mod_FE_FM, "Fasci Mafia Geo" = fasci_mod_FE_FMG ) ) ``` ```{webr} # Format-Definitionen für FE-label gm <- tribble( ~raw, ~clean, ~fmt, "FE: provincia1853", "FE: Provinz", ~fmt ) ``` ```{webr} # Tabellarische Zusammenfassung # Modelle ohne Fixed Effects modelsummary( models = panels[["Panel A"]], coef_omit = "^(?!(sp3m1893\\_n30*)$).*", stars = T, coef_map = c( "sp3m1893_n30" = "Rel. Niederschlag 1893" ), gof_map = gm, notes = c( "Abhängige Variable: Peasants Fasci", "Cluster-robuste Standardfehler: Naheste Wetterstation + Distrikt" ), title = "Acemoglu et al. (2020) -- Panel A: 1st-Stage-Schätzungen ohne Fixed Effects", output = "gt" ) %>% tabopts() # Modelle mit Fixed Effects modelsummary( models = panels[["Panel B"]], coef_omit = "^(?!(sp3m1893\\_n30*)$).*", stars = T, gof_omit = "^(?!(FE.*)$).*", coef_map = c( "sp3m1893_n30" = "Rel. Niederschlag 1893" ), gof_map = gm, notes = c( "Abhängige Variable: Peasants Fasci", "Cluster-robuste Standardfehler: Naheste Wetterstation + Distrikt" ), title = "Acemoglu et al. (2020) -- Panel B: 1st-Stage-Schätzungen mit Fixed Effects", output = "gt" ) %>% tabopts() ``` Die Ergebnisse bestätigen die Vermutung der Autoren, dass der Effekt des relativen Niederschlags auf den Organisationsgrad der Bauern ohne Kontrolle für die betrachteten Kovariablen tendenziell überschätzt wird. Tatsächlich führt multiple Regression mit Fixed Effects und sämtlichen Kovariablen (letzte Spalte in Panel B) zu einem etwa 30\% kleineren geschätzten Effekt als für die einfache Regression in Panel A. #### Second-Stage-Regressionen Mit `fml = Mafia1900 ~ 1 | peasants_fasci ~ sp3m1893_n30` legen wir fest, dass `Mafia1900` auf eine Konstante sowie die angepassten Werte des endogenen Regressors `peasants_fasci` aus der ersten Stufe (instrumentiert durch `sp3m1893_n30`) regressiert werden soll. Wie für die First-Stage-Regressionen schätzen wir jeweils vier verschiedene Spezfikationen, unterscheiden zwischen Modellen mit und ohne Provinz-Fixed-Effects und berechnen cluster-robuste Standardfehler. ```{webr} # Second-Stage-Regression: # Keine Kontrollvariablen # (1) fasciIV_mod <- feols( fml = Mafia1900 ~ 1 # IV-Spezifikation | peasants_fasci ~ sp3m1893_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasciIV_mod, n = 1) ``` ```{webr} # Second-Stage-Regression: # + Fasci-Determinanten fasciIV_mod_F <- feols( fml = Mafia1900 ~ # Fasci predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel # IV-Spezifikation | peasants_fasci ~ sp3m1893_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasciIV_mod_F, n = 1) ``` ```{webr} # Second-Stage-Regression: # + Fasci-Determinanten # + Mafia-Determinanten fasciIV_mod_FM <- feols( fml = Mafia1900 ~ # Fasci predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel # Mafia + sulfurproduction1868_70 + Citrus_groves + Olives_groves + Vineyards + Mafia1885 # IV-Spezifikation | peasants_fasci ~ sp3m1893_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasciIV_mod_FM, n = 1) ``` ```{webr} # Second-Stage-Regression: # + Fasci-Determinanten # + Mafia-Determinanten # + Geographie fasciIV_mod_FMG <- feols( fml = Mafia1900 ~ # Fasci predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel # Mafia + sulfurproduction1868_70 + Citrus_groves + Olives_groves + Vineyards + Mafia1885 # Geographie + lnpop1861 + lnsurface + centreheight + maxheight + slope2 + pa_pdist1856 + port2_pdist1856 + roads1799 + ave_temp + sp3m_ave_n30 + var_sp3m_n30 # IV-Spezifikation | peasants_fasci ~ sp3m1893_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasciIV_mod_FMG, n = 1) ``` ```{webr} # Second-Stage-Regression: # + Provinz-FE fasciIV_mod_FE <- feols( fml = Mafia1900 ~ 1 # Fixed Effecs | provincia1853 # IV-Spezifikation | peasants_fasci ~ sp3m1893_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasciIV_mod_FE, n = 1) ``` ```{webr} # Second-Stage-Regression: # + Fasci-Determinanten # + Provinz-FE fasciIV_mod_FE_F <- feols( fml = Mafia1900 ~ # Fasci predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel # Fixed Effects | provincia1853 # IV-Spezifikation | peasants_fasci ~ sp3m1893_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasciIV_mod_FE_F, n = 1) ``` ```{webr} # Second-Stage-Regression: # + Fasci-Determinanten # + Mafia-Determinanten # + Provinz-FE fasciIV_mod_FE_FM <- feols( fml = Mafia1900 ~ # Fasci predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel # Mafia + sulfurproduction1868_70 + Citrus_groves + Olives_groves + Vineyards + Mafia1885 # Fixed Effects | provincia1853 # IV-Spezifikation | peasants_fasci ~ sp3m1893_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasciIV_mod_FE_FM, n = 1) ``` ```{webr} # Second-Stage-Regression: # + Fasci-Determinanten # + Mafia-Determinanten # + Geographische Faktoren # + Provinz-FE fasciIV_mod_FE_FMG <- feols( fml = Mafia1900 ~ # Fasci predr_peas_fasci + ruralcentre1861 + Rural_rent + Urban_rent + agricola_rel + seminatoritot_rel # Mafia + sulfurproduction1868_70 + Citrus_groves + Olives_groves + Vineyards + Mafia1885 # Geo + lnpop1861 + lnsurface + centreheight + maxheight + slope2 + pa_pdist1856 + port2_pdist1856 + roads1799 + ave_temp + sp3m_ave_n30 + var_sp3m_n30 # Fixed Effects | provincia1853 # IV-Spezifikation | peasants_fasci ~ sp3m1893_n30, vcov = ~ distretto1853 + cl1_stn_sp1893_n30, data = ADD_dat ) summary(fasciIV_mod_FE_FMG, n = 1) ``` ```{webr} # Modelle in benannten Listen sammeln panels <- list( "Panel A" = list( "Keine Controls" = fasciIV_mod, "Fasci" = fasciIV_mod_F, "Fasci + Mafia" = fasciIV_mod_FM, "Fasci + Mafia + Geo" = fasciIV_mod_FMG ), "Panel B" = list( "Keine Controls" = fasciIV_mod_FE, "Fasci" = fasciIV_mod_FE_F, "Fasci + Mafia" = fasciIV_mod_FE_FM, "Fasci + Mafia + Geo" = fasciIV_mod_FE_FMG ) ) ``` ```{webr} # Tabellarische Zusammenfassungen # keine Fixed Effects modelsummary( models = panels[["Panel A"]], coef_omit = "^(?!(fit\\_peasants\\_fasci*)$).*", stars = T, gof_omit = "^(?!(FE.*)$).*", coef_map = c( "fit_peasants_fasci" = "Fasci" ), gof_map = gm, notes = c( "Abhängige Variable: Mafia-Präsenz in 1900", "Instrument: Rel. Niederschlag 1893", "Cluster-robuste Standardfehler: Naheste Wetterstation + Distrikt" ), title = "Acemoglu et al. (2020) -- Panel A: 2SLS-Regressionen", output = "gt" ) %>% tabopts() # mit Fixed Effects modelsummary( models = panels[["Panel B"]], coef_omit = "^(?!(fit\\_peasants\\_fasci*)$).*", stars = T, gof_omit = "^(?!(FE.*)$).*", coef_map = c( "fit_peasants_fasci" = "Fasci" ), gof_map = gm, notes = c( "Abhängige Variable: Mafia-Präsenz in 1900", "Instrument: Rel. Niederschlag 1893", "Cluster-robuste Standardfehler: Naheste Wetterstation + Distrikt" ), title = "Acemoglu et al. (2020) -- Panel B: 2SLS-Regressionen mit Fixed Effects", output = "gt" ) %>% tabopts() ``` Die Ergebnisse für die IV-Schätzung des kausalen Effekts von Bauernorganisationen auf die Präsenz der Mafia zu Beginn des 20. Jahrhundert sind in der obigen Tabelle dargestellt. Die LATE-Schätzung in Spalte 4 von Panel B hat folgende Interpretation: Die Präsenz der Fasci in einer Gemeinde erhöht die Intensität der Mafiaaktivität um etwa 1.56 Punkte.^[Nach @Cutrera1900 entspricht eine Erhöhung des Mafia-Index von 1 auf 2 einem Sprung von einer geringen Präsenz der Mafia zu einer signifikanten Präsenz.] Wie ist dieser Effekt einzuordnen? Wir berechnen zunächst die *Anzahl* an Gemeinden, in denen es im Jahr 1900 Fasci gibt. ```{webr} library(tidyr) # Anz. Gemeinden mit Fasci im Jahr 1900 ADD_dat %>% drop_na() %>% filter( peasants_fasci == 1 ) %>% count() ``` Für diese 84 Gemeinden erwarten wir jeweils einen Effekt von 1.56 Punkten auf die im Jahr 1900 beobachtete *gesamte* Intensität des Einflusses der Mafia in Sizilien, gemessen als gewichtete Summe. ```{webr} # Beitrag der Gemeinden mit Fasci im Jahr 1900 84 * 1.56 ``` Die gesamte Mafia-Intensität für Sizilien zu Beginn des 20. Jahrhunderts berechnen wir analog. ```{webr} # Mafia-Intensität in Sizilien im Jahr 1900 ADD_dat %>% drop_na() %>% count(Mafia1900) %>% mutate(Beitrag = Mafia1900 * n) %>% summarise( Intensität = sum(Beitrag) ) ``` Der Anteil der Gemeinden mit Bauernorganisationen im Jahr 1900 an dieser Summe erlaubt eine Einordnung des Effekts der Fasci auf die Verbreitung der Mafia. ```{webr} # Anteil des geschätzten Effekts 84 * 1.56 / 342 ``` Die Interpretation lautet, dass bis zu 38\% der Stärke der Mafia im Jahr 1900 in Sizilien auf ihren Einsatz gegen die Fasci zurückzuführen sein könnte. ## Case Study: Ökonomische Schocks und Bürgerkriege Ein Kernproblem bei der Schätzung der kausalen Beziehung zwischen wirtschaftlichen Schocks und dem Auftreten von kriegerischen Auseinandersetzungen ist die Endogenität von Variablen, die wirtschaftliche Stabilität messen [vgl. @FearonLaitin2003]. Beispielsweise kann eine Verschlechterung der wirtschaftlichen Bedingungen eines Landes die Wahrscheinlichkeit eines Konflikts im Inland erhöhen. Gleichermaßen beeinflussen kriegerische Handlung im Inland die wirtschaftliche Situation einer Volkswirtschaft negativ. ```{dot} //| fig-width: 6 //| fig-height: 4 //| fig-align: 'center' //| fig-cap: "Simultane Kausalität zw. Kriegsrisiko und Wirtschaftslage" //| label: "fig-IVDAGrain1" digraph IV_DAG { layout=neato fontname="Helvetica,Arial,sans-serif" node [fontname="Helvetica,Arial,sans-serif"] edge [fontname="Helvetica,Arial,sans-serif"] K [pos="3,5"]; BIP [pos="0,0!"]; Konflikt [pos="6,0!"]; K -> BIP; K -> Konflikt; BIP -> Konflikt; Konflikt -> BIP; } ``` In der "Konfliktgleichung" ```{=tex} \begin{align} \textup{Konflikt-Wsk.} = \beta_0 + \beta_1\, \Delta\text{BIP} + \text{Kontrollv.} + \epsilon\label{eq:rainconf1} \end{align} ``` liegt dann Endogenität des Regressors für "wirtschaftliche Stabilität", hier $\Delta\text{BIP}$, aufgrund von simultanter Kausalität vor. Die multiple Regression \eqref{eq:rainconf1} liefert dann eine verzerrte Schätzungen des interessierenden Effekts $\beta_1$. @Migueletal2004 untersuchen die kausale Beziehung zwischen wirtschaftlichen Schocks und dem Auftreten von Bürgerkriegen in 41 afrikanischen Sub-Sahara-Ländern für den Zeitraum von 1981 bis 1999.^[Füge `rainconflict$country %>% unique()` im nächsten WebR-Chunk für eine Liste der Länder aus.] Hierbei adressieren die Autoren das Endogenitätsproblem in einem IV-Ansatz mit Instrumenten, die mit wirtschaftlichen Schocks korreliert sind, aber nicht unmittelbar mit der Wahrscheinlichkeit eines Bürgerkriegs zusammenhängen. Die Identifikationsstrategie basiert auf metereologischen Faktoren: In den betrachteten Agrarökonomien sind Schocks in den Niederschlägen entscheidend für die Entwicklung des Bruttoinlandsprodukts. Ungünstige Wetterbedingungen, die landwirtschaftliche Erträge negativ beeinflussen, haben einen direkten (negativen) Effekt auf die wirtschaftliche Situation. Variation in der Regenmenge ist jedoch keine unmittelbare Determinante der Konfliktwahrscheinlichkeit und folglich exogen in Modellgleichung \eqref{eq:rainconf1}. Maße für die Veränderungen in der Niederschlagsmenge sind somit plausible Instrumente für endogene Variablen, die Veränderungen im Bruttoinlandsprodukt messen. ```{dot} //| fig-width: 6 //| fig-height: 4 //| fig-align: 'center' //| fig-cap: "IV-Design bei simultaner Kausalität zw. Kriegsrisiko und Wirtschaftslage" //| label: "fig-IVDAGrain2" digraph IV_DAG { layout=neato fontname="Helvetica,Arial,sans-serif" node [fontname="Helvetica,Arial,sans-serif"] edge [fontname="Helvetica,Arial,sans-serif"] C [pos="1.5,1.5!"]; Niederschlag [pos="-2,0!"]; BIP [pos="0,0!"]; Konflikt [pos="3,0!"]; Niederschlag -> BIP; C -> BIP; C -> Konflikt; BIP -> Konflikt; Konflikt -> BIP; } ``` Anhand von Länder-Paneldaten mit Informationen über wirtschaftliche Indikatoren, Konfliktereignisse und Niederschläge, isolieren @Migueletal2004 in einem IV-Ansatz die exogene Variation der wirtschaftlichen Entwicklung, d.h. Variation in $\Delta\text{BIP}$ die nicht mit dem Fehlerterm in der Konfliktgleichung korreliert ist. Die Ergebnisse ihrer 2SLS-Schätzungen zeigen, dass negative ökonomische Schocks die Wahrscheinlichkeit eines Bürgerkriegs signifikant erhöhen. @Migueletal2004 zeigen weiterhin, dass ihre Ergebnisse robust gegenüber verschiedenen Modell-Spezifikationen unter Berücksichtigung diverser Kontrollvariablen sind. Wir reproduzieren nachfolgend die zentralen Ergebnisse von @Migueletal2004 anhand eines Auszugs (`rainconflict.csv`) aus dem der Studie zugrundeliegenden Datensatz. Der vollständige Datensatz ist, nach Registrierung, [hier](https://dataverse.harvard.edu/dataset.xhtml?persistentId=doi:10.7910/DVN/27324) verfügbar. ```{webr} library(readr) # Datensatz 'rainconflict' einlesen rainconflict <- read_csv2( file = "datasets/rainconflict.csv" ) ``` ```{webr} # Überblick verschaffen glimpse(rainconflict) ``` ```{r, echo = F} #| tbl-cap: "`rainconflict` -- Polit-ökonomische Charakteristika afrikanischer Länder v. 1981 bis 1999" #| label: tbl-rainconfdata tribble( ~Variable, ~Beschreibung, "any_prior", "Dummyvariable: Konflikt mit mind. 25 Toten im Jahr t", "war_prio", "Dummyvariable: Konflikt mit mind. 1000 Toten im Jahr t", "ccode", "Länderkennung (numerisch)", "country", "Name des Landes", "gdp_g", "(BIP_t - BIP_t-1) / BIP_t-1", "gdp_g_l", "(BIP_t-1 - BIP_t-2) / BIP_t-2", "GPCP_g", "(Niederschlag_t - Niederschlag_t-1) / Niederschlag_t-1", "GPCP_g_l", "GPCP_g in t-1", "GPCP_g_fl", "GPCP_g in t+1", "gdp_1979", "Log Pro-Kopf-BIP in 1979", "polity2_IV", "Diff. zw. Demokratie- und Autokrarie-Score im Jahr t (Skala v. -10 bis 10)", "polity2l", "polity2_IV im Jahr t-1 (Skala v. -10 bis 10)", "ethfrac", "Ethno-linguistische Fragmentierung (Anteil)", "relfrac", "Religiöse Fragmentierung (Anteil)", "oil", "Dummy für Öl-exportierendes Land", "lmtnest", "Log Anteil Gebirgsregionen an Landesfläche", "lpopl1", "Log Bevölkerung im Jahr t-1", "tot_100_g", "Handelsbedingungen (Index, 1995 = 100)", "year", "Jahr der Beobachtung" ) %>% gt() %>% tabopts() ``` ```{webr} # ccode zu Typ 'factor' transformieren rainconflict <- rainconflict %>% mutate( ccode = as.factor(ccode) ) ``` ```{webr} library(modelsummary) # Statistische Zusammenfassung datasummary( formula = All( rainconflict %>% select(-ccode, -year) ) ~ (mean + sd) * Arguments(na.rm = TRUE) + N, fmt = 4, data = rainconflict, title = "`rainconflict.csv` -- Statistische Zusammenfassung", output = "gt" ) %>% tabopts() ``` ```{webr} library(fixest) # (1) ( rncnf_mod1 <- feols( fml = gdp_g ~ GPCP_g + GPCP_g_l, data = rainconflict, vcov = ~ ccode ) ) ``` ```{webr} # (2) rncnf_mod2 <- feols( fml = gdp_g ~ GPCP_g + GPCP_g_l + gdp_1979 + polity2l + ethfrac + relfrac + oil + lmtnest + lpopl1 + year:ccode, data = rainconflict, vcov = ~ ccode ) rncnf_mod2 %>% print(n = 10) ``` ```{webr} # (3) rncnf_mod3 <- feols( fml = gdp_g ~ GPCP_g + GPCP_g_l + year:ccode | ccode, data = rainconflict, vcov = ~ ccode ) rncnf_mod3 %>% print(n = 2) ``` ```{webr} # (4) rncnf_mod4 <- feols( fml = gdp_g ~ GPCP_g + GPCP_g_l + GPCP_g_fl + year:ccode | ccode, data = rainconflict, vcov = ~ ccode ) rncnf_mod4 %>% print(n = 3) ``` ```{webr} # (5) rncnf_mod5 <- feols( fml = gdp_g ~ GPCP_g + GPCP_g_l + tot_100_g + year:ccode | ccode, data = rainconflict, vcov = ~ ccode ) rncnf_mod5 %>% print(n = 3) ``` ```{webr} # Tabellarischer Vergleich mit modelsummary() # (Tabelle 3 in Ditella und Schargrodsky, 2004) modelsummary( models = list( "(1)" = rncnf_mod1, "(2)" = rncnf_mod2, "(3)" = rncnf_mod3, "(4)" = rncnf_mod4, "(5)" = rncnf_mod5 ), stars = T, coef_omit = "^year.*$", gof_omit = "^(?!(R2|Num.Obs.|FE.*)$).*", output = "gt" ) %>% tab_caption( caption = "`rainconflict.csv` -- First-Stage-Regressionen für `gdp_g`" ) %>% tabopts() ``` ```{webr} S1_res <- tibble( x = residuals( feols( fml = GPCP_g ~ GPCP_g_l + year:ccode | ccode, data = rainconflict ) ), y = residuals( feols( fml = gdp_g ~ GPCP_g_l + year:ccode | ccode, data = rainconflict ) ) ) ``` ```{webr} library(ggplot2) library(cowplot) ggplot( data = S1_res, mapping = aes(x = x, y = y)) + geom_point( size = .75, alpha = .5 ) + geom_smooth(method = "loess", span = .5) + coord_cartesian( xlim = c(-.4, .4), ylim = c(-.1, .1) ) + theme_cowplot() ``` ```{webr} feols( fml = any_prio ~ GPCP_g + GPCP_g_l + year:ccode | ccode, data = rainconflict, vcov = ~ ccode ) %>% print(n = 2) ``` ```{webr} feols( fml = war_prio ~ GPCP_g + GPCP_g_l + year:ccode | ccode, data = rainconflict, vcov = ~ ccode ) %>% print(n = 2) ``` ```{webr} rf_res <- tibble( x = residuals( feols( fml = GPCP_g_l ~ GPCP_g + year:ccode | ccode, data = rainconflict ) ), y = residuals( feols( fml = any_prio ~ GPCP_g + year:ccode | ccode, data = rainconflict ) ) ) ``` ```{webr} ggplot( data = rf_res , mapping = aes(x = x, y = y)) + geom_point(pch = 19) + geom_smooth(method = "loess") + coord_cartesian( xlim = c(-.4, .4), ylim = c(-.2, .4) ) + theme_cowplot() ``` ```{webr} # (1) mod_conf_probit <- feglm( fml = any_prio ~ gdp_g + gdp_g_l + gdp_1979 + polity2l + ethfrac + relfrac + oil + lmtnest + lpopl1 + year, data = rainconflict, vcov = ~ ccode, family = binomial(link = "probit") ) #library(marginaleffects) #mod_conf_probit_avge <- mod_conf_probit %>% # avg_slopes() ``` ```{r, echo=F} library(readr) library(dplyr) library(fixest) # Datensatz 'rainconflict' einlesen rainconflict <- read_csv2( file = "datasets/rainconflict.csv" ) # ccode zu Typ 'factor' transformieren rainconflict <- rainconflict %>% mutate( ccode = as.factor(ccode) ) mod_conf_probit <- feglm( fml = any_prio ~ gdp_g + gdp_g_l + gdp_1979 + polity2l + ethfrac + relfrac + oil + lmtnest + lpopl1 + year, data = rainconflict, vcov = ~ ccode, family = binomial(link = "probit") ) ``` Mit dem Paket `marginaleffects` können wir durchschnittliche marginale Einflüsse von Regressoren in einem nicht-linearen Modell (hier ein Probit-Modell) berechnen. Für ein Probit-Modell lässt sich der marginale Effekt eines Regressors $X_j$ auf die Wahrscheinlichkeit, dass das Ereignis eintritt, wie folgt darstellen: Die Wahrscheinlichkeitsfunktion des Probit-Modells ist definiert als \begin{align*} P(Y = 1 | \boldsymbol{x}) = \Phi(\boldsymbol{x}'\boldsymbol\beta), \end{align*} wobei $\Phi$ die kumulative Verteilungsfunktion der Standardnormalverteilung ist, $\boldsymbol{x}$ der Vektor der Regressoren und $\boldsymbol{\beta}$ der Vektor der Koeffizienten. Der marginale Effekt für den Regressor $j$ ist dann die Ableitung von $P(Y = 1 | \boldsymbol{x})$ bezüglich $x_j$, \begin{align*} \frac{\partial P(Y = 1 | \boldsymbol{x})}{\partial x_j} = \phi(\boldsymbol{x}'\boldsymbol{\beta}) \cdot \beta_j, \end{align*} siehe auch @sec-probitreg. Hierbei ist $\phi(\cdot)$ die Dichtefunktion der Standardnormalverteilung und $\beta_j$ der Koeffizient des Regressors $x_j$. Die Funktion `marginaleffects::avg_slopes()` berechnet diesen Effekt anhand der geschätzten Koeffizienten $\widehat{\boldsymbol{\beta}}$ für jede Beobachtung $i=1,\dots,N$ und bestimmt den durchschnittlichen Effekt als das arithmetische Mittel über alle marginalen Effekte: \begin{align*} \textup{Average Marginal Effect}(j) = \frac{1}{N} \sum_{i=1}^{N} \phi(\boldsymbol{x}_i' \widehat{\boldsymbol{\beta}}) \cdot \widehat\beta_j, \end{align*} wobei $N$ die Anzahl der Beobachtungen ist und $\boldsymbol{x}_i$ den Wert von $\boldsymbol{x}$ für die $i$-te Beobachtung darstellt. ```{r} # Marginale Effekte berechnen library(marginaleffects) avg_slopes(mod_conf_probit) ``` ```{webr} # (2) ( mod_conf_ols <- feols( fml = any_prio ~ gdp_g + gdp_g_l + gdp_1979 + polity2l + ethfrac + relfrac + oil + lmtnest + lpopl1 + year, data = rainconflict, vcov = ~ ccode ) ) # (3) feols( fml = any_prio ~ -1 + gdp_g + gdp_g_l + gdp_1979 + polity2l + ethfrac + relfrac + oil + lmtnest + lpopl1 + i(ccode, year), data = rainconflict, vcov = ~ ccode ) %>% summary(n = 10) # (4) feols( fml = any_prio ~ -1 + gdp_g + gdp_g_l + i(ccode, year) | ccode, data = rainconflict, vcov = ~ ccode ) %>% summary(n = 10) ``` ```{webr} #### IV models #### # (5) feols( fml = any_prio ~ + gdp_1979 + polity2l + ethfrac + relfrac + oil + lmtnest + lpopl1 + i(ccode, year) | gdp_g + gdp_g_l ~ GPCP_g + GPCP_g_l, data = rainconflict, vcov = ~ ccode ) %>% summary(n = 10) # (6) feols( fml = any_prio ~ + i(ccode, year) | ccode | gdp_g + gdp_g_l ~ GPCP_g + GPCP_g_l, data = rainconflict, vcov = ~ ccode ) %>% summary(n = 10) # (7) feols( fml = war_prio ~ + i(ccode, year) | ccode | gdp_g + gdp_g_l ~ GPCP_g + GPCP_g_l, data = rainconflict, vcov = ~ ccode ) %>% summary(n = 10) ``` ```{webr} modelsummary( models = list( mod_conf_ols ), output = "gt" ) %>% tabopts() ```