A critical study on the experimental determination of stiffness and viscosity of the human triceps surae by free vibration methods

Depósito de Investigación de la Universidad de Sevilla https://idus.us.es/ This is an Accepted Manuscript of an article published by Sage in Proceedings of the Institution of Mechanical Engineers, Part H: Journal of Engineering in Medicine, Vol. 227, on September 2013, available at: https://doi.org/10.1177/09544119134878 1 A critical study on the experimental determination of stiffness and viscosity of the human triceps surae by free vibration methods F. París-García1, A. Barroso2, J. Cañas2, J. Ribas3, F. París2 1 Faculty of Sport Science, University Pablo de Olavide 2 School of Engineering, University of Seville 3 Faculty of Medicine, University of Seville [email protected], [email protected], [email protected], j[email protected], [email protected] Abstract: Muscles and tendons play an important role in human performance. Their mechanical behaviour can be described by analytical/numerical models including springs and dampers. Free vibration techniques are a widely used approach to the in vivo determination of stiffness and viscosity of muscle-tendon complexes (MTC) involved in sport movements. By looking at data reported in literature it appears that visco-elastic properties of the triceps surae MTC are independent of the modality in which free vibration is induced as well as they do not depend on the composition of the population of subjects submitted to the experiments. This research will critically discuss this important aspect focussing in particular on two studies documented in literature. For that purpose, two equipments will be developed to reproduce literature experiments under the assumption that the oscillating part of the body behaves as a single degree of freedom system: the governing degree of freedom is associated to the vertical displacement of the lower leg or to the rotation of the foot around the ankle articulation. Unlike literature, measurements are now conducted on the same population of subjects in order to draw more general conclusions on the real equivalence of results and validity of the mechanical properties determined experimentally. Free vibration tests are accurately simulated by analytical models describing the response of each vibrating system. It is found that if the two measurement protocols are applied to the same population of individuals as it is done in this study, values of visco-elastic properties of MTC extracted from experimental data are significantly different, the differences presenting a convincing consistency. This result is in contrast with literature and confirms the need to evaluate results of free vibration techniques by taking homogeneous bases of comparison. Keywords: Visco-elastic properties of human triceps surae, in vivo measurements, free vibration techniques, data fitting. 1.- Introduction. Muscular and tendinous fibres play an important role in improving human performance. Since mechanical behaviour of muscles and tendons resembles a 2 visco-elastic spring, their response can be described by simple elastic models including spring and damper elements (see, for example, [1]). Stretch-shortening cycles (SSC) appear in many movements performed by the human body. The muscle-tendon complex (MTC) can work in three different ways [2,3]: (i) excentrically (the muscle opposes elongation); (ii) isometrically (there is neither stretching nor shortening of the muscle); (iii) concentrically (there is a shortening of the muscle). The global mechanical response of MTC depends essentially on two factors: (i) transformation of chemical energy into mechanical energy (i.e. muscular shortening); (ii) efficiency of MTC elastic components in terms of their capacity to store elastic energy at the macroscopic level in the tendons and at the microscopic level in the muscles by means of titin, which is the elastic protein incorporated in the sarcomere that is in turn the smallest functional part of a muscle [4]. In order to establish a direct relationship between the elastic parameters of an MTC (typically stiffness and damping) and its corresponding efficiency in a particular sport discipline, MTC must play an important role in the development of the movements entailed by the discipline in question. For example, the triceps surae (TS), formed by the gastrocnemius (the muscle that generates movement), the soleus (the postural muscle) and the Achilles tendon, is the MTC that permits plantar flexion, a basic phase in the biomechanics of walking [4]. Therefore, TS is of fundamental importance in all sport disciplines involving general displacements. Experimental determination of MTC viscoelastic characteristics is widely documented in literature. There are two possible approaches to solve this problem: (i) protocols based on the use of cadavers such as, for example, the quick-release method [5], the -method [6], and the null-point method [7]; (ii) protocols based 3 on in vivo measurements such as, for example, the quick release method [8,9], application of ultrasounds [10], and free vibration techniques [1,11-14]. The basic idea behind free-vibration techniques that appear to be the approach most commonly used in literature is to induce damped free vibrations in the region of investigation by perturbing through impact another part of the human body. Some characteristic parameter describing the free vibration movement (typically a force) is experimentally measured and then related by a theoretical/analytical/numerical model to visco-elastic properties of tendons and muscles involved in the movement. If free vibration techniques are utilized, it is however necessary to check that experimental measurements are in good agreement with theoretical predictions provided by models within the limits of validity of each assumption made in the modelling stage. From a careful survey of literature it might appear that visco-elastic properties of the human triceps surae MTC be independent of the modality in which free vibration occurs as well as they do not depend on the composition of the population of subjects submitted to experiments. For example, Fukashiro et al. [12] and Babic and Lenarcic [13] considered different vibration movements: the displacement of the lower leg [12] and the rotation of the foot around the ankle articulation [13], respectively. Although Refs. [12,13] followed the same general principle (i.e. free vibration of the MTC under investigation), the movements monitored in the experimental measurements and, consequently, the associated dynamic equations, are completely different and involve different auxiliary muscular resorts. Furthermore, experiments were carried out on unrelated sets of people. However, visco-elastic properties were found to be quite similar. In view of this, the present study will critically review the determination of triceps surae visco-elastic parameters documented in [12,13]. The present 4 investigation is motivated by the fact that although the procedures developed in Refs. [12,13] derive from strong theoretical foundations, the similarity between the results obtained does not seem to have a physical explanation. Therefore, it may be very useful to assess in much deeper detail the real significance of the measured values. The main objective of this paper hence is to assess the true generality and representativity of each procedure in itself as well as to critically compare results. For that purpose, measurements described in [12,13] are reproduced in this research by implementing two experimental set-ups able to carry out the tests. Unlike literature, experimental protocols are carried out on the same population of individuals to avoid further uncertainties in the interpretation of the results. It is found that visco-elastic properties of the human triceps surae MTC extracted from experimental data now are significantly different. The differences found present a convincing consistency. This result is in contrast with literature data and confirms the need to evaluate results of free vibration techniques by considering homogeneous bases of comparison. The manuscript is structured as follows. Section 2 describes the freevibration experimental set-ups, the corresponding theoretical models to simulate dynamic behaviour and the least-square approach used to extract triceps surae visco-elastic properties. Section 3 presents and discusses experimental results in terms of stiffness and viscosity. Finally, Section 4 summarizes the main findings of the study. 5 2.- Determination of visco-elastic properties of the human triceps surae MTC In this research, visco-elastic properties of the human triceps surae were determined by conducting in vivo measurements. For that purpose, two experimental setups implementing the free vibration approach were specifically designed and realized. As mentioned previously, the rationale behind free vibration techniques is to originate a free damped vibration of an isolated part of the body. In the present case, the vibration response was controlled by the viscoelastic properties of the TS (Achilles tendon + soleus + gastrocnemius) complex which is assumed to behave as a 1-degree of freedom dynamic system. The global reaction force developed in the oscillation period was monitored by a load cell and the force values measured over the vibration period were compared with their counterpart predicted by a theoretical model reproducing the 1-degree of freedom dynamic system assumed. In this way, it is possible to build an error function depending on the unknown visco-elastic properties of TS. Finally, stiffness and damping constants can be determined via least-square fitting of experimental data. The experimental equipments, the corresponding theoretical models, and the least-square fitting procedure are described in detail in the rest of this section. 2.1 TS identification procedure based on the vertical displacement of the lower leg 6 The experimental apparatus was designed following the measurement approach described in Ref. [12]. The basic idea is to generate a configuration of the human body by allowing in the lower leg a vertical displacement, which will be used as the reference equation of motion. This vertical displacement will be accompanied by rotations of the foot and of the upper leg. The overview schematic of the testing apparatus is shown in Fig. 1. Since the frontal part of the foot is simply supported by a plate, the vertical displacement considered in the measurement involves the rotation of the foot around the articulation of the ankle. The testing apparatus includes the following main parts (see nomenclature in Fig. 1): 1. Load cell (Interface 1500ASK-300) capturing the evolution of the reaction force as vibration starts. The load cell data acquisition frequency of 1 KHz granted a satisfactory sampling of experimental data. The load capacity of 1500 N is adequate for the typical force values to be measured in the experiments. 2. The impact on the concentrated mass Mw causes the vertical movement of the lower leg as a 1-degree of freedom system. It should be noted that the total oscillating mass M is actually equal to the sum of Mw and the one degree of freedom equivalent mass of the lower limb. 3,4. Elements to accommodate the position of the person currently tested. Both the horizontal and vertical positions can be accommodated in several different locations, which are properly numbered. This allows the chair to be maintained in the same position, for the same person, if the test needs to be repeated. For any morphology of the subject, the three right angles indicated in Fig. 1 must be granted in order to ensure the repeatability of measurements. To this end, the person is asked to maintain the back in contact with the backseat described in the next paragraph. 7 5. The presence of the vertical backseat ensures, additionally to the former function, that the pressure on the load cell exercised by the foot be constant during measurements. The tested person must keep the back and the shoulders in contact with the backseat in order not to induce pressure variations that might be originated even by very simple movements of the subject's body during the test. Since the type of movement expected in the experiment involves the rotation of the thigh around the hip articulation in the plane of Fig. 1, part of the thigh comes in contact with the seat (of 31 cm depth) as movements develop. However, this did not significantly interfere with the above-mentioned rotation of the thigh because the displacement values observed at the knee level were of the order of mm, involving rotations of the thigh of an order of miliradians (≈0.1 degrees). This small angle, together with the size of the seat (only a very short part of the thigh is in contact with the seat, see Fig.2), make the upper leg to behave freely during the vibration. Furthermore, the thigh-seat contact zone recedes when the leg must hold the equilibrium position just before impact. 6. The simply supported beam allows the loading system to be easily adapted to the anthropometrical features of the subject currently tested. 7. The supporting frame included a triangular structure in order to be enough rigid in the plane of Fig. 1 thus not altering the free vibration behaviour of the system. Figure 2 shows the testing apparatus and a person ready for the test. It can be seen that the arms of the subject were held in neutral position in order to avoid any interference with the measurement. The setup can be utilized to test indifferently both legs (right and left) without introducing any variation either in the position or in the interpretation of experimental data. 8 When an impulse is applied to the concentrated mass Mw, the total mass M (i.e. the one-degree of freedom equivalent mass of the lower limb) starts to oscillate vertically as it is indicated in Fig. 1. Figure 3a shows an arbitrary intermediate position for this movement. Whilst distances L (length of the lower leg), R (forefoot distance from the articulation of the ankle) and r (rearfoot distance from the articulation of the ankle) remain constant during the vibration process the length associated to the triceps surae changes. In this study, the TS complex was schematized according to the classical Hill's model [15] by means of a spring of stiffness k and a damper of constant c. It is important to remark that the parallel elastic component (PEC) present in Hill’s model has not been considered in the present study (following [12] and the references therein) as is not likely that the PEC significantly contributes to the induced oscillation movement. Since for the position represented in Fig. 1 (i.e. 90° angle at the knee level) the gastrocnemius muscle does not work, the apparent values of stiffness k and damping c of the TS complex are relative to the soleus muscle and the Achilles tendon. For clarity of representation, the equivalent TS spring-damper system sketched in the figure is inclined by a certain angle with respect to the direction corresponding to the lower leg. However, in the derivation of theoretical model, TS was always considered parallel to the lower leg. The configurations shown in Fig. 3 are referred to a fixed point corresponding to the contact region between the subject's foot and the load cell measuring the reaction force F. The mass represented in Fig.3 is the total mass M involved in the movement, which is the sum of the concentrated mass Mw and the one degree of freedom equivalent mass representing in the model the distributed mass of the 15 3-4. Elements to accommodate the position of the person currently tested. For any morphology of the subject, the three right angles indicated in Fig. 6 must be granted and the rotation of the foot around the ankle articulation must be the only movement allowed. This ensured the validity of measurements. 5. The system transferring the load from the weight to the foot was comprised of a steel wire and two pulleys. The wire was very stiff with respect to tensile forces and flexible enough with respect to bending in order to easily adapt to the local shape of pulley channel. After several trials, the cross section of the wire was selected as a thin rectangular steel lamina of width 17mm and thickness 0.5 mm. The position of the pulley placed at the level of the foot (i.e. the pulley closest to the foot in the direction of transferred load) was tuneable in height thus allowing the final part of the wire to be held horizontal near the subject's foot. The disposition presented in Section 2.1 can be used in a natural form for both legs. However, the disposition associated to present procedure, requires either to conceive a unique very complicated device to test both legs, or to design a simpler one, but applicable only to one leg. The second option was selected and a testing apparatus for the right leg, Fig.7, was designed and built. To try to use the apparatus of Fig. 7 to test also left leg would imply the subject to take a forced posture where the right leg would be placed unnaturally in the right side of the structure, what can significantly affect measurements. In order to derive the theoretical model describing the dynamic behaviour of the system under free vibration the reference scheme of Fig. 8 was considered. The distance R between the origin of the reference system and the line of action of the equivalent force acting on the total mass M, the distance RG between the origin of the reference system and the centre of gravity of the foot, and the distance r 16 between the origin of the reference system and action line of the Achilles tendon can be expressed as follows: ff GGG FF rr rR rR    sin sin sin    (18) Under the hypothesis of infinite stiffness of the wire, the condition of momentum equilibrium written with respect to the equilibrium position (i.e. after applying the mass Mw) can be stated as: frMRRmI Gf    )( 22 (19) where: I has been previously defined,  is the angle that defines the position of the foot during the vibration movement and f is the reaction force developed in the MTC. The reaction force f developed in the MTC is the sum of the contribution of a spring of stiffness k and a damper of damping coefficient c. It can be written: rcrkuckuf           (20) where the relationship u=  r between  and u is derived from Fig. 8. By substituting the last equality of Eq. (20) in the RHS of Eq. (19) and rearranging terms, the dynamic equilibrium equation describing the movement of the system can be obtained: 0)( 2222   krcrMRRmI Gf  (21) Similarly to the previous loading condition, the terms relative to the inertia of the foot are negligible with respect to MR2. Therefore, Eq. (21) reduces to: 0 222   krcrMR  (22) The simplifying assumption leading to Eq. (22) makes it unnecessary to measure the position of the centre of gravity of the foot. 17 The general solution of the differential equation (22) is: )cossin()( tBtAet DD t     (23) hence analogous to the general displacement solution (6) derived in Section 2.1. The values of  and  D are respectively defined as: 2 2 2 MR cr   (24a) 2 2 2 2   MR kr D (24b) By combining Eqs. (24), the visco-elastic parameters k and c of the TS complex can be obtained: )( 22 2 2   D r MR k (25a)  2 2 2 r MR c (25b) In order to determine the values of k and c, one has to determine  and  D first by measuring the force developed in the wire as the impulsive load is applied. This force, which is registered by the load cell, can be expressed as: )()( RgMtFm     (26) The assumption that MsM reasoned in Section 2.1 will also be employed in this Section 2.2, as done in (26) with reference to the term Msg. If the displacement solution of Eq. (23) is derived twice with respect to time the following expression is obtained: )cossin()( tBtAet DaDa t      (27) Where Aa and Ba can be expressed in terms of A, B,  and  D as:   )(2 22 DDF ABA   and   )(2 22 DDF BAB   . 18 By substituting Eq. (27) in the RHS of Eq. (26), the final expression of the force to be measured is: MgtBtAetFF DFDF t m )cossin()(   (28) Where AF and BF can be expressed in terms of M, R, A, B,  and  D as:   )(2 22 DDF ABMRA   and   )(2 22 DDF BAMRB   . Expression (28) includes the unknown parameters  ,  D, AF, BF and M that must be determined with the least-squares fitting procedure described in Section 2.4. Equations (25) must hence be used to evaluate MTC visco-elastic parameters k and c. Similarly to the loading case previously considered, the values determined for parameters k and c may be strongly sensitive to the ratio R/r. Therefore, once forefoot and rearfoot distances R and r were determined for a subject, special cares must be taken in placing the foot always in the correct position (i.e. the second metatarsal head must be in contact with the transmission device). This position must be the same for both loading cases in order to have a homogenous basis of comparison for the visco-elastic parameters extracted from the two in vivo procedures considered in the study. 2.3. Forefoot and rearfoot distances measurement A noninvasive procedure based on pedigraphies, to obtain the lever arms of the forefoot (R) and rearfoot (r), has been developed. The procedure has two phases: the first is to obtain a footprint, by means of a pedigraphy, indicating on it a certain number of anatomical references. In a second phase, a procedure for locating the necessary 3 points of reference is followed. The three points, 19 indicated as a, b and c in Fig. 9 are placed along the longitudinal axis of the foot and their location allows us to obtain the desired lever arms. Point a is located at the position of the resultant of the distribution of pressure in the forefoot. Point b represents the intersection point between the projection of the ankle rotation axis and the longitudinal axis of the foot. Finally, c represents the projection of the Achilles tendon axis on the footprint plane. In this way, the positions a, b and c determine the distances R and r, as indicated in Fig. 9. A more detailed description of the procedure can be found in [20]. 2.4.- Processing of experimental data and extraction of TS viscoelastic parameters Reaction force values were recorded by the load cell at regular time intervals. A typical force-time plot is shown in Fig. 10 where the steady-state value finally reached by F is denoted as Fss. The plot in Fig. 10 is relative to the measurement protocol described in Section 2.1 and was obtained for a nominal weight of 25 kg. The nominal value of applied load was recovered after the vibration period originated by the impulsive load. The steady-state force value indicated in the figure is about 35 kg, hence 10 kg higher than the nominal weight of 25 kg. This is because of the presence of the weight of the lower leg involved in the vibration process. Notice that all values of the vertical axis are negative due to the fact that the load cell is in compression. The equivalent force measured by the load cell for the two experimental procedures described in Sections 2.1 and 2.2 always has the same expression (i.e. Eqs. (12) and (28) for, respectively, the first and the second loading case). The instantaneous force value always depends on five parameters: the damping coefficient  , the vibration frequency of the damped system  D, the total mass M 20 involved in the movement (after the identification of Ms and M), and the coefficients AF and BF. In order to determine the unknown parameters listed above, one can minimize the error functional  corresponding to the difference between the values of force measured by the cell and their counterpart predicted by numerical models:      n i FFD MBAFtFΩ 1 2 expexp ),,,,()(  (29) where: n is the number of control points (i.e. experimentally recorded data) to be included in the fitting procedure; Fexp is the force value experimentally measured at the sampling time texp; F is the force value predicted by the theoretical model based on the assumption of a 1-degree of freedom system. The error functional  was minimized in the least-square sense. Generally speaking, the fitting procedure should always converge to the global minimum. However, there are some aspects to be considered. In the first place, the probability of reaching the global optimum may depend significantly on the ranges of variation chosen for the five unknown parameters involved in fitting procedure. Furthermore, optimization behaviour may be sensitive to the initial values set for the fitting parameters that should always be chosen in view of their physical meaning or on the basis of previous experience acquired for similar problems. In particular, very accurate estimations of the initial values of M,  and  D can be obtained from the experimental data, the first as the value before (and after) the perturbation, the second as the decay of consecutive peaks of the oscillating signal, and the third, measuring the time interval between consecutive peaks of the signal. Other critical aspects are how to select the time interval in which experimental data must be recorded either in terms of the size and position 21 of the interval and the number of control points where the error functional  must be evaluated. A sensitivity study was carried out in this research in order to optimize the performance of the fitting procedure. All cases considered in the parametric analysis are summarized by Fig. 11. The experimental curve taken as an example in Fig.11 corresponds to the first procedure (Section 2.1), reason by which the reaction force (data in the vertical axis) presents negative values, the load cell being always in compression. The following general nomenclature was adopted: case.subcase (tinitial, tfinal, density). Density=1 means that all data points were used to build the error functional while density=0.1 means that only 1 data point out of every 10 data points was used. Cases 1) For fixed sizes (i.e. between 1.4 and 1.7 s) of the data time interval, different time origins were chosen for the experimental force-time plots. Subcase "a" corresponds to tinitial=1.4 s; subcase "b" corresponds to tinitial=1.5 s; subcase "c" corresponds to tinitial=1.6 s; subcase "d" corresponds to tinitial=1.7 s. Cases 2) For a fixed time origin (i.e. 1.4 s), different sizes of the time interval were considered. Subcase "a" corresponds to a time interval of 0.5 s; subcase "b" corresponds to a time interval of 0.6 s; subcase "c" corresponds to a time interval of 0.7 s; subcase "d" corresponds to a time interval of 0.8 s. Cases 3) For a fixed origin (i.e. 1.4 s) and size (i.e. 2.2 s) of the time interval, different densities of the experimental data points were chosen in order to build the error functional  . Subcase "a" corresponds to density=1; subcase "b" corresponds to density=0.5; subcase "c" corresponds to density=0.2; subcase "d" corresponds to density=0.1. The results of sensitivity analysis obtained for the stiffness coefficient k and the damping coefficient c are summarized in Table 1. Remarkably, the result 22 of fitting process does not significantly depend on the size and the position of time data interval as well as it is independent of sampling. The coefficient of variation defined as the ratio between the standard deviation and the average value always was lower than 3%. The same sensitivity study was carried out for the experimental data obtained from the second procedure. Analogous qualitative and quantitative results were obtained, but they have not been included here for the sake of brevity. The results found showing the insensitivity of the procedure to size and position of the interval under study, give confidence on the robustness of the procedure itself. Moreover, the results suggested a criterion for choosing the data interval where to carry out the fitting process. The data interval was set to be the largest in size and with the maximum density by including all possible data points. In general, the largest size of the data time interval is determined by two transient effects at the beginning and end of the vibration period. Experimental data may initially appear irregular just after the instant at which impact occurs because of the transient effect of the impact itself (this is consistent with other study reported in literature: see for example, Fig. 2 of Fukashiro et al [12]). Data relative to the final part of the vibration process when the oscillation is almost completely damped also should be avoided as the signal to noise ratio may get close to 1. By cutting properly the initial and final parts of the data time interval it is possible not to introduce numerical disturbances in the fitting process. Once the values of the five unknown parameters were determined by minimizing the error functional (29), the values of the stiffness K and damping C of the equivalent vibrating system were calculated with Eqs. (13) and (14) for the first load case. The reaction force f developed in the TS was calculated with Eq. 23 (16). Equations (17) finally provided the apparent values of stiffness k and damping c of the TS. For the second loading case, Eq. (16) was again used to evaluate the reaction force f developed in the TS. Finally, Eqs. (25) allowed visco-elastic parameters k and c to be determined. 3.- Results and discussion The apparent values of stiffness k and damping c were determined for a population including seven males and three females. Each subject was submitted to both loading conditions in order to have a homogeneous basis of comparison between the two measurement protocols. The right leg of each subject was tested as this was the case for which both equipments developed were fully usable. The characteristics of the population in terms of sex, age, height, weight and R/r ratio are summarized in Table 2. For comparison purposes, representative values of populations of Ref [12] and [13] have also been included in Table 2. As can be seen, the average age of the populations of the three studies is quite similar as well as the size of the population (6 in Ref. [12], 10 in Ref. [13], 10 in the present study). However, whereas height and weight are quite similar in what respect to the present study and reference [13], they both are significantly higher than those corresponding to the population of reference [12]. For each load level, experimental tests were replicated four times to check the level of repeatability of measurements. For that purpose, tests were performed sequentially without altering the general configuration of the testing apparatus and taking a proper time interval lasting several seconds between each replication. In particular, the foot position was unaltered (no repositioning) during these four tests. In this way, although there is a device to reproduce the position of the foot 24 in the load cell in all tests, uncertainties in the values of R/r are completely avoided. Figures 12 and 13 show the results obtained in the case of a 30 kg applied weight respectively for the first loading case described by the vertical translation of the lower leg and the second loading case described by the rotation of the foot around the ankle articulation. The plots present the evolution of the reaction force (F) sensed by the load cell during the tests. Experimental data appear superimposed on the corresponding theoretical values predicted by the 1-degree of freedom models to which were given in input the parameter values obtained from the fitting procedure described in Section 2.3. Notice that while in Fig. 12 the value of the signal (F) oscillates around a value in between 41 and 42 kg, in Fig. 13 it oscillates around a value in between 32 and 33 kg. In both figures the applied weight is the same (30 kg), but in the first procedure (Fig. 12) there is an additional oscillating mass of the lower part of the leg that produces an additional weight of around 11 kg, whereas in the second procedure (Fig. 13) the additional weight involved in the oscillation is only around 2.5 kg The figures are relative to Subject 1. Since data trends obtained for the other subjects are very similar they will not be reported in the paper for the sake of brevity. Results below each plot in Figures 12 and 13 include the steady-state value of the reaction force fss developed in the triceps surae after the stabilization of the vibration signal. This quantity was obtained from the steady-state force value Fss measured by the load cell (see Fig. 10) and using Eq. (3). The values of apparent stiffness k and viscosity c obtained via least-square fitting are also indicated in the figure. The maximum and minimum values, range of variation, average value, standard deviation and the coefficient of variation for fss, k and c are reported for both testing protocols. The time interval utilized in the fitting process for each test 31 [11] Lafortune, M. A., Hennig, E. M. and Lake, M. J. (1996). Dominant role of interface over knee angle for cushioning impact loading and regulating initial leg stiffness. J Biomech, 29(12), 1523- 1529. [12] Fukashiro, S., Noda, M. and Shibayama, A. (2001). In vivo determination of muscle viscoelasticity in the human leg. Acta Physiol Scand, 172(4), 241-248. [13] Babic, J. and Lenarcic, J. (2004). In vivo determination of triceps surae muscle-tendon complex viscoelastic properties. Eur J Appl Physiol, 92(4-5), 477-484. [14] Faria, A., Gabriel, R., Abrantes, J., Brás, R. and Moreira, H. (2010) Triceps-surae musculotendinous stiffness: Relative differences between obese and non-obese postmenopausal women. Clin Biomech, 24, 866-871. [15] Hill, A. V. (1938). The heat of shortening and the dynamic constants of muscle. Proceedings of the Royal Society of London. Biological Sciences, Series B. 126(843), 136-195. [16] Zatsiorski, V. and Seluyanov, V. (1983). The mass and inertial characteristics of the main segments of the human body. In: Matsui, H. and Kobayahi K (eds). Biomechanics VIII-B, pp 1152-1159. Human kinetics. [17] Zatsiorsky, V.M., Seluyanov, V.N. and Chugunova, L., (1990). In vivo body segment inertial parameters determination using a gammascanner method. In: Berme, N., Cappozzo, A. (Eds.), Biomechanics of Human Movement: Application in rehabilitation, sports, and ergonomics. Bertec Corporation, Worthington, OH, pp. 187–202. [18] Donahue, T.L.H., Hull, M.L., Rashid, M.M. and Jacobs, R.C., (2000). A Finite element model of the human knee joint for the study of tibio-femoral contact. J. Biomech. Eng. 124, 273-280. [19] Strickland, M.A., and Taylor, M. (2009). In-silico wear prediction for knee replacements - methodology and corroboration, Journal of Biomechanics, 42, 1469- 1474. [20] París-García, F. (2010). In-vivo determination of the viscoelastic properties of the Triceps Surae by means of the free vibration technique. Ph.D. Thesis (in Spanish). University of Seville. [21] Clough, R.W. and Penzien, J. (1975). Dynamics of Structures, Prentice-Hall. 32 Figures Fig. 1. Schematic of the experimental setup for determining TS visco-elastic properties based on the vertical displacement of the lower leg. Fig. 2. Global view of the experimental apparatus based on the measurement of lower leg vertical displacements. 33 a) b) Fig. 3. Experimental setup based on the measurement of the vertical displacement of the lower leg: a) Simplified scheme of an arbitrary position taken in the vibration process; b) Forces involved in the equilibrium around the ankle. Fig. 4. a) Equivalent dynamic system including one degree of freedom; b) Force diagram. U c k R r L M F u U   F F tf f Tibia-Fibula in compression Triceps Surae in tension Ankle articulation M K C F U MU KU CU . KU CU . F a) b) M K C F U MU KU CU . KU CU . F a) b) 34 a) b) Fig. 5 Force diagram for the equivalent dynamic system of TS muscle-tendon complex. Fig. 6. Schematic of the experimental setup for determining TS visco-elastic properties based on the rotation of the foot around the ankle articulation. ku cu . f M kc u ku cu . 1 3 2 4 5 x y 90º 90º foot rotation 90º 35 Fig. 7. Global view of the experimental apparatus based on the measurement of the rotation of the foot around the ankle articulation. Fig. 8. Experimental setup based on the measurement of the rotation of the foot around the ankle articulation: simplified scheme of the foot functioning during the free vibration. . M r F ,  F k cf u  R r r G ,  G r f ,  f R G F . M r F ,  F k cf u  R r r G ,  G r f ,  f R G F 36 Fig. 9. Forefoot and rearfoot distances measurements by means of a pedigraphy. Fig. 10. Typical force-time diagram registered by the load cell. F ss time (s) load (kg) -45 -40 -35 -30 -25 1.5 2 2.5 3 3.5 37 Fig. 11. Sensitivity parameters analyzed in order to optimize the fitting process. Fig. 12. Reaction force F vs time results obtained for Subject 1 for the loading case based on the measurement of the vertical displacement of the right leg. 1.4 1.6 1.8 2.0 2.2  44  42  40  38  36  34  32 1.a (1.4/1.9/1) 1.b (1.5/2.0/1) 1.c (1.6/2.1/1) 1.d (1.7/2.2/1) 2.a (1.4/1.9/1) 2.b (1.4/2.0/1) 2.c (1.4/2.1/1) 2.d (1.4/2.2/1) 3.a (1.4/2.2/1) 3.b (1.4/2.2/0.5) 3.c (1.4/2.2/0.2) 3.d (1.4/2.2/0.1) time (s) Load (Kg) Case.Subcase (t initial /t final /density) N Range Maximun Minimun S.D. Coefficient of Variation (%) Force fss N (Mw=30kg) 4 47 1101 1054 20.5 1.9 Stiffness k kN/m (Mw=30kg) 4 13 189 176 6.1 3.4 Viscosity c N·s/m (Mw=30kg) 4 399 1786 1387 170.2 11.0 fss= 1101.49 (N) k= 189.046 (kN/m) c= 1387.146 (N·s/m) fss= 1053.69 (N) k= 175.602 (kN/m) c= 1526.787 (N·s/m) SUBJECT: 1 RIGHT LEG WEIGHT APPLIED: 30 KG PROCEDURE BASED ON THE VERTICAL DISPLACEMENT OF THE LEG fss= 1072.27 (N) k= 177.702 (kN/m) c= 1485.3253 (N·s/m) fss= 1087.43 (N) k= 184.009 (kN/m) c= 1785.964 (N·s/m) STATISTICAL DATA Average 1078.7 1546.3 181.6 1.4 1.6 1.8 2.0 2.2  46  44  42  40  38  36 1.4 1.6 1.8 2.0 2.2  46  44  42  40  38  36 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2.0  46  44  42  40  38  36 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2.0  46  44  42  40  38  36 Reaction Force (kg) Reaction Force (kg) Reaction Force (kg) Reaction Force (kg) Experimental data Fitted evolution 38 Fig. 13. Reaction force F vs time results obtained for Subject 1 in the loading case based on the measurement of the rotation of the foot around the ankle. N Range Maximun Minimun S.D. Coefficient of Variation (%) Force f ss N (M w =30kg) 4 1 839 838 0.4 0.0 Stiffness k kN/m (M w =30kg) 4 12 106 94 5.6 5.6 Viscosity c N·s/m (M w =30kg) 4 51 152 101 23.8 19.9 SUBJECT: 1 RIGHT LEG WEIGHT APPLIED: 30 KG PROCEDURE BASED ON THE ROTATION OF THE FOOT AROUND THE ANKLE f ss = 838.05 (N) k= 93.9518 (kN/m) c= 100.5 (N·s/m) f ss = 838.618 (N) k= 98.7411 (kN/m) c= 103.2(N·s/m) f ss = 838.764 (N) k= 104.477 (kN/m) c= 152.1 (N·s/m) f ss = 838.837 (N) k= 106.315 (kN/m) c= 123.3 (N·s/m) 119.8 STATISTICAL DATA Average 838.6 100.9 2.0 2.2 2.4 2.6 2.8 30 31 32 33 34 35 2.0 2.2 2.4 2.6 2.8 3.0 30 31 32 33 34 35 2.0 2.2 2.4 2.6 2.8 30 31 32 33 34 35 2.0 2.2 2.4 2.6 2.8 30 31 32 33 34 35 Reaction Force (kg) Reaction Force (kg) Reaction Force (kg) Reaction Force (kg) Experimental data Fitted evolution 39 Fig. 14. Variation of the apparent stiffness of TS complex with respect to the nominal applied load for the two measurement protocols and for the ten tested subjects.. 0 50 100 150 200 250 300 0 10 20 30 40 50 Stiffness values k(kN/m) Weights (kg) (subject 1) Vertical displacement Rotation 0 50 100 150 200 250 300 0 10 20 30 40 50 Stiffness values k(kN/m) Weights (kg) (subject 2) Vertical displacement Rotation 0 50 100 150 200 250 300 0 10 20 30 40 50 Stiffness values k(kN/m) Weight (kg) (subject 3) Vertical displacement Rotation 0 50 100 150 200 250 300 0 10 20 30 40 50 Stiffness values k(kN/m) Weight (kg) (subject 4) Vertical displacement Rotation 0 50 100 150 200 250 300 0 10 20 30 40 50 Stiffness values k(kN/m) Weight (kg) (subject 5) Vertical displacement Rotation 0 50 100 150 200 250 300 350 0 10 20 30 40 50 Stiffness values k(kN/m) Weight (kg) (subject 6) Vertical displacement Rotation 0 50 100 150 200 250 300 350 0 10 20 30 40 50 Stiffness values k(kN/m) Weight (kg) (subject 7) Vertical displacement Rotation 0 50 100 150 200 250 300 350 0 10 20 30 40 50 Stiffness values k(kN/m) Weight (kg) (subject 8) Vertical displacement Rotation 0 50 100 150 200 250 300 350 0 10 20 30 40 50 Stiffness values k(kN/m) Weight (kg) (subject 9) Vertical displacement Rotation 0 50 100 150 200 250 300 350 0 10 20 30 40 50 Stiffness values k(kN/m) Weight (kg) (subject 10) Vertical displacement Rotation 40 Fig. 15. Variation of the apparent viscosity of TS complex with respect to the nominal applied load for the two measurement protocols and for the ten tested subjects. 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 1) Vertical displacement Rotation 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 2) Vertical displacement Rotation 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 3) Vertical displacement Rotation 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 4) Vertical displacement Rotation 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 5) Vertical displacement Rotation 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 6) Vertical displacement Rotation 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 7) Vertical displacement Rotation 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 8) Vertical displacement Rotation 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 9) Vertical displacement Rotation 0 500 1000 1500 2000 2500 3000 3500 0 10 20 30 40 50 Viscosity value c( N·sg/m) Weights (kg) (subject 10) Vertical displacement Rotation