Distribution for percentages in proc genmod

[Season]

Thank you for your correction! Yes, I did not remember your message precisely. Sorry for the inconvenience it might cause. I have edited my original post accordingly.

[palolix]

Thank you for your suggestion. Yes, the AIC and BIC values are correct (negative) . The estimates and SE from the glimmix model using beta dist are much lower than the ones using genmod and normal dist.

I added to the glimmix model 

output out=xxx resid=r;

 but nothing happened.  Maybe I am missing something?

proc glimmix data=one;
where Season=2021;
PercDMp=PercDM/100;
class Harvest Variety;
model PercDMp=Harvest*Variety/ dist=beta ddfm=kr;
output out=xxx resid=r;
run;

Thank you very much Season!

[Season]

---

@palolix wrote:

Thank you for your suggestion. Yes, the AIC and BIC values are correct (negative) .

---

That is interesting. You might have noticed a message read "(smaller is better)" following "AIC" and "BIC" in SAS output. This indicates that models with smaller AIC or BIC values are better-fit models than those with larger values of them. Given that the AIC's and BIC's of beta regression model are both negative while their counterparts for the normal model are positive, it is demonstrated by the two statistics that the beta model fits the data dramatically better than the normal model does.

---

@palolix wrote:

The estimates and SE from the glimmix model using beta dist are much lower than the ones using genmod and normal dist.

---

That is yet another interesting discovery. So it is the case that model misspecification leads to enlarged standard errors in your data. Given that the resultant P values are similar, I speculate that the regression coefficient estimates of the normal model are larger than those of the beta model.

---

@palolix wrote:

I added to the glimmix model 

output out=xxx resid=r;

 but nothing happened.  Maybe I am missing something?

proc glimmix data=one;
where Season=2021;
PercDMp=PercDM/100;
class Harvest Variety;
model PercDMp=Harvest*Variety/ dist=beta ddfm=kr;
output out=xxx resid=r;
run;

Thank you very much Season!

---

The OUTPUT statement writes a dataset to a SAS library where datasets are stored. In response to your code, the dataset named "xxx" is stored in the Work library. You can visit that library to find it.

[palolix]

Ok, thank you Season!

[SteveDenham]

Just for fun, consider a modeling approach that doesn't assume homogeneity of variance. Working from your GLIMMIX code, try something like:

/* Changes are in red  */
proc glimmix data=one;
where Season=2021;
PercDMp=PercDM/100;
class Harvest Variety;
model PercDMp=Harvest*Variety/ dist=beta ddfm=kr2;
random _residual_ / group=Variety;
lsmeans Harvest*Variety/slicediff=Harvest adjust=simulate(seed=1);
covtest 'common variance' homogeneity; 
run;

If it turns out that the likelihood ratio test for variance homogeneity for Variety is not significant, try it again grouping by Harvest. I really don't know if either will affect your conclusions, but at least you have dealt with a common assumption (homogeneity) that may not be true for your data. 

The other thing to look at is to try a different method than the default RSPL. Consider method=laplace, so that the error variance component is included in the optimization. That may yield standard errors that are more in line with expectations. 

SteveDenham

[palolix]

Thank you so much Steve for keep helping me. Ok, I tried that code and the likelihood ratio test for variance homogeneity for Variety was significant so I didn't change the group. I am still getting the same significant p-values. I tried laplace but same thing. 

I am still surprised that I am getting these significant p-values with whatever model I try (see asteriks in graph, it doesn't look to me that those differences in perc dry matter between the varieties are significant, specially when adjusting the p-values). 1,2 and 4 in 2024 are the harvest months.