DID回帰モデル
UPDATE: 2024-12-07 15:51:22.995503
はじめに
今回は、差分の差分法の仮定をおさらいし、DID回帰モデルを実際に実行するまでをまとめておく。
DID回帰
9.6.6 Castle-doctrine statutes and homicidesのCheng and Hoekstra (2013)の例を参考にDIDの回帰を実行してみる。書籍のRコードはこちら。
データセットの背景
「Castle-doctrine statutes and homicides」の章のデータセットは、異なるタイミングで行われた銃の規制に関する処置と暴力の関係を分析しているデータセットで実際の研究データを再現しているものになる。Castle Doctrineは、自宅や自分が法的に占有する場所にいる際に、侵入者や脅威に対する自己防衛を正当化する法的な原則のことで、これが犯罪抑止効果を持っているという説と、正当防衛という名の殺傷をしやすくなったという点で、犯罪率を増加させている説がある。Cheng and Hoekstra (2013)の研究では、Castle Doctrineはが殺人率を増加させる可能性を示唆している。殺人件数が8%増加することで、法改正が行われた21州で年間約600件の殺人事件が増加することを示した。下記はアブストラクトを和訳したもの。
2000年から2010年にかけて、20以上の州が、正当防衛のための殺傷力の行使を容易にする法律を可決した。これらの法律には、自宅以外の場所での退却義務の撤廃、危害が迫っているという合理的確信の推定規定の追加、法律に基づいて行動する者の民事責任の撤廃などが含まれる。 本稿では、このような方法で正当防衛を助けることが犯罪抑止につながるのか、あるいは逆に殺人を増加させるのかを検証する。 そのために、法の採用における州内のばらつきを利用した差分研究デザインを適用する。 強盗、窃盗、加重暴行は法律の影響を受けていない。 一方、殺人事件は約8％増加し、これらの殺人事件は警察によって大部分が殺人事件として分類されていることがわかった。このことは、自衛権強化の主な結果が殺人の純増であることを示唆している。最後に、報告された正当化される殺人の相対的増加に関する証拠と、過少報告の程度と性質に関する仮定を用いて、この増加全体が法的に正当化されたものであるかどうかを評価するために、封筒裏計算を提示する。
DIDモデル
Cheng and Hoekstra (2013)では下記のモデルが利用されている。\(D_{it}\)はは本来01のことが多いものの、いくつかの州が年の途中で法改正を行うため、0から1までの範囲を取る変数になっている。7月に法改正した場合、\(D_{it}\)は法改正の採用年よりも前は0、採用年は0.5、採用年よりも後は1となる。\(X_{it}\)はReigin-by-Yearという地域別年別固定とのこと。
\[ Y_{it} = \alpha + \delta D_{it} + \gamma X_{it} + \sigma_{i} + \tau_{t} + \epsilon_{it} \]
パネルデータは下記の通りダウンロードして使用できる。
library(bacondecomp)
library(tidyverse)
library(haven)
library(lfe)
read_data <- function(df){
full_path <- paste("https://github.com/scunning1975/mixtape/raw/master/", df, sep = "")
df <- read_dta(full_path)
return(df)
}
castle <- read_data("castle.dta")
castle## # A tibble: 550 × 185
## state year sid cdl pre2_cdl caselaw anywhere assumption civil homicide_c
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Alab… 2000 1 0 0 0 0 0 0 329
## 2 Alab… 2001 1 0 0 0 0 0 0 379
## 3 Alab… 2002 1 0 0 0 0 0 0 303
## 4 Alab… 2003 1 0 0 0 0 0 0 299
## 5 Alab… 2004 1 0 1 0 0 0 0 254
## 6 Alab… 2005 1 0 1 0 0 0 0 374
## 7 Alab… 2006 1 0.581 0 0 1 0 1 368
## 8 Alab… 2007 1 1 0 0 1 0 1 392
## 9 Alab… 2008 1 1 0 0 1 0 1 357
## 10 Alab… 2009 1 1 0 0 1 0 1 322
## # ℹ 540 more rows
## # ℹ 175 more variables: robbery_gun_r <dbl>, jhcitizen_c <dbl>,
## # jhpolice_c <dbl>, homicide <dbl>, robbery <dbl>, assault <dbl>,
## # burglary <dbl>, larceny <dbl>, motor <dbl>, murder <dbl>,
## # hc_felonywsus <dbl>, jhcitizen <dbl>, jhpolice <dbl>, population <dbl>,
## # police <dbl>, unemployrt <dbl>, income <dbl>, blackm_15_24 <dbl>,
## # whitem_15_24 <dbl>, blackm_25_44 <dbl>, whitem_25_44 <dbl>, …全部はいらないので、カラムを絞っておく。
castle <- castle %>%
select(state, year, sid, cdl, post, l_homicide, popwt) %>%
mutate(year_f = as.factor(year))2006年に法改正した州と法改正しなかった州の並行トレンドを確認したいが、これは難しい問題であれる。
# 2006年に法改正した州
treated_states <-c(
"Alabama", "Alaska", "Arizona", "Florida", "Georgia",
"Indiana", "Kansas", "Kentucky", "Louisiana", "Michigan",
"Mississippi", "Oklahoma", "South Carolina", "South Dakota"
)
# 法改正した州
#treated_states <- c(
# "Alabama", "Alaska", "Arizona", "Florida", "Georgia",
# "Indiana", "Kansas", "Kentucky", "Louisiana", "Michigan",
# "Mississippi", "Missouri", "Montana", "North Dakota", "Ohio",
# "Oklahoma", "South Carolina", "South Dakota",
# "Tennessee", "Texas", "West Virginia"
#)
nevertreated_states <- c(
"Arkansas", "California", "Colorado", "Connecticut", "Delaware",
"Hawaii", "Idaho", "Illinois", "Iowa", "Maine",
"Maryland", "Massachusetts", "Minnesota", "Nebraska", "Nevada",
"New Hampshire", "New Jersey", "New Mexico", "New York", "North Carolina",
"Oregon", "Pennsylvania", "Rhode Island", "Utah", "Vermont",
"Virginia", "Washington", "Wisconsin", "Wyoming"
)
# データフレームにフラグ列を追加
tmp <- castle %>%
mutate(
state_flag = case_when(
state %in% treated_states ~ "Treated", # 2006 Treated states
state %in% nevertreated_states ~ "Never Treated", # Never treated states
TRUE ~ "Other" # 2006以降
)
)
ggplot(tmp, aes(x = year, y = l_homicide, color = state_flag, group = state)) +
geom_line(size = 1.2) +
geom_point(size = 3) +
geom_vline(xintercept = 2005.2, linetype = "solid", color = "gray", size = 1) +
labs(
title = "States That Enacted Castle Doctrine in 2006 \nCompared to States That Did Not Enact Castle Doctrine from 2000 to 2010",
x = "Year",
y = "Log Homicides per 100K",
color = "State Type",
linetype = "State Type"
) +
scale_x_continuous(
breaks = seq(2000, 2010, by = 1),
labels = seq(2000, 2010, by = 1)
) +
scale_y_continuous(
breaks = seq(1.0, 2.5, by = 0.1)
) +
theme_minimal()
そのため、ここでは研究で示されている通り、2006年に法改正した州と法改正しなかった州で並行トレンドを確認する。これを見て並行トレンドがあることを決定することはどうなのかわからないが、法改正した州では雑人率が上がっていそうである。
filtered_data <- tmp %>%
filter(state_flag != "Other") %>%
group_by(state_flag, year) %>%
summarise(l_homicide_m = mean(l_homicide))
# 折れ線グラフを描画
ggplot(filtered_data, aes(x = year, y = l_homicide_m, color = state_flag, group = state_flag)) +
geom_line(size = 1.2) +
geom_point(size = 3) +
geom_vline(xintercept = 2005.2, linetype = "solid", color = "gray", size = 1) +
labs(
title = "Figure 4: \nLog Homicide Rate for the 13 States That Enacted Castle Doctrine in 2006 \nCompared to States That Did Not Enact Castle Doctrine from 2000 to 2010",
x = "Year",
y = "Log Homicides per 100K",
color = "State Type",
linetype = "State Type"
) +
scale_x_continuous(
breaks = seq(2000, 2010, by = 1),
labels = seq(2000, 2010, by = 1)
) +
scale_y_continuous(
breaks = seq(1.0, 2.5, by = 0.1)
) +
theme_minimal()
最後にこれをDIDが回帰モデルで法改正の効果を確認する。研究だと下記のようになっている。
固定効果として、州と年を入れてDID回帰を実行した結果がこちら。研究の数値ともだいたい一致していることがわかる。また、法改正によって、殺人件数が8%増加することも核にできる。
library(estimatr)
lm_robust(
l_homicide ~ cdl,
fixed_effects = state + year_f,
clusters = state,
weights = popwt,
data = castle,
se_type = "stata"
) ## Estimate Std. Error t value Pr(>|t|) CI Lower CI Upper DF
## cdl 0.08006351 0.03589849 2.230275 0.03034154 0.007922786 0.1522042 49ちなみにこの書き方でもOK。
# こちらでもOK
lm_robust(
l_homicide ~ 0 + cdl + state + year_f,
data = castle,
clusters = state,
weights = popwt,
se_type = "stata") ## Estimate Std. Error t value Pr(>|t|)
## cdl 0.080063510 0.03589849 2.2302749 3.034154e-02
## stateAlabama 1.980730466 0.02185005 90.6510657 3.098423e-56
## stateAlaska 1.609874052 0.02153730 74.7481835 3.689040e-52
## stateArizona 1.995793337 0.02197303 90.8292127 2.815813e-56
## stateArkansas 1.858420619 0.02165309 85.8270498 4.443086e-55
## stateCalifornia 1.868509041 0.02165309 86.2929613 3.413649e-55
## stateColorado 1.280527866 0.02165309 59.1383500 3.145575e-47
## stateConnecticut 1.192163178 0.02165309 55.0574221 9.935227e-46
## stateDelaware 1.429889016 0.02165309 66.0362646 1.508435e-49
## stateFlorida 1.732397345 0.02271280 76.2740592 1.381884e-52
## stateGeorgia 1.933580696 0.02175776 88.8685638 8.151706e-56
## stateHawaii 0.751905417 0.02165309 34.7250903 3.609013e-36
## stateIdaho 0.724781538 0.02165309 33.4724338 2.033402e-35
## stateIllinois 1.906231304 0.02165309 88.0350806 1.289889e-55
## stateIndiana 1.707184482 0.02175776 78.4632539 3.492074e-53
## stateIowa 0.499509336 0.02165309 23.0687350 5.720204e-28
## stateKansas 1.417874361 0.02187236 64.8249324 3.699440e-49
## stateKentucky 1.544674409 0.02172254 71.1092892 4.163317e-51
## stateLouisiana 2.520124751 0.02162094 116.5594243 1.465726e-61
## stateMaine 0.459213079 0.02165309 21.2077413 2.440286e-26
## stateMaryland 2.223758034 0.02165309 102.6993511 7.067440e-59
## stateMassachusetts 0.997547505 0.02165309 46.0695273 5.236701e-42
## stateMichigan 1.855818281 0.02148679 86.3701838 3.268178e-55
## stateMinnesota 0.833157919 0.02165309 38.4775575 2.815795e-38
## stateMississippi 2.094474007 0.02175776 96.2633197 1.658955e-57
## stateMissouri 1.872740088 0.02075953 90.2111039 3.927074e-56
## stateMontana 1.053380456 0.02049712 51.3916400 2.747484e-44
## stateNebraska 1.125072168 0.02165309 51.9589721 1.619673e-44
## stateNevada 2.046783348 0.02165309 94.5261662 4.030234e-57
## stateNew Hampshire 0.177530001 0.02165309 8.1988308 9.502878e-11
## stateNew Jersey 1.481415519 0.02165309 68.4159023 2.710811e-50
## stateNew Mexico 2.064171366 0.02165309 95.3291934 2.668486e-57
## stateNew York 1.565356935 0.02165309 72.2925512 1.868517e-51
## stateNorth Carolina 1.867258926 0.02165309 86.2352275 3.526711e-55
## stateNorth Dakota 0.240449420 0.02080404 11.5578231 1.328675e-15
## stateOhio 1.517390169 0.02042797 74.2800110 5.005785e-52
## stateOklahoma 1.727066536 0.02140237 80.6951178 8.927123e-54
## stateOregon 0.835226425 0.02165309 38.5730869 2.503147e-38
## statePennsylvania 1.737502020 0.02165309 80.2426915 1.173484e-53
## stateRhode Island 1.107763109 0.02165309 51.1595914 3.415961e-44
## stateSouth Carolina 2.001204125 0.02182475 91.6942672 1.774375e-56
## stateSouth Dakota 0.804471457 0.02175776 36.9740053 1.862704e-37
## stateTennessee 1.935070423 0.02093393 92.4370541 1.197626e-56
## stateTexas 1.791632973 0.02075317 86.3305848 3.341967e-55
## stateUtah 0.748416850 0.02165309 34.5639785 4.493098e-36
## stateVermont 0.613546350 0.02165309 28.3352825 4.812843e-32
## stateVirginia 1.700505189 0.02165309 78.5340746 3.342204e-53
## stateWashington 1.121204027 0.02165309 51.7803305 1.911758e-44
## stateWest Virginia 1.286389427 0.02053025 62.6582429 1.919152e-48
## stateWisconsin 1.150693078 0.02165309 53.1422172 5.475736e-45
## stateWyoming 0.919208660 0.02165309 42.4516208 2.616748e-40
## year_f2001 0.021185065 0.01718593 1.2326983 2.235701e-01
## year_f2002 0.010314126 0.02349636 0.4389671 6.626150e-01
## year_f2003 0.030194835 0.02265211 1.3329814 1.887028e-01
## year_f2004 0.005179136 0.02537146 0.2041323 8.390956e-01
## year_f2005 0.032462732 0.03351395 0.9686335 3.374851e-01
## year_f2006 0.018511399 0.02932479 0.6312543 5.308080e-01
## year_f2007 -0.012939745 0.03254003 -0.3976562 6.926107e-01
## year_f2008 -0.059908107 0.03233934 -1.8524837 6.998473e-02
## year_f2009 -0.133784945 0.03394592 -3.9411198 2.572925e-04
## year_f2010 -0.182036401 0.03863099 -4.7121862 2.059905e-05
## CI Lower CI Upper DF
## cdl 0.007922786 0.152204235 49
## stateAlabama 1.936821143 2.024639788 49
## stateAlaska 1.566593227 1.653154876 49
## stateArizona 1.951636873 2.039949801 49
## stateArkansas 1.814907111 1.901934127 49
## stateCalifornia 1.824995533 1.912022550 49
## stateColorado 1.237014357 1.324041374 49
## stateConnecticut 1.148649670 1.235676686 49
## stateDelaware 1.386375507 1.473402524 49
## stateFlorida 1.686754267 1.778040423 49
## stateGeorgia 1.889856845 1.977304548 49
## stateHawaii 0.708391909 0.795418926 49
## stateIdaho 0.681268030 0.768295047 49
## stateIllinois 1.862717795 1.949744812 49
## stateIndiana 1.663460631 1.750908334 49
## stateIowa 0.455995828 0.543022845 49
## stateKansas 1.373920205 1.461828518 49
## stateKentucky 1.501021330 1.588327488 49
## stateLouisiana 2.476675838 2.563573664 49
## stateMaine 0.415699570 0.502726587 49
## stateMaryland 2.180244526 2.267271543 49
## stateMassachusetts 0.954033997 1.041061014 49
## stateMichigan 1.812638954 1.898997609 49
## stateMinnesota 0.789644411 0.876671428 49
## stateMississippi 2.050750155 2.138197858 49
## stateMissouri 1.831022251 1.914457925 49
## stateMontana 1.012189957 1.094570954 49
## stateNebraska 1.081558660 1.168585677 49
## stateNevada 2.003269839 2.090296856 49
## stateNew Hampshire 0.134016493 0.221043510 49
## stateNew Jersey 1.437902010 1.524929027 49
## stateNew Mexico 2.020657858 2.107684875 49
## stateNew York 1.521843427 1.608870444 49
## stateNorth Carolina 1.823745417 1.910772434 49
## stateNorth Dakota 0.198642135 0.282256705 49
## stateOhio 1.476338617 1.558441721 49
## stateOklahoma 1.684056869 1.770076203 49
## stateOregon 0.791712917 0.878739934 49
## statePennsylvania 1.693988512 1.781015528 49
## stateRhode Island 1.064249601 1.151276618 49
## stateSouth Carolina 1.957345655 2.045062594 49
## stateSouth Dakota 0.760747606 0.848195309 49
## stateTennessee 1.893002124 1.977138722 49
## stateTexas 1.749927924 1.833338021 49
## stateUtah 0.704903341 0.791930358 49
## stateVermont 0.570032842 0.657059859 49
## stateVirginia 1.656991681 1.744018698 49
## stateWashington 1.077690519 1.164717536 49
## stateWest Virginia 1.245132344 1.327646511 49
## stateWisconsin 1.107179570 1.194206587 49
## stateWyoming 0.875695151 0.962722168 49
## year_f2001 -0.013351352 0.055721483 49
## year_f2002 -0.036903568 0.057531821 49
## year_f2003 -0.015326276 0.075715946 49
## year_f2004 -0.045806732 0.056165003 49
## year_f2005 -0.034886069 0.099811533 49
## year_f2006 -0.040418966 0.077441763 49
## year_f2007 -0.078331384 0.052451894 49
## year_f2008 -0.124896451 0.005080237 49
## year_f2009 -0.202001829 -0.065568060 49
## year_f2010 -0.259668275 -0.104404527 49# こちらでもOK
#library(fixest)
#feols(
# l_homicide ~ cdl | state + year,
# weights = castle$popwt,
# data = castle) %>%
# fixef()参考文献および参考資料
Cheng and Hoekstra (2013)の研究を再現するためのスクリプト
library(bacondecomp)
library(tidyverse)
library(haven)
library(lfe)
read_data <- function(df)
{
full_path <- paste("https://github.com/scunning1975/mixtape/raw/master/",
df, sep = "")
df <- read_dta(full_path)
return(df)
}
castle <- read_data("castle.dta")
#--- global variables
crime1 <- c("jhcitizen_c", "jhpolice_c",
"murder", "homicide",
"robbery", "assault", "burglary",
"larceny", "motor", "robbery_gun_r")
demo <- c("emo", "blackm_15_24", "whitem_15_24",
"blackm_25_44", "whitem_25_44")
# variables dropped to prevent colinearity
dropped_vars <- c("r20004", "r20014",
"r20024", "r20034",
"r20044", "r20054",
"r20064", "r20074",
"r20084", "r20094",
"r20101", "r20102", "r20103",
"r20104", "trend_9", "trend_46",
"trend_49", "trend_50", "trend_51"
)
lintrend <- castle %>%
select(starts_with("trend")) %>%
colnames %>%
# remove due to colinearity
subset(.,! . %in% dropped_vars)
region <- castle %>%
select(starts_with("r20")) %>%
colnames %>%
# remove due to colinearity
subset(.,! . %in% dropped_vars)
exocrime <- c("l_lacerny", "l_motor")
spending <- c("l_exp_subsidy", "l_exp_pubwelfare")
xvar <- c(
"blackm_15_24", "whitem_15_24", "blackm_25_44", "whitem_25_44",
"l_exp_subsidy", "l_exp_pubwelfare",
"l_police", "unemployrt", "poverty",
"l_income", "l_prisoner", "l_lagprisoner"
)
law <- c("cdl")
dd_formula <- as.formula(
paste("l_homicide ~ ",
paste(
paste(xvar, collapse = " + "),
paste(region, collapse = " + "),
paste(lintrend, collapse = " + "),
paste("post", collapse = " + "), sep = " + "),
"| year + sid | 0 | sid"
)
)
#Fixed effect regression using post as treatment variable
dd_reg <- felm(dd_formula, weights = castle$popwt, data = castle)
summary(dd_reg)