本文首发于“生信补给站”公众号 https://mp.weixin.qq.com/s/WJQkN_jtYg19yGz0uq8tyg
前面在scRNA分析|多样本merge 和 harmony去批次提到单细胞的数据形式有很多种。如果提供的是标准的10X的三个文件就可以直接read10X读取，那如果只有矩阵文件如何进行下游分析呢？
如果额外给了细胞水平的meta文件，如何利用呢？
本文以2021年12月发表在nature cancer上的文献数据为例，读取提供的GSE179994中的矩阵和meta数据。
一 载入数据，R包
文章https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE179994只提供了矩阵文件，不能使用Read10X函数的形式，但是也可以很简单的读取。
1.1 读取下载的数据
1) 读取矩阵文件
注意区分rds 和 RData文件的读取方式 。如果是文本文件的话直接read.table 或者 readMM 读取。
library(tidyverse) library(Seurat) data <- readRDS("GSE179994_all.Tcell.rawCounts.rds") data[1:4,1:4] 4 x 4 sparse Matrix of class "dgCMatrix" P1.ut.AAACCTGAGGCACATG-1 P1.ut.AAACCTGGTGGTACAG-1 P1.ut.AAACCTGTCACGGTTA-1 P1.ut.AAACCTGTCCGCTGTT-1 OR4F5 . . . . OR4F29 . . . . OR4F16 . . . . SAMD11
2）读取meta文件
meta <- read.table("GSE179994_Tcell.metadata.tsv",header = T ) head(meta) cellid patient sample celltype cluster 1 P1.ut.AAACCTGAGGCACATG-1 P1 P1.pre <NA> <NA> 2 P1.ut.AAACCTGGTGGTACAG-1 P1 P1.pre CD4 CD4_C8-Treg 3 P1.ut.AAACCTGTCACGGTTA-1 P1 P1.pre CD4 CD4_C8-Treg 4 P1.ut.AAACCTGTCCGCTGTT-1 P1 P1.pre CD4 CD4_C9-Prolif. 5 P1.ut.AAACCTGTCTGGTATG-1 P1 P1.pre <NA> <NA> 6 P1.ut.AAACCTGTCTTGAGAC-1 P1 P1.pre CD4 CD4_C8-Treg
可以看到提供了cell level的注释信息。
1.2 创建seurat对象
依然使用CreateSeuratObject 函数，此处count 为读取的矩阵文件。
sce0 <- CreateSeuratObject(counts = data) sce0 head(sce0@meta.data)
An object of class Seurat 19790 features across 150849 samples within 1 assay Active assay: RNA (19790 features, 0 variable features)
orig.ident nCount_RNA nFeature_RNA P1.ut.AAACCTGAGGCACATG-1 SeuratProject 3265 1241 P1.ut.AAACCTGGTGGTACAG-1 SeuratProject 3871 1783 P1.ut.AAACCTGTCACGGTTA-1 SeuratProject 3945 1816 P1.ut.AAACCTGTCCGCTGTT-1 SeuratProject 21848 4426 P1.ut.AAACCTGTCTGGTATG-1 SeuratProject 2959 1376 P1.ut.AAACCTGTCTTGAGAC-1 SeuratProject 1809 1024
之后就可以进行标准的seurat流程，然后进行个性化挖掘分析了。
二 添加meta信息
文献提供的metadata文件中含有每个细胞的patient ，sample ，celltype 以及 cluster ，是很有利用价值的 。
之前在scRNA分析|Marker gene 可视化 以及 细胞亚群注释--你是如何人工注释的？中提到了如何添加亚群注释（cluster level）结果到metadata的方式，这里介绍下如何添加每个细胞（cell level）的metadata。
1，CreateSeuratObject中的meta.data参数
CreateSeuratObject函数除了简单的过滤条件外 ，还有一个重要的meta.data参数，可以输入提供的meta信息。
sce1 <- CreateSeuratObject(counts = data, meta.data =meta, min.cells = 3, # 可以注释掉 min.features = 200, # 可以注释掉 project = "GSE179994")c head(sce1@meta.data) sce1 <- CreateSeuratObject(counts = data, meta.data =meta, min.cells = 3, # 可以注释掉 min.features = 200, # 可以注释掉 project = "GSE179994")c head(sce1@meta.data)
可以看到虽然添加上了meta的信息列，但是全是NA，这是什么原因呢？
检索之后https://github.com/satijalab/seurat/issues/2715 发现 ，是因为CreateSeuratObject要求meta文件中rownames是count文件的colnames 。
其实 ??CreateSeuratObject函数的帮助文档中也很明确的提到了该点要求。
发现问题后，只需要将meta文件的cellid列转为rownames即可。
meta2 <- meta %>% column_to_rownames("cellid") sce2 <- CreateSeuratObject(counts = data, meta.data =meta2, min.cells = 3, min.features = 200, project = "GSE179994") head(sce2@meta.data)
2，AddMetaData函数
也可以在sce0的基础上，使用AddMetaData函数进行添加
sce3 <- AddMetaData( object = sce0, metadata = meta2 ) head(sce3@meta.data)
3，矩阵添加
建议优先使用上述两种单细胞特有的方式进行添加。
sce4 <- sce0 sce4@meta.data <- sce4@meta.data %>% rownames_to_column("cellid") %>% inner_join(meta) %>% column_to_rownames("cellid") #一定要转回来 head(sce4@meta.data)
注意一定要把cell barcode的列转为rowname ！
OK，搞定！以上三种方式结果是一致的。
建议可以先自行注释，然后参照文献注释结果修正自己的marker库，或者找到一些手动注释的思路。
