A Mathematical Exploration of SDH-b Loss in Chromaffin Cells
Elías Vera-Sigüenza, Himani Rana, Ramin Nashebi, Ielyaas Cloete, Katarína Kl’uvčková, Fabian Spill, Daniel A. Tennant

TL;DR
This paper explores how chromaffin cells survive when a key enzyme (SDH-b) is lost, using a mathematical model to explain their unique metabolic adaptation.
Contribution
The study introduces a mathematical model to explain how chromaffin cells retain Complex I function despite SDH-b loss, distinguishing them from other cell types.
Findings
Retention of Complex I function in SDH-b deficient chromaffin cells is linked to cofactor oxidation.
Mitochondrial proton leaks and control of the proton gradient are critical for cell fitness in SDH-b loss.
The model supports the idea that chromaffin cells avoid lysis by managing mitochondrial swelling and ATP synthase reversal.
Abstract
The succinate dehydrogenase (SDH) is a four-subunit enzyme complex (SDH-a, SDH-b, SDH-c, and SDH-d) central to cell carbon metabolism. The SDH bridges the tricarboxylic acid cycle to the electron transport chain. A pathological loss of the SDH-b subunit leads to a cell-wide signalling cascade that shifts the cell’s metabolism into a pseudo-hypoxic state akin to the so-called Warburg effect (or aerobic glycolysis). This trait is a hallmark of phaeochromocytomas, a rare tumour arising from chromaffin cells; a type of cell that lies in the medulla of the adrenal gland. In this study, we leverage the insights from a mathematical model constructed to underpin the metabolic implications of SDH-b dysfunction in phaeochromocytomas. We specifically investigate why chromaffin cells seemingly have the ability to maintain electron transport chain’s Complex I function when confronted with the loss…
Genes, proteins, chemicals, diseases, species, mutations and cell lines named across the full text — each resolved to its canonical identifier and authoritative record.
Click any figure to enlarge with its caption.
Figure 10
Figure 11
Figure 12
Figure 13
Figure 14
Figure 15
Figure 16
Figure 1
Figure 2
Figure 3
Figure 4
Figure 5
Figure 6
Figure 7
Figure 8
Figure 9- —http://dx.doi.org/10.13039/100016778Paradifference foundation
- —http://dx.doi.org/10.13039/501100000289Cancer Research UK
- —http://dx.doi.org/10.13039/100014013UK Research and Innovation
- —Maria de Maeztu
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Taxonomy
TopicsCancer, Hypoxia, and Metabolism · Adrenal and Paraganglionic Tumors · Neuroblastoma Research and Treatments
Introduction
The succinate dehydrogenase (SDH) is a four-subunit complex (SDH-a, SDH-b, SDH-c, or SDH-d) and an essential component of cell central carbon metabolism. It has two main functions: first, it accepts an electron from succinate to produce fumarate within the tricarboxylic acid (or TCA), and second, it catalyses the reduction of ubiquinone to ubiquinol as the mitochondrial electron transport chain (ETC) complex II (CII) (Paredes et al. 2021; Baysal et al. 2000). As such, it functions as a critical link between the TCA cycle and ETC (Rutter et al. 2010; Paredes et al. 2021).
Many diseases and age related cell pathologies present SDH subunit mutations that cause partial to total loss of its functionality (Baysal et al. 2000). Needless to say, these metabolic disruptions result in catastrophic metabolic consequences (Kim et al. 2015; Moosavi et al. 2019; Huang and Millar 2013). However, the loss of its b subunit (SDH-b), in particular, is of special interest as it sets off a unique cascade of metabolic and signalling activities not seen in other mutations of this enzyme (Kim et al. 2015; Kluckova and Tennant 2018; Goncalves et al. 2021). Dysfunction of the SDH-b subunit begins with a significant increase in intracellular succinate concentration Kluckova and Tennant (2018); Goncalves et al. (2021); Li and Ye (2014); Letouzé et al. (2013); Kl’učková et al. (2020). Relatively high succinate levels lead to inhibition of the 2-oxoglutarate (2-OG)-dependent hypoxia-inducible factor (HIF) prolyl-hydroxylases, stabilising HIF-2 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\alpha $$\end{document} in normoxic (or normal oxygen tensions) conditions leading to a pseudo-hypoxic state with metabolic consequences resembling those seen in cancer described as the ‘Warburg effect’ (or the so-called aerobic glycolysis) a hallmark of phaeochromocytomas Kluckova and Tennant (2018); Goncalves et al. (2021); Li and Ye (2014); Letouzé et al. (2013); Kl’učková et al. (2020). These are rare tumours of the peripheral nervous system and adrenal glands, respectively (Goncalves et al. 2021; Li and Ye 2014; Letouzé et al. 2013; Selak et al. 2005; Brière et al. 2005).
Previous research has shown that in phaeochromocytomas, chromaffin cells (a type of endocrine cell located in the adrenal medulla in charge of producing and secreting catecholamines) preserve ETC complex I (CI) function despite the loss of SDH-b, unlike other cell types (Kluckova and Tennant 2018; Kl’učková et al. 2020; Cardaci et al. 2015; Lorendeau et al. 2017; Hart et al. 2023). SDH-b deficient chromaffin cells appear to have a unique capacity to remodel essential metabolic pathways and retain respiratory fitness under severe oxidative stress (Kl’učková et al. 2020). This adaptation is accompanied by considerable mitochondrial swelling, which is thought to be caused by the cell’s efforts to maintain a high membrane potential and ionic balance (Kolupaev et al. 2017). Furthermore, in SDH-b deficient chromaffin cells, the mitochondrial adenosine triphosphate (ATP) synthase (also known as ETC’s complex V) has been shown to function in reverse, hydrolysing mitochondrial ATP rather than synthesising it to preserve membrane potential at the expense of cellular ATP (Kl’učková et al. 2020; Acin-Perez et al. 2023). This is consistent with other mitochondrial diseases, in which it has been observed that blocking the mitochondrial ATP synthase can reverse its activity and restore cellular energy equilibrium (Acin-Perez et al. 2023; Valdebenito et al. 2023). Similar findings in phaeochromocytomas have shed light on the link between SDH-b loss, mitochondrial enlargement, and CI functionality in SDH-b deficient chromaffin cells. However, the exact mechanisms underlying these events remain unclear (Szarek et al. 2015; Wang et al. 2023; Dahia 2014; Pasini and Stratakis 2009; Amar et al. 2021).
In this study, we employ an in-silico dynamical model based on an immortalised mouse chromaffin cell (imCC) line model to investigate the metabolic consequences of SDH-b dysregulation in phaeochromocytoma (Letouzé et al. 2013; Kl’učková et al. 2020). Our research builds upon the findings of Kl’učková et al. (2020), who conducted extensive metabolic screenings on both wild type (WT) and imCC SDH-b knockout (K.O.) variants to underpin the biology behind loss of SDH-b. Here we demonstrate how we develop a robust computational framework able to predict metabolic shifts and their cell-wide consequences. We then characterise the causal mechanisms underlying the conclusions obtained by Kl’učková et al. (2020). Our main goal is to provide a first step towards a complete computational model able to explain why the loss of SDH-b activity, despite the enzyme’s widespread role in metabolism, predominantly impacts specific cell types such as chromaffin cells (Kl’učková et al. 2020; Brière et al. 2005; Fischer 2008; Fishbein and Nathanson 2012).
Materials and Methods
Cell Culture and Chemicals
Previously characterised immortalised mouse chromaffin cell lines deficient in SDH-b (SDH- \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {b}^{-/-}$$\end{document} CL6 and CL8) as well as their SDH- \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {b}^{+/+}$$\end{document} counterparts were maintained in Dulbecco’s Modified Eagle Medium (DMEM) supplemented with 10 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} fetal bovine serum (FBS) and 1 mM pyruvate. All chemicals, including DMEM and FBS, were obtained from Sigma-Aldrich unless stated otherwise.
Total Cell Protein and Cell Growth Evaluation
1 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\times 10^6$$\end{document} trypsinised cells were washed with PBS and lysed in 60 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l of RIPA buffer for 30 min. Protein concentration in the cleared supernatant was determined using the BCA protein method (Thermo Fisher Scientific).
For cell growth measurements, cells were suspended in 200 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l of Trypsin, followed by the addition of 400 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l of PBS to achieve a total volume of 600 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l. A 10 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l sample was then pipetted into a cell counting grid chamber (Fast Read 102, Kova International). After loading the chamber with the sample, cells distributed across the chamber squares were counted. The grid consists of 10 squares, each with dimensions of 1 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\times $$\end{document} 1 mm and a depth of 0.1 mm, corresponding to a volume of 0.1 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l per square. The formula for determining the cell concentration (cells/ml) is given by:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \frac{\textrm{N}_{\textrm{cc}}}{\textrm{ml}} = \frac{\sum \textrm{N}_{\textrm{cps}} \times \textrm{D}_\textrm{f} \times 10^4}{\textrm{N}}, \end{aligned}$$\end{document}where, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$N_{cc}$$\end{document} is the number of cells counted, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$N_{cps}$$\end{document} is the number of cells per square, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$D_f$$\end{document} is the dilution factor, and N = 5 the number of squares in the grid. The dilution factor ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$D_f$$\end{document} ) is 60 (Based on 600 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l/10 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l, and the conversion factor = \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$10^4$$\end{document} / 5.
Glucose, Lactate and Sodium Measurements
Media was collected from each well, cells were spun down and supernatant was taken for analysis. Levels were measured using a Contour XT glucometer (Bayer).
Metabolic Tracing Experiments
For tracing experiments, cells were plated to be 70 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} confluent after 48 h with 11 mM glucose and 2 mM glutamine, supplemented either in unlabelled form or as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}\hbox {C}_6$$\end{document} -glucose and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}\hbox {C}_5$$\end{document} -glutamine (CK Isotopes). After 48 h, 100 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l of media was removed for extraction and analysis. Cells were pelleted by centrifugation as described above. The remaining media was aspirated, and the empty wells were washed twice with ice-cold saline, after which 500 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l of MeOH was added. Cells were scraped and transferred to a cold Eppendorf tube. Subsequently, 500 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {D}_6$$\end{document} -glutaric acid in ice-cold water (1 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu $$\end{document} g/mL) was added (CDN Isotopes, D-5227) followed by 500 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l of chloroform (pre-chilled to \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$-20\,^\circ $$\end{document} C). After shaking on ice for 15 min and centrifugation, the polar phase was transferred to another tube for derivatisation, which was dried with centrifugation at \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$45\,^\circ $$\end{document} C.
Derivatisation and Gas Chromatography–Mass Spectrometry
Dried extracts were derivatised using a two-step protocol. Samples were first treated with 2 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} methoxamine in pyridine (40 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l, or 20 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l for primary samples, 1 h at \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$60\,^\circ $$\end{document} C), followed by the addition of N-(tert-butyl-dimethylsilyl)-N-methyl-trifluoroacetamide, with 1 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} tert - butyldimethyl-chlorosilane (60 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l, or 30 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l for primary samples, 1 h at \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$60\,^\circ $$\end{document} C). Samples were transferred to glass vials for GC-MS analysis using an Agilent 8890 GC and 5977B MSD system. One \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu $$\end{document} l of sample was injected in splitless mode with helium carrier gas at a flow rate of 1 mL per minute. The initial GC oven temperature was held at \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$100\,^\circ $$\end{document} C for 1 min before ramping to \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$160\,^\circ $$\end{document} C at a rate of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$10\,^\circ $$\end{document} C per minute, followed by a ramp to \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$200\,^\circ $$\end{document} C at a rate of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$5\,^\circ $$\end{document} C per minute, and a final ramp to \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$320\,^\circ $$\end{document} C at a rate of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$10\,^\circ $$\end{document} C per minute with a 5-min hold. Compound detection was carried out in scan mode. Total ion counts of each metabolite were normalised to the internal standard \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {D}_6$$\end{document} -glutaric acid.
Normalisation and Quantification
GC-MS data were analysed using Agilent Mass Hunter software for real-time analysis of data quality before conversion to.CDF format and analysis with in-house MATLAB scripts. Graphs and statistical analyses were performed using GraphPad Prism 9 and MATLAB.
Model Simulations and Code Availability
To partially quantify and parameterise the metabolic fluxes of our model, we employed \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C metabolic flux analysis ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-MFA). Briefly, this method deduces intracellular flux patterns from mass isotopomer distributions measured via mass spectrometry. The process and established protocols were followed as described in Vera-Siguenza et al. (2023); Antoniewicz (2018); Young (2014); Rahim et al. (2022). We conducted this analysis using the Isotopomer Network Compartmental Analysis (INCA) MATLAB routine, suitable for both steady-state and isotopically non-stationary metabolic flux analysis (Young 2014; Rahim et al. 2022).
All simulations and model code were executed using MATLAB and the ODE15s routine from MATLAB’s ODE suite. Source code for the models and figure generation can be freely accessed and obtained from our GitHub repository under an MIT open-source license (Shampine and Reichelt 1997; Dhooge et al. 2003).
Model Construction
Assumptions
Our model consists of a system of ordinary differential equations. Each equation monitors the rate of change in the concentrations of metabolites, critical ions, membrane potentials, and cellular and mitochondrial volumes with respect to time. Their rate of change is proportional to fluxes across four separate compartments: the extracellular space, cytoplasm, inner mitochondrial membrane space, and mitochondrial matrix, denoted by the subscripts e, c, i, and m, respectively (Fig. 1). Each flux in this model, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$j_{Flux}$$\end{document} , is described by a mathematical sub-model grounded in experimental observations or established mathematical concepts (Vera-Sigüenza et al. 2018, 2020, 2019; Keener and Sneyd 2009; Zhou et al. 2005; Salem 2003).Fig. 1. Schematic representation of the four compartments in the model: the extracellular space (e), cytoplasm (c), inner mitochondrial membrane space (i), and mitochondrial matrix (m). Each compartment is connected through various fluxes ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$j_{Flux}$$\end{document} ) that describe the transport and transformation of metabolites and ions. The extracellular environment (e) is assumed to have constant metabolite and ion concentrations. (Figure created with BioRender.com (Perkel 2020))
Our model is based on mass conservation, with each metabolite concentration or ion species considered to be spatially uniform (Keener and Sneyd 2009). As a result, a concentration change in any of the chemical species occurs simultaneously throughout the compartment. Additionally, all extracellular metabolite and ion concentrations are kept constant - this is evocative of a cell immersed in an infinite ion/metabolite solution. Although we recognise that this is a major simplification; for example, medium glucose and glutamine levels would decline with time, while lactate and other metabolites would increase during culture in-vivo, our modelling assumption is motivated by our need for computational efficiency and may more accurately represent the physiological situation with continuous flow of the interstitial fluid. Previous models, albeit making similar assumptions, have shown satisfactory predictive capacity (Vera-Sigüenza et al. 2018, 2020, 2019). Regardless of our assumptions, the model’s modularity enables easy adaption to alternative configurations and modes in future studies, accounting for these dynamic changes as needed.
Ion Channels and Non-metabolism Related Fluxes
We adopted a model developed by Vera-Sigüenza et al. (2018, 2020, 2019), which is based on the so-called Pump-Leak model (Keener and Sneyd 2009; Mori 2012). The ATP sodium/potassium ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Na}^{+}$$\end{document} / \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{K}^{+}$$\end{document} ) pump (NaK-ATPase), is crucial to this paradigm, as it maintains cellular volume against osmotic pressures at the price of cellular ATP levels (Keener and Sneyd 2009; Mori 2012; Kay and Blaustein 2019). This process also requires the transmembrane movement of ions, specifically chloride \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Cl}^{-}$$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Na}^{+}$$\end{document} , and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{K}^{+}$$\end{document} , which is assisted by separate membrane ion channels and co-transporters or secondary active transport (see Supplementary Material S.1) Keener and Sneyd (2009).
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Cl}^{-}$$\end{document} regulation, is essential to maintain cellular osmotic balance. We equipped our model with a \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Cl}^{-}$$\end{document} - \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Na}^{+}$$\end{document} - \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{K}^{+}$$\end{document} co-transporter (Nkcc) encoded by the Slc12a1/2 gene. In chromaffin cells, this co-transporter is crucial for maintaining elevated intracellular \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Cl}^{-}$$\end{document} levels to activate \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Cl}^{-}$$\end{document} -permeable GABA receptors (Xie et al. 2003). This action leads to a depolarised chloride equilibrium potential. The mathematical model we used was adapted from the original works of Vera-Sigüenza et al. (2018), Benjamin and Johnson (1997), and later corrected by Palk et al. (2010), and Gin et al. (2007). Furthermore, to adhere to the principle of mass balance, we introduced two generic efflux channels, one for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Cl}^{-}$$\end{document} and one for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{K}^{+}$$\end{document} (see Supplementary Material). These are modelled as simplified versions from those in Vera-Sigüenza et al. (2018), resembling those in Keener and Sneyd (2009) and Mori (2012) (Fig. 2-see Supplementary Material S.2).
As stated above, the primary mechanism responsible for osmotic balance, the NaK-ATPase, maintains the necessary \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Na}^{+}$$\end{document} electrochemical gradient and energises all secondary active transports (i.e., ion transporters and co-transporters that rely on concentration gradients) in the model (Fig. 2) (Keener and Sneyd 2009). The submodel we employed (see Supplementary Material) consists of a simplified version of the mathematical construct developed by Crampin et al. (2004) by Gin et al. (2007) and later by Maclaren et al. (2012). Our approach accounts for the ATP dependency of the pump, as detailed in studies by Vera-Sigüenza et al. (2018, 2019) and Palk et al. (2010) (see Supplementary Material S.3).
Our model also incorporates mechanisms for cellular pH regulation (Falkenberg and Jakobsson 2010). It includes a model for carbonic anhydrases-essential for catalysing the interconversion between \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {CO}_2$$\end{document} and water and the dissociated ions of carbonic acid- \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{HCO}^{+}_{3}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} -(see Supplementary Material S.4). The pH regulation system also includes the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Na}^{+}$$\end{document} -proton ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} ) antiporter (see Supplementary Material S.5), and a lactate- \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} symporter (described as part of the glycolytic sub-model of our model-as well as a succinate channel). As part of this pH regulatory system, we included a chloride-bicarbonate ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{HCO}^{+}_{3}$$\end{document} ) antiporter-encoded by the Slc4a2 gene (see Supplementary Material S.6). This ubiquitous molecular machine enhances \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Cl}^{-}$$\end{document} influx while expelling \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{HCO}^{+}_{3}$$\end{document} , following the principles outlined by Vera-Sigüenza et al. (2019, 2018), based on Falkenberg and Jakobsson (2010) (Fig. 2).Fig. 2. Schematic diagram of the cytoplasmic fluxes of the model. (Figure created with BioRender.com (Perkel 2020))
The cytoplasmic volume in our model changes as a direct consequence of the osmotic gradient between its neighbouring compartments: extracellular space and mitochondrial matrix (see Fig. 1). In this model (see Supplementary Material S.7), the osmotic gradient between the cytoplasm and the extracellular space influences the cytoplasmic volume, while the gradient between the cytoplasm and the mitochondrial matrix regulates mitochondrial volume (Vera-Siguenza et al. 2023; Vera-Sigüenza et al. 2018, 2020, 2019; Keener and Sneyd 2009; Mori 2012; Su et al. 2022, 2024). Note, however, that this assumes that the plasma membrane cannot withstand hydrostatic pressure gradients. The mathematical submodel we employed is based on and adapted from Kedem and Katchalsky (Vera-Sigüenza et al. 2019; Keener and Sneyd 2009; Palk et al. 2010; Jarzyńska and Pietruszka 2011).
Finally, to quantify the membrane potential ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V_{m}$$\end{document} ), established by the movements of charged ions across the cellular membrane, we employed Kirchhoff’s law via the so-called electric circuit model of the cell for a simple resistor-capacitor circuit (see Supplementary Material S.8) (Vera-Sigüenza et al. 2018, 2020, 2019; Keener and Sneyd 2009; Mori 2012; Su et al. 2022). This last addition assumes that the net sum of all currents in the circuit is zero. This has profound implications in our model as it essentially transforms it into a system of differential-algebraic equations (i.e., a system of equations that contains both differential and algebraic equations) (Mori 2012). The full mathematical descriptions and equations for each of the mechanisms depicted in Fig. 1 can be found in the Supplementary Material accompanying this article.
Glycolytic and TCA Cycle Model
Our glycolytic and tricarboxylic acid cycle (TCA) models are primarily based on the studies by Zhou et al. (2005) and Salem (2003) (Salem et al. 2002). Briefly, these pathways facilitate the delivery of carbon, derived from glucose, directly into the mitochondrial compartment via the glycolytic pathway (see Supplementary Material S.9–S.11). This pathway encompasses ten enzymatic reactions and is localised in the cytosolic compartment of the cell (Figs. 1 and 2). We included a sub-model that quantifies the reactions modulating the anabolism of pyruvate and lactate, and their subsequent incorporation into the cellular respiratory complex (Fig. 3).
Unlike the models proposed by Zhou et al. (2005) and Salem et al. (Salem 2003; Salem et al. 2002), we devised a net glycolysis reaction. Our simplification is driven primarily by the glucose influx facilitated by Glut1, a membrane glucose transporter encoded by the Slc2a1 gene and widely expressed in chromaffin cells, as well as by the production rate of pyruvate (Jóźwiak and Lipińska 2012; Toledo et al. 2013). Besides computational economy, this decision was motivated by our ability to obtain experimental data to accurately parameterise the entire pathway. This assumption allows us to focus on measurable values that can be directly obtained through experimental procedures, thereby allocating most of our computational resources to the electron transport chain. In other words, by aggregating the glycolytic reactions into a single net reaction, we streamline the model without compromising accuracy, focusing computational resources on solving the equations of the electron transport chain.
We acknowledge that this approach introduces potential model limitations, which we will address in the discussion section. However, our glycolytic model is able to calculate the rate of catalysis of carbon-source species as the difference between the rates of substrate production and utilisation, based on the availability of co-factors in the compartment. For our study purposes, this is sufficient (Rawls et al. 2019). The reaction flux sub-model is based on Michaelis-Menten kinetics and incorporates the reaction flux between metabolic species and their corresponding reaction stoichiometric coefficients (Kl’učková et al. 2020; Keener and Sneyd 2009; Lussey-Lepoutre et al. 2015; Pragallapati and Manyam 2019). Despite this, the model aligns well with experimental observations (Zhou et al. 2005; Kl’učková et al. 2020).
The TCA equations take on the following form:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} P_j - U_j = \sum _{k=1}^n \Big (\beta _{kj}- \phi ^p_{kj}\Big ) - \sum _{k=1}^m \Big (\beta _{jk}- \phi ^p_{jk}\Big ). \end{aligned}$$\end{document}Here,
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}&\phi _{kj}= \lambda _{jk}[{\textrm{met}}]_k,\\&\\&\lambda _{jk} = \frac{\textrm{RS}}{\textrm{RS}_0 + \textrm{RS}}, {\textrm{or}} \ \lambda _{jk} = \frac{\textrm{PS}}{\textrm{PS}_0 + \textrm{PS}}, \\&\\&{\textrm{where}}, \ {\textrm{RS}} = \frac{\textrm{NADH}}{\textrm{NAD}^+} \ {\textrm{and}} \ {\textrm{PS}} = \frac{\textrm{ADP}}{\textrm{ATP}}. \end{aligned}$$\end{document}Equation 1 details the production and utilisation of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$j\textrm{th}$$\end{document} metabolic species, which include all reactions resulting in anabolism or catabolism, are determined by the n reaction fluxes forming species j from species k \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\phi _{kj}$$\end{document} and the m reaction fluxes forming species k from species j \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\phi _{jk}$$\end{document} . The superscript p refers to the reaction processes in the TCA and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\beta _{kj}$$\end{document} are the reaction stoichiometric coefficients. The chemical reactions and their corresponding stoichiometries can be found in the supplementary data, a visualisation is provided in Fig. 3.
The rate coefficients \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\lambda _{jk}$$\end{document} are nonlinear functions of metabolite concentration ratios. These model the energetic state measured by ADP/ATP ratios, and the redox state by NADH/ \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {NAD}^+$$\end{document} ratios-cytosolic and mitochondrial alike (Salem 2003). In this light, a particular ratio is only included in the rate coefficients of reactions where they participate as co-substrate or co-product. While these are Michaelis-Menten models (Keener and Sneyd 2009), the concept of including these ratios has been successfully implemented and validated by Zhou et al. (2005) and Salem et al. (2002). It is important to note that \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\lambda _{jk}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\lambda _{kj}$$\end{document} denote the forward and reverse rate coefficients of a reversible reaction, respectively, and should not be defined separately in the context of reversible reactions.Fig. 3. Schematic diagram of the mitochondrial matrix reaction fluxes of the model associated with the glycolytic and TCA cycle model. (Figure created with BioRender.com (Perkel 2020))
Finally, we equipped the cytoplasm with the malate-aspartate shuttle reactions, which, together with lactate dehydrogenase, facilitate the redox conversion of NADH \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\leftrightarrow $$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {NAD}^+$$\end{document} . This sub-model depends on two generic mitochondrial transporters: one for aspartate and one for malate (see Supplementary Material S.12)(Keener and Sneyd 2009). We term them ‘generic’ for two reasons: first, they are not the primary focus of detailed study within our model; second, we lack detailed knowledge of the specific transporters responsible for these metabolite fluxes, such as their dependence on calcium ions ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {Ca}^{2+}$$\end{document} ) (Vera-Sigüenza et al. 2018, 2019; Palk et al. 2010; Su et al. 2022; Ruprecht and Kunji 2020). Nonetheless, these transporters are essential to this study and to understanding the effects on central carbon metabolism in SDH-b deficient chromaffin cells (Letouzé et al. 2013). Full mathematical descriptions and equations for each of the metabolic submodel mechanisms depicted in Fig. 2 can be found in the Supplementary Material accompanying this article.
Equations of the Electron Transport Chain Model
The electron transport chain (ETC) comprises four enzymatic complexes. Complexes I, III, and IV function as proton ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} ) pumps, moving \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} s from the mitochondrial matrix to the inner membrane space (Fig. 4). These utilise electrons from NADH and ubiquinone reduction, as well as oxygen, to pump \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} s against the established gradient. In contrast, complex II (or SDH) does not pump \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} s but instead reduces ubiquinone for complexes III and IV. These molecular machines are coupled to ATP synthase, which uses the established electrochemical gradient (by the aforementioned complexes) to pump \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} s back into the mitochondrial matrix, driving ATP synthesis (Fig. 4).Fig. 4. Schematic diagram of the inner mitochondrial membrane space reaction fluxes of the model associated with the electron transport chain model
The mathematical construct included in our model is largely based on the models by Beard (2005) and Manhas et al. (2020). However, we have modified some equations and parameters to ensure dimensional coherence with our data and the rest of the model. The conceptual model is briefly depicted in Fig. 4, with full details and equations provided in the Supplementary Material accompanying this article and the original works of Vera-Sigüenza et al. (2020), Beard (2005), and Manhas et al. (2020).
Proton Motive Force
The Proton Motive Force (PMF) represents the potential energy stored across a membrane. This force arises from the differential concentration of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} s on either side of the inner mitochondrial membrane (Fig. 4) (Johnson and Scarpa 1979). The PMF is essential for cellular energy production, as it is utilised by the electron transport chain to actively transport protons from the matrix to the intermembrane space, thereby establishing both a concentration gradient and an electrochemical gradient, also known as the membrane potential difference.
We adapted a model of PMF by Beard (2005) and modified it to include information from Johnson and Scarpa (1979). Briefly, the model relates the electrochemical gradient and the [ \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} ] gradient across the inner mitochondrial membrane proportional to the difference in mitochondrial membrane potentials \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta \Psi _m$$\end{document} (see Supplementary Material S.13). While the concentration gradient component reflects the variation in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} concentrations between the intermembrane space and the mitochondrial matrix, in our model, the combined effect of these gradients creates the energy reservoir used by ATP synthase (Fig. 4).
Electron Transport Chain Complex I
Complex I (C1) is a \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} pump situated across the inner mitochondrial membrane (Fig. 4), initiates the electron transport chain by facilitating the transfer of electrons from NADH to ubiquinone (Q), reducing it to ubiquinol ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$QH_2$$\end{document} ). This action energises the transport of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} against the established chemical gradient from the mitochondrial matrix to the inner membrane space, hence contributing significantly to the mitochondrial PMF (see Supplementary Material S.14).
The enzymatic flux of Complex I is modelled as proportional to the difference in concentrations of NADH and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {NAD}^+$$\end{document} , as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} fluxes drive the reaction. The flux is influenced by the free energy change (Gibbs’ free energy \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta G$$\end{document} ) associated with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} movement from one compartment to another.
Electron Transport Chain Complex III
Complex III mediates electron transfer from ubiquinol in the mitochondrial matrix to cytochrome c in the inner mitochondrial space (Fig. 4). Similar to Complex I, Complex III is a \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} pump and contributes to the mitochondrial PMF (see Supplementary Material S.15). The framework for this model is based on the design by Korzeniewski and Zoladz (2001).
In this model, the flux of Complex III is defined by the influence of phosphate ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$P_i$$\end{document} ) serves as a modulatory factor that drives respiratory activities to meet energy demands, a relationship first explored by Katz et al. (1989). Phosphate’s involvement underscores its significance in respiration or ATP synthesis reactions and highlights its potential impact on mitochondrial regulatory dynamics, including changes in volume.
Electron Transport Chain Complex IV
Similar to complexes I and III, Complex IV contributes to proton pumping (from mitochondrial matrix to inner mitochondrial space), which is instrumental in energising the ATP synthase. The model we use was first constructed by Korzeniewski and Zoladz (2001) (see Supplementary Material S.16).
The flux through Complex IV depends on the oxygen concentration at any given time. This has profound consequences for the model kinetics, which vary non-linearly. Complex IV is also influenced by the proportion of reduced cytochrome c in the matrix relative to its total amount. As the concentration of reduced cytochrome c increases, so does the flux through Complex IV, given that reduced cytochrome c serves as a substrate for the reaction.
Adenosine Triphosphate Synthase
The adenosine triphosphate (ATP) synthase plays a critical role in converting adenosine diphosphate (ADP) into ATP within the mitochondrial matrix. This model, constructed by Korzeniewski and Zoladz (2001), centers around \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} movement across the mitochondrial membrane and how it energetically drives the synthesis of ATP (see Supplementary Material S.17). The process depends on the concentration gradients between ADP and ATP, along with phosphate ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$P_i$$\end{document} ) in the mitochondrial matrix, ensuring the reaction is energetically favourable.
Adenine Nucleotide Translocator (ANT)
The adenine nucleotide translocator (ANT) flux plays a critical role in cellular energy management by facilitating the displacement of one negative charge from the mitochondrial matrix to the mitochondrial inner membrane space. In our model (see Supplementary Material S.18), this process is coupled to the electrostatic membrane potential and modelled as a membrane transporter according to Halestrap and Brenner (2003); Korzeniewski and Zoladz (2001); Korzeniewski (2000) and Keener and Sneyd (2009).
Magnesium-ATP Binding
Magnesium ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {Mg}^{2+}$$\end{document} ) plays a crucial role in the stability and function of ATP, which comprises a ribose sugar, adenine, and three negatively charged phosphate groups. By forming bonds with these phosphate groups, magnesium effectively reduces their inherent repulsion (Pasternak et al. 2010; Lacapère et al. 1990). This interaction (see Supplementary Material S.19) is essential for enzymes such as kinases that rely on ATP- \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {Mg}^{2+}$$\end{document} complexes to enhance catalytic efficiency. The inclusion of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {Mg}^{2+}$$\end{document} in our model coordinates the phosphates during ATP hydrolysis. It also serves as a preamble for a future exploration of the metabolic/signalling intersectionality.
The Mitochondrial Phosphate Carrier (PiC–Slc25a3)
The PiC, encoded by the solute carrier family 25a3 (Slc25a3), is a crucial protein within the mitochondria responsible for transporting phosphate across the inner mitochondrial membrane (Ruprecht and Kunji 2020; Seifert et al. 2015). In our model, it facilitates the movement of inorganic phosphate between the matrix and the mitochondrial intermembrane space, coupled to the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} gradient (see Supplementary Material S.20).
The transport mechanism model, based on that by Korzeniewski and Zoladz (2001), involves a cotransport process where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} and dihydrogen phosphate ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {H}_2$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {PO}_4$$\end{document} ) are moved together in a 1:1 ratio, allowing for an electroneutral exchange across the membrane. The association of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$P_i$$\end{document} is in equilibrium, effectively balancing the phosphate species on both sides of the membrane.
Adenyl Kinase
Adenylate kinase (AK) is an essential enzyme for maintaining cellular adenine nucleotide balance (Noma et al. 2001; Nobumoto et al. 1998). In the mitochondrial inner membrane space, AK catalyses the transfer of high-energy phosphates among ATP, ADP, and AMP, a process vital for cellular energy management. This reaction is modelled to proceed via a general linear equation, applicable across different isozymes (see Supplementary Material S.21). The model incorporates an equilibrium constant, X \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\textrm{AK}}$$\end{document} , and an enzyme activity parameter, X \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\textrm{AK}}$$\end{document} , to quantify the flux of nucleotides mediated by AK, effectively describing its role in the energetic equilibrium within cells.
Proton Leak and Potassium Fluxes
In this study, we assume that the effects of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\textrm{Ca}}^{2+}$$\end{document} concentrations and fluxes on membrane potential are secondary compared to the respiratory chain, adenine nucleotide translocator (ANT) current, and proton leaks. Consequently, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\textrm{Ca}}^{2+}$$\end{document} fluxes are not incorporated at this stage but are planned for inclusion in future iterations of the model. However, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{K}^{+}$$\end{document} and magnesium ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\textrm{Mg}}^{2+}$$\end{document} ) ions are integral to the model due to their roles in buffering matrix pH and facilitating ATP synthesis and ANT flux, respectively. The movements of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{K}^{+}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} across the mitochondrial membrane are modelled using the Goldman-Katz-Hutchkin equation, a solution derived from the one-dimensional Nernst-Planck equation (Keener and Sneyd 2009) (see Supplementary Material S.22).
Potassium/Proton Exchanger
The mitochondrial potassium/proton exchanger is vital for the transport of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{K}^{+}$$\end{document} into the mitochondrial matrix in exchange for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} . This exchanger is critical in maintaining the mitochondrial membrane potential and the pH gradient, both essential for ATP production via oxidative phosphorylation (Garlid 1980). In our model, the dynamics of this exchanger are captured using a linear exchange model as outlined by Keener and Sneyd (2009) (see Supplementary Material S.23).
AMP, ADP, and ATP
Transport of ATP, ADP, AMP, and Pi between the cytosol and the mitochondrial inter-membrane space is modelled using linear transfer between compartments (see Supplementary Material S.24).
Succinate Dehydrogenase/Electron Transport Chain Complex II
Our model is an adaptation from a model first constructed by Manhas et al. (2020). At the core of this model are the binding polynomials that describe the likelihood of various molecules binding to specific sites on the SDH. These polynomials arise as the denominator of a rational function that represents the average number of occupied binding sites as a function of ligand (substrate) activity (Keener and Sneyd 2009). This approach enables the capture of the mechanistic dynamic interaction of substrates and inhibitors with the enzyme, providing a quantitative measure of binding affinities and their impacts on its activity. The specific formulations of these binding polynomials highlight the interactions of ligands such as ubiquinone within the SDH (see Supplementary Material S.25).
Additionally, the model adjusts the redox potentials of various SDH-associated redox centres to account for pH variations, which can significantly impact the enzyme’s electron transfer capabilities (see Supplementary Material S.26). These corrections ensure that the model is able to attain physiological conditions more accurately, allowing for a better understanding of how pH shifts influence the redox state of the enzyme (Manhas et al. 2020).
The fluxes associated with SDH, including the transfer of electrons from succinates through various bound states (or SDH-b subunit) (see Supplementary Material S.27) to the eventual reduction of ubiquinone (SDH-c and SDH-d), are based on the established binding and redox potential models. This comprehensive model ensures that all critical aspects of SDH functionality, from substrate binding to electron transfer, are captured with high fidelity.
Model Equations
The full model equations can be found in the Supplementary Material S.28. Moreover, given the extensive list of parameters involved in this model we have provided, as part of the supplementary data, two files that contain all the parameters and their values - see Supplementary Material.
Experimental Data Integration and Model Parameter Fitting
To elucidate the complex metabolic interplay within chromaffin cells, we conducted an isotopic \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C labelling experiment. This method traces carbon atoms through key metabolic pathways, providing insights into substrate consumption, intermediate metabolite dynamics, and end-product generation (Antoniewicz 2018). Our main purpose was to quantify fluxes, parameterise reactions, refine our model, and validate our findings against experimental work (Kl’učková et al. 2020; Gimenez-Roqueplo et al. 2003; Goncalves et al. 2021). For this purpose, we utilised \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}\hbox {C}_6$$\end{document} -glucose (uniformly labelled glucose) due to its proven ability to comprehensively reveal carbon utilisation in central carbon metabolic processes, as assessed by gas chromatography-mass spectrometry (GC-MS) (Antoniewicz 2018; Young 2014). This technique, combined with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-metabolic flux analysis ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-MFA), enabled us to dissect and compare the metabolic fluxes between wild-type (WT) and SDH-b knockout (K.O.) imCC cells (Antoniewicz 2018).
Isotopic Labelling
Following the protocols outlined by Antoniewicz (2018) and Vera-Siguenza et al. (2023), we began by monitoring the proliferation rates of WT and SDH-b K.O. cell lines over a 48-h period (Fig. 5a). We observed that SDH-b K.O. cells had a significantly lower growth rate (0.021756 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {h}^{-1}$$\end{document} ) compared to WT cells (0.034426 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {h}^{-1}$$\end{document} ). This indicates a substantial reduction in growth for cells lacking SDH-b, aligning with findings from previous studies (Kl’učková et al. 2020; Goncalves et al. 2021). The growth rates were calculated by assuming cells continuously divide (Antoniewicz 2018). Thus, we expect exponential growth according to:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} N_x = N_{x,0} \ {\textrm{exp}}\Big (\mu t \Big ), \end{aligned}$$\end{document}where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$N_x$$\end{document} represents the number of cells and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu $$\end{document} is the growth rate (1/hr). Solving for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu $$\end{document} , we obtain:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \mu = \frac{\ln (N_{x,t_2}) - \ln (N_{x,t_1})}{\Delta t}. \end{aligned}$$\end{document}Fig. 5. Extracellular and growth measurements. a WT vs. SDH-b K.O. growth rate monitored over a 48 h period. Rate parameters were obtained assuming exponential doubling time. b WT vs. SDH-b K.O. glucose consumption comparison. SDH-b K.O. cells show a decrease in consumption due to impaired TCA cycle. c WT vs. SDH-b K.O. lactate efflux comparison. SDH-b K.O. cells show a marked increase in lactate efflux. This is indicative of a shift towards aerobic glycolysis. (For absolute concentration values in figures b and c refer to Supplementary material.)
Subsequently, we assessed glucose and lactate levels in the media to gauge metabolic activity. Over the same period, both WT and SDH-b K.O. cells demonstrated distinct glucose consumption patterns. We noticed that SDH-b K.O. exhibited reduced glucose intake compared to WT cells (Fig. 5b). However, SDH-b K.O. cells exhibited an increased lactate production rate, corroborating the expected shift towards aerobic glycolysis, a phenomenon typically observed in these cells (Fig. 5c) (Brière et al. 2005; Pasini and Stratakis 2009; Letouzé et al. 2013; Kluckova and Tennant 2018; Goncalves et al. 2021). Using these measurements, we were able to calculate the glucose and lactate external rates:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} j_{i} = \pm 1000 \frac{\mu V \Delta C_i}{\Delta N_x}, \end{aligned}$$\end{document}here \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta C_i$$\end{document} (mmol/L) is the change in concentration of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$i^{th}$$\end{document} metabolite between two sampling time points, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta N_x$$\end{document} is the change in cell number during the same time period, V (mL) is the culture volume (see Materials and Methods), and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu $$\end{document} (1/h) is the growth rate (Antoniewicz 2018). Note that in Eq. 4, consumption is defined negative and metabolite secretion is defined postive.
Following preliminary assessments, we undertook comprehensive isotopic labelling using \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}\hbox {C}_6$$\end{document} -glucose (details in Materials and Methods). We then analysed the isotopic labelling patterns of key tricarboxylic acid (TCA) cycle metabolites, including lactate, citrate, succinate, malate, fumarate, and aspartate (Fig. 6a). All simulations in our \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-MFA underwent sum of square residuals statistical test (SSR aims at recovering a sparse weight vector in an underdetermined linear model with a known and fixed dictionary matrix) to assess the goodness-of-fit, parameter confidence intervals, and model identifiability evaluation (Fig. 6b). This enhanced the reliability of our findings (Young 2014). The mass isotopomer distributions (MIDs) for these metabolites revealed distinct differences between SDH-b K.O. and WT cells, suggesting a potential rewiring of the TCA cycle in SDH-b deficient cells (Fig. 6a and Supplementary Materials).Fig. 6. Experimental vs. simulation mass isotopomer distributions from \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}\hbox {C}_6$$\end{document} -Glucose labelling experiment (WT vs. SDH-b K.O.). a \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-MFA predicted mass isotopomer distributions of several metabolites across the central carbon metabolism pathways. The x-axis represents the mass isotopomer states (m + 0, m + 1, m + 2, etc.), where m+n indicates the number of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ ^{13}\hbox {C}_6$$\end{document} -labelled carbon atoms incorporated into the metabolite. b Assessment of good-fit using sum of squares due to error (SSR), or the quantity which the least squares procedure attempts to minimise to attain the observed mass isotopomer distributions. The results obtained provide a good fit to experimental data.(see Supplementary Material)
We then performed \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-metabolic flux analysis ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-MFA) to obtain the fluxes for the central carbon metabolism reactions (see Supplementary Material). Our analysis revealed that fluxes from glucose to pyruvate remained consistent across WT and SDH-b K.O. cells, indicating unaltered glycolytic rates regardless of SDH-b status (Fig. 7a and b). Nevertheless, an increase in lactate production in SDH-b K.O. cells signalled enhanced lactate dehydrogenase activity, which is likely to be compensation compensating for TCA cycle disturbances resultant from SDH-b loss (Kay and Blaustein 2019; Kluckova and Tennant 2018; Goncalves et al. 2021).Fig. 7. Comparison of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}\hbox {C}_6$$\end{document} -glucose labelling experiment-derived fluxes (metabolic flux analysis). a imCC WT flux map. b imCC SDH-b K.O. flux map. The flux alterations in succinate-related reactions suggest a bottleneck, prompting cells to employ compensatory mechanisms such as a reversal of the malate-aspartate shuttle and increased reliance on pyruvate carboxylase over pyruvate dehydrogenase
Our simulations also indicate that SDH-b K.O. cells exhibit significantly altered fluxes in TCA intermediate steps, involving succinate, directly attributable to SDH-b loss (Fig. 7b). The flux alterations in succinate-related reactions indicate a bottleneck, prompting cells to employ compensatory mechanisms such as a decoupling of the malate-aspartate shuttle and increased reliance on pyruvate carboxylase over pyruvate dehydrogenase. A similar outcome was reported by Hart et al. (2023). As a consequence, our analysis suggests, a possible increased reliance on exogenous fatty acid pathways could be responsible for sustaining the loss of pyruvate dehydrogenase activity, which directly affects the citrate synthesis arm of the TCA cycle (e.g., lipid metabolism, fatty acid \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\beta $$\end{document} -oxidation and oxidative phosphorylation) Kluckova and Tennant (2018); Goncalves et al. (2021); Letouzé et al. (2013).
Enrichment and Rate Parameters
After obtaining the fluxes of the model (Fig. 7) through \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-Metabolic Flux Analysis (MFA), we proceeded to carry out an enrichment simulation (Antoniewicz 2018). This approach stems from the realisation that the resulting fluxomic data implicitly provides the rates at which enzymes involved in central carbon metabolism catalyse their respective reactions. By simulating enrichment, we could determine key kinetic parameters from our dataset, thereby enabling the model to recapitulate physiological traits of the chromaffin cell model (imCC) with greater accuracy (Young 2014).
An isotopic enrichment measurement quantifies the relative abundance of different isotopes of an element in a sample. For example, in metabolic studies, a sample might be enriched with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-glucose. By measuring the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C/ \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{12}$$\end{document} C ratio in different metabolites, it is possible to track the incorporation and flux of carbon through metabolic pathways, as well as enzymatic reaction rates (Elahee Doomun et al. 2016).
Given that most reactions in our central carbon metabolism model primarily follow Michaelis-Menten kinetics, we sought to determine the maximum enzymatic reaction rate ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V_{max}$$\end{document} ) and the dissociation constant ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$K_d$$\end{document} ) associated with substrate binding (Keener and Sneyd 2009). However, enrichment percentages derived from \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-labelling experiments do not directly equate to concentration measurements. Since our model deals with changes in concentration over time, a transformation was necessary.
To address this, we normalised the enrichment of the isotopomers by setting \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V_{max} = 1$$\end{document} . Using this approach, we were able to standardise our model equations and then deduce the corresponding \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$K_d$$\end{document} and Hill coefficients to align with enrichment rates. In this frame, we could estimate the rate of catalysis for different enzymes in our model.
Once these values were obtained, we applied a non-linear optimisation routine in MATLAB (Levenberg-Marquardt algorithm, from the optimisation toolbox-see Supplementary Material) to refine the estimation of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V_{max}$$\end{document} . Given that we had information on metabolite concentrations (derived from raw ion counts-see Supplementary Data), \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$K_d$$\end{document} , and flux from our \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-MFA simulation, we optimised the parameter values accordingly.
The reasoning behind equating rate parameters with enrichment rates is as follows: enrichment rates, which indicate the proportion of molecules in a metabolite pool labelled with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C, act as an indirect measure of metabolic flux through a pathway. As flux increases, so does the proportion of labelled molecules, similar to how \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V_{max}$$\end{document} represents the maximal catalytic capacity of enzymes and dictates the rate of metabolic reactions. By aligning the Hill function parameters with enrichment rates, we capture the interplay between substrate availability and enzyme activity, allowing the parameters-particularly \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V_{max}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$K_d$$\end{document} -to serve as proxies for enrichment rates.Fig. 8. Enrichment simulation. a and b display the mean enrichment curve of pyruvate over time, aligned with a glycolysis model. The graphs show the model’s fit to the normalised enrichment data, highlighting the model’s accuracy in representing pyruvate dynamics in glycolysis
Figures 8 and 9 illustrate the results of these simulations and parameter estimations. Here, enrichment percentages reflect the proportion of metabolites containing the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C label, indicating flux through metabolic pathways that incorporate \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C from our tracer into the metabolite of interest.
Figure 8a and b show the enrichment curve for pyruvate over time, with a model fit based on a glycolysis model. The curve depicts normalised enrichment over time, where the x-axis represents time (hours) and the y-axis represents normalised enrichment. Similar to Fig. 8b, the model fit aligns closely with the measured data, demonstrating the accuracy of the glycolysis model in capturing pyruvate enrichment.
We also used the enrichment data to characterise how enzymatic activity influences labelled citrate production over time (Fig. 9a and b). Here, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V_{max}$$\end{document} represents the maximum proportion of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}$$\end{document} C-labelled citrate production when substrate (acetyl-CoA and oxaloacetate) availability is at its highest and most effective for citrate synthase activity. Notably, the incorporation of carbon into the TCA cycle occurs over a significantly longer timescale compared to pyruvate metabolism in glycolysis.Fig. 9. Enrichment simulation. a and b depict how enzymatic activity affects the production of citrate over time - setting the incorporation of carbons into the TCA, with the curve indicating the peak rate of citrate synthesis by the citrate synthase abbreviated as (CS) under optimal substrate availability. c and d Depict the enrichment curve for fumarate, including both the original and adjusted model fits, highlighting the adjustments made to the SDH model parameters for better alignment with experimental data. panel d shows the estimated flux over time using the SDH-a model (imCC WT)
Figure 9c presents the enrichment curve for fumarate over time, alongside the original and optimised model fits. The original model by Manhas et al. (2020) deviated significantly from the measured enrichment data, prompting parameter adjustments that led to a more accurate representation of SDH activity. Figure 9d further depicts our estimation of SDH-a model (imCC WT) flux over time.
We repeated the simulations (see Supplementary Material) for all measured metabolites in our labelling experiment. This approach enabled us to systematically characterise the glycolytic and TCA cycle metabolic sub-models, reinforcing the reliability of our kinetic parameter estimations in capturing metabolic activity.
Model Results
Model Benchmarking
Having attained model parameterisation, we benchmarked the model to ensure its physiological relevance. This was carried out to ensure the model could reliably simulate known biological behaviours, thereby confirming its utility in predicting metabolic responses under well known conditions (Hester et al. 2011). As such, we began by evaluating the model’s response to oxygen deprivation (or hypoxia). Essentially, we simulated conditions where oxygen levels were decreased. The model, in turn, should display known significant metabolic shifts (e.g., changes in electron transport chain activity, alterations in ATP production, and variations in metabolite fluxes) (Buchholz et al. 2002; Yang et al. 2019).Fig. 10. Benchmarking simulation (dimensionless time \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau $$\end{document} ). We subjected the chromaffin cell model to hypoxia (a 30 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} reduction in oxygen flux). a Relative activity of the mitochondrial pyruvate carrier, showing a decrease as hypoxia progresses over time. b Relative activity of lactate dehydrogenase, illustrating an increase due to enhanced anaerobic glycolysis. c Relative glycolysis rate, maintaining steady activity despite hypoxic conditions. d Relative activity of the ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Na}^{+}$$\end{document} / \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} ) NHE antiporter. The antiporter increases its activity to remove \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} s from the cytosol, protecting cells from acidification as a result of anaerobic metabolism
The simulation results depicted in Fig. 10 align closely with established physiological responses to hypoxia (Kluckova and Tennant 2018; Lister 1993; Wheaton and Chandel 2011). In Fig. 10a it is possible to see the model’s response to a decrease in oxygen flux leads to a reduction in the flux through the mitochondrial pyruvate carrier. This is consistent with the known decrease in mitochondrial oxidative metabolism under low oxygen conditions (Bender and Martinou 1863; Ruiz-Iglesias and Mañes 2021). The model indicates that this adaptation helps conserve oxygen by reducing its utilisation in the mitochondria.
Figure 10b depicts an increase in lactate dehydrogenase activity-a well-documented reaction to hypoxia, where cells increase lactate production to regenerate \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {NAD}^+$$\end{document} for continued glycolytic ATP production Li and Ye (2014). In Fig. 10c, a slight upregulation in glycolysis activity further supports this metabolic shift, indicating an increased reliance on glycolytic pathways to meet energy demands when oxidative phosphorylation is compromised. Additionally, Fig. 10d demonstrates an increase in the activity of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {Na}^+$$\end{document} / \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {H}^+$$\end{document} antiporter, indicative of the cellular effort to regulate intracellular pH under acidic stress conditions commonly induced by increased lactate production. This response is critical, as the ability of the NHE to regulate acidity changes in the cytoplasm signifies the model’s correct response to hypoxia (Hulikova et al. 2013; Vera-Siguenza et al. 2023).Fig. 11. Benchmarking simulation (dimensionless time \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau $$\end{document} ). We subjected the chromaffin cell model to hypoxia (a 30 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} reduction in oxygen flux). a Acidification of the cytoplasm as a result of a metabolic shift towards anaerobic glycolysis. b The cell exhibits an increase in NADH concentration, due to the inability of the glycolytic pathway to rid lactate at the rate of normoxic respiration. c Increase in lactate concentration due to a shift to anaerobic glycolysis. d Decrease in pyruvate concentration is expected, most pyruvate is catabolised to generate lactate and thus maintain the cellular cofactor redox balance
The model also highlighted significant changes in intracellular metabolite concentrations and pH, reflecting typical cellular adaptations to hypoxia. Figure 11a illustrates a decrease in intracellular pH, which is a consequence of the accumulation of acidic metabolic by-products such as lactate. This acidification is a hallmark of the hypoxic response, often leading to the activation of pH-regulating mechanisms such as the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{Na}^{+}$$\end{document} / \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\textrm{H}^{+}$$\end{document} antiporter observed in Fig. 10d. Figure 11b shows an increase in intracellular NADH concentration, which is expected due to the reduced activity of the electron transport chain under low oxygen conditions, causing an accumulation of reduced cofactors (Wheaton and Chandel 2011). Figure 11c shows a rise in intracellular lactate concentration, further supporting the increased glycolytic activity and lactate production, while Fig. 11d indicates a reduction in intracellular pyruvate concentration, reflecting its increased conversion to lactate. These shifts in metabolite levels are characteristic of the metabolic shift that occurs during hypoxia, aiming to adapt to reduced oxygen availability while maintaining cellular energy balance.
Based on these results (Figs. 10 and 11), the model has demonstrated a robust ability to replicate well-known metabolic shifts associated with mild hypoxia (Wheaton and Chandel 2011; Lister 1993; Paredes et al. 2021; Zhou et al. 2005; Salem 2003). The observed changes in electron transport chain activity, lactate dehydrogenase activity, glycolysis, and intracellular pH align closely with documented physiological responses to hypoxic conditions (Bender and Martinou 1863; Wheaton and Chandel 2011; Ruiz-Iglesias and Mañes 2021).
SDH-b Knockout
Given the validated performance of our model, we proceeded with simulations involving SDH-b knockouts to investigate the resulting metabolic consequences in chromaffin cells. Our model simulations reveal a comprehensive picture of the metabolic consequences of SDH-b knockout in chromaffin cells, particularly focusing on the activity of electron transport chain (ETC) complexes.
As shown in Fig. 12a, our model demonstrates that upon SDH-b loss, Complex I activity is laregly maintained, suggesting a compensatory mechanism wherein the cell has the capacity to retain Complex I activity to maintain electron flow into the ETC despite the impaired function of Complex II. This finding aligns with previous studies indicating cellular adaptation strategies to sustain mitochondrial function under compromised conditions (Kluckova and Tennant 2018; Kl’učková et al. 2020; Goncalves et al. 2021; Pasini and Stratakis 2009; Letouzé et al. 2013). The activities of Complex III and IV exhibited substantial decreases, with Complex III activity dropping by 60 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} and Complex IV activity by approximately 70 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} (Fig. 12b and c). These reductions indicate that the loss of electron supply from ubiquinol, due to the absence of SDH-b function, severely impacts downstream complexes. Despite the compensatory upregulation of Complex I, it is insufficient to fully sustain the activities of Complex III and IV. This observation is consistent with the known interdependencies within the ETC, where disruption in one complex can significantly affect the overall electron transport efficiency (Liu et al. 2002).Fig. 12SDH-b knockout simulations (dimensionless time \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau $$\end{document} ). a Complex I activity is maintained after SDH-b loss, suggesting a compensatory upregulation to sustain electron flow into the ETC despite Complex II impairment. bComplex III activity decreases by approximately 40 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} , indicating reduced electron supply from ubiquinol due to the loss of SDH-b function. c Complex IV activity decreases by around 40 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} , further highlighting the impact on downstream complexes from the absence of SDH-b function. d ATP synthase activity reverses from 1 to −4, showing a severe stress response where ATP synthase hydrolyses cytoplasmic ATP to maintain mitochondrial membrane potential, depleting cellular ATP reserves
A significant observation was the reversal of ATP synthase activity, which shifted from 1 to −4 (Fig. 12d), as observed by (Kl’učková et al. 2020). This severe metabolic stress response indicates that ATP synthase hydrolyses cytoplasmic ATP to pump protons back into the mitochondrial matrix, a critical action to maintain the compartment’s membrane potential in the face of compromised ETC function (Kl’učková et al. 2020). This reversal comes at the cost of depleting cellular ATP reserves, highlighting the metabolic strain imposed by SDH-b knockout. To sustain this proton flux, the model relies heavily on proton leak mechanisms, which balance the activities of Complexes I through IV, reflecting the cell’s attempt to preserve mitochondrial integrity under stress.Fig. 13SDH-b knockout simulations (dimensionless time \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau $$\end{document} ). a Intracellular lactate concentration increases in SDH-b knockout cells, indicating a shift towards anaerobic glycolysis despite unchanged oxygen flux. b Intracellular pyruvate concentration decreases, suggesting accelerated conversion to lactate, supporting enhanced glycolytic flux. cIntracellular glucose concentration decreases, indicating increased glucose uptake and utilisation through glycolysis, aligning with the Warburg effect observed in cancer cells. d Intracellular ATP concentration maintenance, reflecting impaired ATP synthesis due to dysfunctional ETC activity and consistent with the observed ATP synthase reversal, exacerbating cellular energy stress
Upon investigating cytosolic metabolite results, depicted in Fig. 13, we were able to further elucidate the metabolic consequences of SDH-b knockout in chromaffin cells. These observations highlight the pseudohypoxic phenotype characteristic of phaeochromocytomas, as reported by Kl’učková et al. (2020).
Figure 13a shows a marked increase in intracellular lactate concentration in SDH-b knockout cells compared to wild-type (WT) cells. This substantial rise is indicative of a shift towards anaerobic glycolysis, as seen in our hypoxia benchmarking results. Although this is a common response to hypoxia, where cells rely more heavily on glycolysis for ATP production due to compromised oxidative phosphorylation, here the net respiration of the system is not reduced. In fact, the model’s elevated lactate levels corroborate the increased activity of lactate dehydrogenase observed in previous simulations (Fig. 12b), reinforcing the model’s accuracy in mimicking hypoxic-like metabolic responses.
In Fig. 13b, it is possible to observe a significant reduction in intracellular pyruvate concentration in SDH-b knockout cells. This decrease suggests an accelerated conversion of pyruvate to lactate, further supporting the enhanced glycolytic flux. The reduced pyruvate levels align with the expected metabolic adjustments where pyruvate is rapidly utilised to maintain glycolytic ATP production under conditions where the ETC function is impaired (Seifert et al. 2015).
Figure 13c illustrates a decrease in intracellular glucose concentration in SDH-b knockout cells. Looking at the flux rate, this result indicates an uptick in glucose uptake and utilisation through glycolysis, a compensatory mechanism to counteract the diminished ATP production from oxidative phosphorylation. This model also shows a decreased reliance on extracellular glucose, which aligns with our experimental data (Fig. 5b) (Kl’učková et al. 2020; Li and Ye 2014; Wheaton and Chandel 2011).
Lastly, Fig. 13d shows a maintenance of intracellular ATP concentration in SDH-b knockout cells-relative to WT, reflecting the hydrolysis of mitochondrial ATP due to dysfunctional ETC activity-this occurs as a compensatory mechanism via the glycolytic pathway. The ATP levels are consistent with the observed reversal of ATP synthase activity by Kl’učková et al. (2020) (Fig. 12d), where they observed ATP hydrolysis occurs to maintain mitochondrial membrane potential via blockage of ATPsynthase using oligomycin.
The observed shift to a more glycolytic phenotype despite normoxic conditions has important consequences in mitochondrial metabolism. In Fig. 14a, the model indicates that the concentration of mitochondrial pyruvate upon SDH-b knockout is significantly reduced compared to WT cells. This reduction aligns with the earlier observed decrease in cytosolic pyruvate (Fig. 13b) and suggests a disruption in pyruvate import into the mitochondria, likely due to impaired mitochondrial function and a shift towards cytosolic metabolism (Bender and Martinou 1863; Ruiz-Iglesias and Mañes 2021).
Figure 14b shows the concentration of mitochondrial succinate. In SDH-b knockout cells, succinate levels are markedly elevated compared to WT cells. This accumulation of succinate is a hallmark of SDH deficiency and is consistent with previous findings that link elevated succinate levels to the stabilisation of hypoxia-inducible factors (HIFs), promoting a pseudohypoxic phenotype even under normoxic conditions (Kl’učková et al. 2020; Kluckova and Tennant 2018; Brière et al. 2005; Letouzé et al. 2013; Goncalves et al. 2021). Note that the model indicates that the elevated succinate levels impede the proper functioning of the electron transport chain, further corroborating the observed reductions in Complex III and IV activities (Fig. 12b and c).
The model also shows a decrease in fumarate concentration, as seen in Fig. 14c. The decrease in fumarate is an expected consequence due to the bottleneck created by the loss of the SDH-b subunit, which compromises succinate dehydrogenase (SDH) activity. However, there is a slight accumulation of downstream metabolites like malate and fumarate. The buildup of these intermediates results from pyruvate carboxylases redirecting the flux of carbons towards the malate/aspartate shuttle. The model indicates this supports the increased oxidation of NADH into \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {NAD}^+$$\end{document} . This disrupts the TCA cycle and contributes to metabolic shifts associated with the pseudohypoxic response observed in chromaffin cells upon SDH-b loss (Kl’učková et al. 2020).
The mitochondrial membrane potential ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta \Psi $$\end{document} ) significantly decreases in SDH-b knockout cells, as shown in Fig. 14d. This reduction in membrane potential indicates impaired electron transport chain (ETC) function. The diminished ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta \Psi $$\end{document} ) reflects the compromised ability of the ETC to maintain proton gradients across the mitochondrial membrane, crucial for ATP synthesis. This loss of potential, the model indicates, is a direct consequence of the disrupted function of Complex II. The impaired ETC function exacerbates cellular energy stress, forcing the cell to rely more heavily on anaerobic glycolysis to meet its energy demands. The decrease in ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta \Psi $$\end{document} ) further supports the observation of reversed ATP synthase activity (Fig. 12d), where ATP is hydrolysed to pump protons back into the mitochondrial matrix, a desperate measure to preserve mitochondrial integrity and function under stress.Fig. 14SDH-b knockout simulations (dimensionless time \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau $$\end{document} ). a Mitochondrial pyruvate concentration significantly decreases in SDH-b knockout cells, indicating disrupted pyruvate import and a shift towards cytosolic metabolism. b mitochondrial succinate concentration increases markedly in SDH-b knockout cells, a hallmark of SDH deficiency, promoting a pseudohypoxic phenotype by stabilising hypoxia-inducible factors. c Mitochondrial fumarate concentration decreases due to the loss of SDH-b, creating a bottleneck in the TCA cycle, with slight upstream accumulation of intermediates like malate and fumarate. dMitochondrial membrane potential ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta \Psi $$\end{document} ) decreases, reflecting impaired electron transport chain function and energy stress
Mitochondrial Swelling
In addition to the metabolic changes highlighted by our simulations, the model also highlights significant alterations in mitochondrial volume and the proton motive force (PMF) in SDH-b knockout cells.Fig. 15SDH-b knockout simulations (dimensionless time \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau $$\end{document} ). a Mitochondrial NADH concentration decreases in SDH-b knockout cells, indicating impaired oxidative phosphorylation and a bottleneck in the TCA cycle. bMitochondrial volume increases in SDH-b knockout cells, reflecting altered osmotic balance and ionic homeostasis due to disrupted ETC function.c Mitochondrial pH slightly decreases, consistent with increased proton concentration in the intermembrane space caused by reversed ATP synthase activity and impaired ETC function.d Proton motive force (PMF) decreases in SDH-b knockout cells. Despite compensatory reversal of ATP synthase activity, the PMF is maintained but not at optimal levels, indicating energetically costly maintenance of the electrochemical gradient
The results above do not reflect the full consequence of a reversal of the ATP synthase. The model indicates that this activity has significant implications, evident upon examining the model’s proton motive force (PMF) across the mitochondrial membrane. Briefly, the PMF is the electrochemical gradient generated by the ETC, consisting of two components: a pH gradient (difference in proton concentration) and an electrical potential gradient (voltage difference) across the inner mitochondrial membrane (Hollinshead and Tennant 2016; Valdebenito et al. 2023). In normal conditions, ATP synthase utilises the PMF to drive the synthesis of ATP from ADP and inorganic phosphate. Protons flow back into the mitochondrial matrix through ATP synthase, which drives the phosphorylation of ADP. When ATP synthase operates in reverse (due to severe metabolic stress or impaired ETC function), it hydrolyses ATP to pump protons from the mitochondrial matrix to the intermembrane space. This reversal is a compensatory mechanism to maintain the PMF, particularly the electrical potential gradient (Acin-Perez et al. 2023).
The reversal of ATP synthase activity, in this model, helps to maintain the electrical component of the PMF. By hydrolysing ATP and pumping protons out of the matrix, ATP synthase helps sustain the voltage difference across the inner mitochondrial membrane. However, this comes at the cost of depleting the available mitochondrial ATP pool, as ATP is consumed to sustain the proton gradient instead of being produced. The pH gradient (difference in proton concentration) is also affected. Proton pumping by reversed ATP synthase increases the proton concentration in the intermembrane space, which can help maintain the pH gradient (Acin-Perez et al. 2023). Despite this, the overall efficiency of the proton gradient might be reduced due to the impaired ETC function, which normally contributes significantly to proton pumping. The combined effect of ATP synthase reversal and reduced ETC efficiency results in a complex balance. The PMF is partially maintained by the reversed ATP synthase activity, but the system is under stress, leading to an overall less stable and energetically costly maintenance of the PMF.
Figure 15a shows the model’s reaction to SDH-b loss by demarcating a decrease in mitochondrial NADH concentration compared to WT imCCs. This decrease indicates impaired oxidative phosphorylation, as NADH is a crucial substrate for the ETC function, and its reduced levels reflect, in this model, a diminished electron transfer through the chain-as seen in our results above (Figs. 13 and 14). However, the result is also indicative of a bottleneck in the TCA cycle, as the uncoupling of the malate/aspartate shuttle prevents the model from oxidising \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {NAD}^+$$\end{document} for its use in the electron transport chain. These results reflect the metabolic interplay seen in our own experimental data (Fig. 7)
In Fig. 15b, it is possible to see how the loss of SDH-b leads to an increase in mitochondrial volume. This increase, in our model, is a direct response to the altered osmotic balance and ionic homeostasis resulting from disrupted ETC function. Mitochondrial swelling is a well-documented response to metabolic stress caused by SDH-b loss (Kl’učková et al. 2020; Goncalves et al. 2021). Unlike observations by Kl’učková et al. (2020), our model did not achieve the expected fold-change. This could be due to diverse factors, we suspect our model does not account for other metabolic alterations that may affect mitochondrial swelling (i.e., the mitochondrial transition pore mPTP, which is regulated by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {Ca}^{2+}$$\end{document} ).
Concomitantly, the model shows a slight decrease in mitochondrial pH in SDH-b knockout cells (Fig. 15c). The decrease in pH is consistent with the increased proton concentration in the mitochondrial intermembrane space due to the reversed ATP synthase activity and impaired ETC function. However, this acidification is not too drastic, as might be naively expected. This is due, the model indicates, to its capacity of sustaining the proton motive force under metabolic stress (i.e., ATP synthase reversal).
Indeed, Fig. 15d depicts a decrease in the proton motive force (PMF) in SDH-b knockout cells. Despite the compensatory reversal of ATP synthase activity, the PMF is maintained, albeit not at optimal levels. The reduced PMF indicates an overall less efficient and energetically costly maintenance of the electrochemical gradient, crucial for ATP production and mitochondrial function.
Complex 1 Retention
As the SDH-b’s ability to transport electrons through its sulfur-iron clusters (denoted in the model as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$k_{[^{3}Fe-^{4}S]}$$\end{document} ) decreases, initially, there is an observed increase in Complex I (C1) activity. This increase in C1 activity occurs as a compensatory response by the cell to maintain electron flow into the ETC, even as the function of Complex II is impaired due to SDH-b loss-this observation is in accordance to what Kl’učková et al. (2020) report. Furthermore, this adaptive mechanism, the model reveals, helps in sustaining the proton motive force (PMF), ensuring that ATP synthesis can continue, albeit at a potentially reduced efficiency (Acin-Perez et al. 2023).Fig. 16SDH-b knockout simulation. a Surface plot illustrating the relationship between Proton Motive Force (PMF) and ETC Complex I activity as SDH-b electron transport capacity ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$k_{3Fe-4S}$$\end{document} ) decreases. The x-axis represents ETC Complex I activity, the y-axis represents PMF, and the z-axis (or third independent variable) is ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$k_{3Fe-4S}$$\end{document} ). The color of the surface is not indicative of any specific variable but rather represents the model-generated surface formed from individual simulation solutions
However, as the capacity of SDH-b continues to decline beyond a critical point, the compensatory mechanism of upregulating Complex I becomes insufficient. This is depicted in Fig. 16 where the PMF begins to plateau and eventually declines as the ability of the enzyme to transport electrons (all the way to SDH-b K.O.) diminishes. As a result, the cell can no longer sustain the necessary proton gradient solely through Complex I activity, leading to a total collapse of the ETC functionality.
This result highlights the initial adaptive increase in Complex I activity in response to SDH-b loss, followed by a collapse in the system as the compensatory mechanisms are overwhelmed. The maintenance of PMF through proton leaks and the reversal of ATP synthase activity underscores the cell’s efforts to adapt to metabolic stress and sustain mitochondrial function under compromised conditions.
Our previous study (Kl’učková et al. (2020)) reported that despite SDH-b deficiency, chromaffin cells retain Complex I function, which is critical for maintaining metabolic fitness. Their study shows that the retention of Complex I activity helps to sustain mitochondrial NADH oxidoreductase function, thereby preventing a complete metabolic collapse. This supports our model’s prediction of an initial compensatory upregulation of Complex I activity, followed by a critical threshold beyond which the system fails to maintain ETC functionality.
In-depth discussion of the surface plot reveals key insights into the cellular adaptations and limitations in response to SDH-b deficiency. Initially, as SDH-b function decreases, the electron transport chain compensates by upregulating Complex I activity to maintain the PMF. This adaptive response is reflected in the transient increase in Complex I activity, which helps to sustain ATP production and cellular energy balance. However, this compensatory mechanism has its limits. As SDH-b capacity further declines, the system reaches a tipping point where Complex I alone cannot sustain the electron flow and PMF, leading to a gradual decline in ETC efficiency.
Figure 16 highlights the dynamic interplay between SDH-b activity and Complex I function. The observed fold in the solution surface indicates a critical threshold beyond which the compensatory mechanisms can no longer cope with the metabolic stress. At this point, the system experiences a rapid decline in ETC function, leading to impaired mitochondrial respiration and energy production-all of this whilst the proton leak salvages the PMF. This fold represents a ‘sink’ in the mathematical sense, where the system’s state is drawn towards a stable, but suboptimal, condition of impaired Complex I activity and reduced PMF. This proton leak seemingly, is able to dissipate the proton gradient, reducing the efficiency of ATP synthesis but preventing a complete collapse of the mitochondrial membrane potential. This leak mechanism, along with the reversal of ATP synthase activity, highlights the cell’s efforts to adapt to the severe metabolic stress induced by SDH-b loss. The reliance on proton leaks and ATP synthase reversal indicates that this compensation comes at a significant energetic cost, the model predicts this could limit the cells’ ability to cope with additional stressors-despite their perceived metabolic robustness.
Discussion
The primary objective of this study was to develop a comprehensive mathematical model to elucidate the metabolic consequences of SDH-b loss in chromaffin cells, particularly focusing on the electron transport chain (ETC) and associated metabolic pathways. Our findings, derived from detailed simulations and experimental data integration, provide significant insights into the adaptive mechanisms employed by chromaffin cells in response to SDH-b deficiency (Kl’učková et al. 2020).
To validate and extend previous observations, we revisited \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}\hbox {C}_6$$\end{document} -glucose labelling experiments, identifying any discrepancies or inconsistencies between our current and past findings. This approach facilitated a deeper understanding of the underlying biological processes driving our mathematical model (Figs. 6 and 7). By subjecting the cells to fully labelled glucose ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{13}\hbox {C}_6$$\end{document} -glucose), we aimed to gather comprehensive data on metabolite concentrations, flux rates, and enzyme kinetics (Figs. 8, 9 and 10). These data serve as crucial inputs for refining model parameters, ensuring their congruence with experimental observations.
Our model simulations support earlier observations that SDH-b knockout leads to notable alterations in the activity of ETC complexes (Kl’učková et al. 2020; Goncalves et al. 2021; Letouzé et al. 2013). While Complex I activity is maintained through compensatory upregulation, Complexes III and IV exhibit substantial decreases in activity (Fig. 13). This imbalance underscores the critical interdependencies within the ETC, where impairment in one complex can propagate through the chain, affecting overall electron transport efficiency. The maintenance of Complex I activity suggests that cells may activate alternative pathways or regulatory mechanisms to uphold mitochondrial function under stress conditions, aligning with previous studies on cellular adaptation strategies (Kluckova and Tennant 2018; Kl’učková et al. 2020).
A significant finding in our simulations was the reversal of ATP synthase activity, where the enzyme hydrolyses ATP to pump protons back into the mitochondrial matrix (Fig. 12d). This mechanism, the model indicates, serves to maintain the mitochondrial membrane potential in the face of compromised ETC function but at the cost of cellular ATP reserves (Fig. 14d). This reversal is indicative of severe metabolic stress and highlights the cell’s prioritization of maintaining membrane potential over efficient ATP production (Kl’učková et al. 2020).
SDH-b knockout simulations also elucidated a shift towards a pseudohypoxic phenotype, characterised by increased lactate production and reduced pyruvate and glucose levels (Fig. 13). These changes reflect an enhanced reliance on glycolysis for ATP production, mimicking the Warburg effect commonly observed in cancer cells (Goncalves et al. 2021; Wheaton and Chandel 2011). The elevated intracellular lactate and reduced pyruvate concentrations align with a metabolic shift towards anaerobic glycolysis, even under normoxic conditions. This shift is a well-documented response to hypoxia, but in our model, it occurs as a direct consequence of SDH-b loss, independent of oxygen availability (Kluckova and Tennant 2018; Goncalves et al. 2021; Letouzé et al. 2013).
Our results also observed a reproduction of experimental data where upon SDH-b knockout, an increase in mitochondrial volume led to irreversible changes in the cell’s proton motive force (PMF) (Fig. 15). The model indicated that mitochondrial swelling occurs as a response to altered osmotic balance and ionic homeostasis resulting from disrupted ETC function. The maintenance of PMF, despite the reversal of ATP synthase activity, underscores the cell’s effort to preserve mitochondrial integrity under stress. However, this comes at the cost of increased energetic demands, highlighting the metabolic strain imposed by SDH-b knockout (Kl’učková et al. 2020; Goncalves et al. 2021).
Despite the strengths of this model in capturing key metabolic alterations following SDH-b loss, several limitations should be noted. Firstly, the model simplifies the glycolytic pathway by aggregating it into a net reaction to optimize computational efficiency. While this approach provides valuable insights into overall metabolic shifts, it does not account for enzyme-specific regulation within glycolysis. Additionally, the model assumes a constant extracellular environment, meaning it does not simulate the depletion of glucose or the accumulation of lactate over time, which could influence long-term metabolic adaptations.
Furthermore, the model accounts for mitochondrial swelling but does not differentiate between swelling in the intermembrane space (IMS) and matrix, instead treating it as a proportional expansion relative to the cytosol. This simplification means that specific mechanisms driving IMS or matrix expansion cannot be independently assessed. Another limitation is the exclusion of calcium ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {Ca}^{2+}$$\end{document} ) regulation, which plays a key role in mitochondrial metabolism and signaling. The mitochondrial permeability transition pore (mPTP), a known contributor to mitochondrial swelling and apoptosis, is also not explicitly modeled, despite its relevance in SDH-b deficient cells.
In terms of electron transport chain dynamics, while the model captures the compensatory retention of Complex I activity, it does not include detailed feedback mechanisms that regulate long-term mitochondrial adaptations. Additionally, proton leak mechanisms, which balance ATP synthase activity and electron transport, are represented in a simplified manner. Lastly, although the model successfully predicts ATP synthase reversal, it does not incorporate potential alternative ATP salvage pathways that cells may use to counteract energy depletion.
Future iterations of the model should incorporate calcium dynamics, mitochondrial transition pore behavior, and explicit differentiation of IMS and matrix volume changes to improve predictive accuracy. Including time-dependent changes in extracellular metabolite concentrations would also enhance the model’s ability to capture longer-term metabolic shifts in SDH-b deficient cells. Despite these limitations, the current framework provides a robust first step toward a deeper understanding of mitochondrial adaptations to SDH-b dysfunction.
The mathematical model presented in this study provides a first-step towards the construction of a robust framework for understanding the metabolic consequences of SDH-b loss in chromaffin cells. The findings highlight the cell’s adaptive mechanisms, including the compensatory upregulation of Complex I, the reversal of ATP synthase activity, and the metabolic shift towards a pseudohypoxic phenotype, a poorly understood phenomena seen in phaeochromocytomas. These insights not only advance our understanding of the metabolic adaptations in SDH-b deficient cells but also offer potential avenues for new insights in other conditions where SDH plays a centre role. In subsequent studies we aim to refine the model by incorporating additional regulatory mechanisms, like incorporating the intersection of cell-signalling and metabolism ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {Ca}^{2+}$$\end{document} ) and validating the findings with further experimental data.
Supplementary Information
Below is the link to the electronic supplementary material.Supplementary file 1 (pdf 293 KB)
