[R] nls problem

Troels Ring tr|ng @end|ng |rom gvdnet@dk
Fri Apr 3 06:43:47 CEST 2020


Thank you – yes indeed, on suggestion from Rui I noticed that the problem was mainly that I took the parameters from “m” instead of coef(m) i.e

a        b        d 

-1630.71   457.68   -32.11

Instead of

          a           b           d 
-1630.70993   457.67555   -32.11469 

 

And that alone saves me. 

BW
Troels

 

Fra: William Dunlap <wdunlap using tibco.com> 
Sendt: 2. april 2020 18:50
Til: Troels Ring <tring using gvdnet.dk>
Cc: r-help mailing list <r-help using r-project.org>
Emne: Re: [R] nls problem

 

Roundoff/cancelation error: compare the following.  The first is equivalent to your function, the last to fitted().  

 

> with(aedf, t(cbind(1, pH, pH^2) %*% round(coef(m), digits=2)))
            [,1]      [,2]       [,3]         [,4]       [,5]       [,6]       [,7]       [,8]       [,9]     [,10]
[1,] -0.01635965 0.1076024 0.07477337 -0.007166685 -0.1337865 -0.2851877 -0.4804638 -0.6865135 -0.9870137 -1.398027
> with(aedf, t(cbind(1, pH, pH^2) %*% round(coef(m), digits=4)))
           [,1]        [,2]       [,3]       [,4]       [,5]       [,6]       [,7]      [,8]      [,9]     [,10]
[1,] -0.2196286 -0.09864071 -0.1345213 -0.2180792 -0.3463237 -0.4991861 -0.6959856 -0.903394 -1.205598 -1.618608
> with(aedf, t(cbind(1, pH, pH^2) %*% round(coef(m), digits=16)))
           [,1]        [,2]       [,3]       [,4]       [,5]      [,6]      [,7]       [,8]      [,9]     [,10]
[1,] -0.2208705 -0.09989926 -0.1357969 -0.2193638 -0.3476174 -0.500488 -0.697296 -0.9047119 -1.206926 -1.619947

 

Note that your model is linear and could be fitted with

   lm(data=aedf, Flux ~ pH + I(pH^2))

   lm(data=aedf, Flux ~ poly(pH, 2))
The latter uses a more stable parameterization.




Bill Dunlap
TIBCO Software
wdunlap tibco.com <http://tibco.com> 

 

 

On Thu, Apr 2, 2020 at 4:15 AM Troels Ring <tring using gvdnet.dk <mailto:tring using gvdnet.dk> > wrote:

Dear friends - I'm on Win10 with R 6.3.1 and have a very simple problem with
nls which apparently gives a decent fit to the parable below, even without
starting values. But when I then think I know the meaning of the three
parameters a, b, and d it goes very wrong. I guess I am again overlooking
something easy but cannot spot it.

BW
Troels Ring,

Aalborg, Denmark





aedf <- structure(list(Flux = c(-0.141256, -0.154709, -0.215247, -0.302691, 

            -0.32287, -0.511211, -0.605381, -0.813901, -1.11659, -1.76906

), pH = c(7.06273, 7.11182, 7.16182, 7.18818, 7.21455, 7.23818, 

          7.26273, 7.28455, 7.31182, 7.34364)), class = "data.frame",
row.names = c(NA, 


-10L))

m <- with(aedf,nls(Flux~a + b*pH + d*pH^2))



with(aedf,plot(pH,Flux))

with(aedf,lines(pH,fitted(m),lty=1,col="red",lwd=3))



m

# a        b        d 

# -1630.70   457.67   -32.11 



fitted(m)

# 1] -0.22087053 -0.09989926 -0.13579690 -0.21936385 -0.34761742 -0.50048799

# [7] -0.69729602 -0.90471194 -1.20692552 -1.61994657



FPG <- function(pH) -1630.70 + 457.67*pH -32.11*pH^2



FPG(aedf$pH)

# [1] -0.016359649  0.107602395  0.074773375 -0.007166685 -0.133786467

# [6] -0.285187665 -0.480463769 -0.686513537 -0.987013685 -1.398026917



# So why aren't fitted(m) and FPG(aedf$pH) not closer ("equal")?




This email has been scanned by BullGuard antivirus protection.
For more info visit www.bullguard.com <http://www.bullguard.com> 
<http://www.bullguard.com/tracking.aspx?affiliate=bullguard <http://www.bullguard.com/tracking.aspx?affiliate=bullguard&buyaffiliate=smtp&url=/> &buyaffiliate=smt
p&url=/> 

        [[alternative HTML version deleted]]

______________________________________________
R-help using r-project.org <mailto:R-help using r-project.org>  mailing list -- To UNSUBSCRIBE and more, see
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.


This email has been scanned by BullGuard antivirus protection.
For more info visit www.bullguard.com <http://www.bullguard.com/tracking.aspx?affiliate=bullguard&buyaffiliate=smtp&url=/> 

	[[alternative HTML version deleted]]



More information about the R-help mailing list