具有特定状态趋势的固定效应回归答案

【问题标题】：Fixed effects regression with state specific trends具有特定状态趋势的固定效应回归
【发布时间】：2015-12-11 21:27:04
【问题描述】：

我正在处理面板数据，我想估计具有特定状态趋势的固定效应回归。

在 Stata 中，我可以通过以下方式完成此操作，

xi i.state i.year i.state*time
reg y x _I*

上面将创建状态虚拟变量、年份虚拟变量和 50 个（状态 x 时间）虚拟变量，其中时间以数字方式识别趋势（即 1、2、3...）

在 R 中，我可以使用 plm 或 lm 运行固定效果模型，例如，

plm(y ~ x, index = c("state", "year"), effect = "twoways", data = df)
 lm(y ~ x + factor(state) + factor(year), data = df)

我如何像 xi 在 Stata 中那样包含 50 个（状态 x 时间）假人？

我知道interaction() 不是我想要的，因为这会创建一个具有 n 个级别的新因子变量，其中 n =（状态数）x（时间段数）。我想要做的是创建 50 个（状态 x 时间）变量，使得 state1xtime 为 1,2,3 ...当 state == 1 和零时，重复 state2xtime，其中 state == 2，等等。

【问题讨论】：

我不确定你为什么使用factor(year)。我认为只是lm(y ~ x + factor(state) + year + factor(state):year, data = df) - 我假设您将您的年份变量居中。

标签： r stata

【解决方案1】：

您只需将state 与year 交互即可完成此操作。正确的运算符是:，它只包括交互项。

请注意lm 和plm 之间存在细微差别：

对于lm，只需state:year
对于 plm，年份已隐式转换为因子，state:as.integer(year) 也是如此（执行state:year 将为您提供州和年份的所有组合）。

检查：

library(plm)

data("Produc", package = "plm")
produc_plm <- pdata.frame(Produc, index = c("state","year"))

### simple state-specific time trend ###
fe1_Ttrend_lm <- lm(log(gsp) ~ -1 + log(pcap) + log(pc) + log(emp) + unemp + state +state:year, data = Produc)
fe1_Ttrend_plm <- plm(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp + state : as.integer(year), data = produc_plm)

summary(fe1_Ttrend_lm)$coef[1:4,]
summary(fe1_Ttrend_plm)$coef[1:4,]

【讨论】：

【解决方案2】：

要完成 Matifou 的回答，您也可以使用包 fixest：

library(fixest)
data("Produc", package = "plm")
fe_fixest = feols(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp + year::state | state, data = Produc)
# Notice the double colon (on the left the numeric variable, on the right the factor). The alternative also works,
# fe_fixest = feols(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp + interact(year, state) | state, data = Produc)
# fe_fixest = feols(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp + i(year, state) | state, data = Produc)

# Requesting ("standard" standard-errors, otherwise, clustered at state level by default)
coeftable(fe_fixest, se = "standard")[1:4, ]

请注意，如果您不关心交互系数的标准误差，您可以使用以下语法将它们投影出来：

fe_fixest_bis = feols(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp | state[year], data = Produc)
# state[year] means 'state' fixed-effects and 'state x year' interactions
# The interacted terms are projected out, and the estimation is faster

coeftable(fe_fixest_bis, se = "s")

【讨论】：

【解决方案3】：

这就是你要找的东西吗：

dta <- data.frame(state = rep(LETTERS[1:3], 10), 
              time = rep(1:3, each = 10))
dta <- cbind(dta, model.matrix( ~ state - 1, data = dta) * dta$time)

head(dta, 1)
#     state time stateA stateB stateC
# 1     A    1      1      0      0

tail(dta, 1)
#      state time stateA stateB stateC
# 30     C    3      0      0      3

【讨论】：

可能有更好的方法，但此解决方案适用于以下修改，
dta <- data.frame(state = rep(LETTERS[1:3], 10), time = rep(seq(1, 10, 1), each = 3))