# plmパッケージを読み込む
library(plm)Warning: package 'plm' was built under R version 4.4.3# Grunfeldデータを読み込む
data("Grunfeld", package = "plm")
df <- Grunfeld # dfとしていったん保存する
# 中身を見てみる
head(df) firm year inv value capital
1 1 1935 317.6 3078.5 2.8
2 1 1936 391.8 4661.7 52.6
3 1 1937 410.6 5387.1 156.9
4 1 1938 257.7 2792.2 209.2
5 1 1939 330.8 4313.2 203.4
6 1 1940 461.2 4643.9 207.2Grundfeld（1950）による米国の10の製造業のパネルデータ。中身は以下の通り
firm: 企業ID
year: 年
inv: 投資額
value: 企業価値
capital: 企業設備
plmパッケージを用いるためには、どの変数が識別用のデータとなっているのかをindex引数で定義する必要がある。ここでは、firmとyearが企業と年を識別するインデックスになっているので、そちらを指定する。中身を見てみると、一番左列にインデックスが付されていることが分かる。
# plm用のデータとしてデータフレームを再定義する
dfp <- pdata.frame(df, index = c("firm","year"))
# 中身を見てみる
head(dfp) firm year inv value capital
1-1935 1 1935 317.6 3078.5 2.8
1-1936 1 1936 391.8 4661.7 52.6
1-1937 1 1937 410.6 5387.1 156.9
1-1938 1 1938 257.7 2792.2 209.2
1-1939 1 1939 330.8 4313.2 203.4
1-1940 1 1940 461.2 4643.9 207.2以下、企業価値（value）と資産状況（capital）が投資額（inv）をどの程度説明するのかをモデル化する。以下のような単純なデータの場合、推定値はすべて同じになる。
model_pooled <- plm(inv ~ value + capital, data = dfp, mdoel = "pooling")model_fixed <- plm(inv ~ value + capital, data = dfp, mdoel = "within")model_random <- plm(inv ~ value + capital, data = dfp, mdoel = "random")結果を比べてみる。
library(texreg)Warning: package 'texreg' was built under R version 4.4.3Version: 1.39.4
Date: 2024-07-23
Author: Philip Leifeld (University of Manchester)
Consider submitting praise using the praise or praise_interactive functions.
Please cite the JSS article in your publications -- see citation("texreg").screenreg(list(pooled = model_pooled,
fixed = model_fixed,
random = model_random))
=============================================
pooled fixed random
---------------------------------------------
value 0.11 *** 0.11 *** 0.11 ***
(0.01) (0.01) (0.01)
capital 0.31 *** 0.31 *** 0.31 ***
(0.02) (0.02) (0.02)
---------------------------------------------
R^2 0.77 0.77 0.77
Adj. R^2 0.75 0.75 0.75
Num. obs. 200 200 200
=============================================
*** p < 0.001; ** p < 0.01; * p < 0.05ケースだけではなく、時間ダミーも入れて、ケースと時間の両方効果から統制した固定効果モデルを行う場合には、effectを”two-ways”に指定する。
model_fixed_twoways <- plm(inv ~ value + capital, model = "within",
effect = "twoways", data = dfp)
# 結果を比べてみる
screenreg(list(`within case` = model_fixed,
`twoways` = model_fixed_twoways))
==================================
within case twoways
----------------------------------
value 0.11 *** 0.12 ***
(0.01) (0.01)
capital 0.31 *** 0.36 ***
(0.02) (0.02)
----------------------------------
R^2 0.77 0.72
Adj. R^2 0.75 0.67
Num. obs. 200 200
==================================
*** p < 0.001; ** p < 0.01; * p < 0.05基本モデルの場合には、ある時点X_tの値が、当該時点のy_tに作用することをモデル化している。もしも一時点前のXの値（X_t-1）が現在の時点の従属変数y_tに影響することをモデル化したい場合には、独立変数の対応をt時間分をずらしたうえで、モデル化すればよい。
なお、plmで使われているlag関数は、名前がベースパッケージに入っているものと競合するため、plm::と明確にplmの方のものを用いることを指定したほうが安全。以下はk = 1時間分のラグを持たせる場合の例。
library(Formula)
# 一時点前の企業価値が、現在の投資状況に影響をあたえることを想定する場合
model_fixed_lag <- plm(inv ~ plm::lag(value, 1) + plm::lag(capital, 1), model = "within", effect = "twoways", data = dfp)
# 結果を比べる
screenreg(list(`twoways` = model_fixed_twoways,
`twoways-lag` = model_fixed_lag))
=============================================
twoways twoways-lag
---------------------------------------------
value 0.12 ***
(0.01)
capital 0.36 ***
(0.02)
plm::lag(value, 1) 0.10 ***
(0.02)
plm::lag(capital, 1) 0.42 ***
(0.03)
---------------------------------------------
R^2 0.72 0.59
Adj. R^2 0.67 0.51
Num. obs. 200 190
=============================================
*** p < 0.001; ** p < 0.01; * p < 0.05plmパッケージには、Betweenという平均値をケース分コピーする関数が存在する。以下のスクリプトで挙動を確かめることができる。
# firm = 1の平均値 --> 608.02
mean(dfp$inv[dfp$firm==1])[1] 608.02Between(dfp$inv)[1:10] 1 1 1 1 1 1 1 1 1 1
608.02 608.02 608.02 608.02 608.02 608.02 608.02 608.02 608.02 608.02 これを利用して、各企業の平均ダミーをランダム効果モデルに投入する。
# ハイブリッドモデル
model_hybrid <- plm(inv ~ value + capital + Between(value) + Between(capital), model = "random", data = dfp)
# 結果を比較してみてみる
screenreg(list(`random model` = model_random,
`hybrid model` = model_hybrid))
============================================
random model hybrid model
--------------------------------------------
value 0.11 *** 0.11 ***
(0.01) (0.01)
capital 0.31 *** 0.31 ***
(0.02) (0.02)
(Intercept) -8.53
(47.52)
Between(value) 0.02
(0.03)
Between(capital) -0.28
(0.19)
--------------------------------------------
R^2 0.77 0.77
Adj. R^2 0.75 0.77
Num. obs. 200 200
s_idios 52.77
s_id 84.20
============================================
*** p < 0.001; ** p < 0.01; * p < 0.05