Osteocytes as a record of bone formation dynamics: A mathematical model of osteocyte generation in bone matrix

How osteons form: A quantitative hypothesis-testing analysis of cortical pore filling and wall asymmetry

Osteon morphology provides valuable information about the interplay between different processes involved in bone remodelling. The correct quantitative interpretation of these morphological features is challenging due to the complexity of interactions between osteoblast behaviour, and the evolving geometry of cortical pores during pore closing. We present a combined experimental and mathematical modelling study to provide insights into bone formation mechanisms during cortical bone remodelling based on histological cross-sections of quiescent human osteons and hypothesis-testing analyses. We introduce wall thickness asymmetry as a measure of the local asymmetry of bone formation within an osteon and examine the frequency distribution of wall thickness asymmetry in cortical osteons from human iliac crest bone samples from women 16-78 years old. Our measurements show that most osteons possess some degree of asymmetry, and that the average degree of osteon asymmetry in cortical bone evolves with age. We then propose a comprehensive mathematical model of cortical pore filling that includes osteoblast secretory activity, osteoblast elimination, osteoblast embedment as osteocytes, and osteoblast crowding and redistribution along the bone surface. The mathematical model is first calibrated to symmetric osteon data, and then used to test three mechanisms of asymmetric wall formation against osteon data: (i) delays in the onset of infilling around the cement line; (ii) heterogeneous osteoblastogenesis around the bone perimeter; and (iii) heterogeneous osteoblast secretory rate around the bone perimeter. Our results suggest that wall thickness asymmetry due to off-centred Haversian pores within osteons, and that nonuniform lamellar thicknesses within osteons are important morphological features that can indicate the prevalence of specific asymmetry-generating mechanisms. This has significant implications for the study of disruptions of bone formation as it could indicate what biological bone formation processes may become disrupted with age or disease.

Hypoxia inducible factor-1α related mechanism and TCM intervention in process of early fracture healing

As a common clinical disease, fracture is often accompanied by pain, swelling, bleeding as well as other symptoms and has a high disability rate, even threatening life, seriously endangering patients' physical and psychological health and quality of life. Medical practitioners take many strategies for the treatment of fracture healing, including Traditional Chinese Medicine (TCM). In the early stage of fracture healing, the local fracture is often in a state of hypoxia, accompanied by the expression of hypoxia inducible factor-1α (HIF-1α), which is beneficial to wound healing. Through literature mining, we thought that hypoxia, HIF-1α and downstream factors affected the mechanism of fracture healing, as well as dominated this process. Therefore, we reviewed the local characteristics and related signaling pathways involved in the fracture healing process and summarized the intervention of TCM on these mechanisms, in order to inspirit the new strategy for fracture healing, as well as elaborate on the possible principles of TCM in treating fractures based on the HIF molecular mechanism.

Theoretical Evidence of Osteoblast Self-Inhibition after Activation of the Genetic Regulatory Network Controlling Mineralization

Bone is a hard-soft biomaterial built through a self-assembly process under genetic regulatory network (GRN) monitoring. This paper aims to capture the behavior of the bone GRN part that controls mineralization by using a mathematical model. Here, we provide an advanced review of empirical evidence about interactions between gene coding (i) transcription factors and (ii) bone proteins. These interactions are modeled with nonlinear differential equations using Michaelis-Menten and Hill functions. Compared to empirical evidence, the two best systems (among 12⁶=2,985,984 possibilities) use factors of inhibition from the start of the activation of each gene. It reveals negative indirect interactions coming from either negative feedback loops or the recently depicted micro-RNAs. The difference between the two systems also lies in the BSP equation and two ways for activating and reducing its production. Thus, it highlights the critical role of BSP in the bone GRN that acts on bone mineralization. Our study provides the first theoretical evidence of a necessary genetic inhibition for bone mineralization with this work.

Osteon: Structure, Turnover, and Regeneration

Bone is composed of dense and solid cortical bone and honeycomb-like trabecular bone. Although cortical bone provides the majority of mechanical strength for a bone, there are few studies focusing on cortical bone repair or regeneration. Osteons (the Haversian system) form structural and functional units of cortical bone. In recent years, emerging evidences have shown that the osteon structure (including osteocytes, lamellae, lacunocanalicular network, and Haversian canals) plays critical roles in bone mechanics and turnover. Therefore, reconstruction of the osteon structure is crucial for cortical bone regeneration. This article provides a systematic summary of recent advances in osteons, including the structure, function, turnover, and regenerative strategies. First, the hierarchical structure of osteons is illustrated and the critical functions of osteons in bone dynamics are introduced. Next, the modeling and remodeling processes of osteons at a cellular level and the turnover of osteons in response to mechanical loading and aging are emphasized. Furthermore, several bioengineering approaches that were recently developed to recapitulate the osteon structure are highlighted. Impact statement This review provides a comprehensive summary of recent advances in osteons, especially the roles in bone formation, remodeling, and regeneration. Besides introducing the hierarchical structure and critical functions of osteons, we elucidate the modeling and remodeling of osteons at a cellular level. Specifically, we highlight the bioengineering approaches that were recently developed to mimic the hierarchical structure of osteons. We expect that this review will provide informative insights and attract increasing attentions in orthopedic community, shedding light on cortical bone regeneration in the future.

Modelling cell guidance and curvature control in evolving biological tissues

Tissue geometry is an important influence on the evolution of many biological tissues. The local curvature of an evolving tissue induces tissue crowding or spreading, which leads to differential tissue growth rates, and to changes in cellular tension, which can influence cell behaviour. Here, we investigate how directed cell motion interacts with curvature control in evolving biological tissues. Directed cell motion is involved in the generation of angled tissue growth and anisotropic tissue material properties, such as tissue fibre orientation. We develop a new cell-based mathematical model of tissue growth that includes both curvature control and cell guidance mechanisms to investigate their interplay. The model is based on conservation principles applied to the density of tissue synthesising cells at or near the tissue's moving boundary. The resulting mathematical model is a partial differential equation for cell density on a moving boundary, which is solved numerically using a hybrid front-tracking method called the cell-based particle method. The inclusion of directed cell motion allows us to model new types of biological growth, where tangential cell motion is important for the evolution of the interface, or for the generation of anisotropic tissue properties. We illustrate such situations by applying the model to simulate both the resorption and infilling components of the bone remodelling process, and to simulate root hair growth. We also provide user-friendly MATLAB code to implement the algorithms.

The mechanoresponse of bone is closely related to the osteocyte lacunocanalicular network architecture

The explanation of how bone senses and adapts to mechanical stimulation still relies on hypotheses. The fluid flow hypothesis claims that a load-induced fluid flow through the lacunocanalicular network can be sensed by osteocytes, which reside within the network structure. We show that considering the network architecture results in a better prediction of bone remodeling than mechanical strain alone. This was done by calculating the fluid flow through the lacunocanalicular network in bone volumes covering the complete cross-sections of mouse tibiae, which underwent controlled in vivo loading. The established relationship between mechanosensitivity and network architecture in individual animals implies possibilities for patient-specific therapies. A new connectomics approach to analyze lacunocanalicular network properties is necessary to understand skeletal mechanobiology.

Organisms rely on mechanosensing mechanisms to adapt to changes in their mechanical environment. Fluid-filled network structures not only ensure efficient transport but can also be employed for mechanosensation. The lacunocanalicular network (LCN) is a fluid-filled network structure, which pervades our bones and accommodates a cell network of osteocytes. For the mechanism of mechanosensation, it was hypothesized that load-induced fluid flow results in forces that can be sensed by the cells. We use a controlled in vivo loading experiment on murine tibiae to test this hypothesis, whereby the mechanoresponse was quantified experimentally by in vivo micro-computed tomography (µCT) in terms of formed and resorbed bone volume. By imaging the LCN using confocal microscopy in bone volumes covering the entire cross-section of mouse tibiae and by calculating the fluid flow in the three-dimensional (3D) network, we could perform a direct comparison between predictions based on fluid flow velocity and the experimentally measured mechanoresponse. While local strain distributions estimated by finite-element analysis incorrectly predicts preferred bone formation on the periosteal surface, we demonstrate that additional consideration of the LCN architecture not only corrects this erroneous bias in the prediction but also explains observed differences in the mechanosensitivity between the three investigated mice. We also identified the presence of vascular channels as an important mechanism to locally reduce fluid flow. Flow velocities increased for a convergent network structure where all of the flow is channeled into fewer canaliculi. We conclude that, besides mechanical loading, LCN architecture should be considered as a key determinant of bone adaptation.

Heterogeneity of the osteocyte lacuno-canalicular network architecture and material characteristics across different tissue types in healing bone

Various tissue types, including fibrous connective tissue, bone marrow, cartilage, woven and lamellar bone, coexist in healing bone. Similar to most bone tissue type, healing bone contains a lacuno-canalicular network (LCN) housing osteocytes. These cells are known to orchestrate bone remodeling in healthy bone by sensing mechanical strains and translating them into biochemical signals. The structure of the LCN is hypothesized to influence mineralization processes. Hence, the aim of the present study was to visualize and match spatial variations in the LCN topology with mineral characteristics, within and at the interfaces of the different tissue types that comprise healing bone. We applied a correlative multi-method approach to visualize the LCN architecture and quantify mineral particle size and orientation within healing femoral bone in a mouse osteotomy model (26 weeks old C57BL/6 mice). This approach revealed structural differences across several length scales during endochondral ossification within the following regions: calcified cartilage, bony callus, cortical bone and a transition zone between the cortical and callus region analyzed 21 days after the osteotomy. In this transition zone, we observed a continuous convergence of mineral characteristics and osteocyte lacunae shape as well as discontinuities in the lacunae volume and LCN connectivity. The bony callus exhibits a 34% higher lacunae number density and 40% larger lacunar volume compared to cortical bone. The presented correlations between LCN architecture and mineral characteristics improves our understanding of how bone develops during healing and may indicate a contribution of osteocytes to bone (re)modeling.

Modeling Osteocyte Network Formation: Healthy and Cancerous Environments

Advanced cancers, such as prostate and breast cancers, commonly metastasize to bone. In the bone matrix, dendritic osteocytes form a spatial network allowing communication between osteocytes and the osteoblasts located on the bone surface. This communication network facilitates coordinated bone remodeling. In the presence of a cancerous microenvironment, the topology of this network changes. In those situations, osteocytes often appear to be either overdifferentiated (i.e., there are more dendrites than healthy bone) or underdeveloped (i.e., dendrites do not fully form). In addition to structural changes, histological sections from metastatic breast cancer xenografted mice show that number of osteocytes per unit area is different between healthy bone and cancerous bone. We present a stochastic agent-based model for bone formation incorporating osteoblasts and osteocytes that allows us to probe both network structure and density of osteocytes in bone. Our model both allows for the simulation of our spatial network model and analysis of mean-field equations in the form of integro-partial differential equations. We considered variations of our model to study specific physiological hypotheses related to osteoblast differentiation; for example predicting how changing biological parameters, such as rates of bone secretion, rates of cancer formation, and rates of osteoblast differentiation can allow for qualitatively different network topologies. We then used our model to explore how commonly applied therapies such as bisphosphonates (e.g., zoledronic acid) impact osteocyte network formation.

A level-set method for the evolution of cells and tissue during curvature-controlled growth

Most biological tissues grow by the synthesis of new material close to the tissue's interface, where spatial interactions can exert strong geometric influences on the local rate of growth. These geometric influences may be mechanistic or cell behavioural in nature. The control of geometry on tissue growth has been evidenced in many in vivo and in vitro experiments, including bone remodelling, wound healing, and tissue engineering scaffolds. In this paper, we propose a generalisation of a mathematical model that captures the mechanistic influence of curvature on the joint evolution of cell density and tissue shape during tissue growth. This generalisation allows us to simulate abrupt topological changes such as tissue fragmentation and tissue fusion, as well as three dimensional cases, through a level-set-based method. The level-set method developed introduces another Eulerian field than the level-set function. This additional field represents the surface density of tissue-synthesising cells, anticipated at future locations of the interface. Numerical tests performed with this level-set-based method show that numerical conservation of cells is a good indicator of simulation accuracy, particularly when cusps develop in the tissue's interface. We apply this new model to several situations of curvature-controlled tissue evolutions that include fragmentation and fusion.

3D analysis of the osteonal and interstitial tissue in human radii cortical bone

Human cortical bone has a complex hierarchical structure that is periodically remodelled throughout a lifetime. This microstructure dictates the mechanical response of the tissue under a critical load. If only some structural features, such as the different porosities observed in bone, are primarily studied, then investigations may not fully consider the osteonal systems in three-dimensions (3D). Currently, it is difficult to differentiate osteons from interstitial tissue using standard 3D characterization methods. Synchrotron radiation micro-computed tomography (SR-μCT) in the phase contrast mode is a promising method for the investigation of osteons. In the current study, SR-μCT imaging was performed on cortical bone samples harvested from eight human radii (female, 50-91 y.o.). The images were segmented to identify Haversian canals, osteocyte lacunae, micro-cracks, as well as osteons. The significant correlation between osteonal and Haversian canal volume fraction highlights the role of the canals as sites where bone remodelling is initiated. The results showed that osteocyte lacunae morphometric parameters depend on their distance to cement lines, strongly suggesting the evolution of biological activity from the beginning to the end of the remodelling process. Thus, the current study provides new data on 3D osteonal morphometric parameters and their relationships with other structural features in humans.

Melorheostotic Bone Lesions Caused by Somatic Mutations in MAP2K1 Have Deteriorated Microarchitecture and Periosteal Reaction

Melorheostosis is a rare non-hereditary condition characterized by dense hyperostotic lesions with radiographic "dripping candle wax" appearance. Somatic activating mutations in MAP2K1 have recently been identified as a cause of melorheostosis. However, little is known about the development, composition, structure, and mechanical properties of the bone lesions. We performed a multi-method phenotype characterization of material properties in affected and unaffected bone biopsy samples from six melorheostosis patients with MAP2K1 mutations. On standard histology, lesions show a zone with intensively remodeled osteonal-like structure and prominent osteoid accumulation, covered by a shell formed through bone apposition, consisting of compact multi-layered lamellae oriented parallel to the periosteal surface and devoid of osteoid. Compared with unaffected bone, melorheostotic bone has lower average mineralization density measured by quantitative backscattered electron imaging (CaMean: -4.5%, p = 0.04). The lamellar portion of the lesion is even less mineralized, possibly because the newly deposited material has younger tissue age. Affected bone has higher porosity by micro-CT, due to increased tissue vascularity and elevated 2D-microporosity (osteocyte lacunar porosity: +39%, p = 0.01) determined on quantitative backscattered electron images. Furthermore, nano-indentation modulus characterizing material hardness and stiffness was strictly dependent on tissue mineralization (correlation with typical calcium concentration, CaPeak: r = 0.8984, p = 0.0150, and r = 0.9788, p = 0.0007, respectively) in both affected and unaffected bone, indicating that the surgical hardness of melorheostotic lesions results from their lamellar structure. The results suggest a model for pathophysiology of melorheostosis caused by somatic activating mutations in MAP2K1, in which the genetically induced gradual deterioration of bone microarchitecture triggers a periosteal reaction, similar to the process found to occur after bone infection or local trauma, and leads to an overall cortical outgrowth. The micromechanical properties of the lesions reflect their structural heterogeneity and correlate with local variations in mineral content, tissue age, and remodeling rates, in the same way as normal bone. © 2018 American Society for Bone and Mineral Research.

Chronic stress and body composition disorders: implications for health and disease

Recent studies have suggested that body composition is key to health and disease. First, fat tissue is a complex, essential, and highly active metabolic and endocrine organ that responds to afferent signals from traditional hormone systems and the central nervous system but also expresses and secretes factors with important endocrine, metabolic, and immune functions. Second, skeletal muscle mass is an important predictor of health in adult life, while severe mass loss has been associated with the frailty of old age. Studies have shown that skeletal muscle is also an important endocrine organ that secretes factors with autocrine, paracrine, or endocrine actions, which have been associated with inflammatory processes. Third, the bone is also a systemic endocrine regulator playing a pivotal role in health and disease. Finally, proper hydration in humans has been neglected as a health factor, especially in adults. Chronic stress and stress hormone hypersecretion alone or associated with distinct disorders, such as anxiety, depression, obesity, metabolic syndrome, autoimmune disorders, type 2 diabetes mellitus, and polycystic ovary syndrome (PCOS), have been associated with psychological and somatic manifestations, typically, increased fat mass, osteosarcopenia/frailty, cellular dehydration, and chronic systemic inflammation. This review aims to provide new insights into the newly developed concept of stress-related osteosarcopenic obesity and its prevention.

Mean-field approach to evolving spatial networks, with an application to osteocyte network formation

We consider evolving networks in which each node can have various associated properties (a state) in addition to those that arise from network structure. For example, each node can have a spatial location and a velocity, or it can have some more abstract internal property that describes something like a social trait. Edges between nodes are created and destroyed, and new nodes enter the system. We introduce a "local state degree distribution" (LSDD) as the degree distribution at a particular point in state space. We then make a mean-field assumption and thereby derive an integro-partial differential equation that is satisfied by the LSDD. We perform numerical experiments and find good agreement between solutions of the integro-differential equation and the LSDD from stochastic simulations of the full model. To illustrate our theory, we apply it to a simple model for osteocyte network formation within bones, with a view to understanding changes that may take place during cancer. Our results suggest that increased rates of differentiation lead to higher densities of osteocytes, but with a smaller number of dendrites. To help provide biological context, we also include an introduction to osteocytes, the formation of osteocyte networks, and the role of osteocytes in bone metastasis.

Modeling the Effect of Curvature on the Collective Behavior of Cells Growing New Tissue

The growth of several biological tissues is known to be controlled in part by local geometrical features, such as the curvature of the tissue interface. This control leads to changes in tissue shape that in turn can affect the tissue's evolution. Understanding the cellular basis of this control is highly significant for bioscaffold tissue engineering, the evolution of bone microarchitecture, wound healing, and tumor growth. Although previous models have proposed geometrical relationships between tissue growth and curvature, the role of cell density and cell vigor remains poorly understood. We propose a cell-based mathematical model of tissue growth to investigate the systematic influence of curvature on the collective crowding or spreading of tissue-synthesizing cells induced by changes in local tissue surface area during the motion of the interface. Depending on the strength of diffusive damping, the model exhibits complex growth patterns such as undulating motion, efficient smoothing of irregularities, and the generation of cusps. We compare this model with in vitro experiments of tissue deposition in bioscaffolds of different geometries. By including the depletion of active cells, the model is able to capture both smoothing of initial substrate geometry and tissue deposition slowdown as observed experimentally.

Bone remodeling as a spatial evolutionary game

Bone remodeling is a complex process involving cell-cell interactions, biochemical signaling and mechanical stimuli. Early models of the biological aspects of remodeling were non-spatial and focused on the local dynamics at a fixed location in the bone. Several spatial extensions of these models have been proposed, but they generally suffer from two limitations: first, they are not amenable to analysis and are computationally expensive, and second, they neglect the role played by bone-embedded osteocytes. To address these issues, we developed a novel model of spatial remodeling based on the principles of evolutionary game theory. The analytically tractable framework describes the spatial interactions between zones of bone resorption, bone formation and quiescent bone, and explicitly accounts for regulation of remodeling by bone-embedded, mechanotransducing osteocytes. Using tools from the theory of interacting particle systems we systematically classified the different dynamic regimes of the spatial model and identified regions of parameter space that allow for global coexistence of resorption, formation and quiescence, as observed in physiological remodeling. In coexistence scenarios, three-dimensional simulations revealed the emergence of sponge-like bone clusters. Comparison between spatial and non-spatial dynamics revealed substantial differences and suggested a stabilizing role of space. Our findings emphasize the importance of accounting for spatial structure and bone-embedded osteocytes when modeling the process of bone remodeling. Thanks to the lattice-based framework, the proposed model can easily be coupled to a mechanical model of bone loading.

Governing Equations of Tissue Modelling and Remodelling: A Unified Generalised Description of Surface and Bulk Balance

Several biological tissues undergo changes in their geometry and in their bulk material properties by modelling and remodelling processes. Modelling synthesises tissue in some regions and removes tissue in others. Remodelling overwrites old tissue material properties with newly formed, immature tissue properties. As a result, tissues are made up of different “patches”, i.e., adjacent tissue regions of different ages and different material properties, within evolving boundaries. In this paper, generalised equations governing the spatio-temporal evolution of such tissues are developed within the continuum model. These equations take into account nonconservative, discontinuous surface mass balance due to creation and destruction of material at moving interfaces, and bulk balance due to tissue maturation. These equations make it possible to model patchy tissue states and their evolution without explicitly maintaining a record of when/where resorption and formation processes occurred. The time evolution of spatially averaged tissue properties is derived systematically by integration. These spatially-averaged equations cannot be written in closed form as they retain traces that tissue destruction is localised at tissue boundaries. The formalism developed in this paper is applied to bone tissues, which exhibit strong material heterogeneities due to their slow mineralisation and remodelling processes. Evolution equations are proposed in particular for osteocyte density and bone mineral density. Effective average equations for bone mineral density (BMD) and tissue mineral density (TMD) are derived using a mean-field approximation. The error made by this approximation when remodelling patchy tissue is investigated. The specific signatures of the time evolution of BMD or TMD during remodelling events are exhibited. These signatures may provide a way to detect remodelling events at lower, unseen spatial resolutions from microCT scans.

Quantifying the osteocyte network in the human skeleton

Osteocytes form an extensive cellular network throughout the hard tissue matrix of the skeleton, which is known to regulate skeletal structure. However due to limitations in imaging techniques, the magnitude and complexity of this network remain undefined. We have used data from recent papers obtained by new imaging techniques, in order to estimate absolute and relative quantities of the human osteocyte network and form a more complete understanding of the extent and nature of this network. We estimate that the total number of osteocytes within the average adult human skeleton is ~42 billion and that the total number of osteocyte dendritic projections from these cells is ~3.7 trillion. Based on prior measurements of canalicular density and a mathematical model of osteocyte dendritic process branching, we calculate that these cells form a total of 23 trillion connections with each other and with bone surface cells. We estimate the total length of all osteocytic processes connected end-to-end to be 175,000 km. Furthermore, we calculate that the total surface area of the lacuno-canalicular system is 215 m(2). However, the residing osteocytes leave only enough space for 24 mL of extracellular fluid. Calculations based on measurements in lactation-induced murine osteocytic osteolysis indicate a potential total loss of ~16,000 mm(3) (16 mL) of bone by this process in the human skeleton. Finally, based on the average speed of remodelling in the adult, we calculate that 9.1 million osteocytes are replenished throughout the skeleton on a daily basis, indicating the dynamic nature of the osteocyte network. We conclude that the osteocyte network is a highly complex communication network, and is much more vast than commonly appreciated. It is at the same order of magnitude as current estimates of the size of the neural network in the brain, even though the formation of the branched network differs between neurons and osteocytes. Furthermore, continual replenishment of large numbers of osteocytes in the process of remodelling allows therapeutic changes to the continually renewed osteoblast population to be rapidly incorporated into the skeleton.