Package 'qrng' reference manual

Title:	(Randomized) Quasi-Random Number Generators
Description:	Functionality for generating (randomized) quasi-random numbers in high dimensions.
Authors:	Marius Hofert [aut, cre], Christiane Lemieux [aut]
Maintainer:	Marius Hofert <[email protected]>
License:	GPL-2 \| GPL-3
Version:	0.0-10
Built:	2024-09-20 05:50:33 UTC
Source:	CRAN

Compute Quasi-Random Sequences

Description

Computing Korobov, generalize Halton and Sobol' quasi-random sequences.

Usage

korobov(n, d = 1, generator, randomize = c("none", "shift"))
ghalton(n, d = 1, method = c("generalized", "halton"))
sobol  (n, d = 1, randomize = c("none", "digital.shift", "Owen", "Faure.Tezuka",
                                "Owen.Faure.Tezuka"), seed, skip = 0, ...)
korobov(n, d = 1, generator, randomize = c("none", "shift"))
ghalton(n, d = 1, method = c("generalized", "halton"))
sobol  (n, d = 1, randomize = c("none", "digital.shift", "Owen", "Faure.Tezuka",
                                "Owen.Faure.Tezuka"), seed, skip = 0, ...)

Arguments

`n`	number $n$ of points to be generated $\ge 2$ .
`d`	dimension $d$ .
`generator`	`numeric` of length $d$ or length 1 (in which case it is appropriately extended to length $d$ ). All numbers must be in $\{1,\dots,n\}$ and must be (coercible to) integers.
`randomize`	`character` string indicating the type of randomization for the point set. `korobov()` one of the following: "none" no randomization. "shift" a uniform random variate modulo 1. `sobol()` one of the following: `"none"` no randomization. "digital.shift" a digital shift. "Owen" calls `generate_sobol_owen_set()` from package spacefillr. "Faure.Tezuka","Owen.Faure.Tezuka" calls `sobol()` from package randtoolbox with `scrambling` being 2 and 3, respectively. If `randomize` is a `logical`, then it is interpreted as `"none"` if `FALSE` and `"digital.shift"` if `TRUE`.
`method`	`character` string indicating which sequence is generated, generalized Halton or (plain) Halton.
`seed`	if provided, an integer used within `set.seed()` for the non-scrambling `randomize` methods (so `"none"` or `"digital.shift"`) or passed to the underlying `generate_sobol_owen_set()` from package spacefillr (for `"Owen"`) and `sobol()` from package randtoolbox for the scrambling methods. If not provided, the non-scrambling methods respect a global `set.seed()` but the scrambling methods do not (they then use a randomly generated one based on time and process identifier).
`skip`	number of initial points in the sequence to be skipped (`skip = 0` means the sequence starts with at the origin). Note that for the scrambling methods this simply computes `n + skip` points and omits the first `skip`-many.
`...`	additional arguments passed to `sobol()` from package randtoolbox for randomization methods `"Faure.Tezuka"` and `"Owen.Faure.Tezuka"`.

Details

For sobol() examples see demo(sobol_examples).

Note that these procedures call fast C code. The following restrictions apply:

korobov(): n,d must be $\le 2^{31}-1$ .
ghalton(): n must be $\le 2^{32}-1$ and d must be $\le 360$ .
sobol(): if randomize = "none" or randomize = "digital.shift", n must be $\le 2^{31}-1$ and d must be $\le 16510$ .

The choice of parameters for korobov() is crucial for the quality of this quasi-random sequence (only basic sanity checks are conducted). For more details, see l'Ecuyer and Lemieux (2000).

The generalized Halton sequence uses the scrambling factors of Faure and Lemieux (2009).

Value

korobov() and ghalton() return an $(n,d)$ -matrix; for $d=1$ an $n$ -vector is returned.

Author(s)

Marius Hofert and Christiane Lemieux

References

Faure, H., Lemieux, C. (2009). Generalized Halton Sequences in 2008: A Comparative Study. ACM-TOMACS 19(4), Article 15.

l'Ecuyer, P., Lemieux, C. (2000). Variance Reduction via Lattice Rules. Stochastic Models and Simulation, 1214–1235.

Lemieux, C., Cieslak, M., Luttmer, K. (2004). RandQMC User's guide. See https://www.math.uwaterloo.ca/~clemieux/randqmc/guide.pdf

Examples

n <- 1021 # prime
d <- 4 # dimension

## Korobov's sequence
generator <- 76 # see l'Ecuyer and Lemieux
u <- korobov(n, d = d, generator = generator)
pairs(u, gap = 0, pch = ".", labels = as.expression(
      sapply(1:d, function(j) bquote(italic(u[.(j)])))))

## Randomized Korobov's sequence
set.seed(271)
u <- korobov(n, d = d, generator = generator, randomize = "shift")
pairs(u, gap = 0, pch = ".", labels = as.expression(
      sapply(1:d, function(j) bquote(italic(u[.(j)])))))

## Generalized Halton sequence (randomized by definition)
set.seed(271)
u <- ghalton(n, d)
pairs(u, gap = 0, pch = ".", labels = as.expression(
      sapply(1:d, function(j) bquote(italic(u[.(j)])))))

## For sobol() examples, see demo(sobol_examples)
n <- 1021 # prime
d <- 4 # dimension

## Korobov's sequence
generator <- 76 # see l'Ecuyer and Lemieux
u <- korobov(n, d = d, generator = generator)
pairs(u, gap = 0, pch = ".", labels = as.expression(
      sapply(1:d, function(j) bquote(italic(u[.(j)])))))

## Randomized Korobov's sequence
set.seed(271)
u <- korobov(n, d = d, generator = generator, randomize = "shift")
pairs(u, gap = 0, pch = ".", labels = as.expression(
      sapply(1:d, function(j) bquote(italic(u[.(j)])))))

## Generalized Halton sequence (randomized by definition)
set.seed(271)
u <- ghalton(n, d)
pairs(u, gap = 0, pch = ".", labels = as.expression(
      sapply(1:d, function(j) bquote(italic(u[.(j)])))))

## For sobol() examples, see demo(sobol_examples)

Test Functions

Description

Functions for testing low-discrepancy sequences.

Usage

sum_of_squares(u)
sobol_g(u, copula = copula::indepCopula(dim = ncol(u)), alpha = 1:ncol(u), ...)
exceedance(x, q, p = 0.99, method = c("indicator", "individual.given.sum.exceeds",
                                      "sum.given.sum.exceeds"))
sum_of_squares(u)
sobol_g(u, copula = copula::indepCopula(dim = ncol(u)), alpha = 1:ncol(u), ...)
exceedance(x, q, p = 0.99, method = c("indicator", "individual.given.sum.exceeds",
                                      "sum.given.sum.exceeds"))

Arguments

`u`	$(n, d)$ -matrix containing $n$ $d$ -dimensional realizations (of a potential quasi-random number generator). For `sum_of_squares()` these need to be marginally standard uniform and for `sobol_g()` they need to follow the copula specified by `copula`.
`copula`	`Copula` object for which the inverse Rosenblatt transformation exists.
`alpha`	vector of parameters of Sobol's g test function.
`...`	additional arguments passed to the underlying `cCopula()`.
`x`	$(n, d)$ -matrix containing $n$ $d$ -dimensional realizations.
`q`	"indicator" $d$ -vector containing the componentwise thresholds; if a number it is recycled to a $d$ -vector. "individual.given.sum.exceeds", "sum.given.sum.exceeds" threshold for the sum (row sums of `x`).
`p`	If `q` is not provided, the probability `p` is used to determine `q`. "indicator" $d$ -vector containing the probabilities determining componentwise thresholds via empirical quantiles; if a number, it is recycled to a $d$ -vector. "individual.given.sum.exceeds", "sum.given.sum.exceeds" probability determining the threshold for the sum (row sums of `x`) via the corresponding empirical quantile.
`method`	`character` string indicating the type of exceedance computed (see Section Value below).

Details

For examples see the demo man_test_functions.

See ES_np(<matrix>) from qrmtools for another test function.

Value

sum_of_squares() returns an $n$ -vector (numeric(n)) with the rowwise computed scaled sum of squares (theoretically integrating to 1).

sobol_g() returns an $n$ -vector (numeric(n)) with the rowwise computed Sobol' g functions.

exceedance()'s return value depends on method:

"indicator": returns indicators whether, componentwise, x exceeds the threshold determined by q.
"individual.given.sum.exceeds": returns all rows of x whose sum exceeds the threshold determined by q.
"sum.given.sum.exceeds": returns the row sums of those rows of x whose sum exceeds the threshold determined by q.

Author(s)

Marius Hofert and Christiane Lemieux

References

Radovic, I., Sobol', I. M. and Tichy, R. F. (1996). Quasi-Monte Carlo methods for numerical integration: Comparison of different low discrepancy sequences. Monte Carlo Methods and Applications 2(1), 1–14.

Faure, H., Lemieux, C. (2009). Generalized Halton Sequences in 2008: A Comparative Study. ACM-TOMACS 19(4), Article 15.

Owen, A. B. (2003). The dimension distribution and quadrature test functions. Stat. Sinica 13, 1-–17.

Sobol', I. M. and Asotsky, D. I. (2003). One more experiment on estimating high-dimensional integrals by quasi-Monte Carlo methods. Math. Comput. Simul. 62, 255–-263.

Examples

## Generate some (here: copula, pseudo-random) data
library(copula)
set.seed(271)
cop <- claytonCopula(iTau(claytonCopula(), tau = 0.5)) # Clayton copula
U <- rCopula(1000, copula = cop)

## Compute sum of squares test function
mean(sum_of_squares(U)) # estimate of E(3(sum_{j=1}^d U_j^2)/d)

## Compute the Sobol' g test function
if(packageVersion("copula") >= "0.999-20")
    mean(sobol_g(U)) # estimate of E(<Sobol's g function>)

## Compute an exceedance probability
X <- qnorm(U)
mean(exceedance(X, q = qnorm(0.99))) # fixed threshold q
mean(exceedance(X, p = 0.99)) # empirically estimated marginal p-quantiles as thresholds

## Compute 99% expected shortfall for the sum
mean(exceedance(X, p = 0.99, method = "sum.given.sum.exceeds"))
## Or use ES_np(X, level = 0.99) from 'qrmtools'
## Generate some (here: copula, pseudo-random) data
library(copula)
set.seed(271)
cop <- claytonCopula(iTau(claytonCopula(), tau = 0.5)) # Clayton copula
U <- rCopula(1000, copula = cop)

## Compute sum of squares test function
mean(sum_of_squares(U)) # estimate of E(3(sum_{j=1}^d U_j^2)/d)

## Compute the Sobol' g test function
if(packageVersion("copula") >= "0.999-20")
    mean(sobol_g(U)) # estimate of E(<Sobol's g function>)

## Compute an exceedance probability
X <- qnorm(U)
mean(exceedance(X, q = qnorm(0.99))) # fixed threshold q
mean(exceedance(X, p = 0.99)) # empirically estimated marginal p-quantiles as thresholds

## Compute 99% expected shortfall for the sum
mean(exceedance(X, p = 0.99, method = "sum.given.sum.exceeds"))
## Or use ES_np(X, level = 0.99) from 'qrmtools'

Compute Matrices to Arrays

Description

Converting higher-dimensional matrices of quasi-random numbers to arrays of specific formats.

Usage

to_array(x, f, format = c("(n*f,d)", "(n,f,d)"))
to_array(x, f, format = c("(n*f,d)", "(n,f,d)"))

Arguments

`x`	( $n$ , $f$ $d$ )-matrix of quasi-random numbers to be converted.
`f`	factor $f\ge 1$ dividing `ncol{x}`.
`format`	`character` string indicating the output format to which `x` should be converted.

Details

to_array() is helpful for converting quasi-random numbers to time series paths.

Value

(n * f, d)-matrix or (n, f, d)-array depending on the chosen format.

Author(s)

Marius Hofert

Examples

## Basic call
N <- 4 # replications
n <- 3 # time steps
d <- 2 # dimension
set.seed(271) # note: respected for the choice of 'randomize'
x <- sobol(N, d = n * d, randomize = "digital.shift") # higher-dim. Sobol'
stopifnot(dim(to_array(x, f = n)) == c(N * n, d)) # conversion and check
stopifnot(dim(to_array(x, f = n, format = "(n,f,d)")) == c(N, n, d))

## See how the conversion is done
(x <- matrix(1:(N * n * d), nrow = N, byrow = TRUE))
to_array(x, f = n) # => (n * d)-column x was blocked in n groups of size d each
## Basic call
N <- 4 # replications
n <- 3 # time steps
d <- 2 # dimension
set.seed(271) # note: respected for the choice of 'randomize'
x <- sobol(N, d = n * d, randomize = "digital.shift") # higher-dim. Sobol'
stopifnot(dim(to_array(x, f = n)) == c(N * n, d)) # conversion and check
stopifnot(dim(to_array(x, f = n, format = "(n,f,d)")) == c(N, n, d))

## See how the conversion is done
(x <- matrix(1:(N * n * d), nrow = N, byrow = TRUE))
to_array(x, f = n) # => (n * d)-column x was blocked in n groups of size d each

Package 'qrng'

Help Index

Compute Quasi-Random Sequences

Description

Usage

Arguments

Details

Value

Author(s)

References

Examples

Test Functions

Description

Usage

Arguments

Details

Value

Author(s)

References

Examples

Compute Matrices to Arrays

Description

Usage

Arguments

Details

Value

Author(s)

See Also

Examples