Existence of response space + variance contributions
The analysis done at site level can also be done for the overall dataset (all sites together). As in sitelevel analysis, scatter plot of yield gain due to treatment application versus the average control yield were made in R (see Yield analysis 112012.R; Figure 10). There was high variability in crop response to fertilizer with some fields being highly responsive to fertilizer and even doubling of yields, while others had minimal or no response to the fertilizer. At every level of control yield, there is often a yield response range of upto 6 t/ha. This huge “response space” indicates opportunity for important gains through appropriate fertilizer use and amendments. Some sites such as Pampaida had a limited range of control yields (maximum 1.7 t/ha) while sites such as Kiberashi and Sidindi showed a huge range in control yields. Using these three sites as examples (all seasons can be considered normal), the range of control yield may be dictated by low organic resource uses (Pampaida) and gradients influenced by conversion history (Kiberashi) and varying organic resource management by farmers (Sidindi). Thus the data was subjected to further analysis for assessing variability and response patterns based on the combined responses to the fertilizers and amendments. The individual treatment responses within the response space are similar to those shown earlier for sitelevel analysis which shows evidence of control effect.
Figure 10. Variability in treatment yield response at different levels of the unfertilized control treatment yield as observed in (maize) sentinel sites in subSaharan Africa, 20092012
As yield data is available for both maize and sorghum sites, a single model was fitted which included a term for crop, as well as the interaction between crop and treatment which allowed for different treatment effects to be estimated for the two different crops. The interaction term in this model was significant, indicating that the crops responded differently to the treatments. However the standard errors for the Sorghum treatment estimates were a lot larger than those seen in the individual site analyses, and resulted in no differences being found between the treatments (which does not match the site specific analysis). As a result of this we fitted separate crosssite models for the two crops, and compared the variance components. The variance components for the maize site are much larger than those from the sorghum sites, and therefore dominated the crosssite analysis (forcing the SEs for the sorghum estimates to be much larger than in the individual crop model). Thus we considered that the individual crop models should be used and the following standard models are used in the analysis depending on seasons and the need to adjust for the control variable:
lmer(response~ 1 + TrtDesc + (1Site) + (1sitecluster) + (1siteclustfield), data=cropdata, REML=T) #for crops with only one season
lmer(response~ 1 + TrtDesc + (1Site) + (1sitecluster) + (1siteclustfield) + (1siteseas), data=cropdata, REML=T) #for crops with multiple seasons
lmer(response~ 1 + control + TrtDesc + (1Site) + (1sitecluster) + (1siteclustfield), data=cropdata, REML=T) #for crops with only one season and adjusting for control variable
lmer(response~ 1 + control + TrtDesc + (1Site) + (1sitecluster) + (1siteclustfield) + (1siteseas), data=cropdata, REML=T) #for crops with multiple seasons and adjusting for control
Site, season and cluster show low variance contribution (Figure 11). Because the seasons across sites are not common, it is likely that some sitetosite variability is being attributed to seasons.
Figure 11. Variance contribution by different components in the diagnostic trial design based on maize grain yield.
As also in sitelevel analysis, the R scripts are designed to automatically make plot of treatment effects (see Figure 12) and ready to use tables such as table of treatment estimates and 95% confidence interval limits (back transformed onto Ratio Scale), table of variance components and table of correlation of fixed effects among others. This R code is designed as an R function (see “CrossSiteMixed_resent.r”) that is called/required by the main Rcode referred to as “Yield Analysis 112012.r” and available here.
The fixed effects in Figure 12 show that all of the treatments except PK are significantly different to the control. There are also some significant differences between the treatments: the NPK addition treatments are all significantly greater than the NK and PK treatments, the NPK addition treatments along with NPK and NP are also all significantly greater than PK.
Figure 12. Overall effect of treatments on yield increase over the control in SSA according to AfSIS objective 4 dataset
The term adjusting for a ‘control effect’ was significant, p=0.0291, with an estimate of 0.099, and a 95% confidence interval of 0.188 to 0.010. The interpretation of these results is: the estimated decrease in the difference between treatment and control is 0.099 ton per hectare for each additional ton per hectare that the control plot produces, this difference is estimated to be between a 0.01 and a 0.188 ton per hectare decrease. Note for sorghum, control did not have a significant effect on the model results so control variable was dropped from the analysis.
In the subsequent analysis, we predict the response that can be realized from a given treatment for a specific site where the level of control yield is known. The model used here is of the form: log(Yt)log(Yc)~Treatment*log(Yc) where Yt is yield of a treatment, Yc is yield of the unfertilized control. Two random terms were specified varying slope by the site (log(Yc)Site) and (1Field). Analysis of diagnostic trials has shown that cluster level is of less importance (contributes very low to the variance) and is therefore not included in the current model. Normality checks showed that residuals are slightly negatively skewed. At the present, no further transformations have been done.
The model results shown below for maize are accessed using display or summary functions. The full R scripts for this analysis is named “Yield_Prediction code_10_12.R” and can be downloaed by clicking here.
The results show the estimate and the standard error associated with each treatment and control yield. The estimates indicate the rate of increase of yield over the control. Taking NK treatment for example, the estimate shows that when control yield is at its average, a 49% maize yield increment over control is expected following application of NK (the intercept in this model). Since the model results are for the whole range of control yield, the slope term is contained in the estimate of the fixed control term. For the NK treatment under maize, this means that for every unit increase in control yield, response to NK is reduced by 30%. To aid easy interpretation of the above results, transform matrix was developed in R (see Yield_Prediction plotresults_V2_byCrop.r) basically calculating the expected estimate for each of the crop and treatment combinations and includes the associated slope.
The results for the whole range of control yields are plotted below (plotting scripts are included in the R code). As expected, response rate to nutrients and amendments decreased with increasing level of control yield (Figure 13). Steeper slopes indicate high response to a treatment relative to the control yield, used as an indicator of soil fertility, whereas flatter slopes indicate lesser response of a treatment.
Figure 13. Response rate of maize grain yields at different levels of control yield for treatments implemented in diagnostic trials.
Based on the predicted response rate in the above figure, we calculated the actual predicted response for each treatment at different levels of control yield. For this, we converted predicted responses from log scale to actual scale (by taking their exponents) to reflect the actual change in yields. For each level of control, the yield that can be expected from the application of a given treatment among those used in these diagnostic trials can be predicted (Figure 14). From the figure below, the prediction model may be further refined to avoid over/under predictions at high control yields.
Figure 14. Predicted yield change from the maize with NPK compared to the control treatment at different levels of control yield
The models used in the above crosssite analysis does not include soil parameters. Adding available soil information such as P, pH, Ca, Mg, Zn, S, B exchangeable acidity e.t.c., did not improve the models as yet (residual errors did not reduce and BIC/AIC were even higher in some cases). Also, including principal components (PCs) extracted from soil spectra as covariates in the mixed effects models did not improve the model. It is likely that the control yield factor in the model already takes into account the additional effects of these soil parameters. Soil carbon is still missing and it will be interesting to see if it will improve the prediction models, and further exploration with spectraextracted PCs is still necessary as spectra from additional sentinel sites becomes available.
The Rcode transform matrix that also plots the above graphs is available for download here and is named “Yield_Prediction plotresults_V2_byCrop.R”.
The above results do not discriminate between different sites. Using the crop response model cited above, the response of treatment in reference to the control can be predicted also at the site level. The results of this analysis are presented in Table 4. In short, for each sentinel site, the intercept and the slope associated with each treatment and site were extracted from the model, and predicted response calculated using the overall mean of the log control. We assume that responses of 2 treatments are significant if they are 2 SE away from the mean response of each other. Confidence limits were determined at 2 standard deviations from the mean.
As expected for regression models, the predictions seem to work well for moderate response sites but seem to underestimate the most responsive sites such as Pampaida and overestimate responses in nonresponsive Kiberashi site. We expect that with further refinement such as transformations mentioned above, the present models will be more robust.
The proposed model allows the prediction of response for a given treatment when the control yields are known, and could be applied to generate quick farm level information on yield level that can be achieved from given a management.
Table 4. Site level predicted responses (growth rates) at mean log of control of 0.5 in different sentinel sites. Numbers in bracket indicate the lower and upper confidence limits. Values presented here are in log scale.
Kiberashi 
Koloko 
Kontela 
Mbinga 
Nkhata Bay 
Pampaida 
Sidindi 
Thuchila 

NK 
0.51 (0.34, 0.68) 
0.05 (0.21, 0.31) 
0.06 (0.33, 0.2) 
0.58 (0.43, 0.74) 
0.55 (0.4, 0.69) 
0.48a (0.32, 0.64) 
0.48 (0.33, 0.63) 
0.27 (0.12, 0.42) 
NP 
0.76 (0.6, 0.93) 
0.3 (0.03, 0.56) 
0.18 (0.09, 0.45) 
0.84 (0.68, 0.99) 
0.8 (0.65, 0.94) 
0.73 (0.57, 0.89) 
0.73 (0.59, 0.88) 
0.52 (0.37, 0.67) 
NPK 
0.82 (0.65, 0.99) 
0.26 (0, 0.52) 
0.14 (0.12, 0.41) 
0.89 (0.74, 1.05) 
0.86 (0.71, 1) 
0.79 (0.63, 0.94) 
0.79 (0.65, 0.94) 
0.58 (0.43, 0.72) 
NPK+Lime 
0.86 (0.69, 1.04) 
0.2 (0.06, 0.47) 
0.08 (0.18, 0.36) 
0.94 (0.77, 1.1) 
0.9 (0.74, 1.05) 
0.83 (0.66, 1) 
0.83 (0.68, 0.99) 
0.62 (0.46, 0.78) 
NPK+MN 
0.93 (0.77, 1.1) 
0.22 (0.04, 0.48) 
0.1 (0.17, 0.37) 
1.01 (0.85, 1.16) 
0.97 (0.82, 1.11) 
0.9 (0.74, 1.06) 
0.9 (0.76, 1.05) 
0.69 (0.54, 0.84) 
NPK+Manure 
1.08 (0.92, 1.25) 
0.32 (0.05, 0.58) 
0.2 (0.06, 0.47) 
1.16 (1, 1.31) 
1.12 (0.97, 1.27) 
1.05 (0.89, 1.21) 
1.05 (0.91, 1.2) 
0.84 (0.69, 0.99) 
PK 
0.22 (0.05, 0.39) 
0.01 (0.24, 0.28) 
0.1 (0.37, 0.17) 
0.29 (0.14, 0.45) 
0.26 (0.11, 0.4) 
0.19 (0.03, 0.34) 
0.19 (0.04, 0.34) 
0.01 (0.16, 0.13) 
NK 
0.51 (0.34, 0.68) 
0.05 (0.21, 0.31) 
0.06 (0.33, 0.2) 
0.58 (0.43, 0.74) 
0.55 (0.4, 0.69) 
0.48 (0.32, 0.64) 
0.48 (0.33, 0.63) 
0.27 (0.12, 0.42) 