Gakushukun1’s diary

20代エンジニア, 統計的機械学習勉強中 twitter: @a96665004

PythonでBox-Cox変換を試してみる

今回の目的

正の実数に値をとるサンプルの分析のための前処理を行う方法として Box-Cox 変換を紹介する。



サンプル \{y_i\}正規分布でない分布から生成されている場合に, Box-Cox変換を用いることで, その分布を正規分布に近づけることができる.
Box-Cox変換は, 適当な \lambdaが与えられた上で, 次の式で表される:


y_i^{(\lambda)} = \begin{cases}
\frac{y_i^{\lambda} - 1}{\lambda}  & \lambda \neq 0\\
\log y_i & \lambda = 0
\end{cases}

scipy.stats.boxcoxを使うと, 適切な \lambdaを得た上で, 与えたサンプルをBox-Cox変換してくれる.
例として,


\begin{eqnarray*}
y \sim 2\exp(-2x) \\
y \sim {\rm Beta}(2, 5)
\end{eqnarray*}

の二つの分布から得られたサンプルをBox-Cox変換してみる.

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt


if __name__ == '__main__':

    plt.close('all')
    np.random.seed(0)

    y1 = np.random.exponential(scale=0.5, size=1000000)
    y1_lambda, lmax1 = stats.boxcox(y1)
    y2 = np.random.beta(2.0, 5.0, size=1000000)
    y2_lambda, lmax2 = stats.boxcox(y2)

    bins = 100
    plt.rcParams["font.size"] = 15
    plt.figure(figsize=(10, 10))
    plt.subplot(2, 2, 1)
    plt.hist(y1, bins=bins)
    plt.title(r'$y_1$')
    plt.subplot(2, 2, 2)
    plt.hist(y1_lambda, bins=bins)
    plt.title(r'$y_1^{(\lambda)} (\lambda=%f)$' % lmax1)
    plt.subplot(2, 2, 3)
    plt.hist(y2, bins=bins)
    plt.title(r'$y_2$')
    plt.subplot(2, 2, 4)
    plt.hist(y2_lambda, bins=bins)
    plt.title(r'$y_2^{(\lambda)} (\lambda=%f)$' % lmax2)
    plt.tight_layout()
    plt.savefig('figure.png')
    plt.show()

結果:
f:id:gakushukun1:20190429125931p:plain

内部的には, 対数尤度関数


llf(\lambda) = (\lambda - 1) \sum_i(\log(y_i)) - \frac{N}{2} \log(\frac{\sum_i (y_i ^{(\lambda)}- \bar{y^{(\lambda)}})^2 }{N})

を最大化するような \lambda -2.0 < \lambda < 2.0の範囲で見つけて返している.

まとめ

Box-Cox変換は最尤推定によってデータを正規分布からのサンプルに近づけるものであることが分かった. また, 実際にBox-Coxを使うことで, サンプルの分布を正規分布に近づけることができた.