YCSB的zipfian分布生成原理出自这篇论文：Quickly
Generating Billion-Record Synthetic Databases
但是里面似乎并没有给出推导。我最开始是尝试从定义出发推导出论文里的形式，但是失败了。后来我尝试从论文里的形式出发逆向推导出原始式子，才成功了。为了方便理解，这里从我逆向推导出的原始公式出发，正向推导出论文里的最终形式。
假如随机变量X满足取值范围为1到n，参数为θ的Zipfian分布，即
X ∼ Zipfian(n,θ)
那么X只能取1到n的正整数，而且取正整数i的概率正比于
Riemann zeta function的原本定义是:
，但论文里给这个定义加了个参数n：
这样Zipfian(n,θ)取i的概率就可以表示为：
分布函数
现在，我们从0到1的均匀分布中取一个数u。只要我们算出r，使其满足F(n,θ,r−1) ≤ u < F(n,θ,r)，那么r就是满足Zipfian(n,θ)分布的。
但是r是没有解析解的，因为F并不是一个连续函数。然而，我们观察发现，虽然当i比较小时，F(n,θ,i)的变化幅度比较大，但是当i变大之后，F(n,θ,i)随i的变化幅度越来越小趋近于0。所以我们可以认为当i足够大之后，F(n,θ,i)就近似一个连续函数了。我们的求解思路就是直接算出F的前几项，然后当u比较小时，直接根据这些算出来的项把r求出来。当u超过了预计算的最大值之后，就用跟F近似的连续函数把r估计出来。
论文里预计算了前2项：
如果 uζ(θ,n) < 1，那么就取结果r = 1，如果，就取结果r = 2。
否则，取r > 2。u此时的值域为。在进行接下来的步骤前，我们先把它的值域变换回0到1：
接下来我们使用跟F近似的连续函数来估计r应该是多少。有两个选择：估计，或者估计。显然，i比较小时，前者随着i的增大变化更大，用连续函数近似的话误差较大，而且这个误差随着i的增大并不会消失。因此我们选择估计后者。
我们定义一个反向的分布函数G(n,θ,i)，表示在已知X ≥ 3的情况下，X ≥ i的概率：
其中，
我们假设i和n都足够大，这样就有
i取最小值3时，该值为
那么G的近似连续函数：
接下来，我们只需要求出r，使其满足G̃(n,θ,r+1) < u′ ≤ G̃(n,θ,r)即可。
令G̃(n,θ,i) = u′，有
定义η:
定义，则
i = 1 + n(1−η(1−u))α = 1 + n(ηu−η+1)α
由于G̃是关于i的减函数，所以G̃(n,θ,r+1) < u′ ≤ G̃(n,θ,r)等价于r < = 1 + n(ηu−η+1)α < r + 1，即
r = ⌊1 + n(ηu−η+1)α⌋
这就是论文里的公式。
¶ 评估
我们估计的G̃(n,θ,i)在i比较大的时候比较接近G(n,θ,i)，但当i很小时误差就会比较大。这里我们来测试一下误差究竟有多大。
为了减小浮点误差造成的影响，这里将G̃的形式变换成这样：
#include <cstdio>
#include <cmath>
constexpr size_t n = 20000000;
constexpr double theta = 0.99;
double zeta[n + 1];
double G(size_t i) {
return (zeta[n] - zeta[i - 1]) / (zeta[n] - zeta[2]);
}
double G_tilde(double i) {
return (pow(n, 1 - theta) - pow(i - 1, 1 - theta)) /
(pow(n, 1 - theta) - pow(2, 1 - theta));
}
void space(size_t num) {
while (num) {
num -= 1;
putchar(' ');
}
}
int main() {
zeta[0] = 0;
for (size_t i = 1; i <= n; ++i) {
zeta[i] = zeta[i - 1] + pow((double)i, -theta);
}
size_t v[] = {3, 4, 5, 10, 100, 1000, 10000, 100000, 1000000, 10000000};
printf("i\t\tG(i)\t\tG_tilde(i)\n");
for (size_t i : v) {
size_t ret = printf("%zu", i);
space(16 - ret);
ret = printf("%f", G(i));
space(16 - ret);
ret = printf("%f", G_tilde(i));
putchar('\n');
}
return 0;
}
|
i G(i) G_tilde(i)
3 1.000000 1.000000
4 0.980609 0.976770
5 0.966024 0.960231
10 0.922307 0.913352
100 0.782473 0.772490
1000 0.641863 0.633459
10000 0.498228 0.491684
100000 0.351273 0.346658
1000000 0.200898 0.198258
10000000 0.047020 0.046402
|
可以看到误差其实并不大。
为了进一步验证估计的准确性，我们接下来直接验证F的估计是否准确。
当i ≥ 3时，有
接下来我们对比i ≥ 3时，F和F̃的差异：
#include <cstdio>
#include <cmath>
constexpr size_t n = 20000000;
constexpr double theta = 0.99;
double zeta[n + 1];
double F(size_t i) {
return zeta[i] / zeta[n];
}
double F_tilde(double i) {
return 1 - (pow(n, 1 - theta) - pow(i, 1 - theta)) /
(pow(n, 1 - theta) - pow(2, 1 - theta)) *
(1 - zeta[2] / zeta[n]);
}
void space(size_t num) {
while (num) {
num -= 1;
putchar(' ');
}
}
int main() {
zeta[0] = 0;
for (size_t i = 1; i <= n; ++i) {
zeta[i] = zeta[i - 1] + pow((double)i, -theta);
}
size_t v[] = {3, 4, 5, 10, 100, 1000, 10000, 100000, 1000000, 10000000};
printf("i\t\tF(i)\t\tF_tilde(i)\n");
for (size_t i : v) {
size_t ret = printf("%zu", i);
space(16 - ret);
ret = printf("%f", F(i));
space(16 - ret);
ret = printf("%f", F_tilde(i));
putchar('\n');
}
return 0;
}
|
i F(i) F_tilde(i)
3 0.097466 0.100999
4 0.110890 0.116222
5 0.121653 0.128059
10 0.156545 0.164999
100 0.280381 0.289565
1000 0.409298 0.417032
10000 0.541446 0.547469
100000 0.676695 0.680943
1000000 0.815097 0.817527
10000000 0.956724 0.957292
|
误差也很小。
¶ 使用
#include <iostream>
#include <cmath>
#include <random>
#include <cassert>
class ZipfianGenerator {
public:
ZipfianGenerator(
size_t n, double theta
) : n_(n),
alpha_(1 / (1 - theta)),
zeta2_(1 + pow(2.0, -theta))
{
assert(n > 2);
zetan_ = zeta2_;
for (size_t i = 3; i <= n; ++i) {
zetan_ += pow((double)i, -theta);
}
eta_ = (1 - pow(2.0 / n, 1 - theta)) / (1 - zeta2_ / zetan_);
}
template <typename E>
size_t operator()(E &e) {
std::uniform_real_distribution<> dis(0.0, 1.0);
double u = dis(e);
double uz = u * zetan_;
if (uz < 1)
return 1;
if (uz < zeta2_)
return 2;
return 1 + (size_t)(n_ * pow(eta_ * u - eta_ + 1, alpha_));
}
private:
size_t n_;
double alpha_;
double zeta2_;
double zetan_;
double eta_;
};
int main() {
std::random_device rd;
std::mt19937 e(rd());
size_t n = 20000000;
ZipfianGenerator g(n, 0.99);
for (size_t i = 0; i < n; ++i) {
std::cout << g(e) << std::endl;
}
return 0;
}
|
¶ 注意事项
YCSB里的zipfian似乎默认是scrambled zipfian，它的item
count被指定为了1e10，导致它生成的分布跟zipfian的理论分布不一致。建议不要用YCSB生成zipfian，而是用它的公式自己写一个。