2012/08/31

文系のための「数の可視化」(7)

ヒストグラム箱ヒゲ図は、両方とも分散あるいは標準偏差の状況を
視覚的に解りやすく表現したものである。

ヒストグラムは、全体的な分布の形を視覚化し、
データの偏り方を視覚化する

一方、箱ヒゲ図は、全体的な分布の範囲を四分位で表し
データの偏り方を視覚化する

分散標準偏差というのは、全体のバラツキの程度を指標化しているだけ
それに対して、箱ヒゲ図ヒストグラムは、全体のデータの偏り方を観察できる。

また、標準偏差からの乖離の程度を見たり、
あるいは、対応する分布を考えるためにも用いる。

これらは、統計においては、非常に重要なことである。

さて、今回は、この辺りの話には深入りせず、
もう少し、データのバラツキの可視化方法について考えてみる。

「箱ヒゲ図」と「ヒストグラム」の両方の特性を併せ持った可視化方法がある。

ここで紹介する方法は、かなり魅力的。利用価値も高い
それにも関わらず、あまり、出回っていないから不思議である。
未だに、箱ヒゲ図を使った図しか見ない。

悪しきソフトウェア依存の習慣が原因だろうが。
まぁ、そういったことはどうでも良いことである。

さて、今回は以下の3つの方法を紹介する。
  • ヴァイオリン・プロット(violin plot)
  • ビーンズ・プロット(bean plot)
  • ビースウォーム(beeswarm)
基本的には、箱ヒゲ図であるが、ヒストグラムの特徴も持っている。
慣れると、箱ヒゲ図よりも得られる情報量が多く、可読性も高い

これらは、Rの初期状態では使えないパッケージのインストールから始める
Rコンソール上で、次のようなコマンドを入力する。

ダウンロード先の指定の画面が出てきたら、最も近い場所を選ぶ。
別に、どこを選んでも良いのだが、遠いとダウンロードに時間がかかる。

そうそう。ダウンロードだから、インターネットの接続環境が必要

install.packages("vioplot", dep=TRUE)        # ヴァイオリンプロット
install.packages("beanplot", dep=TRUE)     # ビーンプロット
install.packages("beeswarm", dep=TRUE)   # ビースウォーム

インストール作業は、最初の一回のみ。install.packages()関数を用いる。
引数にdep=TRUEとあるが、これは依存関係の確認の有無

Rでは、様々なパッケージを組み合わせるので、
目的のパッケージだけでは動かない場合がある。

TRUEにしておくと、足りないパッケージは自動的にインストール

さて、インストールができたら、今度は、ライブラリを読み込む
ライブラリの読み込みは、毎回必要
読み込まないと使えない。ヘルプも見れない

library(vioplot)
library(beanplot)
library(beeswarm)

では、さっそく。と言いたいところであるが...。

新しく、パッケージをインストールしたら、最初にヘルプを確認すること。

R使いは、ヘルプを見るのが当たり前。いや、見ないと何もできない。
教科書やホームページに頼ってばかりではいけない。

最初は、パッケージついての情報を確認する。

help(package="vioplot")

ここで重要なのは、「Index:」 と書かれている部分から下の部分
左端に「vioplot」とあり、空白の右側に「violin plot」と書いてある。

実は、左端が関数名で、右側はその説明
ふむ。vioplotの説明は、非常に不親切なようである。説明が少ない。

vioplotパッケージには、一つの関数しか含まれていないようであるが、
パッケージによっては、数十もの関数を含む場合もある。

これで、関数名が解ったので、今度はvioplot()関数ヘルプを確認する。

?vioplot        # あるいは、help(vioplot)

最も重要な部分は、「Usage:」の部分。パラメータの確認はここで。
その次に、「Arguments:」の部分で、各パラメータの説明を確認する。
最後に、「Example:」を確認して、使い方を確認する。

まぁ、理解しにくいものや、無理に理解する必要の無いものもあるのだけれど、
その辺りの判断には、ある程度の慣れが必要。

同様に、他のパッケージについても確認。これは各自で。
ひと通り、確認したところで、いよいよ、実際に使ってみる。

今回のデータは、何するか?迷うが、発電所のデータにしてみるか。

このデータは、平成22年度の「地域別発電電力量」のデータ
総務省統計局からダウンロードしたものを編集したもの。詳細は、棒グラフの話。

これをいつものように、

# Windows と Linux の人は次のコマンド
X <- read.table("clipboard", sep=",") 

# Mac の人は次のコマンド
X <- read.table(pipe("pbpaste") , sep=",")

入力して、その後に以下をコピーして、コマンドを実行する。

pnum,max,hydro_pnum,hydro_max,thermal_pnum,thermal_max,atomic_pnum,atomic_max
Hokkaidou,66,7419,53,1234,11,4065,1,2070
Touhoku,228,17206,209,2423,13,11286,2,3274
Tokyo,192,64988,162,8981,25,38696,3,17308
Chubu,197,32828,183,5219,11,23969,1,3617
Hokuriku,138,8057,127,1904,6,4400,1,1746
Kansai,165,34877,149,8196,12,16907,3,9768
Chugoku,110,11986,97,2906,12,7801,1,1280
Shikoku,65,6963,58,1141,4,3797,1,2022
Kyushu,194,20330,139,3279,45,11577,2,5258
Okinawa,22,1919,0,0,21,1919,0,0

前処理として、色を付けるための色の定義を。
これは、棒グラフの話で出てきた色の定義方法。

col.h <- rgb(50, 100, 255, maxColorValue = 255)  # 水力発電所のイメージ色
col.t <- rgb(255, 100, 50, maxColorValue = 255)  # 火力発電所のイメージ色
col.a <- rgb(200, 255, 100, maxColorValue = 255)  # 原子力発電所のイメージ色??

まずは、通常の箱ヒゲ図から。今回は発電所種別ごとの総出力で見てみる。
boxplot(X[,c(4,6,8)], names=c("Hydro","Thermal","Atomic"), col=c(col.h, col.t, col.a))
title(main="Maximum output of electronic plants in Japan", ylab="Maximum output (Kw)")

これは大丈夫のはず。これで、何の変哲もない、普通の箱ヒゲ図が描けた。
次に、ヴァイオリン・プロット。描き方はかなり異なる。
vec.h <- as.vector(X[,4])    # 元のデータから水力発電の総出力をベクトルで抽出
vec.t <- as.vector(X[,6])    # 元のデータから火力発電の総出力をベクトルで抽出
vec.a <- as.vector(X[,8])    # 元のデータから原子力発電の総出力をベクトルで抽出

vioplot(vec.h, vec.t, vec.a, names=c("Hydro","Thermal","Atomic"), col="lightpink")
title(main="Maximum output of electronic plants in Japan", ylab="Maximum output (Kw)")

violot()関数は、ベクトルデータを並べて指定する
そのため、vec.hvec.tvec.a の3つのベクトルを作成し、
それぞれ、水力発電、火力発電、原子力発電の総出力を代入。

実は、色を個別に付けるのは無理なので、ここでは、薄いピンク色。
ソースコードが公開されているので、修正することは難しくないが、
後々になって、説明することが面倒なので、そのままで。

中心の白い「○」が中央値
真ん中の黒い「帯」がIQR(四分位範囲)
真ん中を突っ切る「線」がIQRの±1.5倍のヒゲ

左右に広がる薄いピンクの「山」は、カーネルと呼ばれる密度推定の山
これが、ヒストグラムに相当する。カーネルの話はいずれ。

vioplot()関数は、ヒストグラムの特性と箱ヒゲ図の特性を合わせた特性を持ち、
重要なこととして、箱ヒゲ図で表現されている全て情報も表現している。
問題は、データの指定方法で、この関数だけデータの指定方法が特殊である。

今度は、ビーンプロットを試してみる。

beanplot(X[,c(4,6,8)], names=c("Hydro","Thermal","Atomic"), col="lightblue", 
    main="Maximum output of electronic plants in Japan", ylab="Maximum output (Kw)", cut=0)

beanplot()関数でのデータの選択は、箱ヒゲ図と同じ。

このグラフにおいて、真ん中の太い線は、中央値を表し、
複数の細い「線」は、データが存在している場所を示している。
横に広がる「山」は、vioplot()関数と同様にカーネルの山

vioplot()関数と比較すると beanplot()関数の山は細かく波を打っている。
これは、カーネルのパラメータの違い。カーネルは、様々な描き方ができる。

beaswarm()関数は、他の二つの方法と比較して、
分布に重点を置いているため、より詳細な密度推定が可能となる。
その一方で、箱ヒゲ図と比較して、外れ値やIQRの表現力は低い

注意するべきは、cut パラメータ。この値は「0」にする。
既定のままだと、データが存在しない所まで山裾が広がってしまう

色に関しては、複数のbean(マメ)の色を変えることはできず、
title()関数を使って、後から、タイトルを入れることもできない。

最後はビースウォーム
beeswarm(X[,c(4,6,8)], labels=c("Hydro","Thermal","Atomic"), col=c(col.h, col.t, col.a), pch=16)
title(main="Maximum output of electronic plants in Japan", ylab="Maximum output (Kw)")

beeswarm()関数は、分布をカーネルの山で表現せず、実際の点を布置する。

今回の場合は、対象数が少ないので、見難いが、
対象数がある程度あれば、この方法が最も誤解無く情報を伝えることができる。

boxplot()関数と異なるのは、ラベルの指定パラメータ。
boxplot()関数names であるのに対し、beeswarm()関数では labels となっている。

beeswarm()関数は、他の二つの方法と比較して、推定を用いていないため、
より実態に即した情報を伝えることができるが、外れ値やIQRの表現はできない
ただし、色の指定は可能なので、色分けをしたい場合に使えるというメリットがある。

さてさて。では、今回紹介した方法のうちで、最も適切なものはどれか?
結論から言うと、ケース・バイ・ケースであるが、無難な方法はビースウォーム

うん?色々と欠点があるのに?

実は、ビースウォーム・プロットの欠点は、
箱ヒゲ図と重ね合わせることで克服できる。

boxplot(X[,c(4,6,8)], names=c("Hydro","Thermal","Atomic"))
beeswarm(X[,c(4,6,8)], labels=c("Hydro","Thermal","Atomic"), col="red", pch=16, add=TRUE)
title(main="Maximum output of electronic plants in Japan", ylab="Maximum output (Kw)")

beeswarm()関数のパラメータの最後に add=TRUE としてある。
このようにすることで、前に描画したプロットに、新しいプロットを追加できる

なお、この方法は、他の方法にも適用できるが、
ヴァイオリン・プロットやビーンプロットを箱ヒゲ図に重ねるメリットは無い。


0 件のコメント:

コメントを投稿