CN115607108B - Multi-target reconstruction method based on surface measurement signal blind source separation - Google Patents

Abstract

A multi-target reconstruction method based on blind source separation of surface measurement signals. Step 1: Based on the light transmission model and finite element theory, the structural information and optical characteristic parameters of the reconstructed fluorescent target are used as prior information, and the physical process of the multi-target measurement signal generation is systematically analyzed to show that it meets the conditions for blind source separation. Step 2: Through multi-excitation point excitation, a surface superposition measurement signal that meets the blind source separation conditions is constructed. Step 3: The separation matrix is constructed using the blind source separation method to separate the superimposed measurement signal into the source measurement signal of each target on the surface of the biological tissue, and a linear relationship between each source measurement signal and the distribution of the light source in the body is established. Step 4: Each target is reconstructed separately using the reconstruction algorithm, and the reconstruction results of each target are combined and displayed in a coordinate system. This method provides accurate positioning and morphological restoration capabilities while ensuring the robustness of the reconstruction problem to parameters and algorithms. This is a breakthrough work in multi-source reconstruction using optical molecular technology in pre-clinical research on cancer staging.

Description

A multi-target reconstruction method based on blind source separation of surface measurement signals

Technical Field

The invention belongs to the field of fluorescent molecular tomography, and in particular relates to a multi-target reconstruction method based on blind source separation of surface measurement signals.

Background Art

Fluorescence molecular tomography (hereinafter referred to as FMT) is a new type of optical molecular imaging technology developed in recent years. It has become an important branch of molecular imaging technology and is widely used in pre-clinical research due to its many advantages such as high sensitivity, high temporal resolution, non-ionization, non-radioactivity, low cost, and fast and convenient imaging. This modality uses an external light source to excite the fluorescent probe to emit photons, and then uses a fluorescence acquisition device (a high-sensitivity CCD camera) to collect the fluorescence signal on the surface of the biological body. Combined with mathematical models, the three-dimensional distribution and fluorescence yield of the fluorescence source inside the target are reconstructed, thereby achieving qualitative and quantitative research on biological processes at the cellular and molecular levels in vivo. Considering the absorption and scattering of light by the target itself, accurately and quickly reconstructing the three-dimensional distribution of fluorescent markers inside the target has always been a key issue in FMT. In the past decade, people have made a lot of efforts in reconstruction algorithms and imaging systems. The l ₂ -norm regularization was proposed to solve the inverse problem and combined with some iterative techniques for image reconstruction. Inspired by the theory of compressed sensing, the sparse regularization method with l ₁ -norm penalty function has also attracted widespread attention due to its excellent performance in FMT. Studies have found that the non-convex l _p (0<p<1)-regularization method can produce sparser and better solutions than the l ₁ -norm method in sparse reconstruction. In addition, prior information such as the optical properties and feasible domain of biological tissues is also incorporated into the reconstruction.

In actual tumor research, due to the spread and metastasis of tumor cells, it is necessary to conduct a comprehensive analysis of multiple tumor areas for staging research. Therefore, it is necessary to obtain high-precision FMT multi-light source reconstruction results for quantitative analysis. However, there are many coping methods that perform well in reconstructing a single fluorescent target, but the serious ill-posedness of the FMT inverse problem leads to large differences in the reconstructed concentrations of each target during multi-target reconstruction, making the processing of reconstruction results and the selection of feasible domains very difficult. For FMT multi-target reconstruction, multi-point excitation and multi-angle projection can alleviate the ill-posedness of the problem to a certain extent by increasing the measurement data, but this will bring new problems, such as increasing data acquisition time and redundant data in the measurement data, which is not conducive to the accurate and efficient reconstruction of fluorescence molecular tomography, which is a major obstacle in the practical application of FMT.

Summary of the invention

Aiming at the problem of three-dimensional detection of multiple tumors in clinical practice, the present invention proposes a multi-target reconstruction method based on blind source separation of surface measurement signals, which solves the problem of difficulty in reconstruction of multiple light sources.

In order to achieve the above object, the technical solution adopted by the present invention is as follows:

A multi-target reconstruction method based on blind source separation of surface measurement signals comprises the following steps:

Step 1: Based on the light transmission model and finite element theory, the structural information and optical characteristic parameters of the reconstructed fluorescent target are used as prior information, and the physical process of the multi-target measurement signal generation is systematically analyzed to show that it meets the conditions for blind source separation.

Step 2: Construct a surface superposition measurement signal that meets the blind source separation conditions through multi-excitation point excitation.

Step 3: Use the blind source separation method to construct a separation matrix, separate the superimposed measurement signals into source measurement signals of each target on the surface of the biological tissue, and establish a linear relationship between each source measurement signal and the distribution of the light source in the body.

Step 4: Use the reconstruction algorithm to reconstruct each target separately, and combine the reconstruction results of each target in a coordinate system for display.

Further, the step 1 specifically includes:

In fluorescence molecular tomography, the transmission of light includes two related processes, namely the excitation process and the emission process. The excitation process refers to irradiating certain specific points or specific areas of the imaging object with external excitation light. This part of light enters the interior through the surface of the imaging object and forms a light intensity distribution inside the object. In the present invention, we consider the steady-state imaging mode, that is, the continuous wave imaging mode. The emission process is that the fluorophore inside the object absorbs a certain amount of excitation light energy, and a part of it is converted into photons and released. The released part of light is called emission light, and its wavelength is larger than the wavelength of the excitation light, and its photon energy is less than the excitation light photon. The transmission process of excitation light and emission light can be described by the following two coupled diffusion equations:

Where Ω represents the three-dimensional space occupied by the imaging object, the subscripts x and m represent the excitation light and the emission light respectively, and η is the quantum yield. is the absorption coefficient of the fluorophore to the excitation light, which is proportional to the concentration of the fluorophore. This is the required light yield distribution. is the Dirac function, Θ represents the intensity of the light source.

The solution is based on the Galerkin variational method and combined with the Robin boundary condition:

in, is the boundary of Ω, v(r) is the normal vector of r, and q is a constant. Assume that the imaging area is discretized into _Np points and _Ne tetrahedral units, and the basis function on the discretized grid is used as the test function. Since any function can be expressed as a linear combination of basis functions, the basis function can be used as the test function. So and It can be expressed as

Among them, φ _xk,mk and (ημ _af ) _k represent the function value at node k, represents the basis function at node k.

For each unit, substitute equations (4) and (5) into the diffusion equation to form the stiffness matrix of each unit. Superimpose and assemble all unit stiffness matrices to obtain the following equation:

K _x Φ _x ＝T _x (6)

K _m Φ _m ＝Fx (7)

in, make The distribution vector Φ _x of the excitation light in the imaging object can be obtained by directly solving equation (6), and then substituting it into (7), the distribution vector Φ _m of the emission light inside and on the surface of the imaging object is obtained:

Φ _m ＝K _m ^-1 Fx＝Ax＝AΦ _x S (8)

Since any signal can be represented by a combination of impulse signals, we express the above formula as:

Where AΦ _x is an M×N vector, representing the system matrix; S is an N×1 vector, representing the distribution of internal fluorescent probes; Φ is an M×1 vector, representing the measured value of the photon distribution on the surface of the organism.

For FMT multi-target reconstruction, when one excitation point excites multiple light sources, the emitted light is the superposition of the light emitted by multiple light sources. The δ(r-ri) of each light source in the imaging object is different, so the weights of the fluorescence yield distribution of different light sources are different, ω _i =Φ _xi δ(r-ri). The distribution vector of the emitted light of each light source inside and on the surface of the imaging object can be expressed as:

Φ _mi ＝ω _i AS _i (10)

Assuming that the number of light sources and the number of excitation lights are p, the surface light distribution of p aliased signals detected by the detector is:

Where ω＝[ω _ij ] _p×p ∈R ^p×p is the mixing matrix. If the surface light distribution information can be separated, the unmixed Then, the individual reconstruction of each light source can be achieved. This step is the core of the multi-objective reconstruction strategy.

The only available information in the blind source separation algorithm (hereinafter referred to as BSS) comes from the mixed observation matrix. The BSS model can be expressed as follows:

Φ[n]＝Ws[n] (n＝1,2...,M) (12)

Wherein, Φ[n]＝(Φ ₁ [n],...,Φ _p [n]) ^T is a p×M aliasing measurement matrix, W＝[ω _ij ] _p×k ∈R ^p×k is an unknown mixing matrix, and s[n]＝(s ₁ [n],...,s _k [n]) ^T refers to the separated k×M source matrix. It is usually further assumed that the dimensions of Φ and s are equal, that is, p＝k, and we make this assumption in the rest of the present invention. Each row in Φ[n]＝(Φ ₁ [n],...,Φ _p [n]) ^T contains p source information s ₁ ,s ₂ …s _p , which is consistent with the aliasing of FMT multi-target signals. Therefore, the physical process of generating multi-target measurement signals satisfies the conditions for blind source separation.

Further, the step 2 specifically includes:

In FMT, the superposition of light emitted by multiple light sources can be used as the mixed observation value of BSS, that is, the surface light distribution information _Φm is used as the mixed observation value Φ. Where Φ is a p×M matrix, we use it as the input data matrix of the blind source separation algorithm to satisfy the conditions of blind source separation:

Further, the step three specifically includes:

For mixed observations, our goal is to obtain the surface light distribution information of each light source through the separation algorithm, and then reconstruct each light source separately. The assumption of p = k here means that the number of separated light sources should be consistent with the number of excitation points. The surface light distribution of p independent signals after blind source separation is:

In the specific separation process, we design the separation matrix B = {b _ij } _p×p ∈R ^p×p to satisfy

z[n]＝BΦ[n]＝BWs[n]＝Ps[n] (15)

Among them, z[n]＝(z ₁ [n], ..., z _p [n]) ^T represents the separated or extracted source signal, and P is a permutation matrix, which means that the extracted source point z[n] can be equivalent to the actual source signal s[n].

The algorithm used in the present invention decomposes the FMT data by minimizing the joint correlation function between the estimated non-negative sources. The separation matrix B can be obtained by maximizing the volume of the simplex cov{0,z ₁ ,...z _p } (after separation), which is related to the volume of cov{0,Φ ₁ ,...,Φ _p } (before separation):

vol(cov{0,z ₁ ,...z _p })＝|det(B)|vol(cov{0,Φ ₁ ,...,Φ _p }) (16)

It can be seen from relevant literature that the optimal separation matrix can be obtained by maximizing the volume of the solid region formed by the demixing source vector. The matrix is full rank and the sum of the elements in each row of the matrix is 1. Therefore, solving the volume maximization problem can obtain the optimal solution:

{z _i ,i＝1,...,p}={s _i ,i＝1,...,p} (17)

Therefore, the blind source separation problem can be transformed into the following non-convex optimization problem:

(1) p = k = 2

In this case, we can use the analytical method to find the closed-form solution and integrate the equality constraint in formula (18) into the objective function to obtain:

Assuming det(B)≥0, formula (18) is transformed into

max(B ₁₁ -B ₂₁ ) (20)

stB ₁₁ Φ ₁ [n]+(1-B ₁₁ Φ ₂ [n])≥0

B ₂₁ Φ ₁ [n]+(1-B ₂₁ Φ ₂ [n])≥0

The above constraints can be equivalent to: β≤B ₁₁ ,B ₂₁ ≤α

in:

Therefore, formula (20) is finally transformed into the following form:

The optimal solution is B ₁₁ ^* = α and B ₂₁ ^* = β, so the optimal value is det(B ^* ) = α - β > 0.

Similarly, it can be proved that when det(B)≤0, the optimal value is det(B ^* )＝β-α<0.

(2) p = k > 2

Expand any row of B by algebraic co-factors. Take the i-th row as an example, let B _i ^T = [B _i1 ,B _i2 ,...,B _iN ], then the determinant is:

Where _Bij is the submatrix with the i-th row and j-th column of B deleted. Obviously, when _Bij is fixed with respect to j=1,2,...,N, the determinant det(B) becomes a linear function with respect to _Bi . Therefore, by fixing the other row vectors of B, formula (18) is transformed into the following non-convex maximization problem:

The objective function in formula (23) is still non-convex. By solving the following two linear programming problems, the global optimal solution of the local maximization problem can be obtained:

During the algorithm implementation, we update the matrix B in each iteration until convergence, thereby selecting the optimal solution.

After separation, Φ in the linear model is redefined as z, so the linear relationship between the surface measurement value and the internal fluorescent target distribution after separation is expressed by the following formula:

AΦ _x S＝Ax＝z (26)

Among them, A is the system weight matrix, z is the measured value of the photon distribution of each light source on the surface of the organism after blind source separation, and S is the three-dimensional distribution of the reconstructed target.

Since the surface measurement signal of each light source after blind source separation is still the result of the combined action of multiple excitation lights, for the jth light source, the FMT model can be expressed as follows:

So far, we have successfully transformed the multi-light source reconstruction problem into multiple single-light source reconstruction problems.

Further, the step 4 specifically includes:

In FMT, target reconstruction is equivalent to finding an approximate solution to the equation Ax = z, so we use The matrix formed by the mixture and the matrix z obtained by BSS separation are used to reconstruct each light source separately. In order to obtain a sparser solution, inspired by the theory of compressed sensing, the FMT model can be transformed into the following minimization problem with an l ₁ regularization term:

Using the reconstruction algorithm to reconstruct each light source separately, we get:

S _j =(S _j [1],...,S _j [N]) ^T j＝1,2...,p (29)

Since blind source separation establishes a general multi-light source reconstruction model and has no algorithm dependency, the reconstruction algorithm is not the focus of the present invention. In the following embodiments, we use the classic fast iterative shrinkage threshold algorithm (hereinafter referred to as FISTA) as a specific reconstruction algorithm. FISTA is a very popular gradient-based algorithm, which is improved by the gradient descent method. For the classic LASSO problem, the FISTA algorithm based on the soft threshold function can quickly iterate and solve it. Time complexity, using different approximate function starting points can greatly improve the convergence speed.

Finally, we combine the reconstruction results of each target in one coordinate system for display.

Compared with the prior art, the present invention has the following advantages:

First, from a clinical perspective, the present invention specifically addresses the problem of comprehensive analysis and three-dimensional detection of multiple tumor regions due to the spread and metastasis of tumor cells in actual tumor research.

Second, the present invention introduces the idea of blind source separation based on surface measurement signals into multi-light source reconstruction, transforming the multi-target reconstruction problem into multiple single-target reconstruction problems, which to a certain extent solves the problem of low multi-target reconstruction accuracy caused by the serious ill-posedness of the FMT inverse problem.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG1 is a flow chart of a multi-target reconstruction method based on blind source separation of surface measurement signals provided by an embodiment of the present invention;

Figure 2 shows the USC Digital Mouse simulation model, with all organs and tissues labeled;

FIG3 is a .grid.am file of a digital mouse simulation model without light source segmented by Amira software;

FIG4 is a forward simulation result of the surface measurement signal of dual-target aliasing displayed by tecplot software. The size and shape of the light source is a point light source with a radius of 1.2 mm, and the positions are (23.5, 8, 16.5) mm and (16.1, 8, 16.5) mm respectively; (a) is the surface spot distribution diagram obtained by excitation point 1, and (b) is the surface spot distribution diagram obtained by excitation point 2;

FIG5 is a forward simulation result of the surface measurement signal of the above dual targets after blind source separation displayed by tecplot software, wherein (a) is the surface spot distribution diagram of target 1, and (b) is the surface spot distribution diagram of target 2;

FIG6 is an experimental result of reconstruction using the FISTA algorithm after blind source separation, where (a) is a three-dimensional display of the result, and the real light source position is marked with an arrow; (b) is a cross-sectional view of (a) along the XOY plane at Z=16.5 mm, and the circle represents the real light source position, which is marked with an arrow;

Specific implementation examples

The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. It should be noted that the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention.

The reconstruction method of the present invention is based on the fluorescence molecular tomography technology, which is based on the physical model of light propagation in biological tissues. An external light source is used to excite the fluorescent probe to emit photons, and then a fluorescence acquisition device (a high-sensitivity CCD camera) is used to collect the fluorescence signal on the biological surface. Then, the spatial position of the internal light source target and the signal intensity distribution are reversely calculated by using a reconstruction algorithm in combination with the previously known anatomical structure inside the biological tissue and the corresponding optical parameters of each tissue (such as absorption coefficient, scattering coefficient, etc.), thereby quantitatively displaying the changes in the early lesions inside the biological tissue. Specifically, it involves a multi-target reconstruction method based on blind source separation of surface measurement signals in the field of fluorescence molecular tomography technology. The method can be used to accurately find the locations of multiple tumors in the lesion area, solving the problem of poor multi-target reconstruction results due to the serious ill-posedness of the FMT inverse problem.

The application principle of the present invention is described in detail below in conjunction with the accompanying drawings.

Referring to FIG. 1 , the present invention takes the case where the number of targets is 2 as an example, and the process of the reconstruction method is as follows:

S101: Acquire measurement data of a limited number of angles;

S102: Obtain anatomical structure information and optical characteristic parameters of the reconstructed target, and obtain a forward result of a surface measurement signal of dual-target superposition;

S103: using a blind source separation algorithm, respectively obtaining forward results of source measurement signals of dual targets on the surface of biological tissue after blind source separation;

S104: establishing a linear relationship between each source measurement signal and the distribution of the light source in the body, and reconstructing each target separately using the FISTA algorithm to obtain a separated reconstruction result;

The implementation process of step S101 is as follows:

The FMT imaging system (laser and cooled CCD) and the XCT imaging system (X-ray tube and X-ray detector) are placed in a straight line with the center of the mouse, and the two systems are perpendicular to each other. The excitation light generated by the laser illuminates the surface of the mouse, stimulating the fluorescent targets in the mouse's body to generate emitted fluorescence and propagate to the surface of the mouse, and then these fluorescent signals are collected by the cooled CCD camera. At the same time, the XCT system obtains the structural image of the mouse.

Step S102 specifically includes:

(1) Perform three-dimensional reconstruction of the computed tomography data of the imaging target using 3Dmed software;

(2) Using Amira to segment the volume data, we can obtain the anatomical structure information of the imaging target;

(3) Obtaining optical characteristic parameters of the imaging target using the light flux density and the anatomical structure information of the imaging target;

As mentioned above, the distribution vector _Φm of the emitted light inside and on the surface of the imaging object can be expressed as:

Φ _m = AΦ _x S

For FMT multi-target reconstruction, when one excitation point excites multiple light sources, the emitted light is the superposition of the light emitted by multiple light sources. The δ(r-ri) of each light source in the imaging object is different, so the weights of the fluorescence yield distribution of different light sources are different, ω _i =Φ _xi δ(r-ri). The distribution vector of the emitted light of each light source inside and on the surface of the imaging object can be expressed as:

Φ _mi ＝ω _i AS _i

The BSS model can be expressed as follows:

Φ[n]＝Ws[n] (n＝1,2...,M)

Among them, Φ[n]＝(Φ ₁ [n],...,Φ _p [n]) ^T is a p×M aliasing measurement matrix, W＝[ω _ij ] _p×k ∈R ^p×p is an unknown mixing matrix, and s[n]＝(s ₁ [n],...,s _p [n]) ^T refers to the separated p×M source matrix.

In FMT, the superposition of light emitted by multiple light sources can be used as the mixed observation value of BSS, that is, the surface light distribution information _Φm is used as the mixed observation value Φ. Where Φ is a p×M matrix, we use it as the input data matrix of the blind source separation algorithm to satisfy the conditions of blind source separation:

The specific implementation of step S103 is:

The surface light distribution of p independent signals after blind source separation is:

In the specific separation process, we design the separation matrix B = {b _ij } _p×p ∈R ^p×p to satisfy

z[n]＝BΦ[n]＝BWs[n]＝Ps[n]

Among them, z[n]＝(z ₁ [n], ..., z _p [n]) ^T represents the separated or extracted source signal, and P is a permutation matrix, which means that the extracted source point z[n] can be equivalent to the actual source signal s[n].

The algorithm used in the present invention decomposes the FMT data by minimizing the estimated joint correlation function between non-negative sources. The separation matrix B can be obtained by maximizing the volume of the simplex cov{0,z ₁ ,...z _p } (after separation), which is related to the volume of cov{0,Φ ₁ ,...,Φ _p } (before separation), that is,

vol(cov{0,z ₁ ,...z _p })＝|det(B)|vol(cov{0,Φ ₁ ,...,Φ _p })

Solving the volume maximization problem gives us the optimal solution:

{z _i ,i＝1,...,p}={s _i ,i＝1,...,p}

After separation, Φ in the linear model is redefined as z, so the linear relationship between the surface measurement value and the internal fluorescent target distribution after separation is expressed by the following formula:

AΦ _x S＝Ax＝z

Among them, A is the system weight matrix, z is the measured value of the photon distribution of each light source on the surface of the organism after blind source separation, and S is the three-dimensional distribution of the reconstructed target.

Since the surface measurement signal of each light source after blind source separation is still the result of the combined action of multiple excitation lights, for the jth light source, the FMT model can be expressed as follows:

So far, we have successfully obtained the light distribution information of each independent light source on the surface of the organism.

The step S104 is:

In FMT, target reconstruction is equivalent to finding an approximate solution to the equation Ax = z, so we use The matrix formed by the mixture and the matrix z obtained by BSS separation are used to reconstruct each light source separately. In order to obtain a sparser solution, inspired by the theory of compressed sensing, the FMT model can be transformed into the following minimization problem with an l ₁ regularization term:

Using the inverse algorithm to reconstruct each light source separately, we get:

S _j =(S _j [1],...,S _j [N]) ^T j＝1,2...,p

We use the classic FISTA algorithm as the specific reconstruction algorithm to reconstruct the two separated targets separately, and finally combine the reconstruction results of each target in a coordinate system for display.

The following is a detailed description of the application effect of the present invention in combination with simulation. The steps are as follows:

1. Get the .raw file containing the light source.

Figure 2 is a digital mouse model used in simulation experiments. The digital mouse model is obtained by tomographic reconstruction using tissue information extracted from imaging slice data obtained from real experimental animals through CT technology. The model is a .raw file with a size of 38.0mm×20.8mm×35.0mm, and does not contain a light source. In order to reduce computational complexity and save system resources, the head and tail regions of the mouse are removed, leaving only the trunk. According to its internal tissue structure, it can be divided into six organs and the remaining tissues of the body, from top to bottom, namely, heart, lungs, liver, stomach, kidneys and muscles.

2. Obtain the .grid.am file of the entire simulation model and light source.

FIG3 is a grid structure of the entire simulation model obtained by segmenting the .raw file obtained in FIG1 using Amira software. The file format is .grid.am file.

3. Obtain simulation results of dual-target aliased surface measurement signals.

Figure 4 shows the simulation results of the light distribution information of the dual light sources on the surface of biological tissue. The probe inside the biological tissue emits photons by stimulated radiation. The photons undergo reflection, absorption, scattering and other physical optical processes during the transmission from the inside to the outside, and are transmitted to the surface of the organism to form light distribution information. There is a one-to-one nonlinear mapping relationship between the light distribution information on the biological surface and the light source information. The .grid.am file of the entire simulation model obtained in II is used to simulate the light distribution information on the surface of biological tissue.

4. Use the aliasing surface to measure the signal and generate pre-run data.

The light distribution information on the biological tissue surface is extracted from the forward simulation result file obtained in step III using the Matlab programming language, and the .grid.am grid file generated in step II is read in accordingly. Based on the grid file and the simulation result file, the surface energy distribution vector, the grid node coordinate (X, Y and Z) matrix, the grid internal tetrahedron matrix, the grid internal tetrahedron index matrix and the system matrix are generated, and these results are saved as the corresponding .mat file of Matlab.

5. Perform blind source separation on the .mat file obtained in step IV to obtain a new .mat file after separation.

For FMT multi-target reconstruction, when one excitation point excites multiple light sources, the emitted light is the superposition of the light emitted by multiple light sources. Select p excitation points with the same number of light sources for multi-angle excitation. If the surface light distribution information of the p aliased signals detected by the detector can be separated, that is, the surface light distribution of the unmixed source signal can be obtained, then each light source can be reconstructed separately. We use the BSS algorithm to perform blind source separation on the .mat file obtained in the previous step.

FIG5 shows the simulation results of the measurement signals of the two target surfaces after separation.

6. Use the separated .mat file obtained in step V as the input of the reconstruction algorithm for reconstruction.

The digital mouse model was discretized into a tetrahedral mesh consisting of 13064 nodes and 67649 tetrahedral elements.

FIG6 is the experimental result reconstructed using FISTA, where (a) is a three-dimensional result display diagram. The real light source is two spherical fluorescent sources with a radius of 1.2 mm, which are placed in the liver. The centers are (16.1, 8.0, 16.5) and (23.5, 8.0, 16.5), respectively. The distance between the two light sources is 5 mm, and the positions are marked with arrows; (b) is a cross-sectional view of (a) along the XOY plane at Z = 16.5 mm. The circle represents the position of the real light source, which is marked with an arrow;

To further quantitatively evaluate the ability of the BBS strategy in source localization and shape recovery, we use the location error (LE) and Dice coefficient as quantitative indicators.

LE measures the distance between the center of the true light source position and the center of the reconstructed area. LE is defined as:

LE＝||L _r -L ₀ || ₂

Where L _r is the center of the reconstructed area, L ₀ is the center of gravity of the real fluorescence area. ||.|| ₂ is the operator of Euclidean distance. The lower the LE index, the better the reconstruction effect.

In the embodiment of the present invention, the two light sources LE are 0.56 mm and 0.53 mm respectively, and their position errors are small, which shows that the strategy proposed in the present invention has good position restoration capability for multiple targets.

The Dice coefficient is used to evaluate the morphological similarity between the reconstructed area and the real fluorescence area, which is defined as:

Among them, S _r and S ₀ represent the area of the reconstructed region and the area of the actual fluorescence region, respectively. The higher the Dice coefficient, the higher the similarity between the two regions in position and morphology.

In the embodiment of the present invention, the Dice coefficient is 72%, and the similarity is relatively high, which indicates that the strategy proposed in the present invention has good morphological restoration capability for multiple targets.

Claims (2)

1. A multi-target reconstruction method based on blind source separation of surface measurement signals, characterized in that it comprises the following steps: Step 1: Based on the light transmission model and finite element theory, the structural information and optical characteristic parameters of the reconstructed fluorescent target are used as prior information to systematically analyze the physical process of the multi-target measurement signal generation to show that it meets the conditions for blind source separation; In fluorescence molecular tomography, the transmission of light includes two related processes, namely the excitation process and the emission process. Continuous wave imaging is used. In the emission process, the fluorophore inside the object absorbs a certain amount of excitation light energy, and part of it is converted into photons and released. The released part of light is called emission light, whose wavelength is larger than the wavelength of the excitation light, and its photon energy is less than the excitation light photon. The transmission process of excitation light and emission light can be described by the following two coupled diffusion equations: Where Ω represents the three-dimensional space occupied by the imaging object, the subscripts x and m represent the excitation light and the emission light respectively, and η is the quantum yield. is the absorption coefficient of the fluorophore to the excitation light, which is proportional to the concentration of the fluorophore. is the required distribution of light yield, is the Dirac function, Θ represents the intensity of the light source; The solution is based on the Galerkin variational method and combined with the Robin boundary condition: in, is the boundary of Ω, v(r) is the normal vector of r, and q is a constant; suppose the imaging area is discretized into _Np points and _Ne tetrahedral units, and the basis function on the discretized grid is used as the test function. Since any function can be expressed as a linear combination of basis functions, the basis function can be used as the test function, so and It can be expressed as Among them, Φ _xk,mk and (ημ _af ) _k represent the function value at node k, represents the basis function at node k; For each unit, substitute equations (4) and (5) into the diffusion equation to form the stiffness matrix of each unit. Superimpose and assemble all unit stiffness matrices to obtain the following equation: K _x Φ _x ＝T _x (6) K _m Φ _m ＝Fx (7) in, make The distribution vector Φ _x of the excitation light in the imaging object can be obtained by directly solving equation (6), and then substituting it into (7), the distribution vector Φ _m of the emission light inside and on the surface of the imaging object is obtained: Φ _m ＝K _m ^-1 Fx＝Ax＝AΦ _x S (8) Since any signal can be represented by a combination of impulse signals, the above formula can be expressed as: Where AΦ _x is an M×N vector, representing the system matrix; S is an N×1 vector, representing the distribution of internal fluorescent probes; Φ is an M×1 vector, representing the measured value of the photon distribution on the surface of the organism; For FMT multi-target reconstruction, when one excitation point excites multiple light sources, the emitted light is the superposition of the light emitted by multiple light sources. The δ(r-ri) of each light source in the imaging object is different, so the weights of the fluorescence yield distribution of different light sources are different, ω _i =Φ _xi δ(r-ri). The distribution vector of the emitted light of each light source inside and on the surface of the imaging object can be expressed as: Φ _mi ＝ω _i AS _i (10) Assuming that the number of light sources and the number of excitation lights are p, the surface light distribution of p aliased signals detected by the detector is: Where ω＝[ω _ij ] _p×p ∈R ^p×p is the mixing matrix. If the surface light distribution information can be separated, the unmixed It is then possible to achieve individual reconstruction of each light source, which is the core of the multi-objective reconstruction strategy. In the blind source separation algorithm, hereinafter referred to as BSS, the only available information comes from the mixed observation matrix. The BSS model can be expressed as follows: Φ[n]＝Ws[n] (n＝1,2...,M) (12) Among them, Φ[n]＝(Φ ₁ [n],...,Φ _p [n]) ^T is a p×M aliasing observation matrix, W＝[ω _ij ] _p×k ∈R ^p×k is an unknown mixing matrix, s[n]＝(s ₁ [n],...,s _k [n]) ^T refers to the separated k×M source matrix; it is usually further assumed that the dimensions of Φ and s are equal, that is, p＝k, and it is assumed that: each row in Φ[n]＝(Φ ₁ [n],...,Φ _p [n]) ^T contains p source information s ₁ ,s ₂ …s _p , which is consistent with the aliasing of FMT multi-target signals, so the physical process generated by the multi-target measurement signals meets the conditions for blind source separation; Step 2: Construct a surface superposition measurement signal that meets the blind source separation conditions through multi-excitation point excitation; In FMT, the superposition of light emitted by multiple light sources can be used as the mixed observation value of BSS, that is, the surface light distribution information _Φm is used as the mixed observation value Φ; where Φ is a p×M matrix, which is used as the input data matrix of the blind source separation algorithm to meet the conditions of blind source separation: Step 3: Use the blind source separation method to construct a separation matrix, separate the superimposed measurement signals into source measurement signals of each target on the surface of the biological tissue, and establish a linear relationship between each source measurement signal and the distribution of the light source in the body; For mixed observations, the surface light distribution information of each light source is obtained through the separation algorithm, and then each light source is reconstructed separately. The assumption of p = k here means that the number of separated light sources should be consistent with the number of excitation points. The surface light distribution of p independent signals after blind source separation is: In the specific separation process, the separation matrix B = {b _ij } _p×p ∈R ^p×p is designed to satisfy z[n]＝BΦ[n]＝BWs[n]＝Ps[n] (15) Wherein, z[n]＝(z ₁ [n],...,z _p [n]) ^T represents the separated or extracted source signal, and P is a permutation matrix, which means that the extracted source point z[n] can be equivalent to the actual source signal s[n]; The FMT data is decomposed by minimizing the estimated joint correlation function between non-negative sources. The separation matrix B can be obtained by maximizing the volume of the simplex cov{0,z ₁ ,...z _p }, after separation, which is related to the volume of cov{0,Φ ₁ ,...,Φ _p }, before separation: vol(cov{0,z ₁ ,...z _p })＝|det(B)|vol(cov{0,Φ ₁ ,...,Φ _p }) (16) The optimal separation matrix can be obtained by maximizing the volume of the solid region formed by the demixing source vector. The matrix has full rank and the sum of the elements in each row of the matrix is 1. Therefore, the optimal solution can be obtained by solving the volume maximization problem: {z _i ,i＝1,...,p}={s _i ,i＝1,...,p} (17) Therefore, the blind source separation problem can be transformed into the following non-convex optimization problem: During the algorithm implementation, the matrix B is updated in each iteration until convergence, thereby selecting the optimal solution; After separation, Φ in the linear model is redefined as z, so the linear relationship between the surface measurement value and the internal fluorescent target distribution after separation is expressed by the following formula: AΦ _x S＝Ax＝z (19)where A is the system weight matrix, z is the measured value of the photon distribution of each light source on the surface of the organism after blind source separation, and S is the three-dimensional distribution of the reconstructed target; Since the surface measurement signal of each light source after blind source separation is still the result of the combined action of multiple excitation lights, for the jth light source, the FMT model can be expressed as follows: So far, the multi-light source reconstruction problem is transformed into multiple single-light source reconstruction problems; Step 4: Use the reconstruction algorithm to reconstruct each target separately; in FMT, target reconstruction is equivalent to finding an approximate solution to the equation Ax = z, so use The matrix formed by the mixture and the matrix z obtained by BSS separation are used to reconstruct each light source separately. In order to obtain a sparser solution, the FMT model can be transformed into the following minimization problem with an l ₁ regularization term: Use the reconstruction algorithm to reconstruct each light source separately and get: S _j =(S _j [1],...,S _j [N]) ^T j＝1,2...,p (22) Since blind source separation establishes a universal multi-light source reconstruction model, it has no dependence on the reconstruction algorithm; therefore, in the following embodiments, the classic fast iterative shrinkage threshold algorithm, hereinafter referred to as FISTA, is used as a specific reconstruction algorithm. The FISTA algorithm based on the soft threshold function can be quickly iterated and solved. Time complexity, using different approximate function starting points to improve the convergence speed; finally, the reconstruction results of each target are combined and displayed in a coordinate system. 2. The multi-objective reconstruction method based on blind source separation of surface measurement signals according to claim 1, characterized in that the blind source separation problem in step 3 is transformed into the following non-convex optimization problem: (1) p = k = 2 In this case, we can use the analytical method to find the closed-form solution and integrate the equality constraint in formula (18) into the objective function to obtain: Assuming det(B)≥0, formula (18) is transformed into max(B ₁₁ -B ₂₁ )(25) stB ₁₁ Φ ₁ [n]+(1-B ₁₁ Φ ₂ [n])≥0 B ₂₁ Φ ₁ [n]+(1-B ₂₁ Φ ₂ [n])≥0 The above constraints can be equivalent to: β≤B ₁₁ ,B ₂₁ ≤α in: Therefore, formula (20) is finally transformed into the following form: The optimal solution is B ₁₁ ^* = α and B ₂₁ ^* = β, so the optimal value is det(B ^* ) = α-β>0; Similarly, it can be proved that when det(B)≤0, the optimal value is det(B ^* )＝β-α<0; (2) p = k > 2 Expand any row of B by algebraic co-factors. Take the i-th row as an example, let B _i ^T = [B _i1 ,B _i2 ,...,B _iN ], then the determinant is: Where _Bij is the submatrix of B with the i-th row and j-th column deleted. Obviously, when _Bij is fixed with respect to j=1,2,...,N, the determinant det(B) becomes a linear function with respect to _Bi . Therefore, by fixing the other row vectors of B, formula (20) is transformed into the following non-convex maximization problem: The objective function in formula (23) is still non-convex. By solving the following two linear programming problems, we can obtain the global optimal solution to the local maximization problem: