Skip to main content
Log in

Bone refilling in cortical basic multicellular units: insights into tetracycline double labelling from a computational model

  • Original Paper
  • Published:
Biomechanics and Modeling in Mechanobiology Aims and scope Submit manuscript

Abstract

Bone remodelling is carried out by ‘bone multicellular units’ (\(\text{ BMU }\)s) in which active osteoclasts and active osteoblasts are spatially and temporally coupled. The refilling of new bone by osteoblasts towards the back of the \(\text{ BMU }\) occurs at a rate that depends both on the number of osteoblasts and on their secretory activity. In cortical bone, a linear phenomenological relationship between matrix apposition rate and \(\text{ BMU }\) cavity radius is found experimentally. How this relationship emerges from the combination of complex, nonlinear regulations of osteoblast number and secretory activity is unknown. Here, we extend our previous mathematical model of cell development within a single cortical \(\text{ BMU }\) to investigate how osteoblast number and osteoblast secretory activity vary along the \(\text{ BMU }\)’s closing cone. The mathematical model is based on biochemical coupling between osteoclasts and osteoblasts of various maturity and includes the differentiation of osteoblasts into osteocytes and bone lining cells, as well as the influence of \(\text{ BMU }\) cavity shrinkage on osteoblast development and activity. Matrix apposition rates predicted by the model are compared with data from tetracycline double labelling experiments. We find that the linear phenomenological relationship observed in these experiments between matrix apposition rate and \(\text{ BMU }\) cavity radius holds for most of the refilling phase simulated by our model, but not near the start and end of refilling. This suggests that at a particular bone site undergoing remodelling, bone formation starts and ends rapidly, supporting the hypothesis that osteoblasts behave synchronously. Our model also suggests that part of the observed cross-sectional variability in tetracycline data may be due to different bone sites being refilled by \(\text{ BMU }\)s at different stages of their lifetime. The different stages of a \(\text{ BMU }\)’s lifetime (such as initiation stage, progression stage, and termination stage) depend on whether the cell populations within the \(\text{ BMU }\) are still developing or have reached a quasi-steady state whilst travelling through bone. We find that due to their longer lifespan, active osteoblasts reach a quasi-steady distribution more slowly than active osteoclasts. We suggest that this fact may locally enlarge the Haversian canal diameter (due to a local lack of osteoblasts compared to osteoclasts) near the \(\text{ BMU }\)’s point of origin.

This is a preview of subscription content, log in via an institution to check access.

Access this article

Price excludes VAT (USA)
Tax calculation will be finalised during checkout.

Instant access to the full article PDF.

Institutional subscriptions

Fig. 1
Fig. 2
Fig. 3
Fig. 4
Fig. 5

Similar content being viewed by others

Notes

  1. Over-refilling occurs mainly on trabecular and endosteal surfaces. However, a net bone gain is possible in cortical \(\text{ BMU }\)s that follow a pre-existing Haversian canal whilst reducing the canal’s size (‘type II osteon’) (Parfitt 1983; Robling and Stout 1998)

  2. When \(R=R_H\), osteoblast secretory activity is assumed to stop (see the ‘Osteoblasts’ section below), but refilling may stop before \(R\) reaches \(R_H\) if the population of osteoblasts goes extinct.

References

  • Aubin JE (2008) Mesenchymal stem cells and osteoblast differentiation. In: Bilezikian JP, Raisz LG, Martin TJ (eds) Principles of bone biology, vol 1, 3rd edn. Academic Press, San Diego, pp 85–107

    Chapter  Google Scholar 

  • Buenzli PR, Jeon J, Pivonka P, Smith DW, Cummings PT (2012a) Investigation of bone resorption within a cortical basic multicellular unit using a lattice-based computational model. Bone 50:378–389

    Article  Google Scholar 

  • Buenzli PR, Pivonka P, Gardiner BS, Smith DW (2012b) Modelling the anabolic response of bone using a cell population model. J Theor Biol 307:42–52

    Article  MathSciNet  Google Scholar 

  • Buenzli PR, Pivonka P, Smith DW (2011) Spatio-temporal dynamics of cell distribution in bone multicellular units. Bone 48:918–926

    Article  Google Scholar 

  • Buenzli PR, Thomas CDL, Clement JG, Pivonka P (2012c) Endocortical bone loss in osteoporosis: the role of bone surface availability. preprint: arXiv:1206.6071

  • Cohen J, Harris WH (1958) The three-dimensional anatomy of haversian systems. J Bone Joint Surg Am 40:419–434

    Google Scholar 

  • Cooper DM, Thomas CDL, Clement JG, Hallgrimsson B (2006) Three-dimensional microcomputed tomography imaging of basic multicellular unit-related resorption spaces in human cortical bone. Anat Rec A Discov Mol Cell Evol Biol 288A:806–816

    Article  Google Scholar 

  • Dallas SL, Bonewald LF (2010) Dynamics of the transition from osteoblast to osteocyte. Ann N Y Acad Sci 1192:437–443

    Article  Google Scholar 

  • Eriksen EF, Axelrod DW, Melsen F (1994) Bone histomorphometry. Raven Press, New York

    Google Scholar 

  • Franz-Odendaal TA, Hall BK, Witten PE (2006) Buried alive: how osteoblasts become osteocytes. Dev Dyn 235:176–190

    Article  Google Scholar 

  • Frost HM (1969) Tetracycline-based histological analysis of bone remodeling. Calcif Tissue Res 3:211–237

    Article  Google Scholar 

  • Gori F, Hofbauer LC, Dunstan CR, Spelsberg TC, Khosla S, Riggs BL (2000) The expression of osteoprotegerin and RANK ligand and the support of osteoclast formation by stromal-osteoblast lineage cells is developmentally regulated. Endocrinology 141:4768–4776

    Article  Google Scholar 

  • Harada S-I, Rodan GA (2003) Control of osteoblast function and regulation of bone mass. Nature 423:349–355

    Article  Google Scholar 

  • Hood SJ, Morrin LA, Parks NJ, Bugreeff SE (1985) Bone histomorphometry: analysis of results of juvenile cortical bone study. Annual Report Laboratory for Energy-Related Health Research, UC Davis, pp 70–75

  • Heaney RP (1994) The bone-remodeling transient: implications for the interpretation of clinical studies of bone mass change. J Bone Miner Res 9:1515–1523

    Article  Google Scholar 

  • Hernandez CJ, Beaupré GS, Carter DR (2003) A theoretical analysis of the changes in basic multicellular unit activity at menopause. Bone 32:357–363

    Google Scholar 

  • Ilnicki L, Epker BN, Frost HM, Hattner R (1966) The radial rate of osteon closure evaluated by means of in vivo tetracycline labeling in beagle dog rib. J Lab Clin Med 67:447–454

    Google Scholar 

  • Iqbal J, Sun L, Zaidi M (2009) Coupling bone degradation to formation. Nat Med 15:729–731

    Article  Google Scholar 

  • Jaworski ZFG, Hooper C (1980) Study of cell kinetics within evolving secondary haversian systems. J Anat Lond 131:91–102

    Google Scholar 

  • Jaworski ZFG, Duck B, Sekaly G (1981) Kinetics of osteoclasts and their nuclei in evolving secondary Haversian systems. J Anat 133:397–405

    Google Scholar 

  • Ji B, Genever PG, Patton RJ, Putra D, Fagan MJ (2012) A novel mathematical model of bone remodelling cycles for trabecular bone at the cellular level. Biomech Model Mechanobiol 11:973–982

    Article  Google Scholar 

  • Jilka RL (2002) Osteoblast progenitor fate and age-related bone loss. J Musculoskel Neuronal Interact 2:581–583

    Google Scholar 

  • Jilka RL (2003) Biology of the basic multicellular unit and the pathophysiology of osteoporosis. Med Pediatr Oncol 41:182–185

    Article  Google Scholar 

  • Jones SH (1974) Secretory territories and rates of matrix production of osteoblasts. Calcif Tissue Res 14:309–315

    Article  Google Scholar 

  • Jowsey J (1966) Studies of haversian systems in man and some animals. J Anat 100:857–864

    Google Scholar 

  • Komarova SV, Smith RJ, Dixon SJ, Sims SM, Wahl LM (2003) Mathematical model predicts a critical role for osteoclast autocrine regulation in the control of bone remodeling. J Theor Biol 229:293–309

    Google Scholar 

  • Lee WR (1964) Appositional bone formation in canine bone: a quantitative microscopic study using tetracycline markers. J Anat Lond 98:665–677

    Google Scholar 

  • Ma YL, Cain RL, Halladay DL, Yang X, Zeng Q, Miles R, Chandrasekhar S, Martin TJ, Onyia JE (2001) Catabolic effects of continuous human PTH (1–38) in vivo is associated with sustained stimulation of RANKL and inhibition of osteoprotegerin and gene-associated bone formation. Endocrinology 142:4047–4054

    Article  Google Scholar 

  • Manolagas SC (2000) Birth and death of bone cells: basic regulatory mechanisms and implications for the pathogenesis and treatment of osteoporosis. Endocr Rev 21:115–137

    Article  Google Scholar 

  • Manson JD, Waters NE (1965) Observations on the rate of maturation of the cat osteon. J Anat Lond 99:539–549

    Google Scholar 

  • Marotti G (2000) The osteocyte as a wiring transmission system. J Muskuloskel Neuronal Interact 1:133–136

    Google Scholar 

  • Marotti G, Zallone AZ (1976) Number, size and arrangement of osteoblasts in osteons at different stages of formation. Calcif Tissue Int 21(suppl):96–101

    Google Scholar 

  • Martin RB (2000) Does osteocyte formation cause the nonlinear refilling of osteons? Bone 26:71–78

    Article  Google Scholar 

  • Martin RB, Burr DB, Sharkey NA (1998) Skeletal tissue mechanics. Springer, New York

    Book  Google Scholar 

  • Martin RB, Dannucci GA, Hood SJ (1987) Bone apposition rate differences in osteonal and trabecular bone. Trans Orth Res Soc 12:178

    Google Scholar 

  • Martin TJ (2004) Paracrine regulation of osteoclast formation and activity: milestones in discovery. J Musculoskel Neuronal Ineract 4: 243–253

    Google Scholar 

  • Metz LN, Martin RB, Turner AS (2003) Histomorphometric analysis of the effects of osteocyte density on osteonal morphology and remodeling. Bone 33:753–759

    Article  Google Scholar 

  • Parfitt AM (1994) Osteonal and hemi-osteonal remodeling: the spatial and temporal framework for signal traffic in adult human bone. J Cell Biochem 55:273–286

    Article  Google Scholar 

  • Parfitt AM (1998) Osteoclast precursors as leukocytes: importance of the area code. Bone 23:491–494

    Article  Google Scholar 

  • Parfitt MA (1983) The physiological and clinical significance of bone histomorphometric data. In: Recker RR (ed) Bone histomorphometry: techniques and interpretation. CRC Press, Boca Raton, pp 143–223

    Google Scholar 

  • Pazzaglia UE, Congiu T, Franzetti E, Marchese M, Spagnuolo F, Di Mascio L, Zarattini G (2012a) A model of osteoblast-osteocyte kinetics in the development of secondary osteons in rabbits. J Anat 220:372–383

    Article  Google Scholar 

  • Pazzaglia UE, Congiu T, Marchese M, Dell’Orbo C (2010) The shape modulation of osteoblast-osteocyte transformation and its correlation with the fibrillar organization in secondary osteons. Cell Tissue Res 340:533–540

    Article  Google Scholar 

  • Pazzaglia UE, Congiu T, Marchese M, Spagnuolo F, Quacci D (2012b) Morphometry and patterns of lamellar bone in human Haversian systems. Anat Rec 295:1421–1429

    Article  Google Scholar 

  • Pazzaglia UE, Congiu T, Zarattini G, Marchese M, Quacci D (2011) The fibrillar organisation of the osteon and cellular aspects of its development: a morphological study using the SEM fractured cortex technique. Anat Sci Int 86:128–134

    Article  Google Scholar 

  • Pivonka P, Buenzli PR, Scheiner S, Hellmich R, Dunstan CR (2012) The influence of bone surface availability in bone remodelling—a mathematical model including coupled geometrical and biomechanical regulations of bone cells. Eng Struct 47:134–147

    Article  Google Scholar 

  • Pivonka P, Zimak J, Smith DW, Gardiner BS, Dunstan CR, Sims NA, Martin TJ, Mundy GR (2008) Model structure and control of bone remodeling: a theoretical study. Bone 43:249–263

    Article  Google Scholar 

  • Polig E, Jee WSS (1990) A model of osteon closure in cortical bone. Calcif Tissue Int 47:261–269

    Article  Google Scholar 

  • Qiu S, Rao DS, Palnitkar S, Parfitt AM (2010) Dependence of bone yield (volume of bone formed per unit of cement surface area) on resorption cavity size during osteonal remodeling in human rib: implications for osteoblast function and the pathogenesis of age-related bone loss. J Bone Miner Res 25:423–430

    Article  Google Scholar 

  • Robling AG, Stout SD (1998) Morphology of the drifting osteon. Cells Tissues Organs 164:192–204

    Article  Google Scholar 

  • Roodman GD (1999) Cell biology of the osteoclast. Exp Hematol 27:1229–1241

    Article  Google Scholar 

  • Rumpler M, Woesz A, Dunlop JWC, van Dongen JT, Fratzl P (2008) The effect of geometry on three-dimensional tissue growth. J R Soc Interface 5:1173–1180

    Article  Google Scholar 

  • Ryser MD, Komarova SV, Nigam N (2010) The cellular dynamics of bone remodelling: a mathematical model. SIAM J Appl Math 70:1899–1921

    Article  MATH  MathSciNet  Google Scholar 

  • Ryser MD, Nigma N, Komarova SV (2009) Mathematical modeling of spatio-temporal dynamics of a single bone multicellular unit. J Bone Miner Res 24:860–870

    Article  Google Scholar 

  • Seeman E (2008) Modeling and remodeling: the cellular machinery responsible for the gain and loss of bone’s material and structural strength. In: Bilezikian JP, Raisz LG, Martin TJ (eds) Principles of bone biology, vol 1, 3rd edn. Academic Press, San Diego, pp 3–28

    Google Scholar 

  • Tang Y et al (2009) TGF \(\upbeta \)1-induced migration of bone mesenchymal stem cells couples bone resorption with formation. Nat Med 15:757–766

    Article  Google Scholar 

  • Thomas GP, Baker SU, Eisman JA, Gardiner EM (2001) Changing RANKL/OPG mRNA expression in differentiating murine primary osteoblasts. J Endocrinol 170:451–460

    Article  Google Scholar 

  • van Bezooijen RL, Papapoulos SE, Hamdy NAT, Löwik CWGM (2008) SOST/Sclerostin: an osteocyte-derived inhibitor of bone formation that antagonizes canonical Wnt signaling. In: Bilezikian JP, Raisz LG, Martin TJ (eds) Principles of bone biology, vol 1, 3rd edn. Academic Press, San Diego, pp 139–152

    Chapter  Google Scholar 

  • van Kampen NG (2007) Stochastic processes in physics and chemistry, 3rd edn. Elsevier, Amsterdam

    Google Scholar 

  • van Oers RFM, Ruimerman R, Tanck E, Hilbers PAJ, Huiskes R (2008) A unified theory for osteonal and hemi-osteonal remodeling. Bone 42(2):250–259

    Article  Google Scholar 

  • van Oers RFM, van Rietbergen B, Ito K, Huiskes R, Hilbers PAJ (2011) Simulations of trabecular remodeling and fatigue: Is remodeling helpful or harmful? Bone 48:1210–1215

    Article  Google Scholar 

  • Volpi G, Palazzini S, Cané V, Remaggi F, Muglia MA (1981) Morphometric analysis of osteoblast dynamics in the chick embryo tibia. Anat Embryol 162:393–401

    Article  Google Scholar 

  • Wolfram Research, Inc (2008) Mathematica Edition: Version 7.0. Wolfram Research, Inc., Champaign

  • Zallone AZ (1977) Relationships between shape and size of the osteoblasts and the accretion rate of trabecular bone surfaces. Anat Embryol 152:65–72

    Article  Google Scholar 

Download references

Acknowledgments

Financial support by the Australian Research Council (Project number DP0988427 (PP, DWS) and Project number DE130101191 (PRB)) is gratefully acknowledged.

Author information

Authors and Affiliations

Authors

Corresponding author

Correspondence to Pascal R. Buenzli.

Appendix: Model description

Appendix: Model description

In this appendix, we detail some aspects of the model presented in the Methods section (Sect. 2). The governing equations of all cells and signalling molecules considered in the model are presented, as well as a table listing all the model parameters.

Osteoblast secretory activity. In cortical \(\text{ BMU }\)s, the secretory activity of single osteoblasts varies as refilling proceeds (Marotti et al. 1976; Volpi et al. 1981; Parfitt 1983; Martin et al. 1998; Zallone 1977). The secretion rate depends on the amount of organelles involved in glycoprotein and protein synthesis, and so is likely to depend on the cell’s protoplasmic volume. The precise factors that influence osteoid secretion rate are currently not well known, but may include nutrients, hormones, and regulatory molecules (Parfitt 1994; Manolagas 2000; Aubin 2008), and the local curvature of the cavity (Rumpler et al. 2008; Qiu et al. 2010).

To account for the variation in osteoid secretion rate during refilling in our model, we assume that \({k_\text{ form }}\) depends on the current radius of the cavity [Eq. (8)]: \({k_\text{ form }}(x,t)=\widetilde{k}_\text{ form }\big (R(x,t)\big )\), with \(\widetilde{k}_\text{ form }({R_\text{ H }})=0\), where \({R_\text{ H }}\) is the Haversian canal radius at which osteoblasts have become (nonsynthesising) bone lining cells (Fig. 1b). Marotti et al. (1976) reported osteoid secretion rates at three different radii in canine cortical \(\text{ BMU }\)s. We add to these data a zero secretion rate at \({R_\text{ H }}^\text{ dog }\approx 15\ \upmu \text{ m }\), the average Haversian canal radius in dogs (Parfitt 1983) (see Table 1). The data of Table 1 are used to generate pseudo-human data by appropriate scalings as explained in the paragraph ‘Scaling data from animal models to humans’ below. The pseudo-human data are then interpolated by a continuous piecewise polynomial curve \(\widetilde{k}_\text{ form }(R)\) that has a sharp stepwise increase between \(R={R_\text{ H }}=20\ \upmu \text{ m }\) and \(R=22~\upmu \text{ m }\), and a linear slope for \(R\ge 22~\upmu \text{ m }\), as depicted in Fig. 6. Note that all the active osteoblasts in the cross section at \(x\) secrete at the same rate \({k_\text{ form }}(x,t)\). As long as the cavity radius \(R\) is larger than the ‘target’ Haversian canal radius \(R_H\), the activity of osteoblasts is therefore assumed continuous and uninterrupted, which is probably the most common pattern of bone refilling in cortical \(\text{ BMU }\), at least in an average sense (see also comments in Sect. 4) (Parfitt 1983; Martin et al. 1998; Pazzaglia et al. 2011, 2012b).

Table 1 Osteoid secretion rate versus \(\text{ BMU }\) cavity radius (in dogs)
Fig. 6
figure 6

Assumed functional relationship \(\widetilde{k}_\text{ form }(R)\) (8) between osteoid secretion rate \({k_\text{ form }}\) and \(\text{ BMU }\) cavity radius \(R\) (solid line). Data points (crosses) are taken from Table 1 and scaled according to Eqs. (20), (21)

Osteoblast-to-bone lining cell transition. At the end of the formation phase at the back of a \(\text{ BMU }\), active osteoblasts terminally differentiate into bone lining cells, thus forming a cellular layer on the newly formed bone surface (Parfitt 1983; Martin et al. 1998). To account for such a transition in the model whilst ensuring the spatial exclusion of active osteoblasts and bone lining cells, we assume that active osteoblasts whose osteoid secretion rate \({k_\text{ form }}\) is rapidly dropping become bone lining cells. This occurs between cavity radii 20–22 \(\upmu \text{ m }\) (Fig. 6).

Osteoblast-to-osteocyte transition. A fraction of active osteoblasts is entrapped in the new bone matrix during bone formation and becomes osteocytes (Marotti 2000; Franz-Odendaal et al. 2006; Dallas and Bonewald 2010). This is represented by the sink term \(-\sigma _\text{ OCY }^\text{ prod. }(x,t)\) in Eq. (10). Here, we determine what rate of osteocyte generation \(\sigma _\text{ OCY }^\text{ prod. }(x,t)\) (density produced per unit time) reproduces a given osteocyte density in cortical \(\text{ BMU }\)s, \(\text{ OCY }_\text{ exp. }(x,r)\) (e.g. one determined experimentally). We allow for a possible dependence upon the radial coordinate \(r\) in the cross section (Fig. 1b).

The correspondence between \(\text{ OCY }_\text{ exp. }\) and \(\sigma _\text{ OCY }^\text{ prod. }\) is found by calculating the number of osteocytes in the slice of thickness \(\varDelta x\) at time \(t,\,\varDelta N_\text{ OCY }(x,t)\). This number increases during bone deposition from zero at the onset of bone refilling (\(t=t_\text{ F }\)) (Fig. 1b, stage 1) to a constant final number at the end of bone refilling (Fig. 1b, stage 4). At any time \(t\) in between, \(\varDelta N_\text{ OCY }(x,t)\) is given by (i) integrating spatially the given density \(\text{ OCY }_\text{ exp. }\) from \({R_\text{ c }}\) to \(R(x,t)\) (Fig. 1b) or (ii) integrating temporally the production rate \(\sigma _\text{ OCY }^\text{ prod. }\) from \(t_\text{ F }\) to t. Equating these integrations gives:

$$\begin{aligned} 2\pi \!\!\int \limits _{R(x,t)}^{R_\text{ c }}\mathrm d r\, r\, \text{ OCY }_\text{ exp. }(x,r) = \int \limits _{t_\text{ F }}^t \mathrm d t^{\prime } (\sigma _\text{ OCY }^\text{ prod. }\,S_{{\text{ OB }}_{\text{ a }}})(x,t^{\prime }). \end{aligned}$$

The dependence upon \(t\) occurs in the integral bounds. From standard calculus, the differentiation of the above equation with respect to t is given by the integrant evaluated at the time-dependent bound times the derivative of the bound (with a negative sign for the lower bound). Thus:

$$\begin{aligned}&-2\pi R(x,t) \text{ OCY }_\text{ exp. }\left( x, R(x,t)\right) \tfrac{\partial }{\partial t} R(x,t) \\&\qquad = \sigma _\text{ OCY }^\text{ prod. }(x,t)\,S_{{\text{ OB }}_{\text{ a }}}(x,t). \end{aligned}$$

Using Eq. (6) to substitute \(\frac{\partial }{\partial t}R\) and then Eq. (11) to substitute \(\rho _{{\text{ OB }}_{\text{ a }}}/S_{{\text{ OB }}_{\text{ a }}}\), one finds:

$$\begin{aligned} \sigma _\text{ OCY }^\text{ prod. }(x,t)&= 2\pi \,R\ \text{ OCY }_\text{ exp. }\,{k_\text{ form }}\frac{\rho _{{\text{ OB }}_{\text{ a }}}}{S_{{\text{ OB }}_{\text{ a }}}}\nonumber \\&= {k_\text{ form }}(x,t)\, {{\text{ OB }}_{\text{ a }}}(x,t)\ \text{ OCY }_\text{ exp. }\left( x,R(x,t)\right) \nonumber \\ \end{aligned}$$
(17)

Geometric influence of BMU cavity shrinkage. The tendency of cavity shrinkage to concentrate the density of osteoblasts is represented by the source term \(G\,{{\text{ OB }}_{\text{ a }}}\) in Eq. (10). To determine the form of \(G(x,t)\), we consider the number of active osteoblasts in the slice at \(x\) of thickness \(\varDelta x,\,\varDelta N_{{\text{ OB }}_{\text{ a }}}(x,t) = {{\text{ OB }}_{\text{ a }}}*S_{{\text{ OB }}_{\text{ a }}}\ \varDelta x\). The rate of change of this number is \(\tfrac{\partial }{\partial t}\varDelta N_{{\text{ OB }}_{\text{ a }}}(x,t) = \left[ (\tfrac{\partial }{\partial t}{{\text{ OB }}_{\text{ a }}}) S_{{\text{ OB }}_{\text{ a }}}+ {{\text{ OB }}_{\text{ a }}}\, \tfrac{\partial }{\partial t}S_{{\text{ OB }}_{\text{ a }}}\right] \varDelta x\). Within the slice, \(\varDelta N_{{\text{ OB }}_{\text{ a }}}\) can only be modified by biological creation or elimination, which is described by the source and sink terms \(\mathcal D _{\text{ OB }_\text{ p }}\,{\text{ OB }_\text{ p }}- A_{{\text{ OB }}_{\text{ a }}}{{\text{ OB }}_{\text{ a }}}- \sigma _\text{ OCY }^\text{ prod. }\equiv \tfrac{\partial }{\partial t}{{\text{ OB }}_{\text{ a }}}\!-\! G\ {{\text{ OB }}_{\text{ a }}}\) in Eq. (10). Hence:

$$\begin{aligned} \tfrac{\partial }{\partial t}\varDelta N_{{\text{ OB }}_{\text{ a }}}&= \left[ (\tfrac{\partial }{\partial t}{{\text{ OB }}_{\text{ a }}}) S_{{\text{ OB }}_{\text{ a }}}+ {{\text{ OB }}_{\text{ a }}}\, \tfrac{\partial }{\partial t}S_{{\text{ OB }}_{\text{ a }}}\right] \varDelta x \nonumber \\&= \left[ \tfrac{\partial }{\partial t}{{\text{ OB }}_{\text{ a }}}- G\ {{\text{ OB }}_{\text{ a }}}\right] S_{{\text{ OB }}_{\text{ a }}}\ \varDelta x. \end{aligned}$$
(18)

Dividing by \({{\text{ OB }}_{\text{ a }}}*S_{{\text{ OB }}_{\text{ a }}}\ \varDelta x\), this gives:

$$\begin{aligned} G(x,t)&= - \frac{\tfrac{\partial }{\partial t} S_{{\text{ OB }}_{\text{ a }}}(x,t)}{S_{{\text{ OB }}_{\text{ a }}}(x,t)} = - \frac{\tfrac{\partial }{\partial t} R(x,t)}{R(x,t)\left( 1-\frac{h_{{\text{ OB }}_{\text{ a }}}}{2R(x,t)}\right) }\nonumber \\&= -{k_\text{ form }}(x,t)\,{{\text{ OB }}_{\text{ a }}}(x,t)\,\frac{h_{{\text{ OB }}_{\text{ a }}}}{R(x,t)}, \end{aligned}$$
(19)

where Eqs. (12) and (14) were used for the last equalities.

Scaling data from animal models to humans. To calibrate our model to a human cortical \(\text{ BMU }\), we scale \(\text{ BMU }\)-related data obtained from animal models according to known (or suspected) cross-species differences.

The \(\text{ BMU }\) cavity radius data \(R^\text{ dog }\) from Table 1 was scaled to (pseudo-)human values \(R\) by assuming a linear relationship \(R=R(R^\text{ dog })\) such that \(R({R_\text{ H }}^\text{ dog }) = {R_\text{ H }}\) and \(R({R_\text{ c }}^\text{ dog }) = {R_\text{ c }}\). This gives:

$$\begin{aligned} R(R^\text{ dog }) = {R_\text{ H }}+ \frac{{R_\text{ c }}-{R_\text{ H }}}{{R_\text{ c }}^\text{ dog }-{R_\text{ H }}^\text{ dog }}(R^\text{ dog } -{R_\text{ H }}^\text{ dog }). \end{aligned}$$
(20)

We used the typical values \({R_\text{ H }}^\text{ dog } \approx 15\,\upmu \text{ m }\) and \({R_\text{ H }}\approx 20\,\upmu \text{ m }\) for canine and human Haversian canal radii. In Parfitt (1983), the average cement line radius in dogs is \({R_\text{ c }}^\text{ dog }\approx 60~\upmu \text{ m }\). In Jowsey (1966), it is \({R_\text{ c }}^\text{ dog }\approx 77~\upmu \text{ m }\). Since one of the three \(\text{ BMU }\) cavities in Marotti et al. (1976) has a radius \(65\,\upmu \text{ m }\) in the formation phase, we take here \(70~\upmu \text{ m }\) as the representative canine cement line radius. We used the value \({R_\text{ c }}\approx 100\,\upmu \text{ m }\) for the human cement line radius (Jowsey 1966; Parfitt 1983; Martin et al. 1998).

The osteoid secretion rates \({k_\text{ form }^\text{ dog }}\) of Table 1 estimated by Marotti et al. (1976) were scaled by a factor \(\alpha _{k_\text{ form }}=1.25\) to account for higher secretion rates in humans (Polig and Jee 1990), i.e.:

$$\begin{aligned} {k_\text{ form }}= \alpha _{k_\text{ form }}{k_\text{ form }^\text{ dog }}. \end{aligned}$$
(21)

The model of Polig and Jee (Polig and Jee 1990) suggests that \(\alpha _{k_\text{ form }}\) may take values from 1.27 to 1.7. Here, the specific numerical value of the factor \(\alpha _{k_\text{ form }}\) is chosen so as to reduce the total number of active osteoblasts in the \(\text{ BMU }\) (see ‘Results and Discussion’, Sects. 34). The scaled radius and osteoid secretion rate data of Table 1 are shown in Fig. 6 (crosses).

Scalings are also applied to the data of Table 2 reporting canine osteoblast surface densities \(\rho _{{\text{ OB }}_{\text{ a }}}^\text{ dog }\) at three \(\text{ BMU }\) cavity radii. The radii \(R^\text{ dog }\) are scaled to correspond to human average \(\text{ BMU }\)s as in Eq. (20). The osteoblast surface densities \(\rho _{{\text{ OB }}_{\text{ a }}}^\text{ dog }\) are scaled to pseudo-human data as \(\rho _{{\text{ OB }}_{\text{ a }}}= \alpha _{k_\text{ form }}^{-1}\rho _{{\text{ OB }}_{\text{ a }}}^\text{ dog }\). The choice of the scaling factor \(\alpha _{k_\text{ form }}^{-1}=0.8\) is motivated by the corresponding scaling factor \(\alpha _{k_\text{ form }}=1.25\) performed on \({k_\text{ form }^\text{ dog }}\). Indeed, increased osteoid secretion rate is associated with increased cell volume and thus with a corresponding decrease in cell density (Volpi et al. 1981; Marotti et al. 1976).

Table 2 Osteoblast surface density versus \(\text{ BMU }\) cavity radius (in dogs)

Scaling of the cell density distributions in the BMU. The cell distribution profiles along \(x\) resulting from the numerical simulation are calibrated against experimental values of osteoblast surface densities and total cell numbers. To this end, we have modified the values of a number of parameters compared to the parameters used in Buenzli et al. (2011). We used a scaling scheme of the model parameters such that \({{\text{ OB }}_{\text{ a }}},\,{\text{ OB }_\text{ p }},\,{\text{ OB }_\text{ u }},\,{\text{ OC }_\text{ a }}\), and \({\text{ OB }_\text{ p }}\) would scale uniformly without affecting much their own spatial distribution and without affecting their spatial relation with other cell density profiles. We defined scaling factors \(\alpha _{{\text{ OB }}_{\text{ a }}},\,\alpha _{\text{ OB }_\text{ p }},\,\alpha _{\text{ OB }_\text{ u }},\,\alpha _{\text{ OC }_\text{ a }}\), and \(\alpha _{\text{ OC }_\text{ p }}\), such that if \(\alpha _{{\text{ OB }}_{\text{ a }}}=500\), the density of \({{\text{ OB }}_{\text{ a }}}\)s is multiplied by a factor 500 etc. Based on Eqs. (23)–(26), the new parameter values (denoted below by a prime) ensuring such scaled densities were determined by multiplying the previous parameter values (nonprimed) with an adequate combination of the above scaling factors, according to:

$$\begin{aligned} \begin{array}{ll} {\text{ OB }_\text{ u }^\text{ max }}^{\prime } = \alpha _{\text{ OB }_\text{ u }}\ {\text{ OB }_\text{ u }^\text{ max }}, &{}{\text{ OC }_\text{ p }^\text{ max }}^{\prime } = \alpha _{\text{ OC }_\text{ p }}\ {\text{ OC }_\text{ p }^\text{ max }},\\ D_{\text{ OB }_\text{ u }}^{\prime } = \frac{\alpha _{\text{ OB }_\text{ p }}}{\alpha _{\text{ OB }_\text{ u }}} D_{\text{ OB }_\text{ u }}, &{}D_{\text{ OC }_\text{ p }}^{\prime } = \frac{\alpha _{\text{ OC }_\text{ a }}}{\alpha _{\text{ OC }_\text{ p }}} D_{\text{ OC }_\text{ p }}, \\ n_{{\text{ TGF }\upbeta }}^{\text{ bone }^{\prime }} = n_{\text{ TGF }\upbeta }^\text{ bone }/\alpha _{\text{ OC }_\text{ a }}, &{}N^{\text{ RANK }^{\prime }}_{\text{ OC }_\text{ p }}= N^\text{ RANK }_{\text{ OC }_\text{ p }}/ \alpha _{\text{ OC }_\text{ p }}, \\ \beta ^{\text{ RANKL }^{\prime }}_{\text{ OB }_\text{ p }}= \beta ^\text{ RANKL }_{\text{ OB }_\text{ p }}/ \alpha _{\text{ OB }_\text{ p }}, &{}\beta ^{\text{ OPG }^{\prime }}_{{{\text{ OB }}_{\text{ a }}}} = \beta ^{\text{ OPG }}_{{{\text{ OB }}_{\text{ a }}}} / \alpha {_{{\text{ OB }}_{\text{ a }}}}, \end{array} \end{aligned}$$
(22)

where \({\text{ OB }_\text{ u }^\text{ max }}\) and \({\text{ OC }_\text{ p }^\text{ max }}\) denote parameters scaling the given distributions of \({\text{ OB }_\text{ u }}\)s and \({\text{ OC }_\text{ p }}\)s (Buenzli et al. 2011). We note that it is not possible to control the scaling of \({\text{ OB }_\text{ p }}\) independently from that of \({{\text{ OB }}_{\text{ a }}}\) without affecting their spatial relationship, and so we have always taken \(\alpha _{\text{ OB }_\text{ p }}= \alpha _{{\text{ OB }}_{\text{ a }}}\) when calibrating the model.

Governing equations of cell densities and signalling molecules concentrations. Below we summarise all the governing equations of the model. The material balance equation governing the dynamics of active osteoblasts, Eq. (10), is rewritten using the explicit expressions in Eqs. (17) and (19):

$$\begin{aligned} \tfrac{\partial }{\partial t} {\text{ OC }_\text{ a }}&= {\mathcal{D }_{\text{ OC }_\text{ p }}}(\text{ RANKL })\,{\text{ OC }_\text{ p }}- {\mathcal{A }_{\text{ OC }_\text{ a }}}({\text{ TGF }\upbeta })\,{\text{ OC }_\text{ a }}\nonumber \\&\quad - \tfrac{\partial }{\partial x} \left( {\text{ OC }_\text{ a }}v_{\text{ OC }_\text{ a }}\right) ,\nonumber \\ \tfrac{\partial }{\partial t} {\text{ OB }_\text{ p }}&= {\mathcal{D }_{\text{ OB }_\text{ u }}}({\text{ TGF }\upbeta })\,{\text{ OB }_\text{ u }}- {\mathcal{D }_{\text{ OB }_\text{ p }}}({\text{ TGF }\upbeta })\,{\text{ OB }_\text{ p }},\nonumber \\ \tfrac{\partial }{\partial t}{{\text{ OB }}_{\text{ a }}}&= \mathcal D _{\text{ OB }_\text{ p }}({\text{ TGF }\upbeta }) {\text{ OB }_\text{ p }}- A_{{\text{ OB }}_{\text{ a }}}{{\text{ OB }}_{\text{ a }}}\nonumber \\&\quad + {k_\text{ form }}\, {{\text{ OB }}_{\text{ a }}}\left[ \!\tfrac{h_{{\text{ OB }}_{\text{ a }}}}{R}{{\text{ OB }}_{\text{ a }}}\!-\! \text{ OCY }_\text{ exp. }\left( x,R\right) \!\right] ,\nonumber \\ \tfrac{\partial }{\partial t} {\text{ TGF }\upbeta }&= n_{\text{ TGF }\upbeta }^\text{ bone }\,{k_\text{ res }}\, {\text{ OC }_\text{ a }}- D_{\text{ TGF }\upbeta }\,{\text{ TGF }\upbeta }, \end{aligned}$$
(23)

where

$$\begin{aligned} {\mathcal{D }_{\text{ OC }_\text{ p }}}(\text{ RANKL })&= D_{\text{ OC }_\text{ p }}\pi ^\text{ act }\left( \tfrac{\text{ RANKL }}{k^\text{ RANKL }_{\text{ OC }_\text{ p }}}\right) ,\nonumber \\ {\mathcal{A }_{\text{ OC }_\text{ a }}}({\text{ TGF }\upbeta })&= A_{\text{ OC }_\text{ a }}\pi ^\text{ act }\left( \tfrac{{\text{ TGF }\upbeta }}{k^{\text{ TGF }\upbeta }_{\text{ OC }_\text{ a }}}\right) ,\nonumber \\ {\mathcal{D }_{\text{ OB }_\text{ u }}}({\text{ TGF }\upbeta })&= D_{\text{ OB }_\text{ u }}\pi ^\text{ act }\left( \tfrac{{\text{ TGF }\upbeta }}{k^{\text{ TGF }\upbeta }_{\text{ OB }_\text{ u }}}\right) ,\nonumber \\ {\mathcal{D }_{\text{ OB }_\text{ p }}}({\text{ TGF }\upbeta })&= D_{\text{ OB }_\text{ p }}\pi ^\text{ rep }\left( \tfrac{{\text{ TGF }\upbeta }}{k^{\text{ TGF }\upbeta }_{\text{ OB }_\text{ p }}}\right) , \end{aligned}$$
(24)

and

$$\begin{aligned} \text{ PTH }(x,t)&= \beta _\text{ PTH }/ D_\text{ PTH }, \nonumber \\ \text{ RANK }(x,t)&= N^\text{ RANK }_{\text{ OC }_\text{ p }}\, {\text{ OC }_\text{ p }}, \nonumber \\ \text{ OPG }(x,t)&= \frac{\beta ^\text{ OPG }_{{\text{ OB }}_{\text{ a }}}\, {{\text{ OB }}_{\text{ a }}}\, \pi ^\text{ rep }\! \left( \frac{\text{ PTH }}{k^\text{ PTH }_{\text{ OB },\text{ rep }}}\right) }{\beta ^\text{ OPG }_{{\text{ OB }}_{\text{ a }}}\,{{\text{ OB }}_{\text{ a }}}\,\pi ^\text{ rep }\! \left( \frac{\text{ PTH }}{k^\text{ PTH }_{\text{ OB },\text{ rep }}}\right) /\text{ OPG }_\text{ sat } + D_\text{ OPG }} \nonumber \\ \text{ RANKL }(x,t)&= \frac{\beta ^\text{ RANKL }_{\text{ OB }_\text{ p }}\,{\text{ OB }_\text{ p }}}{1+k^\text{ RANKL }_\text{ RANK }\,\text{ RANK }+ k^\text{ RANKL }_\text{ OPG }\,\text{ OPG }} \nonumber \\&\quad \times \left\{ D_\text{ RANKL }+ \frac{\beta ^\text{ RANKL }_{\text{ OB }_\text{ p }}\,{\text{ OB }_\text{ p }}}{N^\text{ RANKL }_{\text{ OB }_\text{ p }}{\text{ OB }_\text{ p }}\ \pi ^\text{ act }\left( \frac{\text{ PTH }}{k^\text{ PTH }_{\text{ OB },\text{ act }}}\right) }\right\} ^{-1}. \nonumber \\ \end{aligned}$$
(25)

The density distributions \({\text{ OB }_\text{ u }}(x,t)\) and \({\text{ OC }_\text{ p }}(x,t)\) are assumed to be given functions (see Table 3). The dimensionless functions \(\pi ^\text{ act }\) and \(\pi ^\text{ rep }\) express the activation and repression of cell differentiation, apoptosis, or ligand expression by regulatory signalling molecules. These functions are determined by the fraction of receptors on the cell occupied by the concerned signalling molecules [see Pivonka et al. (2008); Buenzli et al. (2011)]. Mathematically, they are so-called Hill functions given by:

$$\begin{aligned} \pi ^\text{ act }(\xi ) = \frac{\xi }{1+\xi }, \qquad \pi ^\text{ rep }(\xi ) = 1-\pi ^\text{ act }(\xi ) = \frac{1}{1+\xi }. \end{aligned}$$
(26)

A slight change in the expression for \(\text{ RANKL }\) in Eq. (25) has been made compared to Buenzli et al. (2011). The production of \(\text{ RANKL }\) is now correctly proportional to the number of cells that express \(\text{ RANKL }\), i.e., we have replaced \(\beta _\text{ RANKL }\) in Ref. [Buenzli et al. 2011, Eq. (19)] by \(\beta _{\text{ OB }_\text{ p }}^\text{ RANKL }\,{\text{ OB }_\text{ p }}\) in Eq. (25). The same inconsistency of having a production rate of \(\text{ RANKL }\) not scaled by the number of osteoblasts is present in the temporal model of Pivonka et al. (2008), which was corrected in Buenzli et al. (2012b). The behaviour of the \(\text{ BMU }\) model is not changed significantly by this correction. Some inconsistent behaviours of the temporal model (Pivonka et al. 2008) were corrected by this change, see Buenzli et al. (2012b) for more details.

Table 3 \(\text{ BMU }\) model parameters and given functions

BMU model parameters and functions. All the parameters used in the model are listed in Table 3. Several parameters were taken over from the temporal model of bone remodelling of Pivonka et al. (2008); Buenzli et al. (2012b) and our previous model of the \(\text{ BMU }\) (Buenzli et al. 2011) (see last column of Table 3). We emphasise that rate parameters for cell differentiation/apoptosis are weighted by Hill functions of the concentration of signalling molecules, Eq. (24), which vary themselves in space and time. These rate parameters were therefore calibrated to known \(\text{ BMU }\) characteristics that they determine and should not be viewed as the actual value of the differentiation/apoptosis rates. For example, \(A_{\text{ OC }_\text{ a }}\) was chosen such that the \({\text{ OC }_\text{ a }}\) population reaches a steady state within 10 days (corresponding to the lifespan of osteoclast nuclei in active osteoclasts Jaworski et al. 1981), but \(1/A_{\text{ OC }_\text{ a }}\ne 10~\text{ days }\).

Rights and permissions

Reprints and permissions

About this article

Cite this article

Buenzli, P.R., Pivonka, P. & Smith, D.W. Bone refilling in cortical basic multicellular units: insights into tetracycline double labelling from a computational model. Biomech Model Mechanobiol 13, 185–203 (2014). https://doi.org/10.1007/s10237-013-0495-y

Download citation

  • Received:

  • Accepted:

  • Published:

  • Issue Date:

  • DOI: https://doi.org/10.1007/s10237-013-0495-y

Keywords

Navigation