我提出了一种新的低偏差准随机序列，它比其他流行的序列(如Sobol和Halton序列)提供了许多实质性的改进。
这篇博客文章前段时间曾出现在“黑客新闻”的头版。
在图1b中，可以看到单位正方形内的点的简单均匀随机采样表现出点的聚集，并且也存在根本不包含点的区域(“白噪声”)。低差异序列和准随机序列是一种以确定性方式构建(无限)序列点的方法，该方法降低了群集(差异)的可能性，同时仍然确保整个空间被均匀覆盖(“蓝色噪声样本”)。
也就是说，准随机序列可用于创建点分布，这些点分布看起来不如晶格规则，但比随机采样更规则(参见图1b)。众所周知的准随机序列包括Halton序列和Sobol序列。它们在许多数值计算应用中扮演着重要的角色，包括物理、金融和近几十年来的计算机图形学。
在一维上创建完全确定的低偏差准随机序列的方法被研究得非常好，并且得到了普遍的解决。在这篇文章中，我几乎只关注开放(无限)序列，首先在一维，然后扩展到更高的维。开序列(即，在$n$中可扩展)的基本优势在于，如果基于有限数量项的结果误差太大，则可以扩展序列，而无需丢弃所有先前计算的点。有很多方法可以构造开放序列。对不同类型进行分类的一种方法是通过构造它们的基础(超)参数的方法：
为简明起见，这篇文章一般集中于比较这种新的加性递归$R$-序列，它属于第一类，因为它是一种基于秩为1的无理数(通常称为Kronecker、Weyl或Richtmyer序列)的递归方法，与Halton序列进行比较，Halton序列基于规范的一维van der Corput序列。规范的Kronecker递推序列定义为：$$R_1(\α)：\；\；t_n=\{s_0+n\α\}，\quad n=1，2，3，…。$$其中$\alpha$是任何无理数。请注意，符号${x\}$表示$x$的小数部分。在计算中，此函数更常见地表示为$$R_1(\α)：\；\；t_n=s_0+n\α\；(\tExtrm{mod}\；1)；\quad n=1，2，3，…。$$对于$s_0=0$，序列$R(\φ)$的前几项是：$$t_n=0.618，\；0.236，\；0.854，\；0.472，\；0.090，\；0.708，\；0.327，\；0.944，\；0.562，\；0.180，\；798，\；416，\；0.034，\；0.652，\；0.271，\；0.888，…。$$。
需要注意的是，$s_0$的值不会影响序列的总体特征，并且在几乎所有情况下都设置为零。但是，特别是在计算中，$s\neq 0$选项提供了一个通常有用的额外自由度。如果$s\neq 0$，则该序列通常被称为‘移位格序列。
如果\(\α=1/\φ\)，则获得给出最小可能差异的\(\α\)值，其中\(\φ\)是黄金比例。即，$$\phi\EQUIV\FRAC{\sqrt{5}+1}{2}\simeq 1.61803398875…。；$$有趣的是，也有无穷多个其他的$\alpha$值也实现了最优差异，并且它们都通过Moebius变换$$\alpha‘=\frac{p\alpha+q}{r\alpha+s}\quad\tExtrm{对于所有整数}\；p，q，r，s\quad\texm{从而}|ps-qr|=1$$进行比较。我们现在将这种递推方法与著名的van der Corput反基序列[van der Corput，1935]进行比较。Van der Corput序列实际上是一族序列，每个序列都由一个唯一的超参数定义。b=2的序列的前几个项是：$$t_n^{[2]}={\frac{1}{2}，{\frac{1}{4}，\frac{3}{4}，\frac{1}{8}，\frac{5}{8}，\frac{3}{8}，\frac{7}{8}，\frac{1}{16}，\frac{1}{4}，\frac{1}{8}，\frac{3}{8}，\frac{1}{16}，\FRAC{9}{16}，\FRAC{5}{16}，\FRAC{13}{16}，\FRAC{3}{16}，\FRAC{11}{16}，\FRAC{7}{16}，\FRAC{15}{16}，…。$$以下部分比较了每个序列的一般特征和有效性。考虑求定整数$$A=\int_0^1 f(X)\tExtrm{d}x$$的任务我们可以用以下公式近似：$$A=\frac{1}{n}\sum_{i=1}^{n}f(X_I)，[0，1]$$中的{quad x_i
如果({x_i\})是随机选择的，这是蒙特卡罗(Monte Carlo)方法；
如果(Xi)是低偏差序列的元素，这就是拟蒙特卡罗方法。
下图显示了逼近与函数有关的定积分的典型误差曲线\(s_n=|A-A_n|\)，\(f(X)=\tExtrm{exp}(\frac{-x^2}{2})，\；[0，1]中的x\)具有：(I)基于加性递归的准随机点，其中S\(\α=1/\φ\)，(蓝色)；(Ii)基于van der Corput序列的准随机点，(橙色)；(Iii)随机选择的点，(绿色)；(Iv)Sobol序列(红色)。结果表明，对于n=10^6的点，随机抽样方法的误差为10^-4，van der Corput序列的误差为10-6，而R(Phi)-序列的误差为10-7，比van der Corput误差和均匀随机抽样的误差分别提高10倍和1000倍，说明随机抽样的误差是随机抽样的10倍，比均匀随机抽样的误差高1000倍，而对于n=10^6的点，随机抽样的误差为10-4，范德科普特序列的误差为5-(simeq 10^-6})，而R(Phi)-序列的误差为5-(simeq 10^-7})，比均匀随机抽样的误差高1 0 0 0倍。
这与基于均匀随机抽样的误差渐近减小$1/\sqrt{n}$的知识是一致的，而基于两个准随机序列的误差曲线趋于$1/n$。
它表明，van der Corput序列为积分问题提供了令人难以置信的一致结果！
结果表明，对于所有值$n$，$R_1(\φ)$-序列比van der Corput序列产生更好的结果。
新的$R_1$序列是使用黄金比的Kronecker序列，是一维拟随机蒙特卡罗(QMC)积分方法的最佳选择之一。
还应该注意的是，尽管理论上\(\α=\φ\)提供了可证明的最优选项，但\(\sqrt{2}\)非常接近最优，并且\(\α\)的几乎任何其他无理性值都为一维积分提供了极好的误差曲线。这就是为什么对于任何素数来说，\(\α=\sqrt{p}\)是非常常用的。此外，从计算的角度来看，在[0，1]$中的$\α\区间中选择一个随机值几乎肯定会是一个无理数(在机器精度范围内)，因此对于低偏差序列来说是一个可靠的选择。为直观起见，上图没有显示Niederreiter序列的结果，因为其结果与Sobol和$R$序列的结果几乎没有区别。本文中使用的Neiderreiter和Sobol序列(及其优化的参数选择)是使用“英特尔MKL库中提供的封闭的、专有的和完全优化的生成器”通过MATHEMICA进行计算的。
目前大多数构造高维低偏差的方法都是简单地(以分量方式)将\(d\)个一维序列组合在一起。为简明起见，本文主要描述Halton序列[Halton，1960]，Sobol序列和\(d-\)维Kronecker序列。
Halton序列是通过使用\(d\)不同的一维van de Corput序列来简单地构造的，每个序列的碱基相对于所有其他序列都是相对质数的。也就是说，两两共同质数。到目前为止，由于其明显的简单性和敏感性，最频繁的选择是选择第一(D)个素数。由(2，3)-Halton序列定义的前625个点的分布如图1所示。虽然许多二维Halton序列是低偏差序列的极好来源，但众所周知，许多是高度有问题的，并且不表现出低偏差。例如，图3显示(11，13)-Halton序列产生高度可见的线。在选择哪些(\(p_1，p_2)\)对是样本，哪些是有问题的对的方法上已经做了很大的努力。这个问题在更高的维度上就更成问题了。
当推广到更高维时，Kronecker递推方法通常面临更大的挑战。也就是说，虽然使用\(\alpha=\sqrt{p}x\)可以产生优秀的一维序列，但即使是找到两个素数对来作为二维情况的基础，也是非常具有挑战性的，而且这些素数对几乎没有问题！作为解决此问题的一种方法，有些人建议使用其他众所周知的无理数，例如\(\φ，\pi，e，…。\)。这些方法会产生中等程度可接受的解决方案，但通常不会使用，因为它们通常不如精心选择的Halton序列好。大量的努力都集中在这些退化问题上。
建议的解决方案包括跳跃/燃烧、跳跃/稀疏。而对于有限序列，加扰是另一种经常用于克服该问题的技术。加扰不能用于创建开放(无限)低差异序列。
类似地，尽管Sobol序列的性能通常较好，但其复杂性以及更重要的是要求非常谨慎地选择其超参数，使得它不那么诱人。
新的$R_d$序列是唯一的不需要任何基参数选择的$d$维低偏差准随机序列。
在这一部分中，我展示了如何构造一类新的\(d-\)维开(无限)低偏差序列，它不需要选择任何基本参数，并且具有突出的低偏差性质。
有许多可能的方法来推广斐波纳契序列和/或黄金比率。下面提出的推广黄金比例的方法并不新鲜[Krcadinac，2005]。特征多项式与许多代数领域有关，包括Perron数和Pisot-Vijayaraghavan数。然而，最新的是这种推广形式与高维低差序列的构造之间的明确联系。我们将黄金比的广义形式\(\φ_d\)定义为唯一的正根$x^{d+1}=x+1$。那是,
对于\(d=2\)，\(\φ_2=1.324717957244746…。\)，它通常被称为塑料常数，并且具有一些漂亮的性质(另请参阅此处)。这个值被推测为最有可能是相关二维问题的最佳值[Hensley，2002]。
对于$d>；3$，虽然该方程的根没有封闭的代数形式，但我们可以很容易地通过标准方法(如牛顿-拉普森)或通过注意对于以下序列\(R_d(\φ_d)\)：$$t_0=t_1=…来获得数值近似。=t_{d}=1；$t_{n+d+1}\；=\；t_{n+1}+t_n，\quad\texm{for}\；n=1，2，3，.$$。
1928年，建筑师兼和尚汉斯·范德兰(Hans Van De Laan)将这一特殊的常数序列$\φ_d称为“和谐数”。这些特殊的值可以非常优雅地表达如下：
我们还具有如下非常好的性质：$$FRAC{t_(n+1)}{t_n}$$这个序列，有时被称为广义的或延迟的Fibonacci序列，[Kak 2004，Wilson 1993]，\(d=2)的序列通常被称为Padovan序列[Stewart，1996，OEIS A000931]，而\(d=3)序列在[OEIS A079398]中被列出，我们还具有以下非常好的性质：$$FRAC=\Lim_{n\to\Ininfty}；\；\frac{t_(n+1)}}{t_n}.$$这个序列，有时称为广义或延迟Fibonacci序列，[Kak 2004，Wilson 1993]，并且\(d=2)的序列通常被称为Padovan序列[Stewart，1996，OEIS A000931]。如前所述，这篇文章的主要贡献是描述了这个广义序列与\(d-\)维低偏差序列的构造之间的显式联系。
主要结果：无参数的(d-无限)维开(φd)序列(R_d(φ_d))与其他已有方法相比，具有更好的低偏差特性。$$\mathbf{t}_n=\{n\pmb{\alpha}\}，n=1，2，3，…。$\texm{其中}\quad\pmb{\alpha}=(\frac{1}{\φ_d}，\frac{1}{\φ_d^2}，\frac{1}{\φ_d^3}，…。\frac{1}{\phi_d^d})，$\texm{and}\；\phi_d\\texm{是}的唯一正根：{x^{d+1}=x+1.$$。
对于二维，图1显示了这个广义序列\(n=150)。与(2，3)-Halton序列、基于\(\sqrt{3}，\sqrt{7})的Kronecker序列、Niederreiter和Sobol序列相比，\(R_2\)-序列的点明显分布得更均匀。(由于Niederreiter和Sobol序列的复杂性，它们是使用英特尔提供的专有代码通过MATHEMICA进行计算的。)。这类序列通常称为Korobov序列[Korobov 1959]，其中基向量$\pmb{\alpha}$是单个实值的函数。
总而言之，在1维中，第$n$项($n$=1，2，3，…)的伪代码。。)。是否定义为。
在2维中，第$n$项($n$=1，2，3，…)的$x$和$y$坐标的伪代码。。)。定义为
第$n$项($n$=1，2，3，…)的三维坐标$x$，$y$和$z$的伪代码。。)。定义为。
G=1.22074408460575947536a1=1.0/ga2=1.0/(g*g)a3=1.0/(g*g*g)x[n]=(0.5+a1*n)%1Y[n]=(0.5+a2*n)%1Z[n]=(0.5+a3*n)%1。
使用上面嵌套的根式公式g=phi_d#将numpy导入为np#，或者您可以直接对其进行硬编码。#φ(1)=1.6180339887498948482#φ(2)=1.32471795724474602596 defφ(D)：对于范围(10)内的i，x=2.0000：x=power(1+x，1/(d+1))返回x。
#维度个数。D=2#所需点数n=50g=φ(D)α=np.zeros(D)，j在范围(D)内：α[j]=power(1/g，j+1)%1z=np.zeros((n，d))#这个数字可以是任何实数。#常见默认设置通常为SEED=0#，但SEED=0.5通常更好。对于范围(N)内的i：Z[i]=(种子+Alpha*(i+1))%1打印(Z)。
我已经编写了这样的代码，以便与整篇文章中使用的数学约定保持一致。但是，出于编程约定和/或效率的原因，有一些修改值得考虑。首先，由于$R_2$是加性递归序列，因此对于非常大的$n$，不需要浮点乘法并保持较高精度的$z$的替换公式是。
其次，对于那些允许矢量化的语言，分数函数代码可以被矢量化如下：
最后，您可以将所有常量乘以$2^{32}$，然后修改frac(.)，从而将此浮点加法替换为整数加法。相应地发挥作用。以下是一些包含代码的演示，由其他人根据此序列进行演示：
新的$R_2$-序列是唯一一个最小填充距离仅下降$1/\sqrt{n}$的二维低偏差拟随机序列。
虽然量化差异的标准技术分析是通过分析\(d^*)-差异来实现的，但我们首先要提到其他几种几何方法(可能还有更直观的方法！)。这种新的序列如何比其他标准方法更可取。如果我们表示点\(i\)和\(j\)之间的距离为\(d_{ij}\)和\(d_0=tExtrm{inf}\；d_{ij}\)，则下图显示了\(d_0(N)\)在R序列、(2，3)-Halton、Sobol、Niederrieter和随机序列中的变化情况，如下图6所示。
与上图一样，每个最小距离度量都使用因子\(1/\sqrt{n}\)进行标准化。可以看到，在n=300个点之后，几乎可以肯定，对于随机分布(绿色)，将会有两个点彼此非常接近。还可以看出，虽然(2，3)-Halton序列(橙色)比随机抽样要好得多，但不幸的是它也渐近趋于零。对于Sobol序列，归一化$d_0$变为零的原因是因为Sobol自己证明了Sobol序列以$1/n$的速率下降-这是好的，但明显比$R_2$更糟糕，后者仅下降$1/\sqrt{n}$。
对于\(R(\φ_2)\)序列(蓝色)，两点之间的最小距离始终在$0.549/\sqrt{n}$和$0.868/\sqrt{n}$之间。请注意，最佳直径为0.868相当于59.2%的填充分数。将其与其他圆形填料进行比较。
还请注意，Bridson Poisson圆盘采样在\(n\)中不可扩展，在典型的推荐设置下，仍然只产生49.4%的包装分数。值得注意的是，这个d_0的概念将\(φ_d)低偏差序列与在\(d\)维上的数/向量的极差逼近紧密地联系在一起[Hensley，2001]。虽然人们对二维的不可逼近数知之甚少，但\(φ_d)的构造可能为研究高维的更高的不可逼近数提供了一个新的窗口。
另一种可视化点分布均匀性的方法是基于二维序列的前\(n\)点创建Voronoi图，然后基于其面积为每个区域上色。在下图中，给出了(I)$R_2$-序列；(Ii)(2，3)-Halton序列；(Iii)基于素数的递归；(Iv)简单随机抽样的基于颜色的Voronoi图。所有的数字都使用相同的色标。再次清楚的是，$R_2$-序列提供了比Halton或简单随机抽样更均匀的分布。这与上面相同，但是根据每个Voronoi单元格的顶点数进行着色。不仅很清楚，与Halton或简单的随机抽样相比，R-序列提供的分布要均匀得多，而且更引人注目的是，对于$n$的键值，它只由六边形组成！如果我们考虑广义Fibonacci数列，则$A_1=A_2=A_3=1；\quad A_{n+3}=A_{n+1}+A{n}$。即，$A_n：$，\Begin{array}{r}1&；1&；1&；2&；3&；4&；5&；7\\9&；\textbf{12}&；16&；21&；28&；37&；\textbf{49}&；65&；86\\114&；\textbf{151}&；200&；265&；351&；465&；\textbf{1897}&；2513&；3329&；4410&；5842&\textbf{7739}&；10252&；13581\\17991&；41824&textbf{23833}&；31572&；41824&；55405&；73396&；\textbf{97229}&；128801&；170625\\226030&；\textbf{299426}&；6655&；525456&；696081&；922111&；\textbf{1221537}&；1618192&；2143648\end{array}$n=A_{9k-2}$或$n=A_{9k+2}$的所有值将仅由六边形组成。
对于$n$的特定值，$R_2$序列的Voronoi网格仅由六边形组成。
$R$-序列是唯一一个可用于通过其Delaunay网格产生$d$维准周期平铺的低偏差准随机序列。
Delaunay三角剖分(有时拼写为“Delone Triangulation”)是Voronoi图的对偶，它提供了查看这些分布的另一种方式。然而，更重要的是，Delaunay三角剖分提供了一种创建平面准周期平铺的新方法。$R_2$序列的Delaunay三角剖分提供了比Halton或Random序列规则得多的模式。更具体地说，点分布的Delaunay三角剖分，其中$n$等于任意广义斐波那契数列$：A_N=$1，1，2，3，4，5，7，9，12，16，21，28，37，…。那么Delaunay三角剖分只由3个相同配对的三角形组成，即平行四边形(菱形)！(除了那些顶点与凸壳相同的三角形。)。此外，
对于$n=A_k$，$R_2$-序列的Delaunay三角剖分形成准周期瓦片，每个准周期瓦片只由三个基三角形(红、黄、蓝)组成，它们总是通过三个平行四边形(菱形)配对形成明确定义的平面准周期瓦片。
请注意这一点。
.