Stata FAQ
Part 2: How can I understand a 3-way continuous interaction? (Stata 11)

If you have reached this page without reading Part 1, you need to go back and review it because it has a lot of explanation that is not duplicated on this page. Go back now.

Here is the regression model from the previous page with two covariates V1 and V2 added. 

Y = b0 + b1X + b2Z + b3W + b4XZ + b5XW + b6ZW + b7XZW + b8V1 + b9V2

Adding covariates to your model does not change the formulas for the simple slopes as long as the covariates are not interacted with the independent variable or the moderator variables. The values of the slopes themselves will change because there are additional variables in the model. The formulas for the intercepts (constants), however, will have to change because the covariates are included in those equations. The equation for the intercept for the condition 1 (HzHw) with covariates becomes, 

intercept 1 =  b0 + b2Hz + b3Hw + b6HzHw + b8m1 + b9m2

Where m1 and m2 are the mean values for the two covariates and all of the other values are as before.

We will add two covariates to our previous model, socst which will be renamed as v1 and ses which will be renamed v2. Everything is this section is the same as the previous page with the exception of adding the covariates to the end of the regress command. 

use http://www.ats.ucla.edu/stat/stata/notes/hsb2, clear

rename write y
rename read x
rename math w
rename science z
rename socst v1
rename ses v2

regress y c.x##c.z##c.w v1 v2

      Source |       SS       df       MS              Number of obs =     200
-------------+------------------------------           F(  9,   190) =   25.36
       Model |  9757.27668     9  1084.14185           Prob > F      =  0.0000
    Residual |  8121.59832   190  42.7452543           R-squared     =  0.5457
-------------+------------------------------           Adj R-squared =  0.5242
       Total |   17878.875   199   89.843593           Root MSE      =   6.538

------------------------------------------------------------------------------
           y |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
           x |  -3.463037   1.716525    -2.02   0.045    -6.848932   -.0771426
           z |  -2.902201   1.670106    -1.74   0.084    -6.196533    .3921302
             |
     c.x#c.z |   .0624411     .03093     2.02   0.045     .0014308    .1234514
             |
           w |  -3.696588   1.831366    -2.02   0.045    -7.309009   -.0841671
             |
     c.x#c.w |   .0796921   .0341529     2.33   0.021     .0123245    .1470597
             |
     c.z#c.w |   .0696334   .0332851     2.09   0.038     .0039775    .1352893
             |
 c.x#c.z#c.w |  -.0013838   .0005968    -2.32   0.021    -.0025611   -.0002065
             |
          v1 |   .2789791    .057697     4.84   0.000     .1651703     .392788
          v2 |  -.9056562   .6914623    -1.31   0.192    -2.269585    .4582727
       _cons |   185.3134   88.56438     2.09   0.038      10.6177    360.0092
------------------------------------------------------------------------------
Note that with the inclusion the covariates that the three-way interaction is even "more" significant than before.

We will begin by computing the intercepts (constants) for the simple regression equations. Once again we turn to the margins command which is the same as the previous example with the addition of the atmeans options to hold the covariates constant at their means. As before we use global macros for the range of x and for the high and low values of w and z. 

quietly sum x
global hi=r(max)
global lo=r(min)
quietly sum w
global Hw=r(mean)+r(sd)
global Lw=r(mean)-r(sd)
quietly sum z
global Hz=r(mean)+r(sd)
global Lz=r(mean)-r(sd)


/* compute intercepts */

margins, at( x=(0) w=($Hw $Lw) z=($Hz $Lz)) atmeans ///
         asbalanced noesample noestimcheck force vsquish

Adjusted predictions                              Number of obs   =        200
Model VCE    : OLS

Expression   : Linear prediction, predict()
1._at        : x               =           0
               z               =    61.75089
               w               =    62.01345
               v1              =      52.405 (mean)
               v2              =       2.055 (mean)
2._at        : x               =           0
               z               =    61.75089
               w               =    43.27655
               v1              =      52.405 (mean)
               v2              =       2.055 (mean)
3._at        : x               =           0
               z               =    41.94911
               w               =    62.01345
               v1              =      52.405 (mean)
               v2              =       2.055 (mean)
4._at        : x               =           0
               z               =    41.94911
               w               =    43.27655
               v1              =      52.405 (mean)
               v2              =       2.055 (mean)

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         _at |
          1  |   56.27364   5.809015     9.69   0.000     44.88818     67.6591
          2  |     44.969   7.499049     6.00   0.000     30.27114    59.66687
          3  |   28.23421   10.18945     2.77   0.006     8.263259    48.20515
          4  |   42.76522   5.457104     7.84   0.000     32.06949    53.46095
------------------------------------------------------------------------------

/* put intercepts in matrix i */

mat i=r(b)

mat list i

i[1,4]
            1.         2.         3.         4.
          _at        _at        _at        _at
r1  56.273644  44.969003  28.234207   42.76522
Next, we will compute the simple slopes for the various combinations of w and z. The margins command is exactly the same as for the example without the covariates because the covariances are not involved in the 3-way interaction. After computing the simple slopes we will place the values in matrix named s. 
/* simple slopes */

margins, dydx(x) at(w=($Hw $Lw) z=($Hz $Lz)) post vsquish

Average marginal effects                          Number of obs   =        200
Model VCE    : OLS

Expression   : Linear prediction, predict()
dy/dx w.r.t. : x
1._at        : z               =    61.75089
               w               =    62.01345
2._at        : z               =    61.75089
               w               =    43.27655
3._at        : z               =    41.94911
               w               =    62.01345
4._at        : z               =    41.94911
               w               =    43.27655

------------------------------------------------------------------------------
             |            Delta-method
             |      dy/dx   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
x            |
         _at |
          1  |   .0356648   .0951734     0.37   0.708    -.1508717    .2222012
          2  |   .1435565   .1371306     1.05   0.295    -.1252146    .4123276
          3  |    .498484   .1884961     2.64   0.008     .1290384    .8679296
          4  |   .0929561   .1236172     0.75   0.452    -.1493291    .3352413
------------------------------------------------------------------------------

/* put slopes into matrix s */

mat s=e(b)

mat list s

s[1,4]
            x:         x:         x:         x:
            1.         2.         3.         4.
          _at        _at        _at        _at
y1  .03566477  .14355649  .49848401  .09295611
Now we have enough information to plot the simple slopes. Please note: You must use the names y and x with the twoway function plot command. Do not use variable names from your data. 
twoway  (function y=i[1,1] + s[1,1]*x, range($lo $hi))  ///
        (function y=i[1,2] + s[1,2]*x, range($lo $hi))  ///
        (function y=i[1,3] + s[1,3]*x, range($lo $hi))  ///
        (function y=i[1,4] + s[1,4]*x, range($lo $hi))  ///
        (scatter y x, msym(oh) jitter(3)),      ///
        legend(order(1 "1 HzHw" 2 "2 HzLw" 3 "3 LzHw" 4 "4 LzLw")) /// 
        xtitle(X) ytitle(Y) scheme(lean1)
Once again it looks like slope 3 is the one that is most different from the others.

Now, we are ready compute the tests of differences in simple slopes. 

/* differences in simple slopes */

lincom [x]1bn._at - [x]2bn._at  /* a) HwHz vs HzLw (1 & 2) */

 ( 1)  [x]1bn._at - [x]2._at = 0

------------------------------------------------------------------------------
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         (1) |  -.1078917   .1559956    -0.69   0.489    -.4136375     .197854
------------------------------------------------------------------------------

lincom [x]1bn._at - [x]3bn._at  /* b) HwHz vs LzHw (1 & 3) */

 ( 1)  [x]1bn._at - [x]3._at = 0

------------------------------------------------------------------------------
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         (1) |  -.4628192   .1955057    -2.37   0.018    -.8460033   -.0796352
------------------------------------------------------------------------------

lincom [x]1bn._at - [x]4bn._at  /* c) HwHz vs LwLz (1 & 4) */

 ( 1)  [x]1bn._at - [x]4._at = 0

------------------------------------------------------------------------------
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         (1) |  -.0572913   .1557822    -0.37   0.713    -.3626188    .2480361
------------------------------------------------------------------------------

lincom [x]2bn._at - [x]3bn._at  /* d) HzLw vs LzHw (2 & 3) */

 ( 1)  [x]2._at - [x]3._at = 0

------------------------------------------------------------------------------
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         (1) |  -.3549275   .2449055    -1.45   0.147    -.8349335    .1250785
------------------------------------------------------------------------------

lincom [x]2bn._at - [x]4bn._at  /* e) HzLw vs LwLz (2 & 4) */

 ( 1)  [x]2._at - [x]4._at = 0

------------------------------------------------------------------------------
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         (1) |   .0506004   .1635762     0.31   0.757    -.2700031    .3712038
------------------------------------------------------------------------------

lincom [x]3bn._at - [x]4bn._at  /* f) LzHw vs LwLz (3 & 4) */

 ( 1)  [x]3._at - [x]4._at = 0

------------------------------------------------------------------------------
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         (1) |   .4055279   .2096843     1.93   0.053    -.0054457    .8165015
------------------------------------------------------------------------------
With some cutting and pasting we can put together a table including each of the tests of differences in simple slopes. 
   Slopes      |      Diff.   Std. Err.      t    P>|t|     [95% Conf. Interval]
---------------+----------------------------------------------------------------
(HzHw vs HzLw) |  -.1078917   .1559956    -0.69   0.489    -.4136375     .197854
(HzHw vs LzHw) |  -.4628192   .1955057    -2.37   0.018    -.8460033   -.0796352 
(HzLw vs LzLw) |  -.0572913   .1557822    -0.37   0.713    -.3626188    .2480361
(LzHw vs LzLw) |  -.3549275   .2449055    -1.45   0.147    -.8349335    .1250785
(HzHw vs LzLw) |   .0506004   .1635762     0.31   0.757    -.2700031    .3712038
(HzLw vs LzHw) |   .4055279   .2096843     1.93   0.053    -.0054457    .8165015
Inspecting the table above it appears slope 1 versus slope 3 (HzHw vs LzHw) is the only one of the tests that is significant. But before we reach that conclusion we will need to adjust the p-value to take into account that these are post-hoc tests and there are a total of six of them. With a Bonferroni correction the p-value becomes .019*6 = .114, which is not statistically significant. The Bonferroni correction may not be the best choice because it is overly conservative but it is certainly easy to implement.

As before, we have collected all of the commands into a do-file to make it easier to use. 

use http://www.ats.ucla.edu/stat/stata/notes/hsb2, clear

rename write y
rename read x
rename math w
rename science z
rename socst v1
rename ses v2

regress y c.x##c.z##c.w v1 v2

quietly sum x
global hi=r(max)
global lo=r(min)
quietly sum w
global Hw=r(mean)+r(sd)
global Lw=r(mean)-r(sd)
quietly sum z
global Hz=r(mean)+r(sd)
global Lz=r(mean)-r(sd)

/* intercepts */
margins, at( x=(0) w=($Hw $Lw) z=($Hz $Lz)) atmeans ///
         asbalanced noesample noestimcheck force vsquish
mat i=r(b)

/* simple slopes */
margins, dydx(x) at(w=($Hw $Lw) z=($Hz $Lz)) post vsquish
mat s=e(b)

/* graph the simple slopes */
twoway  (function y=i[1,1] + s[1,1]*x, range($lo $hi))  ///
        (function y=i[1,2] + s[1,2]*x, range($lo $hi))  ///
        (function y=i[1,3] + s[1,3]*x, range($lo $hi))  ///
        (function y=i[1,4] + s[1,4]*x, range($lo $hi))  ///
        (scatter y x, msym(oh) jitter(3)),      ///
        legend(order(1 "1 HzHw" 2 "2 HzLw" 3 "3 LzHw" 4 "4 LzLw")) /// 
        xtitle(X) ytitle(Y) scheme(lean1)

/* differences in simple slopes */
lincom [x]1bn._at - [x]2bn._at
lincom [x]1bn._at - [x]3bn._at
lincom [x]1bn._at - [x]4bn._at
lincom [x]2bn._at - [x]3bn._at
lincom [x]2bn._at - [x]4bn._at
lincom [x]3bn._at - [x]4bn._at
References

Dawson, J.F. & Richter, A.W. (2004) A significance test of slope differences for three-way interactions in moderated multiple regression analysis. RP0422. Aston University, Birmingham, UK

How to cite this page

Report an error on this page or leave a comment

The content of this web site should not be construed as an endorsement of any particular web site, book, or software product by the University of California.