This is presumably obvious to most if not all R programmers, but I became aware today of a hugely (?) delaying tactic in my R codes. I was working with Jean-Michel and Natesh [who are visiting at the moment] and when coding an MCMC run I was telling them that I usually preferred to code Nsim=10000 as Nsim=10^3 for readability reasons. Suddenly, I became worried that this representation involved a computation, as opposed to Nsim=1e3 and ran a little experiment:
[sourcecode language=”r” gutter=”false”]
> system.time(for (t in 1:10^8) x=10^3)
utilisateur système écoulé
30.704 0.032 30.717
> system.time(for (t in 1:1e8) x=10^3)
utilisateur système écoulé
30.338 0.040 30.359
> system.time(for (t in 1:10^8) x=1000)
utilisateur système écoulé
6.548 0.084 6.631
> system.time(for (t in 1:1e8) x=1000)
utilisateur système écoulé
6.088 0.032 6.115
> system.time(for (t in 1:10^8) x=1e3)
utilisateur système écoulé
6.134 0.029 6.157
> system.time(for (t in 1:1e8) x=1e3)
utilisateur système écoulé
6.627 0.032 6.654
> system.time(for (t in 1:10^8) x=exp(3*log(10)))
utilisateur système écoulé
60.571 0.000 57.103
[/sourcecode]
So using the usual scientific notation with powers is taking its toll! While the calculator notation with e is cost free… Weird!
I understand that the R notation 10^6 is an abbreviation for a power function that can be equally applied to pi^pi, say, but still feel aggrieved that a nice scientific notation like 10⁶ ends up as a computing trap! I thus asked the question to the Stack Overflow forum, getting the (predictable) answer that the R code 10^6 meant calling the R power function, while 1e6 was a constant. Since 10⁶ does not differ from ππ, there is no reason 10⁶ should be recognised by R as a million. Except that it makes my coding more coherent.
[sourcecode language=”r” gutter=”false”]
> system.time( for (t in 1:10^8) x=pi^pi)
utilisateur système écoulé
44.518 0.000 43.179
> system.time( for (t in 1:10^8) x=10^6)
utilisateur système écoulé
38.336 0.000 37.860
[/sourcecode]
Another thing I discovered from this answer to my question is that negative integers are also requesting call to a function:
[sourcecode language=”r” gutter=”false”]
> system.time( for (t in 1:10^8) x=1)
utilisateur système écoulé
10.561 0.801 11.062
> system.time( for (t in 1:10^8) x=-1)
utilisateur système écoulé
22.711 0.860 23.098
[/sourcecode]
This sounds even weirder.
The use of Buffon’s needle to approximate π by a (slow) Monte Carlo estimate is a well-known Monte Carlo illustration. But that a similar experiment can be used for approximating e seems less known, if judging from the 08 January riddle from The Riddler. When considering a sequence of length n exchangeable random variables, the probability of a particuliar ordering of the sequence is 1/n!. Thus, counting how many darts need be thrown on a target until the distance to the centre increases produces a random number N≥2 with pmf 1/n!-1/(n+1)! and with expectation equal to e. Which can be checked as follows
