Seemingly Unrelated Regression with Two-way Cluster Robust Errors

[Tianyuan Zhu]

Dear Stata experts,

I would like to follow a paper using seemingly unrelated regression together with two-way cluster robust errors to do my analysis. Here's the description in the paper: "We utilized OLS regression techniques coupled with Seemingly Unrelated Estimation (SUE). ... Therefore, we utilize two-way cluster robust errors (Petersen,2009; Cameron et al., 2011; Thompson, 2011) which cluster around the both firm as well as time."

In Stata, I learned that I can use suest which allows one-way cluster robust errors. I tried to use vce2way (or vcemway) together with suest for seemingly unrelated estimation with two-way cluster robust errors:

reg y1 x1 x2 x3
est store y1

reg y2 x1 x2 x3
est store y2

reg y3 x1 x2 x3
est store y3

vce2way suest y1 y2 y3, cluster(id year)

vcemway suest y1 y2 y3, cluster(id year)

But it always shows the error message "if not allowed r(101); ". So, I am writing to ask how to solve this problem, or whether there is any other way to do it.

Thank you very much!

[Tianyuan Zhu]

Does anyone know how to solve this problem? Thank you.

[Tianyuan Zhu]

I just contacted the author and he said he used suest first, then calculated the two-way clustered errors by hand using the Thompson 2011 approach. https://www.sciencedirect.com/scienc...04405X10001923

But he do think a user written program has been developed since to automate this approach. Does anyone know which program can be used to calculate the two-way clustered errors after suest?

Thank you very much!

[Andrew Musau]

Install reghdfe from SSC. Here is a reproducible example.

sysuse auto, clear
reg price disp gear foreign
reg mpg disp gear foreign
reg weight disp gear foreign
rename (price mpg weight) (dv#), addnumber(1)
reshape long dv, i(make) j(which)
lab def which 1 "price" 2 "mpg" 3 "weight"
label values which which
*ssc install reghdfe
*net install ftools, from("https://github.com/sergiocorreia/ftools/raw/master/src/")
gen cons=1
reghdfe dv i.which#(c.disp c.gear c.foreign c.cons), noabsorb cluster(make trunk)

Res.:

. reg price disp gear foreign

      Source |       SS           df       MS      Number of obs   =        74
-------------+----------------------------------   F(3, 70)        =     18.73
       Model |   282770459         3  94256819.8   Prob > F        =    0.0000
    Residual |   352294937        70  5032784.81   R-squared       =    0.4453
-------------+----------------------------------   Adj R-squared   =    0.4215
       Total |   635065396        73  8699525.97   Root MSE        =    2243.4

------------------------------------------------------------------------------
       price |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
displacement |   25.79084   5.123478     5.03   0.000     15.57238    36.00929
  gear_ratio |  -397.7235   1150.735    -0.35   0.731    -2692.793    1897.346
     foreign |   3749.929   808.4912     4.64   0.000     2137.444    5362.414
       _cons |   1161.034   4175.162     0.28   0.782    -7166.064    9488.131
------------------------------------------------------------------------------

.
. reg mpg disp gear foreign

      Source |       SS           df       MS      Number of obs   =        74
-------------+----------------------------------   F(3, 70)        =     24.11
       Model |  1241.68065         3   413.89355   Prob > F        =    0.0000
    Residual |  1201.77881        70  17.1682687   R-squared       =    0.5082
-------------+----------------------------------   Adj R-squared   =    0.4871
       Total |  2443.45946        73  33.4720474   Root MSE        =    4.1435

------------------------------------------------------------------------------
         mpg |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
displacement |  -.0399061   .0094629    -4.22   0.000    -.0587792   -.0210329
  gear_ratio |   2.225724    2.12537     1.05   0.299    -2.013192     6.46464
     foreign |  -1.501703   1.493256    -1.01   0.318    -4.479909    1.476502
       _cons |   22.90685   7.711385     2.97   0.004     7.526983    38.28673
------------------------------------------------------------------------------

.
. reg weight disp gear foreign

      Source |       SS           df       MS      Number of obs   =        74
-------------+----------------------------------   F(3, 70)        =     95.68
       Model |  35449028.3         3  11816342.8   Prob > F        =    0.0000
    Residual |  8645150.12        70  123502.145   R-squared       =    0.8039
-------------+----------------------------------   Adj R-squared   =    0.7955
       Total |  44094178.4        73  604029.841   Root MSE        =    351.43

------------------------------------------------------------------------------
      weight |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
displacement |   7.131667   .8025978     8.89   0.000     5.530936    8.732398
  gear_ratio |  -24.90265   180.2638    -0.14   0.891    -384.4274    334.6221
     foreign |  -110.2391   126.6509    -0.87   0.387    -362.8364    142.3582
       _cons |   1720.253   654.0431     2.63   0.010     415.8049    3024.701
------------------------------------------------------------------------------

. . reghdfe dv i.which#(c.disp c.gear c.foreign c.cons), noabsorb cluster(make trunk)
(MWFE estimator converged in 1 iterations)
note: 3.which#c.cons omitted because of collinearity

HDFE Linear regression                            Number of obs   =        222
Absorbing 1 HDFE group                            F(  11,     17) =    2175.56
Statistics robust to heteroskedasticity           Prob > F        =     0.0000
                                                  R-squared       =     0.8261
                                                  Adj R-squared   =     0.8170
Number of clusters (make)    =         74         Within R-sq.    =     0.8261
Number of clusters (trunk)   =         18         Root MSE        =  1311.0179

                                    (Std. Err. adjusted for 18 clusters in make trunk)
--------------------------------------------------------------------------------------
                     |               Robust
                  dv |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
---------------------+----------------------------------------------------------------
which#c.displacement |
              price  |   25.79084   3.506061     7.36   0.000     18.39369    33.18798
                mpg  |  -.0399061   .0085348    -4.68   0.000    -.0579128   -.0218993
             weight  |   7.131667   .9023775     7.90   0.000     5.227817    9.035517
                     |
  which#c.gear_ratio |
              price  |  -397.7235   1131.499    -0.35   0.730    -2784.979    1989.532
                mpg  |   2.225724   1.887239     1.18   0.254    -1.756002     6.20745
             weight  |  -24.90265   160.8498    -0.15   0.879    -364.2661    314.4608
                     |
     which#c.foreign |
              price  |   3749.929   1164.626     3.22   0.005     1292.782    6207.076
                mpg  |  -1.501703   1.530279    -0.98   0.340     -4.73031    1.726904
             weight  |  -110.2391   132.0701    -0.83   0.415    -388.8827    168.4045
                     |
        which#c.cons |
              price  |   -559.219   3311.857    -0.17   0.868    -7546.627    6428.189
                mpg  |  -1697.346   613.3852    -2.77   0.013    -2991.475   -403.2163
             weight  |          0  (omitted)
                     |
               _cons |   1720.253   609.3217     2.82   0.012     434.6964    3005.809
--------------------------------------------------------------------------------------

[Tianyuan Zhu]

Install reghdfe from SSC. Here is a reproducible example.

sysuse auto, clear
reg price disp gear foreign
reg mpg disp gear foreign
reg weight disp gear foreign
rename (price mpg weight) (dv#), addnumber(1)
reshape long dv, i(make) j(which)
lab def which 1 "price" 2 "mpg" 3 "weight"
label values which which
*ssc install reghdfe
*net install ftools, from("https://github.com/sergiocorreia/ftools/raw/master/src/")
gen cons=1
reghdfe dv i.which#(c.disp c.gear c.foreign c.cons), noabsorb cluster(make trunk)

Res.:

. reg price disp gear foreign

Source | SS df MS Number of obs = 74
-------------+---------------------------------- F(3, 70) = 18.73
Model | 282770459 3 94256819.8 Prob > F = 0.0000
Residual | 352294937 70 5032784.81 R-squared = 0.4453
-------------+---------------------------------- Adj R-squared = 0.4215
Total | 635065396 73 8699525.97 Root MSE = 2243.4

------------------------------------------------------------------------------
price | Coef. Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
displacement | 25.79084 5.123478 5.03 0.000 15.57238 36.00929
gear_ratio | -397.7235 1150.735 -0.35 0.731 -2692.793 1897.346
foreign | 3749.929 808.4912 4.64 0.000 2137.444 5362.414
_cons | 1161.034 4175.162 0.28 0.782 -7166.064 9488.131
------------------------------------------------------------------------------

.
. reg mpg disp gear foreign

Source | SS df MS Number of obs = 74
-------------+---------------------------------- F(3, 70) = 24.11
Model | 1241.68065 3 413.89355 Prob > F = 0.0000
Residual | 1201.77881 70 17.1682687 R-squared = 0.5082
-------------+---------------------------------- Adj R-squared = 0.4871
Total | 2443.45946 73 33.4720474 Root MSE = 4.1435

------------------------------------------------------------------------------
mpg | Coef. Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
displacement | -.0399061 .0094629 -4.22 0.000 -.0587792 -.0210329
gear_ratio | 2.225724 2.12537 1.05 0.299 -2.013192 6.46464
foreign | -1.501703 1.493256 -1.01 0.318 -4.479909 1.476502
_cons | 22.90685 7.711385 2.97 0.004 7.526983 38.28673
------------------------------------------------------------------------------

.
. reg weight disp gear foreign

Source | SS df MS Number of obs = 74
-------------+---------------------------------- F(3, 70) = 95.68
Model | 35449028.3 3 11816342.8 Prob > F = 0.0000
Residual | 8645150.12 70 123502.145 R-squared = 0.8039
-------------+---------------------------------- Adj R-squared = 0.7955
Total | 44094178.4 73 604029.841 Root MSE = 351.43

------------------------------------------------------------------------------
weight | Coef. Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
displacement | 7.131667 .8025978 8.89 0.000 5.530936 8.732398
gear_ratio | -24.90265 180.2638 -0.14 0.891 -384.4274 334.6221
foreign | -110.2391 126.6509 -0.87 0.387 -362.8364 142.3582
_cons | 1720.253 654.0431 2.63 0.010 415.8049 3024.701
------------------------------------------------------------------------------

. . reghdfe dv i.which#(c.disp c.gear c.foreign c.cons), noabsorb cluster(make trunk)
(MWFE estimator converged in 1 iterations)
note: 3.which#c.cons omitted because of collinearity

HDFE Linear regression Number of obs = 222
Absorbing 1 HDFE group F( 11, 17) = 2175.56
Statistics robust to heteroskedasticity Prob > F = 0.0000
R-squared = 0.8261
Adj R-squared = 0.8170
Number of clusters (make) = 74 Within R-sq. = 0.8261
Number of clusters (trunk) = 18 Root MSE = 1311.0179

(Std. Err. adjusted for 18 clusters in make trunk)
--------------------------------------------------------------------------------------
| Robust
dv | Coef. Std. Err. t P>|t| [95% Conf. Interval]
---------------------+----------------------------------------------------------------
which#c.displacement |
price | 25.79084 3.506061 7.36 0.000 18.39369 33.18798
mpg | -.0399061 .0085348 -4.68 0.000 -.0579128 -.0218993
weight | 7.131667 .9023775 7.90 0.000 5.227817 9.035517
|
which#c.gear_ratio |
price | -397.7235 1131.499 -0.35 0.730 -2784.979 1989.532
mpg | 2.225724 1.887239 1.18 0.254 -1.756002 6.20745
weight | -24.90265 160.8498 -0.15 0.879 -364.2661 314.4608
|
which#c.foreign |
price | 3749.929 1164.626 3.22 0.005 1292.782 6207.076
mpg | -1.501703 1.530279 -0.98 0.340 -4.73031 1.726904
weight | -110.2391 132.0701 -0.83 0.415 -388.8827 168.4045
|
which#c.cons |
price | -559.219 3311.857 -0.17 0.868 -7546.627 6428.189
mpg | -1697.346 613.3852 -2.77 0.013 -2991.475 -403.2163
weight | 0 (omitted)
|
_cons | 1720.253 609.3217 2.82 0.012 434.6964 3005.809
--------------------------------------------------------------------------------------

Hello Andrew, thank you so much for sharing the code! However, I am a bit confused about the following part in the code. I would be very grateful if you could provide a little more detailed explanation about what this part is doing in the SUE analysis. Thank you.
rename (price mpg weight) (dv#), addnumber(1) reshape long dv, i(make) j(which) lab def which 1 "price" 2 "mpg" 3 "weight" label values which which"

[Andrew Musau]

Hello Andrew, thank you so much for sharing the code! However, I am a bit confused about the following part in the code. I would be very grateful if you could provide a little more detailed explanation about what this part is doing in the SUE analysis. Thank you.
rename (price mpg weight) (dv#), addnumber(1) reshape long dv, i(make) j(which) lab def which 1 "price" 2 "mpg" 3 "weight" label values which which"

You have 3 dependent variables to begin with, so you want to have one variable to estimate the 3 models jointly. For this, you rename the variables, reshape and label them so at to identify which coefficient belongs to which model. See

help reshape

[Tianyuan Zhu]

You have 3 dependent variables to begin with, so you want to have one variable to estimate the 3 models jointly. For this, you rename the variables, reshape and label them so at to identify which coefficient belongs to which model. See

help reshape

This is very helpful! Thank you so much!

[Tianyuan Zhu]

I just found out why it always shows the error message "if not allowed r(101)" when using vcemway together with suest. It is simiply because the syntax in suest does not allow [if] whereas vcemway does. So, I just added [if] after syntax [anything] in line 25 of suest.ado (see below) and saved it as a new ado file. Now, the problem is fixed.

 syntax [anything] [if] [,        ///

[Hong Il Yoo]

As we state the -vcemway- help file as well as the background paper, -vcemway- is intended for use with Stata’s “estimation commands”. As you’ll be able to tell from the -suest- help file, -suest- is a “post-estimation” command. Therefore you cannot combine -suest- with -vcemway-.

I am not familiar enough with -suest.ado- to tell whether simply adding -[if]- to the allowed syntax makes it compatible with -vcemway-; the fact that you no longer see a syntax error message does not imply that combining -vcemway- with your hacked version of -suest- delivers suitably clustered standard errors. At the very minimum, I suggest that you check if -vcemway suest y1 y2 y3, cluster(id id)- and -vcemway suest y1 y2 y3, cluster(year year)- deliver the same results as -suest y1 y2 y3, cluster(id)- and -suest y1 y2 y3, cluster(year)-, respectively.

The whole point of using -suest- here is to obtain a covariance matrix that you can use to test differences across OLS regression coefficients. If I had a task like this, I would follow Andrew's data transformation approach. Once the data have been reshaped as described, you can apply -reghdfe, cluster(id year)-, -vcemway reg, cluster(id year)- or any other OLS regression command that supports two-way clustering.

[Tianyuan Zhu]

As we state the -vcemway- help file as well as the background paper, -vcemway- is intended for use with Stata’s “estimation commands”. As you’ll be able to tell from the -suest- help file, -suest- is a “post-estimation” command. Therefore you cannot combine -suest- with -vcemway-.

I am not familiar enough with -suest.ado- to tell whether simply adding -[if]- to the allowed syntax makes it compatible with -vcemway-; the fact that you no longer see a syntax error message does not imply that combining -vcemway- with your hacked version of -suest- delivers suitably clustered standard errors. At the very minimum, I suggest that you check if -vcemway suest y1 y2 y3, cluster(id id)- and -vcemway suest y1 y2 y3, cluster(year year)- deliver the same results as -suest y1 y2 y3, cluster(id)- and -suest y1 y2 y3, cluster(year)-, respectively.

The whole point of using -suest- here is to obtain a covariance matrix that you can use to test differences across OLS regression coefficients. If I had a task like this, I would follow Andrew's data transformation approach. Once the data have been reshaped as described, you can apply -reghdfe, cluster(id year)-, -vcemway reg, cluster(id year)- or any other OLS regression command that supports two-way clustering.

Thank you so much, Dr. Yoo! I will follow your suggestions.

[Tianyuan Zhu]

As we state the -vcemway- help file as well as the background paper, -vcemway- is intended for use with Stata’s “estimation commands”. As you’ll be able to tell from the -suest- help file, -suest- is a “post-estimation” command. Therefore you cannot combine -suest- with -vcemway-.

I am not familiar enough with -suest.ado- to tell whether simply adding -[if]- to the allowed syntax makes it compatible with -vcemway-; the fact that you no longer see a syntax error message does not imply that combining -vcemway- with your hacked version of -suest- delivers suitably clustered standard errors. At the very minimum, I suggest that you check if -vcemway suest y1 y2 y3, cluster(id id)- and -vcemway suest y1 y2 y3, cluster(year year)- deliver the same results as -suest y1 y2 y3, cluster(id)- and -suest y1 y2 y3, cluster(year)-, respectively.

The whole point of using -suest- here is to obtain a covariance matrix that you can use to test differences across OLS regression coefficients. If I had a task like this, I would follow Andrew's data transformation approach. Once the data have been reshaped as described, you can apply -reghdfe, cluster(id year)-, -vcemway reg, cluster(id year)- or any other OLS regression command that supports two-way clustering.

Hello Dr. Hoo, I just finished the checks you suggested. First, -vcemway suest y1 y2 y3, cluster(id id)/cluster(year year) - and -suest y1 y2 y3, cluster(id)/cluster(year) - delivered the same results. Second, the robust standard errors obtained from -vcemway suest y1 y2 y3, cluster(id year)- are different from those from -suest y1 y2 y3, cluster(id) and -suest y1 y2 y3, cluster(year). Third, the robust standard errors from -vcemway suest y1 y2 y3, cluster(id year) are very close to those using -reghdfe- and -vcemway-, though they are not exactly the same. I listed a few values in the following table.

	Robust Std err
	suest+vcemway	reghdfe+vcemway
x1	0.0039227	0.0039372
x2	0.0001816	0.0001823
x3	0.0077337	0.0078306
x4	0.2423195	0.2472624

[Hong Il Yoo]

interesting, thanks for sharing the results.

[Tianyuan Zhu]

Hello Dr. Hoo, I just finished the checks you suggested. First, -vcemway suest y1 y2 y3, cluster(id id)/cluster(year year) - and -suest y1 y2 y3, cluster(id)/cluster(year) - delivered the same results. Second, the robust standard errors obtained from -vcemway suest y1 y2 y3, cluster(id year)- are different from those from -suest y1 y2 y3, cluster(id) and -suest y1 y2 y3, cluster(year). Third, the robust standard errors from -vcemway suest y1 y2 y3, cluster(id year) are very close to those using -reghdfe- and -vcemway-, though they are not exactly the same. I listed a few values in the following table.

	Robust Std err
	suest+vcemway	reghdfe
x1	0.0039227	0.0039372
x2	0.0001816	0.0001823
x3	0.0077337	0.0078306
x4	0.2423195	0.2472624

A correction: For -reghdfe-, I did not use -vcemway- as its -cluster- option allows for two-way/multi-way cluster robust errors.

[Hong Il Yoo]

A correction: For -reghdfe-, I did not use -vcemway- as its -cluster- option allows for two-way/multi-way cluster robust errors.

I see, thanks for the clarification. I thought what you meant was that -reghdfe y x, cluster(id time) noabsorb- gave you the same results as -vcemway regress y x, cluster(id time)- (they should give the same results).

I went back to -vcemway.ado- and I can see -vcemway- with your "if-hacked" version of -suest- must be doing exactly doing what you were looking for. That is, the combination gives you what you would have obtained if you used -suest- three times by clustering on id, time and group(id time) in turn, and then "manually" computed a 2-way clustered covariance matrix. Since -vcemway- cannot apply the efficient route of applying -_robust, custer(.)- when the command in question is -suest-, it will jump to the less efficient but fail-safe route of actually executing -suest, cluster(.)- three times. Not sure when it is going to be, but next time I update -vcemway-, I will revise it a bit so that it works with the unhacked version of -suest- as well as the usual battery of estimation commands. Thanks for bringing this application to my attention.

[Tianyuan Zhu]

I see, thanks for the clarification. I thought what you meant was that -reghdfe y x, cluster(id time) noabsorb- gave you the same results as -vcemway regress y x, cluster(id time)- (they should give the same results).

I went back to -vcemway.ado- and I can see -vcemway- with your "if-hacked" version of -suest- must be doing exactly doing what you were looking for. That is, the combination gives you what you would have obtained if you used -suest- three times by clustering on id, time and group(id time) in turn, and then "manually" computed a 2-way clustered covariance matrix. Since -vcemway- cannot apply the efficient route of applying -_robust, custer(.)- when the command in question is -suest-, it will jump to the less efficient but fail-safe route of actually executing -suest, cluster(.)- three times. Not sure when it is going to be, but next time I update -vcemway-, I will revise it a bit so that it works with the unhacked version of -suest- as well as the usual battery of estimation commands. Thanks for bringing this application to my attention.

Thank you so much, Dr. Hoo. Look forward to the new version of -vcemway.