
==== Front
NPJ Syst Biol Appl
NPJ Syst Biol Appl
NPJ Systems Biology and Applications
2056-7189
Nature Publishing Group UK London

39289347
426
10.1038/s41540-024-00426-5
Article
Understanding flux switching in metabolic networks through an analysis of synthetic lethals
http://orcid.org/0000-0002-1312-2985
Narasimha Sowmya Manojna 125
Malpani Tanisha 12
Mohite Omkar S. 126
Nath J. Saketha 3
http://orcid.org/0000-0002-9311-7093
Raman Karthik kraman@iitm.ac.in

124
1 https://ror.org/03v0r5n49 grid.417969.4 0000 0001 2315 1926 Centre for Integrative Biology and Systems mEdicine (IBSE), Indian Institute of Technology (IIT) Madras, Chennai, 600 036 India
2 grid.417969.4 0000 0001 2315 1926 Department of Biotechnology, Bhupat Jyoti Mehta School of Biosciences, Indian Institute of Technology (IIT) Madras, Chennai, 600 036 India
3 https://ror.org/01j4v3x97 grid.459612.d 0000 0004 1767 065X Department of Computer Science and Engineering, Indian Institute of Technology (IIT) Hyderabad, Hyderabad, 502 284 India
4 https://ror.org/03v0r5n49 grid.417969.4 0000 0001 2315 1926 Department of Data Science and AI, Wadhwani School of Data Science and AI (WSAI), Indian Institute of Technology (IIT) Madras, Chennai, 600 036 India
5 https://ror.org/0168r3w48 grid.266100.3 0000 0001 2107 4242 Present Address: Neuroscience Graduate Program, University of California San Diego, San Diego, CA 92092 USA
6 grid.5170.3 0000 0001 2181 8870 Present Address: Novo Nordisk Foundation Center for Biosustainability, Technical University of Denmark, 2800 Kgs. Lyngby, Denmark
17 9 2024
17 9 2024
2024
10 10419 3 2024
17 8 2024
© The Author(s) 2024
2024
https://creativecommons.org/licenses/by/4.0/ Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
Biological systems are robust and redundant. The redundancy can manifest as alternative metabolic pathways. Synthetic double lethals are pairs of reactions that, when deleted simultaneously, abrogate cell growth. However, removing one reaction allows the rerouting of metabolites through alternative pathways. Little is known about these hidden linkages between pathways. Understanding them in the context of pathogens is useful for therapeutic innovations. We propose a constraint-based optimisation approach to identify inter-dependencies between metabolic pathways. It minimises rerouting between two reaction deletions, corresponding to a synthetic lethal pair, and outputs the set of reactions vital for metabolic rewiring, known as the synthetic lethal cluster. We depict the results for different pathogens and show that the reactions span across metabolic modules, illustrating the complexity of metabolism. Finally, we demonstrate how the two classes of synthetic lethals play a role in metabolic networks and influence the different properties of a synthetic lethal cluster.

Subject terms

Metabolic engineering
Metabolic engineering
issue-copyright-statement© Springer Nature Limited 2024
==== Body
pmcIntroduction

Robustness in the face of environmental perturbations is an essential attribute of microorganisms1–3. This robustness is often achieved by the presence of multiple alternate pathways that achieve similar metabolic functions4,5. The redundancy introduced by alternate pathways comprises a large fraction of most metabolic networks6–8 and shows surprising variance in their distribution9,10. While some of the alternate pathways are very simple, arising due to gene duplication, other alternate pathways could be extremely complex, with compensating reactions spanning different metabolic subsystems9,10.

A straightforward method of studying these alternate pathways involves the identification of synthetic lethals in a metabolic model10. Synthetic lethals are sets of genes/reactions where only the simultaneous loss of all genes/reactions in the set leads to abrogation of cell growth11. When only one of the reactions is deleted, the cell is able to summon alternate pathways to ensure its survival. In many cases, this is made possible through a complex rerouting of fluxes in the metabolic network which exploits the redundancy in metabolism. However, very little is known about how these organisms reroute their fluxes, and how various reactions in the cell can compensate for one another.

Previous methods such as GIMME, iMAT, and RELATCH, for the study of flux distributions, have focused on a given condition, the final steady state of the cell, without considering the prior reference state of the organism12–14. REMI and deltaFBA are algorithms that integrate differential expression of transcriptome and metabolome with the flux distributions between two different states, a wild-type state and a mutant state15,16. These algorithms require gene expression data, which are not always available. Considering cells have a high order of redundancy and synthetic lethal, we require a method that can computationally predict the rewiring of metabolism. FBA alone has also been used17,18 to study redundancy using synthetic lethals in E. coli and other bacteria. In refs. 17,18, FBA was used to directly optimise the single deletions for the biomass objective and find the differences in the flux distributions. However, this can ignore the biological costs associated with altering flux in an organism, which may result in sub-optimal biomass production.

The resistance to antibiotics offered by the flexibility of metabolic rewiring has been shown to come with a fitness cost19. The advantages of alternate optimal solutions for various context-specific constraint-based modelling have also been illustrated before20. In anticipation of adapting to changing environments, bacteria opt for diauxic growth metabolism, which involves replicating at a sub-optimal growth rate so they can invest resources in new metabolic processes21. This sub-optimal growth is also seen in response to perturbations as the metabolic flux is rerouted in small adjustments first and then adaptively mutates to optimise for growth22. For the bacterium B. subtilis23, studies have shown that the relative change in flux distribution remains the same in mutant knock-outs even if absolute flux changes significantly. Thus, MOMA24, a sub-optimal growth rate-based FBA approach, uses a quadratic minimisation between the mutant and wild-type state to predict a sub-optimal growth rate.

FBA presents a challenge in the form of handling multiple flux solutions (even if there is a unique growth rate). To work around this and to identify this set of reactions that come into effect to rescue the cell from a non-lethal deletion, we propose a novel approach termed minRerouting. By solving a minimum p-norm problem, minRerouting can simultaneously solve for flux distributions that satisfy the stoichiometric constraints, maximise the biomass objective, and also minimise the number of reactions with varying metabolic flux values. We build on the wide body of evidence that supports flux balance-based predictions and further enhance the algorithm to account for the multiple solutions possible in any given conditions. The study further helps us understand the redundancies that facilitate robustness in wild-type organisms against perturbations and mutations. In addition, this approach is ideal for exploring and understanding indispensable changes in cellular metabolism between two different conditions, for instance, a healthy state and a diseased state, and has applications in the biotechnology and pharmaceutical industries.

Robustness arising from double lethal has been studied previously. It has been proposed that double lethal pair robustness in an organism can be ascribed to two classes of reaction pairs—plastic synthetic lethal (PSL) and redundant synthetic lethal (RSL)18. PSL pairs are reaction pairs where only one reaction is active, while the other reaction is inactive. The second reaction becomes active only when the first reaction is inactive. RSL pairs are reaction pairs where both the reactions are active simultaneously, yet, the loss of one does not abrogate growth. It has also been shown that these classes are conserved even across different nutrient conditions. The very presence of two distinct reaction pair classes calls for us to analyse the cause behind such selective activation. Are the inactive reactions more “metabolically costly” than the active ones? What kind of reactions make up the RSL pairs, especially when they are both simultaneously active? We employ a novel workflow to answer these questions and explore the structure of metabolic networks and their underlying redundancy. We also analyse the reaction types contributing to the PSL and RSL classes and uncover interesting patterns in their distribution.

Results

The minRerouting pipeline has been carried out on eight genome-scale metabolic models. We have chosen organisms that represent key bacterial pathogens relevant to humans and are present in the BiGG database25. For E. colil, found in the human gut, we have depicted two models, e_coli_core26 and iML151527. The model e_coli_core represents the simplified versions of only the most crucial pathways needed for its survival. The other models studied are for the bacteria Helicobacter pylori, Klebsiella pneumoniae, Mycobacterium tuberculosis, Salmonella Typhimurium, Shigella sonnnei and Yersinia pestis. These models and the number of single and double lethal reactions predicted for them using Fast-SL are listed in Table 1.Table 1 List of all the models analysed as a part of the study

GSMM	Total number of reactions	Single lethals	Double lethals	Ref.	
Escherichia coli (e_coli_core)	95	14	88	26	
Mycobacterium tuberculosis (iEK1008)	1226	351	157	68	
Helicobacter pylori 26695 (iIT341)	554	252	54	69	
Escherichia coli (iML1515)	2712	253	287	27	
Yersinia pestis (iPC815)	1961	212	189	70	
Shigella sonnei (iSSON_1240)	2693	261	267	71	
Klebsiella pneumoniae (iYL1228)	2262	199	144	72	
Salmonella enterica (STM_v1_0)	2545	330	168	73	
The third and fourth columns indicate the number of single lethals and double lethal pairs identified using Fast-SL63,64.

Functional analysis of double lethals

Reaction submodules differ across organisms and within synthetic lethal pairs

Double lethals have been identified in many species previously, such as E. coli, M. pneumoniae, S. Typhimurium, and S. sonnei using FBA17,18. Our predictions are consistent with previous results regarding the number of synthetic double lethals obtained for different species and their composition of reactions. For E. coli, our results for synthetic lethals paralleled observations made by previous studies17,18.

Table 1 has the distribution of the number of single lethals and synthetic double lethals for each model. We see that the fraction of single lethal or essential reactions as compared to the total reactions is similar across models except for models of two species, M. tuberculosis and H. pylori, which have 29% and 45% of essential reactions. The number of double lethals also varies between the models, implying that the occurrence of a reaction as a double lethal is more nuanced.

Figure 1 shows that more than 500 double lethal pairs are present in at least one organism, indicating the high level of redundancy present across metabolic networks. Among these, few are shared across at least five of the models analysed in the study. The tendency of a reaction pair to appear in multiple species with different adaptations implies that these could be potential super targets for drug therapies. The common synthetic lethal pairs examined across the organisms are all from the pentose phosphate pathway. The second most common reactions are from the glycolytic pathway and different amino acid biosynthesis pathways. Similar to observations made by Barve et al.28, the reactions that were essential belonged to linear or anabolic pathways such as ATP and histidine synthesis, while redundancies were present only in more reticulate pathways such as the pyruvate or glycolysis metabolic pathways.Fig. 1 Distribution of common double lethal pairs across models.

500 reaction pairs are unique to specific models, while very few pairs are present across all eight models.

A broader understanding is obtained by looking at each organism’s submodule distribution of double lethal, as shown in Figs. 2 and 3 for the models iML1515 and iEK1008. While it is intuitive to think the reaction pairs would arise from the same submodule, for all models, we see that at least 50% of the reactions are from different submodules, i.e., the synthetic lethals are inter-pathway. This could be because of the ripple effect caused by deleting one reaction, which causes small changes in all other connected pathways. Inter-pathway synthetic lethals also highlight an organism’s need for cross-talk between pathways, such as energy production and nucleotide metabolism29.Fig. 2 Metabolic subsystem analysis for the model iML1515.

The distribution shows that more than half the lethal pairs are from differing submodules. Only the submodules occurring more than the mean of the distribution are depicted for clarity.

Fig. 3 Metabolic subsystem analysis for the model iEK1008.

The distribution shows lethal pairs from unique submodules, such as the mycolic acid pathway and pyruvate metabolism. Unlike iML1515, cell envelope biosynthesis does not appear in any of the top pairs. Only the submodules occurring more than the mean of the distribution are depicted for clarity.

Previously, knock-out studies of reactions in the pentose phosphate pathway and glycolysis pathway have been performed. The flux rates of reactions for a pyk mutant in E. coli were quantified30 and summarised24. The results qualitatively matched the results of the PYK2-NDPK4 synthetic lethal pair’s flux distribution for the PYK2 mutant for 16 of 17 of the reactions in the minRerouting cluster set. The reaction that had a mismatch was the reaction PCK, which, even in wild-type iML1515, was inactive. The above pathways have also been studied in strain-specific studies31,32, and the experiments show that reactions from the above pathways have similar redundancies as seen by minRerouting. Transaminases are another set of promiscuous enzymes catalysing the amino acid synthesis reactions. Their inter-dependency, which gives E. coli the ability to switch metabolic fluxes, has previously been studied33–35. These reactions for the iML1515 model include VALTA, ALATA_L, VPAMTr, ASNS1, and ASNS2. These enzymes can perform underground metabolism and could be important biotechnological targets or drug targets. For the model iEK1008, two synthetic lethal pairs, namely, DHPS2-FOLD3 and TREP6PP-TREY, present in M. tuberculosis have been illustrated and shown in Supplementary Figs. 28 and 29.

In simpler and closely related bacterial systems, such as E. coli and S. sonnei, more than 50% of synthetic lethal reactions involve cell envelope biosynthesis, where membrane lipid metabolism reactions act as their backups. These systems also exhibit reaction pairs from the cofactor and prosthetic group biosynthesis submodule. Y. pestis, K. pneumoniae, and S. Typhimurium, phylogenetically different from the above two species36, additionally had reaction pairs from submodules associated with amino acid metabolism and glycerophospholipid metabolism. Finally, the specialist species M. tuberculosis and H. pylori have distinct dominant submodules, such as mycolic acid production and haem transport, respectively. The distributions for the other six models are given in the Supplementary Information. Thus, the redundancies are organism-specific, with double lethal reactions emerging from submodules required for biomass synthesis, directly or indirectly. Some submodules, like the pentose phosphate pathway, are shared in all organisms. These findings align with those from earlier experiments10,17,18.

Tendency of reactions in forming double lethal

We next define the redundancy index (RI) for a reaction as the fraction of synthetic lethal pairs in which it occurs within a species (see subsection “Flux redistribution by synthetic lethals”). Several studies indicate that pathways lacking redundancy tend to have more critical functions for survival than those with redundant pathways29,37. However, there is evidence contradicting this notion10,38. Specifically, when a reaction is engaged in multiple pairs as a synthetic lethal and consequently possesses multiple backups, its function likely plays a significant role in ensuring the organism’s survival. Irrespective of the pathway, it has been shown that the essentiality of the metabolite determines the robustness of the related reactions, with the metabolite concentrations remaining unchanged to perturbations despite changes in incoming and outgoing reaction fluxes39. From Fig. 4, we observed that the models had a mean RI in similar ranges despite having different numbers of reactions and different adaptations. The reactions with a high RI were from different submodules for each organism, ranging from central carbon metabolism to transport and, more specifically, lipid metabolism for the bacterium M. tuberculosis.Fig. 4 Mean redundancy index and reaction compensation index of reactions forming synthetic lethal pairs, across organisms.

Reactions from e_coli_core have a higher RCI and RI compared to other models.

A reaction is essential if it is necessary for growth in the given media40. Growth requires the production of bio-molecules, such as fatty acids, amino acids, and purine/ pyrimidine. The reactions that do not stop growth, or can be bypassed when eliminated turn out to be redundant. Thus, we saw that the reactions from submodules related to the production of the above bio-molecules occur in the single lethal list, and reactions that have redundancy in the form of double lethal were, albeit important, but from pathways that impacted biomass production indirectly. This concurs with the observations made about the reactions with a high RI. The minRerouting approach also reveals an essential characteristic of metabolic networks, i.e., the inter-dependencies of the pathways to produce biomass are revealed. These inter-dependencies are due to molecules with a high RI acting as connecting points for other pathways in a network. Suthers et al.29 investigated the topological structure of redundancy and observed comparable findings regarding the essentiality of various reaction modules.

Similarly, the tendency of a reaction to be part of a synthetic cluster, in the form of the fraction of synthetic clusters it occurs in, was calculated as mentioned in the subsection “Flux redistribution by synthetic lethals”. We have previously defined this as the reaction compensation index (RCI)10. The RCI shows similar trends as RI in terms of value, as seen in Fig. 4. However, the reactions that have a high RI do not necessarily have a high RCI. A high RCI indicates reactions that are important for all the key pathways that the synthetic clusters form, but these are not necessarily essential at least at the order of a double lethal.

To further understand how these reactions compensate for each other by rerouting flux, and ensure their viability, we further analysed the results from the minRerouting algorithm.

Effect of media/environment on single and double lethal

From previous results, we observed that flux rerouting occurs through diverse hidden pathways. To ascertain the influence of the environmental conditions on metabolic rewiring and associated costs, we changed the media concentrations in two ways. First, we changed the concentration of up to 10 carbon sources as described in the “Methods” section. Second, we looked at the effect of unblocking all exchange reactions in the model to varying degrees, i.e., −100, −500, and −1000 mmol g DW−1 h−1. In total, the results for four new models (Model_carbon, Model_exchange_1, Model_exchange_2, and Model_exchange_3) were compared to the original for all eight species and the results are given in the Supplementary Information. From Supplementary Table 3, we see that the number of single lethals decreased for all the models where exchange reactions were unblocked. The carbon sources had little to no effect on single lethals. This is intuitive since the single lethal reactions are directly linked to biomass production (amino acid synthesis and purine and pyrimidine metabolism). From Supplementary Table 4, we see that the effect of environment on double lethals was once again more varied for each species, and the number depended on whether the carbon source was present in the model and the number of exchange reactions in the model.

However, there is little overlap between the single and double lethal reactions between the different environments. As seen from Supplementary Fig. 22 for E. coli, the model where carbon sources are changed in the media has the most fraction of reactions in common with the original model. However, by adding exchange reactions to the media, the need for synthesis of amino acid and purine/ pyrimidine synthesis was removed. Hence, the submodules for single lethal in the new models consist of reactions mainly from the cofactor and prosthetic group biosynthesis or the lipopolysaccharide biosynthesis/ recycling submodules. Similarly, the need for central carbon metabolic pathways was removed, and we see from Supplementary Fig. 25 that the double lethals related to the pentose phosphate pathway and lipid metabolism disappeared. Instead, the transport reactions become redundant and form synthetic lethal pairs.

In M. tuberculosis, as seen in Supplementary Figs. 23 and 26, where the number of DLs increases, the trend remains the same with the decrease in dependency and, hence, a decrease in redundancy on the pentose phosphate pathway and the lipid metabolism pathways. A unique example of change in the double lethal is for S. enterica when the exchange reactions are switched on, and specifically, the reaction ENLIPAtex, present only in S. enterica is switched on. From Supplementary Fig. 27, we can see that this reaction enables the circulation of the metabolites from the Lipopolysaccharide Biosynthesis Recycling pathway. Thus, a major portion of the double lethal is from this pathway. In all the species, we see an increase in the double lethals related to transport pathways due to the increase in richness of the environment media. The Supplementary Information depicts the observations for three models for simplicity. The reader can repeat the analysis for the other models using the code available.

Flux redistribution analysis

For each of the p-norms, the resultant minRerouting set is analysed. The properties described in the subsection “Flux redistribution by synthetic lethals” are studied in the following sections.

Stricter optimisation constraints expectedly yield smaller minRerouting sets

The size of the minimal rerouting set, as discussed in the “Introduction” section, is the number of reactions through which flux is rerouted. minRerouting allows the user to input their preferred norm to minimise the rerouted flux since differences in the calculation of the norms result in slightly different results. The L2 norm is called the least squares error norm as it minimises the sum of the squares of the differences. As a result, it tends to be influenced by outliers which can lead to unexpected solutions. However, its solution is unique and stable. On the other hand, the L1 norm gives a sparse solution, but with the possibility of multiple solutions for the minimisation problem. The L0 norm, computationally difficult to calculate, also results in a sparse solution. The comparison of the minRerouting cluster size across three norms, shown in Supplementary Fig. 4, thus reveals that the strictest L0 norm optimisation results in the smallest SL Cluster Size, followed by the L1 and L2 norms, respectively.

Interestingly, the sizes of the clusters calculated from minRerouting for L1 norm are lesser than observations made by Massucci et al.17 as seen in Supplementary Fig. 12. In stressful environments, restructuring metabolism incurs functional and structural costs. Our approach minimises the alterations in flux required, even if it means a partial reduction in growth rate, such that the expenses are reduced. Thus, a smaller cluster size is expected. The cluster size for each double lethal pair of an organism is shown in Fig. 5, along with other properties of the cluster for easy comparison across models.Fig. 5 Cluster size, synthetic accessibility, and the net flux difference between the two metabolic states of reaction pairs.

Mostly, the pairs that have a high synthetic accessibility have a low net flux difference and small cluster size. Outliers showing different properties exist in the top left corner for all the models except e_coli_core.

Size of the common minRerouting Set indicates the efficiency of rerouting

The size of the common minimal rerouting set or synthetic lethal cluster, as discussed in the “Methods” section, is the number of common reactions through which flux is rerouted. The comparison of the common minRerouting Set size across norms is shown in Supplementary Fig. 4.

Once again, the results obtained show that the L0 norm optimisation results in the smallest Common SL Cluster Size, followed by the L1 and L2 norms, respectively. In certain models, the median Common SL Cluster Size is 0. This indicates that the L0 norm, in addition to ensuring minimum SL Cluster size, forces the two flux vectors to take completely exclusive reaction pathways, with no reaction (with modified flux) overlap. Since the organism undergoes complete rewiring, the number of distinct reactions flux is rerouted through between the two flux vectors, normalised over the total number of reactions undergoing a flux change has been termed Synthetic Accessibility. The organism has to divert its energy into activating/deactivating the reactions not common between the two reaction pair deletions. A Synthetic Accessibility of 1 would imply no common reactions between the metabolic states of a synthetic lethal pair. While the range of Synthetic Accessibility is vast, on average, at least 25% of reactions in a cluster are turned on/off when switching states. A low Synthetic Accessibility could mean that the function of the reaction is not being completely replaced but is being compensated for by a collective change in a group of reactions mediated through its pair. A high Synthetic Accessibility could mean that most of the reactions that get activated in the mutant are not needed in the wild-type state.

In literature, there are multiple schools of thought on the role of redundancy in metabolic networks. True redundancy is believed to not be possible as it is evolutionarily unstable41,42. A genuinely redundant gene coding for the redundant reaction would be lost completely through genetic drift since its mutation would incur no fitness costs to the organism. Hence, we see that synthetic accessibility indeed has a vast range. However, the clusters with a small cluster size tend to have a high Synthetic Accessibility. Thus, the backing is efficient, and reactions are not unnecessarily active. However, one group of clusters is unique, as seen in the top left corner of every model in Fig. 5. They are outliers that have a high cluster size and low synthetic accessibility. Thus, Synthetic Accessibility gives us an idea of the efficiency of the reaction pair in the metabolic network and the cost of replacing that function that the organism is willing to bear, even with a decrease in its growth rate.

Net flux difference is indicative of the extent of rerouting

The flux difference is the total flux difference between the two flux vectors representing the deletion of individual reactions that make up the double lethal. The net flux difference across norms is compared in Supplementary Fig. 5. Flux difference shows the same trend for L0, L1, and L2 norms. As seen in Supplementary Fig. 12, the net flux difference between reaction pairs obtained using our L1 norm formulation is significantly lesser than that reported by Massucci et al17. The flux difference gives an indication of the cost of maintaining the redundancy between the reaction pairs. minRerouting thus reveals the flux redistribution that the organisms can undertake to decrease this cost further. Since the net flux difference is comparable between the norms while the cluster size is more for L2, the average change in a reaction flux is lower for L2 norm than L1 and L0 norm. This may be because the L2 norm penalises differences more as it squares them to obtain the error.

The average flux difference of clusters is high in the e_coli_core model, due to its limited ability to adjust to perturbations. Thus, even though organisms do not need many reactions to exist, more reactions provide robustness to the network.

For each pair of a subset of synthetic lethals, we also generated the corresponding knockout models—modelDel1 and modelDel2. From the flux space generated by flux sampling of these two metabolic states (2000 flux distributions each) using Flux Sampling Analysis (FSA)43, we identified the flux distributions between the two that produced the minimal change in the p-norm for L0, L1, and L2-norms. This difference, generated from the flux sampling space, was compared to our observed flux difference values generated from minRerouting. As can be seen from Supplementary Figs. 30 and 31, flux sampling fails to perform as well as minRerouting in terms of minimising the change between two mutant states but does better for iIT341 than for iEK1008. This is probably because we only took 2000 flux solutions while iEK1008 has more reactions than iIT341 norm, and thus, a sample size of 2000 captures its solution space better.

The outlier groups for all models barring e_coli_core, mentioned previously as the top left bubbles, consisting of clusters with a low Synthetic Accessibility and high cluster size, also have a high net flux difference which can be seen in Fig. 6. From the bubble plot, we see that the outliers exist as two subgroups in each species. The outliers that still have a lesser flux difference all have higher Synthetic Accessibility and lesser cluster size, as seen in the bottom left of each panel. The ones with higher flux difference are the second group with lower Synthetic Accessibility and bigger cluster size, as seen in the top right corner of each panel. The pentose phosphate pathway module and reactions from various transport submodules show up as outliers in many of the organisms, along with species-specific reactions. In K. pneumoniae, we see histidine metabolism as a replacement for the TCA cycle as a major subgroup of the outliers. For M. tuberculosis we see the reaction space from mycolic acid metabolism, and for H. pylori, reactions from the urea metabolism are present in the outliers. These outliers are important reactions of interest since they are not bypassed even though their deletion causes a major disruption in the flux distribution of the organism.Fig. 6 Cluster size, synthetic accessibility and the net flux difference between the two metabolic states of outliers across species.

Here, we can see that RSLs have the highest flux difference owing to their cluster size. Pairs with higher synthetic accessibility have lesser flux differences despite large cluster sizes.

Effect of media/environment on flux rerouting in synthetic lethals

The characteristics of minRerouting studied for the models also change when the media concentrations change. As seen in Supplementary Fig. 22, due to a change in environment that led to more reactions being activated, the cluster size almost always increased without an increase in synthetic accessibility. This may indicate that organisms prefer rerouting flux through already active reactions instead of switching on a new reaction. This preference (small synthetic accessibility) may have to do with the need for impromptu transcription and translation of the mRNA for the new reaction. These results are consistent across the models. On the other hand, though minimised, the flux difference shows extremely large values mainly due to a general increase in the flux rates of reactions, including the growth rate. Thus, the difference in flux between alternate routes carrying a high flux will result in a high flux difference.

The rerouting also occurs with two strategies: PSLs and RSLs. What difference does the plasticity or redundancy make to the observations we have seen? We next check how they impact the robustness of metabolism.

Metabolic efficiency analysis

We next analysed the metabolic efficiency of the reactions that comprise double lethal pairs, based on pFBA reaction types40. The relation between these types and the class of synthetic lethal pair has also been studied and the results obtained have been elucidated below.

High occurrence of PSLs over RSLs

Using flux variability analysis and the methodology described in the subsection “Plasticity and redundancy in synthetic lethals”, we have classified the clusters into redundant synthetic lethals (RSLs) and plastic synthetic lethals (PSLs). This classification of synthetic pairs is based on the two strategies taken by the synthetic pairs to maintain robustness. It helps explain different properties seen till now concerning cluster size, flux difference, and composition. Previous analysis17,18 has revealed that PSLs are the more complicated way of acquiring redundancy, needing sophisticated functional organisation with fewer resources, while RSLs represent a more rudimentary strategy. For simple species like E. coli, their results showed that RSLs had more intra-pathway submodule pairs, PSLs had more inter-pathway submodule pairs, and for M. pneumoniae, no pattern was observed for RSLs or PSLs. From Fig. 7, we see that the number of RSLs is higher than PSLs only in the e_coli_core model, which is the simplified metabolic network of E. coli without any functional or structural complexities of metabolic networks. Thus, organisms tend to promote the plasticity of networks. The clusters do not prefer inter-pathway or intra-pathway reaction pairs but seem to depend more on the organism and which pathways are important to it. In E. coli, the most common submodules for pairs were inter-pathway due to their dependency on the cell envelope synthesis pathway and the membrane lipid biosynthesis pathway.Fig. 7 PSL-RSL distribution across models.

While in most models, the fraction of PSLs is much greater than that of RSLs, the trend is reversed in the case of the model e_coli_core. This could be expected because e_coli_core only comprises the metabolic core of E. coli, and hence, most of the double lethal pairs comprise reactions that require to be simultaneously active.

From Supplementary Fig. 13, we see that PSLs exhibit a smaller cluster size for each species than RSLs. This smaller cluster size indicates that transitioning between two states is structurally simpler, with fewer changes needed. This ease of transition is also reflected in the functional costs, as PSLs have a lower flux difference between the two states, as seen in Supplementary Fig. 15. The difference in costs and efficiency is highlighted by the higher Synthetic Accessibility of PSLs as compared to RSLs in Supplementary Fig. 14. Again, these trends are not seen in e_coli_core, which is unique due to limited flexibility stemming from its lesser reactions and connectivity.

For E. coli, we compared the expression of RSL and PSL genes across contexts for 10 synthetic lethal pairs each, as seen in the iModulonDB44. The gene co-expression patterns are represented as heatmaps in Supplementary Figs. 32 and 33. For the subset of lethals analysed, RSL genes (given in Supplementary Table 2) were co-expressed, though not always to the same extent. Specifically, genes coding for the reactions RPI (b2914, b4090), RPE (b3386), TKT1 (b2935, b2465), and TKT2 (b2935, b2465) co-expressed together as these form synthetic lethals with each other. Conversely, PSL genes (Supplementary Table 1) exhibited differential expression for at least two genes coding for the two reactions in a pair. From the heatmap, we can also see that some reactions, such as AACPS3, were expressed in a few environmental conditions, and the corresponding lethal pair reaction 3HAD160 was expressed in mutually exclusive other conditions. Thus, the classification of synthetic lethal as RSLs or PSLs may be dependent on environmental conditions as well.

While the RSLs and PSLs help us understand the above trends, why does the network have differing strategies? Why do not the RSL reactions in a pair compete against each other, with only one reaction becoming evolutionarily stable? In PSLs, why does the inactive reaction stay in the network when it is not needed in wild-type scenarios? To understand why this happens and dive deeper into the types of reactions that contribute to RSL and PSL pairs, the metabolic efficiency analysis of the reactions was performed using the reaction classes returned by pFBA40.

PSLs and RSLs constitute distinct classes of reactions

Using the approach from subsection “Parsimonious FBA (pFBA)”, we see how the type of reaction, pFBA Optimal, blocked, enzymatically less efficient (ELE), metabolically less efficient (MLE), zero flux reaction, or essential reactions, influenced double lethal pair formation. From Supplementary Fig. 2, we see that the Redundancy Index does not depend on the type of reaction. In fact, we see a high composition of ELE and MLE reactions to be redundant which is unexpected yet concurrent with the study by Wang et al.38. We also see that the percentage composition of RSL pairs is lower than PSLs, possibly because of the more evolved nature of PSLs. Yet, the RSLs are found to be composed of mainly pFBA Optima reaction pairs (Fig. 8). In PSLs, while one of the reactions was pFBA Optima (possibly the active reaction in wild type state), the second reaction was seen to be part of any of the other types of reactions. The reason reactions other than pFBA Optima form redundant pairs in nature is puzzling but can be explained by various evolutionary forces acting simultaneously over the organism while under the influence of its environment38. For E. coli, we found the same PSLs to repeat in constraint-based computational experiments9. We found one PSL pair to have been experimentally validated—ASNS1 and ASNS2, with the former being the efficient active enzyme under normal growth conditions45. In literature, we found two enzyme studies that described the efficiency of PSL pairs in M. tuberculosis. One of these pairs is depicted in Supplementary Fig. 28. Previous studies show that DHPS2 is the main enzyme catalysing the formation of folate, an important co-factor, but FOLD3 is present as a less efficient enzyme for the same function46,47. Similarly, there exist two thymidylate synthases in M. tuberculosis that catalyse the biosynthesis of deoxythymidine monophosphate (dTMP or thymidylate). The product of ThyA (TMDS) is the efficient enzyme that is backed up by ThyX (TMDS3)48,49. The role of these redundant processes in the resistance to the drug para-aminosalicylic acid highlights the need for minRerouting to predict alternate pathways of drug targets that may lead to drug resistance50.Fig. 8 Schematic of reaction pair distribution between the RSL and PSL classes for e_coli_core and iYL1228.

a, c Majority of the reaction pairs that are classified as RSL pairs are (pFBA optimal, pFBA optimal) pairs. b, d At least one reaction is pFBA optimal in PSLs. The distribution for the rest of the models can be accessed in the supplementary results provided in the Supplementary Information.

From the metabolic efficiency analysis, it is clear that the type of reaction influences the formation of different classes of synthetic lethals and, in turn, its properties. Comparisons between e_coli_core and iML1515 also reveal that the extent of connections and density of the network also impact the rerouting of metabolism. The usage of our algorithm and these results indicate how the robustness of metabolic networks is more complex than anticipated.

Discussion

In this study, we have looked at how organisms, especially pathogenic bacteria, rely on the presence of redundancy in the form of synthetic lethals to make themselves robust against genetic perturbations. With the help of synthetic lethals and our proposed algorithm, minRerouting, we have analysed how flux is redirected in metabolic networks upon perturbations. Mahadevan and Lovley4 had previously established a relationship between the mode of metabolism and their redundancy. Typically, specialists do not rely on alternative pathways as compared to generalist species. From our results, we see that the bacteria have differing numbers of single and double lethals, implying differing strategies for ensuring robustness. However, more analysis needs to be carried out to further strengthen this hypothesis. These synthetic lethals differ in functionality as well. Yet, there is a set of synthetic lethals that are common across species. These are related to the central carbon metabolism, like the pentose phosphate pathway. Reactions that are from these common pathways are important as their functionalities could be potential targets for drug therapy and combating antimicrobial resistance. They represent nodes in the metabolic network that interact with multiple other nodes simultaneously.

Our optimisation formulation is used to obtain the minimal set of reactions through which fluxes are rerouted. We used three different approaches based on three norms—sparse (L0), linear (L1) and quadratic (L2) approaches—and obtained different minRerouting sets for each of the norms. It is likely that the L0 captures the final steady state of the cell post-adaptation to the knock-out, while the other norms capture the transient response of the cell to the perturbation, i.e., gene deletion, similar to results observed by Shlomi et al.51.

We also showed that the size of the minRerouting set is the smallest when the sparse formulation is used and learned that sparsity forces the use of exclusive reaction rerouting, which results in a null common rerouted reaction set (see Supplementary Fig. 4).

We next proposed a systematic, conditional FVA approach to classify double lethal pairs into PSL (Back-Up Reactions) or RSL (Parallel-Use Reactions) reaction pairs. Parallel to remarks made previously18, the Synthetic Accessibility of PSLs highlights their sophisticated nature in contrast to RSLs, and in all models barring e_coli_core, the PSLs were greater in number. These RSL and PSL pairs have differing properties in terms of the number of reactions through which flux is rerouted, the new set of reactions that become active when the flux is rerouted, as well as the costs of rerouting. We explored the reasons for such a disparity and hypothesised that the RSL pairs should have higher metabolic efficiency as both the reactions are simultaneously active, while reactions in the PSL pairs would have lower metabolic efficiency. The results of our study have proved that this was indeed true. The RSL reaction pairs of most of the organisms (excluding Shigella sonnei), comprise completely (or majorly) of (pFBAOptimal, pFBAOptimal) reaction pairs. As pFBAOptimal reactions are considered crucial for the growth of an organism, our hypothesis was validated. Reactions comprising PSLs must have come about by various evolutionary processes such as horizontal gene transfer, and pleiotropy and may not have been a result of being a back-up for adaptation or robustness. These observations were made possible by the results from the minRerouting algorithm.

An extension of minRerouting can be used to understand the complex metabolic reroutings that occur in several diseases. Particularly, in the case of cancer where the cells re-programme their metabolic activities, rerouting fluxes in such a way that they can continue to proliferate and maintain their malignant properties, minRerouting can help us understand these reroutings and perhaps help in finding better therapeutic cures. Potential synergistic effects of drugs can be unearthed from applications of this algorithm. For some cancers, synthetic lethal drug targets have been tested experimentally and computationally52. minRerouting has previously been used to identify how metabolic flux changes between two synthetic lethal pathways53. Appropriately, bacterial studies of synthetic lethals are also of importance. A previous study used constrain-based modelling to predict the drug targets of multiple M. tuberculosis drugs54. Some of these targets (inha and fas) are redundant, as seen by our analysis, and thus, this may lead to antibiotic resistance of M. tuberculosis. Gene expression experiments for drug treatment of M. tuberculosis with Isoniazid reveal many genes to be induced in response to the subsequent inactivation of the fatty acid synthases caused by the drug55–57. Isoniazid induces a shift from the TCA cycle to the fatty acid metabolism and mycolic acid synthesis. The increase in gene products of cmaA2, mmaA2, mmaA3, mmaA4 (MYC1M1, MYC1M2, MYC1CYC4, and MYC1CYC5, MYC2CYC1, MYC2CYC2, MYC2CYC3) and fbpC (FBPA) seen experimentally match the overall reactions seen in the synthetic lethal clusters or rerouting sets of fatty acid synthase, FAS160, synthetic lethal pairs55,57. Rc3140 (FADE233), and ICL (ICL), genes playing a role in the shift in the pathways above, are also seen in most of the synthetic lethal clusters in the iEK1008 model and have a high RCI of 0.35 and 0.43, respectively56,57.

The submodule distribution and the metabolic efficiency results take us one step closer to exploring the world of synthetic lethals present across metabolic networks. They reveal hidden dependencies between reactions and the influence they have on flux rerouting in a network. This approach of interpreting flux switching is crucial to understanding redundancies in metabolic networks and their differing roles.

Methods

Flux balance analysis

Flux balance analysis (FBA)58,59 is a constraint-based approach that is used to predict the steady-state flux distribution in a given organism’s metabolic network. FBA employs a linear programming (LP) formulation, with an objective to maximise the biomass flux under certain flux and stoichiometric constraints. The formulation of FBA is as follows:1 maxc⊤v≡maxvbio

2 s.t.Sv=0;vLB≼v≼vUB;

Here, v represents the flux vector, the jth entry corresponds to the flux through the jth reaction, and c represents the objective function. Typically cTv = vbio, where vbio is the biomass flux, S represents the stoichiometric matrix of dimensions m × r, with m being the number of metabolites and r being the number of reactions, and vLB and vUB represent the permissible lower and upper bounds of the reaction fluxes. FBA has been experimentally validated in many scenarios60,61, and has widespread applications62. An important extension of FBA is the Minimisation of Metabolic Adjustment (MoMA)24, which seeks to identify a minimally different flux from the wild-type flux (by minimising the L2 norm of this difference). In scenarios where a reaction flux has been altered due to perturbations or additional constraints, MoMA finds the minimal changed flux values compared to the wild-type.

Identification of synthetic lethals

Fast-SL63,64 is an efficient algorithm that identifies synthetic lethals by systematic pruning of the search space and exhaustive enumeration from the remaining reactions. Fast-SL rapidly identifies synthetic lethals and scales well for higher-order lethals. In this paper, we used Fast-SL to identify synthetic lethal reaction pairs in a given genome-scale metabolic model.

minRerouting formulation

Figure 9 represents a toy metabolic network comprising ten reactions and seven metabolites. The metabolites A and B are the ‘input’ metabolites and G is the ‘output’ metabolite. In such a case, we can see that the network consists of three double lethal pairs: {(R1, R2), (R4, R6) and (R5, R6)}. Taking the reaction pair (R4, R6) into consideration, we can see that when reaction R4 is active, and R6 is deleted or inactive, all the fluxes will be routed through reactions R4 and R5, as shown in Fig. 9c. Similarly, when reaction R6 is active, and R4 is deleted or inactive, all the fluxes will be routed through reaction R6, as shown in Fig. 9b. In addition to these changes, the fluxes routed through the remaining reactions could vary based on which pathway is chosen. For instance, the flux through reaction R3 could be significantly higher when the R4 pathway is used than when the R6 pathway is used. We define the ‘minRerouting set’ as the minimal reaction set comprising all the reactions that have a modified flux following the individual deletion of synthetic double lethal reactions. Hence, from Fig. 9, when reaction R4 is active, and reaction R6 is inactive or deleted, the first rerouting set becomes R3, R4, and R5. When reaction R6 is active, and reaction R4 is inactive or deleted, the second rerouting set becomes R3 and R6. The common rerouting set for the double lethal pair (R4, R6) consists of R3 and the complete ‘minRerouting set’ for the double lethal pair becomes R3, R4, R5, and R6.Fig. 9 A sample metabolic network used to illustrate the concept of minRerouting.

The nodes and edges are depicted as metabolites and reactions, respectively. This metabolic network comprises 10 reactions, 7 metabolites, and 3 double lethal pairs. The minRerouting Sets 1 and 2 are highlighted using orange and blue rectangular boxes, and the common reaction between the rerouting sets is R3.

In order to determine the minRerouting set, we first obtain the wild-type flux distribution, vWT, and the set of all lethal pairs in the model. Then, for each lethal pair Ri and Rj, the optimal flux distributions vΔRi and vΔRj, that minimise the distance between the two flux distributions is obtained, by an extension of the MOMA formulation. The generalised p-norm formulation for obtaining the minRerouting of a metabolic network, for a given lethal pair, is as follows:

Step 1: An adaptation of MOMA is performed to obtain the optimal flux distributions vΔRi and vΔRj with minimal flux distance between them:3 min∥vΔRi−vΔRj∥p

4 s.t.SvΔRi=0;SvΔRj=0

5 vLB≼vΔRi≼vUB;vLB≼vΔRj≼vUB

6 vΔRi,Ri=0;vΔRj,Rj=0

7 vΔRi,bio≥(1−γ)vΔRi,bio*

8 vΔRj,bio≥(1−γ)vΔRj,bio*

where vΔRi and vΔRj represent the flux distribution when reactions Ri and Rj are deleted, respectively. vΔRi,bio* and vΔRj,bio* represent the optimal biomass flux in the models when reaction Ri and Rj are deleted, respectively. vΔRi,Ri and vΔRj,Rj represent the flux through reactions Ri and Rj in the models where Ri and Rj are deleted, respectively. γ is the growth rate slack provided for the new flux distributions vΔRi and vΔRj from the optimal biomass vΔRi,bio* and vΔRj,bio*.

Step 2: The flux distributions vΔRi and vΔRj, obtained from Eqs. (3)–(8) are analysed. The reactions that have different flux values in vΔRi and vΔRj are identified as the rerouting set. The size of the rerouting set, the individual reaction flux difference and the total flux difference are subsequently analysed.

For obtaining the L0-norm solution, we used the LP formulation and IBM ILOG CPLEX v12.8 solver as it was one of the few solvers which supported L0-norm optimisation. We used the Gurobi solver for the L1-norm and L2-norm optimisations.

Plasticity and redundancy in synthetic lethals

Previously18, it has been suggested that synthetic lethal reaction pairs can be classified into two categories: plastic synthetic lethal (PSL) and redundant synthetic lethal (RSL). PSL comprises of reaction pairs where one reaction acts as a backup for the other, i.e., the second reaction becomes active when the first reaction is deleted. RSL comprises reaction pairs where both reactions are active simultaneously.

The classification approach proposed by a previous study18 is based on flux vectors which are predicted using FBA. While an FBA solution satisfies all the flux and stoichiometric constraints for a given model, it only represents one possible flux instance from the permissible flux space. Here, we propose a classification approach that is more systematic and thorough, taking into consideration the allowable flux space for each reaction that is part of a double lethal pair. Thus, if solutions where both reactions are not active at the same time can be found, the pair is classified as a PSL.

For instance, using the above approach, if the absolute fluxes of two synthetic lethal reactions, R1 and R2, obtained from FBA, are greater than 0, then, they are classified as RSL reactions. However, there is also a chance that R2 can accommodate zero flux without any change in the optimal biomass flux, while R1 is active. In this case, the reaction pair would have to be classified as a PSL pair. As FBA only picks a single flux instance from the permissible space, we would not be able to correctly classify these reaction pairs. This necessitates a more thorough and systematic manner of classifying the reaction classes.

In order to classify the lethal pairs as PSL or RSL, we performed a flux variability analysis (FVA)65 on the model. FVA is used to obtain the maximum and minimum flux values that a reaction can carry in a model. It solves two LP problems (maximisation and minimisation) for each reaction in the model while constraining the objective function (or) biomass growth rate value. FVA is formulated as follows:9 min/maxvj

10 s.t.Sv=0;vLB≼v≼vUB;vbio=vWT,bio;

The product of the minimum and maximum flux ranges is used to determine the category of the reaction pair. In case a reaction Ri is always active, the product of its minimum and maximum fluxes will always be positive as 0 is not in the range of permissible flux values. When both the reactions are simultaneously active (with a positive or negative flux), while satisfying the biomass constraint, it is considered an RSL pair. The product of the sign of the minimum and maximum fluxes for RSL pairs given as [Ri, Rj] includes the combinations [(<0, <0), (<0, <0)], [(>0, >0), (>0, >0)], [(<0, <0), (>0, >0)], and [(>0, >0), (<0, <0)]. In cases where this is not satisfied, a conditional FVA is performed before the double lethal pair is classified as PSL or RSL.

For each of the ambiguous conditions, two conditional FVAs are performed. In the first FVA, the flux of Ri is constrained to be >0 and in the second, the flux of Ri is constrained to be <0. In this manner, the maximum and minimum fluxes of reaction Rj are obtained when reaction Ri is active. If the reaction Rj can carry a flux value of zero, in either of the two constraint conditions, the reaction pair is considered to be a PSL reaction pair, as Rj can be inactive when Ri is active. However, if Rj is always active under both constraint conditions, the reaction pair is considered to be an RSL. We used this process to determine the classification of PSL and RSL classes instead of relying on a simple FVA because, in FVA, we obtain the maximum and minimum flux values of one reaction, independent of the activity of the other. The whole process is explained pictorially in Fig. 10.Fig. 10 Flowchart depicting the classification of reaction pairs into PSL or RSL pairs.

After the initial FVA, conditional FVAs are performed to classify the ambiguous reaction pairs. The reaction pairs are classified as RSL only when both reactions are simultaneously active.

Parsimonious FBA (pFBA)

pFBA40 employs a bi-level optimisation problem, where first, an optimal flux distribution that maximises the biomass production is identified, followed by the minimisation of total flux through all reactions. In addition to obtaining this flux distribution, pFBA also classifies all the reactions based on their enzymatic/metabolic efficiency. pFBA classifies reactions into a total of six classes: (i) Essential, (ii) enzymatically less efficient (ELE), (iii) metabolically less efficient (MLE), (iv) pFBA optimal, (v) blocked, and (vi) zero flux reactions. We used these reaction classes further, to classify the reactions in the lethal pairs and derive insights into the categorical distribution of these classes across lethal pairs.

Flux redistribution by synthetic lethals

For each synthetic lethal pair, the rerouting set would comprise two subsets—one for each of the double lethal reactions, seen in Supplementary Fig. 1. The reaction sets were analysed, and the properties studied, along with their definitions, are provided below:Size of the synthetic lethal (SL) cluster, number of reactions with a modified flux given by the union of minRerouting subsets 1 and 2 or the complete minRerouting set

Size of the common synthetic lethal cluster, the intersection of the two minRerouting subsets representing the reactions that remain active with a modified flux

Net flux difference, the absolute sum of the net difference between the final flux vectors obtained using minRerouting

Synthetic Accessibility (SA), the fraction of distinct reactions flux is rerouted through for a given SL pair:SA=SymmetricdifferenceofthesubsetsCompleteminReroutingset

Redundancy index (RI), the fraction of synthetic lethal pairs in which a reaction (Ri) occurs within a species:RI=NumberofsyntheticlethalpairsRiispartofTotalnumberofsyntheticlethalpairs

Reaction compensation index (RCI), the fraction of synthetic clusters in which a reaction (Ri) occurs within a species:RCI=NumberofsyntheticlethalclustersRiispartofTotalnumberofsyntheticlethalclusters

Changing the environment of organisms

The media was changed to depict the change in the environment faced by an organism by fixing the lower bound concentrations of the exchange reactions. First, we changed the concentration of 10 carbon source exchange reactions in each model to—1000 mmol g DW−1 h−1. The carbon sources are sucrose, glucose, glycerol, melibiose, cellobiose, alpha-galactose, beta-galactose, arabinose, trehalose, and fructose. These models were hence known as Model_carbon for each species. Second, we looked at the effect of unblocking all exchange reactions in the model to varying degrees, i.e., 100 mmol g DW−1 h−1, 500 mmol g DW−1 h−1, and 1000 mmol g DW−1 h−1. These models were named Model_exchange_1, Model_exchange_2, and Model_exchange_3. In total, four new models were generated for each species with different environments, and the minRerouting analysis was repeated for them.

Metabolic subsystem analysis

A set of reactions that share a similar metabolic function is referred to as a metabolic subsystem66. All the reactions in a genome-scale metabolic model (GSMM) are categorised under different metabolic subsystems. To identify the metabolic subsystem of the reactions that comprise a double lethal pair, the set of all distinct reactions to be analysed is obtained. Then, the subsystem of each of these reactions is obtained by systematically querying the BiGG database25.

Implementation

The implementation of minRerouting and initial analysis were done using MATLAB. The metabolic cost analysis and flux rerouting analysis were done using Python and R. The figures were generated using R. The COBRA Toolbox67 for MATLAB was used for all metabolic network analysis. All code written as a part of this project is open-sourced and can be accessed at https://github.com/RamanLab/minRerouting/.

Supplementary information

Supplementary Information

Supplementary information

The online version contains supplementary material available at 10.1038/s41540-024-00426-5.

Acknowledgements

T.M. thanks IBSE for the post-baccalaureate fellowship. The other authors received no specific funding for this work.

Author contributions

K.R. conceptualised the study. S.M.N., O.S.M., J.S.N., and K.R. designed the mathematical formulation/computational experiments. T.M., O.S.M., and S.M.N. performed the computational experiments. T.M., O.S.M., S.M.N., and K.R. analysed the data and interpreted the results. K.R. supervised the project. S.M.N., O.S.M., and T.M. drafted the manuscript, with inputs from K.R. All authors read and approved the final manuscript.

Data availability

All models analysed during the study are taken from the BiGG database25.

Code availability

The ‘minRerouting’ algorithm is available from https://github.com/RamanLab/minRerouting.

Competing interests

The authors declare no competing interests.

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

These authors contributed equally: Sowmya Manojna Narasimha, Tanisha Malpani.
==== Refs
References

1. Kitano H Biological robustness Nat. Rev. Genet. 2004 5 826 837 10.1038/nrg1471 15520792
Kitano, H. Biological robustness. Nat. Rev. Genet. 5, 826–837 (2004).15520792
2. Wagner, A. Robustness and Evolvability in Living Systems (Princeton University Press, 2005).
3. Wagner A Robustness, evolvability, and neutrality FEBS Lett. 2005 579 1772 1778 10.1016/j.febslet.2005.01.063 15763550
Wagner, A. Robustness, evolvability, and neutrality. FEBS Lett. 579, 1772–1778 (2005).15763550
4. Mahadevan R Lovley DR The degree of redundancy in metabolic genes is linked to mode of metabolism Biophys. J. 2008 94 1216 1220 10.1529/biophysj.107.118414 17981891
Mahadevan, R. & Lovley, D. R. The degree of redundancy in metabolic genes is linked to mode of metabolism. Biophys. J. 94, 1216–1220 (2008).17981891
5. Sambamoorthy G Sinha H Raman K Evolutionary design principles in metabolism Proc. Biol. Sci. 2019 286 20190098 30836874
Sambamoorthy, G., Sinha, H. & Raman, K. Evolutionary design principles in metabolism. Proc. Biol. Sci. 286, 20190098 (2019).30836874
6. Baba T Construction of Escherichia coli K-12 in-frame, single-gene knockout mutants: the Keio collection Mol. Syst. Biol. 2006 2 2006 0008 10.1038/msb4100050
Baba, T. et al. Construction of Escherichia coli K-12 in-frame, single-gene knockout mutants: the Keio collection. Mol. Syst. Biol. 2, 2006–0008 (2006).
7. Goodall EC The essential genome of Escherichia coli k-12 MBio 2018 9 e02096 17 10.1128/mBio.02096-17 29463657
Goodall, E. C. et al. The essential genome of Escherichia coli k-12. MBio 9, e02096–17 (2018).29463657
8. Gerdes S Experimental determination and system level analysis of essential genes in Escherichia coli mg1655 J. Bacteriol. 2003 185 5673 5684 10.1128/JB.185.19.5673-5684.2003 13129938
Gerdes, S. et al. Experimental determination and system level analysis of essential genes in Escherichia coli mg1655. J. Bacteriol. 185, 5673–5684 (2003).13129938
9. Ghim C-M Goh K-I Kahng B Lethality and synthetic lethality in the genome-wide metabolic network of Escherichia coli J. Theor. Biol. 2005 237 401 411 10.1016/j.jtbi.2005.04.025 15975601
Ghim, C.-M., Goh, K.-I. & Kahng, B. Lethality and synthetic lethality in the genome-wide metabolic network of Escherichia coli. J. Theor. Biol. 237, 401–411 (2005).15975601
10. Sambamoorthy G Raman K Understanding the evolution of functional redundancy in metabolic networks Bioinformatics 2018 34 i981 i987 10.1093/bioinformatics/bty604 30423058
Sambamoorthy, G. & Raman, K. Understanding the evolution of functional redundancy in metabolic networks. Bioinformatics 34, i981–i987 (2018).30423058
11. Hartman JL Garvik B Hartwell L Principles for the buffering of genetic variation Science 2001 291 1001 1004 10.1126/science.1056072 11232561
Hartman, J. L., Garvik, B. & Hartwell, L. Principles for the buffering of genetic variation. Science 291, 1001–1004 (2001).11232561
12. Becker SA Palsson BO Context-specific metabolic networks are consistent with experiments PLoS Comput. Biol. 2008 4 1 10 10.1371/journal.pcbi.1000082
Becker, S. A. & Palsson, B. O. Context-specific metabolic networks are consistent with experiments. PLoS Comput. Biol. 4, 1–10 (2008).
13. Zur H Ruppin E Shlomi T iMAT: an integrative metabolic analysis tool Bioinformatics 2010 26 3140 3142 10.1093/bioinformatics/btq602 21081510
Zur, H., Ruppin, E. & Shlomi, T. iMAT: an integrative metabolic analysis tool. Bioinformatics 26, 3140–3142 (2010).21081510
14. Kim J Reed JL RELATCH: relative optimality in metabolic networks explains robust metabolic and regulatory responses to perturbations Genome Biol. 2012 13 R78 10.1186/gb-2012-13-9-r78 23013597
Kim, J. & Reed, J. L. RELATCH: relative optimality in metabolic networks explains robust metabolic and regulatory responses to perturbations. Genome Biol. 13, R78 (2012).23013597
15. Pandey V Hadadi N Hatzimanikatis V Enhanced flux prediction by integrating relative expression and relative metabolite abundance into thermodynamically consistent metabolic models PLoS Comput. Biol. 2019 15 1 23 10.1371/journal.pcbi.1007036
Pandey, V., Hadadi, N. & Hatzimanikatis, V. Enhanced flux prediction by integrating relative expression and relative metabolite abundance into thermodynamically consistent metabolic models. PLoS Comput. Biol. 15, 1–23 (2019).
16. Ravi S Gunawan R δfba-predicting metabolic flux alterations using genome-scale metabolic models and differential transcriptomic data PLoS Comput. Biol. 2021 17 1 18 10.1371/journal.pcbi.1009589
Ravi, S. & Gunawan, R. δfba-predicting metabolic flux alterations using genome-scale metabolic models and differential transcriptomic data. PLoS Comput. Biol. 17, 1–18 (2021).
17. Massucci FA Sagués F Serrano MA Metabolic plasticity in synthetic lethal mutants: viability at higher cost PLoS Comput. Biol. 2018 14 1 20 10.1371/journal.pcbi.1005949
Massucci, F. A., Sagués, F. & Serrano, M. A. Metabolic plasticity in synthetic lethal mutants: viability at higher cost. PLoS Comput. Biol. 14, 1–20 (2018).
18. Güell O Sagués F Serrano MÁ Essential plasticity and redundancy of metabolism unveiled by synthetic lethality analysis PLoS Comput. Biol. 2014 10 e1003637 10.1371/journal.pcbi.1003637 24854166
Güell, O., Sagués, F. & Serrano, M. Á. Essential plasticity and redundancy of metabolism unveiled by synthetic lethality analysis. PLoS Comput. Biol. 10, e1003637 (2014).24854166
19. Andersson DI Levin BR The biological cost of antibiotic resistance Curr. Opin. Microbiol. 1999 2 489 493 10.1016/S1369-5274(99)00005-3 10508723
Andersson, D. I. & Levin, B. R. The biological cost of antibiotic resistance. Curr. Opin. Microbiol. 2, 489–493 (1999).10508723
20. Mahadevan R Schilling C The effects of alternate optimal solutions in constraint-based genome-scale metabolic models Metab. Eng. 2003 5 264 276 10.1016/j.ymben.2003.09.002 14642354
Mahadevan, R. & Schilling, C. The effects of alternate optimal solutions in constraint-based genome-scale metabolic models. Metab. Eng. 5, 264–276 (2003).14642354
21. Fondi M Bosi E Presta L Natoli D Fani R Modelling microbial metabolic rewiring during growth in a complex medium BMC Genom. 2016 17 970 10.1186/s12864-016-3311-0
Fondi, M., Bosi, E., Presta, L., Natoli, D. & Fani, R. Modelling microbial metabolic rewiring during growth in a complex medium. BMC Genom. 17, 970 (2016).
22. Edwards JS Palsson BO Metabolic flux balance analysis and the in silico analysis of Escherichia coli K-12 gene deletions BMC Bioinform. 2000 1 1 10.1186/1471-2105-1-1
Edwards, J. S. & Palsson, B. O. Metabolic flux balance analysis and the in silico analysis of Escherichia coli K-12 gene deletions. BMC Bioinform. 1, 1 (2000).
23. Fischer E Sauer U Large-scale in vivo flux analysis shows rigidity and suboptimal performance of Bacillus subtilis metabolism Nat. Genet. 2005 37 636 640 10.1038/ng1555 15880104
Fischer, E. & Sauer, U. Large-scale in vivo flux analysis shows rigidity and suboptimal performance of Bacillus subtilis metabolism. Nat. Genet. 37, 636–640 (2005).15880104
24. Segrè D Vitkup D Church GM Analysis of optimality in natural and perturbed metabolic networks Proc. Natl Acad. Sci. USA 2002 99 15112 15117 10.1073/pnas.232349399 12415116
Segrè, D., Vitkup, D. & Church, G. M. Analysis of optimality in natural and perturbed metabolic networks. Proc. Natl Acad. Sci. USA 99, 15112–15117 (2002).12415116
25. Schellenberger J Park JO Conrad TM Palsson BØ BiGG: a Biochemical Genetic and Genomic knowledgebase of large scale metabolic reconstructions BMC Bioinform. 2010 11 1 10 10.1186/1471-2105-11-213
Schellenberger, J., Park, J. O., Conrad, T. M. & Palsson, B. Ø. BiGG: a Biochemical Genetic and Genomic knowledgebase of large scale metabolic reconstructions. BMC Bioinform. 11, 1–10 (2010).
26. Orth JD Fleming RM Palsson BØ Reconstruction and use of microbial metabolic networks: the core Escherichia coli metabolic model as an educational guide EcoSal 2010 4 10 1128
Orth, J. D., Fleming, R. M. & Palsson, B. Ø. Reconstruction and use of microbial metabolic networks: the core Escherichia coli metabolic model as an educational guide. EcoSal 4, 10–1128 (2010).
27. Monk JM iML1515, a knowledgebase that computes Escherichia coli traits Nat. Biotechnol. 2017 35 904 908 10.1038/nbt.3956 29020004
Monk, J. M. et al. iML1515, a knowledgebase that computes Escherichia coli traits. Nat. Biotechnol. 35, 904–908 (2017).29020004
28. Barve A Rodrigues JFM Wagner A Superessential reactions in metabolic networks Proc. Natl Acad. Sci. USA 2012 109 E1121 30 10.1073/pnas.1113065109 22509034
Barve, A., Rodrigues, J. F. M. & Wagner, A. Superessential reactions in metabolic networks. Proc. Natl Acad. Sci. USA 109, E1121–30 (2012).22509034
29. Suthers PF Zomorrodi A Maranas CD Genome scale gene/reaction essentiality and synthetic lethality analysis Mol. Syst. Biol. 2009 5 301 10.1038/msb.2009.56 19690570
Suthers, P. F., Zomorrodi, A. & Maranas, C. D. Genome scale gene/reaction essentiality and synthetic lethality analysis. Mol. Syst. Biol. 5, 301 (2009).19690570
30. Marcel E Metabolic flux responses to pyruvate kinase knockout in Escherichia coli J. Bacteriol. 2002 184 152 164 10.1128/JB.184.1.152-164.2002 11741855
Marcel, E. et al. Metabolic flux responses to pyruvate kinase knockout in Escherichia coli. J. Bacteriol. 184, 152–164 (2002).11741855
31. Ishii N Multiple high-throughput analyses monitor the response of E. coli to perturbations Science 2007 316 593 597 10.1126/science.1132067 17379776
Ishii, N. et al. Multiple high-throughput analyses monitor the response of E. coli to perturbations. Science 316, 593–597 (2007).17379776
32. Monk J Multi-omics quantification of species variation of Escherichia coli links molecular features with strain phenotypes Cell Syst. 2016 3 238–251.e12 27667363
Monk, J. et al. Multi-omics quantification of species variation of Escherichia coli links molecular features with strain phenotypes. Cell Syst. 3, 238–251.e12 (2016).27667363
33. Iwasaki T Escherichia coli amino acid auxotrophic expression host strains for investigating protein structure–function relationships J. Biochem. 2020 169 387 394 10.1093/jb/mvaa140
Iwasaki, T. et al. Escherichia coli amino acid auxotrophic expression host strains for investigating protein structure–function relationships. J. Biochem. 169, 387–394 (2020).
34. Schulz-Mirbach H On the flexibility of the cellular amination network in E coli eLife 2022 11 e77492 10.7554/eLife.77492 35876664
Schulz-Mirbach, H. et al. On the flexibility of the cellular amination network in E coli. eLife 11, e77492 (2022).35876664
35. Cotton CA Underground isoleucine biosynthesis pathways in E. coli eLife 2020 9 e54207 10.7554/eLife.54207 32831171
Cotton, C. A. et al. Underground isoleucine biosynthesis pathways in E. coli. eLife 9, e54207 (2020).32831171
36. Fukushima M Kakinuma K Kawaguchi R Phylogenetic analysis of Salmonella, Shigella, and Escherichia coli strains on the basis of the gyrB gene sequence J. Clin. Microbiol. 2002 40 2779 2785 10.1128/JCM.40.8.2779-2785.2002 12149329
Fukushima, M., Kakinuma, K. & Kawaguchi, R. Phylogenetic analysis of Salmonella, Shigella, and Escherichia coli strains on the basis of the gyrB gene sequence. J. Clin. Microbiol. 40, 2779–2785 (2002).12149329
37. He X Zhang J Higher duplicability of less important genes in yeast genomes Mol. Biol. Evol. 2005 23 144 151 10.1093/molbev/msj015 16151181
He, X. & Zhang, J. Higher duplicability of less important genes in yeast genomes. Mol. Biol. Evol. 23, 144–151 (2005).16151181
38. Wang Z Zhang J Abundant indispensable redundancies in cellular metabolic networks Genome Biol. Evol. 2009 1 23 33 10.1093/gbe/evp002 20333174
Wang, Z. & Zhang, J. Abundant indispensable redundancies in cellular metabolic networks. Genome Biol. Evol. 1, 23–33 (2009).20333174
39. Kim P-J Metabolite essentiality elucidates robustness of Escherichia coli metabolism Proc. Natl Acad. Sci. USA 2007 104 13638 13642 10.1073/pnas.0703262104 17698812
Kim, P.-J. et al. Metabolite essentiality elucidates robustness of Escherichia coli metabolism. Proc. Natl Acad. Sci. USA 104, 13638–13642 (2007).17698812
40. Lewis NE Omic data from evolved E. coli are consistent with computed optimal growth from genome-scale models Mol. Syst. Biol. 2010 6 390 10.1038/msb.2010.47 20664636
Lewis, N. E. et al. Omic data from evolved E. coli are consistent with computed optimal growth from genome-scale models. Mol. Syst. Biol. 6, 390 (2010).20664636
41. Ihmels J Collins SR Schuldiner M Krogan NJ Weissman JS Backup without redundancy: genetic interactions reveal the cost of duplicate gene loss Mol. Syst. Biol. 2007 3 86 10.1038/msb4100127 17389874
Ihmels, J., Collins, S. R., Schuldiner, M., Krogan, N. J. & Weissman, J. S. Backup without redundancy: genetic interactions reveal the cost of duplicate gene loss. Mol. Syst. Biol. 3, 86 (2007).17389874
42. Brookfield J Can genes be truly redundant? Curr. Biol. 1992 2 553 554 10.1016/0960-9822(92)90036-A 15336052
Brookfield, J. Can genes be truly redundant? Curr. Biol. 2, 553–554 (1992).15336052
43. Megchelenbrink W Huynen M Marchiori E optGpSampler: An Improved Tool for Uniformly Sampling the Solution-Space of Genome-Scale Metabolic Networks PLoS One 2014 9 e86587 10.1371/journal.pone.0086587 24551039
Megchelenbrink, W., Huynen, M. & Marchiori, E. optGpSampler: An Improved Tool for Uniformly Sampling the Solution-Space of Genome-Scale Metabolic Networks. PLoS One 9, e86587 (2014).24551039
44. Rychel K iModulonDB: a knowledgebase of microbial transcriptional regulation derived from machine learning Nucleic Acids Res. 2020 49 D112 D120 10.1093/nar/gkaa810
Rychel, K. et al. iModulonDB: a knowledgebase of microbial transcriptional regulation derived from machine learning. Nucleic Acids Res. 49, D112–D120 (2020).
45. Humbert R Simoni RD Genetic and biomedical studies demonstrating a second gene coding for asparagine synthetase in Escherichia coli J. Bacteriol. 1980 142 212 220 10.1128/jb.142.1.212-220.1980 6102982
Humbert, R. & Simoni, R. D. Genetic and biomedical studies demonstrating a second gene coding for asparagine synthetase in Escherichia coli. J. Bacteriol. 142, 212–220 (1980).6102982
46. Gengenbacher M Xu T Niyomrattanakit P Spraggon G Dick T Biochemical and structural characterization of the putative dihydropteroate synthase ortholog Rv1207 of Mycobacterium tuberculosis FEMS Microbiol. Lett. 2008 287 128 135 10.1111/j.1574-6968.2008.01302.x 18680522
Gengenbacher, M., Xu, T., Niyomrattanakit, P., Spraggon, G. & Dick, T. Biochemical and structural characterization of the putative dihydropteroate synthase ortholog Rv1207 of Mycobacterium tuberculosis. FEMS Microbiol. Lett. 287, 128–135 (2008).18680522
47. Gibson SER Harrison J Molloy A Cox JAG Cholesterol-dependent activity of dapsone against non-replicating persistent mycobacteria Microbiology 2022 168 001279 10.1099/mic.0.001279
Gibson, S. E. R., Harrison, J., Molloy, A. & Cox, J. A. G. Cholesterol-dependent activity of dapsone against non-replicating persistent mycobacteria. Microbiology 168, 001279 (2022).
48. Hunter JH Gujjar R Pang CKT Rathod PK Kinetics and ligand-binding preferences of mycobacterium tuberculosis thymidylate synthases, thya and thyx PLoS ONE 2008 3 1 10 10.1371/journal.pone.0002237
Hunter, J. H., Gujjar, R., Pang, C. K. T. & Rathod, P. K. Kinetics and ligand-binding preferences of mycobacterium tuberculosis thymidylate synthases, thya and thyx. PLoS ONE 3, 1–10 (2008).
49. Fivian-Hughes AS Houghton J Davis EO Mycobacterium tuberculosis thymidylate synthase gene thyx is essential and potentially bifunctional, while thya deletion confers resistance to p-aminosalicylic acid Microbiology 2012 158 308 318 10.1099/mic.0.053983-0 22034487
Fivian-Hughes, A. S., Houghton, J. & Davis, E. O. Mycobacterium tuberculosis thymidylate synthase gene thyx is essential and potentially bifunctional, while thya deletion confers resistance to p-aminosalicylic acid. Microbiology 158, 308–318 (2012).22034487
50. Mathys V Molecular genetics of para-aminosalicylic acid resistance in clinical isolates and spontaneous mutants of mycobacterium tuberculosis Antimicrob. Agents Chemother. 2009 53 2100 2109 10.1128/AAC.01197-08 19237648
Mathys, V. et al. Molecular genetics of para-aminosalicylic acid resistance in clinical isolates and spontaneous mutants of mycobacterium tuberculosis. Antimicrob. Agents Chemother. 53, 2100–2109 (2009).19237648
51. Shlomi T Berkman O Ruppin E Regulatory on/off minimization of metabolic flux changes after genetic perturbations Proc. Natl Acad. Sci. USA 2005 102 7695 7700 10.1073/pnas.0406346102 15897462
Shlomi, T., Berkman, O. & Ruppin, E. Regulatory on/off minimization of metabolic flux changes after genetic perturbations. Proc. Natl Acad. Sci. USA 102, 7695–7700 (2005).15897462
52. O’Neil NJ Bailey ML Hieter P Synthetic lethality and cancer Nat. Rev. Genet. 2017 18 613 623 10.1038/nrg.2017.47 28649135
O’Neil, N. J., Bailey, M. L. & Hieter, P. Synthetic lethality and cancer. Nat. Rev. Genet. 18, 613–623 (2017).28649135
53. Sahoo S Metabolite systems profiling identifies exploitable weaknesses in retinoblastoma FEBS Lett. 2019 593 23 41 10.1002/1873-3468.13294 30417337
Sahoo, S. et al. Metabolite systems profiling identifies exploitable weaknesses in retinoblastoma. FEBS Lett. 593, 23–41 (2019).30417337
54. Chung BK-S Dick T Lee D-Y In silico analyses for the discovery of tuberculosis drug targets J. Antimicrob. Chemother. 2013 68 2701 2709 10.1093/jac/dkt273 23838951
Chung, B. K.-S., Dick, T. & Lee, D.-Y. In silico analyses for the discovery of tuberculosis drug targets. J. Antimicrob. Chemother. 68, 2701–2709 (2013).23838951
55. Goossens SN Sampson SL Rie AV Mechanisms of drug-induced tolerance in Mycobacterium tuberculosis Clin. Microbiol. Rev. 2020 34 10 1128 10.1128/CMR.00141-20
Goossens, S. N., Sampson, S. L. & Rie, A. V. Mechanisms of drug-induced tolerance in Mycobacterium tuberculosis. Clin. Microbiol. Rev. 34, 10–1128 (2020).
56. Wilson M Exploring drug-induced alterations in gene expression in Mycobacterium tuberculosis by microarray hybridization Proc. Natl Acad. Sci. USA 1999 96 12833 12838 10.1073/pnas.96.22.12833 10536008
Wilson, M. et al. Exploring drug-induced alterations in gene expression in Mycobacterium tuberculosis by microarray hybridization. Proc. Natl Acad. Sci. USA 96, 12833–12838 (1999).10536008
57. Karakousis PC Williams EP Bishai WR Altered expression of isoniazid-regulated genes in drug-treated dormant Mycobacterium tuberculosis J. Antimicrob. Chemother. 2007 61 323 331 10.1093/jac/dkm485 18156607
Karakousis, P. C., Williams, E. P. & Bishai, W. R. Altered expression of isoniazid-regulated genes in drug-treated dormant Mycobacterium tuberculosis. J. Antimicrob. Chemother. 61, 323–331 (2007).18156607
58. Varma A Palsson BO Metabolic flux balancing: basic concepts, scientific and practical use Bio/Technol. 1994 12 994 998 10.1038/nbt1094-994
Varma, A. & Palsson, B. O. Metabolic flux balancing: basic concepts, scientific and practical use. Bio/Technol. 12, 994–998 (1994).
59. Kauffman KJ Prakash P Edwards JS Advances in flux balance analysis Curr. Opin. Biotechnol. 2003 14 491 496 10.1016/j.copbio.2003.08.001 14580578
Kauffman, K. J., Prakash, P. & Edwards, J. S. Advances in flux balance analysis. Curr. Opin. Biotechnol. 14, 491–496 (2003).14580578
60. Varma A Palsson BO Stoichiometric flux balance models quantitatively predict growth and metabolic by-product secretion in wild-type Escherichia coli w3110 Appl. Environ. Microbiol. 1994 60 3724 3731 10.1128/aem.60.10.3724-3731.1994 7986045
Varma, A. & Palsson, B. O. Stoichiometric flux balance models quantitatively predict growth and metabolic by-product secretion in wild-type Escherichia coli w3110. Appl. Environ. Microbiol. 60, 3724–3731 (1994).7986045
61. Edwards JS Ibarra RU Palsson BO In silico predictions of Escherichia coli metabolic capabilities are consistent with experimental data Nat. Biotechnol. 2001 19 125 130 10.1038/84379 11175725
Edwards, J. S., Ibarra, R. U. & Palsson, B. O. In silico predictions of Escherichia coli metabolic capabilities are consistent with experimental data. Nat. Biotechnol. 19, 125–130 (2001).11175725
62. McCloskey D Palsson BO Feist AM Basic and applied uses of genome-scale metabolic network reconstructions of Escherichia coli Mol. Syst. Biol. 2013 9 661 10.1038/msb.2013.18 23632383
McCloskey, D., Palsson, B. O. & Feist, A. M. Basic and applied uses of genome-scale metabolic network reconstructions of Escherichia coli. Mol. Syst. Biol. 9, 661 (2013).23632383
63. Pratapa A Balachandran S Raman K Fast-SL: an efficient algorithm to identify synthetic lethal sets in metabolic networks Bioinformatics 2015 31 3299 3305 10.1093/bioinformatics/btv352 26085504
Pratapa, A., Balachandran, S. & Raman, K. Fast-SL: an efficient algorithm to identify synthetic lethal sets in metabolic networks. Bioinformatics 31, 3299–3305 (2015).26085504
64. Raman K Pratapa A Mohite O Balachandran S Computational prediction of synthetic lethals in genome-scale metabolic models using Fast-SL Methods Mol. Biol. 2018 1716 315 336 10.1007/978-1-4939-7528-0_14 29222760
Raman, K., Pratapa, A., Mohite, O. & Balachandran, S. Computational prediction of synthetic lethals in genome-scale metabolic models using Fast-SL. Methods Mol. Biol. 1716, 315–336 (2018).29222760
65. Mahadevan R Schilling CH The effects of alternate optimal solutions in constraint-based genome-scale metabolic models Metab. Eng. 2003 5 264 276 10.1016/j.ymben.2003.09.002 14642354
Mahadevan, R. & Schilling, C. H. The effects of alternate optimal solutions in constraint-based genome-scale metabolic models. Metab. Eng. 5, 264–276 (2003).14642354
66. Wang H Genome-scale metabolic network reconstruction of model animals as a platform for translational research Proc. Natl Acad. Sci. USA 2021 118 e2102344118 10.1073/pnas.2102344118 34282017
Wang, H. et al. Genome-scale metabolic network reconstruction of model animals as a platform for translational research. Proc. Natl Acad. Sci. USA 118, e2102344118 (2021).34282017
67. Heirendt L Creation and analysis of biochemical constraint-based models using the COBRA Toolbox v.3.0 Nat. Protoc. 2019 14 639 702 10.1038/s41596-018-0098-2 30787451
Heirendt, L. et al. Creation and analysis of biochemical constraint-based models using the COBRA Toolbox v.3.0. Nat. Protoc. 14, 639–702 (2019).30787451
68. Shichun L Synthetic lethality reveals mechanisms of Mycobacterium tuberculosis resistance to β-lactams mBio 2014 5 10 1128
Shichun, L. et al. Synthetic lethality reveals mechanisms of Mycobacterium tuberculosis resistance to β-lactams. mBio 5, 10–1128 (2014).
69. Thiele I Vo TD Price ND Palsson BØ Expanded metabolic reconstruction of helicobacter pylori (i it341 gsm/gpr): an in silico genome-scale characterization of single-and double-deletion mutants J. Bacteriol. 2005 187 5818 5830 10.1128/JB.187.16.5818-5830.2005 16077130
Thiele, I., Vo, T. D., Price, N. D. & Palsson, B. Ø. Expanded metabolic reconstruction of helicobacter pylori (i it341 gsm/gpr): an in silico genome-scale characterization of single-and double-deletion mutants. J. Bacteriol. 187, 5818–5830 (2005).16077130
70. Charusanti P An experimentally-supported genome-scale metabolic network reconstruction for Yersinia pestis CO92 BMC Syst. Biol. 2011 5 1 13 10.1186/1752-0509-5-163 21194489
Charusanti, P. et al. An experimentally-supported genome-scale metabolic network reconstruction for Yersinia pestis CO92. BMC Syst. Biol. 5, 1–13 (2011).21194489
71. Monk JM Genome-scale metabolic reconstructions of multiple Escherichia coli strains highlight strain-specific adaptations to nutritional environments Proc. Natl Acad. Sci. USA 2013 110 20338 20343 10.1073/pnas.1307797110 24277855
Monk, J. M. et al. Genome-scale metabolic reconstructions of multiple Escherichia coli strains highlight strain-specific adaptations to nutritional environments. Proc. Natl Acad. Sci. USA 110, 20338–20343 (2013).24277855
72. Liao Y-C An experimentally validated genome-scale metabolic reconstruction of Klebsiella pneumoniae MGH 78578, iYL1228 J. Bacteriol. 2011 193 1710 1717 10.1128/JB.01218-10 21296962
Liao, Y.-C. et al. An experimentally validated genome-scale metabolic reconstruction of Klebsiella pneumoniae MGH 78578, iYL1228. J. Bacteriol. 193, 1710–1717 (2011).21296962
73. Thiele I A community effort towards a knowledge-base and mathematical model of the human pathogen Salmonella Typhimurium LT2 BMC Syst. Biol. 2011 5 1 9 10.1186/1752-0509-5-8 21194489
Thiele, I. et al. A community effort towards a knowledge-base and mathematical model of the human pathogen Salmonella Typhimurium LT2. BMC Syst. Biol. 5, 1–9 (2011).21194489
