Determination of the Mechanical Properties of Lumbar Porcine Vertebrae with 2D Digital Image Correlation

Abstract

Purpose

To evaluate the strain fields and to calculate the modulus of elasticity and Poisson's ratio of trabecular bone of the 6 lumbar vertebrae of the porcine spine by a 2-dimensional digital image correlation (2D DIC).

Methods

This study was performed through a 2D DIC technique and the specimens were tested under compression. The resulting images were analyzed numerically by 2D DIC. Then, representative regions of interest were examined. The strain fields were determined and stress-strain curves were obtained.

Results

The full field measurement of the strain in the lumbar bone spine was evaluated and with this data, the Young's modulus and Poisson's ratio were determined.

Conclusions

This research highlights the potential applications of noninvasive acquisition techniques in biomechanical analysis. This is useful in the mechanical characterization of bony structures and in the design of prostheses.

Keywords

Digital image correlation Trabecular bone Elastic properties Field measurements

Introduction

Different experimental techniques for exploring mechanical features exhibited by animals in vitro and in vivo have been developed. In certain cases, the intention is to understand the biomechanical behavior of the human spine (1-2). The knowledge of normal conditions of the human spine is an essential starting point for a better understanding of numerous pathologies and injuries. To achieve this goal, the anatomical analogy between quadrupeds and bipeds has been studied by different authors (3-4-5-6). There are also mechanical tests that compare composite materials and natural bone, which can solve orthopedic disorders (7). The relationships between the 2 specimens allow results to be extrapolated from 1 kind to another. For instance, the pork spine exhibits a large morphological and biomechanical similarity with the human spine (8-9-10). One of these similarities is the orientation on the arrangement of the trabecular morphology, which indicates according to Wolff's law (11) that the applied loads over the vertebral body are similar. This is in accordance with densitometry studies that validate this fact (12). Therefore, porcine models have been a good choice for predicting results related to human cases in several reported studies (13-14).

The mechanical behavior of trabecular bone primarily depends on 3 parameters: the structure formed by its cells, the solid volume fraction, and the material properties in the cell wall (15). It is a porous body that displays viscoelastic and anisotropic behavior, and these features complicate its study, making it advisable to assume some simplifications. To achieve this objective, there are diverse mechanical tests to determine the mechanical properties of trabecular bone (16-17-18). Some experimental works give punctual results but they cannot describe the behavior of the whole body. Because most of the numerical analysis of the spine visualizes it as a continuous body and frequently as a homogeneous media, many researchers have chosen to obtain the properties at a macro level. It is reported in the literature that there are strong variations in the mechanical properties of porcine and human specimens (19-20-21-22-23). Nevertheless, most of the authors agree regarding the inherent properties of the specimens; moreover, experimental procedures greatly influence the fluctuations in these results. Characteristics such as race, age, weight, size, morphology, possible diseases, injuries, time of post-mortem testing, and storage in different media before experimentation could have an effect on the resulting evaluations.

On the other hand, the 2D DIC technique has proven to be an efficient method for calculating deformations on flat specimens. Multiple applications have been developed in experimental mechanics fields for this purpose. Scientific research associated with mechanical structures and material characterization has been performed using this technique because the displacement and strain fields can be analyzed without any contact with the specimen. With these objectives in mind, this work focuses on the identification of the mechanical properties for describing some global and local aspects of the entire surface of a sample. It was decided to use 2D DIC to determine the Young's modulus and Poisson's ratio. Of note, his noncontact optical method seems to be an ideal candidate to describe deformation in the porcine trabecular bone of the vertebral column.

Theoretical basis of 2D DIC

2D DIC can be defined as the search for a point p(x, y) in an undeformed image, from a point $p' (x', y')$ in a deformed image. This optical method compares images to obtain deformation measurements without any contact with the body. This advantage eliminates invasiveness from the experiment and it requires low cost equipment in comparison with interferometric techniques like electronic speckle pattern interferometry (ESPI) or digital holographic interferometry (DHI).

Since the discretized value of p (x, y) can be repeated consistently in the distorted image, a small square subset with center in p (x, y) and dimensions of $(2 M + 1) * (2 M + 1)$ is chosen to find a solution. The reference subset is known as f(x, y) and the deformed subset is named $g (x_{i}^{'}, y_{i}^{'})$ . In order to evaluate the similarity between the 2 subsets, a criterion or correlation functions are used. For this latter purpose, the zero-normalized sum of squared differences (ZNSSD) was considered. It is defined as

C = \sum_{i = - M}^{M} \sum_{j = - M}^{M} [\frac{f (x_{i}, y_{j}) - f_{m}}{Δ f} - \frac{g (x_{i}^{'}, y_{j}^{'}) - g_{m}}{Δ g}], Eq. [1]

Eq. [1]

where $f_{m} = \frac{1}{{(2 M + 1)}^{2}} \sum_{i = - M}^{M} \sum_{j = - M}^{M} f (x_{i}, y_{i})$ , $Δ f = \sqrt{\sum_{i = - M}^{M} \sum_{j = - M}^{M} {[f (x_{i}, y_{i}) - f_{m}]}^{2}}$ , and so on. When C ≈ 0, it can be said that f(x_i, y_i) ≈ $g (x_{i}^{'}, y_{i}^{'})$ , and the difference between p(x, y) and $p^{'} (x^{'}, y^{'})$ will define the displacement vector in the deformed images. The problem is then reduced to finding the minimum of equation 1 for obtaining the components of displacement u and v. In the literature, a large number of possibilities for determining the minimum of this function are shown and an extensive review of these methods can be found in (24). For this work, the Newton Raphson method was chosen and implemented in a MATLAB software, as proposed by (25) and improved by (26), to determine the deformation in the vertebral bodies analyzed. The calculation time to reach a solution is shorter than other correlation algorithms. Moreover, it can be considered that the deformation of the subset $g (x_{i}^{'}, y_{i}^{'})$ allows the measurement of larger deformations.

Assuming that the deformed body behaves under the principles of a continuum solid, it is possible to say that the relationship between x to x′ and y to y′ is

\begin{array}{l} x^{'} = x + ξ_{d} \\ y^{'} = y + η_{d}, Eq. [2] \end{array}

Eq. [2]

where ξ_d and η_d are polynomials of degree d. Since complex deformations are likely to be found in trabecular bone, it is convenient to use second-order polynomials. In this way, it is possible to describe a more convoluted pattern of deformations (27). When d = 2, equation 2 becomes

\begin{array}{r} \begin{matrix} x^{'} = x + u + u_{x} Δ x + u_{y} Δ y + 1 / 2 u_{x x} Δ x^{2} + 1 / 2 u_{y y} Δ y^{2} + u_{x y} Δ x Δ y \\ y^{'} = y + v + v_{x} Δ x + v_{y} Δ y + 1 / 2 v_{x x} Δ x^{2} + 1 / 2 v_{y y} Δ y^{2} + v_{x y} Δ x Δ y, \end{matrix} \\ Eq. [3] \end{array}

Eq. [3]

where, $u_{x} = \partial u / \partial x$ , $u_{y} = \partial u / \partial y$ and $Δ x = x_{i} - x_{o}$ , etc. These values are the argument of the correlation function (1), and is defined by

P = {u, v, u_{x}, v_{x}, u_{x x}, v_{x x}, u_{y y}, v_{y y}, u_{x y}, v_{x y}} \cdot Eq. [4]

Eq. [4]

In a next step, a P that satisfies C(P) ≈ 0 is sought. The iteration process of equation 5 is performed until the limit of convergence is reached and an optimal P has been found.

P^{k + 1} = P^{k} - \nabla \nabla C {(P^{k})}^{- 1} \nabla C (P^{k}) . Eq. [5]

Eq. [5]

To reach a subpixel resolution, it is necessary to reconstruct the deformed images from pixel integer values. In this work, bicubic splines were used for the interpolation of the deformed images through equation 6:

I (x, y) = \sum_{i = 0}^{3} \sum_{j = 0}^{3} a_{i j} x^{i} y^{i}, Eq. [6]

Eq. [6]

To ensure the performance of the programmed routines in MATLAB, a comparative experimental test was carried out in order to calculate the deformation of a well-known sample. This validates the proposed algorithm that is used in calibration of the acquired data. A tension test was conducted with a universal Instron-8500 machine (Norwood, MA, USA). The specimens were manufactured with 304L stainless steel and an axial clip on an extensometer of 50 mm to 75 mm, adhering to a resolution of 0.0001 mm. This transducer acquired the deformation data. The test was executed with a crosshead velocity of 2.5 mm/min. Images were taken with an 8-bit Sony CCD camera (model GN258DNF; Sony, Tokyo, Japan), with continuous diffusive illumination provided by an LED source installed at 1.25 m from the grips of the machine. A Computar Zoom lens with varifocal length 18 mm to 108 mm was employed (Computar, Commack, NY, USA). The images were detected at 1 Hz and the resulting curve obtained by 2D DIC and the extensometer are compared in Figure 1A; additionally the resultant strain field is plotted in Figure 1B.

Figure 1A shows that the different techniques obtain results that are in perfect agreement, giving an average strain error of 0.0016, a standard deviation (SD) of 0.0034, and an average relative error of 3.78%. The stress-strain curves of porcine specimens were obtained with these confidence limits (CL). In addition, the Young's modulus E and the compressive Poisson's ratio were determined.

Fig. 1

Strain in a 304 L stainless steel target. (A) Comparison between 2D DIC and extensometer data; (B) strain field in the specimen subjected to 31 kN.

Material and methods

Three spines of different pigs were chosen. Their race was Pietrain and they weighed around 95 kg. They were approximately 7 months old. After the animals were sacrificed, they were refrigerated for 24 h at a temperature below 2°C. Having undertaken this process, 18 prismatic specimens of 12 x 12 x 15 mm were taken from the L1 to L6 vertebrae. It is necessary to chose prismatic specimens to implement 2D DIC, because the surface of the object to study must be flat and parallel to the camera sensor. One side of the specimen was painted with a layer of white paint. Then, random points of black paint were sprayed to attain a unique pattern to correlate the surface. At the end of this preparation, the specimens were tested in a universal machine (MTS 858 Table Top System; MTS, Eden Prairie, MN, USA), with a capacity of 15 kN. The crosshead speed was 2 mm/min, and an approximately 30-N preload was applied. The period between obtaining the specimens and testing them was not greater than 8 h. To record the deformation of the compression test, the optical setup shown in Figure 2 was employed. A diffuse, white light illumination was used together with an 8-bit CCD camera, equipped with an 18 mm-to-108 mm zoom lens and a resolution of 754 x 480 pixels, at a distance of 60 cm from the specimen. This distance was optimized to minimize the errors introduced by out-of-plane deformation. The sampling frequency of the load, the crosshead distance, and the image acquisition were set at 2 Hz.

Fig. 2

Experimental setup.

Determination of mechanical properties

Different regions of interest (ROI) were chosen in order to determine the stress-strain curves. Figure 3 shows a typical ROI that was selected. Within this area, measurement points were generated at an interval of 5 pixels to find the displacement fields u, v. The size of the square subset to evaluate the correlation coefficient was 31 pixels per side.

Fig. 3

ROI of the compression test.

It has been previously indicated that changes in the intensity and poor interpolation algorithms (28) as well as the out-of-plane displacement (29), the convergence of the method, or the heating of the electronic systems introduce errors or noise into the displacement field as shown in Figure 4A. A simple but effective way to reduce this noise is proposed by Pan et al (30): the ROI is approximated by small linear plane subsection functions with equation 7. Figure 4B illustrates the smoothed displacement field within a square window size of 17 elements per side.

\begin{matrix} u_{s} (i, j) = a_{0} + a_{1} x + a_{2} y \\ v_{s} (i, j) = b_{0} + b_{1} x + b_{2} y \end{matrix}, Eq. [7]

Eq. [7]

Fig. 4

Displacement u in a lumbar specimen. (A) Raw displacement; (B) smooth displacement.

Another important advantage of this method is that the displacement field is derived without the need for subsequent numerical derivation. This process is advantageous because it avoids the numerical derivation that has been considered an inconvenient, noisy process in the calculated displacement field (30).

Strains were calculated with the Lagrangian formulation, because the analysis starts with the undeformed body and it progresses to the deformed state, leaving the tensor strain as

ε = 1 / 2 [F^{T} F - 1], Eq. [8]

Eq. [8]

where F stands for the deformation gradient tensor and it can be expressed in terms of the derivatives of displacement as

F = [\begin{matrix} u_{x} + 1 & u_{y} \\ v_{x} & v_{y} + 1 \end{matrix}] . Eq. [9]

Eq. [9]

To compare the engineering strain obtained by our MTS system and the strain field measured by the 2D DIC method, it was applied the equation (10) which allows a direct comparison between both methods,

ε_{i n g} = \sqrt{1 + 2 ε_{y y}} - 1. Eq. [10]

Eq. [10]

By using the average strain of the ROI (Fig. 3) the stress-strain curves of the 18 specimens were obtained. Within the elastic region, the points were fitted to a straight line by the method of ordinary least squares (OLS). To calculate the longitudinal Poisson's ratio, the same procedure was followed, by averaging the strain transverse to longitudinal strain.

Results and discussion

Eighteen strain-stress curves were calculated for the porcine specimens. Figure 5 shows a typical case in which these curves show the strain obtained by 2D DIC, the linear approximation for the OLS, and its counterpart delivered by the MTS-858. The difference by the method of 2D DIC and strain was estimated by the MTS system. In this case a higher divergence than in the case of steel was observed; however, it is still within reasonable parameters. The strain average error in this analysis was 0.0018 with a standard deviation of 0.0122 and a relative error of 16.36%.

Fig. 5

Typical strain-stress curve from trabecular compression test.

Poisson's ratio was calculated with the longitudinal ε_yy (Fig. 6A) and transversal ε_xx (Fig. 6B) strains. The evaluations were performed within the elastic region to obtain the average Poisson's ratio. This procedure was repeated in all cases to calculate the average of E, υ and the standard deviation among all. The results can be seen in Table I.

Fig. 6

Transversal and longitudinal strain for a load of 810 N. (A) ε_xxstrain; (B) ε_yystrain.

Table I

Average modulus E, ν, and standard deviation

Porcine vertebrae	L1	L2	L3	L4	L5	L6
E [MPa]	234.67	212.29	273.50	276.51	183.119	322.119
ν	0.294	0.389	0.326	0.251	0.177	0.402
STD	112.35	85.54	146.45	107.19	65.09	219.45

Important differences in the resulting parameters related to mechanical properties of porcine trabecular bone have been reported (19-20-21-22). Since experimental techniques used for their determination are diverse, different agents could be responsible for these significant variations. For instance, some systems for measuring strain require being in direct physical contact with the system, and it includes an additional parameter that must be taken into account regarding the physical properties of the surroundings for each experiment. Additionally, from punctual measurements, it is contemplated that local results could be extended in order to give information about the whole sample. The 2D DIC method shows an adequate convergence for the calculation of Young's modulus, not only for particular measurements but for all the selected areas that were studied. The consulted literature reports Young's modulus values within a range from 229 MPa to 1080 MPa, which is in good agreement with our reported average value of about 252 MPa. The variation in our reported study samples during the tests can be easily explained by inherent features in the specimen, such as age, pathology, nutrition, its conservation and storage, among others.

From the biomechanical point of view, an important topic is the identification of the mechanical parameters of trabecular bone in vertebral bodies. This would be useful in further research and would help to understand the behavior of the spine under different types of load. From the obtained results, it can be observed that at similar levels of strain in the trabecular bone, a large heterogeneity is exhibited in comparison with an isotropic material such as steel. When trabecular bone is considered as homogeneous in a numerical analysis, reasonable approximations to reality are obtained. It is now possible to develop finite element models with a great level of geometrical detail. However, its weak point is that unlike interferometric techniques, the exact mechanical properties of the material cannot be introduced in small regions (31). Therefore, usually it is assumed that the mechanical properties are homogeneous (32). This situation is more complex when ceramic materials used as rehabilitation materials, for instance, are introduced into the required numerical model (33).

Conclusions

Results of a sensitive, noncontact technique for a highly accurate determination of biomechanical properties are presented in this work. An important topic in biomechanics is the identification of the mechanical parameters of trabecular bone for further investigation into the behavior of the spine under different types of load. In this work it is demonstrated that the 2D DIC method can be used to calculate the Young's modulus and the Poisson's ratio for the trabecular bone of the spine column without being in contact with the body. It was also shown that this method can provide the behavior of the entire specimen surface in compression. From the presented results it can be observed that at similar levels of strain in the trabecular bone a large heterogeneity is exhibited in comparison to an almost ideal material such as stainless steel. However, this resulting feature is not too far from the ideal model, indicating that numerical analysis considered as homogeneous trabecular bone can provide good approximations of reality.

Footnotes

Financial support: The authors would like to acknowledge the financial support they received from the National Council for Science and Technology, the National Polytechnic Institute, the Center for Research in Optics, the Polytechnic University of the Valley of Mexico, and the Institute of Science and Technology of the Federal District.

Conflict of interest: None.

References

Alini

Eisenstein

Ito

et al.

Are animal models useful for studying human disc disorders/degeneration?

Eur Spine J 20081712–19

Jiang

Moreau

Raso

Russell

Bagnall

A comparison of spinal ligaments—differences between bipeds and quadrupeds.

J Anat 1995187Pt 185–91

Cotterill

Kostuik

D'Angelo

Fernie

Maki

An anatomical comparison of the human and bovine thoracolumbar spine.

J Orthop Res 198643298–303

Kumar

Kukreti

Ishaque

Mulholland

Anatomy of deer spine and its comparison to the human spine.

Anat Rec 20002602189–203

Boszczyk

Putz

Comparative and functional anatomy of the mammalian lumbar spine.

Anat Rec 20012642157–168

McLain

F.R.

Yerby

and Moseley

Comparative Morphometry of L4 Vertebrae: Comparison of Large Animal Models for the Human Lumbar Spine.

Spine (Phila Pa 1976) 2002278E200–E206

Fuentes

Sáenz

Igartua

Martinetti

Dolcini

Barandika

Structural characterization and mechanical performance of calcium phosphate scaffolds and natural bones: a comparative study.

J Appl Biomater Biomech 201083159–165

Bozkus

Crawford

Chamberlain

et al.

Comparative anatomy of the porcine and human thoracic spines with reference to thoracoscopic surgical techniques.

Surg Endosc 200512191652–1665

Busscher

Ploegmakers

Verkerke

Veldhuizen

Comparative anatomical dimensions of the complete human and porcine spine.

Eur Spine J 20101971104–1114

10.

Yingling

Callaghan

McGill

The Porcine Cervical Spine as a Model of the Human Lumbar Spine: An Anatomical, Geometric, and Functional Comparison.

J Spinal Disord 1999125415–423

11.

Goldstein

Matthews

Kuhn

Hollister

Trabecular bone remodeling: an experimental model.

J Biomech 199124Suppl 1135–150

12.

Smit

The use of a quadruped as an in vivo model for the study of the spine - biomechanical considerations.

Eur Spine J 2002112137–144

13.

Schmidt

Richter

Claes

Puhl

Wilke

Limitations of the cervical porcine spine in evaluating spinal implants in comparison with human cervical spinal segments: a biomechanical in vitro comparison of porcine and human cervical spine specimens with different instrumentation techniques.

Spine (Phila Pa 1976) 2005 30 11 1275

14.

Tai

C-L

Hsieh

P-H

Chen

W-P

Chen

L-H

Chen

W-J

Lai

P-L

Biomechanical comparison of lumbar spine instability between laminectomy and bilateral laminotomy for spinal stenosis syndrome - an experimental study in porcine model.

BMC Musculoskelet Disord 2008 9 1 84

15.

Ramírez-Díaz

Ortiz-Prado

Schouwenaars-Franssens

Ruiz-Cervantes

Metodología para el modelado del comportamiento mecánico de hueso esponjoso a partir de sus microestructuras.

Ingeniería, investigación y tecnología 2010112199–216 http://www.redalyc.org/articulo.oa?id=40412515007. Accessed July 22, 2015

16.

Borah

Dufresne

Cockman

et al.

Evaluation of Changes in Trabecular Bone Architecture and Mechanical Properties of Minipig Vertebrae by Three-Dimensional Magnetic Resonance Microimaging and Finite Element Modeling.

J Bone Miner Res 20001591786–1797

17.

Goulet

Goldstein

Ciarelli

Kuhn

Brown

Feldkamp

The relationship between the structural and orthogonal compressive properties of trabecular bone.

J Biomech 1994274375–389

18.

Jaasma

Bayraktar

Niebur

Keaveny

Biomechanical effects of intraspecimen variations in tissue modulus for trabecular bone.

J Biomech 2002352237–246

19.

Lin

R-M

Tsai

K-H

Chang

G-L

Distribution and regional strength of trabecular bone in the porcine lumbar spine.

Clin Biomech (Bristol, Avon) 1997125331–336

20.

Mosekilde

Kragstrup

Richards

Compressive strength, ash weight, and volume of vertebral trabecular bone in experimental fluorosis in pigs.

Calcif Tissue Int 1987406318–322

21.

Ryan

Pandit

Apatsidis

Stress distribution in the intervertebral disc correlates with strength distribution in subdiscal trabecular bone in the porcine lumbar spine.

Clin Biomech (Bristol, Avon) 2008237859–869

22.

Teo

JCM

Si-Hoe

Keh

JEL

Teoh

Relationship between CT intensity, micro-architecture and mechanical properties of porcine vertebral cancellous bone.

Clin Biomech (Bristol, Avon) 2006213235–244

23.

Di Mascio

Bellini

Galbusera

Raimondi

Brayda-Bruno

Assietti

Lumbar total disc replacement: a numerical study.

J Appl Biomater Biomech 20108297–101

24.

Pan

Qian

Xie

Asundi

Two-dimensional digital image correlation for in-plane displacement and strain measurement: a review.

Meas Sci Technol 2009 20 6 062001

25.

Bruck

McNeill

Sutton

Peters

III

Digital image correlation using Newton-Raphson method of partial differential correction.

Exp Mech 1989293261–267

26.

Vendroux

Knauss

Submicron deformation field measurements: Part 2.

Improved digital image correlation. Exp Mech 199838286–92

27.

Cary

Deformation measurements by digital image correlation: Implementation of a second-order displacement gradient.

Exp Mech 2000404393–400

28.

Schreier

Braasch

Sutton

Systematic errors in digital image correlation caused by intensity interpolation.

Opt Eng 2000 39 11 2915

29.

Sutton

Yan

Tiwari

Schreier

Orteu

The effect of out-of-plane motion on 2D and 3D digital image correlation measurements.

Opt Lasers Eng 20084610746–757

30.

Pan

Asundi

Xie

Gao

Digital image correlation using iterative least squares and pointwise least squares for displacement field and strain field measurements.

Opt Lasers Eng 2009477-8865–874

31.

Morales-Bonilla

Torres-Torres

Urriolagoitia-Sosa

Hernández-Gómez

Urriolagoitia-Calderón

Interferometric laser detection of nanomechanical perturbations in biological media under ablation conditions.

J Phys Conf Ser 2011305012124–012127

32.

Beltrán-Fernández

Hernández-Gómez

Urriolagoitia-Calderón

González-Rebatú

Urriolagoitia-Sosa

Biomechanics and Numerical Evaluation of Cervical Porcine Models Considering Compressive Loads Using 2-D Classic Computer Tomography CT, 3-D Scanner and 3-D Computed Tomography.

Appl Mech Mater 201024-25287–295 http://www.scientific.net/AMM.24-25.287. Accessed July 22, 2015

33.

Hernández-Gómez

Beltrán-Fernández

Urriolagoitia-Calderón

González-Rebatú

Galán Vera

Urriolagoitia-Sosa

Biomechanical Characterization of a Cervical Corporectomy Using Porcine Specimens, Following an Experimental Approach.

Key Eng Mater 2011478103–111