读入csv文件的表达矩阵成为Seurat对象（CNS图表复现01）

Posted on by

上个月我们组建了：《单细胞CNS图表复现交流群》，见：你要的rmarkdown文献图表复现全套代码来了（单细胞），也分享了单细胞转录组数据分析的流程：

交流群里大家讨论的热火朝天，而且也都开始了图表复现之旅，在这里我还是带大家一步步学习CNS图表吧。如果你也想加入交流群，自己去：你要的rmarkdown文献图表复现全套代码来了（单细胞）找到我们的拉群小助手哈。

今天讲解第一步：读入csv文件的表达矩阵成为Seurat对象。

交流群里共享了该文章的谷歌云文件的百度云版本，方便在中国大陆的朋友们下载。

9G的文件，谷歌云文件的百度云版本

读取两个csv文件

如果是10x数据集，一般来说是每个样本都有3个文件，参考代码是：

但是我们这个CNS文章是直接把表达矩阵给出来了csv文件，所以直接读取，代码如下：

rm(list=ls())
options(stringsAsFactors = F)
library(Seurat)
# Load data 
dir='./'
Sys.time()
raw.data <- read.csv(paste(dir,"Data_input/csv_files/S01_datafinal.csv", sep=""), header=T, row.names = 1)
Sys.time()
dim(raw.data)
raw.data[1:4,1:4]
head(colnames(raw.data)) 
# Load metadata 
metadata <- read.csv(paste(dir,"Data_input/csv_files/S01_metacells.csv", sep=""), row.names=1, header=T)
head(metadata)

# Find ERCC's, compute the percent ERCC, and drop them from the raw data.
erccs <- grep(pattern = "^ERCC-", x = rownames(x = raw.data), value = TRUE)
percent.ercc <- Matrix::colSums(raw.data[erccs, ])/Matrix::colSums(raw.data)
fivenum(percent.ercc)
ercc.index <- grep(pattern = "^ERCC-", x = rownames(x = raw.data), value = FALSE)
raw.data <- raw.data[-ercc.index,]
dim(raw.data)

构建Seurat对象

有了表达矩阵，直接使用 CreateSeuratObject 函数即可，然后慢慢添加这个表达矩阵的一些其它外部属性，全部代码如下：


# Create the Seurat object with all the data (unfiltered)
main_tiss <- CreateSeuratObject(counts = raw.data)
# add rownames to metadta 
row.names(metadata) <- metadata$cell_id
# add metadata to Seurat object 
main_tiss <- AddMetaData(object = main_tiss, metadata = metadata)
main_tiss <- AddMetaData(object = main_tiss, percent.ercc, col.name = "percent.ercc")
# Head to check
head(main_tiss@meta.data)

# Calculate percent ribosomal genes and add to metadata
ribo.genes <- grep(pattern = "^RP[SL][[:digit:]]", x = rownames(x = main_tiss@assays$RNA@data), value = TRUE)
percent.ribo <- Matrix::colSums(main_tiss@assays$RNA@counts[ribo.genes, ])/Matrix::colSums(main_tiss@assays$RNA@data)
fivenum(percent.ribo)
main_tiss <- AddMetaData(object = main_tiss, metadata = percent.ribo, col.name = "percent.ribo")
main_tiss

# Filter cells so that remaining cells have nGenes >= 500 and nReads >= 50000
main_tiss_filtered <- subset(x=main_tiss, subset = nCount_RNA > 50000 & nFeature_RNA > 500)
main_tiss_filtered

我们得到了main_tiss_filtered这个变量，是一个Seurat对象，就可以follow我们的教程后后续分析流程啦。参考代码仍然是是：

文末友情推荐

要想真正入门生物信息学建议务必购买全套书籍，一点一滴攻克计算机基础知识，书单在：什么，生信入门全套书籍仅需160 。
如果大家没有时间自行慢慢摸索着学习，可以考虑我们生信技能树官方举办的学习班：

如果你课题涉及到转录组，欢迎添加一对一客服：详见：你还在花三五万做一个单细胞转录组吗？

号外：生信技能树知识整理实习生招募，长期招募，也可以简单参与软件测评笔记撰写，开启你的分享人生！另外，：绝大部分生信技能树粉丝都没有机会加我微信，已经多次满了5000好友，所以我开通了一个微信好友，前100名添加我，仅需150元即可，3折优惠期机会不容错过哈。我的微信小号二维码在：0元，10小时教学视频直播《跟着百度李彦宏学习肿瘤基因组测序数据分析》

​

Comments are closed.