The two models enable you to test different things, so it boils down to what your experimental question is. If you are only interested in the effect of season, model 2 is appropriate. If you are interested in what the lactation curve looks like and whether there are differences between seasons at various sampling times, then the first model is appropriate. Note that the first model is a split-plot in time, with the whole plot being exactly what model 1 tests (at least for balanced data).  I would recommend finding a copy of SAS for Mixed Models (any edition, but the third edition is best), and particularly Chapter 8: Analysis of Repeated Measures Data for more information.   SteveDenham ... View more

Consider the following: Suppose you did some sort of permutation test, where each observation has an equal probability of being assigned to control, treatment 1 or treatment 2. The specified comparisons make sense under this approach. Now suppose you did the same permutation test, but with restricted sampling so that only control and one of the treatment arms are used. My brain says to me - "Hey, these are not the same thing. Maybe I could expect something different from the two analyses."  And in the end the asymptotic tests, such as those in PHREG, should approach in probability the results with an infinite sample size for the permutation test. (I don't remember where this comes from, but Fisher might have something to say about it).  Now it is very possible that my brain is missing a key element in this argument, in which case you should just ignore any of my speculations.   SteveDenham ... View more

And to add to @ballardw 's comments, calculation of power at this point has some issues as well. Unless you are going to use this calculation to power and set a sample size for an upcoming study, I would recommend against calculating power (1 minus the probability of a type II error). Why? Because you have already collected and analyzed the data. Either you rejected the null hypothesis or you didn't. Working through the probabilities, based on this piece of information, your power is either one or zero for this dataset and analysis. I refer you to Zhang et al. (2019) http://dx.doi.org/10.1136/gpsych-2019-100069 , who finish their paper with this statement: "In general, post hoc power analyses do not provide sensible results."   SteveDenham ... View more

Perhaps you could use the EFFECT statement to fit a spline to the start dates to deal with the suspected change in the treatment effect, and use @PaigeMiller 's suggestion of fitting elapsed days as the independent variable. You would probably need an interaction term. This would fit a response surface for your DiD data. Does that make any sense?   SteveDenham ... View more

Thanks, @Rick_SAS . That property was what I was thinking about when I said that if a smeared estimate is used for the expected value (= eta in your post), then it should also be applied to the confidence bounds. The Swiss Tech paper implies that a direct exponentiation should be adequate.  Thoughts?   SteveDenham  ... View more

A PROC SUMMARY or MEANS involving all covariates would be a good check for empty cells due to cross tabulating, so maybe run that before trying to get the mixed model running.   SteveDenham ... View more

Thanks, @jiltao. That is likely the case, especially as I left out a step. After calculating the sum of squared deviations and before taking a square root, there needs to be a division by an appropriate degrees of freedom value. So we would have the deviations of the values predicted by the model from the observed values, each then squared, summed over all observations, divided by an appropriate, model-based, degrees of freedom value. That looks, at least to me, a lot like a variance due to things not in the model like unmodeled fixed or random effects. If the model was as simple as a mean, it would be the variance, wouldn't it? Taking the square root gives a standard deviation, or in the case of a general linear model, the RMSE.   Since the predicted values are empirical BLUPs, this calculated value is a measure of how closely the empirical BLUPs represent the variability in the raw data.  I really need my copy of Graybill's Theory and Application of the Linear Model to refresh my BLUP knowledge, but it is in a box in the basement somewhere...   SteveDenham     ... View more

This is probably more for @StatsMan , @sbxkoenk  and @jiltao  than the OP.  Could an RMSE-like value be calculated from the residuals obtained from using the OUTPRED= option from the MODEL statement? I can conceive of squaring the raw (unscaled) residuals (observed - predicted), summing them up and taking the square root.  This seems logical, but if it were, there would be a lot more of it done, so I am obviously missing something.  Anyone care to take a stab at this?   SteveDenham  ... View more

You are probably going to need to expand on this, but as a first step why not drop back to      random intercept/subject=blk; and see if that gets the algorithm going. I am a bit confused by having as many random effects as you specify in the current RANDOM statement - you may be exhausting the available data.  Additionally, specifying Date as a repeated measure and then including all the Date related factors, including Date itself, as being nested within blk is probably causing singularity issues.   SteveDenham   ... View more

This is a case where the MIXED approach can be generalized (see what I did there) by using GLIMMIX with a binomial distribution.   SteveDenham ... View more

Remember, this model is giving you the expected sales, conditional on the monthly sales up to the previous month. It is therefore not unlikely to get a lower value than the current sales, especially if there has been any sort of level change in an unmonitored (and thus not included in the model) variable. I guess you could have some dichotomous check in the output that if the predicted value is less than the current value that the output would default to the current value.   SteveDenham ... View more

What sort of measure do you want to use for performance? Machine performance can be obtained from the log, in terms of cpu time and elapsed time. Model performance for PHREG is harder as there isn't a real good equivalent to the RMSE (smaller is better in terms of model performance). I suppose you could look at the relative sizes of the standard errors for various parameters, but that overlooks the danger of overfitting. With only two independent variables that may not be an issue, however.   SteveDenham ... View more