Performance assessment of high-density diffuse optical topography regarding source-detector array topology.

Recent advances in optical neuroimaging systems as a functional interface enhance our understanding of neuronal activity in the brain. High density diffuse optical topography (HD-DOT) uses multi-distance overlapped channels to improve the spatial resolution of images comparable to functional magnetic resonance imaging (fMRI). The topology of the source and detector (SD) array directly impacts the quality of the hemodynamic reconstruction in HD-DOT imaging modality. In this work, the effect of different SD configurations on the quality of cerebral hemodynamic recovery is investigated by presenting a simulation setup based on the analytical approach. Given that the SD arrangement determines the elements of the Jacobian matrix, we conclude that the more individual components in this matrix, the better the retrieval quality. The results demonstrate that the multi-distance multi-directional (MDMD) arrangement produces more unique elements in the Jacobian array. Consequently, the inverse problem can accurately retrieve the brain activity of diffuse optical topography data.

Introduction

Recent developments in functional neuroimaging systems quantitively enhance our understanding of spatially and temporally distributed neural activity in the brain [1,2]. Functional near-infrared spectroscopy (fNIRS) is a new, emerging, and growing technology for monitoring neurological activity in which red and near-infrared light is used to measure changes in Oxy- and Deoxyhemoglobin in brain tissue [3-6]. The fNIRS is an optical neuroimaging technology that is radiation-free, relatively inexpensive, compatible with implanted electronic devices, and portable and can wirelessly record brain activity [7-10]. Although the fNIRS systems are mobile and compact, the resolution and depth of the imaging are less than those obtained by the fMRI [11]. HD-DOT allows brain activity to be mapped in 3D by creating overlaps between the fNIRS channels. This method uses high-density SD arrays to improve the spatial resolution that is comparable to fMRI [10].

SD configuration and number of overlapped channels alongside the inverse problem significantly influences the spatial resolution [10,12]. In this paper, the impact of SD topology on the spatial resolution and hemodynamic reconstruction of HD-DOT has been investigated. We have developed an analytical simulation setup to evaluate the performance of the different combinations of SD on hemodynamic regeneration. Also, this simulation setup can be employed to optimize the arrangement of SD and the number of multi-distance channels and performance of the inverse problem on hemodynamic reconstruction. The SD arrangement and the number of channels determine the elements of the Jacobian matrix. We observed that the MDMD arrangement produces more unique components in the Jacobian matrix. The outcomes of this work indicate that the more individual ingredients in this matrix, the better the reconstruction quality. The Jacobian is the sensitivity matrix, which is computed by the forward model. In this work, the solution of Diffusion Equation inside inhomogeneous media similar to properties of the brain tissue constitutes the base of the forward model.

The forward model is a part of the simulation setup, which plays a fundamental role in confirming the results of this research. Therefore, in the discussion section, the performance of the forward model [13] used in this study is compared with the statistical model on Colin 27 brain template [14].

The rest of the study has continued as follows. The second section describes the simulation setup, which is applied to four different topologies of SD on the forward model. In the following, the synthetic fNIRS data are modeled and generated to simulate hemodynamic changes in the brain. Finally, based on the calculated diffuse reflectance, the inverse method was employed for hemodynamic reconstruction. Part 3 represents the simulation of the forward model and hemodynamic reconstruction. The discussion and conclusion of this study are outlined in sections 4 and 5.

Materials and methods

Modeling approach

The block diagram of the proposed simulation setup illustrated in Fig 1 contains all the steps taken in this article. This scheme employs an analytical forward model that has less computation time. As the number of channels increases, the computational volume in the modeling increases, so the use of analytical models are preferred to numerical approaches [15,16]. Since the perturbative Diffusion Equation (pDE) equations do not have an analytical solution in complex geometries, the simple geometry, that can approximate a semi-infinite medium for thick slabs, is utilized [16].

Fig 1

Illustrates the modeling approach in this work and the different arrangements of source and detector in XY-plane.

(a), (b), (c), and (d) represent the configuration of the SD in four different modes (HD-DOT, Arrangement-1, Arrangement-2, and MDMD topology). (e) Block diagram of analytical simulation setup.

The sources and detectors are aligned in XY-plane on the top of one layer slab geometry. Synthetic cerebral hemodynamic is generated to model more realistic reflectance. The cost function, which is used to modify the regularization parameter, is the correlation between reconstructed and synthetic data. We will prove that the accuracy of reconstruction depends on the SD arrangement, the number of channels, and the energy regularization parameter. Here the forward model plays an essential role since the validation of a model to realistic results mainly depends on it.

Also, the simulation setup is expandable for many quantities of perturbation inside the medium and can be used to evaluate the performance of HD-DOT. The purpose of this scheme is to recover the synthetic hemodynamic in location S1-S9, particularly the recovery of S5. The different arrangement of SD is applied to the simulation setup, and the potential of each combination is analyzed in hemodynamic reconstruction. The first arrangement (Fig 1(a)) represents one repeatable part of HD-DOT [10]; the other remaining provisions are proposed to view the effect of different SD topology on hemodynamic reconstruction (Fig 1(b), 1(c) and 1(d)). The block diagram of Fig 1(e) illustrates the analytical simulation setup. The following sections will describe this diagram in detail.

Arrangement of SD

There is a direct relationship between SD topology and the accuracy of hemodynamic signal reconstruction. The higher the number of channels with different overlapped directions and distances, the higher the efficiency of the recovery. For the same amount of SD, MDMD arrangement (shown in Fig 1(d)), has higher multi-distance, multi-directional channels compared to HD-DOT. Arrangement-1 and arrangement-2 both have the same number of SD, but the lowest quantity of channels belongs to arrangement-1. These four arrangements are illustrated in Fig 2 and compared in terms of multi-distance, multi-directional, and the number of channels.

Fig 2

(a), (b), (c), and (d) indicate the SD location, number of channels, number of distinct directions, and distance in XY-plane for HD-DOT, arrangement-1, arrangement-2, and MDMD respectively.

The simulation setup was created based on the analytical solution of perturbation theory to verify the accuracy of the SD arrangement in the reconstruction of the hemodynamic response, which will be described in the next section.

Theory for reflectance perturbation

The analytical solution for perturbative DE has been obtained for the geometry of Fig 3. The geometry of the boundary for the analytical solution of the perturbative DE is a slab. The slab geometry is widely used for calculation of the reflectance in brain functional imaging [16-20]. It is better to note that the boundary is not limited in the direction of axis X and Y [16].

Fig 3

Location of SD in Cartesian coordinates on the surface of slab geometry and nine perturbed inclusion inside it.

For the sake of simplicity, the geometry has been supposed to have only one layer. The absorption and scattering coefficient of the slab respectively considered to be μ = 0.01mm−1 and and the thickness of slab is equal to40mm and refractive index n = 1.4 [16].

The result obtained by this approach is accurate when the defect causes small perturbation on photon migration. Consequently, the volume of the inhomogeneity (inclusion) is regarded to be small concerning baseline optical properties of the homogeneous medium [16].

The reflectance of each channel has been calculated for several inclusions inside the medium. Reflectance in inhomogeneous media is the superposition of the reflectance inside homogeneous media R0(ρ, t), plus the absorption (δR(ρ, t)) and scattering (δR (ρ, t)) effect of inclusion [16].

Inside the slab medium, nine dynamic perturbations are inserted to simulate the real function of the brain in the cerebral cortex. Each perturbation is located at the center of a voxel at 15mm depth.

The final expression of R(ρ) For each channel derived based on Born approximation [21]. Finally, the diffuse reflectance for the channel between a source located in Sx and a detector at Dx would be: Where R0(ρ) is the reflectance for homogeneous media and the δR(ρ, T, i) is the absorption perturbation of ith inclusion and T is the sampling time of dynamic perturbation.

The channel definition in Eq (2) is based on the source (Sx) and detector (Dx), and the corresponding distance between them (ρ).

The position of SD in Cartesian coordinates on the surface of Slab geometry (XY-plane) and nine perturbed inclusion inside it has displayed in Fig 3. For example, the dynamic perturbation located at the center of S9 voxel has shown in this figure. The distance between adjacent voxels is considered to be 10mm. Next section will describe how this dynamic perturbation is generated.

Synthetic fNIRS data

This section describes the scheme for simulating the perturbation inside the medium. So R(ρ, T) is modulated using the synthetic Hemodynamic Response Function (HRF). The ΔHbO2 and ΔHb concentration in this medium were perturbed in nine regions, to mimic hemodynamic response concerning the duration of the task. The event or task duration was considered to be random to examine the performance of the inverse problem in all possible states. The perturbation is generated by the convolution of the boxcar function with synthetic hemodynamic. Boxcar function (s(t)) is regularly repeated as a rectangular pulse waveform with modulated duty cycle (related to the duration of task). The amplitude of 1 indicates the task, and 0 refers to rest [22]. S(t) is the pulse-width modulated signal: Where, i = 1: N; j = 1,2, …, (N − t), and d represents the duration of each task. The HRF(t) was modeled as a linear combination of two different gamma variant time-dependent function [23]:

With: Where α determines the amplitude, ρ1, and ρ2 regulate the starting, end, and duration of HRF, τ1, and τ2 tune the ascending and descending shape of HRF, while β control the undershoot. The value of p coefficient recommended being five [24]. The HRF profile with a peak amplitude of almost 1555 nM was chosen for HbO2 while the Hb profile is the same as HRF for HbO2 but with an inverted magnitude by 33% attenuation and regulated latency [23]. The change of HbO2 corresponding to each perturbation would be the convolution of the HRF(t) and s(t) plus physiological noise:

The physiological noise was modeled as a linear combination of sinusoids [25]:

The ∅(t) for each perturbation is the average of the ten trials of the Eq (7). The value of amplitude and frequency of the sinusoids would be different for each repetition, while phase θ are equally distributed between 0 and 2π for each trial [23].

Inverse problem: Hemodynamic reconstruction

The detectors record the variation in the light intensity, which is formed by the corresponding source. These detectors represent the optical properties of the channel that is called optical density. The Modified Beer-Lambert law (MBLL) is used to relate changes in optical density to changes in concentration of Oxy and Deoxyhemoglobin under the assumption that the scattering losses are constant (S1 Appendix indicates further information about these equations). The variation in Oxy-Hemoglobin (ΔHbO2(t)) and Deoxy-Hemoglobin (ΔHb(t)) according to Beer law [26] determines the change in absorption coefficient (Δμ(λ)). According to Beer’s law [26]: Where, n represent the number of the light absorbing agent (chromophores) in the tissue. However, in near-infrared wavelength (700-900nm), the dominant absorption changes are caused by the concentration of O2Hb and HHb. As a result, the Δμ can be expressed by: Where, and represent the extinction coefficients of Oxy and Deoxyhemoglobin, respectively. The forward model is formed by using the synthetic fNIRS and perturbation theory in two wavelengths.

The solution of the inverse problem to the forward model is required to estimate the synthetic hemodynamic. The forward model can be rewritten as: Where, is the absorption perturbation for each location of the domain under the head surface. is modulated by synthetic hemodynamic and is the Jacobian matrix (it shows the sensitivity of reflectance to each perturbation in specific depth): Where index “j” refers to the number of channels and index “i” refers to the number of perturbations under the SD array. There are several approaches to solve the inverse problem of Eq (10) [27]. One of the commonly employed methods to provide a solution to the inverse problem is energy regularization [28]. Reconstruction of the hemodynamic is obtained by: Where, is the transposition of the Jacobian matrix, “ϵ” is the energy regularization parameter and “I” is the identity matrix. The optimum value for “ϵ” is found empirically based on the simulation results.

Given the Eq (9), and are a function of the wavelength. By calculating the variation of optical density ΔOD at two wavelengths (in this paper 780 and 820 nm), and assuming that the length of the traveling light is identical in both wavelengths, the values of Δμ(λ1) and Δμ(λ2) are obtained. As a result, the corresponding hemodynamic can be reconstructed: Where in this work λ1 = 780nm and λ2 = 820nm. The extinction coefficients of Oxy and Deoxy hemoglobin for both wavelengths are given by [29].

Results

Forward simulation

The elements of the Jacobian matrix for the given arrangement have been calculated by sweeping one inhomogeneity in a 3D position in the medium under study. The contrasts (δR(ρ)⁄R0(ρ)) of channels (Fig 4) have been obtained for 10mm, 20mm, 30mm, 40mm, 50mm, and 60mm distance between SD. The contrast indicates the sensitivity of reflectance for any perturbation inside 3D geometry. This simulation result shows less sensitivity in profound depth. This feature has been described comprehensively and in detail by [27]. According to this figure, penetration depth increases throughout the growing distance among SD.

Fig 4

The contrast (δRa(ρj)⁄R0(ρj)) in XZ-plane for SD separation from 10mm to 60mm.

The synthetic oxyhemoglobin (Δ(HbO2)) and deoxyhemoglobin (Δ(Hb)) are generated with the corresponding Δμ(λ) at the wavelength of the 780nm and 820nm (Fig 8). The solution of DE based on perturbation theory and Born approximation for nine perturbations is employed. The transient brain activity is modeled in a way to modulate absorption. So Δμ(λ) can be time-dependent (Δμ(λ, t)). Then the reflectance due to perturbations inside the medium depends on Δμ(λ, t) of each perturbation. Consequently, R(ρ, T) is modulated using synthetic hemodynamic.

Fig 8

Indicates the similarity of all reconstructed and synthetic hemodynamic in S1-S9 region for HD-DOT, Arrangement-1, Arrangement-2 and MDMD.

According to the forward model, nine synthetic hemodynamic response is generated at 15mm depth and reconstructed by the solution of the inverse problem. Fig 5 shows the change in synthetic oxy and deoxyhemoglobin. This figure represents that each hemodynamic has a distinct pattern compared to others. The different hemodynamic trends in each voxel make it more challenging to recover hemodynamics, and under these conditions, the capability of the inverse algorithm and the SD arrangement can be explored. It can be noted that the nearby hemodynamic activity acts as a systematic noise. Therefore, the hemodynamic recovery of the S5 region will be more difficult because it is surrounded by eight hemodynamic noises. Thus, for hemodynamic retrieval of this area, the number and angle of observations must be much higher than the number of hemodynamic sources under the inspection field.

Fig 5

Changes in synthetic oxyhemoglobin (Δ(HbO2)) and deoxyhemoglobin (Δ(Hb)).

Inverse procedure: Hemodynamic reconstruction

The correlation rate was used to investigate the similarity between the two synthetic and obtained hemodynamic. The accuracy of reconstruction not only depends on the arrangement and number of the SD, but it also relies on the solution of the energy regularization. The optimum value for the energy regularization parameter (ϵ) was achieved by sweeping this parameter from 10−8 to 10−3 and minimizing the cost function; correlation coefficient (CC) for four different SD topology. By using Eqs (13) and (14), and reconstructed in two wavelengths, the calculated Oxy-Deoxy hemoglobin along with synthetic data are compared in the S5 region in Fig 6 (in the existence of SD arrangement of HD-DOT, Arrangement-1, Arrangement-2, and MDMD).

Fig 6

Visual comparison between reconstructed and synthetic hemodynamics for S5.

Fig 6 shows that the topologies of Arrangement-2, HD-DOT, and MDMD have been able to extract the hemodynamics of the S5 region with greater accuracy, but the Arrangement-1 has inferior performance compared to other configurations.

The trend of CC(ϵ) within Fig 7 represents that change in SD topology leads to accurate hemodynamic reconstruction. It also indicates that the MDMD has superior performance compared to other arrangements. The production of each SD arrangement also analyzed in the rebuilding of all dynamic perturbations of S1-S9 using the value of CC. The result of this comparison, as shown in Fig 8, reveals that almost all the combinations have acceptable performance except for the Arrangement-1, which has poor performance in hemodynamic extraction of the S5 region. This figure also confirms that the MDMD operates properly in extracting all the hemodynamic sources because it has many unique elements in the Jacobian matrix compared to other topologies (Fig 9).

Fig 7

Depicts the performance assessment of different SD arrangements in hemodynamic extraction of the S5 region concerning the regularization parameter.

Fig 9

Represents the number of unique elements in the Jacobian matrix regarding the SD topology.

The performance of all the topologies studied in this work is summarized in Table 1 in terms of the topology of the arrangement, the number of SDs, channels, distances, directions, as well as the unique elements of the Jacobian matrix. This table notes that as much as the arrangement of the SDs creates various distance and direction between the overlapped channels, the better the hemodynamic extraction will be observed.

Table 1

Summarized the details and the performance of all SD topologies studied in this investigation.

SD Topology	Number of Distance	Number of Direction	Number of SD	Number of Channels	Number of Unique Elements in Jacobian Matrix	Total Correlation Coefficient
Arrangement-1	2	6	9	12	18	0.9586
Arrangement-2	4	3	9	16	33	0.9764
HD-DOT [10]	4	9	13	36	35	0.9949
MDMD	5	13	13	36	50	0.9986

A singular value analysis of the Jacobian matrix associated with introduced SD arrangements is used as a benchmark [30]. Besides unique elements of the Jacobian matrix, the singular value analysis of the different methods in Fig 10 indicates that both shape of the singular value spectra and the magnitude for MDMD arrangement is higher than other SD combinations.

Fig 10

Singular value spectra for Jacobian matrix of HD-DOT, MDMD, Arrangement-1 and 2.

It is worth noting that, if the depth information is not necessary, and the objective is to achieve topography, Arrangement-2 can be replaced instead of MDMD and HD-DOT because it has fewer SDs (reduces the complexity of the device) and has acceptable performance compared to these arrangements.

Discussion

In this paper, based on a developed simulation setup, the performance of SD arrangement and their quantity alongside inverse problem on hemodynamic reconstruction has been investigated. The simulation approach consists of a forward model, synthetic fNIRS data generation, Inverse problem, and SD arrangement.

The forward model is an analytical method that is implemented by the solution of the pDE in slab medium. Analytical methods have been developed earlier to study light emission inside the simple geometry such as slab medium [13,31-36]. Numerical methods also are used in complex brain models to study light diffusion in tissues [14,37]. In spite of simplicity and approximation, analytical methods take less time calculation compared to statistical approaches, especially when it comes to investigating the effect of several fNIRS channels on depth sensitivity.

The spatial probability pattern of photons penetrating tissue at the source position, scattering within the tissue, and being exposed at a particular detector spot, determines the spatial sensitivity profile for the SD pairs [14]. To verify the forward model, the spatial sensitivity profile is compared with the results of the Monte Carlo method on Colin 27 brain template. The depth sensitivity for analytical pDE inside slab geometry is defined as follow: Where x represents the distance between source and detector.

The depth sensitivity of analytical pDE is compared to the regression of Monte Carlo (MC) on Colin 27 brain template in Fig 11, and the mismatch between these analytical and numerical methods are illustrated in this figure for given SD separations. Toward SD separation around 30mm-50mm, the mismatch is less than 60% for penetration depth from 1mm to 28mm. In 15mm depth, the mismatch is less than 20% for SD separation of 30mm-50mm.

Fig 11

Represents the depth sensitivity of the analytical pDE in slab medium and depth sensitivity of the regression of MC on Colin 27 geometry for SD separation of 10mm-60mm.

The outcomes of the comparison indicate that the analytical approach is not too close to the results of Monte Carlo methods. It is better to emphasize that there is no analytical solution to light- tissue interaction inside complex geometries like the brain. Besides that, due to the high computation time of statistical models are not a desirable candidate to be a forward model in this simulation setup. Although Monte Carlo methods are accurate in estimation of depth sensitivity, they are not well suited to be used in the forward model. Whenever multiple dynamic perturbations exist in geometry, then numerical methods should be run for each sampling time.

Consequently, the computation time will grow significantly by applying statistical approaches. It is suggested to try the solution of analytical techniques on multilayered geometries as a forward model inside the simulation setup. Another alternative to the computation of the forward model and Jacobian matrices is to use a finite element method (FEM) [27,36]. Recently, both NIRSFAST and Neuro-DOT software have been developed to the solution of the forward model based on FEM estimation [38-41].

The results of FEM data when it is applied to Diffusion Equation will be close to reality compared to the analytical techniques, but the computation time will grow considerably [15]. On the other hand, adding dynamics to the perturbations inside the meshed environment will increase the simulation time significantly. Analytical simulation of 36 channels, including a sampling rate of 30 samples per second, takes less than 10 minutes. For the same sampling rate and channels, the FEM simulation will take almost 43–48 hours, depending on system speed. The geometry and mesh of the medium for simulating the FEM are represented in Fig 12(a) and 12(b), respectively.

Fig 12

(a) Illustrates the geometry of the Slab medium which is 120mm × 120mm × 40mm. (b) represents the mesh of the medium.

According to the simulation approach presented in Fig 1, White Gaussian Noise (WGN) is added to each channel () before the reconstruction of the simulated data to avoid inverse-crime. WGN indicates the instrumentational noise, which depends on the Signal to Noise Ratio (SNR) of the signal acquisition device. The effect of 47 dB SNR on the performance of different SD arrangements in the hemodynamic recovery of region S5 has been investigated. Considering the simulation of Fig 13 for Arrangement-1 and 2, there is a significant change in magnitude and shape of CC versus regularization parameter. While little difference with the noiseless condition for MDMD and HD-DOT is observed, it can be concluded that 47 dB of SNR has no significant effect on the performance of these two methods. The SNR has been swept from 32 dB (worst case condition) up to 52 dB (for given ϵ = 10−5) to compare the performance of SD arrangements on hemodynamic reconstruction. The simulation of Fig 14(a) represents that for any SNR, MDMD performs better than other competitors. Concerning Fig 14(b) in the worst-case condition, the optimal point for the regularization parameter has been shifted. Even with the highest noise through determining the appropriate regularization parameter, MDMD still works better than HD-DOT. For SNR = 32 dB, the reconstructed and synthetic hemodynamic of region S5 are presented in Fig 14(c) and 14(d) for both MDMD and HD-DOT, respectively.

Fig 13

Demonstrates the effect of noise on performance assessment of different SD arrangements in hemodynamic extraction of the S5 region concerning the regularization parameter.

Fig 14

(a) Represents the ability of hemodynamic extraction from the S5-region concerning several SNRs. (b) illustrates that in the worst case condition (SNR = 32 dB), the optimum point for ϵ has changed from 10−5 to 10−4. (c) and (d) Represents reconstructed and synthetic ΔHbO2 and ΔHb for MDMD and HD-DOT respectively.

The forward model has been generated with well-defined optical properties using the Jacobian matrix, and the inversion has been performed using the same matrix. Previously, to avoid inverse-crime, WGN was added to each channel. One can also bring the forward model closer to the more realistic model by adding uncertainty to all elements of the Jacobian matrix. For this purpose, Jacobian matrix elements are multiplied by the Gaussian random coefficient. Therefore the forward model is changed as follows:

The Rnd matrix carries random coefficients with Gaussian distribution, the performance of the two MDMD and HD-DOT methods are close together, the effect of the changes on the forward model will only be investigated on the performance of these two methods.

Random coefficients with two different distributions are applied to the forward model. The mean of both data set is one, and the variance (δ) of the first and second random coefficients are 0.07 and 0.2, respectively (Fig 15(a) and 15(c)). Beside ten simulation runs for 0.07 variance, the distribution of random coefficients is plotted in Fig 15(a). The performance of both the MDMD and HD-DOT methods in hemodynamic extraction of the S5 region is compared in Fig 15(b). Similarly, for the variance of 0.2, the above comparison is repeated 18 times. In this case, the distribution of random data for this variance is shown in Fig 15(c) and the results of analyzing of similarity are illustrated in Fig 15(d).

Fig 15

(a) Represents Gaussian distribution (mean = 1 and δ = 0.07) of Rnd matrix elements. (b) Boxplot regarding random coefficients. (c) Represents Gaussian distribution (mean = 1 and δ = 0.2) of Rnd matrix elements. (d) Boxplot regarding random coefficients.

If the variance of random data in the Rnd matrix is 0.07, according to the results of Fig 15(b), the performance of MDMD is better than HD-DOT. Even with a variance of 0.2, MDMD still delivers better results. It should be noted that there is no similarity between the data extracted from the S5 region and the synthetic data under the variance of 0.2 in some states. If the correlation is less than 0.8, there will be no similarity between the reconstructed and synthetic signals.

Although adding uncertainties to the Jacobian matrix elements in the forward model can reduce the gap between the results of this work and the real imaging applications. However, given the limitations of the analytical model used in this comparison, it is not yet possible to claim that MDMD will have better results than other SD configurations in imaging applications.

This proposed simulation setup is expandable for many numbers of perturbation inside the medium and can be used for performance assessment of HD-DOT. A. Eggebrecht and colleagues in 2014 have used HD-DOT for mapping brain function [10], the proposed model can be used for evaluation of the SD separation and arrangement of SD on performance of HD-DOT.

The outcome of this model can be an instrumentation probe with a specific arrangement of SD array for monitoring stimulation induced hemodynamic. Among noninvasive stimulation approaches such as transcranial direct current stimulation and transcranial magnetic stimulation, low-intensity ultrasound stimulation has the spatial resolution in the order of several millimeters [42]. To control the amount of stimulation and study the effect of stimulation on the brain, a simultaneous recording of the hemodynamic activity of the brain is necessary [43]. Therefore, a non-invasive method is required for recording stimulation-induced hemodynamic with a spatial and temporal resolution equivalent to the stimulation approach [44].

Summarily, this simulation setup can be employed for performance assessment of high or low-density DOT, monitoring of stimulation-induced hemodynamic, and SD array design. Finally, the proposed simulation approach, with declared assumptions and simplifications, can be used by researchers who want to arrange an array of sources and detectors for optical topography.

Conclusion

In this paper, an innovative simulation setup proposed for the performance assessment of a variety of sources and detectors toward the rebuilding of cerebral hemodynamics. MDMD arrangement involves more unique elements in the Jacobian matrix and will be able to reconstruct the neural activity accurately. Meanwhile, the performance of several provisions of SD on the reconstruction of brain function is studied. The result of simulation indicates that raising the number of multi-distance and multi-directional overlapped channels increase the accuracy of brain hemodynamic reconstruction. Also, the simulation setup can be employed for performance assessment of high or low-density DOT, monitoring of stimulation-induced hemodynamic, and SD array design. We believe that modeling and simulation of different SD arrays on hemodynamic extraction optimizes the number of SDs required for accurate spatial imaging. Consequently reduces the additional cost and complexity of device fabrication. Based on the modeling approach and simulation results, the achievements of the MDMD looks more beneficial than other methods. But still, there would be a gap between the outcomes of this study and real imaging applications. It is suggested to use the proposed simulation approach with a modified forward model according to the following suggestions. The brain model in this work is one-layer, while the multi-layered medium can be considered to get closer to the real results. Also, we must find the interaction between the perturbations in the solution of DE equations (whenever the perturbations are not small regarding the baseline optical properties), while Born approximation is used in this work. In this study the diffuse reflectance modulated due to change in absorption coefficient. The changes in the optical scattering coefficient, along with the absorption coefficient, should be taken into account for more accurate simulation of diffuse reflectance.

(DOCX)

Click here for additional data file.

(RAR)

Click here for additional data file.

12 Nov 2019

Submitted filename: Rebutal Letter PAMDMD.DOCX

Click here for additional data file.

20 Dec 2019

PONE-D-19-31493

Performance Assessment of High-Density Diffuse Optical Tomography Regarding Source-Detector Array Topology

PLOS ONE

Dear Mr. Setarehdan,

Thank you for submitting your manuscript to PLOS ONE. After careful consideration, we feel that it has merit but does not fully meet PLOS ONE’s publication criteria as it currently stands. Therefore, we invite you to submit a revised version of the manuscript that addresses the points raised during the review process.

We would appreciate receiving your revised manuscript by Feb 02 2020 11:59PM. When you are ready to submit your revision, log on to https://www.editorialmanager.com/pone/ and select the 'Submissions Needing Revision' folder to locate your manuscript file.

If you would like to make changes to your financial disclosure, please include your updated statement in your cover letter.

To enhance the reproducibility of your results, we recommend that if applicable you deposit your laboratory protocols in protocols.io, where a protocol can be assigned its own identifier (DOI) such that it can be cited independently in the future. For instructions see: http://journals.plos.org/plosone/s/submission-guidelines#loc-laboratory-protocols

Please include the following items when submitting your revised manuscript:

A rebuttal letter that responds to each point raised by the academic editor and reviewer(s). This letter should be uploaded as separate file and labeled 'Response to Reviewers'.

A marked-up copy of your manuscript that highlights changes made to the original version. This file should be uploaded as separate file and labeled 'Revised Manuscript with Track Changes'.

An unmarked version of your revised paper without tracked changes. This file should be uploaded as separate file and labeled 'Manuscript'.

Please note while forming your response, if your article is accepted, you may have the opportunity to make the peer review history publicly available. The record will include editor decision letters (with reviews) and your responses to reviewer comments. If eligible, we will contact you to opt in or out.

We look forward to receiving your revised manuscript.

Kind regards,

Alberto Dalla Mora, Ph.D.

Academic Editor

PLOS ONE

Journal Requirements:

When submitting your revision, we need you to address these additional requirements:

1.Please ensure that your manuscript meets PLOS ONE's style requirements, including those for file naming. The PLOS ONE style templates can be found at http://www.plosone.org/attachments/PLOSOne_formatting_sample_main_body.pdf and http://www.plosone.org/attachments/PLOSOne_formatting_sample_title_authors_affiliations.pdf

2. We note that you have stated that you will provide repository information for your data at acceptance. Should your manuscript be accepted for publication, we will hold it until you provide the relevant accession numbers or DOIs necessary to access your data. If you wish to make changes to your Data Availability statement, please describe these changes in your cover letter and we will update your Data Availability statement to reflect the information you provide.

[Note: HTML markup is below. Please do not edit.]

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Yes

Reviewer #2: Yes

Reviewer #3: Yes

Reviewer #4: No

**********

2. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: N/A

Reviewer #2: No

Reviewer #3: Yes

Reviewer #4: N/A

**********

3. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: Yes

Reviewer #2: No

Reviewer #3: Yes

Reviewer #4: Yes

**********

4. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

Reviewer #2: Yes

Reviewer #3: Yes

Reviewer #4: Yes

**********

5. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: The following paper addresses the effect on optical tomography of the source-detector array configurations. The work is done in the Continuous Wave domain and is based on simulated data. Experimental results have not been presented. The forward solver used is based on an analytical solution of the diffusion equation. It has been considered a perturbative solution for the slab in reflectance configuration where an absorption inclusion is inserted inside. The Born approximation is used to obtain the solution of the perturbed reflectance. The results presented show that the multi-distance multi-directional (MDMD) arrangement of sources and detectors produces more unique elements in the Jacobian matrix and consequently the related inverse problem can better retrieve the brain activity of diffuse optical tomography data compared to other arrangements of sources and detectors. The results are corroborated by several retrievals obtained making use of the mentioned forward model applied on synthetic data generated with the same model on which noise is applied.

These results appear scientifically sound and I suggest the publication on PLOS ONE provided the following remarks are addressed.

The results presented in this paper represent a proof of principle since they are obtained by inverting synthetic data generated by using the same forward solver employed in the inversion procedure and with an added noise accordingly to Eq. (7). Thus, they can be considered a first step of study. At this stage arise a question that deserve at least some comments and explanations. Do the improvements obtained by using MDMD hold in all generality when applied to real cases as real experiments on biological tissues, for instance on brain? Compared to the results presented in this paper in real imaging applications the forward solver used in the inversion procedure can show some deficiencies since the diffusion equation and the Born approximation have limitations and the background medium addressed is homogeneous differently from real media. The gap between the case here addressed and the real applications can have some consequences in the conclusions here formulated? I can understand that the MDMD arrangement can be anyway advantageous for the purpose of imaging applications, however some comments should be spent on this fact. It is true that the authors present a comparison of the depth sensitivity of their analytical model with the Monte Carlo data on Colin 27 brain template. However, the results shown in Fig. 10 and 11 emphasizes difference that, although smaller for larger depth, are always present between analytical model and Monte Carlo regression. The authors simply conclude that Monte Carlo methods although accurate in the estimation of depth sensitivity are not well suited to be used in the forward model due to their long computation time. Is it possible to conclude that the differences observed are not decisive in the inversion procedure?

Minor points:

1) At line 60 the single Ref. 15 is not enough, I would also add Ref. 16.

2) At line 62 it is written “… simple slab boundary condition in semi-infinite geometry are utilized [16]”. The sentence is confusing in the sense that is confused the meaning of geometry and boundary condition. Here the term boundary condition appears misleading. To me would make sense to rephrase as “… the simple geometry, that can approximate a semi-infinite medium for thick slabs, is utilized [16]”

3) At line 96 it is written “The analytical solution for perturbative DE has been obtained in the boundary condition of Fig. 3.” So far, I understand it should be “The analytical solution for perturbative DE has been obtained for the geometry of Fig. 3.” Boundary condition and geometry have a different meaning in this context. The actual boundary condition used to solve DE for the forward solver used in this paper is the extrapolated boundary condition that is not a geometry but a condition to make an energy balance at the external interface of the medium.

4) At line 102 the symbol S (The thickness of the slab) is not defined. NB that S is also the symbol used for the different regions S1-S9 (See Fig. 3).

5) At line 104 is written “Consequently, the volume of the inhomogeneity (inclusion) is regarded to be small concerning baseline optical properties of the homogeneous medium [16].” It should be “Consequently, the volume of the inhomogeneity (inclusion) is regarded to be small concerning baseline optical properties of the homogeneous medium [16].”

6) At line 112 it is written “The final expression of Rpert(\\rho) for each channel derived based on Born approximation [21].” It looks like the sentence miss the final part. The sentence recalls the Born approximation; however, Ref. 21 is mainly related to higher order perturbation theory. Do the authors mean that in Ref. 21 the results for the Born approximation are also presented?

7) At line 136 the acronym HRF is used without definition. I understand that implicitly means Hemodynamic Response Function, however why it should be omitted this definition?

8) At line 272, Eq. (15), index of the sum in Eq. (15) is a real number, while has not been used an index numbering the number of Source-Detector pairs?

9) At line 292, it is mentioned that the results of FEM will be closer to reality than analytical techniques. Maybe could be worth to note that to some extent also FEM data when FEM is applied to the Diffusion Equation, since the intrinsic approximations of this theory affects the FEM data.

10) The actual title of Ref. 16 is: “Light Propagation through Biological Tissue and other Diffusive Media: Theory, Solutions and Software” 2009.

11) In the first row of Fig. 4 is missed the info on the y axis (Z (mm)?).

12) In the second row of Fig. 6 is missed the info on the x axis (Time s?).

Reviewer #2: Author report on a simulation platform for diffuse optical tomography (DOT) adapted for reconstruction of hemodynamic responses. The platform is based on analytical solutions, under the Born approximation, of the diffuse equation for a semi-infinite slab. The work is technically sound and, together with the claimed availability of software and data, will be useful for setting up DOT systems.

At my opinion there are the following points to address and clarify:

Reconstruction: what I don’t understand is whether the reconstruction is forced at the depth of 15 mm or it is performed in the whole volume. In the first case, the point has to be emphasized and better specified in the text and, at my opinion, the method can’t be properly called “tomography”. In the second case, as well, it has to be emphasized in the text and, what I expect, is a figure representing slices in the volume at a defined time Ts to see the reconstructed depth that, in DOT, is typically underestimated.

Noise and Inverse-crime: As far as I’ve understood, the only noise added is on the optical properties of the S5 and nearby voxels. This means that the forward model has been generated with well-defined optical properties using the Jacobian matrix, and the inversion has been performed using the same matrix on unnoisy data ΔR. This is typically called an “inverse-crime”. I suggest to add noise (Gaussian or Poisson) to the simulated data to avoid this problem.

Minors:

the acronym HRF is not specified in the text.

Figure 1 is reported with a very poor resolution, writings are too small.

Reviewer #3: The authors present a simulation setup for the performance assessment of different source-detector configurations in high density diffuse optical tomography. The information content in the paper is well organised. While the work presented is interesting and very useful research, a few concerns are raised below, which are to be addressed prior to any publication.

Major concerns:

1) One of the important aspects of the presented work is the significant reduction in computational time due to the use of analytical forward model, when compared to finite element and Monte Carlo methods. Therefore, it is crucial to provide the readers with the comparison of individual computational times for a sample forward model.

2) From Figure 9, and Table-1, the authors want to show that higher the number of unique elements in Jacobian matrix, better will be the recovery. This can be misleading. While it is logical that a Jacobian with more independent information can give a better recovery, the correct way to observe this is to compare the normalised singular values of different Jacobian matrices corresponding to the SD arrangements. I would recommend looking into: Optics Letters Vol. 26, Issue 10, pp. 701-703 (2001).

3) Figure 10 and 11 have same information, therefore figure 10 can be avoided and two additional sub-plots for SD separation 10 and 20mm can be included in Figure 11.

Minor concerns:

1) Add reference for lines 150-152.

2) Use one abbreviation for oxy hemoglobin (either O2Hb or HbO2), and for deoxy hemoglobin (either HHb or Hb) throughout the manuscript.

3) Important future directions such as the use of multi-layered medium, fit better in the conclusion section rather than discussions. Re-edit the conclusion section to incorporate this information seamlessly.

Reviewer #4: In the manuscript, ‘Performance Assessment of High-Density Diffuse Optical Tomography Regarding Source-Detector Array Topology,’ the authors present a simulation study comparing various source-detector separation distances. While the methods are sound, the choice of methods (Born approximation, limited number of source-detector pairs and therefore a small field of view, analytical model instead of anatomy-based FEM using either diffusion or MonteCarlo, the single layer of optical properties) leads to the results and discussion providing limited information to the community for further advances in simulation, system design, or empirical considerations. As such, the authors are encouraged to add complexity to their modeling procedures or better motivate their choices to place their work in context with the current status of the field.

**********

6. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

Reviewer #2: No

Reviewer #3: No

Reviewer #4: No

[NOTE: If reviewer comments were submitted as an attachment file, they will be attached to this email and accessible via the submission site. Please log into your account, locate the manuscript record, and check for the action link "View Attachments". If this link does not appear, there are no attachment files to be viewed.]

While revising your submission, please upload your figure files to the Preflight Analysis and Conversion Engine (PACE) digital diagnostic tool, https://pacev2.apexcovantage.com/. PACE helps ensure that figures meet PLOS requirements. To use PACE, you must first register as a user. Registration is free. Then, login and navigate to the UPLOAD tab, where you will find detailed instructions on how to use the tool. If you encounter any issues or have any questions when using PACE, please email us at figures@plos.org. Please note that Supporting Information files do not need this step.

Submitted filename: report.pdf

Click here for additional data file.

31 Jan 2020

Dear Professor A. Dalla Mora

We would like to thank you and the esteemed reviewers for providing us with insightful comments, which helped us to revise and improve the quality of the manuscript. We have addressed all the comments, as shown in the revised manuscript. We believe that the contents and the clarity of our paper are much improved in the revised version. Below are point-by-point responses for the four reviewer’s comments. Finally, all the responses in the “Response to Reviewers” have been highlighted in blue.

Editor Comments

Point1: Please ensure that your manuscript meets PLOS ONE's style requirements, including those for file naming. The PLOS ONE style templates can be found at http://www.plosone.org/attachments/PLOSOne_formatting_sample_main_body.pdf and http://www.plosone.org/attachments/PLOSOne_formatting_sample_title_authors_affiliations.pdf

Response: Thanks. The revised manuscript has been edited to meet PLOS ONE’s standards.

Point2: We note that you have stated that you will provide repository information for your data at acceptance. Should your manuscript be accepted for publication, we will hold it until you provide the relevant accession numbers or DOIs necessary to access your data. If you wish to make changes to your Data Availability statement, please describe these changes in your cover letter and we will update your Data Availability statement to reflect the information you provide.

Response: The simulation codes would be available according to PLOS ONE’s policies.

Authors hint to Editor and all Reviewers:

In this response, some simulation results are not included in the revised manuscript. They only are added to make the responses more clear. Those figures from the revised manuscript are highlighted inside the green box. To avoid possible confusion, please consider this hint before reading the answers.

Reviewer #1:

Comments

The following paper addresses the effect on optical tomography of the source-detector array configurations. The work is done in the Continuous Wave domain and is based on simulated data. Experimental results have not been presented. The forward solver used is based on an analytical solution of the diffusion equation. It has been considered a perturbative solution for the slab in reflectance configuration where an absorption inclusion is inserted inside. The Born approximation is used to obtain the solution of the perturbed reflectance. The results presented show that the multi-distance multi-directional (MDMD) arrangement of sources and detectors produces more unique elements in the Jacobian matrix, and consequently, the related inverse problem can better retrieve the brain activity of diffuse optical tomography data compared to other arrangements of sources and detectors. The results are corroborated by several retrievals obtained making use of the mentioned forward model applied to synthetic data generated with the same model on which noise is applied.

These results appear scientifically sound and I suggest the publication on PLOS ONE provided the following remarks are addressed.

Response: We appreciate the positive feedback of the esteemed reviewer.

Major Concerns:

Point 1: The results presented in this paper represent a proof of principle since they are obtained by inverting synthetic data generated by using the same forward solver employed in the inversion procedure and with an added noise accordingly to Eq. (7). Thus, they can be considered a first step of study. At this stage arise a question that deserve at least some comments and explanations. Do the improvements obtained by using MDMD hold in all generality when applied to real cases as real experiments on biological tissues, for instance on brain?

Response: Thanks to the reviewer for raising these points. In this manuscript we have tried to develop a simulation model to understand more about the imaging process and to see the strengths and weaknesses of different blocks such as forward model, SD arrangement and inverse problem.

As you mentioned, this manuscript is the first step of the study. There should be a gap between the results of this investigation with the experiment. The only problem is that in real cases, we do not have access to the inside brain; in another word, access to exact hemodynamic change is not possible. It may need another modality such as ECoG, which is invasive. If we apply different arrangements to tissue-like phantom (as a forward model), there may be a chance to produce synthetic hemodynamic, but it is still challenging. The challenge is how to insert dynamic perturbations inside the phantom. In the Analytical solution, we can add multiple dynamic perturbations. Since we know the pattern of perturbations, we can better analyze the performance of SD arrangement and inversion in hemodynamic reconstruction.

This concern is addressed in lines 358-361 and 389-398 of the revised manuscript.

Point 2: Compared to the results presented in this paper in real imaging applications the forward solver used in the inversion procedure can show some deficiencies since the diffusion equation and the Born approximation have limitations and the background medium addressed is homogeneous differently from real media. The gap between the case here addressed and the real applications can have some consequences in the conclusions here formulated? I can understand that the MDMD arrangement can be anyway advantageous for the purpose of imaging applications, however some comments should be spent on this fact.

Response: Thanks. Based on our knowledge in this field and regarding the previous works, the changes in absorption (or scattering) related to brain activity are small with respect to the base-line values (1% to 5% for oxygenation, with the expected changes in scattering being even smaller) therefore the Born approximation can be used when describing the diffusion of light through brain tissue (Chiarelli et al., 2016). Whenever the intention is to measure hemodynamics from tissues other than the brain, then the Born approximation will not be accurate, and the interaction between the perturbations should be considered in the forward model (Sassaroli, Martelli, & Fantini, 2009).

Please refer to lines 389-398. And the results of Fig. 15.

Point 3: It is true that the authors present a comparison of the depth sensitivity of their analytical model with the Monte Carlo data on Colin 27 brain template. However, the results shown in Fig. 10 and 11 emphasizes difference that, although smaller for larger depth, are always present between analytical model and Monte Carlo regression. The authors simply conclude that Monte Carlo methods although accurate in the estimation of depth sensitivity are not well suited to be used in the forward model due to their long computation time. Is it possible to conclude that the differences observed are not decisive in the inversion procedure?

Response: Thanks for the insightful comment. The forward model in this manuscript is simple and one-layer; we have mentioned this in the discussion section. We could use FEM to solve the diffusion equation in a sophisticated and multi-layered Slab medium. Since the Analytical solution for sophisticated and anatomical mediums are not available. Based on the comment of reviwer#3, we have compared the computation time of FEM and Analytical solution when applied to the Diffusion Equation. Analytical simulation of 36 channels, including a sampling rate of 30 samples per second, takes less than 10 minutes. For the same sampling rate and channels, the FEM simulation will take almost 43-48 hours, depending on system speed.

If we use an accurate forward model, then all this concern can be solved. It can be concluded that the results won't be decisive in the inversion procedure. With the help of insightful comments of reviewer #2, we have added instrumentational noise (Gaussian noise) to each channel, and then after inversion, the improvement obtained by MDMD has been preserved (please refer to the response of Point 2 of the second referee).

In order to cover major concerns in points 1-3, we have added several paragraphs along with simulations. We have added uncertainties to elements of the Jacobian matrix with random distribution. Then we have simulated the forward model. The contrast, after considering the random behavior, is also simulated in Fig 1. The changes are visible.

Fig 1: The contrast (depth sensitivity) in forward model by adding random coefficients.

Please refer to lines 331-357. Still the performance of MDMD is superior than others.

The simulation results illustrated in Fig 15 of the revised manuscript:

Fig 15. (a) Represents Gaussian distribution (mean=1 and δ=0.07) of Rnd matrix elements. (b) Boxplot regarding random coefficients. (c) Represents Gaussian distribution (mean=1 and δ=0.2) of Rnd matrix elements. (d) Boxplot regarding random coefficients.

Minor Points:

Point 1: At line 60 the single Ref. 15 is not enough, I would also add Ref. 16.

Response: Thanks. The reference was added to the revised version of the manuscript.

Point 2: At line 62 it is written “… simple slab boundary condition in semi-infinite geometry are utilized [16]”. The sentence is confusing in the sense that is confused the meaning of geometry and boundary condition. Here the term boundary condition appears misleading. To me would make sense to rephrase as “… the simple geometry, that can approximate a semi-infinite medium for thick slabs, is utilized [16]”.

Response: Thanks. Corrected.

Point3: At line 96 it is written “The analytical solution for perturbative DE has been obtained in the boundary condition of Fig. 3.” So far, I understand it should be “The analytical solution for perturbative DE has been obtained for the geometry of Fig. 3.” Boundary condition and geometry have a different meaning in this context. The actual boundary condition used to solve DE for the forward solver used in this paper is the extrapolated boundary condition that is not a geometry but a condition to make an energy balance at the external interface of the medium.

Response: Thanks. Corrected.

Point4: At line 102 the symbol S (The thickness of the slab) is not defined. NB that S is also the symbol used for the different regions S1-S9 (See Fig. 3).

Response: Thanks. Corrected “the thickness of slab is equal to 40mm” We have removed “S” to avoid confusion with region name.

Point5: At line 104 is written “Consequently, the volume of the inhomogeneity (inclusion) is regarded to be small concerning baseline optical properties of the homogeneous medium [16].” It should be “Consequently, the volume of the inhomogeneity (inclusion) is regarded to be small concerning baseline optical properties of the homogeneous medium [16].”

Response: Thanks Corrected.

Point6: At line 112 it is written “The final expression of Rpert(\\rho) for each channel derived based on Born approximation [21].” It looks like the sentence miss the final part. The sentence recalls the Born approximation; however, Ref. 21 is mainly related to higher order perturbation theory. Do the authors mean that in Ref. 21 the results for the Born approximation are also presented?

Response: Thanks. We have cited this ref because it has a comprehensive view of different theories including Born approximation.

Point7: At line 136 the acronym HRF is used without definition. I understand that implicitly means Hemodynamic Response Function, however why it should be omitted this definition?

Response: Thanks. We have defined this acronym several lines before (line 128). In Fig.2 of “Response to Reviewers” we have illustrated HRF signal. The signal looks like hemodynamic response function that is why it is called HRF in (Bonomini et al., 2015) it is also called HRF. This reference is one our main source of study.

〖∆HbO_2〗_i (t)=〖HRF〗_i (t)*s(t)+ ∅_phy (t) and ∆〖Hb〗_i (t)=-1/3×〖∆HbO_2〗_i (t)

Fig.2: The procedure of Oxy and Deoxy hemoglobin generation.

Point8: At line 272, Eq. (15), index of the sum in Eq. (15) is a real number, while has not been used an index numbering the number of Source-Detector pairs?

Response: Thanks. Agreed. That parameter is changed in the revised manuscript to x_SD (line 278).

Point10: The actual title of Ref. 16 is: “Light Propagation through Biological Tissue and other Diffusive Media: Theory, Solutions and Software” 2009.

Response: Thanks. Corrected.

Point11: In the first row of Fig. 4 is missed the info on the y axis (Z (mm)?).

Response: Done.

Point12: In the second row of Fig. 6 is missed the info on the x axis (Time s?).

Response: Done.

Reviewer#2:

Comments:

Author report on a simulation platform for diffuse optical tomography (DOT) adapted for reconstruction of hemodynamic responses. The platform is based on analytical solutions, under the Born approximation, of the diffuse equation for a semi-infinite slab. The work is technically sound and, together with the claimed availability of software and data, will be useful for setting up DOT systems.

Response: We appreciate the positive feedback of the esteemed reviewer. Of course, the simulation code and a simulation guide will be available according to PLOS One policies. Both of them will help the other researcher to set up a system for their purposes.

Major Concerns

Point 1: At my opinion there are the following points to address and clarify: Reconstruction: what I don’t understand is whether the reconstruction is forced at the depth of 15 mm or it is performed in the whole volume. In the first case, the point has to be emphasized and better specified in the text and, at my opinion, the method can’t be properly called “tomography”. In the second case, as well, it has to be emphasized in the text and, what I expect, is a figure representing slices in the volume at a defined time Ts to see the reconstructed depth that, in DOT, is typically underestimated.

Response: Thanks for the insightful comment of the reviewer. The reconstruction is forced at 15mm depth. However, by the simulation approach, it is possible to do tomography if we add other perturbations in different depth. Then it is required to generate a Jacobian matrix for each depth. We have corrected the “Tomography” to “Topography” in the manuscript.

Point 2: Noise and Inverse-crime: As far as I’ve understood, the only noise added is on the optical properties of the S5 and nearby voxels. This means that the forward model has been generated with well-defined optical properties using the Jacobian matrix, and the inversion has been performed using the same matrix on unnoisy data ΔR. This is typically called an “inverse-crime”. I suggest to add noise (Gaussian or Poisson) to the simulated data to avoid this problem.

Response: Thanks for the valuable comment. To avoid inverse-crime, we have added simulation and explanation of how noise affects the results (please refer to lines 309-330, Fig.13 and Fig.14):

Fig 13. Demonstrates the effect of noise on performance assessment of different SD arrangements in hemodynamic extraction of the S5 region concerning the regularization parameter.

Fig 14. (a) Represents the ability of hemodynamic extraction from the S5-region concerning several SNRs. (b) illustrates that in the worst case condition (SNR=32 dB), the optimum point for ϵ has changed from 〖10〗^(-5) to 〖10〗^(-4). (c) and (d) Represents reconstructed and synthetic ∆HbO_2 and ∆Hb for MDMD and HD-DOT respectively.

Minor Points

Point 1: the acronym HRF is not specified in the text.

Response: Corrected. Please refer to line 128.

Point 2: Figure 1 is reported with a very poor resolution, writings are too small.

Response: Thanks. We have increased the writing font. The resolution has been improved.

Reviewer#3:

Comments:

The authors present a simulation setup for the performance assessment of different source-detector configurations in high density diffuse optical tomography. The information content in the paper is well organised. While the work presented is interesting and very useful research, a few concerns are raised below, which are to be addressed prior to any publication.

Response: We appreciate the positive feedback of the esteemed reviewer.

Major Concerns

Point 1: One of the important aspects of the presented work is the significant reduction in computational time due to the use of analytical forward model, when compared to finite element and Monte Carlo methods. Therefore, it is crucial to provide the readers with the comparison of individual computational times for a sample forward model. Compare with the result of FEM

Response: Thanks for the insightful comment. We have provided a quantitative comparison between Analytical and FEM. Analytical simulation of 36 channels, including a sampling rate of 30 samples per second, takes less than 10 minutes. For the same sampling rate and channels, the FEM simulation (Fig 3 of “Response to Reviewers”) will take almost 43-48 hours, depending on system speed) We have excluded the computation time of 4 extra channels(. Please refer to lines 303-306.

Fig. 12: Geometry developed for FEM simulation of 40 channels

Fig. 3: ∆OD for 40 Channels in two wavelength

Point 2: From Figure 9, and Table-1, the authors want to show that higher the number of unique elements in Jacobian matrix, better will be the recovery. This can be misleading. While it is logical that a Jacobian with more independent information can give a better recovery, the correct way to observe this is to compare the normalised singular values of different Jacobian matrices corresponding to the SD arrangements. I would recommend looking into: Optics Letters Vol. 26, Issue 10, pp. 701-703 (2001).

Response: Thanks for the comment. SDV analysis is applied to the Jacobian matrix of different SD arrangement:

A singular value analysis of the Jacobian matrix associated with introduced SD arrangements is used as a benchmark (Culver, Ntziachristos, Holboke, & Yodh, 2001). Besides unique elements of the Jacobian matrix, the singular value analysis of the different methods in Fig 10 indicates that both shape of the singular value spectra and the magnitude for MDMD arrangement is higher than other SD combinations. Please refer to lines 253-256.

Fig 10. Singular value spectra for Jacobian matrix of HD-DOT, MDMD, Arrangement-1 and 2

Point 3: Figure 10 and 11 have same information, therefore figure 10 can be avoided and two additional sub-plots for SD separation 10 and 20mm can be included in Figure 11.

Response: Done. Fig. 10 has been avoided and Fig. 11 in the manuscript represents depth sensitivity of Analytical and MC including SD=10mm and SD=20mm.

Fig 11. Represents the depth sensitivity of the analytical pDE in Slab medium and depth sensitivity of the regression of MC on Colin 27 geometry for SD separation of 10mm-60mm.

Minor Points

Point 1: Add reference for lines 150-152.

Response: Thanks. The Reference is added.

Point 2: Use one abbreviation for oxy hemoglobin (either O2Hb or HbO2), and for deoxy hemoglobin (either HHb or Hb) throughout the manuscript.

Response: Thanks. Corrected.

Point 3: Important future directions such as the use of multi-layered medium, fit better in the conclusion section rather than discussions. Re-edit the conclusion section to incorporate this information seamlessly.

Response: Thanks. The future directions are removed from discussion section and moved to conclusion as follow:

Based on the modeling approach and simulation results, the achievements of the MDMD looks more beneficial than other methods. But still, there would be a gap between the outcomes of this study and real imaging applications. It is suggested to use the proposed simulation approach with a modified forward model according to the following suggestions. The brain model in this work is one-layer, while the multi-layered medium can be considered to get closer to the real results. Also, we must find the interaction between the perturbations in the solution of DE equations (whenever the perturbations are not small regarding the baseline optical properties), while Born approximation is used in this work. In this study the diffuse reflectance modulated due to change in absorption coefficient. The optical scattering coefficient changes along with the absorption coefficient should be taken into account for more accurate simulation of diffuse reflectance. Please refer to lines 389-398.

Reviewer#4:

Comments:

Major Concerns:

Point 1: In the manuscript, ‘Performance Assessment of High-Density Diffuse Optical Tomography Regarding Source-Detector Array Topology,’ the authors present a simulation study comparing various source-detector separation distances. While the methods are sound, the choice of methods (Born approximation, limited number of source-detector pairs and therefore a small field of view, analytical model instead of anatomy-based FEM using either diffusion or Monte Carlo, the single layer of optical properties) leads to the results and discussion providing limited information to the community for further advances in simulation, system design, or empirical considerations. As such, the authors are encouraged to add complexity to their modeling procedures or better motivate their choices to place their work in context with the current status of the field.

Response: Thanks for the positive feedback and concerns. According to the response to comments of esteemed referees. We have added complexity according to the reviewer’s suggestions:

The instrumentational Noise is added to each channel. The effect of Noise on the previous results are illustrated in Fig13 and Fig 14.

We have added uncertainties to elements of the Jacobian matrix with random distribution. The simulation results can be found in Fig 15 of the revised manuscript.

We have compared the computation time of FEM and Analytical methods quantitatively. Please refer to lines 303-306.

We have mentioned that why we could use Born approximation (please refer to the response to respected reviewer#1). Of course, the number of SDs is limited, but it can be expanded. Increasing the number of SDs to cover a large area of the brain will increase the computation time. But it worth to expand it. However, for the first step of the study, it is better to focus on a small area and discuss everything in detail. Imagine each SD arrangement as a single probe; if we sweep this probe spatially on the brain surface, then we will get the results of expansion, and we can have an expanded field of view.

If we want to use sophisticated slab medium or anatomical medium because of complicated geometry, the analytical solutions are not available (Strangman, Li, & Zhang, 2013). However, we are able to use FEM when it is applied to the Diffusion Equation (Okada & Delpy, 2003).

We do not recommend the use of Monte Carlo as a forward model in this simulation setup. MC can be used to study depth sensitivity in a real brain template for one-time research. However, if we add dynamic perturbation inside the brain template, we need to run, for example, 4000 run for each sampling rate. Instead of MC analysis, we recommend doing an experiment in tissue-like phantom. If we apply different arrangements to tissue-like phantom (as a forward model), there may be a chance to produce synthetic hemodynamic, but it is still challenging. The challenge is how to insert dynamic perturbations inside the phantom.

Authors Modification:

Fig 2 has been modified to better match with the simulation codes.

References

Bonomini, V., Zucchelli, L., Re, R., Ieva, F., Spinelli, L., Contini, D., … Torricelli, A. (2015). Linear regression models and k-means clustering for statistical analysis of fNIRS data. Biomedical Optics Express, 6(2), 615. https://doi.org/10.1364/BOE.6.000615

Chiarelli, A. M., Maclin, E. L., Low, K. A., Mathewson, K. E., Fabiani, M., & Gratton, G. (2016). Combining energy and Laplacian regularization to accurately retrieve the depth of brain activity of diffuse optical tomographic data. Journal of Biomedical Optics, 21(3), 036008. https://doi.org/10.1117/1.JBO.21.3.036008

Culver, J. P., Ntziachristos, V., Holboke, M. J., & Yodh, A. G. (2001). Optimization of optode arrangements for diffuse optical tomography: A singular-value analysis. Optics Letters, 26(10), 701–703.

Martelli, F., Del Bianco, S., & Zaccanti, G. (2005). Perturbation model for light propagation through diffusive layered media. Physics in Medicine & Biology, 50(9), 2159.

Okada, E., & Delpy, D. T. (2003). Near-infrared light propagation in an adult head model I Modeling of low-level scattering in the cerebrospinal fluid layer. Applied Optics, 42(16), 2906. https://doi.org/10.1364/AO.42.002906

Sassaroli, A., Martelli, F., & Fantini, S. (2009). Higher-order perturbation theory for the diffusion equation in heterogeneous media: application to layered and slab geometries. Applied Optics, 48(10), D62-73. https://doi.org/10.1364/AO.48.000D62

Strangman, G. E., Li, Z., & Zhang, Q. (2013). Depth Sensitivity and Source-Detector Separations for Near Infrared Spectroscopy Based on the Colin27 Brain Template. 8(8). https://doi.org/10.1371/journal.pone.0066319

Submitted filename: Response to Reviewers.DOCX

Click here for additional data file.

25 Feb 2020

Performance Assessment of High-Density Diffuse Optical Topography Regarding Source-Detector Array Topology

PONE-D-19-31493R1

Dear Dr. Setarehdan,

We are pleased to inform you that your manuscript has been judged scientifically suitable for publication and will be formally accepted for publication once it complies with all outstanding technical requirements.

Within one week, you will receive an e-mail containing information on the amendments required prior to publication. When all required modifications have been addressed, you will receive a formal acceptance letter and your manuscript will proceed to our production department and be scheduled for publication.

Shortly after the formal acceptance letter is sent, an invoice for payment will follow. To ensure an efficient production and billing process, please log into Editorial Manager at https://www.editorialmanager.com/pone/, click the "Update My Information" link at the top of the page, and update your user information. If you have any billing related questions, please contact our Author Billing department directly at authorbilling@plos.org.

If your institution or institutions have a press office, please notify them about your upcoming paper to enable them to help maximize its impact. If they will be preparing press materials for this manuscript, you must inform our press team as soon as possible and no later than 48 hours after receiving the formal acceptance. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information, please contact onepress@plos.org.

With kind regards,

Alberto Dalla Mora, Ph.D.

Academic Editor

PLOS ONE

Additional Editor Comments (optional):

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. If the authors have adequately addressed your comments raised in a previous round of review and you feel that this manuscript is now acceptable for publication, you may indicate that here to bypass the “Comments to the Author” section, enter your conflict of interest statement in the “Confidential to Editor” section, and submit your "Accept" recommendation.

Reviewer #1: All comments have been addressed

Reviewer #2: All comments have been addressed

Reviewer #4: All comments have been addressed

**********

2. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Yes

Reviewer #2: Yes

Reviewer #4: Yes

**********

3. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: Yes

Reviewer #2: Yes

Reviewer #4: Yes

**********

4. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: Yes

Reviewer #2: No

Reviewer #4: Yes

**********

5. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

Reviewer #2: Yes

Reviewer #4: Yes

**********

6. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: I have read carefully the authors response to my reports and to all the reports of the other reviewers together with the changes introduced in the revised version of the manuscript. Indeed, the manuscript has been significantly improved compared the previous version and the main points raised in the reports has been properly addressed. According to this fact, I recommend its publication on Plos One unaltered.

Reviewer #2: (No Response)

Reviewer #4: The authors have greatly improved the manuscript. The manuscript is now appropriate for this journal.

**********

7. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

Reviewer #2: No

Reviewer #4: No

9 Mar 2020

PONE-D-19-31493R1

Performance Assessment of High-Density Diffuse Optical Topography Regarding Source-Detector Array Topology

Dear Dr. Setarehdan:

I am pleased to inform you that your manuscript has been deemed suitable for publication in PLOS ONE. Congratulations! Your manuscript is now with our production department.

If your institution or institutions have a press office, please notify them about your upcoming paper at this point, to enable them to help maximize its impact. If they will be preparing press materials for this manuscript, please inform our press team within the next 48 hours. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information please contact onepress@plos.org.

For any other questions or concerns, please email plosone@plos.org.

Thank you for submitting your work to PLOS ONE.

With kind regards,

PLOS ONE Editorial Office Staff

on behalf of

Professor Alberto Dalla Mora

Academic Editor

PLOS ONE