空港ユーザーのパッケージ選択に関する数学的モデリング

前段時間、ある空港のオーナーが私に連絡を取り、ユーザーがどのパッケージを選ぶかを計算するための数学的な助けを求めました。私は、数学モデルを構築するために必要なデータが必要だと返答しました。しかし、彼らはデータを提供することを拒否しました。おそらく、データを収集していないか、データの匿名化が面倒だったからでしょう。

そこで、協議の結果、私はデータを必要とせず、データに基づいて調整可能なモデルを構築し、コストを削減するためのいくつかの提案も行いました。その結果、彼らは私により高い報酬を支払うことに同意しました。

彼らの空港では PHP を使用していますが、計算の便宜上、私は pytorch を使用しました。彼らはおそらく、このモデルを実行するために Python ができる人を探す必要があるようです。

モデルパラメータ#

ユーザーの予算が従う確率密度関数を $b(x_1)$ 、流量需要が従う確率密度関数を $t(x_2)$ 、ノードの質の需要が従う確率密度関数を $q(x_3)$ と仮定します。

ユーザーが上記の 3 つの要因に対する重視の割合の確率密度関数を $a(y_1,y_2)$ と定義します：

$y_1\in(0,1)$ 、 $y_2\in(0,1)$ 、 $y_1+y_2<1$ 、そうでなければ $a(y_1,y_2)=0$
$y_1$ は予算の重み、 $y_2$ は流量の重み、 $1-y_1-y_2$ は質の重みです。

各パッケージには価格、流量、質の 3 つの属性があります。ユーザーがどのパッケージを購入するかを考えるとき、階層分析法に従って意思決定を行うと仮定します。上記の $b$ 、 $t$ 、 $q$ のパラメータはパラメータ推定が必要です。

$P=(p_{r1},p_{r2},p_{r3},\cdots,p_{rn})'$ 、 $T=(t_1,t_2,t_3,\cdots,t_n)'$ 、 $Q=(q_1,q_2,q_3,\cdots,q_n)'$ を製品配列ベクトルとし、それぞれ製品の価格、流量、質を表します。製品配列ベクトルは製品の性質によって決まり、パラメータ推定は必要ありません。

import torch

# 定義パッケージの価格、流量、質
packageNumber = 3
packagePrice = [1, 2, 3]
packageTrafficLimit = [1, 2, 3]
packageQuality = [1, 2, 3]

階層分析法による定量化#

予算の満足度は $f(p_r,x_1)=10^{-\frac{p_r}{x_1}}$ と定義できます。ここで、 $p_r$ はパッケージの価格、 $x_1$ は予算です。

流量の満足度は $g(t,x_2)=\sqrt{\displaystyle\frac{4t}{x_2}-\frac{3}{4}}$ （根の下が 0 未満の場合は 0 と定義）と定義します。ここで、 $t$ はパッケージの流量、 $x_2$ はパッケージの需要です。

ノードの質の満足度は $h(q,x_3)=\sqrt[3]{\displaystyle\frac{q}{x_3}}$ と定義します。ここで、 $q$ はノードの質の定量指標、 $x_3$ はノードの質の需要です。

ユーザーは満足度を選択します。

y_1\cdot10^{\displaystyle-\frac{p_r}{x_1}}+y_2\sqrt{\displaystyle\frac{4t}{x_2}-\frac{3}{4}} + y_3\sqrt[3]{\displaystyle\frac{q}{x_3}}

最大のパッケージ

また、すべてのパッケージの満足度が 0.1 未満の場合、ユーザーはそのパッケージを選択しません。

ベクトル形式で表すと、次のようになります：

\text{SATISFY} = [f(P,x_1),g(T,x_2),h(Q,x_3)]'Y

SELECT=\begin{cases} \max_{\text{SATISFY}} \text{PACKAGE},\ \text{SATISFY}>0.1\\ \varnothing,\ \text{SATISFY}\leq 0.1 \end{cases}

対応する Python コードは次のとおりです。

# 定義ユーザーのパラメータとパッケージから満足度へのマッピング
def user_satisfaction(vector_x, vector_y, pr, t, q):
    # 確保XとYはtorchテンソル
    vector_x = torch.tensor(vector_x, dtype=torch.float32)
    vector_y = torch.tensor(vector_y, dtype=torch.float32)
    
    term1 = vector_y[:, 0] * 10 ** (-pr / vector_x[:, 0])
    term2 = vector_y[:, 1] * torch.sqrt(4 * t / vector_x[:, 1] - 0.75)
    term3 = vector_y[:, 2] * (q / vector_x[:, 2]) ** (1 / 3)
    
    return term1 + term2 + term3

ユーザーのモデリング#

ユーザーの予算分布#

予算はおおよそパレート分布に従うと考えられ、初期パラメータを $x_1\sim Par(1,0.6)$ と仮定します。すなわち、

b(x_1) = \frac{0.6}{x_1^{1.6}} (x_1>1)

ここで、0.6 は正確な値を得るためにパラメータ推定が必要です。

パレート分布は通常、不平等な分布、例えば富や収入、そしてそれに対応する予算を説明するために使用されます。現実世界では、富の分布はしばしば非常に不平等であり、少数の人々が大部分の富を持っています。したがって、ユーザーの予算を反映するために、パレート分布を使用することは適切であると考えられます。

さらに、パレート分布は「80/20 ルール」と密接に関連しており、80% の効果（消費など）が 20% の原因（消費者など）から生じるというものです。多くの経済モデルにおいて、この現象は一般的であり、少数の消費者が大部分の消費支出を占めています。

パレート分布の 2 つのパラメータのうち、最初のパラメータはカットオフ値であり、 $x_1<1$ のとき、 $b(x_1)=0$ です。カットオフ値はパラメータ推定できませんが、1 と仮定することは適切です。

流量需要分布#

流量需要はおおよそ正規分布に従うと考えられ、パラメータは $x_2\sim N(130,30^2)$ とします。すなわち、

t(x_2)=\frac{1}{40\sqrt{2\pi}}e^{-\displaystyle\frac{1}{2}\left(\frac{x_2-120}{40}\right)^2}

この 2 つのパラメータはともにパラメータ推定が必要です。

質の需要分布#

質を $1\sim10$ と定量化し、1 が完全なリムジン空港、10 が全 IPLC 多入口スマート解析空港（性能はゲーム加速器に近い）とすると、質の需要もパレート分布に従うと考えられ、パラメータは $x_3\sim Par(1,2)$ とします。すなわち、

q(x_3)=\frac{2}{x_3^2}

このパレート分布も同様に、2 番目のパラメータのみがパラメータ推定が必要です。

重み分布#

$y_1$ 、 $y_2$ 、 $y_3$ の分布は、実際には等辺三角形の中に均等に配置されていることに注意してください。この等辺三角形は三角錐の底面です。この三角錐の 3 つの側面は、直角辺の長さが 1 の直角三角形です。

等辺三角形の 3 つの頂点は、それぞれ 1 つの重みが 1 で他の 2 つの重みが 0 の状況を表します。

$(y_1,y_2)$ を使用する場合、3 つの状況はそれぞれ $(0,0)$ 、 $(1,0)$ 、 $(0,1)$ に対応します。 $y_1,y_2,y_3$ をより均等な等辺三角形にするために、この領域を $(-\sqrt{3}+1,0)$ 、 $(\sqrt{3}-1,1)$ 、 $(\sqrt{3}-1,-1)$ に線形変換する必要があります。この線形変換は簡単に求めることができます。

まず、 $(-\sqrt{3}+1,0)$ と $(0,0)$ を揃えますと、他の 2 つの点は $(\sqrt{3},1)$ 、 $(\sqrt{3},-1)$ になります。したがって、変換行列は

\begin{bmatrix} \sqrt{3} & \sqrt{3} \\ 1 & -1 \end{bmatrix} \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix}^{-1}= \begin{bmatrix} \sqrt{3} & \sqrt{3} \\ 1 & -1 \end{bmatrix}

平行移動ベクトルは $(-1,0)$ であり、すなわち

\begin{bmatrix} y_1' \\ y_2' \end{bmatrix}= \begin{bmatrix} \sqrt{3} & \sqrt{3} \\ 1 & -1 \end{bmatrix} \begin{bmatrix} y_1 \\ y_2 \end{bmatrix}+ \begin{bmatrix} -1 \\ 0 \end{bmatrix}

ただし、3 つの要因の重視については、一般的に人々はノードの質を気にせず、価格と流量のみを考慮することができると考えられます。したがって、半分をカットした二元正規分布を選択することができます。カットオフ前のパラメータは

(y_1,y_2)\sim \mathcal{N}\left( \begin{bmatrix} 1/2 \\ 1/2 \end{bmatrix}, \begin{bmatrix} 1/9 && 0 \\ 0 && 1/9 \end{bmatrix} \right)

すなわち、確率密度関数は

a(y_1,y_2)= \frac{9}{\pi}e^{\displaystyle -\frac{9}{2}\left[(x_1-1/2)^2+(x_2-1/2)^2\right]}

この 2 つのパラメータもパラメータ推定は不要であるべきです。

from torch.distributions import Pareto, Normal
from torch.distributions.multivariate_normal import MultivariateNormal

# 定義ユーザーの予算分布
user_budget_distribution = Pareto(1, 0.6)
# 定義ユーザーの流量需要分布
user_traffic_demand_distribution = Normal(130, 30)
# 定義ユーザーの質の需要分布
user_quality_demand_distribution = Pareto(1, 2)
# 定義ユーザーの重視分布
user_weight_mean = torch.tensor([0.5, 0.5])
user_weight_covariance = torch.tensor([[1 / 9, 0], [0, 1 / 9]])
user_weight_distribution = MultivariateNormal(user_weight_mean, user_weight_covariance)

# サンプル生成
sample_size = 10000     # サンプル数
user_budget_sample = user_budget_distribution.sample((sample_size,))
user_traffic_demand_sample = user_traffic_demand_distribution.sample((sample_size,))
user_quality_demand_sample = user_quality_demand_distribution.sample((sample_size,))
user_weight_sample = user_weight_distribution.sample((sample_size,))

# ユーザーのパッケージ選択をシミュレーション
y_1 = user_weight_sample[:, 0]
y_2 = user_weight_sample[:, 1]
y_3 = 1 - y_1 - y_2
user_sample = torch.stack([user_budget_sample, user_traffic_demand_sample, user_quality_demand_sample], dim=1)
user_weight_sample = torch.stack([y_1, y_2, y_3], dim=1)
user_satisfaction_of_packages = []
for i in range(packageNumber):
    pr = packagePrice[i]
    t = packageTrafficLimit[i]
    q = packageQuality[i]
    satisfactions = user_satisfaction(user_sample, user_weight_sample, pr, t, q)
    user_satisfaction_of_packages.append(satisfactions)
    
user_satisfaction_of_packages = torch.stack(user_satisfaction_of_packages, dim=1)
max_values, max_indices = torch.max(user_satisfaction_of_packages, dim=1)
indices = torch.where(max_values < 0.1, torch.tensor(-1), max_indices)

パラメータ推定#

以上より、推定が必要なパラメータは次のとおりです：

ユーザー予算のパレート分布の $k$ パラメータ： $k_b$
流量需要分布の平均と分散： $\mu_t$ 、 $\sigma^2_t$
質の需要のパレート分布の $k$ パラメータ： $k_q$

これらのパラメータは神経ネットワークを通じて推定できるようですが、私は研究していませんし、初期パラメータは使用できないわけではありません。

コストの最小化#

パラメータ定義#

ユーザーのノードの好みを

R = \begin{bmatrix} r_1 \\ r_2 \\ \vdots \\ r_n \end{bmatrix}

と仮定し、 $\sum_{i=1}^n r_i = 1$ と規定します。

ユーザーが各ノードで実際に使用する流量は

U= \begin{bmatrix} u_1 \\ u_2 \\ \vdots \\ u_n \end{bmatrix} = uR

ここで、 $u$ はユーザーが実際に使用する総流量です。

ここでは、ユーザーが流量を超過しても好まないノードを使用することはないと仮定します。

各ノードが受け入れられるユーザー数の制限を考慮せず、ユーザーはノードの流量を使い切る前にリクエストによってノードが過負荷になることはないと仮定します。したがって、各ノードの限界コストは

C=\begin{bmatrix} c_1 \\ c_2 \\ \vdots \\ c_n \end{bmatrix}

であり、コストは $C'U=uC'R$ となります。

ここでは限界コストが定数であると仮定していますが、一般的に空港は規模の経済（限界コストが使用量の減少に伴って低下する）です。

価格を $p_r$ とすると、限界利益は $p_r-uC'R$ となります。

したがって、利益データが得られ、各パッケージの平均利益は

\bar{r} = p_r - \bar{u}C'R

オーバーブッキング問題#

ノードの流量の限界コストを決定するために、オーバーブッキング比を計算する必要があります。 $c_n=c_{tn}/R_{os}$ とし、ここで $c_{tn}$ は実際の流量の限界コスト、 $R_{os}$ はオーバーブッキング比です。

オーバーブッキングが多すぎて可用性が低下しないようにするために、 $Ra_{max}$ をオーバーブッキング可用性と定義し、意味は： $Ra_{max}$ の確率でオーバーブッキングによって不可用にならないことです。

ここで超パラメータ $Ra_{max}=1-10^{-4}=99.99\%$ と定義します。

パッケージ $g_1,g_2,\cdots,g_n$ に対して、ユーザーが実際に使用する総流量の確率密度関数を $f_{n}(u)$ とすると、特定のノードで使用される流量の確率密度関数は $r_nf_n(ur_n)$ となります。

ノードの実際の流量を $t_r$ とすると、可用性を満たすためには

Pr\{u<t_r\} \geq Ra_{max}

すなわち、

\sum_n\int_0^{t_r} r_nf_n(ur_n) dur_n=\sum_n\int_0^{t_r/r_n}f_n(u)du\geq Ra_{max}

例えば、1 つのパッケージのみがあると仮定し、そのパッケージの制限が 150G で、日本のノードの倍率が 1、ユーザーの好みの程度が 0.5、使用する流量が正規分布 $u\sim N(80,100)$ である場合、 $t_r$ は $58.595$ 、オーバーブッキング比は $2.5423$ であるべきです。

ユーザーが使用する流量は、以前のパラメータ推定から得られたユーザーの需要流量と大体同じですが、一定の誤差があります。これは、ユーザーが自分の必要な流量を過大評価することが多いためです。しかし、ユーザーが実際に使用する流量がユーザーの流量需要と等しいと仮定することは可能であり、流量に余裕を持たせることは問題ではありません。

なぜなら

R_{os}=\frac{t_r}{\sum_n r_nt}

であり、 $t_r$ は方程式を解くことで求められますので、 $R_{os}$ はデータ分析によって得られます。

最大期待利益問題#

最大期待利益問題は、すなわち

\sum_n f_n\bar{r_n} =f_n\left( p_{r,n} - \bar{u_n}C'R\right)

を最大化する問題です。（ここで $f_n$ は選択率です）

$C$ と $R$ は基本的に不変の量であり、ユーザーの需要分布のパラメータ推定後、 $f_n$ はパッケージの属性のみに依存し、 $\bar{u_n}$ も安定しているため、最大利益問題はパッケージ属性を決定する問題となります。

もしパッケージの質の属性が定数であると仮定すると、変数は各パッケージの価格と流量のみになります。

したがって、最大期待利益を計算する方法は次のとおりです：

$C'R$ を統計的に得る
各パッケージの質を決定する
各パッケージの価格または流量を決定し、もう一方を変数として設定する。
変数に対して勾配降下法を適用し、上記の和を最大化する
- $f_n$ はユーザーの数学モデルの後にモンテカルロ法で得られる
- $p_r$ は価格であり、入力変数です
- $u_n$ はユーザーのモデリングに関連しています