tidyverse〜group_by() の使い方

バイオインフォマティクスでは、遺伝子発現データや配列データなど、大量のデータを効率よく処理する必要があります。その際、R の tidyverse パッケージ群は非常に有用です。特に、dplyr の group_by() 関数は、データをグループ単位で処理する際に欠かせない機能です。

本記事では、group_by() の基本から応用まで を、バイオインフォマティクスでの活用例を交えながら解説します。

1. group_by() の基本
1. 1.1 group_by() の役割
2. 1.2 group_by() の基本的な使い方
  1. 出力
2. group_by() の応用
3. バイオインフォマティクスでの実践的な活用例
1. 3.1 RNA-Seq データの解析
2. 3.2 タンパク質データの正規化
4. まとめ
1. group_by() を活用するポイント

1. group_by() の基本

1.1 group_by() の役割

group_by() は、データフレームを特定の列ごとにグループ化し、そのグループ単位で集計や操作を行うための関数です。通常、summarize() や mutate() と組み合わせて使われます。

1.2 group_by() の基本的な使い方

まず、簡単なデータフレームを作成し、group_by() を適用してみます。

library(tidyverse)

# サンプルデータ（遺伝子発現データを想定）
df <- tibble(
  gene = c("GeneA", "GeneA", "GeneB", "GeneB", "GeneC", "GeneC"),
  condition = c("Control", "Treatment", "Control", "Treatment", "Control", "Treatment"),
  expression = c(10.5, 15.2, 8.1, 12.3, 20.3, 18.7)
)

# グループ化して平均値を計算
df %>%
  group_by(gene) %>%
  summarize(mean_expression = mean(expression))

出力

# A tibble: 3 × 2
  gene  mean_expression
  <chr>          <dbl>
1 GeneA          12.85
2 GeneB          10.2
3 GeneC          19.5

group_by(gene) によって遺伝子ごとにグループ化し、各遺伝子の発現量の平均を求めています。

2. group_by() の応用

2.1 複数の列でグループ化

バイオインフォマティクスでは、異なる条件（例えば、異なる実験条件やサンプル群）でデータを比較することが重要です。そのため、group_by() を複数の列に適用することがあります。

df %>%
  group_by(gene, condition) %>%
  summarize(mean_expression = mean(expression))

出力

# A tibble: 6 × 3
  gene  condition mean_expression
  <chr> <chr>             <dbl>
1 GeneA Control           10.5
2 GeneA Treatment         15.2
3 GeneB Control            8.1
4 GeneB Treatment         12.3
5 GeneC Control           20.3
6 GeneC Treatment         18.7

このように、遺伝子ごと、かつ条件ごとの平均発現量を計算できます。

2.2 グループ単位での統計量の計算

例えば、発現量の平均値だけでなく、標準偏差 (sd()) も求めることで、発現のばらつきを評価できます。

df %>%
  group_by(gene) %>%
  summarize(
    mean_expression = mean(expression),
    sd_expression = sd(expression)
  )

2.3 mutate() と group_by() の組み合わせ

場合によっては、グループごとに計算した値を元のデータに追加したいことがあります。これは mutate() と group_by() を組み合わせることで実現できます。

df %>%
  group_by(gene) %>%
  mutate(mean_expression = mean(expression))

このコードでは、各 gene ごとに発現量の平均値を計算し、それを mean_expression 列としてデータフレームに追加しています。

3. バイオインフォマティクスでの実践的な活用例

3.1 RNA-Seq データの解析

RNA-Seq データでは、各遺伝子の発現量を比較する際に group_by() が頻繁に使われます。

例えば、異なる条件（コントロール vs 処理群）で発現変化を見たい場合、group_by() と summarize() を組み合わせます。

rna_seq_data <- tibble(
  gene = rep(c("GeneA", "GeneB", "GeneC"), each = 4),
  sample = rep(c("S1", "S2", "S3", "S4"), times = 3),
  condition = rep(c("Control", "Treatment"), each = 2, times = 3),
  expression = c(12, 14, 25, 30, 8, 10, 18, 22, 20, 22, 35, 40)
)

# 条件ごとの平均発現量
rna_seq_data %>%
  group_by(gene, condition) %>%
  summarize(mean_expression = mean(expression), .groups = "drop")

このように、RNA-Seq データの発現量をグループごとに処理することで、発現変化の解析が容易になります。

3.2 タンパク質データの正規化

質量分析によるプロテオミクスデータでは、グループごとに発現量を正規化することがあります。

proteomics_data <- tibble(
  protein = rep(c("ProteinA", "ProteinB", "ProteinC"), each = 3),
  sample = rep(c("S1", "S2", "S3"), times = 3),
  intensity = c(100, 120, 90, 80, 95, 85, 200, 220, 210)
)

# 各タンパク質ごとの Zスコア正規化
proteomics_data %>%
  group_by(protein) %>%
  mutate(z_score = (intensity - mean(intensity)) / sd(intensity))