Estimation of differential pathlength factor from NIRS measurement in skeletal muscle

The utilization of continuous wave (CW) near-infrared spectroscopy (NIRS) device to measure non-invasively muscle oxygenation in healthy and disease states is limited by the uncertainties related to the differential path length factor ( DPF ). DPF value is required to quantify oxygenated and deoxygenated heme groups ’ concentration changes from measurement of optical densities by NIRS. An integrated approach that combines animal and computational models of oxygen transport and utilization was used to estimate the DPF value in situ . The canine model of muscle oxidative metabolism allowed measurement of both venous oxygen content and tissue oxygenation by CW NIRS under different oxygen delivery conditions. The experimental data obtained from the animal model were integrated in a computational model of O 2 transport and utilization and combined with Beer-Lambert law to estimate DPF value in contracting skeletal muscle. A 2.1 value was found for DPF by fitting the mathematical model to the experimental data obtained in contracting muscle (T3) (Med.Sci.Sports.Exerc.48 (10):2013 – 2020,2016). With the estimated value of DPF , model simulations well predicted the optical density measured by NIRS on the same animal model but with different blood flow, arterial oxygen contents and contraction rate (J.Appl.Physiol.108:1169 – 1176, 2010 and 112:9 – 19,2013) and demonstrated the robustness of the approach proposed in estimating DPF value. The approach used can overcome the semi-quantitative nature of the NIRS and estimate non-invasively DPF to obtain an accurate concentration change of oxygenated and deoxygenated hemo groups by CW NIRS measurements in contracting skeletal muscle under different oxygen delivery and contraction rate.


Introduction
Near-Infrared spectroscopy (NIRS) can be used non-invasively to measure muscle oxygenation to evaluate physiological and pathophysiological conditions (Willingham and McCully, 2017;Grassi et al., 2019).The oxygenation of the heme group concentrations can be measured from different types of spectrometers: continuous wave (CW) NIRS provides concentration changes relative to a baseline; frequency domain (FD) or time domain (TD) NIRS provides absolute concentration (Ferrari et al., 2011;Grassi and Quaresima, 2016;Barstow, 2019).
The utilization of CW NIRS devices is limited by the uncertainties related to the differential path length factor (DPF) which is used in the Beer-Lambert law to quantify oxygenated (ΔHbMbO 2 ) and deoxygenated (ΔHHbMb) heme groups' concentration changes from the measurement of the optical density at different wavelengths.Thus, an incorrect DPF can lead to overestimate or underestimate the concentration changes (Endo et al. 2021;Pirovano et al. 2021).The typical range of values for DPF in brain and skeletal muscle tissues are: gastrocnemius (5.8-5.3),forearm (4.4-3.9),adult head (6.5-5.9) and infant head (5.4-4.7)(Ferrari et al. 1992;Duncan et al. 1995).It has been reported that DPF in tissues with a larger adipose tissue layer is higher for forearm and calf muscle (Ferrari et al. 1992;Duncan et al. 1995;Barstow, 2019).Thus, accurate DPF should be determined for each individual muscle group, taking into consideration the thickness of adipose tissue present (Pirovano et al. 2021), water content (Ferreira et al. 2007), and optical properties such as attenuation of light in the muscle due to absorption and scattering (Ferreira et al. 2007;Barstow, 2019;Endo et al. 2021;Pirovano et al. 2021).Besides the challenge related to the estimate of DPF which might in part be overcome with TD and FD NIRS spectrophotometers, the application of this technology is limited by the semi-quantitative nature of the signal obtained either with CW or TD and FD NIRS spectrophotometers.Specifically, the absolute or relative oxygenated (HbMbO 2 ) and deoxygenated (HHbMb) concentrations measured by NIRS devices result from the contribution of the oxygen content within an unknown microvascular and extravascular volume distribution.In addition, this volume distribution changes during exercise and it might be altered under disease states.Thus, the quantitative relationship between the NIRS signal and the heterogeneous oxygen distribution in the muscle region investigated is not known because a quantitative calibration of the signal in vivo is not available yet.
A quantitative approach based on a mathematical model of oxygen transport and metabolism has been proposed to quantify the contribution of Hb and Mb contribution to the NIRS signals (Koirala et al., 2021a(Koirala et al., , c, 2023;;Lai et al., 2009).The analysis was based on concentration changes normalized to the changes from rest to ischemia because DPF was unknown.Nevertheless, this approach has not been applied to analyze absolute and relative concentration change measured by NIRS.We propose to use a physiologically based model of oxygen transport and utilization to relate optical density measured by NIRS to the oxygenated/deoxygenated heme group concentration changes in the microvascular (arteriole, capillary, venules) and extravascular compartments.To quantify this relationship the analysis was based on NIRS data obtained from a canine model of oxidative metabolism to avoid adipose tissue effects on the NIRS signals.
In this study, the computational model was used to: (1) estimate DPF using optical density measurements in pumped-perfused muscle at different blood flow rates; (2) predict muscle oxygenated and deoxygenated changes in independent experiments from those used to estimate DPF, and determine whether the optical density measured in these independent experiments were similar to those predicted by the simulations using the estimated DPF; and (3) analyze the effect of key physiological parameters on DPF estimation.

Animal model
The experimental data obtained from an animal model of oxidative metabolism reported in previous studies (Hernández et al., 2010b;Goodwin et al., 2012;Sun et al., 2016) were used to estimate and validate the DPF.It should be noted that the trial or experiment notation of this work is the same as that used in the experimental works (Hernández et al. 2010b;Goodwin et al. 2012;Sun et al. 2016) from which the NIRS data are used for the analysis proposed here (Table 1).The experimental protocol allows measurement of muscle venous oxygen (C T ven ) concentration and optical density (OD λ1 and OD λ2 ) at two different wavelengths (λ 1 = 760 and λ 2 = 860 nm).All experimental data were averaged from 5 or 6 measurements for each trial condition.A continuous-wave NIRS system (Oxymon Mk III, Artinis Medical Systems BV) was used to measure the optical density (OD λ1 and OD λ2 ) associated with muscle oxygenation changes under different experimental conditions (i.e., blood flow, arterial O 2 content) and during contraction at different blood flows.The light at two different wavelengths (760 and 860 nm) was emitted and received by two fiber-optic bundles.The optodes were placed over the medial head of the left gastrocnemius and held in place with an elastic band.Experimental control of muscle blood flow was made possible by a peristaltic pump connected to the arterial inflow to the gastrocnemius.The right carotid artery or the right femoral artery was cannulated to route blood to the peristaltic pump (Hernández et al. 2010a).

Mathematical model
A mathematical model of skeletal muscle oxygen transport and metabolism developed in another study (Koirala et al., 2021a(Koirala et al., , 2023) ) was used to estimate the DPF of the continuous wave NIRS used to measure the muscle oxygenation under different experimental conditions.The derivation of the equations to calculate the physiological variables of interest are based on previous models (Lai et al. 2007;Spires et al. 2012Spires et al. , 2013)).The muscle volume (V mus ) is composed of extravascular tissue (cells and interstitial space, V t ) and blood (V b ) which is assumed to have vascular (V b,v ) and microvascular (V b,m ) compartments.The microvascular volume consists of arterioles, capillaries, and venules (V b,m = V art + V cap + V ven ) with volume fraction ω art , ω cap and ω ven .
The mathematical model simulates the spatial distribution of the free O 2 concentration in the capillary C F cap (v) and extravascular C F t (v) tissue compartments (Koirala et al. 2023).The concentrations depend on the tissue location as indicated by the muscle volume variable (v) from the arterial input v = 0 to the venous output v = V mus .Also, the model quantifies the relationship between O 2 concentration dissolved (C F ) and bound (C B ) either to Hb (arterioles, capillaries, venules) or Mb (myocytes).These are used to quantify the oxygenated Hb (HbO 2 ) and Mb (MbO 2 ) concentrations in the microvascular and extravascular compartments of the muscle as: (1) where f b,m and f t are the microvascular and extravascular tissue volume fractions in muscle, C B b,m is the bound O 2 concentration in the microvascular and C B cap and C B t are the spatial average of the bound O 2 concentration in capillary and extravascular tissue extravascular tissue (Koirala et al. 2023), respectively.f b,m is a function of the relative heme concentration changes (ΔC exp Heme ) detected by the NIRS signal and used to quantify the microvascular volume changes.The bound concentration for blood and extravascular compartments is related to the dissolved

Table 1
Experimental conditions of each trial.oxygen concentration by local equilibrium (Koirala et al. 2023).Similarly, the deoxygenated Hb (HHb) and Mb (HMb) contributions to the NIRS signal HHbMb are computed as: The derivation of the equations for interpreting the NIRS signal was documented in a previous human NIRS study (Lai et al., 2009).The simulated oxygenated and deoxygenated heme group concentrations of Hb and Mb relative to the rest condition (R) are calculated as: Thus, the total oxygenated (ΔHbMbO 2 ) and deoxygenated (ΔHHbMb) heme group concentration changes relative to the resting conditions are then computed as: The total oxygenated (HbMbO R 2 ) and deoxygenated (HHbMb R ) heme group concentration at rest were computed with the microvascular and extravascular volume distribution at rest (Table 2) by Eqs. ( 1)-( 4).
Microvascular blood volume change.The oxygenated and deoxygenated NIRS measurements are integrated in the mathematical model to quantify microvascular blood volume changes (Koirala et al., 2021a) associated with the heme concentration changes detected by NIRS.
The relative Hb and Mb concentration changes can be related to the microvascular blood volume fraction (Δf b,m ): Combining Eq. ( 9) with Eq. ( 12) yields a relationship depending on Δ C exp Heme measurement: The microvascular volume fraction in muscle at rest f R b,m is calculated by the product of the blood volume fraction in muscle (f b ) and microvascular volume fraction in blood (f m ) at rest.Since f b,m changes are assumed to take place in capillaries, the volume change takes place in the microvascular compartment (V b,m ): The V b,m change from rest (R) to each trial condition investigated is computed as: When blood volume changes occur, the extravascular (V t ), vascular (V b,v ), as well as arteriole (V art ) and venule (V ven ) volume remains constant at the resting value (See Table 2).Eqs. ( 14) and ( 15) can be rearranged as: Once the capillary volume is determined, the total microvascular volume can be calculated (V b,m = V art + V cap + V ven ); thus, the arteriole, capillary, and venule volume fractions are: The blood volume changes are attributed to the capillary compartment (V cap ), whereas arteriole and venule volumes are the same as those at rest and can be calculated from vascular, extravascular, and microvascular volume distribution (Table 2).The microvascular volume distribution was based on measurements from previous studies (Bebout, 1991;Hogan et al. 1993;Poole et al. 1995).f b,m is used to quantify V cap and the arteriole, capillary, and venule volume fractions (ω art , ω cap , and ω ven ).Details to calculate the arteriole, capillary, and venule contributions to the NIRS signals are reported in our previous studies.
The reaction rate coefficient of the ATPase flux, k ATPase , is a model parameter determined to simulate the NIRS signal for each trial (Table 1).The reaction rate coefficient can be computed as where f t is the extravascular tissue volume fraction in muscle and C ATP is the concentration of ATP in muscle.Because the muscle oxygen uptake (VO 2 ) are available for only four muscle blood flow data, a regression line of the VO 2 -Q data was used in Eq.( 18) to calculate k ATPase for each muscle blood flow of the simulations.

Estimate of DPF.
For each physiological condition investigated, the simulated ΔHbMbO 2 and ΔHHbMb concentrations are used to calculate the optical density with the modified Beer-Lambert law for Oxymon Mk III, Artinis Medical Systems BV as follows: where, d is the source-detector distance, OD λ1 and OD λ2 are the optical density at two different wavelengths λ 1 and λ 2 , respectively; ε HbMbO2,λ1 (1.207 mM -1 cm -1 ) and ε HHbMb,λ1 (0.798 mM -1 cm -1 ) are the extinction coefficient for HbMbO 2 and HHbMb at λ 1 and ε HbMbO2,λ2 (0.761 mM -1 cm - 1 ) and ε HHbMb,λ2 (1.017 mM -1 cm -1 ) are the extinction coefficients for HbMbO 2 and HHbMb at λ 2 .The optimal value of DPF is obtained by minimization of the least-square objective function: where n indicates the number of experimental data.The minimization is accomplished by numerical optimization using the 'lsqcurvefit' library in MATLAB and imposing the user-defined constraints like lower and upper bounds.It should be noted that ΔOD mod,λ1 e ΔOD mod,λ2 are  19) and ( 20) whereas ΔHHbMb and ΔHbMbO 2 are computed by the oxygen transport and utilization model Eqs.( 7) and ( 8).The dependency of ΔHHbMb and ΔHbMbO 2 from the oxygen concentration in the microvascular and extravascular compartments has been previously reported (Koirala et al. 2023).
Once the total oxygenated (ΔHbMbO 2 ) and deoxygenated (ΔHHbMb) heme group concentration changes are calculated from Beer-Lambert law for each time point, the absolute HbMbO 2 and HHbMb concentrations can be estimated as follows: 22) The plots are used to visualize the comparison between experimental (ΔHbMbO 2exp , ΔHHbMb exp ) and simulated (ΔHbMbO 2mod , ΔHHbMb mod ) data for ΔHbMbO 2 and ΔHHbMb obtained under different experimental conditions.Also, a quantitative comparison between simulated and measured value is quantified with the mean relative error (e) defined as: where n is the number of experimental data.

DPF estimation
The DPF parameter value was estimated with the minimization of the objective function (Eq.( 21)) based on the optical density measured and simulated by the mathematical model for different blood flows used in T3 trial.By fitting simulated optical density OD λ1 and OD λ2 to experimental data (Fig. 1a) obtained during different blood flow (T3), the optimal estimate was obtained for DPF (2.1) with a 95% confidence interval of (1.95, 2.21).As model input in these simulations, the parameter (k ATPase ) associated with the energy demand was used to simulate the changes in VO 2 (Fig. 1b).With the estimated value of DPF, the oxygenated (ΔHbMbO 2 ), deoxygenated (ΔHHbMb) (Fig. 2a) and total heme (ΔC Heme ) group concentration changes (Fig. 2b) relative to rest condition were calculated.The oxygenated (HbMbO 2 ) and deoxygenated (HHbMb) heme group concentrations at rest were approximatively 240 and 7 μM, respectively.The total heme group concentration change was used to simulate the microvascular blood volume changes (Fig. 2c) by means of Eq. ( 13).
The capability of the mathematical model in reproducing the experimental data was checked by predicting the NIRS data for B1, B2, CT20, EX45, EX70, and T4 trails obtained at different blood flow rates (Q) during contraction (cf.Table 1).Specifically, for each condition the simulation to predict the optical density OD λ1 and OD λ2 (Fig. 3a), was obtained with the same set of parameter values including that of DPF whereas the parameter (k ATPase ) associated with the energy demand was used to simulate the changes in VO 2 (Fig. 3b).The model predicted well the optical density OD λ1 and OD λ2 for each blood flow of all the investigated trials.The comparison between model simulations and experimental data was also reported in terms of oxygenated and deoxygenated NIRS signals as reported in Fig. 4. The concentration changes were obtained by converting the optical densities using the modified Beer-Lambert law.The model predicted well the changes of ΔHbMbO 2 and Δ HHbMb (Fig. 4a) with O 2 delivery (Q) changes.The mean relative error (e ΔHbMbO2 , e ΔHHbMb ) in predicting both the oxygenated and deoxygenated signals of B1, B2, CT20, EX45, EX70 and T4 experiments was approximately 16%.Model simulations tend to underestimate the experimental data for high blood flow.The higher the O 2 delivery, the lower the amplitude of the oxygenated and deoxygenated heme groups relative to the resting conditions.In both cases, the amplitude of OD λ1 decreased and ΔOD λ2 increased slightly in response to an increase of O 2 delivery as shown in Fig. 3a.

Effect of microvascular blood volume on DPF
Because of the microvascular volume fraction f b,m variability and its changes during contraction (cf.Fig. 4c), we analyzed the effects of f b,m on the estimates of DPF as shown in Fig. 5: with a decrease of f b,m from 7% to 3.5%, the estimated value of DPF increased from 1.45 to 2.2.Experimental evidence showed that DPF is affected by blood volume changes.

Discussion
A mathematical model of muscle O 2 transport and metabolism was used to estimate the DPF for a constant wave NIRS instrument used to measure oxygenated and deoxygenated heme groups' concentration changes in skeletal muscle.With the estimated DPF, model simulations well predicted the optical density of NIRS data measured under different experimental conditions.This result demonstrates the capability of the model to quantify the NIRS data and relate it to O 2 transport and utilization conditions and suggests that the approach used can overcome the semi quantitative nature of the NIRS with sufficient information about the microvascular volume distribution at rest.
The estimate of DPF enables us to quantify the NIRS measurement in terms of absolute heme group concentration (HbMbO 2 , HHbMb) in the whole skeletal muscle as a volume average of the extravascular and B. Koirala et al. microvascular.The approach proposed in this work can be applied to NIRS measurement in human skeletal muscle with the knowledge of muscle blood flow.The estimate of DPF can be obtained by experiments designed to stimulate significant oxygenated and deoxygenated heme group concentration changes for blood flow at different workload.

Differential pathlength factor (DPF)
Besides the motion artifact related to the muscle contractions, the DPF is affected by the optical properties of the skeletal muscle including extravascular tissue (myocytes), adipose tissue, and blood volume.The DPF has been reported to be in the range between 4 and 6 and varies among different tissues (Chance et al. 1992).The DPF value (2.1) estimated in our work is lower than the values reported in the literature.In a human NIRS study, DPF for skeletal muscle was reported to be 4.3 (Ferrari et al. 1992).A possible difference between our DPF value and those reported in several studies could be related to the absence of adipose tissue in the skeletal muscle region investigated in our animal model because the NIRS probe was placed directly on the skeletal muscle.In support of this view, several studies reported a relationship between the adipose tissue layer and DPF.In a study on humans' vastus lateralis, DPF decreased linearly with a decrease of the adipose tissue thickness (Pirovano et al. 2021) and in the absence of adipose tissue, the DPF value is expected to be 2.6 which is close to that estimated in this work.Furthermore, another NIRS study on heart muscle reported a DPF  experimental data for muscle oxygen uptake ( VO 2 ) of the trials B1/B2, CT20, EX45, EX70 and T4 (Hernández et al. 2010b;Goodwin et al. 2012;Sun et al. 2016) for different blood flow (Q) during contraction.Uncertainty about the estimation of DPF is related to the microvascular volume distribution which was assumed to simulate the NIRS signals.Because DPF is particularly sensitive to f b,m (Fig. 5), even small variation of f b,m determined important changes in the estimate of DPF.Thus, our findings should be confirmed with more accurate evaluation of f b,m , even though in our f b,m simulations we used values that are consistent and within the range reported in the literature.Existing experimental techniques such as MRI (Hindel et al. 2017) can be used to acquire quantitative information on the microvascular volume in skeletal muscle.The depth penetration of the light was not explicitly considered in our work.Usually, it is approximately half of the distance between the signal source and detector (Chance et al. 1992).Recently a statistical approach was proposed to provide a more accurate estimate of the penetration depth of light in media (Martelli et al. 2016).The error distribution for the data predicted in Fig. 4 presents an underestimation of the data value for high blood flow.Beside an error in the measurement, a possible reason for this discrepancy between simulated and measured data could be related to the uncertainties associated with the microvascular volume distribution.Indeed, the model predictions assume that the microvascular volume distribution is the same for all animal groups investigated.Although the animal groups present similar muscle characteristics, it is possible that among the groups the microvascular volume was different due to the probe position and group variability.Indeed, small changes in the microvascular volume can affect the DPF estimation (Fig. 5) and thus the NIRS signal.
Another important factor that could affect the estimate of DPF is the microvascular hematocrit which has been reported to be 50% (20-30%) of the systemic hematocrit (45-50%).In a NIRS study, this hematocrit has been considered to affect the NIRS signals (Davis and Barstow, 2013).Regardless of these current limitations in estimating DPF, it should be noted that the model can predict different experimental conditions with the same DPF, indicating that a unique DPF value is consistent with all NIRS measurements.The results of our simulations of the NIRS signal based on a volume average of the heme group concentration within the microvascular and extravascular compartments suggest that it is possible to overcome the semi-quantitative nature of the NIRS signal with accurate measurements of the microvascular blood volume fraction and distribution in skeletal muscles.Thus, a practical solution could be to determine the microvascular blood volume fraction by MRI at rest and use this information in the mathematical model to estimate DPF in contracting muscle.
Our NIRS experimental data were obtained assuming constant optical path length within the range of experimental conditions considered.It has been suggested that the optical path length varies during contraction; thus, the current assumption of considering a constant DPF would underestimate both NIRS signals (Endo et al. 2021).Further computational analysis could be performed with our mathematical model in combination with measurements of DPF.

Blood volume effect on DPF
The analysis of the NIRS data under different oxygen delivery conditions (Koirala et al., 2021a) indicated a primary role of capillary oxygen saturation in determining the sensitivity of the oxygenated and deoxygenated NIRS signal to blood volume change.Model simulations showed a decrease in DPF with an increase of microvascular blood volume.This finding is consistent with other experimental studies that reported a reduced optical path length with an increase in blood volume during exercise (Endo et al. 2021) or cuff occlusion (Ferrari et al. 1992;Hammer et al. 2019).In particular, for a blood volume increase of 28 μM detected by a Time-resolved NIRS during contraction, DPF decreased by 5-10% (Endo et al. 2021).Our analysis showed that an increase of microvascular volume corresponding to 8 μM determined a decrease in DPF of approximately 3%.Thus, the quantitative relationship between blood volume and DPF predicted by our model is close to that observed in this experimental study.Also, another human study attributed a DPF decrease of 7% to an increase in blood volume during cuff occlusion (Ferrari et al. 1992).The relative changes observed in these studies are consistent with those quantified by the mathematical model here proposed.
In conclusion, the integrative approach adopted allows us to estimate the DPF by relating the tissue oxygenation measurement by NIRS to the microvascular and extravascular volume distribution.This unique DPF can be used to quantify the oxygenated and deoxygenated heme group concentration changes from rest for different blood flow rates during contraction.This approach was further tested in predicting the oxygenated and deoxygenated heme group concentrations obtained with additional and independent experiments.This finding indicates that it is possible to quantify the heme group concentration measured by the NIRS device with reliable quantitative data on the microvascular and extravascular volume distribution within the muscle region investigated.

Fig. 1 .
Fig. 1.Trial T3: comparison (a) between model fitting (solid line) and experimental data (open circles) for optical density (OD λ1 and OD λ2 ) and between model simulation and experimental data (grey squares) of muscle oxygen uptake ( VO 2 ) of trial T3 at different blood flow (Q).

Fig. 2 .
Fig. 2. Trial T3: comparison (a) between model simulation and experimental data (open and grey circles) for oxygenated (ΔHbMbO 2 ) and deoxygenated (ΔHHbMb) heme group concentrations changes, (b) heme group concentration changes (ΔC Heme , grey squares) and (c) the microvascular volume fraction (f b,m ) changes with Q, f b,m is estimated from ΔC Heme (Eq.(13)).The heme group concentrations changes are relative to the rest condition (R).

Fig. 3 .
Fig. 3. Comparison between model prediction and experimental data for (a) optical density (OD λ1 and OD λ2 ) and comparison between model simulation and

Fig. 4 .
Fig. 4. Comparison between model prediction and experimental data for (a) oxygenated (ΔHbMbO 2 ) and deoxygenated (ΔHHbMb) heme groups concentration changes and (b) heme group concentration changes (ΔC Heme ) of the trials B1/B2, CT20, EX45, EX70 and T4 for different blood flow (Q) during contraction; and (c) the microvascular volume fraction (f b,m ) changes with Q estimated from ΔC Heme (Eq.(13)).The heme group concentrations changes are relative to the rest condition (R).

Table 2
Vascular, extravascular, and microvascular volume distribution at rest.