ixxmu
commented
4 months ago cell2location：结直肠癌肝转移数据分析（一） by 生信菜鸟团
本次根据cell2location教程，使用一篇文献中的数据进行单细胞+空转数据整合分析，目的是对空间转录组数据进行单细胞亚群注释。
大致分析思路如下：
Step1：预处理文献单细胞数据，整理为scanpy包中的anndata对象，为h5ad格式
Step2：对单细胞h5ad格式数据进行 细胞类型signatures评估，即estimate the reference cell type signatures
Step3：预处理文献空间转录组数据，整理为scanpy包中的anndata对象，为h5ad格式
Step4：利用单细胞亚群signature，对空间转录组数据进行spatial mapping
Step5：空间切片中 注释到的 细胞亚群 可视化，也可以与原文的整合结果进行比较看看
数据介绍
本次数据来自文献：
Spatiotemporal Immune Landscape of Colorectal Cancer Liver Metastasis at Single-Cell Level
Cancer Discov（IF 28.2/Q1）. 2022 Jan;12(1):134-153. doi: 10.1158/2159-8290.CD-21-0316. Epub 2021 Aug 20
从这篇文献中，我需要获取的信息有以下两点：
1.文献中的数据下载地址：http://www.cancerdiversity.asia/scCRLM，单细胞为R语言中的rda格式数据，空间为space ranger 软件的标准输出文件
2.文献中单细胞+空间转录组数据整合分析方法：Seurat包进行整合
3.样本情况：
数据情况如A图：总共有来自24个病人的97例配对的单细胞和空间转录组数据，所有肿瘤均为微卫星稳定型
Cohort1：89个单细胞样本
Cohort2：8个空间转录组样本
Cohort3：27个IHC样本
Cohort4：133个芯片样本
Cohort5：430个常规bulk转录组样本
Step1：预处理文献单细胞数据
文献中提供的单细胞数据为两个rda格式：metadata.rda 和 exprmatrix.rda，整理为scanpy包中的anndata对象，为h5ad格式
先读入R对象数据，转为h5ad格式的scanpy对象AnnData
此部分命令在R中运行：
rm(list=ls())
library(Seurat)
opt <- list(input="./data/scCRLM/data/",
od="./data/scCRLM/")
setwd(opt$od)
## 加载metadata
load(paste0(opt$input, "/metadata.rda"))
head(metadata)
table(metadata$main_cell_type)
table(metadata$sub_cell_type)
table(metadata$orig.ident)
table(metadata$patient)
table(metadata$tissue)
## 加载单细胞表达矩阵
load(paste0(opt$input, "/exprmatrix.rda"))
head(colnames(exprmatrix))
dim(exprmatrix)
dim(metadata)
scRNA <- CreateSeuratObject(counts = exprmatrix, meta.data = metadata)
scRNA
## 简单看一下数据的各种指标
scRNA[["percent.mt"]] <- PercentageFeatureSet(scRNA, pattern = "^MT-")
p <- VlnPlot(scRNA, features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol = 3, pt.size = 0)
ggsave(filename = paste0(opt$od,"/01_qc.png"), width = 12, height = 5.5, bg = "white", plot = p)
数据应该是已经进行过滤后的数据了：有178630个细胞，与文献中的是可以对应上的。回去看了一眼文献的质控，做的比较粗糙。
保存一下数据：
推荐大家使用qs进行大数量
qs::qsave(scRNA, "scRNA.qs")
# scRNA = qs::qread("scRNA.qs")
现在进行Seurat转为h5ad格式
#remotes::install_github("mojaveazure/seurat-disk")
### 方式3的详细版本 https://mp.weixin.qq.com/s/eq5XUXXZ8HusXlBlXVFJVA
library(SeuratData)
library(SeuratDisk)
seurat_object <- scRNA
seu <- DietSeurat(
seurat_object,
counts = TRUE, # so, raw counts save to adata.raw.X
data = TRUE, # so, log1p counts save to adata.X
scale.data = FALSE, # set to false, or else will save to adata.X
features = rownames(seurat_object), # export all genes, not just top highly variable genes
assays = "RNA",
dimreducs = c("pca","umap"),
graphs = c("RNA_nn", "RNA_snn"), # to RNA_nn -> distances, RNA_snn -> connectivities
misc = TRUE
)
# step 2: factor to character, or else your factor will be number in adata
i <- sapply(seu@meta.data, is.factor)
seu@meta.data[i] <- lapply(seu@meta.data[i], as.character)
# step 3: convert
SaveH5Seurat(seu, filename = "scRNA.h5seurat", overwrite = TRUE)
Convert("scRNA.h5seurat", "scRNA.h5ad", assay="RNA", overwrite = TRUE)
现在使用python读取查看数据结构
import scanpy as sc
dir = './data/scCRLM/'
adata = sc.read_h5ad(f'{dir}scRNA.h5ad')
print(adata)
使用scanpy包探索一下anndata对象结构
X为矩阵数据，obs为观察值，var为变量；
scanpy软件使用obs存储细胞信息，var存储基因信息。
#查看AnnData
print(adata)
print(adata.X)
print(adata.obs.index)
print(adata.to_df().head())
AnnData object with n_obs × n_vars = 178630 × 20610
X：为矩阵数据
n_obs：细胞数目 178630
n_vars：基因数目 20610
obs: 相当于细胞的 meta.data
var：相当于基因的 meta.data
uns：存储分析方法及参数
obsm：存储降维坐标等多维信息
obsp：存储点与点之间的关系
layers：存储数据，类似于Seurat中的data、counts卡槽
获取基本的信息
# 获取细胞名
adata.obs.index.to_list()
adata.obs.index.to_list()[1:6]
# 获取基因名
adata.var.index.to_list()
adata.var.index.to_list()[1:10]
# 细胞数目
adata.n_obs
# 基因数目
adata.n_vars
# 查看obs某列的值
adata.obs['main_cell_type'].values.tolist()[1:5]
# 可查看每个细胞barcode对应细胞类型
pd.DataFrame({'barcode':adata.obs.index.to_list(),'main_cell_type':adata.obs.main_cell_type.to_list()})
# 可查看细胞类型种类
adata.obs.main_cell_type.value_counts()
# 获取当前object的Variable featurelevels
adata[:, adata.var.highly_variable].var.index
# 存储旧的细胞标签
adata.obs['cluster']=adata.obs.seurat_clusters
对于h5ad与Seurat格式的互换，有一个更全的帖子推荐给大家：
软件测评：百万级单细胞数据的Anndata和Seurat对象互转
下次见~
https://mp.weixin.qq.com/s/n5Kxri1uoEEi3wy6Yway8A