Soybean: Growth of soybean plants

SoybeanR Documentation

Growth of soybean plants

Description

The Soybean data frame has 412 rows and 5 columns.

Format

This data frame contains the following columns:

Plot

a factor giving a unique identifier for each plot.

Variety

a factor indicating the variety; Forrest (F) or Plant Introduction #416937 (P).

Year

a factor indicating the year the plot was planted.

Time

a numeric vector giving the time the sample was taken (days after planting).

weight

a numeric vector giving the average leaf weight per plant (g).

Details

These data are described in Davidian and Giltinan (1995, 1.1.3, p.7) as “Data from an experiment to compare growth patterns of two genotypes of soybeans: Plant Introduction #416937 (P), an experimental strain, and Forrest (F), a commercial variety.” In order to fit the Nonlinear data we sugest to use the three parameter logistic model as in Pinheiro & Bates (1995).

Source

Pinheiro, J. C. and Bates, D. M. (2000), Mixed-Effects Models in S and S-PLUS, Springer, New York. (Appendix A.27)

Davidian, M. and Giltinan, D. M. (1995), Nonlinear Models for Repeated Measurement Data, Chapman and Hall, London.

Examples

## Not run: 
data(Soybean)
attach(Soybean)

#################################
#A full model (no covariate)

y     = weight          #response
x     = Time            #time

#Expression for the three parameter logistic curve

exprNL = expression((fixed[1]+random[1])/(1 + exp(((fixed[2]+random[2])- x)/(fixed[3]+random[3]))))

#Initial values for fixed effects
initial = c(max(y),0.6*max(y),0.73*max(y))

#A median regression (by default)
median_reg = QRNLMM(y,x,Plot,initial,exprNL)

#Assing the fit

fxd     = median_reg$res$beta
nlmodel = median_reg$res$nlmodel
seqc    = seq(min(x),max(x),length.out = 500)

group.plot(x = Time,y = weight,groups = Plot,type="l",
             main="Soybean profiles",xlab="time (days)",
             ylab="mean leaf weight (gr)",col="gray")
             
lines(seqc,nlmodel(x = seqc,fixed = fxd,random = rep(0,3)),
      lwd=2,col="blue")             

#Histogram for residuals
hist(median_reg$res$residuals)

#########################################
#A model for comparing the two genotypes (with covariates)

y     = weight          #response
x     = Time            #time
covar = c(Variety)-1    #factor genotype (0=Forrest, 1=Plan Introduction)

#Expression for the three parameter logistic curve with a covariate

exprNL = expression((fixed[1]+(fixed[4]*covar[1])+random[1])/
                    (1 + exp(((fixed[2]+random[2])- x)/(fixed[3]+random[3]))))

#Initial values for fixed effects
initial = c(max(y),0.6*max(y),0.73*max(y),3)

# A quantile regression for the three quartiles
box_reg = QRNLMM(y,x,Plot,initial,exprNL,covar,p=c(0.25,0.50,0.75))

#Assing the fit for the median (second quartile)

fxd     = box_reg[[2]]$res$beta
nlmodel = box_reg[[2]]$res$nlmodel
seqc    = seq(min(x),max(x),length.out = 500)

group.plot(x = Time[Variety=="P"],y = weight[Variety=="P"],
             groups = Plot[Variety=="P"],type="l",col="light blue",
             main="Soybean profiles by genotype",xlab="time (days)",
             ylab="mean leaf weight (gr)")
             
group.lines(x = Time[Variety=="F"],y = weight[Variety=="F"],
              groups = Plot[Variety=="F"],col="gray")
             
lines(seqc,nlmodel(x = seqc,fixed = fxd,random = rep(0,3),covar=1),
      lwd=2,col="blue")

lines(seqc,nlmodel(x = seqc,fixed = fxd,random = rep(0,3),covar=0),
      lwd=2,col="black")


## End(Not run)

qrNLMM documentation built on Aug. 18, 2022, 5:05 p.m.