topic Re: ANOVA using SAS 9.3 Proc Mixed in Graphics Programming

ANOVA using SAS 9.3 Proc Mixed

TD21 — Tue, 11 Feb 2014 21:17:26 GMT

I am doing ANOVA using SAS 9.3 Proc Mixed on a field data (RCBD with four treatments (One control vs. three different land uses). The response variable is water table depth (WTD) which is collected every 15-minute. I would like to know which of the two sets of codes is more appropriate to use. I would also welcome any suggestion to improve the set of codes. Seas = Seasons (traditional: winter, spring, summer, and fall), Day = Julian Day, WTD = water table depth averaged over a day. Here are the codes:

ODS GRAPHICS ON;
PROC MIXED DATA=WORK.WTDWOR_2014 plots(only)=(ResidualPanel(marginal))
PLOTS(MAXPOINTS= 2000000);
CLASS Block Treat Year Seas Day;
MODEL WTD = Treat|Year|Seas/ddfm=satterthwaite residual;
RANDOM Block Block*Treat Block*Treat*Year Block*Treat*Year*Seas;
LSMEANS Treat Treat*Year Treat*Year*Seas/adjust = tukey;
RUN;
ODS GRAPHICS OFF;

RANDOM Block Block*Treat Block*Treat*Year Block*Treat*Year*Seas;
LSMEANS Treat Treat*Year Treat*Seas(Year)/adjust = tukey;
RUN;
ODS GRAPHICS OFF;

P.S. On set of codes B, I was trying to let SAS know that season(Seas) is embedded within a year. Is this correct?

Thanks.

Re: ANOVA using SAS 9.3 Proc Mixed

Rick_SAS — Tue, 11 Feb 2014 22:52:51 GMT

Not sure if monitors this forum. This question might get more answers from the folks in the Statistical Procedures Community: https://communities.sas.com/community/support-communities/sas_statistical_procedures

Re: ANOVA using SAS 9.3 Proc Mixed

TD21 — Tue, 11 Feb 2014 23:29:13 GMT

Thanks Rick...

Re: ANOVA using SAS 9.3 Proc Mixed

SteveDenham — Wed, 12 Feb 2014 12:20:42 GMT

Some interesting points to consider here:

1. You include DAY as a class variable, but do not use it in the model. As the analysis now stands, you are averaging over all days, even though they are in different years and seasons. This may lead to problems.

2. Of the two models presented, I would think they are identical. Your code for the LSMEANS is certainly appropriate in either case, given how SAS models nested factors (see the documentation on this) SEAS(YEAR) and SEAS*YEAR set up exactly the same matrix.

3. I guess my main concern is that this is a repeated measures design of some complexity, and the correlations within a plot over time are not being modeled. It may be that this just can't be done easily. You may wish to think of season-day, if the measurements are spaced out at the same times within each season for each year. However, I have a strong hunch that this level of complexity may not be amenable to modeling, due to memory and other constraints.

Let us know what you think.

Steve Denham

Re: ANOVA using SAS 9.3 Proc Mixed

TD21 — Wed, 12 Feb 2014 14:44:34 GMT

Thanks Steve. As you may have recognized, I am collecting the data every 15 minutes. I am both interested in main effects and simple effects (of treatments on WTD on annual basis, and on seasonal basis within a given year, if any). I average 240 samples of 15-minute WTD data to come up with an average for a day of a given season per treatment. The raw data are already on in terms of average per day and not the 15-minute interval. Analyzing by season (traditional 4 seasons) means averaging 120 days of WTD into one value per treatment. My point is, I could add a repeated statement in the model; however, the fact that I am averaging the data from 15-minute interval (240) to daily (120) to season makes the autocorrelation issue less of an issue.

You are right with point number 2. I looked at the covariance parameters estimates and fit statistics (e.g. AIC values), the two models were very similar. I just would like to tell the model that DAY is nested within SEAS, and SEAS within Year and I am not sure if I am doing it right. Based on your comment about similarity in matrices, I think I am close.

You mentioned about DAY. Any suggestion as to how to incorporate that into the model (DAY as nested into SEAS) as a code?

Again, thanks for your time.

Re: ANOVA using SAS 9.3 Proc Mixed

SteveDenham — Fri, 14 Feb 2014 18:17:06 GMT

TD,

The nested terms in your model are already included when you 'pipe' the fixed effects, so I am now in favor of Model A (foolish consistency is the hobgoblin of small minds).

If you have the computing power to handle DAY as a repeated measure, I would try it. How many days of data do you have?

Things to consider:

I would shift to HPMIXED, and then port to GLIMMIX for comparisons of lsmeans (see EX. 46.3 in the HPMIXED documentation)

I would index the days within season, rather than using the Julian date. If they don't align, then forget using a repeated measures approach. If they do align, or would with only minor shifts, then code that looks like:

CLASS Treat Year Seas Day Block Plotid;

MODEL WTD = Treat|Year|Seas|Day ;

RANDOM Block Block*Treat Block*Treat*Year Block*Treat*Year*Seas;

REPEATED Day/subject=plotid(block treat year seas) type=AR(1);

Hope this helps.

Steve Denham

Get everything in ODS datasets and then port it all into GLIMMIX.

Steve Denham

Re: ANOVA using SAS 9.3 Proc Mixed

TD21 — Fri, 14 Feb 2014 19:26:02 GMT

I have 92 days in 2009; 365 x 3 = 1095 days in 2010, 2011, and 2013; and 366 in 2012. So in each treatment I have 1553 days of data. By the way, I didn't include Julian date in my dataset but rather Julian Day. In short, I have a column for year, separate column for season, and day. For example in 2010, I have column/variable for year which is 2010, another column for season (winter, spring, summer, and fall) and another season for Day (1 to 365).

Thanks again.

Re: ANOVA using SAS 9.3 Proc Mixed

SteveDenham — Mon, 17 Feb 2014 16:22:56 GMT

Given that data set up, and the already "binned' nature of the data, would yet another binning by week within season make sense? Then you would have 13 repeated measures on each combination of year, season, block and treatment. If you believe that this is something that reflects an interesting approach, I think it would be the best way to incorporate changes over time. The problem I have with both season and day in the model is that the number of observations per season (and per year) isn't consistent because 4 doesn't divide 365 or 366 equally.

An alternative would be to fit day rather than season, and then address season effects through the use of LSMESTIMATE or CONTRAST statements, viewing the season effect as a composite of the day effect. It will require HPMIXED, followed by GLIMMIX to get this to work though.

Steve Denham

Re: ANOVA using SAS 9.3 Proc Mixed

TD21 — Sun, 06 Apr 2014 03:53:51 GMT

Hi Steve:

I followed your latest advice (add week, analyze with HPMIXED then port to GLIMMIX). I seemed to have successfully run HPMIXED, but not with GLIMMIX. The error was about PARMSDATA which I highlighted below. Any suggestion how to fix this error? Pardon the long post. I just wanted to give you as much background information as I can.

Thanks.

55   PROC HPMIXED DATA = WORK.WTDWOR_Final2014; *plots(only)=(ResidualPanel(marginal)
NOTE: Writing HTML Body file: sashtml1.htm
56                       StudentPanel(conditional)) PLOTS(MAXPOINTS= 2000000);
57       CLASS Block Treat Year Seas Week Day;
58       MODEL WTD = Treat|Year|Seas|Week;
59       RANDOM INT Treat/SUBJECT=Block;
60       REPEATED DAY/SUBJECT=Block*Treat*Year*Seas*Week type=ar(1) R;
61       TEST Treat|Year|Seas|Week;
62       ODS OUTPUT COVPARMS=HPMEstimate;
63   RUN;

NOTE: 129 observations are excluded because of: missing response values (n=129).
NOTE: Convergence criterion (GCONV=1E-8) satisfied.
NOTE: The data set WORK.HPMESTIMATE has 4 observations and 3 variables.
NOTE: PROCEDURE HPMIXED used (Total process time):
real time 2:32:48.21
cpu time 1:24:56.25

64
65   ODS GRAPHICS ON;
66   PROC GLIMMIX DATA = WORK.WTDWOR_Final2014;
67       CLASS Block Treat Year Seas Week Day;
68       MODEL WTD = Treat|Year|Seas|Week;
69       RANDOM INT Treat/SUBJECT=Block;
70       PARMS/PDATA=HPMEstimate hold=1,2 noiter;
71       LSMEANS Treat|Year|Seas|Week/ADJUST=Tukey;
72   RUN;

NOTE: Some observations are not used in the analysis because of: missing response values (n=129).
ERROR: The PARMSDATA= data set must have either 2 more or 1 less observations.
NOTE: The SAS System stopped processing this step because of errors.
NOTE: PROCEDURE GLIMMIX used (Total process time):
real time 40.35 seconds
cpu time 37.02 seconds

Re: ANOVA using SAS 9.3 Proc Mixed

SteveDenham — Wed, 09 Apr 2014 17:48:31 GMT

It looks to me like you aren't fitting the same model in GLIMMIX. Try adding a second random statement

RANDOM DAY/residual SUBJECT=Block*Treat*Year*Seas*Week type=ar(1) R;

and changing the PARMS statement to:

PARMS/PDATA=HPMEstimate hold=1,2,3,4 noiter;

Steve Denham