ホーム | 統計 Top | 何が何でもRで分割表を100%積み上げ棒グラフで表示する
目次
- Rのbarplot()関数で100%積み上げ棒グラフを表示できないので
- 分割表形式の行列データを割合に変換する(自作関数)
- barplot()が返す座標について
- 100%積み上げ棒グラフに信頼区間のエラーバーを追加表示する
Rのbarplot()関数で100%積み上げ棒グラフを表示できないので
分割表データの入力から検定、グラフ作成までを一貫してRでやろうとすると少し面倒である。基本的には以下のようなマトリックスの形でデータを作成・保持するのだが、
> example <- matrix( c(67,128, 155,180), nrow=2, byrow=T) ##テストデータ作成 > example ##表示してみう [,1] [,2] [1,] 67 128 [2.] 155 180
これを100%積み上げ棒グラフで表示したいとき、単純に「棒グラフだったらbarplot()関数でいいんじゃろ」と考えて
> barplot(example, beside=F) ##棒グラフを表示 ## なおbeside=F(論理値FALSE)のとき、複数系列が積み上げグラフとして処理される。 ## TRUEならば横に並べられ、下の例だと「黒→灰色→黒→灰色」の順に4本のバーが描かれる。
などとやっても各棒の高さがバラバラで、そんなん考慮しとらんよ……ということになる。これだと100%ではない、ただの積み上げ棒グラフである。各列の度数が違うのだから当たり前だ。
【自作関数】分割表形式の行列データを割合に変換する【しょぼい】
そこで、分割表を構成する行列要素を、生の数字から割合に変換する関数を作ってみた。ひょっとして既に同等の機能が用意されているのかもしれないが、ことRに関しては調べるより自作したほうが早いので気にしない。
↓↓以下の行をRコンソールにコピペしてリターンを押すことで、関数がインストールされる("x2ratio"というオブジェクト名称を既に使っている場合は、適宜変更) ##function "x2ratio" (C) 2012 SUDO Masaaki ##分割表形式の行列の各要素を、横(デフォルト)もしくは縦方向の合計が1になるような割合値に変換する ##(たぶん)行列要素が全て実数の場合のみ有効 x2ratio <- function(x, byrow=FALSE) { if ( !is.matrix(x) ) return(NA) #引数xが行列以外ならば実行しない if ( !is.logical(byrow) ) return(NA) #引数byrowが論理値以外ならば実行しない if (byrow == FALSE) { return( t( t(x) / apply(x , 2, sum) ) ) } else return( x / apply(x , 1, sum) ) } ##ここまで ##使用例 > example ##サンプルデータ [,1] [,2] [1,] 67 128 [2.] 155 180 > x2ratio(example) ##早速実行してみう [,1] [,2] [1,] 0.3018018 0.4155844 [2.] 0.6981982 0.5844156 > barplot ( x2ratio(example), beside=F ) ##実行結果をネストして棒グラフを表示
これで100%積み上げ棒グラフを表示でき、大勝利である。
なお本関数のオプションとして、縦に(列について:デフォルト)集計するか、横に(行について)集計するかを論理値で選べる。
> x2ratio(example, byrow=F) ##縦に集計(省略時はこちら) [,1] [,2] [1,] 0.3018018 0.4155844 [2.] 0.6981982 0.5844156 > x2ratio(example, byrow=T) ##横に集計 [,1] [,2] [1,] 0.3435897 0.6564103 [2.] 0.4626866 0.5373134
barplot()が返す座標について
さて、上記の方法で何とかグラフを表示できたとする。これにエラーバーを追加したい場合、今のところbarplot()関数自体にはそのようなオプションが無いため、低水準作画関数を用いてシコシコ描画することになる。
##2×4分割表の例 > example24 <- matrix( c(67,128,3,22, 155,180,9,49), nrow=2, byrow=T) > example24 [,1] [,2] [,3] [,4] [1,] 67 128 3 22 [2,] 155 180 9 49 > barplot(x2ratio(example24), beside=F )
ここまでの処理により、積み上げ棒グラフが4本横に並んだものが描かれるはずだ。
Y座標について
さて、このグラフにおける座標値を考えよう。y座標は簡単で、グラフの目盛りの値をそのまま読めば良い。
> barplot(x2ratio(example24), beside=F ) ##グラフを描画 > abline(h=0.5, col=rgb(1,0,0,0.8)) ##y=0.5の位置に水平線を赤色で追加する
このabline()は、Rの描画デバイス内に直線を追加描画するための関数であり、上記処理によってy=0.5の位置に水平な直線が追加された。
X座標について
これに対して、x座標の読み取りは少しやっかいである。
barplot(x2ratio(example24), beside=F ) abline(v=0, col=rgb(1,0,1,0.8)) ##h=0の位置に垂線を紫色で追加する abline(v=1, col=rgb(0,0,1,0.8)) ##h=1の位置に垂線を青色で追加する
つまりx=0となる原点は、1本目のプロットの左端ではなく、y軸線上である。そしてy軸と1本目プロットの左端との間にマージンがあり、棒が増えるごとにマージンの幅が積み増しされてゆく。
これはどういうことなの......とbarplot()関数のヘルプファイルを読むと、"space"という設定項目があり、"the amount of space (as a fraction of the average bar width) left before each bar. May be given as a single number or one number per bar."と記されている。さらにこの値はbeside=Fの場合、デフォルトでは0.2になるとも書かれている。
つまり、1本目のプロットより左にマージン分の0.2が存在し、そこから系列ごとに棒の幅である1と、マージンの幅である0.2が追加されていくってことか。
barplot(x2ratio(example24), beside=F ) abline(v=0.7, col=rgb(1,0,1,0.8), lwd=2) ##h=(0.2+0.5)の位置に垂線を紫色で追加する abline(v=1.2, col=rgb(0,0,1,0.8), lwd=2) ##h=(0.2+1.0)の位置に垂線を青色で追加する abline(v=1.4, col=rgb(1,0,0,0.8), lwd=2) ##h=(0.2+1.0+0.2)の位置に垂線を赤色で追加する text(0.7,0.1,"0.7") ##(x,y)=(0.7,0.1)にラベルを追加。位置は標準で中央揃えとなる text(1.2,0.2,"1.2") text(1.4,0.3,"1.4")
これで、x,y座標を自在に指定して図形ないしテキストを追加表示できる。科学の勝利なのであった。
barplot()関数が返すx座標の再利用
実際には、さらに簡単なx座標の取得方法がある。それは、barplot()関数が返す値の再利用である。
hogehoge <- barplot(x2ratio(example24), beside=F ) ##barplot()関数の結果を、適当な名前で新規オブジェクトに代入すると...... hogehoge [1] 0.7 1.9 3.1 4.3
このベクトルhogehogeが表すものは、barplot()関数によって描画された棒グラフにおける、4本のプロットの各中心に他ならない。
100%積み上げ棒グラフに信頼区間のエラーバーを追加表示する
上記を踏まえて、barplot()関数で描いた積み上げ棒グラフに、二項信頼区間のエラーバーを追加してみよう。
> example24 [,1] [,2] [,3] [,4] [1,] 67 128 3 22 [2,] 155 180 9 49 > binom.test( c(example24[1,1],example24[2,1]) ) Exact binomial test data: c(example24[1, 1], example24[2, 1]) number of successes = 67, number of trials = 222, p-value = 3.299e-09 alternative hypothesis: true probability of success is not equal to 0.5 95 percent confidence interval: 0.2421815 0.3668047 sample estimates: probability of success 0.3018018
たとえば左端の系列である(67, 155)について割合は67/(67+155) = 0.3018、二項分布に基づく95%信頼区間は0.242-0.367と計算される。
エラーバーとして、今回は両端の「ひげ」を付けない線を描く。Rで線分を描くための低水準作画関数である、segments()を用いる。←segment()ではないので注意
barplot(x2ratio(example24), beside=F ) segments(0.7, 0.242, 0.7, 0.367, lwd=2) ##信頼区間の線分を追加する ##なお、座標は左から順に(始点のx)(始点のy)(終点のx)(終点のy) ##lwd=2というのは、線幅を標準の2倍にするという意味のオプション