来源：pAo: A Consensus Genome-Scale Metabolic Model for Aspergillus oryzae Capturing Intraspecies Diversity.（Comput Struct Biotechnol J）| 编译：DNA Lab Space
导读
米曲霉（Aspergillus oryzae）是酱油、清酒、味噌等传统发酵食品的核心菌种，也是现代生物技术中的重要细胞工厂。然而，此前唯一的米曲霉基因组规模代谢模型（GEM）发表于2008年，且仅基于单一菌株RIB40。本研究整合187株米曲霉的基因组数据，构建了泛基因组规模代谢模型pAo，包含2,025个反应，较旧模型代谢覆盖度提升52%，并新发现了细胞色素P450介导的外源物代谢等此前被忽视的途径。基于pAo推导出187个菌株特异模型，并通过285种底物的高通量表型筛选验证了其中8个模型。研究为米曲霉的菌株理性筛选和代谢工程改造奠定了重要资源基础。
研究背景
米曲霉是传统食品发酵的基石微生物，长期安全使用的历史使其获得美国FDA和WHO的GRAS（Generally Recognized as Safe）地位。这一地位推动了其向现代生物技术的拓展——米曲霉以高分泌淀粉酶、蛋白酶等天然酶类著称，具备优异的重组蛋白生产能力，广泛用于苹果酸、曲酸、葡糖酸等有机酸的生产。近年来，以米曲霉菌丝体为原料的可持续肉类替代品等新型食品应用也相继涌现，其作为益生元和益生菌的潜力亦得到证实。
这种多功能性使米曲霉成为食品科学与工业生物技术的核心微生物。然而，要进一步优化其生物技术应用，必须深入理解调控其代谢的分子组件——基因、蛋白质和代谢途径。在这一背景下，基因组规模代谢模型（GEM）展现出强大价值。GEM是基于生物体基因组数据构建的知识库，能够系统性地描述基因-蛋白质-反应之间的关联，为代谢预测和工程改造提供计算平台。但此前的米曲霉GEM（oryzae2008）仅限于单一实验室菌株RIB40，无法反映物种层面的代谢多样性，亟需一个涵盖种内多样性的新一代模型。
研究方法
菌株选择
本研究所用的综合基因组数据集详见van der Luijt等人（待发表）的报道。简言之，我们将新测序的米曲霉NRRL菌株和A1560的基因组组装与NCBI Genome数据库[81]及覆盖各进化枝的米曲霉综合基因组数据库[82]中的公开基因组相结合。此外，我们还纳入了从NCBI Sequence Read Archive[83]获得的原始测序数据经de novo组装得到的基因组。所有基因组均使用BUSCO（v5.8.2）[84,85]评估完整性，基因预测使用AUGUSTUS（v3.4.0）[86]完成。为确保数据质量，仅保留BUSCO完整性评分>94%的基因组。进化枝分配同样详见van der Luijt等人（待发表）的报道。
【数据】菌株数：187株米曲霉+3种曲霉属物种；BUSCO阈值：>94%；软件：BUSCO v5.8.2、AUGUSTUS v3.4.0
构建米曲霉泛GEM
建立GPR（基因-蛋白质-反应规则）以将每个蛋白质序列与其催化的代谢反应相关联。这一过程通过直系同源关系将新蛋白质序列与已知GPR的蛋白质序列进行关联。具体而言，我们使用了3个经过精心整理的GEM的蛋白质序列：米曲霉RIB40的oryzae2008 GEM[49]，以及近缘物种黑曲霉（A. niger）[58]和烟曲霉（A. fumigatus）[46]的2个GEM。
为将同源序列分组为直系同源基因簇（确保重复或高度相似的基因以单一基因而非多次计数的方式表示），使用BPGA[87]进行蛋白质序列聚类，该工具调用USEARCH算法[88]中的clust_fast函数。序列以0.5的同一性阈值分组为直系同源基因簇。重要的是，这一聚类步骤本身并不确定基因功能；GPR随后通过与既有GEM的" guilt by association"原则分配给新序列。
选择0.5阈值的原因在于：(a) 该阈值足够宽松，能将较大比例的新序列与3个既有GEM中的序列关联起来；(b) 同时足够严格，确保既有GEM中的大多数序列能够分开聚类。例如，oryzae2008中注释了GPR的1,346个蛋白质中，使用这些设置鉴定出1,160个序列簇。该方法总共鉴定出11,949个蛋白质序列簇。
对于泛GEM，我们采用与oryzae2008 GEM相同的4个细胞区室：细胞质、线粒体、过氧化物酶体和胞外空间。相比之下，黑曲霉和烟曲霉GEM还将反应分配至其他区室，包括细胞核、高尔基体、内质网、液泡和脂质颗粒。为确保模板模型间的一致性，来自这些区室的反应在我们的泛GEM中均被分配至细胞质。区室化在GEM中很重要，但对表型预测影响最大的反应属于能量代谢和（支链）氨基酸代谢[89]，这两者均定位于保留的线粒体中。被合并的细胞器以脂质为主，且已通过代表囊泡穿梭的转运反应相连，因此合并这些少数区室不应对预测产生不利影响。此外，我们注意到3个既有GEM使用不同的代谢物命名体系。为解决这一问题，我们使用KEGG[60]、ChEBI[90]和PubChem[91]数据库整合了原始代谢物名称。代谢物标识符优先标准化为KEGG标识符，若无法获得KEGG标识符则使用ChEBI标识符，若两者均未找到则保留代谢物名称。
对于在3个既有GEM中均无直系同源物的序列，我们使用RAVEN工具箱[92]的getKEGGModelForOrganism函数对KEGG[60]反应数据库进行BLAST[59]同源搜索。该函数旨在基于蛋白质与KEGG中收录的所有真核生物直系同源物的同源性来重建GEM（KEGG release 105），使我们能够为缺乏直系同源物的序列获取GPR。随后，我们移除getKEGGModelForOrganism引入的所有额外反应，以报告一个独立的GEM。剩余反应的细胞定位使用DeepLoc[93]预测。最后，由于KEGG反应数据库通常缺乏准确的反应可逆性信息，我们使用在线equilibrator工具[94]（https://equilibrator.weizmann.ac.il）更新反应的可逆性。无法通过上述任何方法建立GPR的序列被排除在泛GEM之外。
生物量方程直接采用oryzae2008（表S5）[49]。根据Fritzemeier等人[95]描述的方法，通过测试以下能量载体的耗散反应是否为零来评估错误能量产生循环的存在：ATP、三磷酸胞苷、三磷酸鸟苷、三磷酸尿苷、三磷酸肌苷、还原型烟酰胺腺嘌呤二核苷酸（NADH）、还原型烟酰胺腺嘌呤二核苷酸磷酸（NADPH）、还原型黄素腺嘌呤二核苷酸（FADH2）、还原型黄素单核苷酸（FMNH2）、乙酰辅酶A、L-谷氨酸和泛醌-8。泛GEM中最初存在的3个此类循环（源自其他曲霉属物种的反应）通过将4个反应（OtherAsp_R04962、OtherAsp_R01708、r1736和r1737）约束至其热力学可行方向而消除。
模型质量通过memote[96]评估。共识pAo模型总体得分为57.7%，所有反应均电荷平衡，约80%质量平衡（受限于约430个含通用/聚合伪代谢物的反应，这些反应无法承载唯一分子式）。完整报告见https://doi.org/10.5281/zenodo.18598684。
【数据】蛋白质序列簇总数：11,949；oryzae2008中GPR注释蛋白：1,346；聚类阈值：0.5；区室数：4；memote总体得分：57.7%；电荷平衡：100%；质量平衡：约80%
菌株特异GEM的模型相似性分析
我们使用RAVEN工具箱[92]的compareMultipleModels函数在反应水平上比较菌株特异GEM。简言之，基于反应存在与否构建二元反应存在矩阵，据此计算菌株间的两两Hamming距离。这些距离作为t分布随机邻域嵌入（t-SNE）的输入，以在降维空间中可视化菌株特异GEM。
【数据】分析方法：compareMultipleModels函数、Hamming距离、t-SNE降维
菌株与培养基
8株米曲霉（RIB40、NRRL 471、NRRL 2217、NRRL 3483、NRRL 3488、NRRL 5589、NRRL 5592和NRRL 35890）用于96孔和24孔板中的表型表征。
使用马铃薯葡萄糖琼脂（葡萄糖[20 g/l]、马铃薯浸出液[4 g/l]和琼脂[15 g/l]）平板和斜面，从保存于20%甘油、−70°C的孢子悬液冷冻保藏管中初步复苏菌株并收获分生孢子。
对于PM（表型微阵列）测定，使用丝状真菌接种液（Biolog Inc.，Hayward，CA，USA）作为基础培养基，根据PM板类型进行不同补充：PM1和PM2A板无补充物，PM3B板补充100 mM D-葡萄糖、5 mM无水磷酸二氢钾和2 mM硫酸钠。
对于24孔板生长研究，制备Czapek–Dox培养基（NaNO3 [3 g/l]、KCl [2 g/l]、KH2PO4 [1 g/l]、MgSO4·7H2O [0.5 g/l]和FeSO4 [0.02 g/l]，pH 5.5）。该培养基补充4种碳源（葡萄糖、麦芽糖、甘油或木糖）之一（20 g/l）、0.01%（v/v）Tween 80作为表面活性剂，以及氯霉素[25 μg/ml]以防止细菌污染。
【数据】菌株数：8株；保存条件：20%甘油、−70°C；PDA培养基：葡萄糖20 g/l、马铃薯浸出液4 g/l、琼脂15 g/l；PM3B补充物：100 mM葡萄糖、5 mM KH2PO4、2 mM Na2SO4；Czapek–Dox培养基：NaNO3 3 g/l、KCl 2 g/l、KH2PO4 1 g/l、MgSO4·7H2O 0.5 g/l、FeSO4 0.02 g/l、pH 5.5；碳源浓度：20 g/l；Tween 80：0.01%（v/v）；氯霉素：25 μg/ml
表型微阵列
PM实验使用PM1、PM2A和PM3B板（Biolog Inc.）进行。每个孔接种100 μl新鲜配制的分生孢子悬液（以丝状真菌接种液配制并含适当补充物），使用Opentrons OT-2自动移液工作站（Opentrons Labworks Inc.，Long Island City，NY，USA）分装。向每孔加入氧化还原染料D（Biolog Inc.）——一种基于四唑的染料，通过NAD(P)H依赖性还原指示代谢活性。板在30°C下于OmniLog多板读取仪中孵育。
【数据】接种体积：100 μl/孔；染料：氧化还原染料D（四唑类）；孵育温度：30°C；仪器：OmniLog多板读取仪
研究结果
以泛GEM绘制米曲霉的代谢潜力
我们构建了米曲霉泛GEM——pAo，以捕获该物种的代谢潜力。该模型在先前仅局限于单一实验室菌株RIB40的米曲霉GEM（oryzae2008）基础上大幅扩展[49]。我们的重建流程整合了187株米曲霉的基因组数据（图1A和表S1）。我们预测蛋白质序列，并通过与3个经过整理的曲霉属GEM（来自米曲霉RIB40、黑曲霉和烟曲霉）[46,49,58]的同源搜索推断其代谢功能。
该分析鉴定出11,949个直系同源蛋白质序列簇，分类为严格核心（存在于所有187株中）、软核心（>95%）、外壳（5%至95%）、云（<5%但>1株）或菌株特异（图S1），分别包含6,959、1,882、1,687、611和810个簇。对于在既有GEM中具有直系同源物的簇，我们分配相应的代谢反应。对于其余簇，我们对京都基因与基因组百科全书（KEGG）[60]反应数据库进行BLAST[59]搜索，以鉴定潜在的新代谢功能。
所得模型pAo包含2,025个反应，其中1,426个来自oryzae2008，86个来自其他曲霉属GEM，454个来自KEGG（图1B）。另有59个反应基于实验验证作为缺口填充纳入，以确保在特定底物上的生长（表S2）。这一经过整理的网络由1,397个基因催化，涵盖1,824种代谢物，作为pAo模型，在oryzae2008模型基础上大幅扩展（图1C）[49]。最显著的是，生化反应数量（不包括转运和交换反应）从oryzae2008中的1,069个增加到pAo中的1,652个。
为评估pAo相对于先前米曲霉GEM的覆盖度，我们检查了oryzae2008基因在pAo核心中的代表程度。oryzae2008的1,160个基因中，69个在pAo中为非核心（既非严格核心也非软核心）；这些基因在模板模型中催化124个反应，其中30个与柄曲霉素/黄曲霉毒素生物合成相关——该途径已知具有菌株变异性，因此不应构成pAo核心的一部分[3,61,62]。此外，尽管pAo是oryzae2008的超集，oryzae2008缺少构成pAo核心的49个基因，这些基因共同催化124个反应，其中26个与脂肪酸生物合成相关（图S2）。
为便于下游途径分析，每个反应根据KEGG数据库[60]的子系统信息分配单一途径注释。将这些注释与oryzae2008比较，揭示了pAo中在早期模型中可能被忽视的额外代谢途径（图1D）。最显著的扩展在脂肪酸生物合成（包含38个反应，由4个基因编码）和细胞色素P450介导的外源物（异源化合物）代谢——P450是含血红素的单加氧酶超家族[63]——包含34个新反应，由7个基因催化。尽管这类酶在米曲霉中已有充分文献记载[64]，但oryzae2008中完全缺失。我们还观察到4条参与氨基酸代谢的途径。虽然构成氨基酸代谢核心的所有酶已包含在oryzae2008中，但我们的方法鉴定出54个额外基因扩展了米曲霉的氨基酸代谢，例如开辟了将精氨酸代谢导向丁酸酯或乙醛酸代谢的额外途径（图S3）。
【数据】直系同源簇总数：11,949（严格核心6,959、软核心1,882、外壳1,687、云611、菌株特异810）；pAo反应总数：2,025（来自oryzae2008：1,426；其他曲霉GEM：86；KEGG：454；缺口填充：59）；基因数：1,397；代谢物数：1,824；生化反应数：1,652（oryzae2008为1,069）；脂肪酸生物合成新增：38反应/4基因；P450介导外源物代谢新增：34反应/7基因；氨基酸代谢新增基因：54
菌株特异基因组规模模型揭示显著的种内代谢多样性
我们使用泛GEM为全部187株米曲霉构建了菌株特异GEM。每个模型代表泛GEM中该菌株所含基因对应的反应子集。自发反应和遗传基础未知的反应包含在所有模型中。分析揭示了高度保守的代谢核心，1,152个基因（82.5%）和1,118个GPR相关反应（71.1%）存在于超过95%的菌株中（图2A）。尽管如此，我们观察到显著的菌株特异（90个GPR相关反应）代谢内容。基于各菌株特异GEM中存在的反应进行降维表示，显示米曲霉菌株主要按其分类学进化枝聚类（图2B），表明进化关系反映在代谢网络拓扑中。这一结论得到进化枝内反应内容相似性高于进化枝间相似性以及反应内容的成对Jaccard相似性的支持（图S4和S5）。
【数据】核心基因：1,152个（82.5%）；核心GPR相关反应：1,118个（71.1%）；菌株特异GPR相关反应：90个
高通量表型分型驱动模型校正并揭示一致的底物利用模式
为评估表型多样性，我们对8株系统发育多样的米曲霉菌株在285种底物（190种碳源和95种氮源）上进行了表型微阵列（PM）高通量表型分型。PM上样本的聚类表明生物学重复聚集在一起，证实菌株间差异超过实验变异（图3和图S6）。对所有285种底物的聚类（图S6）还揭示了几个连贯的底物簇，尤其是在氮源中。例如，一个簇包含精氨酸代谢中的关键代谢物（腐胺、胍丁胺和鸟氨酸），另一个组包含大多数测试的二肽。
【数据】测试底物：285种（190种碳源+95种氮源）；菌株数：8株；生物学重复：3个
不同碳源上的生长比较
对8株米曲霉在4种碳源（葡萄糖、麦芽糖、甘油和木糖）上的生长进行了比较。生物量测定（图4A）和最大比生长速率测定（图4B）的结果以各菌株获得的最大值归一化，以平均值±SD（3个生物学重复）表示。生长速率的柱状图比较以任意单位（AU）表示（图5）。
【数据】碳源：葡萄糖、麦芽糖、甘油、木糖；生物学重复：3个；结果表示：平均值±SD，归一化至各菌株最大值
4种碳源下8个菌株特异GEM的通量比较分析
对8个菌株特异GEM在4种碳源下进行了比较通量分析（图6）。通量平衡分析（FBA）预测的最大生长速率以碳摄取归一化至等碳摩尔供应量（图6A）。简约FBA（pFBA）成本以单位生长的总通量表示（图6B）。单反应敲除后保留的生长显示核心反应模拟在RIB40 GEM中进行，附属反应模拟在相关菌株特异GEM中进行（图6C）。相对于葡萄糖的子系统通量份额变化显示受影响最大的7个子系统，每个点代表一个菌株，按底物着色（图6D）。
【数据】碳源：4种；菌株：8株；FBA预测最大生长速率：碳摄取归一化至等碳摩尔供应量；pFBA成本：总通量/单位生长
讨论与解读
本研究为推进米曲霉这一工业重要霉菌的代谢建模提供了3项关键资源：第一，涵盖187株米曲霉全部代谢潜力的综合泛GEM（pAo）；第二，由此推导的187个菌株特异GEM集合；第三，8个经高通量PM数据验证的高精度模型，覆盖系统发育多样的菌株。pAo模型作为最新一代的代谢目录，捕获的反应和代谢物数量远超先前的oryzae2008模型。其核心用途在于作为模板——通过自动化缺口填充和基于同源搜索的快速新菌株注释，简化菌株特异模型的生成。此外，泛GEM可用于评估米曲霉针对特定生物技术应用的整体代谢能力，进而识别最适菌株。
为验证方法，我们选择了来自不同系统发育进化枝的8个菌株特异GEM进行实验校正。这些菌株覆盖米曲霉的主要系统发育多样性，而本研究显示这种多样性反映在代谢网络的显著差异中。生长实验表明，尽管网络结构层面的通量平衡分析显示变异源于调控层面而非网络结构层面，不同菌株间仍存在显著的代谢多样性。
编译者解读
从合成生物学应用角度看，pAo模型的价值在于将"菌株筛选"从经验试错提升为计算驱动的理性选择。研究团队以187株基因组为基底构建泛模型，再反推菌株特异模型，这一"先泛后特"的策略既保证了代谢覆盖度，又保留了菌株间差异的可解析性。尤其值得肯定的是，他们将计算预测与285种底物的高通量表型数据闭环互校——既用实验数据修正模型缺口（59个缺口填充反应），又用模型预测指导菌株选择。局限也很明显：memote评分仅57.7%，约20%反应未达质量平衡；P450介导的外源物代谢虽被纳入，但其底物特异性在模型中仍较粗糙。未来若能将调控层信息（如转录因子网络）整合进泛模型，有望进一步解释"网络结构相同而表型不同"的调控根源。
参考来源
Gilis J, van der Luijt CRB, Feller M, Barba CSG, Sommer MOA. pAo: A Consensus Genome-Scale Metabolic Model for Aspergillus oryzae Capturing Intraspecies Diversity. Comput Struct Biotechnol J. PMCID: PMC13486731.
PMCID: PMC13486731