CN103961103A - Method for performing phase correction on ICA estimation components of plural fMRI data - Google Patents

Abstract

The invention provides a method for performing phase correction on ICA estimation components of plural fMRI data and belongs to the field of analysis of plural fMRI data. Based on the time process component (please see the formula in the specification) estimated by ICA, the phase angle theta[k] is estimated, preliminary phase correction is performed on the ICA estimation components (please see the formulae in the specification) with phase ambiguity, and primary phase correction signals (please see the formulae in the specification) are obtained; the primary phase correction signals (please see the formulae in the specification) are subjected to symbol ambiguity detection and removal through prior information which can be easily acquired and a correlation coefficient method. Because one ICA estimation component (please see the formula in the specification) has larger non-annular degree than the other ICA estimation component (please see the formula in the specification), error influences caused by high-amplitude noise voxel can be avoided. Due to the utilization of the prior information, the symbol ambiguity can be accurately detected and eliminated. When ICA estimation components, acquired under movement stimulation, of 16 pieces of tested plural fMRI data are subjected to phase correction, the accuracy of the method based on the ICA estimation component (please see the formula in the specification) is only 81.25 percent, but the accuracy rate of the method is 100%, and a guarantee is provided for ICA of multiple tested plural fMRI data.

Description

A kind of ICA to plural fMRI data is estimated the method that line phase is proofreaied and correct that is divided into

Technical field

The present invention relates to a kind of ICA analytical method of plural fMRI data, particularly relate to a kind of ICA to plural fMRI data and be estimated the method that line phase is proofreaied and correct that is divided into.

Background technology

Independent component analysis (independent component analysis, ICA) is a kind of analytical method of data-driven, claims again blind source to separate (blind source separation, BSS).ICA, without any prior information, just can estimate at utmost separate source signal and hybrid parameter thereof from mixed signal, obtains extensive use in fields such as voice, image, processing of biomedical signals.Because people are limited to the degree of awareness of brain, from 1998, ICA is at Functional MRI (functional magnetic resonance imaging, fMRI) in the analysis of data, obtained paying much attention to and effectively applying, have and from fMRI data (being mixed signal), estimate at utmost separate space activation brain district composition (spatial activations, be source signal) and time course composition (time courses, be hybrid parameter) ability, and then provide foundation for brain function analysis and clinical diagnosis.

But because do not utilize prior information, the estimated signal of ICA has two kinds of ambiguities, the one, order ambiguity (permutation ambiguity), the order of ICA estimated signal and the order of source signal may be different.In fMRI data analysis, people are often indifferent to ICA and estimate the order of composition, but are concerned about the content of estimating composition, and therefore, order ambiguity can not analyzed and cause adverse effect fMRI.The 2nd, amplitude ambiguity (scaling ambiguity), the amplitude of ICA estimated signal and the amplitude of source signal may be different.For plural fMRI data, the amplitude ambiguity of ICA comprises again amplitude ambiguity and phase ambiguity.Taking relevant (task-related) composition of task as example, ICA activates brain district composition in estimated space and can be expressed as in formula, s _*for source signal; be a random complex constant, d _*for amplitude, θ _*for phase angle.The amplitude ambiguity of ICA is embodied in estimated signal with source signal s _*between complex constant difference, has comprised that respectively amplitude ambiguity is the amplitude d between the two _*difference, and phase ambiguity is the phase angle θ between the two _*difference.In single tested plural fMRI data analysis, amplitude d _*ambiguity can eliminate by method for normalizing (even d _*=1), phase angle θ _*ambiguity show as real part-imaginary part Joint Distribution with respect to s _*real part-imaginary part Joint Distribution occur θ _*rotation, so also do not affect.But, in the time carrying out how tested plural fMRI data analysis (also claiming group analysis, group analysis), phase angle θ _*ambiguity group analysis result is had a serious impact.Still, taking task Related Component as example, suppose that total P is individual tested, P ICA of process and amplitude normalization, the space obtaining is activated brain district composition and is respectively: i=1 ..., P.For different tested, phase place there is the random difference opposite sex.That is to say, with respect to source signal real part-imaginary part Joint Distribution, real part-imaginary part Joint Distribution may produce the different anglecs of rotation when by how tested while being added and then asking for group average signal, there is gross error in the possibility of result.For example,, when how tested in 360 °, present approximately while distributing uniformly, P essence is other than ring type signal average signal can become mistakenly annular signal.

For above-mentioned phase ambiguity problem, have at present a kind of based on estimated activation brain district, the space composition of ICA estimate the anglec of rotation and then proofread and correct method (the Rodriguez of phase place, P.A., Calhoun, V.D., Adali, T., 2012.De-noising, phase ambiguity correction and visualization techniques forcomplex-valued ICA of group fMRI data.Pattern Recognition45,2050-2063).The method is to maximize real part energy be that principle is asked for being intended to will high amplitude voxel (highmagnitude voxels) major part concentrate on the positive axis of complex field real part.In addition, owing to being maximization real part energy, may have 180 ° of rotation errors, also claim symbol ambiguity, high amplitude voxel possibility major part concentrates on the negative semiaxis of complex field real part.For this reason, the method has been checked real part sum.If be less than zero with value, will be multiplied by (1), be equivalent to by rotate 180 °, then will as final phase correction signal; Otherwise, if correlation coefficient is not less than zero, directly output when signal to noise ratio when higher, the method can obtain good correction result.

But because the signal to noise ratio of plural fMRI data is very low, the signal to noise ratio of ICA estimated signal is also lower.Be embodied in, the estimated space of ICA activate brain district composition (as ) in contain high amplitude noise voxel a large amount of, that cannot eliminate with threshold method (removing the voxel of amplitude lower than a certain threshold value).At this moment, if adopt above-mentioned maximization the principle of real part energy, concentrates in the high amplitude voxel of complex field real part positive axis and may contain a large amount of noise voxels, and then cause the anglec of rotation misjudgment and symbol ambiguity inspection mistake.Therefore, high amplitude noise voxel has had a strong impact on the performance that activates the method for correcting phase of brain district composition based on space.

Summary of the invention

The object of the invention is to, provide a kind of ICA to plural fMRI data of more robust to be estimated the method that is divided into line phase correction, effectively avoid the negative effect of high amplitude noise voxel.

Technical scheme of the present invention is, based on the estimated time course composition of ICA instead of brain district composition is activated in space carry out phase angle θ _kestimate, k=1 ..., K, the number that K is independent element, estimates that principle is to maximize real part energy:

${\mathrm{\θ}}_k=\arg \underset{{\mathrm{\θ}}_k}{\max }\underset{t=1}{\overset{N}{\mathrm{\Σ}}}{\left(\mathrm{Re}\right\{{\hat{a}}_{k,t}{e}^{-j{\mathrm{\θ}}_k}\left\}\right)}^2$

In formula, Re represents real part, and t represents time point, N is count total time (i.e. the scanning times to tested full brain).This estimation principle will activate the positive axis that concentrates on complex field real part time response of voxel (activated voxels), namely activation voxel major part is concentrated on to the positive axis of complex field real part.

Phase angle θ _kto ask for process as follows: application real part and imaginary part build matrix calculate A _kcovariance matrix and to cov (A _k) carry out feature decomposition, obtain cov (A _k)=V Λ V ^t, Λ is the diagonal matrix being made up of eigenvalue, V is the orthogonal matrix being made up of characteristic vector:

$V=\left[\begin{array}{cc}\cos \left({\mathrm{\θ}}_k\right)& \sin \left({\mathrm{\θ}}_k\right)\\ -\sin \left({\mathrm{\θ}}_k\right)& \cos \left({\mathrm{\θ}}_k\right)\end{array}\right]$

Try to achieve phase angle θ from matrix V _k=arccos (v ₁₁), v ₁₁for the first row first row element of matrix V, and then obtain with preliminary phase correction signal with

Next, the present invention adopts the prior information and the correlation coefficient process that are easy to obtain, to preliminary phase correction signal with carry out the detection of symbol ambiguity and removal.Prior information comprises the prior information a about time course composition _ref, as the task stimulus signal (stimuli) of fMRI data; Or activate the prior information s of brain district composition about space _ref, as DMN (default mode network) network; The correlation coefficient process that detects and remove symbol ambiguity is specific as follows: calculate with a _refor with s _refcorrelation coefficient, if correlation coefficient is less than zero, there is symbol ambiguity in explanation, will with be multiplied by (1), obtain final phase calibration signal with otherwise, if correlation coefficient is not less than zero, with be final phase calibration signal.

Effect and benefit that the present invention reaches are that, compared with the existing method for correcting phase that activates brain district composition based on space, the method for correcting phase based on time course composition can be avoided the erroneous effects of high amplitude noise voxel.Its main cause is, the estimated time course composition of ICA activates brain district composition and has larger other than ring type degree than space.For example, for the plural fMRI data that gather under Motor stimulation, the other than ring type degree of estimated task procedure component correlation time of ICA is up to 0.89, and space is activated the other than ring type degree of brain district composition and only had 0.29 (maximum other than ring type degree=1, other than ring type degree=0 of annular signal).Larger other than ring type degree means that real part-imaginary component Butut of complex signal presents more obvious deflection θ _keven under the interference of high amplitude noise voxel, also can clearly be detected.And, during due to elimination symbol ambiguity, utilize prior information, correlation coefficient process can detect exactly 180 ° of rotation errors that may exist.The ICA of the 16 tested plural fMRI data that gather under to Motor stimulation is estimated and is divided into line phase timing, and accuracy of the present invention is 100%, and accuracy based on brain district composition is activated in space only has 81.25%.Therefore, the present invention can be in the ICA of how tested plural fMRI data analyzes, and brain district composition is activated in all tested time course compositions and space and be corrected into consistent phase place, and then carry out the group analysis such as correct addition is average.

Brief description of the drawings

Accompanying drawing is the concrete steps that the present invention proofreaies and correct the task Related Component of single tested plural fMRI data.

Detailed description of the invention

Below in conjunction with technical scheme and accompanying drawing, describe a specific embodiment of the present invention in detail.

Supposing the tested plural fMRI data of the list gathering under existing Motor stimulation, is 165 to the scanning times of tested full brain.This plural number fMRI data are carried out to ICA, the task Related Component of obtaining go forward side by side line phase proofread and correct concrete steps as shown in drawings.

The first step, carries out PCA (principle component analysis) dimensionality reduction by tested list plural fMRI data.Estimate independent element number according to theory of information criterion, result is 30, therefore adopt PCA that plural fMRI data are down to 30 dimensions from 165 dimensions.

Second step, the plural fMRI data after adopting plural EBM (entropy bound minimization) algorithm to PCA dimensionality reduction are carried out ICA, obtain 30 plural ICA and estimate compositions, therefrom select task Related Component, comprise time course composition and active region composition between Naokong

The 3rd step, chooses task relevant carry out phase angle θ _*estimate: first build matrix then calculate covariance matrix and feature decomposition cov (A _*)=V Λ V ^t, obtain orthogonal matrix V, and then try to achieve phase angle θ _*=arccos (v ₁₁).

The 4th step, estimates composition to the ICA with phase ambiguity with carry out preliminary phasing, obtain with

The 5th step, carries out the detection of symbol ambiguity and removal to preliminary phase correction signal.First by the task stimulus signal of fMRI data and hemodynamics receptance function (hemodynamic response function) convolution and obtain a _ref, then calculate with a _refcorrelation coefficient, if correlation coefficient <0, there is symbol ambiguity in explanation, will with be multiplied by (1), obtain final phase calibration signal with otherwise, if correlation coefficient>=0, with be final phase calibration signal.

Claims (1)

1. the ICA of plural fMRI data is estimated and is divided into the method that line phase is proofreaied and correct, it is characterized in that, based on the estimated time course composition of ICA instead of brain district composition is activated in space carry out phase angle θ _kestimate, estimate that principle is to maximize real part energy:

${\mathrm{\&theta;}}_k=\arg \underset{{\mathrm{\&theta;}}_k}{\max }\underset{t=1}{\overset{N}{\mathrm{\&Sigma;}}}{\left(\mathrm{Re}\right\{{\hat{a}}_{k,t}{e}^{-j{\mathrm{\&theta;}}_k}\left\}\right)}^2$

In formula, Re represents real part, and t represents time point, N counts total time; Phase angle θ _kto ask for process as follows: application real part and imaginary part build matrix calculate A _kcovariance matrix and to cov (A _k) carry out feature decomposition, obtain cov (A _k)=V Λ V ^t, Λ is the diagonal matrix being made up of eigenvalue, V is the orthogonal matrix being made up of characteristic vector:

$V=\left[\begin{array}{cc}\cos \left({\mathrm{\&theta;}}_k\right)& \sin \left({\mathrm{\&theta;}}_k\right)\\ -\sin \left({\mathrm{\&theta;}}_k\right)& \cos \left({\mathrm{\&theta;}}_k\right)\end{array}\right]$

Try to achieve phase angle θ from matrix V _k=arccos (v ₁₁), v ₁₁for the first row first row element of matrix V, obtain with preliminary phase correction signal with

The prior information that employing is easy to obtain and correlation coefficient process, to preliminary phase correction signal with carry out the detection of symbol ambiguity and removal; Prior information comprises the prior information a about time course composition _refor activate the prior information s of brain district composition about space _ref; The correlation coefficient process that detects and remove symbol ambiguity is specific as follows: calculate with a _refor with s _refcorrelation coefficient, if correlation coefficient is less than zero, there is symbol ambiguity in explanation, will with be multiplied by (1), obtain final phase calibration signal with otherwise, if correlation coefficient is not less than zero, with be final phase calibration signal.