**文書の過去の版を表示しています。**
Visium解析引継ぎ資料
解析フローとしては大きく以下の流れで行った。
- Space Ranger(というよりかはSTAR)と互換性を保つようにGTFファイルを編集
- Space Rangerのmkrefコマンドでインデックスを作成
- Space Rangerのcountコマンドでリファレンスへのマッピングとカウント、および組織の検出
- Web Summaryによる概要の確認
- Loupe Browserによる解析
- Seuratによる解析
- その他下流解析
GTFファイルの編集
10X Genomicsによって提供されているヒトやマウス以外の非モデル生物のカスタムリファレンスを使用したい場合は、GTFファイルを編集し、STARとの互換性を保つようにする必要がある。
互換性を保つためには以下の2つの処理が必要(もしかしたら他にも必要かも)
- GTFファイルの9列目のAttribute列に
gene_id
とtranscript_id
が入っている必要がある。 - コメント行を削除
Before
# Predicted genes for sequence number 1 on both strands # start gene g00001 scaffold10|size694159 AUGUSTUS gene 19 26628 0.11 + . ID=pfu_aug170726_10_00001 scaffold10|size694159 AUGUSTUS transcript 19 26628 0.11 + . ID=pfu_aug170726_10_00001.t1;Parent=pfu_aug170726_10_00001 scaffold10|size694159 AUGUSTUS transcription_start_site 19 19 . + . Parent=pfu_aug170726_10_00001.t1 scaffold10|size694159 AUGUSTUS exon 19 284 . + . Parent=pfu_aug170726_10_00001.t1 scaffold10|size694159 AUGUSTUS start_codon 224 226 . + 0 Parent=pfu_aug170726_10_00001.t1 scaffold10|size694159 AUGUSTUS intron 285 6779 1 + . Parent=pfu_aug170726_10_00001.t1 scaffold10|size694159 AUGUSTUS intron 6929 8036 1 + . Parent=pfu_aug170726_10_00001.t1 scaffold10|size694159 AUGUSTUS intron 8154 8308 1 + . Parent=pfu_aug170726_10_00001.t1
After
scaffold10|size694159 AUGUSTUS gene 19 26628 0.11 + . ID=pfu_aug170726_10_00001; transcript_id "pfu_aug170726_10_00001"; gene_id "pfu_aug170726_10_00001"; scaffold10|size694159 AUGUSTUS transcript 19 26628 0.11 + . ID=pfu_aug170726_10_00001.t1;Parent=pfu_aug170726_10_00001; transcript_id "pfu_aug170726_10_00001.t1"; gene_id "pfu_aug170726_10_00001.t1"; scaffold10|size694159 AUGUSTUS transcription_start_site 19 19 . + . Parent=pfu_aug170726_10_00001.t1; transcript_id "pfu_aug170726_10_00001.t1"; gene_id "pfu_aug170726_10_00001.t1"; scaffold10|size694159 AUGUSTUS exon 19 284 . + . Parent=pfu_aug170726_10_00001.t1; transcript_id "pfu_aug170726_10_00001.t1"; gene_id "pfu_aug170726_10_00001.t1"; scaffold10|size694159 AUGUSTUS start_codon 224 226 . + 0 Parent=pfu_aug170726_10_00001.t1; transcript_id "pfu_aug170726_10_00001.t1"; gene_id "pfu_aug170726_10_00001.t1"; scaffold10|size694159 AUGUSTUS intron 285 6779 1 + . Parent=pfu_aug170726_10_00001.t1; transcript_id "pfu_aug170726_10_00001.t1"; gene_id "pfu_aug170726_10_00001.t1"; scaffold10|size694159 AUGUSTUS intron 6929 8036 1 + . Parent=pfu_aug170726_10_00001.t1; transcript_id "pfu_aug170726_10_00001.t1"; gene_id "pfu_aug170726_10_00001.t1"; scaffold10|size694159 AUGUSTUS intron 8154 8308 1 + . Parent=pfu_aug170726_10_00001.t1; transcript_id "pfu_aug170726_10_00001.t1"; gene_id "pfu_aug170726_10_00001.t1";
インデックスの作成
Space Rangerのmkrefコマンドを使用して、リファレンスのインデックスを作成する。
spaceranger mkref --genome=(インデックス名) --fasta=(リファレンスの のFASTAファイル) --genes=(リファレンスのGTFファイル)
実際に使用したコマンドは以下の通り
spaceranger mkref --genome=akoya_ref_20211015 --fasta=raw/pfu_170726_scaffold.gapclosed.fasta --genes=raw/pfu_aug170726.modID_modified_20201020_filtered.gtf
ちなみにSpace Rangerのダウンロードは以下のリンク先から可能である。
Space Ranger Software Download
また最新版でない可能性があるが、伊藤がインストールしたSpace Rangerは/suikou/files/m208/ito.takumi/work/visium/space_ranger/spaceranger-1.0.0/bin/spaceranger にあるのでそれを使用してもよい。
組織検出およびマッピングとリードのカウント
Space Rangerのcountコマンドを用いて、組織検出、リードのマッピングとカウントを行った。
spaceranger count --id=(作成するディレクトリの名前) --transcriptome=(リファレンスインデックスのPATH) --fastqs=(シーケンスデータが格納されたディレクトリのPATH) --sample=(サンプル名、FASTQファイルの接頭辞?おそらくなくてもOK) --image=(組織切片の画像) --slide=(Visium SlideのID、Slideに記載されている) --area=(Slide上の切片を置いた位置) --localcores=(コア数) --localmem=(メモリ)
実際に使用したコマンドは以下のとおり
spaceranger count --id=mantle_20211015 --transcriptome=/home/ito.takumi/work/pinctada/visium_spaceranger3/ref/akoya_ref_20211015/ --fastqs=/home/ito.takumi/work/pinctada/visium_spaceranger3/mantle/input/ --sample=akoya_mantle --image=/home/ito.takumi/work/pinctada/visium_spaceranger3/image/20210718155118_pinctada_mantle.tif --slide=V10T17-119 --area=A1 --localcores=8 --localmem=64
countコマンドの出力としてouts
ディレクトリが作成され、以下のファイルおよびディレクトリがouts
ディレクトリ内に作成される。
analysis → 下流解析の結果 cloupe.cloupe → Loupe Browserの入力ファイルとして使用 filtered_feature_bc_matrix → filtered_feature_bc_matrix.h5 metrics_summary.csv molecule_info.h5 possorted_genome_bam.bam possorted_genome_bam.bam.bai raw_feature_bc_matrix raw_feature_bc_matrix.h5 spatial web_summary.html → Visium実験のサマリーを表示するHTMLファイル
Web Summaryによる実験概要の確認
Space Rangerのcountコマンドの出力であるWeb Summaryをもとに実験概要の確認を行う。
②シーケンス情報の確認(スポット当たりのリード数が50000リード程度あれば十分なシーケンス量を確保できているといえる)
またこれらをグラフで可視化したものや、簡易的なクラスタリング図などを「Analysis」タブから確認することもできる。
Loupe Browserによる下流解析
Loupe Browserを用いることでGUIベースで空間遺伝子発現解析の解析を行うことができる。各遺伝子の空間遺伝子発現分布の取得、組織のクラスタリング、UMAPの描画などが主な機能となる。ただしクラスタリングに関しては細かいパラメータの調整などが不可能であるため、CUIベースのSeuratで行うことを推奨する。一方各遺伝子の空間遺伝子発現分布の取得においてはどちらでもそこまで変わりはないためLoupe Browserで行った方が楽に思われる(伊藤の感想としては、色調のパラメータ調整がLegendの表示がLoupe Browserの方が楽と感じた)。
入力データとしてはSpace Rangerの出力ファイルである、cloupe.cloupe
を使用する。
①「Categories」を選択することでクラスタリングを行う。(Spatialで組織クラスタリング結果図を、UMAPでUMAP図を表示する)
Seuratによる下流解析
Seuratを用いて、各遺伝子の空間遺伝子発現分布や組織のクラスタリング、UMAPの描画、クラスターのマーカー遺伝子の取得、他の組織の空間遺伝子発現解析データとのマージなどを行うことができる。
Seurat解析手順
1.環境構築
①R 4.1.2をインストール (https://cran.r-project.org/index.html)
②R studioをインストール (https://www.rstudio.com/products/rstudio/download/#download)
③Seuratをインストール
> install.packages('Seurat')
④hdf5rをインストール
> install.packages(“hdf5r”)
2.解析
①必要なライブラリをロードする
library(Seurat)\\library(ggplot2)\\library(patchwork)\\library(dplyr)\\
②Space Rangerの解析結果を読み込む
mantle <- Load10X_Spatial("C:/kougaku/pinctada/visium2/Visium2/Result/mantle/outs_20211015", filename = "filtered_feature_bc_matrix.h5", assay = "mantle", slice = "C:/kougaku/pinctada/visium2/Visium2/Result/Image/Visium_20210718/20210718155118_pinctada_mantle.tif", filter.matrix = TRUE,to.upper = FALSE)hinge <- Load10X_Spatial("C:/kougaku/pinctada/visium2/Visium2/Result/hinge/outs_20211022", filename = "filtered_feature_bc_matrix.h5", assay = "hinge", slice = "C:/kougaku/pinctada/visium2/Visium2/Result/Image/20210902105928_pinctada_hinge.tif", filter.matrix = TRUE,to.upper = FALSE)
③バイオリンプロットの作成
plot1 <- VlnPlot(mantle, features = "nCount_mantle", pt.size = 0.1) + NoLegend()plot2 <- VlnPlot(hinge, features = "nCount_hinge", pt.size = 0.1) + NoLegend()plot1 + plot2
④組織内のリードカウント数の分布を表示
plot3 <- SpatialFeaturePlot(mantle, features = "nCount_mantle") + theme(legend.position = "right")plot4 <- SpatialFeaturePlot(hinge, features = "nCount_hinge") + theme(legend.position = "right")plot3 + plot4
⑤カウント数が0のスポットのデータを除去する。
mantle_filtered <- subset(mantle, subset = nCount_mantle > 1)hinge_filtered <- subset(hinge, subset = nCount_hinge > 1)
⑥データの正規化を行う。
mantle_scaled <- SCTransform(mantle_filtered, assay = "mantle", verbose = FALSE)hinge_scaled <- SCTransform(hinge_filtered, assay = "hinge", verbose = FALSE)
⑦mantle, hingeのデータをマージ
akoya.merge <- merge(mantle_scaled, hinge_scaled)
⑧RNA発現データ上での次元削減とクラスタリング
DefaultAssay(akoya.merge) <- "SCT"VariableFeatures(akoya.merge) <- c(VariableFeatures(mantle_scaled), VariableFeatures(hinge_scaled))akoya.merge <- RunPCA(akoya.merge, verbose = FALSE)akoya.merge <- FindNeighbors(akoya.merge, dims = 1:15)akoya.merge <- FindClusters(akoya.merge, verbose = FALSE)akoya.merge <- RunUMAP(akoya.merge, dims = 1:30)
⑨mergeデータのUMAPを作成
DimPlot(akoya.merge, reduction = "umap")
⑩merge後の組織のクラスタリング図を作成(2組織とも表示)
SpatialDimPlot(akoya.merge)
⑪merge後の組織のクラスタリング図を作成(1組織ずつラベルと一緒に表示)
SpatialDimPlot(akoya.merge,label = TRUE, combine = FALSE)
⑫各クラスターのポジティブマーカーを抽出(閾値は適宜設定)
pbmc.markers <- FindAllMarkers(akoya.merge, only.pos = TRUE, min.pct = 0.25, logfc.threshold = 0.25)View(pbmc.markers)write.csv(pbmc.markers, "cluster_markers.csv")
⑬各遺伝子の空間遺伝子発現を表示(MSI60=pfu-aug170726-71-03468.t1の場合)(_が-に変換されていることに注意)、本日は行いませんでしたが、Loupe Browserで行っているような解析をSeuratでも行う事ができます。
SpatialFeaturePlot(akoya.merge, features = "pfu-aug170726-71-03468.t1")
⑭各遺伝子が発現しているスポットをUMAP上で表示(以下は複数遺伝子を同時に表示したい場合)
FeaturePlot(akoya.merge, features = c("pfu-aug170726-71-03468.t1", "pfu-aug170726-869-31402.t1", "pfu-aug170726-1189-31713.t1", "pfu-aug170726-1547-25368.t1", "pfu-aug170726-3143-12855.t1", "pfu-aug170726-102-06870.t1"), pt.size = 0.2, ncol = 3)
※2021/12/10追記:特定のクラスタのみを抽出し、そのデータで再度PCA解析を行う。 ⑮クラスター7,11のデータを抽出する。
hinge_selected <- subset(hinge_scaled, idents = c(7,11))
再度PCA解析などを行う
hinge_selected <- RunPCA(hinge_selected, verbose = FALSE)hinge_selected <- FindNeighbors(hinge_selected, dims = 1:15)hinge_selected <- FindClusters(hinge_selected, verbose = FALSE)hinge_selected <- RunUMAP(hinge_selected, dims = 1:30)
グラフを描画(crop=TRUEにしないと組織の一部分しか表示されないので注意)
SpatialDimPlot(hinge_selected,label = TRUE,crop = FALSE)