Python在生物信息学中的统计应用

Python凭借其强大的科学计算库和易用性，已成为生物信息学统计分析的核心工具。以下是关键应用领域及示例：

---

1. 基因表达数据分析

利用统计学方法分析RNA-seq数据，检测差异表达基因：

import numpy as np
import pandas as pd
from scipy import stats

# 模拟基因表达矩阵（行：基因，列：样本）
data = pd.DataFrame({
    'GeneA': np.random.normal(10, 2, 20),  # 对照组
    'GeneB': np.random.normal(15, 3, 20)   # 实验组
})

# 独立样本t检验
t_stat, p_value = stats.ttest_ind(data['GeneA'], data['GeneB'])
print(f"基因B差异表达显著性: p值 = {p_value:.4f}")

---

2. 基因组关联研究（GWAS）

通过统计检验识别疾病相关基因位点：

import statsmodels.api as sm

# 模拟基因型数据 (0,1,2 表示等位基因数量)
genotypes = np.random.choice([0,1,2], size=100)
phenotypes = np.where(genotypes > 1, 1, 0)  # 疾病表型

# 逻辑回归分析
model = sm.Logit(phenotypes, sm.add_constant(genotypes))
result = model.fit()
print(result.summary())  # 输出OR值及p值

---

3. 蛋白质结构统计分析

分析蛋白质构象的几何参数分布：

import matplotlib.pyplot as plt
from bio.PDB import PDBParser

# 解析PDB文件获取二面角
parser = PDBParser()
structure = parser.get_structure("protein", "protein.pdb")
phi_angles = []  # 存储Φ角

# 拉氏图绘制
plt.figure(figsize=(8,6))
plt.hist(phi_angles, bins=50, density=True)
plt.xlabel('Φ角 (弧度)')
plt.ylabel('概率密度')
plt.title('蛋白质主链二面角分布')
plt.show()

---

4. 进化树构建与检验

基于遗传距离构建系统发育树：

from Bio.Phylo.TreeConstruction import DistanceTreeConstructor
from Bio.Phylo.Consensus import bootstrap_trees

# 计算遗传距离矩阵
distance_matrix = calculate_evolutionary_distances(sequences)

# 自举法检验树节点稳定性
bootstrapped_trees = bootstrap_trees(distance_matrix, 1000)  
support_values = get_support_values(main_tree, bootstrapped_trees)

---

5. 富集分析

识别基因功能富集通路：

from gseapy import enrichr

# 差异基因列表
diff_genes = ['TP53', 'BRCA1', 'EGFR', 'MYC']

# KEGG通路富集分析
enrichment_results = enrichr(
    gene_list=diff_genes,
    gene_sets='KEGG_2021_Human'
)
enrichment_results.results.head()  # 展示Top5富集通路

---

核心统计库

库名称	应用场景
SciPy	假设检验、概率分布
StatsModels	回归分析、方差分析
scikit-learn	机器学习分类/聚类
Biopython	生物序列统计分析
lifelines	生存分析（如癌症预后）

---

典型案例

1. 癌症亚型分类

   • 使用PCA降维 + K-means聚类
   • 通过轮廓系数$$ s(i) = \frac{b(i) - a(i)}{\max[a(i),b(i)]} $$ 评估聚类质量

2. 药物剂量响应分析

   • 拟合四参数逻辑曲线：
     $$ y = \frac{A-D}{1+(x/C)^B} + D $$ $A$: 上渐近线, $D$: 下渐近线, $C$: IC₅₀值

---

通过Python实现的统计方法，研究人员能够高效处理高通量生物数据，揭示基因-表型关联、进化关系和疾病机制等关键生物学规律。