グラフィックス参考実例集:棒グラフ、ヒストグラムを描く barplot, hist
(グラフィックス参考実例集に戻る。Rのグラフィックスパラメータを参照する。)
barplot 関数は棒グラフを描きます。棒グラフは、棒の長さにより観測値の大きさや度数を表すものです。数値ベクトルを引数として与えるのが、関数barplotもっとも簡単な使い方です。こうするとベクトルのそれぞれの値が棒の高さになります。 引数が行列ならば各列が一本の棒となり、列の各要素の値により棒をブロックに分けます。
似たものにヒストグラムを描く hist 関数がある。見かけは似ているが、棒グラフが類別・尺度(カテゴリー化された)データを扱うのにたいし、hist は間隔・比率尺度(本来の数値)データを適当な区間毎に集計したものを棒グラフ状に表示する。
注意。base パッケージ中の barplot 関数 (hist も?) は対数y軸をサポートしていない。また par(ylog=TRUE) による指定も無効にされる。アドオンパッケージ gregmisc 中には両対数軸をサポートする barplot 関数があるらしい。
barplot01 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot01.png") d <- ceiling(runif(5, 1, 10)) # サンプルデータ barplot(d) } barplot01()
barplot03 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot03.png") d <- ceiling(runif(5, 1, 10)) d.names <- paste("gene", 1:5) # x軸に表示される名前を作る barplot(d, names.arg=d.names, col=c("blue", "green", "red", "yellow", "cyan")) } barplot03()
棒グラフの棒の幅は1になっています(barplotのwidth引数)。棒と棒の間隔は標準で0.2になっています(barplotのspace引数)。この状態で左からn番目の棒の中心は(0.2 + 1) * n - 1/2になっています。ここからarrows関数で矢印を引っぱると誤差範囲が表示できます。
arrowsのangleは矢印の傾きを指定します。今回は直角になるように90にしてあります。
barplot04 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot04.png") d <- ceiling(runif(5, 1, 10)) d.names <- paste("gene", 1:5) d.error <- runif(5, 0, 1) # 誤差範囲のデータ d.error.x <- (0.2 + 1) * 1:5 - 1/2 # エラーバーを表示する位置(x座標) barplot(d, names.arg=d.names, col=c("blue", "green", "red", "yellow", "cyan")) arrows(d.error.x, d, d.error.x, d + d.error, angle=90) arrows(d.error.x, d, d.error.x, d - d.error, angle=90) } barplot04()
別の方法:関数barplotは各棒グラフの中点座標の配列を返すので,誤差範囲を表示するためにこれを利用できます.
barplot04b <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot04b.png") d <- ceiling(runif(5, 1, 10)) d.names <- paste("gene", 1:5) d.error <- runif(5, 0, 1) # 誤差範囲のデータ # 各棒グラフの中点座標を取得 d.error.x <- barplot(d, names.arg=d.names, col=c("blue", "green", "red", "yellow", "cyan")) arrows(d.error.x, d, d.error.x, d + d.error, angle=90) arrows(d.error.x, d, d.error.x, d - d.error, angle=90) } barplot04b()
barplot05 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot05.png") dm <- matrix(ceiling(runif(12, 1, 10)), ncol=3) # サンプルデータ barplot(dm) } barplot05()
棒グラフの間の線はsegments関数で引きます。segmentsには線分の始点と終点を与えないといけないのでまずそれを計算します。
barplot06 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot06.png") dm <- matrix(ceiling(runif(12, 1, 10)), ncol=3) d.seg.xs <- rep((0.2 + 1) * 1:(ncol(dm) - 1), each=nrow(dm)) # 線分の始点のx座標 d.seg.xe <- rep((0.2 + 1) * 2:ncol(dm) - 1, each=nrow(dm)) # 線分の終点のx座標 d.seg.ys <- apply(dm[,1:(ncol(dm) - 1)], 2, cumsum) # 線分の始点のy座標 d.seg.ye <- apply(dm[,2:ncol(dm)], 2, cumsum) # 線分の終点のy座標 barplot(dm) segments(d.seg.xs, d.seg.ys, d.seg.xe, d.seg.ye) # 線分を描く } barplot06()
barplot07 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot07.png") dm <- matrix(ceiling(runif(12, 1, 10)), ncol=3) d.seg.xs <- rep((0.2 + 1) * 1:(ncol(dm) - 1), each=nrow(dm)) d.seg.xe <- rep((0.2 + 1) * 2:ncol(dm) - 1, each=nrow(dm)) d.seg.ys <- apply(dm[,1:(ncol(dm) - 1)], 2, cumsum) d.seg.ye <- apply(dm[,2:ncol(dm)], 2, cumsum) d.names.r <- paste("exp", 1:nrow(dm)) # 凡例に表示する文字列 d.names.c<-paste("gene",1:ncol(dm)) # 棒グラフのグループごとの名前 barplot(dm, names.arg=d.names.c, legend.text=d.names.r) segments(d.seg.xs, d.seg.ys, d.seg.xe, d.seg.ye) } barplot07()
barplot08 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot08.png") dm <- matrix(ceiling(runif(12, 1, 10)), ncol=3) # サンプルデータ barplot(dm, beside=T) } barplot08()
barplot09 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot09.png") dm <- matrix(ceiling(runif(12, 1, 10)), ncol=3) d.names.r <- paste("exp", 1:nrow(dm)) d.names.c<-paste("gene",1:ncol(dm)) d.angle <- c(20,-20) * 1:nrow(dm) barplot(dm, beside=T, angle=d.angle, density=10, col="black", names.arg=d.names.c, legend.text=d.names.r) } barplot09()
標準では棒の幅1、グループ内の棒の間隔0、グループの間隔1で作成されています。このときm番目のグループのn番目の棒(d[n,m]に相当)の中心は(1+1*N)*(m-1)+1*n+1/2になります(ああ、ややこしい)。
arrowsのlengthは矢印の矢じりの部分の長さを指定します。棒の数が多い場合はこの値を小さくしないときれいな図になりません(標準で0.25)。
barplot10 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot10.png") dm <- matrix(ceiling(runif(12, 1, 10)), ncol=3) d.names.r <- paste("exp", 1:nrow(dm)) d.names.c<-paste("gene",1:ncol(dm)) dm.error <- runif(length(dm), 0,1) dm.error.x <- rep((nrow(dm)+1) * (1:ncol(dm) - 1), each = nrow(dm)) # (1+1*N)*(m-1)まで dm.error.x <- dm.error.x + rep(1:nrow(dm) + 1/2, ncol(dm)) # +1*n+1/2 の部分 barplot(dm, ylim=c(0, max(dm + dm.error)), beside=T, names.arg=d.names.c, legend.text=d.names.r) arrows(d.error.x, dm, d.error.x, dm + dm.error, angle=90, length=0.1) arrows(d.error.x, dm, d.error.x, dm - dm.error, angle=90, length=0.1) } barplot10()
別の方法:関数barplotは各棒グラフの中点座標の配列を返すので,誤差範囲を表示するためにこれを利用できます.
barplot10 <- function(){ on.exit(par(old.par <- par(no.readonly = TRUE))) on.exit(dev.off()) png("barplot10.png") dm <- matrix(ceiling(runif(12, 1, 10)), ncol=3) d.names.r <- paste("exp", 1:nrow(dm)) d.names.c<-paste("gene",1:ncol(dm)) dm.error <- runif(length(dm), 0,1) # 各棒グラフの中点座標を取得 dm.error.x <- barplot(dm, ylim=c(0, max(dm + dm.error)), beside=T, names.arg=d.names.c, legend.text=d.names.r) arrows(dm.error.x, dm, dm.error.x, dm + dm.error, angle=90, length=0.1) arrows(dm.error.x, dm, dm.error.x, dm - dm.error, angle=90, length=0.1) } barplot10()
ヒストグラムにより、データの分布の仕方を大雑把に知ることができます。このグラフは関数histで描きます。引数には数値ベクトルを指定します。他に細かなことを設定する引数がいくつかありますが、基本的にはこの引数を指定するだけでヒストグラムを描くことができます。例えば、
> hist(rnorm(100))
で、ヒストグラムが描かれます。 特に何も指定せずデータのみを与えた場合は、データの範囲を log2n + 1 ( n はデータの個数)個の階級に分割し、各階級に属すデータの数を棒グラフとして作図します。通常は絶対度数を使いますが、probability=Tと指定することによって相対度数で作図することもできます。また、引数nclassで階級数を明示的に指定することもできます。
histmisc1 <- function() { data(islands) op <- par(mfrow=c(2, 2)) hist(islands) str(hist(islands, col="gray", labels = TRUE)) hist(sqrt(islands), br = 12, col="lightblue", border="pink") ##-- For non-equidistant breaks, counts should NOT be graphed unscaled: r <- hist(sqrt(islands), br = c(4*0:5, 10*3:5, 70, 100, 140), col='blue1') text(r$mids, r$density, r$counts, adj=c(.5, -.5), col='blue3') sapply(r[2:3], sum) sum(r$density * diff(r$breaks)) # == 1 lines(r, lty = 3, border = "purple") # -> lines.histogram(*) par(op) }
既定の Sturges の方法は単に伝統があり、有名と言うだけで実際は使わない方が良さそう。 Q&Aコーナーの関連質問参照。
> x <- rnorm(1000)^2 > hist(x) > hist(x, breaks="Sturges") # 既定の hist(x) と同一 (古典的な Sturges の方法による分割点決定) > hist(x, breaks="scott") # Scott の方法による分割点決定 > hist(x, breaks="FD") # Frieman-Daiconis の方法による分割点決定
ついでに棒の間の隙間を色々変えてみる。
> par(mfrow=c(2,2)) > barplot(1:5, space=0, col="red", border="red") > barplot(1:5, space=0.1, col="red", border="red") > barplot(1:5, space=0.5, col="red", border="red") > barplot(1:5, space=1, col="red", border="red") > pp=recordPlot(); png("histnoborder.png"); replayPlot(pp); dev.off()
この例では塗りつぶしと枠を同じ色にすることで解決しているが、 そもそも枠を描きたくないという場合はオプションで「border=NA」を指定する。hist()でも同じ。
truehist は総面積が1となるヒストグラム(度数を区間幅で割る)を書く
library(MASS) data(islands) par(mfrow=c(1,2)) hist(islands) truehist(islands)
quantile 関数を用い各棒の面積(標本確率)が等しくなるように分割点を定める
x <- rnorm(400) nbin <- 10 # 区間の数を10にする hist(x,breaks=quantile(x,prob=seq(0,1,length=nbin+1)))