In [4]:
%load_ext autoreload
%autoreload 2
%matplotlib inline
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
import random

In [5]:
class GaussianSampler(object):
	def __init__(self, mean, var):
		self.mean = mean
		self.var = var

	def sample(self):
		return random.normalvariate(self.mean, self.var)

For each product, we'll model its usage time using a Gaussian distribution in days.


In [6]:
products = {"food" : GaussianSampler(30.0, 5.0),
            "poop bags" : GaussianSampler(14.0, 2.0)}

The customer will have an inventory with a remaining time (given in days) for each item.


In [6]:

Probability of a purchase given the number of days remaining will be modeled using an exponential distribution.


In [7]:
class ExponentialPDF(object):
    def __init__(self, avg_rate):
        self.lambd = 1.0 / avg_rate
    
    def probability(self, t):
        return self.lambd * np.exp(-self.lambd * t)

In [8]:
customer_inventory = {"food" : 0,
                      "poop bags" : 0}
item_purchases = {"food": [],
                  "poop bags" : []}
pdf = ExponentialPDF(1.0)
for t in xrange(1, 101):
    items = customer_inventory.keys()
    print "Day", t
    for item in items:
        customer_inventory[item] = max(0, customer_inventory[item] - 1)
        prob = pdf.probability(customer_inventory[item])
        print item, customer_inventory[item], prob
        prob = pdf.probability(customer_inventory[item])
        r = random.random()
        if r < prob:
            additional_days = products[item].sample()
            additional_days = max(0, additional_days)
            print "Purchased", additional_days
            customer_inventory[item] += additional_days
            item_purchases[item].append(1)
        else:
            item_purchases[item].append(0)


Day 1
food 0 1.0
Purchased 32.243591031
poop bags 0 1.0
Purchased 11.7001029395
Day 2
food 31.243591031 2.69824151729e-14
poop bags 10.7001029395 2.2542617268e-05
Day 3
food 30.243591031 7.33458088524e-14
poop bags 9.7001029395 6.12771868856e-05
Day 4
food 29.243591031 1.99374579397e-13
poop bags 8.7001029395 0.00016656866361
Day 5
food 28.243591031 5.41956296232e-13
poop bags 7.7001029395 0.000452780571482
Day 6
food 27.243591031 1.47318995187e-12
poop bags 6.7001029395 0.00123078519974
Day 7
food 26.243591031 4.00454547603e-12
poop bags 5.7001029395 0.00334562104319
Day 8
food 25.243591031 1.08854831987e-11
poop bags 4.7001029395 0.00909434088661
Day 9
food 24.243591031 2.95898111731e-11
poop bags 3.7001029395 0.0247209815739
Day 10
food 23.243591031 8.04334460193e-11
poop bags 2.7001029395 0.067198594994
Day 11
food 22.243591031 2.18640774715e-10
poop bags 1.7001029395 0.18266471967
Purchased 14.1731832483
Day 12
food 21.243591031 5.94327244867e-10
poop bags 14.8732861878 3.47227314043e-07
Day 13
food 20.243591031 1.61554894988e-09
poop bags 13.8732861878 9.43861698108e-07
Day 14
food 19.243591031 4.39151735345e-09
poop bags 12.8732861878 2.56568210255e-06
Day 15
food 18.243591031 1.19373818212e-08
poop bags 11.8732861878 6.97424703695e-06
Day 16
food 17.243591031 3.2449168084e-08
poop bags 10.8732861878 1.89579689877e-05
Day 17
food 16.243591031 8.82059839515e-08
poop bags 9.87328618778 5.15331026039e-05
Day 18
food 15.243591031 2.39768723337e-07
poop bags 8.87328618778 0.000140081496372
Day 19
food 14.243591031 6.51758963679e-07
poop bags 7.87328618778 0.000380780986092
Day 20
food 13.243591031 1.7716645475e-06
poop bags 6.87328618778 0.00103507003512
Day 21
food 12.243591031 4.8158835456e-06
poop bags 5.87328618778 0.00281361206764
Day 22
food 11.243591031 1.309092873e-05
poop bags 4.87328618778 0.0076481905558
Day 23
food 10.243591031 3.55848336844e-05
poop bags 3.87328618778 0.0207899374084
Day 24
food 9.24359103098 9.6729606773e-05
poop bags 2.87328618778 0.0565129090721
Day 25
food 8.24359103098 0.000262938332365
poop bags 1.87328618778 0.153618013804
Day 26
food 7.24359103098 0.000714740490873
poop bags 0.873286187775 0.417577055447
Purchased 14.5560301377
Day 27
food 6.24359103098 0.0019428660884
poop bags 14.4293163255 5.41286931911e-07
Day 28
food 5.24359103098 0.00528125758324
poop bags 13.4293163255 1.471370431e-06
Day 29
food 4.24359103098 0.0143559465199
poop bags 12.4293163255 3.99959950551e-06
Day 30
food 3.24359103098 0.0390235085555
poop bags 11.4293163255 1.08720386569e-05
Day 31
food 2.24359103098 0.106076894189
poop bags 10.4293163255 2.95532651195e-05
Day 32
food 1.24359103098 0.288346893893
poop bags 9.42931632548 8.03341035459e-05
Day 33
food 0.24359103098 0.783808121963
poop bags 8.42931632548 0.000218370733874
Day 34
food 0 1.0
Purchased 36.5529582129
poop bags 7.42931632548 0.000593593197758
Day 35
food 35.5529582129 3.6269906232e-16
poop bags 6.42931632548 0.00161355360296
Day 36
food 34.5529582129 9.85918270305e-16
poop bags 5.42931632548 0.00438609343818
Day 37
food 33.5529582129 2.68000371851e-15
poop bags 4.42931632548 0.0119226380909
Day 38
food 32.5529582129 7.28500540824e-15
poop bags 3.42931632548 0.0324090904698
Day 39
food 31.5529582129 1.98026978214e-14
poop bags 2.42931632548 0.0880970417011
Day 40
food 30.5529582129 5.38293136425e-14
poop bags 1.42931632548 0.239472587597
Day 41
food 29.5529582129 1.46323245113e-13
poop bags 0.429316325479 0.650953983279
Purchased 15.7805632008
Day 42
food 28.5529582129 3.97747818271e-13
poop bags 15.2098795263 2.47989476222e-07
Day 43
food 27.5529582129 1.08119066672e-12
poop bags 14.2098795263 6.74105286863e-07
Day 44
food 26.5529582129 2.93898094244e-12
poop bags 13.2098795263 1.83240815175e-06
Day 45
food 25.5529582129 7.98897849001e-12
poop bags 12.2098795263 4.98100178122e-06
Day 46
food 24.5529582129 2.17162950573e-11
poop bags 11.2098795263 1.35397666294e-05
Day 47
food 23.5529582129 5.90310102358e-11
poop bags 10.2098795263 3.68049015903e-05
Day 48
food 22.5529582129 1.6046292244e-10
poop bags 9.20987952627 0.000100046095191
Day 49
food 21.5529582129 4.36183446209e-10
poop bags 8.20987952627 0.000271953482566
Day 50
food 20.5529582129 1.1856695357e-09
poop bags 7.20987952627 0.000739246209846
Day 51
food 19.5529582129 3.22298395346e-09
poop bags 6.20987952627 0.00200947953898
Day 52
food 18.5529582129 8.76097871412e-09
poop bags 5.20987952627 0.00546233171548
Day 53
food 17.5529582129 2.38148092381e-08
poop bags 4.20987952627 0.0148481570432
Day 54
food 16.5529582129 6.47353632001e-08
poop bags 3.20987952627 0.0403614754766
Day 55
food 15.5529582129 1.75968961446e-07
poop bags 2.20987952627 0.109713865358
Day 56
food 14.5529582129 4.7833323027e-07
poop bags 1.20987952627 0.298233206532
Day 57
food 13.5529582129 1.30024452779e-06
poop bags 0.209879526271 0.81068190596
Purchased 14.2720394249
Day 58
food 12.5529582129 3.53443107245e-06
poop bags 13.4819189512 1.39597292122e-06
Day 59
food 11.5529582129 9.60757975818e-06
poop bags 12.4819189512 3.79464782477e-06
Day 60
food 10.5529582129 2.61161094721e-05
poop bags 11.4819189512 1.03149222275e-05
Day 61
food 9.55295821288 7.09909458082e-05
poop bags 10.4819189512 2.80388656529e-05
Day 62
food 8.55295821288 0.000192973397975
poop bags 9.48191895122 7.62175389949e-05
Day 63
food 7.55295821288 0.000524556081093
poop bags 8.48191895122 0.00020718075126
Day 64
food 6.55295821288 0.00142589126324
poop bags 7.48191895122 0.000563175671356
Day 65
food 5.55295821288 0.00387597431023
poop bags 6.48191895122 0.00153087019368
Day 66
food 4.55295821288 0.0105359905351
poop bags 5.48191895122 0.0041613366292
Day 67
food 3.55295821288 0.0286397916163
poop bags 4.48191895122 0.0113116857413
Day 68
food 2.55295821288 0.0778510251214
poop bags 3.48191895122 0.0307483497997
Day 69
food 1.55295821288 0.211621026915
poop bags 2.48191895122 0.0835826805156
Day 70
food 0.552958212883 0.575245591982
poop bags 1.48191895122 0.22720128162
Day 71
food 0 1.0
Purchased 36.3403500401
poop bags 0.48191895122 0.617597115229
Day 72
food 35.3403500401 4.48622435801e-16
poop bags 0 1.0
Purchased 11.0886620455
Day 73
food 34.3403500401 1.21948221508e-15
poop bags 10.0886620455 4.15479637339e-05
Day 74
food 33.3403500401 3.31489634537e-15
poop bags 9.08866204547 0.000112939074827
Day 75
food 32.3403500401 9.01082249884e-15
poop bags 8.08866204547 0.000307000234826
Day 76
food 31.3403500401 2.44939550581e-14
poop bags 7.08866204547 0.00083451315966
Day 77
food 30.3403500401 6.65814729415e-14
poop bags 6.08866204547 0.00226844195751
Day 78
food 29.3403500401 1.80987208009e-13
poop bags 5.08866204547 0.00616626455202
Day 79
food 28.3403500401 4.91974238714e-13
poop bags 4.08866204547 0.0167616448812
Day 80
food 27.3403500401 1.33732463317e-12
poop bags 3.08866204547 0.0455628746957
Day 81
food 26.3403500401 3.63522524909e-12
poop bags 2.08866204547 0.123852734338
Day 82
food 25.3403500401 9.88156673694e-12
poop bags 1.08866204547 0.336666637155
Day 83
food 24.3403500401 2.68608832977e-11
poop bags 0.0886620454721 0.915154802028
Purchased 12.0726256171
Day 84
food 23.3403500401 7.30154509646e-11
poop bags 11.1612876625 1.42139360516e-05
Day 85
food 22.3403500401 1.98476573554e-10
poop bags 10.1612876625 3.863748408e-05
Day 86
food 21.3403500401 5.39515263266e-10
poop bags 9.16128766254 0.000105027570872
Day 87
food 20.3403500401 1.46655453631e-09
poop bags 8.16128766254 0.000285494537388
Day 88
food 19.3403500401 3.9865085465e-09
poop bags 7.16128766254 0.000776054613107
Day 89
food 18.3403500401 1.0836453741e-08
poop bags 6.16128766254 0.0021095351527
Day 90
food 17.3403500401 2.9456535289e-08
poop bags 5.16128766254 0.00573431107208
Day 91
food 16.3403500401 8.00711646054e-08
poop bags 4.16128766254 0.015587473586
Day 92
food 15.3403500401 2.1765599173e-07
poop bags 3.16128766254 0.0423711462004
Day 93
food 14.3403500401 5.91650327176e-07
poop bags 2.16128766254 0.115176716767
Day 94
food 13.3403500401 1.60827233316e-06
poop bags 1.16128766254 0.31308277625
Purchased 17.3842393453
Day 95
food 12.3403500401 4.37173745845e-06
poop bags 17.5455270079 2.39924411603e-08
Day 96
food 11.3403500401 1.18836144921e-05
poop bags 16.5455270079 6.52182168263e-08
Day 97
food 10.3403500401 3.23030133303e-05
poop bags 15.5455270079 1.77281493683e-07
Day 98
food 9.34035004013 8.78086941403e-05
poop bags 14.5455270079 4.81901062802e-07
Day 99
food 8.34035004013 0.000238688777662
poop bags 13.5455270079 1.30994290213e-06
Day 100
food 7.34035004013 0.000648823366976
poop bags 12.5455270079 3.56079398718e-06

In [9]:
plt.hold(True)
plt.plot(item_purchases["food"], "r-")
plt.plot(item_purchases["poop bags"], "g-")


Out[9]:
[<matplotlib.lines.Line2D at 0x104ee6050>]

Model 2

In the second version of the model, we'll have 4 structures:

  • Item definition -- gives amount, category
  • Item usage model -- functions for modeling item when an item will run out and probabilities
  • Customer inventory state -- for each category, gives an item usage model

In [136]:
class ExhaustableItem(object):
    def __init__(self, usage_rate=None, amount_used_average=None, amount_used_variance=None, average_time_for_purchase_before_exhaustion=None):
        """
        usage_rate is given in times/day -- used to determine when an item is used
        
        amount_used_average and amount_used_variance are given in units/day -- used to determine how
        much is used per usage.
        """
        
        self.usage_rate = usage_rate
        self.amount_used_average = amount_used_average
        self.amount_used_variance = amount_used_variance
        
        self.current_amount = 0.0
        self.last_amount_increase_time = 0.0
        self.exhaustion_time = 0.0
        
        self.epsilon = 1e-5
        
        self.average_time_for_purchase_before_exhaustion = average_time_for_purchase_before_exhaustion
        
    def increase_amount_and_simulate(self, purchase_time, purchased_amount):
        """
        Increase current amount, from a purchase.
        
        purchase_time -- given in seconds since start of model
        
        purchased_amount -- given in units
        """
        
        # simulate forward to current time
        # record amount left
        # add new amount
        # simulate forward until we run out of product
        # record run out time
        # (discard results of second simulation)
        
        # if we don't have any, don't bother simulating
        previous_amount = 0.0
        previous_time = self.last_amount_increase_time
        if self.current_amount != 0.0:
            # we want the last step that does not go past time
            new_time = self.last_amount_increase_time
            new_amount = self.current_amount
            while new_time < purchase_time:
                previous_time = new_time
                previous_amount = new_amount
                new_time, new_amount = self._simulate_usage_step(previous_time, previous_amount)
        
        self.current_amount = previous_amount + purchased_amount
        self.last_amount_increase_time = purchase_time
        
        # start next simulation from last accepted simulation step
        new_time = previous_time
        new_amount = self.current_amount
        while new_amount > self.epsilon:
            previous_time = new_time
            previous_amount = new_amount
            new_time, new_amount = self._simulate_usage_step(previous_time, previous_amount)
        
        self.exhaustion_time = new_time
        
    
    def _simulate_usage_step(self, start_time, initial_amount):
        """
        Simulate 1 step of usage.
        
        da/dt = -R(amount_used_average, amount_used_variance)
        
        timestep is determined from exponential distribution:
        f = \lambda \exp (-\lambda t) where \lambda = 1.0 / usage_rate
        \Delta t sampled from f
        
        Returns time after 1 step and remaining amount
        """
        
        # given in days since last usage
        timestep = random.expovariate(1.0 / self.usage_rate)
        
        # Might be more realistic to model usage
        # in units/time used since items are used
        # in discrete rather than continuous
        # quantities
        
        # given in units/day
        usage_rate = random.normalvariate(self.amount_used_average, self.amount_used_variance)
        
        # can't use a negative amount :)
        if usage_rate < 0.0:
            usage_rate = 0.0
        
        # given in units
        usage_amount = usage_rate * np.sqrt(timestep)
        
        amount = initial_amount - min(usage_amount, initial_amount)
        
        time = start_time + timestep
        
        return time, amount
        
    def probability_purchase(self, time):
        """
        Probability of a purchase at given time.
        
        time is seconds since model start.
        
        Probability is given by an exponential distribution:
        
        f = \lambda \exp (-\lambda t) where \lambda = 1.0 / average_time_for_purchase_before_exhaustion
        
        t is the amount of time from the curent time until the item
        is exhausted
        """
        
        time_until_exhaustion = max(self.exhaustion_time - time, 0.0)
        
        lambd = 1.0 / self.average_time_for_purchase_before_exhaustion
        prob = lambd * np.exp(-lambd * time_until_exhaustion)
        
        return prob

In [137]:
model = ExhaustableItem(usage_rate=1.0, amount_used_average=0.5, amount_used_variance=0.2, average_time_for_purchase_before_exhaustion=3.0)

In [138]:
previous_time = 0.0
previous_amount = 30.0 #lbs
times = [0.0]
amounts = [30.0]
while previous_amount > model.epsilon:
    new_time, new_amount = model._simulate_usage_step(previous_time, previous_amount)
    print new_time, new_amount
    times.append(new_time)
    amounts.append(new_amount)
    previous_time = new_time
    previous_amount = new_amount


0.626878287749 29.8168308891
1.56705742335 29.2447620842
2.39518957935 28.9321363204
3.75359710358 28.0290567534
3.92201625002 27.6856868558
6.06366986417 26.8380050901
6.72668934234 26.6166431218
6.89448198508 26.3065004094
8.30005214079 25.7504684377
10.2370631866 24.6669629332
10.6632313421 24.2062916104
12.0845569966 24.0278605594
13.6235708409 23.1317790939
13.9648412891 22.8262337662
16.8279596189 22.1710921672
18.2023718667 21.6381005816
19.6136615795 21.1776473341
21.7790766644 20.2028591676
23.030503414 20.0163840614
23.1727963577 19.7848670195
23.26190144 19.619108441
24.1930940802 19.545414062
25.6048730506 18.9939460198
25.8961239801 18.607928513
26.5606189317 18.1271970346
26.6122193097 18.0125411531
27.2699773581 17.6375590988
27.3197903426 17.6283093293
27.3802304865 17.4608446683
27.934783539 16.7493324081
29.2875091713 16.3344907227
29.6274366382 16.20113079
29.8079998085 15.8202400107
31.1563830016 15.2320688864
32.0513788702 14.4771392186
33.2720712628 14.0298955338
35.5020287581 13.9377285531
35.9609062079 13.6908064612
36.0788683119 13.6229602464
36.6638674863 13.3412709886
38.3248874413 12.8875706372
40.3566951438 12.1311396632
41.8160364646 11.660182914
43.8093119443 10.9437385783
45.8867807529 10.2419721673
45.8939201638 10.2034217963
46.341297134 9.58055972714
47.4368628708 9.17178482498
47.5537494897 8.80324010853
48.8991982051 8.09073691746
51.0143621705 7.19469650099
51.0981962368 6.95095174078
51.6010923206 6.66523407481
52.162315021 5.96803034534
52.5859567271 5.59012158885
54.7662270234 4.8175634607
55.6423602197 4.33172210163
56.3594972948 4.01446269195
56.9471171841 3.68843966752
57.87725546 3.03316336344
58.9348337409 2.26935248884
61.332033861 1.41413280961
61.8753739941 1.19471742985
65.281858093 0.0

In [139]:
model.increase_amount_and_simulate(0.0, 30.0)

In [140]:
print model.exhaustion_time


77.8627941012

In [141]:
model.exhaustion_time = 0.0
model.current_amount = 0.0
model.increase_amount_and_simulate(0.0, 30.0)
print model.exhaustion_time


74.9227887094

In [142]:
times = range(0, 100)
probs = [model.probability_purchase(float(t)) for t in times]
plt.plot(times, probs)


Out[142]:
[<matplotlib.lines.Line2D at 0x1055fda10>]

In [145]:
poopbag_model = ExhaustableItem(usage_rate=1.0, amount_used_average=2.5, amount_used_variance=0.5, average_time_for_purchase_before_exhaustion=3.0)
poopbag_model.exhaustion_time = 0.0
poopbag_model.current_amount = 0.0
poopbag_model.increase_amount_and_simulate(0.0, 120.0)
print poopbag_model.exhaustion_time


56.2508398654

Model Transaction Occurences

Transaction occurences will be modeled as:

$$p(t) = \max_{i \in items} p_i(t)$$

In [93]:
class TransactionSimulator(object):
    def __init__(self, item_models=None, timestep=None):
        self.item_models = item_models
        self.timestep = timestep
        
        self.time = 0.0 # days?
        
    def probability(self, time):
        prob = 0.0
        for model in self.item_models:
            prob = max(prob, model.probability_purchase(time))
        return prob
    
    def accept(self, time):
        prob = self.probability(time)
        r = random.random()
        return r < prob
    
    def simulate(self, max_steps=None):
        for i in xrange(max_steps):
            accepted = self.accept(self.time)
            self.time += self.timestep
            if accepted:
                return self.time
        return None

In [147]:
poopbag_model = ExhaustableItem(usage_rate=1.0, amount_used_average=2.5, amount_used_variance=0.5, average_time_for_purchase_before_exhaustion=3.0)
food_model = ExhaustableItem(usage_rate=1.0, amount_used_average=0.5, amount_used_variance=0.2, average_time_for_purchase_before_exhaustion=3.0)
poopbag_model.increase_amount_and_simulate(0.0, 120.0)
food_model.increase_amount_and_simulate(0.0, 30.0)
simulator = TransactionSimulator(item_models=[food_model, poopbag_model], timestep=1.0)

In [148]:
transactions = []
steps = []
for i in xrange(100):
    trans_time = simulator.simulate(max_steps=1)
    steps.append(i)
    if trans_time is not None:
        transactions.append(1.0)
    else:
        transactions.append(0.0)

In [149]:
plt.plot(steps, transactions)


Out[149]:
[<matplotlib.lines.Line2D at 0x107fc2250>]

In [101]:



---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
<ipython-input-101-fe9f3e1bca32> in <module>()
----> 1 plt.plot(np.transaction_times)

AttributeError: 'module' object has no attribute 'transaction_times'

In [ ]: