Changes between Initial Version and Version 1 of PanelLinearModel01


Ignore:
Timestamp:
Nov 14, 2016, 10:15:31 AM (5 years ago)
Author:
admin
Comment:

--

Legend:

Unmodified
Added
Removed
Modified
  • PanelLinearModel01

    v1 v1  
     1R言語による多変量時系列分析。
     2
     3複数グループ・複数項目の時系列データで時系列間の関係性・影響を明らかにする。
     4
     5
     6以下、実行可能なR言語ソースコードも用い紹介する。
     7
     8
     9例:複数時系列間の関係性・影響
     10
     11農家ごとの農地の肥沃度・人の各仕事の労働量・各肥料の投入量の各時系列、および、生産量の時系列から、肥沃度・各労働量・各肥料量が生産量に与える影響を明らかにする。
     12国ごとの各政策種類への投資額の時系列、および、GDP から、各政策投資がGDPに与える影響を明らかにする
     13地域ごとの各キャンペーン種類の投入量の各時系列、および、商品売上の時系列から、各地域の各キャンペーンが売上に与える影響を明らかにする。等。
     14
     15R Library
     16
     17複数時系列の時系列分析として、R言語では plm 、 ccgarch 等のパッケージがあるが、今回は plm (Panel Linear Model) を紹介する。
     18
     19※参考:Panel Data Econometrics in R: The plm Package (PDF)
     20
     21
     22Panel Data(交差時系列)
     23
     24計量経済学の領域では、交差系列 (Cross section Data) で時系列 (Time series Data) である系列を
     25
     26『交差時系列(Panel Data)』と呼ぶ。 
     27
     28{{{
     29例)47都道府県の人口を時間を追って調べたもの。
     30
     31交差系列 (Cross section Data):同一時点での様々なData。例) ある時点で47都道府県の人口、人口密度、男女比などを調べたもの。
     32時系列 (Time series Data):同一種類のDataを様々な時点で取ったもの。例) ある都道府県の人口を時間を追って調べたもの。
     33※参考:Wikipedia: 計量経済学 - Wikipedia、パネルデータの意義とその活用(PDF)
     34
     35
     36関数
     37
     38plm package には複数時系列の評価関数として次の4つが用意されている。
     39
     40
     41Functions:
     42
     43plm: estimation of the basic panel models, i.e. within, between and random e
     44ffect models. Models are estimated using the lm function to transformed data,
     45pvcm: estimation of models with variable coecients,
     46pgmm: estimation of generalized method of moments models,
     47pggls: estimation of general feasible generalized least squares models
     48これらはRの線形回帰関数 lm() と同じ引数の形をとり、1引数目に formula、2引数目にdata をとる。また、Option として、次の3つを指定できる。
     49
     50
     51Options:
     52
     53index: this argument enables the estimation functions to identify the structure of the data, i.e. the individual and the time period for each observation,
     54effect: the kind of e
     55ffects to include in the model, i.e. individual e
     56ffects, time e
     57ffects or both,
     58model: the kind of model to be estimated, most of the time a model with
     59xed e
     60ffects or a model with random e
     61ffects.
     62
     63Results:
     64
     65上記各関数の結果は 関数名と同じクラスのオブジェクトに保存され、これらのクラスは panelmodel オブジェクトを継承している。これらのオブジェクトは以下を含んでいる。
     66
     67coefficients
     68residuals
     69fitted.values
     70vcov
     71df.residual
     72call
     73
     74データ
     75
     76Rの組み込みデータ「Produc: 米国の48州の1970-1986のパネルデータ」を用いる
     77
     78library("plm")
     79data("Produc", package = "plm")
     80head(Produc, 10)
     81Produc データ (10行目まで抜粋):
     82
     83> head(Produc, 10)
     84     state year     pcap     hwy   water    util       pc   gsp    emp unemp
     851  ALABAMA 1970 15032.67 7325.80 1655.68 6051.20 35793.80 28418 1010.5   4.7
     862  ALABAMA 1971 15501.94 7525.94 1721.02 6254.98 37299.91 29375 1021.9   5.2
     873  ALABAMA 1972 15972.41 7765.42 1764.75 6442.23 38670.30 31303 1072.3   4.7
     884  ALABAMA 1973 16406.26 7907.66 1742.41 6756.19 40084.01 33430 1135.5   3.9
     895  ALABAMA 1974 16762.67 8025.52 1734.85 7002.29 42057.31 33749 1169.8   5.5
     906  ALABAMA 1975 17316.26 8158.23 1752.27 7405.76 43971.71 33604 1155.4   7.7
     917  ALABAMA 1976 17732.86 8228.19 1799.74 7704.93 50221.57 35764 1207.0   6.8
     928  ALABAMA 1977 18111.93 8365.67 1845.11 7901.15 51084.99 37463 1269.2   7.4
     939  ALABAMA 1978 18479.74 8510.64 1960.51 8008.59 52604.05 39964 1336.5   6.3
     9410 ALABAMA 1979 18881.49 8640.61 2081.91 8158.97 54525.86 40979 1362.0   7.1
     95state(州)、year(年)、gsp (州内総生産)、pcap (民間資本ストック)、pc(公的資本)、emp (労働投入量)
     96
     97※参考:Rパッケージガイドブック 蓮見亮さん記事。
     98
     99
     100Rでの実装例
     101
     102library("plm")
     103
     104# Producデータ
     105data("Produc", package = "plm")
     106data <- Produc
     107
     108# index に用いる グループ列名、時間列名
     109grpCName <- "state"
     110timeCName <- "year"
     111
     112# pdata.frame (panel data frame) 作成。グループ列, 時間列で index 作成 
     113idx <- c(grpCName, timeCName)
     114pdata <- pdata.frame(data, index = idx, drop.index = TRUE, row.names = TRUE)
     115
     116# 目的変数列 指定
     117tidx <- 6
     118
     119###################################
     120# Panel Linear Model 作成
     121###################################
     122# plm(formula, data, subset, na.action, effect = c("individual","time","twoways"),
     123#    model = c("within","random","ht","between","pooling","fd"),
     124#    random.method = c("swar","walhus","amemiya","nerlove"),
     125#    inst.method = c("bvk","baltagi"), index = NULL, ...)
     126#
     127#[model 指定]:
     128# within: Fixed Effect Model [Default]
     129#   (固定効果モデル。個別効果を完全除去。
     130#    目的変数、説明変数、残差のそれぞれをグループごとの平均を回帰式から引き算)
     131# pooling: Pooling Model     (グループ、時間を気にせず回帰)
     132# fd: First-Difference Model (グループごとに回帰。変数の時間差分を用いる)
     133# between: Between Model    (グループごとに回帰。変数の時間平均を用いる)
     134# random: error components model
     135#   (変量効果モデル。残差をグループ間のものとグループ内のものに分割)。
     136#    ramdom.method 変動性の指定ができ swar [Default], walhus, amemiya, nerlove を指定可能。
     137#[effect 指定]:
     138# effect = "time" : 上記処理をグループ毎の処理に置き換える (グループと時間を入れ替え)
     139
     140# Formula作成:目的変数以外を説明変数とする
     141cn <- colnames(pdata)
     142fmla <- as.formula(paste(cn[tidx], " ~ ", paste(cn[-tidx], collapse = " + ")))
     143
     144# 固定効果モデル
     145res.fe <-  plm(formula = fmla, data = pdata, model = "within")
     146#(res.fe <-  plm(formula = fmla, data = data, index = idx, model = "within")) と同じ
     147summary(res.fe)
     148
     149# 変動効果モデル
     150res.re <-  plm(formula = fmla, data = pdata, model = "random")
     151summary(res.re)
     152
     153###################################
     154# 各種検定
     155###################################
     156
     157# 各種豊富な検定関数が用意されている。
     158
     159library(car)
     160# Tests for individual and time effect
     161plmtest(fmla, data = pdata, effect = "twoways", type = "ghm")
     162# Tests of poolability
     163pooltest(fmla, data = pdata, model = "within")
     164# F tests of e^Kffects based on the comparison of the within and the pooling models
     165pFtest(fmla, data = pdata, effect = "twoways")
     166# Hausman test based on the comparison of two sets of estimates
     167phtest(res.fe, res.re)
     168# Tests of serial correlatio
     169pwtest(fmla, data = pdata)
     170# Locally robust tests for serial correlation or random effects
     171pbsytest(fmla, data = pdata, test = "j")
     172#Conditional LM test for AR(1) or MA(1) errors under random effect
     173pbltest(fmla, data = data, alternative = "onesided") #pdata.frame ではなく元 data.frame を用いる
     174
     175library(lmtest)
     176#General serial correlation test
     177pbgtest(res.fe, order = 2)
     178#Wooldridge's test for serial correlation in short FE panels
     179pwartest(fmla, data = pdata)
     180#Wooldridge's first-difference-based test
     181pwfdtest(fmla, data = pdata)
     182
     183# 他にも Cross-sectional dependence の検定も用意されている。
     184
     185実行結果
     186
     187残差、回帰係数、R^2、等、各種回帰結果が得られている。
     188
     189> # 固定効果モデル
     190> res.fe <-  plm(formula = fmla, data = pdata, model = "within")
     191> summary(res.fe)
     192Oneway (individual) effect Within Model
     193
     194Call:
     195plm(formula = fmla, data = pdata, model = "within")
     196
     197Balanced Panel: n=48, T=17, N=816
     198
     199Residuals :
     200     Min.   1st Qu.    Median      Mean   3rd Qu.      Max.
     201-1.54e+04 -7.66e+02 -7.04e+01  1.64e-09  7.71e+02  1.68e+04
     202
     203Coefficients :
     204         Estimate  Std. Error t-value  Pr(>|t|)   
     205pcap   1.1832e+04  1.5029e+04  0.7873   0.43136   
     206hwy   -1.1833e+04  1.5029e+04 -0.7873   0.43134   
     207water -1.1831e+04  1.5029e+04 -0.7872   0.43141   
     208util  -1.1832e+04  1.5029e+04 -0.7873   0.43134   
     209pc     1.2001e-01  2.1647e-02  5.5437 4.086e-08 ***
     210emp    3.4932e+01  8.8757e-01 39.3571 < 2.2e-16 ***
     211unemp -1.6813e+02  6.0568e+01 -2.7759   0.00564 **
     212---
     213Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
     214
     215Total Sum of Squares:    1.3289e+11
     216Residual Sum of Squares: 4.978e+09
     217R-Squared      :  0.96254
     218      Adj. R-Squared :  0.89766
     219F-statistic: 2793.44 on 7 and 761 DF, p-value: < 2.22e-16
     220
     221> # 変動効果モデル
     222> res.re <-  plm(formula = fmla, data = pdata, model = "random")
     223summary(res.re)
     224> Oneway (individual) effect Random Effect Model
     225   (Swamy-Arora's transformation)
     226
     227Call:
     228plm(formula = fmla, data = pdata, model = "random")
     229
     230Balanced Panel: n=48, T=17, N=816
     231
     232Effects:
     233                   var  std.dev share
     234idiosyncratic  6541433     2558 0.185
     235individual    28839985     5370 0.815
     236theta:  0.8853 
     237
     238Residuals :
     239    Min.  1st Qu.   Median  3rd Qu.     Max.
     240-11500.0   -847.0     61.7    849.0  19800.0
     241
     242Coefficients :
     243               Estimate  Std. Error t-value  Pr(>|t|)   
     244(Intercept) -1.7738e+03  1.1556e+03 -1.5350    0.1252   
     245pcap         1.4200e+04  1.5615e+04  0.9094    0.3634   
     246hwy         -1.4201e+04  1.5615e+04 -0.9094    0.3634   
     247water       -1.4199e+04  1.5615e+04 -0.9093    0.3634   
     248util        -1.4201e+04  1.5615e+04 -0.9094    0.3634   
     249pc           1.3118e-01  1.7636e-02  7.4380 2.613e-13 ***
     250emp          3.4150e+01  7.6819e-01 44.4552 < 2.2e-16 ***
     251unemp       -2.7338e+02  5.9394e+01 -4.6028 4.840e-06 ***
     252---
     253Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
     254
     255Total Sum of Squares:    1.8368e+11
     256Residual Sum of Squares: 5743500000
     257R-Squared      :  0.96873
     258      Adj. R-Squared :  0.95923
     259F-statistic: 3576.08 on 7 and 808 DF, p-value: < 2.22e-16
     260
     261検定結果例
     262
     263> # Tests for individual and time effect
     264> plmtest(fmla, data = pdata, effect = "twoways", type = "ghm")
     265
     266        Lagrange Multiplier Test - two-ways effects (Gourieroux, Holly and
     267        Monfort)
     268
     269data:  fmla
     270chisq = 3156.201, df = 2, p-value < 2.2e-16
     271alternative hypothesis: significant effects
     272
     273参考資料:時系列分析入門資料
     274
     275Tokyo.R で講師をした時系列分析の入門資料です。[R][データマイニング] R言語による時系列分析
     276
     277}}}