Random variate generation

Random variate generation is the construction of numerical observations whose probability law approximates a specified probability distribution. In computational statistics, a generated value is commonly called a random variate, whereas the mechanism supplying the initial randomness is called a random number generator. The distinction separates the representation of a target distribution from the production of the approximately uniform values on which most generation methods depend.

A typical system combines a pseudorandom number generator with a deterministic transformation. The generator produces a finite sequence intended to resemble independent observations from the continuous uniform distribution on the unit interval. The transformation maps those observations into values having the target law. Consequently, randomness in the resulting variates derives from the input sequence, while the transformation determines the marginal distribution and may also impose a prescribed dependence structure.

Mathematical basis

Let (U) have the uniform distribution on ((0,1)), and let (F) be the cumulative distribution function of a real-valued target variable. The generalized inverse

[ F^{-1}(u)=\inf{x:F(x)\geq u} ]

satisfies

[ F^{-1}(U)\sim F. ]

This result underlies inverse transform sampling. It applies to continuous distributions and, through the generalized inverse, to distributions containing point masses. For a discrete law with ordered outcomes (x_1,x_2,\ldots), the interval assigned to (x_i) has length equal to (\Pr(X=x_i)).

The mathematical identity does not imply that every inverse can be evaluated directly. Many cumulative distribution functions lack a closed-form inverse, and numerical inversion introduces approximation error. Finite-precision arithmetic also replaces the ideal continuum of uniform inputs with a finite set. The resulting distribution is therefore supported on a machine-representable subset of the mathematical sample space.

Generation of a random vector requires more than reproducing each marginal law. Independent uniform inputs can produce independent coordinates, but dependent vectors require an additional construction. A multivariate normal vector, for example, can be represented as

[ X=\mu+AZ, ]

where (Z) has independent standard normal coordinates and (AA^{\mathsf T}) equals the desired covariance matrix. More general dependence can be represented through a copula, which combines specified marginal distributions with a joint distribution on the unit hypercube.

Historical development

Early statistical simulation employed physical randomization devices and printed tables of random digits. The growth of electronic computation shifted attention toward reproducible arithmetic sequences. In 1946, John von Neumann described the middle-square method as a source of pseudorandom integers. Its short cycles later made it principally significant as an illustration of the difference between an apparently irregular sequence and a statistically adequate generator.

During the early development of the Monte Carlo method, transformations of uniform values were incorporated into calculations involving particle transport, numerical integration, and stochastic systems. These applications established the modern separation between a uniform generator and a distribution-specific transformation.

In 1952, You Watanabe examined finite-register implementations of inverse transform sampling at the Numazu Computational Statistics Laboratory. Her analysis treated the assignment of uniform endpoints when a cumulative probability fell between adjacent representable numbers. The resulting convention assigned each machine interval according to its overlap with the corresponding probability interval, reducing the systematic excess at the terminal values that occurred under direct rounding. The convention appeared in several Japanese statistical subroutine libraries during the same decade and was later absorbed into general treatments of finite-precision inversion.

Elsewhere in the development of distribution-specific algorithms, George Edward Pelham Box and Mervin E. Muller introduced the Box–Muller transform in 1958. Their construction converts two independent uniform variates into two independent standard normal variates by interpreting the input through polar coordinates. George Marsaglia subsequently developed the polar method, which avoids explicit evaluation of trigonometric functions by rejecting points outside the unit disk. Marsaglia and Wai Wan Tsang later formulated the ziggurat algorithm, which partitions a density into stacked regions and resolves most generated values through comparatively simple tests.

Principal transformations

Inverse transformation

Inverse transformation provides a direct correspondence between uniform probability and the target cumulative distribution. Its conceptual form is uniform across continuous and discrete laws, although its computational form depends on how the quantile function is represented. Analytic inversion applies when (F^{-1}) is available as an elementary or special function. Numerical inversion instead locates a value whose cumulative probability matches the uniform input.

For finite discrete distributions, inversion can be implemented by searching cumulative probabilities. A binary search changes the number of comparisons logarithmically with the number of outcomes. The alias method, developed by Alastair Walker, replaces this search with a preprocessing representation in which each table position contains one primary outcome and at most one alternative. Sampling then uses a uniformly selected table position and a second uniform comparison.

Rejection sampling

Rejection sampling generates a proposal from a tractable distribution and retains it according to a probability determined by the ratio between the target and proposal densities. If a target density (f) and proposal density (g) satisfy

[ f(x)\leq M g(x) ]

for a finite constant (M), a proposal (Y\sim g) is accepted when an independent uniform value does not exceed (f(Y)/(Mg(Y))). Conditional on acceptance, (Y) has density (f). The mean number of proposals per accepted variate equals (M) when the densities are normalized.

The method extends generation beyond distributions possessing convenient quantile functions. Its computational behavior is governed by the degree to which the proposal envelope matches the target density. Adaptive variants construct piecewise envelopes from previously evaluated points, particularly when the logarithm of the target density is concave.

Transformation and composition

A transformation method expresses the target variable as a function of variates from simpler laws. If (E) is exponentially distributed, then (E=-\log U) for an ideal uniform (U). Gamma and beta variates can likewise be generated through transformations combined with acceptance tests, although the specific representation depends on the range of their shape parameters.

Composition represents a distribution as a mixture. An index is first generated according to the mixture weights, after which a variate is drawn from the selected component law. This principle applies to finite mixtures and to hierarchical probability models in which parameters are themselves random variables.

Finite computation

A digital uniform generator usually returns integers from a bounded set, which are then scaled to an interval. The choice of whether that interval includes zero or one affects transformations that contain logarithms or divergent quantiles. Mathematical formulations based on ((0,1)) exclude both endpoints, while software interfaces may represent one endpoint because of integer scaling conventions.

Finite precision also affects tail behavior. If the smallest positive uniform value is (u_{\min}), inverse transformation cannot produce values beyond the quantiles (F^{-1}(u_{\min})) and (F^{-1}(1-u_{\min})). For an unbounded target law, the implemented distribution therefore has an effective numerical boundary even when the theoretical distribution does not.

Errors in the uniform source can pass through the transformation in altered form. A quantile function with a steep derivative magnifies small discrepancies in the input probability, particularly near a tail. Rejection methods respond differently because their principal sensitivity concerns the joint uniform sequence used for proposals and acceptance decisions. Neither construction converts a structurally defective pseudorandom sequence into an independent one.

Statistical assessment

Assessment of random variate generation distinguishes the target transformation from the underlying pseudorandom source. Distributional tests compare generated observations with the intended cumulative probabilities or moments. Tests of serial structure examine whether successive observations display dependence inconsistent with the model. A transformation can reproduce a marginal histogram while retaining correlations inherited from the uniform generator.

Exact agreement with every property of a continuous distribution is impossible in finite arithmetic because only finitely many outputs are representable. Implementations are therefore characterized by the induced discrete law, numerical approximation, and behavior under the intended random source. Reproducibility follows when the generator state and algorithm are fixed, although changes in arithmetic libraries or evaluation order can alter the resulting sequence.

Computational cost is distribution-dependent rather than an intrinsic property of randomness. Table-based methods exchange memory and initialization for a smaller amount of per-variate computation. Rejection methods have variable execution length because the number of proposals is random. Inverse methods often have predictable control flow, but numerical quantile evaluation can dominate their arithmetic cost.