计算药物靶点发现与验证¶

# ========== R脚本：差异基因筛选 + 网络分析 ==========

# 1. 差异表达分析（DESeq2）
library(DESeq2)                                      # 加载DESeq2包

# 读取基因表达矩阵
counts <- read.csv("gene_counts.csv", row.names=1)  # 读取count矩阵
coldata <- data.frame(
    condition=c(rep("disease",3), rep("control",3))  # 分组信息
)

# DESeq2 分析
dds <- DESeqDataSetFromMatrix(counts, coldata,       # 创建DESeq对象
                               design=~condition)
dds <- DESeq(dds)                                     # 运行差异分析
res <- results(dds, alpha=0.05)                       # 提取结果

# 筛选差异基因（|log2FC| > 1, padj < 0.05）
degs <- subset(res, abs(log2FoldChange) > 1 &        # 两倍变化
                    padj < 0.05)                      # 校正p值<0.05
deg_genes <- rownames(degs)                           # 差异基因列表
write.csv(degs, "DEGs.csv")                           # 保存结果
cat("差异基因数量:", length(deg_genes), "\n")

# ========== Python脚本：数据库交叉筛选 ==========
import pandas as pd                                   # 数据处理

# 读取疾病相关基因（从GeneCards/DisGeNET下载）
disease_genes = pd.read_csv("disease_genes.txt",     # 疾病基因列表
                            sep="\t")
disease_set = set(disease_genes["gene_symbol"])        # 转为集合

# 读取药物靶点（从DrugBank/TTD下载）
drug_targets = pd.read_csv("drugbank_targets.txt",   # DrugBank靶点
                            sep="\t")
drug_set = set(drug_targets["gene_symbol"])            # 转为集合

# 读取差异基因
deg_set = set(pd.read_csv("DEGs.csv")["gene"])        # DEG集合

# 取交集：疾病基因 ∩ 药物靶点 ∩ 差异基因
potential_targets = disease_set & drug_set & deg_set   # 三者交集
print(f"潜在靶点数量: {len(potential_targets)}")
print("靶点列表:", potential_targets)

# 保存韦恩图数据
# 用 matplotlib-venn 画韦恩图
from matplotlib_venn import venn3                      # 导入韦恩图工具
import matplotlib.pyplot as plt

plt.figure(figsize=(8, 8))
venn3([disease_set, drug_set, deg_set],                # 三个集合
      set_labels=('Disease Genes', 'Drug Targets', 'DEGs'))
plt.title("Potential Drug Targets")
plt.savefig("venn_targets.pdf")                        # 保存图片

# 用 STRING 数据库构建 PPI 网络
# STRING API 查询
import requests                                        # HTTP请求

genes = list(potential_targets)                        # 潜在靶点列表
string_api = "https://string-db.org/api/tsv/network"  # STRING API

params = {
    "identifiers": "%0d".join(genes),                  # 基因列表
    "species": 9606,                                    # 人类
    "caller_identity": "bioinformatics_script"
}

response = requests.get(string_api, params=params)     # 发送请求
ppi_data = pd.read_csv(pd.io.common.StringIO(response.text),
                        sep="\t")                      # 解析结果

# 保存网络文件（Cytoscape可导入）
ppi_data.to_csv("ppi_network.tsv", sep="\t", index=False)

# 用 networkx 计算网络拓扑
import networkx as nx                                  # 导入网络分析库

G = nx.from_pandas_edgelist(ppi_data,                 # 创建网络
                             source="preferredName_A",
                             target="preferredName_B")

# 计算度中心性（找hub基因/核心靶点）
degree = dict(G.degree())                              # 每个节点的度
sorted_degree = sorted(degree.items(),                 # 按度排序
                       key=lambda x: x[1], reverse=True)

print("Top 10 Hub基因（核心靶点候选）:")
for gene, d in sorted_degree[:10]:
    print(f"  {gene}: degree = {d}")

# 用 AutoDock Vina 做分子对接
# 安装
conda install -c bioconda autodock-vina  # 安装AutoDock Vina

# 准备受体蛋白（从PDB下载）
# 1. 去PDB数据库（https://www.rcsb.org/）下载靶蛋白3D结构
# 2. 用PyMOL去除水分子和配体
# 3. 转换为PDBQT格式

# 准备配体分子（药物小分子）
# 从PubChem/ZINC下载，转为PDBQT格式

# 运行对接
vina --receptor protein.pdbqt \    # 受体蛋白文件
     --ligand drug.pdbqt \          # 配体药物文件
     --center_x 10.0 \              # 对接盒子中心X坐标
     --center_y 20.0 \              # 对接盒子中心Y坐标
     --center_z 30.0 \              # 对接盒子中心Z坐标
     --size_x 20 \                   # 盒子大小X
     --size_y 20 \                   # 盒子大小Y
     --size_z 20 \                   # 盒子大小Z
     --out docking_result.pdbqt \   # 输出结果
     --log docking_log.txt          # 日志文件

# 结合能（affinity）< -7 kcal/mol → 结合力强
# 结合能 < -5 kcal/mol → 中等结合力

# === 药物靶点发现速查 ===

# 常用数据库：
# DrugBank — https://go.drugbank.com/ （药物-靶点）
# TTD — https://db.idrblab.net/ttd/ （治疗靶点）
# GeneCards — https://www.genecards.org/ （基因-疾病）
# DisGeNET — https://www.disgenet.org/ （基因-疾病关联）
# STRING — https://string-db.org/ （PPI网络）
# PDB — https://www.rcsb.org/ （蛋白结构）
# PubChem — https://pubchem.ncbi.nlm.nih.gov/ （药物化合物）

# 常用工具：
# Cytoscape — PPI网络可视化
# AutoDock Vina — 分子对接
# SwissTargetPrediction — 反向靶点预测
# clusterProfiler — 通路富集分析