DOI: 10.5061/dryad.3j9kd51ss
title: Data and code from: Thermodynamics underpinning the microbial community-level nitrogen networks

description: Nitrogen species often serve as crucial electron donors or acceptors in
 microbial catabolism, enabling the synthesis of adenosine triphosphate
 (ATP). Although theoretically any nitrogen redox reactions could be an
 energy source, it remains unclear why specific reactions are predominantly
 utilized. This study evaluates energetically superior reactions from 988
 theoretically plausible combinations involving 11 nitrogen species, oxygen
 gas, hydrogen ion, and water. Our analysis of the similarity between this
 model-based energetically superior network and the actual microbial
 community-level nitrogen network, reconstructed as a combination of
 enzymatic reactions, showed increased link overlap rates with
 thermodynamic weighting on reaction rates. In particular, existing
 microbial reactions involving solely nitrogen species and additionally
 oxygen, such as anaerobic ammonia oxidation (ANAMMOX) and complete and
 partial nitrification, were frequently identified as energetically
 superior among the examined reactions. The alignment of these reactions
 with thermodynamically favorable outcomes underscores the critical role of
 thermodynamics not only in individual metabolic processes but also in
 shaping the broader network interactions within ecosystems, consequently
 affecting biodiversity and ecological functions.
description: All codes were written in the Wolfram Language platform using
 Mathematica 12.
description: # Hierarchical Structure within the Nnetwork Directory ```
 yourdirectory/Nnetwork │ ├── Data │ ├── DfG0.csv │ ├── N_reactions.csv │
 ├── N_reactions_not_half.csv │ └── REDOXlist │ ├── REDOXlist.csv │ ├──
 REDOXlist_Text.csv │ ├── REDOXlist_kdef.csv │ ├── divdata_all.csv │ └──
 edges_all.csv │ ├── Programs │ ├── 0_Make_Redox_reaction_List.nb │ ├──
 1_Preparation.nb │ ├── 2_err_response.nb │ ├──
 3_time_constrained_response.nb │ ├── 4_SS_export.nb │ ├──
 5_SS_visualize.nb │ ├── Fig_2a.nb │ ├── Fig_3b.nb │ ├── Fig_S2.nb │ ├──
 Fig_S3.nb │ ├── Functions_X.nb │ ├── Functions_1.nb │ ├── Functions_2.nb │
 └── Functions_3.nb │ └── Results └── kdef_XX └── Iin_YY └── p_ZZ ``` --- #
 Description of Each Directory and File --- ## Data --- All CSV files in
 this directory are used as input data. These files are processed by
 `../Programs/0_Make_Redox_reaction_List.nb`. --- #### N\_reactions.csv *
 Each line contains a string of half-reactions. It is acceptable to include
 overall redox reactions and blank lines. * If the first character of a
 line is "#", it is treated as a comment. * Insert a space
 between stoichiometric numbers and chemical formulas; spaces should also
 be inserted before and after the "+" sign and the
 "&lt;-&gt;" sign, which denotes reversible reactions. Even if
 the stoichiometric number is 1, it should not be omitted but explicitly
 noted. Failure to adhere to these rules will result in errors during the
 execution of `0_Make_Redox_reaction_List.nb`. *
 `0_Make_Redox_reaction_List.nb` extracts only those half-reactions that
 include the electron denoted by "e" from this file and returns
 all possible overall redox reactions. Reactions that do not include
 "e" are ignored. * The script first splits the reaction formula
 at the "&lt;-&gt;" position. Each split segment is further
 divided at the "+" positions, creating a list of sets of
 stoichiometric numbers and chemical species. * Any string can be used as a
 symbol for chemical species; however, using symbols like "+",
 "-", "=" will cause an error. If it is necessary to
 include valence information in the string, using notations like
 "NH4p1" is recommended. * The symbols for chemical species are
 used to reference the standard Gibbs energy change of formation for that
 species from the "DfG0.csv" file. Therefore, the strings for
 chemical species in this file and "DfG0.csv" must be identical.
 #### N\_reactions\_not\_half.csv * Add reactions that are not derived from
 the half-reactions in `N_reactions.csv` but involve nitrogen species
 dynamics to this file. Leave blank if there are no such reactions. * The
 notation in the file should follow the same guidelines as
 `N_reactions.csv`. #### DfG.csv * This table summarizes the chemical
 species appearing in both `N_reactions.csv` and
 `N_reactions_not_half.csv`, along with their standard Gibbs energy change
 of formation. Enter the string of the chemical species in the first
 column, their state in the second column, and the numerical value of the
 standard Gibbs energy of formation in the third column. Columns beyond the
 third are not used in processing and can be used freely for remarks and
 other notes. --- #### REDOXlist --- This directory stores CSV files
 exported from `../Programs/0_Make_Redox_reaction_List.nb`. The files
 `REDOXlist_kdef.csv` or `REDOXlist.csv` within this directory are imported
 and processed by `../Programs/1_Preparation.nb` for further analysis. ---
 * ##### REDOXlist.csv * This file contains detailed reaction information
 formatted for processing in the Wolfram Language, as outlined in the table
 below. Columns 4, 12, and 13 contain additional data that is not utilized
 in this study. This list contains partial blanks because some reactions do
 not involve electron transfer and therefore lack information related to
 electron exchange. | Column | Description | | ------ |
 ------------------------------------------------------------ | | 1 | Redox
 or acid-base reaction | | 2 | Number of electrons transferred (if
 applicable) | | 3 | ΔrGº (kJ/mol-N) | | 4 | ΔrGºe (kJ/mol-e) (if
 applicable) | | 5 | Forward reaction rate | | 6 | Reverse reaction rate |
 | 7 | Forward reaction rate constant | | 8 | Reverse reaction rate
 constant | | 9 | ΔrG | | 10 | Edge of forward reaction | | 11 | Edge of
 reverse reaction | | 12 | Proportion of each nitrogen species in forward
 reaction flow | | 13 | Proportion of each nitrogen species in reverse
 reaction flow | * #### REDOXlist\_kdef.csv * This file modifies
 `REDOXlist.csv` by replacing columns 7 and 8 with `kdef` (a constant) to
 simulate scenarios where reaction rates are unaffected by thermodynamic
 factors, using identical rate constants for all reactions. This list
 contains partial blanks because some reactions do not involve electron
 transfer and therefore lack information related to electron exchange. *
 #### REDOXlist\_Text.csv * Converts `REDOXlist.csv` into more standard
 notation for chemical equations and functions. The table below details the
 format and content. This list contains partial blanks because some
 reactions do not involve electron transfer and therefore lack information
 related to electron exchange. | Column | Description | | ------ |
 ----------------------------------------- | | 1 | Redox and acid-base
 reaction | | 2 | Oxidation reaction (if applicable) | | 3 | Reduction
 reaction (if applicable) | | 4 | Number of electrons moved (if applicable)
 | | 5 | ΔrGº (kJ/mol-N) | | 6 | ΔrGºe (kJ/mol-e) (if applicable) | | 7 |
 Forward reaction rate | | 8 | Reverse reaction rate constant | | 9 | ΔrG |
 * #### divdata\_all.csv * Contains `divdata` for 988 reactions, which
 includes detailed reaction dynamics data. For function specifics of
 `divdata`, see `../Programs/Functions_0.nb`. * #### edges\_all.csv *
 Details the proportion of each nitrogen species involved in the flow of
 988 reactions. Although this file is created, this file is not utilized in
 the current study. --- ## Programs --- This directory contains Mathematica
 notebooks written in the Wolfram Language. These notebooks were created
 using Mathematica 12. Each notebook is designed for specific steps or
 processes, and should be executed in numerical order as indicated by the
 numbers at the beginning of the file names. --- ###
 0\_Make\_Redox\_reaction\_List.nb * This notebook creates each file in
 `../Data/REDOXlist`. 1. Specify the absolute path of the Nnetwork
 directory as the string `dir` and execute the `Functions_0.nb` where
 functions are defined. 2. Import `N_reactions.csv`, combine half-reactions
 comprehensively to obtain a list of all reactions (`listall`), and
 calculate and associate the standard Gibbs energy change for each
 reaction. 3. Import `N_reactions_not_half.csv` and add all reactions and
 their standard Gibbs energy changes to `listall`. 4. Retain only one of
 any duplicate overall reactions that have the same number of transferred
 electrons and delete the others. 5. Determine the reaction rates for each
 reaction. Refer to the Methods section in the manuscript. 6. Export the
 resulting files to the `Data/REDOXlist` directory. ### 1\_Preparation.nb
 Perform numerical differentiation under specified conditions `kdef`,
 `Iin`, and `p`. Set the reaction rate constant k to kdef when `p =
 "def"`. In this case, import
 `../Data/REDOXlist/REDOXlist_kdef.csv`; otherwise, import
 `../Data/REDOXlist/REDOXlist.csv`. The numerical differentiation uses
 Mathematica's `NDSolve` function. 1. Create the output directory for
 generated files: `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin]
 value]/p_[p value]/Model` Here, the placeholders represent the Log10
 transformed values of each parameter, or the value of `p`. 2. Export
 `model.csv` (differential equation model), `init.csv` (initial values),
 and `cond.csv` (parameter conditions). 3. Using these files, perform
 parallel computations of numerical differentiation with the NDSolve
 function under specified O2 conditions from ωlist and H+ conditions from
 hlist. If no errors occur, the results are exported to:
 `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin] value]/p_[p
 value]/O2_[O2 value at calculation time in log10]_H_[H+ value at
 calculation time in log10].csv` If an error occurs, create an empty file
 named: `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin] value]/p_[p
 value]/err/O2_[O2 value at calculation time in log10]_H_[H+ value at
 calculation time in log10].csv` If the execution of NDSolve exceeds the
 `timelimit` seconds, the calculation is interrupted, and no results are
 exported. ### 2\_err\_response.nb Execute this notebook when CSV files
 exist in the directory `../Results/kdef_[Log10[kdef]
 value]/Iin_[Log10[Iin] value]/p_[p value]/err` to retry numerical
 differentiation. If no CSV files are present, this notebook can be
 ignored. 1. Retrieve the O2 and H+ conditions from the error CSV files in
 the directory. 2. Retry the numerical differentiation under these
 conditions with a modified AccuracyGoal. 3. If no errors occur, export the
 results to: `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin]
 value]/p_[p value]/O2_[O2 value at calculation time in log10]_H_[H+ value
 at calculation time in log10].csv` and delete the corresponding error CSV
 file. 4. Continue retrying until all error CSV files are processed and
 deleted. ### 3\_time\_constrained\_response.nb This notebook searches for
 and retries numerical differentiation for any missing O2 and H+ condition
 files where the execution of the NDSolve function was previously
 interrupted. 1. Retry the numerical differentiation under the missing
 conditions with a modified AccuracyGoal. 2. If no errors occur, export the
 results to: `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin]
 value]/p_[p value]/O2_[O2 value at calculation time in log10]_H_[H+ value
 at calculation time in log10].csv` Continue the process until there are no
 more missing condition files. ### 4\_SS\_export.nb Compile the
 steady-state concentrations from the files located at
 `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin] value]/p_[p
 value]/O2_[O2 value at calculation time in log10]_H_[H+ value at
 calculation time in log10].csv`. Calculate the steady-state Gibbs energy
 changes of reactions (-ΔrG), reaction rates (r), and work (W = -ΔrG * r).
 1\. **Specify H+ Condition**: * Select any H+ condition from the hlist
 (list of H+ conditions used for calculations). * Access the file
 `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin] value]/p_[p
 value]/O2_[log value in ωlist]_H_[log value in hlist].csv`. * Retrieve the
 last line and sort it by {O2 condition, tmax time nitrogen species
 concentration (e.g., N2, N2H4, etc.)}. * Export the steady-state
 concentrations for the specified pH to `../Results/kdef_[Log10[kdef]
 value]/Iin_[Log10[Iin] value]/p_[p value]/0Summary/SS/H_[specified H+ log
 value].csv`. 2\. **Calculate reaction rate and Gibbs energy change at
 steady state**: * Calculate the reaction rate (ri) and the Gibbs energy
 change (∆rGi) for the i-th reaction in steady state under each O2 and H+
 condition: {{r1, ∆rG1}, {r2, ∆rG2}, ..., {r988, ∆rG988}}. * Export these
 calculations to `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin]
 value]/p_[p value]/0Summary/SS/RRGibbs/O2_[O2 log value]_H_[H+ log
 value].csv`. * Ensure alignment of rows 2 to 495 with the rightward
 reactions and rows 496 to 990 with the leftward reactions as per the
 `Data/REDOXlist/REDOXlist.csv`. 3\. **ri Selected Edges**: * Arrange the
 data by reaction rate ri and its associated chemical reaction formula i,
 along with the directed edge and stoichiometric contributions for each. *
 Sort the list in descending order of reaction rate and extract the top θ
 reactions to form a sorted list of {chemical reaction formula',
 r', edge', edge_weight'}. * Output the sorted reaction data
 to `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin] value]/p_[p
 value]/0Summary/RateSelectedEdges/[θ value]/Reactions/O2_[O2 log
 value]_H_[H+ log value].csv` and the edges data to a corresponding edges
 file. 4\. **-∆rGi Selected Edges**: * Sort by Gibbs energy change ∆rGi in
 ascending order and select the top θ where the Gibbs energy change is
 negative. * Export as `../Results/kdef_[Log10[kdef] value]/Iin_[Log10[Iin]
 value]/p_[p value]/0Summary/GibbsSelectedEdges/[θ value]/Reactions/O2_[O2
 log value]_H_[H+ log value].csv` and the corresponding edges data. 5\. **W
 Selected Edges (when θ ≥ 1)**: * Calculate the work done Wi = ri * (-∆rGi)
 and arrange by work, sorting in descending order. * Extract and export the
 top θ elements as detailed reaction and edges files. 6\. **Power Selected
 Edges (when θ &lt; 1)**: * Calculate total work Wtot = ∑_i max(Wi, 0). *
 Obtain top reactions contributing to θ * Wtot and export the detailed
 reaction and edges files to the respective directories. ###
 5\_SS\_visualize.nb This notebook calculates edge overlap rates and
 generates various figures (e.g., Figures 1, 2b-d, 3a, 4). Note that the
 figures produced in this notebook may differ slightly from those in the
 published paper, as the eps files used in the paper were processed further
 after being generated by this notebook. 1. **Settings**: * Specify the
 directory containing the numerical results as `tmpdir` by setting
 parameters `p`, `Iin`, `kdef`, `θ`. "Other settings" include the
 definition of several variables and functions essential for calculations
 and the visualization of graphs. This setup enables the creation of both
 enzyme-level and community-level network graphs. 2. **Steady State
 Concentrations**: * Visualize the concentrations of each chemical species
 at steady state along with the combinations of major chemical species
 compositions. This visualization helps in understanding the distribution
 and dominance of chemical species under different conditions. 3. **Power
 Selected Edges**: * Retrieve and rank top reactions by work done, convert
 these reactions to edges, and calculate and visualize edge overlap rates.
 Additionally, visualize key reactions and the pH and O2 conditions under
 which these reactions were selected as key. 4. **Rate Selected Edges**: *
 Identify and rank top reactions by reaction rate, convert these reactions
 to edges for visualization, and calculate edge overlap rates. 5. **Gibbs
 Selected Edges**: * Select top reactions based on negative Gibbs energy
 changes, convert to edges, and visualize and calculate edge overlap rates.
 6. **Summary of Edge Overlap Rates**: * Arrange and present the figures
 showing edge overlap rates. ### Functions\_\*.nb This file compiles all
 the functions utilized across the various notebooks. For detailed
 explanations of specific functions and their applications, refer to the
 descriptions within the file itself. ### Fig\_\*\*.nb This notebook is
 dedicated to visualizing specific figures, denoted as Figure **, used in
 the paper. --- ## Results --- This directory contains the results
 generated from the programs whose filenames start with numbers from 0 to
 4. --- The structure of the Results directory is hierarchically organized
 based on the parameter conditions used in the calculations. It consists of
 directories named after the log10 values of the kdef and Iin parameters,
 as well as directories named after the p parameter. Within these
 directories, the results of numerical differentiation performed under
 various O2 and H+ concentration conditions are automatically saved as
 `O2_xx._H_yy.csv` files. Additionally, within the `p_ZZ` directories,
 other subdirectories serve specific functions (see `Programs` for
 details): * `0Summary` contains the summary of steady-state analysis
 results. * `err` temporarily stores files in cases where numerical
 computation errors occur. * `Model` stores the calculation conditions and
 the system of differential equations in Wolfram format, as described
 above. ---