Visualizing dynamic three-dimensional changes of human reticular dermal collagen under mechanical strain

In clinical practice, plastic surgeons are often faced with large skin defects that are difficult to close primarily. Management of large skin wounds e.g. burns or traumatic lacerations requires knowledge of skin biomechanic properties. Research into skin microstructural adaptation to mechanical deformation has only been performed using static regimes due to technical limitations. Here, we combine uniaxial stretch tests with fast second harmonic generation imaging and we apply this for the first time to investigate dynamic collagen rearrangement in reticular human dermis. Ex vivo human skin from the abdomen and upper thigh was simultaneously uniaxially stretched while either periodically visualizing 3D reorganization, or visualizing 2D changes in real time. We determined collagen alignment via orientation indices and found pronounced variability across samples. Comparing mean orientation indices at the different stages of the stress strain curves (toe, heel, linear) showed a significant increase in collagen alignment during the linear part of the mechanical response. We conclude that fast SHG imaging during uni-axial extension is a promising research tool for future studies on skin biomechanic properties.


Introduction
Management of large skin wounds e.g. burns or traumatic lacerations requires knowledge of skin biomechanic properties. Lack of proper attention for tensionfree closure will result in an increased risk of impaired wound healing and scar hypertrophy. Human skin is characterised by its nonlinear, viscoelastic and anistropic behavior, with vast local and biological variations that are clinically relevant. For example, intraoperative rapid expansion is a surgical technique in which a constant load is applied to skin and stress relaxation takes place, increasing its expandibility and surface area [1]. The underlying mechanism of this phenomenon called creep is thought to be based on alignment of collagen fibers in the direction of the applied load. Furthermore, physiological (aging, Sundamage) or pathological changes (collagen disorders, scarring) result in clinically relevant differences in these biomechanical properties, increasing the importance of insight into the underlying microstructural differences [2][3][4].
Human skin is comprised of two main layers: the epidermis, comprised of compact layers of cells, which acts as a barrier against outside hazards, and the dermis as the central load bearing structure. The dermis is comprised of a dense network of collagen and elastin fibers in a ground substance made up of hyaluronic acid and glycosaminglycans among others [5]. The organization of collagen fibers differs in the deep and superficial layers and divides the dermis in a papillary and reticular layer. The thinner papillary layer is made up of a network of randomly oriented thin bundles. In contrast, the reticular dermis consists of thick haphazardly structured bundles [6]. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
Collagen is a family of protein biopolymers of which type I is the most abundantly present in skin. It is comprised of triple helical fibrils that are assembled into fibers and then into fiber bundles [7].
Research into deformation-induced microstrucural changes of skin has previously been carried out under static conditions, without the possibility to observe the underlying dynamics [8], because a suitable microscopy technique that is sufficiently fast, specific for relevant structures, and of high enough resolution to visualize the biomechanically relevant structures, has been lacking. Furthermore, no technique was available that can be used on unsectioned and unprocessed tissue, which is vital in the dynamic asssessment of microstructural deformation. However, non-linear optics provides the solution here. Specifically, Second Harmonic Generation (SHG) microscopy fulfills all requirements for the dynamic assessment of skin biomechanics. SHG relies on the non-linear susceptibility of dermal proteins to provide a sufficient signal. Proteins like collagen type I and II are not only abundantly present in the human body but are structured on multiple length scales in parallel fashion, increasing the intensity of the generated signal. This allows for imaging of unprocessed, and uncut tissue [7,9].
Research into skin biomechanics has previously focused mainly on porcine and murine skin [7,[9][10][11], as it constitutes an accesible model. However these models fundamentally differ from human skin, biologically and mechanically. Although porcine skin, the most appropriate substitute, has similarities in collagen and elastin content, large difference in tangent moduli exist with human skin [12]. Gallagher et al assessed mechanical properties of human and porcine ex vivo burn models and found porcine skin to be an inadequate substitute [13].
Furthermore, aside from a recent study by Ueda et al [14] in which, after chemical clearing, the entire dermis was imaged during bi-axial deformation, other studies have focused on the papillary dermal layer as its random basket weave pattern of fine collagen fibers provides an accessible reference point for subsequent comparison. However, this layer is not the main load bearing structure of skin and research into the reticular dermal layer has been lacking.
This study is original because for the first time, we perform dynamic stretch measurements combined with fast SHG imaging to investigate dynamic collagen rearrangement in untreated reticular dermis of human ex-vivo samples. Our hypothesis is that with this combination, increased insight into the main load bearing structure of human skin when subjected to controlled mechanical tension,will be generated, which is of significance for reconstructive surgery.

Human skin samples
Ten ex vivo skin samples were obtained from either healthy patients (n = 5) or cadavers (n = 5). Healthy skin involved surgical waste tissue from patients undergoing abdominal dermolipectomy surgery at the Red Cross Hospital (Beverwijk, the Netherlands). Under national guidelines, surgical waste tissue can be used for scientific purposes when patients have not objected after oral informed consent (The Act on the Medical Treatment agreement, article 7:467 BW). Due to the scarcity of elective surgeries taking place during the COVID-19 pandemic, an alternative source of ex vivo skin tissue was found in the form of cadaver skin from individuals donating their body to science. This was approved separately by the medical ethical committee of the Amsterdam UMC -location VUmc (Biobank reference number 2017.098).
Split thickness skin (0.8 mm thick) was taken from full thickness abdominoplasty samples or directly from cadaver upper thighs, using a commercial D80 cordless dermatome (Humeca, Borne, the Netherlands). The use of split thickness skin by means of a dermatome had several advantages: equal thickness across the sample resulting in homogeneous in-plane tension. The split skin graft consists of the epidermis and the reticular dermal layer, and allowed for examination of the reticular dermal layer at a well determined depth. Ex-vivo skin was then stored in RPMI 1640 medium (Sigma Aldrich) and stored at 4°C. Samples were used within 48 h after graft taking. Using a scalpel, split-thickness skin was manually cut into 4 × 2 cm dog-bone shaped [15] samples before each experiment.

Biomechanical assessment
To apply uniaxial mechanical strain, a stretching device was developed using two DSZY-1 linear actuators (Drive-Systems Europe Ltd). The actuators were bolted on a custom-made detachable frame, with a close fit relative to the motorized microscopy stage. Fixed to the actuator arms were two in line load cells (210 series, Richmond Industries) with a maximum capacity of 1000 Newtons. Their force readout was sampled at a rate of 2 Hz. Connected to the load cells were two stainless steel gripping arms designed to firmly hold the skin samples in place. The gripping arms contained six vertical puncture holes through which the outer ends of the samples were additionally fixed by piercing them with six steel needles. Proprietary software (Electronics workshop, VU Amsterdam) running on a laptop computer connected to the stretching device was used to impose stretching regimens on the samples and for simultaneous data acquisition.

Microscopic assessment of collagen rearrangement
Second Harmonic Generation (SHG) images were acquired using a custom-built multiphoton microscope (figure 1). The laser source was a femtosecond pulsed laser (FSP-03, Litilit, Vilnius Lithuania) (85fs, 13.3 MHz) with a center wavelength of 1050 nm. An acousto-optic modulator (AOM) (MT250-A0.5-1064, AA Opto Electronic, Orsay, France) adjusted the power on the sample, resulting in an average power of 50mW. A quarter-wave plate (Standa, Italy) placed after the AOM ensured circular polarization. In plane imaging involved scanning the laser beam in the x and y direction using galvanometer mirrors. An oil-immersion objective with a high NA (Nikon S Fluor 40x/1.30, Nikon) was used to focus the incident light on the sample. Epi-detected SHG signal was filtered using ausing a dichroic mirror dichroic mirror (FF872, DM1), a bandpass filter at 525 nm and was detected using a photomultiplier tube (H16201-40, Hamamatsu, Japan). To facilitate the use of an oil immersion objective, an aluminum objective cover, adjustable in height, separated the objective and 0.17 mm coverslip at the working distance of the objective. To achieve three-dimensional imaging (through 2D stacking), a tunable lens (Optotune Switzerland AG) was used, allowing for 30 μm of axial imaging at 1 μm increments. Following ROIs during stretching involved moving the skin stretching device using the motorized microscopy stage.
Microscopy data was recorded using in-house developed LabView software. Single images were 1000 ×1000 pixels in size, with a field of view of 500 × 500 μm. Image stacks were 30 images in depth with 1 μm slice spacing. Raw image stacks were saved as 32-bit tiff files and were falsed-colored in red. Lower resultion images obtained during two-dimensional imaging were 500 × 500 pixels in size with a field of 500 × 500 μm.

Assessment of collagen alignment
To measure the three-dimensional orientation of collagen fibers, we used the weighted vector summation algorithm created by Lui et al [11], which has previously been applied to the murine dermis. To define fiber orientation in three-dimensional space, the azimuthal (j) and polar (θ) angles are determined. These angles lie in the plane of optical sectioning for θ, and perpendicular to this plane for j ( figure 1(d)). Briefly, a voxel of predefined size is first projected onto a two-dimensional plane and the orientation of the center pixel in this plane is determined by the addition of all weighted vectors. The azimuthal angle θ is determined by projecting fibers onto the xy plane. The azimuthal angle j is determined by first establishing β and γ, which are obtained by projecting the fibers to xz and yz plane. We used a voxel window of 11 × 11 × 11 pixels, and a sampling ratio between the xy-and z-dimension of 2, corresponding with a lateral and axial resolution of 0.5 and 1 μm/pixel, respectively.
To quantify overall fiber orientation from these two angles, the three-dimensional orientation index is calculated based on directional statistics [9,16]. The formulas used to calculate the orientation index have been described before and are provided as supplementary material S1. The value of this three-dimensional orientation index ranges from zero to one, with zero meaning fibers are strongly disorganized and one meaning fibers are perfectly parallel.

Stretching protocol
Samples were placed with the reticular dermis facing the objective in the skin stretching device. Ambient temperature and humidity were kept stable at 20°C and 45% respectively during the experiment in the enclosure of the microscope. The humidity of the sample was constant as the duration of the experiment was kept short and the protective objective cap and cover slip were in contact with the dermal side of the sample, reducing water evaporation.
Without preconditioning, samples were subjected to uniaxial deformation to a maximum of 160 percent of their original length at a rate of 0.5 millimeter per second to obtain stress strain relationships. The maximum of 160% was chosen to avoid failure which could damage the microscope.
Subsequently, this stretching process was repeated at the same strain rate, but with intermittent pauzes to allow for three-dimensional imaging.

Data analysis
Applied stretch ratios (λ) were calculated by dividing sample length (L) by the intitial sample length (L0).
Force at the relevant moments of deformation was transformed into nominal stress by dividing by the original cross-sectional area of the samples.
For every sample, stress strain curves were made by plotting nominal stress against the applied strain ratio. Toe, heel, and linear parts of the stress-strain curve were identified manually for all ten curves. The corresponding microscopy datapoints (orientation index) were binned to these three stages of the stress strain curve, as a direct correlation between mechanical data and microscopy data was not possible because data for some samples crossed multiple separate ROIs. Differences in mean orienation index were assessed non-parametrically using Dunn's multiple comparison test.
Statistical analysis and graph preparation were done using PRISM 9 for Mac.

Collagen microstructural assessment
To assess the real time imaging capabilities of the experimental setup described above, representative regions of interest (ROI) were followed throughout the stretching process by imaging at a low resolution, at a rate of 3 Hz. Time series image stacks of reticular dermal tissue were acquired, depicting stretch-and relaxation-based reorganization of collagen on the microstructural scale (figure 2). Reticular dermal collagen was found to be in its characteristic bundle organization, often arranged in thick bundles of varying diameter. The supramolecular structure of the collagen matrix, which seemed quasi random at lower strain rates, swiftly oriented in the direction of uniaxial stretch at higher strain rates ( figure 2(a), (b)). Braidlike structures unraveled with increasing stretch (figure 2(c), (d)), with individual fibers going from crimped to a straight structure. The uncrimped bundles returned to their crimped state with the sample's return to a more relaxed position (figure 2(e), (f), Supplementary video 1), suggesting that we are below the skin's plastic threshold, its yield strength, i.e. within the physiological boundaries of reversible (and time-dependent) deformation.
For analysis of the full collagen response to stretching, we considered three-dimensional images at higher resolution necessary, therefore in a consecutive set of experiments stretching was halted at regular intervals, for high resolution 3D imaging for approximately 5 min, during which some relaxation of the tissue inevitably took place. No preconditioning was performed, to more accurately mimick the physiological response. Thus for every sample, a series of threedimensional images was obtained at discrete intervals.
For all ten samples, three-dimensional images were obtained at regularly spaced intervals in the stretching regimen which, contrary to the two-dimensional data, allowed for quantitative synthesis. Representative results for a single ROI are shown in figure 3. Representative sections from three-dimensional z-stacks are plotted, showing gradual alignment in the direction of stretch (horizontal axis) and uncrimping of individual fibers.
Orientation analysis using Liu's method is also shown in figure 3. For both θ and j, fiber orientation is visually depicted for the two-dimensional sections using a color-coded angle distribution. For θ, zero and 180 degrees correspond to the axis of deformation, while 90 degrees corresponds to in-plane orientation for j. Quantitative analysis of fiber re-organization at increasing stretch ratios is shown in the bottom panel, in the form of frequency plots for all possible angles of j and θ. In θ direction, at low stretch ratio a partial orientation is present already, as evidenced by the two peaks in the θ angle distribution. The fibers progressively align along the stretch direction with increasing stretch ratio and start to form a single peak around zero and 180 degrees. From the j orientation analysis it appears that the main orientation was in plane already but is further enhanced at high stretch ratios.

Correlation with tissue mechanics
Stress-strain curves were obtained for all ten samples. In general, these were comparable to previously described results of human skin biomechanics [17], meaning they showed the typical J shape characteristic of the non-linear response, as can be seen in figure 4 for an illustrative example.
The stress strain curve is plotted alongside a curve with the orientation indices of the same sample at given intervals. For this sample, the exact same ROI was tracked during the course of stretching, allowing for a direct correlation of the skin's collagen rearrangement and its mechanical response to stretch. In this case, collagen reorganization accurately followed the stress strain curve. This implies that stretching of the skin is enabled by reorientation of the collagen fibers up to an orientation index of 0.3, after which a small further stretch occurs without increased alignment. The correlation between stress-strain and orientation index is in agreement with the results of Lynch et al on murine tissue, but is in contradiction with the classical notion that associates the heel region to an alignment of the fibers and the linear one to the stretching of the fibers [10].
However, collagen orientation indices as a function of strain showed pronounced variability across samples and ROIs, with only few cases where the orientation index followed the stress strain curve as in figure 4. Furthermore, following individual ROIs proved challenging and thus for some samples, the calculated orientation index datapoints covered multiple ROIs ( figure 5). Nonetheless, comparing mean orientation indices at the different stages of the stress strain curves (toe, heel, linear) showed that a significant increase in collagen alignment happens only during the linear part of the mechanical response.
Concurrently, samples also showed vast differences in extensibility as shown by the large variance in onset of the heel phase, a phenomenon also seen in the literature [17] (figure 5). Variability in skin biomechanics has been attributed to areas of the body, age, and sample orientation (anisotropy).
Assessing these curves, a dichotomy between the two skin sources was apparent; skin surgically excised for cosmetic purposes showed greater extensibility than skin obtained from cadavers, providing at least a partial explanation of the identified variability.
The highly anisotropic character of skin is clinically relevant in that extensibility and skin tension have preferred directionality, depending on its place on the body. This anisotropy manifests itself as skin tension lines, caused by either preferential arrangment of collagen or tension effect of underlying musculature. Over the years, starting with Langer, different terms like static, dynamic, or resting skin tension lines were coined to try and encapsulate the effect that, depending on pronounced directional variation, human skin allows for optimal placement of surgical incisions [18]. The microstructural mechanism behind these skin tension lines was investigated for one abdominal skin sample by cutting two separete samples perpendicular to each other. One sample, along resting skin tension lines, and one perpendicular to these lines, were evaluated both biomechanically and microscopically ( figure 6). The stress strain curves show very different behaviour. The sample cut parallel to the skin lines shows a high degree of pre-orientation at low stretch ratio ( figure 6(B), notably in the direction of the applied stretch, see figure 6(D)) and a smaller extensibility (figure 6(B)) compared to the sample cut perpendicular to the skin lines. Also, the parallel sample exhibited reorganization to higher degree of parallelism than the perpendicular sample.

Discussion
In this article we demonstrate the feasibility of real time dynamic visualization of stretched human dermal collagen, using a robust experimental setup and a fast image processing algorithm by Liu et al [9]. We provide a first overview of our work on simultaneous mechanical testing of human skin and three-dimensional microstructural assessment of collagen rearrangement in the reticular dermal layer.
The tissue samples show a large variation in preexistent orientation, in onset of the heel phase in the stress strain curve, and in the extensibility. A first observation is that in-plane (θ) reorganization in reticular dermal collagen occurs in the direction of deformation. This has been previously described for the papillary dermis, but not for the reticular dermis. While these fibers also show a clear reorganization in the plane perpendicular to optical sectioning (j) towards the axes of deformation, the difference was not as strong as the reorganization in the plane of optical sectioning. This can, in part, be explained by preexisting organization in this plane. However limited z-stack depth of the threedimensional images may also play a role.
The statistical analysis of the collagen reorganization and stress strain curves in all 10 samples allowed us to conclude that collagen reorganization happens mainly in the linear region of the stress strain curve. Although this is in accordance with previous results by Lynch et al on murine tissue, they are incongruent with previously held beliefs of an initial aligment in the heel region, and subsequent stretching of fibers in the linear region [10].
Mechanical responses of all samples showed the well-reported J-shaped curve and three-dimensional reorganization was shown to take place mainly in the linear part of the deformation response, for both surgically excised skin and cadaver skin. However absolute differences in the stress strain curves of the two groups were large and did not consistently follow the neat correlation seen in figure 4. Also in figure 6, the correlation is less pronounced than in figure 4, for instance. However we consider the variation in itself an important result, warranting publication. This large variability is also reported in literature and can be attributed to, for example, age disparity or differences in abdominal skin compared to upper thigh skin, or to local heterogeneity. Although we tested a potential grouping in two sets for the orientation indices as a function of stretch for the abdominal and cadaver groups, this turned out not to be the case. Therefore, we combined the orientation index data into a single figure, to enable on one hand an overall conclusion, but also the large variation in the data can be appreciated. Due to the origin of the skin tissue, demographic information was limited to the cadaver skin tissues which exhibited a stark uniformity in age. Although to be expected, a more heterogenous demographic distribution in skin samples would allow for a more thorough assessment of the effect of different characteristics on the skin's microstructural response to deformation. Persons undergoing dermolipectomy are typically multiple decades younger than individuals donating their bodies to science, however no exact demographic information could be obtained due to regulatory constraints. Not surprisingly, we found that surgically excised skin showed more extensibility than cadaver skin.
Our results are in concurrence with recent work by Ueda et al [14]. The authors performed biaxial  stretching experiments on full thickness skin, after chemical clearing, finding multiple modes in fiber orientation distribution. We similarly found two modes in the three-dimensional fiber orientation distribution (figure 3), underscoring the hypothesis by Ueda et al concerning a major distribution interspersed with a smaller number of differently orientated fibers. Although Ueda et al describe impressive penetration depth in visualizing full thickness skin, imaging tissue dynamics in untreated samples will increase clinically relevant insight into the relation between collagen micro-orientation and skin macroscopic features such as extensibility and elasticity. Imaging collagen dynamics in real-time (as depicted in figure 2) and in high resolution can potentially increase our understanding of stretch-induced changes, as well as the supramolecular changes involved with relaxtion related processes like creep. The relevance of sample orientation was examined in the context of relaxed skin tension lines, which are crucial in planning the orientation of surgical incisions. Here, a clear difference is exhibited between the mechanical respons from skin parallel to these lines compared to perpendicularly orientated skin. Although only shown in one experiment, our data suggest that this marked anisotropy, a well-known characteristic of human skin, has an underlying origin in three-dimensional collagen orientation and reorganization response. Therefore, it will be of interest to study if the observed variations in onset of the heel phase in the stress strain curve and in the extensibility can be related to preexistent collagen fiber orientation, enabling a deeper insight into the relation between collagen microstructure and extensibility and elasticty of skin.
We therefore aim in the future to extend our experiments to investigate the correlation of the skin's viscoelastic properties to collagen pre-orientation, but also to factors such as fiber density and thickness. A larger, more homogenous range of samples will allow for a more in depth analysis of storage and loss moduli [19], and a subsequent direct correlation with collagen dynamics. We aim to further increase our imaging speed, to minimize the waiting time and possible relaxation during the stress strain measurements. Furthermore, we will extend this framework to other extracellular matrix molecules such as elastin to further elucidate the microstructural basis of human skin's unique mechanical properties. Not limited to skin tissue, the combination of SHG imaging and mechanical deformation can also be applied to a range of biological tissues and applications including tissue engineering.

Conclusion
When used in conjunction with mechanical testing, SHG imaging is a powerful tool that allows for visualization of dynamic three-dimensional changes to collagen in biological tissues. Figure 6. Analysis of orientation and its effect on skin microstructural changes while undergoing uniaxial deformation. Figure 6(a): schematic representation of the experiment. Surgically excised abdominal skin was cut into identical samples. The first sample was cut along skin tension lines, the second perpendicular to these lines. Figure 6(b) shows the difference in stress strain response for both samples. Orientation indices at discrete intervals are plotted in paired colors, showing differential responses to uniaxial deformation. Finally, we thank Humeca (Borne, the Netherlands) for kindly providing us with a D80 cordless dermatome.

Data availability statement
The data is not being shared publicly because another paper is being written using this data. The data that support the findings of this study are available upon reasonable request from the authors.