PROC NLMIXED with categorical covariates in nominal scale

[issac]

Hi everybody

I am working with 9.3 version and want to apply NLMIXED with class effects measured in nominal scale. The problem is that unlike other procedures like PROC MIXED, GLIMMIX, PHREG, etc. there is no CLASS statement and I am a little bit confused how to define these variables. following in my submitted code:

proc nlmixed data = test ;

lambda=exp(b0+bdisage*disage+blos*los+bmar*mar+bemp*emp+bagents*agents+badmsource*admsource+bseq*seq+e);

ll=-lambda*rtime**(alpha+1)+rstatus*(LOG(alpha+1)+alpha*LOG(rtime)+LOG(lambda));

MODEL rtime~GENERAL(ll);

RANDOM e~NORMAL(0,s2) SUBJECT=id;

PARMS b0=1 bdisage=0 blos=0 bmar=0 bemp=0 bagents=0 badmsource=0 bseq=0 s2=1 alpha=0;

run;

Apart from "disage" and "los" which are continuous, all other effects (from "mar" to "admsource") are nominal scaled effects, and "seq" is ordinal scaled variable.

Part of the results from log window is as below:

NOTE: Character value converted to numeric for argument 2 of '*' operation.

NOTE: Character value converted to numeric for argument 2 of '*' operation.

NOTE: Character value converted to numeric for argument 2 of '*' operation.

NOTE: Character value converted to numeric for argument 2 of '*' operation.

408  ll=-lambda*rtime**(alpha+1)+rstatus*(LOG(alpha+1)+alpha*LOG(rtime)+LOG(lambda));

409  MODEL rtime~GENERAL(ll);

410  RANDOM e~NORMAL(0,s2) SUBJECT=id;

411  PARMS b0=1 bdisage=0 blos=0 bmar=0 bemp=0 bagents=0 badmsource=0 bseq=0 s2=1 alpha=0;

412  run;

NOTE: Execution error for observation 1.

NOTE: PROCEDURE NLMIXED used (Total process time):

      real time           0.24 seconds

      cpu time            0.09 seconds

Any helpful comments would be totally appreciative!

Thanks!

Issac

[SteveDenham]

You will probably have to use PROC GLMMOD to generate the design matrix.  That would be the singular form.  You might also choose to use PROC LOGISTIC to get a non-singular form using the OUTDESIGN= and OUTDESIGNONLY oprtions on the PROC LOGISTIC statement.  This will give a design matrix that you can use to construct lambda=.

Good luck.

Steve Denham

[issac]

Steve;

I did this in proc logistic but when I push the new SAS-Data-Set into the NLMIXED, I got an error. I think this is caused since the new data set have new variables labeled differently, say marMARRIED, marNEVER_MARRIED, but I think I need the mar as an categorical covariate in my model, not decoded levels of them. I want to get an estimation for mar as a single explanatory variable, not for its levels, since it is meaningless. Can we do that with the help of design statements?

Thanks so much!

Issac   

[SteveDenham]

Issac,

Mapping the old variables to the new ones is one task, and it will have to be addressed.  It is much easier to address this in the MODEL statement, than to use a data step and rename variables.

As far as an effect for mar as a single test, in NLMIXED the only way I know to address it is through the use of appropriate contrast statements.  So for a variable with four levels, coded under the GLMMOD singular method as var_1 to var_4, you might have:

CONTRAST 'Four level variable' var_1 - var_2=0, var_2 - var_3=0, var_3 - var_4=0;

Does this help?

Steve Denham

[issac]

Steve;

I am afraid it does not work. Since I have to create a likelihood function with the help of lambda and ll, and then assign it to the response variable 'rtime' (MODEL rtime~GENERAL(ll);). And in the ll and lambda definition, I need to stick with 'mar', 'emp', 'agents', 'admsource' as whole covariates, not level-specific, such as below:

lambda=exp(b0+bdisage*disage+blos*los+bmar*mar+bemp*emp+bagents*agents+badmsource*admsource+bseq*seq+e);

[SteveDenham · Accepted solution]

Issac,

I guess I am missing something here.  Suppose disage has four levels.  If you use GLMMOD to generate the design matrix, it will create five columns that are indicator variables.  So if a subject belongs to the first age class, it would have a vector that looks something like (1 1 0 0 0), the second class (1 0 1 0 0), the third (1 0 0 1 0) and the fourth (1 0 0 0 1).  Each column has a name, named cleverly col1 through col5.  The first column (col1) is the intercept, and will not be used in the NLMIXED code.

Your lambda statement will then look like (for ONLY disage):

lambda = exp(b0+b12*col2+b13*col3+b14*col4+b15*col5);

I use the nomenclature for the betas as bij, where i denotes the group variable, and j denotes the column in the design matrix.  This is only to keep things straight.  By the time you get through coding all of the betas, this may seem tedious, but it is going to help in the contrast statements that give the equivalent of testing of groups.

The contrast statement to give a 3 degree of freedom test of whether any of the betas are different from one another within age is then:

CONTRAST 'Disage' b12-b13=0, b13-b14=0, b14-b15=0;

You should be able to generalize from this.  As far as I know, this is the only way to deal with class variables in any current version of NLMIXED.

Steve Denham

[issac]

Steve

Really appreciate for your helpful comments. Given that we complete theses steps, I think, we may reach an analysis of likelihood parameter estimates table which has 4 lines for disage categories (with one line as an omitted, which can be though as a reference contrast). But what about overall test for disage? I mean Type III Analysis of Effects Table? As far as I know, this test also does not depend on which category of disage is the omitted category. Can we have something as well in the procedure, or we need to get it done manually?

Thanks!

Issac

[issac]

Steve;

I think I found something helpful. For the overall test, we should fit the model with and without the set of indicators.  A likelihood ratio chi-square statistic is then found by taking twice the positive difference in the log-likelihoods for the two models and the DF is the difference in the number of parameters estimated.

Anyway, I face a silly problem in running the code. I submit the following (even without any categorical variables)

proc nlmixed data=test;

lambda=exp(b0+bdisage*disage+blos*los+bseq*seq+e);

ll=-lambda*rtime**(alpha+1)+rstatus*(LOG(alpha+1)+alpha*LOG(rtime)+LOG(lambda));

MODEL rtime~GENERAL(ll);

RANDOM e~NORMAL(0,s2) SUBJECT=id;

PARMS b0=1 bdisage=0 blos=0 bseq=0 s2=1 alpha=0;

run;

and get a prompt result with specification table as below

Data Set	TMP1.READMISSION
Dependent Variable	rtime
Distribution for Dependent Variable	General
Random Effects	e
Distribution for Random Effects	Normal
Subject Variable	id
Optimization Technique	Dual Quasi-Newton
Integration Method	Adaptive Gaussian Quadrature

and this print in the log :

...

NOTE: Execution error for observation 1.

NOTE: PROCEDURE NLMIXED used (Total process time):

      real time           0.04 seconds

      cpu time            0.04 seconds

What is this about?

[SteveDenham]

Issac,

Two points:  Regarding the likelihood ratio test, this is an excellent idea.  It does require multiple fits of the nested models, where using the CONTRAST statement properly should require only a single run.  I have to admit that writing good orthogonal contrast statements is not always easy.

As to the execution error for observation 1, my best guess is that for the data in hand and the starting values proposed for the parameters, something is uncalculable.  My first guess would be log(rtime).  If rtime is zero in observation 1, then this statement will return a missing value, and execution will stop.

Steve Denham

[issac]

Steve;

Totally correct. We can also get some other options with CONTRAST in NLMIXED. It is a good point. Besides, the problem is exactly what you guessed, then I change it to a very little number (near zero) and now it works.

Thanks a lot!

Issac