Construction of pathway models of steroidogenesis using the mEPN notation
A model of human steroidogenesis is presented in Fig. 1, representing a framework of the reactions that produce biologically active steroids under normal conditions. Diagrams were compiled using the modified Edinburgh pathway notation (mEPN) using the network editing software yED (http://www.yworks.com). The editable version of this diagram is available as a ‘.graphml’ file that can be opened in yED (Additional file 1). Steroid and enzyme interaction information was obtained from Miller et al. [1] and associated references. Active pathways in different tissues (namely the adrenal reticularis, glomerulosa and fasciculata, the testis, ovary, prostate and placenta, Additional file 2) are highlighted with red boxes to easily visualise the steroidogenic reactions that take place in each tissue.
mEPN is based on the principles of ‘process diagrams’ and is designed to be unambiguous yet concise [3]. A detailed protocol of how to edit mEPN diagrams has been recently published [4]. Both the biological entities (such as proteins or steroids) and the way they interact with each other (such as phosphorylation or dimerisation) are represented as components in the pathway. Small biochemicals such as steroids are represented by hexagons, all proteins (enzymes) by rounded rectangles and all genes by parallelograms. Black rectangles allow for parameterisation of models whereby initial token input on nodes feeding into the pathway can be defined (Fig. 1b). It can be expanded to produce large, clear and informative pathway models [5, 6].
We have chosen to label steroids and their intermediates based on names in common use in the biological community but have also linked nodes to the Chemspider database (http://www.chemspider.com). This provides a reference to the exact biochemical structure a molecule represents and gives alternate names. The link can be opened by pressing the F8 key when the node is highlighted in yED (Fig. 1c). There are number of naming conventions for biochemical molecules. The comprehensive diagram of the major steroidogenic pathways (Fig. 1a) contains nodes to represent both the gene and the protein produced for each enzyme isoform, and therefore each gene node is linked to Ensembl (http://www.ensembl.org) [7].
Construction and parameterisation for dynamic flow of a cell-specific model of rat Leydig cell steroidogenesis using previously-published experimental data
Whilst the formalised framework diagram has utility as a resource in its own right, the ability to parameterise the pathway and run and test simulations immeasurably adds to its overall value. We constructed and parametrised a specific model of rat Leydig cell steroidogenesis during postnatal development and adulthood (Fig. 2) focussing on the specific reactions that Leydig cells use to produce their main steroid product: androgens. The simplified Leydig cell version (Fig. 2a) uses a single protein node to represent all isoforms of a particular enzyme and so an Ensembl link is not included. The editable version of the simplified diagram is also available as a ‘.graphml’ file that can be opened in yED, (Additional file 3).
Parameterisation of pathway diagrams constructed in yED to run as a signalling Petri nets (SPNs) in the software Graphia Professional is a logical process and no formal training in mathematical modelling is necessary for the user. A detailed description of how to parameterise yED diagrams so that they can be run as SPNs in Graphia Professional (Kajeka, Edinburgh, UK, formerly BioLayout Express3D) [2] can be found in Livigni et al. [4]. By representing a complex network as a Petri net, the SPN method models signal flow as the pattern of token accumulation at protein nodes over time.
Parameterisation of the cell-specific diagram was achieved using previously published enzyme activity data measured at three different stages of rat Leydig cell maturation [8]. In yED, tokens representing the enzyme activity in pmol/minute/million cells were added to the arrows connecting the black input nodes of each of the steroidogenic enzyme nodes using the notation a–b,c;d–e,f where ‘a–b’ are the first and last time blocks that you would like the number of tokens ‘c’ to be added to the model and ‘d–e’ are the first and last time blocks that you would like the number of tokens ‘f’ to be added to the model (and so on for the number of variable inputs required), as shown in Fig. 2a. The first 20 time blocks represent the ‘progenitor’ Leydig cell stage present at around postnatal day (pnd) 21 in the rat. Time blocks 21–50 represent the ‘immature’ Leydig cell stage present at around pnd 35 and time blocks 51–100 represent the mature adult Leydig cell that constitute all of the Leydig cells in the testis from pnd 90 onwards. The variation in enzyme activity at the three stages is illustrated by the black bar graphs next to the enzyme input nodes in Fig. 2a. The editable version of the parameterised diagram is available as a ‘.graphml’ file that can be opened in yED (Additional file 4).
The parameterised diagram was run as a SPN in Graphia Professional over 100 time blocks, 500 runs and with standard normal stochastic distribution and consumptive transitions [2]. The token flow was visualised as an animation seen in Additional file 5. Figure 2b shows a screenshot of the animated output of the token flow at each of the three stages (screenshots taken at time block 20 representing progenitor, 50 representing immature and 100 representing adult Leydig cell stages). This provides an overview of all of the nodes in the diagram and helps visualise the pathways of token flow. Two nodes were selected for specific visualisation in Fig. 2c: testosterone and 3α-androstanediol (‘3α-diol’). If these outputs are taken as a prediction of the relative production rate of these two steroids at immature and adult Leydig cell stages, we would predict that 3α-diol is more abundant than testosterone in immature Leydig cells and that testosterone is more abundant than 3α-diol in adult Leydig cells. This prediction of the model is consistent with previous experimental measurement [8], showing that our model appropriately recapitulates the in vivo situation and thus has utility for hypothesis testing.
Using the model to predict steroid production in Hsd17b3 deficiency
Enzymes of the 17-beta hydroxysteroid dehydrogenase class catalyse the conversion between 17-keto and 17-hydroxy-steroids. Different isoforms of the enzyme are expressed in different steroidogenic tissues. 17-beta hydroxysteroid dehydrogenase type 3 is the isoform expressed by Leydig cells in humans (HSD17B3), mice and rats (Hsd17b3) [9] and preferentially catalyses the conversion of androstenedione to testosterone and androstanedione to DHT. Male humans with mutations in HSD17B3 present with varying degrees of physiological undervirilisation and plasma androstenedione levels at the time of puberty are usually ten times normal levels [10]. To demonstrate the predictive power of our in silico model, we mimicked a loss of function mutation in Hsd17b3 by removing the token input from the HSD17B node of the rat Leydig cell model (Fig. 3a, b), and re-ran the simulation. This time tokens accumulated at the androstenedione node, with no testosterone produced, consistent with the circulating androgen profile observed in patients with a loss of function of HSD17B3 (Fig. 3c). When visualised as a graph (Fig. 3d), androstenedione production is shown to increase as Leydig cells mature during postnatal life, whereas no testosterone is produced. This simple demonstration shows the utility of the model for predicting outcomes of genetic or pharmacological manipulations before beginning any laboratory or in vivo work, and has the potential to be scaled with multiple knockouts modelled simultaneously.
