Author's Accepted Manuscript

Rheological behaviour of reconstructed skin C. Pailler-Mattei, L. Laquièze, R. Debret, S. Tupin, G. Aimond, P. Sommer, H. Zahouani

www.elsevier.com/locate/jmbbm

PII: DOI: Reference:

S1751-6161(14)00163-5 http://dx.doi.org/10.1016/j.jmbbm.2014.05.030 JMBBM1172

To appear in: Journal of the Mechanical Behavior of Biomedical Materials

Received date:7 February 2014 Revised date: 16 May 2014 Accepted date: 27 May 2014 Cite this article as: C. Pailler-Mattei, L. Laquièze, R. Debret, S. Tupin, G. Aimond, P. Sommer, H. Zahouani, Rheological behaviour of reconstructed skin, Journal of the Mechanical Behavior of Biomedical Materials, http://dx.doi.org/10.1016/j.jmbbm.2014.05.030 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting galley proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

Rheological behaviour of reconstructed skin C. Pailler-Mattei

a,b*

, L. Laquièze a, R. Debret c, S. Tupin a, G. Aimond c, P. Sommer c, H.

Zahouani a

a

Ecole Centrale de Lyon, Laboratoire de Tribologie et Dynamiques des Systèmes, Université de

Lyon, UMR-CNRS 5513, Ecully, France b

Faculté de Pharmacie de Lyon (ISPB), Laboratoire de Biophysique, Université de Lyon, Lyon,

France c

Institut de Biologie et Chimie des Protéines, UMR-CNRS 5305, Université de Lyon, Lyon,

France.

*Author for correspondence: Cyril Pailler-Mattei E-mail: [email protected]

Figures: 8; references: 59; words without refs, captions and Appendices section: 5711

1. Introduction The histological and the biochemical methods employed to characterise reconstructed tissue models are particularly well-controlled. The main limitation of existing models concerns the determination of their mechanical properties (Monteiro-Riviere et al., 1997; Ponec et al., 2002; Netzlaff et al., 2005). In vitro skin models are based on a three-dimensional (3D) scaffold colonised by dermal cells (fibroblasts) covered by epidermal cells (keratinocytes) to form the 2 outer layers of the skin, the dermis and the epidermis. The hypodermis layer is generally not represented in artificial skins. The 3D scaffold consists of natural or synthetic polymers, or a mixture of both types, assembled to mimic the fibrillar structures of the dermis that support and anchor dermal cells (Drury and Mooney, 2003; Mano et al., 2007; Cheung et al., 2007). Hence the 3D scaffold is used as a template to guide cell growth and tissue development. Several parameters can be modulated to influence the formation of the tissue such as the nature of the polymer, its chemical composition, porosity and pore interconnectivity, 3D organisation and mechanical properties (Shahabeddin et al., 1990; Peppas et al., 2006; Harley et al., 2007; Kanungo and Gibson, 2010). Although the skin structure is highly complex, the mechanical behaviour of whole skin is mainly linked to the dermis structure (Echinard, 1998). The dermis contains a fibrillar network mainly composed of collagen and elastic fibres that contribute to the tissue’s resistance and elasticity, respectively (Silver et al., 2001; Silver et al., 2003). In vivo, the dermis also comprises skin appendages, nerves, and blood and lymph vessels. Under physiological conditions, the dermal fibrillar network is stabilised and consolidated by the enzymatic crosslinking of free lysine residues contained in collagen and elastic fibres. All these fibres are immersed in a liquid medium composed of proteoglycans, ions and water. The surrounding liquid allows the displacement of the proteins that make up the fibres and the fibres themselves, to slide against each other and thus confer to the tissue part of its viscosity. The upper layer of the skin, the epidermis, has a stratified and cohesive structure resulting from the differentiation of the keratinocytes from the basal to the horny layer. In vivo, for indentation tests, the epidermis layer has less influence on the mechanical properties of whole skin than that of the dermis layer. In contrast, engineered in vitro models have been used to demonstrate that the epidermal layer contributes to the mechanical property of whole skin when tensile tests are performed (Ebersole et al., 2010). Soft tissues like skin present anisotropic and viscoelastic mechanical properties (linear or nonlinear depending on the level of deformation). Many techniques, such as stress relaxation tests (Jachowicz et al., 2007), static indentation (Pailler-Mattei et al., 2013), dynamic indentation (Boyer et al., 2009) and oscillating shear tests (Gerhardt et al., 2012; Lamers et al., 2013) have been developed to characterise the viscoelastic behaviour of skin. When the experimental conditions exceed the small strain regime, the quasi-linear viscoelasticity (QLV) method is applied to characterise the viscoelastic behaviour of soft tissues, and corresponds to an adaptation of linear viscoelasticity appropriate for non-linear materials (Fung, 1993). On the basis of classical tensile tests, biaxial traction, applied multi-axial loading and rich set of deformations tests, soft biological tissues have been reported to present anisotropic behaviour (Flynn et al., 2011; Groves et al., 2013; Jor et al., 2011; Limbert, 2011). Generally, the experimental conditions required for these kind of tests are very severe and the mechanical behaviour of skin is non-linear (Flynn et al., 2013; Ni Annaidh et al., 2012). The aim of this

paper is to characterise the viscoelastic mechanical properties of reconstructed skin using an indentation load-relaxation test. Indentation load-relaxation tests are often encountered in the literature for characterising the time-dependent behaviour of biological materials (Ledoux and Blevins, 2007; Hidaka et al., 2011) and cutaneous tissues (Jachowicz et al., 2007). In this study, the generalised Maxwell and Kelvin-Voigt rheological models are first proposed to model the time-dependent behaviour of reconstructed skin. Furthermore, the analytical solutions of the two previous models used to obtain the time-dependent load for a spherical indentation test are also suggested. Relaxation tests are then performed on each structure composing the reconstructed skin: 3D-scaffold, 3D-scaffold+dermis (dermal equivalent), 3D-scaffold+dermis+epidermis (reconstructed skin). The experimental relaxation curves are analysed using the generalised Maxwell and Kelvin-Voigt models. The results permit understanding and modelling the effect of each biological structure on the global mechanical response of the reconstructed skin. Finally, whenever possible, the viscoelastic parameters of each structure composing the reconstructed skin are compared with the values obtained in the literature.

2. Experimental details This study was performed with reconstructed skin samples licensed by BASF® Beauty Care Solutions France. All the human tissues were obtained after informed consent was received. The principle underlying the preparation of reconstructed skin is summarised in Fig. 1a. 2.1. Isolation of human dermal fibroblasts and keratinocytes Normal dermal fibroblasts and keratinocytes were isolated from human adult skin biopsies. The fibroblasts were amplified in medium consisting of Dulbecco’s modified Eagle medium (DMEM) (Invitrogen, Carlsbad, CA, USA) supplemented with 10% foetal bovine serum (FBS) (HyClone, Logan, UT, USA), 25mg/l gentamycin (Sigma-Aldrich, Saint-Quentin Fallavier, France), 100,000 UI/l penicillin (Sigma) and 1mg/l amphotericin B (Sigma-Aldrich). The keratinocytes were amplified in K-FSM (Invitrogen) medium with 25mg/l gentamycin, 100,000 UI/l penicillin and 1mg/l amphotericin B. Cells were used at early passages (below 5) to construct the skin equivalent. 2.2. Scaffold preparation The 3D-scaffold used in this study is a biomaterial with a base of collagen, chitosan and glycosaminoglycan. Eighteen 3D-scaffold substrates were prepared. Briefly, 1.5% bovine collagens (95% type I, 5% type III) were solubilised in an acetic acid solution at a concentration of 0.025M. Chitosan from crab-deacetylated chitin was dissolved at 0.7% in an acetic acid solution at a concentration of 0.084M. Chondroitin-4-sulfate from bovine cartilage was prepared at 2.5% in deionised water. These compounds were mixed to obtain the following final composition in dry matter: 72% collagen, 20% chitosan and 8% chrondroitin-4-sulfate. The final preparation was poured into Snapwell inserts (Corning Life Sciences, New York, NY, USA) at 400 mg/cm2, then freeze-dried overnight (-40°C/+37°C, 0.1mbar) and sterilized by γ-irradiation (10 kGy) (Fig. 1b, c). Six 3D-scaffold substrates were rehydrated with DMEM for 21 days before the mechanical measurements were performed to characterise the acellular dermal equivalents.

2.3. Engineered dermis Dermal fibroblasts were inoculated into twelve 3D-scaffold substrates at a density of 2 ×105 cells/cm2. The fibroblasts were seeded on the top of each 3D-scaffold substrate rehydrated previously with DMEM. The dermal equivalents were incubated at 37°C, 5% CO2 in DMEM with 10% FBS and antibiotics for 21 days. The culture medium was refreshed every day. Six dermal equivalents were set aside after 21 days to perform the mechanical characterisation.

2.4. Reconstructed skin As the interaction between the fibroblasts and keratinocytes is essential for the development of skin tissue, Keratinocytes were seeded at a density of 5 ×106 cells/cm² and cultured for 3 days in the keratinocyte medium composed of 3:1 mixture of DMEM and Ham’s F12 (Invitrogen), supplemented respectively with 10% FBS, 10 ng/ml epidermal growth factor (Sigma-Aldrich), 0.12 IU/ml insulin (Lilly, Saint-Cloud, France), 0.4 µg/ml hydrocortisone (Sigma-Aldrich), 5 µg/ml triiodo-L-thyronine (Sigma-Aldrich), 24.3 µg/ml adenine (Sigma-Aldrich) and antibiotics, as described for the fibroblast culture medium. Then, the culture was elevated to the air–liquid interface and cultured for an additional 14 days in a simplified keratinocyte medium containing DMEM supplemented with bovine serum albumin 1.6 g/l (Sigma-Aldrich), 0.12 IU/ml insulin, 0.4 µg/ml hydrocortisone, 10 µg/ml L-ascorbic acid and antibiotics. The medium was changed every 2 days. Finally, six reconstructed skin samples were obtained to perform the mechanical characterisation. All the 3D-scaffold, dermal equivalent and reconstructed skin samples had a surface area of 0.8 cm² and a thickness of about 2 mm (Fig. 1b). 2.5. Histological procedures The samples were fixed with Bouin’s fixative (Sigma-Aldrich) and embedded in paraffin. Sections (8µm) were de-paraffinized and rehydrated prior to counterstaining with Meyer’s hematoxylin (Sigma-Aldrich) for 20 minutes and eosin Y (Sigma-Aldrich) for 2 minutes and images were acquired using a Leica DM 750 optical microscope. 2.6. Experimental device An original light load indentation device, based on the technique developed previously for cutaneous tissues in vivo (Pailler-Mattei et al., 2009; Pailler-Mattei et al., 2011), has been designed to study the mechanical properties of reconstructed tissues. Traditionally, the indentation test consists in recording the penetration depth, δ , of a rigid indenter as a function of applied normal load, Fz , during a loading/unloading experiment (Pailler-Mattei et al., 2008). In the present study, the indentation device was used in loadrelaxation mode. The relaxation test is similar to the indentation test, with a holding step at the maximum penetration depth. The relaxation tests are performed in controlled displacement mode. The z-displacement is obtained from the National Instrument displacement table and controlled by a displacement sensor (Fig. 2). The maximum displacement during the

loading/unloading cycle can reach about 500 μm with a resolution of 1 μm. The experimental device offers a wide range of indenting velocities from 1 to 100 μm.s-1. 2.7. Indentation load-relaxation test The load-relaxation tests were carried out on the 3D-scaffold, the dermal equivalent (dermis+3Dscaffold) and the reconstructed skin (epidermis+dermis+3D-scaffold) samples. The loadrelaxation tests were repeated five times for each sample, and there were six samples for each category of material. The tests were carried out at a controlled constant penetration depth, δ 0 = 50 μ m , and the displacement was applied for a constant indentation speed δ = 50 μ m.s-1 . These parameters (i.e. indentation speed and depth) were fixed for all the experimental measures to avoid complicating the rheological analysis of the three samples. Load–relaxation data were collected at a sampling rate of 1 kHz. The indenter used was a spherical PTFE indenter, with a radius of curvature R=1.6 mm. The tests were performed in an air-conditioned room with controlled temperature (22–24°C) and relative humidity (30–40%). Before the tests, the samples were kept in an incubator at a constant temperature of 37°C and 5% CO2. The tissue samples were removed from their culture-plate and deposited on a glass slide for testing. In order to facilitate reading the results, the outcome of a single sample for each category of material will be presented to illustrate the analysis.

3. Analytical solutions

Most of the approaches used to resolve the viscoelastic indentation problem are based on the functional equation developed by Lee and Radok (Lee and Radok, 1960). The method of functional equations consists in using the elastic solution to solve the viscoelastic problem by replacing the elastic parameters with their corresponding viscoelastic operators. The method of functional equations can be seen as an extension of the Laplace transform method (Lee, 1955; Ting, 1966). The Laplace transform method consists in eliminating the explicit time dependence of the viscoelastic problem by applying the Laplace transform to the time dependent moduli and solving the corresponding equivalent elasticity problem in the Laplace domain. It has been shown that the viscoelastic solution always satisfies the same boundary conditions as the elastic solution from which it was obtained (Radok, 1957). The viscoelastic behaviour of biomaterials is a central focus of bioengineering research (Patel et al., 2005; Geerligs et al., 2008; Darling et al., 2008; Pan and Xiong, 2010). Many biomaterials and biological tissues exhibit non-linear mechanical time-dependent behaviour often characterised by indentation load-relaxation tests (Gefen and Haberman, 2007; Screen, 2008; Anssari-Benam et al., 2011; Chin et al., 2012). In our case, because the value of the penetration depth applied to the samples tested during the relaxation tests was very low ( δ 0 =50 µm) compared to the total thickness of the biomaterials (2mm), we could assume that the tissue strain level remains in the linear regime and within the hypothesis of small deformation. This hypothesis allows neglecting the effect of the substrate on which the samples were placed. Usually, phenomenological rheological models based on springs and dashpots in series or in parallel are used to model the mechanical behaviour of viscoelastic materials. It is often proposed in the literature to model the viscoelastic behaviour of the biological tissues with the Maxwell model and/or the Kelvin-Voigt model (Forgacs et al., 1998; Elias et al., 2012; Gras et

al., 2013). These two models will be used to analyse the mechanical and rheological behaviours of our different samples (3D-scaffold, dermal equivalent, reconstructed skin). In this part, the viscoelastic solutions of the force-penetration relationship for the spherical indentation load-relaxation tests are indicated for the generalised Maxwell and Kelvin-Voigt models. All the mathematical details concerning the analytical solutions for the spherical indentation load-relaxation tests, based on the generalised Maxwell model and the generalised Kelvin-Voigt model, had been developed and provided in the appendices section. 3.1. The generalised Maxwell model The generalised Maxwell model consists of a spring connected in parallel with n-branches containing a spring and a dashpot connected in series (Fig. 3). Based on the generalised Maxwell model, an analytical solution using a functional equation for the spherical indentation loadrelaxation test is proposed as: G n − it ⎞ ⎛ ηi 3K + ⎜ G0 + ∑ Gi e ⎟ G ⎜ ⎟ n − it ⎞ = 1 i 16 12 3 2 ⎛ ⎝ ⎠ Fz ( t ) = R δ 0 ⎜ G0 + ∑ Gi e ηi ⎟ Gi ⎜ ⎟ n 3 i =1 ⎝ ⎠ 3K + 4 ⎜⎛ G + G e− ηi t ⎟⎞ i ⎜ 0 ∑ ⎟ i =1 ⎝ ⎠

(3.1)

3.2. The generalised Kelvin-Voigt model The generalised Kelvin-Voigt model consists of a spring connected in series with n-branches, containing a spring and a dashpot connected in parallel (Fig. 4). Based on the generalised Kelvin-Voigt model, an analytical solution using a functional equation is proposed for the spherical indentation load-relaxation test as: −1

G n ⎛ 1 − i t ⎞⎞ 1 ⎛ ηi + ∑ ⎜1 − e ⎟⎟ −1 3K + ⎜ Gi ⎜ ⎟⎟ ⎜ n G G ⎛ ⎞ t − ⎛ ⎞ i = 1 3 1 i 0 16 2 2 1 1 ⎝ ⎠⎠ ηi ⎝ ⎜ ⎟ = + − Fz ( t ) R δ0 ⎜1 e ⎟ ∑ −1 ⎜ ⎟ G ⎜ G0 i =1 Gi ⎟ 3 n ⎛ 1 − i t ⎞⎞ ⎛ ⎝ ⎠⎠ ⎝ 1 + 3K + 4 ⎜ ⎜ 1 − e ηi ⎟ ⎟ ⎜ ⎟⎟ ⎜ G0 ∑ i =1 Gi ⎝ ⎠⎠ ⎝

(3.2)

4. Results and discussion

In this part, we propose to use the previous analytical models (equations (3.1) and (3.2)) to analyse the mechanical behaviour and the contribution of the different components by increasing

the complexity of the skin model (3D-scaffold, scaffold+dermis+epidermis) during the load relaxation-tests.

3D-scaffold+dermis,

3D-

The experimental data obtained from the relaxation curves (Fig. 5) were fitted to the analytical models (equations (3.1) and (3.2)) by the simplex method (Dantzig, 1963; Nelder and Mead, 1965) to obtain G0 , Gi , K and ηi . From the optimised values of G0 , Gi and K obtained using the simplex method, it is possible to obtain the values of the corresponding Young moduli E0 and

Ei , as: 9 KG0 3K + G0 9 KGi Ei = 3K + Gi

E0 =

(4.1) (4.2)

Moreover, in the following when choosing between the different rheological models, or to quantify their similarities compared to the experimental curves, an error criterion, noted ξ , has been suggested. It calculates the smallest error, in percentage, between the experimental relaxation curve and the theoretical relaxation curves obtained from the different rheological models. It is defined as: th 1 n Fz ( t )i − Fz ( t )i × 100 ξ= ∑ n i =1 Fz ( t )i

(4.3)

where Fz ( t ) is the experimental value of the normal load measured during the load-relaxation tests and Fzth ( t ) is the theoretical value of the normal load for the load-relaxation tests obtained from the rheological models. n is the total number of points in the curve and i is the ith point on the curve. Finally, in the next part of this paper, all data are indicated as mean ± SD (SD: standard deviation). The standard deviation values were calculated from the whole of the results obtained on the same tissue sample group (3D-scaffold, dermal equivalent, reconstructed skin). 4.1. Histological analysis Tissue organisation is usually described by histological analyses. This descriptive approach provides information on cell density and cohesion, material deposition and the size of anatomical compartments that are important data directly implicated in the mechanical behaviour of tissue. Although the participation of each entity in mechanical response is presented in specific paragraphs of the results section, a concise histological interpretation of the reconstructed skin model is presented below.

The different samples were embedded in paraffin and then stained with hematoxylin and eosin dyes (Fig. 1d, e and f). The 3D-scaffold alone appears very lacunar due to highly interconnected pores (Fig. 1d). Pink stain mainly corresponds to chitosan deposits while collagen and proteoglycans cannot be distinguished from the polymer surface. The dermal equivalent seems more cohesive and dense (Fig. 1e). Fibroblast nuclei appear deep blue. After 21 days in cultures, cells were organised in 2 to 3 layers over the 3D-scaffold, and had colonised the scaffold. Fibrillar and amorphous structures are also observable within the pores indicating extracellular matrix deposition by fibroblasts. The reconstructed skin shows classical features of a stratified and fully differentiated epidermis (Fig. 1f). The basal layer is characterised by columnar cells perpendicularly anchored to the dermis surface. The cells then become more polygonal and the nuclei flatten and disappear progressively as the keratinocytes differentiate and move to the skin surface. The pink coloration of the outer layer indicates cornification corresponding to the terminal differentiation state of keratinocytes. The epidermis is highly cohesive due to strong cell/cell interactions and very little extracellular space. These histological features are in agreement with the standard characterisation of this reconstructed skin model (Black et al., 2005).

4.2. Analysis of the experimental data Load relaxation tests for each biological sample were performed and variation of the normal load, Fz ( t ) , as a function of time is reported in figure 5. The values obtained regarding the initial normal load, Fz ( t = 0 ) , were clearly distinguishable for the three samples. For the 3Dscaffold, the initial normal load measured was around 0.3mN, whereas for the dermal equivalent and the reconstructed skin the initial normal loads measured were about 0.6mN and 0.9mN, respectively. Moreover, the relaxation curves decreased from their initial load value to 0.5mN within the first 10 seconds for the dermal equivalent and the reconstructed skin, while the value remained unchanged during this period for the 3D-scaffold alone (Fig. 5). We assigned these differences observed during the short observation periods to the incremental insertion of biological material inside the sample. The addition of fibroblasts in the 3D-scaffold increased its initial stiffness, which may explain the higher value of the initial normal loads measured. The increase of the initial normal loads measured for the reconstructed skin compared to the dermal equivalent can be attributed to the presence of keratinocytes that form a rigid layer on the surface of the dermal equivalent, thus increasing the initial stiffness of the structure as a whole. However, after a longer period of time (t >100s), the three types of sample exhibited similar relaxation curves converging toward the same value of the normal load measured, i.e. 0.2mN < Fz ( t > 100s ) < 0.3mN (Fig. 5). The fibroblasts and the keratinocytes no longer seemed to influence the mechanical response of the samples in which they were present. The mechanical response of the 3D-scaffold structure appeared to correspond to the predominant response for every sample tested. 4.3. Mechanical behaviour and rheological model for the 3D-scaffold

The mechanical properties and the structure of the 3D-scaffold have already been studied in the literature (Freyman et al., 2001; Harley et al., 2007). The mechanical behaviour of the 3Dscaffold is very difficult to interpret because it relies on many parameters, such as the percentage of collagen, pore size, volume and distribution, and the state of hydration (Harley et al., 2007). This part of the paper provides an evaluation of the mechanical and the rheological behaviour of the 3D-scaffold in order to determine its role in the behaviour of reconstructed skin. Fig. 6 represents the average mechanical behaviour of the 3D-scaffold during the load-relaxation test. The relaxation profile of the 3D-scaffold is due to the coupled interaction between the fluid displacement within the pores and the intrinsic properties of collagen. The results indicate that the one-branch Maxwell model (equation. (3.1), with i=1) seems to fit well with the experimental data curves obtained for the 3D-scaffold (Fig. 6). The error criterion, calculated from equation (4.3), for the one-branch Maxwell model is ξ =3.3 ± 1.3%. The rheological parameters extracted from the rheological model are indicated in table 1. The results show that the Young’s modulus E0 is around 3.7 ± 1.2 kPa while the Young’s modulus E1 is around 0.4 ± 0.1 kPa. The viscosity η1 is about 10.8 ± 1.9 kPa.s, which is in the same order of magnitude as the viscosities found in the literature for this kind of biomaterial (Elias et al., 2012). This value of viscosity is mainly related to the behaviour of the collagen fibres inside the 3D-scaffold matrix (Newman et al., 1997; Knapp et al., 1997; Wang and Spector, 2009). The relaxation function, Ψ ( t ) , for the one-branch generalised Maxwell model is given by:

Ψ ( t ) = E0 + E1 e

−

t

τ1

(4.4)

where E0 and E1 are the Young’s moduli obtained from equation (4.1) and equation (4.2). The relaxation time is defined as τ 1 =

η1 E1

and η1 is directly extracted from the equation (3.1).

The values of the relaxation function, Ψ ( t ) , for t → 0 and t → ∞ are indicated in table 1. As t → 0 , the value of the relaxation function obtained is equal to 4.1 ± 1.1 kPa, which corresponds to the elastic mechanical response of the overall 3D-scaffold matrix. As t → ∞ , the mechanical response observed corresponds to the behaviour of the 3D-scaffold after the rearrangement of its internal structure under the effect of the external force. The value of the relaxation function is then equal to E0 . These values of the relaxation function, or Young’s modulus, are in good agreement with the scaffold Young’s modulus values previously published in the literature (Webba et al., 2006; Harley et al., 2007; Elias et al., 2012). We observed almost no difference between Ψ ( t ) and Ψ ( t ) . The two values of the relaxation function are very close to each other. t →0

t →∞

Therefore, the shape of the relaxation curve is essentially related to the viscosity of the 3Dscaffold and not to its elastic properties.

4.4. Mechanical behaviour and rheological model for the dermal equivalent

The dermal fibroblasts were inoculated into the 3D-scaffold matrix to form the dermal equivalent. In this part, we study the role of fibroblasts in the rheological behaviour of the dermal equivalent. Fig. 7 represents the average mechanical behaviour of the dermal equivalent during the loadrelaxation test. It can be seen that adding fibroblasts inside the 3D-scaffold matrix has a significant effect on the rheological behaviour of the dermal equivalent compared to that of the matrix alone, especially for the short period (Fig. 5). Two different rheological models were tested to analyse the mechanical behaviour of the dermal equivalent: the one-branch Maxwell model (equation (3.1), with i=1) and the one-branch KelvinVoigt model (equation (3.2), with i=1). The two theoretical curves obtained from the one-branch Maxwell and Kelvin-Voigt models are reported in Fig. 7. The one-branch Kelvin-Voigt model fits the experimental dermal equivalent curve better than the one-branch Maxwell model, especially at the beginning. The error criterion, ξ , is 13.7 ± 2.7% for the Maxwell model and 4.6 ± 1.1% for the Kelvin-Voigt model. The fibroblasts seem to introduce an instantaneous elastic component to the mechanical behaviour of the dermal equivalent, modelled by the Young’s modulus of the spring E0 in the Kelvin-Voigt model (Fig. 7). The rheological parameters extracted from the one-branch Kelvin-Voigt model are indicated in Table 2. The results show that the Young’s modulus E0 is around 6.7 ± 1.1 kPa, the Young’s modulus E1 is around 6.4 ± 1.4 kPa and the viscosity η1 is about 160 ± 4.6 kPa.s. The relaxation function, Ψ ( t ) , for the one-branch generalised Kelvin-Voigt model is given by: t − * E0 ⎛ τ1 ⎜ E1 + E0 e Ψ (t ) = E0 + E1 ⎜⎝

⎞ ⎟ ⎟ ⎠

(4.5)

where E0 and E1 are the Young’s moduli obtained from equations (4.1) and (4.2). The relaxation time is defined as τ 1* =

η1 E1

and η1 is extracted from equation (3.2).

The values of the relaxation function, Ψ ( t ) , for t → 0 and t → ∞ are indicated in Table 2. As t → 0 , the value of the relaxation function obtained is equal to 6.7 ± 1.1 kPa, which corresponds to the instantaneous elastic response of the dermal equivalent. This Young’s modulus of the dermal equivalent is of the same order of magnitude as the values already observed in the literature (Zahouani et al., 2009). We note that the initial relaxation function of the dermal equivalent is higher than that of the 3D-scaffold. This implies that fibroblasts contribute to the overall mechanical response of the reconstructed tissue, and reflects the fact that the mechanical stress applied during the relaxation test is conducted into the cells. Regarding the dermal equivalent, the cells are attached to the 3D-scaffold through focal adhesion plaques involving integrin subunits, thereby connecting the extracellular matrix fibres to the intracellular cytoskeleton formed by the actin filament network (Silver et al., 2003). Since the actin network provides stiffness to the cell, many studies have modelled the cytoskeleton to better understand

its mechanical behaviour. Interestingly, when cross-linked with rigid components actin networks exhibit exceptional elastic behaviour that reflects the mechanical properties of individual filaments (Gardel et al., 2004). More recently, a theory for predicting the relaxation behaviour of actin networks cross-linked with the flexible cross-linker α-actinin has been proposed (Fallqvist and Kroon, 2013). Under these conditions, the deformation was resolved by a viscoelastic model. This is consistent with our results showing an increase of the dermal Young’s modulus and also a significant increase of viscosity and consequently a decrease of dermal relaxation time. The dermal relaxation time, τ1* , which corresponds to the time necessary for the sample to return to an equilibrium state after external solicitation (Table 2), is shorter than the 3D-scaffold relaxation time (Table 1). As t → ∞ , the dermal equivalent relaxation function and the 3D-scaffold relaxation function are similar. The fibroblasts do not seem to affect the mechanical relaxation response of the dermal equivalent. The influence of the fibroblasts on the whole dermal structure appears to be negligible compared to the 3D-scaffold matrix response (Fig. 5). 4.5. Mechanical behaviour and rheological model for the reconstructed skin Keratinocytes were seeded onto the dermal equivalent to form the reconstructed skin. Reconstructed skin is a stratified structure composed of two different layers, which are the dermal matrix and the epidermis layer. In this section, the effect of the keratinocyte layer on the mechanical behaviour of the reconstructed skin is studied. To analyse the mechanical behaviour of the reconstructed skin the one-branch Kelvin-Voigt model was tested first. The experimental reconstructed skin load-relaxation curve and the theoretical one-branch Kelvin-Voigt curve are shown in Fig. 8. The rheological parameters extracted from the one-branch Kelvin-Voigt model are as follow: E0 = 6.8 ± 1.4 kPa, E1 = 2.5 ± 1.5 kPa, η1 = 43 ± 3.1 kPa.s. These values are close to those observed for the dermal equivalent, except for viscosity. The analysis of the experimental curve and the theoretical curve shows that both curves fit correctly ( ξ =14.1 ± 3.2%), except for the very short time values (Fig. 8). As a consequence, the one-branch Kelvin-Voigt model does not seem sufficient to accurately model the relaxation behaviour of reconstructed skin, especially for the very short time. A more complex rheological model, the two-branch Kelvin-Voigt model, is then proposed to analyse the beginning of the relaxation curve of the reconstructed skin. The behaviour of the reconstructed skin at the beginning of the relaxation curve probably corresponds to the effect of the epidermis. When trying to understand the specific effect of the epidermis on the whole reconstructed skin, the previous average parameters from the dermal equivalent obtained with the one-branch Kelvin-Voigt model were imposed on the two-branch Kelvin-Voigt model. The error criterion, ξ , is 5.7 ± 2.4 % for the two-branch Kelvin-Voigt model. The rheological parameters extracted from the two-branch Kelvin-Voigt model are indicated in Table 3. The results show that the Young’s modulus E2 is 7.6 ± 1.1 kPa and the viscosity η2 is around 18.6 ± 5.9 kPa.s. The relaxation function, Ψ ( t ) , for the two-branch generalised Kelvin-Voigt model is given by:

t − * E0 ⎛ τ1 ⎜ E1 + E0 e Ψ (t ) = E0 + E1 ⎜⎝

t − * ⎞ E0 ⎛ τ2 ⎟+ ⎜ E + E0 e ⎟ E0 + E2 ⎜ 2 ⎠ ⎝

⎞ ⎟ ⎟ ⎠

(4.6)

where E0 , E1 and E2 are the Young’s moduli obtained from equations (4.1) and (4.2). τ1* and

τ 2* are the relaxation times defined as: τ 1* =

η1 E1

and τ 2* =

η2 E2

. Viscosity, η1 , is imposed from the

rheological model of the dermal equivalent and viscosity, η2 , is extracted from the equation (3.2). The values of the reconstructed skin relaxation function, Ψ ( t ) , for t → 0 and t → ∞ are indicated in Table 3. As t → 0 , the value of the relaxation function is equal to 13.4 kPa. This Young’s modulus value is of the same order of magnitude as those observed in the literature for human skin and range from 10 kPa to 18 kPa (Pailler-Mattei et al., 2008). This refers to in vivo parameters where the apparent Young’s modulus is associated with skin tension. Such tensions also exist in vitro as demonstrated through the contraction of floating collagen gels by dermal fibroblasts (Brown, 2013). To avoid gel contraction, the collagen structures are anchored to the cell culture device to provide an opposing force to cell contraction which is a physiological cell response generally observed in the skin wound healing processes close to wounds before reepithelialization, and in pathological conditions such as fibrosis. In the present reconstructed skin model, contraction is limited by the chitosan polymer inside the 3D-architecture that provides a sufficient opposing resistance force (Vaissiere et al., 2000). The addition of stratified and fully differentiated keratinocytes significantly increases the instantaneous Young’s modulus, Ψ ( t ) , of the reconstructed skin. Terminal differentiation of t →0

keratinocytes forms a rigid layer (horny layer) on the skin surface, influencing the instantaneous elastic response of the reconstructed tissue and increasing its global rigidity. As described for fibroblasts in the dermal equivalent, the cytoskeleton of keratinocytes should strongly influence the global behaviour of the reconstructed skin. This assumption is even more realistic for differentiated keratinocytes whose cytoplasm is filled with cross-linked keratins. The keratin network is responsible for the extreme stiffness and resilience of this cell type (Lulevich et al., 2010). Moreover, we cannot exclude that the seeding of keratinocytes on the dermis also stimulates fibroblast metabolism, which leads to more collagen fibres, more tensile stress in the matrix and thus a higher Young’s modulus. The two-branch model introduces two relaxation times, τ1* and τ 2* . The relaxation time τ 2* (Table 3) is much faster than the relaxation time τ1* imposed by the dermal equivalent parameters obtained (Table 2). The rigid keratinocyte layer on the reconstructed skin transferred the normal load applied completely throughout the entire volume of the sample during the relaxation test. Therefore the displacement of the fluid contained in the reconstructed tissue was greatly influenced by the pressure applied on the skin, affecting its rheological behaviour. The shape of the beginning of the relaxation curve of the reconstructed skin results, in part, from the displacement of fluid inside the reconstructed tissue (Fig. 8). Therefore the relaxation time τ 2* reflects the fluid displacement inside the reconstructed skin.

As t → ∞ , the relaxation function of the reconstructed skin is 5.1 ± 0.8 kPa (Table 3). The value of the relaxation function of the reconstructed skin is slightly higher than those of the dermal equivalent and the 3D-scaffold. The epidermis formed by keratinocytes appeared to increase the global stiffness of the reconstructed skin and affect its mechanical response to solicitation over long periods of time. However, the increase of the relaxation function of the reconstructed skin compared to the relaxation functions of the dermal equivalent and the 3D-scaffold was slight. Therefore, the epidermis effect for solicitation over long periods of time was not predominantly involved in the global mechanical response of the reconstructed skin, compared to the collagen effect, as observed in Fig. 5. As mentioned previously, the reconstructed skin is a stratified or anisotropic material composed of a rigid layer (epidermis) set on a fibrillar system (dermal equivalent). In this study, an axisymmetric indenter (spherical indenter) was used for the relaxation tests, which has the effect of partially masking the anisotropic behaviour of the material tested. However, part of the anisotropic response of the reconstructed skin was taken into account by the rheological models. Indeed, the two-branch Kelvin-Voigt model was assumed to model the specific effect of the dermal layer and that of the epidermis layer. It allowed observing the decoupled mechanical response from each of the main structures of the reconstructed skin and thus partially taking into account the anisotropy of the material.

5. Conclusions

Quantitative rheological analysis of samples of three biological tissues (3D-scaffold, dermal equivalent and reconstructed skin) was performed and load-relaxation tests carried out to evaluate their viscoelastic behaviour. These tissues exhibited specific viscoelastic behaviour as a function of their own structures. Different viscoelastic models were proposed to interpret the experimental data and quantify time-dependent parameters. The results indicated that the rheological behaviour of the 3D-scaffold, which is the minimal element common to all the samples, was correctly modelled with a simple one-branch Maxwell model. However, the mechanical behaviour of the dermal equivalent and the reconstructed skin were modelled by a one-branch Kelvin-Voigt model and a two-branch Kelvin-Voigt model, respectively. A rheological transition was then observed between a Maxwell model and a Kelvin-Voigt model, which differed according to the biological structure of the samples tested. Our results showed that fibroblasts and keratinocytes increased the stiffness and the viscosity of the dermal equivalent and the reconstructed skin and affected their instantaneous elastic response. In contrast, the fibroblasts and the keratinocytes did not significantly influence the dermal equivalent and the mechanical response of the reconstructed skin after a long period of solicitation. During the load-relaxation tests we observed that over a long period of solicitation the mechanical response of the dermal equivalent and the reconstructed skin converged toward the mechanical response of the 3D-scaffold. The viscous properties of the 3D-scaffold thus greatly influence the long-term mechanical responses of the dermal equivalent and the reconstructed skin. It is thus tempting to speculate that a stress-strain applied on the reconstructed tissue would cause an initial response of cellular cytoskeletons while a prolonged effort would activate a larger network of fibres of the extracellular compartment. Finally, reconstructed skin is mainly used as a skin substitute for grafts. In this regard, although already very efficient, it could be improved still further by more faithful mimicking of skin in

vivo, especially with respect to mechanical similarities for a high level of deformation. Consequently, the non-linear viscoelastic and anisotropic behaviour of reconstructed skin could be investigated. The relaxation tests were not adapted for studying the non-linear behaviour of reconstructed skin. In general, studies of the non-linear behaviour of biological tissues are performed using tensile tests. However, because of the current techniques used for manufacturing reconstructed skin samples (circular shape, small diameter, small thickness, etc.) it is difficult to plan tests with large deformation regimes and thus obtain information on their non-linear behaviour. The indentation load-relaxation tests did not allow measuring the anisotropic behaviour of the reconstructed skin. Nevertheless, the indentation load-relaxation tests coupled with rheological models provided partial information on the stratified structure of the biomaterial. Although the tensile tests appeared to be one of the best way of characterising the anisotropic behaviour of the materials (Groves et al., 2013; Ni Annaidh et al., 2012), further information could be obtained by carrying out indentation load-relaxation tests. Indeed a nonaxisymmetric indenter such as a dihedral punch should be used to improve the anisotropic analysis of reconstructed skin based on indentation load-relaxation tests and relaxation tests should also be performed in different directions.

Conflict of interest The authors state no conflict of interests

Acknowledgments We would like to thank BASF® Beauty Care Solution France for providing us with the 3Dscaffold, dermal equivalent and reconstructed skin samples.

Appendices

In this part, the analytical solutions for the spherical indentation load-relaxation tests based on the generalised Maxwell model and the generalised Kelvin-Voigt model were demonstrated. A. The generalised Maxwell model The generalised Maxwell model consists of a spring connected in parallel with n-branches containing a spring and a dashpot connected in series (Fig. 3). Based on the generalised Maxwell model, an analytical solution using a functional equation for the spherical indentation loadrelaxation test is proposed. For the generalised Maxwell model, the stress, σ 0 , in the upper branch is given by:

σ 0 = G0 ε 0

(A.1)

where G0 is the shear modulus and ε 0 the strain (Fig. 3). In all the lower n-branches containing the spring and the dashpot in series, the stress in the ith branch, σ i , is the same in the spring and in the dashpot, therefore:

σ i = Gi ε ia = ηi εib (A.2) th a where Gi is the shear modulus of the spring in the i branch, ε i the spring strain in the ith branch, ηi the dashpot viscosity in the ith branch and εib the dashpot strain rate in the ith branch. The total strain, ε i , in the ith branch is the sum of the spring strain, ε ia , and the dashpot strain, ε ib , as: ε i = ε ia + ε ib (A.3) The time derivative of equation (A.3) is written as: σ σ εi = εia + εib = i + i (A.4) Gi ηi where εi is the total strain rate, εia the spring’s strain rate, σ i the total stress in the ith branch and σi the total stress rate in the ith branch. The total strain, ε , in the generalised Maxwell model is the same in all the branches, thus the total strain is: ε = ε 0 = ε i Using the Laplace transform, equation (A.1) and equation (A.4) can be rewritten in the Laplace domain as in (Gefen and Haberman 2007):

εˆ ( s ) =

σˆ 0 ( s ) G0

(A.5)

⎛ 1 1 ⎞ + ⎟ σˆ i ( s ) ⎝ Gi ηi s ⎠

εˆ ( s ) = ⎜

(A.6)

where εˆ ( s ) , σˆ 0 ( s ) and σˆ i ( s ) are the Laplace transform of ε , σ 0 and σ i , respectively. The Laplace transform variable is noted s. The total stress, σ , in the n-branch generalised Maxwell model is the sum of the stress in the n

upper branch and in the lower n-branches and thus the total stress is: σ = σ 0 + ∑ σ i . In the i =1

Laplace domain the total stress is written as: n

σˆ ( s ) = σˆ 0 ( s ) + ∑ σˆ i ( s )

(A.7)

i =1

Using equations (A.5) and (A.6) in equation (A.7), we obtain: −1

⎛ 1 1 ⎞ σˆ ( s ) = G0εˆ ( s ) + ∑ ⎜ + ⎟ εˆ ( s ) ηi s ⎠ i =1 ⎝ Gi n

Equation (A.8) can still be written as follows: σˆ ( s ) = Gˆ ( s ) εˆ ( s ) n ⎛ 1 1 ⎞ with Gˆ ( s ) = G0 + ∑ ⎜ + ⎟ ηi s ⎠ i =1 ⎝ Gi

(A.8)

(A.9)

−1

being the shear modulus in the Laplace domain.

The effect of the contact geometry is now considered. For indentation with a spherical indenter, the normal load, Fz ( t ) , is related to the penetration depth, δ ( t ) , by (Johnson 2001): 3 4 12 R δ ( t ) 2 E* ( t ) (A.10) 3 Where R is the radius of curvature of the indenter and E * ( t ) is the reduced Young’s modulus of the material tested. Subsequently, to simplify writing, the reduced Young’s modulus will be 2 1 −ν m2 1 1 −ν ind + with Eind , ν ind and Em , ν m being the respective noted E * . E * is defined as: * = E Eind Em Young’s modulus and Poisson’s ratio of the indenter and the indented biomaterial. In this study, since Eind (PTFE indenter)  Em (biomaterials: 3D-scaffold, 3D-Scaffold+dermal equivalent,

Fz ( t ) =

reconstructed skin), it can be assumed that E * is the reduced Young’s modulus of the indented 1 1 −ν m2 biomaterial, * ≈ . E Em The force-penetration depth relationship for the spherical indentation of an elastic material can be converted to a Laplace transformed viscoelastic force-penetration depth relationship as in (Nixon, 1965):

3 4 1 Fˆz ( s ) = R 2 δˆ ( s ) 2 Eˆ * ( s ) 3

(A.11)

When a sample is subjected to a suddenly applied constant displacement or step displacement, the displacement boundary condition in the contact area can be expressed by: δ (t ) = δ0 H (t ) (A.12) where δ 0 is a constant and H ( t ) is the Heaviside unit step function. The Laplace Transform of equation (A.12) is: 1 δˆ ( s ) = δ 0 Hˆ ( s ) = δ 0 s (A.13) By substituting equation (A.13) in equation (A.11), the force-penetration relationship in Laplace domain becomes: 3

4 1 δ 2 Fˆ ( s ) = R 2 0 Eˆ * ( s ) 3 s

(A.14)

The reduced Young’s modulus E * is defined in terms of shear modulus, G, and bulk modulus, K, as: 3K + G E * = 4G (A.15) 3K + 4G A frequent assumption is that the time-dependent deformation is restricted to the shear components of the deformation, such that the bulk modulus is assumed to be time-independent (Cheng and Cheng, 2005; Oyen, 2006). The simplification is equivalent to assuming a timedependent Young’s modulus and a time-independent Poisson’s ratio (Johnson, 2001). Thus, in the Laplace domain, the reduced Young’s modulus is given by: 3K + Gˆ ( s ) Eˆ * ( s ) = 4 Gˆ ( s ) 3K + 4 Gˆ ( s ) By replacing equation (A.16) in equation (A.14), we obtain: 3 3K + Gˆ ( s ) 4 12 δ 0 2 ˆ ˆ Fz ( s ) = R 4G (s) s 3 3K + 4Gˆ ( s ) Equation (A.17) can still be written using equation (A.8) as follows:

(A.16)

(A.17)

−1 n ⎛ ⎛ 1 1 ⎞ ⎞ 3K + ⎜ G0 + ∑ ⎜ + ⎟ ⎟ 3 −1 ⎜ n ηi s ⎠ ⎟ ⎞ i =1 ⎝ Gi 1 δ 2 ⎛ ⎛ ⎞ 16 1 1 ⎝ ⎠ Fˆz ( s ) = R 2 0 ⎜ G0 + ∑ ⎜ + ⎟ ⎟⎟ −1 n 3 η s ⎜ G s ⎛ i =1 ⎝ i i ⎠ ⎠ ⎛ 1 1 ⎞ ⎞ ⎝ 3K + 4 ⎜ G0 + ∑ ⎜ + ⎟ ⎟ ⎜ ηi s ⎠ ⎟ i =1 ⎝ Gi ⎝ ⎠

(A.18)

The viscoelastic solution of the force-penetration relationship can be obtained by implementing the inverse Laplace transform of Fˆ ( s ) as in (Nixon, 1965): G n − it ⎞ ⎛ ηi 3K + ⎜ G0 + ∑ Gi e ⎟ G ⎜ ⎟ n − it ⎞ 1 = i 16 12 3 2 ⎛ ⎝ ⎠ Fz ( t ) = R δ 0 ⎜ G0 + ∑ Gi e ηi ⎟ Gi ⎜ ⎟ n − t 3 ⎛ i =1 ⎝ ⎠ 3K + 4 ⎜ G + G e ηi ⎞⎟ i ⎜ 0 ∑ ⎟ i =1 ⎝ ⎠

(A.19)

B. The generalised Kelvin-Voigt model The generalised Kelvin-Voigt model consists of a spring connected in series with n-branches, containing a spring and a dashpot connected in parallel (Fig. 4). Based on the generalised Kelvin-Voigt model, an analytical solution using a functional equation is proposed for the spherical indentation load-relaxation test as: For a generalised Kelvin-Voigt model, the stress, σ 0 , in the first spring is given by:

σ 0 = G0 ε 0

(B.1)

where G0 is the shear modulus and ε 0 the strain (Fig. 4) In all the n-branches containing the spring and the dashpot in parallel, the stress in the ith branch, σ i , is the sum of the stresses of the spring and the dashpot; therefore: σ i = Gi ε i + ηi εi (B.2) where Gi and ηi are the shear modulus of the spring and the dashpot viscosity, respectively, in the ith branch, while ε i and εi are the strains of the spring and the dashpot strain rate, respectively, in the ith branch. The total stress, σ , in the generalised Kelvin-Voigt model is the same in all the branches, thus the total stress is: σ = σ 0 = σ i Using the Laplace transform, equation (B.1) and equation (B.2) can be rewritten in the Laplace domain as in (Gefen and Haberman, 2007):

σˆ 0 ( s ) = G0 εˆ ( s )

(B.3)

σˆ i ( s ) = ( Gi + ηi s ) εˆ ( s )

(B.4)

where εˆ ( s ) , σˆ 0 ( s ) and σˆ i ( s ) are the Laplace transform of ε , σ 0 and σ i , respectively. The Laplace transform variable is noted s. The total strain, ε , in the n-branch generalised Kelvin-Voigt model is the sum of the strain in n

each branch, thus the total strain is: ε = ε 0 + ∑ ε i . In the Laplace domain the total strain is i =1

written as: n

εˆ ( s ) = εˆ0 ( s ) + ∑ εˆi ( s )

(B.5)

i =1

Using equations (B.3) and (B.4) in equation (B.5), we obtain: σˆ ( s ) n σˆ ( s ) +∑ εˆ ( s ) = G0 i =1 Gi + ηi s n ⎛ 1 ⎞ 1 ˆ ε (s) = ⎜ + ∑ ⎟ σˆ ( s ) G G η s + = 1 i i i ⎝ 0 ⎠

(B.6)

Equation (B.6) can still be written as follows: σˆ ( s ) = Gˆ ( s ) εˆ ( s )

(B.7)

−1

n ⎛ 1 ⎞ 1 with Gˆ ( s ) = ⎜ +∑ ⎟ the shear modulus in the Laplace domain. ⎝ G0 i =1 Gi + ηi s ⎠

Using the previous results and assumptions obtained from the generalised Maxwell model, equation (A.17) is written, in the case of the Kelvin-Voigt model, as: −1

3

n ⎞ 16 1 δ 2 ⎛ 1 1 +∑ Fˆz ( s ) = R 2 0 ⎜ ⎟ s ⎝ G0 i =1 Gi + ηi s ⎠ 3

−1

n ⎛ 1 ⎞ 1 3K + ⎜ +∑ ⎟ ⎝ G0 i =1 Gi + ηi s ⎠ −1 n ⎛ 1 ⎞ 1 3K + 4 ⎜ +∑ ⎟ ⎝ G0 i =1 Gi + ηi s ⎠

(B.8)

The viscoelastic solution of the force-penetration relationship can be obtained by implementing the inverse Laplace transform of Fˆ ( s ) (Nixon, 1965): −1

G n ⎛ 1 − i t ⎞⎞ 1 ⎛ ηi + ∑ ⎜1 − e ⎟⎟ −1 3K + ⎜ Gi ⎜ ⎟⎟ ⎜ n G G ⎛ ⎞ t − ⎛ ⎞ i = 1 3 1 i 0 16 1 1 ⎝ ⎠⎠ ⎝ + ∑ ⎜ 1 − e ηi ⎟ ⎟ Fz ( t ) = R 2 δ 0 2 ⎜ −1 ⎜ ⎟ G ⎜ ⎟ 3 G G n ⎛ 1 − i t ⎞⎞ ⎠⎠ ⎝ 0 i =1 i ⎝ 1 ⎛ + 3K + 4 ⎜ ⎜ 1 − e ηi ⎟ ⎟ ⎜ ⎟⎟ ⎜ G0 ∑ i =1 Gi ⎝ ⎠⎠ ⎝

(B.9)

References

Anssari-Benam, A., Bader, D.L., Screen, H.R.C., 2011. Anisotropic time-dependant behaviour of the aortic valve. J. Mech. Behav. Biomed. Mater. 4, 1603-1610. Black, A.F., Bouez, C., Perrier, E., Schlotmann, K., Chapuis, F., Damour, O., 2005. Optimization and characterization of an engineered human skin equivalent. Tissue Eng. 11, 723733. Boyer, G., Laquieze, L., Le Bot, A., Laquieze, S., Zahouani, H., 2009. Dynamic indentation on human skin in vivo: ageing effects. Skin Res. Technol. 15, 55-67. Brown, R.A., 2013. In the beginning there were soft collagen-cell gels: towards better 3D connective tissue models? Exp. Cell. Res. 319, 2460-2469. Cheng, Y.T., Cheng, C.M., 2005. Relationships between initial unloading slope, contact depth, and mechanical properties for conical indentation in linear viscoelastic solids. J. Mater. Res. 20, 1046-1053. Cheung, H.Y., Lau, K.T., Lu, T.P., Hui, D., 2007. A critical review on polymer-based bioengineered materials for scaffold development. Compos. Part B-Eng. 38, 291-300. Chin, L.K., Calabro, A., Walker, E., Derwin, K.A., 2012. Mechanical properties of tyramine substituted-hyaluronan enriched fascia extracellular matrix. J. Biomed. Mater. Res. A 100, 786793. Dantzig, G.B., 1963. Linear programming and extensions, Princeton University Press, Princeton. Darling, E.M., Topel, M., Zauscher, S., Vail, T.P., Guilak, F., 2008. Viscoelastic properties of human mesenchymally-derived stem cells and primary osteoblasts, chondrocytes and adipocytes. J. Biomech. 41, 454-464. Drury, J.L., Mooney, D.J., 2003. Hydrogels for tissue engineering: scaffold design variables and applications. Biomaterials 24, 4337-4351. Ebersole, G.C., Anderson, P.M., Powell, H.M., 2010. Epidermal differentiation governs engineered skin biomechanics. J. Biomech. 43, 3183-3190. Echinard, C., 1998. Under the epidermis, the dermis, or how to interpret skin cultures. Ann. Chir. Plast. Esth. 43, 197-205.

Elias, P.Z., Spector, M., 2012. Viscoelastic characterization of rat cerebral cortex and type I collagen scaffolds for central nervous system tissue engineering. J. Mech. Behav. Biomed. Mater. 12, 63-73. Fallqvist, B., Kroon, M., 2013. A chemo-mechanical constitutive model for transiently crosslinked actin networks and a theoretical assessment of their viscoelastic behaviour. Biomech. Model. Mechanobiol. 12, 373-382. Flynn, C., Taberner, A., Nielsen, P., 2011. Estimating material parameters of structurally based constitutive relation for skin mechanics. Ann. Biomed. Eng. 39, 1935-1946. Flynn, C., Taberner, A.J., Nielsen, P.M.F., Fels, S., 2013. Simulating the three-dimensional deformation of in vivo facial skin. J. Mech. Behav. Biomed. Mater. 28, 484-94. Forgacs, G., Foty, R.A., Shafrir, Y., Steinberg, M.S., 1998. Viscoelastic properties of living embryonic tissues: a quantitative study. Biophys. J. 74,2227–2234 Freyman, T.M., Yannas, I.V., Yokoo, R., Gibson, L.J., 2001. Fibroblast contraction of a collagen-GAG matrix. Biomaterials 22, 2883-2891. Fung, Y.C., 1993. Biomechanics: Mechanical Properties of Living Tissues, Springer-Verlag New York. Gardel, M.L., Shin, J.H., MacKintosh, F.C., Mahadevan, L., Matsudaira, P., Weitz, D.A., 2004. Elastic Behavior of Cross-Linked and Bundled Actin Networks. Science 304, 1301-1305. Geerligs, M., Peters, G.W.M., Ackermans, P.A.J., Oomens, C.W.J., Baaijens, F.P.T., 2008. Linear viscoelastic behavior of subcutaneous adipose tissue. Biorheol. 45, 677-688. Gefen, A., Haberman, E., 2007. Viscoelastic properties of ovine adipose tissue covering the gluteus muscles. J. Biomech. Eng. 129, 924-930. Gerhardt, L.C., Schmidt, J., Sanz-Herrera, J.A., Baaijens, F.P.T., Ansari, T., Peters, G.W.M., Oomens, C.W.J., 2012. A novel method for visualising and quantifying through-plane skin layer deformations. J. Mech Behav. Biomed. Mater. 14, 199-207. Gras, L.L., Mitton, D., Viot, P., Laporte, S., 2013. Viscoelastic properties of the human sternocleidomastoideus muscle of aged women in relaxation. J. Mech. Behav. Biomed. Mater. 27, 77-83 Groves, R.B., Coulman, S.A., Birchall, J.C., Evans, S.L., 2013. An anisotropic, hyperelastic model for skin: Experimental measurements, finite element modelling and identification of parameters for human and murine skin. J. Mech. Behav. Mater. 18, 167-180. Harley, B.A., Leung, J.H., Silva, E.C., Gibson, L.J., 2007. Mechanical characterization of collagen-glycosaminoglycan scaffolds. Acta Biomater. 3, 463-474.

Hidaka, K., Moine, L., Collin, G., Labarre, D., Grossiord, J.L., Huang, N., Osuga, K., Wada, S., Laurent, A., 2011. Elasticity and viscoelasticity of embolization microspheres. J. Mech. Behav. Biomed. Mater. 4, 2161-2167. Jachowicz, J., McMullen, R., Prettypaul, D., 2007. Indentometric analysis of in vivo skin and comparison with artificial skin models. Skin Res. Technol. 13, 299-309. Johnson, K.L., 2001. Contact mechanics, Cambridge University Press, Cambridge. Jor, J., Nash, M., Nielsen, P., Hunter, P., 2011. Estimating material parameters of a structurally based constitutive relation for skin mechanics. Biomech. Model. Mechanobiol. 10, 767-778. Kanungo, B.P., Gibson, L.J., 2010. Density–property glycosaminoglycan scaffolds. Acta Biomater. 6, 344-353.

relationships

in

collagen–

Knapp, D.M., Barocas, V.H., Moon, A.G., Yoo, K., Petzold, L.R., Tranquillo, R.T., 1997. Rheology of reconstituted type I collagen gel in confined compression. J. Rheol. 41, 971-993. Lamers, E., Van Kempen, T.H.S., Baaijens, F.P.T., Peters, G.W.M., Oomens, C.W.J., 2013. Large amplitude oscillatory shear properties of human skin. J. Mech. Behav. Biomed. Mater. 28, 462-470. Ledoux, W.R., Blevins, J.J., 2007. The compressive material properties of the plantar soft tissue. J. Biomech. 40, 2975-2981. Lee, E.H., Radok, J.R., 1960. The contact problem for viscoelastic bodies. J. Appl. Mech. 27, 438-444. Lee, E.H., 1955. Stress analysis in viscoelastic bodies. Q. Appl. Math. 13, 183-190. Limbert, G., 2011. A mesostructurally based anisotropic continuum model for biological soft tissues decoupled invariant formulation. J. Mech. Behav. Biomed. 4, 1637-1657. Lulevich, V., Yang, H.Y., Isseroff, R.R., Liu, G.Y., 2010. Single cell mechanics of keratinocyte cells. Ultramicroscopy 110, 1435-1442. Mano, J.F., Silva, G.A., Azevedo, H.S., Malafaya, P.B., Sousa, R.A., Silva, S.S., Boesel, L.F., Oliveira, J.M., Santos, T.C., Marques, A.P., Neves, N.M., Reis, R.L., 2007. Natural origin biodegradable systems in tissue engineering and regenerative medicine: present status and some moving trends. J. R. Soc. Interface 4, 999-1030. Monteiro-Riviere, N.A., Inman, A.O., Snider, T.H., Blank, J.A., Hobson, D.W., 1997. Comparison of an in vitro skin model to normal human skin for dermatological research. Microsc. Res. Techniq. 37, 172-179.

Netzlaff, F., Lehr, C.M., Wertz, P.W., Schaefer, U.F., 2005. The human epidermis models EpiSkin, SkinEthic and EpiDerm: an evaluation of morphology and their suitability for testing phototoxicity, irritancy, corrosivity, and substance transport. Eur. J. Pharm. Biopharm. 60, 167178. Nelder, J.A., Mead, R., 1965. A Simplex Method for Function Minimization. The Comp. J. 7, 308-313. Newman, S., Cloître, M., Allain, C., Forgacs, G., Beysens, D., 1997. Viscosity and elasticity during collagen assembly in vitro: relevance to matrix-driven translocation. Biopolymers 41, 337-347. Ni Annaidh, A.N., Bruyere, K., Destrade, M., Gilchrist, M.D., Maurini, C., Ottenio, M., Saccomandi, G., 2012. Automated estimation of collagen fibre dispersion in the dermis and its contribution to the anisotropic behaviour of skin. Ann. Biomed. Eng. 40, 1666-1678. Nixon, F.E., 1965. Handbook of Laplace transformation, fundamentals, applications, tables, and examples, 2nd eds Prentice Hall, New Jersey. Oyen, M.L., 2006. Analytical techniques for indentation of viscoelastic materials. Philos. Mag. 86, 5625-5641. Pailler-Mattei, C., Bec, S., Zahouani, H., 2008. In vivo measurements of the elastic mechanical properties of human skin by indentation tests. Med. Eng. Phys. 30, 599-606. Pailler-Mattei, C., Nicoli, S., Pirot, F., Vargiolu, R., Zahouani, H., 2009. A new approach to describe the skin surface physical properties in vivo. Colloids Surf. B 68, 200-206. Pailler-Mattei, C., Guerret-Piécourt, C., Zahouani, H., Nicoli, S., 2011. Interpretation of the human skin biotribological behaviour after tape stripping. J. R. Soc. Interface 8, 934-941. Pailler-Mattei, C., Debret, R., Vargiolu, R., Sommer, P., Zahouani, H., 2013. In vivo skin biophysical behaviour and surface topography as a function of ageing. J. Mech. Behav. Biomed. 28, 474-483. Pan, Y., Xiong, D., 2010. Stress relaxation behavior of nano-hydroxyapatite reinforced poly(vinyl alcohol) gel composites as biomaterial. J. Mater. Sci. 45, 5495-5501. Patel, P.N., Smith, C.K., Patrick, C.W., 2005. Rheological and recovery properties of poly(ethylene glycol) diacrylate hydrogels and human adipose tissue. J. Biomed. Mater. Res. A 73, 313-319. Peppas, N.A., Hilt, J.Z., Khademhosseini, A., Langer, R., 2006. Hydrogels in biology and medicine: from molecular principles to bionanotechnology. Adv. Mater. 18, 1345-1360. Ponec, M., Boelsma, E., Gibbs, S., Mommaas, M., 2002. Characterization of reconstructed skin models. Skin Pharmacol. Appl. 15, 4-17.

Radok, J.R., 1957. Viscoelastic stress analysis. Q. Appl. Math. 15, 198-202. Screen, H.R.C., 2008. Investigating load relaxation mechanics in tendon. J. Mech. Behav. Biomed. Mater. 1, 51-58. Shahabeddin, L., Berthod, F., Damour, O., Collombel, C., 1990. Characterization of skin reconstructed on a chitosan-cross-linked collagen-glycosaminoglycan matrix. Skin Pharmacol. 3, 107-114. Silver, F.H., Freeman, J.W., DeVore, D., 2001. Viscoelastic properties of human skin and processed dermis. Skin Res. Technol. 7, 18-23. Silver, F.H., Siperko, L.M., Seehra, G.P., 2003. Mechanobiology of force transduction in dermal tissue. Skin Res. Technol. 9, 3-23 Ting, T.C., 1966. The contact stresses between a rigid indenter and a viscoelastic half-space. J. Appl. Mech. 88, 845-854. Vaissiere, G., Chevallay, B., Herbage, D., Damour, O., 2000. Comparative analysis of different collagen-based biomaterials as scaffolds for long-term culture of human fibroblasts. Med. Biol. Eng. Comput. 38, 205-210. Wang, T.W., Spector, M., 2009. Development of hyaluronic acid-based scaffolds for brain tissue engineering. Acta Biomater. 5, 2371-2384. Webba, K., Hitchcocka, R.W., Smeala, R.M., Lib, W., Grayb, S.D., Trescoa, P.A., 2006. Cyclic strain increases fibroblast proliferation, matrix accumulation, and elastic modulus of fibroblastseeded polyurethane constructs. J. Biomech. 39, 1136–1144. Zahouani, H., Pailler-Mattei, C., Sohm, B., Vargiolu, R., Cenizo, V., Debret, R., 2009. Characterization of the mechanical properties of a dermal equivalent compared with human skin in vivo by indentation and static friction tests. Skin Res. Technol. 15, 68-76.

Table 1- Mechanical parameters of the 3D-scaffold obtained with the one-branch

generalised Maxwell model. τ 1 is the relaxation time of the 3D-scaffold defined as τ 1 = All data are indicated as mean ± SD (SD: standard deviation). E1 η1 E0 (kPa) τ 1 (s) (kPa) (kPa.s) 3Dscaffold

3.7±1.2

0.4±0. 10.8±1. 1 9

Ψ (t )

η1

E1

Ψ (t )

t →0

t →∞

(kPa)

(kPa)

4.1±1.1

3.7±1.2

28.2± 2.2

.

Table 2- Mechanical parameters of the dermal equivalent obtained with the one-branch generalised Kelvin-Voigt model. τ1* is the relaxation time of the dermal equivalent defined as

τ 1* =

η1 E1

. All data are indicated as mean ± SD (SD: standard deviation).

E0 (kPa) Dermal equivalent

6.7±1.1

E1 (kPa) 6.4±1.4

η1 (kPa.s) 160±4. 6

τ1* (s) 20.8± 3.8

Ψ (t )

Ψ (t )

t →0

t →∞

(kPa)

(kPa)

6.7±1.1

3.3±0.6

Table 3- Mechanical parameters of the reconstructed skin obtained with the two-branches generalised Kelvin-Voigt model. The parameters E0 , E1 and η1 , have been imposed from the

one-branch Kelvin-Voigt model used for the dermal equivalent: E0 = 6.7 kPa, E1 = 6.4 kPa,

η1 = 160 kPa.s. τ 2* is one of the two relaxation time of the reconstructed skin, defined as η τ 2* = 2 . All data are indicated as mean ± SD (SD: standard deviation). E2

Reconstructe d skin

η2

E2 (kPa)

(kPa.s)

7.6±1.5

18.6±5.9

τ 2* (s) 1.9±0.4

Ψ (t )

t →0

(kPa) 13.4

Ψ ( t ) (kPa)

t →∞

5.1±0.8

(a)

(b)

(c)

(d)

(f)

(e)

50 µm

50 µm

50 µm

Fig. 1- (a) Principle of reconstructed skin preparation. The 3D-scaffold used in this study is a biomaterial with a base of collagen, chitosan and glycosaminoglycan. Dermal fibroblasts were inoculated into the 3D-scaffold substrates during 3 weeks. Finally, keratinocytes were seeded onto the dermal equivalent to obtain the reconstructed skin (after 17 days). (b) and (c) Macroscopic view of 3D-scaffold and inserts in culture-plate. (d), (e) and (f) Hematoxylin-eosin staining of 3D-scaffold alone, dermal equivalent and reconstructed skin. Scale bars = 50 µm.

translation screw for macroscopic approach

displacement table displacement sensor

spherical indenter

Fig. 2 - Schematic representation of the micro-indentation device

0

G0

Gi

i b i

a i

n-branch

Fig. 3- Representation of a generalised n-branch Maxwell model. The generalised Maxwell model consists of a spring (shear modulus G0 ) connected in parallel with n-branch containing a

spring (shear modulus Gi ) and a dashpot (viscosity ηi ) connected in series.

G0

Gi

i 0

i

n-branch Fig. 4- Representation of a generalised n-branch Kelvin-Voigt model. The generalised Kelvin-Voigt model consists of a spring (shear modulus G0 ) connected in serie with n-branch

containing a spring (shear modulus Gi ) and a dashpot (viscosity ηi ) connected in parallel

Normal load (mN)



1

3D-scaffold Dermal equivalent Reconstructed skin

0.8

0.6

0.4

0.2

0

0

50

100

150

200

Time (s)

Fig. 5- Representative load-relaxation curves obtained on the three different tissues (3Dscaffold, dermal equivalent, reconstructed skin). The tests were performed in controlled penetration depth, δ 0 = 50 μ m , with a constant indentation speed δ = 50 μ m.s-1 . The used indenter was a spherical PTFE indenter, with a radius of curvature R=1.6 mm.

Normal load (mN)

0.3

0.2

E0

0.1

3D-scaffold One-branch Maxwell E1

0

0

50

100

1

150

200

Time (s)

Fig. 6- Comparison between the 3D-scaffold experimental load-relaxation curve and the theoretical load-relaxation curve obtained with a one-branch Maxwell model. The phenomenological rheological model corresponding to the generalised one-branch Maxwell model is reported on the graph.

0.6

Dermal equivalent One branch Maxwell One branch Kelvin-Voigt E1

Normal load (mN)

0.5 0.4

E0 0.3

E0

0.2

1

0.1

E1 0

0

50

1

100

150

200

Time (s) Fig. 7- Comparison between the dermal equivalent experimental load-relaxation curve and the theoretical load-relaxation curves obtained with a one-branch Maxwell model and the onebranch Kelvin-Voigt model. The phenomenological rheological models corresponding to the generalised one-branch Maxwell model and the generalised one-branch Kelvin-Voigt model are respectively reported on the graph, under and above the experimental curve. The error criterion, ξ , for the Maxwell model and for the Kelvin-Voigt are respectively 13.7 ± 2.7% and 4.6 ± 1.1%.

1

Reconstructed skin One branch Kelvin-Voigt Two branches Kelvin-Voigt

Normal load (mN)

0.8

0.6

E1

E2

1

2

E0

0.4 E1

0.2

0

E0

1

0

50

100

150

200

Time (s)

Fig. 8- Comparison between the reconstructed skin experimental load-relaxation curve and the theoretical load-relaxation curves obtained with a one-branch Kelvin-Voigt model and a twobranches Kelvin-Voigt model. The phenomenological rheological models corresponding to the generalised one-branch Kelvin-Voigt model and the generalised two-branches Kelvin-Voigt model are reported on the graph under and above the experimental curve respectively.