yield.saemix | R Documentation |

The`yield.saemix`

contains data from winter wheat experiments.

yield.saemix

This data frame contains the following columns:

- site
the site number

- dose
dose of nitrogen fertiliser (kg/ha)

- yield
grain yield (kg/ha)

- soil.nitrogen
end-of-winter mineral soil nitrogen (NO3- plus NH4+) in the 0 to 90 cm layer was measured on each site/year (kg/ha)

The data in the `yield.saemix`

comes from 37 winter wheat experiments carried out between 1990 and 1996
on commercial farms near Paris, France. Each experiment was from a different site.
Two soil types were represented, a loam soil and a chalky soil. Common winter wheat varieties were used.
Each experiment consisted of five to eight different nitrogen fertiliser rates, for a total of 224 nitrogen treatments.
Nitrogen fertilizer was applied in two applications during the growing season. For each nitrogen treatment,
grain yield (adjusted to 150 g.kg-1 grain moisture content) was measured. In addition,
end-of-winter mineral soil nitrogen (NO3- plus NH4+) in the 0 to 90 cm layer was measured on each site-year
during February when the crops were tillering. Yield and end-of-winter mineral soil nitrogen measurements
were in the ranges 3.44-11.54 t.ha-1 , and 40-180 kg.ha-1 respectively.

Makowski, D., Wallach, D., and Meynard, J.-M (1999). Models of yield, grain protein, and residual mineral nitrogen responses to applied nitrogen for winter wheat. Agronomy Journal 91: 377-385.

data(yield.saemix) saemix.data<-saemixData(name.data=yield.saemix,header=TRUE,name.group=c("site"), name.predictors=c("dose"),name.response=c("yield"), name.covariates=c("soil.nitrogen"),units=list(x="kg/ha",y="t/ha",covariates=c("kg/ha"))) # Model: linear + plateau yield.LP<-function(psi,id,xidep) { x<-xidep[,1] ymax<-psi[id,1] xmax<-psi[id,2] slope<-psi[id,3] f<-ymax+slope*(x-xmax) #' cat(length(f)," ",length(ymax),"\n") f[x>xmax]<-ymax[x>xmax] return(f) } saemix.model<-saemixModel(model=yield.LP,description="Linear plus plateau model", psi0=matrix(c(8,100,0.2,0,0,0),ncol=3,byrow=TRUE,dimnames=list(NULL, c("Ymax","Xmax","slope"))),covariate.model=matrix(c(0,0,0),ncol=3,byrow=TRUE), transform.par=c(0,0,0),covariance.model=matrix(c(1,0,0,0,1,0,0,0,1),ncol=3, byrow=TRUE),error.model="constant") saemix.options<-list(algorithms=c(1,1,1),nb.chains=1,seed=666, save=FALSE,save.graphs=FALSE,displayProgress=FALSE) # Plotting the data plot(saemix.data,xlab="Fertiliser dose (kg/ha)", ylab="Wheat yield (t/ha)") saemix.fit<-saemix(saemix.model,saemix.data,saemix.options) # Comparing the likelihoods obtained by linearisation and importance sampling # to the likelihood obtained by Gaussian Quadrature saemix.fit<-llgq.saemix(saemix.fit) { cat("LL by Importance sampling, LL_IS=",saemix.fit["results"]["ll.is"],"\n") cat("LL by linearisation, LL_lin=",saemix.fit["results"]["ll.lin"],"\n") cat("LL by Gaussian Quadrature, LL_GQ=",saemix.fit["results"]["ll.gq"],"\n") } # Testing for an effect of covariate soil.nitrogen on Xmax saemix.model2<-saemixModel(model=yield.LP,description="Linear plus plateau model", psi0=matrix(c(8,100,0.2,0,0,0),ncol=3,byrow=TRUE,dimnames=list(NULL, c("Ymax","Xmax","slope"))),covariate.model=matrix(c(0,1,0),ncol=3,byrow=TRUE), transform.par=c(0,0,0),covariance.model=matrix(c(1,0,0,0,1,0,0,0,1),ncol=3, byrow=TRUE),error.model="constant") saemix.fit2<-saemix(saemix.model2,saemix.data,saemix.options) # BIC for the two models { cat("Model without covariate, BIC=",saemix.fit["results"]["bic.is"],"\n") cat("Model with covariate, BIC=",saemix.fit2["results"]["bic.is"],"\n") pval<-1-pchisq(-2*saemix.fit["results"]["ll.is"]+2*saemix.fit2["results"]["ll.is"],1) cat(" LRT: p=",pval,"\n") } #' @keywords datasets

Embedding an R snippet on your website

Add the following code to your website.

For more information on customizing the embed code, read Embedding Snippets.