An alternating rank-r and rank-1 shift-invariant cpd method for spatially compressed multi-subject fMRI

Abstract

The application discloses a spatial compression multi-subject fMRI alternating rank R and rank 1 shift invariant CPD method, and belongs to the field of medical signal processing. Firstly, the spatial dimension of three-dimensional multi-subject fMRI data is subjected to average down-sampling processing, so that the relative position relationship between feature points is retained, and the spatial dimension is effectively reduced. Secondly, an alternating rank R and rank 1 least square method is adopted, so that the shift invariant CPD model is effectively relaxed. Finally, spatial orthogonal constraints are added to shared brain spatial activation components. The application can quickly and effectively extract shared brain function information of multi-subjects, and has a good development prospect in future large-scale high-dimensional fMRI data research and intelligent medical treatment and the like.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of medical signal processing, in particular to an alternating rank R and rank 1 shift-invariant CPD method for spatially compressed multi-subject fMRI. BACKGROUND

[0002] fMRI is a neuroimaging technique that measures the changes in blood dynamics (blood oxygen level dependent signal) caused by neuronal activity using magnetic resonance imaging. Due to its non-invasiveness, repeatability, and high resolution, it is widely used in brain research and clinical treatment. Multi-subject fMRI data can be considered as a five-dimensional tensor, which contains three-dimensional brain spatial information, one-dimensional time information, and one-dimensional subject intensity. Typically, the three-dimensional brain space is processed in one dimension, so that the multi-subject fMRI data becomes a three-dimensional tensor. Canonical polyadic decomposition (CPD) decomposes the three-dimensional multi-subject fMRI data into shared brain spatial activation components, shared time process components, and subject-specific intensity difference information. These information can provide reliable population characteristics for brain cognitive and brain disease research.

[0003] However, the spatial dimension of three-dimensional multi-subject fMRI data is much larger than the time dimension and the subject dimension, which also leads to slow updating and easy falling into local optimization of the CPD method using alternating least squares (ALS). Moreover, the CPD based on ALS requires the data to strictly comply with the CPD model. In addition, there are differences in response speed and blood dynamics between different subjects, which leads to the problem of time difference in multi-subject fMRI data. To solve this problem, shift-invariant CPD allows the time process components of different subjects to have time differences based on CPD, and estimates the time delay of each subject. However, shift-invariant CPD is also based on ALS, so shift-invariant CPD still has the problem of slow updating and easy falling into local optimization. SUMMARY

[0004] The purpose of the present application is to provide a shift-invariant CPD method of alternating rank R and rank 1 least squares method with spatial dimension average down-sampling compression. First, the spatial dimension of three-dimensional multi-subject fMRI data is processed by average down-sampling, which not only preserves the relative position relationship between feature points, but also effectively reduces the spatial dimension. Second, the alternating rank R and rank 1 least squares method is used to effectively relax the shift-invariant CPD model. Finally, the shared brain spatial activation component is added with spatial orthogonal constraint, which effectively improves the robustness of the algorithm to noise and reduces the influence of spatial difference between subjects.

[0005] The technical scheme of the present application is based on the time shift invariant CPD method, using average downsampling to reduce the spatial dimension size of three-dimensional multi-subject fMRI data; then using the shift invariant CPD based on alternating rank R and rank 1 least square method to update the shared brain spatial activation component, shared time process component, subject intensity difference information and subject time delay information; meanwhile, in each iteration, the shared brain spatial activation component is updated twice, and the spatial orthogonal constraint is added. The specific steps are as follows:

[0006] Step 1: input three-dimensional multi-subject fMRI data , component number and average downsampling width , wherein , , and represent the number of voxels (spatial dimension) in the , and directions of the three-dimensional brain space respectively, represents the number of scans (time dimension), represents the number of subjects (subject dimension); is the component number, which is an integer and satisfies .

[0007] Step 2: average down-sampling of the input three-dimensional multi-subject fMRI data in the one-dimensional brain space dimension with a specified width . The three-dimensional multi-subject fMRI data after sampling is :

[0008]

[0009] wherein , , , .

[0010] Step 3: randomly initialize the compressed brain spatial activation component , shared time process component and subject intensity difference information , initialize the subject time delay information as a zero matrix, initialize the iteration number , maximum iteration number and convergence condition . Wherein ; ; ; ; is an integer and satisfies .

[0011] Fourth step: First update of the joint mixing matrix . According to the shared time process component , the subject time delay information and the subject intensity difference information , generate the joint mixing matrix :

[0012]

[0013] wherein, is the shared time process component containing the time delay information of the first subject. is the shared time process component containing the time delay information of the second subject. is the shared time process component containing the time delay information of the third subject. is the shared time process component containing the time delay information of the fourth subject. is the shared time process component containing the time delay information of the fifth subject. is the shared time process component containing the time delay information of the sixth subject.

[0014] Fifth step: First update of the spatial compressed brain spatial activation component :

[0015]

[0016] wherein, “ ” is the pseudo-inverse operation, “ ” is the matrix transpose, is the mod 1 expansion (tensor matrix) of the tensor .

[0017] Sixth step: Second update of the compressed brain spatial activation component . Apply the orthogonal constraint to so that the components are mutually independent:

[0018]

[0019] wherein, is the left singular matrix, “ ” is the reduced singular value decomposition.

[0020] Seventh step: Second update of the joint mixing matrix :

[0021]

[0022] At the same time, set , the number of iterations , the maximum number of iterations , convergence condition , preparing for cyclic iteration updating shared time process component , column vector data of subject time delay information and subject intensity difference information .

[0023] Eighth step: select the joint mixing matrix , and convert it into matrix :

[0024]

[0025] Ninth step: update the column vector of shared time component .

[0026]

[0027] wherein, , , is the exponential function, is the dot product operation, is the frequency domain form (discrete Fourier transform) of , is the th row vector of . Tenth step: update the column vector of subject time delay information

[0028] according to the literature “KUANG L D, LIN Q H, GONG X F, et al. Shift-invariant canonical polyadic decomposition of complex-valued multi-subject fMRI data with a phase sparsity constraint[J]. IEEE transactions on medical imaging, 2020, 39(4): 844-853.”.

[0029]

[0030] wherein, is the th row vector of ​​​​​column vector.

[0031] Eleventh step: updating subject intensity difference information column vector .

[0032]

[0033] Twelfth step: calculating mean square error:

[0034]

[0035] If and , and , the inner iteration number , jump to the ninth step; if or , and , the inner iteration number , jump to the eighth step; if , jump to the thirteenth step.

[0036] Thirteenth step: calculating mean square error:

[0037]

[0038] Wherein, “⊙” is Khatri-Rao product, “ ” represents the Frobenius norm of matrix, that is, the square root of the sum of squares of all elements. If or , jump to the fourteenth step; otherwise, , jump to the fourth step.

[0039] Fourteenth step: restoring the shared brain space components before compression . According to formula (2), the joint mixing matrix is obtained, and then according to formula (12) (13), it is restored:

[0040]

[0041]

[0042] Fifteenth step: outputting the shared brain space components , shared time process components , subject intensity difference information and subject time delay information .

[0043] The effect achieved by the present application is that the task-related components of the multi-subject fMRI data of the task state can be quickly and effectively estimated. In the fMRI data analysis of the 16-subject tapping finger task, compared with the shift-invariant CPD method, the performance of the task-related shared brain spatial information is improved by 20%, the performance of the shared time information is improved by 8%, and the running time required for iteration is reduced by 91% under the condition of setting the same initial value. Compared with the independent component analysis and shift-invariant canonical polyadic decomposition (ICASCP), the performance of the task-related shared brain spatial information and shared time information is improved by 9% and 4% respectively. Therefore, the present application can quickly and effectively extract the brain function information shared by multiple subjects, and has good development prospects in future large-scale high-dimensional fMRI data research and intelligent medical treatment and the like. BRIEF DESCRIPTION OF DRAWINGS

[0044] Figure 1 is the workflow diagram of the alternating rank R and rank 1 shift-invariant CPD method for spatially compressed multi-subject fMRI provided by the present application. DETAILED DESCRIPTION

[0045] One specific embodiment of the present application will be described in detail below in combination with the technical solutions and the drawings.

[0046] First step: input the three-dimensional multi-subject task state fMRI data of the tapping finger task . Among them, the number of subjects , the number of scans , the number of brain spatial voxels 59610.

[0047] Second step: according to formula (1), the three-dimensional multi-subject fMRI data is subjected to average down-sampling with a specified width in the one-dimensional brain spatial dimension to obtain the sampled multi-subject fMRI data .

[0048] Third step: randomly initialize the compressed shared brain spatial activation components , the shared time process components and the subject intensity difference information , initialize the subject time delay information , initialize the number of iterations , the maximum number of iterations and the convergence condition .

[0049] Step 4: Update the joint mixing matrix for the first time according to equation (2). .

[0050] Step 5: According to formula (3), update the compressed brain spatial activation components for the first time. .

[0051] Step 6: According to formula (4), activate the compressed brain spatial components. Apply orthogonal constraints to update the compressed spatial activation components of the brain for the second time. .

[0052] Step 7: Update the joint mixing matrix for the second time according to equation (5). At the same time, set Number of iterations Maximum number of iterations Convergence condition We are preparing to perform an inner loop iteration.

[0053] Step 8: Select the joint mixing matrix according to equation (6). The List Transform it into a matrix .

[0054] Step 9: Update the shared time process components according to equation (7). column vectors .

[0055] Step 10: Update the subject's time delay information according to equation (8). column vectors .

[0056] Step 11: Update the subject intensity difference information according to equation (9). column vectors .

[0057] Step 12: Calculate the mean square error according to equation (10). ,like and ,at the same time number of inner iterations Jump to step nine; if or ,at the same time Then the number of inner iterations Jump to step eight; if Proceed to step thirteen.

[0058] Step 13: Calculate the mean square error according to equation (11). ,like or If yes, go to step 14; otherwise, go to step 4.

[0059] Step 14: Reduce the shared brain space components before compression according to formula (2) (12) (13) .

[0060] Step 15: Output the shared brain space components , shared time course components , subject intensity difference information and subject time delay information .

Claims

1. A spatially compressed multi-subject fMRI alternating rank R and rank-1 shift-invariant CPD method, characterized in that, Based on the time-shift invariant CPD method, spatial dimension averaging downsampling is used to reduce the spatial dimension of three-dimensional multi-subject fMRI data; then, time-shift invariant CPD based on alternating rank R and rank-1 least squares is used to update the shared brain spatial activation components, shared temporal process components, subject intensity difference information, and subject time delay information; at the same time, in each iteration, the shared brain spatial activation components are updated twice and spatial orthogonal constraints are added. Includes the following steps: Step 1: Input three-dimensional multi-subject fMRI data , number of components and average downsampling width ,in, , , and Representing the three-dimensional brain space , and Number of voxels in direction Belongs to the spatial dimension. Indicates the number of scans. Belongs to the time dimension. Indicates the number of subjects. Belongs to the subject dimension, are integers and satisfy ; Step 2: Process the input three-dimensional multi-subject fMRI data Specify width in one-dimensional brain space Average downsampling, sampled three-dimensional multi-subject fMRI data : in, , , , ; Step 3: Randomly initialize the compressed brain spatial activation components Shared time process components Information on differences in intensity between subjects Initialize subject time delay information The zero matrix is ​​initialized with the number of iterations. Maximum number of iterations and convergence conditions ,in, , , , , are integers and satisfy ; Step 4: First update of the joint mixing matrix According to the shared time process components Subject time delay information Information on differences in intensity between subjects According to equation (2), a joint mixing matrix is ​​generated. : in, For including the first The shared time process component of time delay information among subjects. for Time shifted One point, It is by Move to the left , Or move to the right , get; Step 5: First update of brain spatial activation components after spatial compression : in, This is a pseudo-inverse operation. For matrix transpose, For tensor The modulo 1 expansion; Step 6: Second update of compressed brain space activation components ,right Apply orthogonal constraints to ensure that the components are uncorrelated: in, It is a singular matrix with left and right sides. To reduce singular value decomposition; Step 7: Second update of the joint mixing matrix : At the same time, set Number of iterations Maximum number of iterations Convergence conditions Prepare for iterative updates of shared time process components. Subject time delay information Information on differences in intensity between subjects column vector data ; Step 8: Select the joint mixing matrix The List Transform it into a matrix : Step 9: Update the shared time component column vectors : in, , , It is an exponential function. It's a dot product operation. for , yes The frequency domain form after Fourier transform, yes The frequency domain form after Fourier transform, for The Row vectors; Step 10: Update the subject's time delay information column vectors : in, for The Column vectors; Step 11: Update the subject intensity difference information column vectors : Step 12: Calculate the mean square error: like and ,at the same time number of inner iterations Jump to step nine; if or ,at the same time Then the number of inner iterations Jump to step eight; if Proceed to step thirteen; Step 13: Calculate the mean square error: Where ⊙ represents the Khatri-Rao product, denoted by Frobenius norm, which is the square root of the sum of the squares of all elements; if or If yes, then proceed to step fourteen; otherwise, Jump to step four; Step 14: Restore the shared brain space components before compression According to equation (2), the joint mixing matrix is ​​obtained. Then, restore according to equations (12) and (13): Step 15: Output shared brain space components Shared time process components Information on differences in subject intensity Time delay information of the subjects .

Images