Mô phỏng Monte Carlo với Python: xây dựng một biểu đồ trên bay

Question

Tôi có một câu hỏi khái niệm về việc xây dựng một biểu đồ khi đang bay với Python. Tôi đang cố gắng tìm ra nếu có một thuật toán tốt hoặc có thể là một gói hiện có.

Tôi đã viết một hàm, chạy mô phỏng Monte Carlo, được gọi là 1.000.000.000 lần và trả về số nổi 64 bit ở cuối mỗi lần chạy. Dưới đây là chức năng cho biết:

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 là một đối tượng dataframe Pandas có cột dán nhãn là 'length' và 'gvw', đó là chiều dài xe tải và trọng lượng, tương ứng. head là khoảng cách giữa hai xe tải liên tiếp, span là chiều dài cầu. Chức năng này đặt ngẫu nhiên các xe tải trên cầu miễn là tổng chiều dài của xe tải nhỏ hơn chiều dài cầu. Cuối cùng, tính toán trọng lượng trung bình của các xe tải có sẵn trên cầu/foot (tổng trọng lượng hiện có trên cầu chia cho chiều dài cầu).

Kết quả là tôi muốn xây dựng một biểu đồ dạng bảng cho thấy sự phân bố của các giá trị trả về, có thể được vẽ sau này. Tôi đã có một số ý tưởng trong tâm trí:

1. Giữ thu thập các giá trị trả lại trong một vector NumPy, sau đó sử dụng chức năng histogram hiện khi phân tích Montecarlo được hoàn tất. Điều này sẽ không khả thi, vì nếu tính toán của tôi là chính xác, tôi sẽ cần 7,5 GB bộ nhớ cho vectơ đó (1.000.000.000 64 bit nổi ~ 7,5 GB)

2. Khởi tạo một mảng có nhiều ô và phạm vi nhất định . Tăng số lượng các mục trong thùng phù hợp bằng một ở cuối mỗi lần chạy. Vấn đề là, tôi không biết phạm vi giá trị tôi sẽ nhận được. Thiết lập một biểu đồ với phạm vi và kích thước thùng thích hợp là không xác định. Tôi cũng phải tìm ra cách gán các giá trị cho các thùng đúng, nhưng tôi nghĩ nó có thể thực hiện được.

3. Thực hiện điều đó một cách nhanh chóng. Sửa đổi phạm vi và kích thước thùng mỗi khi hàm trả về một số. Điều này sẽ là quá khó khăn để viết từ đầu tôi nghĩ.

Vâng, tôi đặt cược có thể có cách tốt hơn để xử lý vấn đề này. Bất cứ ý tưởng sẽ được hoan nghênh!

Trên lưu ý thứ hai, tôi đã thử nghiệm chạy hàm trên chỉ với 1.000.000.000 lần để nhận giá trị lớn nhất được tính toán (đoạn mã bên dưới). Và điều này mất khoảng một giờ khi span = 200. Thời gian tính toán sẽ tăng nếu tôi chạy nó trong thời gian dài hơn (vòng lặp while chạy lâu hơn để lấp đầy cầu bằng xe tải). Có cách nào để tối ưu hóa điều này bạn nghĩ không?

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 

Cảm ơn!

Answer 1

Dưới đây là giải pháp khả thi, với kích thước thùng cố định và thùng có dạng [k * size, (k + 1) * size [. Hàm finalizebins trả về hai danh sách: một với số đếm bin (a) và giá trị còn lại (b) có giới hạn bin thấp hơn (giới hạn trên được suy ra bằng cách thêm binsize).

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()