Dear R community,
I have currently some problems with non linear regression analysis in R.
My data correspond to the degradation kinetic of a pollutant in two
different soil A and B, x data are time in day and y data are pollutant
concentration in soil.
In a first time, I want to fit the data for the soil A by using the Ct =
C0*exp(-k*Tpst) with Ct the concentration of pollutant at time t, C0 is the
initial concentration (i.e. t=0), Tps is the time in days.
The table containing data is called «tabMika»
Here, you can find my R program:
> tabMika<-read.delim("RMika.txt")
> tabMika
Tps SolA Solb
1 0 32.97 35.92
2 0 32.01 31.35
3 1 21.73 22.03
4 1 23.73 18.53
5 2 19.68 18.28
6 2 18.56 16.79
# and continue like that until 29 days.
> library(nls)
> attach(tabMika)
> Expon<-function(Tps,parm){
+ C0<-parm[1]
+ k<-parm[2]
+ }
> DegSA.nls<-nls(SolA~C0*exp(-k*Tps),start=c(C0=35, k=1),tabMika)
> summary(DegSA.nls)
Formula: SolA ~ C0 * exp(-k * Tps)
Parameters:
Estimate Std. Error t value Pr(>|t|)
C0 25.682104 1.092113 23.52 < 2e-16 ***
k 0.087356 0.007582 11.52 6.36e-13 ***
---
Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` '
1
Residual standard error: 2.584 on 32 degrees of freedom
Correlation of Parameter Estimates:
C0
k 0.7101
# until this step, I think (or may be hope) that it is OK.
#Then, I would like to use bootstrap function to obtain 10 estimations
(more of course but 10 is retained for this example) of the parameters C0
and t. I use this way:
> library(bootstrap)
> theta<-function(tabMika){coef(eval(DegSA.nls$call))}
> bootSolA.nls<-bootstrap(tabMika,10,theta)
Warning message:
multi-argument returns are deprecated in: return(thetastar, func.thetastar,
jack.boot.val, jack.boot.se,
> bootSolA.nls
$thetastar
[,1] [,2] [,3] [,4] [,5] [,6]
C0 25.68210358 25.68210358 25.68210358 25.68210358 25.68210358 25.68210358
k 0.08735615 0.08735615 0.08735615 0.08735615 0.08735615 0.08735615
[,7] [,8] [,9] [,10]
C0 25.68210358 25.68210358 25.68210358 25.68210358
k 0.08735615 0.08735615 0.08735615 0.08735615
$func.thetastar
NULL
$jack.boot.val
NULL
$jack.boot.se
NULL
$call
bootstrap(x = tabMika, nboot = 10, theta = theta)
# as you can notify, the 10 estimations of C0 and k are the same value
(????!!!!!!), so, my program is wrong but I cannot find where is the problem.
# In a second time, I would like to fit the same data with a new function
Ct = C0*Tpsb. I did as follow and obtained the following results:
> Pat<-function(Tps,parm){
+ a<-parm[1]
+ b<-parm[2]
+ }
> PatSA.nls<-nls(SolA~a*(Tps^b),start=c(a=35, b=-1),tabMika)
Error in numericDeriv(form[[3]], names(ind), env) : Missing value or an
Infinity produced when evaluating the model
> summary(PatSA.nls)
Error in summary(PatSA.nls) : Object "PatSA.nls" not found
# I cannot understand the mean of «Error in numericDeriv(form[[3]],
names(ind), env) : Missing value or an Infinity produced when evaluating
the model», please, can you help me.
Thank you very much in advance.
Michael C
Michaël COEURDASSIER, PhD
Department of Environmental Biology
UC INRA EA3184MRT
Institute for Environmental Sciences and Technology
University of Franche-Comte
Place Leclerc
25030 Besançon cedex
FRANCE
Tel : +33 (0)381 665 741
Fax : +33 (0)381 665 797
E-m@il: michael.coeurdassier@univ-fcomte.fr
http://lbe.univ-fcomte.fr/
[[alternative HTML version deleted]]
