[R] nls & fitting
Gabor Grothendieck
ggrothendieck at gmail.com
Sun May 21 23:07:22 CEST 2006
You have nearly as many parameters as data points which may
cause fundamental problem singularity problems but one thing to
try, just in case, is to transform it to unconstrained. For example,
let
A1 = 0.1 + A1x^2
and then substitute A1 with the right hand side so that it becomes a function
of the parameter A1x instead of A1 and similarly for the other parameters
so that it becomes unconstrained.
On 5/21/06, Lorenzo Isella <lorenzo.isella at gmail.com> wrote:
> Dear All,
> I may look ridiculous, but I am puzzled at the behavior of the nls with
> a fitting I am currently dealing with.
> My data are:
>
> x N
> 1 346.4102 145.428256
> 2 447.2136 169.530634
> 3 570.0877 144.081627
> 4 721.1103 106.363316
> 5 894.4272 130.390552
> 6 1264.9111 36.727069
> 7 1788.8544 52.848587
> 8 2449.4897 25.128742
> 9 3464.1016 7.531766
> 10 4472.1360 8.827367
> 11 6123.7244 6.600603
> 12 8660.2540 4.083339
>
> I would like to fit N as a function of x according to a function
> depending on 9 parameters (A1,A2,A3,mu1,mu2,mu3,myvar1,myvar2,myvar3),
> namely
> N ~
> (log(10)*A1/sqrt(2*pi)/log(myvar1)*exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1))
>
> +log(10)*A2/sqrt(2*pi)/log(myvar2)*exp(-((log(x/mu2))^2)/2/log(myvar2)/log(myvar2))
>
> +log(10)*A3/sqrt(2*pi)/log(myvar3)*exp(-((log(x/mu3))^2)/2/log(myvar3)/log(myvar3)))
>
> (i.e. N is to be seen as a sum of three "bells" whose parameters I need
> to determine).
>
>
> So I tried:
> out<-nls(N ~
> (log(10)*A1/sqrt(2*pi)/log(myvar1)*exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1))
>
> +log(10)*A2/sqrt(2*pi)/log(myvar2)*exp(-((log(x/mu2))^2)/2/log(myvar2)/log(myvar2))
>
> +log(10)*A3/sqrt(2*pi)/log(myvar3)*exp(-((log(x/mu3))^2)/2/log(myvar3)/log(myvar3)))
> ,start=list(A1 = 85,
> A2=23,A3=4,mu1=430,mu2=1670,mu3=4900,myvar1=1.59,myvar2=1.5,myvar3=1.5 )
> ,algorithm = "port"
> ,control=list(maxiter=20000,tol=10000)
> ,lower=c(A1=0.1,A2=0.1,A3=0.1,mu1=0.1,mu2=0.1,mu3=0.1,myvar1=0.1,myvar2=0.1,myvar3=0.1)
> )
>
> getting the error message:
> Error in nls(N ~ (log(10) * A1/sqrt(2 * pi)/log(myvar1) *
> exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1)) + :
> Convergence failure: singular convergence (7)
>
>
> I tried to adjust tol & maxiter, but unsuccessfully.
> If I try fitting N with only two "bells", then nls works:
>
> out<-nls(N ~
> (log(10)*A1/sqrt(2*pi)/log(myvar1)*exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1))
>
> +log(10)*A2/sqrt(2*pi)/log(myvar2)*exp(-((log(x/mu2))^2)/2/log(myvar2)/log(myvar2))
> )
> ,start=list(A1 = 85, A2=23,mu1=430,mu2=1670,myvar1=1.59,myvar2=1.5 )
> ,algorithm = "port"
> ,control=list(maxiter=20000,tol=10000)
> ,lower=c(A1=0.1,A2=0.1,mu1=0.1,mu2=0.1,myvar1=0.1,myvar2=0.1)
> )
>
> out
> Nonlinear regression model
> model: N ~ (log(10) * A1/sqrt(2 * pi)/log(myvar1) *
> exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1)) + log(10) *
> A2/sqrt(2 * pi)/log(myvar2) *
> exp(-((log(x/mu2))^2)/2/log(myvar2)/log(myvar2)))
> data: parent.frame()
> A1 A2 mu1 mu2 myvar1 myvar2
> 84.920085 40.889968 409.656404 933.081936 1.811560 2.389215
> residual sum-of-squares: 2394.876
>
> Any idea about how to get nls working with the whole model?
> I had better luck with the nls.lm package, but it does not allow to
> introduce any constrain on my fitting parameters.
> I was also suggested to try other packages like optim to do the same
> fitting, but I am a bit unsure about how to set up the problem.
> Any suggestions? BTW, I am working with R Version 2.2.1
>
> Lorenzo
>
> ______________________________________________
> R-help at 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
>
More information about the R-help
mailing list