作者,Evil Genius
生活总是无情,牛马也很无奈,除了接受现实,别无其他办法。
公司大调整的日子,也进入倒计时。
不过还是推荐大家走好每一步,生物信息这个行业不算什么暴利行业,但是有个相对安稳的生活还是能做到的,将来大家毕业参加工作,比如进了诺禾、华大等,生活还是不成问题的。
但是随着现在越来越卷,公司的进入门槛也在提高,18年我参加工作基本上是个硕士就行了,不要求有什么代码语言基础,来了公司会培训,25年不行了,不是科班出身,没有代码基础,没有相关的一些基础,很难进去。
以后的路怎么走,谁也不知道,最好自己再说吧。
visium分析思路
今天我们更新脚本
这个图片在课上讲过,强调的是在不同条件下,配受体共定位分布差异,那么今天我们需要更新的是,借助单细胞 + 空间的多组学数据,分析不同条件下的共定位配受体,并且分析条件之间差异的共定位配受体,23年培训的共定位绘图大家有空复习一下。
所有的分析都建立在一定的基础之上,linux 、R、python生信必备,虽然让大家好好学习这个,11月打算开番外,2025番外--linux、R、python培训, 不过实际情况也能想到,没多少人把基础当回事,用的时候才知道着急了。
分两步进行
第一步,单细胞分析不同组别的差异配受体对。
第二步,首先验证同组间配受体的空间共定位情况,去掉单细胞分析得到共定位较差的配受体对,然后分析不同组间共定位差异的配受体对,进行空间分析与展示,这个时候基本上就找到了靶点。
首先实现第一步,实现单细胞不同条件下的差异配受体对,做好细胞类型,放在adata.obs['celltype']下面,不同组别放在adata.obs['batch']下面。
import scanpy as sc
import squidpy as sq
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from scipy.stats import mannwhitneyu
from statsmodels.stats.multitest import multipletests
import warnings
import sys
import os
warnings.filterwarnings('ignore')
# 设置参数
input_file = sys.argv[1] if len(sys.argv) > 1 else "data.h5ad"
output_dir = "./lr_analysis_results"
os.makedirs(output_dir, exist_ok=True)
print("=== 单细胞差异配受体对分析 ===")
# 1. 加载数据
print(f"\n1. 加载数据: {input_file}")
adata = sc.read(input_file)
print(f"数据维度: {adata.n_obs} 个细胞, {adata.n_vars} 个基因")
print(f"批次分布: {adata.obs['batch'].value_counts().to_dict()}")
print(f"细胞类型分布: {adata.obs['celltype'].value_counts().to_dict()}")
# 2. 使用squidpy的ligrec分析配受体相互作用
print("\n2. 使用squidpy ligrec分析配受体相互作用...")
# 第一步:分析每个批次内的配受体相互作用
batches = adata.obs['batch'].unique()
print(f"分析批次: {batches}")
# 存储每个批次的结果
batch_results = {}
for batch in batches:
print(f"\n分析批次: {batch}")
batch_adata = adata[adata.obs['batch'] == batch].copy()
# 使用squidpy的ligrec分析配受体相互作用
sq.gr.ligrec(
batch_adata,
cluster_key="celltype",
copy=True,
use_raw=False,
n_perms=1000
)
# 获取结果
if 'ligrec' in batch_adata.uns:
lr_results = batch_adata.uns['ligrec']
batch_results[batch] = lr_results
# 筛选显著的结果 (pvalue < 0.05)
significant_pairs = lr_results[lr_results['pvalue'] < 0.05]
print(f"批次 {batch} 中发现 {len(significant_pairs)} 个显著配受体对")
# 保存结果
significant_pairs.to_csv(f"{output_dir}/significant_lr_pairs_{batch}.csv", index=False)
# 显示top 10显著对
if len(significant_pairs) > 0:
top_pairs = significant_pairs.nsmallest(10, 'pvalue')
print(f"批次 {batch} Top 10 显著配受体对:")
for idx, row in top_pairs.iterrows():
print(f" {row['ligand']} - {row['receptor']}: pvalue={row['pvalue']:.4f}")
# 3. 比较不同批次间的差异配受体
print("\n3. 比较不同批次间的差异配受体...")
# 收集所有批次中出现的配受体对
all_lr_pairs = set()
for batch, results in batch_results.items():
if len(results) > 0:
pairs = set(zip(results['ligand'], results['receptor']))
all_lr_pairs.update(pairs)
print(f"总共发现 {len(all_lr_pairs)} 个独特的配受体对")
# 分析配受体对在不同批次中的出现情况
lr_presence = {}
for ligand, receptor in all_lr_pairs:
pair_name = f"{ligand}_{receptor}"
presence = []
for batch in batches:
if batch in batch_results and len(batch_results[batch]) > 0:
batch_pairs = set(zip(batch_results[batch]['ligand'], batch_results[batch]['receptor']))
present = 1 if (ligand, receptor) in batch_pairs else 0
else:
present = 0
presence.append(present)
lr_presence[pair_name] = presence
# 创建出现矩阵
presence_df = pd.DataFrame(lr_presence, index=batches).T
print("\n配受体对在不同批次中的出现情况:")
print(presence_df)
# 4. 识别批次特异的配受体对
print("\n4. 识别批次特异的配受体对...")
batch_specific_pairs = {}
for batch in batches:
batch_specific = []
for pair_name in presence_df.index:
row = presence_df.loc[pair_name]
# 只在当前批次出现,其他批次不出现
if row[batch] == 1 and row.sum() == 1:
batch_specific.append(pair_name)
batch_specific_pairs[batch] = batch_specific
print(f"批次 {batch} 特异的配受体对: {len(batch_specific_pairs[batch])} 个")
# 5. 识别共同出现的配受体对
print("\n5. 识别共同出现的配受体对...")
common_pairs = []
for pair_name in presence_df.index:
row = presence_df.loc[pair_name]
if row.sum() == len(batches): # 在所有批次中都出现
common_pairs.append(pair_name)
print(f"在所有批次中都出现的配受体对: {len(common_pairs)} 个")
# 6. 可视化结果
print("\n6. 可视化分析结果...")
# 创建热图显示配受体对的出现情况
plt.figure(figsize=(12, 8))
sns.heatmap(presence_df, cmap="YlGnBu", cbar_kws={'label': '出现 (1) / 不出现 (0)'})
plt.title("配受体对在不同批次中的出现情况")
plt.tight_layout()
plt.savefig(f"{output_dir}/lr_presence_heatmap.png", dpi=300, bbox_inches='tight')
plt.show()
# 创建批次特异性配受体对的可视化
plt.figure(figsize=(10, 6))
batch_counts = [len(pairs) for pairs in batch_specific_pairs.values()]
plt.bar(batch_specific_pairs.keys(), batch_counts)
plt.title("各批次特异的配受体对数量")
plt.ylabel("特异配受体对数量")
plt.xticks(rotation=45)
plt.tight_layout()
plt.savefig(f"{output_dir}/batch_specific_lr_pairs.png", dpi=300, bbox_inches='tight')
plt.show()
# 7. 保存汇总结果
print("\n7. 保存汇总结果...")
# 创建汇总表格
summary_data = []
for batch in batches:
if batch in batch_results and len(batch_results[batch]) > 0:
total_pairs = len(batch_results[batch])
significant_pairs = len(batch_results[batch][batch_results[batch]['pvalue'] < 0.05])
specific_pairs = len(batch_specific_pairs[batch])
else:
total_pairs = 0
significant_pairs = 0
specific_pairs = 0
summary_data.append({
'batch': batch,
'total_lr_pairs': total_pairs,
'significant_pairs': significant_pairs,
'batch_specific_pairs': specific_pairs
})
summary_df = pd.DataFrame(summary_data)
summary_df.to_csv(f"{output_dir}/lr_analysis_summary.csv", index=False)
print("汇总结果:")
print(summary_df)
# 保存共同出现的配受体对
if common_pairs:
common_df = pd.DataFrame({'common_lr_pairs': common_pairs})
common_df.to_csv(f"{output_dir}/common_lr_pairs.csv", index=False)
# 保存批次特异的配受体对
for batch, pairs in batch_specific_pairs.items():
if pairs:
specific_df = pd.DataFrame({f'{batch}_specific_pairs': pairs})
specific_df.to_csv(f"{output_dir}/{batch}_specific_lr_pairs.csv", index=False)
这样我们就得到了单细胞不同条件下,差异的配受体对。
接下来第二步,分析差异性配受体对的空间模式,是否空间存在共定位,空间分布是否有差异(临近包括一个spot临近的6个spot)