topic Re: Marginal Effects with Proc glimmix in Statistical Procedures

Marginal Effects with Proc glimmix

Jerem — Sun, 16 Jun 2013 14:39:09 GMT

Hi everybody!

I'm quite new with SAS procedures, so my question might be sounds easy for a lot of you but I couldn't find the answer so far....

Here is my problem, I'm trying to calculate the marginal effects at the means of my independent variables. For this research, it was indeed decided not to use odds ratio but the marginal effects.

I'm using the proc glimmix to integrate the multi-level dimension of my data. Here is the code I wrote:

proc glimmix data=candidates;

class commune ;

model binouverture (event='1') = EurostatA EurostatC sizelog ENP

/ solution dist=binary link=logit ddfm=satterth ;

random intercept / subject= commune solution;

run;

Based on this, what would be the easiest procedure to calculated the marginal effects?

Many thanks for your insights!

Jerem

Re: Marginal Effects with Proc glimmix

SteveDenham — Mon, 17 Jun 2013 13:10:32 GMT

The code you have gives the conditional effects. I assume that there are multiple measurements for each level of the variable 'commune'. If so, then the marginal values may be obtained by treating this as a repeated measures design using the following code:

proc glimmix data=candidates;

class commune ;

model binouverture (event='1') = EurostatA EurostatC sizelog ENP

/ solution dist=binary link=logit ddfm=satterth ;

random intercept / subject= commune solution residual; /* Converts to a repeated measures design */

run;

Do you wish to get marginal estimates for each level of commune? If so, then I would suggest:

proc glimmix data=candidates;

class commune ;

model binouverture (event='1') = EurostatA EurostatC sizelog ENP commune commune*EurostatA commune*EurostatC commune*sizelog commune*ENP

/ solution dist=binary link=logit ddfm=kenwardrogers ; /* Addresses individual commune levels, and adds the Kenward-Rogers adjustment for repeated measures designs */

random intercept / subject= commune solution residual; /* Converts to a repeated measures design */

lsmeans commune/ilink ci; /* Outputs the estimates for each level of commune, on both the logit and original scale */

run;

I hope this is helpful.

Steve Denham

Re: Marginal Effects with Proc glimmix

Jerem — Mon, 17 Jun 2013 14:14:02 GMT

Dear Steve,

Thanks very much for your quick response! It already helps although I'm not familiar enough with multi-level regression to be sure that I'm on the good tracks.

So, just to be sure, here is my intention:

-I've got 1012 observations (=1012 electoral lists). Those 1012 were present in 262 different 'communes'. In other words, in each municipality(=communes), you have a certain amount of electoral lists: sometimes only a couple, sometimes more. I'm seeking to use a logistic regression explaining whether or not these electoral lists include a certain kind of candidates (If yes=1; no=0). My independent variables are : EurostatA EurostatC sizelog ENP as written in the model.

Because lists are embedded in distinct communes. I want to include the multi-level structure of my data. For that purpose, I was told to use a fixed effects model.

Furthermore, to report my results, I'm asked to present the marginal effects of the independent variables on the dependent variable.

Is that possible with the solution you proposed?

Many thanks in advance!

Jeremy

Re: Marginal Effects with Proc glimmix

SteveDenham — Mon, 17 Jun 2013 14:42:55 GMT

So an individual record would look something like:

commune list binouverture eurostata eurostatc sizelog enp

The list variable ranges from one or two to some larger number, commune is indexed so that there are 262 different municipalities.

Now comes the question, with several possiblities: Are the municipalities a sample that you wish to use to infer to an entire population where you know how many sampling units exist? Are the municipalities a sample that represent a population that could be considered infinite in some sense (i.e., a population of all POSSIBLE units? Are the municipalities the ONLY ones for which you wish to draw conclusions?

One of these is where your research question is directed, and each would be addressed by a different kind of analysis in SAS. I await your answer. If the second scenario is the most likely, then the model I presented earlier is appropriate for obtaining marginal estimates. If either of the other scenarios applies, then the analysis should be changed.

Steve Denham

Re: Marginal Effects with Proc glimmix

Jerem — Mon, 17 Jun 2013 15:48:09 GMT

Actually, the observation is the electoral list itself (I did not code individual candidates). So for each of 1012 lists, I coded the information for :

-commune (one of the 262 municipalities),

-Eurostat (the municipality is either eurostatA, eurostatB, eurostatC. Reference used is EurostatB),

-sizelog (log of the population of the municipality),

-enp (effective number of parties in the municipality).

-binouverture of the list (Yes/no),

In other words, for each list from the same municipality, the variables 'commune' (that identifies the municipality), 'EurostatA', EurostatC', 'ENP', 'sizelog' is always the same while the dependent variable 'binouverture' varies for each list.

Regarding the key questions you asked: the 262 municipalities is the population. Actually, we coded information for all the 1012 lists presented in all the 262 municipalities at the last local elections. I don't want to infer my results to other case studies but evaluate the impact of independent variables on 'binouverture' for that particular election. To some extent, it's an electoral report, not theory-building. Thus statistical significance is not important because my sample is my population of interest in that paper.

The goal is thus to assess the marginal effect of each of the independent variables 'Eurostata' 'Eurostatc' 'Sizelog' 'Enp' on the dependent variable 'binouverture' but taking into account the fact that lists are presented in 262 municipalities.

Does it make clearer? Thanks for your help Steve!

Jeremy

Re: Marginal Effects with Proc glimmix

SteveDenham — Mon, 17 Jun 2013 16:05:31 GMT

Excellent answers. OK, everything is a fixed effect, and you have multiple measurements on the 262 municipalities.

proc glimmix data=candidates;

class commune ;

model binouverture (event='1') = EurostatA EurostatC sizelog ENP commune commune*EurostatA commune*EurostatC commune*sizelog commune*ENP

/ solution dist=binary link=logit ddfm=kenwardrogers ; /* Addresses individual commune levels, and adds the Kenward-Rogers adjustment for repeated measures designs */

random intercept / subject= commune solution residual; /* Converts to a repeated measures design. Even though this is a RANDOM statement, the residual option means treat the values within a commune as correlated. Note also the inclusion of commune and commune by covariate interactions in the model statement. */

lsmeans commune/ilink ci; /* Outputs the estimates for each level of commune, on both the logit and original scale */

run;

Be careful about the interpretation of the solution. The values under commune*eurostatA for instance are deviations from the slope estimated by the value under EurostatA, for each commune. The standard reference level is the last level of commune, so these values are the differences in slopes between each commune and the reference.

The least squares means here are the expected marginal means, at the mean level of EurostatA, EurostatC, sizelog, and ENP. You may want to estimate marginal means for specific values of these four continuous covariates. For that, see the documentation of the LSMEANS statement and the AT= option.

Steve Denham

Re: Marginal Effects with Proc glimmix

Jerem — Mon, 17 Jun 2013 19:11:36 GMT

I've got my answers regarding the model now, thanks!

I submitted the model but it is still charging. Is proc glimmix a time-consuming procedure?

By the way, for the line lsmeans commune/ilink ci SAS proposed to change ci by cl. Is that correct in your opinion?

I'll keep you updated!

Jeremy

Re: Marginal Effects with Proc glimmix

Jerem — Mon, 17 Jun 2013 19:57:38 GMT

Actually, just got the answer in the meantime: after the Iteration History it is mentioned "Did not converge"...

I've to mention that I used only 691 observations out of the 1012 (for reasons specific to the subject under investigation).

When I integrate all observations, then there is the results which 'stop' after the descriptions of the dimensions of the model.

Number of Observations Read	1013
Number of Observations Used	1012

1	0	614
2	1	398

R-side Cov. Parameters	1
Columns in X	2104
Columns in Z per Subject	0
Subjects (Blocks in V)	262
Max Obs per Subject	14

Re: Marginal Effects with Proc glimmix

SteveDenham — Tue, 18 Jun 2013 15:08:23 GMT

Welcome to the fine tuning part of GLIMMIX. It really is not much fun. It involves the use of the NLOPTIONS statement.

If it looks like the objective function is closing in on some value, it may be that all you need to do is add:

nloptions maxiter=100;

to the code. GLIMMIX has a default of 20 iterations, which is very often not enough for a complex binary model to reach convergence. Try inserting this line, and see if results are obtained. If not, we will need to fine tune the convergence criteria.

Good luck.

Steve Denham

Re: Marginal Effects with Proc glimmix

Jerem — Tue, 18 Jun 2013 15:51:54 GMT

Well, at least I could handle GLIMMIX in the future thanks to your help 😉

The only thing that I do not understand is that it does not show the results. Here is the only result given by SAS.

Results of the The GLIMMIX Procedure

Number of Observations Read 1013

Number of Observations Used 1012

Response Profile

Ordered Total

Value binouverture Frequency

1 0 614

2 1 398

The GLIMMIX procedure is modeling the probability that binouverture='1'.

Dimensions

R-side Cov. Parameters 1

Columns in X 1315

Columns in Z per Subject 0

Subjects (Blocks in V) 262

Max Obs per Subject 14

proc glimmix data=candidats_ouverture;

nloptions maxiter=100;

class commune ;

model binouverture (event='1') = EurostatA EurostatC sizelog ENP commune commune*EurostatA commune*EurostatC commune*sizelog commune*ENP

/ solution dist=binary link=logit ddfm=kenwardroger ; /* Addresses individual commune levels, and adds the Kenward-Rogers adjustment for repeated measures designs */

random intercept / subject= commune solution residual; /* Converts to a repeated measures design. Even though this is a RANDOM statement, the residual option means treat the values within a commune as correlated. Note also the inclusion of commune and commune by covariate interactions in the model statement. */

lsmeans commune/ilink cl ; /* Outputs the estimates for each level of commune, on both the logit and original scale */

run;

Re: Marginal Effects with Proc glimmix

SteveDenham — Tue, 18 Jun 2013 16:42:29 GMT

If there is no Iteration History in the output, the log must have something indicating a problem. What NOTEs and WARNINGs might be showing up?

If there is an Iteration History, could you please post it?

Steve Denham

Re: Marginal Effects with Proc glimmix

Jerem — Tue, 18 Jun 2013 18:13:44 GMT

Yes, sorry Steve. There is indeed a note:

AVERTISSEMENT: Obtaining minimum variance quadratic unbiased estimates as starting values for the covariance parameters failed.

Re: Marginal Effects with Proc glimmix

SteveDenham — Wed, 19 Jun 2013 17:33:51 GMT

And so it goes...

With the full dataset, it can't get good starting values. It looks like the reduced set (691 records) did something but didn't converge. Maybe we can use the parameters from that as starting points. There should be something in the output that looks like Parameter Values at the Last Iteration. Making sure that the order is maintained, try adding a PARMS statement. Check the documentation for its use, but it will look something like:

parms (value1) (value2); /* Where value1 and value2 are from the output parameter values at the last iteration */

If this doesn't help, then simplifying the model is the next option, or going to a conditional model (G side repeated effects rather than R side). Estimates will be conditional (subject specific) instead of marginal (population averaged). I hope that the PARMS statement will help.

Steve Denham

Re: Marginal Effects with Proc glimmix

1zmm — Thu, 20 Jun 2013 12:57:30 GMT

One problem could be that the number of columns in X [=fixed effects], 1315, exceeds the number of observations, 1012. Reduce the number of fixed-effect parameters in your MODEL statement to be less than the number of observations.

Re: Marginal Effects with Proc glimmix

SteveDenham — Thu, 20 Jun 2013 17:04:59 GMT

I am afraid that is my fault. To include the separate slopes models for communes, I suggested adding the interactions. Since this is all about marginal effects, those interactions should be deleted. Try instead:

proc glimmix data=candidates;

nloptions maxiter=100;

class commune ;

model binouverture (event='1') = EurostatA EurostatC sizelog ENP commune

/ solution dist=binary link=logit ddfm=kenwardrogers ; /* Addresses individual commune levels, and adds the Kenward-Rogers adjustment for repeated measures designs. This also applies the averaged effects of the continuous covariates across all communes */

lsmeans commune/ilink ci; /* Outputs the estimates for each level of commune, on both the logit and original scale */

run;

I hope this addresses the column issue, which should in turn address the "can't get started" problem. All thanks should go to for pointing out the (what should have been) obvious source of the problem.

Steve Denham

Re: Marginal Effects with Proc glimmix

Jerem — Mon, 24 Jun 2013 16:17:18 GMT

Many thanks to both of you!

I've been extremely busy over the last days so I couldn't work that much on the model. I've thus need to have a closer look at this later.

Now the model do not converge....I tried different dependent variable and this is always the same problem: "did not converge".

Is it a problem in the modelling or the data?

Best,

Jérémy

Re: Marginal Effects with Proc glimmix

SteveDenham — Mon, 24 Jun 2013 16:23:58 GMT

"Did not converge" isn't the same as "Will not noverge". We can address the first with the NLOPTIONS statement, and by adjusting the convergence criterion in the PROC GLIMMIX statement. In order to be of assistance, though, I need to see code and output for that code, especially the iteration history

Steve Denham

Re: Marginal Effects with Proc glimmix

Jerem — Mon, 24 Jun 2013 16:42:04 GMT

Hi Steve,

I included the code as well as the SAS output. Is that all you need?

I'm really grateful for your help here!

proc glimmix data=candidates;

nloptions maxiter=100;

class commune ;

model binouverture (event='1') = EurostatA EurostatC sizelog ENP commune

/ solution dist=binary link=logit ddfm=kenwardrogers ; /* Addresses individual commune levels, and adds the Kenward-Rogers adjustment for repeated measures designs. This also applies the averaged effects of the continuous covariates across all communes */

lsmeans commune/ilink cl; /* Outputs the estimates for each level of commune, on both the logit and original scale */

run;

Jérémy

Le Système SAS

The GLIMMIX Procedure

Data Set	WORK.CANDIDATES
Response Variable	binouverture
Response Distribution	Binary
Link Function	Logit
Variance Function	Default
Variance Matrix Blocked By	Commune
Estimation Technique	Residual PL
Degrees of Freedom Method	Kenward-Roger
Fixed Effects SE Adjustment	Prasad-Rao-Kackar-Harville-Jeske

Commune	252	Aiseau-Presles Amay Amel Andenne Anderlues Anhée Ans Anthisnes Antoing Arlon Assesse Ath Attert Aubange Aubel Awans Aywaille Baelen Bassange Bastogne Beaumont Beauraing Beauvechain Beloeil Berloz Bernissart Bertogne Bertrix Beyne-Heusay Binche Bièvre Blegny Bouillon Boussu Braine-l'Alleud Braine-le-Château Braine-le-Comte Braives Brugelette Brunehaut Burdinne Burg-Reuland Büllingen Bütgenbach Celles Cerfontaine Chapelle-Lez-Herlaimont Charleroi Chastre Chaudfontaine Chaumont-Gistoux Chimay Chiny Chièvres Châtelet Ciney Colfontaine Comblain-au-Pont Comines-Warneton Courcelles Court-St-Etienne Couvin Crisnée Dalhem Dinant Dison Doische Donceel Dour Durbuy Ecaussinnes Eghezée Ellezelles Enghien Engis Erezée Erquelinnes Esneux Estaimpuis Estinnes Etalle Eupen Faimes Farciennes Fauvillers Fernelmont Ferrières Fexhe-le-Haut-Clocher Flemalle Fleurus Flobecq Floreffe Florennes Florenville Fléron Fontaine-l'Evêque Fosses-la-Ville Frameries Frasnes-lez-Anvaing Froidchapelle Gedinne Geer Gembloux Genappe Gerpinnes Gesves Gouvy Grace-Hollogne Grez-Doiceau Habay Ham-sur-Heure-Nalinnes Hamois Hannut Hastière Havelange Hensies Herbeumont Herstal Herve Honnelles Hotton Houffalize Houyet Huy Hélecine Héron Incourt Ittre Jalhay Jemeppe-sur-Sambre Jodoigne Juprelle Jurbise Kelmis La Bruyère La Hulpe La Louvière La Roche-en-Ardenne Lasne Le Roeulx Lens Les-Bons-Villers Lessines Leuze-en-Hainaut Libramont-Chevigny Lierneux Limbourg Lincent Liège Lobbes Lontzen Malmedy Manage Manhay Marche-en-Famenne Marchin Martelange Meix-devant-Virton Merbes-le-Château Messancy Mettet Modave Momignies Mons Mont-St-Guibert Mont-de-l'Enclus Montigny-le-Tilleul Morlanwelz Mouscron Musson Namur Nandrin Nassogne Neufchâteau Neupré Nivelles Ohey Olne Onhaye Oreye Orp-Jauche Ottignies-LLN Ouffet Oupeye Paliseul Pecq Pepinster Perwez Philippeville Plombières Pont-à-Celles Profondeville Péruwelz Quaregnon Quiévrain Quévy Raeren Ramilies Rebecq Remicourt Rixensart Rochefort Rouvroy Rumes Saint Nicolas Saint-Georges-sur-Meuse Saint-Ghislain Saint-Hubert Sainte-Ode Sambreville Sankt Vith Seneffe Seraing Silly Sivry-Rance Soignies Sombreffe Somme-Leuze Soumagne Spa Sprimont Stavelot Stoumont Tellin Theux Thimister-Clermont Thuin Tinlot Tintigny Tournai Trooz Tubize Vaux-sur-Sûre Verlaine Verviers Vielsalm Villers-la-Ville Villers-le-Bouillet Viroinval Virton Visé Vresse-sur-Semois Walcourt Walhain Wanze Waremme Wasseiges Waterloo Wavre Welkenraedt Wellin Yvoir

Number of Observations Read	804
Number of Observations Used	803

1	0	476
2	1	327

R-side Cov. Parameters	1
Columns in X	257
Columns in Z per Subject	0
Subjects (Blocks in V)	252
Max Obs per Subject	14

Optimization Technique	None
Parameters	0
Lower Boundaries	0
Upper Boundaries	0
Fixed Effects	Profiled
Residual Variance	Profiled
Starting From	Data

0	1.340781E154	2.00000000	.
1	1.340781E154	1.13896683	.
2	1.340781E154	2.00000000	.
3	1.340781E154	2.00000000	.
4	1.340781E154	0.91900698	.
5	3460.7243613	2.00000000	.
6	1.340781E154	2.00000000	.
7	3728.5751363	2.00000000	.
8	1.340781E154	2.00000000	.
9	3993.8654206	2.00000000	.
10	4127.2354391	2.00000000	.
11	1.340781E154	2.00000000	.
12	4394.1129597	2.00000000	.
13	4525.2340423	2.00000000	.
14	4658.3500927	2.00000000	.
15	4791.7156	2.00000000	.
16	4924.1913727	2.00000000	.
17	5059.8010277	2.00000000	.
18	5191.3961115	2.00000000	.
19	5324.1823992	2.00000000	.
20	1.340781E154	2.00000000	.
21	1.340781E154	2.00000000	.
22	1.340781E154	2.00000000	.
23	1.340781E154	2.00000000	.
24	1.340781E154	2.00000000	.
25	1.340781E154	0.99810273	.
26	1.340781E154	0.81393176	.
27	1.340781E154	2.00000000	.
28	3148.7088811	2.00000000	.
29	3281.2624666	2.00000000	.
30	1.340781E154	2.00000000	.
31	3547.6660221	2.00000000	.
32	3677.9846413	2.00000000	.
33	3814.0718304	2.00000000	.
34	1.340781E154	2.00000000	.
35	1.340781E154	0.14868120	.
36	4210.5072254	2.00000000	.
37	4346.0003696	2.00000000	.
38	1.340781E154	2.00000000	.
39	1.340781E154	0.09425822	.
40	1.340781E154	0.09072392	.
41	4877.3197427	2.00000000	.
42	5009.6182614	2.00000000	.
43	1.340781E154	2.00000000	.
44	1.340781E154	0.06120650	.
45	5409.7912614	2.00000000	.
46	1.340781E154	2.00000000	.
47	1.340781E154	2.00000000	.
48	16935.289247	2.00000000	.
49	1.340781E154	2.00000000	.
50	1.340781E154	2.00000000	.
51	1.340781E154	0.99810273	.
52	1.340781E154	0.81393176	.
53	1.340781E154	2.00000000	.
54	3148.7088811	2.00000000	.
55	3281.2624666	2.00000000	.
56	1.340781E154	2.00000000	.
57	3547.6660221	2.00000000	.
58	3677.9846413	2.00000000	.
59	3814.0718304	2.00000000	.
60	1.340781E154	2.00000000	.
61	1.340781E154	0.14868120	.
62	4210.5072254	2.00000000	.
63	4346.0003696	2.00000000	.
64	1.340781E154	2.00000000	.
65	1.340781E154	0.09425822	.
66	1.340781E154	0.09072392	.
67	4877.3197427	2.00000000	.
68	5009.6182614	2.00000000	.
69	1.340781E154	2.00000000	.
70	1.340781E154	0.06120650	.
71	5409.7912614	2.00000000	.
72	1.340781E154	2.00000000	.
73	1.340781E154	2.00000000	.
74	16935.289247	2.00000000	.
75	1.340781E154	2.00000000	.
76	1.340781E154	2.00000000	.
77	1.340781E154	0.99810273	.
78	1.340781E154	0.81393176	.
79	1.340781E154	2.00000000	.
80	3148.7088811	2.00000000	.
81	3281.2624666	2.00000000	.
82	1.340781E154	2.00000000	.
83	3547.6660221	2.00000000	.
84	3677.9846413	2.00000000	.
85	3814.0718304	2.00000000	.
86	1.340781E154	2.00000000	.
87	1.340781E154	0.14868120	.
88	4210.5072254	2.00000000	.
89	4346.0003696	2.00000000	.
90	1.340781E154	2.00000000	.
91	1.340781E154	0.09425822	.
92	1.340781E154	0.09072392	.
93	4877.3197427	2.00000000	.
94	5009.6182614	2.00000000	.
95	1.340781E154	2.00000000	.
96	1.340781E154	0.06120650	.
97	5409.7912614	2.00000000	.
98	1.340781E154	2.00000000	.
99	1.340781E154	2.00000000	.

Did not converge.

Residual (VC)	1.0328	.

Re: Marginal Effects with Proc glimmix

SteveDenham — Mon, 24 Jun 2013 17:02:35 GMT

I may have spoken too hastily. This has all the symptoms of quasi-separation--for at least some of the communes, all of the values are either 1 or zero, and all of the continuous covariates are identical. I think you will either have to relinquish the hope of finding marginal estimates (i.e., consider commune to be a random effect), or find a way to collapse the communes into some larger categorization that eliminates the quasi-separation. There is a major problem in trying to estimate 258 parameters from only 803 observations. A rule of thumb for getting decent estimates is at least 5, and preferably 10, observations per parameter. I think you may be able to build a respectable model if you could collapse the commune variable into 70 or 80 "super-commune" levels.

Steve Denham