從讀取數(shù)據(jù)開始。
Read10X ()函數(shù) 讀取 cellranger pipeline輸出的10X單細胞數(shù)據(jù)，返回一個(UMI)count矩陣。此矩陣中行（row）代表基因，列（col）代表細胞
library(dplyr)
library(Seurat)
library(patchwork)
# 加載 PBMC dataset
pbmc.data <- Read10X(data.dir = "pbmc3k_filtered_gene_bc_matrices/filtered_gene_bc_matrices/hg19/")
接下來我們使用 count 矩陣來創(chuàng)建一個 Seurat 對象。
該對象充當一個容器荷并，它包含單細胞數(shù)據(jù)(如count矩陣)和分析結(jié)果(如 PCA 或聚類)。
比如,count matrix儲存在pbmc[["RNA"]]@counts.
#使用raw data(non-normalized data)創(chuàng)建Seurat 對象
pbmc <- CreateSeuratObject(counts = pbmc.data, project = "pbmc3k", min.cells = 3, min.features = 200)
pbmc
Seurat 的 scRNA-seq 數(shù)據(jù)的標準預(yù)處理流程青扔。
包括 QC 源织、數(shù)據(jù)標準化、縮放和高度可變基因的檢測微猖。
1. 質(zhì)量控制（QC）及選擇細胞進行下游分析
Seurat允許你簡單的過濾一下你的細胞谈息，質(zhì)控的一般指標包括：
(1)在每個細胞里檢測到的基因
-->低質(zhì)量細胞或者空的droplets只有非常少量的基因
-->細胞doublets 或者 multiplets 有非常高的基因count數(shù)
(2)在一個細胞內(nèi)檢測到的分子總數(shù)
(3)讀取到線粒體基因組的比例
-->低質(zhì)量/死細胞通常有很高的線粒體污染
-->我們使用PercentageFeatureSet()函數(shù)計算線粒體QC指標 （計算來自一系列特征的計數(shù)百分比）
-->使用所有以MT-開頭的基因作為線粒體基因
# The "[[ ]]"operator can add columns to object metadata. This is a great place to stash QC stats
pbmc[["percent.mt"]] <- PercentageFeatureSet(pbmc, pattern = "^MT-")
head(pbmc@meta.data, 5)
我們還可以可視化QC指標，并用它們來過濾細胞：
（1）將unique基因count數(shù)超過2500凛剥，或者小于200的細胞過濾掉侠仇。
（2）把線粒體count數(shù)占5%以上的細胞過濾掉
# Visualize QC metrics as a violin plot
VlnPlot(pbmc, features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol = 3)
pbmc <- subset(pbmc, subset = nFeature_RNA > 200 & nFeature_RNA < 2500 & percent.mt < 5)
（畫小提琴一直報錯，暫未解決）
我們還可以用FeatureScatter函數(shù)來可視化特征-特征之間的關(guān)系，可以使用Seurat對象里的任何東西逻炊，如對象中的列互亮、PC分數(shù)等。
plot1 <- FeatureScatter(pbmc, feature1 = "nCount_RNA", feature2 = "percent.mt")
plot2 <- FeatureScatter(pbmc, feature1 = "nCount_RNA", feature2 = "nFeature_RNA")
plot1 + plot2
2. 數(shù)據(jù)標準化
在去除掉不想要的細胞后余素，就可以標準化數(shù)據(jù)了豹休。在默認情況下，我們使用global-scaling標準化方法溺森，稱為“LogNormalize”慕爬，這種方法是利用總的表達量對每個細胞里的基因表達值進行標準化，乘以一個scale factor（默認值是10000）屏积，再用log轉(zhuǎn)換一下。標準化后的數(shù)據(jù)存放在pbmc[["RNA"]]@data里磅甩。
pbmc <- NormalizeData(pbmc, normalization.method = "LogNormalize", scale.factor = 10000)
3. 鑒定高度變化的基因
接下來計算數(shù)據(jù)集中表現(xiàn)出細胞間高度差異的特征子集(即炊林，它們在某些細胞中高表達，而在其他細胞中低表達)卷要。在下游分析中關(guān)注這些基因有助于突出單細胞數(shù)據(jù)集中的生物信號渣聚。
pbmc <- FindVariableFeatures(pbmc, selection.method = "vst", nfeatures = 2000)
#識別出10個差異最大的基因
top10 <- head(VariableFeatures(pbmc), 10)
# 帶標簽或不帶標簽作圖展示差異基因
plot1 <- VariableFeaturePlot(pbmc)
plot2 <- LabelPoints(plot = plot1, points = top10, repel = TRUE)
plot1 + plot2
