Моделирование Монте-Карло с помощью Python: построение гистограммы на лету

У меня есть концептуальный вопрос о построении гистограммы на лету с помощью Python. Я пытаюсь выяснить, есть ли хороший алгоритм или, возможно, уже существующий пакет.

Я написал функцию, которая запускает моделирование Монте-Карло, вызывается 1 000 000 000 раз и возвращает 64-битное плавающее число в конце каждого запуска. Ниже указанная функция:

def MonteCarlo(df,head,span):
    # Pick initial truck
    rnd_truck = np.random.randint(0,len(df))
    full_length = df['length'][rnd_truck]
    full_weight = df['gvw'][rnd_truck]

    # Loop using other random trucks until the bridge is full
    while True:
        rnd_truck = np.random.randint(0,len(df))
        full_length += head + df['length'][rnd_truck]
        if full_length > span:
            break
        else:
            full_weight += df['gvw'][rnd_truck]

    # Return average weight per feet on the bridge
    return(full_weight/span)

df - это объект фрейма данных Pandas, имеющий столбцы, помеченные как 'length' и 'gvw', которые представляют собой длину и вес грузовика, соответственно. head - расстояние между двумя идущими подряд грузовиками, span - длина моста. Функция случайным образом размещает грузовики на мосту, если общая длина автопоезда меньше длины моста. Наконец, вычисляется средний вес грузовиков на мосту на фут (общий вес на мосту, деленный на длину моста).

В результате я хотел бы построить табличную гистограмму, показывающую распределение возвращаемых значений, которую можно построить позже. У меня было несколько идей:

  1. Продолжайте собирать возвращаемые значения в виде вектора numpy, а затем используйте существующие функции гистограммы после завершения анализа Монте-Карло. Это невозможно, поскольку, если мой расчет верен, мне потребуется 7,5 ГБ памяти только для этого вектора (1000000000 64-битных чисел с плавающей запятой ~ 7,5 ГБ)

  2. Инициализируйте массив numpy с заданным диапазоном и количеством ячеек. Увеличивайте количество предметов в соответствующей корзине на единицу в конце каждого прогона. Проблема в том, что я не знаю, какой диапазон значений я бы получил. Настройка гистограммы с диапазоном и подходящим размером ячейки неизвестна. Мне также нужно выяснить, как присвоить значения правильным ячейкам, но я думаю, что это выполнимо.

  3. Сделайте это как-нибудь на лету. Изменяйте диапазоны и размеры ячеек каждый раз, когда функция возвращает число. Думаю, это было бы слишком сложно написать с нуля.

Что ж, я уверен, что может быть лучший способ справиться с этой проблемой. Любые идеи приветствуются!

Во-вторых, я протестировал выполнение указанной выше функции 1000000000 раз только для получения наибольшего вычисленного значения (фрагмент кода ниже). И это займет около часа, когда span = 200. Время вычислений увеличилось бы, если бы я использовал его для более длинных пролетов (цикл while работает дольше, чтобы заполнить мост грузовиками). Как вы думаете, есть ли способ оптимизировать это?

max_w = 0
i = 1
    while i < 1000000000:
        if max_w < MonteCarlo(df_basic, 15., 200.):
            max_w = MonteCarlo(df_basic, 15., 200.)
    i += 1
print max_w

Спасибо!


person marillion    schedule 06.08.2013    source источник
comment
Присвоение значения ячейке - это просто двоичный поиск. Однако вы не можете изменить диапазон на лету, а это значит, что вы должны знать его заранее или хранить все. Или, по крайней мере, сделайте некоторые предположения: например, объедините данные в небольшие бункеры заданного размера (таким образом, вам не нужно хранить слишком много данных) и расширяйте список ящиков всякий раз, когда данные переполняют их.   -  person    schedule 07.08.2013
comment
@arbautjc спасибо за ответ. В конце я немного отредактировал пост, связанный с проблемами производительности, однако он имеет более низкий приоритет по сравнению с проблемой гистограммы, которая у меня есть. Я немного надеялся, что может быть научный пакет, способный на это.   -  person marillion    schedule 07.08.2013
comment
Я предлагаю вам быструю и грязную реализацию с использованием хеш-таблицы вместо отсортированных списков (что намного проще).   -  person    schedule 07.08.2013


Ответы (1)


Вот возможное решение с фиксированным размером ячейки и ячейками вида [k * size, (k + 1) * size [. Функция finalizebins возвращает два списка: один с количеством бункеров (a), а другой (b) с нижними границами бункера (верхняя граница выводится путем добавления размера бункера).

import math, random

def updatebins(bins, binsize, x):
    i = math.floor(x / binsize)
    if i in bins:
        bins[i] += 1
    else:
        bins[i] = 1

def finalizebins(bins, binsize):
    imin = min(bins.keys())
    imax = max(bins.keys())
    a = [0] * (imax - imin + 1)
    b = [binsize * k for k in range(imin, imax + 1)]
    for i in range(imin, imax + 1):
        if i in bins:
            a[i - imin] = bins[i]
    return a, b

# A test with a mixture of gaussian distributions

def check(n):
    bins = {}
    binsize = 5.0
    for i in range(n):
        if random.random() > 0.5:
            x = random.gauss(100, 50)
        else:
            x = random.gauss(-200, 150)
        updatebins(bins, binsize, x)
    return finalizebins(bins, binsize)

a, b = check(10000)

# This must be 10000
sum(a)

# Plot the data
from matplotlib.pyplot import *
bar(b,a)
show()

введите описание изображения здесь

person Community    schedule 06.08.2013