Hi,
I've noticed I get different results fitting a function to some data on
my laptop to when I do it on my computer at work.
Here's a code snippet of what I do:
##------------------------------------------------------------------
require(circular) ## for Bessel function I.0
## Data:
dd <- c(0.9975948929787, 0.9093316197395, 0.7838819026947,
0.9096108675003, 0.8901804089546, 0.2995955049992, 0.9461286067963,
0.8248071670532, 0.2442084848881, 0.2836948633194, 0.7353935241699,
0.5812761187553, 0.8705610632896, 0.8744471669197, 0.7490273118019,
0.9947383403778, 0.9154829382896, 0.8659985661507, 0.6448246836662,
0.8588128685951, 0.7347437739372, -0.1645197421312, 0.970999121666,
0.8038327097893, 0.9558997154236, 0.6846113204956, 0.6286814808846,
0.9201356172562, 0.9422197341919, 0.3470877110958, 0.4154576957226,
0.0721184238791, 0.14151956141, -0.6142936348915, -0.4688512086868,
0.6805665493011, 0.3594025671482, 0.8991097211838, 0.7656877636909,
0.9282909035683, 0.9454715847969, 0.9766132831573, 0.4316343963146,
0.62679708004, 0.2093886137009, 0.3937581181526, 0.4254160523415,
0.8684504628181, 0.3844584524632, 0.9578431844711, 0.956972181797,
0.4456568360329, 0.9793710708618, 0.5825698971748, 0.929228246212,
0.9211971759796, 0.9407976865768, 0.821156680584, 0.2048042863607,
0.6473184227943, 0.9456319212914, 0.7021154165268, 0.9761978387833,
0.1485801786184, 0.2195029109716, 0.5378785729408, 0.8304615020752,
0.8596342802048, 0.950027525425, 0.9102076888084, 0.5108731985092,
0.7200184464455, 0.3571084141731, 0.9765330553055, -0.143017962575,
0.8576183915138, 0.1283493340015, -0.3226098418236, 0.7031792402267,
0.8708582520485, 0.56754809618, 0.060470353812, 0.8015220761299,
0.7363410592079, 0.671902179718, 0.8082517385483, 0.9468197822571,
0.9729647636414, 0.7919752597809, 0.9539568424225, 0.4840737581253,
0.850653231144, 0.5909016132355, 0.8414449691772, 0.9699150323868)
xlims <- c(-1,1)
bw <- 0.05
b <- seq(xlims[1],xlims[2],by=bw) ; nb <- length(b)
h <- hist( dd, breaks=b, plot=FALSE)
FisherAvgdPdf <- function(theta,theta0,kappa){
A <- kappa/(2*sinh(kappa))
A * I.0( kappa*sin(theta)*sin(theta0) ) * exp(
kappa*cos(theta)*cos(theta0) )
}
nls(dens ~ FisherAvgdPdf(theta,theta0,kappa),
data = data.frame( theta=acos(h$mids), dens=h$density ),
start=c( theta0=0.5, kappa=4.0 ),
algorithm="port", lower=c(0,0), upper=c(acos(xlims[1]),500),
control=list(warnOnly=TRUE) )
##------------------------------------------------------------------
On one machine, nls converges, and on the other it doesn't. Any ideas
why, and which is "right"? I can't see anything in R News that
could
be relevant.
The different R versions (and computers) are:
> R.version
_
platform i686-pc-linux-gnu
arch i686
os linux-gnu
system i686, linux-gnu
status
major 2
minor 11.1
year 2010
month 05
day 31
svn rev 52157
language R
version.string R version 2.11.1 (2010-05-31)
and
> R.version
_
platform x86_64-pc-linux-gnu
arch x86_64
os linux-gnu
system x86_64, linux-gnu
status
major 2
minor 10.1
year 2009
month 12
day 14
svn rev 50720
language R
version.string R version 2.10.1 (2009-12-14)
Thanks,
Phil Bett
--
------------------------------------------------------------------------
Dr Philip Bett (Room 2.009) http://www.astro.uni-bonn.de/~pbett
Argelander-Institut f?r Astronomie Tel : +49 (0)228-73-5084
Auf dem H?gel 71, D-53121 Bonn Email: p.e.bett at dunelm.org.uk
Rheinische Friedrich-Wilhelms-Universit?t Bonn
Since one is a 32-bit and the other one a 64-bit, and therefore the compiler is also different, you can get different numerical results easily, even with identical versions of R (and most of us do not have outdated R installations around). I just tried your example on 3 different systems with R-2.13.0 and all told me "singular convergence"... Uwe Ligges On 18.06.2011 15:44, Philip Bett wrote:> Hi, > I've noticed I get different results fitting a function to some data on > my laptop to when I do it on my computer at work. > > Here's a code snippet of what I do: > > ##------------------------------------------------------------------ > require(circular) ## for Bessel function I.0 > > ## Data: > dd <- c(0.9975948929787, 0.9093316197395, 0.7838819026947, > 0.9096108675003, 0.8901804089546, 0.2995955049992, 0.9461286067963, > 0.8248071670532, 0.2442084848881, 0.2836948633194, 0.7353935241699, > 0.5812761187553, 0.8705610632896, 0.8744471669197, 0.7490273118019, > 0.9947383403778, 0.9154829382896, 0.8659985661507, 0.6448246836662, > 0.8588128685951, 0.7347437739372, -0.1645197421312, 0.970999121666, > 0.8038327097893, 0.9558997154236, 0.6846113204956, 0.6286814808846, > 0.9201356172562, 0.9422197341919, 0.3470877110958, 0.4154576957226, > 0.0721184238791, 0.14151956141, -0.6142936348915, -0.4688512086868, > 0.6805665493011, 0.3594025671482, 0.8991097211838, 0.7656877636909, > 0.9282909035683, 0.9454715847969, 0.9766132831573, 0.4316343963146, > 0.62679708004, 0.2093886137009, 0.3937581181526, 0.4254160523415, > 0.8684504628181, 0.3844584524632, 0.9578431844711, 0.956972181797, > 0.4456568360329, 0.9793710708618, 0.5825698971748, 0.929228246212, > 0.9211971759796, 0.9407976865768, 0.821156680584, 0.2048042863607, > 0.6473184227943, 0.9456319212914, 0.7021154165268, 0.9761978387833, > 0.1485801786184, 0.2195029109716, 0.5378785729408, 0.8304615020752, > 0.8596342802048, 0.950027525425, 0.9102076888084, 0.5108731985092, > 0.7200184464455, 0.3571084141731, 0.9765330553055, -0.143017962575, > 0.8576183915138, 0.1283493340015, -0.3226098418236, 0.7031792402267, > 0.8708582520485, 0.56754809618, 0.060470353812, 0.8015220761299, > 0.7363410592079, 0.671902179718, 0.8082517385483, 0.9468197822571, > 0.9729647636414, 0.7919752597809, 0.9539568424225, 0.4840737581253, > 0.850653231144, 0.5909016132355, 0.8414449691772, 0.9699150323868) > > xlims <- c(-1,1) > bw <- 0.05 > b <- seq(xlims[1],xlims[2],by=bw) ; nb <- length(b) > h <- hist( dd, breaks=b, plot=FALSE) > > FisherAvgdPdf <- function(theta,theta0,kappa){ > A <- kappa/(2*sinh(kappa)) > A * I.0( kappa*sin(theta)*sin(theta0) ) * exp( > kappa*cos(theta)*cos(theta0) ) > } > > > nls(dens ~ FisherAvgdPdf(theta,theta0,kappa), > data = data.frame( theta=acos(h$mids), dens=h$density ), > start=c( theta0=0.5, kappa=4.0 ), > algorithm="port", lower=c(0,0), upper=c(acos(xlims[1]),500), > control=list(warnOnly=TRUE) ) > > ##------------------------------------------------------------------ > > On one machine, nls converges, and on the other it doesn't. Any ideas > why, and which is "right"? I can't see anything in R News that could be > relevant. > > > The different R versions (and computers) are: > > R.version > _ > platform i686-pc-linux-gnu > arch i686 > os linux-gnu > system i686, linux-gnu > status > major 2 > minor 11.1 > year 2010 > month 05 > day 31 > svn rev 52157 > language R > version.string R version 2.11.1 (2010-05-31) > > > and > > > R.version > _ > platform x86_64-pc-linux-gnu > arch x86_64 > os linux-gnu > system x86_64, linux-gnu > status > major 2 > minor 10.1 > year 2009 > month 12 > day 14 > svn rev 50720 > language R > version.string R version 2.10.1 (2009-12-14) > > > > > Thanks, > > Phil Bett > > >
Prof. John C Nash
2011-Jun-20 17:55 UTC
[R] different results from nls in 2.10.1 and 2.11.1
Interesting! I get nice convergence in both 32 and 64 bit systems on 2.13.0. I agree the older versions are a bit of a distraction. The inconsistent behaviour on current R is a concern. Maybe Philip, Uwe, and I (and others who might be interested) should take this off line and see what is going on. I'll offer to act as tabulator of results and have set up a temporary wiki on http://macnash.telfer.uottawa.ca/PBettProb to allow tries to be saved for comparison. Needs a login to "write". Please use username=Ruser password=statistics Hopefully there won't be wikispam, or I'll have to close it off. JN On 06/20/2011 06:00 AM, r-help-request at r-project.org wrote:> Date: Mon, 20 Jun 2011 11:25:54 +0200 > From: Uwe Ligges <ligges at statistik.tu-dortmund.de> > To: p.e.bett at dunelm.org.uk > Cc: r-help at r-project.org, Philip Bett <pbett at astro.uni-bonn.de> > Subject: Re: [R] different results from nls in 2.10.1 and 2.11.1 > Message-ID: <4DFF1222.7040202 at statistik.tu-dortmund.de> > Content-Type: text/plain; charset=ISO-8859-1; format=flowed > > Since one is a 32-bit and the other one a 64-bit, and therefore the > compiler is also different, you can get different numerical results > easily, even with identical versions of R (and most of us do not have > outdated R installations around). > > I just tried your example on 3 different systems with R-2.13.0 and all > told me "singular convergence"... > > Uwe Ligges > >