## Abstract

We introduce, and provide examples of, the application of the bond graph formalism to explicitly represent biophysical processes between and within modular biological compartments in ApiNATOMY. In particular, we focus on modelling scenarios from acid–base physiology to link distinct process modalities as bond graphs over an ApiNATOMY circuit of multiscale compartments. The embedding of bond graphs onto ApiNATOMY compartments provides a semantically and mathematically explicit basis for the coherent representation, integration and visualisation of multiscale physiology processes together with the compartmental topology of those biological structures that convey these processes.

## 1. Introduction

Two significant challenges for the VPH and Physiome projects are (i) to build better tools for combining biophysically based models involving multiple types of biological process (e.g. cell signalling, ion channel electrophysiology, metabolism and soft tissue mechanics) and (ii) to find ways of linking the biophysical models with the anatomical context in which they are applied. In this paper, we propose a new framework that addresses both of these issues.

### 1.1. Background

The VPH community^{1} and the IUPS Physiome^{2} community [1] are together developing the standards, tools and databases required to support quantitative multiscale modelling of the human body for a wide range of pharmaceutical and clinical applications. The primary goal is to link molecular level information with cell, tissue, organ and whole body models that are based on anatomical data and biophysical principles. Such an approach is needed (e.g. for medical diagnostics) if we are to exploit combined measurements in one individual of molecular biomarkers, physiological measurements and clinical imaging.

At the centre of this, effort is a set of markup language standards for encoding models and data. CellML [2] is used for encoding models based on ordinary differential and algebraic equations; FieldML [3] is used for models with spatially varying fields (typically governed by partial differential equations); SED-ML [4] codes the simulation protocol (including, for example, the numerical algorithms for solving the models and the display formats for specified outputs from solution of the models) and BiosignalML [5] is a data standard for encoding the input and output time varying signals.

A range of open source desktop and web-browser tools have been developed for building and solving Physiome models, including, for example, OpenCOR [6]. This tool provides a desktop environment for creating or reading CellML and SED-ML files, displaying and changing model and solution parameters, solving the models, graphing the output variables and uploading the models back into the user's workspace on PMR.

A database of models, the Physiome Model Repository [7], has also been developed. It currently contains about 600 *exposures* (models, based on peer-reviewed publications, that have been exposed for public access) and another 400 collaborative *workspaces*, used for models currently under development and therefore with restricted access.

### 1.2. Challenge

The current challenge for the VPH/Physiome community is to leverage representational standards in engineering and biology to improve the interoperability (and therefore the re-usability) and to clarify the meaning and readability (and therefore the verifiability) of biophysically based models. This representational challenge has to be met at two levels, namely (i) at a purely process level, to encapsulate biophysical processes and (ii) at a location level in which the anatomical, cellular, subcellular and molecular structures that participate in, or convey, a process are described. Furthermore, overcoming this challenge entails the coherent combination of process and location representation in describing the physiological meaning, or semantics, of a model.

In this paper, we introduce a solution that combines the standardised representation of process and location. In particular, we propose a new framework for developing multiscale biophysical models in physiology by combining (i) *ApiNATOMY*, a representation of functional anatomy, with (ii) *bond graphs*, an engineering methodology that represents mass and energy-conserving processes consistently and powerfully. We also demonstrate the application of bond graphs in ApiNATOMY by generating and linking three biophysical scenarios from different scales, namely:

A. the fluid mechanics of blood and urine flow in the kidney;

B. the biochemical and diffusive processes between blood and urine of the proximal tubule of the nephron and

C. the multistate dynamics of a membrane transporter in the renal tubular epithelium.

## 2. Methods and results

### 2.1. Process knowledge: bond graphs

The theory of bond graphs was developed by Paynter [8] for electro-mechanical systems, with the theory of network thermodynamics, developed by Katchalsky and co-workers [9], providing additional thermodynamic concepts. The application of bond graphs to biological systems has been the focus of a number of studies (e.g. [10–17]).

Bond graphs provide a systematic way of constructing models of processes that obey mass conservation and energy conservation. (*N.B.* For readers not familiar with bond graphs, we provide a tutorial at this link: https://doi.org/10.17608/k6.auckland.5425879.v1.) The starting point is that power (or energy) is the common currency across all physical systems, and that power is always the product of a driving *potential* (** μ**), expressed as energy (in joules) per unit of some quantity and a

*flow*(

**), expressed as that quantity per second. For example, in electrical systems, the quantity being conveyed is charge in coulombs (C) and electrical power is J C**

*q*^{−1}× C s

^{−1}. (expressing J C

^{−1}as ‘volts’ and C s

^{−1}as ‘amps’ unnecessarily hides the underlying concepts of energy and power.) Similarly, J m

^{−1}(rather than newtons) and m s

^{−1}, or J m

^{−3}(rather than pascals) and m s

^{−3}, express the equivalent concepts for mechanics. For biochemical reactions, the flow of a substance expressed in mol s

^{−1}is driven by a chemical potential expressed as J mol

^{−1}. (

*N.B*. The seven units of the SI system under the newly proposed definitions (http://www.bipm.org/en/measurement-units/rev-si/) are now based on constants that are consistent with the use of joules and seconds (together covering energy and power), metres, moles, entropy, coulombs and candela.) These seven units cover all of continuum physics (table 1).

Another key feature of bond graphs is the clear separation of processes that *transmit* energy (conveyed by a *bond* with its associated pair of variables *potential* and *flow*) from processes that *store* energy (either statically like a capacitor or spring or dynamically like an inductor or inertial mass) or *dissipate* energy (like a resistor or viscous damper). A bond graph depicts the flow of energy around a system with the topology of that system giving rise to the conservation of mass and conservation of energy, but with the constitutive laws (table 2) defined empirically for the storage and dissipative mechanisms.

The key elements of the bond graph approach are shown in figure 1. The black arrows are bonds that carry a potential ** μ_{i}** and a flow

**. The central row in figure 1**

*υ*_{i}*a*illustrates the topology of a simple system (e.g. an electrical network of capacitors, resistors and inductors; a mechanical system of springs, dashpots and masses or a biochemical system of solute concentrations and membrane-bound reactions). A junction or node, where bonds meet, is labelled with either a potential

**or a flow**

*μ*_{i}**, and every bond impinging on that node has the value of the potential or flow at that node. The flows in the bonds at a potential node (also called a 0-node) sum to zero (mass conservation) and the potentials in the bonds at a flow node (also called a 1-node) sum to zero (energy conservation). In both cases, transmission of power through nodes is conserved () because with**

*υ*_{i}**constant for 0-nodes and with**

*μ***constant for 1-nodes. (Note that labelling a 0-node with**

*υ***and a 1-node with**

*μ*_{i}**, and the convention that the bonds impinging on the node adopt the value at that node, is a convenient shorthand—the full bond graph diagram is shown in figure 1**

*υ*_{i}*b*where every bond has been annotated with its potential-flow conjugate pair.)

The top band of figure 1*a* shows the constitutive relations for the three types of mechanism—static energy storage (a capacitor for which potential ** μ** is a function of

**, the integral of flow**

*q***), energy dissipation (a resistance for which flow**

*υ***is a function of potential**

*υ***) and dynamic energy storage (an inductor for which**

*μ***is a function of**

*μ***, the derivative of flow**

*a***). The equivalent quantities for a range of physical systems are shown in table 2. The conservation equations generated by the bond graph topology are shown in the lower band of figure 1**

*υ**a*. For the potential nodes, these equations are (i.e. conservation of current for an electrical system or mass for a mechanical or biochemical system), and for the flow node, the equation is (representing voltage drops around an electrical circuit, or force balance for a mechanical system, or the stoichiometry of a chemical reaction for a biochemical system). Also, shown below the conservation equations at the two potential nodes are the differential equations that express

**(the stored charge in an electrical capacitor, displacement of a mechanical spring or solute amount dissolved in a solvent) or**

*q***(the rate of change of current in an electrical system, or acceleration in a mechanical system) with the flow**

*a***.**

*υ*In summary, bond graphs provide a graphical notation for the set of linear constraint equations (the conservation laws), but note that constitutive laws can be nonlinear (and usually are in biological systems).

### 2.2. Location knowledge: ApiNATOMY

The management of structural knowledge of (i) the scales at which compartments are defined as well as (ii) the parthood and connectivity relationships between compartments are key to provide a mechanistic explanation to a correlation between biophysical measurements (or model variables) from two distinct compartments. The consistent management of such knowledge provides the means to determine the anatomical route by which physical influence (that mechanistically accounts for the correlation) is transferred from one body region to another (also discussed further in [18]).

Given the above knowledge management goals, a methodology such as ApiNATOMY^{3} [19] is required to provide formal knowledge articulation for multiscale compartments and, in particular, to articulate topological information about key anatomical structures and processes of physiological relevance (figure 2) such as

a.

*Conduits*. geometric layered compartments from any scale, representing conduits that convey the movement of solutes. Conduits fall exclusively into one of three topological categories of tube, bag or cyst. Tubes, bags and cysts are all conduits of a kind that are topologically distinguished as follows:i. a tube is a conduit with both terminal borders open;

ii. a bag is a conduit that has one open and one closed terminal border and

iii. a cyst is a conduit that has both terminal borders closed.

These geometric material objects are dependent upon a variety of features, in particular

— the dimensional scales of their thickness and length;

— a layered structure; and

— an intrinsic polarity to support the relative orientation of geometric material objects within other geometric material objects.

b.

*Trees*. the longitudinal linkage of contiguous conduits of different scales to create canonical arborisation models of mesothelial, endothelial, neural and epithelial trees. A tree is an important kind of material object in the ApiNATOMY model, because many circuits of conduits have tree-like structures or tree-like parts. Furthermore, the number of levels in trees can quickly become very high and a discrete reconstruction of tree structures from simpler segmental objects becomes impractical. Thus, the category of tree as a material object in its own right is necessary for any practical treatment of the topology of physiological circuitry. Trees are characterised by the following features:i. a number of levels,

ii. the degree of branching between consecutive levels,

iii. an association with a conduit type for each level,

iv. a distribution of the length for the conduit associated with each level and

v. a distribution of the thickness for the conduit associated with each level.

c.

*Coalescences*. the radial linkage of conduits that overlap over their outermost layer—a coalescence between two conduits occurs side to side, such that the outermost layer of the conduits is shared;d.

*Housed multiscale assemblies*. the embedding of assemblies of conduits, trees and coalescences within layers of a housing conduit of a scale that is larger than that of the assembled constituent components.e.

*Bond graphs*. the application of the bond graph formalism to formally represent the transfer of energy or mass along and across conduits in terms of physical laws.

The ApiNATOMY method thus outlined takes into account the above requirements to manage topological knowledge about key anatomical structures.

### 2.3. Multiscale scenarios: combining knowledge about process and location

In the examples below, we show how equations generated by the bond graph approach can be expressed in CellML and solved with OpenCOR. We also show how the quantities and processes defined by the bond graph model can be expressed within the anatomical framework established by ApiNATOMY.

#### 2.3.1. Long-range advection scenario: renal blood supply

Here, we focus on the flow of blood over renal arterial, capillary and venous systems. The key anatomical components of this system are depicted in ApiNATOMY in figure 4.

We give an example in figure 3 of the bond graph formulation of blood flow over three modules of branching blood vessels: arterial (figure 4, tube segments B–G), capillaries (figure 4, tube segments H–J) and venous (figure 4, tube segments K–P).

Arterial vasculature is represented by a single module that includes a compliance represented by the parameter *E*_{1}, two viscous resistance terms *R*_{1} and *R*_{2}, and two inertial terms *I*_{1} and *I*_{2}.

Elastance is a function of the arterial wall, while the viscous and inertial terms are a function of the material properties of blood (viscosity and density, respectively) and the dimensions (radius and length) of the blood vessels.

Viscous resistance *R* can be evaluated from Poiseuille's equation, which describes the relation between pressure drop and steady blood flow through a uniform vessel:

where *μ* is blood viscosity, *l* is the length of the vessel and *r* is the radius of the vessel.

Vessel elastance, *E*, can be calculated from blood vessel properties:
where *h* and *E _{m}* are the thickness and elasticity of the vessel wall, respectively.

The inertial term *I* can be obtained numerically using the blood vessel properties:
where *ρ* is blood density.

The capillary bed module contains *n*_{2} identical capillaries. Blood in these vessels flows into a single representative module for the venous side. Both the capillary module and the venous module contain compliant, inertial and resistive components similar to the arterial one, but with different parameter values.

Compared to figure 3, a more granular depiction of bond graph modules (corresponding to the anatomical segments in figure 4) is shown in figure 5. This depiction sets the stage to locate a bond graph module onto a two-layer ApiNATOMY tube segment representing a cardiovascular conduit. The details of this overlay are elaborated in figure 6.

To illustrate how the pressure in one vascular segment can influence the wall elastance and hence resistance to flow of another vessel, we include the bond graph representation of the vasa vasorum for the renal artery in figure 7. The elastance of the renal artery wall (represented by *E*_{1}) is a function of renal artery pressure *μ*_{1}. The elastance of the vasa vasorum, which is embedded in the renal artery wall, is represented by *E*_{v}, and this can be described empirically as a function of *μ*_{1}.

#### 2.3.2. Subcellular transepithelial transport scenario

reabsorption is one of the main cellular mechanisms in the kidney. filtered by the glomerulus is combined with H^{+} excreted from the epithelial cell (via an ATP-ase proton-pump and the passive sodium–hydrogen exchanger 3 (NHE3) transporter) to form CO_{2} and water (via the carbonic anhydrase reaction in the tubule). Then, CO_{2} returns into the epithelial cell (via the diffusion channel in the brush border) to regenerate in the epithelial cell (via the carbonic anhydrase reaction in the epithelial cell). Finally, diffuses into the blood through leak channels in the basolateral membrane of the epithelial cell.

Figure 8 gives the bond graph corresponding to the reabsorption model in the kidney. We assume that every that is transported out of the epithelial cell to the blood is immediately replaced with a H^{+} from the hydration reaction. We are also assuming the cell to be a well-mixed system. In this diagram, symbols representing ionic concentration measurements are located in the precise ApiNATOMY layer representing the corresponding anatomical compartment from where such measurements were taken.

The six reactions are modelled by the six Re components. The reactions for this model are (c, cytosol; b, basolateral extracellular; a, apical extracellular) as follows:

The state derivatives are as follows: and the reaction flows are as follows:

#### 2.3.3. Sodium–hydrogen exchanger 3 transporter mechanism scenario

The NHE3 (also known as the antiporter) was studied experimentally and explained by a biophysical model by Aronson [20]; further experiments and modelling were conducted by Weinstein [21]. Chang & Fujita [22] examined the kinetics of the model and calculated the kinetic rates. This section looks at the bond graph-based model of NHE3 based on the model of Chang & Fujita [22].

The model illustrated in figure 9 is based on the eight-state biomolecular cycle presented in the Chang & Fujita model. This figure shows an electroneutral antiporter in which a single transporter site, *E*, may be oriented towards the outside, *E*_{o}, or internal, *E*_{i}, membrane face. It also includes three solute species: Na^{+}, H^{+} and . This transporter mediates the 1 : 1 exchange of Na^{+} for H^{+} with competitive binding and transport of .

The nine reactions are modelled by the nine **Re** components, and there is no slippage reaction component which the enzyme undergoes conformational changes without transporting species:

The nine sets of reaction kinetic parameters are given in the Chang & Fujita model and listed in the first two columns of table 3, and the third column gives the corresponding equilibrium constants.

The vector of nine equilibrium constants is converted into the vector of 10 thermodynamic constants *K* of table 4 using the formula of [11].

The corresponding rated constants are then computed as discussed by Gawthrop *et al.* [11] and listed in the final column of table 3.

The eight states of the transporter form a conserved moiety: where is the total number of transporters corresponding to the experimental situation. Here, a normalised model was used where, in effect, a value of was used. The corresponding membrane current can then be scaled accordingly using the physiological value of .

The biomolecular system of the is described in stoichiometric form as follows: where is the system state, is the stoichiometric matrix and is the vector of reaction flows. In the case of the , state derivatives are as follows:

and the reaction flows are as follows

The bond graph model of is implemented in CellML, simulated using OpenCOR and illustrated an ApiNATOMY context using bond graphs in figure 10.

## 3. Discussion

Bond graphs provide a highly modular and systematic way of expressing the topology of a process network in graphical form. This topological structure is amenable to express both (i) the equations governing conservation of mass and conservation of energy for the network and (ii) the constitutive relations that capture the empirically determined properties of the components of the network. These two sets of equations map precisely to a declarative framework such as CellML and can be solved with freely available open source software such as OpenCOR. Furthermore, because biophysical models are generated from an underlying physics framework and are given a biological context through the embedding with ApiNATOMY, annotations needed to describe both the physical and biological meaning of the model parameters and variables are automatically available as CellML metadata.

From the perspective of knowledge management about biological structure and location, ApiNATOMY provides a flexible modularisation of compartments that correspond to the field of ‘0’ nodes in bond graphs. Specifically, based on our experience with the above three biophysical scenarios, an ApiNATOMY compartment is semantically equivalent to a well-mixed, equipotential field tracked by a variable associated with a 0-node in a bond graph. In the case of a chemical reaction, this 0-node corresponds to the chemical potential of a single chemical species in a well-mixed compartment.

The explicit mapping of bond graph nodes to ApiNATOMY compartmental location sets the stage to coherently link bond graphs for the (i) NHE3 mechanism, (ii) subcellular transepithelial transport and (iii) long-range renal advection. In particular, such linkage would require (i) the bond graph and compartments of the NHE3 mechanism in figure 10 to plug into the Re4 component of the model in figure 8, and the resulting solution to be plugged into (ii) the compartmental coalescence of the peritubular capillary (figure 4, J) and proximal convoluted tubule (figure 4, 15) to bridge to the higher-scale bond graph model of long-range advective flow of blood and urine. To this end, future work in support for this strategy in multiscale model integration includes the following:

a) developing new interfaces in OpenCOR to create models based on bond graphs and to associate the parameters, variables and compartments of the CellML model with the corresponding components in ApiNATOMY. These interfaces will be tailored to the context appropriate for the model or part of the model (e.g. electrophysiology, metabolism, signalling pathways and mechanics).

b) building a modular catalogue of multiscale compartments in ApiNATOMY to create a canonical map of anatomical components across scales (see mockup of this map in figure 11).

c) generating a catalogue of structure–function modules combining ApiNATOMY compartments and associated bond graphs for repeating canonical units in physiology such as vessel segments, primary functional tissue units [23] and membrane transporters.

## 4. Conclusion

Bond graphs support cross-talk between processes from different domains of physics and encourage unit harmonisation across physical domains. The power of the bond graph approach is particularly relevant when different physical mechanisms interact with one another. In such integrative scenarios, the 1-node concept, in particular, is a powerful mechanism for ensuring consistency between, for example, molar flow rates, associated charge transfer and metabolic energy considerations. Conversely, the 0-node concept is essential to match an equipotential field with anatomical location—a key feature that informs the interoperation between bond graphs and ApiNATOMY.

Implementing a combination between ApiNATOMY and bond graphs is anticipated to foster significant knowledge modularisation in physiology by

i. encouraging the communal development of catalogues of modular structure–function units at different scales;

ii. supporting the coherent visualisation and curation of multiscale physiology knowledge by clinicians, physiologists, systems biologists and modellers and

iii. enhancing publication services in physiology model building, sharing and verifiability, such as the upcoming

*Physiome*journal that is to be launched by the IUPS and VPH communities to support reproducible, reusable biophysically based modelling, and associated biomedical knowledge management frameworks.

## Data accessibility

This article has no additional data.

## Competing interests

We declare we have no competing interests.

## Funding

We received no funding for this study.

## Footnotes

One contribution of 9 to a theme issue ‘The virtual physiological human: translating the VPH to the clinic’.

- © 2018 The Author(s)

Published by the Royal Society. All rights reserved.