往期系列

亨貝格實驗室單細胞轉錄組數據分析（一）

亨貝格實驗室單細胞轉錄組數據分析（二）

亨貝格實驗室單細胞轉錄組數據分析（三）

亨貝格實驗室單細胞轉錄組數據分析（四）

亨貝格實驗室單細胞轉錄組數據分析（五）

亨貝格實驗室單細胞轉錄組數據分析（六）

Hemberg-lab單細胞轉錄組數據分析（七） - 導入10X和SmartSeq2數據Tabula Muris

收藏|北大生信平臺 “單細胞分析，染色質分析” 視頻和PPT分享

該如何自學入門生物信息學

生物信息之程序學習

收藏|你想要的生信學習系列教程 - 寶典在手，生信無憂

表達矩陣質控

UMI表達定量（UMI）

UMI表達定量簡介

基因定量后會整理成一個行為基因（或轉錄本）細胞列為的表達矩陣。雖然前面做了原始數據質控和測序數據質控移除了一部分從讀取數層面就不合格的細胞，還需要進一步根據表達矩陣移除其它類型低質量細胞。如果未能識別并移除低質量細胞會混淆下游分析中的有意義的生物信息。

鑒于scRNASeq還沒有標準方法，后續用到的質控質量值的標準會因為實驗方法不同而有很大變化。因此，執行質控時，我們是通過數據集內部比較找到異常細胞，而不是依賴于其它獨立的質量標準。因此比較不同的建庫方法獲得的不同數據集時需要格外注意。

董數據集

Yoav Gilad在Yoav Gilad?我們使用芝加哥大學實驗室的3個不同來源的誘導多能性干細胞（iPSC）的數據集（http://jdblischak.github.io/singleCellSeq/analysis/）。細胞分選采用Fluidigm C1微流控臺，使用同時UMIs狀語從句：ERCC spike in進行質控為了保證柯林斯重復性，數據是2016年3月15生成的原始數據的拷貝，于存儲tung文件夾數下。

library(SingleCellExperiment) library(scater) options(stringsAsFactors = FALSE)

讀入數據和注釋：

molecules <- read.table("tung/molecules.txt", sep = "\t") anno <- read.table("tung/annotation.txt", sep = "\t", header = TRUE)

查看讀入的數據：

head(molecules[ , 1:3])

結果

## NA19098.r1.A01 NA19098.r1.A02 NA19098.r1.A03 ## ENSG00000237683 0 0 0 ## ENSG00000187634 0 0 0 ## ENSG00000188976 3 6 1 ## ENSG00000187961 0 0 0 ## ENSG00000187583 0 0 0 ## ENSG00000187642 0 0 0

查看注釋數據

head(anno)

結果

## individual replicate well batch sample_id ## 1 NA19098 r1 A01 NA19098.r1 NA19098.r1.A01 ## 2 NA19098 r1 A02 NA19098.r1 NA19098.r1.A02 ## 3 NA19098 r1 A03 NA19098.r1 NA19098.r1.A03 ## 4 NA19098 r1 A04 NA19098.r1 NA19098.r1.A04 ## 5 NA19098 r1 A05 NA19098.r1 NA19098.r1.A05 ## 6 NA19098 r1 A06 NA19098.r1 NA19098.r1.A06

數據包含3個個體，3次重復和9個批次。

通過使用SingleCellExperiment（SCE）和scater包標準化分析過程。首先構建SCE對象：

umi <- SingleCellExperiment(assays = list(counts = as.matrix(molecules)), colData = anno )

移除不在任何細胞表達的基因：

keep_feature <- rowSums(counts(umi) > 0) > 0 umi <- umi[keep_feature, ]

定義對照特征基因集：ERCC spike-ins和線粒體基因（作者提供，見http://jdblischak.github.io/singleCellSeq/analysis/qc-filter-ipsc.html）：

isSpike(umi, "ERCC") <- grepl("^ERCC-", rownames(umi)) isSpike(umi, "MT") <- rownames(umi) %in% c("ENSG00000198899", "ENSG00000198727", "ENSG00000198888","ENSG00000198886", "ENSG00000212907", "ENSG00000198786","ENSG00000198695", "ENSG00000198712", "ENSG00000198804","ENSG00000198763", "ENSG00000228253", "ENSG00000198938","ENSG00000198840")

計算質量矩陣：

umi <- calculateQCMetrics(umi,feature_controls = list(ERCC = isSpike(umi, "ERCC"), MT = isSpike(umi, "MT")) ) str(umi)

SingleCellExperiment對象結構。str（結構）是一個很好的工具，可以用來查看數據的結構組成。（RStudio中的View對于較大的對象會給出更好的展現方式。）

如下面的展示SingleCellExperiment有10個數據槽（slots），使用@索引。比如想獲得降維后的結果，使用umi@reduceDims，使用umi@colData可以獲取質控信息，都有哪些質控變量，進一步的使用umi@colData@listData$total_features_by_counts可以獲得每個細胞的基因數。umi@rowRanges@elementMetadata@listData基因的屬性信息。

（在Rstudio中有個便利的地方，輸入變量名再輸@，$可以彈出對應的子屬性，方便快速輸入。）

# 獲取原始讀入的表達矩陣head(umi@assays$data@listData$counts)[,1:3]NA19098.r1.A01 NA19098.r1.A02 NA19098.r1.A03 ENSG00000237683 0 0 0 ENSG00000187634 0 0 0 ENSG00000188976 3 6 1 ENSG00000187961 0 0 0 ENSG00000187583 0 0 0 ENSG00000187642 0 0 0str(umi)## Formal class 'SingleCellExperiment' [package "SingleCellExperiment"] with 10 slots## ..@ int_elementMetadata:Formal class 'DataFrame' [package "S4Vectors"] with 6 slots ## .. .. ..@ rownames : NULL ## .. .. ..@ nrows : int 18726 ## .. .. ..@ listData :List of 3 ## .. .. .. ..$ is_spike_ERCC: logi [1:18726] FALSE FALSE FALSE FALSE FALSE FALSE ... ## .. .. .. ..$ is_spike : logi [1:18726] FALSE FALSE FALSE FALSE FALSE FALSE ... ## .. .. .. ..$ is_spike_MT : logi [1:18726] FALSE FALSE FALSE FALSE FALSE FALSE ... ## .. .. ..@ elementType : chr "ANY" ## .. .. ..@ elementMetadata: NULL ## .. .. ..@ metadata : list() ## ..@ int_colData :Formal class 'DataFrame' [package "S4Vectors"] with 6 slots ## .. .. ..@ rownames : NULL ## .. .. ..@ nrows : int 864 ## .. .. ..@ listData : Named list() ## .. .. ..@ elementType : chr "ANY" ## .. .. ..@ elementMetadata: NULL ## .. .. ..@ metadata : list() ## ..@ int_metadata :List of 3 ## .. ..$ version :Classes 'package_version', 'numeric_version' hidden list of 1 ## .. .. ..$ : int [1:3] 1 4 1 ## .. ..$ spike_names : chr [1:2] "ERCC" "MT" ## .. ..$ size_factor_names: chr(0) ## ..@ reducedDims :Formal class 'SimpleList' [package "S4Vectors"] with 4 slots ## .. .. ..@ listData : Named list() ## .. .. ..@ elementType : chr "ANY" ## .. .. ..@ elementMetadata: NULL ## .. .. ..@ metadata : list() ## ..@ rowRanges :Formal class 'CompressedGRangesList' [package "GenomicRanges"] with 5 slots ## .. .. ..@ unlistData :Formal class 'GRanges' [package "GenomicRanges"] with 7 slots ## .. .. .. .. ..@ seqnames :Formal class 'Rle' [package "S4Vectors"] with 4 slots ## .. .. .. .. .. .. ..@ values : Factor w/ 0 levels: ## .. .. .. .. .. .. ..@ lengths : int(0) ## .. .. .. .. .. .. ..@ elementMetadata: NULL ## .. .. .. .. .. .. ..@ metadata : list() ## .. .. .. .. ..@ ranges :Formal class 'IRanges' [package "IRanges"] with 6 slots ## .. .. .. .. .. .. ..@ start : int(0) ## .. .. .. .. .. .. ..@ width : int(0) ## .. .. .. .. .. .. ..@ NAMES : NULL ## .. .. .. .. .. .. ..@ elementType : chr "ANY" ## .. .. .. .. .. .. ..@ elementMetadata: NULL ## .. .. .. .. .. .. ..@ metadata : list() ## .. .. .. .. ..@ strand :Formal class 'Rle' [package "S4Vectors"] with 4 slots ## .. .. .. .. .. .. ..@ values : Factor w/ 3 levels "+","-","*": ## .. .. .. .. .. .. ..@ lengths : int(0) ## .. .. .. .. .. .. ..@ elementMetadata: NULL ## .. .. .. .. .. .. ..@ metadata : list() ## .. .. .. .. ..@ seqinfo :Formal class 'Seqinfo' [package "GenomeInfoDb"] with 4 slots ## .. .. .. .. .. .. ..@ seqnames : chr(0) ## .. .. .. .. .. .. ..@ seqlengths : int(0) ## .. .. .. .. .. .. ..@ is_circular: logi(0) ## .. .. .. .. .. .. ..@ genome : chr(0) ## .. .. .. .. ..@ elementMetadata:Formal class 'DataFrame' [package "S4Vectors"] with 6 slots ## .. .. .. .. .. .. ..@ rownames : NULL ## .. .. .. .. .. .. ..@ nrows : int 0 ## .. .. .. .. .. .. ..@ listData : Named list() ## .. .. .. .. .. .. ..@ elementType : chr "ANY" ## .. .. .. .. .. .. ..@ elementMetadata: NULL ## .. .. .. .. .. .. ..@ metadata : list() ## .. .. .. .. ..@ elementType : chr "ANY" ## .. .. .. .. ..@ metadata : list() ## .. .. ..@ elementMetadata:Formal class 'DataFrame' [package "S4Vectors"] with 6 slots ## .. .. .. .. ..@ rownames : NULL ## .. .. .. .. ..@ nrows : int 18726 ## .. .. .. .. ..@ listData :List of 9 ## .. .. .. .. .. ..$ is_feature_control : logi [1:18726] FALSE FALSE FALSE FALSE FALSE FALSE ... ## .. .. .. .. .. ..$ is_feature_control_ERCC: logi [1:18726] FALSE FALSE FALSE FALSE FALSE FALSE ... ## .. .. .. .. .. ..$ is_feature_control_MT : logi [1:18726] FALSE FALSE FALSE FALSE FALSE FALSE ... ## .. .. .. .. .. ..$ mean_counts : num [1:18726] 0.2824 0.0301 2.6389 0.2384 0.0116 ... ## .. .. .. .. .. ..$ log10_mean_counts : num [1:18726] 0.108 0.0129 0.561 0.0929 0.005 ... ## .. .. .. .. .. ..$ n_cells_by_counts : int [1:18726] 201 24 728 178 10 11 20 632 798 19 ... ## .. .. .. .. .. ..$ pct_dropout_by_counts : num [1:18726] 76.7 97.2 15.7 79.4 98.8 ... ## .. .. .. .. .. ..$ total_counts : int [1:18726] 244 26 2280 206 10 11 21 1638 2873 19 ... ## .. .. .. .. .. ..$ log10_total_counts : num [1:18726] 2.39 1.43 3.36 2.32 1.04 ... ## .. .. .. .. ..@ elementType : chr "ANY" ## .. .. .. .. ..@ elementMetadata: NULL ## .. .. .. .. ..@ metadata : list() ## .. .. ..@ elementType : chr "GRanges" ## .. .. ..@ metadata : list() ## .. .. ..@ partitioning :Formal class 'PartitioningByEnd' [package "IRanges"] with 5 slots ## .. .. .. .. ..@ end : int [1:18726] 0 0 0 0 0 0 0 0 0 0 ... ## .. .. .. .. ..@ NAMES : chr [1:18726] "ENSG00000237683" "ENSG00000187634" "ENSG00000188976" "ENSG00000187961" ... ## .. .. .. .. ..@ elementType : chr "ANY" ## .. .. .. .. ..@ elementMetadata: NULL ## .. .. .. .. ..@ metadata : list() ## ..@ colData :Formal class 'DataFrame' [package "S4Vectors"] with 6 slots ## .. .. ..@ rownames : chr [1:864] "NA19098.r1.A01" "NA19098.r1.A02" "NA19098.r1.A03" "NA19098.r1.A04" ... ## .. .. ..@ nrows : int 864 ## .. .. ..@ listData :List of 50 ## .. .. .. ..$ individual : chr [1:864] "NA19098" "NA19098" "NA19098" "NA19098" ... ## .. .. .. ..$ replicate : chr [1:864] "r1" "r1" "r1" "r1" ... ## .. .. .. ..$ well : chr [1:864] "A01" "A02" "A03" "A04" ... ## .. .. .. ..$ batch : chr [1:864] "NA19098.r1" "NA19098.r1" "NA19098.r1" "NA19098.r1" ... ## .. .. .. ..$ sample_id : chr [1:864] "NA19098.r1.A01" "NA19098.r1.A02" "NA19098.r1.A03" "NA19098.r1.A04" ... ## .. .. .. ..$ is_cell_control : logi [1:864] FALSE FALSE FALSE FALSE FALSE FALSE ... ## .. .. .. ..$ total_features_by_counts : int [1:864] 8368 8234 7289 7985 8619 8659 8054 9429 8243 9407 ... ## .. .. .. ..$ log10_total_features_by_counts : num [1:864] 3.92 3.92 3.86 3.9 3.94 ... ## .. .. .. ..$ total_counts : int [1:864] 63322 63976 43630 53922 70887 68051 58069 89082 60053 87489 ... ## .. .. .. ..$ log10_total_counts : num [1:864] 4.8 4.81 4.64 4.73 4.85 ... ## .. .. .. ..$ pct_counts_in_top_50_features : num [1:864] 17.6 16.5 18.2 16.8 16.2 ... ## .. .. .. ..$ pct_counts_in_top_100_features : num [1:864] 25.2 24.4 26.1 24.4 24.1 ... ## .. .. .. ..$ pct_counts_in_top_200_features : num [1:864] 34.9 34.5 36 34.1 34 ... ## .. .. .. ..$ pct_counts_in_top_500_features : num [1:864] 50.2 50.5 51.5 49.8 50 ... ## .. .. .. ..$ total_features_by_counts_endogenous : int [1:864] 8324 8190 7248 7942 8573 8616 8009 9384 8201 9361 ... ## .. .. .. ..$ log10_total_features_by_counts_endogenous : num [1:864] 3.92 3.91 3.86 3.9 3.93 ... ## .. .. .. ..$ total_counts_endogenous : int [1:864] 57252 58967 39420 49076 65244 63508 53456 83920 56280 82287 ... ## .. .. .. ..$ log10_total_counts_endogenous : num [1:864] 4.76 4.77 4.6 4.69 4.81 ... ## .. .. .. ..$ pct_counts_endogenous : num [1:864] 90.4 92.2 90.4 91 92 ... ## .. .. .. ..$ pct_counts_in_top_50_features_endogenous : num [1:864] 12.5 13 13.1 12.3 12.5 ... ## .. .. .. ..$ pct_counts_in_top_100_features_endogenous : num [1:864] 20 20.7 21 19.8 20.2 ... ## .. .. .. ..$ pct_counts_in_top_200_features_endogenous : num [1:864] 29.8 30.8 31.1 29.5 30.1 ... ## .. .. .. ..$ pct_counts_in_top_500_features_endogenous : num [1:864] 45.9 47.3 47.4 45.9 46.7 ... ## .. .. .. ..$ total_features_by_counts_feature_control : int [1:864] 44 44 41 43 46 43 45 45 42 46 ... ## .. .. .. ..$ log10_total_features_by_counts_feature_control: num [1:864] 1.65 1.65 1.62 1.64 1.67 ... ## .. .. .. ..$ total_counts_feature_control : int [1:864] 6070 5009 4210 4846 5643 4543 4613 5162 3773 5202 ... ## .. .. .. ..$ log10_total_counts_feature_control : num [1:864] 3.78 3.7 3.62 3.69 3.75 ... ## .. .. .. ..$ pct_counts_feature_control : num [1:864] 9.59 7.83 9.65 8.99 7.96 ... ## .. .. .. ..$ pct_counts_in_top_50_features_feature_control : num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ pct_counts_in_top_100_features_feature_control: num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ pct_counts_in_top_200_features_feature_control: num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ pct_counts_in_top_500_features_feature_control: num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ total_features_by_counts_ERCC : int [1:864] 31 31 28 30 33 30 32 32 29 33 ... ## .. .. .. ..$ log10_total_features_by_counts_ERCC : num [1:864] 1.51 1.51 1.46 1.49 1.53 ... ## .. .. .. ..$ total_counts_ERCC : int [1:864] 1187 1277 1121 1240 1262 1308 1203 1217 1217 1339 ... ## .. .. .. ..$ log10_total_counts_ERCC : num [1:864] 3.07 3.11 3.05 3.09 3.1 ... ## .. .. .. ..$ pct_counts_ERCC : num [1:864] 1.87 2 2.57 2.3 1.78 ... ## .. .. .. ..$ pct_counts_in_top_50_features_ERCC : num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ pct_counts_in_top_100_features_ERCC : num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ pct_counts_in_top_200_features_ERCC : num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ pct_counts_in_top_500_features_ERCC : num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ total_features_by_counts_MT : int [1:864] 13 13 13 13 13 13 13 13 13 13 ... ## .. .. .. ..$ log10_total_features_by_counts_MT : num [1:864] 1.15 1.15 1.15 1.15 1.15 ... ## .. .. .. ..$ total_counts_MT : int [1:864] 4883 3732 3089 3606 4381 3235 3410 3945 2556 3863 ... ## .. .. .. ..$ log10_total_counts_MT : num [1:864] 3.69 3.57 3.49 3.56 3.64 ... ## .. .. .. ..$ pct_counts_MT : num [1:864] 7.71 5.83 7.08 6.69 6.18 ... ## .. .. .. ..$ pct_counts_in_top_50_features_MT : num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ pct_counts_in_top_100_features_MT : num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ pct_counts_in_top_200_features_MT : num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. .. ..$ pct_counts_in_top_500_features_MT : num [1:864] 100 100 100 100 100 100 100 100 100 100 ... ## .. .. ..@ elementType : chr "ANY" ## .. .. ..@ elementMetadata: NULL ## .. .. ..@ metadata : list() ## ..@ assays :Reference class 'ShallowSimpleListAssays' [package "SummarizedExperiment"] with 1 field ## .. ..$ data: NULL ## .. ..and 14 methods. ## ..@ NAMES : NULL ## ..@ elementMetadata :Formal class 'DataFrame' [package "S4Vectors"] with 6 slots ## .. .. ..@ rownames : NULL ## .. .. ..@ nrows : int 18726 ## .. .. ..@ listData : Named list() ## .. .. ..@ elementType : chr "ANY" ## .. .. ..@ elementMetadata: NULL ## .. .. ..@ metadata : list() ## ..@ metadata : list()

總結

以上是生活随笔為你收集整理的Hemberg-lab单细胞转录组数据分析（八） - Scater包输入导入和存储的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。