立命館大学の柴田悠氏への私信なのですが、パネルデータ分析の固定効果モデルに年次などタイムダミーが入っているときに、ある瞬間に個々が同じような行動をとって、説明変数の値を一斉に変化させたときに、それがタイムダミーの係数として観測されない事を、シミュレーションで示したいと思います*1。
1. モデル
計量モデルは以下になります。
Yが被説明変数、Xが説明変数、aが固定効果、βとγは推定する係数、Iが個々の数、iが個々を表す添字、Tが期間の長さ、tが期を表す添字とします。εは誤差項です。
2. 係数の真の値
Iは10、Tは5として、真の値はβが0.015、が{0.00 0.05 0.04 0.06 0.02}、
が{0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0}を設定します。
3. 観測値の生成
観測値になるXは一様分布から乱数で0から1を生成し、t=3のときだけ+1するようにしましょう。ある瞬間に個々が同じような行動をとった例とします。誤差項εは標準偏差0.1の正規分布から生成します。Xとεが決まれば、真の値から観測値Yも定まります。
Rで生成するためのコードは以下になります。
set.seed(201408)
I <- 10 # 個体の数
ai <- (1:i)/10 # 固定効果
T <- 5 # 期間の長さ
gamma <- c(0, (5^(1:(T-1)) %% 7)/100) # タイムダミーの効果
epsilin <- rnorm(I*T, sd=0.1)
df <- data.frame(
i = numeric(I*T), t = numeric(I*T), Y = numeric(I*T), X = numeric(I*T)
)
for(i in 1:I){
for(t in 1:T){
l <- (i-1)*T + t
df$i[l] <- i
df$t[l] <- t
df$X[l] <- runif(1) + ifelse(t==3, 1, 0)
df$Y[l] <- ai[m] + 0.15*df$X[l] + gamma[t] + epsilin[l]
}
}
生成した観測値の一部を見ると、こうなります。t=3のときに、明らかにXが大きくなっていることを確認してください。
| i | t | Y | X |
|---|---|---|---|
| 1 | 1 | 1.0242424 | 0.90365560 |
| 1 | 2 | 1.2447202 | 0.74466442 |
| 1 | 3 | 1.2130497 | 1.03975109 |
| 1 | 4 | 1.3497088 | 0.29219997 |
| 1 | 5 | 1.0248503 | 0.02165759 |
| 2 | 1 | 1.3007814 | 0.63130932 |
| 2 | 2 | 1.3448006 | 0.86396715 |
| 2 | 3 | 1.2383986 | 1.17330402 |
| 2 | 4 | 1.2283682 | 0.66754810 |
| 2 | 5 | 1.1475519 | 0.75067122 |
4. 推定結果
ちゃっちゃと固定効果モデルで推定してみましょう。plmパッケージを使います。
library(plm)
plmdata <- plm.data(df, index=c("i", "t"))
td <- as.factor(df$t)
rfe <- plm(Y ~ X + td, model="within", data=plmdata)
summary(rfe)
推定結果は以下になります。
Oneway (individual) effect Within Model
Call:
plm(formula = Y ~ X + td, data = plmdata, model = "within")Balanced Panel: n=10, T=5, N=50
Residuals :
Min. 1st Qu. Median 3rd Qu. Max.
- 0.16200 -0.05120 -0.00315 0.04550 0.16500
Coefficients :
Estimate Std. Error t-value Pr(>|t|)
X 0.188583 0.049956 3.7750 0.0005953 ***
td2 0.067277 0.041042 1.6392 0.1101244
td3 -0.015109 0.059223 -0.2551 0.8001284
td4 0.113756 0.040291 2.8233 0.0077890 **
td5 0.043061 0.040308 1.0683 0.2926920
- -
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Total Sum of Squares: 0.54784
Residual Sum of Squares: 0.284
R-Squared : 0.4816
Adj. R-Squared : 0.33712
F-statistic: 6.50321 on 5 and 35 DF, p-value: 0.00022576
Xが0.188583と推定されている一方で、t=3のタイムダミーの係数は-0.015109となっています。観測数が十分ではないのか精度に欠けますが、タイムダミーの値は過少に推定されていますから、タイムダミーの係数として観測されない事が分かります。
また、t=3のインパクトを大きくとっても、推定結果が大きく変わることはありません。観測数を十分増やすと、推定値は真の値に近づきます。この例ではI=300ぐらいにすると、Xの係数βが0.1541469になりました*2。
