Unscented Transformation（U変換）とは
Unscented Transformation（U変換）は、確率変数 $x$ が非線形関数 $f$ によって変換された確率変数 $y = f(x)$ の平均 $\bar{y}$ と共分散 $P_y$ を推定する手法です。
$x$ の平均 $\bar{x}$ と共分散 $P_x$ が既知である場合、単純にモンテカルロ法で多数のサンプルを生成して変換し、その平均と共分散を計算することも考えられます。
$$ \bar{y} \approx \frac{1}{N} \sum_{i=1}^N f(x_i) $$ $$ P_y \approx \frac{1}{N} \sum_{i=1}^N (f(x_i) - \bar{y})(f(x_i) - \bar{y})^T $$
しかし、この方法では高精度な推定を得るために膨大な数のサンプル $N$ が必要となり、計算コストが非常に高くなります。
U変換は、この問題を解決するために、線形近似を行うことなく、少数の代表点（シグマ点）を用いることで、変換後の確率変数の統計的性質を効率的かつ高精度に推定します。これにより、モンテカルロ法の利点（非線形性への対応）を保ちつつ、計算量を大幅に削減できます。
参考: UKF （Unscented Kalman Filter）って何？
Pythonによる実装
ここでは、2次元の入力 $X=(X_1, X_2)$ を1次元の出力 $Y=f(X)=X_1 \cdot X_2$ に変換する例を考えます。
1. 入力確率変数の設定
入力 $X$ の平均ベクトル $\bar{x}$ と共分散行列 $P_x$ を定義します。
$$ \bar{x} = [E[X_1], E[X_2]] $$ $$ P_x = \begin{pmatrix} \text{Var}[X_1] & \text{Cov}[X_1, X_2] \ \text{Cov}[X_2, X_1] & \text{Var}[X_2] \end{pmatrix} $$
例として、$X_1$ が平均0、分散1、$X_2$ が平均1、分散4、共分散が2であると仮定します。
import matplotlib.pyplot as plt
import numpy as np
import random
import math
import scipy.linalg
# 入力xの次元数
n = 2
# 出力yの次元数
m = 1
# xの平均ベクトルと共分散行列
x_mean = np.array([0.0, 1.0])
x_P = np.array([[1.0, 2.0], [2.0, 4.0]])
print(f"入力xの平均: {x_mean}")
print(f"入力xの共分散行列:\n{x_P}")
# 非線形変換関数 f(x) = x[0] * x[1]
def f(x):
return np.array([x[0] * x[1]])
2. シグマ点の計算
シグマ点は、入力確率変数 $x$ の平均と共分散を正確に表現するように選ばれた代表点です。これらの点は、以下の式で計算されます。
$$ \sigma_0 = \bar{x} \tag{1} $$ $$ \sigma_i = \bar{x} + (\sqrt{(n+\lambda)P_x})_i \quad \text{for } i=1, \dots, n \tag{2} $$ $$ \sigma_i = \bar{x} - (\sqrt{(n+\lambda)P_x})_i \quad \text{for } i=n+1, \dots, 2n \tag{3} $$
ここで、$ (\sqrt{(n+\lambda)P_x})_i $ は、行列の平方根（コレスキー分解など）の $i$ 列目を表します。
パラメータ $\lambda$ は、以下の式で計算されます。 $$ \lambda = \alpha^2 (n + \kappa) - n \tag{4} $$
- $\alpha$: シグマ点の平均からの広がりを決定するスカラー値（通常0から1の間の小さい正の値）。
- $\kappa$: 通常は0に設定される調整パラメータ。シグマ点の広がりをさらに調整します。
# パラメータ設定
alpha = 0.5
kappa = 0.0
# 式4: lambdaの計算
lam = alpha**2 * (n + kappa) - n
print(f"lambda: {lam}")
# シグマ点を格納する配列 (2n+1個のシグマ点)
sigma_points = np.zeros((n, 2 * n + 1))
# 式1: 最初のシグマ点は平均自身
sigma_points[:, 0] = x_mean
# (n+lambda) * Px の平方根を計算 (コレスキー分解を使用)
sqrt_term = scipy.linalg.cholesky((n + lam) * x_P, lower=True)
# 式2, 3: 残りのシグマ点を計算
for i in range(n):
sigma_points[:, i + 1] = x_mean + sqrt_term[:, i]
sigma_points[:, i + n + 1] = x_mean - sqrt_term[:, i]
print("\n生成されたシグマ点:")
for i in range(2 * n + 1):
print(f" シグマ点 {i}: {sigma_points[:, i]}")
3. 変換と統計量の計算
計算されたシグマ点それぞれを非線形関数 $f$ で変換し、$y_{\sigma}$ を求めます。 $$ y_{\sigma, i} = f(\sigma_i) \tag{5} $$
次に、各シグマ点に対応する重み $w_i$ を計算します。 $$ w_0 = \frac{\lambda}{n + \lambda} \tag{6} $$ $$ w_i = \frac{1}{2(n + \lambda)} \quad \text{for } i=1, \dots, 2n \tag{7} $$
最後に、これらの重みと変換されたシグマ点 $y_{\sigma}$ を用いて、$y$ の平均 $\bar{y}$ と共分散 $P_y$ を計算します。 $$ \bar{y} \approx \sum_{i=0}^{2n} w_i y_{\sigma, i} \tag{8} $$ $$ P_y \approx \sum_{i=0}^{2n} w_i (y_{\sigma, i} - \bar{y})(y_{\sigma, i} - \bar{y})^T \tag{9} $$
# 変換されたシグマ点を格納する配列
sigma_y = np.zeros((m, 2 * n + 1))
# 式5: 非線形変換
for i in range(2 * n + 1):
sigma_y[:, i] = f(sigma_points[:, i])
# 重み関数を格納する配列
w = np.zeros(2 * n + 1)
# 式6, 7: 重み関数の計算
w[0] = lam / (n + lam)
for i in range(1, 2 * n + 1):
w[i] = 1 / (2 * (n + lam))
print(f"\n計算された重み: {w}")
# 式8: yの平均の計算
y_mean = np.sum(w * sigma_y, axis=1)
# 式9: yの共分散の計算
y_P = np.zeros((m, m))
for i in range(2 * n + 1):
diff = sigma_y[:, i] - y_mean
y_P += w[i] * np.outer(diff, diff) # outer積で(m,m)行列を生成
print(f"\n出力yの推定平均: {y_mean}")
print(f"出力yの推定共分散行列:\n{y_P}")
この実装により、非線形変換後の確率変数の平均と共分散を、少数のシグマ点から効率的に推定できることがわかります。このU変換は、Unscented Kalman Filter (UKF) の基盤となる技術です。