
ゆうき( @BASEBALLY15 )です^^
(読者からの質問)カーネル密度推定で両側3%の値を取得したいです・・・
分かりました^^
それでは,今回は,その質問にお答えします!

カーネル密度推定についてのご質問
Question
Pythonでカーネル密度推定を用いて、分布を推定し、その分布に対して正規分布の平均値±3σの様な、上限・下限を取りたいと思っているのですが、カーネル密度推定で出た分布の面積の両側3%にあたるX軸の値を 求めるという事は可能でしょうか?
Answer
自分自身、勉強不足で解決することが難しかったため、テラテイルで質問しました.
その回答が、上のURLとなります.
この回答を参考にし,以下に両側3%のx軸値を求めました.
カーネル密度推定の両側3%のx軸値について
今回は,アイリス(iris)のデータを使いました.
import numpy as np
import pandas as pd
from pandas import Series, DataFrame
import matplotlib.pyplot as plt
import seaborn as sns
sns.set_style("whitegrid")
#ガウシアン分布のためのライブラリ
from scipy.stats import gaussian_kde
#pdfを全区間で積分するためのライブラリ
from scipy.integrate import cumtrapz
%matplotlib inline
#irisのデータを使用
iris_data = sns.load_dataset("iris")
#表示
iris_data
sepal_length | sepal_width | petal_length | petal_width | species | |
---|---|---|---|---|---|
0 | 5.1 | 3.5 | 1.4 | 0.2 | setosa |
1 | 4.9 | 3.0 | 1.4 | 0.2 | setosa |
2 | 4.7 | 3.2 | 1.3 | 0.2 | setosa |
3 | 4.6 | 3.1 | 1.5 | 0.2 | setosa |
4 | 5.0 | 3.6 | 1.4 | 0.2 | setosa |
#sepal_lengthのデータを抽出
iris_sepal_length = iris_data["sepal_length"]
#表示
iris_sepal_length
0 5.1 1 4.9 2 4.7 3 4.6 4 5.0 ... 145 6.7 146 6.3 147 6.5 148 6.2 149 5.9
#sepal_lengthのデータから,kdeで推定する
iris_kernel = gaussian_kde(iris_sepal_length)
#積分範囲
iris_xs = np.linspace(0, 10, num=1000)
#yの値を取得
iris_ys = iris_kernel(iris_xs)
#積分して累積分布関数を求める
iris_cdf = cumtrapz(iris_ys, iris_xs)
以下のコードを記述し,カーネル密度推定で求めた両側3%のx軸値を求めたいと思います.
#cdf(x) = 0.03となるxを求める(両側)
idx3= np.searchsorted(iris_cdf, 0.03)
idx97 = np.searchsorted(iris_cdf, 0.97)
#実行
print(iris_xs[idx3])
print(iris_xs[idx97])
3%の時のxの値 = 4.374374374374375
97%の時のxの値 = 7.607607607607608
最後に,この結果を可視化したいと思います.
#可視化
plt.plot(iris_xs, iris_ys, label="KDE")
#3%の範囲を表示
plt.fill_between(iris_xs[:idx3], 0, iris_ys[:idx3], facecolor="r", alpha=0.5)
#97%の範囲を表示
plt.fill_between(iris_xs[idx97:], 0, iris_ys[idx97:], facecolor="r", alpha=0.5)
#凡例を右上に表示
plt.legend(loc="upper right")

終わりに
今回は,カーネル密度推定に関して,ご質問を頂いたので回答しました.
まだまだ,力不足な部分があるので,今後,力をつけていきたいと思います^^
それでは・・・
早急な回答ありがとうございます。
知りたかった結果です。
ありがとうございました。
こちらこそ、勉強になりました。
ありがとうございます!
今後も、よろしくお願い致します。