Re: [R] different results from nls in 2.10.1 and 2.11.1

2011-06-20 Thread Uwe Ligges 
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





__
R-help@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.

Re: [R] different results from nls in 2.10.1 and 2.11.1

2011-06-20 Thread Mike Marchywka 

 Date: Mon, 20 Jun 2011 11:25:54 +0200
 From: lig...@statistik.tu-dortmund.de
 To: p.e.b...@dunelm.org.uk
 CC: r-help@r-project.org; pb...@astro.uni-bonn.de
 Subject: Re: [R] different results from nls in 2.10.1 and 2.11.1
 
 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.

I guess the best all around approach is to dump the results from your
FisherAvgPdf function and get some idea what trajectory it takes
in the different cases. This is presumably not just an issue with R
and could have something to tell you about your data vs the model.
Even if it converged without incident, you'd probably want to 
look at more details.  You apparently are trying to fit a histogram to
a pdf and it is not too hard to just plot fit over histo and
spot check parameter space. 
 
Consider for example, 
 
plot(h$mids,h$density,type=b)
points(h$mids,FisherAvgdPdf(acos(h$mids),3.527e-8,3.198),pch=20)
points(h$mids,FisherAvgdPdf(acos(h$mids),3.527e-8,3.298),pch=13)
points(h$mids,FisherAvgdPdf(acos(h$mids),4.527e-8,3.198),pch=14)

and the find various measures for how good/bad these are etc. 
 
e-h$density-FisherAvgdPdf(acos(h$mids),3.527e-8,3.198)
sum(e*e)
e-h$density-FisherAvgdPdf(acos(h$mids),3.527e-8,3.298)
sum(e*e)

the error surface as function of parameters seems smooth etc,
 
foo - function (x,y) { e-h$density-FisherAvgdPdf(acos(h$mids),x,y) ; sum(e*e) 
}

x=(2+(1:100)*.02)/1e8;
y=3+((1:100)/100);


df-data.frame(x=0,y=0,z=0);
for ( i in x ) { for ( j in y ) { gg-data.frame(i,j,foo(i,j));str(gg); 
colnames(gg)-colnames(df); df-rbind(df,gg); }}
sel=which(df$x!=0)
scatterplot3d(df$x[sel],df$y[sel],df$z[sel])

 
 
 
 
  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

Re: [R] different results from nls in 2.10.1 and 2.11.1

2011-06-20 Thread Prof. John C Nash 
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-requ...@r-project.org wrote:
 Date: Mon, 20 Jun 2011 11:25:54 +0200
 From: Uwe Ligges lig...@statistik.tu-dortmund.de
 To: p.e.b...@dunelm.org.uk
 Cc: r-help@r-project.org, Philip Bett pb...@astro.uni-bonn.de
 Subject: Re: [R] different results from nls in 2.10.1 and 2.11.1
 Message-ID: 4dff1222.7040...@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
 


__
R-help@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.

[R] different results from nls in 2.10.1 and 2.11.1

2011-06-18 Thread Philip Bett 

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
day31
svn rev52157
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
day14
svn rev50720
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.b...@dunelm.org.uk
Rheinische Friedrich-Wilhelms-Universität Bonn

__
R-help@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.