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)
______________________________________________
[email protected] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.