The Robust Randomized Quasi Monte Carlo method, applications to integrating singular functions

This repository contains the code of the paper The Robust Randomized Quasi Monte Carlo method, applications to integrating singular functions by Emmanuel Gobet, Matthieu Lerasle and David Métivier.

The paper contains multiple figures. In all, the code $M$ is the number of realizations used to show the estimators' density. Since we are interested in the tails of this distribution, we typically need $M$ very large, which makes the simulation quite heavy. Some simulations are lightweight and can run in less than an hour, they are in a Julia Pluto notebook. Other heavy simulations are given a script meant to be run in parallel on a Slurm or SGE cluster.

For all the Randomized Quasi Monte Carlo simulations, we use $M$ Nested Uniform Scramble of the Sobol' sequence¹.

Lightweight simulations

These simulations take less than an hour and can then be put into a Pluto Notebook. The best way to visualize the notebooks is to open the html version of the notebook in a browser. From there, you can "Run or Edit" with Julia (if installed on your machine) or with the free² Binder cloud service. The user is free to change parameters to explore beyond what is shown in the paper.

• Figure 1 (Figure 1.a and Figure 1.b) - Exactly reproducible in the self-contained Pluto notebook. It runs for approximately 10 minutes with the paper settings of $M = 10^7$. This notebook only showcases robust estimators without any Quasi Monte Carlo. The two following "versions" are exactly the same, the first one is just rendered as an html page.
  • The friendly version of Figure_1.jl in a web browser to see directly the resulting notebook.
  • Figure_1.jl
• Figure 5 - Exactly reproducible in the Pluto notebook. It runs for approximately 40 minutes with the paper settings of $M = 10^4$. This notebook shows the whole workflow to do multiple Robust Randomized Quasi Monte Carlo. Thanks to Pluto's capabilities, we only show the interesting and interactive pieces and hide³ the background functions. The two following "versions" are exactly the same, the first one is just rendered as an html page.
  • The friendly version of Figure_5.jl.
  • The script version of the notebook Figure_5.jl.

Heavy simulations

These simulations use high-dimensional Randomized Quasi Monte Carlo simulations with $M = 10^4$. In order to run, one must typically use HPC/cluster computing. Unfortunately, the authors were not able to find a simple and efficient way to do distributed computations on a cluster, i.e., getting the same result regardless of the number of cores used without affecting the performance. It is not as simple as fixing a seed in serial computation. Hence, here we provide the script that can be used to produce similar figures. Since $M\times n$ is quite high, the result should be very similar to the paper's figures.

• Figure 2
  • Figure 2.a
  • Figure 2.b
• Figure 3
  • Figure 3.a
  • Figure 3.b
• Figure 4
• Figure 6

Programming language and tooling used

Julia

The Julia programming language is used in this project. One advantage of Julia is that it is both readable and fast (compiled). There is no C/C++/Fortan wrapper. All the code is Julia (and the packages used too). Hence, it is very easy to take a look at the randomization methods used in the paper, see here for example for Nested Uniform Scrambling⁴.

Main packages used

• QuasiMonteCarlo.jl developed by multiple authors. The documentation can be found here. @dmetivie contributed to the Quasi Monte randomization method and docs of the package.
• RobustMeans developed by @dmetivie.

Pluto notebooks

Most readers are probably familiar with Jupyter notebooks, Pluto's notebooks have plenty of noteworthy differences, see the official website. Here are some:

• Works only for Julia.
• Reactive, i.e., if you change one parameter somewhere in the notebook, all other parameters depending on function, plots, text, etc. will be automatically reevaluated accordingly. This allows one to quickly explore the effects of changing parameters, methods, etc.
• They are reproducible by nature, i.e., all the exact information on versioning used to generate the result is contained in the notebook.
• They are Julia files, i.e., one can open them as a simple and readable Julia script or as the notebook version.

Clusters

Part of this work used the École Polytechnique IDCS mesocentre (Cholesky). The final part of the work has been done using the MESO@LR cluster.

Footnotes

1. The "true" Sobol' sequence where the first point starts at 0. Note that the Sobol.jl@v.1.4 package used in QuasiMonteCarlo.jl implements a truncated version. Instead, one can use the implementation in the (depreacted) RandomizedQuasiMonteCarlo.jl package or in QMCGenerators.jl. ↩

2. This cloud service is (currently) free without registration. Hence, in one click, it loads the notebook and runs it. Setup might take some time. ↩

3. Hidden cells can be shown again by clicking on the eye symbol next to them. All the code is visible in the .jl script. ↩

4. Obviously reading this piece of code without context is not very helpful. ↩