CN120014155A

CN120014155A - A method for constructing a high-resolution 3D atlas of cerebral hemooxygen dynamic response

Info

Publication number: CN120014155A
Application number: CN202510018632.9A
Authority: CN
Inventors: 苑立波; 张珈豪; 谢耀庆; 杨世泰
Original assignee: Guilin University of Electronic Technology
Current assignee: Guilin University of Electronic Technology
Priority date: 2025-01-07
Filing date: 2025-01-07
Publication date: 2025-05-16

Abstract

The present invention discloses a method for constructing a high-resolution cerebral hemodynamic response 3D atlas, including: selecting a rectangular network or a honeycomb network to detect several brain regions, arranging several optical fiber probes; dividing the area below the basic unit into several local areas, calculating the time-resolved partial path length; executing a signal inspection strategy at wavelengths 1 and 2 to obtain the first-measured diffuse light signal of the two wavelengths; after a short interval, obtaining the re-measured diffuse light signal of wavelength 1; executing a basic unit with a basic unit signal joint analysis strategy to obtain a change in the absorption coefficient of the local area of wavelength 1; obtaining the re-measured diffuse light signal of wavelength 2; executing a basic unit with a basic unit signal joint analysis strategy to obtain a change in the absorption coefficient of the local area of wavelength 2; solving a set of equations for the same wavelength and local area to obtain a change in blood oxygen in the local area; and obtaining a 3D atlas of cerebral hemodynamic response. The present invention improves the spatial resolution of the 3D atlas of cerebral hemodynamic response, and provides a technical reserve for a full spectrum of 3D brain depth information.

Description

Construction method of high-resolution cerebral blood oxygen kinetics response 3D map

Technical Field

The invention belongs to the technical field of near infrared spectrum, and particularly relates to a construction method of a high-resolution cerebral blood oxygen dynamic response 3D map.

Background

With the continued maturation of near infrared spectroscopy, various targeted near infrared instruments have emerged, with functional near infrared spectroscopy (fNIRS) having great potential as a technology. fNIRS is a noninvasive technique which can be used for monitoring the change conditions of oxyhemoglobin and deoxyhemoglobin during brain activities, is currently in a continuously perfected and developed stage, and is mainly applied to the fields of advanced cognition, brain science, neurology, psychology and the like in natural environments. Brain tissue activity can lead to changes in blood oxygen content, which in turn affects the optical properties of blood oxygen, primarily in terms of absorption and scattering effects. fNIRS is based on the principle that changes in blood oxygen and blood volume in the brain are reflected by monitoring changes in the concentration of oxyhemoglobin and deoxyhemoglobin in the blood, to infer brain activity. Since red blood cells have good absorption properties for light having a wavelength below 650nm, whereas water absorbs light having a wavelength above 900nm, fNIRS employs near infrared light having a wavelength in the range of 650nm to 900nm [ Guo Jianghui, high temporal resolution fNIRS and its pulmonary disease application research [ D ]. University of electronic technology, 2023.doi:10.27005/d.cnki.gdzku.2023.003016 ]. Within this band, there is good detection capability for oxyhemoglobin and deoxyhemoglobin. The fNIRS system is mainly divided into three types, namely a first Continuous Wave (CW) spectrum system, a second Frequency Domain spectrum system (FD), and a third Time Domain spectrum system (TD).

The time domain fNIRS technique employs incident light that is a pulse with a short pulse width, typically on the order of picoseconds (1 picosecond = 10 ^-12 seconds). Typically the measurement of the light intensity in the time domain fNIRS is by recording the time of the photons by a high sensitivity detector or a single photon counter. Since the distance light travels in the tissue is proportional to the flow time and diffuse reflection in the tissue is random, there may be different light intensities at different times (typically hundreds of picoseconds) after the light originates from the light source. When the detector receives a large number of photons, the transmission distance of the photons generates a probability distribution, and the more the number of photons is, the more obvious the distribution is, and accordingly, the distribution of the light intensity of the outgoing light with time can be obtained, and the distribution is called a time Point spread function (Temporal Point SpreadFunction, TPSF) in optics. By analyzing the photon time-of-flight coding information contained in the TPSF, the optical parameters of the tissue can be obtained, so that the time-domain measurement method for the HbO and HbR concentration variation [Torricelli A,Contini D,Pidderi A,et al.Time domain functional NIRS imaging for humanbrain mapping[J/OL].NeuroImage,2014,85:28-50.]. in the tissue is suitable for each acquisition channel, namely each source-detector pair, and does not need to perform multi-channel detection analysis like a continuous wave near infrared spectrum technology, and currently common measurement methods include a fitting method, a moment calculation method and a gating method (also called a time gating or time window method). The principle of the gating method is to evaluate the hemodynamic changes of different areas of the brain by using the separation of TPSF signals in different time gates to represent different photon arrival times and the influence of multi-layer tissue characteristics on the flight path of the probe light with time [Lange,Frédéric,and Ilias Tachtsidis."Clinical brain monitoring with time domain NIRS:a review and future perspectives."Applied Sciences 9.8(2019):1612.].

The average partial path length (MEAN PARTIAL PATHLENGTHS, MPP) is the average flight path length of diffuse light received in the fNIRS measurement over different media, representing the sensitivity of the measured optical density to absorption changes in the layer, and Steinbrink et al extended the average partial path length to time domain diffuse light measurements in 2000, suggesting that the time dependent average partial path length (Time Dependent Mean Partial Pathlengths,TMPP)[Steinbrink,J.,et al."Determining changes in NIR absorption using a layered model of the human head."Physics in Medicine&Biology 46.3(2001):879.].TMPP describes the effect of absorption changes of different layers on time resolved reflectivity, useful for quantifying the effect of multi-layer tissue properties on the flight path of probe light over time.

In the study of diffusion of light in a multilayer turbid medium, the theoretical expression of each TMPP is typically calculated by resolving the target medium using a solution to the radiation transmission equation approximated by the diffusion equation. However, as the composition of the turbid medium and the optical parameters change more complex, the above solution method is no longer robust. A large number of documents calculate TMMP[Vera,Demián A.,et al."Retrieval of chromophore concentration changes in a digital human head model using analytical mean partial pathlengths ofphotons."Journal ofBiomedical Optics 29.2(2024):025004-025004.]. monte carlo photon transmission simulation by using monte carlo photon transmission simulation (Monte Carlo photon transport simulation), which is a computer simulation technology and is widely applied in the fields of biological optical application, medical imaging, photoelectric equipment design and the like. The basic idea is that a large number of photons are randomly generated, motion trail simulation is carried out according to a certain physical rule, and parameters such as reflection, transmission, absorption and the like under different wavelengths on the surface or inside an object are calculated through statistical analysis of the trail, so that information such as the tissue structure, the composition or the morphology of the object is deduced. Monte Carlo simulation is used as a stochastic method of solving radiation transfer equations, which in the general case govern photon migration, where light can not only be absorbed, but also scattered in an anisotropic manner. Thus, MC simulations are gold standards for simulating light transmission experiments in multi-layered turbid media, as they enable accurate control of all parameters involved in photon migration processes in turbid media. In time-domain based MC simulation, TMPP of different layers can be obtained by time-resolved measurement and formula derivation.

A key issue in near infrared photon migration applications is how to increase the penetration depth of the measurement and reduce the measurement area simultaneously, especially for probing the brain. Scalp, skull and cerebrospinal membranous layers attenuate the optical changes produced by the cerebral cortex and when probing some human tissue (e.g., breast), local localization is also required to detect and characterize small deep focal spots. Time-resolved diffuse reflectance measurements of zero source detector distance solve the above problems to a certain extent compared to longer source detector distance time-resolved diffuse reflectance measurements, which are advantageous in that (1) the number of photons collected at any arrival time is increased, (2) higher contrast is generated in local areas, (3) better spatial resolution is provided, (4) brain tissue regions of interest are easily located [Pifferi,Antonio,et al."Time-resolved functional near-infrared spectroscopy at null source-detector separation."Biomedical Optics.OpticaPublishing Group,2008.].

The self-detection-inter-detection combined time domain fNIRS is based on a bidirectional functional diffusion time-resolved optical fiber brain-computer interface probe and a single photon avalanche diode with extremely high time resolution, has the self-detection function of a single optical fiber probe and the inter-detection function of a double optical fiber probe, can effectively improve the transverse and longitudinal spatial resolutions of detection through time resolution, can enhance the detection effect of brain tissue depth, greatly improves the resolution of a traditional functional near-infrared brain hemodynamic response diagram, and enables the acquisition of a high-resolution dynamic brain blood flow and brain blood oxygen 3D diagram to be possible.

The fNIRS optical fiber probes are arranged in a certain space mode and a certain distance from the source detector, so that a wide brain area is detected, and the method can be applied to single-layer or fault detection imaging of brain tissues. fNIRS of tomographic imaging, also known as diffuse optical tomography (Diffuse Optical Tomography, DOT), can provide higher spatial resolution 3D images than fNIRS of conventional 2D tomographic imaging. DOT obtains time and space distribution information of optical signals by measuring diffusion light emitted from the boundary of biological tissues, establishes a forward model based on an optical transmission process, and performs inversion reconstruction on original data to obtain optical parameters (such as absorption coefficients and scattering coefficients) of different spatial positions inside the tissues, so that tomography is performed on the tissues to obtain information of hemoglobin content, blood oxygen saturation, blood volume and the like. Because diffuse light is used, the spatial resolution of DOT is relatively low (typically 1cm or less), and imaging quality is to be improved. The DOT of the continuous wave measurement mode requires fNIRS equipment supporting flexible alignment of a large number of optical fibers and is often difficult to apply and analyze [ Song Bowen, zhao Yanyu ] diffuse optical imaging methods and applications (special solicitation) [ J ]. Laser and optoelectronics advances 2024,61 (08): 11-36 ]. The time domain fNIRS combined with self-detection and mutual detection is applied to DOT, so that the number and spacing limit of optical fiber probes can be reduced, the transverse and depth resolution and detection depth of DOT can be improved, and a high-resolution brain nerve activity and blood oxygen dynamics response map can be obtained.

The geometrical arrangement of probes of the fNIRS system is mainly square arrangement and diamond arrangement, and the arrangement can increase the utilization rate of a transmitting probe and a receiving probe [ Ma Pei, shen Moshuang, shen Huijuan, and the like ]. A functional near-infrared brain imaging system is reviewed [ J ]. An optical instrument, 2022,44 (05): 1-13 ]. However, the probe geometry of most commercial multi-channel fNIRS systems is not suitable for exploring or locating functional activation of the human cerebral cortex, which provides a grid of points spaced about 30mm apart, whereas the human cerebral cortex is composed of about 10-15 mm wide gyres that can be used to demonstrate independent functions. The detection signals of each pair of source and detector with the interval of 30mm have too low discrimination of the functional brain areas of the cerebral cortex, and the depth specificity is still to be improved. Therefore, such multi-channel measurements may ignore the activity of partial cortical gyros (the detection depth of the false negative error )[Yamada,Toru,et al."Wearable fNIRS device of higher spatial resolution realized by triangular arrangement of dual-purpose optodes."Neural Imaging and Sensing 2023.Vol.12365.SPIE,2023.]. time domain fNIRS is not limited by the fiber probe spacing, the zero source detector distance and the small source detector distance are integrated into the geometric arrangement of the probe, which can be used to increase channel density and spatial resolution, providing accurate local brain region localization.

In order to overcome the defects of the prior art, a more accurate high-resolution brain blood oxygen dynamic response 3D map is constructed, a detection network is selected according to channel density, a plurality of brain areas are detected, a plurality of optical fiber probes of an optical fiber brain blood oxygen detection system are deployed in the detection network, a proper detection network is selected according to the calculated channel density of different detection networks, brain area areas corresponding to basic units are divided into local areas, the path length of a time resolution part is solved, a signal inspection strategy is adopted to obtain time resolution reflectivity signals, a primitive signal joint analysis strategy is executed for each basic unit, curve inversion and absorption change equation set solving are carried out on the sequentially obtained time resolution reflectivity signals, local area brain blood oxygen change values are obtained, a brain blood oxygen dynamic response three-dimensional RGB map is constructed, a brain blood oxygen dynamic response 3D map is obtained, a brain tissue blood oxygen dynamic response 3D map is obtained, and a brain depth information full-target reserve technology is provided.

Disclosure of Invention

The invention aims to provide a construction method of a high-resolution cerebral blood oxygen dynamic response 3D map, which aims to improve the spatial resolution and the detection depth of a time domain fNIRS by utilizing a zero-spacing source-detector optical fiber probe of an optical fiber cerebral blood oxygen detection system and multi-channel time resolution measurement from the surface to the domain so as to obtain the high-resolution cerebral blood oxygen dynamic response 3D map.

In order to achieve the aim of the invention, the invention adopts the following technical scheme:

a method for constructing a high-resolution cerebral hemodynamic response 3D map, comprising the steps of:

1) Selecting N _brain brain regions to be detected, surrounding the brain regions corresponding to the regions of brain scalp by adopting a detection network, wherein the detection network is a rectangular network with the size of MxN taking a rectangle as a basic unit or a hexagonal honeycomb network with the circumference of C taking a regular triangle as a basic unit, for the rectangular network, setting the rectangle with the basic unit of M ₁xN₂, M ₁ being less than or equal to M, N ₂ being less than or equal to N, for the hexagonal honeycomb network, setting the regular triangle with the basic unit of C ₁, C ₁ being less than or equal to C, each intersection point in the detection network being a probe arrangement point, for placing an optical fiber probe of an optical fiber cerebral blood oxygen detection system, defining that absorption scattering occurs when the optical fiber probe emits into the brain regions, the optical path received by another optical fiber probe is a mutual detection channel and the absorption scattering occurs when the optical fiber probe emits into the brain regions, the optical path received by the same optical fiber probe is a self-detection channel, and the interval ρ of the channel is the interval between an emitting optical fiber and a receiving optical fiber;

2) The channel densities of rectangular and hexagonal cellular networks are calculated as:

wherein: Is the number of channels of the ith basic unit, S _i (j) is the surface area of the ith basic unit covering the jth brain region, and when S _i (j) is zero, let S _i (j) be 1;

Setting a plurality of optical fiber probes on the rectangular network and the hexagonal honeycomb network respectively, if the channel density of the rectangular network is greater than or equal to that of the hexagonal honeycomb network, adopting the rectangular network as a detection network, otherwise adopting the hexagonal honeycomb network as the detection network;

3) If a rectangular network is adopted, taking a brain region below each basic unit as a rectangular region, dividing the rectangular region into n partial regions, dividing the partial regions into n-1 measuring regions and other regions, wherein the thickness of each measuring region is not less than 1mm, not more than 1cm, the rest regions except the measuring regions are other regions, if a hexagonal honeycomb network is adopted, taking a square prism region below each basic unit, the square prism region is not more than three times the side length of the basic unit, not less than two times the side length of the basic unit, dividing the square prism region into n partial regions, dividing the partial regions into n-1 measuring regions and other regions, and dividing the thickness of each measuring region into n-1 measuring regions and other regions, wherein the thickness of each measuring region is not less than 1cm, and the rest regions are removed;

4) Setting seven parameters including length, width, thickness, refractive index, absorption coefficient, scattering coefficient and anisotropic factor in each local area, on the basis, calculating time-dependent average partial path length L _i,j (rho, lambda, t) of the jth local area of different channel spacing rho, wavelength lambda and detection time t of different basic units i by adopting Sobber sampling through a time-resolved Monte Carlo simulation method based on a voxel structure, wherein the number of simulated photons is 10 ⁸;

5) The signal inspection strategy is as follows, one of the optical fiber probes is selected as a core probe, a laser diode coupled with the core probe is driven to emit near infrared pulse light with the light intensity S and the dominant wavelength lambda, meanwhile, the optical fiber probes in all basic units comprising the core probe are enabled to receive brain area diffuse light signals, the time for receiving the signals is set to be the integral multiple of T _r,T_r, after the signals are received, all optical fiber probe numbers for receiving brain area diffuse light signals are put into a waiting queue in a clockwise sequence, the queue head is taken out, the optical fiber probe corresponding to the queue head number is enabled to serve as the core probe, the operation is repeated until no optical fiber probe number exists in the queue, if the deployed optical fiber probe is not used as the core probe to emit pulse light, the optical fiber probe is enabled to serve as the core probe, the operation is repeated until all the optical fiber probes are used as the core probe to emit pulse light, the time T _λ＝T_r xR for measuring brain area diffuse light signals is calculated, the pulse light dominant wavelength of the laser diode is set to be lambda ₁, the signal inspection strategy is executed to obtain brain area diffuse light signals measured for the first time, the light probe length lambda ₁ is adjusted to be lambda ₂, and the light signal for measuring the brain area diffuse light signals for the first time is processed to obtain the light signal with the first measurement strategy;

6) After the interval delta T, repeatedly executing a signal inspection strategy at the wavelength lambda ₁ for one time to obtain brain region diffuse light signals of the second measurement of the wavelength lambda ₁, and converting the brain region diffuse light signals of all channels of the wavelength lambda ₁ into time resolution reflectivity signals by using software;

7) Executing a primitive signal joint analysis strategy on each basic unit to obtain the absorption coefficient change value of the local area of each basic unit under the wavelength lambda ₁;

8) Repeatedly executing a signal inspection strategy once at the wavelength lambda ₂ to obtain brain area diffuse light signals measured for the second time at the wavelength lambda ₂, converting the brain area diffuse light signals of all channels at the wavelength lambda ₂ into time resolution reflectivity signals by using software, and executing a primitive signal joint analysis strategy on each basic unit to obtain the absorption coefficient change value of the local area of each basic unit at the wavelength lambda ₂;

9) The relationship between the absorption coefficient change Δu _a,j (λ) and the blood oxygen concentration in the local region is expressed as:

Δu_a,j(λ)＝ε_HbO(λ)Δ[HbO]_j+ε_HbR(λ)Δ[HbR]_j (1)

Wherein ε _HbO (λ) is the molar extinction coefficient of oxyhemoglobin at wavelength λ, ε _HbR (λ) is the molar extinction coefficient of deoxyhemoglobin at wavelength λ, ΔHbO _j is the change in the concentration of oxyhemoglobin in the jth local region, and ΔHbR _j is the change in the concentration of oxyhemoglobin in the jth local region;

The absorption coefficient changes obtained in step 7) and step 8) are expressed by formula (1), and the equations for the same wavelength and local area are combined:

and solving a simultaneous equation set to obtain the change of the concentration of the oxyhemoglobin and the deoxyhemoglobin of all local areas of each basic unit.

10 Constructing a brain three-dimensional map file based on an AAL template or a Brodmann partition template, setting the brain surface, the brain region position and the brain region edge, adjusting the output layout size, the background color, the surface transparency, the brain region color, the brain region edge color and the image resolution of the three-dimensional map, then further dividing the region of the corresponding brain region in the brain three-dimensional map according to the division of the local region of the detected brain in the step 3), enabling each local region to display the absorption coefficient change normalized RGB value, finally obtaining a brain hemodynamic response three-dimensional RGB map, and dividing the brain hemodynamic response depth map layer by layer according to the brain hemodynamic response three-dimensional RGB map, thus obtaining the 3D brain hemodynamic response map.

Further, the delta T should be greater than or equal to the time of solving the absorption coefficient variation value minus the wavelength lambda ₁ to measure the diffuse light signal of the primary brain regionAnd measuring a primary brain region diffuse light signal at a wavelength lambda ₂ A kind of electronic device.

The system comprises a shell, a double-fiber collimator, a reflector and a focusing lens, wherein the system is used for detecting a time-resolved local area from a detection channel, the driving circuit transmits a control signal to the light source selection circuit, and simultaneously transmits a timing signal to the TCSPC, the light source selection circuit enables the pulse laser 1 or the pulse laser 2 to transmit light with a specified wavelength, the light signal reaches a transmitting optical fiber connected with the optical fiber probe through the optical fiber coupler and the optical fiber splitter, then the light signal reaches a reflecting surface of the collimator lens and the reflector through the MPO optical fiber jumper, then the light beam angle adjustment is carried out by the focusing lens, the light signal is injected into the brain cortex, the light signal scattered by brain tissue is absorbed by the brain cortex, and the same focusing surface and the SPC optical fiber jumper are connected with the single-photon counter through the single-photon collimator lens, and the SPC optical fiber optic fiber bridge, and the SPC optical fiber bridge is then the SPC optical fiber bridge is connected with the SPC optical fiber to the SPC optical fiber detector, and the SPC optical fiber bridge is then the SPC optical fiber bridge is amplified, and the SPC signal is transmitted to the PC after the SPC optical fiber signal is amplified, and the SPC signal reaches the single-photon counter and the SPC optical fiber detector through the single-photon counter after the SPC optical fiber bridge, and the SPC optical fiber bridge is connected with the SPC optical fiber.

Further, the elementary signal joint analysis strategy specifically comprises the following steps:

The method comprises the following steps of executing the following operations on channel signals in a basic unit, if the basic unit deploys four optical fiber probes, extracting time resolution reflectivity signals of 12 mutual detection channels and 4 self-detection channels between the optical fiber probes, dividing each time resolution reflectivity signal into m time windows, 16m is more than or equal to n, if the basic unit deploys three optical fiber probes, extracting time resolution reflectivity signals of 6 mutual detection channels and 3 self-detection channels between the optical fiber probes, dividing each time resolution reflectivity signal into m time windows, 9m is more than or equal to n, if the basic unit deploys two optical fiber probes, extracting time resolution reflectivity signals of 2 mutual detection channels and 2 self-detection channels between the optical fiber probes, dividing each time resolution reflectivity signal into m time windows, 4m is more than or equal to n, if the basic unit deploys one optical fiber probe, dividing the time resolution reflectivity signal of 1 self-detection channel of the optical fiber probe into m time windows by using a sliding window method, sliding step length Deltaw, and calculating the average path length of a local path of a part of a [ lambda ] j and a local path of a [ lambda ] t (57 t, a time average area of a [ rho 25 ]) of the time of the basic unit:

Wherein R (K _z,K_l,N_l,γ_z,γ_l, ρ, t) is a time resolution reflectivity multilayer approximate diffusion equation of the channel, and K _z,K_l,γ_z,γ_l,N_l is a equation set solution obtained by solving the approximate diffusion equation by using an eigenfunction method;

And carrying out curve inversion and absorption change equation set solving on the time resolution reflectivity signals measured for the first time and the second time of each channel, and carrying out absorption change value solving on the obtained linear equation sets of all channels.

Further, the curve inversion and absorption change equation set solving method specifically comprises the following steps:

measuring the response function s (t) of the system, taking an initial time-resolved reflectivity signal R ₀ (ρ, λ, t), and obtaining an optimization curve y (t) by convolution, namely Calculating a mean square fitting error:

∑(R_m(ρ,λ,t)-y(t))²(4)

Wherein R _m (ρ, λ, t) is the measured time resolved reflectivity signal;

The square sum of the mean square fitting error is minimized through a least square method, iteration times are set to hundreds of thousands times, a threshold value is set to be T, and when the square sum of the mean square fitting error is smaller than the threshold value T, a fitting time resolution reflectivity signal R (rho, lambda, T) is obtained and is used for theoretical calculation of subsequent absorption coefficient change;

The width of each time window of the fitted time-resolved reflectivity signal R (ρ, λ, t) is set to Δt, the first and last several time windows are truncated, and the number of detected photons per time window is converted using the modified lambert beer law to describe the small independent absorption coefficient variation of the local region:

Where R _g is the number of photons detected in the g-th time window of the fitted time-resolved reflectance signal of the second measurement, R _0,g is the number of photons detected in the g-th time window of the fitted time-resolved reflectance signal of the first measurement, deltau _a,j (lambda) is the change in absorption coefficient of the j-th local area relative to the second measurement of the first measurement, L _g,i,j (ρ, lambda) is the channel spacing ρ of the base unit i, the time-dependent average partial path length of the j-th local area in the g-th time window of the wavelength lambda, and the calculation formula is:

Wherein L _i,j (ρ, λ, t) is the time-dependent average partial path length of the jth partial region of the basic unit i, the channel spacing ρ, the wavelength λ and the detection instant t;

Converting formula (5) into:

Obtaining a linear equation set consisting of a relational expression (7) corresponding to all time windows:

where L _g,i,j (ρ, λ) is the channel spacing ρ of the basic unit i, the time-dependent average partial path length of the jth partial region within the jth time window of wavelength λ.

Further, the absorption change value solving specifically includes the following steps:

Combining linear equation sets of all channels obtained by curve inversion and absorption change equation set solution into a linear equation set:

wherein: For the j-th local area of the g-th time window of the c-th channel of the basic unit i, R _c,g is the number of photons detected in the g-th time window of the second measurement c-th channel, R _c,0,g is the number of photons detected in the g-th time window of the first measurement c-th channel, ρ _c is the spacing of the c-th channel, λ is the dominant wavelength of the near infrared pulse light;

obtaining a coefficient matrix of formula (9):

If the rank of the coefficient matrix is equal to the rank of the augmentation matrix and the rank of the coefficient matrix is equal to the number n of the local areas, the value of the absorption coefficient change Deltau _a,j (lambda) is solved, otherwise, according to the previously solved absorption coefficient change value, deltau _a,j (lambda) is predicted based on a Swin-transducer prediction model.

Compared with the traditional single-channel large-distance source detector distance time domain fNIRS detection, the invention adopts the zero-space source-detector optical fiber probe of the optical fiber cerebral blood oxygen detection system to realize the multi-channel light diffusion signal joint analysis from the surface to the domain, improves the transverse and longitudinal spatial resolution of the local area of the brain tissue, and obtains the high-resolution cerebral blood oxygen dynamic response 3D map.

Drawings

In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings used in the description of the embodiments or the prior art will be briefly described, and it is also possible for a person skilled in the art to obtain other drawings from these drawings without inventive effort.

Fig. 1 is a schematic diagram of a network configuration of a hexagonal honeycomb network using triangles as basic units, and a schematic diagram of a rectangular network configuration using rectangles as basic units.

FIG. 2 is a schematic representation of the process of obtaining a 3D map of the hemodynamic response of the brain according to the present invention.

FIG. 3 is a schematic diagram of the positions of the optical fiber probe corresponding to the brain region in the embodiment.

Fig. 4 is a schematic diagram of a process of obtaining a cerebral hemodynamic response 3D map through a detection network in an embodiment (a) a tiled image of a position of an optical fiber probe placed in the embodiment corresponding to a brain region, (b) a schematic diagram of time-resolved multi-channel local area detection in a basic unit, and (c) a cerebral hemodynamic response 3D map obtained through the detection network.

FIG. 5 is a flow chart for constructing a high resolution cerebral hemodynamic response 3D map.

Detailed Description

The present invention will be further described in detail with reference to the drawings and examples, which are only for the purpose of illustrating the invention and are not to be construed as limiting the scope of the invention.

As shown in the figure, the construction method of the high-resolution cerebral hemodynamic response 3D map of the embodiment comprises the following steps:

1) Selecting a dorsal forehead cortex (Brodmann zone 9) and a forehead zone (Brodmann zone 10) as brain zones to be detected, and surrounding the area of the brain zone corresponding to the brain scalp (the widest part of the dorsal forehead cortex is about 58mm and the widest part of the forehead zone is about 31 mm) by adopting a detection network, wherein the detection network is pre-selected as a rectangular network with the size of 65mm x58mm taking a rectangle as a basic unit, the basic unit is a rectangle with the size of 10mm x10mm, or a hexagonal honeycomb network with the circumference of 25.6mm taking a regular triangle as a basic unit, and the basic unit is a regular triangle with the side length of 1.07 mm;

Setting 29 and 36 optical fiber probes on the rectangular network and the hexagonal honeycomb network respectively, and calculating the channel density of the rectangular network and the hexagonal honeycomb network, wherein the channel density of the rectangular network is larger than that of the hexagonal honeycomb network, so that the rectangular network is used as a detection network;

Defining that the optical path emitted by the optical fiber probe and entering the brain region is a mutual detection channel and the optical path emitted by the optical fiber probe and entering the brain region is an absorption scattering channel, wherein the optical path received by the same optical fiber probe is a self-detection channel, the distance rho of the channel is the distance between the emitting optical fiber and the receiving optical fiber, namely 1.414cm or 1cm, numbering all the optical fiber probe arrangement points, setting the time resolution of a single photon avalanche diode coupled with the optical fiber probe to be 50ps, and arranging 19 optical fiber probes at the probe arrangement points with different numbers;

3) Taking a brain region area below each basic unit as a rectangular body area, wherein the length, width and height of the rectangular body area are 3cm, dividing the rectangular body area into 28 local areas, dividing the local area into 27 measurement areas and other areas, wherein the thickness of each measurement area is 5mm, and the rest areas except the measurement areas are other areas;

5) The signal inspection strategy is as follows, one of the optical fiber probes is selected as a core probe, a laser diode coupled with the core probe is driven to emit near infrared pulse light with the light intensity of 40mW and the dominant wavelength of lambda, meanwhile, the optical fiber probes in all basic units comprising the core probe are enabled to receive brain area diffuse light signals, the time for receiving the signals is set to be 5ns, after the signals are received, all optical fiber probe numbers for receiving the brain area diffuse light signals are put into a waiting queue in a clockwise sequence, the head of the queue is taken out, the optical fiber probes corresponding to the head of the queue are enabled to be used as the core probe, the operation is repeated until no optical fiber probe numbers are in the queue, if the deployed optical fiber probes are not used as the core probe to emit pulse light, the optical fiber probes are enabled to be used as the core probe to emit pulse light, the time T _λ for measuring the brain area diffuse light signals is set to be 735 ns, the primary wavelength of the pulse light of the laser diode is set to be 735nm, the signal inspection strategy is performed once, the first measured brain area diffuse light signals with the wavelength of 735nm are obtained, the pulse light of the laser diode is adjusted to be 905nm, and the primary signal inspection strategy is performed once, and the first measurement strategy of the brain area diffuse light signals with the wavelength of 905nm is obtained once;

6) Repeatedly executing a signal inspection strategy once at the wavelength of 735nm after the interval of delta T to obtain brain region diffuse light signals measured for the second time at the wavelength of 735nm, and converting the brain region diffuse light signals of all channels at the wavelength of 735nm into time resolution reflectivity signals by using software;

7) Executing a primitive signal joint analysis strategy on each basic unit to obtain the absorption coefficient change value of the local area of each basic unit under the wavelength of 735 nm;

8) Repeatedly executing a signal inspection strategy once at the wavelength of 905nm to obtain brain area diffuse light signals measured for the second time at the wavelength of 905nm, converting the brain area diffuse light signals of all channels at the wavelength of 905nm into time resolution reflectivity signals by using software, and executing a primitive signal joint analysis strategy on each basic unit to obtain the absorption coefficient change value of the local area of each basic unit at the wavelength of 905 nm;

8) The relationship between the absorption coefficient change Δu _a,j (λ) and the blood oxygen concentration in the local region is expressed as:

Δu_a,j(λ)＝ε_HbO(λ)Δ[HbO]_j+ε_HbR(λ)Δ[HbR]_j (1)

9) Constructing a brain three-dimensional map file of a Brodmann partition template, setting the brain surface, the brain region position and the brain region edge, adjusting the output layout size, background color, surface transparency, brain region color, brain region edge color and image resolution of the three-dimensional map, then further dividing the region of the corresponding brain region in the brain three-dimensional map according to the division of the local region of the detected brain region in the step 3), enabling each local region to display an absorption coefficient change normalized RGB value, finally obtaining a brain blood oxygen dynamics response three-dimensional RGB map, and dividing a brain blood oxygen dynamics response depth map layer by layer according to the brain blood oxygen dynamics response three-dimensional RGB map to obtain a brain blood oxygen dynamics response 3D map.

Specifically, ΔT should be greater than or equal to the time of solving the absorption coefficient variation value minus the wavelength λ ₁ to measure the primary brain region diffuse light signalAnd measuring a primary brain region diffuse light signal at a wavelength lambda ₂ And (2) solving time for estimating the change value of the absorption coefficient according to the operation environment of the software is 2-5ms, and let DeltaT be 6ms.

Specifically, the optical fiber cerebral blood oxygen detection system comprises a driving circuit, a light source selection circuit, a transmitting part consisting of a pulse laser 1 and a pulse laser 2, an optical fiber coupler, an optical fiber branching device, an MPO optical fiber jumper wire, an optical fiber probe, an optical fiber combiner, a single photon counter, a driving circuit 2, a receiving part consisting of a time-dependent single photon counter (TCSPC) and a SPAD and a PC signal processing end, wherein the optical fiber probe consists of a shell, a double optical fiber collimator, a reflector and a focusing lens; the system performs zero source-detector distance time resolution local area detection through a self-detection channel by the following process that a driving circuit transmits a control signal to a light source selection circuit and simultaneously transmits a timing signal to a TCSPC, the light source selection circuit enables a pulse laser 1 or a pulse laser 2 to emit light with specified wavelength, the light signal reaches an emitting optical fiber connected with an optical fiber probe through an optical fiber coupler and an optical fiber branching device and then reaches an emitting optical fiber connected with the optical fiber probe through an MPO optical fiber jumper, and then after passing through a collimating lens and a reflecting surface of a reflector, a focusing lens performs light beam angle adjustment, the method comprises the steps of injecting into cerebral cortex, absorbing scattered partial light signals by brain tissue, returning to an MPO optical fiber jumper wire through a receiving optical fiber connected with an optical fiber probe after passing through the same focusing lens, reflecting surface and collimating lens, reaching a single photon counter through an optical fiber combiner, sending a stop signal to TCSPC after detecting the signal of the single photon counter by a drive circuit 2, sending the single photon signal to SPAD by TCSPC, and sending the signal to a PC signal processing end for detection data processing after SPAD conversion and amplification.

Specifically, the primitive signal joint analysis strategy specifically comprises the following steps:

The method comprises the following operations on channel signals in a basic unit, if the basic unit deploys four optical fiber probes, extracting time resolution reflectivity signals of 12 mutual detection channels and 4 self-detection channels between the optical fiber probes, dividing each time resolution reflectivity signal into 3 time windows, if the basic unit deploys three optical fiber probes, extracting time resolution reflectivity signals of 6 mutual detection channels and 3 self-detection channels between the optical fiber probes, dividing each time resolution reflectivity signal into 4 time windows, if the basic unit deploys two optical fiber probes, extracting time resolution reflectivity signals of 2 mutual detection channels and 2 self-detection channels between the optical fiber probes, dividing each time resolution reflectivity signal into 7 time windows, if the basic unit deploys one optical fiber probe, extracting time resolution reflectivity signals of 1 self-detection channel of the optical fiber probes, dividing the time resolution reflectivity signals into 28 time windows by using a sliding window method, sliding step length of 100ps, and recalculating a local path length of a j (57L, 57 m) of a j-th local path of a region of a channel distance ρ, a wavelength λ and a detection moment t of the basic unit by using a sliding window method:

Specifically, the curve inversion and absorption change equation set solving method specifically comprises the following steps:

∑(R_m(ρ,λ,t)-y(t))²(4)

Wherein R _m (ρ, λ, t) is the measured time resolved reflectivity signal;

The square sum of the mean square fitting error is minimized through a least square method, iteration times are set to hundreds of thousands times, the threshold value is 0.001, and when the square sum of the mean square fitting error is smaller than the threshold value 0.001, a fitting time resolution reflectivity signal R (rho, lambda, t) is obtained and is used for theoretical calculation of subsequent absorption coefficient change;

The width of each time window of the fitted time-resolved reflectance signal R (ρ, λ, t) is set to 500ps, the beginning and the last one are truncated, and the number of detected photons for each time window is converted using the modified lambert beer law to describe the small independent absorption coefficient variation of the local region:

Converting formula (5) into:

Specifically, the solution of the absorption change value specifically includes the following steps:

obtaining a coefficient matrix of formula (9):

If the rank of the coefficient matrix is equal to the rank of the augmentation matrix and the two are equal to the number of the local areas 28, the value of the absorption coefficient change Deltau _a,j (lambda) is solved, otherwise, according to the previously solved absorption coefficient change value, deltau _a,j (lambda) is predicted based on a Swin-transducer prediction model.

The foregoing has outlined and described the basic principles, features, and advantages of the present invention. It will be understood by those skilled in the art that the present invention is not limited to the embodiments described above, and that the above embodiments and descriptions are merely illustrative of the principles of the present invention, and various changes and modifications may be made without departing from the spirit and scope of the invention, which is defined in the appended claims. The scope of the invention is defined by the appended claims and equivalents thereof.

Claims

1. A method for constructing a high-resolution 3D atlas of cerebral hemooxygen dynamics response, characterized in that the method steps include:

1) Select n _brain regions to be detected, and use a detection network to surround the area of the cerebral scalp corresponding to the brain region. The detection network is a rectangular network with a size of MxN using a rectangle as a basic unit or a hexagonal honeycomb network with a perimeter of C using an equilateral triangle as a basic unit. For the rectangular network, the basic _unit is set to be a rectangle of _M1xN2 , _M1 is less than or equal to M, and _N2 is less than or equal to N. For the hexagonal honeycomb network, the basic unit is set to be an equilateral triangle with a perimeter of _C1 , _C1 is less than or equal to C. Each intersection in the detection network is a probe arrangement point for placing the optical fiber probe of the optical fiber brain blood oxygen detection system. The optical path emitted by the optical fiber probe into the brain region, which undergoes absorption and scattering, and is received by another optical fiber probe is defined as a mutual detection channel, and the optical path emitted by the optical fiber probe into the brain region, which undergoes absorption and scattering, and is received by the same optical fiber probe is defined as a self-detection channel. The channel spacing ρ is the spacing between the transmitting optical fiber and the receiving optical fiber;

2) The calculation formula for the channel density of the rectangular network and the hexagonal honeycomb network is:

Where: is the number of channels of the ith basic unit, S _i ( j ) is the surface area of the ith basic unit covering the jth brain region, and when S _i ( j ) is zero, let S _i ( j ) be 1;

A number of optical fiber probes are respectively set in the rectangular network and the hexagonal honeycomb network. If the channel density of the rectangular network is greater than or equal to the channel density of the hexagonal honeycomb network, the rectangular network is used as the detection network, otherwise the hexagonal honeycomb network is used as the detection network;

3) All fiber optic probe arrangement points are numbered, the time resolution of the single-photon avalanche diode coupled to the fiber optic probe is set to T, and R fiber optic probes are deployed at probe arrangement points with different numbers; if a rectangular network is used, the brain area below each basic unit is taken as a rectangular body area, the length, width and height of the rectangular body area are not greater than three times the side length of the basic unit and not less than twice the side length of the basic unit, and it is divided into n local areas, and the local areas are divided into n-1 measurement areas and other areas, the thickness of each measurement area is not less than 1mm and not more than 1cm, and the remaining area excluding the measurement area is the other area; if a hexagonal honeycomb network is used, the brain area below each basic unit is taken as a regular triangular prism area, the edge length of the regular triangular prism area is not greater than three times the side length of the basic unit and not less than twice the side length of the basic unit, and it is divided into n local areas, and the local areas are divided into n-1 measurement areas and other areas, the thickness of each measurement area is not less than 1mm and not more than 1cm, and the remaining area excluding the measurement area is the other area;

4) Seven parameters are set for each local area, including length, width, thickness, refractive index, absorption coefficient, scattering coefficient and anisotropy factor. On this basis, a voxel-based time-resolved Monte Carlo simulation method is used to simulate 10 ⁸ photons. Sobel sampling is used to calculate the time-dependent average partial path length L _i,j (ρ,λ,t) of the jth local area with different channel spacing ρ, wavelength λ and detection time t for different basic units i;

5) The signal inspection strategy is as follows: select one of the fiber optic probes as the core probe, drive the laser diode coupled to the core probe to emit near-infrared pulse light with a light intensity of S and a main wavelength of λ, and at the same time, let the fiber optic probes in all basic units including the core probe receive the brain area diffuse light signal. The time for receiving the signal is set to _Tr , where _Tr is an integer multiple of T. After the signal is received, put all the fiber optic probe numbers that receive the brain area diffuse light signal into the waiting queue in clockwise order, take out the head of the queue, and let the fiber optic probe corresponding to the head number be the core probe. Repeat the above operation until there is no fiber optic probe number in the queue. If there is a deployed fiber optic probe that does not emit pulse light as a core probe, let it be the core probe and repeat the above operation until all fiber optic probes have emitted pulse light as core probes; calculate the time T _λ = _Tr xR for measuring a brain area diffuse light signal once, set the main wavelength of the laser diode pulse light to λ ₁ , execute the signal inspection strategy once, obtain the first measured brain area diffuse light signal with wavelength λ ₁ , and then adjust the main wavelength of the laser diode pulse light to λ ₂ , execute the signal inspection strategy once, and obtain the first measured brain area diffuse light signal of wavelength λ ₂ ;

6) After an interval of ΔT, the signal inspection strategy is repeatedly executed once at wavelength λ ₁ to obtain the second measured brain diffuse light signal of wavelength λ ₁ , and the brain diffuse light signals of all channels of wavelength λ ₁ are converted into time-resolved reflectance signals by using software;

7) Execute the basic unit signal joint analysis strategy for each basic unit to obtain the absorption coefficient change value of the local area of each basic unit under the wavelength λ ₁ ;

8) Repeat the signal inspection strategy once at wavelength λ ₂ to obtain the second measured brain diffuse light signal at wavelength λ ₂ , use software to convert the brain diffuse light signals of all channels at wavelength λ ₂ into time-resolved reflectance signals, perform the basic unit signal joint analysis strategy on each basic unit, and obtain the absorption coefficient change value of the local area of each basic unit at wavelength λ ₂ ;

9) The relationship between the absorption coefficient change Δu _a,j (λ) in the local area and the blood oxygen concentration is expressed as:

Δu _a,j (λ)＝ε _HbO (λ)Δ[HbO] _j +ε _HbR (λ)Δ[HbR] _j (1)

Wherein: ε _HbO (λ) is the molar extinction coefficient of oxygenated hemoglobin at wavelength λ, ε _HbR (λ) is the molar extinction coefficient of deoxygenated hemoglobin at wavelength λ, Δ[HbO] _j is the change in oxygenated hemoglobin concentration in the jth local area, and Δ[HbR] _j is the change in oxygenated hemoglobin concentration in the jth local area;

The change in absorption coefficient obtained in step 7) and step 8) is expressed by equation (1), and the equations for the same wavelength and local area are combined:

The simultaneous equations are solved to obtain the changes in the concentrations of oxyhemoglobin and deoxyhemoglobin in all local areas of each basic unit.

10) constructing a brain three-dimensional image file based on the AAL template or the Brodmann partition template, setting the brain surface, brain region position, and brain region edge, adjusting the three-dimensional image output layout size, background color, surface transparency, brain region color, brain region edge color, and image resolution, and then further dividing the corresponding brain region in the brain three-dimensional image according to the division of the local area of the detected brain region in step 3), so that each local area displays the normalized RGB value of the absorption coefficient change, and finally obtaining a cerebral hemooxygen dynamics response three-dimensional RGB image, according to the cerebral hemooxygen dynamics response three-dimensional RGB image, dividing the cerebral hemooxygen dynamics response depth map layer by layer, and obtaining a cerebral hemooxygen dynamics response 3D atlas.

2. A method for constructing a high-resolution 3D atlas of cerebral hemooxygen dynamics response according to claim 1, characterized in that ΔT should be greater than or equal to the time for solving the absorption coefficient change value minus the time for measuring the diffuse light signal of the brain area at a wavelength of λ ₁ The diffuse light signal of the brain area is measured once at a wavelength of λ ₂ The sum.

3. According to the method for constructing a high-resolution cerebral blood oxygen dynamic response 3D atlas according to claim 1, it is characterized in that the optical fiber cerebral blood oxygen detection system is composed of a driving circuit, a light source selection circuit, and a transmitting part composed of a pulsed laser 1 and a pulsed laser 2, an optical fiber coupler, an optical fiber splitter, an MPO optical fiber jumper, an optical fiber probe, an optical fiber combiner, a single photon counter and a driving circuit 2, and a receiving part composed of a time-correlated single photon counter (TCSPC) and a SPAD, and a PC signal processing end, wherein the optical fiber probe is composed of a housing, a dual-fiber collimator, a reflector and a focusing lens; the process of the system performing zero-source-detector distance time-resolved local area detection through a self-detection channel is as follows: the driving circuit transmits a control signal to the light source selection circuit, and simultaneously transmits a timing signal to the TCSPC The light source selection circuit makes the pulse laser 1 or the pulse laser 2 emit light of a specified wavelength. The optical signal passes through the optical fiber coupler and the optical fiber splitter, and then reaches the transmitting optical fiber connected to the optical fiber probe via the MPO optical fiber jumper. Then, after passing through the collimating lens and the reflecting surface of the reflector, the beam angle is adjusted by the focusing lens and injected into the cerebral cortex. Part of the light signal scattered by the brain tissue is absorbed, passes through the same focusing lens, reflecting surface and collimating lens, and returns to the MPO optical fiber jumper via the receiving optical fiber connected to the optical fiber probe, and then reaches the single photon counter through the optical fiber combiner. After the driving circuit 2 detects the signal of the single photon counter, it sends a stop signal to the TCSPC, and then the TCSPC sends a single photon signal to the SPAD. After the signal is converted and amplified by the SPAD, it is sent to the PC signal processing end for detection data processing.

4. The method for constructing a high-resolution 3D atlas of cerebral hemodynamic response according to claim 1, characterized in that the basic signal joint analysis strategy specifically comprises the following steps:

The following operations are performed on the channel signals in the basic unit: if the basic unit deploys four fiber optic probes, the time-resolved reflectivity signals of 12 mutual detection channels and 4 self-detection channels between the fiber optic probes are extracted, and each time-resolved reflectivity signal is divided into m time windows, 16m≥n; if the basic unit deploys three fiber optic probes, the time-resolved reflectivity signals of 6 mutual detection channels and 3 self-detection channels between the fiber optic probes are extracted, and each time-resolved reflectivity signal is divided into m time windows, 9m≥n; if the basic unit deploys two fiber optic probes, the time-resolved reflectivity signals of 2 mutual detection channels and 2 self-detection channels between the fiber optic probes are extracted, and each time-resolved reflectivity signal is divided into m time windows, 4m≥n; if the basic unit deploys one fiber optic probe, the time-resolved reflectivity signal of 1 self-detection channel of the fiber optic probe is extracted, and the sliding window method is used to divide the time-resolved reflectivity signal into m time windows, with a sliding step size of Δw, and m≥n; the channel spacing ρ, wavelength λ and the time-dependent average partial path length L of the jth local area at the detection time t of the basic unit are recalculated using the formula _i,j (ρ,λ,t):

Where: R( _Kz , _Kl , _Nl , _γz , _γl , ρ, t) is the time-resolved reflectivity multilayer approximate diffusion equation of the channel, _Kz , _Kl , _γz , _γl , _Nl are the solutions of the system of equations obtained by solving the approximate diffusion equation using the eigenfunction method;

The time-resolved reflectivity signals of the first and second measurements of each channel are subjected to curve inversion and the absorption change equations are solved, and the absorption change values of the obtained linear equations of all channels are solved.

5. The method for constructing a high-resolution cerebral hemodynamic response 3D atlas according to claim 2, characterized in that the curve inversion and absorption change equation group solution specifically include the following steps:

The system response function s(t) is measured, an initial time-resolved reflectivity signal R ₀ (ρ,λ,t) is taken, and the optimization curve y(t) is obtained by convolution, that is, Calculate the mean squared fitting error:

∑(R _m (ρ,λ,t)-y(t)) ² (4)where: R _m (ρ,λ,t) is the measured time-resolved reflectivity signal;

The sum of squares of the mean square fitting error is minimized by the least squares method, and the number of iterations is set to hundreds of thousands, and the threshold is T. When the sum of squares of the mean square fitting error is less than the threshold T, the fitting time-resolved reflectivity signal R(ρ,λ,t) is obtained for the subsequent theoretical calculation of the absorption coefficient change.

The width of each time window of the fitted time-resolved reflectivity signal R(ρ,λ,t) is set to Δt, and the first and last time windows are discarded. In order to describe the small independent absorption coefficient changes in the local area, the modified Lambert-Beer law is used to convert the number of detected photons in each time window into:

Where: _Rg is the number of photons detected in the g-th time window of the fitted time-resolved reflectivity signal of the second measurement, R0 _, g is the number of photons detected in the g-th time window of the fitted time-resolved reflectivity signal of the first measurement, Δu _a,j (λ) is the change in the absorption coefficient of the j-th local area in the second measurement relative to the first measurement, Lg _,i,j (ρ,λ) is the time-dependent average partial path length of the j-th local area in the g-th time window of the channel spacing ρ of the basic unit i and the wavelength λ, and the calculation formula is:

Where: Li _,j (ρ,λ,t) is the channel spacing ρ, wavelength λ of basic unit i and the time-dependent average partial path length of the jth local area at detection time t;

Transform formula (5) into:

The linear equation system consisting of equation (7) corresponding to all time windows is obtained:

Where Lg _,i,j (ρ,λ) is the time-dependent average partial path length of the jth local area in the gth time window with the channel spacing ρ and wavelength λ of the basic unit i.

6. The method for constructing a high-resolution cerebral hemooxygen dynamics response 3D atlas according to claim 2, characterized in that the absorption change value is solved, specifically comprising the following steps:

The linear equations of all channels obtained by solving the curve inversion and absorption variation equations are combined into one linear equation:

Where: is the time-dependent average partial path length of the jth local area in the gth time window of the cth channel of the basic unit i, R _c,g is the number of photons detected in the gth time window of the cth channel in the second measurement, R _c,0,g is the number of photons detected in the gth time window of the cth channel in the first measurement, ρ _c is the spacing of the cth channel, and λ is the main wavelength of the near-infrared pulse light;

The coefficient matrix of formula (9) is obtained:

If the rank of the coefficient matrix is equal to the rank of the augmented matrix, and both are equal to the number of local regions n, then the value of the absorption coefficient change Δu _a,j (λ) is solved, otherwise Δu _a,j (λ) is predicted based on the Swin-Transformer prediction model according to the previously obtained absorption coefficient change value.