Hierarchical structure of bone and anisotropy
This section briefly reviews the hierarchical structure of bone across the scales, and its consequences on the elastic anisotropy, mainly based on Weiner and Wagner (1998) and Rho et al. (1998). The bone structure is illustrated on Fig. 1.2.
The nanoscale is the length scale of the basic constituents of bone: min-eral crystals (mainly hydroxyapatite), collagen molecules (type I), non-collagenous proteins, and water. The mineral crystals form plate-like structures of about 50 25 1:5 4 nm called platelets. Collagen molecules are made of three long amino acids chains organized in a triple helix and are approximately 1 1:5 nm in diameter and 300 nm long. To form the so-called mineralized collagen fibrils, col-lagen molecules align along a given direction, and mineral crystal grow in between the collagen molecules with their crystallographic c-axis – the stiﬀest direction – aligned with the collagen molecules. Mineralized collagen fibrils then assemble into arrays or bundles called fibers. The mineral part provides stiﬀness and strength to the more compliant collagen part.
At the microscale, fibers arrange parallel to each other in layers of about 3 7 m thick called lamellae. In cortical bone, the lamellae are arranged concentrically around the Haversian canals ( 50 m) to form cylindrical layered structures called osteons ( 100 – 300 m in diameter). The Harvesian canals are aligned with the long axis of bone and contains a blood vessel and nerve fibers. Inside a lamella, the main direction of the fibers is tilted. Several patterns of tilting angle along the radius of the osteons have been described in the literature, with diﬀerent expected consequences on anisotropy (Reisinger et al., 2011a). However, Reisinger et al. (2011b) observed using nanoindentation that osteons are generally stiﬀer in the axial direction. The space between osteons is filed with the interstitial matrix, which is made of remnants of osteons.
At the mesoscale (millimeter scale), which is the intermediate scale between the highest level of microstructure (lamellae and osteons) and the whole organ, bone tissue can be seen as a two phase composite: long cylindrical soft fibers (Haversian canals) are aligned in a hard tissue matrix (osteons and interstitial tissue). Bone is stiﬀer in the axial direction at this length scale, due to a cumulative eﬀect of the preferred orientation of the Harversian porosity and the preferred direction of alignment of the mineralized collagen fibrils inside the matrix.
Measurement of cortical bone elasticity at the millimeter-scale
This section reviews state-of-the-art methods for the measurement of cortical bone elasticity on millimeter-sized specimens and discusses their respective advantages and shortcomings. These methods can be classified in two categories: 1) mechanical testing methods, in which elasticity is obtained from the stress-strain relationship in low frequency or quasi-static regime and 2) ultrasonic methods, based on the fact that, under some particular conditions, the velocity of a wave propagating in an elastic medium can be directly related to a stiﬀness coeﬃcient.
For both approaches, basic equations are recalled. An extensive discussion of the elasticity of anisotropic media, and particularly of wave propagation, can be found in textbooks, such as Auld (1990) and Royer and Dieulesaint (2000).
Resonant Ultrasound Spectroscopy
Summarizing the particularities of bone elasticity assessment and the limitations of the conventional methods described above, the three main requirements for an enhanced method are:
The method should be able to completely characterize the anisotropic elasticity of a material from a single specimen, without requiring 45 -oriented faces. This is linked to the spatial variability of bone elasticity.
It should be fully applicable on small specimens (of the order of a few mm3), since the thickness of the cortical shell of bones limits the specimen size.
It should be as accurate and reproducible as possible. Hence, it should not assume idealized state of stress and strain but rather take into account the complete stress-strain relationship. Moreover, delicate contact or bounding between the sample and the apparatus should be avoided.
Resonant Ultrasound spectroscopy (RUS) is a method that has been developed to satisfies these requirements, originally for the characterization of geological materials (see the historical review in section 1.4.2), and that could therefore be beneficial for the study of bone mechanics. The following section introduces the basic principles of RUS, followed by a brief historical review. Finally, the diﬃculties of the application of RUS to bone tissue samples are introduced.
The basic principle of RUS is to use the frequencies of mechanical resonance of a sample of material to infer its stiﬀness tensor. Indeed, the resonant frequencies of an elastic body depend only on its geometry, mass density and stiﬀness tensor. If the two first parameters are directly measured, elasticity can then be estimated by matching model-predicted frequencies to the experimental results (inverse problem approach). This requires a setup for the measurement of the resonant frequencies, a method to compute the vibration modes of a sample given its elasticity and finally a procedure to solve the inverse problem. These elements are described in a generic way in the next paragraphs, mostly based on the book by Migliori and Sarrao (1997).
RUS and attenuative materials such as bone
In this section, we discuss the resonant behavior of materials that are highly attenu-ative for mechanical waves at ultrasonic frequencies, resulting in important damping of the resonant modes. The few studies in which RUS was applied to highly attenu-ative materials, including bone, were in purpose excluded from the above historical review and will be considered here.
Many RUS applications, and all that have been cited in the historical review above, were done on materials in which the damping of the mechanical vibrations at ultrasonic frequencies is very low, such as crystals. The quality factor Q, usually introduced to characterize the level of damping, is defined as inverse of the relative width of a resonant peak at half amplitude Q = f = f (Lakes, 2009). For a metal or for crystal samples, Q is usually well above 1000 (Table 1.1 for some values), producing very sharp peaks in the spectrum at positions of the resonant frequencies (Fig. 1.4, right). In these cases, it is not a problem to find the resonant frequencies at the positions of the local maxima of the spectrum, and close peaks are well resolved (Fig. 1.6). As the damping in the material increases, the Q factor decreases, and the peaks broaden. Then, close peaks may overlap and the frequencies may not correspond anymore to maxima of the spectrum (Fig. 1.6). For bone, the Q factor is of the order of 20 (Table 1.1), which corresponds to the worst case in Fig. 1.6.
In geophysics, RUS has been applied to rock samples with a Q of about 150 by Ulrich et al. (2002). In this last study, authors investigated the side-length ratios of a rectangular parallelepiped that minimizes peak overlapping. In the context of bio-materials, Lee et al. (2002) applied RUS to cortical bone but could not measure all the elastic moduli, due to the low Q factor in the range 10-30. They were able to find the three first resonant modes and used them to estimate the three shear moduli (as the first resonant modes in RUS are generally mostly dependent on the shear moduli). Lebedev (2002) proposed the use of a signal processing method to retrieve resonant frequencies with Q = 50, based on the fact that the frequency response of a sample has a known functional form: it is a sum of Lorentzian line-shapes. The work of Lebedev, as a possible path to overcome the diﬃculties reported by Lee et al. in the application of RUS to bone, was the starting point of the present work.
Resonant Frequencies Calculation
RUS is based on a comparison of calculated and measured resonant frequencies. We used a well-documented method to calculate the frequencies of mechanical resonance for a solid of given elasticity and rectangular parallelepiped geometry (Migliori and Sarrao, 1997). Briefly, the resonant angular frequencies ! were found by searching the stationary points of the Lagrangian L = 2 ZV !2ui2 Cijkl @xj @xl dV; (2.1) 1 @ui @uk.
where V and are respectively the specimen’s volume and mass density, Cijkl are the stiﬀness constants and u is the displacement field. By expanding the displace-ment field in a set of polynomial functions (Rayleigh-Ritz method), the stationary Receiver equation @L = 0 was written as a generalized eigenvalue problem !2Ma = Ka; (2.2).
where M and K denote respectively the mass and stiﬀness matrices of the vibra-tion problem. Eq. (2.2) was solved numerically, giving the eigenvalues !2 and the eigenvectors a, which contain the coeﬃcients of the polynomial expansion.
Resonant Frequencies Measurement
Transducers were built from small piezoceramic elements (diameter of 3 mm and length of 2 mm) polarized in compression and suspended to copper thin films. This construction reduces the transmission of energy to the mechanical holder and then prevents any unwanted system resonance (Migliori and Sarrao, 1997). Support films were glued to steel cylinders connected to the ground, which provided isolation from electromagnetic perturbation and cross-talking between transducers. With these transducers the measurement system was non-resonant in the frequency range of interest.
The specimen was placed between the two transducers, held by two opposite corners to be as close as possible to stress-free boundaries conditions (Fig. 2.2). A Vectorial Network Analyzer (Bode 100, Omicron electronics GmbH, Klaus, Austria) was used to record the frequency response of the specimen between 100 kHz and 280 kHz with 50 Hz resolution. The receiver output signal was preconditioned before recording by a broadband charge amplifier (HQA-15M-10T, Femto Messtechnik GmbH, Berlin, Germany). Phase and magnitude of the frequency response were saved for further signal processing.
Preliminary measurements showed that the number of resonant frequencies which can be retrieved may vary from one measurement to another on a same speci-men. This is expected due to the diﬃculty to excite the vibrational modes which are associated to a motion of the specimen’s corner nearly parallel to the surface of the transducers. This diﬃculty is well known for low damping material (Stekel et al., 1992) and is even more critical in the case of high damping, where the signal is weak and the modal overlap is strong. Although in previous preliminary tests we ob-served that it was possible to determine the full stiﬀness tensor of a bone-mimicking material based on a single positioning of the specimen (Bernard et al., 2011), the quality of the measurement can be increased by combining a few measurements of the specimen. Indeed, intermediate repositioning introduces small variations of the transducer-specimen coupling which in practice increases the probability of exciting and detecting weakly excited modes. Accordingly, the measurements were repeated, removing and replacing the specimen six times. This limited the number of missed modes in the frequency range of interest. Repeated measurements also provided a way to estimate uncertainty on the resonant frequencies.
Uncertainty on the elastic constants
Following Migliori et al. (1993), we used the curvature of the cost function F around the minimum to estimate the uncertainty on each particular constant. The method consists in finding the largest possible changes on the constants which lead to an increase of 2% in F .
The uncertainty on the measured resonant frequencies also provide a way to es-timate uncertainty on the elastic constants by the mean of Monte-Carlo simulation. We generated 100 random sets of frequencies from independent normal distributions centered on the experimental mean values and using experimental standard devia-tions. This 100 sets were used to estimate the distribution of the elastic constants values. The mean and standard deviation on each constant was then calculated.
From the measured elastic constants the 6 6 stiﬀness matrix was constructed and numerically inverted to obtain the compliance matrix (Auld, 1990). From the latter, engineering moduli (Young’s moduli, shear moduli and Poisson’s ratios) were calcu-lated (Bower, 2009). This was repeated for the 100 sets of constants obtained from Monte-Carlo simulation, allowing to estimate the uncertainties (standard deviation) on the engineering moduli.
Ultrasonic Velocities Measurement
For comparison purpose, we repeated BWV measurements as described by (Granke et al., 2011). Briefly, the elastic constants on the diagonal of the stiﬀness tensor were deduced from velocities of longitudinal and shear ultrasonic waves propagating in the principal directions of the specimen. Velocities were obtained from the time-of-flight of an ultrasonic pulse propagating trough the specimen between a pair of transducers contacting two opposite faces. Central frequencies of 2.25 MHz and 1 MHz were used for longitudinal and shear wave velocities measurements respectively (V105RM and V152RM, Panametrics Inc., Waltham, MA). Signals were acquired using an oscilloscope (TDS 2012, Tektronix Inc., Beaverton, OR) and processed with MatLab.
Repeatability of the method, evaluated to be 3:2% and 4:7% for longitudinal and shear elastic constants respectively (Granke et al., 2011), was used to evaluate uncertainty on the elastic constants measured from BWV.
Table of contents :
1.1 Context and motivation
1.2 Hierarchical structure of bone and anisotropy
1.3 Assessment of cortical bone elasticity at the millimeter-scale
1.3.1 Mechanical testing
1.3.2 Ultrasonic waves methods
1.4 Resonant Ultrasound Spectroscopy
1.4.1 Basic principles
1.4.2 Historical review
1.4.3 RUS and attenuative materials such as bone
1.5 Outline of the thesis
2 RUS measurement of cortical bone elasticity: a feasibility study
2.2 Materials and Methods
2.2.2 Resonant Frequencies Calculation
2.2.3 Resonant Frequencies Measurement
2.2.4 Signal Processing
2.2.5 Elastic Constants Estimation
2.2.6 Uncertainty on the elastic constants
2.2.7 Engineering moduli
2.2.8 Ultrasonic Velocities Measurement
3 RUS for viscoelastic characterization of anisotropic attenuative solid materials
3.2 Resonant frequencies computation
3.2.1 Rectangular Parallelepiped
3.3 Samples and measurement setup
3.3.2 Setup for the measurement of the frequency responses
3.3.3 Additional elasticity measurements
3.4 Processing of the measured spectra
3.4.1 Estimation of the resonant frequencies in time domain
3.4.2 Non-linear fitting in frequency domain
3.5 Elastic constants estimation
3.5.1 Bayesian formulation of the RUS inverse problem
3.5.2 Automated pairing of the resonant modes
3.6 Damping factors estimation
4 Bayesian mode identification and stiffness estimation using Markov Chain Monte Carlo
4.2.1 Forward problem
4.2.2 Bayesian formulation of the inverse problem
4.2.3 Gibbs sampling
4.3 Application 1 – Data from Ogi et al., 2002
4.3.1 Prior distributions
4.3.2 Results and discussion
4.3.3 Additional results
4.4 Application 2 – Data from Bernard al., 2013
4.4.1 Prior distributions
4.4.2 Results and discussion
5 Application to a collection of human tibial cortical bone specimens
5.2 Materials and methods
5.2.2 Measurement setup and signal processing
5.2.3 Estimation of the elastic properties
5.4 Statistical analysis
Summary and conclusion
A Appendix: Transformed elastic parameters
A.1 Isotropic symmetry
A.2 Cubic symmetry
A.3 Transversely Isotropic symmetry
A.4 Orthotropic symmetry
B Appendix: Cutting protocol for human tibia specimens
B.2 Specimen preparation
B.3 Transversal cuts
B.4 Radial cuts
B.5 Axial cuts
B.6 C+2/3 specimen