挖掘公共單細胞數(shù)據(jù)集時衩匣,會遇到常見各種單細胞測序數(shù)據(jù)格式。現(xiàn)總結如下粥航,方便自己日后調(diào)用琅捏,以創(chuàng)建Seurat對象
(1)barcodes.tsv.gz
、features.tsv.gz
递雀、matrix.mtx.gz
(2)表達矩陣
(3)h5
(4)h5ad
格式一:barcodes.tsv.gz
柄延、features.tsv.gz
、matrix.mtx.gz
【☆】
- 這是cellranger上游比對分析產(chǎn)生的3個文件,分別代表細胞標簽(barcode)搜吧、基因ID(feature)市俊、表達數(shù)據(jù)(matrix)
- 一般先使用
read10X()
對這三個文件進行整合,得到行為基因滤奈、列為細胞的表達矩陣(為稀疏矩陣dgCMatrix格式摆昧,節(jié)約內(nèi)存);然后再配合CreateSeuratObject()
函數(shù)創(chuàng)建Seurat對象 - 示例數(shù)據(jù)集:GSE166635蜒程,創(chuàng)建代碼如下----
dir="./data/HCC2/filtered_feature_bc_matrix/"
list.files(dir)
#[1] "barcodes.tsv.gz" "features.tsv.gz" "matrix.mtx.gz"
counts <- Read10X(data.dir = dir)
class(counts)
#[1] "dgCMatrix"
#attr(,"package")
#[1] "Matrix"
scRNA <- CreateSeuratObject(counts = counts)
scRNA
#An object of class Seurat
#33694 features across 9112 samples within 1 assay
#Active assay: RNA (33694 features, 0 variable features)
- 如上
Read10X()
函數(shù)接受的參數(shù)為目錄名绅你,該目錄包含了所需的三個配套文件;值得注意的是三個文件名只能分別是barcodes.tsv.gz
昭躺、features.tsv.gz
忌锯、matrix.mtx.gz
,然后read10X
函數(shù)可以自動加載领炫。如上截圖那樣就是需要修改的~
關于barcodes.tsv.gz
偶垮、features.tsv.gz
、matrix.mtx.gz
三個文件的格式與內(nèi)容
- 一般來說直接使用
read10X()
不會出現(xiàn)什么問題帝洪,但今天遇到GSE148192數(shù)據(jù)集時似舵,出現(xiàn)了報錯~~
dir = "./GSE148192_RAW/GSM4462451/"
list.files(dir)
#[1] "barcodes.tsv.gz" "features.tsv.gz" "matrix.mtx.gz"
counts = Read10X(dir)
#Error in dimnamesGets(x, value) :
# invalid dimnames given for “dgTMatrix” object
- 所以這個GSE ID提供的數(shù)據(jù)格式可能是有點問題,接下來就通過對比GSE166635的GSM5076750(可以正常讀入)與GSE148192的GSM4462451(讀入失敗)葱峡,探索下這三個文件的格式
(1)barcodes.tsv.gz
-
GSM5076750的格式:如下看出就簡單的一列啄枕,為細胞的barcode標簽信息
-
GSM4462451的格式:如下看出,區(qū)別在于多了行名族沃,以及三列細胞注釋信息
(2)features.tsv.gz
-
GSM5076750的格式:如下可以看出均為基因的注釋信息,前兩列為基因ID
-
GSM4462451的格式:如下看出泌参,區(qū)別在于同樣多了行名脆淹,以及額外兩列信息
(3)matrix.mtx.gz
- GSM5076750的格式:如下(前三行為注釋信息,其中第三行為total number genes沽一、cells盖溺、counts),結合上述細胞標簽與基因名信息铣缠,知道了前兩列分別為基因和細胞的索引烘嘱,第三列為表達信息。
利用這種方式實現(xiàn)了高效的儲存數(shù)據(jù)(值得借鑒學習)蝗蛙。以第四行為例:表示barcodes.tsv.gz
文件里第一個細胞的features.tsv.gz
第33665個基因的counts數(shù)為22蝇庭。
-
GSM4462451的格式:如下看出,區(qū)別有兩點:第一列為細胞索引捡硅、第二列為基因索引哮内,并且第3列是非整型數(shù)據(jù)。
經(jīng)過一番探索壮韭,將GSM4462451的
barcodes.tsv.gz
北发、features.tsv.gz
行名刪除纹因;matrix.mtx.gz
的第一列與第二列調(diào)換,第三列改為整型后琳拨,read10X()
便可以順利都成功瞭恰。我認為GSM4462451這幾個文件應該是作者自己制作的,吐槽一下~~狱庇。不過了解了一番這三個文件的格式也是有所收獲惊畏。
格式二:直接提供表達矩陣
- 這種是最方便的,直接創(chuàng)建Seurat即可
- 示例數(shù)據(jù):GSE144320
scRNA <- CreateSeuratObject(counts = counts)
scRNA
格式三:h5格式文件
- 使用
Read10X_h5()
函數(shù)僵井,讀入表達矩陣陕截,在創(chuàng)建Seurat對象 - 示例數(shù)據(jù):GSE138433
sce <- Read10X_h5(filename = GSM4107899_LH16.3814_raw_gene_bc_matrices_h5.h5")
sce <- CreateSeuratObject(counts = sce)
格式四:h5ad格式
- 需要安裝,使用
SeuratDisk
包的兩個函數(shù)批什; - 先將后
h5ad
格式轉(zhuǎn)換為h5seurat
格式农曲,再使用LoadH5Seurat()
函數(shù)讀取Seurat對象。 - 示例數(shù)據(jù)集:GSE153643
#remotes::install_github("mojaveazure/seurat-disk")
library(SeuratDisk)
Convert("GSE153643_RAW/GSM4648565_liver_raw_counts.h5ad", "h5seurat",
overwrite = TRUE,assay = "RNA")
scRNA <- LoadH5Seurat("GSE153643_RAW/GSM4648565_liver_raw_counts.h5seurat")
#注意一下驻债,我之前載入時乳规,表達矩陣被轉(zhuǎn)置了,需要處理一下~
以上是我目前了解到的針對不同數(shù)據(jù)來源合呐,創(chuàng)建Seurat對象的幾種方式暮的。如遇新的方法,會繼續(xù)補充~~