You could use the LSMEANS statement in PROC GLM. For example,
lsmeans x1*x2 / slice=x1 diff;
You could get the test for the x2 effect for each level of x1 for y1score to y4score respectively.
The problem here is with your data. Many of the LSMEANS are nonestimable, so the LSMEANS statement above will probably not produce anything useful.
The reasons for nonestimable LSMEANS in your case might be --
1. There is an empty cell for x1*x2 which makes some lsmeans non-estimable. (I did not look at all other two way interactions so there might be other interaction terms with empty cells too).
2. You have confoundings in your data/model, so some effects have 0 df and therefore cannot be tested. This is indicated in the Log --
NOTE: H Matrix for x1*x2 has zero d.f.
NOTE: H Matrix for x2*x4 has zero d.f.
As a matter of fact, if you add the E option in the MODEL statement you will see the General Form of Estimable Functions table, which confirms the confounding in your data/model.
You might want to remove the effects that are confounded with other effects in the model. For example,
PROC GLM DATA = data1; CLASS x1 x2 x3 x4; MODEL y1scores y2scores y3scores y4scores = x1 x2 x3 x4 x1*x3 x1*x4 x2*x3 x3*x4 ; MANOVA H = _ALL_ ; lsmeans x1*x3 / slice=x1 diff; lsmeans x1*x4 / slice=x1 diff; run; quit;
... View more