[R] different forms of nls recommendations

Bernardo Rangel Tura tura at centroin.com.br
Sun Mar 21 22:49:08 CET 2010


On Sat, 2010-03-20 at 14:55 -0800, emorway wrote:
> Hello, 
> 
> Using this data:
> http://n4.nabble.com/file/n1676330/US_Final_Values.txt US_Final_Values.txt 
> 
> and the following code i got the image at the end of this message:
> 
> US.final.values<-read.table("c:/tmp/US_Final_Values.txt",header=T,sep=" ")
> US.nls.1<-nls(US.final.values$ECe~a*US.final.values$WTD^b+c,data=US.final.values,start=list(a=2.75,b=-0.95,c=0.731),trace=TRUE)
> f.US1<-function(x){coef(US.nls.1)["a"]*x^coef(US.nls.1)["b"]+coef(US.nls.1)["c"]}
> xvals.US1<-seq(min(US.final.values$WTD),max(US.final.values$WTD),length.out=75)
> yvals.US1<-f.US1(xvals.US1)
> Rsq.nls.1<-sum((predict(US.nls.1)-mean(US.final.values$ECe))^2/sum((US.final.values$ECe-mean(US.final.values$ECe))^2))
> plot(US.final.values$WTD,US.final.values$ECe,col="red",pch=19,cex=.75)
> lines(xvals.US1,yvals.US1,col="blue")
> 
> but the r^2 wasn't so hot.  
> Rsq.nls.1
> [1] 0.2377306
> 
> So I wanted to try a different equation of the general form a/(b+c*x^d)
> 
> US.nls.2<-nls(US.final.values$ECe~(a/(b+c*US.final.values$WTD^d)),data=US.final.values,start=list(a=100.81,b=73.7299,c=0.0565,d=-6.043),trace=TRUE,algorithm="port")
> 
> but that ended with "Convergence failure: false convergence (8)".  I tried


Hi emorway,

Do you have 657 obs and 4 parameters to fit.
In my opinion you have few obs...
I think do you  fit in steps:

US.nls.2<-nls(ECe~(a/(b + c *
WTD^d)),data=US.final.values,start=list(a=100.81,b=73.7299,c=0.0565,d=-6.043),trace=TRUE,algorithm="port")
temp_nls1 <- nls(ECe~(100/(73 + .05 *
WTD^d)),data=US.final.values,start=list(d=-6.043),trace=TRUE,algorithm="port")
temp_nls2 <- nls(ECe~(100/(73 + .05 *
WTD^d)),data=US.final.values,start=list(d=-1.01613),trace=TRUE,algorithm="port")
temp_nls3 <- nls(ECe~(100/(73 + c *
WTD^(-1.01613))),data=US.final.values,start=list(c=0.05),trace=TRUE,algorithm="port")
temp_nls4 <- nls(ECe~(100/(73 + c *
WTD^(-1.01613))),data=US.final.values,start=list(c=-14.7127),trace=TRUE,algorithm="port")
temp_nls5 <- nls(ECe~(100/(b-14.7127 *
WTD^(-1.01613))),data=US.final.values,start=list(b=73),trace=TRUE,algorithm="port")
temp_nls6 <- nls(ECe~(100/(b-14.7127 *
WTD^(-1.01613))),data=US.final.values,start=list(b=70.4936),trace=TRUE,algorithm="port")
temp_nls7 <- nls(ECe~(a/(70.4936-14.7127 *
WTD^(-1.01613))),data=US.final.values,start=list(a=100),trace=TRUE,algorithm="port")
  0:     2243.9898:  100.000
  1:     2122.8218:  106.219
  2:     1359.8819:  187.530
  3:     1359.8819:  187.530



-- 
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil



More information about the R-help mailing list