2

I'd like to estimate $e$ using monte carlo methods. What is the monte carlo method that converges faster when $n \to \infty$?

Can you show some simulations using different $n$'s for this method?

Leo Ribeiro
  • 173
  • 1
  • 6
  • See https://stats.stackexchange.com/questions/193990/approximate-e-using-monte-carlo-simulation for multiple approaches. Does it answer your question? – Tim Apr 16 '18 at 09:44
  • @Tim It doesn't answer my question because some responses use $e$ itself to estimate it and others aren't monte carlo methods, at least, not unbiased, because the expected value is changing. I think that the Aksakal's answer is good, but I'd like other approaches. – Leo Ribeiro Apr 16 '18 at 15:08

1 Answers1

3

Given that there are some very efficient deterministic methods to calculate $e$ (up to a given number of decimal places), I would be surprised if there is any stochastic Monte-Carlo method that can hold a candle to these. Nevertheless, I'll try to get the ball rolling by giving one possible estimator. To be clear, I make absolutely no claim to efficiency here, but I'll give an estimator, and hopefully others will be able to offer better methods.

I will assume for the purposes of this question that you are able to generate and use a sequence of $n$ uniform pseudo-random variables $U_1, \cdots , U_n \sim \text{IID U}(0,1)$ and you then need to estimate $e$ by some method using basic arithmetic operations.$^\dagger$ The present method is motivated by a simple result involving uniform random variables:

$$\mathbb{E} \Bigg( \frac{\mathbb{I}(U_i \geqslant 1 / e) }{U_i} \Bigg) = \int \limits_{1/e}^1 \frac{du}{u} = 1.$$

Estimating $e$ using this result: We first order the sample values into descending order to obtain the order statistics $u_{(1)} \geqslant \cdots \geqslant u_{(n)}$ and then we define the partial sums:

$$S_n(k) \equiv \frac{1}{n} \sum_{i=1}^k \frac{1}{u_{(i)}} \quad \text{for all } k = 1, .., n.$$

Now, let $m \equiv \min \{ k | S(k) \geqslant 1 \}$ and then estimate $1/e$ by interpolation of the ordered uniform variables. This gives an estimator for $e$ given by:

$$\hat{e} \equiv \frac{2}{u_{(m)} + u_{(m+1)}}.$$

This method has some slight bias (owing to the linear interpolation of the cut-off point for $1/e$) but it is a consistent estimator for $e$. The method can be implemented fairly easily but it requires sorting of values, which is more computationally intensive than deterministic calculation of $e$. This method is slow, since it involves sorting of values.

Implementation in R: The method can be implemented in R using runif to generate uniform values. The code is as follows:

EST_EULER <- function(n) { U <- sort(runif(n), decreasing = TRUE);
                           S <- cumsum(1/U)/n;
                           m <- min(which(S >= 1));
                           2/(U[m-1]+U[m]); }

Implementing this code gives convergence to the true value of $e$, but it is very slow compared to deterministic methods.

set.seed(1234);

EST_EULER(10^3);
[1] 2.715426

EST_EULER(10^4);
[1] 2.678373

EST_EULER(10^5);
[1] 2.722868

EST_EULER(10^6); 
[1] 2.722207

EST_EULER(10^7);
[1] 2.718775

EST_EULER(10^8);
[1] 2.718434

> exp(1)
[1] 2.718282

$^\dagger$ Obviously we want to avoid any method that makes use of any transformation that involves an exponential or logarithm. The only way we could employ these would be if we already have $e$.

Ben
  • 91,027
  • 3
  • 150
  • 376