The PLM Procedure

(View the complete code for this example.)

This example conducts an analysis similar to Example 15 in Examples: GLIMMIX Procedure: Examples: GLIMMIX Procedure. It uses simulated data to perform multiple comparisons among predicted values in a model with group-specific trends that are modeled through regression splines. The estimable functions are formed using nonpositional syntax with constructed effects. Consider the data in the following DATA step. Each of the 100 observations for the continuous response variable y is associated with one of two groups.

data spline;
   input group y @@;
   x = _n_;
   datalines;
 1    -.020 1    0.199 2    -1.36 1    -.026
 2    -.397 1    0.065 2    -.861 1    0.251
 1    0.253 2    -.460 2    0.195 2    -.108
 1    0.379 1    0.971 1    0.712 2    0.811
 2    0.574 2    0.755 1    0.316 2    0.961
 2    1.088 2    0.607 2    0.959 1    0.653
 1    0.629 2    1.237 2    0.734 2    0.299
 2    1.002 2    1.201 1    1.520 1    1.105
 1    1.329 1    1.580 2    1.098 1    1.613
 2    1.052 2    1.108 2    1.257 2    2.005
 2    1.726 2    1.179 2    1.338 1    1.707
 2    2.105 2    1.828 2    1.368 1    2.252
 1    1.984 2    1.867 1    2.771 1    2.052
 2    1.522 2    2.200 1    2.562 1    2.517
 1    2.769 1    2.534 2    1.969 1    2.460
 1    2.873 1    2.678 1    3.135 2    1.705
 1    2.893 1    3.023 1    3.050 2    2.273
 2    2.549 1    2.836 2    2.375 2    1.841
 1    3.727 1    3.806 1    3.269 1    3.533
 1    2.948 2    1.954 2    2.326 2    2.017
 1    3.744 2    2.431 2    2.040 1    3.995
 2    1.996 2    2.028 2    2.321 2    2.479
 2    2.337 1    4.516 2    2.326 2    2.144
 2    2.474 2    2.221 1    4.867 2    2.453
 1    5.253 2    3.024 2    2.403 1    5.498
;

The following statements fit a model with separate trends for the two groups; the trends are modeled as B-splines.

proc orthoreg data=spline;
   class group;
   effect spl = spline(x);
   model y = group spl*group / noint;
   store ortho_spline;
   title 'B-splines Comparisons';
run;

Results from this analysis are shown in Output 88.6.1. The "Parameter Estimates" table shows the estimates for the spline coefficients in the two groups.

By default, the ORTHOREG procedure constructs B-splines with seven knots. Since B-spline coefficients satisfy a sum-to-one constraint and since the model contains group-specific intercepts, the last spline coefficient for each group is redundant and is set to 0.

The following statements make a prediction for the input data set by using the SCORE statement with PROC PLM and graph the observed and predicted values in the two groups:

proc plm restore=ortho_spline;
   score data=spline out=ortho_pred predicted=p;
run;

proc sgplot data=ortho_pred;
   series  y=p x=x / group=group name="fit";
   scatter y=y x=x / group=group;
   keylegend "fit" / title="Group";
run;

The prediction plot in Output 88.6.2 suggests that there is some separation of the group trends for small values of x and for values that exceed about x = 40.

Output 88.6.2: Observed Data and Predicted Values by Group

In order to determine the range on which the trends separate significantly, the PLM procedure is executed in the following statements with an ESTIMATE statement that applies group comparisons at a number of values for the spline variable x:

%macro GroupDiff;
   %do x=0 %to 75 %by 5;
      "Diff at x=&x" group 1 -1 group*spl [1,1  &x] [-1,2  &x],
   %end;
   'Diff at x=80' group 1 -1 group*spl [1,1 80] [-1,2 80]
%mend;

proc plm restore=ortho_spline;
   show effects;
   estimate %GroupDiff / adjust=simulate seed=1 stepdown;
run;

For example, the following ESTIMATE statement compares the trends between the two groups at x = 25:

estimate 'Diff at x=25' group 1 -1 group*spl [1,1 25] [-1,2 25];

The nonpositional syntax is used for the group*spl effect. For example, the specification requests that the spline be computed at x = 25 for the second level of variable group. The resulting coefficients are added to the vector for the estimate after being multiplied with –1.

Because comparisons are made at a large number of values for x, a multiplicity correction is in order to adjust the p-values to reflect familywise error control. Simulated p-values with step-down adjustment are used here.

Output 88.6.3 displays the "Store Information" for the item store and information about the spline effect (the result of the SHOW statement).

Output 88.6.4 displays the results from the ESTIMATE statement.

Notice that the "Store Information" in Output 88.6.3 displays the classification variables (from the CLASS statement in PROC ORTHOREG), the constructed effects (from the EFFECT statement in PROC ORTHOREG), and the model effects (from the MODEL statement in PROC ORTHOREG). Output 88.6.4 shows that at the 5% significance level the trends are significantly different for and for . Between those values you cannot reject the hypothesis of trend congruity.

To see this effect more clearly, you can filter the results by adding the a filtering statement to the previous PROC PLM step:

proc plm restore=ortho_spline;
   estimate %GroupDiff / adjust=simulate seed=1 stepdown;
   filter adjp > 0.05;
run;

This produces Output 88.6.5, which displays the subset of the results in Output 88.6.4 that meets the condition in the FILTER expression.