2.7 Two-Way Analysis of Variance
Let us start by creating a copy of program effort and
grouping it into categories 0-4, 5-14 and 15+. This time
we will make the copying, recoding and labeling in just
one command, using effort_g for effort in groups.
. recode effort (0/4=1 "Weak") (5/14=2 "Moderate") /// > (15/max=3 "Strong"), gen(effort_g) label(effort_g) (20 differences between effort and effort_g)
Here's a table showing steeper declines in countries with strong programs, with a smaller difference between weak and moderate:
. tabulate effort_g, summarize(change)
rECODE of |
effort | Summary of % Change in CBR between
(Program | 1965 and 1975
Effort) | Mean Std. Dev. Freq.
------------+------------------------------------
Weak | 5 4 7
Moderate | 9.3333333 7.393691 6
Strong | 27.857143 6.3358391 7
------------+------------------------------------
Total | 14.3 11.810343 20
An Additive Model
Let us create dumy variables for moderate and strong programs.
. gen effort_mod = effort_g==2 // or effort >=5 & effort < 15 . gen effort_str = effort_g==3 // or effort >=15 & !missing(effort)
The comments show how we could create these dummies directly from the original variable, bypassing the creation of the grouped factor.
We are now ready to fit the additive model with weak programs in low setting countries as the reference cell:
. regress change setting_med setting_high effort_mod effort_str
Source | SS df MS Number of obs = 20
-------------+------------------------------ F( 4, 15) = 13.55
Model | 2075.80829 4 518.952073 Prob > F = 0.0001
Residual | 574.39171 15 38.2927807 R-squared = 0.7833
-------------+------------------------------ Adj R-squared = 0.7255
Total | 2650.2 19 139.484211 Root MSE = 6.1881
------------------------------------------------------------------------------
change | Coef. Std. Err. T P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
setting_med | -1.680829 3.854967 -0.44 0.669 -9.897497 6.535839
setting_high | 2.387565 4.45653 0.54 0.600 -7.111304 11.88643
effort_mod | 3.836269 3.574561 1.07 0.300 -3.782727 11.45527
effort_str | 20.6715 4.339232 4.76 0.000 11.42265 29.92036
_cons | 5.379275 3.10526 1.73 0.104 -1.23943 11.99798
------------------------------------------------------------------------------
Compare these estimates with the results in Table 2.15 and the anova with Table 2.16.
Countries with strong family planning programs show steeper declines than countries with weak programs at the same level of social setting, on average 21 percentage points more.
This statement is based on the assumption of additivity, namely that the difference in outcomes across categories of program effort is the same at every level of setting. We will test this assumption below.
Factor Variables
Before we do that, let us show how we can get exactly the
same results, and admitedly with less effort, using
factor variables. We simply use the categorical
versions of our predictors with the i. prefix
to instruct Stata to treat them as factors:
. regress change i.setting_g i.effort_g
Source | SS df MS Number of obs = 20
-------------+------------------------------ F( 4, 15) = 13.55
Model | 2075.80829 4 518.952073 Prob > F = 0.0001
Residual | 574.39171 15 38.2927807 R-squared = 0.7833
-------------+------------------------------ Adj R-squared = 0.7255
Total | 2650.2 19 139.484211 Root MSE = 6.1881
------------------------------------------------------------------------------
change | Coef. Std. Err. T P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
setting_g |
2 | -1.680829 3.854967 -0.44 0.669 -9.897497 6.535839
3 | 2.387565 4.45653 0.54 0.600 -7.111304 11.88643
|
effort_g |
2 | 3.836269 3.574561 1.07 0.300 -3.782727 11.45527
3 | 20.6715 4.339232 4.76 0.000 11.42265 29.92036
|
_cons | 5.379275 3.10526 1.73 0.104 -1.23943 11.99798
------------------------------------------------------------------------------
You just have to remember that setting 2 and 3 are medium and high, and effort 2 and 3 are moderate and strong programs.
Fitted Values
Let us reproduce Table 2.17 in the notes, showing fitted means
by setting and effort. I will use predict to
generate fitted values, and then simply tabulate them by the
two relevant factors:
. predict anovafit
(option xb assumed; fitted values)
. label var anovafit "Two-way anova fit"
. tabulate setting_g effort_g , summarize(anovafit) means
Means of Two-way anova fit
Social | RECODE of effort (Program
Setting | Effort)
(Grouped) | Weak Moderate Strong | Total
-----------+---------------------------------+----------
Low | 5.3792748 9.2155437 . | 7.5714285
Medium | 3.6984456 7.5347152 24.369947 | 8.5999999
High | 7.7668395 11.603108 28.438341 | 23.749999
-----------+---------------------------------+----------
Total | 5.0000001 9.3333331 27.857142 | 14.3
Can you get the missing cell in the upper right corner? How about the margins?
A Two-Factor Interaction
Let us now consider a model with an interaction between social setting and program effort, so differences in fertility decline by effort will vary by setting. To do this "by hand" we need to create four dummy variables. At this point it becomes hard to create reasonable names in 12 characters or less. I'll use one character for each variable and three for each category. An alternative is to use camelCasing, which saves the spaces used by the underscores.
. gen se_med_mod = setting_med * effort_mod
. gen se_med_str = setting_med * effort_str
. gen se_hi_mod = setting_high * effort_mod
. gen se_hi_str = setting_high * effort_str
.
. regress change setting_med setting_high effort_mod effort_str ///
> se_med_mod se_med_str se_hi_mod se_hi_str
note: se_hi_str omitted because of collinearity
Source | SS df MS Number of obs = 20
-------------+------------------------------ F( 7, 12) = 8.15
Model | 2189.45 7 312.778571 Prob > F = 0.0009
Residual | 460.75 12 38.3958333 R-squared = 0.8261
-------------+------------------------------ Adj R-squared = 0.7247
Total | 2650.2 19 139.484211 Root MSE = 6.1964
------------------------------------------------------------------------------
change | Coef. Std. Err. T P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
setting_med | 3.333333 5.05937 0.66 0.522 -7.690086 14.35675
setting_high | 6.333333 7.155029 0.89 0.393 -9.256136 21.9228
effort_mod | 8.583333 4.732607 1.81 0.095 -1.728132 18.8948
effort_str | 19.33333 6.692917 2.89 0.014 4.75072 33.91595
se_med_mod | -14.58333 8.578579 -1.70 0.115 -33.27445 4.107784
se_med_str | -.3333333 9.797427 -0.03 0.973 -21.68009 21.01343
se_hi_mod | -6.583333 9.959379 -0.66 0.521 -28.28296 15.11629
se_hi_str | (omitted)
_cons | 2.666667 3.577515 0.75 0.470 -5.128068 10.4614
------------------------------------------------------------------------------
Oops, Stata dropped a term. Why? Because there are no countries with strong programs in low settings, so we have only eight groups, but are trying to represent their means using nine parameters, which is obviously one too many.
Fortunately this doesn't affect testing. We can use the sum of squares to build the hierarchical anova in Table 2.21, and Stata can test the interaction for us, dropping automatically the redundant term:
. test se_med_mod se_med_str se_hi_mod se_hi_str
( 1) se_med_mod = 0
( 2) se_med_str = 0
( 3) se_hi_mod = 0
( 4) o.se_hi_str = 0
Constraint 4 dropped
F( 3, 12) = 0.99
Prob > F = 0.4318
We have no evidence that the differences by effort vary with social setting.
This makes the issue of interpreting parameters moot, but it may be worth noting briefly the problems caused by the empty cell. As things stand, the coefficient of moderate effort compares moderate with weak at low setting, but the coefficient of strong effort compares strong with weak at high setting. (Table 2.20 in the notes may help see this point. When the term for high and strong is dropped, the only difference between weak and strong programs at high setting is the coefficient of strong.)
The parametrization I like best for this problem combines the main effects of effort with the interactions, so we obtain differences between strong and weak, and between moderate and weak programs, at each level of setting. This allows us to omit the difference between strong and weak programs at low setting, which is the one we can't identify. I will not pursue this tack at this point, but will return to the idea when we have to interpret a significant interaction.
Factor Interactions
We can fit the same model using factor variables, we
just need to use the i. prefix to indicate that
grouped setting and effort should be treated as factors, not
covariates, combined with a hash # to indicate
interactions.
Stata understands a single or double hash, typed without spaces after a factor specification:
A single hash,
as in i.setting#i.effort, takes the first cell
in the cross-tabulation of setting and effort groups
as the baseline and creates eight dummies, one for each
of the other cells.
In this parametrization each combination of levels is compared
directly to the baseline.
A double hash,
as in i.setting##i.effort,
creaes main effects and interactions, and corresponds to the
parametrization we used above. Recall that in this case
"main" effects are really differences between categories of
one factor when the other factor is at the baseline, and the
"interactions" are additional differences when the
other factor is not at the baseline.
In our example we have the added complication that one of the cells, (1,3) is empty. Stata detects this and drops the last term:
. regress change i.setting_g##i.effort_g
note: 1b.setting_g#3.effort_g identifies no observations in the sample
note: 3.setting_g#3.effort_g omitted because of collinearity
Source | SS df MS Number of obs = 20
-------------+------------------------------ F( 7, 12) = 8.15
Model | 2189.45 7 312.778571 Prob > F = 0.0009
Residual | 460.75 12 38.3958333 R-squared = 0.8261
-------------+------------------------------ Adj R-squared = 0.7247
Total | 2650.2 19 139.484211 Root MSE = 6.1964
------------------------------------------------------------------------------
change | Coef. Std. Err. T P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
setting_g |
2 | 3.333333 5.05937 0.66 0.522 -7.690086 14.35675
3 | 6.333333 7.155029 0.89 0.393 -9.256136 21.9228
|
effort_g |
2 | 8.583333 4.732607 1.81 0.095 -1.728132 18.8948
3 | 19.33333 6.692917 2.89 0.014 4.75072 33.91595
|
setting_g#|
effort_g |
1 3 | (empty)
2 2 | -14.58333 8.578579 -1.70 0.115 -33.27445 4.107784
2 3 | -.3333333 9.797427 -0.03 0.973 -21.68009 21.01343
3 2 | -6.583333 9.959379 -0.66 0.521 -28.28296 15.11629
3 3 | (omitted)
|
_cons | 2.666667 3.577515 0.75 0.470 -5.128068 10.4614
------------------------------------------------------------------------------
Unfortunately we can't control which of the dummies is dropped. All we can do is change the baseline, or combine categories, but that is not particularly helpful.
Continue with 2.8 Analysis of Covariance Models