A single model may contain a number of linear equations. In such a model, it is often unrealistic to expect that the equation errors would be uncorrelated. A set of equations that has contemporaneous cross-equation error correlation (i.e. the error terms in the regression equations are corrlated) is called a seemingly unrelated regression (SUR) system. At first look, the equations seem unrelated, but the equations are related through the correlation in the errors. The systemfit R package allows a user to specify multiple equations and fit them in an SUR. After doing so, one can perform tests on coefficients across the equations.

We will illustrate SUR using the **hsb2** dataset, predicting **
read** and **math** with the overlapping sets of coefficients and then
comparing some coefficients across the two equations. We will first define our
model equations as R formulas.

First let’s install the systemfit R package which requires installation of three other packages: nnet, mgcv, and quantreg.

install.packages("nnet") install.packages("mgcv") install.packages("quantreg") install.packages("systemfit") install.packages("foreign") install.packages("car") library(foreign) library(systemfit) hsb2 <- read.dta("https://stats.idre.ucla.edu/stat/stata/notes/hsb2.dta") r1 <- read~female + as.numeric(ses) + socst r2 <- math~female + as.numeric(ses) + science

Once the equations have been defined, they can be passed in a list to the **
systemfit** command. A summary of the **systemfit** first shows a summary
of the system (where N = 400), then the separate equations, and then how the
residuals of the two equations are related. These are followed by the OLS fits
of the separate equations.

fitsur <- systemfit(list(readreg = r1, mathreg = r2), data=hsb2) summary(fitsur)systemfit results method: OLS N DF SSR detRCov OLS-R2 McElroy-R2 system 400 392 22835.2 3227.86 0.405103 0.342707 N DF SSR MSE RMSE R2 Adj R2 readreg 200 196 12550.9 64.0351 8.0022 0.400037 0.390854 mathreg 200 196 10284.4 52.4712 7.2437 0.411171 0.402158 The covariance matrix of the residuals readreg mathreg readreg 64.0351 11.4952 mathreg 11.4952 52.4712 The correlations of the residuals readreg mathreg readreg 1.00000 0.19831 mathreg 0.19831 1.00000 OLS estimates for 'readreg' (equation 1) Model Formula: read ~ female + as.numeric(ses) + socst Estimate Std. Error t value Pr(>|t|) (Intercept) 20.6824980 2.9789550 6.94287 5.5019e-11 *** femalefemale -1.5111280 1.1510793 -1.31279 0.19079 as.numeric(ses) 1.2183658 0.8399004 1.45061 0.14849 socst 0.5699327 0.0562967 10.12373 < 2.22e-16 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 8.002195 on 196 degrees of freedom Number of observations: 200 Degrees of Freedom: 196 SSR: 12550.883066 MSE: 64.035118 Root MSE: 8.002195 Multiple R-Squared: 0.400037 Adjusted R-Squared: 0.390854 OLS estimates for 'mathreg' (equation 2) Model Formula: math ~ female + as.numeric(ses) + science Estimate Std. Error t value Pr(>|t|) (Intercept) 19.305181 2.998047 6.43925 9.0557e-10 *** femalefemale 1.160903 1.041641 1.11449 0.266432 as.numeric(ses) 1.399639 0.742390 1.88531 0.060867 . science 0.575330 0.054328 10.58993 < 2.22e-16 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 7.243704 on 196 degrees of freedom Number of observations: 200 Degrees of Freedom: 196 SSR: 10284.364144 MSE: 52.471246 Root MSE: 7.243704 Multiple R-Squared: 0.411171 Adjusted R-Squared: 0.402158

We may be interested in comparing the effect of **female** on **read**,
controlling for **ses** and **socst**, to the effect of **female** on
**math**, controlling for **ses** and **science**. For this, we will
use the** linear.hypothesis** command from the **car** package. To do
this, we create a “restriction” on the model system. We will force the
coefficient of **female** to be the same in both equations and then compare
such a system fit to the one seen when the coefficients are not equal.

library(car) restriction <- "readreg_femalefemale- mathreg_femalefemale" linearHypothesis(fitsur, restriction, test = "Chisq")Linear hypothesis test (Chi^2 statistic of a Wald test) Hypothesis: readreg_femalefemale - mathreg_femalefemale = 0 Model 1: restricted model Model 2: fitsur Res.Df Df Chisq Pr(>Chisq) 1 393 2 392 1 2.9626 0.08521 . --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1