Strategies for Enhancing in vitro Degradation of Linuron by Variovorax sp. Strain SRS 16 Under the Guidance of Metabolic Modeling.

Phenyl urea herbicides are being extensively used for weed control in both agricultural and non-agricultural applications. Linuron is one of the key herbicides in this family and is in wide use. Like other phenyl urea herbicides, it is known to have toxic effects as a result of its persistence in the environment. The natural removal of linuron from the environment is mainly carried through microbial biodegradation. Some microorganisms have been reported to mineralize linuron completely and utilize it as a carbon and nitrogen source. Variovorax sp. strain SRS 16 is one of the known efficient degraders with a recently sequenced genome. The genomic data provide an opportunity to use a genome-scale model for improving biodegradation. The aim of our study is the construction of a genome-scale metabolic model following automatic and manual protocols and its application for improving its metabolic potential through iterative simulations. Applying flux balance analysis (FBA), growth and degradation performances of SRS 16 in different media considering the influence of selected supplements (potential carbon and nitrogen sources) were simulated. Outcomes are predictions for the suitable media modification, allowing faster degradation of linuron by SRS 16. Seven metabolites were selected for in vitro validation of the predictions through laboratory experiments confirming the degradation-promoting effect of specific amino acids (glutamine and asparagine) on linuron degradation and SRS 16 growth. Overall, simulations are shown to be efficient in predicting the degradation potential of SRS 16 in the presence of specific supplements. The generated information contributes to the understanding of the biochemistry of linuron degradation and can be further utilized for the development of new cleanup solutions without any genetic manipulation.

Introduction

Phenyl urea herbicides are among the most widely used herbicides for weed control in several crops (mostly cereals) through their pre- or post-emergence applications. These substances interrupt electron transfer in photosystem II, leading to the formation of reactive oxygen species and resulting in cell damage (Liu, 2010). The increased rate of application of xenobiotics such as the phenyl urea herbicides in recent years enhanced their burden to the environment due to their persistence in the surroundings (Hasanuzzaman et al., 2020). These compounds are consistently found to have negative effects on the ecosystem, including hazards to human health (de Souza et al., 2020; Garcês et al., 2020). Linuron has been one of the most widely applied phenyl urea herbicides in agriculture practice, also reported as an environmental pollutant (Dejonghe et al., 2003; Horemans et al., 2016). Remediating the environment from the accumulated toxic substances is the focus of a growing number of researches efforts. The degradation or the removal of such pollutants is reported to be conducted through physical, (photo-)chemical, and chemical processes (Katsumata et al., 2011; Reddy and Kim, 2015; Kovács et al., 2016; Hao et al., 2019; Bhat et al., 2020). Microorganisms play a major role in the biological removal of linuron. So far, different strategies have been approached for enhancing biodegradation and bioremediation (Raman and Chandra, 2009; Pimviriyakul et al., 2020). By promoting the growth of soil microbial degraders, soil amendments are being used in order to accelerate the removal rate of pollutants (Bao et al., 2020). However, the selection of the amendments is mostly based on trial and error.

Biodegradation of linuron is generally initiated by amidase hydrolases, which leads to the formation of a more toxic intermediate, 3,4-dichloroaniline (DCA). DCA further degrades to metabolites that can be consumed in the central metabolism of the microorganisms. However, partial degradation of linuron produces chloroanilines that are more toxic than linuron itself. Bacterial genera such as Arthrobacter, Bacillus, Comamonas, Pseudomonas, Sphingobium, and Variovorax are able to degrade/transform linuron in various environments either in isolation or as part of a consortium (Turnbull et al., 2001; Dejonghe et al., 2003; Lerner et al., 2020; Zhang et al., 2020). Variovorax SRS 16 is widely studied for its ability to utilize linuron as a sole carbon and nitrogen source. Variovorax SRS 16 also possesses a modified chlorocatechol ortho-cleavage pathway, allowing the utilization of linuron as a substrate for growth. The presence of the amidase (specific linuron hydrolase) and gene clusters which are responsible for the degradation of linuron and dichloroaniline in SRS 16 has been described in detail with the support of proteomic studies (Bers et al., 2011). Based on the relatively extensive phenotypic information, together with the publication of its genome sequence (Sørensen et al., 2005; Bers et al., 2011), Variovorax SRS 16 is recognized as a good model for exploring biodegradation of phenyl urea herbicides.

Processing genomic information into a metabolic model is increasingly used as a route for generating a predictive tool to elucidate and manipulate cellular biochemical activity. Genome-scale metabolic modeling has been proven as an efficient approach to decode the genomic and functional information for a specific phenotype by investigating the gene–protein interactions on a cellular level (Thiele and Palsson, 2010; Gu et al., 2019). To construct metabolic models, some preliminary information on physiological requirements (mainly of growth) of the organism is needed for better curation and to ensure the appropriate functioning of the model (Covert et al., 2001). This approach generally follows some defined steps such as collecting basic information related to the conditions required for growth of cells together with description of the metabolism of the organism based on its genome sequence. Conversion of this information into a mathematical framework, as a model, is the next step. Further, the behavior of the model can be predicted under specific growth conditions along with biomass generation and exchange (uptake or release) of relevant compounds (metabolites). In the final step, experiments are carried out to validate the predictions (Devoid et al., 2013). Genome-scale modeling is considered advantageous as it allows screening of multiple conditions in a short time and lowering of cost by providing a limited set of solutions that can be further tested. Solutions that enhance a desired or improved behavior of the organism can be predicted based on media supplementations without any requirement of genetic modification (García-Jiménez et al., 2018).

The importance of these models as a tool for rapid biodegradation and for the development of effective bioremediation strategies is recognized (Scheibe et al., 2009) and demonstrated in a few studies. For example, a metabolic model of Ralstonia eutropha H16, a microorganism with multiple potential applications in environmental biotechnology, was used to analyze its metabolic characteristics (Park et al., 2011). The genome-scale model of Pseudomonas veronii 1YdBTEX2, a degrader of monoaromatic compounds, was used together with ‘omics data (transcriptomics and exometabolomics) to characterize active metabolic processes at different growth stages. The matching of predictions with bacterial response indicates the importance of modeling to enhance the success rate of bio-inoculants in complex systems (Hadadi et al., 2020). Such studies help to predict the changes in cellular processes under a specific given condition (e.g., xenobiotic biodegradation). In addition, such integrative “omics modeling” approaches have also been used for optimization of predefined processes (Calmels et al., 2019). Xu et al. (2019) used models constructed for indigenous species whose abundance in soil was affected by exposure to the herbicide atrazine. The community modeling revealed interspecies metabolic interactions that support enhanced growth and degradation. The simulations pointed at a contribution of non-degrader species to the in situ process achieved thorough trophic exchanges with the degrader (Faust, 2019). The effect of specific metabolites on the degradation rate of atrazine by Paenarthrobacter aurescens TC1 was predicted by a genome-scale model and was further validated in vitro (Ofaim et al., 2020). The predictive biology facilitated the identification of optimal conditions for degradation that can promote the development of bioaugmentation or biostimulation approaches for better controlling the degradation processes. One bottleneck of bioaugmentation is that in many cases, the exogenous degrading microorganisms fail to establish in a new environment and do not produce enough biomass to eliminate contaminants in situ. Stimulating the growth of native degraders through adding specific compounds to the contaminated site is a practical alternative to bioaugmentation. The present study aims to elucidate the effect of nutritional supplements on linuron degradation by Variovorax SRS 16 (as a model for native degraders) and estimate their potential usage as biostimulants in order to harness the full potential of indigenous degraders. SRS 16 was selected due to the availability of a fully sequenced genome and extensive biochemical knowledge of its full pathway in linuron degradation – a pollutant that belongs to a large contaminant group of the phenyl urea herbicides. A genome-scale metabolic model of SRS 16 was constructed following automatic and manual protocols and was applied to explore the stimulating effect of different supplements through iterative simulations Compounds pointed by simulations as strong enhancers of degradation were further validated through wet laboratory experiments.

Materials and Methods

Constraint-Based Reconstruction

As a first step, a draft model of Variovorax sp. SRS 16 was reconstructed by following a bottom-up approach as elaborated by Ofaim et al. (2020). Briefly, the whole-genome sequence of SRS 16 was retrieved from the genome database of NCBI and annotated using RAST (Overbeek et al., 2014). An initial metabolic reconstruction was obtained by analyzing the annotated genome sequence through Model SEED (Faria et al., 2018). The metabolic reconstruction has three main components: a complete list of reactions and their associated genes, information of interactions between the gene and protein associated with the metabolic reaction, and the components involved in biomass generation (Devoid et al., 2013). The list of reactions contains the biomass (which includes all the biomass constituents and their fractional status), cytosolic, transport, and exchange reactions. The model was improved by adding potentially missing reactions through an automated gap-filling process (Henry et al., 2010). The improved draft model was exported to Systems Biology Markup Language (SBML) format from Model SEED.

The reactions in the original model were complemented by reactions described for Variovorax SRS 16 in additional genomic annotation databases including BiGG (Norsigian et al., 2020), KEGG (Kanehisa et al., 2016), UniProt (The UniProt Consortium, 2019), IMG (Chen et al., 2019), and MetaCyc (Karp, 2002). Reactions were included according to EC accessions. Reaction stoichiometries were verified for all reactions, and non-balanced reactions were fixed. Directionality and stoichiometry of reactions were determined according to the KEGG scheme. Specific reactions involved in linuron degradation were added based on a literature survey (Figure 1). The reconstruction was curated to verify that it allows the generation of all biomass components under physiologically feasible minimal conditions – minimal salt solution and linuron (Supplementary File 1). The Variovorax sp. SRS 16 model in SBML format is provided in Supplementary File 2.

FIGURE 1

Possible fate of linuron and 3,4-DCA in SRS 16. Reactions in black were included in the automatic reconstruction. Reactions in red were manually added to the model. Pathway 1 (modified ortho-cleavage pathway) was retrieved from BioCyc based on direct evidence for Variovorax SRS 16 strain (https://biocyc.org/META/NEW-IMAGE?type=PATHWAY&object=PWY-7496). Pathways 2 and 3 (ortho- and meta-cleavage pathways, respectively) were completed by following the generic KEGG pathways.

Simulations of Growth and Degradation of Linuron Using Flux Balance Analysis

Simulations were carried out using flux balance analysis (FBA), allowing us to depict cellular processes based on cellular reconstructions (Rana et al., 2020). Briefly, FBA follows the law of mass conservation and consider the metabolic framework as a static stoichiometric matrix (metabolites × reactions). FBA describes the predicted flux distribution through linear optimization by targeting specific cellular objectives (mainly growth). In the present study, the growth of the model was chosen as the objective function through the maximization of the biomass reaction. Flux variability analysis (FVA) was carried out to account for the possible flow of fluxes involved in secretion and uptake of all metabolites (Mahadevan and Schilling, 2003). The simulations were applied to predict biomass generation and linuron degradation over a time period under a range of defined conditions – in silico growth media. All simulations were carried out under definitions that follow experimentally verified viable conditions in minimal media with linuron (Supplementary File 3) and 120 exchange metabolites (one at a time) representing an alternative carbon/nitrogen source or other supplements (Supplementary File 4).

Dynamic modeling was used for the prediction of the profile of consuming metabolites typical to the biomass increase and linuron degradation across time. To this end, we simulated the behavior of our metabolic model across time. The simulation process follows (Xu et al., 2019; Ofaim et al., 2020) and is illustrated in Supplementary File 5. Briefly, the model works under the following assumptions: (1) a finite start amount of media components is available; (2) a maximal amount of uptake a single cell can acquire from the media in a given time point is defined (the lower bound of the exchange reaction value); (3) new substrate concentrations in each time point are determined by the predicted substrate concentration for the previous step augmented with any additional substrates provided or consumed in the current iteration. The maximum uptake was set to a ratio of up to 1 unit of each metabolite available in the media; (4) after each time tick, the biomass amount was updated according to the flux amount of the biomass reaction in the model at this time tick. As the biomass production rate serves as a proxy for the size of the population in the simulated environment and substrate uptake/secretion is mainly affected by population size, the model was used to evaluate the actual substrate uptake and growth rate given the supplied media across time.

Simulations were carried out until a state where time cycles did not lead to an increase in biomass was reached. Initial concentration values for all metabolites were set to a fixed amount of 50 units (represented as the initial lower bounds, LB, of the exchange reactions). The upper bound of the exchange reactions was set to 1,000 units to allow the secretion of all the exchange metabolites. Reversible non-exchange reactions’ lower bound was −1,000 units, and for non-exchange non-reversible reactions, it is 0 units; the upper bound for non-exchange reactions was 1,000 units. At each growth cycle, the generation of biomass is updated on the basis of flux flow in the biomass reaction. The algorithm assumes that media components are available to all growing cells with equal probability. The growth rate and the flux of substrate (consumed or secreted) associated with the model in the given media conditions are calculated at every time point, where, in each cycle, we try to maximize the biomass of each member of the community and then we fixate the biomass and minimize the uptake of metabolites of each member while maintaining the max biomass found. The cycles of predictions were carried out until the saturation point appeared in the biomass production (no growth recorded reflecting the exhaustion of given resources). Starting with one bacterial cell, the flux balance model was used to predict the uptake of carbon and nitrogen sources (including linuron) across time.

All model simulations were done on an Intel i7 quad-core server with 32 GB of memory, running Linux. The development programming language of our simulators was JAVA, and our linear programming software was IBM CPLEX.

In vitro Experiments of Variovorax sp. SRS 16 Strain: Growth and Linuron Degradation

Variovorax sp. strain SRS 16 (NCBI: txid282217, kindly provided by Dirk Springael, KU Leuven) was revived from stock cultures stored in glycerol at −80°C. The purity and authenticity of the strain were checked by 16S rRNA gene sequencing. The minimal medium (Supplementary File 3) described by Sørensen and Aamand (2003) was prepared, and bacterial cells were grown on agar plates at 25°C. For the quantitative analysis, 250 ml flasks were used to prepare the minimal medium (50 ml of the medium per flask) and autoclaved. The effects of the seven selected substrates – representing strong, moderate, and weak enhancers of linuron degradation – on growth and linuron degradation were tested. The supplements (filter sterilized) were added to a final concentration of 0.12 mM in media (equivalent to 30 ppm of linuron). The linuron concentration was selected on the basis of a previous study. MS, MS + C, and MS + N consist of only minimal medium, minimal media added with only carbon (glucose), and minimal media added with only nitrogen (ammonium salt), respectively. MS was treated as a negative control. The autoclaved minimal medium was supplemented with the substrates (separately) and linuron (30 ppm) followed by bacterial inoculation to a final OD (at 600 nm) of 0.05–0.1. The mother culture was raised in MSCN (succinic acid and ammonium salt) medium with linuron, by inoculating fresh agar plate-grown bacterial cells followed by incubation at 25°C (120 rpm) for 24 h.

All the inoculated flasks were incubated at 25°C (120 rpm) for 7 days. Bacterial growth and linuron degradation were monitored at definite intervals (zeroth day, third day, fifth day, and seventh day) from the day of inoculation. For bacterial growth, 200 μl was taken, and OD was measured at 600 nm by using Infinite® 200 PRO (Tecan Trading AG, Switzerland). Linuron degradation was measured through HPLC with standard procedures. Briefly, 1 ml of sample was taken and centrifuged at 10,000 × g for 5 min. The supernatant was filtered through a 0.22 μm PTFE syringe filter and transferred to the HPLC vials for the detection of residual linuron. Linuron and DCA were analyzed by using Agilent 1100 HPLC (Waldbronn, Germany). Detection of linuron was done at 240 nm using the external calibration method (sensitivity 0.5 mg L–1). The mobile phase of 70% methanol at a flow rate of 1 ml min–1 with a reverse-phased (Phenomenex, Torrance, CA) of 250 mm length × 4.6 mm inner diameter with particle size 5 μm was used for separation.

Linuron biodegradation and biomass buildup experiments were carried out in biological triplicates. The effect of supplements was statistically analyzed by performing repeated-measures ANOVA at p < 0.05 in SPSS v19.

Results

Reconstruction of Genome-Scale Metabolic Network for Variovorax sp. Strain SRS 16

The genome of Variovorax sp. strain SRS 16 is about 7.7 Mb in size, including its chromosome (5.7 Mb) and four plasmids with an approximate size of <1 Mb each (Öztürk et al., 2020). The complete sequence was retrieved from NCBI and annotated for reconstruction using the Model SEED pipeline. This initial reconstruction contained a list of 2,150 gene–protein–reaction associations that were classified as exchange, transport, and cytosolic as well as a list of all relevant metabolites and a biomass reaction. Based on the taxonomic classification, the biomass reaction was defined as Gram-negative bacteria in the Model SEED reconstruction pipeline. The composition of the biomass reaction summarizes the fractional contribution of generalized microbial biomass precursors (e.g., amino acids and lipids) to the synthesis of a new cell and is similar to the previously published genome-scale reconstruction of Escherichia coli strain K-12 (Monk et al., 2017). Initial simulations were carried out for debugging and removing futile or erroneously energy-generating loops. To this end, all external fluxes were blocked (upper and lower bounds set to zero). Next, a minimal medium was used, verifying that growth requires the supply of both carbon and nitrogen sources. After establishing no growth under infeasible conditions, we tested growth (biomass production) under experimentally verified conditions (succinic acid and ammonium salt as a carbon and nitrogen source, respectively) as described by Sørensen et al. (2009). Fine tuning of growth simulations to correctly represent the bacteria’s biology was done by manual curation. This included manual gap filling (addition of spontaneous and literature-supported reactions) and curation of reaction directionality.

The reconstructed metabolic network presented here covers 19.2% of the open reading frames (ORFs) present in the genome (Table 1). The coverage is lower than gold standard models of other Gram-negative models such as iAF1260 for Escherichia (Feist et al., 2010) (29%), iPC1209 for Pectobacterium (Wang et al., 2015), and iPAO1 for Pseudomonas (25.8%) (Zhu et al., 2018). It is, however, similar to coverage obtained for the genome-scale reconstruction of the closely related genus Rhodoferax ferrireducens (from the same family Comamonadaceae) – 15.6% (Risso et al., 2009). Total number of reactions is similar to gold standard reconstruction including the updated model of E. coli str. K-12 (Orth et al., 2011). A detailed description of the network including the reactions, metabolites, genes, and compartments that comprise the network is provided in Table 1. The model is also available as an SBML file (Hucka et al., 2003) in Supplementary File 2. The SBML file can be used with tools such as MATLAB or other SBML-compliant software. The minimal media used for simulations are available in Supplementary File 1. The set of metabolites, reactions, and exchanges are provided in Supplementary Files 1, 4.

TABLE 1

Functional details of the Variovorax sp. SRS 16 model.

Serial no.	Category	Total number
1	Metabolites	2,185
2	Exchange reactions	127
3	Reversible reactions	1,432
4	Transport reactions	277
5	Biochemical reactions	1,746
6	Total reactions	2,150
7	Genes associated with a reaction*	1,425 (total 7,411) (Öztürk et al., 2020)

Adapting and Curating iRZ1425 for Modeling Linuron Degradation

Model reactions involved in linuron degradation and their complete link to the core metabolism were added manually based on the detailed reports of the relevant pathways in the specific strain as well as in other bacterial species, forming three possible pathways (Figure 1). In Pathway 1, degradation occurs via modified ortho-cleavage, where in the first step, there is a ring cleavage in DCA; then, chloro intermediates transform into maleylacetate and finally enter the citric acid cycle. The pathway was described by Bers et al. (2011) based on the study of linuron degradation in Variovorax SRS 16 and is supported in part by molecular evidence. Notably, none of the reactions in the modified ortho-cleavage pathway was identified by automatic reconstruction. Alternatively, we also included two hypothetical pathways (not specifically detected for SRS 16), Pathways 2 and 3 (Figure 1), for the degradation of the herbicide hydrolysis product (Arora, 2015; Hussain et al., 2015). Pathway 2 follows the removal of successive dehalogenation and the formation of aniline from DCA. The product aniline enters the ortho-cleavage pathway through catechol to the core metabolism. Pathway 3 begins with the dehalogenation of DCA with the formation of 4-chloroaniline. 4-Chloroaniline is further converted to catechol derivatives, which leads to the meta-cleavage pathway. Some of the enzymes involved in these pathways were detected as part of the automatic reconstruction, and missing reactions were manually added (Figure 1). In the present study, simulations were carried out considering all three pathways, with Pathway 1 only and with Pathways 2 and 3 together; the same results were obtained when considering all the pathways together and with Pathway 1 only. The Variovorax sp. SRS 16 strain is documented to grow on linuron alone as it utilizes linuron as a sole carbon and nitrogen source in minimal media (Sørensen and Aamand, 2003). Sørensen et al. (2009) reported a strong enhancing effect of nitrogen (ammonium sulfate) on the degradation rate but not of carbon. These experimental observations were successfully captured by model simulations, indicating a stronger impact of nitrogen supplement (ammonium sulfate) over a carbon supplement (glucose) as enhancers of degradation. Simulation outcomes are consistent while using different pathway alternatives (Figure 2). In accordance with the documented mineralization of linuron (Sørensen et al., 2009), simulations indicate that linuron is converted into central-metabolism metabolites and biomass, where none of the degradation by-products (Figure 1) accumulate in the medium. The predicted production of biomass is fully correlated with linuron degradation (Supplementary File 6), providing additional evidence of full degradation.

FIGURE 2

Simulated linuron degradation in the modified media SRS 16 model considering Pathway 1 (A) and Pathways 2 and 3 (B). Simulations were carried out in minimal mineral media (MS) with linuron supplemented by nitrogen (ammonium sulfate) or carbon (succinate) sources (MS_N and MS_C, respectively).

Simulation-Based Predictions for Potential Enhancers of Linuron Degradation

The SRS 16 linuron degradation behavior under the influence of specific carbon (succinic acid) and nitrogen (ammonium salt) sources was extensively studied and reported (Sørensen et al., 2009). Here, we used simulations to further explore the effect of media supplements that can serve as potential degradation enhancers. We simulated growth in 120 different media combinations, each supplementing the linuron-containing minimal mineral media with a single exchange metabolite. The list of metabolites, following omission of 20 toxic and other non-relevant substances (considering future application in soil) such as diuron and dipeptides, is provided in Supplementary File 4. The effect of 28 exchange supplements selected across the full scale of degradation efficiency is shown in Figure 3. Selected metabolites were chosen to represent the biochemical diversity of carbon (simple sugars, organic acids, and biopolymers) and nitrogen (mainly amino acids, amine, and ammonium) sources that can act as future biostimulants in terms of regulation, costs, and accessibility. Simulations for all compounds were carried out considering both the SRS 16-specific pathway (Pathway 1, Supplementary File 7) and generic degradation pathways (Pathways 2 and 3, Figure 3). Different supplements were predicted to have variable impact on the growth of Variovorax SRS 16 and the degradation rate of linuron. In a reference medium containing linuron only (MS), 40% of the linuron was degraded at the 11th simulation round. The predictions stratify a group of enhancers with variable degrees of linuron degradation in the 11th iteration (>40–100%), in comparison to non-enhancers (40% linuron degradation, as in MS). All the carbon sources and several amino acids (arginine, lysine, methionine, cysteine, and tryptophan) were classified as non-enhancers vs. nitrogen sources that are predicted to expedite degradation at various degrees, with some variations depending on the simulation pathways used. Overall, considering the different pathway options, the rate of degradation broadly follows the same trend where most amino acids act as moderate enhancers, where on the 11th round, 50–75% of linuron is degraded. Fast enhancers include glutamine and asparagine (∼100% degraded in the final, 14th, round in comparison to 60% in MS). Though groups of enhancers vs. non-enhancers are similar when considering both simulation pathways (Figure 3 and Supplementary File 7), differences are observed in internal ranking as presented in Figure 4. The major difference found is that in Pathway 1, methionine and cysteine influence the degradation rate positively, whereas in Pathways 2 and 3, they are categorized as non-enhancers.

FIGURE 3

Predictions for linuron degradation by SRS 16 in minimal media supplemented with linuron and 28 compounds (selected carbon and nitrogen sources). Compounds are ranked according the linuron amount in the ninth simulation rounds. All media contain linuron. Reference medium: MS = minimal solution + linuron + no supplement, MP = 2-methyl propanamine. A total of seven metabolites (green font) were selected for experimental validation. Simulations assuming Pathways 2 and 3; simulations assuming alternative pathways are provided in Supplementary File 7.

FIGURE 4

Linuron degradation status (at round 10) under the influence of supplements with the ortho- and meta-cleavage pathway (Pathways 2 and 3) vs. the modified ortho-cleavage pathway (Pathway 1). Substrates in red were selected for further validation in the lab experiment. The x and y scales represent the amount of linuron (%).

In vitro Validation of Potential Enhancers of Linuron Degradation by Variovorax sp. SRS 16

Seven metabolites representing a range of predicted degradation enhancement potential and biochemical characteristics (Figure 3) were selected for validation. These metabolites include glutamine and asparagine as strong enhancers; aniline as a moderate enhancer; arginine, oxalic acid, and catechol as weak enhancers (Figure 3); and methionine whose enhancement potential depends on the choice of simulation pathway (non-enhancer in Pathway 1 and strong enhancer in Pathways 2 and 3, Figure 4). The simulations (growth and degradation) for the seven selected substrates assuming degradation Pathways 2 and 3 are shown in Figure 5. Simulations considering Pathways 1 and degradation Pathways 1–3 are shown in Supplementary File 8. Overall, the laboratory experiments supported predictions in the majority of the cases. Growth experiments (Figure 5B) are fully consistent with predictions, dividing supplement as non-enhancers vs. enhancers. In agreement with predictions, a slow (negligible) growth rate was found in MS and in MS supplemented by arginine, catechol, aniline, and methionine. A significantly higher growth rate was recorded in the presence of glutamine and asparagine (p < 0.05). Degradation results (Figure 5D) are overall consistent with growth results (Figure 5B), with significant stratification of enhancers (glutamine and asparagine) vs non-enhancers (p < 0.05). Aniline and methionine had a non-significantly stronger enhancement effect in comparison to MS and were classified into the non-enhancers group (p < 0.05 repeated-measures ANOVA), consistent with simulations in Pathways 2 and 3.

FIGURE 5

Predictions (based on Pathways 2 and 3) vs. observations for biomass and linuron degradation by Variovorax SRS 16. (A,B) Biomass predicted and observed. (C,D) Linuron degradation predicted and observed.

The results provide an insight that some supplements (associated with nitrogen metabolism) have a strong enhancing effect that helps in the conversion of linuron into potential cellular building blocks.

Discussion

Linuron is considered a stable phenyl urea herbicide. Under field conditions, the half-life of linuron varies from 30 to 150 days in different soil types, with an estimated average half-life of 60 days (Swarcewicz et al., 2013). DCA, its primary degradation intermediate, is reported to be more toxic. Therefore, complete mineralization of this compound is desired. Variovorax SRS 16, a member of a genus reported as an efficient degrader of phenyl urea herbicide, is widely known for its ability to mineralize linuron (Sørensen et al., 2009). Stimulating the growth of this degrader through adding compounds to a contaminated site is a potential strategy for enhancing its activity while bypassing the need of introducing exogenous microbes. Here, we describe the construction of a genome-scale metabolic model of Variovorax sp. strain SRS 16 aiming at finding model-based solutions for accelerating its linuron degradation activity. Construction followed automatic and manual protocols, and iterative simulations were applied for predicting potential compounds that promote the degradation of linuron. Predictions for significant degradation enhancers were confirmed by wet laboratory experiments.

Although, overall, the study demonstrates that model predictions can guide metabolic exploration toward predicting biochemical routes leading to enhanced degradation activity, several limitations of this work should be acknowledged. First, further adaptations are required for the current version of iRZ1425 to meet the requirements of gold standard models (Thiele and Palsson, 2010). Mainly, the biomass composition is an approximation that relies on automatic protocols and genomic information (Oh et al., 2007; Henry et al., 2010), similar to other recent works (Bordel et al., 2019; Naizabekov and Lee, 2020; Ofaim et al., 2020). The construction of a biomass reaction specific to the SRS 16 reaction requires layers of experimental data that are currently missing, including the precursors of general biomass constructions (for example, proteins, lipids, and species-specific components). Collection of currently lacking experimental data will be incorporated in future versions of the model to define an accurate biomass objective function using, for example, the BOFdat computational platform (Lachance et al., 2019). Second, though the current analysis demonstrates the use of model-based solutions for accelerating linuron degradation activity, the mechanism explaining the different effects of the metabolites tested remains unaddressed here. Simulations, supported by experimental evidence from the current study as well as previous reports (Sørensen et al., 2009), clearly point to nitrogenous compounds being advantageous in comparison to carbon ones. However, nitrogen content per se is insufficient for predicting the stimulating effect of compounds, and degradation efficiency cannot be ranked according to N content or to the C:N ratio. For example, the efficient stimulants glutamine and asparagine have C:N ratios of 5:2 and 4:2, respectively, in comparison to 6:4 in the non-efficient supplement arginine (Supplementary File 4). Simulation outcomes cannot be predetermined based solely on the C and/or N content of media supplements but rather by considering the complementary stoichiometric effect of multiple chemical reaction cascades, or pathways, that construct the genome-scale metabolic network. Similarly, simulation outcomes cannot be predicted based on compound structural properties, for example, different biostimulation potential of the two peripheral aromatic compounds, aniline and catechol.

Thus, the simulation predictions reflect the complexity of cellular processes and the balance between different metabolic pathways contributing to growth. The simulations provide multiple optimal solutions for the internal fluxes involved in central metabolism and hence are not sufficient for deciphering the mechanism by which environmental (media dependent) conditions affect linuron degradation. However, the integration of simulation data with ‘omics data (e.g., metabolomics, proteomics, and transcriptomics) can reduce the solution space and shed light on the fine-tuning of cellular activity, considering the need to balance between myriads of constraints. The model provides a platform for future sequential studies through the integration of metabolomic and transcriptomic data (similar to Hadadi et al., 2020), allowing the exploration of the mechanism of linuron degradation as well as other aspects of Variovorax sp. strain SRS 16 biochemistry. Such additional experiments are required in order to complement model predictions regarding the mechanism behind accelerating the phenomenon. For example, transcriptomics/proteomics studies targeting the effect of different media supplements (e.g., amino acids) would be beneficial to track the molecular processes considering different degradation scenarios.

Though current analysis requires additional layers on information in order to provide a description of the degradation mechanism, the analysis provides several insights that can be further explored in future studies. Model construction suggests three alternative routes for the degradation of linuron. Pathway 1 (modified ortho-cleavage pathway) relies on a detailed description of linuron degradation at strain level (Öztürk et al., 2020). Variovorax SRS 16 was reported to contain the gene coding for the amidase LibA, which hydrolyzes linuron in N,O-dimethyl hydroxylamine and DCA. The catabolic pathway for DCA degradation into components of central metabolism is suggested to involve conversion of 4,5-dichlorocatechol to oxoadipate (modified ortho-cleavage pathway) (Bers et al., 2011). Alternative degradation pathways of DCA were reported in other bacterial genera (Dejonghe et al., 2002; Breugelmans et al., 2010). In the current study, automatic construction based on the SRS 16 genome sequence suggested the partial presence of enzymes in chloroaniline or chlorocatechol metabolism (Pathways 2 and 3), where enzymes from the modified ortho-cleavage pathway (Pathway 1) were not detected. The missing/partial degradation pathways were completed through a gap-filling procedure (Figure 1). Here, the possibility of the presence of the ortho-cleavage pathway and meta-cleavage pathway along with the existing modified ortho-cleavage pathway is proposed based on genome data. The presence of partial Pathways 2 and 3 indicates that SRS 16 might be using these pathways for linuron mineralization under different conditions. Simulations in the SRS 16 model with complete Pathways 2 and 3 provided similar results to those obtained with Pathway 1 for most of the supplements (such as glutamine, asparagine, and aniline), with some exceptions such as methionine.

In the present study, it is shown that reconstruction is important for the study of bacterial physiology, allowing the simulation of degradation outcomes considering multiple conditions. It is well known that the degradation process is influenced by various physicochemical and nutritional factors (Kanissery and Sims, 2011). The availability of carbon and nitrogen sources (nutrients) significantly affects the efficiency of degradation (Yassir et al., 1998; Wu et al., 2011). Here, simulations suggest the importance of specific supplements to increase the degradation rate. Simulations predicted variations in the enhancement potential of different supplements, further supported by experimental validation. As a single cell is a complex network of reactions and metabolites, supplements can have a variable effect associated with the activation of different sets of reactions under each media condition. Integration of model predictions with omics data is useful for revealing the metabolic shifts in the microbial cell (Alam et al., 2010). To date, several strategies have been adopted in order to accelerate the biodegradation of herbicides in the environment, ranging from optimizing media to constructing transgenics (Azab et al., 2018; Li et al., 2020). Genome-scale modeling approaches along with omics technology have established their importance in biotechnology and medicine (Zhang and Hua, 2016). They are continuously appreciated for their role in the degradation of pollutants (Faust, 2019; Xu et al., 2019; Cardozo et al., 2020) and in other environmental problems related to agricultural soils (Mazzola and Freilich, 2017). Recently, Ofaim et al. (2020) demonstrated the use of genome-scale modeling for optimizing atrazine degradation by P. aurescens TC1. The current study demonstrated that genome-scale modeling allows the optimization of microbial activity, leading to herbicide biodegradation. The integration of predictive biology helped to screen the effect of the higher number of supplements on the degradation by SRS 16. Such supplement-derived amendment might reflect exchanges in the indigenous community where the degrader can be supported by non-degrader species (Xu et al., 2019) and accelerate the degradation rate.

In this way, computational biology can reduce cost, time, and effort (mainly related to gene manipulation) and save resources to identify the optimal solutions for a specific phenomenon. The association of genome-scale metabolic models with the omics approach can help to unravel the hidden facts of metabolism (Massaiu et al., 2019). Thus, there is a strong possibility of improvement in degradation processes through the contribution of modeling as an environmental biotechnology approach.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Materials, further inquiries can be directed to the corresponding author.

Author Contributions

KD, RZ, ZR, and SF designed the study and drafted the manuscript. RZ constructed the metabolic model. KD participated in model construction, led the design and analysis of computational and experimental results, and wrote the manuscript. SM participated in computational data analysis. NK and DB performed the wet laboratory experiments. BÖ, RA, and HE analyzed and reviewed the results. All authors carried out writing and improvement of manuscript.

Conflict of Interest

BÖ was employed by the Leibniz Institute DSMZ, German Collection of Microorganisms and Cell Cultures, Braunschweig, Germany. The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.