カーネル密度推定グラフの描画 with Excel
イントロダクション
Step 0
下表はある調査により得た1変数の標本です(n=30)。
DL
この標本をヒストグラムに通してみると,二峰性の分布であることがわかります。
こうした特徴を持つ標本で,カーネルカーブを描画したいと思います。以下,Excelでこれをおこなうためのstepです。
作図にあたり参考にした資料
- [PDF]鈴木大慈『データ解析 第十回「ノンパラメトリック密度推定法」』
- [Web]カーネル密度推定 ―"Wikipedia"
- [Web]データの分布とヒストグラム・密度推定 ―"R-Tips"
工程
Step 1前提
カーネル密度推定量f̂(x)は下式にて求めます。
ここで前提として,「R」のdensity()関数デフォルトの設定,すなわち引数bwをnrd0でプロットしたときと同様の成果物を得ることを目的に据えたいと思います。すなわちバンド幅hは,Silvermanの方法
(ただしsは標準偏差)で求め,カーネル関数は
すなわちガウス関数を利用するものとします。
以下これらの式は,順に❶式,❷式,❸式と呼称します。
Step 2「2」式の計算(1)
❷式をシートの上方で求めます。
とはいえ現状では当該部分に利用可能なスペースが残されていないので,シートの最上部にあらたに3つの行を追加したうえで,
ここ(セルA1, A2)に見出し「n」「h」を用意することにします。
Step 3「2」式の計算(2)
「n」については関数で求めるか,直接タイプしておきます(下式は“関数”の場合)。
つづいて「h」を求めます。これは,❷式をExcel上の計算式にて配置します。
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) |
Step 4「3」式の計算
この作例でいう変数名「x」の直右の列から,❶式にいういくらかのxを用意します。換言すれば,成果物としたいグラフ横軸のレンジをイメージして,この数字を手繰っておく必要があります。
具体的に,ここでは0-100の範囲でグラフを描きたいと考えています。ということで0から100までを任意のスパンで埋めていく必要があるのですが,当然,より短いスパンのほうがより忠実なグラフを描画することができます。もっとも,その分処理すべき計算の量がかさむこととなるので,どの程度の“忠実”を目指すかで適宜斟酌すればいいかと思います。
ともあれ,この作例では下図のように0-100を増分2で埋めることに決めました。
下図彩色セルで❸式を計算します。
具体的に,これは下式のとおりです。
B5 | =( 1 /SQRT(2*PI()) ) * EXP( -(((B$4-$A5) /$B$2)^2 /2) ) |
---|
入力後,この式を表中すべての領域にコピーしておきます。
Step 5「1」式の計算
表の最下行の直下のセル(下図彩色部分)で,❶式を計算します。
B35 | =( 1 /($B$1*$B$2) ) * SUM(B5:B34) |
---|
これを右方の残りのセルにすべてコピーします。
Step 6グラフの作成
下図強調部分を選択した状態から,
平滑線散布図を挿入します。これにより,下図のグラフが出力されます。
Step 7完成
任意の書式設定を重ね,カーネルカーブの完成です。
拡張
累積確率線を追加する
参考までに,上のグラフを確率を計算する目的で使用する意図があるとして,ここに累積確率線を加えてみたいと思います。ただしxの各区間の確率は,簡便な台形公式を利用して求めます。
C36 | =(B35+C35)*(C4-B4) /2 |
---|---|
C37 | =SUM($C$36:C36) |
上のシートにいう「累積Δ」を既存のグラフにあたらしい系列として加えると,グラフは下図のように変化します。オレンジの線が目的の累積確率線です。
これによりたとえばxが70弱のとき,下側確率は0.8であることがざっくりと確認できます(=下側80%点は,xで70弱が対応)。
その他の参照
このサイトの関連How-toです。
メインサイト「ひとりマーケティングのためのデータ分析」の度数分布に関するHow-toです。