turn on suggestions

Auto-suggest helps you quickly narrow down your search results by suggesting possible matches as you type.

Showing results for

Find a Community

- Home
- /
- Analytics
- /
- Stat Procs
- /
- How do I visualize fitted values from a mixed-effe...

Topic Options

- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

08-25-2016 02:32 AM

Dear All,

For a logistic mixed-effects model, I am trying to visualize the fitted response on top of the original data points. What confuses me is that the fitted curve seems to run much too low as compared to the (eye-balled) mean of the data points. Could anybody advice me on why this may be happening (or simply point out what I am misunderstanding)?

The syntax is very simple (as implemnted in v.9.4):

proc glimmix data= data1;

class Site;

model Success/Trials = Explanatory /solution link=logit err=binomial;

random Site;

run;

I then extract the intercept (a) and coefficient (b) from solutions and plot the fitted curve at the original (ilink) scale as

y=exp(a+b*Explanatory)/(1+(a+b*Explanatory))

Here, I expect the resultant prediction y to apply to the mean value of the random effect, and hence to – on average – hit the middle of the original data points (expressed as the ratio Success/Trials). But as said, when visually examined, the fitted function seems to run much too low at the ilink scale.

I am now wondering whether I am misunderstanding something basic about the default parameterization of the random effect in glimmix? Or is there some other reason why this should not work as outlined above?

If I am completely off the mark then how would you create the plot described above, with the idea of showing the original data and the fitted function at the original scale, for immediate and intuitive interpretation by the non-expert viewer?

With very many thanks for any advice you may offer,

Tomas

Accepted Solutions

Solution

08-26-2016
09:31 AM

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

08-25-2016 08:26 AM

You can use the OUTPUT statement to compute the predicted values on the inverse link scale. You can create the predicted values accounting for the random effect, but probably you are asking for the marginal predicted values.

For a binomial response, the observed values are proportions. You can compute that quantity in PROC LMIMIX itself, or in a separate DATA Step. In the following, I compute in in the PROC, then use PROC SGPLOT to overlay a scatter plot of the data and a curve for the predicted values.

```
data data1;
input Site Success Trials Explanatory;
datalines;
1 20 30 10
2 10 18 8
3 15 30 5
4 20 55 5
5 5 60 3
11 20 28 9
12 11 16 8
13 14 25 6
14 18 45 4
15 7 54 2
;
proc glimmix data= data1;
class Site;
model Success/Trials = Explanatory /solution link=logit err=binomial;
random Site;
proportion = Success / Trials; /* create new variable for graphing */
output out=PredOut predicted(noblup ilink)=PredMarg /* marginal mean */
/ symbols; /* write proportion variable to output data */
run;
proc sort data=PredOut;
by Explanatory;
run;
proc sgplot data=PredOut;
scatter x=Explanatory y=proportion / datalabel=Site;
series x=Explanatory y=PredMarg;
run;
```

All Replies

Solution

08-26-2016
09:31 AM

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

08-25-2016 08:26 AM

You can use the OUTPUT statement to compute the predicted values on the inverse link scale. You can create the predicted values accounting for the random effect, but probably you are asking for the marginal predicted values.

For a binomial response, the observed values are proportions. You can compute that quantity in PROC LMIMIX itself, or in a separate DATA Step. In the following, I compute in in the PROC, then use PROC SGPLOT to overlay a scatter plot of the data and a curve for the predicted values.

```
data data1;
input Site Success Trials Explanatory;
datalines;
1 20 30 10
2 10 18 8
3 15 30 5
4 20 55 5
5 5 60 3
11 20 28 9
12 11 16 8
13 14 25 6
14 18 45 4
15 7 54 2
;
proc glimmix data= data1;
class Site;
model Success/Trials = Explanatory /solution link=logit err=binomial;
random Site;
proportion = Success / Trials; /* create new variable for graphing */
output out=PredOut predicted(noblup ilink)=PredMarg /* marginal mean */
/ symbols; /* write proportion variable to output data */
run;
proc sort data=PredOut;
by Explanatory;
run;
proc sgplot data=PredOut;
scatter x=Explanatory y=proportion / datalabel=Site;
series x=Explanatory y=PredMarg;
run;
```

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

08-25-2016 01:18 PM

Adding to Rick's valuable response, note that the marginal predicted value corresponds to predicted linear predictor (.e.g, logit) where the random effects are at their theoretical average value (=0). The inverse-link of this predicted value is **not** the predicted probability averaged over the random effects. This is because the inverse link is a non-linear function. The inverse-link does correspond to the predicted probability when the random effects take the value of 0. Thus, the predicted p from the inverse-link function is not the mean of the proportions over the sites.

Read SAS Global Forum articles by Walt Stroup or his excellent 2015 article in Agronomy Journal (open access) for details on this.

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

08-25-2016 01:22 PM

If you want to investigate the alternative that LVM mentions, you can get the predictions that include the random effects with

`predicted(ilink)=Pred`

However, you won't get a nice smooth curve in the fit plot.

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

08-26-2016 09:32 AM

Very many thanks for your quick and valuable pointers!