ポアソン分布~賢い在庫管理、宝くじまで

ポアソン分布

世の中にある事象が一定の頻度で、ある期間中にいつでも起きうるが、正確な発生時刻が知ることができない。例えば自動車が交差点を通過、赤ちゃんが病院で誕生、スーパーで商品が販売、機械が故障との事象があるとして、機械は月に1回故障、スーパーは1日に牛乳50本を販売、病院に1時間に3人の赤ちゃんが誕生、15分の間に40台の車が交差点を通過という具体的事象の発生頻度に対して、特定の期間例えば1日における独立事象の発生数の確率を求めるには以下ポアソン分布が使用可能となる。ポアソン分布とは、単位時間あたりに平均\(λ\)回起こる事象が期間\(t\)に\(k\)回起きる確率を表すのに使われる確率分布のこと。
$$P[X(t)=k]=\frac{(λt)^k e^{-λt}}{k!}$$
ただし、tは連続時間の長さ(期間)、λは単位時間に事象発生数(期待値)、eは自然対数の底(ネイピア数)。事象が期間t内に発生数がkの確率を計算する。ポアソン分布の期待値、分散ともλである。

応用例

在庫管理への応用例

【在庫管理の使用例】某ECサイトでは、1日に平均でPCを2台販売している。在庫切れの確率を5%以下に抑えたいので、どれだけの量の在庫が必要なのか。ただし、\(λ=2, t=1\)、単位はそれぞれ台、天とする。
1日にPCを0台販売の確率:
$$P[X(1)=0]=\frac{(2*1)^0 e^{-2*1}}{0!}=e^{-2}$$
1日にPCを1台販売の確率:
$$P[X(1)=1]=\frac{(2*1)^1 e^{-2*1}}{1!}=2e^{-2}$$
1日にPCを2台販売の確率:
$$P[X(1)=2]=\frac{(2*1)^2 e^{-2*1}}{2!}=2e^{-2}$$
1日にPCを3台販売の確率:
$$P[X(1)=3]=\frac{(2*1)^3 e^{-2*1}}{3!}=\frac{4}{3}e^{-2}$$
1日にPCを4台販売の確率:
$$P[X(1)=4]=\frac{(2*1)^4 e^{-2*1}}{4!}=\frac{2}{3}e^{-2}$$
1日にPCを5台販売の確率:
$$P[X(1)=5]=\frac{(2*1)^5 e^{-2*1}}{5!}=\frac{4}{15}e^{-2}$$
なので、1日にPCを3台以上販売の確率:
$$ P[X(1)>2]=1-\sum_{i=0}^{i=2}P[X(1)=i]\\\scriptsize=1-e^{-2}-2e^{-2}-2e^{-2}=0.3233$$
1日にPCを4台以上販売の確率:
$$P[X(1)>3]=1-\sum_{i=0}^{i=3}P[X(1)=i]\\\scriptsize=1-e^{-2}-2e^{-2}-2e^{-2}-\frac{4}{3}e^{-2}=0.1429$$
1日にPCを5台以上販売の確率:
$$P[X(1)>4]=1-\sum_{i=0}^{i=4}P[X(1)=i]\\\scriptsize=1-e^{-2}-2e^{-2}-2e^{-2}-\frac{4}{3}e^{-2}-\frac{2}{3}e^{-2}=0.0527$$
1日にPCを6台以上販売の確率:
$$P[X(1)>5]=1-\sum_{i=0}^{i=5}P[X(1)=i]\\\scriptsize=1-e^{-2}-2e^{-2}-2e^{-2}-\frac{4}{3}e^{-2}-\frac{2}{3}e^{-2}-\frac{4}{15}e^{-2}=0.0166$$
つまり1日にPCを6台以上販売の確率が1.66%なので、PCを6台以上(最低6台)の在庫が確保できれば在庫切れの確率を5%以下に抑えられる。

ナンバーズ予測への応用例

ナンバーズ数字0~9(ロト6数字1-43、ロト7数字1-37)の出た回数を1000回分もしくは2000回分で平均して、\(λ_i\)にしておく。最近の9回、ナンバーズ数字0~9(ロト6数字1-43、ロト7数字1-37)の出現回数+今度も出る(+1)の確率をポアソン分布で計算して大きさが大から小の順で並べて、確率が大きい数字が予測の結果とする。注意!:あくまで確率論からの予測(可能性)なので、ポアソン分布から宝くじが当たる意味ではない。

発明者が残した言葉

「私がランダム現象を記述する1つの確率分布を確立した」 by フランス数学者・物理学者ポアソン。

1+

ロボット・ドローン部品お探しなら
ROBOT翔・電子部品ストア

粒子フィルタ

はじめに

ガウス分布に従う時系列確率変数の非線形関数は必ずしもガウス分布にならない、また非ガウス分布の確率変数が存在するから、優雅なるカルマンフィルタファミリーには限界あり、無理な場合がある。非線形、非ガウス時系列空間モデルに対して、1990年代から未知の状態密度関数p(x)に対して、既知の密度関数q(x)、p(x)⊂q(x) から逐次モンテカルロ法(逐次重要性サンプリングSequential Importance Samplingに、リサンプリングResampling)による粒子フィルタ(パーティクル・フィルタ、PF)が提案されて以来、研究・応用が活発になっている。粒子フィルタがROSのSLAM Gmappingアルゴリズムに組まれている。

粒子フィルタ

状態空間モデル

・システムモデル(状態方程式)
$$x_t = f_t(x_{t-1},v_t)$$

・観測モデル(観測方程式)
$$y_t = h_t(x_t,w_t)$$

ただし、\( v_t、w_t\)はフィルタリングをかける事前に既知で、非ガウス分布可能、\( y_t\)は2項分布またポアソン分布に従うcount dataのモデル化なども可能になる。また、粒子フィルタには、\( v_t、w_t\)はガウス分布でも構わない。
以下に粒子\( x_t^{[i]} \)~\( f_t (x_t^{[i]} | x_{t−1}^{[i]}, v_t) \)として、粒子フィルタリングの手順を示す。

粒子フィルタリングの手順

1、初期化\( (t=0) \)
初期分布\( f_t (x_0^{[i]}, v_0) \)に従って、n個の粒子\( \{x^{[i]}_0 | i = 1,2, ⋯ ,n\} \)を無作為に発生させる。

2、一期先予測\( (t=t+1) \)
粒子\( x_t^{[i]} \)を\( f_t (x_t^{[i]} | x_{t−1}^{[i]}, v_t) \)に従って状態推移させ、それぞれの類似値が\(x_t^{[i]}\)である粒子の集合\( \{\hat{x}^{[i]}_t | i = 1,2, ⋯ ,n\} \)を発生させる。

3、フィルタリング
・重み計算
粒子\( \hat{x}_t^{[i]} \) の重み\( w_t^{[i]} = p(y_k | \hat{x}_k^{[i]}) \)を計算する。
ただし、\( p(y_k | \hat{x}_k^{[i]}) \)は、尤度であり、観測値または類似値ではないことに注意する。

・重みの正規化
$$\hat{w}_t^{[i]} = \frac{w_t^{[i]}}{ \sum_{i=1}^n w_t^{[i]} }$$

・リサンプリング
粒子\( \hat{x}^{[i]}_0 \)を\(\hat{w}_t^{[i]}\)に従った確率でリサンプリングし、粒子集合\( \{x^{[i]}_t | i = 1,2, ⋯ ,n\} \)を発生させる。

・2に戻る

参考文献

粒子フィルタ、著者:樋口知之 先生

1+

ロボット・ドローン部品お探しなら
ROBOT翔・電子部品ストア

SLAM~Unscentedカルマンフィルタ

はじめに

UKF(Unscented kalman filter)では無損変換(Unscented Transform)の線形化手法を利用し、つまり(2*n次元+1)個Sigma Pointの線形回帰に基づき、確率変数の非線形関数を線形化する。確率変数の拡張を考慮するので、この線形化は拡張カルマンフィルタEKFで使用するテイラー級数線形化より比較的正確に状態推定できる他、EKFと違い、状態・観測ヤコビアン行列を求める必要でなく、但しEKFと同様にUKFも予測と更新の手順から時系列の状態推定を行う。状態変数、観測変数ともガウス分布に当てるのがカルマンフィルタKFの共通特徴という。またガウス分布に当てないパーティクル(和文粒子、英文Particle)フィルタPFがある。そもそもガウス分布に従う確率変数の非線形関数はガウス分布にならないから、EKF、UKFよりPFの方が精度よく、ROSのSLAM Gmappingアルゴリズムに組まれている。

Unscented Kalman Filter

状態空間モデル

状態モデル(状態方程式)
$$x_t = F \cdot x_{t-1} + B \cdot u_{t-1} + q_{t-1} \\\scriptsize=\begin{bmatrix}1&0&0&0\\0&1&0&0\\0&0&1&0\\0&0&0&0\end{bmatrix}\begin{bmatrix}x_{t-1}\\y_{t-1}\\ψ_{t-1}\\v_{t-1}\end{bmatrix}+\begin{bmatrix}Δt\cdot cosψ_{t-1}&0\\Δt\cdot sinψ_{t-1}&0\\0&Δt\\1&0\end{bmatrix}\begin{bmatrix}v_{t-1}\\ψ_{t-1}\end{bmatrix}+q_{t-1}$$

観測モデル(観測方程式)
$$z_t = H \cdot x_t + r_t \\=\begin{bmatrix}1&0&0&0\\0&1&0&0\end{bmatrix}\begin{bmatrix}x_t\\y_t\\ψ_t\\v_t\end{bmatrix}+r_t$$

予測ステップ

$$χ_{t-1}=(x_{t-1} , x_{t-1}+γ\sqrt{{\sum}_{t-1}} , x_{t-1}-γ\sqrt{{\sum}_{t-1}})\\
\bar{χ}_t^{*}=g(u_t,χ_{t-1})\\
\bar{x}_t=\sum_{i=0}^{2n} w_m^{[i]} \bar{χ}_t^{*[i]}\\
\bar{\sum}_t=\sum_{i=0}^{2n} w_c^{[i]}(\bar{χ}_t^{*[i]}-\bar{x}_t)(\bar{χ}_t^{*[i]}-\bar{x}_t)^T +R_t\\
\bar{χ}_{t}=(x_{t} , x_{t}+γ\sqrt{{\sum}_{t}} , x_{t}-γ\sqrt{{\sum}_{t}})\\
\bar{Z}_{t}=h(\bar{χ}_{t})\\
\hat{z}_{t}=\sum_{i=0}^{2n} w_m^{[i]} \bar{Z}_{t}^{[i]}\\
S_t=\sum_{i=0}^{2n} w_c^{[i]}(\bar{Z}_t^{[i]}-\hat{z}_t)(\bar{Z}_t^{[i]}-\hat{z}_t)^T +Q_t\\
\bar{\sum}^{x,z}_t=\sum_{i=0}^{2n} w_c^{[i]}(\bar{χ}_t^{[i]}-\bar{x}_t)(\bar{Z}_t^{[i]}-\hat{z}_t)^T
$$

更新ステップ

$$K_t=\bar{\sum}^{x,z}_t S^{-1}_t\\
x_t = \bar{x}_t+K_t(z_t-\hat{z}_t)\\
{\sum}_t = \bar{\sum}_t-K_t S_t K^T_t$$

参考文献

Probabilistic Robotics/Sebastian Thrun, Wolfram Burgard and Dieter Fox.

1+

ロボット・ドローン部品お探しなら
ROBOT翔・電子部品ストア

SLAM~拡張カルマンフィルタ

概要

SLAM(ロボットの自己位置推定、マッピング同時に行うこと)に使用可能なセンサは、色々とあるが、それぞれ単一のセンサでは完璧ではないため、複数のセンサ情報を、統計的・確率的に組み合わせて、より精密で安定した自己位置を得る手法が一般的だ。状態値(位置、姿勢またはポーズ)、観測値とも正規分布(ガウス分布)に近似可能の場合、カルマンフィルタKFが適用される。また、非線形問題に対して、拡張カルマンフィルタEKFが適用される。勿論、UKF他のフィルタが存在する。なぜかROSに公式に採用されたのはEKFのようだ。ここでは、タイヤパルスセンサ(ホール式や、光学式エンコーダ)からのodom情報と距離センサ(レーザや、カメラなど)からの位置情報を組み合わせて、マッピング地図がある前提でロボットの自己位置、姿勢に拡張カルマンフィルタをかけて推定してみる。

状態空間モデル

$$x_t = f_{t-1}(x_{t-1}^{true}+p_{t-1}) + q_{t-1}\\=\scriptsize\begin{bmatrix}x_{t-1}+v_{t-1}\cdot Δt\cdot cosψ_{t-1} & 0 & 0 \\0 & y_{t-1}+v_{t-1}\cdot Δt\cdot sinψ_{t-1} & 0\\0 & 0 & ψ_{t-1}+ω_{t-1}\cdot Δt\end{bmatrix}+q_{t-1}$$
$$y_t = H_t(x_t^{true}) + r_t\\=\begin{bmatrix}x_t\\y_t\end{bmatrix}+ r_t$$

拡張カルマンフィルタの適用

予測ステップ

事前状態推定値

$$\bar{x}_t=\scriptsize\begin{bmatrix}\hat{x}_{t-1}+v_{t-1}\cdot Δt\cdot cos\hat{ψ}_{t-1} & 0 & 0 \\0 & \hat{y}_{t-1}+v_{t-1}\cdot Δt\cdot sin\hat{ψ}_{t-1} & 0\\0 & 0 & \hat{ψ}_{t-1}+ω_{t-1}\cdot Δt\end{bmatrix}$$

線形化近似(ヤコビアン行列)

$$\hat{F}_t =\frac{\partial f_t(x)}{\partial x}|_{x=\hat{x}_{t-1}}\\=\begin{bmatrix}1 & 0 & -v_{t-1}\cdot Δt\cdot sin\hat{ψ}_{t-1} \\0 & 1 & v_{t-1}\cdot Δt\cdot cos\hat{ψ}_{t-1}\\0 & 0 & 1\end{bmatrix}$$

事前誤差共分散行列

$$\bar{P_t} = \hat{F}_{t-1}\hat{P}_{t-1}\hat{F}_{t-1}^T + Q_{t-1}$$

フィルタリングステップ

カルマンゲイン行列

$$K_t = \bar{P_t}H_t^T(H_t\bar{P_t}H_t^T + R_t)^{-1}$$

状態推定値

$$\hat{x_t} = \bar{x_t} + K_t(y_t – H_t\bar{x_t})$$

事後誤差共分散行列

$$\hat{P_t} = (I-K_tH_t)\bar{P_t}$$

関連記事

9軸IMUセンサ 6軸/9軸フュージョン 低遅延 USB出力 補正済み ROS対応
研究開発用 台車型ロボット キット

2+

ロボット・ドローン部品お探しなら
ROBOT翔・電子部品ストア

研究開発用 台車型ロボット キット

低コスト台車ロボットプラットフォーム

外形寸法は400x395x212.5mm、可搬重量は約40kg、最高速度は約0.5m/s、アルミ合金により軽量で高剛性なロボットキットとなります。SLAM、ROS、AI等の機能が実装する台車型ロボット、自律移動ロボット、自律走行ロボット、無人搬送ロボット、AGVの研究開発用プラットフォームとしてご利用いただけます。

主な規格

外形寸法 L400mm×W395mm×H212.5mm
台車重量 約10kg
積載重量 40kg(台車重量を除く)
駆動方式 駆動輪2個(二輪駆動)、前部キャスター1個
最高速度 約1.0m/s
駆動輪直径 125mm
駆動輪素材 アルミホイル、ゴム
ブラシ付きモータ 定格電圧12V、定格電流3A、定格回転数77RPM、定格トルク2.79Nm、減速比51、ストール電流7A
光学式エンコーダ  500PPR(モータ・シャフト)/13500PPR(ギヤボックス・シャフト)
モータドライバ エンコーダAB相入力2系統、モータ制御入力2系統、最大电流30A

販売元

研究開発用 台車ロボット キット
カスタマイズ、特注もお承りいたします。

関連記事

9軸IMUセンサ 6軸/9軸フュージョン 低遅延 USB出力 補正済み ROS対応
点検向け自律移動ロボットRED(薄型・小型)
SLAM~拡張カルマンフィルタ
SLAM~Unscentedカルマンフィルタ
粒子フィルタ
オイラー角~ジンバルロック~クォータニオン

2+

ロボット・ドローン部品お探しなら
ROBOT翔・電子部品ストア