- 学习资料来源：
- scanpy主页：https://scanpy.readthedocs.io/en/stable/
- 官网：https://scanpy-tutorials.readthedocs.io/en/latest/spatial/integration-scanorama.html【注意教程有两个版本，这里是latest版本的学习笔记】
这篇教程主要介绍怎么使用scanppy进行多个Visium空间转录组数据分析，以及与单细胞数据的整合分析。
教程将使用Scanorama paper - code进行数据整合与label transfer。
加载函数库
# 安装包
pip install scanorama
# conda 版本：conda install -c bioconda scanorama -y
import scanpy as sc
import anndata as an
import pandas as pd
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
import seaborn as sns
import scanorama
sc.logging.print_versions()
sc.set_figure_params(facecolor="white", figsize=(8, 8))
sc.settings.verbosity = 3
outdir = '/Pub/Users/project/scanpy/stRNA/'
读取数据
教程将使用两个小鼠大脑(矢状)空间数据，来源：10x genomics website
adata_spatial_anterior = sc.datasets.visium_sge(sample_id="V1_Mouse_Brain_Sagittal_Anterior")
adata_spatial_posterior = sc.datasets.visium_sge(sample_id="V1_Mouse_Brain_Sagittal_Posterior")
adata_spatial_anterior.var_names_make_unique()
adata_spatial_posterior.var_names_make_unique()
两个数据的情况：
adata_spatial_anterior
AnnData object with n_obs × n_vars = 2695 × 32285
obs: 'in_tissue', 'array_row', 'array_col'
var: 'gene_ids', 'feature_types', 'genome'
uns: 'spatial'
obsm: 'spatial'
adata_spatial_posterior
AnnData object with n_obs × n_vars = 3355 × 32285
obs: 'in_tissue', 'array_row', 'array_col'
var: 'gene_ids', 'feature_types', 'genome'
uns: 'spatial'
obsm: 'spatial'
QC并可视化
sc.pp.calculate_qc_metrics(adata_spatial_anterior, inplace=True)
sc.pp.calculate_qc_metrics(adata_spatial_posterior, inplace=True)
for name, adata in [
("anterior", adata_spatial_anterior),
("posterior", adata_spatial_posterior),
]:
fig, axs = plt.subplots(1, 4, figsize=(12, 3))
fig.suptitle(f"Covariates for filtering: {name}")
sns.distplot(adata.obs["total_counts"], kde=False, ax=axs[0])
sns.distplot(
adata.obs["total_counts"][adata.obs["total_counts"] < 20000],
kde=False,
bins=40,
ax=axs[1],
)
sns.distplot(adata.obs["n_genes_by_counts"], kde=False, bins=60, ax=axs[2])
sns.distplot(
adata.obs["n_genes_by_counts"][adata.obs["n_genes_by_counts"] < 4000],
kde=False,
bins=60,
ax=axs[3],
)
plt.savefig(outdir + "01-sns.distplot_"+name+".png")
结果图：
anterior：
posterior：
标准化
使用normalize_total函数对数据进行标准化，然后高变基因选择。
还可以选择SCTransform or GLM-PCA进行数据标准化。
for adata in [
adata_spatial_anterior,
adata_spatial_posterior,
]:
sc.pp.normalize_total(adata, inplace=True)
sc.pp.log1p(adata)
sc.pp.highly_variable_genes(adata, flavor="seurat", n_top_genes=2000, inplace=True)
数据整合
现在使用Scanorama对两个数据进行整合，Scanorama对每个数据集会返回两个list对象，一个是整合后的embeddings，一个是corrected counts。
Note：空间数据也可以使用前面教程分享的 BBKNN or Ingest进行数据整合。
adatas = [adata_spatial_anterior, adata_spatial_posterior]
adatas_cor = scanorama.correct_scanpy(adatas, return_dimred=True)
两个数据连接在一起，并保存在整合后的embeddings：adata_spatial.obsm['scanorama_embedding']中。
此外，我们将计算UMAP，以可视化的结果和定性评估数据整合的结果。
Note：我们使用 un _ merge = “ only”策略连接这两个数据集，以便在连接的 anndata 对象中保持两个图像来自 visium 数据集。
adata_spatial = sc.concat(
adatas_cor,
label="library_id",
uns_merge="unique",
keys=[
k
for d in [
adatas_cor[0].uns["spatial"],
adatas_cor[1].uns["spatial"],
]
for k, v in d.items()
],
index_unique="-",
)
sc.pp.neighbors(adata_spatial, use_rep="X_scanorama")
sc.tl.umap(adata_spatial)
sc.tl.leiden(adata_spatial, key_added="clusters")
sc.pl.umap(
adata_spatial, color=["clusters", "library_id"], palette=sc.pl.palettes.default_20
)
plt.savefig(outdir + "02-spatial.png")
结果图：
我们还可以在空间坐标系中可视化聚类结果。为此，我们首先需要将cluster颜色保存在字典中。然后我们可以绘制Anterior and Posterior的Visium组织的矢状图，并排在一起。
clusters_colors = dict(
zip([str(i) for i in range(18)], adata_spatial.uns["clusters_colors"])
)
fig, axs = plt.subplots(1, 2, figsize=(15, 10))
for i, library in enumerate(
["V1_Mouse_Brain_Sagittal_Anterior", "V1_Mouse_Brain_Sagittal_Posterior"]
):
ad = adata_spatial[adata_spatial.obs.library_id == library, :].copy()
sc.pl.spatial(
ad,
img_key="hires",
library_id=library,
color="clusters",
size=1.5,
palette=[
v
for k, v in clusters_colors.items()
if k in ad.obs.clusters.unique().tolist()
],
legend_loc=None,
show=False,
ax=axs[i],
)
plt.tight_layout()
plt.savefig(outdir + "02-spatial-1.png")
结果图：从聚类图中，我们可以清楚地看到两个组织的皮层分层(参见 Allen 脑图谱)。此外，由于两个组织中的簇之间存在明显的连续性，因此数据集合整合似乎工作得很好。
与单细胞数据进行整合
scanorama还可以对 scRNA-seq 数据集和空间转录组学数据集之间进行数据整合。这样的整合特别有用，因为它允许我们将从 scRNA-seq 数据集中确定的细胞类型标签转移到 Visium 数据集。
本次将使用来自Tasic et al.的数据， 这个数据利用 smart-seq 技术对小鼠皮层进行了分析。
数据编号为：GSE115746，可以直接下载处理好的h5ad格式：https://hmgubox.helmholtz-muenchen.de/f/4ef254675e2a41f89835/?dl=1
由于数据集是从小鼠皮层cortex产生的，我们将提取Visium数据集的子集，以便只选择皮层的spots部分。注意，整合也可以在整个大脑切片上进行，但是它会产生假阳性的细胞类型分配，因此应该更加谨慎地进行解释。
dir = '/Pub/Users/zhangjuan/project/scanpy/data/GSE115746/'
adata_cortex = sc.read(dir+"./adata_processed.h5ad")
提取子集：cortex
adata_anterior_subset = adata_spatial_anterior[
adata_spatial_anterior.obsm["spatial"][:, 1] < 6000, :
]
adata_posterior_subset = adata_spatial_posterior[
(adata_spatial_posterior.obsm["spatial"][:, 1] < 4000)
& (adata_spatial_posterior.obsm["spatial"][:, 0] < 6000),
:,
]
使用Scanorama进行整合：
adatas_anterior = [adata_cortex, adata_anterior_subset]
adatas_posterior = [adata_cortex, adata_posterior_subset]
# Integration.
adatas_cor_anterior = scanorama.correct_scanpy(adatas_anterior, return_dimred=True)
adatas_cor_posterior = scanorama.correct_scanpy(adatas_posterior, return_dimred=True)
连接数据集并将整合的embeddings分配到 anndata 对象中：
Note：这里采用join="outer"与uns_merge="first"的策略将两个数据连接起来，这是因为我们希望保留visium数据的 obsm['coords']以及 images。
adata_cortex_anterior = sc.concat(
adatas_cor_anterior,
label="dataset",
keys=["smart-seq", "visium"],
join="outer",
uns_merge="first",
)
adata_cortex_posterior = sc.concat(
adatas_cor_posterior,
label="dataset",
keys=["smart-seq", "visium"],
join="outer",
uns_merge="first",
)
在上一步，我们已经将每个visium数据集与scRNA-seq数据集集成在一个共同的embedding中。在这样的embedding空间中，我们可以计算样本之间的距离，并使用这些距离作为权重，用于将细胞类型标签从scRNA-seq数据集转移到Visium数据集。
这种方法与Seurat中的TransferData函数非常相似(see paper)。这里，我们用一个简单的python函数重新实现了标签传递函数。
首先，让我们在相同的embedding空间中计算visium数据集和scRNA-seq数据集之间的余弦距离：
from sklearn.metrics.pairwise import cosine_distances
distances_anterior = 1 - cosine_distances(
adata_cortex_anterior[adata_cortex_anterior.obs.dataset == "smart-seq"].obsm[
"X_scanorama"
],
adata_cortex_anterior[adata_cortex_anterior.obs.dataset == "visium"].obsm[
"X_scanorama"
],
)
distances_posterior = 1 - cosine_distances(
adata_cortex_posterior[adata_cortex_posterior.obs.dataset == "smart-seq"].obsm[
"X_scanorama"
],
adata_cortex_posterior[adata_cortex_posterior.obs.dataset == "visium"].obsm[
"X_scanorama"
],
)
然后，让我们将标签从scRNA-seq数据集传播到visium数据集：
def label_transfer(dist, labels):
lab = pd.get_dummies(labels).to_numpy().T
class_prob = lab @ dist
norm = np.linalg.norm(class_prob, 2, axis=0)
class_prob = class_prob / norm
class_prob = (class_prob.T - class_prob.min(1)) / class_prob.ptp(1)
return class_prob
class_prob_anterior = label_transfer(distances_anterior, adata_cortex.obs.cell_subclass)
class_prob_posterior = label_transfer(distances_posterior, adata_cortex.obs.cell_subclass)
class_prob_[anterior-posterior]对象是一个numpy形状数组(cell_type, visium_spots)，包含每个spot分配给每个细胞类型的权重。这个值本质上告诉我们，从表达值的角度来看，这些spots与来自scRNA-seq数据集的所有其他注释细胞类型有多相似。
将class_prob_[anterior-posterior]对象转换成dataframe，并分配到对应的anndata数据中。
cp_anterior_df = pd.DataFrame(
class_prob_anterior, columns=np.sort(adata_cortex.obs.cell_subclass.unique())
)
cp_posterior_df = pd.DataFrame(
class_prob_posterior, columns=np.sort(adata_cortex.obs.cell_subclass.unique())
)
cp_anterior_df.index = adata_anterior_subset.obs.index
cp_posterior_df.index = adata_posterior_subset.obs.index
adata_anterior_subset_transfer = adata_anterior_subset.copy()
adata_anterior_subset_transfer.obs = pd.concat(
[adata_anterior_subset.obs, cp_anterior_df], axis=1
)
adata_posterior_subset_transfer = adata_posterior_subset.copy()
adata_posterior_subset_transfer.obs = pd.concat(
[adata_posterior_subset.obs, cp_posterior_df], axis=1
)
然后，我们能够探索细胞类型如何从scRNA-seq数据集转移到visium数据集。让我们先看看神经元皮层。
sc.pl.spatial(
adata_anterior_subset_transfer,
img_key="hires",
color=["L2/3 IT", "L4", "L5 PT", "L6 CT"],
size=1.5,
)
plt.savefig(outdir + "03-adata_anterior_subset_transfer.png")
sc.pl.spatial(
adata_posterior_subset_transfer,
img_key="hires",
color=["L2/3 IT", "L4", "L5 PT", "L6 CT"],
size=1.5,
)
plt.savefig(outdir + "03-adata_posterior_subset_transfer.png")
结果图如下：
anterior：
posterior：
有趣的是，这种方法似乎是有效的，因为在前矢状slide和后矢状slide中，连续的皮层神经元层都可以被正确识别。
同样，我们可视化星形胶质细胞astrocytes和少突胶质细胞oligodendrocytes。
sc.pl.spatial(
adata_anterior_subset_transfer, img_key="hires", color=["Oligo", "Astro"], size=1.5
)
plt.savefig(outdir + "03-adata_anterior_subset_transfer_Oligo_Astro.png")
sc.pl.spatial(
adata_posterior_subset_transfer, img_key="hires", color=["Oligo", "Astro"], size=1.5
)
plt.savefig(outdir + "03-adata_posterior_subset_transfer_Oligo_Astro.png")
结果图：
anterior：
posterior：
在本次教程中，我们展示了如何使用scanpy整合分析多张slices切片数据，展示了注释后的单细胞数据在空间数据上的映射，我们展示了这种利用Scanorama的数据集成性能的方法是有用的，并且为探索性分析提供了一个直接的工具。
然而，对于 label transfer 分析，我们建议分析人员探索更有原则的方法，基于细胞类型反褶积，这可能会提供更准确和可解释的结果。参见最近的方法，例如：
