Broadband Time Domain Di ﬀ use Optical Reﬂectance Spectroscopy: A Review of Systems, Methods, and Applications

: This review presents recent developments and a wide overview of broadband time domain di ﬀ use optical spectroscopy (TD-DOS). Various topics including physics of photon migration, advanced instrumentation, methods of analysis, applications covering multiple domains (tissue chromophore, in vivo studies, food, wood, pharmaceutical industry) are elaborated. The key role of standardization and recent studies in that direction are discussed. Towards the end, a brief outlook is presented on the current status and future trends in broadband TD-DOS. monomer uptake in softwood samples. The clear increase in the absorption spectrum ( a ) of the natural wood around 1100–1200 nm following soaking in the monomer bath (treated wood) is clearly reﬂecting the absorption spectrum of the monomers (MIX) [151]. The concentration of the monomer at di ﬀ erent positions along the wood sample ( b ) demonstrates a precision in monomer assessment better than 1% in mass [152].


Introduction
The study of the propagation of light in highly scattering media has great importance in a wide range of applications (biological media, pharmaceutical, food, wood) [1][2][3][4][5][6]. Over the last decade, advancements in instrumentation and novel methods of analysis and methodologies have fostered innovation and breakthroughs in diffuse optical spectroscopy (DOS) [7][8][9]. In particular, probing human tissue with light brings in a slew of benefits thanks to non-invasiveness, the potential for compact bedside or wearable monitoring, and cost-effective solutions [10][11][12]. That makes diffuse optics an attractive tool to explore the human body in vivo [12][13][14][15]. To harvest maximum benefit, various practical implementation approaches have been developed by the scientific community. Typically, the measurement scheme includes a light source and a detector placed at a known distance, and the acquired data are processed using diffusion models to extract the optical properties (absorption, scattering) of the medium under investigation. The continuous wave (CW) approach uses CW light sources to estimate optical properties [16][17][18][19]. Thus, only light attenuation is measured. This approach is easier to implement and cost-effective. However, it suffers from the coupling between absorption and scattering properties that cannot be independently evaluated from a single measurement. The solution to the coupling could be found by adding another dimension to the measurement. That can be achieved performing CW measurements at multiple source-detector distances or a different approach: frequency domain (FD) and time domain (TD): frequency domain (FD) and time domain (TD) acquire photon 1.
Future prospects.
The first part, the physics of photon migration, provides the fundamentals of temporal propagation of photons in highly scattering media. The gamut of broadband TD-DOS instruments and advancements is discussed in the second section. The next section elaborates on various analysis methods and methodologies to extract optical properties and estimate tissue constituents. A wide range of broadband TD-DOS applications is covered across various domains (tissue constituent characterization, in vivo studies, food, wood, standardization techniques) are deeply elaborated in this section. A brief summary on the future prospects of the broadband TD-DOS is presented at the end. This review will provide a broad overview of the history and recent advancements of broadband TD-DOS as shown in the timeline Figure 1, covering the following aspects: 1.
Future prospects.
The first part, the physics of photon migration, provides the fundamentals of temporal propagation of photons in highly scattering media. The gamut of broadband TD-DOS instruments and advancements is discussed in the second section. The next section elaborates on various analysis methods and methodologies to extract optical properties and estimate tissue constituents. A wide range of broadband TD-DOS applications is covered across various domains (tissue constituent characterization, in vivo studies, food, wood, standardization techniques) are deeply elaborated in this section. A brief summary on the future prospects of the broadband TD-DOS is presented at the end.

Physics of Time Domain Diffuse Optical Spectroscopy (TD-DOS)
The photon propagation in a diffusive medium can be described by the radiative transfer equation (RTE), originally developed to describe the transport of non-interacting neutrons in nuclear science [54].
The assumption of non-interacting photons is still valid because, inside a diffusive medium, a multiple scattering regime is present, and consequently the light coherence vanishes after a few scattering events. Hereafter the RTE equation: 1 v ∂L ∂t = −s · ∇L(r, s, t) − (µ a + µ s )L(r, s, t) + µ s 4π p(s, s ) L(r, s, t)ds + Q(r, s, t), where L(r, s, t) is the radiance at the position r, which is the power flowing in the direction of the unit vector s, per unit of solid angle and area normal to s. v is the photon velocity in the medium, µ a and µ s is are the absorption and scattering coefficients, p(s, s ) is the scattering phase function, which gives the probability that a photon running in the direction s is scattered to the direction s , and Q(r, s, t) is the spatial, temporal, and angular distribution of the source [55].
When a short pulse is injected in a diffusive medium and detected on the external surface, Equation (1), together with proper boundary conditions, rules the distribution of the photon time-of-flight (DTOF) after injection of a short pulse. The DTOF is related both to the volume probed by the light in the medium and to its optical properties (µ a and µ s ). In fact, photons having travelled long paths give information on a larger volume than photon having travelled short paths ( Figure 2). This fact is key for detecting information from deep regions of a tissue when the source and detector are placed on the surface on the same side, which is the most used configuration for non-invasive monitoring of biological tissues [56].
where , , is the radiance at the position , which is the power flowing in the direction of the unit vector , per unit of solid angle and area normal to . is the photon velocity in the medium, and is are the absorption and scattering coefficients, , is the scattering phase function, which gives the probability that a photon running in the direction is scattered to the direction ′, and , , is the spatial, temporal, and angular distribution of the source [55].
When a short pulse is injected in a diffusive medium and detected on the external surface, Equation (1), together with proper boundary conditions, rules the distribution of the photon time-offlight (DTOF) after injection of a short pulse. The DTOF is related both to the volume probed by the light in the medium and to its optical properties ( and ). In fact, photons having travelled long paths give information on a larger volume than photon having travelled short paths ( Figure 2). This fact is key for detecting information from deep regions of a tissue when the source and detector are placed on the surface on the same side, which is the most used configuration for non-invasive monitoring of biological tissues [56]. Sensitivity profiles (log 10), which represent the relative change in the detected signal when a localized defect is placed at every position of the volume, in a homogeneous diffusive medium. At three time intervals (0.5, 1, 1.5 ns) for: (a) reflectance geometry (3 cm source-detector distance) and (b) collinear transmittance geometry. Each white line represents the contour edge of the contrast every decade. Figure 2. Sensitivity profiles (log 10), which represent the relative change in the detected signal when a localized defect is placed at every position of the volume, in a homogeneous diffusive medium. At three time intervals (0.5, 1, 1.5 ns) for: (a) reflectance geometry (3 cm source-detector distance) and (b) collinear transmittance geometry. Each white line represents the contour edge of the contrast every decade.
Another interesting advantage of the time domain approach is the fact that the shape of the DTOF is strongly modified by the absorption and scattering properties of the probed portion of the medium. Photons running a long path in the medium carry information about the probability of absorption, thus affecting the tail of the DTOF (Figure 3a). Instead, photons running a short path are more sensitive to the scattering properties of the medium, thus modifying the temporal position of the barycenter of the DTOF (Figure 3b).
These two effects open the possibility to estimate the spatial distribution of both optical properties (µ a and µ s ) simultaneously from the shape of time-resolved DTOF using analytical, numerical, or stochastic models. Another advantage of such an approach relies on the fact that the estimation is unaffected by the amplitude of the DTOF, thus rejecting every source of error related to the optical coupling between tissue and delivery and/or detection optics.
Another interesting advantage of the time domain approach is the fact that the shape of the DTOF is strongly modified by the absorption and scattering properties of the probed portion of the medium. Photons running a long path in the medium carry information about the probability of absorption, thus affecting the tail of the DTOF (Figure 3a). Instead, photons running a short path are more sensitive to the scattering properties of the medium, thus modifying the temporal position of the barycenter of the DTOF (Figure 3b). These two effects open the possibility to estimate the spatial distribution of both optical properties ( and simultaneously from the shape of time-resolved DTOF using analytical, numerical, or stochastic models. Another advantage of such an approach relies on the fact that the estimation is unaffected by the amplitude of the DTOF, thus rejecting every source of error related to the optical coupling between tissue and delivery and/or detection optics.

Instrumentation/Technology for Broadband TD-DOS
The key technical challenges for time domain diffuse optical spectrometry reside in the need to cover a wide spectral range with high sensitivity, down to the single-photon level, and with picosecond temporal resolution. Two different approaches can be envisaged to perform spectroscopy, which is either a sequential scheme, where the laser emitting wavelength within a relatively narrow (<10 nm) band is scanned progressively over the entire spectral range, or a parallel scheme, where the whole spectral range is acquired simultaneously. On the source side, a laser covering a wide spectral range is needed together with-for sequential wavelength scanning-a proper mechanism for filtering the desired wavelength out of a broad emission or alternatively a pulsed tunable laser. Furthermore, when the application involves biological media, a reasonably short total measurement time is also a must. Different realizations have been proposed in the course of the last decades, based on technical advancements mainly in the generation of supercontinuum laser sources and of timecorrelated single-photon counting (TCSPC) detection [57].
The very first demonstration of broadband TD-DOS was realized by exploiting the white light generated from a water cell by a focused 10 Hz amplified laser pulse and a streak-camera for detection [58]. This setup was used to acquire the first in vivo absorption and reduced scattering spectra of the female breast in the 650-820 nm range [59]. In the first realization, a table-top amplified Ti:Sapphire laser was used as an injection source based on a 76 MHz femtosecond self-locking oscillator operated at 790 nm, followed by a pulse stretcher, a first-stage amplifier with intra-cavity

Instrumentation/Technology for Broadband TD-DOS
The key technical challenges for time domain diffuse optical spectrometry reside in the need to cover a wide spectral range with high sensitivity, down to the single-photon level, and with picosecond temporal resolution. Two different approaches can be envisaged to perform spectroscopy, which is either a sequential scheme, where the laser emitting wavelength within a relatively narrow (<10 nm) band is scanned progressively over the entire spectral range, or a parallel scheme, where the whole spectral range is acquired simultaneously. On the source side, a laser covering a wide spectral range is needed together with-for sequential wavelength scanning-a proper mechanism for filtering the desired wavelength out of a broad emission or alternatively a pulsed tunable laser. Furthermore, when the application involves biological media, a reasonably short total measurement time is also a must. Different realizations have been proposed in the course of the last decades, based on technical advancements mainly in the generation of supercontinuum laser sources and of time-correlated single-photon counting (TCSPC) detection [57].
The very first demonstration of broadband TD-DOS was realized by exploiting the white light generated from a water cell by a focused 10 Hz amplified laser pulse and a streak-camera for detection [58]. This setup was used to acquire the first in vivo absorption and reduced scattering spectra of the female breast in the 650-820 nm range [59]. In the first realization, a table-top amplified Ti:Sapphire laser was used as an injection source based on a 76 MHz femtosecond self-locking oscillator operated at 790 nm, followed by a pulse stretcher, a first-stage amplifier with intra-cavity dumping for reduction of the repetition rate to 10 Hz, a second stage multi-pass amplifier, and finally a pulse compressor. The outcome was a 10 Hz train of 200 fs pulses with 450 mJ energy, which were then focused onto a water cuvette to produce a white light pulsed emission, representing a considerable achievement at that time.
In a more recent realization ( Figure 4) [60], the supercontinuum was generated into a 100 cm long photonic crystal fiber starting from a self-locked Ti:Sapphire laser. Light was injected into and collected from the sample under study using a couple of fibers for reflectance measurements or a free beam cell for transmittance measurements. Detection was achieved using an imaging spectrometer for spectral dispersion (Chromex, Model 250 IS) followed by a streak camera for temporal dispersion (Hamamatsu, Model C5680), resulting in 2D (wavelength-time) images projected onto a camera sensor.
The key advantages of such an approach were the parallel spectral acquisition and the high temporal Appl. Sci. 2019, 9, 5465 6 of 38 resolution (ps) offered by the streak-camera. Also, the almost continuous acquisition over time and wavelength permitted adjustable optimization by slicing temporal and spectral photon distributions at post-processing. As limitations, we can quote the far-from-single-photon sensitivity of the detection, the non-linearity of both the temporal and spectral axes, the limited spectral range imposed by the poor quantum efficiency of the S1 photocathode in the streak-camera under use, and in general the cost/complexity of the workstation.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 6 of 40 dumping for reduction of the repetition rate to 10 Hz, a second stage multi-pass amplifier, and finally a pulse compressor. The outcome was a 10 Hz train of 200 fs pulses with 450 mJ energy, which were then focused onto a water cuvette to produce a white light pulsed emission, representing a considerable achievement at that time. In a more recent realization ( Figure 4) [60], the supercontinuum was generated into a 100 cm long photonic crystal fiber starting from a self-locked Ti:Sapphire laser. Light was injected into and collected from the sample under study using a couple of fibers for reflectance measurements or a free beam cell for transmittance measurements. Detection was achieved using an imaging spectrometer for spectral dispersion (Chromex, Model 250 IS) followed by a streak camera for temporal dispersion (Hamamatsu, Model C5680), resulting in 2D (wavelength-time) images projected onto a camera sensor. The key advantages of such an approach were the parallel spectral acquisition and the high temporal resolution (ps) offered by the streak-camera. Also, the almost continuous acquisition over time and wavelength permitted adjustable optimization by slicing temporal and spectral photon distributions at post-processing. As limitations, we can quote the far-from-single-photon sensitivity of the detection, the non-linearity of both the temporal and spectral axes, the limited spectral range imposed by the poor quantum efficiency of the S1 photocathode in the streak-camera under use, and in general the cost/complexity of the workstation.
A different concept was based on mode-locked lasers, which were sequentially tuned and permitted the use of single-channel TCSPC detection [57], as shown in Figure 5 [22]. As sources, two tunable mode-locked lasers were used, namely, a synchronously pumped dye laser in the 610-700 nm range and an actively mode-locked Ti:Sapphire laser in the 700-1050 nm range. For detection, a microchannel plate photomultiplier (MCP) with S1 surface (R1564U, Hamamatsu, Japan) was used, followed by a chain for TCSPC (SPC130, Becker & Hickl, Germany). The key complexity of this system resided in the handling of the laser sources. Both of them were tunable lasers whose spectral scanning required delicate alignment to provide optimal lasing over the entire range but also automatic control-operated via mechanical actuators-for rotating the intracavity Lyot's filters. In addition, the Ti:Sapphire laser demanded automatic switching of the output mirror to avoid spurious laser emission as well as a fine adjustment of the cavity length to compensate for chromatic dispersion. A replica of the instrument response function (IRF) was implemented to compensate for variations in pulse shape and detector response during the wavelength tuning and to provide real-time feedback for automatic optimization of the cavity length. This reference IRF was achieved by splitting a small portion of the laser that was directly fed to the MCP. Despite the inherent complexity of the A different concept was based on mode-locked lasers, which were sequentially tuned and permitted the use of single-channel TCSPC detection [57], as shown in Figure 5 [22]. As sources, two tunable mode-locked lasers were used, namely, a synchronously pumped dye laser in the 610-700 nm range and an actively mode-locked Ti:Sapphire laser in the 700-1050 nm range. For detection, a microchannel plate photomultiplier (MCP) with S1 surface (R1564U, Hamamatsu, Japan) was used, followed by a chain for TCSPC (SPC130, Becker & Hickl, Germany). The key complexity of this system resided in the handling of the laser sources. Both of them were tunable lasers whose spectral scanning required delicate alignment to provide optimal lasing over the entire range but also automatic control-operated via mechanical actuators-for rotating the intracavity Lyot's filters. In addition, the Ti:Sapphire laser demanded automatic switching of the output mirror to avoid spurious laser emission as well as a fine adjustment of the cavity length to compensate for chromatic dispersion. A replica of the instrument response function (IRF) was implemented to compensate for variations in pulse shape and detector response during the wavelength tuning and to provide real-time feedback for automatic optimization of the cavity length. This reference IRF was achieved by splitting a small portion of the laser that was directly fed to the MCP. Despite the inherent complexity of the system, the fully automated operation was achieved, permitting the first acquisition of in vivo absorption and scattering spectra of various biological tissues by time domain diffuse optics over a wide spectral range (600-1000 nm) [61]. Surely, the adoption of TCSPC permitted to reach high sensitivity and dynamic range, and cover the critical region of high water absorption (above 900 nm), where signal attenuation is dramatically increased as compared to the 700-800 nm range. system, the fully automated operation was achieved, permitting the first acquisition of in vivo absorption and scattering spectra of various biological tissues by time domain diffuse optics over a wide spectral range (600-1000 nm) [61]. Surely, the adoption of TCSPC permitted to reach high sensitivity and dynamic range, and cover the critical region of high water absorption (above 900 nm), where signal attenuation is dramatically increased as compared to the 700-800 nm range. The introduction of compact commercial supercontinuum sources based on fiber amplified lasers coupled to photonic crystal fibers [62][63][64] and of rugged silicon photomultipliers [65] usable in diffuse optics [43,66] made it possible to realize movable rack-based systems suitable for use in clinical settings, as shown in Figure 6 [23,67]. Automatic wavelength selection of the supercontinuum source (SC450, Fianium, UK) was achieved by dispersing the light using a rotating prism and selecting a narrow band with a slit. The system was mounted on a cart and permitted non-invasive clinical studies on the application of broadband diffuse optics to bone tissues [29] and thyroid [68]. The introduction of compact commercial supercontinuum sources based on fiber amplified lasers coupled to photonic crystal fibers [62][63][64] and of rugged silicon photomultipliers [65] usable in diffuse optics [43,66] made it possible to realize movable rack-based systems suitable for use in clinical settings, as shown in Figure 6 [23,67]. Automatic wavelength selection of the supercontinuum source (SC450, Fianium, UK) was achieved by dispersing the light using a rotating prism and selecting a narrow band with a slit. The system was mounted on a cart and permitted non-invasive clinical studies on the application of broadband diffuse optics to bone tissues [29] and thyroid [68].
The implementation of a TCSPC detection scheme-granting single-photon sensitivity, but operated in a parallel configuration-required the use of multianode detectors and multichannel TCSPC boards. A realization is shown in Figure 7 [69]. The source was a supercontinuum generated in a photonic crystal fiber (NL-2.4-800, Blaze Photonics, UK) by a self-locking mode-locked Ti:Sapphire oscillator (Tissa 50, CDP Systems Corp., Russian Federa-tion). Diffusely re-emitted light was dispersed by an imaging spectrometer (SP-2150, Acton Research, USA) and detected by a 16-channel multianode PMT (R5900U-01-L16, Hamamatsu, Japan) connected to a TCSPC board (SPC-630, Becker & Hickl, Germany). Routing bits, associated with the PMT channels hit by the impinging photons were processed by the TCSPC router electronics (PML-16, Becker & Hickl, Germany) to separate signals arising from different wavelengths. This configuration enabled the fully parallel acquisition and permitted to track dynamic variations in absorption and scattering spectra following hemodynamic changes. The key constrains were a limited spectral range (520-850 nm) imposed by the available multichannel PMTs as well as reduced dynamics due to the need to accommodate all wavelengths with different signal attenuation in the same acquisition. Appl. Sci. 2019, 9, x FOR PEER REVIEW 8 of 40 The implementation of a TCSPC detection scheme-granting single-photon sensitivity, but operated in a parallel configuration-required the use of multianode detectors and multichannel TCSPC boards. A realization is shown in Figure 7 [70]. The source was a supercontinuum generated in a photonic crystal fiber (NL-2.4-800, Blaze Photonics, UK) by a self-locking mode-locked Ti:Sapphire oscillator (Tissa 50, CDP Systems Corp., Russian Federa-tion). Diffusely re-emitted light was dispersed by an imaging spectrometer (SP-2150, Acton Research, USA) and detected by a 16channel multianode PMT (R5900U-01-L16, Hamamatsu, Japan) connected to a TCSPC board (SPC-630, Becker & Hickl, Germany). Routing bits, associated with the PMT channels hit by the impinging photons were processed by the TCSPC router electronics (PML-16, Becker & Hickl, Germany) to separate signals arising from different wavelengths. This configuration enabled the fully parallel acquisition and permitted to track dynamic variations in absorption and scattering spectra following hemodynamic changes. The key constrains were a limited spectral range (520-850 nm) imposed by the available multichannel PMTs as well as reduced dynamics due to the need to accommodate all wavelengths with different signal attenuation in the same acquisition. Using a similar scheme and a compact detection stage with an integrated spectrometer and multichannel detection (PML Spec, Becker & Hickl, Germany), it was possible to construct a clinically compatible parallel broadband diffuse optical spectrometer as shown in Figure 8 [70]. The parallel acquisition permitted to follow the in vivo dynamics of Indiocynanine Green, a contrast agent approved for clinical use, which is used for tracking blood content, as well as extravasation in leaking vasculature as for brain injury [71][72][73][74][75] or tumors [76,77]. As a pilot study, the disentanglement of intracerebral from extracerebral absorption changes on phantoms and healthy volunteers was demonstrated [70]. Using a similar scheme and a compact detection stage with an integrated spectrometer and multichannel detection (PML Spec, Becker & Hickl, Germany), it was possible to construct a clinically compatible parallel broadband diffuse optical spectrometer as shown in Figure 8 [70]. The parallel acquisition permitted to follow the in vivo dynamics of Indiocynanine Green, a contrast agent approved for clinical use, which is used for tracking blood content, as well as extravasation in leaking vasculature as for brain injury [71][72][73][74][75] or tumors [76,77]. As a pilot study, the disentanglement of intracerebral from extracerebral absorption changes on phantoms and healthy volunteers was demonstrated [70].

Analysis
A variety of approaches have been applied for solving and approximating the solution of the RTE. These can be divided into three main categories: i) Analytical, ii) Numerical, iii) Stochastic solutions. The recent advancements in these categories fostered by the power given by parallel computing (graphical processing unit (GPU), multicore, cluster, etc.) allow the solution of the RTE in

Analysis
A variety of approaches have been applied for solving and approximating the solution of the RTE. These can be divided into three main categories: i) Analytical, ii) Numerical, iii) Stochastic solutions. The recent advancements in these categories fostered by the power given by parallel computing (graphical processing unit (GPU), multicore, cluster, etc.) allow the solution of the RTE in almost every geometry in reasonable computation time. Hereafter, a resume of these methods is presented.

Analytical Models
One of the first approximations derived from the RTE is the diffusion equation (DE) [78]. Starting from the RTE, the radiance, the source, and the phase function are decomposed into a spherical harmonic series. Depending on the number N of terms retained in the calculus, the approximation is called PN. As an example, for P1, the functions are approximated as the sum of isotropic and a slightly anisotropic term: where g is the anisotropy factor, defined as the average value of (s · s ), which is the cosine of the scattering angle. This approximation, in turns, brings to a coupled equation between Φ(r, t) (the fluence) and J(r, t) (the flux), whose solution allows one to derive the radiance, which in turn help the derivation of the signal transmitted or reflected. Solutions to the P3 approximation have been recently derived in a semi-infinite medium [79]. Moreover, by introducing the approximations of (i) µ a µ s = µ s (1 − g), (ii) an isotropic source term, and (iii) a smooth temporal variation of J(r, s, t), Equation (2) becomes: where D = 1 3µ s is the diffusion coefficient. The first Equation of (3) is the well-known diffusion equation (DE) and, together with the Fick's Law (second Equation (3)), has been widely used for analytically solving diffusion problems in a variety of geometries, such as: layered laterally infinite media, layered cylinder, sphere, parallelepiped [80][81][82][83][84][85]. The diffusion approximation was shown to work quite well in situations in which the diffusion is higher than the absorption and where the source-detector distance allows the detection of photons having undergone a number of scattering events such that the radiance L becomes almost isotropic [86]. Furthermore, the RTE equation has been recently solved in the time domain for a layered medium by Liemert et al. [80], showing a perfect agreement with Monte Carlo simulations, which is considered the gold-standard for model testing.

Numerical Models
When the complexity of the tissue to model is high, either because of a heterogeneous distribution of optical properties inside or because of the complex shape (e.g., mouse, hand, finger, brain) or both, no explicit solutions can be found to the RTE or to the DE. In this case, the typical approach consists in the numerical calculation of the solution using a finite-element approach (FEM), which is well suited for the DE [87,88]. To get the time domain solution, the typical approach is to apply the FEM to the spatial coordinates and a finite-difference scheme for the evolution in time.
The advantage of this approach is that, once the optical properties of the medium have been defined in the mesh, the solution relies on solving a series of sparse linear systems with the same matrix, which can be efficiently implemented on GPUs [88]. At each temporal step, the fluence distribution in the medium is computed for each source, subsequently, a projection operator (which does not change with time) is used to compute the signal on the detectors. The computation time depends on the finesse of the mesh, which in turn affects the matrix dimensions, and on the number of light sources. The number of detectors does not affect the computational time because it is a simple projection operation once the fluence is determined.
To further speed up the finite-difference scheme, which is the main bottleneck of this approach, a superposition of solutions obtained in the temporal frequency domain has been recently proposed [89].

Stochastic Models
Another approach, which unifies both the accuracy of the RTE and the flexibility of modeling a complex tissue, is based on Monte Carlo (MC) simulations. Provided the distribution of optical properties, the phase function, and the boundary conditions, the photon propagation is managed in a stochastic way. Basically, a probability distribution is generated for the occurrence of both the scattering and absorption events in the medium, whilst the phase function is used for sampling the direction at every scattering event [78,90]. A variety of techniques can be implemented for reducing the variance of an MC simulation. A comparison of the different approaches can be found in [91].
An interesting advantage of MC in time domain is the possibility of generating all the combinations of absorptions by using only one simulation with fixed scattering. This requires one to save, for each detected photon, its pathlength in every region of the sample. Once this data is stored, the absorption is introduced by weighting the photon contribution with the formula shown in Equation (4): where i is the pathlength of the photon inside the region i out of the N regions. The MC computation is well suited for taking advantage of massively parallel computation. In fact, each photon is independent of the others and, using a parallel approach, a multitude of photons can be launched in parallel [92][93][94][95][96][97]. The introduction of GPUs has dramatically accelerated the computation time of MC simulations, paving the way to the possibility of using them as a forward model in inverse problems [97,98].

Perturbative Approaches
A typical case-of-study in TD-DOS is the localization and characterization of a small localized perturbation of the optical properties over an otherwise homogeneous background tissue. This can be representative, as an example, of a lesion in the breast, activation in the brain, etc. Moreover, in diffuse optical tomography (DOT), which consists in the 3D reconstruction of the optical properties by means of measurements on the boundary of the tissue, a key issue is the knowledge of the sensitivity of each measurement for a given localized change of optical properties [99,100], which in turn is the signal variation due to a small localized perturbation in every internal point of the sample.
The most common approach to the perturbative problem relies on the expansion of the fluence Φ(r, t) as: where Φ 0 is the unperturbed fluence, and δΦ = δΦ abs + δΦ sca is the perturbation to the fluence due to a variation of absorption and scattering. Equation (5) leads to a DE for δΦ involving the locally perturbed optical properties δµ a and δµ s . If |δΦ| Φ 0 , the resulting approximation is named Born approximation, under whom analytical expressions can be computed in simple geometries, such as layered slab, cylinder, sphere, and parallelepiped. Nevertheless, more accurate solutions can be even found iteratively [101].
If the perturbation is strong, the Born approximation is weak, and other strategies were proposed based on higher order approximation [102][103][104][105].
Exploiting MC methods, Sassaroli et al. [106,107] introduced a fast method for computing the perturbed fluence, which is based on: (i) the storage of the pathlength of each photon both in the background medium and in the perturbation volume, and (ii) the storage of the number of scattering interactions inside the perturbation volume [106,108]. This method is based on one single simulation and a re-weighting rule for each trajectory based on the aforementioned information. The approach is exact, but the result is strongly limited by the noise of the MC simulation. Finally, an approach based on the angular decomposition of the radiance derived by MC for retrieving the sensitivity map for both absorption and scattering was proposed [47].

Multivariate Approach
The temporal curves can also be analyzed by considering multivariate tools. In their simple formulation, the above mentioned methods do not take the spectral dimension of the data into account: they evaluate each wavelength separately. A procedure involving multivariate methods could be used to process TD-DOS data across the entire wavelength range. MADSTRESS is a multivariate tool based on linear regression and on a two-dimensional (2D) interpolation procedure. This method was used to calculate the absorption and scattering coefficients of apples and fructose powder [109,110].

Final Considerations
As mentioned at the beginning of the section, nowadays the computational power has strongly accelerated the possibility to use almost all methods to solve an inverse problem, ranging from a simple homogeneous medium to a 3D tomography. The DE can be easily overcome by MC, RTE, PN solutions.
TD-DOS aims at the determination of the optical properties of homogeneous or layered tissues (e.g., brain tissue, muscle, breast, etc.) in a wide range of wavelengths. From this, by using the Beer-Lambert law and a semi-empirical scattering power law [29,111,112], it is possible to retrieve the chromophore concentrations and structural information (see Section 5.2).
To make the analysis more robust, spectrally-constrained approaches were proposed, in which all the DTOFs collected from tissue at different wavelengths were analyzed as a whole with the constrain given by the optical properties of the expected constituents (e.g., water, lipid, collagen, hemoglobin) [113,114].
When the shape of the sample is complex or the accuracy of MC is needed, it is possible to generate a look-up table (LUT) of simulations and apply the inverse problem by searching and/or interpolating results in the LUT [115,116]. This avoids the repetitive computation of MC, which, however, nowadays efficient, could not reach the computation speed of an analytical formula.
Finally, an advanced fitting method was proposed to be used when there is an a priori knowledge of the probability distribution of the unknowns (typically absorption and scattering) or uncertainties in the experimental measurement (e.g., uncertainty on the temporal calibration). It is possible to include these uncertainties in the analysis process using a Bayesian approach and thus improve the estimation of the optical properties [117,118].

Characterization of Absorbers/Tissue Constituents
Absorption spectra of various absorbers are essential for the quantitative analysis of highly scattering media. In tissue optics, the characterization of highly scattering tissue constituents like collagen, elastin, lipid, thyroglobulin using traditional spectrophotometers pose several challenges [119]. An integration sphere could uncouple the effects of scattering, but falls short to address uncoupling from other effects like fluorescence, finite boundary effects, physical structure of tissue constituents, and so on. TD-DOS can be used to characterize highly scattering absorbers [120,121]. By careful design of the measurement protocol, TD-DOS can effectively address the above-mentioned challenges to extract the absolute absorption spectra of various absorbers [119].
Collagen is a structural protein and an abundant tissue constituent present in various soft and hard tissues. Recent works have shown the correlation of collagen density with breast cancer [122,123], other studies aimed at the non-invasive characterization of collagen in bone [26,52]. The first time domain characterization of collagen tissue was presented by Taroni et al. over the range of 600-1100 nm [124,125]. Konugolu Venkata Sekar et al. extended this range to 500-1700 nm by addressing challenges related to fluorescence contamination and decrease in reduced scattering beyond 1100 nm [121]. Figure 9 shows the broadband collagen spectrum with its characteristic peaks around 910, 1030, 1200, and 1500 nm. This opens up the opportunity to explore collagen content in shortwave near infrared (SWIR) range, which has great potential for penetration depth due to low tissue scattering. however, nowadays efficient, could not reach the computation speed of an analytical formula.
Finally, an advanced fitting method was proposed to be used when there is an a priori knowledge of the probability distribution of the unknowns (typically absorption and scattering) or uncertainties in the experimental measurement (e.g., uncertainty on the temporal calibration). It is possible to include these uncertainties in the analysis process using a Bayesian approach and thus improve the estimation of the optical properties [117,118].

Characterization of Absorbers/Tissue Constituents
Absorption spectra of various absorbers are essential for the quantitative analysis of highly scattering media. In tissue optics, the characterization of highly scattering tissue constituents like collagen, elastin, lipid, thyroglobulin using traditional spectrophotometers pose several challenges [119]. An integration sphere could uncouple the effects of scattering, but falls short to address uncoupling from other effects like fluorescence, finite boundary effects, physical structure of tissue constituents, and so on. TD-DOS can be used to characterize highly scattering absorbers [120,121]. By careful design of the measurement protocol, TD-DOS can effectively address the above-mentioned challenges to extract the absolute absorption spectra of various absorbers [119]. Collagen is a structural protein and an abundant tissue constituent present in various soft and hard tissues. Recent works have shown the correlation of collagen density with breast cancer [122,123], other studies aimed at the non-invasive characterization of collagen in bone [26,52]. The first time domain characterization of collagen tissue was presented by Taroni et al. over the range of 600-1100 nm [124,125]. Konugolu Venkata Sekar et al. extended this range to 500-1700 nm by addressing challenges related to fluorescence contamination and decrease in reduced scattering Figure 9. Absorption (black solid squares) and reduced scattering (red triangles) spectra of Type I collagen (bovine Achilles tendon) [122].
Elastin is a protein responsible for the elasticity in the tissues such as skin arteries, ligaments, and lungs. The ratio of elastin to collagen content is a cursor for dermal skin aging [126] and reduced elasticity of the tissue was correlated to breast cancer [127,128]. Figure 10 shows the broadband elastin spectrum measured using elastin extracted from bovine neck ligament [129]. The key features of the elastin spectrum are peaks at 910, 1025, 1185, and 1275 nm. Elastin has similar spectral features to collagen, thus making it difficult to uncouple the two constituents in human tissue, where they often coexist. However, the peak at 1275 nm (and the relative weight of other peaks) could be used to distinguish elastin from collagen, as shown by a close comparison of the two absorption spectra [129].
Appl. Sci. 2019, 9, x FOR PEER REVIEW 14 of 40 beyond 1100 nm [121]. Figure 9 shows the broadband collagen spectrum with its characteristic peaks around 910, 1030, 1200, and 1500 nm. This opens up the opportunity to explore collagen content in shortwave near infrared (SWIR) range, which has great potential for penetration depth due to low tissue scattering. Elastin is a protein responsible for the elasticity in the tissues such as skin arteries, ligaments, and lungs. The ratio of elastin to collagen content is a cursor for dermal skin aging [126] and reduced elasticity of the tissue was correlated to breast cancer [127,128]. Figure 10 shows the broadband elastin spectrum measured using elastin extracted from bovine neck ligament [129]. The key features of the elastin spectrum are peaks at 910, 1025, 1185, and 1275 nm. Elastin has similar spectral features to collagen, thus making it difficult to uncouple the two constituents in human tissue, where they often coexist. However, the peak at 1275 nm (and the relative weight of other peaks) could be used to distinguish elastin from collagen, as shown by a close comparison of the two absorption spectra [129]. Lipids play a crucial role in the onset and progression of various fat-related diseases like obesity, metabolic disorder, lifestyle-related morbidity [130,131]. Hence, the optical characterization of lipids could enable their non-invasive quantification for in vivo studies with various diagnostic aims. The time domain broadband (600-1100 nm) diffuse optical characterization of lipid was initially performed on purified pork fat, as reported in Figure 10a [125]. The extension of the characterization Lipids play a crucial role in the onset and progression of various fat-related diseases like obesity, metabolic disorder, lifestyle-related morbidity [130,131]. Hence, the optical characterization of lipids could enable their non-invasive quantification for in vivo studies with various diagnostic aims. The time domain broadband (600-1100 nm) diffuse optical characterization of lipid was initially performed on purified pork fat, as reported in Figure 10a [125]. The extension of the characterization to the SWIR window is shown in Figure 11b, where the spectrum measured using TD-DOS [132] is overlaid on the one obtained with CW DOS performed at a small source-detector distance on beef fat. Combining the results in the two spectral ranges, one could obtain a broadband spectrum across the 600-1700 nm range. In the NIR range, an intense and sharp 930 nm peak is unique to lipid. Features of lipids in the SWIR range slightly overlap with the collagen spectrum. However, on careful notice, the lipid peak at 1400 nm is less intense and red-shifted as compared to the collagen peak. Recent studies on the near infrared (NIR) characterization of human thyroid in vivo revealed some unique features [53,135,136]. The quantification of tissue constituents using typical chromophores (lipid, water, collagen, oxy, and de-oxy hemoglobin) resulted in an overestimation of some tissue constituents. This led to the hypothesis of missing constituents that are unique to thyroid tissue. Figure 12 shows the characterization of thyroid-specific tissue constituents: thyroglobulin, tyrosine, and iodine [119]. Thyroglobulin, the most abundant protein in thyroid tissue, has a monotonous decreasing absorption spectrum with a shoulder and a peak, respectively, at 1000 and 1190 nm. Multiple features are found in the tyrosine spectrum with 6 peaks (875, 915, 1090, 1145, 1185, 1290 nm) spread across the broadband 600-1350 nm range. Iodine was found to have an exponential decaying absorption spectrum with minimum to no absorption at long NIR wavelengths. The inset in Figure 12 describes the anatomy of the human thyroid and the breakdown of tissue constituents. Though the potential of the in vivo characterization of the thyroid so far has only been limitedly explored, a few examples are presented in the next section. Recent studies on the near infrared (NIR) characterization of human thyroid in vivo revealed some unique features [53,68,133]. The quantification of tissue constituents using typical chromophores (lipid, water, collagen, oxy, and de-oxy hemoglobin) resulted in an overestimation of some tissue constituents. This led to the hypothesis of missing constituents that are unique to thyroid tissue. Figure 12 shows the characterization of thyroid-specific tissue constituents: thyroglobulin, tyrosine, and iodine [119]. Thyroglobulin, the most abundant protein in thyroid tissue, has a monotonous decreasing absorption spectrum with a shoulder and a peak, respectively, at 1000 and 1190 nm. Multiple features are found in the tyrosine spectrum with 6 peaks (875, 915, 1090, 1145, 1185, 1290 nm) spread across the broadband 600-1350 nm range. Iodine was found to have an exponential decaying absorption spectrum with minimum to no absorption at long NIR wavelengths. The inset in Figure 12 describes the anatomy of the human thyroid and the breakdown of tissue constituents. Though the potential of the in vivo characterization of the thyroid so far has only been limitedly explored, a few examples are presented in the next section.

In Vivo Characterization of Biological Tissues
TD-DOS combines distinctive features, such as non-invasiveness, deep tissue penetration, and capacity to provide functional information, which makes it a potentially powerful means for the characterization of biological tissues in vivo either for physiology studies or for clinical diagnostics.
Since the initial developments of DOS, researchers have shown interest in breast cancer-related applications. They include diagnosis (detection and discrimination between malignant and benign breast lesions), as well as therapy monitoring and cancer risk assessment. Each of these goals can take advantage of a thorough optical characterization of breast tissue, as can be achieved through broadband TD-DOS spectroscopy.
Breast tissue is inherently heterogeneous with two major components: fibro-glandular and adipose tissue. Inter-subject variability is very strong, spanning from fully adipose (BI-RADS category 1) to extremely dense breasts (BI-RADS category 4), where BI-RADS category is the clinical metrics to quantify "breast density", that is the fraction of radio-opaque fibro-glandular tissue in the breast, based on visual inspection of X-ray mammographic images [125]. Figure 13a shows typical examples of absorption spectra obtained in vivo from breasts of different types (i.e., composition), performing TD-DOS measurements in transmittance geometry, and interpreting the data with the solution of the diffusion equation for a homogeneous slab model. Notwithstanding the simple modeling, which is most often used to interpret breast data, the measured absorption properties are really sensitive to tissue composition. At short wavelengths, the major features are always the trailing edge of hemoglobin absorption, increasing upon decreasing wavelength, and the minor peak of deoxygenated hemoglobin (deoxyHb) around 760 nm. Instead, strong variability is observed above 900 nm, as a result of different tissue composition. The spectrum

In Vivo Characterization of Biological Tissues
TD-DOS combines distinctive features, such as non-invasiveness, deep tissue penetration, and capacity to provide functional information, which makes it a potentially powerful means for the characterization of biological tissues in vivo either for physiology studies or for clinical diagnostics.
Since the initial developments of DOS, researchers have shown interest in breast cancer-related applications. They include diagnosis (detection and discrimination between malignant and benign breast lesions), as well as therapy monitoring and cancer risk assessment. Each of these goals can take advantage of a thorough optical characterization of breast tissue, as can be achieved through broadband TD-DOS spectroscopy.
Breast tissue is inherently heterogeneous with two major components: fibro-glandular and adipose tissue. Inter-subject variability is very strong, spanning from fully adipose (BI-RADS category 1) to extremely dense breasts (BI-RADS category 4), where BI-RADS category is the clinical metrics to quantify "breast density", that is the fraction of radio-opaque fibro-glandular tissue in the breast, based on visual inspection of X-ray mammographic images [125]. Figure 13a shows typical examples of absorption spectra obtained in vivo from breasts of different types (i.e., composition), performing TD-DOS measurements in transmittance geometry, and interpreting the data with the solution of the diffusion equation for a homogeneous slab model. Notwithstanding the simple modeling, which is most often used to interpret breast data, the measured absorption properties are really sensitive to tissue composition. At short wavelengths, the major features are always the trailing edge of hemoglobin absorption, increasing upon decreasing wavelength, and the minor peak of deoxygenated hemoglobin (deoxyHb) around 760 nm. Instead, strong variability is observed above 900 nm, as a result of different tissue composition. The spectrum of the adipose breast (red symbols) is dominated by the absorption peaks of lipids around 930 and 1040 nm, while the highly fibro-glandular breast (blue symbols) is characterized by very strong water absorption around 980 nm, with collagen weakly contributing as a shoulder on the trailing edge of the water peak.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 17 of 40 of the adipose breast (red symbols) is dominated by the absorption peaks of lipids around 930 and 1040 nm, while the highly fibro-glandular breast (blue symbols) is characterized by very strong water absorption around 980 nm, with collagen weakly contributing as a shoulder on the trailing edge of the water peak. Interpreting the absorption spectra with the Beer law enables the estimate of tissue composition: oxy-and deoxyhemoglobin, water lipid, and collagen (or a subset) are thus quantified. From the two forms of hemoglobin, more commonly used physiological parameters are then derived: total hemoglobin content (tHb) and oxygenation level (SO2). Figure 13b reports the reduced scattering spectra measured from the same three subjects. The spectrum of the dense (fibrous) breast typically shows a steeper slope than that of adipose or mixed breasts. That reflects the different composition: lower lipid and higher collagen content. Actually, the reduced scattering spectrum is generally interpreted with a simple empirical approximation to Mie theory [112,113] and provides information on the microscopic structure of tissue: μs' = aλ -b , where the amplitude a is related to the number of scattering centers (cell membranes, subcellular organelles, etc.) and the slope b depends on their size, being higher for smaller scattering centers.  Table 1 displays the results obtained interpreting the absorption and reduced scattering spectra shown in Figure 13. Both tissue composition and scattering parameters related to the microscopic tissue structure can vary markedly with breast type.
Clinical applications often impose limitations on overall measurement time. Thus, it is important to determine whether tissue quantification can be achieved even operating at a limited number of wavelengths. The scattering spectra can be modeled with just two parameters (a and b), while the absorption should generally provide information on (at least) 5 constituents for thorough tissue characterization. Thus, measurements need to be performed at least at (at least) 7 distinct wavelengths, spread over a broad spectral range, including marked absorption features of all the constituents of interest (i.e., 600 to 800 nm for blood, and including the range of 900 to 1100 nm for water, lipids, and collagen). The acquisition of the full spectrum provides more robust results, but satisfactory outcomes can be obtained even from data collected at a small set of wavelengths (7), provided that global data fitting with spectral constrains is performed [113]. Actually, using that spectral approach, optically derived tissue parameters were effectively used in the discrimination Interpreting the absorption spectra with the Beer law enables the estimate of tissue composition: oxy-and deoxyhemoglobin, water lipid, and collagen (or a subset) are thus quantified. From the two forms of hemoglobin, more commonly used physiological parameters are then derived: total hemoglobin content (tHb) and oxygenation level (SO 2 ). Figure 13b reports the reduced scattering spectra measured from the same three subjects. The spectrum of the dense (fibrous) breast typically shows a steeper slope than that of adipose or mixed breasts. That reflects the different composition: lower lipid and higher collagen content. Actually, the reduced scattering spectrum is generally interpreted with a simple empirical approximation to Mie theory [112,113] and provides information on the microscopic structure of tissue: µ s ' = aλ -b , where the amplitude a is related to the number of scattering centers (cell membranes, subcellular organelles, etc.) and the slope b depends on their size, being higher for smaller scattering centers. Table 1 displays the results obtained interpreting the absorption and reduced scattering spectra shown in Figure 13. Both tissue composition and scattering parameters related to the microscopic tissue structure can vary markedly with breast type. Clinical applications often impose limitations on overall measurement time. Thus, it is important to determine whether tissue quantification can be achieved even operating at a limited number of wavelengths. The scattering spectra can be modeled with just two parameters (a and b), while the absorption should generally provide information on (at least) 5 constituents for thorough tissue characterization. Thus, measurements need to be performed at least at (at least) 7 distinct wavelengths, spread over a broad spectral range, including marked absorption features of all the constituents of interest (i.e., 600 to 800 nm for blood, and including the range of 900 to 1100 nm for water, lipids, and collagen). The acquisition of the full spectrum provides more robust results, but satisfactory outcomes can be obtained even from data collected at a small set of wavelengths (7), provided that global data fitting with spectral constrains is performed [113]. Actually, using that spectral approach, optically derived tissue parameters were effectively used in the discrimination between malignant and benign breast lesions [13], in the optical measurement of breast density (important risk factor for developing cancer [134]), and in the identification of women at high risk for their high breast density [123].
The knowledge of the optical properties of the human head is quite valuable for the development of non-invasive optical technologies that aim at imaging or monitoring brain functions or injuries. Understanding the biological variability in the optical properties of head layers (scalp-skull, gray matter, cerebrospinal fluid, white matter) is beneficial for setting up simulation scenarios, optimizing existing technologies, analyzing actual measurements with a priori information (e.g., on the optical pathlength). Figure 14 shows the results of a multilaboratory exercise where 3 different approaches-namely a multi-distance CW system, a 5-wavelength time-resolved instrument, and a broadband TD-DOS workstation-were challenged in retrieving the optical properties of the forehead of 9 subjects [135]. In particular, the broadband TD-DOS system (blue squares) was the only one apt to cover a wide spectral range (600-1100 nm), thus revealing the clear contribution of lipids and water to tissue absorption. The broadband TD-DOS data were analyzed using a homogeneous solution of the diffusion equation. Simulations highlighted that the spatial region determining the retrieved parameters is the brain for what concerns the absorption, and the skull for the scattering. In the case of the few-wavelength TD-DOS system, exploiting a multi-distance approach and longer time-per-wavelength acquisition, it was also possible to apply a two-layer model, to disentangle the scalp-skull properties (black circles) from the brain properties (gray circles). The feasibility to apply heterogeneous models is certainly of great interest also for broadband spectroscopy, with the aim to describe biological tissues more accurately and with higher informative content. Key challenges are surely the strong light attenuation above 900 nm, which limits, in particular, the availability of long-lived (>1 ns) photons carrying information from deeper structures (brain). Also, fast and accurate models of photon migration are needed to overcome the key limitations of the diffusion approximation when dealing with short (<200 ps) propagation times or thin (<2 mm) layers. In that direction, especially attractive is the proposal of analytical multi-layered models comparable with Monte Carlo in terms of accuracy, but quick enough to be used in multi-parameter inversion procedures [136,137].
Appl. Sci. 2019, 9, x FOR PEER REVIEW 18 of 40 between malignant and benign breast lesions [13], in the optical measurement of breast density (important risk factor for developing cancer [134]), and in the identification of women at high risk for their high breast density [123]. The knowledge of the optical properties of the human head is quite valuable for the development of non-invasive optical technologies that aim at imaging or monitoring brain functions or injuries. Understanding the biological variability in the optical properties of head layers (scalp-skull, gray matter, cerebrospinal fluid, white matter) is beneficial for setting up simulation scenarios, optimizing existing technologies, analyzing actual measurements with a priori information (e.g., on the optical pathlength). Figure 14 shows the results of a multilaboratory exercise where 3 different approachesnamely a multi-distance CW system, a 5-wavelength time-resolved instrument, and a broadband TD-DOS workstation-were challenged in retrieving the optical properties of the forehead of 9 subjects [135]. In particular, the broadband TD-DOS system (blue squares) was the only one apt to cover a wide spectral range (600-1100 nm), thus revealing the clear contribution of lipids and water to tissue absorption. The broadband TD-DOS data were analyzed using a homogeneous solution of the diffusion equation. Simulations highlighted that the spatial region determining the retrieved parameters is the brain for what concerns the absorption, and the skull for the scattering. In the case of the few-wavelength TD-DOS system, exploiting a multi-distance approach and longer time-perwavelength acquisition, it was also possible to apply a two-layer model, to disentangle the scalp-skull properties (black circles) from the brain properties (gray circles). The feasibility to apply heterogeneous models is certainly of great interest also for broadband spectroscopy, with the aim to describe biological tissues more accurately and with higher informative content. Key challenges are surely the strong light attenuation above 900 nm, which limits, in particular, the availability of longlived (>1 ns) photons carrying information from deeper structures (brain). Also, fast and accurate models of photon migration are needed to overcome the key limitations of the diffusion approximation when dealing with short (<200 ps) propagation times or thin (<2 mm) layers. In that direction, especially attractive is the proposal of analytical multi-layered models comparable with Monte Carlo in terms of accuracy, but quick enough to be used in multi-parameter inversion procedures [136,137]. reduced scattering spectra of the human forehead. Median values over 9 subjects are shown for data collected with different systems (CW, TD at 5 wavelengths, and broadband TD). A homogeneous model was applied to data from all systems. For the TD system at discrete wavelengths results are displayed also when the data are interpreted with a two-layer model [138].
In the last years, several other TD-DOS applications started to be considered, stemming from presently unmet clinical needs. One of them is the optical characterization of bone. The aging population suffers more and more from bone-related pathologies, but in several widespread cases (e.g., osteoporosis) effective non-invasive means are currently not available for early diagnosis, which instead is key to efficient treatment. An initial attempt towards the non-invasive optical Figure 14. (a) Absorption and (b) reduced scattering spectra of the human forehead. Median values over 9 subjects are shown for data collected with different systems (CW, TD at 5 wavelengths, and broadband TD). A homogeneous model was applied to data from all systems. For the TD system at discrete wavelengths results are displayed also when the data are interpreted with a two-layer model [138].
In the last years, several other TD-DOS applications started to be considered, stemming from presently unmet clinical needs. One of them is the optical characterization of bone. The aging population suffers more and more from bone-related pathologies, but in several widespread cases (e.g., osteoporosis) effective non-invasive means are currently not available for early diagnosis, which instead is key to efficient treatment. An initial attempt towards the non-invasive optical characterization of bone tissue was made more than a decade ago, performing measurements in the 650-1000 nm range on the calcaneus of 7 female volunteers spanning over a broad age range (26-82 y): Measurement feasibility was shown and possible age dependence of the optical properties was suggested to reflect age-related changes in composition and structure of bone tissue [26].
More recently, data were collected from 17 healthy volunteers at 6 locations where the bone is rather superficial and thus accessible to optical probing, with the aim of identifying the best location for further systematic optical studies on bone [29,138]. From Figure 15, it appears that, notwithstanding significant inter-subject variability (measured by the standard deviation), on average different locations are characterized by different optical properties. In particular, lower absorption, with smaller water contribution at 980 nm (calcaneus and trochanter) corresponds to higher scattering.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 19 of 40 characterization of bone tissue was made more than a decade ago, performing measurements in the 650-1000 nm range on the calcaneus of 7 female volunteers spanning over a broad age range (26-82 y): Measurement feasibility was shown and possible age dependence of the optical properties was suggested to reflect age-related changes in composition and structure of bone tissue [26]. More recently, data were collected from 17 healthy volunteers at 6 locations where the bone is rather superficial and thus accessible to optical probing, with the aim of identifying the best location for further systematic optical studies on bone [29,138]. From Figure 15, it appears that, notwithstanding significant inter-subject variability (measured by the standard deviation), on average different locations are characterized by different optical properties. In particular, lower absorption, with smaller water contribution at 980 nm (calcaneus and trochanter) corresponds to higher scattering.
Looking for an ideal location for clinical diagnostics of bone, three parameters were considered: the thickness of the superficial skin layer, the penetration depth of photons at the wavelength of operation, and the bone volume at the probed location. At the distal radius location, superficial tissue (skin) thickness is small, the high light penetration depth is foreseen based on the measured optical properties, and the bone volume is high, thus providing a strong contribution to the overall collected signal. All in all, these observations suggest the distal radius location as the best one for future diffuse optical studies of bone-related pathologies. Thyroid nodules are common, but presently their nature cannot be identified non-invasively. Fine needle aspiration followed by cytological analysis is the established standard to diagnose thyroid pathologies. However, besides being an invasive procedure, it may miss up to one-third of malignancies in palpable nodules [139].
A non-invasive optical diagnostics capable to estimate tissue composition and blood parameters may be of real interest. Thus, after the optical characterization of its specific constituents (tyrosine and thyroglobulin), described in Section 5.1, the in vivo characterization of the thyroid has recently started [53]. TD-DOS data were collected in reflectance geometry from 6 healthy volunteers on the thyroid area and at a reference lateral position on the neck. At the thyroid location, two sourcedetector distances d (d = 1.5 cm and 2.5 cm) were tested to identify best conditions for probing the gland (as a trade-off between signal level and minimum contribution from surrounding tissues).
These initial measurements reported in Figure 16 show that in the range of 600-1100 nm, the absorption properties at the thyroid location are dominated (around 600 nm) by the tail of hemoglobin and thyroglobulin absorption peaks at shorter wavelengths. Strong water absorption at 980 nm is the other dominant spectral feature. Finally, over the entire spectral range, the absorption is very high, spanning from 0.2 to 1.2 cm −1 , which may at least in part be attributed to thyroglobulin. Similar absorption properties are obtained at the two source-detector separations, with just slightly Looking for an ideal location for clinical diagnostics of bone, three parameters were considered: the thickness of the superficial skin layer, the penetration depth of photons at the wavelength of operation, and the bone volume at the probed location. At the distal radius location, superficial tissue (skin) thickness is small, the high light penetration depth is foreseen based on the measured optical properties, and the bone volume is high, thus providing a strong contribution to the overall collected signal. All in all, these observations suggest the distal radius location as the best one for future diffuse optical studies of bone-related pathologies.
Thyroid nodules are common, but presently their nature cannot be identified non-invasively. Fine needle aspiration followed by cytological analysis is the established standard to diagnose thyroid pathologies. However, besides being an invasive procedure, it may miss up to one-third of malignancies in palpable nodules [139].
A non-invasive optical diagnostics capable to estimate tissue composition and blood parameters may be of real interest. Thus, after the optical characterization of its specific constituents (tyrosine and thyroglobulin), described in Section 5.1, the in vivo characterization of the thyroid has recently started [53]. TD-DOS data were collected in reflectance geometry from 6 healthy volunteers on the thyroid area and at a reference lateral position on the neck. At the thyroid location, two source-detector distances d (d = 1.5 cm and 2.5 cm) were tested to identify best conditions for probing the gland (as a trade-off between signal level and minimum contribution from surrounding tissues).
These initial measurements reported in Figure 16 show that in the range of 600-1100 nm, the absorption properties at the thyroid location are dominated (around 600 nm) by the tail of hemoglobin and thyroglobulin absorption peaks at shorter wavelengths. Strong water absorption at 980 nm is the other dominant spectral feature. Finally, over the entire spectral range, the absorption is very high, spanning from 0.2 to 1.2 cm −1 , which may at least in part be attributed to thyroglobulin. Similar absorption properties are obtained at the two source-detector separations, with just slightly lower absorption below 800 nm for the shorter source-detector distance. In agreement with what expected based on thyroid physiology, the absorption at the reference location is on average lower over the entire spectral range, reflecting lower blood content and no thyroid-specific contributions.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 20 of 40 lower absorption below 800 nm for the shorter source-detector distance. In agreement with what expected based on thyroid physiology, the absorption at the reference location is on average lower over the entire spectral range, reflecting lower blood content and no thyroid-specific contributions. The measurements on the thyroid with d = 1.5 cm resulted in higher reduced scattering than other measurements.
It is important to highlight that just a few subjects were analyzed and that the inter-subject variability is marked. Still, this preliminary study showed that the non-invasive characterization of the thyroid by optical means is feasible and could effectively be further pursued. Finally, it is interesting to consider adipose tissue. It has long been regarded just as an energy depot, while it is now becoming evident that it is an endocrine organ, deeply involved in the development of many pathologies, like obesity, diabetes, and metabolic diseases that are growing in impact on the quality of life and in burden on the healthcare systems [131,140]. Visceral fat, located mainly within the abdominal cavity has an important link with pathologies. Thus, the non-invasive in vivo characterization of abdominal fat and the capacity to monitor its pathologic changes may certainly be of interest for pathophysiology studies.
An interesting work on the optical assessment of abdominal adipose tissue was recently performed with a different approach to diffuse optics (a combination of frequency domain measurements at selected wavelengths and broadband CW measurements interpreted with a homogeneous model of the abdomen). That study suggested that DOS may be of help in monitoring the effects of a low-calorie diet on the obese subjects, and more generally in investigating the role of adipose tissue in metabolic disorders [141].
However, TD-DOS data collected at 610-1010 nm from on the abdomen of 3 healthy volunteers and interpreted with the diffusion approximation for a homogeneous semi-infinite medium showed that, depending on the subject, the absorption spectrum can be dominated by lipid absorption (as expected when the adipose tissue layer is properly probed) or by water and blood, when the subcutaneous adipose tissue is thin and the underlying muscle is reached [142]. That clearly evidenced the need to model the abdomen as a multi-layer structure, when aiming at the quantification of the optical properties and of the composition of adipose tissue. Further data in support of this need were recently obtained performing TD-DOS on 4 healthy volunteers at 3 different source-detector separations d (d = 1, 2, 3 cm), probing increasing depths. As evident from Figure 17, the measured absorption properties markedly depend on depth [143]. With d = 3 cm, when the adipose tissue layer is thick (>2 cm), it dominates the estimated absorption properties, which clearly resemble the spectral features of lipid. Instead, for thinner subcutaneous adipose tissue, the The measurements on the thyroid with d = 1.5 cm resulted in higher reduced scattering than other measurements.
It is important to highlight that just a few subjects were analyzed and that the inter-subject variability is marked. Still, this preliminary study showed that the non-invasive characterization of the thyroid by optical means is feasible and could effectively be further pursued.
Finally, it is interesting to consider adipose tissue. It has long been regarded just as an energy depot, while it is now becoming evident that it is an endocrine organ, deeply involved in the development of many pathologies, like obesity, diabetes, and metabolic diseases that are growing in impact on the quality of life and in burden on the healthcare systems [131,140]. Visceral fat, located mainly within the abdominal cavity has an important link with pathologies. Thus, the non-invasive in vivo characterization of abdominal fat and the capacity to monitor its pathologic changes may certainly be of interest for pathophysiology studies.
An interesting work on the optical assessment of abdominal adipose tissue was recently performed with a different approach to diffuse optics (a combination of frequency domain measurements at selected wavelengths and broadband CW measurements interpreted with a homogeneous model of the abdomen). That study suggested that DOS may be of help in monitoring the effects of a low-calorie diet on the obese subjects, and more generally in investigating the role of adipose tissue in metabolic disorders [141].
However, TD-DOS data collected at 610-1010 nm from on the abdomen of 3 healthy volunteers and interpreted with the diffusion approximation for a homogeneous semi-infinite medium showed that, depending on the subject, the absorption spectrum can be dominated by lipid absorption (as expected when the adipose tissue layer is properly probed) or by water and blood, when the subcutaneous adipose tissue is thin and the underlying muscle is reached [142]. That clearly evidenced the need to model the abdomen as a multi-layer structure, when aiming at the quantification of the optical properties and of the composition of adipose tissue. Further data in support of this need were recently obtained performing TD-DOS on 4 healthy volunteers at 3 different source-detector separations d (d = 1, 2, 3 cm), probing increasing depths. As evident from Figure 17, the measured absorption properties markedly depend on depth [143]. With d = 3 cm, when the adipose tissue layer is thick (>2 cm), it dominates the estimated absorption properties, which clearly resemble the spectral features of lipid. Instead, for thinner subcutaneous adipose tissue, the absorption properties also depend on water (at 980 nm) and blood (at short wavelengths) of the underlying muscle. On the other end, the contribution of the skin also needs to be considered, as it also shows marked intersubject variability, and possibly age dependence. Thus, a layered model, accounting for skin, adipose tissue layer and underlying muscle, is needed to allow a proper investigation of the adipose tissue. absorption properties also depend on water (at 980 nm) and blood (at short wavelengths) of the underlying muscle. On the other end, the contribution of the skin also needs to be considered, as it also shows marked intersubject variability, and possibly age dependence. Thus, a layered model, accounting for skin, adipose tissue layer and underlying muscle, is needed to allow a proper investigation of the adipose tissue.

Wood
An interesting field for diffuse optics with potential industrial applications is the nondestructive spectroscopy of wood. This field has been explored at large by different NIR techniques directly focused on non-destructive retrieval of wood composition or conservation status [144][145][146][147][148].
The key advantage of broadband TD-DOS in this field is the intrinsic capability of disentangling absorption from scattering contributions, resulting in a more accurate estimate of the optical properties. An example of an application is the monitoring of archeological wood and in general of age-related degradation of wood artifacts, as shown in Figure 18 [149]. A clear reduction in the scattering coefficient is observed upon wetting-as for waterlogged wood-and even further by artificially procured degradation, due to the decrease of refractive index mismatch among cells and disruption of the cellular structure. The absorption spectra clearly display the water contribution in wet samples (around 980 nm) as well as the increased absorption of oxidative products in the region <800 nm. The strong optical anisotropy of wood materials is evident from the large difference in the scattering coefficient using a source-detector pair parallel (up) or perpendicular (down) with respect to the wood fibers.
The capability to detect water content within the diffusive medium can be exploited for the monitoring of wood moisture-a key parameter with relevance to timber industrial applications [150]. Figure 19 shows a sample specimen kept in a chamber under controlled humidity and continuously monitored using broadband TD-DOS. Upon increasing the humidity level-and therefore the moisture content in wood-the water absorption peak around 980 nm increases correspondingly, demonstrating a sensitivity of the technique down to a few percents of moisture content.

Wood
An interesting field for diffuse optics with potential industrial applications is the non-destructive spectroscopy of wood. This field has been explored at large by different NIR techniques directly focused on non-destructive retrieval of wood composition or conservation status [144][145][146][147][148]. The key advantage of broadband TD-DOS in this field is the intrinsic capability of disentangling absorption from scattering contributions, resulting in a more accurate estimate of the optical properties. An example of an application is the monitoring of archeological wood and in general of age-related degradation of wood artifacts, as shown in Figure 18 [149]. A clear reduction in the scattering coefficient is observed upon wetting-as for waterlogged wood-and even further by artificially procured degradation, due to the decrease of refractive index mismatch among cells and disruption of the cellular structure. The absorption spectra clearly display the water contribution in wet samples (around 980 nm) as well as the increased absorption of oxidative products in the region <800 nm. The strong optical anisotropy of wood materials is evident from the large difference in the scattering coefficient using a source-detector pair parallel (up) or perpendicular (down) with respect to the wood fibers.
The capability to detect water content within the diffusive medium can be exploited for the monitoring of wood moisture-a key parameter with relevance to timber industrial applications [150]. Figure 19 shows a sample specimen kept in a chamber under controlled humidity and continuously monitored using broadband TD-DOS. Upon increasing the humidity level-and therefore the moisture content in wood-the water absorption peak around 980 nm increases correspondingly, demonstrating a sensitivity of the technique down to a few percents of moisture content. , and finally to degraded sample (black diamonds). Likewise, the absorption spectrum (left) is increased for the wet sample in the water-absorbing region (around 980 nm) and for the degraded sample in the red region (<800 nm) due to oxidation products. The strong optical anisotropy of wood materials is evident from the large difference in the scattering coefficient in the configuration parallel (up) or perpendicular (down) to the wood fibers [149].  , and finally to degraded sample (black diamonds). Likewise, the absorption spectrum (left) is increased for the wet sample in the water-absorbing region (around 980 nm) and for the degraded sample in the red region (<800 nm) due to oxidation products. The strong optical anisotropy of wood materials is evident from the large difference in the scattering coefficient in the configuration parallel (up) or perpendicular (down) to the wood fibers [149].
Appl. Sci. 2019, 9, x FOR PEER REVIEW 22 of 40 Figure 18. Monitoring the effect of degradation on wood samples measured along the fibers (up), perpendicular to them (down). The scattering spectrum (right) is progressively decreased moving from the dry sample (open diamonds), to wet sample (gray diamonds), and finally to degraded sample (black diamonds). Likewise, the absorption spectrum (left) is increased for the wet sample in the water-absorbing region (around 980 nm) and for the degraded sample in the red region (<800 nm) due to oxidation products. The strong optical anisotropy of wood materials is evident from the large difference in the scattering coefficient in the configuration parallel (up) or perpendicular (down) to the wood fibers [149].  A further attractive industrial application was pursued in the quantification of monomer uptake in wood polymer composites [151]. The key interest derives from recent attempts to generate structurally strong wood materials starting from low-quality softwood structures soaked in a (possibly biodegradable) monomer undergoing polymerization. Figure 20a shows the alteration in the absorption spectrum of wood caused by the uptake of a mix of monomers. The wood specimen was immersed in a solution containing the monomer mix (45% glycidyl methacrylate, 45% methylmethacrylate, and 10% ethylene glycol dimethacrylate) and uptake was forced applying mild vacuum to the chamber, while the sample was continuously monitored by broadband TD-DOS. The difference between the treated and untreated wood corresponds to the absorption spectrum of the monomer mix, thus permitting absolute quantitation. Figure 20b demonstrates capabilities to monitor the monomer uptake non-destructively across a 1 cm thick wood sample at different positions on the surface along the specimen with a precision better than 1% [152].
Appl. Sci. 2019, 9, x FOR PEER REVIEW 23 of 40 A further attractive industrial application was pursued in the quantification of monomer uptake in wood polymer composites [151]. The key interest derives from recent attempts to generate structurally strong wood materials starting from low-quality softwood structures soaked in a (possibly biodegradable) monomer undergoing polymerization. Figure 20a shows the alteration in the absorption spectrum of wood caused by the uptake of a mix of monomers. The wood specimen was immersed in a solution containing the monomer mix (45% glycidyl methacrylate, 45% methylmethacrylate, and 10% ethylene glycol dimethacrylate) and uptake was forced applying mild vacuum to the chamber, while the sample was continuously monitored by broadband TD-DOS. The difference between the treated and untreated wood corresponds to the absorption spectrum of the monomer mix, thus permitting absolute quantitation. Figure 20b demonstrates capabilities to monitor the monomer uptake non-destructively across a 1 cm thick wood sample at different positions on the surface along the specimen with a precision better than 1% [152].  [151]. The concentration of the monomer at different positions along the wood sample (b) demonstrates a precision in monomer assessment better than 1% in mass [152].

Food
Besides wood, also food products, and in particular fruits, have been investigated using optical techniques (mainly NIR spectroscopy) [2,153,154], and in particular broadband TD-DOS. Figure 21 displays the absorption and reduced scattering spectra of different fruits and vegetables [155]. The absorption peak of chlorophyll around 670 nm and of water around 980 nm are clearly visible. In the case of kiwifruit, chlorophyll content is so high that only the falling edge of the 670 nm chlorophyll peak is measurable. That initial broadband spectroscopy of food products by TD-DOS permitted to identify the best operating wavelengths for compact off-the-field systems applied to the nondestructive quality assessment of food [156].
Differently from standard NIR spectroscopy, TD-DOS is largely independent of the optical contact or surface pigmentation. Figure 22 clearly shows that the absorption (left) and reduced scattering (right) spectra of different varieties of apples are not affected by skin removal [157]. This property is due to the physics of TD-DOS, since the photon temporal distribution is determined by the combined effect of scattering (creating a spreading of photon paths) and absorption (ruling the survival probability of each random path), as discussed and shown in Section 2. Therefore, thin superficial layers or simply source or contact instabilities will simply change the signal amplitude without affecting the temporal shape, therefore with no impact whatsoever on the retrieved optical properties.
Another field of food quality control and food safety is in the area of dairy products. Several studies have been conducted in this area using photon time-of-flight NIR spectroscopy, enabling independent measurements of the NIR absorption and scattering spectra [158][159][160][161].  [151]. The concentration of the monomer at different positions along the wood sample (b) demonstrates a precision in monomer assessment better than 1% in mass [152].

Food
Besides wood, also food products, and in particular fruits, have been investigated using optical techniques (mainly NIR spectroscopy) [2,153,154], and in particular broadband TD-DOS. Figure 21 displays the absorption and reduced scattering spectra of different fruits and vegetables [155]. The absorption peak of chlorophyll around 670 nm and of water around 980 nm are clearly visible. In the case of kiwifruit, chlorophyll content is so high that only the falling edge of the 670 nm chlorophyll peak is measurable. That initial broadband spectroscopy of food products by TD-DOS permitted to identify the best operating wavelengths for compact off-the-field systems applied to the non-destructive quality assessment of food [156].
Differently from standard NIR spectroscopy, TD-DOS is largely independent of the optical contact or surface pigmentation. Figure 22 clearly shows that the absorption (left) and reduced scattering (right) spectra of different varieties of apples are not affected by skin removal [157]. This property is due to the physics of TD-DOS, since the photon temporal distribution is determined by the combined effect of scattering (creating a spreading of photon paths) and absorption (ruling the survival probability of each random path), as discussed and shown in Section 2. Therefore, thin superficial layers or simply source or contact instabilities will simply change the signal amplitude without affecting the temporal shape, therefore with no impact whatsoever on the retrieved optical properties.     Another field of food quality control and food safety is in the area of dairy products. Several studies have been conducted in this area using photon time-of-flight NIR spectroscopy, enabling independent measurements of the NIR absorption and scattering spectra [158][159][160][161].

Pharmaceuticals
The quality control in the production of pharmaceutical tablets is strict, with tolerances of the amount of active pharmaceutical ingredient (API) of 1%. To secure that the production is providing pharmaceutical tablets with an API content within this narrow interval, typically random tablets are taken out in the production line for wet chemistry analysis. This analysis is time-consuming and labor intense. With results of the analysis available first after a long delay. If the tablet does not fulfill the quality control, it means that large batches of tablets have to be discarded due to the long analysis delay. Direct in-line optical spectroscopic analysis is considered a better technique, provided that the accuracy of the measurements can be sufficiently good. DOS, as well as Raman transmission spectroscopy [162,163], are studied for such analysis. A complication is the highly scattering character of the tablets, providing a very long and spread path of the light as it propagates through a tablet. The scattering is caused by the porous character of the tablet, with changes in the refractive index as light propagates through the tablet. Any alterations in the size distribution of the filling material or the compression force in the production of the tablet will alter the light path as it passes through the tablet. This will cause alterations in measured DOS or Raman spectra, making precise absolute concentration measurements challenging.
The biophotonics team at Lund University had a long-lasting collaboration with AstraZeneca in Molndal, Sweden to assess the potential of near-infrared photon time-of-flight spectroscopy (pTOF or TD-DOS) as a more sophisticated version of near-infrared spectroscopy (NIRS) to allow absorption spectra to be measured independent of the scattering properties of the tablet. Figure 23 shows the extraction of absorption and reduced scattering spectra of tablets with different scattering, showing that the scattering contribution is effectively eliminated by TD-DOS. In this way, the uncertainties caused by variations in scattering properties could be eliminated in the assessment of API concentration, and the scattering properties could be linked to the matrix properties of the tablet, related to releasing properties of the API, including the delay and duration of the release after swallowing the tablet. These measurements could be complemented by NIR gas spectroscopy to further asses the matrix properties of the tablets [64,158,159,[162][163][164][165][166][167][168][169][170][171][172][173].

Pharmaceuticals
The quality control in the production of pharmaceutical tablets is strict, with tolerances of the amount of active pharmaceutical ingredient (API) of 1%. To secure that the production is providing pharmaceutical tablets with an API content within this narrow interval, typically random tablets are taken out in the production line for wet chemistry analysis. This analysis is time-consuming and labor intense. With results of the analysis available first after a long delay. If the tablet does not fulfill the quality control, it means that large batches of tablets have to be discarded due to the long analysis delay. Direct in-line optical spectroscopic analysis is considered a better technique, provided that the accuracy of the measurements can be sufficiently good. DOS, as well as Raman transmission spectroscopy [162,163], are studied for such analysis. A complication is the highly scattering character of the tablets, providing a very long and spread path of the light as it propagates through a tablet. The scattering is caused by the porous character of the tablet, with changes in the refractive index as light propagates through the tablet. Any alterations in the size distribution of the filling material or the compression force in the production of the tablet will alter the light path as it passes through the tablet. This will cause alterations in measured DOS or Raman spectra, making precise absolute concentration measurements challenging. The biophotonics team at Lund University had a long-lasting collaboration with AstraZeneca in Molndal, Sweden to assess the potential of near-infrared photon time-of-flight spectroscopy (pTOF or TD-DOS) as a more sophisticated version of near-infrared spectroscopy (NIRS) to allow absorption spectra to be measured independent of the scattering properties of the tablet. Figure 23 shows the

Standardization and Performances Assessment
In the biophotonics field, the need is recognized for shared procedures for performance assessment of newly developed instruments to ensure the reproducibility of results and reliable comparisons between different instruments [170]. Diffuse optics is pursuing this approach [49][50][51] as it is undergoing great technology progresses, like miniaturization of pulsed laser sources [171], an increase of detector light-harvesting capability [172], high throughput timing electronics [173][174][175], and innovative measurement technique [176][177][178]. In the medical imaging field, in particular, it is crucial to enable a solid validation of new or upgraded instruments in laboratory settings, performing measurements on tissue-mimicking phantoms instead of on humans and anticipating possible technical issues that could occur during expensive clinical trials. Consequently, important human and financial investments are being spread, leading to the use of many different materials for the fabrication of phantoms. Phantoms must rely on well-defined recipes, devised thanks to a proper knowledge of the spectral characteristics of the different phantom constituents, and on precise final characterizations to verify their targeted optical properties. Therefore, the role of TD-DOS facilities is crucial. Our purpose here is not to provide an exhaustive review of all the different materials involved in this process, but to make examples where broadband TD-DOS allowed one to solve specific issues.
Ref. [179] presents a phantom kit composed of a homogeneous epoxy resin bulk phantom with a hole, in which different epoxy resin rods can move along one direction (Figure 24a). Each rod contains a different absorption perturbation [180], thus allowing to move the inhomogeneity with respect to the position of the source-detector pair. The kit has been fabricated using titanium dioxide particles to obtain the targeted µ s and black toner powder to obtain the desired µ a . The optical properties of the bulk phantom are shown in Figure 24b. The absorption coefficient is quite flat up to 1100 nm, with a maximum difference in the range of about 0.05 cm −1 due to faint absorption peaks proper of epoxy resin. At longer wavelengths, larger absorption peaks of the epoxy resin matrix are clearly visible.   In some cases, even if perishable within a few weeks, the high flexibility of liquid phantoms is required to finely change the optical properties during measurement or to freely move various perturbations. For this reason, phantoms based on water dilutions of Intralipid-20% and India ink have been widely employed. In a multi-center study, the optical properties of these phantoms were characterized in a broad wavelength range (between 633 and 916 nm) [181]. Reference values were identified for the µ a of India ink and the µ s of Intralipid-20%, with uncertainty <2%, providing a solid base for the preparation of liquid phantoms.
Multimodal imaging is emerging as a novel approach to improve the diagnostic performance of a single medical examination by combining the strengths of different measurement techniques. In the field of diffuse optics, there is quite a strong effort in combining, for instance, ultrasound imaging with diffuse optics [182,183] or different diffuse optics techniques [184]. An example is the combination of DOS and diffuse correlation spectroscopy (DCS) [1]. DCS exploits the autocorrelation function of light speckles intensity, which depends on the speed of moving scattering particles (like red blood cells) to measure tissue microvascular blood flow. Hence, it can provide a complimentary information to DOS. A joint phantom for multimodal optical measurements has to be liquid to allow the proper Brownian motion of the scattering particles. With respect to a standard liquid phantom, in this case, there is also the need to tune the Brownian diffusion coefficient D b independent of the optical properties. For this reason, the possibility was tested to tune the phantom viscosity by adding glycerol, which has high viscosity and negligible absorption [184]. Unfortunately, its refractive index is very close to that of fat droplets, reducing their scattering power [185]. A TD-DOS system was employed to characterize the µ s variation with the relative content of glycerol in water between 600 and 1100 nm (Figure 25b), making it possible to reliably compensate for the decreased scattering power with a calibrated increase of fat droplets concentration, and thus enabling the independent tuning of optical and dynamic properties (Figure 25a).  To evaluate the performance of a system in retrieving the relative concentrations of different tissue constituents, it is possible to use water-based phantoms with the addition of other chromophores. Water and lipids are major tissue constituents, thus water-based phantoms are the best starting point for mimicking the absorption spectrum of tissues. In this case, lipid dilutions must be replaced with a more realistic fat, like the purified pork fat. Being water and lipid immiscible, different recipes have been developed, based on the simple emulsion of the two using a professional homogenizer (nature emulsion), or with the addition of a thickener, like Agar or Triton-X [186]. The optical properties of these phantoms were characterized on a broad spectral range (600-1100 nm), as shown in Figure 26 for phantoms containing Triton-X. This allowed one to identify the range of tissue optical properties that can be covered in the 3 cases. While nature emulsions are more limited due to the need to mix them in a similar ratio, Agar and Triton-X allowed one to simulate a much broader range of tissue compositions. To evaluate the performance of a system in retrieving the relative concentrations of different tissue constituents, it is possible to use water-based phantoms with the addition of other chromophores. Water and lipids are major tissue constituents, thus water-based phantoms are the best starting point for mimicking the absorption spectrum of tissues. In this case, lipid dilutions must be replaced with a more realistic fat, like the purified pork fat. Being water and lipid immiscible, different recipes have been developed, based on the simple emulsion of the two using a professional homogenizer (nature emulsion), or with the addition of a thickener, like Agar or Triton-X [186]. The optical properties of these phantoms were characterized on a broad spectral range (600-1100 nm), as shown in Figure 26 for phantoms containing Triton-X. This allowed one to identify the range of tissue optical properties that can be covered in the 3 cases. While nature emulsions are more limited due to the need to mix them in a similar ratio, Agar and Triton-X allowed one to simulate a much broader range of tissue compositions. In this paper, we have discussed broadband TD-DOS exclusively. An alternative method to measure the same tissue properties is based on frequency domain (FD) recordings of the diffuse light propagating through turbid media such as e.g., human tissue. Frequency domain measurements can be limited to one modulation frequency (as compared to the equivalent of all possible modulation frequencies in TD-DOS), and can thus be conducted with few photons for the same signal-to-noise ratio. This benefit always comes to the expense of less information provided. For some applications, a smaller amount of information might suffice, and then FD measurements may be a good way to go. Frequency domain measurements may be of particular interest if a dynamic event is to be recorded in a homogeneous medium and thus the recording time is limited by the rate of changes in tissue. Another example could be the immediate response of tissue to rapid changes in tissue blood perfusion. However, in most other cases, TD-DOS may be the preferred technique for its much higher informative content. In particular, the depth information is retrieved in a straight forward way by time gating or tagging late photons in time of flight curves. With improving data analysis methods for heterogeneous structures, and with further advances in the development of experimental set-up components, we foresee TD-DOS to address pressing needs like depth sensitivity across various domains of highly scattering media (e.g., optical tomography of deep tissue).

Future Outlook
In the last decade, the advancement in source and detector technologies brought in the benefits of considerable cost reduction and miniaturization of broadband time domain instruments as portable devices [41]. Recent studies on small footprint sources and detectors could further promote the miniaturization of TD-DOS towards wearable devices [187,188]. Some initiatives were taken in this direction by EU funded projects. Another key challenge was found in the analysis of heterogeneous media, where recent advances in multi-layer tissue analysis provide hope to further increase the accuracy of analyzed data, thus enabling reliable quantification of tissue constituents in heterogeneous tissue structures [48,152,189]. There is increasing interest in TD-DOS to explore new human tissue types and locations [52,53,190], a further continuation in this direction could add to the exploration and invention of new methodologies to address pressing clinical needs using broadband TD-DOS. In this paper, we have discussed broadband TD-DOS exclusively. An alternative method to measure the same tissue properties is based on frequency domain (FD) recordings of the diffuse light propagating through turbid media such as e.g., human tissue. Frequency domain measurements can be limited to one modulation frequency (as compared to the equivalent of all possible modulation frequencies in TD-DOS), and can thus be conducted with few photons for the same signal-to-noise ratio. This benefit always comes to the expense of less information provided. For some applications, a smaller amount of information might suffice, and then FD measurements may be a good way to go. Frequency domain measurements may be of particular interest if a dynamic event is to be recorded in a homogeneous medium and thus the recording time is limited by the rate of changes in tissue. Another example could be the immediate response of tissue to rapid changes in tissue blood perfusion. However, in most other cases, TD-DOS may be the preferred technique for its much higher informative content. In particular, the depth information is retrieved in a straight forward way by time gating or tagging late photons in time of flight curves. With improving data analysis methods for heterogeneous structures, and with further advances in the development of experimental set-up components, we foresee TD-DOS to address pressing needs like depth sensitivity across various domains of highly scattering media (e.g., optical tomography of deep tissue).

Future Outlook
In the last decade, the advancement in source and detector technologies brought in the benefits of considerable cost reduction and miniaturization of broadband time domain instruments as portable devices [41]. Recent studies on small footprint sources and detectors could further promote the miniaturization of TD-DOS towards wearable devices [187,188]. Some initiatives were taken in this direction by EU funded projects. Another key challenge was found in the analysis of heterogeneous media, where recent advances in multi-layer tissue analysis provide hope to further increase the accuracy of analyzed data, thus enabling reliable quantification of tissue constituents in heterogeneous tissue structures [48,152,189]. There is increasing interest in TD-DOS to explore new human tissue types and locations [52,53,190], a further continuation in this direction could add to the exploration and invention of new methodologies to address pressing clinical needs using broadband TD-DOS. Figure 27 summarizes the current status and future trends in broadband TD-DOS. One trend of growth is foreseen in the extension of broadband TD-DOS to the SWIR region (1000-1900 nm) [121,132]. Though few attempts were made till date, a dedicated approach in this direction would open up many more opportunities as the depth sensitivity could be higher due to relatively low scattering in that wavelength range. The measurement of tissue chromophore performed by broadband TD-DOS could complement other technologies like an ultrasound for morphological information, elastography to estimate tissue stiffness, diffuse correlation spectroscopy to monitor blood flow. Some exploration is already under progress in this direction [36,187,188,191]. A further push is foreseen in the near future to bring in a multimodal broadband TD-DOS, where the specificity of TD-DOS will be combined with complementary information from other modalities. A preliminary attempt to standardize diffuse optics has been performed in the past [49][50][51]. Recently, renewed attempts to further strengthen it could bring in the much-needed standards that would harmonize TD-DOS measurements across laboratories [192].
Also, an open data of broadband, well-characterized, stable, reproducibly recipe to make phantom kits, and anatomically correct anthropomorphic tissue phantoms for clinical exploration could further support the cause of standardization [193].  Figure 27 summarizes the current status and future trends in broadband TD-DOS. One trend of growth is foreseen in the extension of broadband TD-DOS to the SWIR region (1000-1900 nm) [121,132]. Though few attempts were made till date, a dedicated approach in this direction would open up many more opportunities as the depth sensitivity could be higher due to relatively low scattering in that wavelength range. The measurement of tissue chromophore performed by broadband TD-DOS could complement other technologies like an ultrasound for morphological information, elastography to estimate tissue stiffness, diffuse correlation spectroscopy to monitor blood flow. Some exploration is already under progress in this direction [36,187,188,191]. A further push is foreseen in the near future to bring in a multimodal broadband TD-DOS, where the specificity of TD-DOS will be combined with complementary information from other modalities. A preliminary attempt to standardize diffuse optics has been performed in the past [49][50][51]. Recently, renewed attempts to further strengthen it could bring in the much-needed standards that would harmonize TD-DOS measurements across laboratories [192]. Also, an open data of broadband, wellcharacterized, stable, reproducibly recipe to make phantom kits, and anatomically correct anthropomorphic tissue phantoms for clinical exploration could further support the cause of standardization [193].

Summary and Conclusions
We have presented an up-to-date review in broadband TD-DOS covering various topics of physics, instrumentations, methods of analysis, application to tissue chromophore characterization, in vivo studies, food, wood, pharmaceuticals, and standardization. A renewed focus across the world in biophotonics along with recent advancements in technology could bring further growth in the field. To catalyze this growth, an open data approach and a standardization attempt need to be

Summary and Conclusions
We have presented an up-to-date review in broadband TD-DOS covering various topics of physics, instrumentations, methods of analysis, application to tissue chromophore characterization, in vivo studies, food, wood, pharmaceuticals, and standardization. A renewed focus across the world in biophotonics along with recent advancements in technology could bring further growth in the field. To catalyze this growth, an open data approach and a standardization attempt need to be initiated across laboratories to harvest maximum benefits by harmonizing the measurement protocols. Taking maximum benefits of technological advancements, the vision for the next decade is to push the limits of broadband TD-DOS to bring in small footprint, standardized, multimodal, application-specific TD-DOS devices [41].