基礎医学生物学研究者のブログ

ggplot2による棒グラフの作成方法

日付；2022/07/15（土）
Rによるデータ解析（その４）

最近では、どうもNatureやCellやScienceだけではなく、他のジャーナルでも棒グラフに個々のデータをいプロットする論文が増えてきた。ジャーナル側からの要求なのか、それともただイキってNature、Cell、Scienceを真似しているだけなのかは知らない。

尚、基礎生物学研究必要なのは平均値と標準偏差であると確信している。そして、基礎生物学研究で標準誤差を使っている論文があれば、それはエセ研究者によるデータのバラツキを視覚的に小さく見せかけるためのチート技術なので、要注意である。そういったエセ統計解析を見破るためには、現時点ではこの棒グラフにドットを乗せるのはけっこう有用である。ちなみに、標準誤差は、母集団の平均値の分布というか範囲を推定するものであり、それを適切に使用するためにはかなりのN数が必要だし、それを計算したら母集団の情報、すなわち、サプライヤーで飼育している同じ週齢のマウス全体の情報を意味することになり、そんなものに興味はないと思う。標準誤差を使っている論文は、暗にその様なデータを示していることになる。一方、標準偏差はデータのバラツキなので、個々の実験で必要なのはやっぱり標準偏差であると考えている。突き詰めればキリがないような気もするが。標準誤差は疫学のデータでは必須であると思う。

そういうことなので、ggplot2を使って棒グラフの上に各データポイントをプロットする方法について記す。

目次

必要なライブラリ

library(tidyverse)

重要なこと

元のデータセットOriginalDataSetから平均値と標準偏差を計算して、Mean_SDというデータセットをつくり、次にそのデータセットにDrugやTimeなどの列を新しく作っていくが、この新しく作ったDrugやTimeなどは元のデータセットのそれとと同じ値、同じ文字列、同じ列名にしたほうが、ggplot2を使った作図上都合が良いように思う。言うても、data = OriginalDataSetとか、使うデータセットをいちいち指定しているから、問題はないのだが、個人的な整理整頓のつもりで。

1. 単純なドットあり棒グラフを作る。

これは最も単純なドットありの2群の棒グラフを作るためのコードである。エンドポイントだけの2群の場合などが、これに相当するだろうと思う。

1.1. 平均値と標準偏差を計算する。

まずは平均値と標準偏差を計算する。これは、OriginalDataSetというデータセットを使って、Drug毎の平均値と標準偏差を計算して、新しく Mean_SDというデータセットを作っている。

Mean_SD <- OriginalDataSet %>% group_by(Drug) %>%
  summarise(
    Mean = mean(Weight),
    SD = sd(Weight))

1.2. 棒グラフを書く。

次に、棒グラフを書く。研究ではこれはよく使う。動物実験をしていると、やっぱり個体毎のデータが把握できるとすごく便利だと思う。

ggplot(Mean_SD, aes(Drug, Mean)) +
  geom_col(color = "black", fill = c("blue", "orange")) + 
  geom_errorbar(aes(ymin=Mean-SD, ymax=Mean+SD), width = 0.2) +
  geom_jitter(
    data = OriginalDataSet,
    mapping = aes(x = Drug, y = Weight),
    position = position_jitter(width = 0.1, height = 0, seed = NA),
    size = 1.5)+
  theme_classic()

2. 経時的なドットあり棒グラフを作る。

次に、すこし複雑な棒グラフについて書く。上記は単純にエンドポイントでの2群の棒グラフだったが、これは24、48、72時間での2群の棒グラフである。

2.1. 平均値と標準偏差を計算する。

これは以前に書いた方法と同じである。これは、OriginalDataSetというデータセットを使って、Drug_Day毎の平均値と標準偏差を計算して、新しく Mean_SDというデータセットを作っている。

Mean_SD <- OriginalDataSet %>% group_by(Drug_Day) %>%
summarise(
mean = mean(Percentage),
sd = sd(Percentage))

次に、後々グラフを作りやすいようにMean_SDというデータセットにいちいちラベルしていく。というのも、どうやらこのときは投与薬剤（Drug）と投与日（Day）をくっつけたデータセットを作ってしまったらしいので、それをちまちまとラベルし直していく。separete()とか使えばよいのだが、ここではわかりやすくする必要があったので、新しく列を作ることにした。mutate()を使ってDrugとTimeという列をMean_SDに追加する。

Mean_SD <- mutate(
Mean_SD,
Drug =
case_when(
str_detect(Drug_Day, "DMSO") == "TRUE" ~ "DMSO",
str_detect(Drug_Day, "XX") == "TRUE" ~ "XX"
))

Mean_SD <- mutate(
Mean_SD,
Time =
case_when(
str_detect(Drug_Day, "D1") == "TRUE" ~ "24",
str_detect(Drug_Day, "D2") == "TRUE" ~ "48",
str_detect(Drug_Day, "D3") == "TRUE" ~ "72"
))

2.2. 棒グラフを書く。

以下はDMSOとXXという試薬の2群の棒グラフを、24、48、72時間毎に並べて書くというコードである。

ggplot(Mean_SD, aes(Time, mean, fill = Drug)) +
geom_bar(data = Mean_SD, mapping = aes(Time, mean, fill = Drug), stat="identity", position="dodge", color = "black", width = 0.95) +
geom_errorbar(aes(ymin=mean-sd, ymax=mean+sd), width = 0.6, position=position_dodge(0.95)) +
geom_jitter(
data = OriginalDataSet,
mapping = aes(x = Day, y = Percentage, group = Drug),
position = position_jitterdodge(jitter.width = 0.61, jitter.height = 0, dodge.width = 0.95, seed = NA), size = 1.5) +
coord_cartesian(ylim = c(0, 70), clip = "on") +
scale_fill_grey(start=1.0, end=0.5) +
theme_classic() +
theme(plot.margin=margin(0.1, 0.1, 0.1, 0.1, "cm"), axis.ticks.length=unit(0.4, "cm"), legend.position = "none")

2.3. 別バージョンの棒グラフをかく。

これは、上記よりも自然なggplotのドットありの棒グラフ（geom_dotplot()）だと思うが、これ、近いデータのドットがどうやっても重なってしまったので、どうも好きになれなかった。重なってしまったがために、例えば「Authers mentioned that the data were obtained at least fro 3 independent experiments, but the dots of Figure XXXX seem just 2 points. Authors must check the data again.」とかいうレビューをこんなところで食らいたくない。 

ggplot(Mean_SD, aes(Time, mean, fill = Drug)) +
geom_bar(data = Mean_SD, mapping = aes(Time, mean, fill = Drug), stat="identity", position="dodge", color = "black") +
geom_errorbar(aes(ymin=mean-sd, ymax=mean+sd), width = 0.6, position=position_dodge(0.95)) +
geom_dotplot(
mapping = aes(x = Day, y = Percentage, fill = Drug),
data = OriginalDataSet,
binaxis = "y",
stackdir = "centerwhole",
position=position_dodge(0.95),
dotsize = 2.0) +
scale_fill_manual(values = c("white", "gray"))+
coord_cartesian(ylim = c(0, 70), clip = "on") +
theme_classic() +
theme(plot.margin=margin(0.1, 0.1, 0.1, 0.1, "cm"), axis.ticks.length=unit(0.3, "cm"), legend.position = 'none')

3. 経時的な多群のドットあり棒グラフ

次は、多群、この場合は4群の棒グラフを経時的に並べる方法である。

3.1. 平均値と標準偏差を計算する。

まずは例に漏れずに平均値と標準偏差を計算する。

Mean_SD <- OriginalDataSet %>% group_by(Time_Drug) %>%
summarise(
mean = mean(Percent_of_gated_cells),
sd = sd(Percent_of_gated_cells))

次に、時点や薬剤濃度をわかりやすくラベルする。

Mean_SD <- mutate(
Mean_SD,
Drug =
case_when(
str_detect(Time_Drug, "DMSO") == "TRUE" ~ "DMSO",
str_detect(Time_Drug, "0.125uM") == "TRUE" ~ "0.125uM",
str_detect(Time_Drug, "0.25uM") == "TRUE" ~ "0.25uM",
str_detect(Time_Drug, "0.5uM") == "TRUE" ~ "0.5uM"
))

Mean_SD <- mutate(
Mean_SD,
Time =
case_when(
str_detect(Time_Drug, "12_") == "TRUE" ~ "12",
str_detect(Time_Drug, "24_") == "TRUE" ~ "24",
str_detect(Time_Drug, "36_") == "TRUE" ~ "36"
))

3.2. 棒グラフを書く。

ggplot(Mean_SD, aes(Time, mean, fill = Drug)) +
geom_bar(data = Mean_SD, mapping = aes(Time, mean, fill = Drug), stat="identity", position="dodge", color = "black", width = 0.95) +
geom_errorbar(aes(ymin=mean-sd, ymax=mean+sd), width = 0.6, position=position_dodge(0.95)) +
geom_jitter(
data = OriginalDataSet,
mapping = aes(x = Time, y = Percent_of_gated_cells, group = Drug),
position = position_jitterdodge(jitter.width = 0.41, jitter.height = 0, dodge.width = 0.95, seed = NA), size = 1.5) +
coord_cartesian(ylim = c(0, 70), clip = "on") +
scale_fill_grey(start=1.0, end=0.25) +
theme_classic() +
theme(plot.margin=margin(0.1, 0.1, 0.1, 0.1, "cm"), axis.ticks.length=unit(0.4, "cm"), legend.position = "none")