PROC MIXED AR(1) with RANDOM and REPEATED fails to ESTIMATE Random Effect coefficients

[TBeans]

Hi All,

I'm running multiple versions of a model through the %MACRO below.  I want to specify an AR(1) (or ARMA 1,1) for the R matrix and produce estimates of the in-group mixed effect estimates (in other words, the group-specific coefficients).

The program runs fine when the REPEATED statement is not used, but the marginal Var-Covar matrix does not appear to have the correct form.  I've seen many examples of using RANDOM and REPEATED together to produce an AR(1) Covariance structure, but none with the ESTIMATE function being run -- Am I doing something wrong, misunderstanding how RANDOM and/or REPEATED interact??

     %macro ME_Estimate(InDat, FEtab, REtable, groupvar, FE, RE);
          %local i;
          %let nRE = %sysfunc(countw(&RE.));

          DATA temp; set &InDat.; RUN;
          PROC SQL noprint;
          SELECT DISTINCT &GROUPVAR. INTO :GRP_LIST separated by " " FROM temp;
          SELECT COUNT(DISTINCT CITY) INTO :N_grp FROM temp;
          QUIT;

         
          PROC MIXED data= &InDat.;
          class &groupvar. qtr;
          MODEL response =   &FE.    / DDFM=KR Solution;
          RANDOM &RE. / TYPE=ARMA(1,1) subject=&groupvar. Solution;

          REPEATED qtr  / TYPE=ARMA(1,1) subject=&groupvar.; * HERE IS THE PROBLEM ;

          %DO i = 1 %TO &N_grp;
                    %LET CURR_GRP =%SCAN(&GRP_LIST,&I.);
                    %DO tests = 1 %TO &nRE.;
                             %LET CURR_RE = %SCAN(&RE., &tests.);
                              ESTIMATE "&CURR_RE._&CURR_GRP." &CURR_RE. 1 | &CURR_RE. 1 / subject %est_Zeroes(&i.);
                    %END;
          %END;
          ods output SolutionF = &FEtab.  SolutionR = REdeviate  estimates = &REtable.;
          RUN;
          PROC DATASETS noprint; DELETE TEMP; RUN;
      %mend;

I then call the macro with something like:

       %let FEvars1 = Var1 Var2 Var3;

        %let REvars1_1 = Var1;

       %ME_Estimate(InDat=V&V._ALL(where=(qtr >= '1JAN1987'd)),
                          FEtab = FEout1_1, REtable= REout1_1, groupvar=GRP1, 

                           FE = &FEvars1., RE = &REvars1_1.);

And the Est_Zeroes macro is simply:

     %macro est_Zeroes(numzeroes);
               %local i;
               %do i = 1 %to %eval(&numzeroes-1);
                    0
               %end;
                    1;
     %mend;

As I said, it works fine without the REPEATED statement.  When I include REPEATED the Std Errors and/or the estimates from the Random Effects Solutions are zero!!  I'm not sure what's going on.  

[SteveDenham]

One thing that helps me is to step back from the macro version and see if I can get PROC MIXED to behave with regular input.

PROC MIXED data= InDat.;
          class groupvar qtr fixed_effects;
          MODEL response =  fixed_effects qtr fixed_effect*qtr   / DDFM=KR Solution;
          RANDOM intercept /  subject=groupvar Solution;
 
          REPEATED qtr  / TYPE=ARMA(1,1) subject=groupvar;
run;

Here I used groupvar to mean the subject on which the repeated measures are observed. I introduce fixed_effects to cover all fixed effects in the model, qtr is the time indicator, and I make sure I include the time by fixed effect interaction.

In the RANDOM statement, I do not specify any TYPE=, as I only want to model variability across subjects.

The REPEATED statement covers the ARMA structure.

Try fitting this to see if it gives sensible results.

I see the following as possible causes of your problem--failure to include the repeated effect variable in the model statement, or modeling both the RANDOM and REPEATED effect identically.  The latter will lead to zeroes (as you have found), the former fails to apply the repeated nature of observations to the subjects.

Steve Denham

[TBeans]

Hi Steve!

I certainly agree with the general principle of wokring in regular, non-macro version.  I've beening trying your suggestions but I think I'm running into the limitations of my dataset and available processing power.

1) When I include my fixed effects variables in the CLASS statement then I get the Error: The SAS System stopped processing this step because of insufficient memory.  I can avoid this by reducing the model to a single covariate (which I really don't want to do) but even then the model fails:

2) In the greatly simplified model I get the Warning: Stopped because of infinite likelihood  (actually, I got this after I upped the maxiter to over 1000).

So the issue really seems to be inclusion of the Fixed Effect covariates in the CLASS statement.  I don't fully understand why this is necessary -- but whenever I don't do this, I return to my original problem: zero estimated random effects.

Could you explain why it is necessary to include the FE-Vars in CLASS?

Given these results, would you recommend trying to shift into PROC HPMIXED?  I've already started to try this, but I'm getting similar issues with infinite likelihood values (i.e., No valid objective function could be found) or a Non-PD G matrix....

Thanks very much for your help!

[Rick_SAS]

Most likely you are running out of memory because you put a continuous variable on the CLASS statement. Only put the categorical (discrete) fixed effects on the CLASS statement. Include the continuous effects in the MODEL statement, but not in the CLASS statement.

[SteveDenham]

@Rick_SAS has pointed this out, and it is where you were at the start. It isn't necessary to include them as CLASS variables, provided it makes sense to consider each of them as a continuous regressor, as you are not going to be able to get point estimates without extensive use of the AT= option in the ESTIMATE and LSMESTIMATE statement.  If you are going to go this route, consider centering all of them (but don't rescale). 

Now as far as the infinite likelihood problem.  That is most often due to multiple observations per subject at a time point.  If this occurs, it will show up right away, and never get started.  However, if you hit an infinite likelihood after 1000 iterations, it is something else, most likely a near floating point overflow.  Recentering your covariates will help.  

Also, are the covariates time-varying?  If they are not constant across all time points, then this becomes a very different problem.  You may want to look into PROC PANEL in SAS/ETS to get around this.

Finally, some basic info--how many levels of qtr do you hope to fit, how many fixed effects do you hope to include in the model, and how many subjects do you have?  I probably should have asked that to begin with.

Steve Denham

[TBeans]

@Rick_SAS That makes sense -- I figured it was the inclusion of a continuous variable in the CLASS statement that was eating up the memory.  I do have some non-continuous regressors that I guess I should include the CLASS statement.  But the failed model estimtations I referred to in my earlier reply did not include these discrete variables.... the final model woud/should, so thanks for clarifying!

@SteveDenham My observations -- 11 cities which are my SUBJECT, and each contain (in this particular version) 57 quarters of data.  I'm trying to pull out Mixed Effects specific to each city.  

By way of background: Before attempting an explicit AR(1) residual covariance specification, I implemented a Crossed Effects model in which I independently grouped by City and Quarter.  I was doing this in R (for the the lme4 package makes such implementation easy), I specifically changed to SAS upon learning of its wider arrange of covariance structure otpions.  In the crossed effects model, the idea was that controlling for cross-sectional variations (i.e., the qurater groups) helped control for some of the serial/AR correlation.

 Regarding your other questions:

A) I don't have multiple observations: Each observation is 1 city in 1 quarter.

B) All of the covariates are time-varying.

C) I would like to include six (6) or seven (7) Fixed Effects, and then permute which covariates I treat as random effects.

I'm fairly new to SAS, but I haven't looking into PROC PANEL because I'm specifically intersted in extracting the CITY-level effects, and -- AFAIK -- Panel methodologies will either be "fixed effects" in which I can only pull out a constant term (and not a "slope" estiamte) for each group, and "random effects" using the FGLS approach (and is rather similar to the LMM approach) doesn't give you a way to pull out group-specific estimates.  Why do you say it becomes a very different problem?

I'll try recentering the variables!  That's a great idea -- and now I remember having to do this rebasing in R sometime ago.

Please let me know what you think is reasonable or illogical in what I'm thinking now!

Many thanks to you both!

[SteveDenham · Accepted solution]

I'm really at a loss then.  Is this something like what you want to do:

PROC MIXED data= InDat.;
          class city qtr ;
          MODEL response =  fixed_effects qtr fixed_effect*qtr   / DDFM=KR Solution;
          RANDOM intercept /  subject=city Solution;
 
          REPEATED qtr  / TYPE=ARMA(1,1) subject=city;
run;

You may need to move to HPMIXED or HPLMIXED to get anything to work.

Steve Denham

[TBeans]

Sorry for the long delayed reply.

Yes in the end moving to PROC HPMIXED overcame the limitations I was facing with the MIXED procedure.

Thanks.