BDAstyle

ビジネスデータ分析ツールの作成 with Excel

カーネル密度推定グラフの描画 with Excel

イントロダクション

Step 0

下表はある調査により得た1変数の標本です(n=30)。

DL

カーネル密度推定・元データ

この標本をヒストグラムに通してみると,二峰性の分布であることがわかります。

ヒストグラムでデータの分布を見る→二峰型

こうした特徴を持つ標本で,カーネルカーブを描画したいと思います。以下,エクセルでこれをおこなうためのstepです。

Excelで描いたKernelCurve

作図にあたり参考にした資料

工程

Step 1前提

カーネル密度推定量f̂(x)は下式にて求めます。

[LaTeX] \displaystyle \hat{f}_h(x) = \frac{1}{nh} \sum_{i=1}^{n} K \left( \frac{x-x_i}{h} \right)

ここで前提として,「R」のdensity()関数デフォルトの設定,すなわち引数bwnrd0でプロットしたときと同様の成果物を得ることを目的に据えたいと思います。すなわちバンド幅hは,Silverman の方法

[LaTex] \displaystyle h = \frac{0.9 \times min \Bigl\{ s, \frac{IQR}{1.34} \Bigr\} }{n^{\frac{1}{5}}}

(ただしsは標準偏差)で求め,カーネル関数は

[LaTex] \displaystyle K(x) = \frac{1}{\sqrt{2\pi}}e^{-\frac{x^2}{2}}

すなわちガウス関数を利用するものとします。

以下これらの式は,順に式,式,式と呼称します。

Step 2「2」式の計算(1)

式をシートの上方で求めます。

とはいえ現状では当該部分に利用可能なスペースが残されていないので,シートの最上部にあらたに3つの行を追加したうえで,

ここ(セルA1, A2)に見出し「n」「h」を用意することにします。

見出しn,hの追加

Step 3「2」式の計算(2)

「n」については関数で求めるか,直接タイプしておきます(下式は“関数”の場合)。

つづいて「h」を求めます。これは,式をエクセル上の計算式にて配置します。

B1 =COUNT(A5:A34)
B2 =( 0.9 * MIN( STDEV.S(A5:A34), (QUARTILE.INC(A5:A34, 3)-QUARTILE.INC(A5:A34, 1)) /1.34 ) ) /B1^(1/5)

n,hを求める

Step 4「3」式の計算

この作例でいう変数名「x」の直右の列から,式にいういくらかのxを用意します。換言すれば,成果物としたいグラフ横軸のレンジをイメージして,この数字を手繰っておく必要があります。

具体的に,ここでは0-100の範囲でグラフを描きたいと考えています。ということで0から100までを任意のスパンで埋めていく必要があるのですが,当然,より短いスパンのほうがより忠実なグラフを描画することができます。もっとも,その分処理すべき計算の量がかさむこととなるので,どの程度の“忠実”を目指すかで適宜斟酌すればいいかと思います。

ともあれ,この作例では下図のように0-100を増分2で埋めることに決めました。

0-100まで増分2でxを作成

下図彩色セルで式を計算します。

具体的に,これは下式のとおりです。

B5 =( 1 /SQRT(2*PI()) ) * EXP( -(((B$4-$A5) /$B$2)^2 /2) )

入力後,この式を表中すべての領域にコピーしておきます。

セルB5に計算式を入力し,表を埋める

Step 5「1」式の計算

表の最下行の直下のセル(下図彩色部分)で,式を計算します。

B35 =( 1 /($B$1*$B$2) ) * SUM(B5:B34)

これを右方の残りのセルにすべてコピーします。

セルB35に計算式を入力し,表の最も右の列まで埋める

Step 6グラフの作成

下図強調部分を選択した状態から,

xとxについてのカーネル関数の合計を選択

平滑線散布図を挿入します。これにより,下図のグラフが出力されます。

先のデータ範囲をソースにして,散布図を挿入する

Step 7完成

任意の書式設定を重ね,カーネルカーブの完成です。

Excelで作成したKernelCurve

その他の参照