## Preparation

In textbooks, calculations are mainly performed by Excel functions. Although Excel has an excellent GUI, it does not have enough API libraries to connect to external web systems or data analysis tools. Therefore, we will use Python to perform the same calculations as in the textbook. Here are the preparations for this.

### github

• The jupyter notebook file on github is here .

• If you want to run it on google colaboratory here

### Author’s environment

This is the author’s environment.

The author's environment.

ProductName: Mac OS X
ProductVersion: 10.14.6
BuildVersion: 18G2022

Python -V

Python 3.7.3


import numpy as np
import scipy
from scipy.stats import binom

%matplotlib inline
%config InlineBackend.figure_format = 'svg'

import matplotlib
import matplotlib.pyplot as plt
import seaborn as sns

print("numpy version :", np.__version__)
print("matplotlib version :", matplotlib.__version__)
print("sns version :",sns.__version__)

numpy version : 1.16.2
matplotlib version : 3.0.3
sns version : 0.9.0


## Overview

In the end-of-book commentary2, the authors explain six tools that they use frequently, with specific examples. You can refer to this and apply it when you need it in a real business situation.

1. Gamma-Poisson Least Sensitive Model
• Based on the data of “When did you buy recently” and “When did you visit recently” (Recent Purchase Period: Recency), it tells us which brand, which facility, and which period we should focus our resources relatively. 2.
1. negative binomial distribution
• This is a tool that corrects for the difference between the consumer household panel’s data on your brand and the actual sales figures. This means that it is very useful for benchmarking trial rates, repeat rates, and purchase frequency during forecasting.
• It will tell you how much market share you can get in a newly created category. It also allows you to simulate the market share based on your marketing plan.
1. trial model/repeat model
• Using data from concept tests, concept use tests, and household panel data, you can predict the sales of a new product in the first year of its launch.

VPP Model (Volume per Purchase)

• Helps you determine the size of your product. 6.
1. delishley NBD model
• It provides a concrete example of how Delishley S is calculated for NBD category K. It predicts the quarterly purchase rate, the number of quarterly purchases, and the percentage of 100% loyal customers for Colgate in Table 1-4 of the textbook, as explained in Explanation 1, 1-6.

## 2-1. Gamma-Poisson Least Sensitive Model

We can build an NBD model by calculating $m$ and $k$ from the data of “when did you buy recently” and “when did you visit recently”. the formula describing the NBD model can be calculated as follows.

$$P\left(r \right) = \frac{\left(1 + \frac{M}{K} \right)^{-K} \cdot \Gamma\left(K + r \right)}{\Gamma\left(r + 1 \right)\cdot \Gamma\left(K \right)} \cdot \left(\frac{M}{M+K} \right)^r$$

Let $mt$ denote the corresponding average value of $M$ over the period $t$ for a given product, and let $K$ denote $k$. Since the penetration rate can be calculated by subtracting the number of people who never buy this product from 100%, we have

\begin{aligned} P(t) &=1-P_0\left(r=0 \right) \\ &= 1 - \frac{\left(1 + \frac{mt}{k} \right)^{-k} \cdot \Gamma\left(k + 0 \right)}{\Gamma\left(0 + 1 \right)\cdot \Gamma\left(k \right)} \cdot \left(\frac{mt}{mt+k} \right)^0 \\ &= 1 - \left(1 + \frac{mt}{k} \right)^{-k} \end{aligned}

It becomes Thus, the penetration rate in a period $t$ and $t-1$ is

$$P\left(t \right) - P\left(t-1 \right) = \left(1+\frac{m\times t}{k}\right)^{-k} - \left(1+\frac{m\times \left(t-1 \right)}{k}\right)^{-k}$$

This becomes

To apply this to an arbitrary period, we use two variables $t_1$ and $t_2$ and define $f \left(x \right)$ as follows.

$$f\left(t_1,t_2,m,k \right) = \left(1+\frac{m\times t_1}{k} \right)^{-k} - \left(1+\frac{m\times t_2}{k} \right)^{-k}$$

Table 10-1 in the textbook can be expressed using the common function $f\left(x \right)$ as follows

Gamma distributionActual values
$\displaystyle f\left(t_1=0,t_2= \frac{14}{31}\right)$43.9%
$\displaystyle f\left(t_1=\frac{14}{31},t_2=1 \right)$25.6%
$\displaystyle f\left(t_1=1,t_2=2 \right)$19.1%
$\displaystyle f\left(t_1=2,t_2=3 \right)$5.1%
$\displaystyle f\left(t_1=3,t_2=4 \right)$1.5%
$\displaystyle f\left(t_1=4,t_2=5 \right)$0.7%
$\displaystyle f\left(t_1=5,t_2=6 \right)$1.4%
$\displaystyle f\left(t_1=6,t_2=\infty \right)$2.7%

### Derivation of m,k by scipy’s curve_fig

In general, to perform least-squares fitting of nonlinear functions, we can use the curve_fit module of scipy . According to the scipy website

  scipy.optimize.curve_fit(f, xdata, ydata, p0=None, sigma=None, absolute_sigma=False, check_finite=True, bounds=(-inf, inf), method=None, jac =None, **kwargs)[source].


In addition, $xdata$ and $ydata$ are

  xdata : array_like
The independent variable where the data is measured. Must be an M-length sequence or an (k,M)-shaped array for functions with k predictors.

ydata : array_like
The dependent data, a length M array - nominally f(xdata, ...)


It is defined as To solve the fitting problem in the textbook.

$$f\left(t_1,t_2,m,k \right) =\left(1+frac{m\times t_1}{k} \right)^{-k} -\left(1+frac{m\times t_2}{k} \right)^{-k}$$

For the function $\displaystyle f\left(x \right)$ defined above, we will be solving a fitting problem for a two-variable function, where the two variables specify the time period (to find the number of purchases in two weeks to a month, we use $\displaystyle t_1=\frac{ 14}{31}, t_2=1$) and are defined as a two-dimensional array as follows

x = np.array([.
[0.0 ,14/31 ,1.0 ,2.0 ,3.0 ,4.0 ,5.0 ,6.0 ],
[14/31 ,1.0 ,2.0 ,3.0 ,4.0 ,5.0 ,6.0 , 10000.0].
])


x[0] is the array of $t_1$, and x[1] is the array of $t_2$. x[1,7]=10000.0 is originally $\infty$, but infinity is not acceptable in actual numerical calculations, so 10000 is used, which is effectively infinity. This value can be as small as 100.

The actual code for fitting is as follows

import json
import numpy as np

from scipy.optimize import curve_fit
from scipy.special import gamma

def _get_delta_nbd(x, m, k):
return (1 + m * x[0] / k )**(-k) - (1 + m * x[1] / k )**(-k)

x = np.array([.
[0.0 ,14/31 ,1.0 ,2.0 ,3.0 ,4.0 ,5.0 ,6.0 ],
[14/31 ,1.0 ,2.0 ,3.0 ,4.0 ,5.0 ,6.0 , 10000.0].
])

y = [0.439, 0.256, 0.191, 0.051, 0.015, 0.007, 0.014, 0.027])

parameters, covariances = curve_fit(_get_delta_nbd, x, y)
print('parameters : ', parameters)
print('covariances : ', covariances)

parameters : [1.37824241 4.14429889]
covariances : [[ 0.00284656 -0.03699629]]
[-0.03699629 1.57449471]]


The resulting $m$ and $k$ are

\begin{aligned} m&= 1.378 \\ k&= 4.144 \end{aligned}

and $m$ and $k$ used by the author.

\begin{aligned} m&= 1.37552 \\ k&= 4.061 \end{aligned}

The value is almost equal to.

## 2-2. Negative binomial distribution

This section explains how to correct the panel data by using the difference between the actual sales data and the data obtained from the panel data.

First of all, the panel data gives us the following information

• (A) : Number of households
• (B) : Penetration rate
• (C) : Average number of purchases
• (D) : Average number of items purchased
• (E) : Average purchase price

Please refer to Table 10-2 below for specific values. It is P281 in the textbook. From here, we can create a

Sales by Panel Data Sales = number of households x penetration rate x average number of purchases x average number of units purchased x average unit purchase price
The panel data sales can be obtained as follows.

In addition, we know the following results.

• Sales amount

Using this ratio of panel data sales to actual sales, we can correct the panel data and various parameters.

To do this, the textbook makes a number of important assumptions. The following are some of the assumptions made in the textbook, which will help you to understand the subsequent calculations smoothly.

### Assumptions

• Actual current sales: 5.89 billion yen
• Sales based on panel data: 4.12 billion yen (actual sales ratio: 70%, calculated by AxBxCxDxE)
• Average number of items purchased per purchase is the same in reality as in the panel data
• Average unit price per purchase is the same in reality as in the panel data
• $K$ is the same in reality as in the panel data

Again, only “sales” are known as actual results. In the textbook example, we only know that the sales are 5.89 billion yen.

Table 10-2 Corrected panel data for a brandAdjusted household panel data for one year
ItemBefore correctionAfter correction
(A)Total number of households in 2008 (thousands)4997349973
(B)Penetration rate15.0%17.4%
(C)Average number of purchases2.503.07
(D)Average number of items purchased per transaction1.101.10
(E)Average unit price per purchase200 yen200 yen
(F)Percentage of customers who purchased two or more times50%55%
(G)Annual Sales (AxBxCxDxE)4.12 billion yen5.89 billion yen
(H)Ratio of G to actual70%100%
Table 10-3 Calculating the correctionOne year of household panel data
ItemBefore correctionAfter correction
(I)Brand $m$:(BxCxD)0.41250.5893
(J)Brand $k$0.098990.09899
(K)$P_0$(probability of never buying)85.00%82.53%
(L)$P_1$(Probability of buying once)6.79%7.00%
(M)$P_{+2}=100\%-P_0-P_1$8.21%10.47%
(N)Percentage of buyers who purchased two or more times by model:$\left(\frac{M}{B}\right)$54.76%59.95%

### Steps in the correction

1. brand’s $m=$penetration rate x average number of purchases x average number of units purchased
2. $k$ of brand

$$P\left(r \right) = \frac{\left(1 + \frac{M}{K} \right)^{-K} \cdot \Gamma\left(K + r \right)}{\Gamma\left(r + 1 \right)\cdot \Gamma\left(K \right)} \cdot \left(\frac{M}{M+K} \right)^r$$

by substituting $\displaystyle K=k, M=m=0.4125, r=0$ into

$$P_0=\frac{\left(1+\frac{m}{k} \right)^{-k}\cdot \Gamma\left(k+0 \right)}{\Gamma\left(0+1 \right)\cdot \Gamma\left(k \right)}=\left(1+\frac{0.4125}{k} \right)^{-k} =0.85$$

We obtain the nonlinear equation Where $\displaystyle P_0$ is the probability of never having made a purchase, so it can be calculated from (1 - penetration) and

$$P_0=1 - 015 = 0.85$$

We use the fact that $$P_0=1 - 015 = 0.85$$

#### Solving nonlinear equations numerically

Equation to find $k$

$$\left(1+\frac{0.4125}{k} \right)^{-k} =0.85$$

is nonlinear and cannot be solved analytically. We will use numerical methods to solve it using a computer. Here, we will use Python’s newton method to obtain the solution. In the textbook, we use Excel to obtain the value of $k$, but either method is fine.

The result is $$k=0.09899$$ as a result.

#### python code

The python code to get $k$ is as follows

from scipy.optimize import newton

MIN_k = 0
MAX_k = 1.0

def check_k(k):
if MIN_k < k and k < MAX_k:
return True
else:
return False

def get_k(m, P0):

def func(k, m=m, P0=P0):
return (1 + m / k) ** (-1 * k) - P0

k = None
try:
for initial_k in [(i + 1) * 0.01 for i in range(100)]:
k = newton(func, initial_k)
if check_k(k):
return k
else:
if not check_k(k):
return None
except:
return None

m = 0.4125
P0 = 0.85

print("k = {:,.5f}".format(get_k(m, P0)))

k = 0.09893


and the value is almost equal to the textbook even using python.

#### 3. P_1, the probability of buying once

$P_1$ is the same as $P_0$.

$$P\left(r \right) = \frac{\left(1 + \frac{m}{k} \right)^{-k} \cdot \Gamma\left(k + r \right)}{\Gamma\left(r + 1 \right)\cdot \Gamma\left(k \right)} \cdot \left(\frac{m}{m+k} \right)^r$$

Just substitute $\displaystyle k=0.09899, m=0.4125, r=1$ into However, since it contains a gamma function, calculations by python or excel are required.

\begin{aligned} P_1&=P\left(1 \right) \\ &= \frac{\left(1 + \frac{0.4125}{0.09899} \right)^{-0.09899} \cdot \Gamma\left(0.09899 + 1 \right)}{\Gamma\left(1 + 1 \right)\cdot \Gamma\left(0.09899 \right)} \\ & \quad \quad \quad \times \left(\frac{0.4125}{0.4125+0.09899} \right)^1 \\ &= 0.0679 \end{aligned}

The following example shows how to do this

### Percentage of 100% loyal customers of Colgate

The numbers on the diagonal in Table 10-10 are the percentage of buyers who decide to buy only Colgate toothpaste, so we can divide them by the percentage of toothpaste purchases to get the percentage of Colgate loyal customers.

### Average number of purchases of Colgate

By multiplying the probability of purchase of Colgate by the number of purchases and summing them, we can calculate the expected value of the number of purchases (average number of purchases).

### python code

Here is the python code we used to calculate in Table 10-8, 9, and 10. The derivation of the negative binomial and Delishley distributions is difficult, but the calculation of the results themselves is not very complicated.

import numpy as np
import math
from scipy.special import gamma

def get_nbd(M, T, K, R):
return ((1 + M * T / K)**(-1 * K)) * \
(gamma(K + R) / math.factorial(R) / gamma(K)) * \
((M * T / (M * T + K)) ** R)

def get_p_rj_0(r, a, S, R):
return (math.factorial(R)/ math.factorial(r) / math.factorial(R - r)) * \
(gamma(S) / gamma(a) / gamma(S - a)) * \
(gamma(a + r) * gamma(S - a + R - r) / gamma(S + R)))

def print01():
for R in range(0,11):
print('R={} | '.format(R), end='')
for r in range(R + 1):
print('{:.3f} | '.format(round(get_p_rj_0(r=r, a=1.2 * 0.25, S=1.2, R=R), 3)), end='')
print()

def print02():
for R in range(0,11):
print('R={} | '.format(R), end='')
for r in range(R + 1):
print('{:.1f} % | '.format(round(100 * get_nbd(M=1.46, T=1, K=0.78,R=R) * get_p_rj_0(r=r, a=1.2 * 0.25, S=1.2, R=R), 3)), end='')
print()

print('Table 10-9 Percentage of category purchases by number of purchases when S=0.12')
print01()
print()
print()
print('Table 10-10 Percentage of category and brand purchases by number of purchases when S=0.12')
print02()

Table 10-9 Proportion of category by number of purchases when S=0.12
R=0 | 1.000 |
R=1 | 0.750 | 0.250 |
R=2 | 0.648 | 0.205 | 0.148 |
R=3 | 0.587 | 0.182 | 0.125 | 0.106 | R=4
R=4 | 0.545 | 0.168 | 0.113 | 0.091 | 0.083 |
R=5 | 0.514 | 0.157 | 0.105 | 0.083 | 0.072 | 0.069 |
R=6 | 0.489 | 0.149 | 0.099 | 0.078 | 0.066 | 0.060 | 0.059 |
R=7 | 0.468 | 0.143 | 0.094 | 0.074 | 0.062 | 0.055 | 0.052 | 0.052 |
R=8 | 0.451 | 0.137 | 0.090 | 0.070 | 0.059 | 0.052 | 0.048 | 0.045 | 0.046 |
R=9 | 0.437 | 0.132 | 0.087 | 0.068 | 0.057 | 0.050 | 0.045 | 0.042 | 0.040 | 0.041 |
R=10 | 0.424 | 0.128 | 0.084 | 0.066 | 0.055 | 0.048 | 0.043 | 0.040 | 0.038 | 0.037 | 0.038 |

Table 10-10 Ratios of category and brand by number of purchases when S=0.12
R=0 | 43.9 % |
R=1 | 16.7 % | 5.6 % |
R=2 | 8.4 % | 2.6 % | 1.9 % |
R=3 | 4.6 % | 1.4 % | 1.0 % | 0.8 % |
R=4 | 2.6 % | 0.8 % | 0.5 % | 0.4 % | 0.4 % |
R=5 | 1.5 % | 0.5 % | 0.3 % | 0.2 % | 0.2 % | 0.2 % | R=6
R=6 | 0.9 % | 0.3 % | 0.2 % | 0.1 % | 0.1 % | 0.1 % | 0.1
R=7 | 0.6 % | 0.2 % | 0.1 % | 0.1 % | 0.1 % | 0.1 % | 0.1 % | 0.1
R=8 | 0.3 % | 0.1 % | 0.1 % | 0.1 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0
R=9 | 0.2 % | 0.1 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 %
R=10 | 0.1 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 % | 0.0 %


## Summary

The above is my attempt to break down the explanations at the end of “Strategic Thinking in Probability” in my own way. I am not a marketing expert, nor do I have any practical experience. I usually work in IT-related fields, such as web system development and machine learning model development. In that context, I usually use Poisson distribution and Gamma distribution, but I never thought that they are applied to marketing in this way.

At first, my friend who specializes in marketing told me “What is the negative binomial distribution? According to my friend, overseas companies such as P&G usually use mathematics in their marketing, but in Japan, it seems that there is still a long way to go. Dr. Arenberg, a great authority on marketing, published a paper that is the basis of this book many decades ago. However, I believe that probability and statistics will be applied more and more to marketing in Japan in the future as well, as a result of “Strategy Theory of Probability Thinking”.