name: cobrapy description: “约束基代谢建模(COBRA)。FBA、FVA、基因敲除、通量采样、SBML模型,用于系统生物学和代谢工程分析。”
COBRApy - 约束基重构与分析
概述
COBRApy 是一个用于约束基重构与分析(COBRA)的Python库,对系统生物学研究至关重要。使用基因组规模代谢模型,进行细胞代谢的计算仿真,开展代谢工程分析,并预测表型行为。
核心能力
COBRApy 提供了全面工具,组织成几个关键领域:
1. 模型管理
加载现有模型从仓库或文件:
from cobra.io import load_model
# 加载捆绑的测试模型
model = load_model("textbook") # E. coli 核心模型
model = load_model("ecoli") # 完整 E. coli 模型
model = load_model("salmonella")
# 从文件加载
from cobra.io import read_sbml_model, load_json_model, load_yaml_model
model = read_sbml_model("path/to/model.xml")
model = load_json_model("path/to/model.json")
model = load_yaml_model("path/to/model.yml")
以多种格式保存模型:
from cobra.io import write_sbml_model, save_json_model, save_yaml_model
write_sbml_model(model, "output.xml") # 首选格式
save_json_model(model, "output.json") # 用于 Escher 兼容性
save_yaml_model(model, "output.yml") # 人类可读
2. 模型结构与组件
访问和检查模型组件:
# 访问组件
model.reactions # 所有反应的 DictList
model.metabolites # 所有代谢物的 DictList
model.genes # 所有基因的 DictList
# 通过 ID 或索引获取特定项
reaction = model.reactions.get_by_id("PFK")
metabolite = model.metabolites[0]
# 检查属性
print(reaction.reaction) # 化学计量方程
print(reaction.bounds) # 通量约束
print(reaction.gene_reaction_rule) # GPR 逻辑
print(metabolite.formula) # 化学式
print(metabolite.compartment) # 细胞位置
3. 通量平衡分析 (FBA)
执行标准 FBA 仿真:
# 基本优化
solution = model.optimize()
print(f"目标值: {solution.objective_value}")
print(f"状态: {solution.status}")
# 访问通量
print(solution.fluxes["PFK"])
print(solution.fluxes.head())
# 快速优化(仅目标值)
objective_value = model.slim_optimize()
# 更改目标
model.objective = "ATPM"
solution = model.optimize()
节俭 FBA(最小化总通量):
from cobra.flux_analysis import pfba
solution = pfba(model)
几何 FBA(找中心解):
from cobra.flux_analysis import geometric_fba
solution = geometric_fba(model)
4. 通量变异分析 (FVA)
确定所有反应的通量范围:
from cobra.flux_analysis import flux_variability_analysis
# 标准 FVA
fva_result = flux_variability_analysis(model)
# 在 90% 最优性下的 FVA
fva_result = flux_variability_analysis(model, fraction_of_optimum=0.9)
# 无循环 FVA(消除热力学不可行循环)
fva_result = flux_variability_analysis(model, loopless=True)
# 针对特定反应的 FVA
fva_result = flux_variability_analysis(
model,
reaction_list=["PFK", "FBA", "PGI"]
)
5. 基因和反应敲除研究
执行敲除分析:
from cobra.flux_analysis import (
single_gene_deletion,
single_reaction_deletion,
double_gene_deletion,
double_reaction_deletion
)
# 单敲除
gene_results = single_gene_deletion(model)
reaction_results = single_reaction_deletion(model)
# 双敲除(使用多进程)
double_gene_results = double_gene_deletion(
model,
processes=4 # CPU 核心数
)
# 使用上下文管理器手动敲除
with model:
model.genes.get_by_id("b0008").knock_out()
solution = model.optimize()
print(f"敲除后生长率: {solution.objective_value}")
# 上下文退出后模型自动恢复
6. 生长培养基与最小培养基
管理生长培养基:
# 查看当前培养基
print(model.medium)
# 修改培养基(必须重新分配整个字典)
medium = model.medium
medium["EX_glc__D_e"] = 10.0 # 设置葡萄糖摄取
medium["EX_o2_e"] = 0.0 # 厌氧条件
model.medium = medium
# 计算最小培养基
from cobra.medium import minimal_medium
# 最小化总导入通量
min_medium = minimal_medium(model, minimize_components=False)
# 最小化组分数量(使用 MILP,较慢)
min_medium = minimal_medium(
model,
minimize_components=True,
open_exchanges=True
)
7. 通量采样
采样可行通量空间:
from cobra.sampling import sample
# 使用 OptGP 采样(默认,支持并行处理)
samples = sample(model, n=1000, method="optgp", processes=4)
# 使用 ACHR 采样
samples = sample(model, n=1000, method="achr")
# 验证采样
from cobra.sampling import OptGPSampler
sampler = OptGPSampler(model, processes=4)
sampler.sample(1000)
validation = sampler.validate(sampler.samples)
print(validation.value_counts()) # 应全为 'v' 表示有效
8. 生产包络
计算表型相平面:
from cobra.flux_analysis import production_envelope
# 标准生产包络
envelope = production_envelope(
model,
reactions=["EX_glc__D_e", "EX_o2_e"],
objective="EX_ac_e" # 醋酸生产
)
# 带碳产率
envelope = production_envelope(
model,
reactions=["EX_glc__D_e", "EX_o2_e"],
carbon_sources="EX_glc__D_e"
)
# 可视化(使用 matplotlib 或 pandas 绘图)
import matplotlib.pyplot as plt
envelope.plot(x="EX_glc__D_e", y="EX_o2_e", kind="scatter")
plt.show()
9. 缺口填充
添加反应以使模型可行:
from cobra.flux_analysis import gapfill
# 准备含候选反应的通用模型
universal = load_model("universal")
# 执行缺口填充
with model:
# 移除反应以创建缺口用于演示
model.remove_reactions([model.reactions.PGI])
# 查找需添加的反应
solution = gapfill(model, universal)
print(f"需添加的反应: {solution}")
10. 模型构建
从头构建模型:
from cobra import Model, Reaction, Metabolite
# 创建模型
model = Model("my_model")
# 创建代谢物
atp_c = Metabolite("atp_c", formula="C10H12N5O13P3",
name="ATP", compartment="c")
adp_c = Metabolite("adp_c", formula="C10H12N5O10P2",
name="ADP", compartment="c")
pi_c = Metabolite("pi_c", formula="HO4P",
name="磷酸盐", compartment="c")
# 创建反应
reaction = Reaction("ATPASE")
reaction.name = "ATP 水解"
reaction.subsystem = "能量"
reaction.lower_bound = 0.0
reaction.upper_bound = 1000.0
# 添加带化学计量的代谢物
reaction.add_metabolites({
atp_c: -1.0,
adp_c: 1.0,
pi_c: 1.0
})
# 添加基因-反应规则
reaction.gene_reaction_rule = "(gene1 and gene2) or gene3"
# 添加到模型
model.add_reactions([reaction])
# 添加边界反应
model.add_boundary(atp_c, type="exchange")
model.add_boundary(adp_c, type="demand")
# 设置目标
model.objective = "ATPASE"
常见工作流
工作流1:加载模型并预测生长
from cobra.io import load_model
# 加载模型
model = load_model("ecoli")
# 运行 FBA
solution = model.optimize()
print(f"生长率: {solution.objective_value:.3f} /h")
# 显示活跃通路
print(solution.fluxes[solution.fluxes.abs() > 1e-6])
工作流2:基因敲除筛选
from cobra.io import load_model
from cobra.flux_analysis import single_gene_deletion
# 加载模型
model = load_model("ecoli")
# 执行单基因敲除
results = single_gene_deletion(model)
# 查找必需基因(生长率 < 阈值)
essential_genes = results[results["growth"] < 0.01]
print(f"找到 {len(essential_genes)} 个必需基因")
# 查找影响最小的基因
neutral_genes = results[results["growth"] > 0.9 * solution.objective_value]
工作流3:培养基优化
from cobra.io import load_model
from cobra.medium import minimal_medium
# 加载模型
model = load_model("ecoli")
# 计算最大生长率 50% 的最小培养基
target_growth = model.slim_optimize() * 0.5
min_medium = minimal_medium(
model,
target_growth,
minimize_components=True
)
print(f"最小培养基组分数: {len(min_medium)}")
print(min_medium)
工作流4:通量不确定性分析
from cobra.io import load_model
from cobra.flux_analysis import flux_variability_analysis
from cobra.sampling import sample
# 加载模型
model = load_model("ecoli")
# 首先检查最优性下的通量范围
fva = flux_variability_analysis(model, fraction_of_optimum=1.0)
# 对于范围大的反应,采样以理解分布
samples = sample(model, n=1000)
# 分析特定反应
reaction_id = "PFK"
import matplotlib.pyplot as plt
samples[reaction_id].hist(bins=50)
plt.xlabel(f"通过 {reaction_id} 的通量")
plt.ylabel("频率")
plt.show()
工作流5:使用上下文管理器进行临时更改
使用上下文管理器进行临时修改:
# 在上下文外模型保持不变
with model:
# 临时更改目标
model.objective = "ATPM"
# 临时修改边界
model.reactions.EX_glc__D_e.lower_bound = -5.0
# 临时敲除基因
model.genes.b0008.knock_out()
# 优化带更改
solution = model.optimize()
print(f"修改后生长率: {solution.objective_value}")
# 所有更改自动恢复
solution = model.optimize()
print(f"原始生长率: {solution.objective_value}")
关键概念
DictList 对象
模型使用 DictList 对象表示反应、代谢物和基因 - 同时像列表和字典一样行为:
# 通过索引访问
first_reaction = model.reactions[0]
# 通过 ID 访问
pfk = model.reactions.get_by_id("PFK")
# 查询方法
atp_reactions = model.reactions.query("atp")
通量约束
反应边界定义可行通量范围:
- 不可逆:
lower_bound = 0, upper_bound > 0 - 可逆:
lower_bound < 0, upper_bound > 0 - 使用
.bounds同时设置两个边界以避免不一致
基因-反应规则 (GPR)
链接基因到反应的布尔逻辑:
# AND 逻辑(两者都需)
reaction.gene_reaction_rule = "gene1 and gene2"
# OR 逻辑(任一足够)
reaction.gene_reaction_rule = "gene1 or gene2"
# 复杂逻辑
reaction.gene_reaction_rule = "(gene1 and gene2) or (gene3 and gene4)"
交换反应
表示代谢物导入/导出的特殊反应:
- 命名约定带前缀
EX_ - 正通量 = 分泌,负通量 = 摄取
- 通过
model.medium字典管理
最佳实践
- 使用上下文管理器 进行临时修改以避免状态管理问题
- 验证模型 使用
model.slim_optimize()确保可行性 - 检查解状态 优化后 -
optimal表示成功求解 - 使用无循环 FVA 当热力学可行性重要时
- 适当设置 fraction_of_optimum 在 FVA 中探索次优空间
- 并行化 计算密集型操作(采样、双敲除)
- 首选 SBML 格式 用于模型交换和长期存储
- 使用 slim_optimize() 当仅需目标值时用于性能
- 验证通量采样 以确保数值稳定性
故障排除
不可行解:检查培养基约束、反应边界和模型一致性
优化慢:尝试不同求解器(GLPK、CPLEX、Gurobi)通过 model.solver
无界解:验证交换反应有适当的上限边界
导入错误:确保正确的文件格式和有效的 SBML 标识符
参考资料
详细工作流和 API 模式,参考:
references/workflows.md- 全面逐步工作流示例references/api_quick_reference.md- 常见函数签名和模式