## 1.

## Introduction

Forensic examiners are often asked to assess the age of traumatic bruises, which may result from an accident or physical abuse. The current protocol involves visual inspection of the bruise and relies on differences between the absorption spectra of the intact skin, extravasated hemoglobin, and the products of its subsequent biochemical decomposition. However, the perceived color of the bruise may vary strongly with the subject’s age, natural skin tone, ambient lighting, and so on. Consequently, the present approach—in addition to requiring the involvement of a highly trained expert—produces bruise age estimates that are subjective and of limited accuracy.^{1}2.^{–}^{3}

Objective analysis of bruises by measurements of diffuse reflectance spectra (DRS) in the visible spectral range was recently proposed as an approach to overcome the above limitations.^{4}5.^{–}^{6} A simple mathematical model of the processes governing the bruise evolution was introduced, i.e., mass diffusion and biochemical transformations of the extravasated hemoglobin.^{6} However, the inverse analysis of the DRS from test subjects, aiming to assess the (assumingly universal) model parameters, requires knowledge of the lesion geometry (i.e., the epidermal thickness and depth of vascular injury), which was unavailable. Thus, a presumed depth distribution of relevant chromophores was used, which may have contributed to the fact that the obtained parameter values were rather scattered and varied between different reports.^{6}^{,}^{7} (An account of additional limitations of the used approaches is left for Sec. 4.)

Pulsed photothermal radiometry (PPTR) presents a viable technique to assess the missing information. It is based on measurements of transient change of mid-infrared (IR) emission upon pulsed laser irradiation. From such radiometric signals, temperature depth profiles immediately after laser exposure can be reconstructed.^{8} This provides valuable information on depth distribution and concentrations of selected absorbers, as was recently demonstrated in various applications *in vivo*.^{9}10.^{–}^{11}

In parallel with providing basic characterization of the affected site, PPTR profiling also allows quantitative assessment of the hemoglobin dynamics and thus represents a complementary approach to DRS for analysis of bruise evolution. As we describe in present report, this is achieved by combining PPTR measurements in human subjects with numerical simulations of light transport and energy deposition in bruised skin. The high sensitivity of such an approach to parameters of the analytical bruise evolution model was demonstrated in our preliminary tests, where temperature profiles reconstructed from experimental PPTR signals were compared with those obtained from simulated ones.^{12}^{,}^{13} The results indicated that the hemoglobin diffusivity in skin is significantly lower than the values assessed earlier from studies based on DRS. This warranted further development and application of our approach in order to improve our knowledge of the bruise evolution processes.

Reconstruction of laser-induced temperature profiles from PPTR signals represents a severely ill-posed inverse problem. This leads to inevitable loss of accuracy and uniqueness in the iteratively obtained solutions, especially in the presence of experimental noise. However, our recent analysis indicated that this effect does not reflect in the considerable uncertainties of the assessed values for key bruise evolution parameters.^{14}

Nevertheless, we apply here a more streamlined methodology, where simulated PPTR signals are fitted directly to the experimentally obtained signals. This approach avoids the above controversy by eliminating the reconstructions of either depth profile. At the same time, it significantly reduces the computational complexity. This allows us to introduce objective multivariate fitting to extract reliable values of the relevant skin characteristics and parameters of the bruise evolution model.

Furthermore, we also introduce simultaneous fitting of multiple PPTR signals obtained at various times after incidental injury. This enables even more accurate and reliable assessment of the bruise evolution parameters as compared with analysis of individual PPTR signals. At the same time, this approach provides a unique insight into temporal variations of these parameters over the course of bruise resolution.

## 2.

## Methodology

## 2.1.

### Mathematical Model of Bruise Evolution

Bruising is initiated by the rupture of one or more blood vessels, which results in the rapid formation of a subcutaneous blood pool. The subsequent evolution of the bruise is governed by mass diffusion of extravasated hemoglobin toward the skin surface. The inflammatory response of the tissue initiates the biochemical transformation of hemoglobin into biliverdin and bilirubin. At the same time, bruise products are being removed through the lymphatic system. The process is depicted in Fig. 1. The epidermal thickness is denoted by ${d}_{\text{epi}}$, while ${d}_{\text{der}}$ marks the depth location of the blood pool formation measured from the epidermal–dermal (ED) junction.

In a laterally homogenous part of the bruise, the time evolution of the hemoglobin volume fraction, ${f}_{\mathrm{h}}(z,t)$, can be modeled using the one-dimensional (1-D) diffusion equation:^{6}

## (1)

$$\frac{\partial {f}_{\mathrm{h}}}{\partial t}=D{\nabla}^{2}{f}_{\mathrm{h}}-\frac{{f}_{\mathrm{h}}}{\tau}.$$Here, $D$ denotes the hemoglobin mass diffusivity (in units of ${\mathrm{mm}}^{2}/\mathrm{h}$) and $\tau $ is the characteristic time of its decomposition due to lymphatic system drainage and biochemical transformation to biliverdin in the heme oxygenase process. Biliverdin is, in turn, converted to bilirubin in the process of biliverdin reductase.^{6} Spectral properties of these compounds account for the green/yellow hue of the bruise several days after injury. However, because a green laser ($\lambda =532\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{nm}$) is used in our PPTR setup, their contribution to the induced temperature rise will be negligible. The dynamics of bilirubin generation and decomposition can thus be omitted from our analysis. (A quantitative verification of this claim is provided in Sec. 4.) The convection of hemoglobin due to pressure gradients is also neglected since its contribution is believed to be small.^{6}

The blood reservoir is assumed to form instantaneously at the time of injury. In our version of the model, the hemoglobin source function is box-like with the amplitude fixed at 100 vol% and a variable duration $T$. The latter controls the amount of extravasated blood and thus accounts for the severity of a specific injury. Because the epidermis is impenetrable to hemoglobin, a zero-flux boundary condition is applied at the ED junction (depth ${d}_{\text{epi}}$).

An example of a bruise evolution with time is presented in Fig. 2, indicating the characteristic variation of color from the initial reddish hue to darker, blue/indigo tones, and ultimately to green/yellow/brown discoloration.

## 2.2.

### Pulsed Photothermal Depth Profiling

PPTR involves time-resolved measurements of mid-IR emission following pulsed laser irradiation of the sample. From the radiometric signal transient, $\mathrm{\Delta}S(t)$, the laser-induced temperature depth profile, $\mathrm{\Delta}T(z,t=0)$, can be obtained by solving the integral equation

The underlying details and functional form of the kernel function $K(z,t)$ were presented earlier.^{8} The theory was later augmented to account for the strong variation of the tissue absorption coefficient within the acquisition spectral band.^{15}16.^{–}^{17}

By representing the PPTR signal and laser-induced temperature profile as vectors $\mathit{S}$ and $\mathit{T}$, respectively, they become related by a matrix equation,

## (3)

$$\mathit{S}=\mathit{K}\mathit{T};\phantom{\rule[-0.0ex]{1em}{0.0ex}}{K}_{i,j}=K({z}_{j},{t}_{i})\mathrm{\Delta}z.$$Due to the defect rank of the kernel matrix $\mathit{K}$, solving Eq. (3) for the initial temperature profile $T$ presents a severely ill-posed inverse problem.^{8} A unique solution does not exist in general and the most probable result is commonly obtained by iterative minimization of the residual norm squared, ${\Vert \mathit{S}-\mathit{K}\mathit{T}\Vert}^{2}$. Our custom reconstruction code involves the projected $v$-method with a non-negativity constraint.^{17} The accuracy and robustness of the process were further improved by nonuniform binning of the PPTR signals.^{18}

## 2.3.

### Experimental Setup and Procedure

The study protocol was approved by the Medical Ethics Committee of the Republic of Slovenia. The study involves healthy volunteers with fair skin (Fitzpatrick types I–II), 20 to 60 years of age, with fresh bruises resulting from sports-related accidents or similar minor injuries with a known time of occurrence.

In each subject, PPTR measurements were performed two to six times within 2 weeks after the injury. The same measurement sites on the bruise and nearby healthy skin were used for the entire measurement sequence. Prior to measurement, the superficial layer of dehydrated skin cells was removed by tape stripping. The site was cleaned with medical-grade ethanol, rehydrated using physiological solution, and air dried.

All test sites were irradiated with single 1-ms pulses at 532 nm, emitted from a medical-grade laser (DualisVP, Fotona, Ljubljana, Slovenia). The applied wavelength is absorbed well by epidermal melanin and intravascular as well as extravasated hemoglobin within the dermis. At a spot size of $\sim 5\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}$, the radiant exposure was around $0.5\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{J}/{\mathrm{cm}}^{2}$. Note that for the involved skin types and comparable irradiation conditions, the threshold for skin injury is around 7 to $10\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{J}/{\mathrm{cm}}^{2}$,^{19}^{,}^{20} so our subjects did not incur any discomfort or adverse effects.

Mid-IR emission from the tissue surface was recorded within the acquisition spectral band of a fast IR camera ($\lambda =3.5$ to $5.1\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mu \mathrm{m}$; FLIR SC7500) over a period of 5 s at a rate of 1000 frames per second. PPTR signals $\mathrm{\Delta}S(t)$ were produced by lateral averaging of the data over a subwindow corresponding to an area of $1.5{\times 1.5\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}}^{2}$ on the skin surface. The manufacturer provided calibration system (Hypercal™) was used for nonlinear conversion of the signal amplitudes to radiometric temperature values.

## 2.4.

### Simulation of Pulsed Photothermal Radiometry Measurements

Our simulation begins by computing the hemoglobin depth distribution ${f}_{\mathrm{h}}(z,t)$ at the selected time point, using Eq. (1) and the assumed model parameter values [Fig. 3(a)]. In order to correlate such a concentration profile with the PPTR data, we simulate laser energy deposition in bruised skin using the 1-D multilayer weighted-photon Monte Carlo modeling of photon transport in multi-layered tissue (MCML) technique.^{21}

Our optical skin model is composed of an epidermal layer and 50 layers within the dermis to account for the hemoglobin depth distribution ${f}_{\mathrm{h}}(z,t)$. The bottom layer is semi-infinite and starts at the depth where the blood reservoir is formed (${d}_{\text{der}}$). Optical properties are computed separately for each layer.

The absorption coefficient of the epidermal layer, ${\mu}_{a,\mathrm{epi}}$, is calculated using the well-known equations^{22}

## (4)

$${\mu}_{a,\mathrm{epi}}={f}_{\mathrm{mel}}\text{\hspace{0.17em}}\hspace{0.17em}{\mu}_{a,\mathrm{mel}}+(1-{f}_{\mathrm{mel}})\text{\hspace{0.17em}}{\mu}_{a,\mathrm{base}},$$## (5)

$${\mu}_{a,\mathrm{mel}}=6.6\times {10}^{10}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}{\left(\frac{\lambda}{\mathrm{nm}}\right)}^{-3.33},$$## (6)

$${\mu}_{a,\mathrm{base}}=0.0244\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}+8.53\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}\text{\hspace{0.17em}}\mathrm{exp}(-\frac{\lambda -154\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{nm}}{66.2\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{nm}}),$$Absorption in the dermis involves contributions from normally present absorbers and additional hemoglobin extravasated upon the injury. For the intravascular present hemoglobin, we assume a hematocrit of 41% and oxygen saturation level ${S}_{\mathrm{bl}}=80\%$. The absorption of blood is calculated as a linear combination of the values for oxygenated and deoxygenated whole blood.^{23} We also take into account the optical screening within the vessels using the customary correction factor and average vessel size distribution.^{24}^{,}^{25} By combining the obtained effective blood absorption coefficient, ${\mu}_{\text{bl}}$, with the baseline absorption (${\mu}_{a,\text{base}}$), we derive the absorption coefficient for intact dermis.

For each dermal sublayer, we must consider also the contribution of extravasated hemoglobin:

## (7)

$${\mu}_{a,i}={f}_{\mathrm{bl}}\hspace{0.17em}\hspace{0.17em}{\mu}_{\mathrm{bl}}+(1-{f}_{\mathrm{bl}}){\mu}_{a,\mathrm{base}}+{f}_{i}{\mu}_{\mathrm{h}}.$$Here, ${f}_{\mathrm{bl}}$ and ${f}_{i}$ are the fractional volumes of intravascular hemoglobin and extravasated hemoglobin concentration in the specific dermal sublayer, respectively. Oxygen saturation of the latter is believed to be lower than for intravascular hemoglobin due to the altered homeostatic environment. Thus, we apply an oxygenation level of ${S}_{\mathrm{h}}=40\%$.

The reduced scattering coefficient for epidermis and dermis is calculated as

## (8)

$${\mu}_{s,\mathrm{derm}}^{\prime}={f}_{\mathrm{bl}}{\mu}_{s,\mathrm{bl}}^{\prime}+(1-{f}_{\mathrm{bl}}){\mu}_{s,\mathrm{base}}^{\prime},$$^{26}The reduced scattering coefficient for epidermis thus equals ${\mu}_{s,\mathrm{base}}^{\prime}$. The epidermal and dermal layers are assigned the same refractive index value, $n=1.4$.

By also considering the thermal properties of the involved tissues, the simulated energy deposition profiles can be converted into the corresponding temperature profiles [Fig. 3(b)].^{17} Temporal changes of the latter are evidently not as pronounced as the variations of the hemoglobin profile in Fig. 3(a). This reflects the baseline contribution to optical absorption by normal dermis (containing 2% of vascular blood) and the complex dependence of light fluence at any given depth on both absorption and scattering properties of the shallower and even deeper skin layers.

From these data, the corresponding PPTR signals are easily predicted using Eq. (3). The results [Fig. 3(c)] essentially represent the time evolution of temperature in the superficial part of the skin, convolved with a function accounting for partial absorption of the IR emission contributions from subsurface tissue layers.^{8}^{,}^{15} The observed similarity between the initial signal jumps (from the value of 0, corresponding to the tissue baseline temperature) thus follows directly from the fact that the initial temperature rise within the epidermis (due to melanin absorption) is only marginally affected by the dynamics of extravasated hemoglobin deeper inside the skin [see Fig. 3(b)].

Similarly, the subsequent shape of the PPTR signals in Fig. 3(c) is dominated by thermal relaxation of the epidermis, which is again very similar for all three presented examples. Contributions from heat diffusing toward the skin surface from deeper heated regions are superimposed onto such backgrounds due to linearity of Eqs. (2) and (3).^{27} Owing to the diffusive nature of heat transport, however, the amplitudes of individual signal contributions inevitably decrease (nearly inversely proportional) with the absorber’s depth and also become increasingly similar to one another (see, e.g., Fig. 3 in Milanič et al.^{17}). More specifically, higher spatial frequency components of the initial temperature profile have a lesser effect on the resulting PPTR signal, and their contribution is progressively reduced with increasing depth. Both findings follow rigorously from the singular value decomposition of the kernel matrix $\mathit{K}$, which is a mathematical theory beyond the scope of the present manuscript.^{8}^{,}^{28}

Based on the same theory, the related inverse problem—i.e., assessment of the laser-induced temperature profile from PPTR signals—is characterized as severely ill-posed. In practical applications, this leads to the inherent nonuniqueness of the solutions and an inevitable loss of accuracy, which is accentuated at larger depths and lower signal-to-noise ratios.^{8}^{,}^{17}^{,}^{18} Nevertheless, by controlling the latter, applying a sufficiently detailed physical model and a carefully selected inversion algorithm, one can often extract valuable information from PPTR signals which appear very similar to one another.^{11}^{,}^{29}30.^{–}^{31}

## 2.5.

### Fitting of Measured Data

The experimentally obtained PPTR signals are matched with simulated signals by objective multidimensional fitting using the nonlinear least-squares method (function *lsqnonlin* in MATLAB Optimization Toolbox by Mathworks, USA). The initial 2 s of the transient PPTR signal, which corresponds to 2000 uniformly distributed data points, are used in the fitting process. As already pointed out above, this approach eliminates the reconstruction of either depth profile with several beneficial consequences.

In the first step of our analysis, we fit the simulated PPTR signals to those measured on intact skin. In our numerical model, healthy skin is represented by two homogeneous layers representing the epidermis and semi-infinite dermis. The best match is sought by varying the epidermal thickness ${d}_{\text{epi}}$, epidermal melanin content ${f}_{\text{mel}}$, and blood content in normal dermis ${f}_{\mathrm{bl}}$.

In the subsequent analysis of bruised skin, the parameter values assessed from the nearby healthy site are used as constants. The difference between the PPTR signals measured in the bruise and healthy skin is attributed exclusively to the presence of extravasated hemoglobin. In the fitting process, we thus vary the bruise evolution parameters, i.e., hemoglobin mass diffusivity $D$, decomposition time $\tau $, dermal thickness ${d}_{\text{der}}$, and source duration $T$, to obtain the optimal match between the simulated and measured PPTR signals.

The information on the parameter values of the presented bruise evolution model in literature is limited and somewhat controversial. Randeberg et al.^{6} used a hemoglobin mass diffusivity of $D=0.01\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{2}/\mathrm{h}$ based on the earlier reported diffusivity of myoglobin in skeletal muscles and observed agreement of their model predictions with experimental data. In a subsequent study, Stam et al.,^{7} who also incorporated lateral diffusion of chromophores in the bruise, reported a significantly lower optimal value, $D=0.0015\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{2}/\mathrm{h}$.

Similarly, the hemoglobin decomposition times, as assessed by matching the model predictions with measurements *in vivo*, were reported in a range of 72 to 120 h.^{6} However, considerably shorter decomposition times ($\tau =6$ to 12 h) were observed for some bruises with a distinct whitish central region.^{6}

In both studies referenced above,^{6}^{,}^{7} the depth of blood pool formation was assumed to match the dermal thickness. Moreover, as the latter value was also not available for a particular person and anatomic location, average anatomical values were applied in their analyses. On the other hand, both groups acknowledged the strong influence of the depth of blood spill on the spatial distribution of extravasated hemoglobin and consequently the evolution of the bruise color (and DRS).^{6}^{,}^{7} Accordingly, they pointed out the importance of having such information provided by an independent measurement on an individual patient basis.

## 3.

## Results

## 3.1.

### Healthy Skin

Relevant characteristics of intact skin are assessed from PPTR measurements on a test site near the bruise as described above. The obtained information, specific to the particular person and anatomic location, is used in the subsequent analysis of the bruise.

The result for a bruise located on the arm of a 27-year-old subject (15_AK; see Fig. 2) is presented in Fig. 4(a). An excellent match was obtained between the measured and simulated PPTR signals (solid and dashed lines, respectively). The corresponding model parameters for healthy skin are listed in Table 1, complemented with the algorithm-provided confidence intervals. The numbers in parentheses are the ratios of the latter and the former and serve as rough estimates of the relative errors.

## Table 1

Skin model parameters for the best-fitting simulated pulsed photothermal radiometry (PPTR) signal in healthy skin [see Fig. 4(a)]. The values are complemented with confidence intervals. The numbers in parenthesis mark the ratios of the latter to the former.

Parameter | Fitted value |
---|---|

depi (mm) | 0.161±0.001 (0.006) |

fmel (%) | 1.100±0.001 (0.001) |

fbl (%) | 2.41±0.01 (0.004) |

The matrix of correlations between the obtained skin parameter values is presented in Table 2. For example, the value of ${d}_{\mathrm{epi}}$ is very strongly and negatively correlated with both ${f}_{\text{mel}}$ and ${f}_{\mathrm{bl}}$. This means that, e.g., a small increase of the epidermal thickness can be largely compensated by a slightly smaller (relative) reduction of the melanin and/or blood volume fractions to produce similar PPTR signals. The correlation between ${f}_{\text{mel}}$ and ${f}_{\mathrm{bl}}$ is also rather strong but positive. At higher melanin content, the energy transmitted into dermis is reduced, but a similar temperature rise can still be obtained by increasing the blood volume fraction.

## Table 2

Matrix of correlations between the three healthy skin parameters listed in Table 1.

depi | fs | fbl | |
---|---|---|---|

depi | 1 | −0.86 | −0.91 |

fmel | −0.86 | 1 | 0.62 |

fbl | −0.91 | 0.62 | 1 |

As a visual representation of the implied structure and chromophore distribution in the healthy skin site, we present in Fig. 4(b) the temperature depth profiles as reconstructed from the measured and simulated PPTR signals in Fig. 4(a). The prominent subsurface temperature peak results from absorption of laser light in the epidermal melanin. Its width and amplitude are indicative of the epidermal thickness and melanin concentration, respectively. The smaller temperature rise in deeper layers is primarily due to absorption in intravascular hemoglobin.

## 3.2.

### Bruised Skin

We proceed to analysis of PPTR signals acquired from the incidental bruise (see Fig. 2). The measurements were performed at 16, 43, 92, and 115 h after the injury. The individual characteristics of skin as assessed from the nearby healthy site (Table 1) are applied here as constants. Only the four parameters of the bruise evolution model are varied to obtain the best match between the measured and simulated PPTR signals.

Figure 5(a) shows the measured and best-fitting simulated PPTR signal in a fully developed bruise, 92 h after the injury. The corresponding model parameter values are listed in Table 3. Figure 5(b) presents the laser-induced temperature depth profiles as reconstructed from the signals in Fig. 5(a). The somewhat erratic behavior of both solutions below $z=0.55\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}$ reflects the loss of accuracy with depth, inherent to PPTR depth profiling and primarily due to the ill-posedness of the involved inverse problem accentuated by experimental noise and rounding errors. The corresponding depth profile of extravasated hemoglobin, ${f}_{\mathrm{h}}(z)$, calculated using Eq. (1), is also plotted for illustration (dash-dotted line).

## Table 3

Bruise model parameters as assessed from PPTR measurement 92 h after injury [see Fig. 5(a)]. The values are complemented with confidence intervals, and the numbers in parentheses mark the ratios of the latter to the former.

Parameter | Fitted value |
---|---|

D (×10−4 mm2/h) | 7.4±0.2 (0.03) |

τ (h) | 87±2 (0.02) |

dder (mm) | 0.35±0.01 (0.03) |

T (sh) | 38±1 (0.03) |

The matrix of correlations between the four bruise model parameters is presented in Table 4. All parameters are strongly correlated as evidenced by the fact that the absolute values of all matrix elements exceed 0.5. The strongest (negative) correlation is indicated between $\tau $ and $T$. At this particular time point and bruise structure, the effect of prolonged source duration $T$ on the amount and distribution of extravasated hemoglobin could thus be approximately compensated, e.g., by a suitably reduced value of $\tau $.

## Table 4

Matrix of correlations between the bruise evolution parameters listed in Table 3.

D | τ | dder | T | |
---|---|---|---|---|

D | 1 | −0.57 | 0.59 | 0.56 |

τ | −0.57 | 1 | −0.72 | −0.99 |

dder | 0.59 | −0.72 | 1 | 0.63 |

T | 0.56 | −0.99 | 0.63 | 1 |

## 3.3.

### Simultaneous Fitting of Multiple Pulsed Photothermal Radiometry Measurements

Our models of skin structure and bruise evolution produce simulated PPTR signals which can be fitted very closely to the measured ones, both in healthy and in bruised skins, and the obtained parameters’ values appear to be assessed with tight confidence margins. However, due to strong correlations between individual parameters, the indicated confidence intervals do not reflect the actual accuracy of the reported values. Other combinations of parameter values may yield a similarly good match with the measured signal.

In order to achieve a more robust assessment of the bruise model parameters, we have implemented simultaneous fitting of multiple PPTR signals obtained at various stages of bruise evolution. The three parameters of intact skin, determined separately from the measurement at a healthy site (Table 1), are again used as constants.

In our first approach, we thus vary the parameters of the bruise evolution model ($D$, ${d}_{\text{der}}$, $\tau $, and $T$) to obtain the best overall match between the simulated and measured PPTR signals at all four time points (approach I). The resulting parameter values are presented in Fig. 6 (open circles). The hemoglobin diffusivity ($D=4.0\times {\mathrm{mm}}^{2}/\mathrm{h}$) is considerably lower than the value obtained from a single PPTR signal (Table 3). The other parameters have also changed in turn, in accordance with the correlations seen in Table 4. For example, the lower value of $D$ results in a lesser amount of extravasated hemoglobin diffusing from the blood pool into the dermis. However, a similar hemoglobin concentration can be obtained at delayed time points by prolonging the decomposition time (to $\tau =171\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{h}$ versus the earlier 87 h). Similarly, the blood pool persists somewhat longer ($T=46\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{h}$) and is located slightly deeper inside the skin (${d}_{\text{der}}=0.5\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}$).

These changes likely primarily result from the ill-posedness of the attempted inverse analysis (reflected in strong correlations between the assessed parameters), combined with relatively high-noise levels in the measured PPTR signals. On the other hand, they might also be indicative of a certain variability of the parameter values over the course of the bruise healing process. To evaluate this hypothesis, we test the following several alternative approaches to the PPTR signal analysis.

First, we allow the decomposition time $\tau $ to assume independent values at each time point, while the remaining model parameters are varied jointly, like before, to obtain the best match with all measured PPTR signals (approach II, see Table 5). The resulting values of $D$ and ${d}_{\text{der}}$ (Fig. 6, orange squares; column II in Table 6) are fairly similar to those assessed with the previous approach (column I). Meanwhile, the hemoglobin decomposition time exhibits a considerable variation over time, increasing from $\tau \sim 30\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{h}$ at the first two time points to over 202 h at 5 days post injury [Fig. 6(b)]. Such an effect is physiologically plausible, because the tissue inflammation response, driving the decomposition and removal of extravasated hemoglobin, would have subsided over time.

## Table 5

Summary of the fitting approaches (denoted I–IV) applied in our analysis, indicating the model parameters whose values are assumed to stay constant over the course of bruise development.

I | D, τ, dder, T |

II | D, dder, T |

III | D, T |

IV | T |

## Table 6

Bruise evolution parameters as assessed from PPTR measurements at different times after injury using the four fitting approaches (marked I–IV).

Approach # | I | II | III | IV |
---|---|---|---|---|

D (×10−4 mm2/h) | 4.0±0.1 | 4.6±0.3 | 10.3±0.5 | 7.6±0.4 |

9.1±0.2 | ||||

6.2±0.2 | ||||

7.2±1.1 | ||||

τ (h) | 171±16 | 36±4 | 33±7 | 65±19 |

26±1 | 66±7 | 104±14 | ||

131±4 | 56±1 | 67±1 | ||

202±4 | 101±1 | 110±2 | ||

dder (mm) | 0.50±0.01 | 0.46±0.01 | 0.56±0.04 | 0.52±0.13 |

0.73±0.11 | 0.72±0.03 | |||

0.43±0.01 | 0.38±0.01 | |||

0.42±0.01 | 0.40±0.10 | |||

T (h) | 46±2 | 37±1 | 45±4 | 42±4 |

The notion that the above effect is significant is supported by the analysis of the quality of the fits. As evidenced by the residual norms presented in Fig. 7, approach II results in considerably lower residual values (orange columns) as compared with approach I (white). Moreover, this is valid not only overall but also separately for each time point, which strongly suggests that $\tau $ does vary during the bruise healing process. This should thus be accounted for in the analysis of PPTR signals, which will then likely yield more accurate estimates for the remaining model parameters, too.

Encouraged by the above experience, we present next an even more ambitious analysis, where ${d}_{\text{der}}$ is also allowed to vary independently, in addition to $\tau $ (approach III). The additional decrease of the residual norms at each time point indicates that introducing the additional degrees of freedom in the analysis was justified (Fig. 7, red columns). The obtained parameter values are presented in Fig. 6 (red diamonds) and Table 6 (column III). Despite the rather large confidence intervals, the indicated increase of ${d}_{\text{der}}$ at the first two time points appears convincing [Fig. 6(c)]. The weighted average (${d}_{\text{der}}=0.44\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}$) closely resembles the constant value assessed with the previous approach. The variability of $\tau $ is now markedly smaller than that seen in approach II [Fig. 6(b)]. This suggests that the latter might have been a result of overcompensation; i.e., an attempt of the algorithm to account for variations of two or more model parameters by varying only $\tau $. Meanwhile, the value of $D=10.3\times {\mathrm{mm}}^{2}/\mathrm{h}$ is now considerably higher than assessed with the former approaches [Fig. 6(a)].

In the last fitting approach, $\tau $, ${d}_{\text{der}}$, and $D$ are all allowed to vary independently at each time point (approach IV). The resulting values of $\tau $ and ${d}_{\text{der}}$ assume very similar values and trends as seen in the previous approach (blue triangles in Fig. 6; column IV in Table 6). The obtained values of $D$ lie between the (constant) values assessed using the approaches II and III. Their time evolution exhibits a pattern very similar to that seen for ${d}_{\text{der}}$ [Figs. 6(a) and 6(c)]—a somewhat elevated value at the first time point and a further increase at $t=43\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{h}$, if the last two values are taken as the reference. However, the overall residual norm has not dropped significantly in comparison with approach III (Fig. 7).

The source duration $T$ is an inherent property of each specific bruise in our model, reflecting the amount of blood extravasated upon the initial injury. Consequently, the same value of $T$ must be used to describe all PPTR measurements in the time sequence.

Figure 8(a) presents the entire time sequence of measured (solid lines) and best-matching simulated PPTR signals (dashed line) obtained using fitting approach IV. The laser-induced temperature profiles reconstructed from these signals are plotted in Fig. 8(b). The corresponding profiles of extravasated hemoglobin, ${f}_{\mathrm{h}}(z,t)$ are presented in Fig. 8(c).

## 4.

## Discussion

In our experience, PPTR depth profiling is a very suitable technique for quantitative analysis of the bruise healing process. By combining the analytical model of bruise evolution with Monte Carlo simulations of pulsed laser energy deposition, PPTR signals measured in volunteers with incidental injuries could be successfully replicated [Figs. 4(a), 5(a), and 8(a)]. By applying an objective multivariate fitting algorithm, we could thus extract relevant quantitative information on the structure and composition of intact skin (i.e., epidermal thickness, blood, and melanin concentrations) as well as key processes governing the bruise evolution, e.g., hemoglobin mass diffusion and biochemical decomposition.

However, despite the above and the small confidence margins indicated for the assessed parameter values (Tables 1 and 3), strong correlations between the model parameters (Tables 2 and 4) limit the accuracy of the assessed values. In other words, similarly good matches could be also obtained with other combinations of the parameter values, with minute uncontrolled effects (such as experimental noise or deficiencies of our simplified models) dictating which combination is reported as optimal in any given example.

Simultaneous fitting of multiple PPTR signals, measured at various times after the injury, is a viable approach for more robust assessment of the bruise model parameters. When all the parameters were considered to be constant with time (i.e., our fitting approach I), the assessed values for $D$, $\tau $, and ${d}_{\text{der}}$ (Table 6, column I) were significantly different from those obtained from analysis of a single PPTR signal (Table 3). The rationale behind this approach is that each model parameter will have the strongest influence on the simulated PPTR signals at different stages of bruise evolution.^{13} For example, the values of $D$ and ${d}_{\text{der}}$ will mostly affect the hemoglobin profile ${f}_{\mathrm{h}}(z)$ at earlier times, while the values of $\tau $ and $T$ will have a stronger effect at later stages of the bruise evolution.

In addition, by allowing selected parameters to assume independent values at each time point, we were able to provide an insight into potential variations of these parameters over the course of the bruise healing process. We have implemented and tested several versions of such a simultaneous fitting, ranging from the most robust (approach II; only $\tau $ can vary with time) to the most ambitious with 13 free parameters (approach IV; $D$, $\tau $, and ${d}_{\text{der}}$ can all vary independently between the four time points).

As evidenced by our analysis of the residual norms (Fig. 7), the introduction of varying $\tau $ (approach II) followed by ${d}_{\text{der}}$ (approach III), each yielded a significant improvement of the match between the simulated and measured PPTR signals at all time points as compared with the approach with constant parameters (approach I). At the same time, the confidence intervals of the assessed values do not increase excessively and typically remain smaller than the differences between the parameter values at various time points (Fig. 6, Table 6). Both observations strongly suggest that the indicated trends are likely significant.

Moreover, we believe that these trends are plausible from the point of view of bruise physiology. Specifically, the trend of increasing $\tau $ [Fig. 6(b)] means a decreasing rate of hemoglobin degradation and lymphatic removal as the inflammatory tissue response, triggered by extravasation of blood, is gradually subsiding. As mentioned earlier in this text, the trend indicated by the approach II (orange squares) is most likely exaggerated. But this is not particularly problematic because solutions III and IV, which also yield significantly lower residual norms, indicate a much more gradual increase of $\tau $. The latter is also more plausible when considering that the assessed values of $\tau $ represent the respective effective values which are valid for the entire time interval from the injury to the PPTR measurement.

The indicated initial increase of ${d}_{\text{der}}$ [Fig. 6(c)] and its subsequent relaxation to a lower (assumingly normal) value can be interpreted as swelling (edema) of the injured site. This is a common result of trauma and typically lasts for a couple of days. As the dermal layer expands, the blood source is (mathematically speaking) temporarily pushed to a deeper location.

Our last fitting approach IV, where $D$ was also allowed to vary over time, resulted in a similar overall residual norm value as the approach III (Fig. 7). Moreover, the resulting values of $\tau $ and ${d}_{\text{der}}$ and their trends are very similar between the two approaches [Figs. 6(b) and 6(c)]. The assessed values of $D$ fall between the values assessed using the approaches II and III [Fig. 6(a)].

Based on the above, we consider this approach as potentially informative, while admittedly also somewhat speculative. We find it particularly intriguing that the assessed time evolution of $D$ is very similar to that of ${d}_{\text{der}}$ [blue triangles in Figs. 6(a) and 6(c), respectively]. If the last two time points are used as a reference (assumingly representing a lesion where the inflammation and edema have already subsided), then both variables exhibit a somewhat elevated value at the first time point and a further increase toward $t=43\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{h}$. Such a correlation could be attributed to the fact that the dermal swelling indicates an excess of extracellular fluid. This should decrease the density of the collagen/elastin matrix, thus enabling faster diffusion of the relatively large hemoglobin molecules.

The obtained values of hemoglobin diffusivity are in the range of $D=6$ to $9\times {10}^{-4}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{2}/\mathrm{h}$ (approach IV), or 4 to $10\times {10}^{-4}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{2}/\mathrm{h}$ if approaches I–III are taken into account. All these values are considerably lower than the last reported value, $15\times {10}^{-4}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{2}/\mathrm{h}$.^{7} We attribute this discrepancy to the fundamental advantage of our experimental approach, which assesses directly the depth distribution of extravasated hemoglobin combined with simultaneous analysis of measurements at several time points after injury.

In comparison, the DRS, used in the former studies of bruise dynamics, has several disadvantages.^{6}^{,}^{7} First, DRS lacks intrinsic sensitivity to chromophore depth, as the light scattered at all tissue layers is lumped together in the recorded spectra. Consequently, inverse analysis of DRS requires accurate knowledge of the lesion geometry, which was lacking (and thus assumed) in the discussed studies.

Second, their analyses were based on analytical solutions for diffuse approximation (DA) of light transport in the tissue, derived earlier for structures with two or three optically homogeneous layers.^{32}^{,}^{33} Clearly, representing the analytically predicted hemoglobin profile ${f}_{h}(z)$ with a two-step function (as the first layer must be assigned to the epidermis) is a very crude approximation, which all but prevents extraction of its specifics in the attempted inverse analysis. Third, the DA is valid only for scattering-dominated tissues. This assumption can be challenged in a bruise, containing unusually large concentrations of strongly absorbing hemoglobin. Finally, the DA is known to be invalid near the irradiated tissue boundaries, which makes its direct application in analysis of DRS controversial. As demonstrated in a recent report, the fitting of simulated DRS spectra from human skin with the DA solutions produces significant systematic errors in the extracted parameter values.^{34} Due to these limitations, we believe that the bruise model parameters extracted from the DRS measurements in the mentioned studies are quite inaccurate and potentially burdened with significant systematic errors.

In this study, we set the oxygenation level of extravasated hemoglobin to ${S}_{\mathrm{h}}=40\%$. This corresponds to the range of oxygenation values reported by Randeberg et al.^{6} Our initial tests confirmed that using the value of 80%, the same as for intravascular hemoglobin, yielded inferior results (i.e., higher residual norms when fitting the experimental PPTR data). Using ${S}_{\mathrm{h}}=0$ also produced results that were just slightly worse than with ${S}_{\mathrm{h}}=40\%$.

The omission of bilirubin dynamics from our analysis is based on its very low absorption coefficient at the 532-nm wavelength used in our PPTR setup (${\mu}_{a}=0.27\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}$) as compared with the values for oxy- and deoxy-hemoglobin (20.8 and $14.8\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}$, respectively).^{23}^{,}^{35} We have verified the validity of this approach by including the corresponding dynamical equation and characteristic times for bilirubin generation and decomposition as provided in literature.^{6} For the example analyzed in present study, we can thus assess that the maximal bilirubin concentration is reached at the location of blood spill (i.e., 0.9 mm below the skin surface), 53 h postinjury. At that point, the bilirubin contribution to skin absorption is only $0.08\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}$. This value is dwarfed by the contribution of extravasated hemoglobin at the same time point, $1.92\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}$.

The relative contribution of bilirubin to light absorption at 532 nm may be larger at later stages of bruise development when most of extravasated hemoglobin has been removed. However, its absolute contribution will be lower than the value given above, while even the baseline absorption of the dermis (containing 2.4% of intravascular blood) at the same wavelength amounts to $0.38\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}$. Our analysis (for the skin and bruise parameters as assessed in present study) thus shows that the maximal contribution of bilirubin to radiometric signals at the involved time points reaches, at most, 0.008 K, which is far below the noise level in our measurements.

In present report, we also present the temperature depth profiles reconstructed from measured and simulated PPTR signals in order to obtain a visual representation of the predicted chromophore distribution inside the bruised skin. Note, however, that the iterative reconstruction process is intrinsically prone to a certain loss of accuracy and uniqueness, which may vary with the level of experimental noise (absent in our simulations), complexity of the target structure, and so on. These issues were bypassed in the present analysis by directly fitting the measured PPTR signals rather than using reconstructed temperature profiles.

In addition, the simplicity of our skin model, especially the assumption of homogeneously distributed epidermal melanin, contributes to the observed differences between the profiles obtained from measured versus simulated signals. These are, however, mostly limited to the superficial skin layers and are thus not detrimental for our analysis of the hemoglobin dynamics in the bruise. In a follow-up study, we are testing slightly more complex skin models, so far with encouraging provisional results.

We would also like to comment on the high-residual norm values observed in the presented example for PPTR measurement at $t=43\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{h}$, regardless of the fitting approach (I–IV; Fig. 7). Note that this time point is very close to the assessed source duration $T$. We hypothesize that the box-like blood source function results in particularly unrealistic predictions at such times. Because the same effect also occurred in other subjects, we are currently exploring more realistic blood source functions and the first results are very promising.

In our continuation of the presented work, we plan to test the robustness of our analysis, in particular the fitting approaches III and IV, and the indicated temporal variations of the bruise model parameters $D$, $\tau $, and ${d}_{\text{der}}$. Within the limits of our current protocol, which allows for 30 analyzed cases, we will also check for any differentiation in the assessed parameter values and trends with respect to the subject’s age, gender, location of the bruise, and so on. Obtaining such information would represent an important step toward future development of a technique for objective determination of age in traumatic bruises.

## 5.

## Conclusions

Photothermal temperature profiling combined with Monte Carlo modeling provides valuable information on the depth distribution of extravasated hemoglobin in bruised skin. This approach allows *in vivo* assessment of the hemoglobin mass diffusion, decomposition rate, and dermal thickness in bruised skin. Simultaneous fitting of a time series of PPTR signals improves the fitting stability and results in lower confidence intervals of the determined parameters as compared with the fitting of individual signals. The obtained results indicate plausible variations in the inflammatory response, dermal thickness (attributed to edema), and possibly the hemoglobin mass diffusivity over the course of the bruise healing process.

## Acknowledgments

This work was supported by The Slovenian Research Agency through Grants P1-0192 and PR-04360. The authors thank Fotona D.D. (Ljubljana, Slovenia) for lending us the laser system for use in this study.

## References

## Biography

**Luka Vidovič** received his BS degree in physics from the University of Ljubljana in 2011. He is employed as a junior research assistant at Jožef Stefan Institute and is currently pursuing his PhD degree in physics. His research efforts focus on applications of pulsed photothermal radiometry and more recently also diffuse reflectance spectroscopy for quantitative characterization of human skin and cutaneous malformations.

**Matija Milanič** received his PhD degree in physics from the University of Ljubljana in 2008. Since 2003, he has been working as a researcher at Jozef Stefan Institute, Slovenia. Currently, he is employed as a postdoctoral researcher at IET, NTNU, Norway. The principal focus of his research is on the interaction of light with biomedical tissues, including photothermal radiometry, spectroscopy imaging, and numerical simulations.

**Boris Majaron** is a senior staff member at Jožef Stefan Institute in Ljubljana, Slovenia. After having obtained his PhD degree in physics from the University of Ljubljana in 1993, his research interests shifted gradually from laser physics and spectroscopy of solid-state laser materials to biomedical optics, specifically laser ablation and characterization of biological tissues using pulsed photothermal radiometry and diffuse reflectance spectroscopy. He is a member of the editorial board for *Lasers in Surgery and Medicine.*