Open Reaction Systems Research Articles

Thermodynamic and stoichiometric laws ruling the fates of growing systems

We delve into growing open chemical reaction systems (CRSs) characterized by autocatalytic reactions within a variable volume, which changes in response to these reactions. Understanding the thermodynamics of such systems is crucial for comprehending biological cells and constructing protocells because it sheds light on the physical conditions necessary for their self-replication. Building on our recent work, where we developed a thermodynamic theory for growing CRSs featuring basic autocatalytic motifs with regular stoichiometric matrices, we now expand this theory to include scenarios where the stoichiometric matrix has a nontrivial left kernel space. This extension introduces conservation laws, which limit the variations in chemical species due to reactions, thereby confining the system's possible states to those compatible with its initial conditions. By considering both thermodynamic and stoichiometric constraints, we clarify the environmental and initial conditions that dictate the CRSs' fate—whether they grow, shrink, or reach equilibrium. We also find that the conserved quantities significantly influence the equilibrium state achieved by a growing CRS. These results are derived independently of specific thermodynamic potentials or reaction kinetics, therefore underscoring the fundamental impact of conservation laws on the growth of the system. Published by the American Physical Society 2024

Multiple sulfur isotopes evidence deep intra-slab transport of sulfate-rich fluids

Sub-arc mantle oxidation is often attributed to sulfates from seawater transported in aqueous fluids derived from slab dehydration, but sulfur speciation and origin(s) in high-pressure fluids remain weakly constrained. Here, we use the four sulfur stable isotopes (32S, 33S, 34S, 36S) to decipher processes, e.g. mixing of two end-member reservoirs and open-system reaction occurring during sulfur devolatilization, and deduce the nature and origin of sulfur species in aqueous slab-derived fluids. We analyzed high pressure metamorphic rocks from the Monviso massif Lower Shear Zone (LSZ; Western Alps), recording evidence of an interface-parallel, channelized fluid flow near peak burial metamorphic conditions (∼550 °C and 2.7 GPa). We sampled variably metasomatized Fe-Ti metagabbros, hybrid (talc) schists and serpentinites containing sulfides supposed to preserve the aqueous fluid composition, as well as other alpine metagabbro, metabasalt and serpentinites from greenschist to eclogite-facies devoid of the studied metasomatic overprint.The millimeter to centimeter-sized sulfide crystals found in abundance within LSZ lithologies are mainly pyrites containing inclusions of the peak burial mineral assemblage (e.g. garnet, omphacite, rutile) indicating their crystallization at peak conditions or during incipient exhumation. The sulfides in the non-metasomatic metagabbros and metabasalts cover a narrow range in δ’34S of 0.89 ± 0.63 ‰ and Δ’33S of 0.017 ± 0.01 (n = 6) that could be inherited from hydrothermal sulfides formed on the seafloor. In contrast, the LSZ sulfides from the Fe-Ti metagabbros, serpentinites and hybrid schists (n = 24) are enriched in 34S and 33S (bulk δ’34S from 4.91 to 21.32 ‰ and Δ’33S from 0.02 to 0.06 ‰) compared to the sulfides from the non-metasomatic rocks. Mixing between reduced sulfur species cannot explain all the data without requiring an end-member composition as high as δ34S ∼ 22 ‰ and Δ33S ∼ 0.06 ‰, i.e. values rarely recorded by sulfides in the literature. Instead, we propose that this 33S-34S-enrichment results from open-system sulfate reduction to sulfides in aqueous fluids at 550°C, with an initial δ34S ∼ 12 ± 3 ‰ and Δ33S ∼ 0.025 ± 0.010 ‰. Because this signature is different from seawater sulfate, we suggest that sulfate was produced by hydrothermal sulfide oxidation during antigorite breakdown at 650°C, in agreement with Mg, Ni and Cr enrichments in the metasomatized HP mafic lithologies. This model is consistent with relatively oxidized 34S-enriched arc lavas and can potentially explain recycled sulfur with positive δ34S signature recorded by ocean island basalts.

Mathematical Model for Chemical Reactions in Electrolytes Applied to Cytochrome c Oxidase: An Electro-Osmotic Approach

This study introduces a mathematical model for electrolytic chemical reactions, employing an energy variation approach grounded in classical thermodynamics. Our model combines electrostatics and chemical reactions within well-defined energetic and dissipative functionals. Extending the energy variation method to open systems consisting of charge, mass, and energy inputs, this model explores energy transformation from one form to another. Electronic devices and biological channels and transporters are open systems. By applying this generalized approach, we investigate the conversion of an electrical current to a proton flow by cytochrome c oxidase, a vital mitochondrial enzyme contributing to ATP production, the ‘energetic currency of life’. This model shows how the enzyme’s structure directs currents and mass flows governed by energetic and dissipative functionals. The interplay between electron and proton flows, guided by Kirchhoff’s current law within the mitochondrial membrane and the mitochondria itself, determines the function of the systems, where electron flows are converted into proton flows and gradients. This important biological system serves as a practical example of the use of energy variation methods to deal with electrochemical reactions in open systems. We combine chemical reactions and Kirchhoff’s law in a model that is much simpler to implement than a full accounting of all the charges in a chemical system.

Self-Organization of Enzyme-Catalyzed Reactions Studied by the Maximum Entropy Production Principle.

The self-organization of open reaction systems is closely related to specific mechanisms that allow the export of internally generated entropy from systems to their environment. According to the second law of thermodynamics, systems with effective entropy export to the environment are better internally organized. Therefore, they are in thermodynamic states with low entropy. In this context, we study how self-organization in enzymatic reactions depends on their kinetic reaction mechanisms. Enzymatic reactions in an open system are considered to operate in a non-equilibrium steady state, which is achieved by satisfying the principle of maximum entropy production (MEPP). The latter is a general theoretical framework for our theoretical analysis. Detailed theoretical studies and comparisons of the linear irreversible kinetic schemes of an enzyme reaction in two and three states are performed. In both cases, in the optimal and statistically most probable thermodynamic steady state, a diffusion-limited flux is predicted by MEPP. Several thermodynamic quantities and enzymatic kinetic parameters, such as the entropy production rate, the Shannon information entropy, reaction stability, sensitivity, and specificity constants, are predicted. Our results show that the optimal enzyme performance may strongly depend on the number of reaction steps when linear reaction mechanisms are considered. Simple reaction mechanisms with a smaller number of intermediate reaction steps could be better organized internally and could allow fast and stable catalysis. These could be features of the evolutionary mechanisms of highly specialized enzymes.

Open system reaction between hydrous melt and gabbroic rock in the Finero Mafic Complex

Open system reaction between hydrous melt and gabbroic rock in the Finero Mafic Complex

Chemical thermodynamics for growing systems

We consider growing open chemical reaction systems (CRSs), in which autocatalytic chemical reactions are encapsulated in a finite volume and its size can change in conjunction with the reactions. The thermodynamics of growing CRSs is indispensable for understanding biological cells and designing protocells by clarifying the physical conditions and costs for their growing states. In this paper, we establish a thermodynamic theory of growing CRSs by extending the Hessian geometric structure of nongrowing CRSs. The theory provides the environmental conditions to determine the fate of the growing CRSs; growth, shrinking, or equilibration. We also identify thermodynamic constraints: one to restrict the possible states of the growing CRSs and the other to further limit the region where a nonequilibrium steady growing state can exist. Moreover, we evaluate the entropy production rate in the steady growing state. The growing nonequilibrium state has its origin in the extensivity of thermodynamics, which is different from the conventional nonequilibrium states with constant volume. These results are derived from general thermodynamic considerations without assuming any specific thermodynamic potentials or reaction kinetics; i.e., they are obtained based solely on the second law of thermodynamics.

Open System Reaction between Hydrous Melt and Gabbroic Rock in the Finero Mafic Complex

Open System Reaction between Hydrous Melt and Gabbroic Rock in the Finero Mafic Complex

Arc-related peridotite blocks exhumed to the Eastern Block of the North China Craton prior to 2.47 Ga

AbstractThe Songling and Majiayu peridotite blocks occur in the Eastern Block, North China Craton (NCC). Geothermobarometry data show that the Songling peridotites began exhumation from a depth of c. 70 km (c. 23 kbar). During exhumation, the Songling peridotites were intruded by an upper-crust-derived, high-δ18O (up to +7.28‰ in zircon) trondhjemitic dyke at 2.47 Ga and experienced granulite-facies metamorphism. The Songling peridotites have hybrid mantle wedge (HMW) -like high SiO2 (> 45 wt%), high FeOt (c. 10 wt%) content, high modal orthopyroxene abundance (> 35%) and high ϵNd(t) (+18.4 to +21.4), which were generated by the reaction between peridotite and eclogite-derived melts. The clinopyroxenes from the Songling peridotites were in equilibrium with a Nb-, Zr- and Ti-depleted arc-like magma. The Majiayu peridotites are characterized by depletion of Nb, Zr and Hf and are highly enriched in FeOt, Th and light rare earth elements (LREEs), which can be interpreted as an open system reaction between hydrous melts and fore-arc mantle peridotites. These two peridotite blocks are considered to be arc-related mantle peridotites that experienced melt extraction and metasomatism in different styles. They were exhumed to the north margin of the North China Craton during the c. 2.47 Ga arc–continent collision along the Zunhua structural belt.

Carbon and hydrogen isotope fractionation during aerobic oxidation of short-chain alkanes in experimental incubations of vent fluids

Aerobic oxidation of short-chain alkanes was observed in gas samples from the Lutao intertidal hydrothermal vents in Taiwan, during storage at 20 °C for up to 29 months without adding bacterial strains and replenishing substrates. The carbon isotope fractionation factors (εC) of methane (C1), ethane (C2), and propane (C3), were calculated using the Rayleigh fractionation equation to be −37.1 ± 7.5‰, −14.8 ± 4.8‰, and −4.7 ± 5.2‰, respectively. The hydrogen isotope fractionation factor (εH) of methane was determined to be −281 ± 187‰. DNA sequencing of the 16S rRNA gene in the vent fluids suggests that aerobic oxidation is dominated by methanotrophs of the genera Methylomicrobium and Methylophaga, which use the ribulose monophosphate pathway (RuMP). The degrees of isotope fractionation (εC and εH values) herein are larger than previously reported values, possibly due to the limited O2 supply and low abundance of aerobic methane-oxidizing bacteria in the experiments. Since the fractionation factor of methane is higher than those of ethane and propane, the aerobic oxidation of thermogenic or microbial alkanes could produce a carbon isotope reversal, which is frequently noted as a trait of abiotic hydrocarbons. This work demonstrates that in addition to anaerobic microbial oxidation, aerobic oxidation with a low cell density can also produce significant isotope fractionation of alkanes in geological closed/semi-closed environments and open flow reaction systems that are characterized by moderate temperatures and a limited supply of substrates and O2; these environments include cold seeps, mud volcanoes, and low-temperature hydrothermal plumes/aquifers/reservoirs.

A fumarole in a one-pot: synthesis, crystal structure and properties of Zn- and Mg-analogs of itelmenite and a synthetic analog of glikinite

Anhydrous sulfate minerals are abundant in the active fumaroles with highly oxidizing conditions on the scoria cones of the Tolbachik volcano. The mineral itelmenite, ideally Na2CuMg2(SO4)4, containing isomorphous admixture of Zn, was described in 2018, whereas glikinite, ideally Zn3O(SO4)2, was described in 2020. Synthetic analogs of both minerals were obtained during studies of phase formation in the Na2SO4–CuSO4–MgSO4–(ZnSO4) systems which lead to essentially different results. Solid-state syntheses resulted in formation of several compounds previously known as minerals only. Both Zn- and Mg-containing analogs of itelmenite were prepared and exhibit slight deviations from the ideal Na2CuM2(SO4)4 stoichiometry. The Mg compound could be prepared single-phase which allowed the study of its thermal expansion and IR spectroscopy. Na2CuMg2(SO4)4 and Na2CuZn2(SO4)4 were evaluated for Na+-ion diffusion. For the Zn compound, several by-products were observed which are synthetic analogs of puninite Na2Cu3O(SO4)2, as well as hermannjahnite CuZn(SO4)2 and glikinite-type (Zn,Cu)3O(SO4)2. All of them were prepared via solid-state reactions in open systems. The Na2CuMg2(SO4)4, Na2CuZn2(SO4)4 and (Zn,Cu)3O(SO4)2 were structurally characterized by the single-crystal XRD. In the Zn-bearing system, the admixture of Cu2+ likely controls the formation of itelmenite-type and glikinite-type phases. The results of the experiments allowed to deduce possible scenarios of the formation processes of itelmenite and some other endemic fumarolic minerals. Our study shows that outstanding mineralogical diversity observed in the fumaroles of the Tolbachik scoria cones is not only due to the formation from the gas enriched by transition metals and involves also intensive exchange with the host basaltic scoria. Similar processes seem also to be responsible for the recrystallization of many other mineral species observed in high-temperature fumaroles resulted from the recent eruptions.

Stoichiometric network analysis of entropy production in chemical reactions.

We use stoichiometric network analysis (SNA) to obtain a flux-based approach for the evaluation and description of the entropy production and exchange for chemical reactions in open systems. For non-equilibrium stationary states (NESS) the production and exchange are expressed as functions over the convex cone of the stationary reaction rates, revealing the reaction pathways and elementary flux modes (EFM) responsible for both entropy production and balance. The analysis of the entropy production of EFMs leads to a unique description of the contribution of the coupling between linear or cyclic reaction network paths, and the fluxes due to the matter exchange of the system with the system's surroundings. Network stoichiometry leads to an independent proof and confirmation of Prigogine's theorem of minimum entropy production for the linear regime of non-equilibrium thermodynamics. Moreover, the non-linear thermodynamic regime allows us to test the validity of the General Evolution Criterion (GEC) for NESS in isothermal chemical networks in mechanical equilibrium.

HISTORY OF CHEMISTRY AS A PART OF ASSESSMENT OF STUDENTS’ UNDERSTANDING OF THE LAW OF CONSERVATION OF MASS

The aim of this research was to explore students’ ideas about chemical reactions and difficulties in understanding the law of conservation of mass in such reactions by using an approach that started from presentations of scientists’ work associated with the law. The developed test items relied on: 1) the historical contents that illustrate the experimental work of three scientists (Lavoisier, Landolt and Lomonosov); 2) the description of school experiments and 3) real life situation. In this way, students would have an opportunity to show understanding of the law of conservation of mass in two contexts, one based on the stories from the history of chemistry and the other contemporary, based on school laboratory experiments and real life situation. Students of different ages were selected for the research: the seventh and the eighth grade of primary school (age 13–14), and the second year of grammar school (age 16). The research involved a total of 301 students. The results indicated that students’ difficulties were mostly associated with the predictions and explanations of mass changes in open reaction systems in which a gas was a reactant than with the reactions in which a gas was a product. Keywords: assessment in chemistry, grammar school, history of chemistry, experimental results, law of conservation of mass.

Wegscheider’s condition and passivity of open chemical reaction systems

The second law imposes a number of constraints on the set of feasible reaction rate coefficients that a chemical network can accept, known in the literature on CRNT (Chemical Reaction Network Theory) as Wegscheider’s condition. One main consequence is that any closed chemical reaction system operating at isothermal conditions, if compatible with thermodynamics, must be complex balanced, otherwise the Wegscheider’s condition is violated. In this contribution we make use of this observation to conclude passivity of open reaction systems with respect to a Gibbs energy function, what in turn allows us to point on mass and energy exchanges with the environment, as the main causes that make chemical reactors unstable. As it will be discussed, control of open reaction systems can take advantage of this fact to design input-output stabilizing feed-back controllers.

Clumped isotope effects during OH and Cl oxidation of methane

A series of experiments were carried out to determine the clumped (13CH3D) methane kinetic isotope effects during oxidation of methane by OH and Cl radicals, the major sink reactions for atmospheric methane. Experiments were performed in a 100L quartz photochemical reactor, in which OH was produced from the reaction of O(1D) (from O3 photolysis) with H2O, and Cl was from photolysis of Cl2. Samples were taken from the reaction cell and analyzed for methane (12CH4, 12CH3D, 13CH4, 13CH3D) isotopologue ratios using tunable infrared laser direct absorption spectroscopy. Measured kinetic isotope effects for singly substituted species were consistent with previous experimental studies. For doubly substituted methane, 13CH3D, the observed kinetic isotope effects closely follow the product of the kinetic isotope effects for the 13C and deuterium substituted species (i.e., 13,2KIE=13KIE×2KIE). The deviation from this relationship is 0.3‰±1.2‰ and 3.5‰±0.7‰ for OH and Cl oxidation, respectively. This is consistent with model calculations performed using quantum chemistry and transition state theory. The OH and Cl reactions enrich the residual methane in the clumped isotopologue in open system reactions. In a closed system, however, this effect is overtaken by the large D/H isotope effect, which causes the residual methane to become anti-clumped relative to the initial methane. Based on these results, we demonstrate that oxidation of methane by OH, the predominant oxidant for tropospheric methane, will only have a minor (∼0.3‰) impact on the clumped isotope signature (Δ13CH3D, measured as a deviation from a stochastic distribution of isotopes) of tropospheric methane. This paper shows that Δ13CH3D will provide constraints on methane source strengths, and predicts that Δ12CH2D2 can provide information on methane sink strengths.

Residual water in hydrous minerals as a kinetic factor for omphacite destabilization into symplectite in the eclogites of Vårdalsneset (WGR, Norway)

Symplectitic intergrowths of sodic plagioclase+diopside±amphibole as replacement of omphacite are commonly found in eclogites. These symplectites are interpreted as the exhumation-related decompression of eclogite into the amphibolite facies. The role of aqueous fluid in symplectite development, which would act as a catalyst and favor open-system reaction, has been suggested but not yet clearly established. In the Vårdalsneset outcrop of the Western Gneiss Region (Norway), eclogites that were not amphibolitized either display a primary eclogitic dry paragenesis (garnet+omphacite+rutile±quartz) or a paragenesis including phengite. In the last case, omphacite is partly transformed into symplectite. The two groups have been further distinguished from a combined petrological, geochemical, and thermochemical study. Group I samples have a fine-grained unaltered microstructure, with medium Al2O3 (14–16wt.%), high Fe2O3 (13–16wt.%) and TiO2 (1.4–2.4wt.%). Group II samples have a coarse-grained microstructure and are characterized by the presence of symplectites and phengite. They display higher Al2O3 (17.5–23wt.%) and lower Fe2O3 (5.5–8wt.%) and TIO2 (0.2–0.5wt.%) contents than Group I. The P–T estimates for samples from both Groups I and II lead to similar conditions for peak eclogite metamorphism: temperatures range from 590 to 720°C and pressures from 15 to 25kbar. Perple_X modeling indicates that for Group I eclogites, temperature range is similar to the temperature range of the water saturation curve, whereas for Group II eclogites, due to slightly different chemical composition, the water saturation curve is located at much higher temperatures (770–900°C), so that OH− remains in residual phengite at peak eclogite temperature. With no residual hydrous phase in the eclogite assemblage, and although some structural water may persist in nominally anhydrous minerals, Group I eclogites were preserved without change during exhumation. In contrast, Group II eclogites retain hydrous minerals, and this hydrous component is available as a kinetic factor for symplectite formation during exhumation. This suggests that symplectite development can be promoted by an internally derived hydrous component, signaled by the presence of a hydrous mineral, and, so, does not necessarily require an influx of external fluid.

Introduction of NAT-testing for infections for screening blood donors in Republic of Kazakhstan

Aim. To evaluate the effectiveness of NAT-screening (nucleic acid amplification technologies) for infections in blood donors in Kazakhstan. Methods. Statistical data of blood donors screening examinations in the Republic of Kazakhstan in 2012-2014 were evaluated. Results. In 2014, the number of examined donors increased by 3.4% compared with 2012. The number of deferrals due to positive screening results for serological markers decreased by 10.9%, while the share of such donors decreased by 13.8% [p 0.01; odds ratio (OR) - 0.86, 95% confidence interval (CI 95%) - 0.83-0.88); χ2=136.76]. In 2014, 100% of donations were screened using NAT-testing (312,510 donors). Most of the NAT-screening in Kazakhstan is performed using closed automated systems. In 2012, 1 Blood Center conducted a polymerase chain reaction screening by open circuit polymerase chain reaction systems, in 5 blood centers polymerase chain reaction was performed with manual sample preparation. In 2014, the number of deferrals due to positive NAT-testing results has increased by 44.3%, the share of such donors - by 38.7% (p 0.01; OR=1.39, 95% CI=1.15-1.67); χ2=11.82). Seronegative NAT-positive samples were discovered according to the results of discriminant test, including human immunodeficiency virus - 2 (0.8%) samples, hepatitis B virus -182 (73.4%), hepatitis C virus - 60 (24.2%), negative result - 4 (1.6%). Conclusion. The introduction of screening NAT-testing of donated blood prevented transfusion of blood infected with: human immunodeficiency virus - 1 in 150,000 donations, hepatitis B virus - 1 in 1.650 donations, hepatitis C virus - 1 in 5,000 donations.

Diffusion-controlled metamorphic reaction textures in an ultrahigh-pressure impure calcite marble from Dabie Shan, China

The overall metamorphic reaction in an impure calcite marble from the ultra-high pressure (UH P ) zone at Changpu, eastern Dabie Shan, China, is partitioned into local reaction domains defined by isolated grains of omphacite, garnet, epidote, quartz and ilmenite within the calcite matrix. These reaction domains witness different stages of the metamorphic evolution. Chemically homogeneous omphacite in some instances has partial or complete coronas of quartz which are interpreted as a prograde phenomenon, i.e . caused by replacement of a plagioclase precursor by omphacite with increasing pressure. Rare relics of talc in epidote indicate that the first breakdown reaction during exhumation was talc + omphacite + garnet = amphibole + epidote ± quartz (still within the eclogite facies). During this early retrogression stage, omphacite was partially or completely rimmed by chemically heterogeneous amphibole. Most garnets developed a thick corona of epidote plus some quartz and calcite during this early retrogression or are completely pseudo-morphed. Where garnet occurs near or next to omphacite, this corona consists partly or entirely of amphibole. The second retrograde stage is commonly characterized by the breakdown of omphacite + quartz to albite + diopside + amphibole, often in the form of symplectite. Localized reactions become evident as omphacite at its site can decompose to three distinct types of symplectite (albite + diopside, albite + amphibole, albite + diopside + amphibole). Matrix quartz grains develop rims of individual diopside crystals at distances up to 2–3 millimetres from symplectites; additional outer rims of albite develop adjacent to garnet or corona epidote. Matrix quartz remote from omphacite sites is locally transformed into albite, where an Al-bearing mineral (mostly epidote) is found nearby. Epidote in contact with quartz develops albite coronas as it decomposes into albite + minute grains of Ce-rich epidote (allanite). Primary epidote and corona epidote remote from any quartz develop a complex compositional zoning with enrichment of mainly Ce 3+ at the outer rim, indicating exposure to a fluid undersaturated in alumina. Finally, large ilmenite grains in the matrix develop thick overgrowths of titanite, but in several instances an intermediate step of rutile formation is preserved. All these features indicate a direct correlation of textures with transport distances of dissolved chemical species in a non-pervasively infiltrating pore fluid. Simplified open-system reactions were derived for each site at each metamorphic stage and show that Al must have been quite mobile during prograde quartz corona formation, whereas both Si and Al were relatively immobile during the two retrograde stages – the relative sequence of mobility being Na + > Mg 2+ , Ca 2+ , Fe 2+ >> Si 4+ > Al 3+ , Ti 4+ . Mass-balance between sites is consistent and no significant gain or loss of cations is indicated.

Structures of chaos in open reaction systems

By numerically simulating the Bray-Liebhafsky (BL) reaction (the hydrogen peroxide decomposition in the presence of hydrogen and iodate ions) in a continuously fed well stirred tank reactor (CSTR), we find "structured" types of chaos emerging in regular order with respect to flow rate as the control parameter. These chaotic "structures" appear between each two successive periodic states, and have forms and evolution resembling to the neighboring periodic dynamics. More precisely, in the transition from period-doubling route to chaos to the arising periodic mixture of different mixed-mode oscillations, we are able to recognize and qualitatively and quantitatively distinguish the sequence of "period-doubling" chaos and chaos consisted of mixed-mode oscillations (the "mixed-mode structured" chaos), both appearing in regular order between succeeding periodic states. Additionally, between these types of chaos, the chaos without such recognizable "structures" ("unstructured" chaos) is also distinguished. Furthermore, all transitions between two successive periodic states are realized through bifurcation of chaotic states. This scenario is a universal feature throughout the whole mixed-mode region, as well as throughout other mixed-mode regions obtained under different initial conditions.

Thermodynamically consistent model calibration in chemical kinetics

BackgroundThe dynamics of biochemical reaction systems are constrained by the fundamental laws of thermodynamics, which impose well-defined relationships among the reaction rate constants characterizing these systems. Constructing biochemical reaction systems from experimental observations often leads to parameter values that do not satisfy the necessary thermodynamic constraints. This can result in models that are not physically realizable and may lead to inaccurate, or even erroneous, descriptions of cellular function.ResultsWe introduce a thermodynamically consistent model calibration (TCMC) method that can be effectively used to provide thermodynamically feasible values for the parameters of an open biochemical reaction system. The proposed method formulates the model calibration problem as a constrained optimization problem that takes thermodynamic constraints (and, if desired, additional non-thermodynamic constraints) into account. By calculating thermodynamically feasible values for the kinetic parameters of a well-known model of the EGF/ERK signaling cascade, we demonstrate the qualitative and quantitative significance of imposing thermodynamic constraints on these parameters and the effectiveness of our method for accomplishing this important task. MATLAB software, using the Systems Biology Toolbox 2.1, can be accessed from http://www.cis.jhu.edu/~goutsias/CSS lab/software.html. An SBML file containing the thermodynamically feasible EGF/ERK signaling cascade model can be found in the BioModels database.ConclusionsTCMC is a simple and flexible method for obtaining physically plausible values for the kinetic parameters of open biochemical reaction systems. It can be effectively used to recalculate a thermodynamically consistent set of parameter values for existing thermodynamically infeasible biochemical reaction models of cellular function as well as to estimate thermodynamically feasible values for the parameters of new models. Furthermore, TCMC can provide dimensionality reduction, better estimation performance, and lower computational complexity, and can help to alleviate the problem of data overfitting.

Thermodynamically consistent Bayesian analysis of closed biochemical reaction systems

BackgroundEstimating the rate constants of a biochemical reaction system with known stoichiometry from noisy time series measurements of molecular concentrations is an important step for building predictive models of cellular function. Inference techniques currently available in the literature may produce rate constant values that defy necessary constraints imposed by the fundamental laws of thermodynamics. As a result, these techniques may lead to biochemical reaction systems whose concentration dynamics could not possibly occur in nature. Therefore, development of a thermodynamically consistent approach for estimating the rate constants of a biochemical reaction system is highly desirable.ResultsWe introduce a Bayesian analysis approach for computing thermodynamically consistent estimates of the rate constants of a closed biochemical reaction system with known stoichiometry given experimental data. Our method employs an appropriately designed prior probability density function that effectively integrates fundamental biophysical and thermodynamic knowledge into the inference problem. Moreover, it takes into account experimental strategies for collecting informative observations of molecular concentrations through perturbations. The proposed method employs a maximization-expectation-maximization algorithm that provides thermodynamically feasible estimates of the rate constant values and computes appropriate measures of estimation accuracy. We demonstrate various aspects of the proposed method on synthetic data obtained by simulating a subset of a well-known model of the EGF/ERK signaling pathway, and examine its robustness under conditions that violate key assumptions. Software, coded in MATLAB®, which implements all Bayesian analysis techniques discussed in this paper, is available free of charge at http://www.cis.jhu.edu/~goutsias/CSS%20lab/software.html.ConclusionsOur approach provides an attractive statistical methodology for estimating thermodynamically feasible values for the rate constants of a biochemical reaction system from noisy time series observations of molecular concentrations obtained through perturbations. The proposed technique is theoretically sound and computationally feasible, but restricted to quantitative data obtained from closed biochemical reaction systems. This necessitates development of similar techniques for estimating the rate constants of open biochemical reaction systems, which are more realistic models of cellular function.