Modelling of Thermal Hyperemia in the Skin of Type 2 Diabetic Patients

The microcirculatory response to thermal stimulation involves both an axon reflex and NO-mediated activation. The analysis of the microcirculatory flow following thermal stimulation may therefore enhance the detection of any impairment of the small unmyelinated fibres that are involved in the axon reflex. The aim of this work is to establish a method of non-invasive measurement of small fibre impairment. The microcirculatory flow in response to local heating is measured by using a laser Doppler instrument, and mathematically modelled to extract a set of quantitative parameters. The results confirm that there is a significant difference in the parameters modelling the axon reflex between diabetic and control subjects, while no significant difference is found in the parameters modelling the NO-mediated activation.


INTRODUCTION
Diabetic neuropathy occurs in about 60% of Type 2 Diabetes Mellitus (T2DM) patients [1].There are several forms of diabetic neuropathy, depending on the type of nerve fibres involved.An autonomic form of diabetic neuropathy is the consequence of small fibre damage, but at present, there are no standard, non-invasive methods for directly investigating them in humans.The most common test used in clinical neurophysiology, i.e., standard electroneurography, is only able to evaluate the state of large myelinated nerve fibres [2].
In the clinical course of autonomic diabetic neuropathy, there is typical impairment of the thermoregulatory functions due to the sudomotor disorder [3] and impaired heat and cold induced vasomotor responses within the cutaneous microcirculation [4,5,6]; in fact, the failure of small peripheral unmyelinated fibres (type C fibres) causes impairment of the arteriolar smooth muscle cells and precapillary sphincter functions.
The cutaneous blood flow is one of the principal components involved in the body thermoregulation processes [7].In physiological conditions, local skin heating provokes a blood flow response characterized by an initial transient peak of vasodilatation, followed by a "plateau" phase [8,9].This pattern of response is the result of the excitation of a retrograde peripheral nervous pathway, called the axon reflex, consisting mainly of C fibres.The initial peak of vasodilatation is closely related to the axon reflex, while the following plateau phase is mediated by the release of Nitric Oxide (NO) from the endothelial cells [10,11].Moreover, Hodges et al. [8] demonstrated the existence of a correlation between the axon reflex and the sympathetic function.A decrease in the amplitude of the axonic flare in diabetic subjects was also observed by Noury et al. [12], who used thermal stimulation, and by Kramer et al. [13], via electrical stimulation of the tissue.
Laser Doppler analysis is a very practical, non-invasive and fast method for estimating the skin blood-flow responses under well-defined conditions.More specifically, by applying a special probe to a limb of the subject, we are able to evaluate the variability of the vasomotor activity of the microcirculation in relation to physiological environmental stimuli in both the temporal and frequency domains.From these evidences, we can assume that it is possible to correlate the blood flow response to a local skin-heating stimulus with the degree of peripheral unmyelinated fibre damage.Therefore, since C fibres control peripheral microcirculation and are widespread in peripheral body tissues [14], the purpose of this study is to develop a non-invasive method for measuring the involvement of small nerve fibres in diabetic neuropathies.

Subjects
The trial involved a group of 28 subjects with type 2 diabetes-related polyneuropathy, in whom it was possible to evoke a sensitive conduction from the sural nerve.The patients showed symptoms suggesting polyneuropathy and dysautonomy (orthostatic hypotension, intestinal peristalsis disorders, hypoesthesia of the lower limbs and chronic symmetric painful syndrome) as described extensively in [15].Patients included 13 females and 15 males aged from 63 to 75 years (mean = 68.5 years, SD , 3.46 years).The study was approved by the local Ethics Committee.Participants agreed to participate in the study after receiving comprehensive information about the purpose of the research and the procedures involved, and all signed the informed consent form.
The study also included a control group of healthy volunteers who only underwent the Laser Doppler flowmetry.The control group consisted of 19 subjects, aged between 22 and 55 years (mean = 29.9years, SD = 10.2 years).

Experimental Settings
The cutaneous perfusion signals were acquired in the central third of the forearm (between the ulna and radius), as illustrated in Figure 1.The acquisitions were mostly performed on the right arm; however in a few cases, it was necessary to perform on the left arm.The measurements were carried out in a sitting position, after acclimation in a controlled-temperature room (at 21°C).Great attention and cooperation were required from the subjects to reduce movement artifacts.However, a few artifacts were present in some samples, and were automatically removed using an algorithm based on amplitude thresholds, as described below.
A Periflux System 5000 (Perimed Inc.) was used for transcutaneous oxygen monitoring in the experiment.The device was connected through the analogue interface to a computer for digital conversion and data storage.The digitizing procedure was performed by means of a National Instrument NI-6009 USB board.Periflux System 5000 is a modular device which includes different units that can be used either alone or in combination.The units employed in this study were the following: • PF 5010 LDPM (Laser Doppler Perfusion Monitoring) unit.This module contains a low-power light source (a diode laser with a wavelength of 633 nm) which emits a laser beam capable of reaching a depth of about 1.5 mm.This value ensures the measurement of perfusion in the microvessels of the superficial plexus.• PF 5020 Temp unit.This unit allows for measuring the skin temperature, as well as supplying thermal stimulation at a constant temperature.
Both units were connected to a probe (Perimed thermostatic laser Doppler probe 457), illustrated in Figure 1, which conveyed the light from the source to the skin (and from the skin to the receiver), by means of 2 optical fibres (fibre separation was 0.25 mm).The probe had a cylindrical shape, 10 mm diameter, and 8 mm height, enabling for positioning on the skin without any patient discomfort.The two measurements taken during the experiment, cutaneous perfusion and skin temperature, were acquired through the same probe.Experimental setting.
While the temperature signal was automatically calibrated and expressed in Celsius degrees, the instrument measured the flowmetry signals in arbitrary units (PU -Perfusion Unit), and required system calibration by measuring known values of perfusion.The correct calibration of the device, which was verified according to the manufacturer's specifications, entailed two steps:

•
Measuring of the perfusion in the motility standard; • Measuring of the perfusion on the "zeroing" disc.
With the device properly calibrated, the measurements were 250 ± 15 PU on the motility standard and 0 ± 1 PU on the zeroing disc, respectively.When these values were within the two ranges, no further calibration was performed.The calibration was checked every two weeks.
During each measurement session, two signals were acquired from each subject: • Cutaneous blood flow in the tissue volume underneath the probe.This signal was obtained through the PF 5010 LDPM unit and digitalised (from the NI-6009 USB board) at a sampling frequency of 100 Hz.

•
Skin temperature, measured by the PF 5020 Temp unit.This signal was also sampled at 100 Hz.In addition, this unit is also able to provide local warming of the skin, heating the tissue underneath the probe to a user-defined temperature (43°C in this study).

Thermal Stimulation and Data Acquisition
A single measurement was performed on each subject and the cutaneous perfusion signal and the skin temperature signal were acquired.The laser-Doppler probe was placed on the forearm skin, and secured with medical tape to reduce any displacements that could cause artifacts.Each acquisition lasted 35 minutes and consisted of two phases: • The first phase (3 minutes) where the baseline flow was recorded without any local warming of the skin; • The second phase (32 minutes) where the skin under the probe was heated to 43°C by the laser Doppler probe.
During the first phase, no stimulation was applied; therefore, the skin temperature and blood flow remained more or less constant, as shown in Figure 2. Subsequently, the second (heating) phase was started by manually pressing the "Heating on/off" button on the front panel of the PF 5020 Temp unit.The time constants of the heating system are very short; therefore, the probe temperature only takes a few seconds to reach the desired level (43°C).
The duration of the thermal stimulus (32 minutes) and the skin temperature value (43°C) were selected to enable visualization of the initial peak of vasodilatation and the following plateau.In fact, as reported in [9], the response of local skin warming consists of an initial peak caused by a local axon reflex, followed by an NO-mediated vasodilatation (the plateau phase).A die-away phase can be expected 60 minutes after the beginning of the local warming [8,10].However, we hypothesize that the most probable effects relating to diabetic neuropathy concern the early phases of the vascular response to the stimulus, directly related to nervous activation; hence the duration was selected to include the initial peak and the onset of the plateau phase.
During the trial, no subjects reported pain or discomfort caused by local heating of the skin.Furthermore, there was no evidence of systemic reactions due to the local warming.A spectral analysis of the acquired signals did not reveal any variations in the cardiac frequency which remained stable in all patients during the trial.

Pre-processing
The blood flow signal showed a significant oscillation synchronized with the cardiac rhythm (Figure 2, top), and a few movement artefacts which appeared as large spikes.The signal was pre-processed to extract the average value of perfusion during each cardiac cycle and for spike detection and removal.
The vasomotion components (endothelial, myogenic, neurogenic, respiratory and cardiac) [16,17] present in the cutaneous microcirculation caused a great variability in the instantaneous values of the blood flow signal.However, this variability did not affect the study, as its aim was to analyse the variations of the mean skin blood flow responsible for the heat transfer, and find a trend in the vasodilatory response caused by a local thermal stimulus.On the other hand, the oscillation amplitude of the pulse wave was directly correlated with the pressure wave and peripheral resistance and elasticity, and this will be analysed in future work.The first pre-processing step consisted of a low-pass filtering (IIR Butterworth, cutoff frequency 2.5 Hz) which reduces acquisition noise and high-frequency disturbances.Spike detection and removal were implemented via non-linear filtering of the signal.A moving average was applied to evaluate local mean and local standard deviation of the signal.Each sample deviating from the local average by more than three times the local variance was labelled as a spike, and was smoothed out from the signal.Subsequently, a median filter was employed to remove the oscillatory components of the vasomotion without altering the trend of the acquired signals.The advantage of this filter with respect to a smoothing filter lies in its ability to achieve considerable attenuation of the heartbeat oscillation with a very low distortion of the trend of the signal.An example of a signal before and after the pre-processing is illustrated in Figure 3.

Modelling
The local heating of the skin induces a vasodilatory response in the cutaneous microcirculation.This response is divided into two phases, which have been identified as an initial peak, derived from a local axon reflex, and a plateau phase, mediated by the releasing of NO.The signal obtained after the pre-processing represents, with good approximation, the trend of the mean blood flow generated by the overlap of these two effects.
In order to identify the analytical functions that best represent the experimental signals, the initial peak and the plateau phase were modelled by using two independent functions consistent with the different physiological sources of the two effects.The mathematical expression of the fitting functions was selected on the basis of a heuristic and physiological hypothesis, after a visual inspection of the mean blood flow signals.
The vasodilatation caused by the axon reflex occurred within the first 5 minutes of the local warming, and its effect produced the initial peak of the flow signal, as Sample blood flow signal before (black) and after (gray) pre-processing.
described by [10].The typical peak had an almost symmetrical "bell" shape, and was thus modelled by a Gaussian function: ( where a 1 is the "bell" amplitude (the maximum amplitude of the initial peak), b 1 represents the mean value of the distribution (the instant when the peak occurs), and c 1 is the Gaussian standard deviation (related to the temporal duration of the peak).
The optimal values of the parameters were determined by fitting the Gaussian curve to the initial part of the flow signal using a least squares method.The result of the fitting procedure on a sample signal is shown in Figure 4. Afterwards, the residual trend was evaluated as the difference between the filtered signal and the Gaussian function.This residual trend constituted the results of the NO-mediated portion of the vasodilatation.As shown by the example in Figure 5, this signal showed a rapid amplitude increase, followed by a sustained vasodilatation phase and, finally, a die-away phase [10].However, the duration of the acquisition window did not always allow for determining the extent of the decreasing phase.Accordingly to the expected shape of the curve, we selected a representation of the plateau phase based on the difference of two exponential functions:  Fitting of the initial axon-related peak with a Gaussian curve.
where a is proportional to the plateau amplitude, b is the time constant of the descent phase, and c represents the time constant of the ascent phase.Since a rapid increase against a slow decline can be observed, it is reasonable to expect that b is much smaller than c.An example of the result of the fitting procedure is shown in Figure 6.Sample blood flow residual signal (dashed line) and optimal exponential fitting (solid line).

Data analysis
Data pre-processing, curve fitting and statistical analysis of the results were implemented using Matlab software (Matlab R2008b, The Mathworks Inc.).The lowpass filtering was performed with a zero-phase digital filtering, implemented with an IIR Butterworth filter with order 20.The median filtering was performed using a moving window of 15.5 s.This value represents a compromise between the reduction of the amplitude fluctuations and the maintenance of the trend of the mean blood flow.The curve-fitting phase was performed with the Levenberg-Marquardt algorithm [18] for solving the non-linear least squares problem.

RESULTS
In order to verify the intra-individual reproducibility, a set of 9 measurements was performed on the same subject (healthy subject, age 26).The experimental settings were the same as those described in section 2.2.The probe was positioned at different locations of the forearm, both on the volar and dorsal sides.The parameters b 1 and c 1 calculated for these signals showed much smaller intra-subject variability (b 1 = 105.86± 7.51 s, c 1 = 95.48 ± 11.00 s) than the inter-subject variability reported in Table 1.This suggests that the axon-reflex parameters are not influenced by the positioning of the probe on the forearm, thus indicating good reproducibility of the method.The fitting procedure was applied independently for each subject to obtain the values of the six parameters describing the flow response.The distribution of each parameter in both diabetic and control groups was visually analysed using box plots, as illustrated in Figure 7.The average value and the standard deviation of the six parameters were also evaluated, as shown in Table 1.A two-sided t-test was used to assess the significance of the differences.A Welch's t-test was selected, as the two groups had different sizes and variances, as suggested by the variance data reported in Table 1.The test indicated that a highly significant difference ( p < 0.01) existed in the values of b 1 , c 1 , and c between the diabetic and the control groups.
In order to assess the presence of a linear dependence between the different parameters, the correlation matrix was calculated.As shown in Table 2, a high correlation exists between the parameters b 1 and c 1 , indicating an almost linear relation between these two parameters.The graphical representation of this relation is shown in  The strong relation between the two parameters is evident.

DISCUSSION
Our results indicate that, in a reproducible condition of thermal stimulation, several parameters show a significant difference between the controls and the diabetic patients.Among those parameters, a 1 , b 1 and c 1 have a direct correlation with unmyelinated fibres leading to the vasomotor response.The parameters b 1 and c 1 are closely linked to the temporal properties of the initial peak of the blood flow, while a 1 is more closely   related to the amount of active microvessels recruited by the stimulus.This latter variable cannot be considered a specific consequence of autonomic nervous activity because it depends on many other factors as well, first and foremost, the individual variability of the micro-vessel anatomy.Moreover, the laser Doppler signal is proportional to the concentration of moving blood cells and the velocity of these erythrocytes.A possible neuro-mediated variability of a 1 could however explain the result of a synchronized spatial recruitment of a very large number of vasomotor units (as happens during intense, massive autonomic activation, for instance due to sudden, widespread body exposure to cold, stress, Valsalva manoeuvre, deep breathing, etc.).Indeed, the very significant difference in b 1 and c 1 which directly express specific proprieties of nerve fibres such as the temporal synchronizing of the reflex response, is extremely interesting.Both parameters are indicative of a delay in the vasomotor response related to peripheral nervous influences, as can be expected in diabetic patients.A large variability of the a and b parameters was observed, and the difference shows a higher p value despite the great distance between mean values.Furthermore, as regards the parameters describing the plateau phase (a and b), the values of p between the controls and the diabetic patients confirm the initial assumption, namely, that this phase is not directly due to a nervous activity but rather, to the NO discharge by the endothelial cells.With this interpretation of the results, the c parameter, which is lower in the diabetic group, is the greatest proof supporting our hypothesis that in the presence of dysautonomic neuropathy, the extent of the micro-vessels activity is not compromised; instead, its temporal dynamics are impaired.There are several clinical implications of the present results.First, the result of thermoregulatory impairment might be able to provide more detailed information about diabetes, as most studies focus on the averaged or maximum blood flow values [5,19] and disregard the temporal dynamics of the response.The scarcity of clinical studies is mainly due to the difficulty in using methods for specific identification of a small fibre neuropathy in a diabetic population and in differentiating this form from the large fibre polyneuropathies where the impairment of tactile conscious sensibility depends on the reliability of anamnestic and/or psychometric data [20].There are only few studies investigating thermoregulatory impairment as a risk factor for the other diabetes-related complications, such as infections and ulceration of the lower limbs [21,22], but the hypothesis in the literature about the microcirculatory pathogenesis is unable to distinguish between the oxidative damage of the endothelium and the impairment of small nerve functionality.Other important perspectives that this work could stimulate include the possibility of using thermoregulatory parameters to estimate other clinical consequences of small fibre failure.There is a long list of pathological conditions that could be investigated, ranging from painful neuropathies to cardiovascular autonomic impairment and gastrointestinal tract disturbances in diabetics.On the other hand, the main limitation of the method is the variability of the amplitude of the perfusion values that are obtained using the laser Doppler instrument.According to our results, this variability is mainly correlated with the positioning of the probe and the total returning light level [23]; however it has a minor impact on the time constants of the flow signal that are relatively independent of any scaling factor applied to the signal and present the most significant variation.
Our next project will aim to provide a more precise definition of the indexes found, by means of extensive measurements on large diabetic populations, and to perfect a non-invasive and practical method of diagnosing small fibre neuropathy.
We set the experimental duration of this study at 35 minutes; however, as the axon reflex is manifested during the first few minutes of local cutaneous heating, this duration could be drastically reduced (to about 5 -10 minutes) in order to ensure greater patient comfort.

CONCLUSIONS
In this experiment, we represent the microcirculatory response due to an axon reflex by means of two objective parameters: b 1 , proportional to the instant of maximum response, and c 1 , related to the temporal duration of this phenomenon.
The main finding of this study is the difference that exists in these two parameters between diabetic patients and control subjects, which could be employed to noninvasively quantify the presence of a small fibre neuropathy.
We conducted a highly reproducible, non-invasive and rapid analysis based on laser Doppler investigation of the time-domain characteristics of the response of the microcirculatory blood flow, triggered by a thermal stimulus.Several quantitative indexes were evaluated, and some of these have been recognized as promising clinical parameters for the assessment of small-fibre functionality in diabetic autonomic neuropathies.
Figure 3.Sample blood flow signal before (black) and after (gray) pre-processing.
signal Gaussian model of the initial peak

Figure 4 .
Figure 4.Fitting of the initial axon-related peak with a Gaussian curve.

Figure 5 .
Figure 5.Sample blood flow signal, Gaussian fitting, and residual of the fit.

Figure 8 ,
Figure 8, which shows the plot of three sample cases with different values of b 1 and c 1 .The strong relation between the two parameters is evident.

Figure 7 .
Figure 7.Box plot of the distribution of the parameters.

Figure 8 .
Figure 8.Three sample measurements with different Gaussian parameter values.