[R] Interactions in a nls model

Anthony Lawrence Nguy-Robertson anthony.robertson at huskers.unl.edu
Tue Feb 8 22:12:20 CET 2011


I am interested in testing two similar nls models to determine if the 
lines are statistically different when fitted with two different data 
sets; one corn, another soybean. I know I can do this in linear models 
by testing for interactions. See Introductory Statistics with R by 
Dallgaard p212-218 for an example. I have two different data sets I am 
comparing to lai. ci.re should have very little difference between corn 
and soybean, while ci.gr, there should be a difference. If I use the 
simplistic form described in Dallgaard (see example in code below) it 
does not work correctly. What do I need to do to test for this 
interaction? Thank you!

My data is located here: ftp://snrftp.unl.edu/Outgoing/example/

Here is my example code:

###load data###
data <- read.csv("example.csv", header=TRUE)
eq <- function(x,a,b,c) {a+b*exp(x^c)}

##create non-linear models##
nls.ci.re <- nls(lai~eq(ci.re,a,b,c),data=data, start=c(a=1,b=1,c=1))
nls.ci.gr <- nls(lai~eq(ci.gr,a,b,c),data=data, start=c(a=1,b=1,c=1))

#create non-linear models for corn#
nls.ci.re.corn <- 
nls(lai~eq(ci.re,a,b,c),data=subset(data,crop=="corn"), 
start=c(a=1,b=1,c=1))
nls.ci.gr.corn <- 
nls(lai~eq(ci.gr,a,b,c),data=subset(data,crop=="corn"), 
start=c(a=1,b=1,c=1))

#create non-linear models for soybean#
nls.ci.re.soybean <- 
nls(lai~eq(ci.re,a,b,c),data=subset(data,crop=="soybean"), 
start=c(a=1,b=1,c=1))
nls.ci.gr.soybean <- 
nls(lai~eq(ci.gr,a,b,c),data=subset(data,crop=="soybean"), 
start=c(a=1,b=1,c=1))

#test interactions according to Introductory Statistics with R by 
Dalgaard p. 213#
corn.soybean.interactions.ci.re.a <-
     abs((summary(nls.ci.re.corn)$coeff[1,1]-
     summary(nls.ci.re.soybean)$coeff[1,1])/
     sqrt(summary(nls.ci.re.corn)$coeff[1,2]^2+
     summary(nls.ci.re.soybean)$coeff[1,2]^2))

corn.soybean.interactions.ci.gr.a <-
     abs((summary(nls.ci.gr.corn)$coeff[1,1]-
     summary(nls.ci.gr.soybean)$coeff[1,1])/
     sqrt(summary(nls.ci.gr.corn)$coeff[1,2]^2+
     summary(nls.ci.gr.soybean)$coeff[1,2]^2))

corn.soybean.interactions.ci.re.a.p.value <- 
pt(corn.soybean.interactions.ci.re.a,df=summary(nls.ci.re.corn)$df[2],lower.tail=FALSE)
corn.soybean.interactions.ci.gr.a.p.value <- 
pt(corn.soybean.interactions.ci.gr.a,df=summary(nls.ci.gr.corn)$df[2],lower.tail=FALSE)



More information about the R-help mailing list