R: [Fwd: R-SIG-Finance Digest, Vol 60, Issue 18]

library("plm")

fdi <- read.csv("C:/data/mydata.csv", na.strings=".")
fdiplm<-plm.data(fdi, index = c("id_code_id", "year"))

fdi_test<-(lfdi_2000~ lagdlfdi+ laglnstock2000+ lagtradegdp +lagdlgdp)

fdi_test_fe <- plm(fdi_test, data=fdiplm, model="within")
fdi_test_re <- plm(fdi_test, data=fdiplm, model="random")

summary (fdi_test_fe)

summary (fdi_test_re)

phtest(fdi_test_re, fdi_test_fe)

Giovani,

Thank you so much for your comments. I am a bit new to R, and to these
mailing lists, so I apologize for being sparse on the details and examples.
I am using Stata 9.2, which might be the answer to my problem, as you
described. I have done quite a bit of internet searching, and did not read
anywhere about the use of a different method for calculating the chi-sq
value, so thanks for that.

 One more issue I have been thinking about. I am assuming your Plm package
knows that the FE is the consistient model, as the same results arrive if
the code is phtest(femod, remod) or phtest(remod, femod). The order does
matter in Stata.

For complteness I am going to post my results using the same Grumfeld
dataset for both stata 9.2 (by hand calculation and canned procedure) and
R.  I am using the Plm package version 1 1-2.

Regards,
Steve



 ## begin Stata9.2 output##
xtreg inv value capital, robust re;

Random-effects GLS regression                   Number of obs      =
200
Group variable (i): firmid                      Number of groups   =
10

R-sq:  within  = 0.7668                         Obs per group: min =
20
       between = 0.8196                                        avg =
20.0
       overall = 0.8061                                        max =
20

Random effects u_i ~ Gaussian                   Wald chi2(3)       =
77.70

corr(u_i, X)       = 0 (assumed)                Prob > chi2        =
0.0000


------------------------------------------------------------------------------
             |               Robust
      invest |      Coef.   Std. Err.      z    P>|z|     [95% Conf.
Interval]

-------------+----------------------------------------------------------------
       value |   .1097811   .0197587     5.56   0.000     .0710547
.1485076
     capital |    .308113   .0418387     7.36   0.000     .2261107
.3901153
       _cons |  -57.83441   24.67795    -2.34   0.019    -106.2023
-9.466507

-------------+----------------------------------------------------------------

     sigma_u |   84.20095
     sigma_e |  52.767964
         rho |  .71800838   (fraction of variance due to u_i)

------------------------------------------------------------------------------

. matrix bfe=e(b);

. matrix vfe=e(V);

. estimates store remod;

. xtreg inv value capital, robust fe;

Fixed-effects (within) regression               Number of obs      =
200
Group variable (i): firmid                      Number of groups   =
10

R-sq:  within  = 0.7668                         Obs per group: min =
20
       between = 0.8194                                        avg =
20.0
       overall = 0.8060                                        max =
20

                                                F(2,188)           =
40.23

corr(u_i, Xb)  = -0.1517                        Prob > F           =
0.0000

------------------------------------------------------------------------------
             |               Robust
      invest |      Coef.   Std. Err.      t    P>|t|     [95% Conf.
Interval]

-------------+----------------------------------------------------------------
       value |   .1101238    .019378     5.68   0.000     .0718975
.1483501
     capital |   .3100653    .042795     7.25   0.000     .2256452
.3944854
       _cons |  -58.74393   23.37422    -2.51   0.013    -104.8534
-12.63449

-------------+----------------------------------------------------------------
     sigma_u |  85.732501
     sigma_e |  52.767964
         rho |  .72525012   (fraction of variance due to u_i)

------------------------------------------------------------------------------

 ###Hausman by hand###

. estimates store femod;

. matrix vre=e(V);

. matrix bre=e(b);

. matrix bdif=bfe-bre;

. matrix list bdif;

bdif[1,3]
         value     capital       _cons
y1  -.00034265  -.00195236   .90952273

. matrix bdifp=bdif';

. matrix dv=vfe-vre;

. matrix dvi=inv(dv);

. matrix list bdif;

bdif[1,3]
         value     capital       _cons
y1  -.00034265  -.00195236   .90952273

. matrix list bdifp;

bdifp[3,1]
                 y1
  value  -.00034265
capital  -.00195236
  _cons   .90952273

. matrix list dvi;

symmetric dvi[3,3]
              value     capital       _cons
  value  -7739.3615
capital   5808.2905   -5305.811
  _cons   3.6641311   .98569198  -.00051157

. matrix chisq=bdif*dvi*bdifp;

. matrix list chisq;

symmetric chisq[1,1]
            y1
y1  -.01956929
###Hausman canned###
.  hausman femod remod;

                 ---- Coefficients ----
             |      (b)          (B)            (b-B)
sqrt(diag(V_b-V_B))
             |     femod        remod        Difference          S.E.

-------------+----------------------------------------------------------------
       value |    .1101238     .1097811        .0003427               .
     capital |    .3100653      .308113        .0019524        .0089965

------------------------------------------------------------------------------
                           b = consistent under Ho and Ha; obtained from
xtreg
            B = inconsistent under Ha, efficient under Ho; obtained from
xtreg

    Test:  Ho:  difference in coefficients not systematic

                  chi2(2) = (b-B)'[(V_b-V_B)^(-1)](b-B)
                          =    -0.01    chi2<0 ==> model fitted on these
                                        data fails to meet the asymptotic
                                        assumptions of the Hausman test;
                                        see suest for a generalized test ##
end Stata9.2 output ##

##begin Output R, using PLM 1.1-2###

test<-data(Grunfeld, package="Ecdat")

fm <- inv~value+capital
femod <- plm(fm, Grunfeld, model="within")
summary(femod)

Oneway (individual) effect Within Model

Call:
plm(formula = fm, data = Grunfeld, model = "within")

Balanced Panel: n=10, T=20, N=200

Residuals :
    Min.  1st Qu.   Median  3rd Qu.     Max.
-184.000  -17.600    0.563   19.200  251.000

Coefficients :
        Estimate Std. Error t-value  Pr(>|t|)
value   0.110124   0.011857  9.2879 < 2.2e-16 ***
capital 0.310065   0.017355 17.8666 < 2.2e-16 ***
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1

Total Sum of Squares:    2244400
Residual Sum of Squares: 523480
F-statistic: 309.014 on 2 and 188 DF, p-value: < 2.22e-16

remod <- plm(fm, Grunfeld, model="random")
summary(remod)

Oneway (individual) effect Random Effect Model
   (Swamy-Arora's transformation)

Call:
plm(formula = fm, data = Grunfeld, model = "random")

Balanced Panel: n=10, T=20, N=200

Effects:
                   var  std.dev share
idiosyncratic 2784.458   52.768 0.282
individual    7089.800   84.201 0.718
theta:  0.86122

Residuals :
   Min. 1st Qu.  Median 3rd Qu.    Max.
-178.00  -19.70    4.69   19.50  253.00

Coefficients :
              Estimate Std. Error t-value Pr(>|t|)
(Intercept) -57.834415  28.898935 -2.0013  0.04536 *
value         0.109781   0.010493 10.4627  < 2e-16 ***
capital       0.308113   0.017180 17.9339  < 2e-16 ***
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1

Total Sum of Squares:    2381400
Residual Sum of Squares: 548900
F-statistic: 328.837 on 2 and 197 DF, p-value: < 2.22e-16

phtest(femod, remod)

        Hausman Test

data:  fm
chisq = 2.3304, df = 2, p-value = 0.3119
alternative hypothesis: one model is inconsistent

###end Plm###





On Mon, May 18, 2009 at 6:01 AM, Millo Giovanni <
Giovanni_Millo at generali.com> wrote:

Dear Steve,

I got your inquiry courtesy of Christian Kleiber, who brought it to our
attention: please next time you post anything re a given package,
include the maintainer's address. We cannot guarantee to parse all the
daily digests of the R system!

Your problem: can you please provide a reproducible example? Else it is
difficult to help, not knowing your data, your results and even the
Stata version you're using.

In the following I replicate what you might have done on a well-known
dataset.

From Stata10, on the usual Grunfeld data taken from package "Ecdat":

## begin Stata10 output ##
. xtreg inv value capital

Random-effects GLS regression                   Number of obs      =
200
Group variable: firm                            Number of groups   =
10

R-sq:  within  = 0.7668                         Obs per group: min =
20
      between = 0.8196                                        avg =
20.0
      overall = 0.8061                                        max =
20

Random effects u_i ~ Gaussian                   Wald chi2(2)       =
657.67
corr(u_i, X)       = 0 (assumed)                Prob > chi2        =
0.0000

------------------------------------------------------------------------
------
        inv |      Coef.   Std. Err.      z    P>|z|     [95% Conf.
Interval]
-------------+----------------------------------------------------------
------
      value |   .1097811   .0104927    10.46   0.000     .0892159
.1303464
    capital |    .308113   .0171805    17.93   0.000     .2744399
.3417861
      _cons |  -57.83441   28.89893    -2.00   0.045    -114.4753
-1.193537
-------------+----------------------------------------------------------
------
    sigma_u |   84.20095
    sigma_e |  52.767964
        rho |  .71800838   (fraction of variance due to u_i)
------------------------------------------------------------------------
------

. estimates store remod

. xtreg inv value capital, fe

Fixed-effects (within) regression               Number of obs      =
200
Group variable: firm                            Number of groups   =
10

R-sq:  within  = 0.7668                         Obs per group: min =
20
      between = 0.8194                                        avg =
20.0
      overall = 0.8060                                        max =
20

                                               F(2,188)           =
309.01
corr(u_i, Xb)  = -0.1517                        Prob > F           =
0.0000

------------------------------------------------------------------------
------
        inv |      Coef.   Std. Err.      t    P>|t|     [95% Conf.
Interval]
-------------+----------------------------------------------------------
------
      value |   .1101238   .0118567     9.29   0.000     .0867345
.1335131
    capital |   .3100653   .0173545    17.87   0.000     .2758308
.3442999
      _cons |  -58.74393   12.45369    -4.72   0.000    -83.31086
-34.177
-------------+----------------------------------------------------------
------
    sigma_u |  85.732501
    sigma_e |  52.767964
        rho |  .72525012   (fraction of variance due to u_i)
------------------------------------------------------------------------
------
F test that all u_i=0:     F(9, 188) =    49.18              Prob > F =
0.0000

. estimates store femod

. hausman femod remod

                ---- Coefficients ----
            |      (b)          (B)            (b-B)
sqrt(diag(V_b-V_B))
            |     femod        remod        Difference          S.E.
-------------+----------------------------------------------------------
------
      value |    .1101238     .1097811        .0003427        .0055213
    capital |    .3100653      .308113        .0019524        .0024516
------------------------------------------------------------------------
------
                          b = consistent under Ho and Ha; obtained from
xtreg
           B = inconsistent under Ha, efficient under Ho; obtained from
xtreg

   Test:  Ho:  difference in coefficients not systematic

                 chi2(2) = (b-B)'[(V_b-V_B)^(-1)](b-B)
                         =        2.33
               Prob>chi2 =      0.3119

.
## end Stata10 output ##

while from plm I get

## begin R putput ##

data(Grunfeld, package="Ecdat")
fm <- inv~value+capital

femod <- plm(fm, Grunfeld)
remod <- plm(fm, Grunfeld, model="random")

phtest(femod, remod)

       Hausman Test

data:  fm
chisq = 2.3304, df = 2, p-value = 0.3119
alternative hypothesis: one model is inconsistent

## end R output ##

which, besides testifying to the goodness and parsimony of an
object-oriented approach as far as screen output is concerned, looks
rather consistent to me.

I cannot but guess that the problem might stem from different RE
estimates: previous versions of Stata used the Wallace-Hussein method by
default for computing the variance of random effects. Now Stata uses
Swamy-Arora, which has been the default of 'plm' since the beginning.
Yet as plm() allows to choose, you can experiment with different values
for the 'random.method' argument in order to see if you get the Stata
result. I suggest you start by comparing the coefficient estimates you
get from Stata and R: FE should be unambiguous, RE might vary as said
above, and for good reason.

You also didn't tell us whether your by-hand calculation agrees with
phtest() output? (I guess it does not)

Please let us know, possibly with a reproducible example and providing
all the above info
Giovanni

PS please also make sure you're not using any VEEEEERY old version of
'plm' (prior to, say, 0.3): these had a bug in the p-value calculation
which made it depend on the order of models compared (so that in the
wrong case you got p.value=1).

Giovanni Millo
Research Dept.,
Assicurazioni Generali SpA
Via Machiavelli 4,
34132 Trieste (Italy)
tel. +39 040 671184
fax  +39 040 671160

----------------------------------------------------------------------
--

Subject:
[R-SIG-Finance] Chi-sq Hausman test---R vs Stata
From:
Steven Archambault <archstevej at gmail.com>
Date:
Sun, 17 May 2009 23:14:13 -0600
To:
r-sig-finance at stat.math.ethz.ch

To:
r-sig-finance at stat.math.ethz.ch


Hi all,

I am running a panel time series regression testing Fixed Effects and
Random Effects. I decided to calculate the chi-sq value for the
Hausman test in both R (Phtest) and Stata. I get different results.
Even within Stata, calculating the Chi-sq value with the canned
procedure or by hand (using
matrices) gives different results. So, the question should come up

there as

well.

Does anybody have any insight on how to pick which results to use? I
guess the one that gives the result I want? Having different programs
give quite different values for the same tests is frustrating me.  I'd

be interested in any feedback folks have!

Thanks,
Steve

      [[alternative HTML version deleted]]