R version 2.10.1 (2009-12-14)
Copyright (C) 2009 The R Foundation for Statistical Computing
ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> ## ----------------------------------------------------------------------------
> ## Greene (2008)
> ##
> ## data
> ## availability: electronic
> ## firms: General Motors, US Steel, General Electric, Chrysler 
> ## errors: none
> ##
> ## analysis
> ## result: OLS point estimates reproducible (minor deviation for GM_value) 
> ##         Standard errors differ (Greene has ML standard errors, minor dev. for GE_const)
> ##         SUR point estimates reproducible (minor deviation for GE_const)
> ##         Standard errors fully reproducible
> ##         ML point estimates reproducible (minor deviation for US_capital)
> ##         Standard errors fully reproducible 
> ## ----------------------------------------------------------------------------
> 
> ## preliminaries
> source("start.R")
Loading required package: kinship
Loading required package: survival
Loading required package: splines
Loading required package: nlme
Loading required package: lattice
[1] "kinship is loaded"
Loading required package: Formula
Loading required package: MASS
Loading required package: sandwich
Loading required package: zoo
Loading required package: Matrix
Loading required package: car
Loading required package: lmtest
> 
> ## data pre-processing
> gr <- subset(Grunfeld, firm %in% c("General Motors", "US Steel", "General Electric", "Chrysler"))
> gr$firm <- factor(gr$firm)
> pgr <- plm.data(gr, c("firm", "year"))
> 
> 
> ## SUR estimation (Example 16.8, pp. 534-536)
> 
> ## Table 16.4 (p. 536)
> ## OLS
> sf_OLS <- systemfit(invest ~ value + capital, data = pgr, method = "OLS")
> gsummary(sf_OLS, digits = 5)

Method: OLS 
                 (Intercept)   value capital     R^2   sigma^2
Chrysler            -6.18996 0.07795 0.31572 0.91358  176.3203
                    13.50648 0.01997 0.02881                  
General.Electric    -9.95631 0.02655 0.15169 0.70531  777.4463
                    31.37425 0.01557 0.02570                  
General.Motors    -149.78245 0.11928 0.37144 0.92135 8423.8751
                   105.84212 0.02583 0.03707                  
US.Steel           -49.19832 0.17486 0.38964 0.47086 9299.6040
                   148.07537 0.07420 0.14237                  
> ## correlation matrix for the 4 OLS residual vectors is
> summary(sf_OLS)$residCor
                    Chrysler General.Electric General.Motors   US.Steel
Chrysler          1.00000000      -0.06792257     -0.2729521  0.3379535
General.Electric -0.06792257       1.00000000      0.2792929  0.4281071
General.Motors   -0.27295213       0.27929288      1.0000000 -0.2614613
US.Steel          0.33795352       0.42810705     -0.2614613  1.0000000
> ## covariance matrix for OLS residual vectors is
> summary(sf_OLS)$residCov * (17/20)
                   Chrysler General.Electric General.Motors   US.Steel
Chrysler          149.87222        -21.37565      -282.7564   367.8402
General.Electric  -21.37565        660.82939       607.5331   978.4503
General.Motors   -282.75642        607.53314      7160.2939 -1967.0464
US.Steel          367.84024        978.45025     -1967.0464  7904.6634
> 
> ## alternatively
> lm_GM <- lm(invest ~ value + capital, data = gr, subset = firm == "General Motors")
> lm_US <- lm(invest ~ value + capital, data = gr, subset = firm == "US Steel")
> lm_GE <- lm(invest ~ value + capital, data = gr, subset = firm == "General Electric")
> lm_CH <- lm(invest ~ value + capital, data = gr, subset = firm == "Chrysler")
> 
> ## FGLS
> sf_FGLS <- systemfit(invest ~ value + capital, data = pgr, method = "SUR",
+   control = systemfit.control(methodResidCov = "noDfCor"))
> gsummary(sf_FGLS, digits = 5)

Method: SUR 
                 (Intercept)   value capital     R^2   sigma^2
Chrysler             0.93655 0.06785 0.31465 0.91206  179.4081
                    11.59044 0.01705 0.02606                  
General.Electric   -19.71195 0.03464 0.13685 0.69594  802.1461
                    26.58279 0.01279 0.02250                  
General.Motors    -160.68640 0.12051 0.38002 0.92098 8463.7099
                    90.40985 0.02187 0.03311                  
US.Steel            21.16398 0.13037 0.44852 0.45621 9557.1076
                   116.17517 0.05737 0.12251                  
>   
> ## ML via iterated FGLS
> sf_ML <- systemfit(invest ~ value + capital, data = pgr, method = "SUR", 
+   control = systemfit.control(methodResidCov = "noDfCor", maxiter = 7, tol = 1e-10))  
> gsummary(sf_ML, digits = 5)

Method: SUR 
                 (Intercept)   value capital     R^2   sigma^2
Chrysler             2.58096 0.06564 0.31373 0.91126  181.0458
                    11.54392 0.01698 0.02617                  
General.Electric   -24.09904 0.03808 0.13109 0.68683  826.1905
                    25.80157 0.01217 0.02223                  
General.Motors    -179.40723 0.12481 0.38018 0.92053 8512.3018
                    86.66046 0.02086 0.03266                  
US.Steel            36.46299 0.12437 0.43675 0.45468 9584.0462
                   106.18099 0.05191 0.11711                  
> ## ML estimate of Sigma is
> summary(sf_ML)$residCov
                    Chrysler General.Electric General.Motors   US.Steel
Chrysler          153.888927         2.516339      -325.4111   427.0102
General.Electric    2.516339       702.261962       615.1667  1288.6454
General.Motors   -325.411067       615.166673      7235.4566 -2455.1368
US.Steel          427.010242      1288.645379     -2455.1368  8146.4393
>