ixxmu
commented
1 year ago 读取不同格式的单细胞转录组数据及遇到问题的解决办法 by 生信菜鸟团
几种不同格式单细胞转录组数据的读取方式，以及最近在读取单细胞转录组数据的时候遇到的问题。
读取h5格式的单细胞文件
以下面这个数据集为例，h5格式的单细胞文件 GSE175687
GEO Accession viewer (nih.gov)
批量导入数据：
rm(list=ls())
options(stringsAsFactors = F)
library(Seurat)
library(dplyr)
library(hdf5r)
library(data.table)
dir='GSE175687_raw/'
samples=list.files( dir )
samples
sceList = lapply(samples,function(pro){
print(pro)
sce =CreateSeuratObject(counts = Read10X_h5( file.path(dir,pro)) ,
project = gsub('_filtered_feature_bc_matrix.h5','',gsub('^GSM[0-9]*_','',pro) ) ,
min.cells = 5,
min.features = 300 )
return(sce)
})
sceList
samples
整合数据
sce.all=merge(x=sceList[[1]],
y=sceList[ -1 ],
add.cell.ids = gsub('_filtered_feature_bc_matrix.h5','',gsub('^GSM[0-9]*_','',samples)))
读取10X格式的单细胞文件
文件输入格式
10X格式文件结构：
├── [ 160] Sample1
│ ├── [2.2M] barcodes.tsv.gz
│ ├── [221K] features.tsv.gz
│ └── [4.5M] matrix.mtx.gz
├── [ 160] Sample2
│ ├── [2.2M] barcodes.tsv.gz
│ ├── [221K] features.tsv.gz
│ └── [4.0M] matrix.mtx.gz
可以看到，每个10X的样品都是一个独立的文件夹，每个文件夹里面都是 matrix.mtx, genes.tsv (or features.tsv), and barcodes.tsv 这样的3个文件，非常有规律。
library(data.table)
sce =CreateSeuratObject(counts = Read10X("./GSE135710_raw/") ,
#project = gsub('.gz','',gsub('^GSM[0-9]*_','') ) ,
min.cells = 5,
min.features = 300 )
library(stringr)
head(rownames(sce@meta.data))
phe=str_split(rownames(sce@meta.data),'-',simplify = T)
head(phe)
tail(phe)
table(phe[,2])
sce@meta.data$orig.ident=paste0('p',phe[,2])
table(sce@meta.data$orig.ident)
fread 读取txt.gz格式文件
数据集：GSE183625
GEO Accession viewer (nih.gov)
导入数据
rm(list=ls())
options(stringsAsFactors = F)
library(Seurat)
library(ggplot2)
library(clustree)
library(cowplot)
library(dplyr)
getwd()
library(data.table)
dir='GSE183625_RAW/'
samples=list.files( dir )
samples
library(data.table)
sceList = lapply(samples,function(pro){
print(pro)
ct=fread(file.path( dir ,pro),data.table = F)
ct[1:4,1:4]
rownames(ct)=ct[,1]
ct=ct[,-1]
sce=CreateSeuratObject(counts = ct ,
project = gsub('.expression_matrix.txt.gz','',gsub('^GSM[0-9]*_','',pro) ) ,
min.cells = 5,
min.features = 300)
return(sce)
})
sceList
samples
整合数据
sce.all=merge(x=sceList[[1]],
y=sceList[ -1 ],
add.cell.ids = gsub('.expression_matrix.txt.gz','',gsub('^GSM[0-9]*_','',samples) ) )
fread读取单细胞测序数据csv格式
数据集：GSE129516
GEO Accession viewer (nih.gov)
导入数据：
rm(list=ls())
options(stringsAsFactors = F)
library(Seurat)
library(ggplot2)
library(clustree)
library(cowplot)
library(dplyr)
library(data.table)
dir='GSE129516_RAW/'
samples=list.files( dir )
samples
sceList = lapply(samples,function(pro){
print(pro)
ct=fread(file.path( dir ,pro),data.table = F)
ct[1:4,1:4]
rownames(ct)=ct[,1]
ct=ct[,-1]
sce =CreateSeuratObject(counts = ct ,
project = gsub('_filtered_gene_bc_matrices.csv.gz','',gsub('^GSM[0-9]*_','',pro) ) ,
min.cells = 5,
min.features = 300 )
return(sce)
})
sceList
samples
整合数据：
sce.all=merge(x=sceList[[1]],
y=sceList[ -1 ],
add.cell.ids = gsub('_filtered_gene_bc_matrices.csv.gz','',gsub('^GSM[0-9]*_','',samples) ) )
最近在读取单细胞转录组数据的时候遇到的问题
同样都是10X数据，在文件读取的时候也是会出现问题的，有的GEO上传上去的数据并不标准，标准的文件如下。
最近遇到的数据集是这样的，需要把第一列和第一行都切割掉。
第一行不切割的话，会在读取的时候报错，文件读不进去。而且由于feature是一列的，所以需要加个参数gene.colum=1。
第一列不切割的话，会在使用marker gene做dotplot的时候报错，因为匹配不到细胞和基因。
去head一下数据的话，发现的确是这样。
把第一列切割后就没问题了。
或者有条件的话可以自己跑一下上游定量流程，得到定量的数据。
https://mp.weixin.qq.com/s/M15kWdH8eDONfakNhY-enA