Multi-level Mixed modeling

[Bhupi]

Hi I have dataset on evaluating the performance of maize and three sorghum species grown on marginal lands where maize is supposed to perform poor. We wanted to find out if any sorghum species perform better and can be recommended to grow for biofuel production. We knew that maize wouldn't grow well even if we supply extra nitrogen fertilizer. Therefore we kept nitrogen fertilizer source as main plot variables and crop species as sub-plot. We wanted to see if any of the crop species would perform better than others including maize and if that consistency can be achieved with nitrogen management and across years. The study was conducted from 2011 to 2015 and at three locations which are little bit different from each other based on their marginality ratings, however, that was not the criteria of their selection. Several variables were measured over years that include yield, structural and non-structural carbohydrates. We are also interested to find out that what is the optimal nitrogen fertilizer rates for each crop species and should I consider only polynomial function or is there any easy way to fit several models and compare (linear, quadratic, linear plus plateau, and quadratic plus plateau). A follow-up question when we fit any of these model, it is possible that one crop species have quadratic response but other have quadratic plus plateau or linear response. Can we fit different functions for different cultivars?

I don't know if what should be the best approach to carry this analysis and answer above mentioned questions.

I look forward to receive your valuable suggestion.

Here is the model that I have written as per my current knowledge.

Proc MIXED data=Final1 namelen=200;
class Year Rep Location N_Appl__kg_ha Species;
model Tot_Biomass__Mg_ha = Year|Location|N|N|N|Species/ddfm=kr;
Random intercept N_Appl__kg_ha/subject=Rep(Location);
Repeated year/type=un subject=N_Appl__kg_ha*Species*Rep(Location);
run;

Thanks,

Bhupinder

[sld]

Hello Bhupinder,

I'd like to clarify some details about your design.

How many REPs were there at each location?

Is the N fertilizer applied to REPs or to whole plots within a REP? 

How many levels of N fertilizer are there?

Did you use the same physical experimental units (REPs, whole plots, subplots) each year?

Did you reapply N fertilizer each year? If so, did each REP or whole plot (whichever is appropriate) receive the same level of N fertilizer each year?

Was each subplot planted with the same species each year?

Susan

[Bhupi]

Hi Susan,

Thanks for taking time and putting thoughts on this issue. Here are the responses to your queries.


How many REPs were there at each location?
Four

Is the N fertilizer applied to REPs or to whole plots within a REP?
N fertilizers rates were applied to whole plots within each rep.

How many levels of N fertilizer are there?
Five (0, 50, 100, 150, and 200 kg N/ha)

Did you use the same physical experimental units (REPs, whole plots, subplots) each year?
Yes experimental units were same each year. Plots were established in 2011 and received same treatment until the end of study (2015).

Did you reapply N fertilizer each year? If so, did each REP or whole plot (whichever is appropriate) receive the same level of N fertilizer each year?
Yes, N fertilizer was applied each year to respective plots. And, yes each plot received same fertilizer rate and same crop each year during the study period.

Was each subplot planted with the same species each year?
Yes. For example, if a plot was assigned corn and n rate 150 in 2011 then same N rate was applied and same crop was planted in subsequent years. Similarly, if another plot was assigned N rate 0 and crop as dual purpose sorghum in 2011 then same treatments were applied each subsequent years.

All above given responses are applicable to each site.

One unique observation: In 2015, at one site, all the crops failed that means we didn't have harvest that year at that site. This has happened only at one site which is supposed to be least marginal among three sites. This shows that there is a potential risk of growing crops at marginal sites and at this moment I am thinking to consider as part of the analysis rather than excluding from the overall analysis. As we were growing these species continuously over years, so this failure is not surprising.

Please let me know if you have more questions.

Thank you,
Bhupinder

[sld · Accepted solution]

My apologies for taking so long to get back to this. I have some ideas, but no clear solution, which slowed me down.

The MIXED code that you propose for the cubic regression on nitrogen level looks fine to me. I'd use a centered value of N to minimize potential problems with collinearity and make interpretation potentially easier. Because the model regresses on N, at some point you could ponder a covariance structure with random slopes.

Another way to invoke a polynomial regression is with the EFFECT statement, which is available for GLIMMIX but not MIXED:

http://blogs.sas.com/content/iml/2017/09/07/polynomial-effects-regression-sas.html

and

https://support.sas.com/resources/papers/proceedings13/437-2013.pdf

To fit plateau models, you'd have to switch to NLMIXED.

This paper https://support.sas.com/resources/papers/proceedings16/SAS4720-2016.pdf provides an illustration of hierarchical models.

A "random parameters" (aka random coefficients) approach is intriguing. In this paper  http://ageconsearch.umn.edu/bitstream/117949/2/jaae417-Tumusiime.pdf  (which illustrates a linear plateau model, among others) there were only 6 treatment levels of N, so maybe there would be some hope with 5 levels. 

This paper illustrates repeated measures in NLMIXED https://link.springer.com/article/10.3758%2Fs13428-013-0397-z

Hopefully this is not redundant to resources that you've already uncovered. It's a challenging analysis, and I'm just speculating fairly wildly here. 

Susan

[Bhupi]

Hi Susan,

Thanks for putting thoughts into my analysis. I'm actually now running the random parameter approach. I did have some hiccups in running it successfully but not bad. I'm getting help on it. i will soon see comparison of nested and non-nested models to pick the best fit.

Regards,

Bhupinder

[sld]

Hi Bhupinder,

Cool! I'd be interested to see what you come up with, if you'd like to share.

Susan