predicting outcome values with lsmeans in quadratic regression

[palolix]

Dear SAS Communities,

I was trying to predict outcome values for a continuous variable (DM) after fitting a quadratic regression in proc glm

I would like to know the dry matter (DM) values when hedonic (categorical predictor variable) is 5, 5.5, and 6.

I was trying this model but I get this error: Wrong number of AT variable values

proc glm data=one;
model DM= hedonic hedonic*hedonic;
lsmeans/at hedonic = (5  5.5  6) pdiff;
run;

I would greatly appreciate your help!

[Ksharp]

@Rick_SAS wrote a blog about this:
https://blogs.sas.com/content/iml/2014/02/19/scoring-a-regression-model-in-sas.html

[Rick_SAS]

If hedonic is a "categorical predictor variable," you should include a CLASS statement? For a classification variable, you can only evaluate the model at levels that are in the model. If you model hedonic as a continuous variable, you can evaluate the model at values that are not in the input data.

From your description, it sounds like you want the predicted value of the model at specific values of the independent variable. As KSharp says, the simplest way to get that information is to "score" (or evaluate) the model at specific values of hedonic by using the STORE statement to save the model and the PLM procedure to evaluate the model at the specified values. Here is an example that uses the Sashelp.Class data:

proc glm data=sashelp.class plots=none;
/* class age; */ /* ??? Treat as categorical?  */
model Weight= age age*age;
store out=GLMModel;
quit;

data AtLevels;
input Age @@;
datalines;
12 12.5 14 16
;

proc plm restore=GLMModel;
  score data=AtLevels out=GLMOut;
run;

proc print; run;

The input data does not contain any subjects with Age=12.5. If you include the CLASS statement, the model cannot form a prediction at Age=12.5. However, if you treat Age as a continuous variable, the model can calculate the predicted value.

[palolix]

Thank you so much for your reply Rick_SAS and Ksharp! With this approach and also with the following one I get the predicted DM values for the 9 levels of the predictor when using means of the reps, but if I use the raw data (N =555), I get predicted DM values for 555 observations (for the individual reps). Is there a way to get the predicted DM values for the 9 levels of the predictor (hedonic) using the raw data instead of the means?

data poly_data;
set one_means;
hedonic2 = hedonic**2;
hedonic3 = hedonic**3;
run;

proc glm data=poly_data;

    model DM = hedonic hedonic2 hedonic3/p solution alpha=.05 clparm;

run;

Thank you very much!

[Rick_SAS]

/* create fake data set that has the same names as the OP's data */
%let Orig_Data = sashelp.class;
data poly_data;
set &Orig_Data;
DM = Weight;
hedonic = Age;
hedonic2 = hedonic**2;
hedonic3 = hedonic**3;
keep DM hedonic hedonic2 hedonic3;
run;

proc glm data=poly_data;
    model DM = hedonic hedonic2 hedonic3/p solution alpha=.05 clparm;
    output out=GLMOut Predicted=Pred;
quit; 

Of course, any obs with the same value of hedonic will have the same predicted value. If you want to see only the unique predictions, you can sort by hedonic and then output only the first observation for each level:

proc sort data=GLMOut;
by hedonic;
run;

data GLMOut_NoDups;
set GLMOut;
by hedonic;
if first.hedonic;
run;

proc print; run;

[palolix]

It worked, thank you so much!! I guess I have to omit the statements 

DM = Weight 

hedonic = Age 

[palolix]

Thank you so much for your support. One more question. Is there a way to get values for hedonic in increments of .05? To include a statement like this in this code?

proc glm data=poly_data;
model DM = hedonic hedonic2/p solution alpha=.05 clparm;
output out=GLMOut Predicted=Pred;
quit;

proc sort data=GLMOut;
by hedonic;
run;

data GLMOut_NoDups;
set GLMOut;
by hedonic;
do hedonic = 1 to 9 by .5;
if first.hedonic;
run;

[Ksharp · Accepted solution]

You could get any value of "hedonic " you want, make a test dataset to include the value you want to predict.

%let Orig_Data = sashelp.class;
data poly_data;
set &Orig_Data;
DM = Weight;
hedonic = Age;
hedonic2 = hedonic**2;
hedonic3 = hedonic**3;
keep DM hedonic hedonic2 hedonic3;
run;

data test;
do hedonic = 1 to 9 by .5;
 hedonic2 = hedonic**2;
 hedonic3 = hedonic**3;
 output;
end;
run;

data poly_data2;
 length dsname dsn $ 40;
 set poly_data test indsname=dsname;
 dsn=dsname;
run;

proc glm data=poly_data2;
    model DM = hedonic hedonic2 hedonic3/p solution alpha=.05 clparm;
    output out=GLMOut Predicted=Pred;
quit; 

proc print data=GLMOut(where=(dsn='WORK.TEST'));run;
run;

[palolix]

Fantastic, thank you so much Ksharp!

Is there a way to also include the observed DM values in the table?

Thanks!

[Ksharp]

What do you mean by "observed DM values" ?
In TEST dataset , you don't have observed DM, you only have predicted value of DM .
If you want to include observed DM from original dataset, you could remove the WHERE condition in PROC PRINT:

proc print data=GLMOut;run;
run;


Or pick up the value of hedonic you need:
proc print data=GLMOut(where=(hedonic in (1 2 3 4)));run;
run;

[palolix]

Oh ok, thank you for the clarification!

[Rick_SAS]

In addition to what KSharp said in his reply, the observed DM values are not important when scoring the model. The predicted values depend only on the explanatory variable(s). The DM is only used during the fitting processes. It is also used to form residuals, which are useful for assessing the goodness of the fit.

[palolix]

Thank you very much for your input Rick_SAS!