使用R语言处理结构方程模型

在回归问题中，我们有时候需要同时估计多个相互依赖的方程，这种就称为联立方程组。联立方程组在减化为一个方程后称为约简型方程。在这种情况下不能使用普通最小二乘来估计方程系数，因为有的自变量与干扰项存在相关，由此会造成估计量的有偏和不一致。这种自变量称为问题变量。

解决这个问题的方法是采用工具变量方法，即找到一个变量，其特征是它与问题变量相关，而与干扰项无关。用工具变量法得到的自变量系数是一致的。与工具变量法等价的是两阶段最小二乘法，首先用最小二乘从模型的约简型方程中构造出所有问题变量的拟合值，之后，用上一步得出的拟合值来代替问题变量，再用最小二乘对结构方程直接进行回归，就可得到一致估计量。在R语言中为联立方程建模的工具包是systemfit。

下面的算例来自Michael P.Murry的现代计量经济学一书下册。研究人员想对鳕鱼市场的供给与需求建立联立方程模型。供给方程形式为：
PRICE = c+QTY+WINDSPD+RAINY+COLD+STORMY+MIXED+WINDSPD2
需求方程形式为：QTY=c+PRICE+DAY1+DAY2+DAY3+DAY4
其中价格和销量为内生变量，其它和天气和日期有关的为外生变量。所有的外生变量均用来作工具变量。

#首先载入扩展包
library(systemfit)
#读入数据再建立方程
demand <- QTY ~ PRICE+DAY1+DAY2+DAY3+DAY4
supply <- PRICE ~ QTY+WINDSPD+RAINY+COLD+STORMY+MIXED+WINDSPD2
system <- list(demand,supply)
#明确工具变量名称
inst <- ~WINDSPD+WINDSPD2+COLD+RAINY+MIXED+STORMY+DAY1+DAY2+DAY3+DAY4
#利用systemfit命令来获得结果
result2sls <- systemfit(method="2SLS",formula=system,inst=inst)
summary(result2sls)

systemfit results
method: 2SLS

         N  DF     SSR  detRCov   OLS-R2 McElroy-R2
system 222 208 62.2896 0.057822 0.186339   0.236081

      N  DF     SSR      MSE     RMSE       R2   Adj R2
eq1 111 105 49.3550 0.470047 0.685600 0.184330 0.145489
eq2 111 103 12.9346 0.125579 0.354371 0.193912 0.139129

The covariance matrix of the residuals
          eq1       eq2
eq1 0.4700473 0.0347262
eq2 0.0347262 0.1255787

The correlations of the residuals
         eq1      eq2
eq1 1.000000 0.142932
eq2 0.142932 1.000000


2SLS estimates for 'eq1' (equation 1)
Model Formula: QTY ~ PRICE + DAY1 + DAY2 + DAY3 + DAY4
Instruments: ~WINDSPD + WINDSPD2 + COLD + RAINY + MIXED + STORMY + DAY1 +
    DAY2 + DAY3 + DAY4

              Estimate Std. Error  t value   Pr(>|t|)  
(Intercept)  8.5395767  0.1560670 54.71736 < 2.22e-16 ***
PRICE       -0.9337678  0.3452685 -2.70447  0.0079832 **
DAY1        -0.0121607  0.2083902 -0.05836  0.9535764  
DAY2        -0.5259261  0.2023729 -2.59880  0.0106999 *
DAY3        -0.5626912  0.2070420 -2.71776  0.0076899 **
DAY4         0.1000506  0.2028977  0.49311  0.6229655  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.6856 on 105 degrees of freedom
Number of observations: 111 Degrees of Freedom: 105
SSR: 49.354962 MSE: 0.470047 Root MSE: 0.6856
Multiple R-Squared: 0.18433 Adjusted R-Squared: 0.145489


2SLS estimates for 'eq2' (equation 2)
Model Formula: PRICE ~ QTY + WINDSPD + RAINY + COLD + STORMY + MIXED + WINDSPD2
Instruments: ~WINDSPD + WINDSPD2 + COLD + RAINY + MIXED + STORMY + DAY1 +
    DAY2 + DAY3 + DAY4

               Estimate  Std. Error  t value  Pr(>|t|)
(Intercept) -3.08018739  5.87202508 -0.52455 0.6010212
QTY          0.05056356  0.12470245  0.40547 0.6859704
WINDSPD      1.44654086  4.06741515  0.35564 0.7228365
RAINY       -0.00366122  0.09328092 -0.03925 0.9687676
COLD         0.04688691  0.07791593  0.60176 0.5486544
STORMY       0.38896484  0.14312697  2.71762 0.0077155 **
MIXED        0.21230097  0.09611835  2.20875 0.0294083 *
WINDSPD2    -0.22889804  0.70527628 -0.32455 0.7461792
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.354371 on 103 degrees of freedom
Number of observations: 111 Degrees of Freedom: 103
SSR: 12.934605 MSE: 0.125579 Root MSE: 0.354371
Multiple R-Squared: 0.193912 Adjusted R-Squared: 0.139129