QMCPy Documentation¶
Discrete Distribution Class¶
Abstract Discrete Distribution Class¶
-
class
qmcpy.discrete_distribution._discrete_distribution.
DiscreteDistribution
¶ Discrete Distribution abstract class. DO NOT INSTANTIATE.
-
gen_samples
(*args)¶ ABSTRACT METHOD to generate samples from this discrete distribution.
- Parameters
args (tuple) – tuple of positional argument. See implementations for details
- Returns
n x d array of samples
- Return type
ndarray
-
pdf
(x)¶ ABSTRACT METHOD to evaluate pdf of distribution the samples mimic at locations of x.
-
plot
(dim_x=0, dim_y=1, n=128, point_size=5, color='c', show=True, out=None)¶ Make a scatter plot from samples. Requires dimension >= 2.
- Parameters
dim_x (int) – index of first dimension to be plotted on horizontal axis.
dim_y (int) – index of the second dimension to be plotted on vertical axis.
n (int) – number of samples to draw as self.gen_samples(n)
point_size (int) – ax.scatter(…,s=point_size)
color (str) – ax.scatter(…,color=color)
show (bool) – show plot or not?
out (str) – file name to output image. If None, the image is not output
- Returns
fig,ax (tuple) from fig,ax = matplotlib.pyplot.subplots(…)
-
Lattice¶
-
class
qmcpy.discrete_distribution.lattice.lattice.
Lattice
(dimension=1, randomize=True, order='natural', seed=None, z_path=None)¶ Quasi-Random Lattice nets in base 2.
>>> l = Lattice(2,seed=7) >>> l.gen_samples(4) array([[0.076, 0.78 ], [0.576, 0.28 ], [0.326, 0.53 ], [0.826, 0.03 ]]) >>> l.gen_samples(1) array([[0.076, 0.78 ]]) >>> l Lattice (DiscreteDistribution Object) d 2^(1) randomize 1 order natural seed 7 mimics StdUniform >>> Lattice(dimension=2,randomize=False,order='natural').gen_samples(4, warn=False) array([[0. , 0. ], [0.5 , 0.5 ], [0.25, 0.75], [0.75, 0.25]]) >>> Lattice(dimension=2,randomize=False,order='linear').gen_samples(4, warn=False) array([[0. , 0. ], [0.25, 0.75], [0.5 , 0.5 ], [0.75, 0.25]]) >>> Lattice(dimension=2,randomize=False,order='mps').gen_samples(4, warn=False) array([[0. , 0. ], [0.5 , 0.5 ], [0.25, 0.75], [0.75, 0.25]])
References
[1] Sou-Cheng T. Choi, Yuhan Ding, Fred J. Hickernell, Lan Jiang, Lluis Antoni Jimenez Rugama, Da Li, Jagadeeswaran Rathinavel, Xin Tong, Kan Zhang, Yizhi Zhang, and Xuan Zhou, GAIL: Guaranteed Automatic Integration Library (Version 2.3) [MATLAB Software], 2019. Available from http://gailgithub.github.io/GAIL_Dev/
[2] F.Y. Kuo & D. Nuyens. Application of quasi-Monte Carlo methods to elliptic PDEs with random diffusion coefficients - a survey of analysis and implementation, Foundations of Computational Mathematics, 16(6):1631-1696, 2016. springer link: https://link.springer.com/article/10.1007/s10208-016-9329-5 arxiv link: https://arxiv.org/abs/1606.06613
[3] D. Nuyens, The Magic Point Shop of QMC point generators and generating vectors. MATLAB and Python software, 2018. Available from https://people.cs.kuleuven.be/~dirk.nuyens/
[4] Constructing embedded lattice rules for multivariate integration R Cools, FY Kuo, D Nuyens - SIAM J. Sci. Comput., 28(6), 2162-2188.
[5] L’Ecuyer, Pierre & Munger, David. (2015). LatticeBuilder: A General Software Tool for Constructing Rank-1 Lattice Rules. ACM Transactions on Mathematical Software. 42. 10.1145/2754929.
-
__init__
(dimension=1, randomize=True, order='natural', seed=None, z_path=None)¶ - Parameters
dimension (int) – dimension of samples
randomize (bool) – If True, apply shift to generated samples. Note: Non-randomized lattice sequence includes the origin.
order (str) – ‘linear’, ‘natural’, or ‘mps’ ordering.
seed (int) – seed the random number generator for reproducibility
z_path (str) – path to generating vector. z_path should be formatted like ‘lattice_vec.3600.20.npy’ where ‘name.d_max.m_max.npy’ and d_max is the maximum dimenion and 2^m_max is the max number samples supported
-
gen_samples
(n=None, n_min=0, n_max=8, warn=True, return_unrandomized=False)¶ Generate lattice samples
- Parameters
n (int) – if n is supplied, generate from n_min=0 to n_max=n samples. Otherwise use the n_min and n_max explicitly supplied as the following 2 arguments
n_min (int) – Starting index of sequence.
n_max (int) – Final index of sequence.
return_unrandomized (bool) – return samples without randomization as 2nd return value. Will not be returned if randomize=False.
- Returns
(n_max-n_min) x d (dimension) array of samples
- Return type
ndarray
Note
Lattice generates in blocks from 2**m to 2**(m+1) so generating n_min=3 to n_max=9 requires necessarily produces samples from n_min=2 to n_max=16 and automatically subsets. May be inefficient for non-powers-of-2 samples sizes.
-
pdf
(x)¶ pdf of a standard uniform
-
set_seed
(seed)¶ See abstract method.
-
Sobol’¶
-
qmcpy.discrete_distribution.sobol.sobol.
DigitalNet
¶
-
class
qmcpy.discrete_distribution.sobol.sobol.
Sobol
(dimension=1, randomize='LMS', graycode=False, seed=None, z_path=None, dim0=0)¶ Quasi-Random Sobol nets in base 2.
>>> s = Sobol(2,seed=7) >>> s.gen_samples(4) array([[0.673, 0.063], [0.398, 0.788], [0.913, 0.564], [0.125, 0.288]]) >>> s.gen_samples(1) array([[0.673, 0.063]]) >>> s Sobol (DiscreteDistribution Object) d 2^(1) randomize 1 graycode 0 seed [61616 58565] mimics StdUniform dim0 0 >>> Sobol(dimension=2,randomize=False,graycode=True).gen_samples(n_min=2,n_max=4) array([[0.75, 0.25], [0.25, 0.75]]) >>> Sobol(dimension=2,randomize=False,graycode=False).gen_samples(n_min=2,n_max=4) array([[0.25, 0.75], [0.75, 0.25]])
References
[1] Marius Hofert and Christiane Lemieux (2019). qrng: (Randomized) Quasi-Random Number Generators. R package version 0.0-7. https://CRAN.R-project.org/package=qrng.
[2] Faure, Henri, and Christiane Lemieux. “Implementation of Irreducible Sobol’ Sequences in Prime Power Bases.” Mathematics and Computers in Simulation 161 (2019): 13–22. Crossref. Web.
[3] F.Y. Kuo & D. Nuyens. Application of quasi-Monte Carlo methods to elliptic PDEs with random diffusion coefficients - a survey of analysis and implementation, Foundations of Computational Mathematics, 16(6):1631-1696, 2016. springer link: https://link.springer.com/article/10.1007/s10208-016-9329-5 arxiv link: https://arxiv.org/abs/1606.06613
[4] D. Nuyens, The Magic Point Shop of QMC point generators and generating vectors. MATLAB and Python software, 2018. Available from https://people.cs.kuleuven.be/~dirk.nuyens/
[5] Paszke, A., Gross, S., Massa, F., Lerer, A., Bradbury, J., Chanan, G., … Chintala, S. (2019). PyTorch: An Imperative Style, High-Performance Deep Learning Library. In H. Wallach, H. Larochelle, A. Beygelzimer, F. d extquotesingle Alch'e-Buc, E. Fox, & R. Garnett (Eds.), Advances in Neural Information Processing Systems 32 (pp. 8024–8035). Curran Associates, Inc. Retrieved from http://papers.neurips.cc/paper/9015-pytorch-an-imperative-style-high-performance-deep-learning-library.pdf
[6] I.M. Sobol’, V.I. Turchaninov, Yu.L. Levitan, B.V. Shukhman: “Quasi-Random Sequence Generators” Keldysh Institute of Applied Mathematics, Russian Acamdey of Sciences, Moscow (1992).
[7] Sobol, Ilya & Asotsky, Danil & Kreinin, Alexander & Kucherenko, Sergei. (2011). Construction and Comparison of High-Dimensional Sobol’ Generators. Wilmott. 2011. 10.1002/wilm.10056.
[8] Paul Bratley and Bennett L. Fox. 1988. Algorithm 659: Implementing Sobol’s quasirandom sequence generator. ACM Trans. Math. Softw. 14, 1 (March 1988), 88–100. DOI:https://doi.org/10.1145/42288.214372
-
__init__
(dimension=1, randomize='LMS', graycode=False, seed=None, z_path=None, dim0=0)¶ - Parameters
dimension (int) – dimension of samples
randomize (bool) – Apply randomization? True defaults to LMS. Can also explicitly pass in ‘LMS’: Linear matrix scramble with DS ‘DS’: Just Digital Shift
graycode (bool) – indicator to use graycode ordering (True) or natural ordering (False)
seeds (list) – int seed of list of seeds, one for each dimension.
z_path (str) – path to generating matricies. z_path sould be formatted like gen_mat.21201.32.msb.npy with name.d_max.m_max.msb_or_lsb.npy
dim0 (int) – first dimension
-
gen_samples
(n=None, n_min=0, n_max=8, warn=True, return_jlms=False)¶ Generate samples
- Parameters
n (int) – if n is supplied, generate from n_min=0 to n_max=n samples. Otherwise use the n_min and n_max explicitly supplied as the following 2 arguments
n_min (int) – Starting index of sequence.
n_max (int) – Final index of sequence.
return_jlms (bool) – return the LMS matrix without digital shift. Only applies when randomize=’LMS’ (the default).
- Returns
(n_max-n_min) x d (dimension) array of samples
- Return type
ndarray
-
pdf
(x)¶ pdf of a standard uniform
-
set_randomize
(randomize)¶ Reset the randomization
- Parameters
randomize (str) – randomization type. Either ‘LMS’: linear matrix scramble with digital shift ‘DS’: just the digital shift
-
set_seed
(seeds)¶ Reset the seeds
- Parameters
seeds (int/list/None) – new seeds
-
Halton¶
-
class
qmcpy.discrete_distribution.halton.halton.
Halton
(dimension=1, generalize=True, randomize=True, seed=None)¶ Quasi-Random Halton nets.
>>> h = Halton(2,seed=7) >>> h.gen_samples(4) array([[0.166, 0.363], [0.666, 0.696], [0.416, 0.03 ], [0.916, 0.474]]) >>> h.gen_samples(1) array([[0.166, 0.363]]) >>> h Halton (DiscreteDistribution Object) d 2^(1) generalize 1 randomize 1 seed 7 mimics StdUniform
References
[1] Marius Hofert and Christiane Lemieux (2019). qrng: (Randomized) Quasi-Random Number Generators. R package version 0.0-7. https://CRAN.R-project.org/package=qrng.
[2] Owen, A. B. “A randomized Halton algorithm in R,” 2017. arXiv:1706.02808 [stat.CO]
-
__init__
(dimension=1, generalize=True, randomize=True, seed=None)¶ - Parameters
dimension (int) – dimension of samples
generalize (bool) – generalize the Halton sequence?
randomize (bool/str) – If False, does not randomize Halton points. If True, will use ‘QRNG’ randomization as in [1]. Supports max dimension 360. You can also set radnomize=’QRNG’ or randomize=’Halton’ to explicitly select a randomization method. Halton OWEN supports up to 1000 dimensions.
seed (int) – seed the random number generator for reproducibility
Note
See References [1] and [2] for specific randomization methods and differences.
-
gen_samples
(n=None, n_min=0, n_max=8, warn=True)¶ Generate samples
- Parameters
- Returns
(n_max-n_min) x d (dimension) array of samples
- Return type
ndarray
-
pdf
(x)¶ ABSTRACT METHOD to evaluate pdf of distribution the samples mimic at locations of x.
-
Korobov¶
-
class
qmcpy.discrete_distribution.korobov.korobov.
Korobov
(dimension=1, generator=[1], randomize=True, seed=None)¶ Quasi-Random Korobov nets.
>>> k = Korobov(2,generator=[1,3],seed=7) >>> k.gen_samples(4) array([[0.982, 0.883], [0.232, 0.633], [0.482, 0.383], [0.732, 0.133]]) >>> k.gen_samples(4) array([[0.982, 0.883], [0.232, 0.633], [0.482, 0.383], [0.732, 0.133]]) >>> k Korobov (DiscreteDistribution Object) d 2^(1) generator [1 3] randomize 1 seed 7 mimics StdUniform >>> Korobov(2,generator=[3,1],seed=7).gen_samples(4) array([[0.982, 0.883], [0.732, 0.133], [0.482, 0.383], [0.232, 0.633]])
References
[1] Marius Hofert and Christiane Lemieux (2019). qrng: (Randomized) Quasi-Random Number Generators. R package version 0.0-7. https://CRAN.R-project.org/package=qrng.
-
__init__
(dimension=1, generator=[1], randomize=True, seed=None)¶ - Parameters
dimension (int) – dimension of samples
generator (ndarray of ints) – generator in {1,..,n-1} either a vector of length d or a single number (which is appropriately extended)
randomize (bool) – randomize the Korobov sequence? Note: Non-randomized Korobov sequence includes origin
seed (int) – seed the random number generator for reproducibility
-
gen_samples
(n=None, n_min=0, n_max=8, warn=True)¶ Generate samples
- Parameters
n (int) – number of samples
- Returns
n x d (dimension) array of samples
- Return type
ndarray
-
pdf
(x)¶ pdf of a standard uniform
-
set_seed
(seed)¶ See abstract method.
-
IID Standard Uniform¶
-
class
qmcpy.discrete_distribution.iid_std_uniform.
IIDStdUniform
(dimension=1, seed=None)¶ A wrapper around NumPy’s IID Standard Uniform generator numpy.random.rand.
>>> dd = IIDStdUniform(dimension=2,seed=7) >>> dd.gen_samples(4) array([[0.076, 0.78 ], [0.438, 0.723], [0.978, 0.538], [0.501, 0.072]]) >>> dd IIDStdUniform (DiscreteDistribution Object) d 2^(1) seed 7 mimics StdUniform
-
__init__
(dimension=1, seed=None)¶
-
gen_samples
(n)¶ Generate samples
- Parameters
n (int) – Number of observations to generate
- Returns
n x self.d array of samples
- Return type
ndarray
-
pdf
(x)¶ ABSTRACT METHOD to evaluate pdf of distribution the samples mimic at locations of x.
-
True Measure Class¶
Abstract Measure Class¶
Uniform¶
-
class
qmcpy.true_measure.uniform.
Uniform
(sampler, lower_bound=0.0, upper_bound=1.0)¶ >>> u = Uniform(Sobol(2,seed=7),lower_bound=[0,.5],upper_bound=[2,3]) >>> u.gen_samples(4) array([[1.346, 0.658], [0.797, 2.471], [1.826, 1.909], [0.25 , 1.22 ]]) >>> u Uniform (TrueMeasure Object) lower_bound [0. 0.5] upper_bound [2 3]
Gaussian¶
-
class
qmcpy.true_measure.gaussian.
Gaussian
(sampler, mean=0.0, covariance=1.0, decomp_type='PCA')¶ Normal Measure.
>>> g = Gaussian(Sobol(2,seed=7),mean=[1,2],covariance=[[9,4],[4,5]]) >>> g.gen_samples(4) array([[-1.567, 3.268], [ 2.411, 1.376], [-2.784, -0.638], [ 3.846, 4.804]]) >>> g Gaussian (TrueMeasure Object) mean [1 2] covariance [[9 4] [4 5]] decomp_type pca
-
__init__
(sampler, mean=0.0, covariance=1.0, decomp_type='PCA')¶ - Parameters
sampler (DiscreteDistribution/TrueMeasure) – A discrete distribution from which to transform samples or a true measure by which to compose a transform
mean (float) – mu for Normal(mu,sigma^2)
covariance (ndarray) – sigma^2 for Normal(mu,sigma^2). A float or d (dimension) vector input will be extended to covariance*eye(d)
decomp_type (str) – method of decomposition either “PCA” for principal component analysis or “Cholesky” for cholesky decomposition.
-
Brownian Motion¶
-
class
qmcpy.true_measure.brownian_motion.
BrownianMotion
(sampler, t_final=1, drift=0, decomp_type='PCA')¶ Geometric Brownian Motion.
>>> bm = BrownianMotion(Sobol(4,seed=7),t_final=2,drift=2) >>> bm.gen_samples(2) array([[0.583, 0.514, 2.726, 3.967], [0.958, 2.973, 3.259, 3.951]]) >>> bm BrownianMotion (TrueMeasure Object) time_vec [0.5 1. 1.5 2. ] drift 2^(1) mean [1. 2. 3. 4.] covariance [[0.5 0.5 0.5 0.5] [0.5 1. 1. 1. ] [0.5 1. 1.5 1.5] [0.5 1. 1.5 2. ]] decomp_type pca
-
__init__
(sampler, t_final=1, drift=0, decomp_type='PCA')¶ - Parameters
sampler (DiscreteDistribution/TrueMeasure) – A discrete distribution from which to transform samples or a true measure by which to compose a transform
t_final (float) – end time for the Brownian Motion.
drift (int) – Gaussian mean is time_vec*drift
decomp_type (str) – method of decomposition either “PCA” for principal component analysis or “Cholesky” for cholesky decomposition.
-
Lebesgue¶
-
class
qmcpy.true_measure.lebesgue.
Lebesgue
(sampler)¶ >>> Lebesgue(Gaussian(Sobol(2,seed=7))) Lebesgue (TrueMeasure Object) transform Gaussian (TrueMeasure Object) mean 0 covariance 1 decomp_type pca >>> Lebesgue(Uniform(Sobol(2,seed=7))) Lebesgue (TrueMeasure Object) transform Uniform (TrueMeasure Object) lower_bound 0 upper_bound 1
-
__init__
(sampler)¶ - Parameters
sampler (TrueMeasure) – A true measure by which to compose a transform.
-
Kumaraswamy¶
-
class
qmcpy.true_measure.kumaraswamy.
Kumaraswamy
(sampler, a=2, b=2)¶ For
we have
PDF:
CDF:
Inverse CDF:
-
__init__
(sampler, a=2, b=2)¶ - Parameters
sampler (DiscreteDistribution/TrueMeasure) – A discrete distribution from which to transform samples or a true measure by which to compose a transform
a (ndarray) – a
b (ndarray) – b
-
Integrand Class¶
Abstract Integrand Class¶
-
class
qmcpy.integrand._integrand.
Integrand
¶ Integrand abstract class. DO NOT INSTANTIATE.
-
f
(x, *args, **kwargs)¶ Evalute transformed integrand based on true measures and discrete distribution
- Parameters
x (ndarray) – n x d array of samples from a discrete distribution
*args – other ordered args to g
**kwargs (dict) – other keyword args to g
- Returns
length n vector of funciton evaluations
- Return type
ndarray
-
f_periodized
(x, ptransform='NONE', *args, **kwargs)¶ Periodized transformed integrand.
-
g
(t, *args, **kwargs)¶ ABSTRACT METHOD for original integrand to be integrated.
- Parameters
t (ndarray) – n x d array of samples to be intput into orignal integrand.
- Returns
n vector of function evaluations
- Return type
ndarray
-
Keister Function¶
-
class
qmcpy.integrand.keister.
Keister
(sampler)¶ .
The standard example integrates the Keister integrand with respect to an IID Gaussian distribution with variance 1./2.
>>> k = Keister(Sobol(2,seed=7)) >>> x = k.discrete_distrib.gen_samples(2**10) >>> y = k.f(x) >>> y.mean() 1.807... >>> k.true_measure Lebesgue (TrueMeasure Object) transform Gaussian (TrueMeasure Object) mean 0 covariance 2^(-1) decomp_type pca >>> k = Keister(Gaussian(Sobol(2,seed=7),mean=0,covariance=2)) >>> x = k.discrete_distrib.gen_samples(2**10) >>> y = k.f(x) >>> y.mean() 1.808... >>> yp = k.f_periodized(x,'c2sin') >>> yp.mean() 1.805...
References
[1] B. D. Keister, Multidimensional Quadrature Algorithms, Computers in Physics, 10, pp. 119-122, 1996.
-
__init__
(sampler)¶ - Parameters
sampler (DiscreteDistribution/TrueMeasure) – A discrete distribution from which to transform samples or a true measure by which to compose a transform
-
g
(x)¶ ABSTRACT METHOD for original integrand to be integrated.
- Parameters
t (ndarray) – n x d array of samples to be intput into orignal integrand.
- Returns
n vector of function evaluations
- Return type
ndarray
-
Custom Function¶
-
class
qmcpy.integrand.custom_fun.
CustomFun
(true_measure, g)¶ Custom user-supplied function handle.
>>> cf = CustomFun( ... true_measure = Gaussian(Sobol(2,seed=7),mean=[1,2]), ... g = lambda x: x[:,0]**2*x[:,1]) >>> x = cf.discrete_distrib.gen_samples(2**10) >>> y = cf.f(x) >>> y.mean() 4.009...
-
__init__
(true_measure, g)¶ - Parameters
true_measure (TrueMeasure) – a TrueMeasure instance.
g (function) – a function handle.
-
European Option¶
-
class
qmcpy.integrand.european_option.
EuropeanOption
(sampler, volatility=0.5, start_price=30, strike_price=35, interest_rate=0, t_final=1, call_put='call')¶ European financial option.
>>> eo = EuropeanOption(Sobol(4,seed=7),call_put='put') >>> eo EuropeanOption (Integrand Object) volatility 2^(-1) call_put put start_price 30 strike_price 35 interest_rate 0 >>> x = eo.discrete_distrib.gen_samples(2**12) >>> y = eo.f(x) >>> y.mean() 9.210... >>> eo = EuropeanOption(BrownianMotion(Sobol(4,seed=7),drift=1),call_put='put') >>> x = eo.discrete_distrib.gen_samples(2**12) >>> y = eo.f(x) >>> y.mean() 9.220...
-
__init__
(sampler, volatility=0.5, start_price=30, strike_price=35, interest_rate=0, t_final=1, call_put='call')¶ - Parameters
sampler (DiscreteDistribution/TrueMeasure) – A discrete distribution from which to transform samples or a true measure by which to compose a transform
volatility (float) – sigma, the volatility of the asset
start_price (float) – S(0), the asset value at t=0
strike_price (float) – strike_price, the call/put offer
interest_rate (float) – r, the annual interest rate
t_final (float) – exercise time
call_put (str) – ‘call’ or ‘put’ option
-
g
(x)¶ See abstract method.
-
Asian Option¶
-
class
qmcpy.integrand.asian_option.
AsianOption
(sampler, volatility=0.5, start_price=30.0, strike_price=35.0, interest_rate=0.0, t_final=1, call_put='call', mean_type='arithmetic', multi_level_dimensions=None)¶ Asian financial option.
>>> ac = AsianOption(Sobol(4,seed=7)) >>> ac AsianOption (Integrand Object) volatility 2^(-1) call_put call start_price 30 strike_price 35 interest_rate 0 mean_type arithmetic dimensions 2^(2) dim_fracs 0 >>> x = ac.discrete_distrib.gen_samples(2**10) >>> y = ac.f(x) >>> y.mean() 1.773... >>> level_dims = [2,4,8] >>> ac2 = AsianOption(Sobol(seed=7),multi_level_dimensions=level_dims) >>> ac2 AsianOption (Integrand Object) volatility 2^(-1) call_put call start_price 30 strike_price 35 interest_rate 0 mean_type arithmetic dimensions [2 4 8] dim_fracs [0. 2. 2.] >>> y2 = 0 >>> for l in range(len(level_dims)): ... new_dim = ac2._dim_at_level(l) ... ac2.true_measure._set_dimension_r(new_dim) ... x2 = ac2.discrete_distrib.gen_samples(2**10) ... level_est = ac2.f(x2,l=l).mean() ... y2 += level_est >>> y2 1.793...
-
__init__
(sampler, volatility=0.5, start_price=30.0, strike_price=35.0, interest_rate=0.0, t_final=1, call_put='call', mean_type='arithmetic', multi_level_dimensions=None)¶ - Parameters
sampler (DiscreteDistribution/TrueMeasure) – A discrete distribution from which to transform samples or a true measure by which to compose a transform
volatility (float) – sigma, the volatility of the asset
start_price (float) – S(0), the asset value at t=0
strike_price (float) – strike_price, the call/put offer
interest_rate (float) – r, the annual interest rate
t_final (float) – exercise time
mean_type (string) – ‘arithmetic’ or ‘geometric’ mean
multi_level_dimensions (list of ints) – list of dimensions at each level. Leave as None for single-level problems
-
g
(x, l=0)¶ See abstract method.
-
Multilevel Call Options with Milstein Discretization¶
-
class
qmcpy.integrand.ml_call_options.
MLCallOptions
(sampler, option='european', volatility=0.2, start_strike_price=100.0, interest_rate=0.05, t_final=1.0)¶ Various call options from finance using Milstein discretization with
timesteps on level
.
>>> mlco = MLCallOptions(Sobol(seed=7)) >>> mlco MLCallOptions (Integrand Object) option european sigma 0.200 k 100 r 0.050 t 1 b 85 >>> y = 0 >>> for level in range(4): ... new_dim = mlco._dim_at_level(level) ... mlco.true_measure._set_dimension_r(new_dim) ... x = mlco.discrete_distrib.gen_samples(2**10) ... y += mlco.f(x,l=level).mean() >>> y 10.406...
References:
[1] M.B. Giles. Improved multilevel Monte Carlo convergence using the Milstein scheme. 343-358, in Monte Carlo and Quasi-Monte Carlo Methods 2006, Springer, 2008. http://people.maths.ox.ac.uk/~gilesm/files/mcqmc06.pdf.
-
__init__
(sampler, option='european', volatility=0.2, start_strike_price=100.0, interest_rate=0.05, t_final=1.0)¶ - Parameters
sampler (DiscreteDistribution/TrueMeasure) – A discrete distribution from which to transform samples or a true measure by which to compose a transform
option_type (str) – type of option in [“European”,”Asian”]
volatility (float) – sigma, the volatility of the asset
start_strike_price (float) – S(0), the asset value at t=0, and K, the strike price. Assume start_price = strike_price
interest_rate (float) – r, the annual interest rate
t_final (float) – exercise time
-
g
(samples, l)¶
-
get_exact_value
()¶ Print exact analytic value, based on s0=k.
-
Linear Function¶
-
class
qmcpy.integrand.linear0.
Linear0
(sampler)¶ >>> l = Linear0(Sobol(100,seed=7)) >>> x = l.discrete_distrib.gen_samples(2**10) >>> y = l.f(x) >>> y.mean() -2.654...e-08 >>> ytf = l.f_periodized(x,'C1SIN') >>> ytf.mean() 8.288...e-10
-
__init__
(sampler)¶ - Parameters
sampler (DiscreteDistribution/TrueMeasure) – A discrete distribution from which to transform samples or a true measure by which to compose a transform
-
g
(x)¶ ABSTRACT METHOD for original integrand to be integrated.
- Parameters
t (ndarray) – n x d array of samples to be intput into orignal integrand.
- Returns
n vector of function evaluations
- Return type
ndarray
-
Stopping Criterion Algorithms¶
Abstract Stopping Criterion Class¶
-
class
qmcpy.stopping_criterion._stopping_criterion.
StoppingCriterion
(allowed_levels, allowed_distribs)¶ Stopping Criterion abstract class. DO NOT INSTANTIATE.
-
__init__
(allowed_levels, allowed_distribs)¶ - Parameters
distribution (DiscreteDistribution) – a DiscreteDistribution
allowed_levels (list) – which integrand types are supported: ‘single’, ‘fixed-multi’, ‘adaptive-multi’
allowed_distribs (list) – list of names (strings) of compatible distributions
-
integrate
()¶ ABSTRACT METHOD to determine the number of samples needed to satisfy the tolerance.
- Returns
- tuple containing:
solution (float): approximation to the integral
data (AccumulateData): an AccumulateData object
- Return type
-
plot
(*args, **kwargs)¶ Create a plot relevant to the stopping criterion object.
-
set_tolerance
(*args, **kwargs)¶ ABSTRACT METHOD to reset the absolute tolerance.
-
Guaranteed Lattice Cubature (qMC)¶
-
class
qmcpy.stopping_criterion.cub_qmc_lattice_g.
CubQMCLatticeG
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=1024.0, n_max=34359738368.0, fudge=<function CubQMCLatticeG.<lambda>>, check_cone=False, ptransform='Baker')¶ Stopping Criterion quasi-Monte Carlo method using rank-1 Lattices cubature over a d-dimensional region to integrate within a specified generalized error tolerance with guarantees under Fourier coefficients cone decay assumptions.
>>> k = Keister(Lattice(2,seed=7)) >>> sc = CubQMCLatticeG(k,abs_tol=.05) >>> solution,data = sc.integrate() >>> solution 1.806... >>> data Solution: 1.8068 Keister (Integrand Object) Lattice (DiscreteDistribution Object) d 2^(1) randomize 1 order natural seed 7 mimics StdUniform Lebesgue (TrueMeasure Object) transform Gaussian (TrueMeasure Object) mean 0 covariance 2^(-1) decomp_type pca CubQMCLatticeG (StoppingCriterion Object) abs_tol 0.050 rel_tol 0 n_init 2^(10) n_max 2^(35) LDTransformData (AccumulateData Object) n_total 2^(10) solution 1.807 error_bound 0.005 time_integrate ...
Original Implementation:
References
[1] Lluis Antoni Jimenez Rugama and Fred J. Hickernell, “Adaptive multidimensional integration based on rank-1 lattices,” Monte Carlo and Quasi-Monte Carlo Methods: MCQMC, Leuven, Belgium, April 2014 (R. Cools and D. Nuyens, eds.), Springer Proceedings in Mathematics and Statistics, vol. 163, Springer-Verlag, Berlin, 2016, arXiv:1411.1966, pp. 407-422.
[2] Sou-Cheng T. Choi, Yuhan Ding, Fred J. Hickernell, Lan Jiang, Lluis Antoni Jimenez Rugama, Da Li, Jagadeeswaran Rathinavel, Xin Tong, Kan Zhang, Yizhi Zhang, and Xuan Zhou, GAIL: Guaranteed Automatic Integration Library (Version 2.3) [MATLAB Software], 2019. Available from http://gailgithub.github.io/GAIL_Dev/
- Guarantee
This algorithm computes the integral of real valued functions in
with a prescribed generalized error tolerance. The Fourier coefficients of the integrand are assumed to be absolutely convergent. If the algorithm terminates without warning messages, the output is given with guarantees under the assumption that the integrand lies inside a cone of functions. The guarantee is based on the decay rate of the Fourier coefficients. For integration over domains other than
, this cone condition applies to
(the composition of the functions) where
is the transformation function for
to the desired region. For more details on how the cone is defined, please refer to the references below.
-
__init__
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=1024.0, n_max=34359738368.0, fudge=<function CubQMCLatticeG.<lambda>>, check_cone=False, ptransform='Baker')¶ - Parameters
integrand (Integrand) – an instance of Integrand
abs_tol (float) – absolute error tolerance
rel_tol (float) – relative error tolerance
n_init (int) – initial number of samples
n_max (int) – maximum number of samples
fudge (function) – positive function multiplying the finite sum of Fast Fourier coefficients specified in the cone of functions
check_cone (boolean) – check if the function falls in the cone
-
integrate
()¶ See abstract method.
Guaranteed Sobol’ Cubature (qMC)¶
-
class
qmcpy.stopping_criterion.cub_qmc_sobol_g.
CubQMCSobolG
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=1024.0, n_max=34359738368.0, fudge=<function CubQMCSobolG.<lambda>>, check_cone=False)¶ Quasi-Monte Carlo method using Sobol’ cubature over the d-dimensional region to integrate within a specified generalized error tolerance with guarantees under Walsh-Fourier coefficients cone decay assumptions.
>>> k = Keister(Sobol(2,seed=7)) >>> sc = CubQMCSobolG(k,abs_tol=.05) >>> solution,data = sc.integrate() >>> solution 1.807... >>> data Solution: 1.8074 Keister (Integrand Object) Sobol (DiscreteDistribution Object) d 2^(1) randomize 1 graycode 0 seed [61616 58565] mimics StdUniform dim0 0 Lebesgue (TrueMeasure Object) transform Gaussian (TrueMeasure Object) mean 0 covariance 2^(-1) decomp_type pca CubQMCSobolG (StoppingCriterion Object) abs_tol 0.050 rel_tol 0 n_init 2^(10) n_max 2^(35) LDTransformData (AccumulateData Object) n_total 2^(10) solution 1.807 error_bound 0.005 time_integrate ...
Original Implementation:
References
[1] Fred J. Hickernell and Lluis Antoni Jimenez Rugama, Reliable adaptive cubature using digital sequences, 2014. Submitted for publication: arXiv:1410.8615.
[2] Sou-Cheng T. Choi, Yuhan Ding, Fred J. Hickernell, Lan Jiang, Lluis Antoni Jimenez Rugama, Da Li, Jagadeeswaran Rathinavel, Xin Tong, Kan Zhang, Yizhi Zhang, and Xuan Zhou, GAIL: Guaranteed Automatic Integration Library (Version 2.3) [MATLAB Software], 2019. Available from http://gailgithub.github.io/GAIL_Dev/
- Guarantee:
This algorithm computes the integral of real valued functions in
with a prescribed generalized error tolerance. The Fourier coefficients of the integrand are assumed to be absolutely convergent. If the algorithm terminates without warning messages, the output is given with guarantees under the assumption that the integrand lies inside a cone of functions. The guarantee is based on the decay rate of the Fourier coefficients. For integration over domains other than
, this cone condition applies to
(the composition of the functions) where
is the transformation function for
to the desired region. For more details on how the cone is defined, please refer to the references below.
-
__init__
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=1024.0, n_max=34359738368.0, fudge=<function CubQMCSobolG.<lambda>>, check_cone=False)¶ - Parameters
integrand (Integrand) – an instance of Integrand
abs_tol (float) – absolute error tolerance
rel_tol (float) – relative error tolerance
n_init (int) – initial number of samples
n_max (int) – maximum number of samples
fudge (function) – positive function multiplying the finite sum of Fast Fourier coefficients specified in the cone of functions
check_cone (boolean) – check if the function falls in the cone
-
integrate
()¶ See abstract method.
Bayesian Lattice Cubature (qMC)¶
-
class
qmcpy.stopping_criterion.cub_qmc_bayes_lattice_g.
CubBayesLatticeG
(integrand, abs_tol=0.01, rel_tol=0, n_init=256, n_max=4194304, alpha=0.01, ptransform='C1sin')¶ Stopping criterion for Bayesian Cubature using rank-1 Lattice sequence with guaranteed accuracy over a d-dimensional region to integrate within a specified generalized error tolerance with guarantees under Bayesian assumptions.
>>> k = Keister(Lattice(2, order='linear', seed=123456789)) >>> sc = CubBayesLatticeG(k,abs_tol=.05) >>> solution,data = sc.integrate() >>> solution 1.808... >>> data Solution: 1.8082 Keister (Integrand Object) Lattice (DiscreteDistribution Object) d 2^(1) randomize 1 order linear seed 123456789 mimics StdUniform Lebesgue (TrueMeasure Object) transform Gaussian (TrueMeasure Object) mean 0 covariance 2^(-1) decomp_type pca CubBayesLatticeG (StoppingCriterion Object) abs_tol 0.050 rel_tol 0 n_init 2^(8) n_max 2^(22) LDTransformBayesData (AccumulateData Object) n_total 256 solution 1.808 error_bound 7.36e-04 time_integrate ...
- Adapted from
- Reference
[1] Sou-Cheng T. Choi, Yuhan Ding, Fred J. Hickernell, Lan Jiang, Lluis Antoni Jimenez Rugama, Da Li, Jagadeeswaran Rathinavel, Xin Tong, Kan Zhang, Yizhi Zhang, and Xuan Zhou, GAIL: Guaranteed Automatic Integration Library (Version 2.3) [MATLAB Software], 2019. Available from http://gailgithub.github.io/GAIL_Dev/
- Guarantee
This algorithm attempts to calculate the integral of function f over the hyperbox [0,1]^d to a prescribed error tolerance tolfun:= max(abstol,reltol*| I |) with guaranteed confidence level, e.g., 99% when alpha=0.5%. If the algorithm terminates without showing any warning messages and provides an answer Q, then the following inequality would be satisfied:
Pr(| Q - I | <= tolfun) = 99%
This Bayesian cubature algorithm guarantees for integrands that are considered to be an instance of a gaussian process that fall in the middle of samples space spanned. Where The sample space is spanned by the covariance kernel parametrized by the scale and shape parameter inferred from the sampled values of the integrand. For more details on how the covariance kernels are defined and the parameters are obtained, please refer to the references below.
Bayesian Digital Net Cubature (qMC)¶
-
class
qmcpy.stopping_criterion.cub_qmc_bayes_net_g.
CubBayesNetG
(integrand, abs_tol=0.01, rel_tol=0, n_init=256, n_max=4194304, alpha=0.01)¶ Stopping criterion for Bayesian Cubature using digital net (Sobol) sequence with guaranteed accuracy over a d-dimensional region to integrate within a specified generalized error tolerance with guarantees under Bayesian assumptions.
>>> k = Keister(Sobol(2, seed=123456789)) >>> sc = CubBayesNetG(k,abs_tol=.05) >>> solution,data = sc.integrate() >>> solution 1.811... >>> data Solution: 1.8110 Keister (Integrand Object) Sobol (DiscreteDistribution Object) d 2^(1) randomize 1 graycode 0 seed [77469 77107] mimics StdUniform dim0 0 Lebesgue (TrueMeasure Object) transform Gaussian (TrueMeasure Object) mean 0 covariance 2^(-1) decomp_type pca CubBayesNetG (StoppingCriterion Object) abs_tol 0.050 rel_tol 0 n_init 2^(8) n_max 2^(22) LDTransformBayesData (AccumulateData Object) n_total 256 solution 1.811 error_bound 0.015 time_integrate ...
- Adapted from
https://github.com/GailGithub/GAIL_Dev/blob/master/Algorithms/IntegrationExpectation/cubBayesNet_g.m
- Reference
[1] Sou-Cheng T. Choi, Yuhan Ding, Fred J. Hickernell, Lan Jiang, Lluis Antoni Jimenez Rugama, Da Li, Jagadeeswaran Rathinavel, Xin Tong, Kan Zhang, Yizhi Zhang, and Xuan Zhou, GAIL: Guaranteed Automatic Integration Library (Version 2.3) [MATLAB Software], 2019. Available from http://gailgithub.github.io/GAIL_Dev/
- Guarantee
This algorithm attempts to calculate the integral of function f over the hyperbox [0,1]^d to a prescribed error tolerance tolfun:= max(abstol,reltol*| I |) with guaranteed confidence level, e.g., 99% when alpha=0.5%. If the algorithm terminates without showing any warning messages and provides an answer Q, then the following inequality would be satisfied:
Pr(| Q - I | <= tolfun) = 99%
This Bayesian cubature algorithm guarantees for integrands that are considered to be an instance of a gaussian process that fall in the middle of samples space spanned. Where The sample space is spanned by the covariance kernel parametrized by the scale and shape parameter inferred from the sampled values of the integrand. For more details on how the covariance kernels are defined and the parameters are obtained, please refer to the references below.
Multilevel qMC Cubature¶
-
class
qmcpy.stopping_criterion.cub_qmc_ml.
CubQMCML
(integrand, abs_tol=0.05, alpha=0.01, rmse_tol=None, n_init=256.0, n_max=10000000000.0, replications=32.0)¶ Stopping criterion based on multi-level quasi-Monte Carlo.
>>> mlco = MLCallOptions(Lattice(seed=7)) >>> sc = CubQMCML(mlco,abs_tol=.05) >>> solution,data = sc.integrate() >>> solution 10.442... >>> data Solution: 10.4424 MLCallOptions (Integrand Object) option european sigma 0.200 k 100 r 0.050 t 1 b 85 Lattice (DiscreteDistribution Object) d 2^(6) randomize 1 order natural seed 985802 mimics StdUniform Gaussian (TrueMeasure Object) mean 0 covariance 1 decomp_type pca CubQMCML (StoppingCriterion Object) rmse_tol 0.019 n_init 2^(8) n_max 10000000000 replications 2^(5) MLQMCData (AccumulateData Object) levels 7 dimensions [ 1. 2. 4. 8. 16. 32. 64.] n_level [4096. 512. 256. 256. 256. 256. 256.] mean_level [1.005e+01 1.807e-01 1.033e-01 5.482e-02 2.823e-02 1.397e-02 7.290e-03] var_level [8.376e-05 2.660e-05 1.911e-05 1.594e-05 3.660e-06 1.478e-06 3.424e-07] bias_estimate 0.007 n_total 188416 time_integrate ...
References
[1] M.B. Giles and B.J. Waterhouse. ‘Multilevel quasi-Monte Carlo path simulation’. pp.165-181 in Advanced Financial Modelling, in Radon Series on Computational and Applied Mathematics, de Gruyter, 2009. http://people.maths.ox.ac.uk/~gilesm/files/radon.pdf
-
__init__
(integrand, abs_tol=0.05, alpha=0.01, rmse_tol=None, n_init=256.0, n_max=10000000000.0, replications=32.0)¶ - Parameters
integrand (Integrand) – integrand with multi-level g method
abs_tol (float) – absolute tolerance
alpha (float) – uncertaintly level. If rmse_tol not supplied, then rmse_tol = abs_tol/norm.ppf(1-alpha/2)
rmse_tol (float) – root mean squared error If supplied (not None), then absolute tolerance and alpha are ignored in favor of the rmse tolerance
n_max (int) – maximum number of samples
replications (int) – number of replications on each level
-
integrate
()¶ See abstract method.
-
set_tolerance
(abs_tol=None, alpha=0.01, rmse_tol=None)¶ See abstract method.
- Parameters
integrand (Integrand) – integrand with multi-level g method
abs_tol (float) – absolute tolerance. Reset if supplied, ignored if not.
alpha (float) – uncertaintly level. If rmse_tol not supplied, then rmse_tol = abs_tol/norm.ppf(1-alpha/2)
rel_tol (float) – relative tolerance. Reset if supplied, ignored if not. Takes priority over aboluste tolerance and alpha if supplied.
-
CLT qMC Cubature (with Replications)¶
-
class
qmcpy.stopping_criterion.cub_qmc_clt.
CubQMCCLT
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=256.0, n_max=1073741824, inflate=1.2, alpha=0.01, replications=16.0)¶ Stopping criterion based on Central Limit Theorem for multiple replications.
>>> k = Keister(Lattice(seed=7)) >>> sc = CubQMCCLT(k,abs_tol=.05) >>> solution,data = sc.integrate() >>> solution 1.380... >>> data Solution: 1.3807 Keister (Integrand Object) Lattice (DiscreteDistribution Object) d 1 randomize 1 order natural seed 1093 mimics StdUniform Lebesgue (TrueMeasure Object) transform Gaussian (TrueMeasure Object) mean 0 covariance 2^(-1) decomp_type pca CubQMCCLT (StoppingCriterion Object) inflate 1.200 alpha 0.010 abs_tol 0.050 rel_tol 0 n_init 2^(8) n_max 2^(30) MeanVarDataRep (AccumulateData Object) replications 2^(4) solution 1.381 sighat 4.61e-04 n_total 2^(12) error_bound 3.56e-04 confid_int [1.38 1.381] time_integrate ...
-
__init__
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=256.0, n_max=1073741824, inflate=1.2, alpha=0.01, replications=16.0)¶ - Parameters
integrand (Integrand) – an instance of Integrand
inflate (float) – inflation factor when estimating variance
alpha (float) – significance level for confidence interval
abs_tol (float) – absolute error tolerance
rel_tol (float) – relative error tolerance
n_max (int) – maximum number of samples
replications (int) – number of replications
-
integrate
()¶ See abstract method.
-
Multilevel MC Cubature¶
-
class
qmcpy.stopping_criterion.cub_mc_ml.
CubMCML
(integrand, abs_tol=0.05, alpha=0.01, rmse_tol=None, n_init=256.0, n_max=10000000000.0, levels_min=2.0, levels_max=10.0, alpha0=- 1.0, beta0=- 1.0, gamma0=- 1.0)¶ Stopping criterion based on multi-level monte carlo.
>>> mlco = MLCallOptions(IIDStdUniform(seed=7)) >>> sc = CubMCML(mlco,abs_tol=.05) >>> solution,data = sc.integrate() >>> solution 10.440... >>> data Solution: 10.4406 MLCallOptions (Integrand Object) option european sigma 0.200 k 100 r 0.050 t 1 b 85 IIDStdUniform (DiscreteDistribution Object) d 2^(6) seed 7 mimics StdUniform Gaussian (TrueMeasure Object) mean 0 covariance 1 decomp_type pca CubMCML (StoppingCriterion Object) rmse_tol 0.019 n_init 2^(8) levels_min 2^(1) levels_max 10 theta 2^(-2) MLMCData (AccumulateData Object) levels 7 dimensions [ 1. 2. 4. 8. 16. 32. 64.] n_level [7.795e+05 1.496e+04 5.916e+03 2.244e+03 7.560e+02 2.770e+02 1.060e+02] mean_level [1.005e+01 1.777e-01 1.071e-01 5.340e-02 2.990e-02 1.167e-02 7.812e-03] var_level [1.959e+02 1.441e-01 4.433e-02 1.093e-02 2.940e-03 7.304e-04 2.261e-04] cost_per_sample [ 1. 2. 4. 8. 16. 32. 64.] n_total 804118 alpha 0.942 beta 1.893 gamma 1.000 time_integrate ...
Original Implementation:
References
[1] M.B. Giles. ‘Multi-level Monte Carlo path simulation’. Operations Research, 56(3):607-617, 2008. http://people.maths.ox.ac.uk/~gilesm/files/OPRE_2008.pdf.
-
__init__
(integrand, abs_tol=0.05, alpha=0.01, rmse_tol=None, n_init=256.0, n_max=10000000000.0, levels_min=2.0, levels_max=10.0, alpha0=- 1.0, beta0=- 1.0, gamma0=- 1.0)¶ - Parameters
integrand (Integrand) – integrand with multi-level g method
abs_tol (float) – absolute tolerance. Reset if supplied, ignored if not.
alpha (float) – uncertaintly level. If rmse_tol not supplied, then rmse_tol = abs_tol/norm.ppf(1-alpha/2)
rel_tol (float) – relative tolerance. Reset if supplied, ignored if not. Takes priority over aboluste tolerance and alpha if supplied.
n_init (int) – initial number of samples
n_max (int) – maximum number of samples
levels_min (int) – minimum level of refinement >= 2
levels_max (int) – maximum level of refinement >= Lmin
alpha0 (float) – weak error is O(2^{-alpha0*level})
beta0 (float) – variance is O(2^{-bet0a*level})
gamma0 (float) – sample cost is O(2^{gamma0*level})
Note
if alpha, beta, gamma are not positive, then they will be estimated
-
integrate
()¶ See abstract method.
-
set_tolerance
(abs_tol=None, alpha=0.01, rmse_tol=None)¶ See abstract method.
- Parameters
abs_tol (float) – absolute tolerance. Reset if supplied, ignored if not.
alpha (float) – uncertaintly level. If rmse_tol not supplied, then rmse_tol = abs_tol/norm.ppf(1-alpha/2)
rel_tol (float) – relative tolerance. Reset if supplied, ignored if not. Takes priority over aboluste tolerance and alpha if supplied.
-
Guaranteed MC Cubature¶
-
class
qmcpy.stopping_criterion.cub_mc_g.
CubMCG
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=1024.0, n_max=10000000000.0, inflate=1.2, alpha=0.01)¶ Stopping criterion with guaranteed accuracy.
>>> k = Keister(IIDStdUniform(2,seed=7)) >>> sc = CubMCG(k,abs_tol=.05) >>> solution,data = sc.integrate() >>> solution 1.803... >>> data Solution: 1.8039 Keister (Integrand Object) IIDStdUniform (DiscreteDistribution Object) d 2^(1) seed 7 mimics StdUniform Lebesgue (TrueMeasure Object) transform Gaussian (TrueMeasure Object) mean 0 covariance 2^(-1) decomp_type pca CubMCG (StoppingCriterion Object) inflate 1.200 alpha 0.010 abs_tol 0.050 rel_tol 0 n_init 2^(10) n_max 10000000000 MeanVarData (AccumulateData Object) levels 1 solution 1.804 n 13473 n_total 14497 error_bound 0.050 confid_int [1.754 1.854] time_integrate ...
Original Implementation:
References
[1] Fred J. Hickernell, Lan Jiang, Yuewei Liu, and Art B. Owen, “Guaranteed conservative fixed width confidence intervals via Monte Carlo sampling,” Monte Carlo and Quasi-Monte Carlo Methods 2012 (J. Dick, F. Y. Kuo, G. W. Peters, and I. H. Sloan, eds.), pp. 105-128, Springer-Verlag, Berlin, 2014. DOI: 10.1007/978-3-642-41095-6_5
[2] Sou-Cheng T. Choi, Yuhan Ding, Fred J. Hickernell, Lan Jiang, Lluis Antoni Jimenez Rugama, Da Li, Jagadeeswaran Rathinavel, Xin Tong, Kan Zhang, Yizhi Zhang, and Xuan Zhou, GAIL: Guaranteed Automatic Integration Library (Version 2.3) [MATLAB Software], 2019. Available from http://gailgithub.github.io/GAIL_Dev/
- Guarantee:
This algorithm attempts to calculate the mean, mu, of a random variable to a prescribed error tolerance, _tol_fun:= max(abstol,reltol*|mu|), with guaranteed confidence level 1-alpha. If the algorithm terminates without showing any warning messages and provides an answer tmu, then the follow inequality would be satisfied:
.
-
__init__
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=1024.0, n_max=10000000000.0, inflate=1.2, alpha=0.01)¶ - Parameters
integrand (Integrand) – an instance of Integrand
inflate – inflation factor when estimating variance
alpha – significance level for confidence interval
abs_tol – absolute error tolerance
rel_tol – relative error tolerance
n_init – initial number of samples
n_max – maximum number of samples
-
integrate
()¶ See abstract method.
CLT MC Cubature¶
-
class
qmcpy.stopping_criterion.cub_mc_clt.
CubMCCLT
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=1024.0, n_max=10000000000.0, inflate=1.2, alpha=0.01)¶ Stopping criterion based on the Central Limit Theorem.
>>> k = Keister(IIDStdUniform(2,seed=7)) >>> sc = CubMCCLT(k,abs_tol=.05) >>> solution,data = sc.integrate() >>> solution 1.801... >>> data Solution: 1.8010 Keister (Integrand Object) IIDStdUniform (DiscreteDistribution Object) d 2^(1) seed 7 mimics StdUniform Lebesgue (TrueMeasure Object) transform Gaussian (TrueMeasure Object) mean 0 covariance 2^(-1) decomp_type pca CubMCCLT (StoppingCriterion Object) inflate 1.200 alpha 0.010 abs_tol 0.050 rel_tol 0 n_init 2^(10) n_max 10000000000 MeanVarData (AccumulateData Object) levels 1 solution 1.801 n 5741 n_total 6765 error_bound 0.051 confid_int [1.75 1.852] time_integrate ... >>> ac = AsianOption(IIDStdUniform(), ... multi_level_dimensions = [2,4,8]) >>> sc = CubMCCLT(ac,abs_tol=.05) >>> solution,data = sc.integrate()
-
__init__
(integrand, abs_tol=0.01, rel_tol=0.0, n_init=1024.0, n_max=10000000000.0, inflate=1.2, alpha=0.01)¶ - Parameters
-
integrate
()¶ See abstract method.
-
Accumulate Data Class¶
Abstract Accumulate Data Class¶
LD Sequence Transform Data (qMC)¶
-
class
qmcpy.accumulate_data.ld_transform_data.
LDTransformData
(stopping_crit, integrand, true_measure, discrete_distrib, basis_transform, m_min, m_max, fudge, check_cone, ptransform)¶ Update and store transformation data based on low-discrepancy sequences. See the stopping criterion that utilize this object for references.
-
__init__
(stopping_crit, integrand, true_measure, discrete_distrib, basis_transform, m_min, m_max, fudge, check_cone, ptransform)¶ - Parameters
stopping_crit (StoppingCriterion) – a StoppingCriterion instance
integrand (Integrand) – an Integrand instance
true_measure (TrueMeasure) – A TrueMeasure instance
discrete_distrib (DiscreteDistribution) – a DiscreteDistribution instance
basis_transform (method) – Transform ynext, combine with y, and then transform all points. For cub_lattice this is Fast Fourier Transform (FFT). For cub_sobol this is Fast Walsh Transform (FWT)
m_min (int) – initial n == 2^m_min
m_max (int) – max n == 2^m_max
fudge (function) – positive function multiplying the finite sum of basis coefficients specified in the cone of functions
check_cone (boolean) – check if the function falls in the cone
-
update_data
()¶ See abstract method.
-
Mean Variance qMC Data (for Replications)¶
-
class
qmcpy.accumulate_data.mean_var_data_rep.
MeanVarDataRep
(stopping_crit, integrand, true_measure, discrete_distrib, n_init, replications)¶ Update and store mean and variance estimates with repliations. See the stopping criterion that utilize this object for references.
-
__init__
(stopping_crit, integrand, true_measure, discrete_distrib, n_init, replications)¶ - Parameters
stopping_crit (StoppingCriterion) – a StoppingCriterion instance
integrand (Integrand) – an Integrand instance
true_measure (TrueMeasure) – A TrueMeasure instance
discrete_distrib (DiscreteDistribution) – a DiscreteDistribution instance
n_init (int) – initial number of samples
replications (int) – number of replications
-
update_data
()¶ See abstract method.
-
Multilevel qMC Data¶
-
class
qmcpy.accumulate_data.mlqmc_data.
MLQMCData
(stopping_crit, integrand, true_measure, discrete_distrib, n_init, replications)¶ Accumulated data for quasi-random sequence calculations, and store multi-level, multi-replications mean, variance, and cost values. See the stopping criterion that utilize this object for references.
-
__init__
(stopping_crit, integrand, true_measure, discrete_distrib, n_init, replications)¶ Initialize data instance
- Parameters
stopping_crit (StoppingCriterion) – a StoppingCriterion instance
integrand (Integrand) – an Integrand instance
true_measure (TrueMeasure) – A TrueMeasure instance
discrete_distrib (DiscreteDistribution) – a DiscreteDistribution instance
replications (int) – number of replications on each level
-
update_data
()¶ See abstract method.
-
Multilevel MC Data¶
-
class
qmcpy.accumulate_data.mlmc_data.
MLMCData
(stopping_crit, integrand, true_measure, discrete_distrib, levels_init, n_init, alpha0, beta0, gamma0)¶ Accumulated data for IIDDistribution calculations, and store multi-level mean, variance, and cost values. See the stopping criterion that utilize this object for references.
-
__init__
(stopping_crit, integrand, true_measure, discrete_distrib, levels_init, n_init, alpha0, beta0, gamma0)¶ Initialize data instance
- Parameters
stopping_crit (StoppingCriterion) – a StoppingCriterion instance
integrand (Integrand) – an Integrand instance
true_measure (TrueMeasure) – A TrueMeasure instance
discrete_distrib (DiscreteDistribution) – a DiscreteDistribution instance
levels_init (int) – initial number of levels
n_init (int) – initial number of samples per level
alpha0 (float) – weak error is O(2^{-alpha0*level})
beta0 (float) – variance is O(2^{-beta0*level})
gamma0 (float) – sample cost is O(2^{gamma0*level})
-
update_data
()¶ See abstract method.
-
Mean Variance MC Data¶
-
class
qmcpy.accumulate_data.mean_var_data.
MeanVarData
(stopping_crit, integrand, true_measure, discrete_distrib, n_init)¶ Update and store mean and variance estimates.
-
__init__
(stopping_crit, integrand, true_measure, discrete_distrib, n_init)¶ - Parameters
stopping_crit (StoppingCriterion) – a StoppingCriterion instance
integrand (Integrand) – an Integrand instance
true_measure (TrueMeasure) – A TrueMeasure instance
discrete_distrib (DiscreteDistribution) – a DiscreteDistribution instance
n_init (int) – initial number of samples
-
update_data
()¶ See abstract method.
-
Utilities¶
-
qmcpy.util.latnetbuilder_linker.
latnetbuilder_linker
(lnb_dir='./', out_dir='./', fout_prefix='lnb4qmcpy')¶ - Parameters
- Returns
- path to file which can be passed into QMCPy’s Lattice or Sobol’ in order to use
the linked latnetbuilder generating vector/matrix e.g. ‘my_poly_lat_vec.10.16.npy’
- Return type
- Adapted from latnetbuilder parser: