topic Re: Proc Glimmix in Statistical Procedures

Proc Glimmix

Ashwini_uci — Wed, 15 Jan 2014 19:28:56 GMT

Hi,

I am working on hierarchical models using PROC GLIMMIX. This is my first time working on this proc.

I tried following code but the log returned an error that says "ERROR: Invalid or missing data."

I am not sure if my code is incorrect ot there is anything missing in the dataset.

Also, what do Insert in place of n in th e=model statement.

proc glimmix data= xxxxx ;

class hospid female dm dmcx htn_c aids alcohol ANEMDEF arth race1 ZIPINC_QRTL hosp_location h_contrl hosp_teach

bldloss chf chrnlung coag depress drug hypothy liver lymph lytes mets neuro obese para perivasc psych pulmcirc renlfail tumor

ulcer valve wghtloss cararrhythmia;

model died/n = female dm dmcx htn_c aids alcohol ANEMDEF arth race1 ZIPINC_QRTL hosp_location h_contrl hosp_teach

bldloss chf chrnlung coag depress drug hypothy liver lymph lytes mets neuro obese para perivasc psych pulmcirc renlfail tumor

ulcer valve wghtloss cararrhythmia / solution;

random intercept / subject=hospid;

run;

Any help will be greatly appreciate!

Thanks

Ashwini

Re: Proc Glimmix

SteveDenham — Wed, 15 Jan 2014 20:16:35 GMT

Hi Ashwini,

That error indicates that something is missing in the dataset, and if you have the highlighted "n" in the model statement, then it is very likely that is what is missing.

So, if we can answer that question, then it will probably solve some of your problems. The events/trials syntax that you are using implies that you know how many events (in this case, died, assuming that is a summed value across the many class variables) and how many trials are available. For this it could be died plus lived, or if it is an intent to treat analysis, it could be number enrolled. Just be sure that the number reflects the cell "total".

With this many categorical variables, I would watch out for quasi-separation. You may need to collapse categories or eliminate some variables if the crosstabs reveal zeroes.

Steve Denham

Re: Proc Glimmix

Ashwini_uci — Thu, 16 Jan 2014 01:41:36 GMT

Hi Steve,

Thanks for your response.

I read a little more and ran another code which is as follows:

proc glimmix data=library.nismicathcabg1_3;

class hospid;

model died (event=last) = female dm dmcx htn_c aids alcohol ANEMDEF arth race1 ZIPINC_QRTL hosp_location h_contrl hosp_teach bldloss chf chrnlung coag depress drug hypothy liver lymph lytes mets neuro obese para perivasc psych pulmcirc renlfail tumor ulcer valve wghtloss cararrhythmia / dist=binary link=logit ddfm=bw solution ;

random intercept / subject=hospid solution;

run;

It returned me incomplete output, where at the end it said "Did not converge".

I am attaching the output for your reference to the original post.

I am not sure why the output is incomplete and showing "did not converge".

Appreciate your advice !

Thanks

Ashwini

Re: Proc Glimmix

SteveDenham — Thu, 16 Jan 2014 17:46:32 GMT

The default number of iterations in GLIMMIX is twenty. I see that the objective function looks to be leveling out, but that there is still some iterating to do.

Add the following line to your GLIMMIX code:

nloptions maxiter=200;

If it still does not converge, try adding pconv=1e-6 to the PROC GLIMMIX statement.

Steve Denham

Re: Proc Glimmix

Ashwini_uci — Fri, 17 Jan 2014 23:21:48 GMT

Hi Steve,

Thanks again for your advice. i tried both but to no avail. I still get the same output that has "Did not converge" in it.

I also tried using your suggestions in a very parsimonious model but that didn't help either.

Not sure what is wrong with the code...

Re: Proc Glimmix

SteveDenham — Mon, 20 Jan 2014 19:26:31 GMT

I really don't think anything is wrong with your code. Could you please share the output, especially the iteration history part? I can be more helpful if I see how that is behaving.

Steve Denham

Re: Proc Glimmix

Ashwini_uci — Mon, 20 Jan 2014 21:35:04 GMT

These are the codes

proc glimmix data=library.nismicathcabg1_3;

class hospid;

model died (event=last) = female dm dmcx htn_c aids alcohol ANEMDEF arth race1 ZIPINC_QRTL hosp_location h_contrl hosp_teach

bldloss chf chrnlung coag depress drug hypothy liver lymph lytes mets neuro obese para perivasc psych pulmcirc renlfail tumor

ulcer valve wghtloss cararrhythmia / dist=binary link=logit ddfm=bw solution ;

random intercept / subject=hospid solution;

nloptions maxiter=200;

run;

proc glimmix data=library.nismicathcabg1_3 pconv=1e-6;

class hospid;

model died (event=last) = female dm dmcx htn_c aids alcohol ANEMDEF arth race1 ZIPINC_QRTL hosp_location h_contrl hosp_teach

bldloss chf chrnlung coag depress drug hypothy liver lymph lytes mets neuro obese para perivasc psych pulmcirc renlfail tumor

ulcer valve wghtloss cararrhythmia / dist=binary link=logit ddfm=bw solution ;

random intercept / subject=hospid solution;

run;

Re: Proc Glimmix

lvm — Wed, 22 Jan 2014 14:57:48 GMT

You have a very large number of predictor variables in your model. I am guessing that the log likelihood function for your data and model is fairly flat in the region of the maximum value, so that the optimization method is bouncing around and not converging at the actual maximum. I suggest you try a much simpler fixed effect part of your model. Start with just a few predictor variables so that you can obtain convergence. Then add more predictor variables to see if you can find the circumstances where problems develop.

Other comments: GLMMs with a binary conditional distribution (your situation) can give biased results. You should try METHOD=LAPLACE in the procedure statement. You should also try fitting the model without a random effect, just to look at the parameter values (to give you a sense of the magnitude and sign of the coefficient estimates). This would be just for exploratory work (you should not ignore the random effect for the final analysis). You should be able to see that many of the predictors will likely not be significant (due to an overfit).

Re: Proc Glimmix

SteveDenham — Fri, 24 Jan 2014 20:05:17 GMT

In addition, presenting the iteration history helps me in interpreting/troubleshooting. Note that from iteration 5 onward, the objective function moves by a nearly fixed quantity (20 and a tiny fraction) and from iteration 13 on the tiny fraction goes away. This is a sign of an over parameterized model, with a collinearity problem. Some of these predictors are nearly completely confounded, such that the effective X matrix has only 17 linearly independent columns (I think that's right. Hopefully, @lvm will correct me if I messed this up), and the 20 redundant columns just keep adding on with each iteration.

So, use subject matter knowledge to boil the independent variables down to a reasonable list. If this is an exploratory study, you could try PROC GLMSELECT and use the LASSO option to get a handle on the independent variables. At least run PROC REG and check the collinearity diagnostics.

Steve Denham

Re: Proc Glimmix

Ashwini_uci — Sat, 25 Jan 2014 00:45:34 GMT

Thanks @Ivm and Steve for your useful suggestions! I tried to minimize the number of predictor variables and I ran the following model. And I could get some output which seems right. But it gives me only estimates and p-values. I am attaching the output for your reference. Is there any way I can also get the odds ratio by using proc glimmix?

proc glimmix data=library.data;

class hospid female dm htn_c ANEMDEF arth race1 ZIPINC_QRTL bldloss chf chrnlung coag depress drug hypothy

liver lymph lytes mets neuro obese para perivasc pulmcirc renlfail wghtloss cararrhythmia;

model died (event='1') = female dm htn_c ANEMDEF arth race1 ZIPINC_QRTL bldloss chf chrnlung coag depress drug hypothy liver lymph lytes mets neuro obese para perivasc pulmcirc renlfail wghtloss cararrhythmia /dist=binary link=logit ddfm=bw solution;

random intercept / subject=hospid;

run;

Re: Proc Glimmix

Ashwini_uci — Sat, 25 Jan 2014 00:46:36 GMT

proc glimmix data=library.data;

class hospid female dm htn_c ANEMDEF arth race1 ZIPINC_QRTL bldloss chf chrnlung coag depress drug hypothy

liver lymph lytes mets neuro obese para perivasc pulmcirc renlfail wghtloss cararrhythmia;

random intercept / subject=hospid;

run;

Re: Proc Glimmix

lvm — Sun, 26 Jan 2014 01:12:15 GMT

Just add the ODDSRATIO option on the model statement. See the documentation for suboptions, but this should give you what you want.

Re: Proc Glimmix

Ashwini_uci — Tue, 18 Mar 2014 20:46:20 GMT

Hello lvm,

As discussed in this post, I finally got the right code for Proc Glimmix and I have been using it thus far. Only today when I started using it on much larger dataset (over 250000 records) . It is the same code as I have posted above. The Proc Glimmix worked fine on the 100000 cases. But on the new dataset witha around 250000 cases, the output I am getting is not complete. Here is the code.I am not sure why it is not yielding a complete. Any help is greatly appreciated. I am also attaching the output for your reference.

proc glimmix data=library.nismi091011_gender5 ; *using oddsratio*/;

class hospid gender racenew ZIPINC_QRTL_1 dm_all cm_htn_c cm_ANEMDEF cm_arth cm_bldloss cm_chf

cm_chrnlung cm_coag cm_depress cm_drug cm_hypothy cm_liver cm_lymph cm_lytes cm_mets cm_neuro cm_obese cm_para

cm_perivasc cm_pulmcirc cm_renlfail cm_wghtloss cararrhythmia;

model died (event='1') = age gender racenew ZIPINC_QRTL_1 dm_all cm_htn_c cm_ANEMDEF cm_arth

cm_bldloss cm_chf cm_chrnlung cm_coag cm_depress cm_drug cm_hypothy cm_liver cm_lymph cm_lytes cm_mets cm_neuro cm_obese cm_para cm_perivasc cm_pulmcirc cm_renlfail cm_wghtloss cararrhythmia /oddsratio dist=binary link=logit ddfm=bw solution;

random intercept / subject=hospid;

where pcionly=1 and stemi=0;

weight discwt;

title 'Logi Reg in-hosp mortality vs gender in POST-PCI MI patients using "where" option with RACE +income for nonstemi with pcionly';

run;

quit;

Re: Proc Glimmix

Ashwini_uci — Tue, 18 Mar 2014 20:47:38 GMT

Hello Steve,

proc glimmix data=library.nismi091011_gender5 ; *using oddsratio*/;

class hospid gender racenew ZIPINC_QRTL_1 dm_all cm_htn_c cm_ANEMDEF cm_arth cm_bldloss cm_chf

cm_chrnlung cm_coag cm_depress cm_drug cm_hypothy cm_liver cm_lymph cm_lytes cm_mets cm_neuro cm_obese cm_para

cm_perivasc cm_pulmcirc cm_renlfail cm_wghtloss cararrhythmia;

model died (event='1') = age gender racenew ZIPINC_QRTL_1 dm_all cm_htn_c cm_ANEMDEF cm_arth

random intercept / subject=hospid;

where pcionly=1 and stemi=0;

weight discwt;

title 'Logi Reg in-hosp mortality vs gender in POST-PCI MI patients using "where" option with RACE +income for nonstemi with pcionly';

run;

quit;

Re: Proc Glimmix

lvm — Tue, 18 Mar 2014 21:00:30 GMT

You may simply be running out of memory. What does the Log say?

Re: Proc Glimmix

lvm — Tue, 18 Mar 2014 21:23:51 GMT

By the way, the output says you have 87,303 observations. You may have a problem with your datafile. It comes back tot he Log: what does it say?

Re: Proc Glimmix

Ashwini_uci — Tue, 18 Mar 2014 21:38:40 GMT

The log looks fine. Doesn't say anything unusual.

Re: Proc Glimmix

Ashwini_uci — Tue, 18 Mar 2014 21:43:14 GMT

The number of observations is right. I am working on a subset (87,303) of the whole sample (250,000) and after weighting this number of sub-sample goes up to almost 400,000.

Re: Proc Glimmix

Ashwini_uci — Tue, 18 Mar 2014 22:21:33 GMT

Strange enough, none of the models that worked successfully on the earlier small dataset is working on this new larger dataset.

Re: Proc Glimmix

Ashwini_uci — Wed, 19 Mar 2014 01:31:50 GMT

The log is showing me the warning as below:

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

Never encountered it before and not sure what it means. Do you have any idea what it exactly means and how to fix it?

Thanks for your time and help!

Ashwini