& Center for Interdisciplinary Scientiﬁc Computation Illinois Institute of Technology [email protected] mypages.iit.edu/~hickernell with Sou-Cheng Choi, Mike McCourt, Jagadeeswaran Rathinavel, Aleksei Sorokin and the rest of the GAIL an QMCPy teams partially supported by SigOpt and NSF-DMS-1522687 Thanks to the organizers, especially during these unusual times Slides at speakerdeck.com/fjhickernell/quasi-monte-carlo-software Google Colaboratory notebook at tinyurl.com/QMCPyTutorial Blog at qmcpy.wordpress.com/ MCQMC 2020, August 2020
about QMC. How do you try it? QMC will give you 100 times the accuracy in the same amount of time as simple MC Just replace your IID random points with low discrepancy (LD) points 2/12
about QMC. How do you try it? QMC will give you 100 times the accuracy in the same amount of time as simple MC Often Just replace your IID random points with low discrepancy (LD) points Sometimes 2/12
about QMC. How do you try it? QMC will give you 100 times the accuracy in the same amount of time as simple MC Often Just replace your IID random points with low discrepancy (LD) points Sometimes Where can I get accurate, eﬃcient, easy to use QMC software to try for my problem? How can I make my great QMC software available for others? 2/12
about QMC. How do you try it? QMC will give you 100 times the accuracy in the same amount of time as simple MC Often Just replace your IID random points with low discrepancy (LD) points Sometimes Where can I get accurate, eﬃcient, easy to use QMC software to try for my problem? How can I make my great QMC software available for others? Let’s try to help 2/12
uses low discrepancy (LD) sequences Ti are random Xi may be deterministic or random T1, T2 · · · IID ∼ F X1, X2 · · · LD ∼ F Ti do not know about one another {Xi}n i=1 represent F well Fn(t1, . . . , tn) = F(t1) · · · F(tn) F{Xi}n i=1 (x) ≈ F(x) 3/12
software? (apologies to those I missed) LD Sequence Generators Multi-Level, Stopping Criteria, Applications MATLAB Statistics Toolbox— Sobol’ and Halton Mike Giles—Multi-Level (Quasi-)Monte Carlo in C++, MATLAB, Python, and R Marius Hofert & Christiane Lemieux—qrng R package, Sobol’ and Halton Guaranteed Automatic Integration Library (GAIL)—Stopping criteria in MATLAB BRODA, Sergei Kucherenko—Sobol’ in C, MATLAB, and Excel UQLab—Framework for Uncertainty Quantiﬁcation in MATLAB Art Owen—randomized Halton in R OpenTURNS—An Open source initiative for the Treatment of Uncertainties, Risks ’N Statistics in Python Pieterjan Robbe—LD sequences in Julia PyTorch—Sobol’ in Python Pierre L’Ecuyer—Lattice Builder and Stochastic Simulation in C/C++ and Java Dirk Nuyens—Magic Point Shop and QMC4PDE in MATLAB, Python, and C++ John Burkhardt—variety in C++, Fortran, MATLAB, & Python QMCPy—Python package incorporating and connecting the work of diﬀerent groups 4/12
software? (apologies to those I missed) LD Sequence Generators Multi-Level, Stopping Criteria, Applications Mike Giles—Multi-Level (Quasi-)Monte Carlo Marius Hofert & Christiane Lemieux—qrng Sobol’ and Halton Guaranteed Automatic Integration Library (GAIL)—Stopping criteria Art Owen—randomized Halton PyTorch—Sobol’ Pierre L’Ecuyer—Lattice Builder Dirk Nuyens—Magic Point Shop QMCPy—Python package incorporating and connecting the work of diﬀerent groups 4/12
better than IID Integration/Expectation Arising in ﬁnance, uncertainty quantiﬁcation, Bayesian inference, ... µ := E[f(X)] = X f(x) (x) dx ≈ 1 n n i=1 f(Xi) =: ^ µn LD points give faster convergence than IID Design of Computer Experiments LD can be more space ﬁlling (even than Latin hypercube sampling) for use in constructing surrogate models and uncertainty quantiﬁcation Global Optimization LD points are more space ﬁlling and ﬁnd good starting points for local methods 5/12
do we need? µ := E[f(X)] expectation = X f(x) (x) dx integration ≈ 1 n n i=1 f(Xi) =: ^ µn LD Generator producing {X1, X2, . . . } that mimics the distribution with PDF , e.g., uniform 6/12
do we need? µ := T g(t) λ(t) dt = · · · = E[f(X)] expectation = X f(x) (x) dx integration ≈ 1 n n i=1 f(Xi) =: ^ µn LD Generator producing {X1, X2, . . . } that mimics the distribution with PDF , e.g., uniform True Measure that deﬁnes the original integral, e.g., Lebesgue 6/12
do we need? µ := T g(t) λ(t) dt = · · · = E[f(X)] expectation = X f(x) (x) dx integration ≈ 1 n n i=1 f(Xi) =: ^ µn LD Generator producing {X1, X2, . . . } that mimics the distribution with PDF , e.g., uniform True Measure that deﬁnes the original integral, e.g., Lebesgue Integrand g, which deﬁnes the original integral, plus the transformed version, f, to ﬁt the LD generator 6/12
do we need? µ := T g(t) λ(t) dt = · · · = E[f(X)] expectation = X f(x) (x) dx integration ≈ 1 n n i=1 f(Xi) =: ^ µn LD Generator producing {X1, X2, . . . } that mimics the distribution with PDF , e.g., uniform True Measure that deﬁnes the original integral, e.g., Lebesgue Integrand g, which deﬁnes the original integral, plus the transformed version, f, to ﬁt the LD generator Stopping Criterion that determines how large n should be to ensure that |µ − ^ µn| ε 6/12
Jupyter Notebooks on Google Colaboratory Prerequisites—No Python, Jupyter, or QMC knowledge assumed Goals Show you how QMC software works Interest you in using/contributing Directions Point your browser to https://tinyurl.com/QMCPyTutorial Open the ﬁle in Google Colaboratory (may need to push a button at the top of your browser) Make a copy of this ﬁle onto your own Google drive account File → Save a copy in Drive Pause for questions 7/12
Coded almost entirely by Aleksei Sorokin (BS/MS ’21) Relies on code from several groups (see above) Documentation and testing by Sou-Cheng Choi and Jagadeeswaran R. Funded and encouraged by Mike McCourt at SigOpt 8/12
community software Our research groups are typically expert at only part of the whole picture: LD sequence generators Increasing eﬃciency, e.g., MLMC, MDM Stopping criteria Use cases A community library lets us all take advantage of the best 10/12
community software Our research groups are typically expert at only part of the whole picture: LD sequence generators Increasing eﬃciency, e.g., MLMC, MDM Stopping criteria Use cases A community library lets us all take advantage of the best Provides a consistent interface for pieces from diﬀerent places 10/12
community software Our research groups are typically expert at only part of the whole picture: LD sequence generators Increasing eﬃciency, e.g., MLMC, MDM Stopping criteria Use cases A community library lets us all take advantage of the best Provides a consistent interface for pieces from diﬀerent places Enables reproducible computational research 10/12
community software Our research groups are typically expert at only part of the whole picture: LD sequence generators Increasing eﬃciency, e.g., MLMC, MDM Stopping criteria Use cases A community library lets us all take advantage of the best Provides a consistent interface for pieces from diﬀerent places Enables reproducible computational research Tedious stuﬀ only done once 10/12
community software Our research groups are typically expert at only part of the whole picture: LD sequence generators Increasing eﬃciency, e.g., MLMC, MDM Stopping criteria Use cases A community library lets us all take advantage of the best Provides a consistent interface for pieces from diﬀerent places Enables reproducible computational research Tedious stuﬀ only done once Many eyes ﬁnd and correct errors MATLAB’s Sobol’ generator’s scrambling corrected in MATLAB 2017a after Tony Jiménez Rugama noticed the problem PyTorch’s Sobol’ generator found to be wrong unless double precision is proactively speciﬁed; also missing the ﬁrst point; reported at https://github.com/pytorch/pytorch/issues/32047 10/12
contribute Try out QMCPy and then Easy Submit your bugs and feature requests as issues to https://github.com/QMCSoftware/QMCSoftware/issues Moderately Diﬃcult Ask your students or collaborators to try QMCPy themselves and submit their bugs and feature requests Email us your blog post to add to qmcpy.wordpress.com/ Heroic Add a feature or use case and make a pull request at https://github.com/QMCSoftware/QMCSoftware/pulls so that we can included it in our next release Questions or suggestions? Email us at [email protected] 11/12
contribute Try out QMCPy and then Easy Submit your bugs and feature requests as issues to https://github.com/QMCSoftware/QMCSoftware/issues Moderately Diﬃcult Ask your students or collaborators to try QMCPy themselves and submit their bugs and feature requests Email us your blog post to add to qmcpy.wordpress.com/ Heroic Add a feature or use case and make a pull request at https://github.com/QMCSoftware/QMCSoftware/pulls so that we can included it in our next release Questions or suggestions? Email us at [email protected] If you are wedded to another language, think of designing your software so that others can add to it easily. 11/12