Zenn
👌

【シングルセル解析】[Seurat] FindAllMarkersの結果をUpSet Plot

2023/02/15に公開
R
#
bioinformatics
#
seurat
#
complexheatmap
#
upset
tech

SeuratのDEG検出機能では、「1つのcluster vs 残りの全て」という比較が行われる。
例えば、「CD4 T cell vs 末梢血のその他の細胞」 という構図だと、その他の細胞にはTregやCD8 TなどT細胞も含まれることになる。
CD4 T cellのDEGとしてT細胞マーカーが得られるだろうが、「Treg vs 残りの全て」の時もT細胞マーカーがDEGとして得られるだろう。
このように各クラスターのDEGはuniqueではない場合がある。

各クラスターのDEGの積集合/差集合を見ることで、共通のDEG/独自のDEGを見つけることができる。
良く使用されるのはベン図であるが、ベン図はクラスター数が多くなると、一気に見づらくなる。

UpSet Plotは集合数が多くても見やすいのでぜひ参考にしてほしい。

データの用意

Seuratのチュートリアルページから取得できるPBMC 2700細胞のデータを使用。
https://satijalab.org/seurat/articles/pbmc3k_tutorial.html

デモデータ処理

チュートリアルと同様に処理を進める。デモ用にcluster番号の情報と細胞名の情報をmeta.dataスロットに保存している。

使用するパッケージの読み込み
library(Seurat)
library(dplyr)

matrix.data <- Read10X(data.dir = "pbmc3k_filtered_gene_bc_matrices/filtered_gene_bc_matrices/hg19")
seurat_object <- CreateSeuratObject(counts = matrix.data, project = "pbmc3k", min.cells = 3, min.features = 200)
seurat_object[["percent.mt"]] <- PercentageFeatureSet(seurat_object, pattern = "^MT-")
seurat_object <- subset(seurat_object, subset = nFeature_RNA > 200 & nFeature_RNA < 2500 & percent.mt < 5)
seurat_object <- seurat_object %>% 
  NormalizeData() %>% 
  FindVariableFeatures() %>% 
  ScaleData() %>% 
  RunPCA() %>% 
  RunUMAP(dims = 1:10) %>% 
  FindNeighbors(dims = 1:10) %>% 
  FindClusters(resolution = 0.3)
seurat_object$Cluster <- paste0("cluster_", seurat_object@active.ident)
Idents(seurat_object) <- "Cluster"
DimPlot(seurat_object, reduction = "umap", label = TRUE, pt.size = 0.5) + NoLegend()

Seurat::FindAllMarkers機能で各クラスタの発現変動遺伝子 (DEG)を取得

deg <- Seurat::FindAllMarkers(object = seurat_object,
                              only.pos = T, 
                              logfc.threshold = 0.5)

degオブジェクトはこのようなdata.frame


各クラスターのDEGの数を確認

table(deg$cluster)

クラスターごとにDEGを分割して、リストにまとめる。

deg_list <- split(deg$gene, deg$cluster)

deg_listの中身。要素が文字列ベクトルのリスト


UpSet plot

UpSet plotはComplexHeatmapパッケージから作成することができる。

1) ライブラリ読み込み

library(ComplexHeatmap)

2) リスト内要素の共通項のパターンを集計

mt <- ComplexHeatmap::make_comb_mat(deg_list)

オブジェクトmtの中身。共通項のパターンを2進数で表現 (code)し、その数 (size)が保存されている。

list_to_matrix()機能を使うと、各遺伝子がどのクラスターのDEGだったのかを01に変換している様子が見て取れる。



make_comb_mat()関数のmode=引数は"distinct"がデフォルト。
modeの考え方。https://jokergoo.github.io/ComplexHeatmap-reference/book/upset-plot.html


3) UpSet plotの描画

関数の引数にmtオブジェクトを指定

ComplexHeatmap::UpSet(mt)

UpSet plotのオプション

size順に並び替え

ComplexHeatmap::UpSet(mt, comb_order = order(ComplexHeatmap::comb_size(mt)))

size閾値を設ける

mt2 <- mt[ComplexHeatmap::comb_size(mt)>=50]
ComplexHeatmap::UpSet(mt2, 
                      comb_order = order(ComplexHeatmap::comb_size(mt2)))

cluster表示順を指定

ComplexHeatmap::UpSet(mt, set_order = sort(names(deg_list)))

組み合わせパターンの数を限定する

1つのクラスターにしか検出されていないもののみ
mt2 <- mt[ComplexHeatmap::comb_degree(mt)==1]
ComplexHeatmap::UpSet(mt2)

2つ以上のクラスターに共通のもののみ
mt2 <- mt[ComplexHeatmap::comb_degree(mt)>=2]
ComplexHeatmap::UpSet(mt2)

組み合わせパターンの数に応じて色を変える。

ComplexHeatmap::UpSet(mt, comb_col = c("red","blue","green","black")[ComplexHeatmap::comb_degree(mt)])

UpSet plotはvenn plotで表現するよりも和集合、積集合が見やすくなる。

venn::venn(deg_list, zcolor = "style")

Discussion