Proc Glimmix - confusion about specifying R-side vs. G-side parameters and output interpretation

[nstdt]

I have repeated (over time) data on 9 subjects, subjected to two types of 'Phases' - Treatment and Control, similar to below table (for just 2 subjects - second columns is outcome Y - in count format). There are a total of 4 phases, alternating in an ABAB design.

A	10
	6
	5
B	20
	15
	7

I converted the outcomes to proportions and am trying to fit a mixed effects logit model. The first model only has 2 random effects - Intercept and Phase :

title "Basic logit model";
proc glimmix data=lambert2;
class case;
model y = phase/solution dist=binomial 	link=logit cl;
random Intercept phase/sub=case;
run;
OUTPUT:
Dimensions 
G-side Cov. Parameters 2 
Columns in X 2 
Columns in Z per Subject 2 
Subjects (Blocks in V) 9 
Max Obs per Subject 59 


Estimated G matrix is not positive definite. 


Fit Statistics 
-2 Res Log Pseudo-Likelihood 1847.10 
Generalized Chi-Square 241.97 
Gener. Chi-Square / DF 0.53 



Covariance Parameter Estimates 
Cov Parm Subject Estimate Standard
Error 
Intercept case 0 . 
phase case 0 . 

Q1: In the Output Dimensions Table, which are the columns in X and which variables go in Z? Are both Intercept and Phase random errors considered G-side parameters? Both covariance parameters are 0 - what does this mean?

Next, I dummy coded the Phase variable to reflect the ABAB design (with three dummies ) and re-ran the Mixed Logit model with what I thought were R-side errors using compound symmetry (I mean that I am modeling the within subject correlation using compound symmetry) as below:

title "Model with only R-side CS errors";
proc glimmix data=model_dat2;
class case phase;
model y = flag_2 flag_3 flag_4/solution dist=binomial 	link=logit cl;
random  Intercept/sub=case type=cs residual;
run;

OUTPUT:
Dimensions 
R-side Cov. Parameters 2 
Columns in X 4 
Columns in Z per Subject 0 
Subjects (Blocks in V) 9 
Max Obs per Subject 59 

Fit Statistics 
-2 Res Log Pseudo-Likelihood 1718.07 
Generalized Chi-Square 141.32 
Gener. Chi-Square / DF 0.31 



Covariance Parameter Estimates 
Cov Parm Subject Estimate Standard
Error 
CS case 0.006232 0.006211 
Residual   0.3099 0.02071 

For the model with the ABAB phase dummies as above, why is number of columns in Z=0? I have specified "Intercept" as a random effect (at least that's my intention)? Is my syntax incorrect? Is the variance of the Intercept not considered a G-side random effect?

Then I tried running the syntax below but this does not converge:

title "Model with only R-side CS errors -with Intercept";
proc glimmix data=model_dat2;
class case phase;
model y = flag_2 flag_3 flag_4/solution dist=binomial 	link=logit cl;
/*random Intercept flag_2 flag_3 flag_4/sub=case;*/
random  phase/sub=case type=cs residual;
random intercept/sub =case;
run;

OUTPUT:
Dimensions 
G-side Cov. Parameters 1 
R-side Cov. Parameters 2 
Columns in X 4 
Columns in Z per Subject 1 
Subjects (Blocks in V) 9 
Max Obs per Subject 59 

Did not converge. 



Covariance Parameter Estimates 
Cov Parm Subject Estimate Standard
Error 
Intercept case 0.2221 . 
CS case -0.03368 . 
Residual   0.3067 . 

How do I interpret the above results?

[SteveDenham]

I'll hit the last one as it is one I have seen many times when fitting an R side effect to a binomial dataset.  The most likely cause of failure to converge is that you have run up against the default maximum iteration limit of 20 iterations. Try adding an NLOPTIONS statement with maxiter=5000 for example. If you still are not converging, your data is probably "pathological" with respect to the residuals.

Backing up one step, if you fit only an R side structure, the lack of Z columns is dependent solely on the syntax chosen, and the number of effects in your random / residual statement.

SteveDenham

[nstdt]

Thanks so much for the reply! I will try increasing the interations. Meanwhile, I want to know how I can get G-side effects for the Intercept in the model shown below(I have specified a random Intercept - but number of columns in Z=0 per the SAS output - what is wrong with my syntax???):

title "Model with only R-side CS errors";
proc glimmix data=model_dat2;
class case phase;
model y = flag_2 flag_3 flag_4/solution dist=binomial 	link=logit cl;
random  Intercept/sub=case type=cs residual;
run;

OUTPUT:
Dimensions 
R-side Cov. Parameters 2 
Columns in X 4 
Columns in Z per Subject 0 
Subjects (Blocks in V) 9 
Max Obs per Subject 59 

[jiltao]

I do not understand your data -- is it a crossover design with repeated measures? Or just a repeated measures data? If you can share your data that would be helpful.

The statement --

random Intercept phase/sub=case;

specifies 2 columns in the Z matrix for each subject -- a column of 1 (intercept) and a column for Phase.

The Estimated G matrix is not positive definite message indicates that the model might not work for your data.  

The statement --

random  Intercept/sub=case type=cs residual;

models the correlations in the (linearized) residuals. It is similar to the REPEATED statement in PROC MIXED. So there is no G-side random effect, therefore no Z matrix. The syntax above specifies INTERCEPT as the repeated effect, which might not be appropriate.

The statements -

random  phase/sub=case type=cs residual;
random intercept/sub =case;

specify a random intercept model as well as a CS correlated "residuals". This is likely an overspecified model and therefore encountered some convergence issue.

Again, if you can share your data (even with the fake y values) I might be able to provide more assistance.

Thanks,

Jill

[nstdt]

I have an ABAB design with repeated time observations within a single phase (A - baseline phase, B- treatment phase. Each subject gets both phases, repeated many times, about 10 to 30 per subject per phase). There are 9 subjects in total. What is the best wyay to model this data using a mixed model structure - I was trying to see if each subject is a 'unit', and it makes sense to have random effects for the Intercept and maybe the Phase for each subject. But there is also the correlation amongst 'outcome' scores for a specific subject within a Phase - surely that has to be something like AR(1) or CS?

Here is data from one single subject 'case A1':

case	phase	time	outcome
case A1	1	1	7
case A1	1	2	9
case A1	1	3	8
case A1	1	4	6
case A1	1	5	7
case A1	1	6	4
case A1	1	7	5
case A1	1	8	10
case A1	1	9	2
case A1	1	10	0
case A1	2	15	3
case A1	2	16	8
case A1	2	17	8
case A1	2	18	6
case A1	2	19	10
case A1	2	20	10
case A1	2	21	10
case A1	2	22	8
case A1	2	23	3
case A1	2	24	4

[jiltao · Accepted solution]

The following paper has an example that is similar to your situation --

Microsoft Word - SAS1919-2015.docx

See Pages 4-6.

You might need to recode the values for time -- right now the values seem to be continuous throughout the data set, which might need to be changed. Also, Phase might need to be in the CLASS statement.

Hope this helps,

Jill

[nstdt]

Thanks so much! The values per se of the 'time' variable aren;t important - I'm just interested in the correlation over time within a Phase for a specific subject.