作者，Evil Genius
今日参考文献，Spatial transcriptomic characterization of pathologic niches in IPF（Science Advances）。
首先搜集一下marker
单细胞 空间分析确定了三个与疾病相关的生态位，具有独特的细胞组成和定位。这种IPF生态位的空间特征将有助于识别破坏疾病驱动生态位的药物靶点，并有助于开发与疾病相关的体外模型。
知识背景
- 细胞的位置、周围邻居和组织结构，如IPF中的蜂窝形态或成纤维细胞灶，对于确定细胞表型、疾病状态以及最终的细胞功能和通讯至关重要
- 了解每个生态位空间内细胞和分子的相互作用对于更好地解剖独特的疾病病理生物学和指导下一代药物的发现至关重要。
- 空间转录组学通过在完整组织的自然组织环境中提供空间解析的基因表达来解决这一问题，而不会诱导细胞应激、死亡或分离偏差。最近建立的Visium用于福尔马林固定和石蜡包埋(FFPE)技术可以对保存的FFPE组织切片进行无偏置的基于探针的近全转录组mRNA分析，而Xenium原位技术可以对多达500个选定基因进行亚细胞分辨率。
分析结果1、空间和scRNA-seq分析定位健康和IPF肺组织中的细胞和基因
- cell2location
- 病理学家对HE染色组织切片的组织病理学注释
- 细胞类型共定位，异常的细胞共定位引起了疾病的发生和变化
- Xenium平台验证Visium的分析结果
分析结果2、三种不同的niche出现在IPF肺组织中
- 非负矩阵分解(NMF)
- 将cell2location逐细胞类型输出矩阵转换为细胞类型和生态位矩阵，以及spot和生态位矩阵。这能够确定每个spot内的主要细胞生态位和表征每个生态位的细胞类型，抽象出生物系统的复杂性，并确定组织内最相关的细胞邻近区。
- 生态位的空间组织，利用位置数据并分析相邻的相互作用。
分析结果3、在气道周围发现异常的基底样细胞
分析结果4、SPP1 巨噬细胞亚群定位于纤维化气道的管腔
进一步研究气道巨噬细胞生态位，我们首先旨在了解IPF中细胞巨噬细胞亚型的异质性，这些亚型都表达spp1 -被认为是许多纤维化疾病中常见的纤维化巨噬细胞标志物。空间Visium数据揭示了巨噬细胞亚型在IPF中的差异定位，并提示气道巨噬细胞生态位优先定位于纤维化气道的管腔，可以在mRNA和蛋白质水平上验证这一点。气道巨噬细胞生态位也出现在纤维化生态位附近，特别是异常的基底细胞，这是由预测的和独特的细胞间通讯支持的。
分析结果5、免疫细胞病灶是由ipf特异性异位内皮细胞募集的
空间转录组学数据以及H&E组织结构显示，免疫生态位是IPF的一个显著组织病理学特征，以多种淋巴样细胞类型的病灶形式存在，可能通过异位内皮支气管血管的不同分泌信号募集。
来看看代码部分，我们看一下重点部分
重点、hotspot（空间打分）
代码语言：javascript复制import scanpy as sc
import squidpy as sq
import numpy as np
import pandas as pd
import os
import sys
import seaborn as sb
import matplotlib.pyplot as plt
from matplotlib import colors
#import scvi
import anndata as ad
import warnings
warnings.filterwarnings("ignore")
from collections import Counter
import ipywidgets as widgets
from ipywidgets import interact, interact_manual
plt.rcParams['figure.figsize'] = (6, 6)
from IPython.core.display import display, HTML
import random
#Define a colour map for gene expression
colors2 = plt.cm.Reds(np.linspace(0, 1, 128))
colors3 = plt.cm.Greys_r(np.linspace(0.7,0.8,20))
#colorsComb = np.vstack([colors3, colors2])
#mymap = colors.LinearSegmentedColormap.from_list('my_colormap', colorsComb)
from matplotlib import colors
colorsComb = np.vstack([plt.cm.Reds(np.linspace(0, 1, 128)), plt.cm.Greys_r(np.linspace(0.7, 0.8, 0))])
mymap = colors.LinearSegmentedColormap.from_list('my_colormap', colorsComb)
# Helper function to split list in chunks
def chunks(lista, n):
for i in range(0, len(lista), n):
yield lista[i:i n]
plt.rcParams['figure.figsize'] = (6, 5)
sc.set_figure_params(dpi=100, vector_friendly=True)
def mysize(w, h, d):
fig, ax = plt.subplots(figsize = (w, h), dpi = d)
return(fig.gca())
plt.rcParams['figure.figsize'] = (6, 5)
sc.set_figure_params(dpi=100, vector_friendly=True)
sc.settings.figdir = "./figures/"
import scvi
## frequently used variables
from matplotlib import colors
import matplotlib.pyplot as plt
colorsComb = np.vstack([plt.cm.Reds(np.linspace(0, 1, 128)), plt.cm.Greys_r(np.linspace(0.7, 0.8, 0))])
mymap = colors.LinearSegmentedColormap.from_list('my_colormap', colorsComb)
## Along these Lines, a colourmap diverging from gray to red
gray_red = colors.LinearSegmentedColormap.from_list("grouping", ["lightgray", "red", "darkred"], N = 128)
## Some more Colour Maps
gray_violet = colors.LinearSegmentedColormap.from_list("grouping", ["lightgray", "mediumvioletred", "indigo"], N = 128)
gray_blue = colors.LinearSegmentedColormap.from_list("grouping", ["lightgray", "cornflowerblue", "darkblue"], N = 128)
def mysize(w, h, d):
fig, ax = plt.subplots(figsize = (w, h), dpi = d)
return(fig.gca())
#plt.rcParams['figure.figsize'] = (6, 5)
#sc.set_figure_params(dpi=120, vector_friendly=True)
import matplotlib.colors as colors
c_low = colors.colorConverter.to_rgba('orange', alpha = 0)
c_high = colors.colorConverter.to_rgba('red',alpha = 1)
cmap_transparent = colors.LinearSegmentedColormap.from_list('rb_cmap',[c_low, c_high], 512)
import matplotlib.colors as colors
c_low2 = colors.colorConverter.to_rgba('green', alpha = 0)
c_high2 = colors.colorConverter.to_rgba('darkblue',alpha = 1)
cmap_transparent2 = colors.LinearSegmentedColormap.from_list('rb_cmap',[c_low2, c_high2], 512)
print(f"squidpy=={sq.__version__}")
print(f"scanpy=={sc.__version__}")
import cell2location as c2l
from cell2location.utils import select_slide
##load
adata_vis = sc.read(f".h5ad")
Niche_NMF_palette = dict(zip(adata_vis.obs.Niche_NMF.cat.categories, adata_vis.uns["Niche_NMF_colors"]))
sc.pl.umap(adata_vis, color="Niche_NMF", palette=Niche_NMF_palette,)
import hotspot
def compute_scores_hotspot(adata,
genes,
layer=None,
n_neighbors = 30,
neighborhood_factor = 3,
gene_symbols = None,
use_rep = "X_scVI",
model = 'danb'
):
if use_rep is None:
if layer is None:
index = pynndescent.NNDescent(adata.X, n_neighbors=n_neighbors 1)
else:
index = pynndescent.NNDescent(adata.layers[layer], n_neighbors=n_neighbors 1)
else:
index = pynndescent.NNDescent(adata.obsm[use_rep], n_neighbors=n_neighbors 1)
ind, dist = index.neighbor_graph
ind, dist = ind[:, 1:], dist[:, 1:]
ind = pd.DataFrame(ind, index=list(range(adata.n_obs)))
neighbors = ind
weights = compute_weights(dist, neighborhood_factor=neighborhood_factor)
weights = pd.DataFrame(weights, index=neighbors.index,
columns=neighbors.columns)
if layer is None:
if scipy.sparse.issparse(adata.X):
umi_counts = adata.X.sum(axis=1).A1
else:
umi_counts = adata.X.sum(axis=1)
else:
if scipy.sparse.issparse(adata.layers[layer]):
umi_counts = adata.layers[layer].sum(axis=1).A1
else:
umi_counts = adata.layers[layer].sum(axis=1)
if gene_symbols is None:
counts_dense = Hotspot._counts_from_anndata(
adata[:, adata.var_names.isin(genes)], layer, dense=True)
else:
counts_dense = Hotspot._counts_from_anndata(
adata[:, adata.var[gene_symbols].isin(genes)], layer, dense=True)
scores = modules.compute_scores(
counts_dense,
model,
umi_counts,
neighbors.values,
weights.values,
)
return scores
import pynndescent
import scipy.sparse
import hotspot
from hotspot import modules
from hotspot.knn import compute_weights
from hotspot.hotspot import Hotspot
senescence =["GLB1","TP53","SERPINE1","CDKN1A","CDKN1B","CDKN2B"]
scores = compute_scores_hotspot(adata_vis, senescence,layer="log1p",n_neighbors = 30,neighborhood_factor = 3,gene_symbols = None,use_rep = "X_scVI",model = 'danb')
adata_vis.obs["senescence_hs_score"] = scores
#for i in fibs:
gene = "senescence_hs_score"
vmins = None
vcenters= None
vmaxs = 2
img_keys='lowres'
cmaps= "bwr"
palettes=Niche_NMF_palette
group= None #"Fibrotic"
import matplotlib.pyplot as plt
from cell2location.utils import select_slide
fig, (ax1, ax2, ax3, ax4, ax5,ax6) = plt.subplots(1, 6, figsize=(26,4), )
plt.suptitle(" control ", y=1.05)
with plt.rc_context({'axes.facecolor': 'white','figure.figsize': [4, 4]}):
slide = select_slide(adata_vis, "90_C1_RO-730_Healthy_processed_CM")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax1, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
slide = select_slide(adata_vis, "91_A1_RO-727_Healthy_processed_CM")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax2, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
slide = select_slide(adata_vis, "91_B1_RO-728_Healthy_processed_CM")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax3, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
slide = select_slide(adata_vis, "1217_0002_processed_aligned")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax4, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
slide = select_slide(adata_vis, "92_A1_RO-3203_Healthy_processed_CM")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax5, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
slide = select_slide(adata_vis, "1217_0004_processed_aligned")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax6, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
plt.savefig("./figures/spatial_plot_senescence_hs_score_healthy_vmax2.pdf")
fig, (ax7, ax8, ax9, ax10, ax11,) = plt.subplots(1, 5, figsize=(22,4), )
plt.suptitle(" IPF ", y=1.05)
with plt.rc_context({'axes.facecolor': 'white','figure.figsize': [4, 4]}):
slide = select_slide(adata_vis, "90_A1_H237762_IPF_processed_CM")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax7, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
slide = select_slide(adata_vis, "1217_0001_processed_aligned")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax8, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
slide = select_slide(adata_vis, "91_D1_24513-17_IPF_processed_CM")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax9, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
slide = select_slide(adata_vis, "1217_0003_processed_aligned")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax10, show=False,vcenter= vcenters,palette=palettes, legend_loc=None)
slide = select_slide(adata_vis, "92_D1_RO-3736_IPF_processed_CM")
sc.pl.spatial(slide, cmap=cmaps,color=gene,size=1.3,img_key=img_keys,vmin=vmins, use_raw=False,layer="log1p",groups=group, vmax=vmaxs, ax=ax11, show=False,vcenter= vcenters,palette=palettes)
plt.savefig("./figures/spatial_plot_senescence_hs_score_IPF_vmax2.pdf") 
