[R] Need suggestions for finding dose response using nls

Ken kef at plantpath.wisc.edu
Fri Mar 4 21:21:10 CET 2005


I am relatively new to R and am looking for advice, ideas or both...

I have a data set that consists of pathogen population sizes on 
individual plant units in an experimental field plot.  However, in 
order to estimate the pathogen population sizes I had to destroy the 
plant unit and could not determine if that plant unit became diseased 
or to what extent it would have become diseased. I collected disease 
data from each plot and have estimates of the disease in an 
experimental plot several days after the pathogen population sizes were 
estimated.  From this data I would like to determine a dose-response 
relationship in the field between the pathogen and host.

i. e.

P(disease) = sum P(disease | pathogen population size) x P(pathogen 
population size)  ... I think


I determined that the pathogen population sizes among plant units are 
well described by the lognormal distribution.  I am also assuming that 
P(disease | pathogen population size) is could be described by the 
normal distribution described by two parameters lam and tau...my dose 
response.

I would like to use nonlinear least squares to estimate lam and tau.  I 
included the R code and some data below.  I'm not sure if this code is 
correct, but it give reasonable estimates of lam and tau for the first 
data set.  However, I included another data set where the estimates of 
lam and tau do not make sense.  I have never done nonlinear regression. 
  Is the code correct?  Is this a sound approach to estimating a dose 
response in the field?  Does anyone have other ideas about how to 
approach this data?  I am currently  doing this analysis ad hoc, but if 
the results are promise I might like to do some planned experiments.

R 2.0.1
mac os 10.2

Thanks in advance for any responses
Ken

Example data...

mn is the mean of the log transformed pathogen population size
ss is the standard deviation of the log transformed pathogen population 
size
dap27 is the disease 7 days after the pathogen population sizes were 
estimated
dap31 is the disease 11 days after the pathogen population sizes were 
estimated

mn	ss	dap27	dap31
6.762	0.492	0.582	0.567
4.382	1.63	0.393	0.394
2.472	2.449	0.181	0.19
6.66	0.495	0.698	0.541
4.282	1.401	0.345	0.435
1.08	2.423	0.16	0.196
6.636	0.362	0.704	0.526
3.638	1.445	0.12	0.509
1.854	2.07	0.148	0.075


fm1 <- nls(dap27 ~ pnorm((mn - lam) / ((ss^2) + (tau^2))) , data = dg, 
start=c(lam = 6, tau = 2), trace = TRUE)

summary(fm1)

plot(dap27 ~ mn, dg)

kfc<-as.numeric(kf<-lm(dg$ss ~ dg$mn)$coefficients)
mm<-mean(dg$mn, na.rm = T)
vv<-var(dg$mn, na.rm = T)

for(i in 1:5){
n<-100
mn<-rnorm(n, mm, sqrt(vv))
xx<- mn * kfc[2] + kfc[1]; lter<-rnorm(n, 0, 1)
ss<-xx + lter
newdat<-data.frame(mn,ss)
lines(lowess(newdat$mn , predict(fm1, newdat), f = 2/3), col = i)
}

The nonlinear regression did not fit for this data set (or at least 
dap27).

mn	ss	dap27	dap31
0.31	0.94	0.01	0.31
0.09	3.33	0.01	.
3.3	2.43	0.07	0.22
5.9	1.78	0.46	0.91
6.03	1.26	0.33	0.96
4.7	1.96	0.09	0.67
2.67	2.04	0.06	0.51
4.7	1.78	0.31	0.93
3.07	2.74	0.02	0.39
5.13	2.6	0.27	0.68
3.38	2.24	0.05	0.44
5.9	1.7	0.49	.
5.06	2.18	0.08	0.54
2.99	2.64	0.09	0.38
5.41	2.73	0.08	0.46
1.65	2.96	0.04	0.39




More information about the R-help mailing list