首页
学习
活动
专区
圈层
工具
发布
社区首页 >问答首页 >LCG是否像我的代码所显示的那样无法通过Kolmogorov-Smirnov测试?

LCG是否像我的代码所显示的那样无法通过Kolmogorov-Smirnov测试?
EN

Stack Overflow用户
提问于 2020-01-12 20:22:05
回答 2查看 140关注 0票数 1

我使用以下Python代码向学生演示随机变量的生成:

代码语言:javascript
复制
import numpy as np
import scipy.stats as stats

def lcg(n, x0, M=2**32, a=1103515245, c=12345):
    result = np.zeros(n)
    for i in range(n):
        result[i] = (a*x0 + c) % M
        x0 = result[i]

    return np.array([x/M for x in result])

x = lcg(10**6, 3)
print(stats.kstest(x, 'uniform'))

根据维基百科的说法,默认参数是glibc使用的参数。打印代码的最后一行

代码语言:javascript
复制
KstestResult(statistic=0.043427751892089805, pvalue=0.0)

P值为0.0表明,如果x的元素真正按照均匀分布分布,则观测基本上永远不会发生。我的问题是:我的代码中是否有bug,或者带有给定参数的LCG没有通过带有10**6副本的Kolmogorov-Smirnov测试?

EN

回答 2

Stack Overflow用户

回答已采纳

发布于 2020-01-13 04:13:56

你的代码有问题,它使得均匀分布像这样

我对你的LCG实现做了一些修改,现在一切都很好了(Python3.7,蟒蛇,Win10 x64)

代码语言:javascript
复制
import numpy as np
import scipy.stats as stats
import matplotlib.pyplot as plt

def lcg(n, x0, M=2**32, a=1103515245, c=12345):
    result = np.zeros(n)
    for i in range(n):
        x0 = (a*x0 + c) % M
        result[i] = x0

    return np.array([x/float(M) for x in result])

#x = np.random.uniform(0.0, 1.0, 1000000)
x = lcg(1000000, 3)
print(stats.kstest(x, 'uniform'))

count, bins, ignored = plt.hist(x, 15, density=True)
plt.plot(bins, np.ones_like(bins), linewidth=2, color='r')
plt.show()

哪种打印

代码语言:javascript
复制
KstestResult(statistic=0.0007238884545415214, pvalue=0.6711878724246786)

和情节

更新

正如@pjs指出的,你最好在循环中除以浮点数(M),不需要对整个数组进行第二次传递

代码语言:javascript
复制
def lcg(n, x0, M=2**32, a=1103515245, c=12345):
    result = np.empty(n)
    for i in range(n):
        x0 = (a*x0 + c) % M
        result[i] = x0 / float(M)

    return result
票数 2
EN

Stack Overflow用户

发布于 2020-01-13 18:01:55

为了补充谢韦林的回答,我的代码不能正常工作的原因是result是一个浮点数的数组。我们可以在第二次迭代中看到两个实现之间的差异。在第一次迭代之后,x0 = 3310558080

代码语言:javascript
复制
In [9]: x0 = 3310558080

In [10]: float_x0 = float(x0)

In [11]: (a*x0 + c) % M
Out[11]: 465823161

In [12]: (a*float_x0 + c) % M
Out[12]: 465823232.0

In [13]: a*x0
Out[13]: 3653251310737929600

In [14]: a*float_x0
Out[14]: 3.6532513107379297e+18

所以这个问题与浮点数的使用有关。

票数 1
EN
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接:

https://stackoverflow.com/questions/59703748

复制
相关文章

相似问题

领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档