※本記事は、元々筆者が大学で学部学生向けに書いたチュートリアルをnote用に書き直したものである。手っ取り早くグラフを作りたいという人は参考にしてほしい。

はじめに

論文には定量的な結果を示すグラフを載せるのが一般的である。グラフを作成する方法はいくつかあり、Excelでもできないことはないが、RやPythonを用いればより洗練されたデザインのグラフを作成できる。
ggplot2はRのパッケージの１つで、シンプルかつスマートにデータを可視化できる。ggplot2に関するチュートリアル記事は多くネット上に存在するが、情報が断片的だったり、既存のデータセットを用いたものがほとんとで、実践的なものが少なかった。

本記事では、実際の研究のようにデータの準備からグラフ描画までの手順を一通り記述し、論文で用いられるようなグラフの作成方法を紹介する。一連のコードはR Studioで用いることを前提に書かれている。

なお、RやR Studioのインストールや基本的な操作については省略するので、ご了承願いたい。

下準備

データの準備

一般的なチュートリアルでは Iris などの既存のデータセットを使うことが多いが、ここでは自分の実験で得たデータを解析する練習として、オリジナルのデータを用いる。
基本的に Excel でデータをまとめる。今回は Control 群、実験群 Aの2群を用意したものとする（下のSampleを参照）。データには 2つの項目があり、group は群の名前を、value は計測値を入力してある。すなわち、1 行目に項目名が書いてあり、各項目の下にデータが記載されたファイルが必要となる。

データの読み込み

普通は csv ファイルで読み込むことが多いが、openxlsx というパッケージを使えば、Excel ファイルをxlsx ファイルでそのまま読み込むことができる。作業ディレクトリをSampleファイルがあるフォルダに設定し、このファイルを読み込んでみよう。

install.packages("openxlsx")
library(openxlsx)
data <- read.xlsx("Sample.xlsx")
data

ggplot2パッケージのインストール

ggplot2 は、tidyverse というデータ分析用パッケージ集の一部である。tidyverse には便利なツールが多いので、まとめてインストールしてしまおう。

install.packages("tidyverse")
library(ggplot2)

箱ひげ図を描く

本節では、Sample のデータをもとに箱ひげ図を描く方法を紹介する。

キャンバスの作成

ggplot関数を用いて、プロットの下地となるキャンバスを用意する。特に何も指定しなければ、地が灰色の、グリッド付きのキャンバスができる。

plot1 <- ggplot(data, aes(x = group, y = value))
plot1

背景の変更

グリッドはなくてもよいので、背景を真っ白なシンプルなものにする。
ggplot2では「＋」を使って指示を追加していくことができる。次のコードは、キャンバスにclassic（＝真っ白な背景）というテーマを追加することを意味している。

plot2 <- plot1 +theme_classic()
plot2

箱ひげ図の描画

geom_boxplot関数を用いて、キャンバスの上に箱ひげ図を描く。

plot3 <- plot2 + geom_boxplot()
plot3

着色

箱ひげ図の内部に色を付ける。エステティック属性のfillにgroupを指定することで、groupごとに色分けがされる。なお、fillはggplot内で指定してもよい。
また枠の色はデフォルトだと薄い黒色なので、純粋な濃い黒色にする。colorを指定して枠の色を変えることができる。
geom_boxplot内に書き加えるため、オブジェクトは上と同じplot3とする。

plot3 <- plot2 + geom_boxplot(aes(fill = group), color = "black")
plot3

色の変更

fillを指定した後は、scale_fill_manual関数を用いて、グラフの内部を好きな色にすることができる。ここでは、Control群を灰色、実験群を赤色にしてみる。デフォルトではアルファベット順にグループが表示されるため、ここでも色の指定はアルファベット順に行っている。

plot4 <- plot3 + scale_fill_manual(values = c("red3", "grey40"))
plot4

表示順の変更

デフォルトでは各グループがアルファベット順に表示される。一般的にControlは最も左に配置されるため、グループの表示を並べ替えて、Control群が左に来るようにする。scale_x_discrete関数を用いる。

plot5 <- plot4 + scale_x_discrete(limit = c("Control", "A"))
plot5

凡例を消去

グラフに着色した際、自動的に凡例が表示される。しかし、グループ名はx軸に表記されているため、わざわざ凡例を表示する必要はないだろう。凡例の位置を決めるlegend.positionをnoneに指定して、凡例を非表示にする。

plot6 <- plot5 + theme(legend.position = "none")
plot6

軸ラベルの変更

x軸とy軸には、そこに書かれている文字や数値が何を意味するのかを表すタイトルが表示される。現在は、ggplot内でｘはgroup、ｙはvalueと指定したため、それぞれこの名前になっている。x軸はグループ名であることが容易にわかるので非表示にし、y軸はvalueが何の値を示すのかわかるように名前を変えよう（ここではvalueは何らかの長さだとする）。x軸はxlab、y軸はylabで名前を指定できる。

plot7 <- plot6 + xlab("") + ylab("Length(μm)")
plot7

タイトル

labs関数を用いてグラフにタイトルをつける。

plot8 <- plot7 + labs(title = "Sample")
plot8

タイトルの調整

タイトルが左上に表示されているため、中央に移動させる。同時に文字サイズも調整する。

plot9 <- plot8 + theme(plot.title = element_text(hjust = 0.5, size = 12))
plot9

文字の調整

軸に表示されている文字や数値のサイズや色を調整する。デフォルトだと文字が小さく、薄い黒色で少し見にくいので、サイズを大きくし色も濃くする。　

plot10 <- plot9 + theme(axis.text.x = element_text(size = 12, color = "black"), axis.text.y = element_text(size = 12, color = "black"))
plot10

Swarm plotの追加（任意）

箱ひげ図などのグラフの上に各データの実測値の分布を表す点をプロットすると、より詳細な情報が得られて理解しやすくなることも多い。そのようなプロットができる関数にはgeom_jitterなどがあるが、ここではgeom_quasirandomを紹介する。quasirandomはjitterに比べてプロットが満遍なく配置される。

install.packages("ggbeeswarm")
library(ggbeeswarm) 

plot11 <- plot10 + geom_quasirandom(size = 1)
plot11

統計学的な検定

論文中のグラフには、標本間に有意差があることを示すアスタリスク（＊）が描かれていることが一般的である。有意差の有無は統計学的検定を行うことで明らかとなる。

ここでは、用いているデータが正規分布に従うと仮定し、t検定を行って2群間の平均値を比較する。
まずはControlとAのデータをそれぞれ抽出する必要がある。抽出方法は色々あるが、たとえば下のコードが使える。

Control <- data$value[data$group == "Control"]
A <- data$value[data$group == "A"]

次に、t検定を行う。

t.test(Control, A)

結果は下のようになる。

Welch Two Sample t-test

data:  Control and A
t = 16.404, df = 31.078, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 20.22835 25.97165
sample estimates:
mean of x mean of y 
    79.65     56.55 

有意水準を0.05とすれば、p<0.05であるからControl群と実験群Aには平均値に有意な差がある。

有意差の表示

この結果をもとに、グラフにアスタリスクを追加する。アスタリスクの数は一般にp<0.05なら１つ、p<0.01なら２つ、p<0.001なら３つ、p<0.0001なら４つである。したがって今回は４つとなる。
ここではannotate関数を用いて手作業で文字と線を追加する方法を紹介する。下のコードで、"text"は文字を、"segment"は線を描画する。それぞれ指定した座標に描画されるが、線は始点と終点の両方を指定する必要がある。

plot12 <- plot11 + annotate("text", x = 1.5, y = 92, label = "****", size = 7) + annotate("segment", x = 1, xend = 2, y = 92, yend = 92)
plot12

なお、統計学的検定とアスタリスクの表示にはggsignifというパッケージを用いる方法もある。それについては他の方の記事などを参照されたい。

また、今回はt検定を用いたが、データによっては他の検定が適切である場合もある。統計学的処理の詳細については本記事で扱う範疇を超えているので、自分で学習を進めてほしい。

ｙ軸の範囲の調整

上のグラフでは、タイトルとアスタリスクの間隔が狭くて少し窮屈に感じるかもしれない。y軸方向の表示範囲を変えることで、余白を設けてみよう。

plot13 <- plot12 + ylim(50,100)
plot13

グラフの保存

ここまでで作ったグラフを画像ファイルとして保存する。プロットウィンドウのExportから出力してもよいが、詳細な設定をするならggsave関数を用いるとよい。画像は、はじめに読み込んだデータファイルが存在するフォルダに出力される（ここではSampleファイルがあるフォルダ）。
今回用いたデータは2群からなるため、画像は縦に細長くすると見た目がよいだろう。下のコードではwidthを2、heightを4として出力した。

ggsave(filename = "boxplot.png", plot = plot13, device = "png", scale = 1, width = 2, height = 4, units = c("in"), dpi = 350) 

これで、下のようなグラフを作ることができた。

棒グラフを描く

データを可視化する際、棒グラフは最もよく用いられる。本節では棒グラフの描き方を紹介する。

上の箱ひげ図の描き方を読んだ上で本節を読んだほうが理解しやすいが、基本的にはコードをコピペすればグラフが作成できる。

箱ひげ図と同様にSampleのデータを用いて棒グラフを描いていこう。

平均値、標準偏差、標準誤差の計算

棒グラフは2つの標本の平均値を示すものとする。そのため、まずは平均値の計算をしなければならない。
そのような計算には、tidyverse内のdplyrパッケージを用いるとよい。

library(dplyr)

下のコードでは、データ内の標本をグループ名によってグループ分けし、valueの値をもとに標準偏差、平均値を計算した上で、さらに標準誤差を計算させている。

data.summary <- data %>% group_by(group) %>% summarize(n = n(), sd = sd(value), value = mean(value), se = sd / sqrt(n))
data.summary

このように、パイプ演算子%>%を用いると複数の関数を連結させることができる。
計算結果は下のようになる。

A tibble: 2 x 5
  group       n    sd value    se
<chr>   <int> <dbl> <dbl> <dbl>
1 A          20  3.24  56.6 0.724
2 Control    20  5.40  79.6 1.21 

グラフ描画

まず、シンプルな棒グラフを作成する。geom_colが棒グラフを描く関数である。平均値をもとに描画するため、データはdata.summaryのほうを指定していることに注意されたい。一方でキャンバスの元となるggplot内は大元のdataを指定している。
なお今回は、グラフの内部は白くし、枠だけに色がついたデザインにする。また棒グラフの幅を0.5とし、スマートな見た目にする。

plot1 <- ggplot(data = data, aes(x = group, y = value, color = group)) + geom_col(data = data.summary, fill = "white", width = 0.5) + theme_classic()
plot1

エラーバー（標準誤差）の表示

棒グラフにはエラーバーが付けられることが一般的である。エラーバーが表すものは標準偏差（standard deviation: sd）や標準誤差(standard error: se)などがあるが、ここでは標準誤差を用いる。geom_errorbar関数を用いて、平均値の上下にエラーバーを追加する。

plot2 <- plot1 + geom_errorbar(aes(ymin = value - se, ymax = value + se), data = data.summary, width = 0.2)
plot2

1つ1つのデータ点を追加

棒グラフに重ねて、１つ１つのデータをSwarm plotでプロットする。ここでもgeom_quasirandomを用いて、内部が白い、枠だけの点を追加する。

plot3 <- plot2 + geom_quasirandom(data = data, width = 0.2, shape = 1, size = 2)
plot3

ここからは箱ひげ図と同様の手順である。

色の変更

plot4 <- plot3 + scale_color_manual(values = c("red3", "grey40"))
plot4

グループの表示順の変更

plot5 <- plot4 + scale_x_discrete(limit = c("Control", "A"))
plot5

凡例を消去

plot6 <- plot5 + theme(legend.position = "none")
plot6

軸ラベルの変更

plot7 <- plot6 + xlab("") + ylab("Length(μm)")
plot7

文字サイズ・色の変更

plot8 <- plot7 + theme(axis.text.x = element_text(size = 12, color = "black"), axis.text.y = element_text(size = 12, color = "black"))
plot8

タイトル

plot9 <- plot8 + labs(title = "Sample") + theme(plot.title = element_text(hjust = 0.5), text = element_text(size = 12))
plot9

x軸とグラフの密着、y軸の表示範囲を変更

x軸と棒グラフの間に空間がないほうが見た目がよい。x軸とグラフを密着させると同時に、y軸方向の範囲も100まで増加させる。

plot10 <- plot9 + scale_y_continuous(expand = c(0,0), limits = c(0,100))
plot10

有意差の表示

plot11 <- plot10 + annotate("text", x = 1.5, y = 92, label = "****", size = 7) + annotate("segment", x = 1, xend = 2, y = 92, yend = 92)
plot11

グラフの保存

ggsave(filename = "barplot.png", plot = plot11, device = "png", scale = 1, width = 2, height = 4, units = c("in"), dpi = 350)

おわりに

本記事では、論文に掲載されるようなデザインの箱ひげ図と棒グラフの作成方法を紹介した。ただ、私自身はggplot2を学び始めたばかりであるし、プログラミングにもあまり詳しくない。そのため、もっと簡単な方法や美しくデザインする方法が他にあると思う。また、もしかしたら記事に書いたコードには不備があるかもしれない。そのあたりはどうかご容赦願いたい。そして、より良い方法があれば、ぜひコメントで教えていただきたい。

本記事がグラフを作成する際の一助となれば幸いである。ありがとうございました。

参考書籍