There are 2 issues here. The first is not related to your solution, but to the reason of the problem. Using your example, with a random statement, if you change the ref from self to team for example, the type 3 results are different, as well as the results from lsmeans (the proportion and standard error for class which appear in both outputs, and the proportions do not add up to 1). However, when running multinomial regression without a random effect, the results are identical. This is the source of the problem. Even fixing covariance parameters or changing estimation method did not work. I find it difficult to interpret, when the results differs depending on the choice of the reference, while it does not affect the results in absence of a random effect. proc glimmix data=school; freq Count; class School Program(ref=first) style; model Style(ref = "self")= program / dist=mult link=glogit ; random school /group = style; store log; run; proc plm restore=log; lsmeans program / ilink e plots=none; run; proc glimmix data=school; freq Count; class School Program(ref=first) style; model Style(ref = "team")= program / dist=mult link=glogit ; random school /group = style; run; It is because of this problem that your solution was interesting, to get the % for the reference group that can add up to 100%. Your solution worked fine with my dataset & model with main effects, but not when I try to get it for an interaction (with an interaction term in plm, even if I had more columns to the cont dataset). It takes a really long time to run and after I get an error. ERROR: Error generating code. However, with your example, I’m able to reproduce the lsmeans results with %nlmeans. proc glimmix data=school; freq Count; class School Program(ref=first) style; model Style(order=data)=School| program / dist=mult link=glogit solution; store log; run; proc plm restore=log; lsmeans school*program / ilink e plots=none; ods output coef=c; run; data cont; length label $40; infile datalines missover; input label k1-k18; set=1; datalines; P(self,S1,regular) 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 P(team,S2,afternoon) 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 ;run; %nlmeans(instore=log, link=glogit, coef=c, contrasts=cont) So it is possible to make it work, I just need to figure out what goes wrong with my own data & model. EDIT: The ERROR: Error generating code. happens when in %NLEstimates in the section " proc nlmixed data = _tParm;". This was resolved by simplifying the model (I had lots of double & triple interactions). If I keep only some of them, the macro runs fine. I suppose that otherwise they are to many parameters to estimate. So, %nlmeans can work nicely to get proportions for the reference category, for main effets as well as interactions, but only for "small" models. Writing down the data "cont" can also become quite long for models with lots of terms.
... View more
