July 16, 2026
We prove fundamental space lower bounds for exact random sampling using an entropy source of i.i.d. uniform bits. A classic result from information theory shows that generating \(n\) discrete random variables \(X_1, \dots, X_n\) requires at least \(H(X_1, \dots, X_n)\) input random bits on average, where \(H\) is the [1] entropy function. How much space must a random sampling algorithm use in order to approach this information-theoretically optimal entropy bound?
We prove that any random sampling algorithm that is exact for arbitrary discrete target distributions and consumes at most \(H(X_1,\ldots,X_n)+\varepsilon n+o(n)\) input bits in expectation for every output process must use \(\Omega(\log(1/\varepsilon))\) bits of space. In fact, i.i.d. sampling from the single distribution \(\mathrm{Bernoulli}(1/3)\) already forces at least \((1/{5.116201}-o(1))\log(1/\varepsilon)\) bits of space. If the sampler handles a family of infinitely many Bernoulli distributions, we show a sharper bound of at least \(\log(1/\varepsilon)\) bits of space. We also prove lower bounds for general i.i.d. sampling: for almost every distribution on \(k\) outcomes, the space is at least \((1/(k+1)-o(1))\log(1/\varepsilon)\) bits.
The proof technique is based on a graph-theoretic analysis of the amount of information that any algorithm can store in its state. Finite state spaces force short cycles around the state-transition graph, and the loss around such cycles reduces to Diophantine lower bounds on fractional parts of integer combinations of log-probabilities. To the best of our knowledge, these results comprise the first known space lower bounds for entropy-efficient random sampling.
We are concerned with the problem of generating a sequence \(X_1, X_2, \dots\) of discrete random variables using an entropy source that emits independent and identically distributed random bits \(U_1, U_2, \dots\). This problem is known in the literature as random variate generation [2] or random sampling, and is a central object of study in theoretical and applied computer science. In online random sampling, the target distribution for each \(X_i\) is dynamically presented at the \(i\)th iteration of sampling, as shown in the following interaction pattern.
Input: Sequence \(U_1, U_2, \dots \overset{\rm iid}{\sim} \mathrm{Uniform}\set{0,1}\) For \(i=1,2,\dots\) Receive next target distribution \(p_i\) Emit next random variable \(X_i \sim p_i\)
Online random sampling arises whenever random variables must be generated sequentially from distributions that are revealed or updated during execution, as each \(p_i\) may itself depend on the previously generated variables or auxiliary sources of randomness. Examples include autoregressive sequence generation in large language models, where the next-token distribution depends on the previously sampled prefix; particle filtering, where resampling distributions are updated after each observation; bandit and reinforcement-learning algorithms, where action distributions depend on the history of actions and rewards; and discrete-event simulations, where transition probabilities depend on the current system state. In all these settings, the random sampler receives target distributions \(p_1, p_2, \dots\) that are presented adaptively over time.
The complexity of an online random sampling algorithm is measured according to three resources that scale with the iteration \(n\): its space usage \(S(n)\), runtime \(R(n)\), and entropy cost \(T(n)\), which is the average number of consumed random input bits from the entropy source. A classic result from [3] gives an “entropy-optimal” sampler whose entropy cost is tightly concentrated around the [1] entropy: \[\begin{align} H(X_1, \dots, X_n) \le T(n) < H(X_1, \dots, X_n) + 2, \end{align}\] where \(H(Z) \mathrel{\vcenter{:}}=-\sum_{z}\mathbb{P}(Z=z) \log(\mathbb{P}(Z=z))\) for any discrete random element \(Z\), and all logarithms are taken with base two. Algorithms such as DDG sampling [3] and the interval algorithm [4] achieve an asymptotically optimal entropy cost of \(T(n) \leq H(X_1, \dots, X_n) + O(1)\), but require unbounded space \(S(n) \to \infty\) as \(n \to \infty\).
In recent work, [5] develop a sampling algorithm for the case that each \(p_i\) is a rational distribution with bounded denominator at most \(d\). Its entropy cost satisfies \[\begin{align} T(n) \le H(X_1, \dots, X_n) + \varepsilon n + \log(d/\varepsilon), \end{align}\] using constant space \(S(n) = O(\log(d/\varepsilon))\). They also offer a conjecture for the general case.
Any online random sampling algorithm that generates exact samples from a sequence of arbitrary discrete distributions, within \(\varepsilon > 0\) of the information-theoretically optimal entropy rate, using a stream of i.i.d. random bits as the entropy source, requires \(\Omega(\log(1/\varepsilon))\) bits of space for an auxiliary state that is carried over between rounds.
In this article, we formally prove this conjecture, establishing the first known space lower bound on entropy-efficient online random sampling.
Theorem 1 (store=theorem:space-lb,restate-keys=note=Restated). Any online random sampler that can generate exact samples from a sequence of arbitrary discrete probability distributions using an input sequence of i.i.d. random bits
requires \(\log(\abs{\mathbb{S}}) \geq \Omega(\log(1/\varepsilon))\) bits of persistent space to achieve an entropy cost of \(T(n) \leq H(X_1,\dots,X_n) + \varepsilon n + o(n)\), where \(\mathbb{S}\) is the state space of the algorithm.
In particular, suppose the sampler can generate samples from a family \(\mathbb{D}\) of probability distributions over the natural numbers \(\mathbb{N}\).
If \(\mathrm{Bernoulli}(1/3) \in \mathbb{D}\), then \(\log(\abs{\mathbb{S}}) \geq (1/5.116201 - o(1))\log(1/\varepsilon)\).
If \(\mathrm{Bernoulli}(1/(2^m+1)) \in \mathbb{D}\) for infinitely many \(m \in \mathbb{N}\), then \(\abs{\mathbb{S}} \geq 1/\varepsilon\).
1 is stronger than the original conjecture: it shows that a space lower bound holds even if the only output distribution handled by the sampler is \(\mathrm{Bernoulli}(1/3)\). By slightly enlarging the set of distributions to include Bernoulli distributions with rational probability \(1/(2^m+1)\), we can obtain a tighter space lower bound of \(\ceil{\log(1/\varepsilon)}\) bits.
We also establish an analogous space lower bound for random sampling with i.i.d. outputs, which is the setting studied by [6].
Theorem 2 (store=theorem:space-lb-iid,restate-keys=note=Restated). Let \(k \geq 2\) be an integer. For almost all discrete distributions \(p\) over \(k\) outcomes, any sampler generating i.i.d. samples from \(p\) with an entropy loss bound of \(\varepsilon\) bits per sample requires at least \((1/(k+1)-o(1))\log(1/\varepsilon)\) bits of persistent space.
In the theorem statement, “almost all” refers to Lebesgue-almost every point \(p \in (0,1)^{k-1}\) with \(\sum_{i=1}^{k-1}p_i < 1\), under the standard parameterization \(p \mapsto (p_1, \dots, p_{k-1}, 1-\sum_{i=1}^{k-1}p_i)\) of the \((k-1)\)-dimensional probability simplex. Specializing 2 to \(k=2\), we obtain a lower bound of \((1/3-o(1))\log(1/\varepsilon)\) bits for almost every Bernoulli distribution, improving the constant \(1/5.116201\) in from 1. Moreover, 2 applies to a less stringent setting than from 1: the former need only generate i.i.d. samples from a fixed distribution \(p\) with entropy loss at most \(\varepsilon\), while the latter must handle infinitely many Bernoulli distributions.
The rest of this paper is organized as follows. 2 reviews existing algorithms for random sampling, illustrates the tradeoff between entropy efficiency and space complexity, and provides intuition for our information-theoretic lower bounds. 3 formally defines online random samplers and entropy loss. We then introduce two new concepts used to prove the \(\Omega(\log(1/\varepsilon))\) lower bound on the size of the state space. In 4, we define an abstract “state information content” function for an arbitrary sampler and show how to construct this function based only on the sampler’s transition function. In 5, we prove a very general bound showing the state-space size needed to ensure an entropy loss of at most \(\varepsilon\) bits per sample when generating i.i.d. samples from a fixed distribution. 6 applies this bound to Bernoulli distributions to prove 1. 7 applies the same bound, together with Diophantine approximation on manifolds, to prove 2. 8 discusses related work in random sampling and coding theory. 9 concludes with open questions.
To motivate our lower bounds, we first present a simple example of a random sampler that maintains persistent space. This example demonstrates both the entropy-space tradeoff and the concept of state information content, which is a main tool in our proof of 1.
Consider the problem of converting an i.i.d. sequence \((U_1, U_2, U_3, \dots)\) of uniform bits into an i.i.d. sequence \((X_1, X_2, X_3, \ldots)\) of uniform trits; i.e., \(U_i \sim \mathrm{Uniform}\set{0,1}\) and \(X_i \sim \mathrm{Uniform}\set{0,1,2}\). The standard inversion method [2] generates \(X_1\) by first partitioning the unit interval \(J_1 \mathrel{\vcenter{:}}=[0,1]\) into three equally sized subintervals: \[\begin{align} [0,1] = \underbrace{\left[0,\textstyle\frac{1}{3}\right]}_{J_{10}} \; \bigcup \; \underbrace{\textstyle\left[\frac{1}{3}, \frac{2}{3}\right]}_{J_{11}} \; \bigcup \; \underbrace{\textstyle\left[\frac{2}{3},1\right]}_{J_{12}}. \label{eq:interval-partition-J1} \end{align}\tag{1}\] It then draws a uniform number \(U \sim \mathrm{Uniform}[0,1]\), and returns \(X_1 \in \set{0,1,2}\) such that \(U \in J_{1X_1}\). To implement this algorithm using an entropy source that provides random bits, we recognize that the expansion \(U = (0.U_1 U_2 U_3 \dots )_2\) satisfies \(U_i \sim \mathrm{Uniform}\set{0,1}\). We keep drawing input bits until we can determine whether \(U \in J_{10}\), \(U \in J_{11}\), or \(U \in J_{12}\), and set \(X_1\) accordingly.
The top row of 1 shows an example of this algorithm in action. The random bits \(U_1U_2U_3 = \texttt{100}\) let us determine that \(U \in J_{11}\), so \(X_1 = 1\). Generating an output in this manner costs 3 input bits on average, whereas the output trit has \(\log(3) \approx 1.585\) bits of entropy and an optimal algorithm [3] uses \({\approx}\,2.667\) bits. To generate \(n\) i.i.d. outputs \(X_1, \dots, X_n\), a stateless sampler repeats this process \(n\) times from scratch, consuming a total of \(3n\) random bits on average.
[4] introduce a random sampling method called the interval algorithm that reduces the entropy loss of naïve inversion sampling by maintaining state across iterations. The key idea is shown in the second row of 1. After generating the first output \(X_1=1\), we have determined that \(U \in I_1 \mathrel{\vcenter{:}}=[4/8, 5/8] \subset [1/3,2/3] = J_{11}\). Rather than resample from scratch, which would discard the three bits used so far, we split \(J_{11} \eqqcolon J_2\) into three new subintervals \[\begin{align} \textstyle\left[\frac{1}{3},\frac{2}{3}\right] = \underbrace{\textstyle\left[\frac{3}{9},\frac{4}{9}\right]}_{J_{2,0}} \; \bigcup \; \underbrace{\textstyle\left[\frac{4}{9}, \frac{5}{9}\right]}_{J_{2,1}} \; \bigcup \; \underbrace{\textstyle\left[\frac{5}{9},\frac{6}{9}\right]}_{J_{2,2}}. \label{eq:interval-partition-J2} \end{align}\tag{2}\] We then draw fresh bits until we can determine that \(U \in J_{2,X_2}\). This approach effectively reuses the first three bits \(\texttt{100}\) drawn in the first iteration, by using the fact that \(U \mid U \in [4/8,5/8]\) follows a uniform distribution over \([4/8, 5/8]\). In the example from 1, when \(U_4 = \texttt{1}\) we can determine that \(X_2 = 2\) using just one extra bit, which is impossible using the stateless sampler. More generally, the interval algorithm requires 5 bits on average to produce two outputs \(X_1, X_2\), instead of the 6 bits used by the stateless sampler.
The entropy loss decreases rapidly with \(n\); for \(n=10\) the stateful algorithm uses \(17.8\) bits on average, which is close to the information-theoretically optimum of roughly \(16.4\) bits using the method of [3] and much smaller than the \(30\) bits incurred by the stateless algorithm. However, this improvement in entropy incurs a cost in terms of space complexity. The expected number of bits needed to represent the state \((I_n, J_n)\) of the interval algorithm grows as \(\Theta(n)\), which is unbounded as \(n \to \infty\) (a detailed analysis of the space and time complexity of the [4] interval method is given in [7]).
Many entropy-efficient algorithms in the literature exhibit a similar entropy-space tradeoff to the one described in 2.2. A key idea for our space lower bounds is to analyze the information contained in the state of an online random sampler. Recall that the information content of a random variable \(X\) drawn from a distribution \(p\) is defined as \(\log(1/p(X))\) [1]. For example, each outcome of a fair coin flip has \(\log(1/0.5) = 1\) bit of information, and an outcome with \(p(X) = 1/3\) has \(\log(3) \approx 1.585\) bits of information. We give a brief survey of the information stored in the state of various sampling algorithms, before giving a formal definition in 4.
For the method described in 2.2, given a state \(I = [I_l, I_r] \subseteq [J_l, J_r] = J\), the total number of bits consumed so far is given by \(\log(1/(I_r - I_l))\), and the information content of the outputs so far is \(\log(1/(J_r - J_l))\) bits (for the case of all output trits in 1, \(\log(1/(J_r - J_l))\) is the same as the number of outputs multiplied by \(\log(3)\)). Because the output information cannot exceed the input information at any time, \(\log(1/(I_r - I_l)) - \log(1/(J_r - J_l)) \geq 0\) provides an upper bound on the information content which the algorithm can extract from the state for future samples. Hence, we call \(\log({J_r-J_l}) - \log({I_r-I_l})\) the state information content.
The online random sampling algorithms of [7] maintain a state consisting of a pair of integers \((Z,M)\), where \(M \geq 1\) and \(Z\) is uniformly distributed over \(\set{0,1,\ldots,M-1}\). Therefore, the information content of the state \((Z,M)\) is given by the information content of \(Z\), namely \(\log M\), because each outcome has probability \(1/M\). These uniform states were used by [8] for the special case of generating a sequence of discrete uniform random variables using i.i.d. random bits. [fig:states-rr] shows an example. A uniform integer \(U\) is drawn between 0 and 15. If the pink cell with label 2 is selected, then \(\mathrm{A}\) is returned. The conditional distribution of \(U\) is uniform over \(\set{0,\ldots,4}\); so the uniform state is \((Z\mathrel{\vcenter{:}}= 2, M\mathrel{\vcenter{:}}= 5)\) with \(Z \mid M \sim \mathrm{Uniform}\set{0,\ldots,M-1}\).
Many sampling algorithms explicitly traverse a data structure called the discrete distribution generating (DDG) tree [3], [5], [9]–[11]. The online sequential version of the algorithm maintains a state consisting of a binary tree, together with the current traversal position (which is a leaf node, after producing a sample). As noted by [12], the usable information in such a DDG state is given by the distribution of leaves in the tree that share the same label as the leaf at the current position (and this information can be reduced to an interval state by fixing an ordering among the leaves sharing the same label). Therefore, the state information content can be defined as the negative log-probability of the leaf, under the distribution of leaves that share the same output label. [fig:states-ddg] shows an example. A random leaf is selected with log-probability proportional to its depth. If leaf \((\mathrm{A},1)\) is selected, then \(\mathrm{A}\) is returned. The conditional distribution over the \(\mathrm{A}\) leaves is \(M \mathrel{\vcenter{:}}=\set{0 \mapsto 4/5, 1 \mapsto 1/5}\); so the nonuniform state is \((Z \mathrel{\vcenter{:}}= 1, M)\) with \(Z \mid M \sim M\).
Lastly, consider the problem where all target distributions \(p_1, p_2, \ldots\) are known in advance, and the samples are generated in batches to improve entropy efficiency (e.g., by applying the interval algorithm to the joint distribution \(p_{1:n}(x_1,\dots,x_n) \mathrel{\vcenter{:}}= p_1(x_1)\dots p_n(x_n)\) of a batch of \(n\) samples), while no state is maintained between batches. In the model where outputs can only be emitted one at a time, the sampler generates \(X_1,\ldots,X_n\) but only emits \(X_1\) in the first round; the state stored for the second round is simply the remaining samples \(X_2, \ldots, X_n\), which have information content \(\log(1/p(X_2, \ldots, X_n))\) bits. In general, the state information content of this batched sampler is the sum of the information content of the samples which have been generated but not yet emitted. Although batched samplers with knowledge of future distributions do not satisfy our definition of online random samplers in 1, for the special case of generating i.i.d. samples (i.e., there is only one possible target distribution), our lower bounds apply.
By maintaining state, an online random sampler can store information that has been consumed from the input source but not yet emitted in an output sample. Suppose a sampler is at state \(s\), reads \(r\) input bits, produces output \(x \sim p\), and moves to a new state \(s'\). The information loss of this transition is the difference between the consumed and produced information: \[0 \le \left(r + h(s)\right) - \left(I_p(x) + h(s')\right),\] where \(h(s)\) and \(h(s')\) denote the information content of states \(s\) and \(s'\), and \(I_p(x) \mathrel{\vcenter{:}}=\log(1/p(x))\) is the information content of the output sample. For the loss to be small, say at most \(\eta\), we need \[h(s') \in [h(s)+r-I_p(x)-\eta, h(s)+r-I_p(x)].\] As the number \(r\) of input bits is an integer, the fractional state information content must satisfy \[\modZ{h(s')} \in \set{\modZ{v} \mid v \in [h(s)-I_p(x)-\eta, h(s)-I_p(x)]},\] where \(\modZ{x}\) is the fractional part of a real number \(x\). A transition that achieves a small loss must therefore choose a next state \(s'\) whose fractional information content lies just before the ideal position \(\modZ{h(s)-I_p(x)}\) that would achieve zero loss; the gap between the two becomes permanently lost information. Therefore, to guarantee a small loss, the possible states of the algorithm should have information content whose fractional parts are densely spaced in \([0,1)\).
One strategy is to use states with information content uniformly spaced over the unit interval. 3 illustrates this idea for \(S=5\) states \(s_0,\ldots,s_4\) with fractional information content \(\modZ{h(s_i)} = i/5\) for \(i=0,\dots,4\). By choosing an appropriate new state after each sample, such a sampler could theoretically ensure that fewer than \(\eta \mathrel{\vcenter{:}}= 1/S=0.2\) bits are lost after each sample.
We construct distributions that force the sampler to lose nearly \(1/S\) bits of information on average. Consider a low-entropy \(\mathrm{Bernoulli}(\delta)\) distribution, for some small \(\delta > 0\). When the sampler generates the likely outcome \(X=0\), the output sample has information content \(\alpha \mathrel{\vcenter{:}}=\log(1/(1-\delta)) = \Theta(\delta)\), which is small but positive. So if the sampler retained the same state \(s' = s\), it would need to consume at least one fresh bit to produce the output, losing \(1-\alpha\) bits. Otherwise, the sampler must transition to a new state \(s' \ne s\), and the gap of \(1/S\) between states means that at least \(1/S-\alpha\) bits will be lost. 3 illustrates the \(\alpha\) bits of output information using blue shading, and the \(1/S-\alpha\) bits of information lost over a transition is shown using orange shading, with a label for the \(s_3 \to s_2\) information loss. Thus, sampling the likely outcome \(X=0\) contributes \((1-\delta)(1/S-\alpha)\) to the average information loss, which approaches \(1/S\) as \(\delta \to 0\).
In the more general case where the state information content is not uniformly spaced, we see that the blue intervals remain the same, while the orange intervals may shrink or expand, so certain transitions may become more efficient and others less efficient. However, the overall loss after \(S\) transitions remains \(1-S\alpha\) because the total length of the orange intervals is invariant. Thus, the loss per sample remains \(1/S\), so we need \(S \geq 1/\varepsilon\) to ensure that the loss per sample is at most \(\varepsilon\).
In the remainder of this paper, we formalize this intuition to establish 1, and then generalize it to prove lower bounds for even weaker samplers with i.i.d. outputs to prove 2.
We now formally define online random samplers and information loss.
The base-two logarithm is denoted by \(\log\). The fractional part of a real number \(x\) is denoted by \(\modZ{x} \mathrel{\vcenter{:}}= x - \floor{x}\). The set of finite binary strings is denoted by \(\set{0,1}^*\), and the length of a binary string \(\mathbf{c} \in \set{0,1}^*\) is \(\abs{\mathbf{c}}\). The set of natural numbers is \(\mathbb{N}\mathrel{\vcenter{:}}=\set{0,1,\dots}\). The support of a distribution \(p\) is \(\mathop{\mathrm{supp}}(p) \mathrel{\vcenter{:}}=\set{x \in \mathbb{N}\mid p(x) > 0}\). The greatest lower bound of a set \(S \subset \mathbb{R}\) is \(\inf S\), which is defined to be \(+\infty\) if \(S\) is empty.
Definition 1. Let \(\mathbb{D}\) be a set of discrete probability distributions over \(\mathbb{N}\mathrel{\vcenter{:}}=\set{0,1,\dots}\). Let \(\mathbb{S}\) be a set of states and \(s_0 \in \mathbb{S}\) a designated initial state. An online random sampler* is a partial function \[f : \mathbb{D}\times \mathbb{S}\times \set{0,1}^* \rightharpoonup \mathbb{N}\times \mathbb{S}\] such that the following conditions hold.*
**Prefix Free: For every \(p \in \mathbb{D}\) and \(s \in \mathbb{S}\), if \((p, s, \mathbf{c}) \in \mathop{\mathrm{dom}}(f)\) then \((p, s, \mathbf{c}') \notin \mathop{\mathrm{dom}}(f)\) for every \(\mathbf{c}, \mathbf{c}' \in \set{0,1}^*\) such that \(\mathbf{c}\) is a proper prefix of \(\mathbf{c}'\).
**Exhaustive: For every \(p \in \mathbb{D}\) and \(s \in \mathbb{S}\), we have \(\sum_{\substack{\mathbf{c} \mid (p,s,\mathbf{c}) \in \mathop{\mathrm{dom}}(f)}} 2^{-\abs{\mathbf{c}}} = 1\).
**Exact: For all distribution sequences \(p_1, p_2, \ldots \in \mathbb{D}\) and output sequences \(x_1, x_2, \ldots \in \mathbb{N}\) and for each finite prefix length \(n \in \mathbb{N}\), we have \[\sum_{\substack{(s_1, \mathbf{c}_1),\;\dots,\;(s_n, \mathbf{c}_n) \\ s.t. \forall i \in \set{1, \ldots, n}, \; f(p_i, s_{i-1}, \mathbf{c}_i) = (x_i, s_i)}} \left( 2^{-\sum_{i=1}^n \abs{\mathbf{c}_i}} \right) = \prod_{i=1}^n p_i(x_i). \qedhere \label{eq:online-random-sampling-exact}\qquad{(1)}\]
To avoid states with undefined behavior, we assume without loss of generality that each state \(s\) is reachable from the initial state \(s_0\), i.e., there exist distributions \(p_1, \ldots, p_k \in \mathbb{D}\), outputs \(x_1, \ldots, x_k \in \mathbb{N}\), input bit sequences \(\mathbf{c}_1, \ldots, \mathbf{c}_k \in \set{0,1}^*\), and states \(s_1, \ldots, s_k = s\) such that \(f(p_i, s_{i-1}, \mathbf{c}_i) = (x_i, s_i)\) for each \(i \in \set{1, \ldots, k}\). This condition is equivalent to restricting \(\mathbb{S}\) to be the subset of states reachable from \(s_0\).
1 generalizes the formulation in [7]. As we are not concerned with the computational complexity at sampling time, we allow noncomputable dependence of the sampler on the distributions. With this definition, our lower bounds apply in full generality, with no assumptions on how irrational probabilities are represented.
Definition 2. The information loss* of an online random sampler \(f\) when sampling \(x_1, \ldots, x_n\) from target distributions \(p_1, \ldots, p_n\) using input bit sequences \(\mathbf{c}_1, \ldots, \mathbf{c}_n\) is \[\sum_{i=1}^n \left( \abs{\mathbf{c}_i} - \log\frac{1}{p_i(x_i)} \right).\] The entropy loss of \(f\) when sampling from \(p_1, \ldots, p_n\) is \[\sum_{\substack{(s_1,x_1,\mathbf{c}_1),\;\ldots,\;(s_n,x_n,\mathbf{c}_n) \\ s.t. \forall i \in \set{1, \ldots, n}, \; f(p_i, s_{i-1}, \mathbf{c}_i) = (x_i, s_i)}} \left( \sum_{i=1}^n \abs{\mathbf{c}_i} - \log{\frac{1}{p_i(x_i)}} \right) 2^{-\sum_{i=1}^{n}\abs{\mathbf{c}_i}},\] which can be understood as the expected information loss given uniformly random input bits. The sampler \(f\) is said to achieve an entropy loss bound of \(\varepsilon\) bits per sample if for all \(p_1,p_2,\ldots \in \mathbb{D}\), the entropy loss of \(f\) when sampling from \(p_1, \ldots, p_n\) is at most \(\varepsilon n + o(n)\) as \(n \to \infty\).*
An online random sampler maintains some state between each iteration. We can analyze the information content of the state given a description of how the sampler operates, as in the interval method from 2. However, to prove general lower bounds on the state space, we need a way to analyze the information content of a state without knowing the inner workings of the sampler. The key fact we need about the state information content is that it gives an upper bound on how much information the sampler can extract from the stored state for future samples, in place of drawing fresh random bits from the entropy source. Therefore, we can define the information content of a state for an arbitrary online random sampler as follows.
Definition 3. Let \(f\) be an online random sampler with state space \(\mathbb{S}\) and target distribution set \(\mathbb{D}\). The state information content* is a function \(h : \mathbb{S}\to \mathbb{R}_{\geq 0}\) such that \(h(s_0) = 0\), and for each pair of states \(s, s' \in \mathbb{S}\), \[h(s') + \log(1/p(x)) \leq h(s) + \abs{\mathbf{c}},\] whenever \(f(p, s, \mathbf{c}) = (x, s')\) for some target distribution \(p \in \mathbb{D}\), some \(x\) in the support of \(p\), and some finite sequence of input bits \(\mathbf{c}\).*
This definition is similar to the data processing inequality, because the input state information is independent of the input bits and the output state information is independent of the output sample. The definition implies that any state information content function must satisfy a system of difference constraints, which we can compactly express in an \(\abs{\mathbb{S}} \times \abs{\mathbb{S}}\) matrix: \[\label{eq:difference-constraints} h(s_j) - h(s_i) \leq A_{i, j},\tag{3}\] where \[\label{eq:difference-constraints-matrix} A_{i, j} \mathrel{\vcenter{:}}= \inf \set*{ \abs{\mathbf{c}} - \log\frac{1}{p(x)} \; \Bigg\vert\; p \in \mathbb{D}, x \in \mathbb{N}, \mathbf{c} \in \set{0,1}^*, f(p, s_i, \mathbf{c}) = (x, s_j) }.\tag{4}\]
Informally, \(A_{i,j}\) is the minimal possible marginal information loss incurred when transitioning from state \(s_i\) to \(s_j\), including any extra information put into the state. When \(s_i\) contains more information than \(s_j\), we may have \(A_{i,j} < 0\), but the quantity \(A_{i,j} - (h(s_j) - h(s_i)) \geq 0\) (i.e., slack in 3 ) will represent permanently lost information that is contained in neither the output \(x \sim p\) nor the next state \(s_j\). To find a solution to this system of difference constraints 3 , we can apply classical techniques used in graph algorithms.
Theorem 3 (Feasibility of difference constraints [13]). Let \(G\) be the directed graph with vertex set \(\mathbb{S}\) and edge weights \(A_{i, j}\) between each pair of vertices \(s_i, s_j\). If \(G\) contains no negative-weight cycles, then a feasible solution to the system 3 is given by \[h(s_i) = \delta(s_0, s_i),\] where \(s_0\) is the initial state and \(\delta(s_0, s_i)\) is the length of the shortest path from \(s_0\) to \(s_i\) in \(G\). If \(G\) contains a negative-weight cycle, then there is no feasible solution for the system.
A shortest path from \(s_0\) to each \(s_i\) in \(G\) can be computed using the Bellman-Ford algorithm. We rule out the existence of negative-weight cycles in \(G\) by analyzing the probability of input and output sequences corresponding to such a cycle. A negative-weight cycle would imply that the sampler generates a certain output sequence with greater probability than it correctly should.
Theorem 4 (Nonexistence of negative cycles). Let \(G\) be the directed graph with vertex set \(\mathbb{S}\) and edge weights \(A_{i, j}\) between each pair of vertices \(s_i, s_j\). There are no negative-weight cycles in \(G\).
Proof. Assume, for the sake of contradiction, that there is a negative-weight cycle \((i_0, i_1, \dots, i_k = i_0)\) in \(G\), i.e., \(A_{i_0, i_1} + A_{i_1, i_2} + \cdots + A_{i_{k-1}, i_k} < 0\). Define \[\epsilon \mathrel{\vcenter{:}}=- (A_{i_0, i_1} + A_{i_1, i_2} + \cdots + A_{i_{k-1}, i_k}) / k. \label{eq:epsilon}\tag{5}\] For each \(j \in \set{1, 2, \ldots, k}\), let \(p_j \in \mathbb{D}\) and \(x_j \in \mathbb{N}\) and \(\mathbf{c}_j \in \set{0,1}^*\) satisfy \[f(p_j, s_{i_{j-1}}, \mathbf{c}_j) = (x_j, s_{i_j}), \qquad \abs{\mathbf{c}_j} - \log(1/p_j(x_j)) < A_{i_{j-1}, i_j} + \epsilon. \label{eq:bad-cycle}\tag{6}\] Such values exist by the definition of \(A_{i, j}\) as an infimum in 4 . Because state \(s_{i_0}\) is reachable from \(s_0\) by 1, we may let \((q_1,\ldots,q_l)\), \((y_1,\ldots,y_l)\), \((\mathbf{b}_1,\ldots,\mathbf{b}_l)\), and \((r_0=s_0, r_1, \ldots, r_{l-1}, r_l=s_{i_0})\) be sequences of distributions, samples, input bit sequences, and states, respectively, such that \(f(q_j, r_{j-1}, \mathbf{b}_j) = (y_j, r_j)\) for each \(j \in \set{1,\ldots,l}\), representing one possible sequence of transitions from \(s_0\) to \(s_{i_0}\).
Then, for any \(n\in\mathbb{N}\), when the sampler is run with an input bit stream whose prefix is \(\mathbf{b}_1\cdots\mathbf{b}_l (\mathbf{c}_1\cdots\mathbf{c}_k)^n\) and a target distribution sequence whose prefix is \(q_1 \cdots q_l (p_1 \cdots p_k)^n\), it produces output symbols \((y_1 \cdots y_l) (x_1 \cdots x_k)^n\). Define \[N \mathrel{\vcenter{:}}= 1 + \left\lfloor \frac{\sum_{j=1}^l \abs*{\mathbf{b}_j}}{\sum_{j=1}^k (\log(1/p_j(x_j)) - \abs{\mathbf{c}_j})} \right\rfloor > \frac{\sum_{j=1}^l \abs{\mathbf{b}_j}}{\sum_{j=1}^k (\log(1/p_j(x_j)) - \abs{\mathbf{c}_j})},\] which is positive by 5 6 .
When running the sampler on the target distribution sequence \(q_1 \cdots q_l (p_1 \cdots p_k)^N\) using uniform random bits, the probability of generating the sequence \((y_1 \cdots y_l) (x_1 \cdots x_k)^N\) is at least \[\begin{align} 2^{-\abs*{\mathbf{b}_1\cdots\mathbf{b}_l (\mathbf{c}_1\cdots\mathbf{c}_k)^N}} &= 2^{-\sum_{j=1}^l \abs*{\mathbf{b}_j} - N \sum_{j=1}^k \abs{\mathbf{c}_j}} \\ &> 2^{-N \sum_{j=1}^k \log(1/p_j(x_j))} \\ &= (p_1(x_1)\cdots p_k(x_k))^N \\ &\geq (q_1(y_1)\cdots q_l(y_l)) (p_1(x_1)\cdots p_k(x_k))^N. \end{align}\] But generating the output \((y_1 \cdots y_l) (x_1 \cdots x_k)^N\) from the distributions \(q_1 \cdots q_l (p_1 \cdots p_k)^N\) with probability greater than \((q_1(y_1)\cdots q_l(y_l)) (p_1(x_1)\cdots p_k(x_k))^N\) violates the exactness condition ?? . ◻
One last property we will need about the state information content is that it is nonnegative. This result also verifies that every \(A_{i,j} > -\infty\). The proof is very similar to that of 4.
Lemma 1 (Nonnegativity of state information content). Let \(G\) be the directed graph with vertex set \(\mathbb{S}\) and edge weights \(A_{i, j}\) between each pair of vertices \(s_i, s_j\). For each \(s_i \in \mathbb{S}\), the shortest paths from \(s_0\) to \(s_i\) have nonnegative weight, i.e., \(\delta(s_0, s_i) \geq 0\).
Proof. Assume, for the sake of contradiction, that there is a state \(s \in \mathbb{S}\) such that \(\delta(s_0, s) < 0\). Let \((s_0, s_{i_1}, s_{i_2}, \ldots, s_{i_k} = s)\) be a shortest path from \(s_0\) to \(s\) in \(G\), so that \[\begin{align} A_{0, i_1} + A_{i_1, i_2} + \cdots + A_{i_{k-1}, i_k} = \delta(s_0, s) < 0. \end{align}\] Define \(\epsilon \mathrel{\vcenter{:}}=- (A_{0, i_1} + A_{i_1, i_2} + \cdots + A_{i_{k-1}, i_k}) / k\). Then, for each \(j \in \set{1, 2, \ldots, k}\), there exist \(p_j \in \mathbb{D}\), \(x_j \in \mathbb{N}\), and \(\mathbf{c}_j \in \set{0,1}^*\) such that \(f(p_j, s_{i_{j-1}}, \mathbf{c}_j) = (x_j, s_{i_j})\) and \(\abs{\mathbf{c}_j} - \log(1/p_j(x_j)) < A_{i_{j-1}, i_j} + \epsilon\), by the definition of \(A_{i, j}\). Then, when running the sampler on the target distribution sequence \(p_1 \cdots p_k\) using uniform random bits, the probability of generating the sequence \(x_1 \cdots x_k\) is at least \[\begin{align} 2^{-\abs*{\mathbf{c}_1\cdots\mathbf{c}_k}} > 2^{-\sum_{j=1}^k \log(1/p_j(x_j))}= p_1(x_1)\cdots p_k(x_k), \end{align}\] by our choice of \(\epsilon\), \(p_j\), \(x_j\), and \(\mathbf{c}_j\). But generating the output \(x_1 \cdots x_k\) from the distributions \(p_1 \cdots p_k\) with probability greater than \(p_1(x_1)\cdots p_k(x_k)\) violates the exactness condition ?? . ◻
[theorem:difference-constraints-feasibility, theorem:no-negative-cycles, lemma:state-information-content-nonnegative] imply the following result.
Corollary 1. For any online random sampler, there exists a state information content function as defined in 3.
The state information content function lets us prove lower bounds on the state space required to achieve a given entropy loss bound \(\varepsilon\). In particular, generating any output sample with probability less than \(1\) requires transitioning from one state to another and/or consuming at least one fresh random bit from the input source. For low-entropy target distributions, each transition must therefore lose information, unless the fractional parts of the possible state information content values are sufficiently dense in \(\mathbb{R}/\mathbb{Z}\). Two lemmas and a key theorem formalize these ideas.
Lemma 2 (Information loss over cycles). Let \((s_{i_0}, s_{i_1}, \ldots, s_{i_k}=s_{i_0})\) be any cycle of states in \(\mathbb{S}\), and let distributions \(p_1, \ldots, p_k\), outputs \(x_1, \ldots, x_k\), and input bit sequences \(\mathbf{c}_1, \ldots, \mathbf{c}_k\) satisfy \(f(p_j, s_{i_{j-1}}, \mathbf{c}_j) = (x_j, s_{i_j})\) for each \(j \in \set{1, 2, \ldots, k}\). Then, the information lost over the course of iterating through the cycle once can be lower bounded as follows: \[\sum_{j=1}^k (\abs{\mathbf{c}_j} - \log(1/p_j(x_j))) \geq \modZ*{- \sum_{j=1}^k \log(1/p_j(x_j))}. \qedhere\]
Proof. The string lengths \(\abs{\mathbf{c}_j}\) are integers and the sum on the left-hand side is nonnegative by 4. The result follows from the fact that if \(x - y \ge 0\) where \(x\) is an integer, then \(x \ge \ceil{y}\), so that \(x - y \ge \ceil{y} - y = -\floor{-y} - y \eqqcolon\modZ{-y}\). ◻
Lemma 3 (State space information loss). When using an online random sampler to generate repeated samples from a target distribution \(p\), the information lost over \(\abs{\mathbb{S}}\) samples is at least the minimum loss incurred from iterating through a cycle of states in \(\mathbb{S}\): \[\sum_{j=1}^{\abs{\mathbb{S}}} (\abs{\mathbf{c}_j} - \log(1/p(x_j))) \geq \begin{align}[t] &h(s_{i_{\abs{\mathbb{S}}}}) - h(s_{i_0}) \\ &+ \inf\set*{ \modZ*{\,-\sum_{\mathclap{x\in \mathop{\mathrm{supp}}(p)}} a_x \log(1/p(x))} \; \middle| \; a_x \in \mathbb{N}, 1 \leq \sum_{\mathclap{x\in \mathop{\mathrm{supp}}(p)}} a_x \leq \abs{\mathbb{S}} } \end{align}\] for every sequence of states \(s_{i_0}, s_{i_1}, \ldots, s_{i_{\abs{\mathbb{S}}}}\), outputs \(x_1, \ldots, x_{\abs{\mathbb{S}}}\), and input bit sequences \(\mathbf{c}_1, \ldots, \mathbf{c}_{\abs{\mathbb{S}}}\), such that \(f(p, s_{i_{j-1}}, \mathbf{c}_j) = (x_j, s_{i_j})\) for each \(j \in \set{1, 2, \ldots, \abs{\mathbb{S}}}\). Here, \(h\) is a state information content function for \(f\).
Proof. By the pigeonhole principle, at least one state must be repeated among \(s_{i_0}, s_{i_1}, \ldots, s_{i_{\abs{\mathbb{S}}}}\), so we may take \(L < R\) satisfying \(s_{i_L} = s_{i_R}\). Then, we can decompose the sequence of states into an initial path from \(s_{i_0}\) to \(s_{i_L}\), followed by a cycle from \(s_{i_L}\) to itself (at index \(i_R\)), followed by a final path from \(s_{i_R}\) to \(s_{i_{\abs{\mathbb{S}}}}\).
Over the initial and final paths, the information difference can be bounded by the state information content function \(h\): \[\sum_{j=1}^L (\abs{\mathbf{c}_j} - \log(1/p(x_j))) \geq \sum_{j=1}^L (h(s_{i_j}) - h(s_{i_{j-1}})) = h(s_{i_L}) - h(s_{i_0}),\] and similarly \[\sum_{j=R+1}^{\abs{\mathbb{S}}} (\abs{\mathbf{c}_j} - \log(1/p(x_j))) \geq h(s_{i_{\abs{\mathbb{S}}}}) - h(s_{i_R}).\] The information difference over the cycle is bounded by 2: \[\begin{align} \sum_{j=L+1}^R (\abs{\mathbf{c}_j} - \log(1/p(x_j))) &\geq \modZ*{- \sum_{j=L+1}^R \log(1/p(x_j))} \\ &\geq \inf\set*{ \modZ*{\,-\sum_{\mathclap{x\in \mathop{\mathrm{supp}}(p)}} a_x \log(1/p(x))} \; \middle| \; a_x \in \mathbb{N}, 1 \leq \sum_{\mathclap{x\in \mathop{\mathrm{supp}}(p)}} a_x \leq \abs{\mathbb{S}}}, \end{align}\] where the second inequality follows because the count vector \(a_x = \abs{\set{j \in \set{L+1, \ldots, R} \mid x_j = x}}\) for each \(x \in \mathop{\mathrm{supp}}(p)\) is feasible for the infimum. Combining these three bounds yields the result: \[\begin{align} ~ & \sum_{j=1}^{\abs{\mathbb{S}}} (\abs{\mathbf{c}_j} - \log(1/p(x_j))) \\ ={} & \sum_{j=1}^L (\abs{\mathbf{c}_j} - \log(1/p(x_j))) + \sum_{j=L+1}^R (\abs{\mathbf{c}_j} - \log(1/p(x_j))) + \sum_{j=R+1}^{\abs{\mathbb{S}}} (\abs{\mathbf{c}_j} - \log(1/p(x_j))) \\ \geq{} & h(s_{i_L}) - h(s_{i_0}) + \modZ*{- \sum_{j=L+1}^R \log(1/p(x_j))} + h(s_{i_{\abs{\mathbb{S}}}}) - h(s_{i_R}) \\ \geq{} & h(s_{i_{\abs{\mathbb{S}}}}) - h(s_{i_0}) + \inf\set*{ \modZ*{\,-\sum_{\mathclap{x\in \mathop{\mathrm{supp}}(p)}} a_x \log(1/p(x))} \; \middle| \; a_x \in \mathbb{N}, 1 \leq \sum_{\mathclap{x\in \mathop{\mathrm{supp}}(p)}} a_x \leq \abs{\mathbb{S}}}. \end{align}\] ◻
Theorem 5. For a discrete probability distribution \(p\) and an integer \(M \ge 1\), define \[\Gamma_p(M) \mathrel{\vcenter{:}}= \inf\set*{ \modZ*{\,-\sum_{\mathclap{x\in \mathop{\mathrm{supp}}(p)}} a_x \log(1/p(x))} \;\middle|\; a_x \in \mathbb{N}, 1 \leq \sum_{\mathclap{x\in \mathop{\mathrm{supp}}(p)}} a_x \leq M}. \label{eq:gamma-p}\qquad{(2)}\] Let \(f\) be an online random sampler with state space \(\mathbb{S}\), and \(p \in \mathbb{D}\). If \(f\) achieves entropy loss at most \(\varepsilon\) bits per sample, then \[\abs{\mathbb{S}} \ge \frac{\Gamma_p(\abs{\mathbb{S}})}{\varepsilon}. \qedhere\]
Proof. Run the sampler to generate \(n\abs{\mathbb{S}}\) i.i.d. samples from \(p\). Partition these samples into \(n\) consecutive blocks of \(\abs{\mathbb{S}}\) samples each. For each block \(b \in \set{0,\dots,n-1}\) starting at state \(s_{b\abs{\mathbb{S}}}\) and ending at \(s_{(b+1)\abs{\mathbb{S}}}\), 3 gives \[\sum_{j=b\abs{\mathbb{S}}+1}^{(b+1)\abs{\mathbb{S}}}\left(\abs{\mathbf{c}_j} - \log(1/p(x_j))\right) \ge h(s_{(b+1)\abs{\mathbb{S}}}) - h(s_{b\abs{\mathbb{S}}}) + \Gamma_p(\abs{\mathbb{S}}).\] Summing over the blocks \(b=0,\dots,n-1\), the \(h\) terms telescope as \[\sum_{j=1}^{n \abs{\mathbb{S}}}\left(\abs{\mathbf{c}_j} - \log(1/p(x_j))\right) \ge h(s_{n\abs{\mathbb{S}}}) - h(s_0) + n\Gamma_p(\abs{\mathbb{S}}) \ge n \Gamma_p(\abs{\mathbb{S}}),\] where the second inequality follows from the fact that \(h(s_0) = 0\) and \(h(s_{n\abs{\mathbb{S}}}) \ge 0\). Taking expectations gives entropy loss at least \(n \Gamma_p(\abs{\mathbb{S}})\) for generating \(n\abs{\mathbb{S}}\) samples. By the assumption on the entropy loss we have \(n \Gamma_p(\abs{\mathbb{S}}) \le \varepsilon n\abs{\mathbb{S}} + o(n\abs{\mathbb{S}})\). Dividing by \(\varepsilon n\) and taking \(n \to \infty\) gives \(\Gamma_p(\abs{\mathbb{S}})/\varepsilon\le \abs{\mathbb{S}}.\) ◻
Proof. To establish , consider the target distribution \(p \mathrel{\vcenter{:}}=\mathrm{Bernoulli}(1/3)\). The fractional part of the negative information content of each outcome is \(\modZ{-\log(3)} = \modZ{-\log(3/2)} = 2-\log(3)\). A useful result from number theory, discussed further in 6.1, is that \[\begin{align} \min_{1 \leq a \leq M} \modZ*{-a \log 3} \geq M^{-4.116201 - o(1)}, \label{eq:apply-diophantine-approx} \end{align}\tag{7}\] where the \(o(1)\) term goes to \(0\) as \(M \to \infty\). Therefore, 5 gives the result: \[\begin{align} \abs{\mathbb{S}} \ge \frac{\Gamma_p(\abs{\mathbb{S}})}{\varepsilon} \ge \frac{\abs{\mathbb{S}}^{-4.116201 - o(1)}}{\varepsilon} \quad\implies\quad \log(\abs{\mathbb{S}}) \ge \left(\frac{1}{5.116201} - o(1)\right) \log\frac{1}{\varepsilon}, \end{align}\] where the \(o(1)\) term now goes to \(0\) as \(\varepsilon\to 0\).
To establish , consider the target distribution \(p_m \mathrel{\vcenter{:}}=\mathrm{Bernoulli}(1/(2^m+1))\) for some \(m \in \mathbb{N}\). The fractional part of the negative information content of each outcome is \(\modZ*{-\log(2^m+1)} = \modZ*{-\log(1+2^{-m})} = 1 - \log(1+2^{-m}) \eqqcolon 1 - \alpha_m\). For any \(M \geq 1\), we have \[\Gamma_{p_m}(M) = \min_{1 \le a \le M}\modZ*{-a\thinspace \alpha_m} \geq \min_{1 \le a \le M}(1 - a\thinspace \alpha_m) = 1 - M\alpha_m.\] Let \(\mathcal{M} \subset \mathbb{N}\) denote the set of all \(m\) for which \(p_m \in \mathbb{D}\). There exist infinitely many such \(m\) by assumption, so taking \(m \to \infty\) along the subsequence ensures that \(\alpha_m \to 0\). Therefore, 5 gives the reuslt: \[\begin{align} \abs{\mathbb{S}} \ge \frac{\Gamma_{p_m}(\abs{\mathbb{S}})}{\varepsilon} \ge \frac{1 - \abs{\mathbb{S}}\thinspace \alpha_m}{\varepsilon} \quad\implies\quad \abs{\mathbb{S}} \ge \sup_{m \in \mathcal{M}} \frac{1 - \abs{\mathbb{S}}\thinspace \alpha_m}{\varepsilon} = \frac{1}{\varepsilon}. \end{align}\] ◻
The bound in of 1 makes use of the fact 7 that integer multiples of \(\log(3)\) cannot be “too close” to integers, or equivalently, that \(\log(3)\) cannot be approximated too well by rational numbers with small denominators. More precisely, [14]’s theorem on linear forms in logarithms [14] shows that there exists a constant \(c > 0\) such that \[\begin{align} \modZ{-a\log(3)} \geq \frac{1}{a^c} && (a = 1, 2, \dots). \end{align}\] Approximation by rational numbers, also called Diophantine approximation, is an important topic studied in number theory. The hardness of approximation of a single number (resp. a set of numbers) is measured by the irrationality exponent (resp. the linear independence measure), which we define following the exposition of [15].
Definition 4. The irrationality exponent* of an irrational real number \(\alpha\), denoted \(\mu(\alpha)\), is the minimum real number such that, for every \(\epsilon > 0\), there exists \(q_0(\epsilon)\) such that \[\abs*{\alpha - \frac{p}{q}} \geq \frac{1}{q^{\mu(\alpha) + \epsilon}}\] for all integers \(p\) and \(q\) with \(q > q_0(\epsilon)\).*
Definition 5. The linear independence measure* of a set of \(\mathbb{Q}\)-linearly independent real numbers \(\alpha_0, \ldots, \alpha_k\), denoted \(\nu(\alpha_0, \ldots, \alpha_k)\), is the minimum real number such that, for every \(\epsilon > 0\), there exists \(H_0(\epsilon)\) such that \[\abs*{p\alpha_0 + q_1\alpha_1 + \cdots + q_k\alpha_k} \geq H^{-\nu(\alpha_0, \ldots, \alpha_k) - \epsilon}\] for all integers \(p, q_1, \ldots, q_k\) satisfying \(H > H_0(\epsilon)\), where \(H \mathrel{\vcenter{:}}=\max\set{\abs{q_1},\ldots,\abs{q_k}}\).*
An immediate consequence of the definition of a linear independence measure is that \[\modZ*{q_1\alpha_1 + \cdots + q_k\alpha_k} \geq \max\set{\abs{q_1},\ldots,\abs{q_k}}^{-\nu(1, \alpha_1, \ldots, \alpha_k) - o(1)}\] where \(o(1)\) describes the limiting behavior as \(\max\set{\abs{q_1},\ldots,\abs{q_k}} \to \infty\).
Although [14] provided the first general proof that sets of logarithms of algebraic numbers have finite linear independence measures, the constants were far from optimal. Extensive work has been devoted to improving the constant bounds for specific numbers, and the current best known bound on the irrationality exponent of \(\log(3)\) is \(\mu(\log(3)) \leq 5.116201\), which is a corollary of the stronger result \(\nu(1, \ln(2), \ln(3)) \leq 4.116201\) [16]. This bound implies 7 , which yields the space lower bound of \((1/5.116201 - o(1))\log(1/\varepsilon)\) bits when sampling from \(\mathrm{Bernoulli}(1/3)\).
We now generalize the Diophantine approximation analysis to prove the space lower bound in 2, which applies to a large class of samplers that emit i.i.d. outputs from a fixed distribution.
When considering i.i.d. samples from a target distribution \(p\) where \(\set{\modZ{\log p(x)} \mid x \in \mathop{\mathrm{supp}}(p)}\) contains only one element, it suffices to analyze the irrationality exponent of the single number \(\modZ{\log p(x)}\) to obtain a state-space lower bound. In general, when \(\abs{\set{\modZ{\log p(x)} \mid x \in \mathop{\mathrm{supp}}(p)}} > 1\), our proof method from 1 requires analyzing the linear independence measure of the set \(\set{1} \cup \set{\modZ{\log p(x)} \mid x \in \mathop{\mathrm{supp}}(p)}\). With a more careful analysis, we can attain lower bounds even if there is a certain \(\mathbb{Q}\)-linear dependence among the \(\modZ{\log p(x)}\) values, as long as their \(\mathbb{Q}\)-linear combinations with nonnegative coefficients do not include any nonzero integers, as shown in the following theorem.
Theorem 6 (Space lower bound from linear independence measure). Let \(p\) be a finite discrete distribution over support \(\set{x_1, \ldots, x_k}\). Let \(\set{1, \alpha_1, \ldots, \alpha_m}\) be a \(\mathbb{Q}\)-linearly independent set of real numbers and \(q \in \mathbb{N}_{>0}\), such that for each \(i \in \set{1,\ldots,k}\) we can write \[\log p(x_i) = \beta_i / q + \lambda_{i,1} \alpha_1 + \ldots + \lambda_{i,m} \alpha_m\] for some coefficients \(\beta_i \in \mathbb{Z}\) and \(\lambda_{i,j} \in \mathbb{N}\) with \(\sum_{j=1}^m \lambda_{i,j} > 0\). If \(\set{1, \alpha_1, \ldots, \alpha_m}\) has finite linear independence measure \(C \mathrel{\vcenter{:}}=\nu(1, \alpha_1, \ldots, \alpha_m)\), then generating i.i.d. samples from \(p\) with an entropy loss bound of \(\varepsilon\) bits per sample requires a state space of size \(\abs{\mathbb{S}} \geq (1/\varepsilon)^{1/(1+C)-o(1)}\).
Proof. By 5, it suffices to show there exists a function \(\delta(M) \to 0\) as \(M\to\infty\) such that \[\Gamma_p(M) \ge M^{-C-\delta(M)}. \label{eq:gamma-p-bound}\tag{8}\] After setting \(M = \abs{\mathbb{S}}\), 5 gives \[\begin{align} \abs{\mathbb{S}} \ge \frac{\Gamma_p(\abs{\mathbb{S}})}{\varepsilon} \ge \frac{\abs{\mathbb{S}}^{-C-\delta(\abs{\mathbb{S}})}}{\varepsilon} \quad\implies\quad \abs{\mathbb{S}} \ge \left(\frac{1}{\varepsilon}\right)^{1/(C+1+\delta(\abs{\mathbb{S}}))} \geq \left(\frac{1}{\varepsilon}\right)^{1/(C+1) - o(1)}. \end{align}\] The linear-independence hypotheses imply that \(\Gamma_p(M)>0\) for every fixed \(M\) so 5 rules out bounded \(\abs{\mathbb{S}}\) as \(\varepsilon\to0\). Therefore \(\abs{\mathbb{S}}\to\infty\) as \(\varepsilon\to0\), so \(\delta(\abs{\mathbb{S}})=o(1)\) in the \(\varepsilon\to0\) limit.
We now prove 8 : \[\begin{align} \Gamma_p(M) ={} & \min\set*{\modZ*{\sum_{i=1}^k a_i \log p(x_i)} \; \middle| \; a_i \in \mathbb{N}, \, 1 \leq \sum_{i=1}^k a_i \leq M} \\ ={} & \min\set*{\modZ*{\sum_{i=1}^k a_i \left(\frac{\beta_i}{q} + \sum_{j=1}^m \lambda_{i,j}\alpha_j \right)} \; \middle| \; a_i \in \mathbb{N}, \, 1 \leq \sum_{i=1}^k a_i \leq M} \\ \geq{} & \frac{1}{q} \min\set*{\modZ*{\sum_{i=1}^k a_i q\sum_{j=1}^m \lambda_{i,j}\alpha_j} \; \middle| \; a_i \in \mathbb{N}, \, 1 \leq \sum_{i=1}^k a_i \leq M} \tag{9} \\ ={} & \frac{1}{q} \min\set*{\modZ*{ \sum_{j=1}^m \left(q\sum_{i=1}^k a_i \lambda_{i,j}\right) \alpha_j} \; \middle| \; a_i \in \mathbb{N}, \, 1 \leq \sum_{i=1}^k a_i \leq M} \\ \geq{} & \frac{1}{q} \min\set*{\modZ*{\sum_{j=1}^m b_j \alpha_j} \; \middle| \; b_j \in \mathbb{N}, \, 1 \leq \sum_{j=1}^m b_j \leq M q \sum_{i=1}^k \sum_{j=1}^m \lambda_{i,j}} \tag{10} \\ \geq{} & M^{-C - \delta(M)}. \tag{11} \end{align}\] We can extract the factor of \(1/q\) in 9 using the identity \(\modZ{z/q+r} \geq \modZ{z+rq}/q = \modZ{rq}/q\) for \(r \in \mathbb{R}\) and \(q \in \mathbb{N}_{>0}\) and \(z \in \mathbb{Z}\). We then reduce to 10 by reparameterizing and bounding the new coefficients \(b_j\). Lastly, 11 follows from [remark:modZ-bound-from-linear-independence-measure]. The fixed constants \(1/q\) and \(q\sum_{i=1}^k\sum_{j=1}^m \lambda_{i,j}\) are absorbed by increasing \(\delta(M)\) by \(O(1/\log M)\), which still leaves \(\delta(M) \to 0\). ◻
The Khintchine-Groshev theorem [17] shows that the linear independence measure of almost all tuples of \(k\) real numbers (together with \(1\)) is \(k\). This result nearly suffices to lower-bound the space complexity for i.i.d. sampling from almost all finite discrete distributions, except that the set of probability distributions over \(k\) outcomes is a \((k-1)\)-dimensional manifold, which has measure zero in \(\mathbb{R}^k\). Therefore, showing that \(\nu(1, \modZ{\log p(x_1)}, \ldots, \modZ{\log p(x_k)}) = k\) for almost all distributions requires a stronger result providing lower bounds for Diophantine approximation on manifolds.
Theorem 7 (Diophantine approximation on manifolds [18]). Let \(U \subseteq \mathbb{R}^d\) be an open set, and let \(\mathbf{f}(U) = \set{\mathbf{f}(\mathbf{u}) \mid \mathbf{u} \in U}\) be a submanifold of \(\mathbb{R}^k\). Assume that for almost all \(\mathbf{u} \in U\), the set of partial derivatives of \(\mathbf{f}\) at \(\mathbf{u}\) up to some fixed finite order spans \(\mathbb{R}^k\). Then, for almost all \(\mathbf{u} \in U\), the numbers \((1, f_1(\mathbf{u}), \ldots, f_k(\mathbf{u}))\) are \(\mathbb{Q}\)-linearly independent with linear independence measure \[\nu(1, f_1(\mathbf{u}), \ldots, f_k(\mathbf{u})) \leq k,\] where \((f_1(\mathbf{u}), \ldots, f_k(\mathbf{u})) = \mathbf{f}(\mathbf{u})\).
We can now prove 2, which is restated below.
Proof. Apply 7 with \[\begin{align} U &= \set*{(p_1,\ldots,p_{k-1}) \mid p_i > 0, \sum_{i=1}^{k-1} p_i < 1}; \\ \mathbf{f}(p_1,\ldots,p_{k-1}) &= \left(\log(p_1),\ldots,\log(p_{k-1}), \log\!\left(1-\sum_{i=1}^{k-1} p_i\right)\right). \end{align}\] At any point, the first-order partial derivatives are \[\frac{\partial}{\partial p_j} \mathbf{f}(p_1,\ldots,p_{k-1}) = \frac{1}{\ln(2)} \left(\ldots,p_j^{-1},\ldots,-\left(1-\sum_{i=1}^{k-1} p_i\right)^{-1}\right),\] where only the \(j\)th and \(k\)th coordinates are nonzero. Furthermore, the second-order partial derivative with respect to \(p_1\) is \[\frac{\partial^2}{\partial p_1^2} \mathbf{f}(p_1,\ldots,p_{k-1}) = \frac{1}{\ln(2)} \left(-p_1^{-2},0,\ldots,0,-\left(1-\sum_{i=1}^{k-1} p_i\right)^{-2}\right),\] which, together with the \(k-1\) first-order partial derivatives, forms a linearly independent set and therefore spans \(\mathbb{R}^k\). Therefore, the theorem applies to this manifold and yields \[\nu(1,\log(p_1),\ldots,\log(p_k)) \leq k\] for almost all finite discrete distributions \((p_1,\ldots,p_k)\).
Finally, 6 with \((\alpha_1,\ldots,\alpha_m)=(\log(p_1),\ldots,\log(p_k))\) yields \[\abs{\mathbb{S}} \geq (1/\varepsilon)^{1/(1+k)-o(1)}\] for almost all distributions, and taking the logarithm gives the desired result. ◻
Previous work on random sampling has mainly focused on either entropy efficiency or computational complexity, with algorithms matching lower bounds on both fronts. [3] derive a DDG tree algorithm that uses the exact minimum possible number of input random bits, whether sampling from a single distribution or from a sequence of distributions; however, this algorithm uses unbounded space and per-sample computation time when sampling an infinite sequence [7]. Computationally optimal sampling algorithms depend on the computational model and format of the target distribution; of particular note is the [19] alias method [19], which represents the target distribution in a linear-size data structure and generates each sample using a constant number of word-level operations, matching a trivial lower bound. Although the original alias method has significant entropy loss, the randomness recycling method of [7] reduces the entropy loss to \(\varepsilon\) bits per sample while increasing the space by \(O(\log(1/\varepsilon))\) bits and adding a constant number of \(O(\log(1/\varepsilon))\)-bit arithmetic operations per sample. Therefore, our lower bound in 1 shows that for rational distributions, the space usage of randomness recycling [7] is optimal up to constant factors, and the runtime is optimal, up to the cost of arithmetic on the state.
Our lower bounds refer to the persistent space \(s \in \mathbb{S}\) that is carried by a sampler across each round of online random sampling. Sampling algorithms may also use temporary space, which is allocated and released within a single round. While temporary space must be unbounded to handle irrational probabilities, its expectation may be finite. For this setting, [20] study the problem of online sampling given CDF oracles that produce dyadic rational approximations to irrational probabilities, up to any precision. They show that over \(n\) rounds of online random sampling, the expected overall space (i.e., maximum of temporary and persistent) must grow as \(\Omega(\log n)\) [20]. They also introduce an algorithm that matches both their lower bound on expected space and our lower bound on persistent space, using deterministically \(\Theta(\log(1/\varepsilon))\) bits of persistent space to achieve an entropy loss bound of \(\varepsilon\) bits per sample. Their algorithm [20] uses an iteration counter \(t\) to scale a parameter \(\Delta_t\) which upper-bounds \(M\), in addition to the uniform state \((Z,M)\), and thus its persistent space is \((3+o(1))\log(1/\varepsilon)\) bits. However, a variant with constant \(\Delta = \widetilde{\Theta}(1/\epsilon)\) requires only \((2+o(1))\log(1/\varepsilon)\) bits of persistent space \((Z,M)\), which is optimal up to a factor of \(2\) by 1, .
Orthogonally to this article’s focus on asymptotic per-sample entropy loss bounds, several works have studied the non-asymptotic efficiency of generating a single sample from a distribution. For this case, the [3] DDG tree algorithm has an entropy loss bounded by \(2\). Explicitly representing the entire DDG tree as in [11] requires exponential space [21], if the target distribution has finite support and arbitrary rational probabilities. This space cost can be bypassed by instead lazily generating the tree during traversal [9] Algorithm 1, [3, p. 384] (at the cost of higher runtime during sampling) or using rejection sampling [5], [10] (at the cost of slightly higher entropy cost).
There exist results in the information theory literature that provide lower bounds on the arithmetic precision or delay (similar to our space complexity) of entropy-efficient codes. [22] shows that for arithmetic codes, achieving a redundancy (similar to our entropy loss per sample) of \(\varepsilon\) requires \(\Omega(\log(1/\varepsilon))\) bits of precision for representing the probabilities; although their model is quite restricted, similar techniques may show that online random samplers based on discretizations of the interval method of [4] require \(\Omega(\log(1/\varepsilon))\) space. [23] show that for almost all distributions, achieving a redundancy of \(\varepsilon\) requires \(\Omega(\log(1/\varepsilon))\) delay; although their lower bound is quite general, and a delay of \(\Omega(\log(1/\varepsilon))\) would imply a space complexity of \(\Omega(\log(1/\varepsilon))\) if each symbol had to be explicitly stored over the entire window of the delay, an efficient coder may compress the symbols at intermediate stages to achieve a lower space complexity. Random sampling and coding theory are closely related but subtly distinct, and to our knowledge the results in our paper do not imply, nor are implied by, the aforementioned results in coding theory. It remains an interesting question to what extent lower bounds or proof techniques can be adapted from one setting to the other.
It is well known that zero entropy loss can be achieved for any dyadic distribution (for which each probability is in \(\set{0}\cup\set{2^{-n} \mid n \in \mathbb{N}}\)) using no persistent state [3]. Our analysis of entropy loss purely in terms of the state information content function does not rule out the possibility of a sampler that achieves an entropy loss of \(\varepsilon= 0\) for particular non-dyadic distributions using a finite state space. However, we conjecture that such a sampler does not exist.
Another open question is what constant factor for the \(\Omega(\log(1/\varepsilon))\) space is necessary, and to what extent it depends on the target distribution set \(\mathbb{D}\). [7] show that \(2\) is an upper bound on the constant factor whenever \(\mathbb{D}\) consists of finitely many rational distributions over finite support, and a variant on the method of [20] achieves this same bound for arbitrary distributions. This constant factor could potentially be reduced to \(1\), according to 1, . On the lower-bound side, the state information content cannot take entirely arbitrary values, because there must be some relationship between states that allows the algorithm to generate correct samples. For example, a state space of size \(\abs{\mathbb{S}}=2\) is no better than one of size \(\abs{\mathbb{S}}=1\), as a consequence of the exactness condition in 1. A more fine-grained analysis may therefore strengthen the lower bounds we have presented.
An orthogonal research direction is sampling from continuous distributions, iteratively outputting new bits in the binary expansion of the generated real number. In particular, given a cumulative distribution function \(F\) satisfying \(F(0)=0\) and \(F(1)=1\), the problem is to generate outputs \(X_1,X_2,\ldots\) following the conditional distributions \[(X_{n+1} \mid X_1,\ldots,X_n) \sim \mathrm{Bernoulli}\left( \frac{F(2^{-n}+\sum_{i=1}^{n}X_i/2^i) - F(2^{-n-1}+\sum_{i=1}^{n}X_i/2^i)}{F(2^{-n}+\sum_{i=1}^{n}X_i/2^i) - F(\sum_{i=1}^{n}X_i/2^i)} \right).\] study finite-state generators for sampling from particular distributions, and pose several related questions about which distributions can be sampled in this way and how large the state space must be in order to come within \(\varepsilon\) of the entropy-optimal algorithm [3]. Question ‘(v)’, about which distributions have finite-state generators, was resolved by [24], who showed that a smooth CDF has a finite-state generator if and only if it is a polynomial with rational coefficients, whose derivative has no irrational roots on \([0,1]\). However, question ‘(iii)’, about which distributions have entropy-optimal finite-state generators; and question ‘(iv)’, about how large of a state space is needed in order to come within \(\varepsilon\) of the entropy-optimal algorithm for \(F(x)=x^3\); remain unanswered to the best of our knowledge. Our lower bounds depend crucially on the fact that an online random sampler must be able to produce correct samples from any of its target distributions at any time step, whereas for continuous samplers the output bits approach the \(\mathrm{Bernoulli}(1/2)\) distribution, which allows significant optimizations. In fact, the \(\varepsilon\) mentioned by [3] is a uniform bound on the entropy loss relative to the optimal algorithm at each iteration count, rather than a per-iteration loss which accumulates as \(\varepsilon n\) in our problem setting. Therefore, entropy-space lower bounds for continuous distributions remain an open question.
This material is based upon work supported by the National Science Foundation under Award No. 2311983. Any opinions, findings and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of the National Science Foundation.