[R] fitting "power model" in nls()

Rolf Turner r.turner at auckland.ac.nz
Sun Dec 2 20:50:26 CET 2007


Rule number 1:  Read the help for the function you are using.
You must supply starting values for the fit --- which the code you gave
doesn't do.

Rule number 2:  Don't use nls()!  Endless grief results.

Instead try:

foo <- function(par,x,y){
  Const <- par[1]
  B <- par[2]
  A <- par[3]
  sum((y-(Const+B*(x^A)))^2)
  }

fit <- optim(c 
(1,2.5,0.2),foo,x=area,y=richness,method="BFGS",control=list 
(maxit=1000))

plot(area,richness)
ccc <- fit$par
curve(ccc[1]+ccc[2]*(x^ccc[3]),from=0,to=2000,add=TRUE,col="red") #  
Fit doesn't look too bad.

HTH.

		cheers,

			Rolf Turner

On 3/12/2007, at 8:08 AM, Milton Cezar Ribeiro wrote:

> Dear all,
> I am still fighting against my "power model".
> I tryed several times to use nls() but I can´t run it.
> I am sending my variables and also the model which I would like to  
> fit.
> As you can see, this "power model" is not the best model to be fit,  
> but I really need also to fit it.
>
> The model which I would like to fit is Richness = B*(Area^A)
>
> richness<-c 
> (44,36,31,39,38,26,37,33,34,48,25,22,44,5,9,13,17,15,21,10,16,22,13,20 
> ,9,15,14,21,23,23,32,29,20,
> 26,31,4,20,25,24,32,23,33,34,23,28,30,10,29,40,10,8,12,13,14,56,47,44, 
> 37,27,17,32,31,26,23,31,34,
> 37,32,26,37,28,38,35,27,34,35,32,27,22,23,13,28,13,22,45,33,46,37,21,2 
> 8,38,21,18,21,18,24,18,23,22,
> 38,40,52,31,38,15,21)
> area<-c 
> (26.22,20.45,128.68,117.24,19.61,295.21,31.83,30.36,13.57,60.47,205.30 
> ,40.21,
> 7.99,1.18,5.40,13.37,4.51,36.61,7.56,10.30,7.29,9.54,6.93,12.60,
> 2.43,18.89,15.03,14.49,28.46,36.03,38.52,45.16,58.27,67.13,92.33,1.17,
> 29.52,84.38,87.57,109.08,72.28,66.15,142.27,76.41,105.76,73.47,1.71,30 
> 5.75,
> 325.78,3.71,6.48,19.26,3.69,6.27,1689.67,95.23,13.47,8.60,96.00,436.97 
> ,
> 472.78,441.01,467.24,1169.11,1309.10,1905.16,135.92,438.25,526.68,88.8 
> 8,31.43,21.22,
> 640.88,14.09,28.91,103.38,178.99,120.76,161.15,137.38,158.31,179.36,21 
> 4.36,187.05,
> 140.92,258.42,85.86,47.70,44.09,18.04,127.84,1694.32,34.27,75.19,54.39 
> ,79.88,
> 63.84,82.24,88.23,202.66,148.93,641.76,20.45,145.31,27.52,30.70)
> plot(richness~area)
>
> I tryed to fit the following model:
>
> m1<-nls(richness ~ Const+B*(area^A))
>
> Thanks a lot,
> miltinho
> Brazil.
>
>
>
>  para armazenamento!
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org 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.

######################################################################
Attention: 
This e-mail message is privileged and confidential. If you are not the 
intended recipient please delete the message and notify the sender. 
Any views or opinions presented are solely those of the author.

This e-mail has been scanned and cleared by MailMarshal 
www.marshalsoftware.com
######################################################################



More information about the R-help mailing list