Influence of indentation test factors on the mechanical response of the skin

This study proposes in vivo tests and design of experiments to determine the influence of experimental factors on the mechanical response of the soft tissue. The experimental factors considered are: room temperature (A), indentation velocity (B), indenter temperature (C), pump pressure (D) and muscle activation (E). An inverse method was developed to obtain the constants for constitutive equations of a multilayer biological model (skin, hypodermis, and muscle) through the use of indentation tests in combination with a finite element method. For each combination of the experimental factors, two groups of constants were established from the inverse method. Sixteen combinations of experimental conditions and their corresponding constants for the Mooney-Rivlin constitutive equations were obtained to be used in further numerical models. The factor D and factor interactions ADE, CDE, and ACDE were statistically significant with respect to skin mechanical response. Therefore, it can be concluded that there is not a current equation able to represent the mechanical properties of the skin under all the experimental conditions considered in this study


Introduction
The skin is a biologic and multilayer material with intrinsic biomechanical properties that are suited to its function to regulate thermal exchange, maintain water equilibrium, and protect against ultraviolet radiation [1].The study of the mechanical properties of soft tissues is valuable to evaluate the progress of diseases such as tumors [2] and treatment effects [3], to simulate medical surgeries [4], to understand the aging process [5,6] and to evaluate the efficacy of cosmetic products [7].Changes in the tissue elasticity is an indicator of the onset and progression of diseases in specific organs [3].Therefore, the palpation technique is a medical practice used in the exploration of penis, hernias, lymph nodes, breast, skin, abdominal regions and others.For example, the evolution of breast cancer [2,8,9] is associated with an increase in stiffness of the pathological tissue.However, a lack of quantitative tests prevents the identification of the evolution of a disease according to specific mechanical properties.
The effects of the following experimental factors on the mechanical behavior of tissues have been studied in the literature: (a) different non-invasive techniques [10].The most commonly used methods report a wide range of variability from one study to another.For example, Young's modulus of skin varies between 0.42 MPa and 0.85 MPa in torsion testing [11], between 4.6 MPa and 20 MPa in tensile testing [12], between 0.05 MPa and 0.15 MPa in suction testing [13], [14] and between 1.1 kPa [15] and 8 kPa [16] in indentation testing.This wide range of reported values shows the influence that techniques and experimental conditions have on the estimation of mechanical properties.(b) bodily region of human skin, which are related to the elastic and anisotropic properties [17,18]; (c) methods of sample attachment and procedures for obtaining samples [17,19] (e.g., precompression stiffens mechanical response) [20]; (d) experimental scale (nano, micro and macroscopic), which measures different tissue properties and enables the determination of the individual responses of each tissue layer [13]; (e) testing method (in vivo, in vitro, or ex vivo), which alter the behaviors of biological tissues [19] (e.g., changes in stiffness values of up to 50 % have been reported for different testing methods) [21]; (f) elastic behavior, which can depend on subject characteristics such as age, gender, color and type of skin pigmentation, dietary habits and use of tobacco, vascularization, and presence of disease [22,23]; (g) room humidity and temperature, which affect skin properties [24] (e.g., the Young's modulus of stratum corneum decreases when hydration increases) [24].All of these studies have used or developed different constitutive equations to represent soft tissue mechanical behavior.A possible explanation for differences in reported mechanical behaviors of skin soft tissue may be because no study has simultaneously analyzed the combinatory effect of different experimental conditions on the mechanical properties of soft tissue.
In this study, it is hypothesized that numerous experimental conditions applied at the same time can affect the estimation of constitutive equations that represent mechanical response of the tissue.It was decided to lock noticeable factors on the best value, such as (a) indentation technique, (b) forearm, (d) millimeters scale, (e) in vivo and (f) specific population, it is known they have significant incidence.Thus, the aim was to establish the experimental factors (such as room temperature, indentation velocity, indenter temperature, and muscle activation) that affect the mechanical response of multilayer soft tissues, specifically of the skin (a layer of skin, hypodermis and muscle).Furthermore, an inverse analysis is proposed to establish constitutive equations of the skin.It compares in vivo experimental measurements with the results of a finite element method analysis (FEA) that reproduces indentation conditions

Materials and methods
An inverse analysis was carried out using the following steps: first, an indentation test was performed to obtain experimental results of indentation force, distance and rate; then, a FEA with identical geometry and indentation protocol and initial a priori mechanical properties was simulated; finally, an optimization algorithm solved the inverse problem to find the constants for the Mooney-Rivlin constitutive equation that fit the initial experimental curve based on the mechanical parameters used in the numerical model (Fig. 1).The method was validated and applied in previous studies on linear and non-linear mechanical behavior and on monolayer and multilayer materials [25][26][27].

In vivo indentation test
The experiment was initially restricted to be performed as an in vivo indentation test on multilayer soft tissues because the test does not alter the composition of skin and allows to obtain reliable force and displacement values.On the other hand, forearm is relatively flat, easily accessible, and less disturbed by the involuntary movements of the body.For these reasons, the indentation tests were carried out on the right ventral forearm zones of right-handed people.The experimental assembly (Fig. 2) includes a device for ensuring that a forearm is maintained perpendicular to the indenter, a pump that exerts pressure on an individual's forearm to minimize involuntary movements, and a mechanism to activate the muscle.Previously, the assembly was used to study in vivo the mechanical properties of human skin [28].
As it is known, mechanical properties of the skin depend on age, gender, habits, and others [22].To reduce the influence of age and gender on the mechanical response, tests were performed on four male individuals between 17 and 19 years old.All of them had a body mass index of 20 and practice football 4 to 6 hours per week.Physical activity should not strengthen the upper body because it can generate variations in mechanical response that are due to muscle toning and not to the composition of the skin.Random procedures were performed according to the experimental design.A spherical tip indenter with a 6.35 mm diameter was selected because rounded edges minimize stress concentration and eliminate any possibility of pain.The indentation depth was 3 mm, which exceeds twice the thickness of the first two layers to guarantee a multilayer tissue response.The experiment was conducted using a TA-XT2i Texture Analyzer, which measures the normal forces related to stainless steel indenter in contact with the skin.During indentation, reaction force and indentation depth are simultaneously measured.Hence, a force-displacement curve can be generated out of these data, which is later used for inverse analysis.
An experimental design was carried out to decide how many replications of each treatment would be required to obtain the maximum information of the test at the lowest cost [29].A split-plot experiment was designed to identify the incidence of the experimental factors in soft tissue mechanical response.The factors considered in this study included the following: -Room temperature (A), which may contribute to skin stiffness.
-Indentation velocity (B), to identify time-dependent skin properties.
-Indenter temperature (C), which may contribute to skin stiffness.
-Pump pressure magnitude (D), to ensure perpendicularity between the sample and the indenter.
-Muscle activation (E), which may change multilayer stiffness.
Before conducting the main experiment, the reproducibility in controlling each experimental factor was studied.Skin temperature is quite difficult to control (between low [30 • C] and high [32 • C] levels [30]) because it depends on an air conditioning system and on the time needed to adjust a room's temperature, the relative humidity is controlled with the air conditioned system [31].To ensure thermal conditioning, subjects were placed in the study room 15 min before testing.Due to this restriction, it was not possible to completely randomize the runs, which made it necessary to execute the experiment in two stages.This situation leads to a split-plot experimental design in which the whole plot corresponds to the temperature factor, and the subplots correspond to combinations of the levels of the other four factors.A fractional design with two complete plot replications with a single block was proposed, leading to 64 runs.Each factor was considered in two levels of experimentation that by convention are called "low" (-1) and "high" (+1).The low and high levels of each factor are indicated in

Muscle activation Relaxed Active
Table 1.Factors and their levels in the design of the experiment.
were selected carefully to avoid any damage over the volunteer but keeping the experimental changes.For instance, difference of 10 • C in environment temperature can reflect 2 • C on the skin [30], the pressure of the pump is needed to fix the forearm and limit all involuntary movements, and the indentation velocity was estimated between 300 um/s and 3 mm/s, considering that the time of tissue response and load application is extended from 1 to 10 seconds.
The epidermis and dermis layers were assumed to form only one skin mechanical layer [13].Thus, the mechanical response was attributed to the cutis and hypodermis layers.The skin thickness of human forearm was not directly measured in this experiment, and it was assumed that the total skin thickness was approximately 2.4 mm, as measured by Hendricks et al. [13].

Numerical simulation
An axisymmetric finite element model of the skin was created.The contact between the indenter and the skin was modeled as a spherical rigid body (the indenter) and a plane deformable hyperelastic body (skin).After the convergence analysis, a model with 896 planar 4-noded elements (CAX4R) was used (Fig. 3).The friction coefficient (µ) between the stainless steel indenter and the skin was determined by using a device capable of simultaneously measuring normal and tangential forces.The friction coefficient obtained was 0.20 ± 0.01.
The mesh was divided into three different sets to represent the skin layer (top elements), the subcutaneous or hypodermis layer (intermediate elements), and the muscle layer (bottom elements); see Fig. 3.The three layers were assumed to be perfectly bonded to each other.The thicknesses of the skin (dermis and epidermis) layer, hypodermis layer and muscle layer were 1.2 mm, 1.2 mm [12], and 4 mm, respectively.
It is globally accepted that most biological soft tissues mechanical behavior can be modeled by hyperelastic material models [18].Some authors [13,32,33] employed Mooney-Rivlin models to characterize human skin.The skin layers were modeled using a compressible hyperelastic Mooney-Rivlin formulation, which is represented by the following equation: Where U is the strain energy function; C 10 , C 01 , and D 1 are empirical temperature-dependent material parameters; I 1 and I 2 are alternative invariants of the left Cauchy Green deformation tensor; and J is the Jacobian of the deformation gradient.
For the muscle layer, the following coefficients were used: C 10 = 4.25 kPa, C 01 = 0 Pa, D 1 = 2.36 MPa-1 [32].For the skin and the hypodermis, the coefficient D1 = 2000 MPa-1 was calculated from a bulk modulus of 1 MPa [19], and the coefficients C 10 and C 01 were obtained from the inverse analysis.The FEA were performed using Abaqus 6.11 (Dassault Systemes).

Inverse algorithm
According to the procedure described previously, an initial equivalent Young's modulus that approximately characterizes the skin and hypodermis layers are considered.An implemented version of a Levenberg-Marquardt optimization algorithm was used in Matlab R to find the differences (error) between the experimental measurement curves and the numerical model results and to propose a new value of equivalent Young's modulus with a minimum error [34].
As stiffness for a spherical indenter does not have a linear relation and there are not available force values of computational and experimental data for exactly the same displacement values, the following type of equation is used to adjust data: Where m and k are initially assumed as constants, line P is force, and h is displacement.
This algorithm interpolates between a Gauss-Newton algorithm and a method of gradient descent, which makes the solution more robust.The algorithm uses either a linear or non-linear equation of adjustment.After running multiple simulations to reach the optimal force-displacement relationship that matches the experimental measurements, the constants of Mooney-Rivlin equation are determined.
The constants C 10 and C 01 are defined using the following equation that approximates a Young's modulus [35]: The equivalent Young's moduli for the cutis E c and the hypodermis E h are iterated in a finite element model using the Mooney-Rivlin Eq. ( 1 ).
To find the differences between the experimental and computational curves, m and k from Eq. ( 2 ) were compared.Relative errors were estimated and defined as a logic disjunction (OR): > 0.005 (5)   Firstly, E c is taken as a fixed value while E h is iterated in steps of 0.05kPa until the error is lower than 5 %.Secondly, the value of E h found previously is fixed.Then, E c is iterated in 0.002KPa steps until an error of 2.5 % is attained.Finally, the found E c and E h values are used to establish the constants of the constitutive equation that correspond to the mechanical properties of each layer.
Additionally, it is known that the skin is stiffer than the hypodermis [16] and so the following relationship must be fulfilled: When this condition is not met, E c is iterated until the error is lower than 5 % and E h is iterated in the contrary case.The iterative process allows the resultant simulation curve to be compared to the experimental curve until the error is less than 2.5

Design of experiment
The results of the experiment show that factor D (pump pressure) and factor interactions ADE (temperature -pump pressure -muscle activation), CDE (indenter temperature -pump pressure -muscle activation), and ACDE (temperature -indenter temperature -pump pressure -muscle activation) are statistically significant.The R-squared value of the experimental design is 84.86 %.Normal distribution, equal variances, and randomization of the residuals were verified to guarantee the validity of the experimental design conclusions.
It is useful to fit a regression model to the experimental data in order to predict the value of mechanical response in different values of the factors studied.The statistical model in coded units [-1,1] of the experiment with interactions to the third order is represented by the following: Where y i j k is a predicted response at the point (x i , x j , x k ), that corresponds to an encoded level between [-1,1] of x 1 =A, x 2 =B, x 3 =C, x 4 =D, x 5 =E, and C i corresponds to an estimated coefficient of individual effects, C i j corresponds to an estimated coefficient of double interaction effects, C i j k corresponds to an estimated coefficient of triple interaction effects.Equation ( 7) and the estimated coefficients that are shown in Table 2 can be used to predict the mechanical response of the skin based on the significant factors.

Numerical model
In the design of experiments proposed, the randomization of the individuals who participated in the experiments guaranteed that the following numerical models did not correspond to a specific individual.
According to the experimental design, an inverse analysis was carried out for each combination of factors.The constitutive equation under a Mooney-Rivlin model (Eq. 1) for the skin and hypodermis were obtained for each combination of factors (Table 3).Because B (indentation velocity) was the only factor that had no significant impact on mechanical response, it was not included in the analysis.This factor and its interactions with other factors could be considered replications of the other measurements.
Examples of displacement versus reaction force curves are shown in Fig. 4; the experimental data interpolation curve (Eq.2) is compared with the computational curve.Although it is not possible to identify the individual behaviors of the layers, a change in the initial slope (between 0 and 1 mm) can be seen when only the upper layer (skin) is changed.On the other hand, when only the lower layer (hypodermis) properties are changed, a change in the curve slope can be observed between 1 mm and 2 mm.The relationship between elasticity moduli allows to determine the change in the curve slope and to improve their estimation.Sudden changes in slope and Ec/Eh≤1 relationships contradict the literature; this indicates that the hypodermis is stiffer than the skin, and that there is an error in the estimation of the thickness of the skin, or muscle and skin constitutive equations selected.This answer corresponds to the factor levels [-1 + 1 + 1 + 1] and [+ 1-1 + 1-1]; Eq. ( 6) was neglected when comparing the curves.
Fig. 5 shows an example of the stress state from a combination of experimental factors at different levels, and it can be seen that there is the reaction of all layers.Although the comparison of the skin reaction force was made at 3mm, for computational savings the model has 2 mm depth of indentation.The stress distribution is directly related to the rigidity of the layers; the middle

Conclusions
Low values of factors A (room temperature), C (indenter temperature), D (pump pressure), and their corresponding interactions result in lower values of the mechanical response.The reason that the two first factors (in their lower levels) directly the mechanical response is related to the temperature delta between the skin and the indenter.The skin temperature is affected by the room temperature.Hence, when the room temperature decreases, the skin temperature decreases as well.Then, the temperature delta between the indenter and the skin is smaller in the lower levels of these two factors than in the higher ones.This leads to the conclusion that the skin is more rigid when the temperature delta is higher.However, the aforementioned conclusion is inferred from the room temperature, since no skin temperature was measured.Thus, future work can be focused on properly validating this statement.Furthermore, as expected, the effect of factor D and its interactions on the mechanical response indicate that the skin is more rigid when the skin is pre-stressed.
One limitation of this study was the inability to measure the thicknesses of the different layers (skin and hypodermis) in each individual, which led to the use of fixed values in the numerical analysis.In general, as shown in Fig. 4, there is a good approximation between experimental and computational curves.However, a lack of precision in the constitutive equations can be observed in the following condition: [-1+1+1+1] (temperature -indenter temperaturepump pressure -muscle activation) in which the results represent a hypodermis layer that is more rigid than the skin layer, indicating a possible error in the estimations of the thickness of each layer of the skin for a given individual, or an error when considering constant coefficients for muscle properties.
In addition, in this preliminary study, the contribution of the important structural components of the dermis: collagen, elastin, and ground substance [1], was neglected.
It is important to take into account that layer thickness is statistically significant with respect to the mechanical response of multilayer tissues and that evidently that the precision of the constitutive equations also depends on this factor [36].However, physiological variations can alter the mechanical behavior of the skin [22] and are specific to each individual.Although it is known that specific morphological conditions [22,37] of volunteers have an incidence on the mechanical response, this study does not pretend to explain it.The randomization of the volunteers allows controlling the effects of extraneous variables.It is assumed that, on average, morphological factors will affect treatment conditions equally; so any significant differences between conditions can fairly be attributed to the independent variables.Inter-subject variability was not studied.
Indentation velocity did not influence skin mechanical properties, which could be due to the low velocities that were used during the experiment.This correlates with a study by Su et al. [38], who carried out indentations at different strain rates to evaluate the viscoelastic response of multi-layer skin and who could not find simulation results corresponding to the experimental data.Future work could include analyses of higher velocities or larger levels of differences.
The elastic modulus of skin is affected by various factors, such as the amount of deformation (due to non-linear stress-strain behavior) and skin thickness.It also varies considerably for different experimental techniques.In this study, the skin elastic modulus varied between 0.372 kPa and 420 kPa, while the hypodermis elastic modulus varied between 0.338 kPa and 59.2 kPa.These values are of the same magnitude order as reported by others [16,38].Pailler-Mattei et al. [16] reported an elastic modulus of 35 kPa for the cutis and 2 kPa for the hypodermis.Their indentation tests were conducted for a constant indentation speed of 400 µm/s, and the indenter used was conical, whereas the mechanical model was a simple assembly of three springs.Su et al. [38] found a Young's modulus between 20.2 kPa and 24.4 kPa for the hypodermis under a constant indentation speed of 0.2 mm/s.These values are within the range found in this study, and the lowest variability encountered by Sue et al. [38] may be because their simulations accounted for real thicknesses.Nevertheless, it is clear that the wide range of values obtained as responses for the skin and hypodermis elastic moduli reflect the great influence that experimental factors have on mechanical properties.
The following was found for each combination of factors when establishing the factors that affect the skin mechanical response: (a) a constitutive equation that characterized the mechanical behavior of two layers (skin and hypodermis) (Table 3) and (b) an equation that predicted the skin reaction force (Eq.7).
The wide range of constitutive equation coefficients that were obtained to model the skin, hypodermis and skin reflect the great influence that skin thickness estimation and experimental factors have on mechanical properties.
The method proposed in this study can be useful to determine the mechanical properties of the skin in a patient-specific manner for future in silico models or to quantify the evolution of skin properties during skin treatment.
Future studies must include measurements of the thickness of each layer of the skin for improving the precision of the constitutive equations.Moreover, the creation and application of a new constitutive equation that consider all experimental condition reduce the wide range of coefficients obtained, and it allows to design models for particular uses or specifications.

Figure 1 .
Figure 1.Inverse analysis diagram to derive skin mechanical properties.

Figure 2 .
Figure 2. Experimental assembly: device for ensuring the position of the forearm, the mechanism of activating the muscle, and the exposed experimentation area are indicated.

k e x p
e r i me nt al −k c om p u t at i onal k e x pe r i me nt al − m c om p u t at i onal e x p e r i me nt al

Table 1 ,
and its values

Table 2 .
Estimated coefficients and effects of reaction force at 3-mm indentation depth.

Table 3 .
Constitutive equation coefficients under a Mooney-Rivlin model that represent the skin under specific experimental conditions.
layer (the hypodermis) is the least rigid, presenting the lowest values for stress distribution.This pattern of stress distribution is consistent with all other numerical models.http://ciencias.javeriana.edu.co/investigacion/universitas-scientiarum