[R] Need suggestions for finding dose response using nls

From: Ken <kef_at_plantpath.wisc.edu>
Date: Sat 05 Mar 2005 - 07:21:10 EST


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

______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Sat Mar 05 07:25:31 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:30:40 EST