时间戳

首页 > 解决方案 > 如何在蒙特卡洛方法中定义准确的样本数

python - 如何在蒙特卡洛方法中定义准确的样本数

问题描述

我正在开发蒙特卡洛（dx）积分的模拟，我发现自己想知道在蒙特卡洛方法中确定样本的确切数量（N）以将解近似为定积分的最佳方法。

这是一个简单的实现代码：

import math
import random

class Montecarlo:

    def __init__ (self):
        print ("Inicializa..")

    def fy(self, ri, a, b):
        res = math.pow((b-a)*ri+a, 2.0)+math.sqrt((b-a)*ri+a)        
        return res

    def integral (self, a, b, N):
        suma = 0.0
        ri = 0.0
        for i in range (N):
            ri = random.random()
            suma+=self.fy(ri,a,b)

        res=((b-a)/N)*suma
        return res
    if __name__ == "__main__":
        monte = Montecarlo()
        res = monte.integral(10.0,27.0,N)
        print("Res: ", res)

标签： pythonmontecarlo

解决方案

赢得 Monte Carlo，您可以计算模拟的统计误差 (stddev)。它下降为 1/sqrt(N)。您可以设定您的目标 - 比如说，使误差低于 2% - 现在轻松计算您需要的许多样本 (N)。

我修改了您的代码并添加了第二动量、西格玛和模拟误差的计算

import math
import random

class Montecarlo:

    def __init__(self):
        print ("Inicializa..")

    def fy(self, ri, a, b):
        res = math.pow((b-a)*ri+a, 2.0) + math.sqrt((b-a)*ri+a)
        return res

    def integral (self, a, b, N):
        sum = 0.0
        var = 0.0
        for i in range(N):
            ri = random.random()
            v  = self.fy(ri, a, b)
            sum += v
            var += v*v

        sum /= float(N)
        var /= float(N)

        sigma = ( var - sum*sum ) * float(N)/float(N-1)
        error = sigma / math.sqrt(N)

        return ((b-a) * sum, (b-a)*error)

if __name__ == "__main__":
    N = 100000
    monte = Montecarlo()
    res, err = monte.integral(10.0, 27.0, N)
    print("Res: {0}, Err: {1}".format(res, err))

推荐阅读