[PYTHON] Simulation of the contents of the wallet

Distribution of bills and coins

If you always pay to minimize the weight of your wallet, what is the distribution of your wallet?

Simulation

Let's check it by simulation.

――10,000 yen is infinite. ――Buy the product 10,000 times and see the distribution. (Ignore the first 100 times) --The price of the product is 100 yen to 9999 yen, [Benford's Law](https://ja.wikipedia.org/wiki/%E3%83%99%E3%83%B3%E3%83%95 % E3% 82% A9% E3% 83% BC% E3% 83% 89% E3% 81% AE% E6% B3% 95% E5% 89% 87). --Use Walker's Alias Method to generate random numbers (rand_from_prob). --Minimize weight after payment is calculated by Combinatorial Optimization (http://qiita.com/Tsutomu-KKE@github/items/bfbf4c185ed7004b5721).

Calculated in Python

First, let's make a weight (wgt) of 100 yen to 9999 yen according to Benford's law.

python


import numpy as np, matplotlib.pyplot as plt
from math import log
#Benford's law
wgt = np.array([log((i+1)/i, 10000) for i in range(100, 10000)])
wgt /= sum(wgt)
plt.plot(wgt)
plt.xlabel('Price')

image

Defines a change that returns the number after minimizing the weight.

python


from pulp import *
money_val = (1, 5, 10, 50, 100, 500, 1000, 5000)
money_wgt = (1.0, 3.7, 4.5, 4.0, 4.8, 7.0, 1.0, 1.0)
def change(price):
    m = LpProblem() #Mathematical model
    x = [LpVariable('x%d'%i, lowBound=0, cat=LpInteger)
         for i in range(len(money_val))] #Quantity after payment
    m += lpDot(money_wgt, x) #Objective function(Weight after payment)
    m += lpDot(money_val, x) == price #Amount after payment
    m.solve()
    return [int(value(i)) for i in x]

Let's simulate. Let's try the distribution of 1000 yen bills.

python


price = 0 #Current money
warm, nrun = 100, 10000
res = []
for i, p in enumerate(rand_from_prob(wgt, warm+nrun)):
    price -= p
    if price < 0:
        price += 10000
    if price:
        res.append(change(price))
a = np.array(res[-nrun:])
plt.hist(a[:,6], bins=5, range=(0, 5)) #Distribution of 1000 yen bills

image

It's an equal probability. The same was true for other coins and 5000 yen bills.

Conclusion

It seems to be equal probability. The distribution of the total amount also has an equal probability of 0 yen to 9999 yen.

python


import pandas as pd
from itertools import product
r2, r5 = range(2), range(5)
ptn = [np.dot(money_val, n) for nn in 
       product(r5, r2, r5, r2, r5, r2, r5, r2)]
plt.hist(ptn)
print(pd.DataFrame(ptn).describe())
>>>
                 0
count  10000.00000
mean    4999.50000
std     2886.89568
min        0.00000
25%     2499.75000
50%     4999.50000
75%     7499.25000
max     9999.00000

image

that's all

Recommended Posts