I think the PARMS statement in GLIMMIX doesn't like parentheses. Try this PARMS statement:   PARMS .,.,., 0, 0,., 0, 0, ., ., 0, 0, ., ., . /hold=3,4,6,7,10,11; If that does not work, then I think you may have to pass the values in a dataset. Non-zero entries can be obtained by fitting a UN covariance structure to separate analyses using a BY year; statement. Note that the matrix must be a full 15x15 matrix. See the example at this link in the PLM documentation: https://documentation.sas.com/doc/en/statug/15.2/statug_plm_examples07.htm .  I realize this is getting pretty complicated.   SteveDenham   ... View more

You could do that, but be aware that the covariance structure will include estimates of covariance between individuals in Year 1 and individuals in Year 2 with this approach. They are obviously independent and there should be no covariance.   Let me check into why the PARMS statement is not behaving the way I thought it would.   SteveDenham ... View more

Least squares means are not arithmetical averages. They are the expected values IF the data are balanced. You have missing values.  Consequently, the lsmeans will not be equal to the means.   SteveDenham ... View more

As I look over the output, I see that you have 83 observations, coded as VID, but not sequentially, based on the PROC PRINT. I am going to assume that GROUP is some sort of pen/pasture grouping of animals that have all of the 8 possible combinations of sex, frame and supplement level represented. The key here is that YEAR and PERIOD are completely confounded, per @mkeintz 's comment--if an observation has PERIOD = 1 or 2, it must have YEAR = 1, and if PERIOD is in 3, 4,or 5, then the record must have YEAR=2.   So you need to remove YEAR from the MODEL.  We will get back to how to test for a year effect after we go through the other factors. Sex, Sire, Period and Trmt are all fixed effects, with GROUP a random effect. That means that we will need to get creative outside the MODEL and RANDOM statement. First, let's look at how to model Period as a REPEATED effect. With only 2 or 3 levels per year, I would try the following structure for the covariance matrix. First, the covariance between periods in different years should be fixed at zero. Next, the number of parameters for YEAR=1 is 3 (2 variances, 1 covariance) and for YEAR=2 is 6 (3 variances, 3 covariances). If I fit something like that I would use the PARMS statement with a HOLD option. So, before we do LSMEANs and SLICES, the code I am suggesting would look like:   PROC MIXED; Class Group Sire Sex TRMT Period VID; Model kgPERdayForageDMI=Sire Sex TRMT|Period; random Group; REPEATED/type= UN r rcorr subject=VID(period); PARMS (.) (.) (.) 0 0 (.) 0 0 (.) (.) 0 0 (.) (.) (.) /hold=3,4,6,7,10,11; RUN; If that runs without errors or warnings, then the following code can be inserted:   lsmeans TRMT TRMT*Period/pdiff=ALL adjust=TUKEY; SLICE TRMT*Period/SLICEBY=period diff adjust=TUKEY; LSMESTIMATE period 'Year 1 averaged over TRMT' 1 1 0 0 0 divisor=2, 'Year 2 averaged over TRMT' 0 0 1 1 1 divisor=3, 'Year 1 minus Year 2, =Year effect' 3 3 -2 -2 -2 divisor=6; The LSMESTIMATE statement calculates the year 1 mean, year 2 mean and the difference, which tests the null hypothesis that the Year means are equal, holding all other parameters at their mean values.   SteveDenham         ... View more

I am not really sure what is causing the difference. It could be a difference in algorithm, such as a difference in whether Period is handled as another repeated measure, with time nested in period.  To check you might try the following:   proc mixed data=Result maxfunc=5000 convf=1e-6; class subject period treatment time; model result = baseline period subject treatment time treatment*time/ddfm=kr2; repeated period time / type=UN@CS subject=subject; lsmeans treatment / pdiff=control; contrast 'linear' treatment -3 -1 1 3; lsmestimate Treatment*time 'Linear treatment at time 1' -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0 0 0 0, 'Linear treatment at time 2' 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0 0 0, 'Linear treatment at time 3' 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0 0, 'Linear treatment at time 4' 0 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0, 'Linear treatment at time 5' 0 0 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0, 'Linear treatment at time 6' 0 0 0 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3/FTEST; run; quit; On my machine, the linear treatment at time 1 comes in with a P value of 0.0019, so I applied a multiplicity adjustment:   proc mixed data=Result maxfunc=5000 convf=1e-6; class subject period treatment time; model result = baseline period subject treatment time treatment*time/ddfm=kr2; repeated period time / type=UN@CS subject=subject; lsmeans treatment / pdiff=control; contrast 'linear' treatment -3 -1 1 3; lsmestimate Treatment*time 'Linear treatment at time 1' -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0 0 0 0, 'Linear treatment at time 2' 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0 0 0, 'Linear treatment at time 3' 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0 0, 'Linear treatment at time 4' 0 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0, 'Linear treatment at time 5' 0 0 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0, 'Linear treatment at time 6' 0 0 0 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3/FTEST adjust=simulate(seed=1) stepdown adjdfe=row; run; quit; This gives an adjusted P value of 0.0109, which is within striking distance of the 0.0105 you are looking for, so it may have a different adjustment method   SteveDenham ... View more

I think you just plug in the value of interest.  Here is code, with the two control rates and the two treated rates you have:   proc power; coxreg sides=2 alpha=.05 hazardratio=.564 stddev=0.477 power = 0.8 eventprob = 0.15 0.18 0.2 0.25 ntotal= .; run; The only thing that bothers me is the hazardratio estimate, but if you are comfortable with that, this should work for you.   SteveDenham   ... View more

I haven't done genetics by environment interaction analyses in years. The best I can come up with to get started in SAS is actually PROC HPLMIXED. The one example is Computing BLUPs for a Large Number of Subjects. If that example isn't really applicable due to its assumption of known variance and covariance values, then HPMIXED would almost certainly be a better choice than MIXED.  And I realize that my answer is pretty vague.   SteveDenham ... View more

I think you have it down for GEE and GENMOD. GLIMMIX allows for a variety of correlated data, including multilevel effects. It also deals with missingness up to missing at random (it does not eliminate records that have missing values for model factors)  Use the first example in the PROC GEE documentation for a good comparison of marginal and random effect models. I have two concerns. The first is that I am not clear on the use of Drug as a response variable when your research question talks about comparing levels of Drug. I would consider Drug as a fixed effect to be included in the model, and the response to be something measured on the patient (cured-not cured, for example). The second is that I cannot tell if there are two levels of clustering here--patient level and hospital level. If each patient is measured one time then there are no patient level clusters - the patient level effect is the "residual" or scale estimate. Since you mention that the design is repeated this would be treated as a patient/student level cluster. However, hospital needs to be specified as either fixed (inference space is then repeated studies at the specified hospitals) or random (inference space is repeated studies at the greater population of hospitals, of which the ones in the data represent a "random" sample).   From this, I can see a PROC GEE approach for the narrow inference space of the sample of hospitals, or a PROC GLIMMIX approach for the broad inference space of "all hospitals".   SteveDenham ... View more

First, I want to support what @PaigeMiller said about residuals being normally distributed. For the non-Gaussian distributions you were examining, the check on normality in the residuals isn't as necessary, but it does let you know if the variance components are normally distributed.   Some other things to note: It appears the response variable is a count, which would imply that a Poisson  or negative binomial distribution is a likely candidate, unless the counts get quite large, where the gamma treats the discrete values as continuous. No need to check binomial or beta, as they only have support on the zero to one interval.   Looking at your code, I don't see an NLOPTIONS statement.  Many times for Poisson or negative binomial, you are going to need more than the default 20 iterations to get convergence. We routinely use NLOPTIONS=5000. The rest of the code looks like it should give you comparisons between the levels of Trmt at each Week, so I am a bit confused that the gamma distribution approach did not give you all of the lsmeans.  Check the output for any messages that may be printed, as well as for WARNINGs or NOTEs in the log.   One more thing to consider is a generalized estimating equation approach using PROC GEE. Your current code is fitting a marginal model in GLIMMIX, so switching to PROC GEE (or GENMOD for a wider selection of possible distributions) may solve your problems. Note that the LSMEANS statement in GEE doesn't support the slicediff= option. LSMESTIMATE or SLICE statements should give you the specific within week comparisons you want.   SteveDenham ... View more

You get eventprob in one of two ways. If you have historical or analytical results available, you can use that value. If you don't, but do have some idea of the range where it might be, you can plug in a "row vector" of possible values, and PROC POWER will calculate ntotal for each of those values, given the other parameters.   SteveDenham ... View more

If WEEK is dichotomous, you might try removing it from the CLASS statement to see if that addresses your issue.  I am not familiar with any algorithms that apply exact tests in generalized estimating equations.   SteveDenham ... View more

    Here you go:   lsmestimate Treatment*time 'Linear treatment at time 1' -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0 0 0 0, 'Linear treatment at time 2' 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0 0 0, 'Linear treatment at time 3' 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0 0, 'Linear treatment at time 4' 0 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0 0, 'Linear treatment at time 5' 0 0 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3 0, 'Linear treatment at time 6' 0 0 0 0 0 -3 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 3/FTEST; From the documentation: "The FTEST option produces an F test that joint test the rows of the LSMESTIMATE against zero." I assume that this is what you are looking for when you ask for the correct overall linear time by treatment interaction. Note that actually this is a time by linear treatment interaction F test. If you want a linear time by linear treatment overall test, it gets a bit more complicated.  We use a macro function to generate the string of coefficients.  In this case, the matrix is the cross product of the linear time vector, which for 6 timepoints is -5 -3 -1 3 5 and the linear treatment vector -3 -1 1 3; to get:   lsmestimate treatment*time 'linear treatment by linear time' 15 9 3 -3 -9 15 5 3 1 -1 -3 -5 -5 -3 -1 1 3 5 -15 -9 -3 3 9 15; I hope you find this useful.   SteveDenham ... View more