Procedure: Generating a Sobol’s Sequence¶
Description and Background¶
The numbers in a Sobol’s sequence are generated directly as binary fractions of length \(w\) bits. These numbers are created from a special set of \(w\) binary fractions \(V_i\) called direction numbers. In Sobol’s original method, the j-th number \(X_j\) is generated by the operation XOR among those direction numbers \(V_i\) such that the i-th bit of \(j\) is not zero, see comments below. See DiscSobol for additional discussion.
Detailed description of the procedure¶
The \(i\)-th number in Sobol’s sequence is
where \(b_j(i)\) is the \(j\)-th digit of the binary expansion of \(i\); \(v_j\) is the \(j\)-th direction number and \(\oplus\) is the binary (bitwise) XOR operation. The number \(b_j(i)\) is a bit and thus it only takes values \(0\) or \(1\), i.e. for an integer \(i\) we have its binary expansion
As for any finite number \(i\), its binary expansion has a finite number of non-zero digits, the Sobol’s number is computed with a XOR sum of a finite number of direction numbers and thus no convergence properties are needed in the above equations.
An important element in Sobol’s sequence are the direction numbers. The direction numbers are binary (dyadic) fractions \(v_j\in(0,1)\) and the resulting Sobol’s number is also a binary fraction \(x_i\in(0,1)\). The direction numbers are represented as
where \(m_j\) takes integer values, \(1<m_j<2^j\), i.e. the direction number \(v_j\) is a shift of the binary digits in \(m_j\). The direction numbers are computed by the recursive relation
where the initial values \(m_1,\ldots,m_d\) are odd integers and the numbers \(a_0,a_1,\ldots,a_d\) are given by the coefficients of \(p_d(x)\), a primitive polynomial of degree \(d\) in \({Z}_2,\) which is written as
Recall that a primitive polynomial is a polynomial that cannot be represented as the product of two polynomials of lower degree. Primitive polynomials are available in tables in the literature, see Press et al. (1994). For a primitive polynomial, the numbers \(a_0\) and \(a_d\) take the value one by definition, and \(a_0\) is not used in the computations.
Running example¶
In what follows, the subindex \(_{\textbf{2}}\) refers to binary numbers, everywhere else usual representation (base ten) is assumed.
Take the primitive polynomial \(p_3(x)=x^3+x+1\), and thus \(d=3\) and \((a_0,a_1,a_2,a_3)=(1,0,1,1)\). The recurrence relation for the direction numbers is
Using a different primitive polynomial will give a different recurrence relation. To iterate the recurrence relation, initial odd values \(m_1,m_2,m_3\) are needed. Take \(m_1=1,m_2=3,m_3=7.\)
Now we compute a few extra numbers \(m_j\). We have
and
The direction numbers are computed by shifting the values obtained, e.g.
and
The first Sobol’s number is obtained with the first direction number
As the binary expansion of two is \(10_{\textbf{2}}\) then
and
The following table summarises the construction performed, with the Sobol’s sequence on the two rightmost columns, the next to last in binary and the last in base ten.
| \(i,j (10)\) | \(i,j\) | \(m_j(10)\) | \(m_j\) | \(v_j\) | \(x_i\) | \(x_i(10)\) |
|---|---|---|---|---|---|---|
| 1 | 1 | 1 | 1 | .1 | .1 | .5 |
| 2 | 10 | 3 | 11 | .11 | .11 | .75 |
| 3 | 11 | 7 | 111 | .111 | .01 | .25 |
| 4 | 100 | 5 | 101 | .0101 | .111 | .875 |
| 5 | 101 | 7 | 111 | .00111 | .011 | .375 |
| 6 | 110 | .001 | .125 | |||
| 7 | 111 | .101 | .625 | |||
| 8 | 1000 | .0101 | .3125 | |||
| 9 | 1001 | .1101 | .8125 | |||
| 10 | 1010 | .1001 | .5625 |
Additional Comments, References and Links¶
It is important to note that for the first Sobol’s number, only the first direction number was needed; then for the following two Sobol’s numbers the second direction number was included; for the following four Sobol’, the third direction number was included. Without iterating again the recursive relation, \(31\) Sobol’s numbers can be constructed using the first five direction numbers in the above table.
In short, to construct the first \(2^k-1\) Sobol’s numbers, we need \(k\) direction numbers. If more Sobol’s numbers are needed, then the recursive equation must be iterated to obtain direction numbers as required.
Note that by selecting odd initial values \(m_1,\ldots,m_d\), all the subsequent \(m_{d+1},m_{d+2},\ldots\) are guaranteed to be odd numbers and thus the \(i\)-th bit of the direction number \(v_i\) is one. This has the important consequence of allowing Sobol’s numbers to lie in consecutive finer binary meshes. In other words, a latin hypercube is constructed with the first \(2^k-1\) Sobol’s points.
The construction of multivariate Sobol’s sequences is achieved by using different primitive polynomials for each dimension. For a table with different primitive polynomials see Press et al. (1994). Sobol’s gave a list of recommended primitive polynomials, to avoid high correlations between different dimensions.
An alternative version of Sobol’s sequence was due to Antonov and Saleev, who proposed taking instead
where \(g_j(i)\) is the \(j\)-th digit of the Gray code representation of \(i\). This different Sobol’s proposal is faster than the original, as it simplifies the computation to \(x_{i+1}=x_i\oplus v_c,\) where \(b_c\) is the rightmost zero bit in the representation of \(i\).
References:
Antonov, Saleev (1979). USSR Comput. Math. Math. Phys. 19, 252-256.
Bratley, Fox (1988), ACM Trans. Math. Soft. 14(1), 88-100.
Press et al. (1994). Numerical Recipes in C, Cambridge.