Abstract
Progressive stiffening of collagen tissue by bioapatite mineral is important physiologically, but the details of this stiffening are uncertain. Unresolved questions about the details of the accommodation of bioapatite within and upon collagen's hierarchical structure have posed a central hurdle, but recent microscopy data resolve several major questions. These data suggest how collagen accommodates bioapatite at the lowest relevant hierarchical level (collagen fibrils), and suggest several possibilities for the progressive accommodation of bioapatite at higher hierarchical length scales (fibres and tissue). We developed approximations for the stiffening of collagen across spatial hierarchies based upon these data, and connected models across hierarchies levels to estimate mineralizationdependent tissuelevel mechanics. In the five possible sequences of mineralization studied, percolation of the bioapatite phase proved to be an important determinant of the degree of stiffening by bioapatite. The models were applied to study one important instance of partially mineralized tissue, which occurs at the attachment of tendon to bone. All sequences of mineralization considered reproduced experimental observations of a region of tissue between tendon and bone that is more compliant than either tendon or bone, but the size and nature of this region depended strongly upon the sequence of mineralization. These models and observations have implications for engineered tissue scaffolds at the attachment of tendon to bone, bone development and graded biomimetic attachment of dissimilar hierarchical materials in general.
1. Introduction
Progressive stiffening of collagen tissue by bioapatite mineral (‘bioapatite’) underlies much of musculoskeletal function and development. The spatial and temporal accommodation of bioapatite by collagen determines this stiffening and thus has important consequences for bone and for partially mineralized tissues that serve mechanical roles. Data for accommodation of bioapatite at the lowest levels of collagen hierarchical structure are now becoming available, although debate persists. Our objective, in this paper, is to apply this new nanostructural information, in conjunction with some recent nanoscale simulations of collagen mechanics, to develop estimates of the mechanical behaviour of partially mineralized collagen tissues. We then apply these models to interpret how partially mineralized tissues provide effective connection of tendon to bone.
Temporally, bioapatite accumulates on a collagen template during development to form mature bone [1]. Spatially, bioapatite increases from 0 vol% in tendon to 50 vol% in bone across the tendontobone insertion site (the ‘insertion’) [2,3]. During development, bioapatite crystals may nucleate and grow in channels that arise in ‘gap regions’ at the terminations of tropocollagen molecules within collagen fibrils before accumulating on extrafibrillar surfaces, and possibly within ‘overlap regions’ (figure 1). However, studies have also demonstrated extrafibrillar mineral in the absence of intrafibrillar mineral [2]. Mineralization patterns of developing tissues may differ considerably from those of mature tissues.
Our focus is prediction of how intrafibrillar and extrafibrillar mineralization stiffens collagen fibrils, fibres and tissues. We connected mechanical moduli across several hierarchical levels using homogenization methods that bound or estimate the degree to which bioapatite can stiffen tissue composed primarily of collagen and water. The coarsest of these methods involves mechanical responses of fictitious tissues obtained by arranging tissue constituents in parallel or in series. Increasingly accurate estimates and tighter bounds can be obtained by considering the arrangements of collagen, water and mineral.
A substantial body of literature focuses on estimates of bone biomechanics based upon inferred patterns of mineralization [4–7]. The science of relating the specific architecture of bone to its macroscale mechanical properties is quite mature, with solid tools and estimates in place for estimating anisotropic material properties [4,8], quantifying bioapatite–collagen interactions [9] and dissecting the key contributors to bone mechanics [10], dissecting failure mechanisms [11]. For the microscale features of bone, a host of technologies exist for relating patientspecific scans to individualized micro mechanical models [12–14].
This study adds to this literature in two ways. First, it quantifies mechanical implications of the most uptodate nanostructural data for collagen/bioapatite morphology. Second, it extends to partially mineralized tissues. Five possible sequences of mineralization were studied to evaluate the range of bioapatite accrual patterns supported by the literature, ranging from nucleation initiating in gap regions, as may occur during bone development, to nucleation initiating outside fibrils, as may occur along the insertion. We quantified how combinations of intrafibrillar and extrafibrillar bioapatite stiffens collagen, and how the specific sequence of mineralization affects tissuelevel mechanical responses. These tissuelevel responses were then studied in the context of attachment of dissimilar materials and engineered grafts to surgically reattach tendon to bone.
1.1. Background: nanoscale physiology of bone
Tendon and bone have nearly identical hierarchical collagen structures ranging from the nano to milliscale. Bone contains approximately 50 vol% bioapatite incorporated into and onto this structure. The insertion presents graded transitions in bioapatite, composition and meso, but not in nanoscale organization of collagen. It also presents gradients in protein content [15–18]. The tendinous end transitions towards bone from aligned, predominantly (approx. 98%) type I collagen fibres to a fibrocartilaginous layer composed of types I, II and III collagen and the proteoglycan aggrecan. This transitions to mineralized fibrocartilage dominated by type II collagen, but with significant amounts of type X collagen and aggrecan. Beyond this is bone proper, with mineralized type I collagen. Although the insertion is classically defined as containing four zones, no sharp boundaries are evident between zones, and the transitional tissue is compositionally graded.
As a firstorder simplification, we focus on type I collagen (‘collagen’), the main organic component of tendon, bone and the insertion. Collagen's nanoscale structure involves triplehelix collagen molecules, 300 × 1.5 nm diameter, that selfassemble into fibrils in a wellaligned but twisted and quasitriclinic molecular packing [19]. Within fibrils, collagen molecules stagger with an axial period of D = 67 nm [20] that includes an overlap region (approx. 27 nm, or 0.4 D), and a gap region (approx. 40 nm, or 0.6 D) where onefifth of triplehelix collagen molecules terminate (figure 1). Intermolecular spaces between triplehelix collagen molecules contain water and a small amount of noncollagenous proteins, including proteoglycans [21]. Collagen fibrils in tendons have a typical diameter of a few hundred nanometres and are decorated with proteoglycans, which form a matrix between fibrils. Fibrils assemble into fibres, fibres into fascicles and fascicles into tendon. We calculated mechanical properties at three levels. Fibrillevel moduli (superscript fibril) represent moduli of fibrils and their internal components, and were based upon data from molecular models. Fibrelevel moduli (superscript fibre) represent moduli of repeating unit cells incorporating fibres and materials including decorins that bind and crosslink fibrils. Tissuelevel properties (no superscript) were derived by averaging over a measured collagen fibre orientation distribution. A complete list of nomenclature can be found in the electronic supplementary material.
Accommodation of bioapatite within collagen is only partially understood [22,23]. A beautiful model system for partially mineralized collagen is avian tendon [24], but unfortunately the mineralization patterns observed there seem to differ from those seen in mammalian bone [25]. In mature bone, bioapatite crystals appear predominately as thin plates with approximate dimensions 2.1 × 30 × 40 nm [26–29]. The 2.1 nm minimum dimension of its Bravais lattice is the most certain of these dimensions [25]. The L_{m} ≈ 30 nm caxes of bioapatite crystals within gap regions align with the tropocollagen molecule axes. Atomic force microscopy studies [30,31] report bioapatite crystals with widths and lengths ranging from 30 to upwards of 200 nm. Recent work suggests that intrafibrillar bioapatite exists predominantly within endtoend gaps between triplehelix collagen molecules, in five defined orientations; three distinct area fractions (α = {0.20, 0.28, 0.58}) of bioapatite are possible within a transverse cross section of a gap region filled to capacity with bioapatite (figure 1).
Volume fractions of bioapatite are central to estimates of stiffening. α = 0.58 corresponds to a tissuelevel volume fraction of bioapatite ϕ_{m} ≈ 21%, based on the relationship , where is the fibrillevel volume fraction of bioapatite, and the area fraction of fibrils in mature tendon [32]. The precise amount of bioapatite in the intrafibrillar spaces of the overlap regions has not been established, but is bounded at 0.6 that of the gap regions; even with this maximum addition (ϕ_{m} ≈ 0.21(1 + 0.6) = 0.33), intrafibrillar bioapatite cannot account for the volume fraction of bioapatite present in fully mineralized bone. Consistent with this, extensive bioapatite is observed exterior to collagen fibrils [3]. The maximum volume fraction of bioapatite that can be accommodated by bone is therefore ϕ_{m} ≈ 0.41, if bioapatite cannot accrue in the overlap region, and ϕ_{m} ≈ 0.53 if it can. Both lie within the range reported for wet bone [4]. We explored the progressive stiffening of collagen by bioapatite within these constraints.
2. Material and methods
We modelled stiffening of collagen by bioapatite within gap regions, on the exterior of fibrils, and possibly within overlap regions. Our focus was prediction and bounding of the ways that bioapatite stiffens collagen. The stiffening was sensitive to the nanoscale structures and interactions of collagen and bioapatite. Although models exist for the structures of fully mineralized and nonmineralized collagen [3], the sequence of bioapatite accumulation during development and the bioapatite distributions within partially mineralized tissues at the insertion are not known [22]. We therefore studied the range of possibilities described below. The nanoscale mechanical interactions between bioapatite and collagen are not known, but are an area of focus by us and others. The interactions likely involve strong adhesion at low stress levels with little effect on tropocollagen mechanics, and sliding at higher stress levels [33]. In the absence of other information and as a first approximation, we model complete adhesion between collagen and bioapatite.
2.1. Models of the sequence of mineralization
Five plausible sequences of mineralization were modelled (figure 2). Models began with unmineralized collagen fibrils (top row, figure 2), followed by prescribed bioapatite accumulation into gap regions, onto the exterior of collagen fibrils and within overlap regions:

— model A (‘gapnucleated’) began with filling of gap regions (row 2, figure 2), and proceeded with extrafibrillar mineralization that initiated at the mineralized gap regions (row 3), then extended the entire length of the fibril (row 4). The first stage of mineralization (0 ≤ ϕ_{m} ≤ 0.21) involved inserting 2.1 nm thick and 30 nm high bioapatite platelets into the 0.4 D (40 nm) spaces between the Cterminus of one triplehelix tropocollagen molecule and the Nterminus of the next. Platelets were assumed to contact the Nterminus of one molecule and to extend 10 nm short of the Cterminus of the next. Bounds and estimates on stiffening by these platelets involved different spatial sequences of filling gaps, ranging from filling the maximum allowable space of one gap region before proceeding to the next (lower stiffness bound), to filling all gaps simultaneously with equal volumes of bioapatite (upper stiffness bound). The second stage (0.21 ≤ ϕ_{m} ≤ 0.41) involved formation of an extrafibrillar sheath from the Nterminus end of each gap region, past the Cterminus end, and then on to the next gap region. Model A is supported by observations of long bioapatite platelets spanning multiple collagen fibrils [3];

— model B (‘nucleationinhibited’) began as in model A with filling of gap regions. Subsequently, nucleation of extrafibrillar bioapatite was treated as rare. Following nucleation, growth of extrafibrillar bioapatite along a nucleation front was promoted. This was modelled by a sheath that began at one point along the fibril then extended to encompass the entire fibril;

— model C (‘nucleationpromoted’) also followed model A for intrafibrillar mineralization. Thereafter, nucleation of small, randomly distributed patches of extrafibrillar bioapatite was promoted, but subsequent growth of these patches was suppressed;

— model D was analogous to model C, except that extrafibrillar bioapatite accrued prior to intrafibrillar bioapatite; and

— model E was analogous to model A, except that bioapatite was also allowed into the intermolecular spaces of the overlap regions. The sequence of bioapatite accrual in model E was (i) as platelets within the gap regions (0 ≤ ϕ_{m} ≤ 0.21); (ii) in an unknown structure within the overlap regions (0.21 ≤ ϕ_{m} ≤ 0.33, rows 2–3); and (iii) as an extrafibrillar ‘gapnucleated’ sheath (0.33 ≤ ϕ_{m} ≤ 0.53, rows 3–5).
2.2. Bounds and estimates for the stiffening of collagen
Bounds and estimates of stiffening by bioapatite are summarized here and presented in complete detail in the electronic supplementary material.
2.2.1. Stiffening of fibrils
For intrafibrillar mineralization, fibrillevel longitudinal modulus bounds and estimates were assembled from structural arrangements of four types of regions: unmineralized collagen fibril within a gap region, modulus ; unmineralized collagen fibril within an overlap region, modulus ; mineralized collagen fibril within a gap region modulus ; and mineralized collagen fibril within an overlap region modulus .
and were estimated from ligament moduli reported by Stabile et al. [34], which, when scaled to the fibril level, yield E^{fibril} ≈ 600 MPa. This is consistent with moduli found for fibre recruitment models of tendon mechanics [35]. The overlap regions (height 0.4 D) and gap regions (height 0.6 D) were taken to resist load in series, with the assumption that , where β ≈ 0.8, because the number of collagen molecules in the gap regions is fourfifths that of overlap regions (figure 1). Then, the total deformation of a fibril subjected to uniaxial loading was , where σ is engineering stress (tensile force divided by initial crosssectional area). Thus 2.1
Bioapatite crystals (E_{ha} = 110 GPa) [36] incorporate into fibrils with the 30 nm caxis in parallel with collagen molecules. Bioapatite is believed to replace water without otherwise altering collagen [2]. The modulus of a mineralized gap region was then a parallel summation: 2.2where α(ϕ_{m}) = ϕ_{m}D/(ρL_{m}) denotes the degree of stiffening by bioapatite.
Stiffening of collagen in overlap regions was modelled with the Hashin–Shtrikman (HS) bounds (see electronic supplementary material, equations (2.1)–(2.3)). HS bounds are the tightest bounds for an unstructured mixture, and were appropriate, because intermolecular spaces within overlap regions cannot accommodate a single Bravais lattice of bioapatite. These bounds involved estimates for moduli and volume fractions of collagen molecules, bioapatite and the largely aqueous matrix between fibrils. Using a density of 1.12 g cm^{−3} for tendon [37], approximately 1.3–1.4 g cm^{−3} for collagen [38], and 1.0 g cm^{−3} for the matrix yielded a volume fraction of collagen in tendon of approximately 30–40%. Adjusting for the volume fraction of fibrils in fibres (approximately the area fraction ρ = 0.80) suggests that fibrils contain 50 vol% of water/matrix and collagen molecules. Using the body temperature bulk modulus of water and Poisson ratio v = 0.499999 yielded E = 13 kPa [39] for this effective matrix. The modulus of the remaining material within an overlap region was estimated from the threedimensional HS upper bound: matching this upper bound to and taking v = 0.3, this was then estimated as . Bioapatite accruing in the overlap regions replaced aqueous matrix and accrued according to , whereas the fibrillevel volume fraction of collagen remained at 50%; the volume fraction of aqueous matrix decreased as . The modulus of mineralized overlap regions was then taken as the threedimensional HS lower bound (bioapatite embedded within a fibril lattice surrounded by aqueous matrix), .
Fibril moduli were constructed from those of the above four regions. The upper bound involved adding equal volumes of bioapatite platelets to all gap regions simultaneously. Three types of regions were then combined in series: overlap regions, partially mineralized gap regions and unmineralized regions of gap regions representing space close to the Ctermini of tropocollagen molecules that never mineralize. The lower bound was obtained by requiring that a single gap region be filled to the maximum level (area fraction α_{max} = 58%) before gap regions in the next gap region begin to accrue bioapatite.
In model E, after bioapatite platelets filled in the maximum allowable volume in gap regions, bioapatite was deposited in overlap regions. The upper bound fibril modulus involved adding bioapatite to all overlap regions simultaneously, and the lower bound involved filling a single overlap region entirely before adding bioapatite to the next.
2.2.2. Longitudinal stiffening of fibres
Linear response of fibres. Adjacent collagen fibrils were modelled as bound together with an extrafibrillar matrix (EFM, volume fraction 1 − ρ = 0.20) containing a relatively small amount of noncollagenous proteins that serve to bind adjacent collagen fibrils [31]. We neglected EFM contributions to longitudinal stiffness because its estimated modulus (E_{EFM} ≈ 0.287 MPa, see electronic supplementary material) was small compared with that of collagen fibrils, but included it below in estimates of transverse moduli.
Bounds and estimates of stiffening involved accounting for extrafibrillar bioapatite in parallel with four types of possibly mineralized fibre segments: unmineralized gap regions, fully mineralized gap regions, unmineralized overlap regions (models A–D) and mineralized overlap regions (model E). Upper bounds on stiffening were estimated by modelling extrafibrillar bioapatite as thin sheaths around fibrils that grew radially increasing mineralization. Lower bounds were estimated by modelling extrafibrillar bioapatite as thick sheaths that grew longitudinally with increasing mineralization. Bounds and estimates for models A, B and E are included in the electronic supplementary material. Models C and D involved random accumulations of extrafibrillar bioapatite whose contributions were estimated using finiteelement methods, as described below. A key feature of each of these models was a percolation threshold representing a critical bioapatite volume fraction at which stiffening by apatite rose rapidly from near the lower bounds to near the upper bounds.
Nonlinear response of fibres. Although our primary interest was linear approximations to collagen behaviour, mechanical measurements show sea cucumber dermal collagen fibres [40–43] and associated molecular dynamics (MD) models for fibrils [33] show significant nonlinearity at large strains. Although these MD simulations have not yet been validated directly, we believe that they represent the best estimates for the mechanics of type I mammalian collagen fibrils. We therefore studied the effects of mineralizing these fibres as a first effort at exploring the effects of nonlinearity given available data. As with the linear estimates, we note that the details of collagen–mineral nanoscale interactions are unclear at present and therefore model as a first approximation perfect adhesion between the two. Validation and modification of these MD results and identification of models of nanoscale collagen–mineral interactions are ongoing areas of research in our laboratories.
We studied the uniaxial relationship between engineering stress σ and engineering strain ε (ratio of fibril distension to initial fibril length) in a partially mineralized mammalian type I collagen fibril that, when unmineralized, followed the stress–strain relationship predicted by the MD simulations of [44], . Mineralization was taken to follow the scheme of model A. Because data for are available for only a single loading rate, viscoelastic effects could not be modelled. An effective nonlinear, continuum constitutive response of hydrated collagen within a fibril, , was estimated from as follows. Collagen in the gap and overlap regions was idealized as identical continua. A uniaxial stress σ on the fibril was elevated to σ/0.8 in the gap regions owing to gaps between tropocollagen molecules. The total strain of a fibril was then the weighted sum of the strain in gap regions and in the overlap regions: 2.3and 2.4 was calculated from equation (2.4) and a fifth degree polynomial fit to data for (dasheddotted curve, figure 3a, blue in online version). The uniaxial constitutive response of partially mineralized fibres was estimated based upon g^{fibril}(ε) (red solid curve, figure 3a) and upon methods analogous to those described above. Details are provided in the electronic supplementary material.
2.2.3. Transverse stiffening of fibres and fibrils
Transverse (x_{2} and x_{3}direction, where the x_{1}direction parallels the fibre axis) properties of fibres and fibrils contribute to the axial moduli of tendon and bone. Fibres have a distribution of orientations, and estimating the axial moduli of tendon and bone involves integrating the full stiffness tensors of fibres over this distribution. Little data exist for rigorous estimates of these transverse properties, and even the simple estimates presented here required a number of assumptions. Fortunately, this ambiguity does not cloud interpretation of predicted moduli of partially mineralized tissues: these estimates are fully appropriate for tangent moduli of tendon and for tissue containing bioapatite volume fractions above the percolation threshold, and the transverse moduli of fibres have little contribution to the structural response of the insertion at bioapatite concentrations below the percolation threshold [45].
Transverse moduli for partially mineralized tissues were studied analytically using HS bounds, and using finiteelement procedures described in §2.2.4. Details may be found in the electronic supplementary material.
2.2.4. Stochastic finiteelement simulation for random bioapatite accumulation
Finiteelement analysis was used for Monte Carlo simulations of random accumulations of extrafibrillar bioapatite (models C–E) and for estimating transverse moduli. The domain studied (figure 4) was a rectangular section of a compliant EFM (red) containing a hexagonal array of fibrils (volume fraction 0.8). Fibrils had a uniform diameter of 155 nm [46]. Fibrils were divided longitudinally into six pairs of overlap and gap regions; gap regions were divided into regions that could accept mineral (height L_{m}) and regions that could not (height 0.6 D − L_{m}). In simulations that allowed addition of randomly accruing bioapatite, fibrils were surrounded by sheaths of matrix elements that could be replaced with bioapatite (mesh not shown).
Periodic unit cell boundary conditions were applied to an eighthspace model. The three surfaces with outward normals in the positive 1, 2 and 3directions were constrained to maintain these normal directions, and to be free of shear traction, whereas the remaining three surfaces were constrained to have no shear tractions and no displacement in the direction of the surface normal. In each instance, one of the three faces with a positive outward normal was perturbed a prescribed amount, whereas the other two such faces were allowed to shrink inwards owing to Poisson contraction. The effective linear longitudinal and transverse moduli of the unit cell were determined from the forces necessary to sustain these isometric levels of stretch.
The EFM was approximated as isotropic, and fibrils were approximated as containing up to four transversely isotropic regions: gap regions with maximal intrafibrillar mineralization, unmineralized gap regions, mineralized overlap regions and unmineralized overlap regions. Poisson's ratios were set equal for each type of region: . The transverse and longitudinal moduli of these regions were as described above. The shear modulus was set to . Finiteelement results showed effective moduli to be relatively insensitive to choices of , and .
Simulations began with estimates of the effects of replacing unmineralized gap and overlap regions with mineralized gap and overlap regions within fibrils. Subsequent effects of extrafibrillar mineralization involved replacement at random of elements having the properties of EFM with elements having the properties of bioapatite. These Monte Carlo simulations were performed using standard techniques using scripts written in the Matlab environment in combination with Comsol finiteelement analysis software. Careful convergence studies were performed. Models with no extrafibrillar mineralization could achieve convergence within 5% with as few as 20 000 quadratic, reduced integration tetrahedral elements (figure 4), whereas models with extrafibrillar mineralization required upwards of 100 000 elements.
2.3. Stiffening of tissues
The spatially varying tissuelevel stiffness tensor C(x) was estimated at points x along the axis of a insertion from fibrelevel moduli C^{fibre}(ϕ_{m}), published values for the spatial variation of ϕ_{m}(x), and published values of the spatial variation of the angular deviation s(x) that defined the circular orientation distribution of fibres according to a von Mises distribution [47]. C^{fibre} was assembled from spatially varying moduli estimates as described in the electronic supplementary material. Assuming a fibre distribution that is symmetric about the axis of the insertion, C(x) was calculated according to: 2.5where θ is the inclination angle from the tendontobone axis, R(θ) is a rotation from a coordinate system with the x_{1}direction oriented θ from the tendontobone axis to the global coordinate system in which C(x) is defined, and , in which b is the root of , I_{0}(b) and I_{1}(b) are the modified Bessel functions of order 0 and 1, respectively, and s is in radians.
3. Results and discussion
3.1. Progressive stiffening of collagen fibres by bioapatite
Bioapatite increased the stiffness of collagen monotonically, but the degree of stiffening was sensitive to the sequence and structure of bioapatite accumulation. The effective nonlinear, fibrillevel constitutive response g^{fibril}(ε) of a hydrated collagen continuum was estimated from equation (2.4) and the data of [44]. The inverse of this function, (g^{fibril})^{−1}(σ), followed these data closely (figure 3a, solid curve, red in online version). At the next hierarchical level, the fibrelevel tangent stiffness of collagen increased linearly with straining, but never approached that of a nonlinear, fully mineralized fibre (figure 3b).
Mineralization according to model A stiffened such fibres monotonically (figure 3b,c). We note once more that while we can be confident of the responses at lower strains, the upper limits of these curves are highly uncertain, because the strain levels associated with permanent changes to collagen are not yet known, and those associated with slippage between collagen and bioapatite are even less certain (cf. [33,48]). Stiffening was less pronounced beyond ϕ_{m} = 0.21, the point at which gap channels became filled to capacity with bioapatite and intrafibrillar mineralization completed (figure 3d). For 0.21 < ϕ_{m} < 0.3, extrafibrillar sheaths extended over gap regions containing intrafibrillar mineral, yielding little additional stiffening (figure 3c,d). Significant stiffening was observed once more for ϕ_{m} > 0.3, as extrafibrillar sheaths extended over the overlap regions.
Trend lines in figure 3d do not converge at the stiffness of bone, because the levels of strain that must be accommodated just prior to percolation in this model are also high for input data to provide trustworthy constitutive estimates; extrapolating out would yield stiffnesses of a mineralized fibril that far exceed the moduli of bone. The high levels of strain in the compliant phase associated with the approach to percolation are evident in the curves associated with 0.21 < ϕ_{m} < 0.3. Here, the curves switch from concave up to concave down, associated with the slight reduction in stiffness of in the last datapoint of the Gautieri data, at a strain level approximately twice the average strain in the fibril. These predictions require validation to be trusted, especially because that single datapoint might represent either yield or numerical error. However, even with this uncertainty they highlight the amplification of strain that occurs in the approach to percolation and the importance of understanding the role of structural effects on the material response of a composite tissue.
Despite these limitations, it is useful to examine the relevance of the predictions to physiological levels of strain. One clear limit is the strain level at which fibrils and fibres begin to disassemble when loaded. The possible range of strains for disassembly is large, with estimates ranging from 0.07 to 0.5 [33]. These seem to be high values compared with the normal physiological range (a few percentage strain in bone, and up to 10% in tendon), but relating these measures to strain levels in vivo is premature because of the many unknowns including the level of prestrain in fibrils and the heterogeneity in fibril strains that underlie the gradual recruitment of fibres during stretching. For load transfer across a region of tissue with a gradation in mineralization, a uniform tensile loading would be expected to induce a spectrum of strains in partially mineralized collagen fibrils based upon their degree of mineralization. Although further study is needed before nonlinear estimates such as these can be applied to model load transfer across a tendontobone insertion site, the results provide an important baseline datapoint for ongoing study.
We next directed attention to linearized moduli to assess the five models of mineralization. The HS bounds (shaded region, figure 5; yellow in online version) formed a broad envelope for stiffening by bioapatite, with the lower bound extending beneath the lower cutoff of the graph. The HS upper bound displayed a slight ‘bump’ at ϕ_{m} = 0.21, where the composite considered shifted from a fourphase (mineralized gap regions, unmineralized gap regions, unmineralized overlap regions and EFM) to a fivephase composite (including extrafibrillar bioapatite). Incorporating microstructure into the bounds led to tighter bounds (darker shaded region, figure 5 orange in online version). For ϕ_{m} < 0.21 in models A–C (figure 5a), the upper bound required adding equal amounts of bioapatite to all gap regions simultaneously, resulting in relatively rapid initial stiffening and little stiffening thereafter, whereas the lower bound required adding bioapatite one gap region at a time, resulting in more gradual stiffening.
For models A–C, subsequent bioapatite accrued on the fibril exterior. The upper bound involved uniform sheaths of bioapatite that increased in thickness with increasing ϕ_{m}. The lower bound required either a sheath of bioapatite growing from one end of an infinitely long fibre (model B, lightshaded region) or periodically repeating belts growing from mineralized gap regions (model A, darkshaded region, orange in online version). The upper bounds associated with models A and B differed, because the latter involved addition of extrafibrillar bioapatite to a homogenized, internally mineralized fibril, whereas the former incorporated the spatial distribution of bioapatite (equations in the electronic supplementary material).
For model E (figure 5c), an intermediate region existed in which bioapatite accrued within the overlap regions. Because bioapatite that exists in overlap regions is likely disorganized, threephase (bioapatite, fibril, water) HStype bounds were used.
The above bounds aided in the interpretation of the estimates of stiffening. For models A–C, the lower bound was the trendline for intrafibrillar mineralization: simulations with gap regions mineralized in random order (circles, figure 5a) followed the lower bound, as occurs with particle reinforced composites [49]. We note that random gapbygap mineralization is supported by the literature [22]. Thereafter, the three models diverged. Model A (extrafibrillar sheaths extending from gap regions) followed the series lower bound, with appreciable stiffening only after the sheaths extended beyond mineralized gap regions. The linear model reduced stiffening slightly at low ϕ_{m} relative to the nonlinear model. Model B (uniform sheath extended a single site) stiffened more rapidly as the sheath encompassed overlap regions at lower ϕ_{m}. Model C (random extrafibrillar additions) followed model B initially, but then stiffened rapidly towards the upper bound as a continuous network of mineralized EFM formed for ϕ_{m} beyond the percolation threshold. In all cases, the final modulus of a fibril exceeded the nominal modulus of E_{bone} = 20 GPa for cortical bone, indicating that the composite need not be fully dense with extrafibrillar bioapatite to reach this stiffness.
Model D (figure 5b, random extrafibrillar mineralization preceding intrafibrillar mineralization) followed models A–C until a percolation threshold, then stiffened beyond the fourphase HS bounds, as expected for a composite in which the stiff phase has a preferred alignment. Subsequent intrafibrillar mineralization yielded little stiffening.
Model E (figure 5c) included bioapatite in the overlap regions (0.21 ≤ ϕ_{m} ≤ 0.33). Stiffening here followed the HS lower bound, continuing the trend from the first region with little change in curvature, suggesting that stiffening is insensitive to whether gap versus overlap mineralize first. Thereafter, model E paralleled model A, except with a discontinuity in stiffening as the extrafibrillar sheaths extended over the different regions. Discontinuities disappeared for randomly accruing extrafibrillar bioapatite as in model C (circles, figure 5a).
The sequence of mineralization determined the nature of the mineralized structures that arose, and, in turn, the degree to which a specific volume fraction of bioapatite could stiffen a fibre. In each case, percolation was a key determinant of the degree of stiffening. Early formation of a sheath on a fibre's exterior (model D), as may occur in development, leads to significant stiffening at approximately onethird of the mineral volume fraction required if intrafibrillar mineralization precedes extrafibrillar mineralization (models A–C).
A number of other sequences of mineralization are possible in addition to those considered. Two important cases merit further discussion here, both relating to the possibility of bioapatite residing in the overlap regions, and uncertainty about its structure. First, how would a fibre behave if bioapatite accrued simultaneously and randomly in both gap and overlap regions? The likely answer is evident from data already presented: the results for model C indicate that the stiffening owing to random intrafibrillar accumulation would roughly follow that of the first two regimes of figure 5c, except that it would rise above that of the lower bounds to result in more rapid stiffening. Second, what if the mineral within the overlap regions is structured? A putative structure for bioapatite in the overlap regions has not yet been pinned down, but strong evidence exists that a connected structure is possible (cf. [29,50–52]). Structured bioapatite, especially if elongated along the axes of fibrils, would lead to a connected network and earlier percolation [53]. Following percolation, bioapatite would contribute almost as if it were added in parallel to the fibrils (cf. [54,55]) and, as observed from model D, would lead to a stiffening that would likely exceed the HS upper bound. Further characterization of the morphology of bipartite within the overlap regions is critically important.
Several important limitations must be borne in mind when interpreting these results. The first is the decomposition of fibrils into elastic matrix, bioapatite and hydrated collagen molecules that are not affected by surface interactions with bioapatite. This simplification would be expected to cause difficulty at very low levels of mineralization, but nevertheless enables a reasonable picture of how accrual of bioapatite in specified patterns stiffens collagen. The key is that uncertainty in the assumed properties and interactions are offset in the gap regions by only moderate levels of mineralization owing to the strong contrast in moduli: as gap regions fill with only a small fraction of bioapatite, the mechanics of the region becomes dominated by the bioapatite. The second is the nature of mineralization of the overlap regions. The structures and bioapatite/collagen interactions in this region are unknown. The third is validation: the absence of experimental data to validate the MD simulations upon which these estimates are based requires that nonlinear estimates be recognized as reliable only at smaller strain levels.
3.1.1. Stiffening transverse to a fibre axis
Transverse modulus estimates (figure 6a) were dominated by the EFM modulus at low levels of mineralization. Intrafibrillar mineralization had minimal (lower bound) stiffening effect, whereas extrafibrillar mineralization could lead to percolation more rapid than that observed for axial moduli. This was sharpest for model A, with a bioapatite network linking mineralized gaps of neighbouring fibres forming shortly after the onset of extrafibrillar mineralization. Formation of such bioapatite bridges is consistent with TEM observations [3]. Subsequent mineralization had a smaller stiffening effect. Model B predicted more gradual stiffening, because the percolated network surrounded both gap and overlap regions at the lowest levels of extrafibrillar mineralization. The random patterns of extrafibrillar bioapatite in model C resulted in further percolation delay. Model D (figure 6b) was analogous to model C but with extrafibrillar mineralization preceding intrafibrillar mineralization. The estimate for model E (figure 6c) was initially lower than the HS bounds, which are valid only for nonstructured solids. Significant stiffening occurred only following extrafibrillar mineralization. We note that these observations might explain AFM measurements reporting a sharp mineralization front at the ligamenttobone insertion: the models predict a sharp rise in transverse stiffness over a narrow range of the gradual rise in mineralization.
3.2. Progressive stiffening of collagenous tissues by bioapatite
We studied effects of gradations in bioapatite at the insertion. A compliant region between tendon and bone has been measured experimentally by tracking local strain between unmineralized and mineralized tissues [15]. This compliant tissue may serve to reduce stress concentrations [56] and is understood to arise from a tradeoff between competing gradients of collagen organization and bioapatite [45]. The details of this zone are of interest, because the natural tissue is not regenerated following healing or surgical repair, and postsurgical failure rates reach 94% [57]. An important factor is surely that type II collagen in the insertion is expected to have mechanical properties that differ from type I collagen in bone and tendon. We hypothesized that, in addition to this, the nanoscale architecture of bioapatite accommodation is important for the spatial variation of stiffness along the insertion, and therefore investigated how the stiffness of a collagen fibril increases as bioapatite accumulates within and upon it, and as bioapatite links neighbouring fibrils. Results showed that the size of the compliant region, underestimated by previous methods, is highly sensitive to the way that mineral accrues on fibres across the insertion site.
The dominant factor in the size of this region is the onset of mineralization: this occurs at about the peak of disorganization in collagen at the insertion, and about halfway along its length (figure 7). The specific pathway by which bioapatite accumulates on the interior and exterior of a fibril can extend the width of the compliant zone. For model D, in which an extrafibrillar sheath of bioapatite forms around fibrils prior to mineralization of the fibril's interior, the rise in stiffness owing to bioapatite is relatively rapid. In all four of the other models, which we consider to be more likely candidates based upon our observations, the pattern of bioapatite accrual delays the rapid rise to the stiffness of bone, in some cases by delaying percolation to nearly 80% of the way from tendon to bone (figure 8). The minimum value of modulus predicted between tendon and bone is higher than that reported in [15], suggesting that other factors not considered, including a spatial gradient change in the collagen content, might play an additional role.
Results suggest that no modification of collagen properties by bioapatite owing to nanoscale interactions are needed to explain the stiffness of bone in terms of the properties of its constituents. Whereas many engineering polymer matrix nanocomposites are dominated in their behaviour by a stiffening and embrittlement that occurs at the polymer/nanoreinforcement interface, this might not occur in bone.
These results suggest that the pattern of mineralization that exists between tendon and bone may be tailored to maximize the size of a compliant band of tissue in healthy adults. Following both healing and repair surgery, the functionally graded tissue region and associated compliant zone are not regenerated, perhaps an important factor in the high rates of retear that occur [57]. Results presented here suggest that graft design can be optimized to account for this. We are currently pursuing this in our own laboratories through nanofibrous networks with gradations in bioapatite content [58] and orientation distribution [59] that mimic those of the insertion. Such nanofibrous electrospun mats are biocompatible, and can be tailored with cues for directed cell growth [59]. Although a wide range of mechanical properties are possible in this system, this study suggests that a fundamental rethinking of this approach is warranted involving stiffening and toughening over the entire range of spatial hierarchies [60,61]. Indeed, MD analysis has shown that microstructural details are key to constitutive and strength response of proteins, including hydrogen bond and crosslinkdriven unfolding mechanisms [62], and dissipative sliding of molecules [63]. While sacrificial bonds of bioapatite between mineralized collagen fibres [31] and the heterogeneity of such links [48,64] have long been known to be central to bone's toughness, the current results suggest an important role for the hierarchically derived stiffening of collagen by bioapatite, and suggest this as an important future direction for the design of scaffolds that guide postsurgical healing.
4. Conclusion
The nanoscale details of the accumulation of bioapatite on partially mineralized collagen fibrils are central to the fibril mechanics. Although the sequence and structure of mineralization are not yet known, reasonable estimates can be made. We presented models for a reasonable range of estimates and showed how these are important in the case of attachment of dissimilar materials. Mineralization sequences cited most frequently in the literature suggest that the nanostructure of collagen is optimized to accumulate bioapatite in a way that delays percolation. This in turn expands the size of the compliant region of partially mineralized tissue between tendon and bone and provides for a slower spatial gradient in tissue modulus, which has important ramifications in terms of macroscale stress transfer, and for cues to cells [65,66]. Nanostructurally based mechanical differences between developing and mature tissues may underlie poor healing outcomes in torn tendontobone attachments [23,67]. In combination with observations that intermediate levels of mineralization are optimal for toughening of collagen [48], the current result that intermediate levels of mineralization might be distributed optimally across hierarchies to minimize stiffening motivate a picture of the insertion as an optimized, tough, compliant band.
Funding statement
This work was supported in part by the NSF (CAREER 844607), the NIH (R01 AR055580), and a joint NSFNIH multiscale modelling grant (U01EB016422).
Acknowledgements
Y.X.L. acknowledges a graduate fellowship from the Fannie Stevens Murphy Foundation.
 Received September 11, 2013.
 Accepted November 26, 2013.
 © 2013 The Author(s) Published by the Royal Society. All rights reserved.