WO2015085889A1 - 磁共振快速参数成像方法和系统 - Google Patents

磁共振快速参数成像方法和系统 Download PDF

Info

Publication number
WO2015085889A1
WO2015085889A1 PCT/CN2014/093168 CN2014093168W WO2015085889A1 WO 2015085889 A1 WO2015085889 A1 WO 2015085889A1 CN 2014093168 W CN2014093168 W CN 2014093168W WO 2015085889 A1 WO2015085889 A1 WO 2015085889A1
Authority
WO
WIPO (PCT)
Prior art keywords
image sequence
target image
matrix
magnetic resonance
parameter
Prior art date
Application number
PCT/CN2014/093168
Other languages
English (en)
French (fr)
Inventor
彭玺
梁栋
刘新
郑海荣
Original Assignee
中国科学院深圳先进技术研究院
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 中国科学院深圳先进技术研究院 filed Critical 中国科学院深圳先进技术研究院
Priority to US15/103,845 priority Critical patent/US10018698B2/en
Publication of WO2015085889A1 publication Critical patent/WO2015085889A1/zh

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • G01R33/56—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution
    • G01R33/561—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution by reduction of the scanning time, i.e. fast acquiring systems, e.g. using echo-planar pulse sequences
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/50—NMR imaging systems based on the determination of relaxation times, e.g. T1 measurement by IR sequences; T2 measurement by multiple-echo sequences
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • G01R33/56—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution
    • G01R33/5608—Data processing and visualization specially adapted for MR, e.g. for feature analysis and pattern recognition on the basis of measured MR data, segmentation of measured MR data, edge contour detection on the basis of measured MR data, for enhancing measured MR data in terms of signal-to-noise ratio by means of noise filtering or apodization, for enhancing measured MR data in terms of resolution by means for deblurring, windowing, zero filling, or generation of gray-scaled images, colour-coded images or images displaying vectors instead of pixels
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • G01R33/56—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution
    • G01R33/561—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution by reduction of the scanning time, i.e. fast acquiring systems, e.g. using echo-planar pulse sequences
    • G01R33/5611—Parallel magnetic resonance imaging, e.g. sensitivity encoding [SENSE], simultaneous acquisition of spatial harmonics [SMASH], unaliasing by Fourier encoding of the overlaps using the temporal dimension [UNFOLD], k-t-broad-use linear acquisition speed-up technique [k-t-BLAST], k-t-SENSE

Definitions

  • the present invention relates to the field of imaging technology, and in particular, to a magnetic resonance rapid parameter imaging method and system.
  • Magnetic resonance imaging is an imaging technique that uses nuclear nucleus to generate signals and reconstruct imaging in a magnetic field. Magnetic resonance imaging is widely used in imaging because of its safety, speed, and accuracy. Among them, magnetic resonance parameter imaging provides a method for quantitative analysis of biochemical characteristics of living tissue. However, due to the need to continuously acquire a series of images, the magnetic resonance parameters are slow to image, which limits their widespread clinical use.
  • Compressed sensing technology is an advanced technology for accelerating magnetic resonance imaging.
  • the magnetic resonance parameter imaging method based on compressed sensing includes two steps of image sequence reconstruction and parameter estimation. Due to the existence of under-sampling and sampling noise, it is inevitable to introduce reconstruction errors in the reconstruction of modern advanced level of compression-based magnetic resonance parameter image sequences, which will affect the accuracy of parameter estimation.
  • a method of rapid magnetic parameter imaging of magnetic resonance comprising:
  • the step of reconstructing the target image sequence according to the undersampled magnetic resonance signal and the prior information, the step of obtaining the target image sequence includes:
  • the target image sequence after initialization is reconstructed by a method based on nonlinear filtering sparse reconstruction, and then the target image sequence is obtained.
  • the step of reconstructing the initialized target image sequence by using a nonlinear filtering sparse reconstruction method to solve the target image includes:
  • the step of obtaining a low rank approximation to the Hankel matrix such that the time signal corresponding to each spatial position of the target image sequence is exponentially attenuated comprises:
  • a magnetic resonance rapid parameter imaging system comprising:
  • a signal acquisition module configured to scan the target to obtain an undersampled magnetic resonance signal
  • a prior information obtaining module for acquiring a priori information of the parameter model
  • An image reconstruction module configured to reconstruct a target image sequence according to the undersampled magnetic resonance signal and the prior information, to obtain a target image sequence
  • a parameter estimation module configured to substitute the target image sequence into a parameter estimation model, obtain a target parameter, and generate a parameter image.
  • the image reconstruction module includes:
  • An initializing unit configured to initialize a target image sequence by a zero-padded Bolly transform method
  • the operation unit is configured to reconstruct the initialized target image sequence by a sparse reconstruction method based on nonlinear filtering, and then obtain a target image sequence.
  • the operation unit includes:
  • a sparse constraint unit configured to perform soft threshold processing on coefficients of the target image sequence in a sparse variation domain based on sparse sampling theory
  • a matrix construction unit configured to construct a time signal corresponding to each spatial position of the target image sequence into a Hankel matrix according to the parameter estimation model
  • a model constraint unit configured to obtain a low rank approximation for the Hankel matrix, such that a time signal corresponding to each spatial position of the target image sequence is exponentially attenuated
  • a signal fidelity unit configured to substitute the undersampled magnetic resonance signal into a K space of the target image sequence such that the target image sequence is consistent with the undersampled magnetic resonance signal, and the operation unit repeatedly performs the sparse constraint unit and the matrix construction unit The operation of the model constraint unit and the signal fidelity unit until the desired target image sequence converges to obtain the final target image sequence.
  • the model constraint unit comprises:
  • a low rank approximation matrix solving unit configured to obtain a low rank approximation matrix for the Hankel matrix by singular value decomposition
  • a Hankel matrix recovery unit configured to perform an average operation enhancement signal along an anti-angle direction of the low rank approximation matrix, and restore a Hankel structure of the low rank approximation matrix
  • a signal extracting unit configured to extract the target image sequence from the low rank approximation Hankel matrix The time signal corresponding to each spatial position of the column.
  • the above magnetic resonance rapid parameter imaging method and system use the undersampled magnetic resonance signal and the prior information to reconstruct the target image sequence to obtain the target image sequence, and substitute the target image sequence into the parameter estimation model to generate the parameter image. Since the a priori information of the parameter is introduced in the image sequence reconstruction step, the reconstruction error generated in the compressed sensing image reconstruction can be corrected, and the accuracy of generating the parameter image by the parameter estimation model is improved.
  • FIG. 1 is a flow chart of a magnetic resonance rapid parameter imaging method in an embodiment
  • FIG. 2 is a flowchart of a method for reconstructing a target image sequence according to an undersampled magnetic resonance signal and a priori information in an embodiment
  • 3 is a flowchart of a method for reconstructing the target image sequence to obtain the target image by a method based on nonlinear filtering sparse reconstruction;
  • FIG. 4 is a flow chart of a method for obtaining a low rank approximation for a Hankel matrix in such an embodiment that the time signal corresponding to each spatial position of the target image sequence is exponentially attenuated;
  • Figure 5 is a schematic structural view of a magnetic resonance parameter imaging system in an embodiment
  • FIG. 6 is a schematic structural diagram of an image reconstruction module in an embodiment
  • FIG. 7 is a schematic structural diagram of an operation unit in an image reconstruction module according to an embodiment
  • FIG. 8 is a schematic structural diagram of a model constraint unit in an arithmetic unit in an example.
  • a magnetic resonance rapid parameter imaging method includes the following steps:
  • step S10 an undersampled magnetic resonance signal of the target is acquired.
  • the target is scanned by a magnetic resonance scanner to acquire an undersampled magnetic resonance signal of the target.
  • Step S20 Acquire a priori information of the parameter model.
  • the preset parameter is a T2 parameter, that is, taking T2-weighted imaging as an example.
  • Magnetic resonance imaging MRI is an imaging technique that uses reconstructed imaging of signals generated by the resonance of a nucleus in a magnetic field. Since the MRI signal is weak, in order to improve the signal-to-noise ratio of the MRI signal, it is required to repeatedly use the same pulse sequence. The interval of this repeated excitation is called the repetition time (TR), and the RF pulse is collected to the echo signal. The time between them is called echo delay time (TE).
  • TR repetition time
  • TE echo delay time
  • an image highlighting a certain tissue characteristic parameter can be obtained by adjusting the repetition time TR and the echo time TE. This image is called a weighted image.
  • T2-weighted imaging refers to a weighted image that highlights the difference in tissue T2 relaxation (transverse relaxation).
  • the T2-weighted magnetic resonance image sequence generally satisfies the following exponential decreasing relationship:
  • r 0,l represents the proton density distribution of the first water component
  • r n represents the T2 weighted image corresponding to the nth echo time
  • r is the spatial coordinate
  • DTE is the echo interval
  • T 2l is the first species Hydrogen proton T2 value of water component.
  • the value of L ranges from 1 to 3.
  • a priori information is all prior knowledge about signal characteristics, such as sparsity, smoothness, exponential decay characteristics, and so on.
  • the exponential decay characteristic of the T2-weighted image sequence is determined by the physical model of the signal, and is a prior knowledge that can be known in advance, that is, the required prior information.
  • Step S30 performing target image sequence reconstruction according to the undersampled magnetic resonance signal and the prior information to obtain a target image sequence.
  • the target image sequence can be obtained by reconstructing the target image sequence using the prior information and the undersampled magnetic resonance signal.
  • the T2 weighted image sequence is obtained, thereby reducing the error generated by the T2 image sequence reconstruction and making the T2 parameter estimation. The result is more accurate.
  • step S40 the target image sequence is substituted into the parameter estimation model to obtain a target parameter, and a parameter image is generated.
  • the T2-weighted image sequence is substituted into (1) under the assumption of Gaussian noise.
  • the T2 parameter can be estimated, and then a parameter image can be generated.
  • the method of the present invention can also be used to image other parameters, such as T2 parameter imaging based on a single exponential, multi-exponential model, T2* parametric imaging, caused by interactions between water molecules and surrounding macromolecules.
  • T1 ⁇ relaxation parameter imaging or myelin fractional imaging and the like can also be used to image other parameters, such as T2 parameter imaging based on a single exponential, multi-exponential model, T2* parametric imaging, caused by interactions between water molecules and surrounding macromolecules.
  • T1 ⁇ relaxation parameter imaging or myelin fractional imaging and the like can also be used to image other parameters, such as T2 parameter imaging based on a single exponential, multi-exponential model, T2* parametric imaging, caused by interactions between water molecules and surrounding macromolecules.
  • step S30 includes:
  • Step S31 the target image sequence is initialized by a method of zero-padding Fourier transform.
  • step S33 the initialized target image sequence is reconstructed by a method based on nonlinear filtering sparse reconstruction, and then the target image sequence is obtained.
  • the time signals corresponding to each spatial position of the target image sequence may be arranged in the following form:
  • the matrix of the above form is called a Hankel matrix, where M represents the total number of echo times, that is, the number of images of the image sequence. As long as the image sequence satisfies the formula (1), the above Hankel matrix is a low rank matrix, and the rank is equal to L.
  • the rank of Hankel matrix is introduced into the target image sequence reconstruction model as a constraint condition.
  • the T2 image sequence is obtained by solving the method based on nonlinear filtering sparse reconstruction.
  • the solution formula is as follows:
  • ⁇ denotes the target image sequence
  • y is the undersampled data of K space
  • F is the Fourier transform matrix
  • W denotes the undersampling operator
  • B denotes the sparse transform matrix of the time domain (eg principal component analysis PCA transform)
  • ⁇ denotes the image domain Sparse transformation matrix (such as wavelet transform).
  • ⁇ controls the relative distance between the reconstructed image and the sampled data, usually proportional to the noise level.
  • the function H[r 1 (r), r 2 (r), ..., r M (r)] represents an operation for forming a Hankel matrix
  • W represents a set of image domain coordinate points.
  • step S33 includes:
  • Step S331 based on the sparse sampling theory, soft threshold processing is performed on coefficients of the target image sequence in a certain sparse variation domain.
  • the target image sequence is sparse in the sparse transform domain.
  • a certain sparse variation domain such as PCA, wavelet transform
  • the target image sequence is iterated using a formula.
  • the formula used is:
  • SoftThresholding() is a soft threshold function with a threshold of ⁇ .
  • Step S332 according to the parameter estimation model, construct a time signal corresponding to each spatial position of the target image sequence into a Hankel matrix.
  • the time signals corresponding to each spatial position of the target image sequence are arranged in the following form:
  • M represents the total number of echo times, that is, the number of images of the image sequence.
  • M represents the total number of echo times, that is, the number of images of the image sequence.
  • the above Hankel matrix is a low rank matrix, and the rank is equal to L.
  • Step S333 obtaining a low rank approximation for the Hankel matrix such that the time signal corresponding to each spatial position of the target image sequence is exponentially attenuated.
  • a low rank approximation Hankel matrix of the Hankel matrix is obtained.
  • the rank of the Hankel matrix is selected according to the parameter estimation model, and the exponential decay signal corresponding to the spatial position of the image sequence is extracted from the low rank approximation matrix of the Hankel matrix.
  • Step S334 substituting the undersampled magnetic resonance signal into the K space of the target image sequence such that the target image sequence is consistent with the undersampled magnetic resonance signal.
  • the undersampled magnetic resonance signal is made to coincide with the K-space data of the target image sequence.
  • y is the undersampled data of the K space
  • F is the Fourier transform matrix
  • W C represents the matrix operation of the unsampled position of the selected K space.
  • step S335 the above steps are repeated until the algorithm converges to obtain a target image sequence.
  • step S331, step S332, step S333, and step S334 are repeatedly performed until the algorithm of the target image sequence obtained converges to obtain a final target image sequence.
  • step S333 includes:
  • Step S3331 Find a low rank approximation matrix for the Hankel matrix by singular value decomposition.
  • U and V are left and right singular value vector matrices, respectively, and diagonal array S is a singular value matrix.
  • u k and v k denote the kth column of U and V, respectively, and ⁇ k is the kth singular value.
  • the low rank approximation matrix can be expressed as
  • step S3333 the enhancement signal is averaged along the diagonal direction of the low rank approximation matrix, and the Hankel structure of the low rank approximation matrix is restored.
  • Step S3335 extracting a time signal corresponding to each spatial position of the target image sequence from the low rank approximation Hankel matrix.
  • the anti-corner direction of the low rank approximation matrix is averaged to enhance the signal, the Hankel structure of the low rank approximation matrix is restored, and the time corresponding to each spatial position of the target image sequence is extracted from the low rank approximation Hankel matrix. signal.
  • the above method introduces the a priori information of the parameter as the prior knowledge into the image sequence reconstruction process, and corrects the error generated in the compressed sensing image reconstruction, and improves the passing parameters. Estimate the accuracy of the model generation parameter image.
  • a magnetic resonance rapid parameter imaging system includes a signal acquisition module 10, an a priori information acquisition module 20, an image reconstruction module 30, and a parameter estimation module 40.
  • the signal acquisition module 10 is configured to scan the target to obtain an undersampled magnetic resonance signal.
  • the target under-sampling magnetic resonance signal can be acquired by scanning the target with a magnetic resonance scanner.
  • the a priori information obtaining module 20 is configured to acquire a priori information of the parameter model.
  • the T2 parameter is taken as an example.
  • the T2-weighted magnetic resonance image sequence generally satisfies the following exponential decreasing relationship:
  • r 0,l represents the proton density distribution of the first water component
  • r n represents the T2 weighted image corresponding to the nth echo time
  • r is the spatial coordinate
  • DTE is the echo interval
  • T 2l is the first species Hydrogen proton T2 value of water component.
  • L usually ranges from 1 to 3.
  • a priori information is all prior knowledge about signal characteristics, such as sparsity, smoothness, exponential decay characteristics, and so on.
  • the exponential decay characteristic of the T2-weighted image sequence is determined by the physical model of the signal, and is a prior knowledge that can be known in advance, that is, the required prior information.
  • the image reconstruction module 30 is configured to perform target image sequence reconstruction according to the undersampled magnetic resonance signal and the prior information to obtain a target image sequence.
  • the target image sequence can be obtained by reconstructing the target image sequence using the prior information and the undersampled magnetic resonance signal.
  • the information that can be known in the T2 parameter estimation is introduced in the process of reconstructing the target image sequence, the error generated by the reconstruction of the target image sequence can be further reduced, so that the parameter estimation is more accurate.
  • the parameter estimation module 40 is configured to substitute the target image sequence into the parameter estimation model, obtain the target parameter, and generate the parameter image.
  • the T2-weighted image sequence is substituted into (1) under the assumption of Gaussian noise.
  • the T2 parameter can be estimated, and then a parameter image can be generated.
  • the method of the present invention can also be used to image other parameters, such as T2 parameter imaging based on a single exponential, multi-exponential model, T2* parametric imaging, caused by interactions between water molecules and surrounding macromolecules.
  • T1 ⁇ relaxation parameter imaging or myelin fractional imaging and the like can also be used to image other parameters, such as T2 parameter imaging based on a single exponential, multi-exponential model, T2* parametric imaging, caused by interactions between water molecules and surrounding macromolecules.
  • T1 ⁇ relaxation parameter imaging or myelin fractional imaging and the like can also be used to image other parameters, such as T2 parameter imaging based on a single exponential, multi-exponential model, T2* parametric imaging, caused by interactions between water molecules and surrounding macromolecules.
  • the image reconstruction module 30 includes:
  • the initializing unit 31 is configured to initialize the target image sequence by a method of zero-padding Fourier transform.
  • F is a Fourier transform matrix
  • W represents an undersampling operator
  • H represents a conjugate transpose.
  • the operation unit 33 is configured to reconstruct the initialized target image sequence by a method based on nonlinear filtering sparse reconstruction, and then obtain a target image sequence.
  • the time signals corresponding to each spatial position of the target image sequence may be arranged in the following form:
  • the matrix of the above form is called a Hankel matrix, where M represents the total number of echo times, that is, the number of images of the image sequence. As long as the image sequence satisfies the formula (1), the above Hankel matrix is a low rank matrix, and the rank is equal to L.
  • the rank of Hankel matrix is introduced into the target image sequence reconstruction model as a constraint condition.
  • the T2 image sequence is obtained by solving the method based on nonlinear filtering sparse reconstruction.
  • the solution formula is as follows:
  • ⁇ denotes the target image sequence
  • y is the undersampled data of K space
  • F is the Fourier transform matrix
  • W denotes the undersampling operator
  • B denotes the sparse transform matrix of the time domain (eg principal component analysis PCA transform)
  • ⁇ denotes the image domain Sparse transformation matrix (such as wavelet transform).
  • ⁇ controls the relative distance between the reconstructed image and the sampled data, usually proportional to the noise level.
  • the function H[r 1 (r), r 2 (r), ..., r M (r)] represents an operation for forming a Hankel matrix
  • ⁇ represents a set of coordinate points of an image domain.
  • the operation unit 33 includes:
  • the sparse constraint unit 331 is configured to perform soft threshold processing on the coefficients of the target image sequence in a certain sparse variation domain based on the sparse sampling theory.
  • the target image sequence is sparse in the sparse transform domain.
  • a certain sparse variation domain such as PCA, wavelet transform
  • the target image sequence is iterated using a formula.
  • the formula used is:
  • SoftThresholding() is a soft threshold function with a threshold of ⁇ .
  • the matrix construction unit 332 is configured to construct a time signal corresponding to each spatial position of the target image sequence into a Hankel matrix according to the parameter estimation model.
  • the time signals corresponding to each spatial position of the target image sequence are arranged in the following form:
  • M represents the total number of echo times, that is, the number of images of the image sequence.
  • M represents the total number of echo times, that is, the number of images of the image sequence.
  • the above Hankel matrix is a low rank matrix, and the rank is equal to L.
  • the model constraint unit 333 is configured to perform a low rank approximation on the Hankel matrix to make the target image sequence
  • the time signal corresponding to each spatial position of the column is exponentially attenuated.
  • a low rank approximation Hankel matrix of the Hankel matrix is obtained.
  • the rank of the Hankel matrix is selected according to the parameter estimation model, and the exponential decay signal corresponding to the spatial position of the image sequence is extracted from the low rank approximation matrix of the Hankel matrix.
  • the signal fidelity unit 334 is configured to substitute the undersampled magnetic resonance signal into the K space of the target image sequence such that the target image sequence is consistent with the undersampled magnetic resonance signal.
  • the formula is: The image sequence is made to coincide with the K-space data of the target image sequence.
  • y is the undersampled data of the K space
  • F is the Fourier transform matrix
  • W C represents the matrix operation of the unsampled position of the selected K space.
  • the arithmetic unit 33 repeatedly performs the operations of the sparse constraint unit 331, the matrix construction unit 332, the model constraint unit 333, and the signal fidelity unit 334 until the desired target image sequence converges to obtain a final target image sequence.
  • the model constraint unit 333 includes:
  • the low rank approximation matrix solving unit 3331 is configured to obtain a low rank approximation matrix for the Hankel matrix by singular value decomposition.
  • U and V are left and right singular value vector matrices, respectively, and diagonal array S is a singular value matrix.
  • u k and v k denote the kth column of U and V, respectively, and ⁇ k is the kth singular value.
  • the low rank approximation matrix can be expressed as:
  • the Hankel matrix recovery unit 3333 is configured to perform an average operation enhancement signal along the diagonal direction of the low rank approximation matrix to restore the Hankel structure of the low rank approximation matrix.
  • the signal extracting unit 3335 is configured to extract, from the low rank approximation Hankel matrix, a time signal corresponding to each spatial position of the target image sequence.
  • the anti-corner direction of the low rank approximation matrix is averaged to enhance the signal, the Hankel structure of the low rank approximation matrix is restored, and the time corresponding to each spatial position of the target image sequence is extracted from the low rank approximation Hankel matrix. signal.
  • the above method introduces the exponential decay signal of the image sequence into the image sequence reconstruction process as a priori knowledge, and corrects the error generated in the compressed sensing image reconstruction and improves the passing parameters. Estimate the accuracy of the model generation parameter image.

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Condensed Matter Physics & Semiconductors (AREA)
  • Radiology & Medical Imaging (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Signal Processing (AREA)
  • General Health & Medical Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Artificial Intelligence (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

一种磁共振快速参数成像的方法和系统。该方法包括：获取目标的欠采样磁共振信号（S10）；获取参数模型的先验信息（S20）；根据欠采样磁共振信号和先验信息进行目标图像序列重构，得到目标图像序列（S30）；将目标图像序列代入参数估计模型，得到目标参数，并生成参数图像（S40）。

Description

磁共振快速参数成像方法和系统 技术领域
本发明涉及成像技术领域，特别是涉及一种磁共振快速参数成像方法和系统。
背景技术
磁共振成像(MRI)是一种利用原子核在磁场内共振产生信号并重建成像的一种成像技术，磁共振成像因具有安全、快速、准确的优点而在成像领域中被广泛使用。其中，磁共振参数成像提供了可以定量分析生命组织生物化学特性的方法。然而，由于需要连续采集一系列图像，磁共振参数成像速度缓慢，从而限制了其临床的广泛使用。
压缩感知技术是一种加速磁共振成像的先进技术。基于压缩感知的磁共振参数成像方法包括图像序列重构和参数估计两个步骤。由于欠采样以及采样噪声的存在，现代先进水平的基于压缩感知的磁共振参数图像序列重构中将难免引入重建误差，进而影响参数估计的准确性。
发明内容
基于此，有必要针对磁共振参数成像不准确的问题，提供一种能提高磁共振参数成像准确性的磁共振快速参数成像方法。
此外，还有必要提供一种能提高磁共振参数成像准确性的磁共振快速参数成像系统。
一种磁共振快速参数成像的方法，所述方法包括：
获取目标的欠采样磁共振信号；
获取参数模型的先验信息；
根据欠采样磁共振信号和先验信息进行目标图像序列重构，得到目标图 像序列；
将所述目标图像序列代入参数估计模型，得到目标参数，并生成参数图像。
在其中一个实施例中，所述根据欠采样磁共振信号和先验信息进行目标图像序列重构，得到目标图像序列的步骤包括：
通过补零傅立叶变换方法初始化目标图像序列；
通过基于非线性滤波稀疏重建的方法对初始化后的目标图像序列进行重构，进而求解得到目标图像序列。
在其中一个实施例中，所述通过基于非线性滤波稀疏重建的方法对初始化后的目标图像序列进行重构，进而求解得到目标图像的步骤包括：
基于稀疏采样理论，对所述目标图像序列在某一稀疏变化域内的系数做软阈值处理；
根据所述参数估计模型，将目标图像序列每一空间位置对应的时间信号构造成Hankel矩阵；
对所述Hankel矩阵求低秩近似，以使得所述目标图像序列每一空间位置所对应的时间信号呈指数衰减；
将所述欠采样磁共振信号代入所述目标图像序列的K空间，以使得目标图像序列与所述欠采样磁共振信号一致；
重复上述步骤直至算法收敛，得到目标图像序列。
在其中一个实施例中，所述对所述Hankel矩阵求低秩近似，以使得所述目标图像序列每一空间位置所对应的时间信号呈指数衰减的步骤包括：
通过奇异值分解对所述Hankel矩阵求低秩近似矩阵；
沿所述低秩近似矩阵的反对角线方向作平均值运算增强信号，恢复所述低秩近似矩阵的Hankel结构；
从所述低秩近似Hankel矩阵中提取所述目标图像序列每一空间位置对应的时间信号。
一种磁共振快速参数成像系统，所述系统包括：
信号获取模块，用于对目标进行扫描得到欠采样磁共振信号；
先验信息获取模块，用于获取参数模型的先验信息；
图像重建模块，用于根据欠采样磁共振信号和先验信息进行目标图像序列重构，得到目标图像序列；
参数估计模块，用于将所述目标图像序列代入参数估计模型，得到目标参数，并生成参数图像。
在其中一个实施例中，所述图像重建模块包括：
初始化单元，用于通过补零博立叶变换方法初始化目标图像序列；
运算单元，用于通过基于非线性滤波的稀疏重建方法对初始化后的目标图像序列进行重构，进而求解得到目标图像序列。
在其中一个实施例中，所述运算单元包括：
稀疏约束单元，用于基于稀疏采样理论，对所述目标图像序列在某一稀疏变化域内的系数做软阈值处理；
矩阵构造单元，用于根据所述参数估计模型，将目标图像序列每一空间位置对应的时间信号构造成Hankel矩阵；
模型约束单元，用于对所述Hankel矩阵求低秩近似，以使得所述目标图像序列每一空间位置所对应的时间信号呈指数衰减；
信号保真单元，用于将所述欠采样磁共振信号代入所述目标图像序列的K空间，以使得目标图像序列与欠采样磁共振信号一致，运算单元通过重复执行稀疏约束单元、矩阵构造单元、模型约束单元和信号保真单元的运算，直到所求目标图像序列收敛，得到最终的目标图像序列。
在其中一个实施例中，所述模型约束单元包括：
低秩近似矩阵求解单元，用于通过奇异值分解对所述Hankel矩阵求低秩近似矩阵；
Hankel矩阵恢复单元，用于沿所述低秩近似矩阵的反对角线方向作平均值运算增强信号，恢复所述低秩近似矩阵的Hankel结构；
信号提取单元，用于从所述低秩近似Hankel矩阵中提取所述目标图像序 列每一空间位置对应的时间信号。
上述磁共振快速参数成像方法和系统，运用欠采样磁共振信号和先验信息进行目标图像序列重构得到目标图像序列，并将目标图像序列代入参数估计模型中生成参数图像。由于在图像序列重构步骤中引入了参数的先验信息，进而可以修正压缩感知图像重建中所产生的重建误差，提高通过参数估计模型生成参数图像的准确性。
附图说明
图1为一个实施例中磁共振快速参数成像方法的流程图；
图2为一个实施例中根据欠采样磁共振信号和先验信息进行目标图像序列重构得到目标图像序列的方法的流程图；
图3为一个实施例中通过基于非线性滤波稀疏重建的方法对所述目标图像序列重构进而求解得到所述目标图像的方法的流程图；
图4为一个实施例中对Hankel矩阵求低秩近似，以使得目标图像序列每一空间位置所对应的时间信号呈指数衰减的方法的流程图；
图5为一个实施例中磁共振参数成像系统的结构示意图；
图6为一个实施例中图像重建模块的结构示意图；
图7为一个实施例中图像重建模块中运算单元的结构示意图；
图8为一个实例中运算单元中模型约束单元的结构示意图。
具体实施方式
如图1所示，在一个实施例中，一种磁共振快速参数成像方法包括如下步骤：
步骤S10，获取目标的欠采样磁共振信号。
具体的，基于压缩感知理论，通过磁共振扫描仪对目标进行扫描，获取目标的欠采样磁共振信号。
步骤S20，获取参数模型的先验信息。
在本实施例中，预设参数为T2参数，即以T2加权成像为例。磁共振成像(magnetic resonance imaging，MRI)是利用原子核在磁场内共振产生的信号经重建成像的成像技术。由于MRI信号很弱，为提高MRI信号的信噪比，要求重复使用同一种脉冲序列，这个重复激发的间隔时间称为重复时间(repetition time，TR)，而射频脉冲放射后到采集回波信号之间的时间称为回波时间(echo delay time，TE)。为了评判被检测组织的各种参数，通过调节重复时间TR、回波时间TE，可以得到突出某种组织特征参数的图像，此图像称为加权像。在射频脉冲的激发下，人体组织内氢质子吸收能量处于激发状态；射频脉冲终止后，处于激发状态的氢质子恢复其原始状态，这个过程称为弛豫。T2加权成像是指突出组织T2弛豫(横向弛豫)差别的加权像。
T2加权磁共振图像序列一般满足如下指数递减关系：
Figure PCTCN2014093168-appb-000001
其中，r0,l表示第l种水成分的质子密度分布，rn表示第n个回波时间所对应的T2加权图像，r为空间坐标，DTE是回波间隔，T2l为第l种水成分的氢质子T2数值。通常情况下，L的取值范围为1～3。
先验信息即关于信号特性的一切先验知识，如稀疏性、平滑性、指数衰减特性等。在本实施例中，T2加权图像序列的指数衰减特性是由信号的物理模型决定的，是可预先获知的先验知识，也即所需的先验信息。
步骤S30，根据欠采样磁共振信号和先验信息进行目标图像序列重构，得到目标图像序列。
利用先验信息和欠采样磁共振信号进行目标图像序列重构，能够得到目标图像序列。在本实施例中，由于在目标图像序列重构的过程中引入了T2参数的先验信息，得到T2加权图像序列，从而可以减少T2图像序列重构所产生的误差，并使得T2参数估计的结果更加准确。
步骤S40，将目标图像序列代入参数估计模型，得到目标参数，并生成参数图像。
本实施例中，在假定高斯噪声的条件下，将T2加权图像序列代入(1) 式，运用最小二乘法就可估计得到T2参数，继而可生成参数图像。
在其他实施例中，本发明的方法还能用于实现其他参数的成像，例如基于单指数、多指数模型的T2参数成像，T2*参数成像、由水分子与周围大分子之间相互作用引起的T1ρ弛豫参数成像或髓鞘分数成像等。
如图2所示，在一个实施例中，步骤S30包括：
步骤S31，通过补零傅立叶变换的方法初始化目标图像序列。
本实施例中，通过补零快速傅立叶变换(zero-padding FFT)方法对目标图像序列进行初始化：ρ(0)＝FHWHy，其中ρ(0)表示初始化的目标图像序列，y为K空间的欠采样数据，F为傅立叶变换矩阵，W表示欠采样算子，H表示共轭转置。
步骤S33，通过基于非线性滤波稀疏重建的方法对初始化后的目标图像序列进行重构，进而求解得到目标图像序列。
本实施例中，具体的，以T2参数成像为例，根据T2参数估计模型，目标图像序列的每一空间位置对应的时间信号可排列成如下形式：
Figure PCTCN2014093168-appb-000002
上述形式的矩阵称为Hankel矩阵，其中，M表示回波时间总数，即图像序列的图像个数。只要图像序列满足式(1)，上述Hankel矩阵就是一个低秩矩阵，且秩等于L。
将Hankel矩阵的秩作为约束条件引入目标图像序列重建模型中通过基于非线性滤波稀疏重建的方法求解，得到T2图像序列。求解公式如下：
Figure PCTCN2014093168-appb-000003
其中ρ表示目标图像序列，y为K空间的欠采样数据，F为傅立叶变换 矩阵，W表示欠采样算子，B表示时间域的稀疏变换矩阵(如主成分分析PCA变换)，Ψ表示图像域的稀疏变换矩阵(如小波变换)。ε控制重建图像与采样数据间的相对距离，通常与噪声等级成正比。函数H[r1(r),r2(r),...,rM(r)]表示形成Hankel矩阵的运算，W表示图像域坐标点集合。
如图3所示，在一个实施例中，步骤S33包括：
步骤S331，基于稀疏采样理论，对目标图像序列在某一稀疏变化域内的系数做软阈值处理。
本实施例中，通过对目标图像序列在某一稀疏变化域(如PCA、小波变换)内的系数进行软阈值处理，使得目标图像序列在该稀疏变换域内是稀疏的。
具体的，利用公式对目标图像序列进行迭代。所用公式为：
(a)α(k)＝ΨBρ(k)，
(b)a(k+1)＝SoftThresholding(a(k),t)，
(c)r(k+1)＝BHΨHa(k+1)，
其中α为变换域系数，ρ表示目标图像序列，B表示时间域的稀疏变换矩阵(如PCA变换)，Ψ表示图像域的稀疏变换矩阵(如小波变换)，上标k表示第k次迭代，SoftThresholding()是软阈值函数，阈值为τ。
步骤S332，根据参数估计模型，将目标图像序列每一空间位置对应的时间信号构造成Hankel矩阵。
本实施例中，根据T2参数估计模型，将目标图像序列的每一空间位置对应的时间信号排列成如下形式：
Figure PCTCN2014093168-appb-000004
其中，M表示回波时间总数，即图像序列的图像个数。只要图像序列满足式(1)，上述Hankel矩阵就是一个低秩矩阵，且秩等于L。
步骤S333，对Hankel矩阵求低秩近似，以使得目标图像序列每一空间位置所对应的时间信号呈指数衰减。
本实施例中，求解出Hankel矩阵的低秩近似Hankel矩阵，其中Hankel矩阵的秩是根据参数估计模型选取的，从Hankel矩阵的低秩近似矩阵中提取图像序列空间位置所对应的指数衰减信号。
步骤S334，将欠采样磁共振信号代入目标图像序列的K空间，以使得目标图像序列跟欠采样磁共振信号一致。
本实施例中，利用公式：
Figure PCTCN2014093168-appb-000005
使得欠采样磁共振信号与目标图像序列的K空间数据一致。其中y为K空间的欠采样数据，F为傅立叶变换矩阵，WC表示取选取K空间未采样位置的矩阵运算。
步骤S335，重复上述步骤直至算法收敛，得到目标图像序列。
本实施例中通过重复执行步骤S331、步骤S332、步骤S333和步骤S334的运算，直到所求目标图像序列的算法收敛，得到最终的目标图像序列。
如图4所示，在一个实施例中，步骤S333包括：
步骤S3331，通过奇异值分解对Hankel矩阵求低秩近似矩阵。
本实施例中，利用公式：[U S V]＝svd(H)对Hankel矩阵做奇异值分解。其中U和V分别是左、右奇异值向量矩阵，对角阵S为奇异值矩阵。令uk和vk分别表示U和V的第k列，σk为第k个奇异值。取前L个最大的奇异值及其对应的奇异向量形成低秩近似矩阵，该低秩近似矩阵可以表示为
Figure PCTCN2014093168-appb-000006
步骤S3333，沿低秩近似矩阵的反对角线方向作平均值运算增强信号，恢复低秩近似矩阵的Hankel结构。
步骤S3335，从低秩近似Hankel矩阵中提取目标图像序列每一空间位置对应的时间信号。
本实施例中，对低秩近似矩阵的反对角线方向作平均值运算增强信号，恢复低秩近似矩阵的Hankel结构，并从低秩近似Hankel矩阵中提取目标图像序列每一空间位置对应的时间信号。
上述方法通过约束图像序列所形成的Hankel矩阵的低秩特性，将参数的先验信息作为先验知识引入图像序列重构过程中，修正了压缩感知图像重建中所产生的误差，提高了通过参数估计模型生成参数图像的准确性。
如图5所示，一种磁共振快速参数成像系统包括：信号获取模块10、先验信息获取模块20、图像重建模块30和参数估计模块40。
信号获取模块10，用于对目标进行扫描得到欠采样磁共振信号。
本实施例中，通过磁共振扫描仪对目标进行扫描就能获取目标的欠采样磁共振信号。
先验信息获取模块20，用于获取参数模型的先验信息。
在本实施例中，以T2参数为例。
T2加权磁共振图像序列一般满足如下指数递减关系：
Figure PCTCN2014093168-appb-000007
其中，r0,l表示第l种水成分的质子密度分布，rn表示第n个回波时间所对应的T2加权图像，r为空间坐标，DTE是回波间隔，T2l为第l种水成分的氢质子T2数值。通常而言，通常情况下L的取值范围为1～3。
先验信息即关于信号特性的一切先验知识，如稀疏性、平滑性、指数衰减特性等。在本实施例中，T2加权图像序列的指数衰减特性是由信号的物理模型决定的，是可预先获知的先验知识，也即所需的先验信息。
图像重建模块30，用于根据欠采样磁共振信号和先验信息进行目标图像序列重构得到目标图像序列。
利用先验信息和欠采样磁共振信号进行目标图像序列重构能够得到目标图像序列。在本实施例中，由于目标图像序列重构的过程中引入了T2参数估计中可预先获知的信息，从而可以进一步减少目标图像序列重构产生的误差，使得参数估计更加准确。
参数估计模块40，用于将目标图像序列代入参数估计模型，得到目标参数，并生成参数图像。
本实施例中，在假定高斯噪声的条件下，将T2加权图像序列代入(1) 式，运用最小二乘法就可估计得到T2参数，继而可生成参数图像。
在其他实施例中，本发明的方法还能用于实现其他参数的成像，例如基于单指数、多指数模型的T2参数成像，T2*参数成像、由水分子与周围大分子之间相互作用引起的T1ρ弛豫参数成像或髓鞘分数成像等。
如图6所示，在一个实施例中，图像重建模块30包括：
初始化单元31，用于通过补零傅立叶变换的方法初始化目标图像序列。
本实施例中，通过补零傅立叶变换方法对目标图像序列进行初始化：ρ(0)＝FHWHy，其中ρ(0)表示初始化的目标图像序列，y为K空间的欠采样数据，F为傅立叶变换矩阵，W表示欠采样算子，H表示共轭转置。
运算单元33，用于通过基于非线性滤波稀疏重建的方法对初始化后的目标图像序列重构，进而求解得到目标图像序列。
具体的，以T2参数成像为例，根据T2参数估计模型，目标图像序列的每一空间位置对应的时间信号可排列成如下形式：
Figure PCTCN2014093168-appb-000008
上述形式的矩阵称为Hankel矩阵，其中，M表示回波时间总数，即图像序列的图像个数。只要图像序列满足式(1)，上述Hankel矩阵就是一个低秩矩阵，且秩等于L。
将Hankel矩阵的秩作为约束条件引入目标图像序列重建模型中通过基于非线性滤波稀疏重建的方法求解，得到T2图像序列。求解公式如下：
Figure PCTCN2014093168-appb-000009
其中ρ表示目标图像序列，y为K空间的欠采样数据，F为傅立叶变换矩阵，W表示欠采样算子，B表示时间域的稀疏变换矩阵(如主成分分析PCA 变换)，Ψ表示图像域的稀疏变换矩阵(如小波变换)。ε控制重建图像与采样数据间的相对距离，通常与噪声等级成正比。函数H[r1(r),r2(r),...,rM(r)]表示形成Hankel矩阵的运算，Ω表示图像域坐标点集合。
如图7所示，在一个实施例中，运算单元33包括：
稀疏约束单元331，用于基于稀疏采样理论，对目标图像序列在某一稀疏变化域内的系数做软阈值处理。
本实施例中，通过对目标图像序列在某一稀疏变化域(如PCA、小波变换)内的系数进行软阈值处理，使得目标图像序列在该稀疏变换域内是稀疏的。
具体的，利用公式对目标图像序列进行迭代。所用公式为：
(a)α(k)＝ΨBρ(k)，
(b)a(k+1)＝SoftThresholding(a(k),t)，
(c)r(k+1)＝BHΨHa(k+1)，
其中α为变换域系数，ρ表示目标图像序列，B表示时间域的稀疏变换矩阵(如PCA变换)，Ψ表示图像域的稀疏变换矩阵(如小波变换)，上标k表示第k次迭代，SoftThresholding()是软阈值函数，阈值为τ。
矩阵构造单元332，用于根据参数估计模型，将目标图像序列每一空间位置对应的时间信号构造成Hankel矩阵。
本实施例中，根据T2参数估计模型，将目标图像序列的每一空间位置对应的时间信号排列成如下形式：
Figure PCTCN2014093168-appb-000010
其中，M表示回波时间总数，即图像序列的图像个数。只要图像序列满足式(1)，上述Hankel矩阵就是一个低秩矩阵，且秩等于L。
模型约束单元333，用于对Hankel矩阵求低秩近似，以使得目标图像序 列每一空间位置所对应的时间信号呈指数衰减。
本实施例中，求解出Hankel矩阵的低秩近似Hankel矩阵，其中Hankel矩阵的秩是根据参数估计模型选取的，从Hankel矩阵的低秩近似矩阵中提取图像序列空间位置所对应的指数衰减信号。
信号保真单元334，用于将欠采样磁共振信号代入目标图像序列的K空间，以使得目标图像序列跟欠采样磁共振信号一致。
本实施例中，通过公式：
Figure PCTCN2014093168-appb-000011
以使得图像序列与目标图像序列的K空间数据一致。其中y为K空间的欠采样数据，F为傅立叶变换矩阵，WC表示取选取K空间未采样位置的矩阵运算。
运算单元33通过重复执行稀疏约束单元331、矩阵构造单元332、模型约束单元333和信号保真单元334的运算，直到所求目标图像序列收敛，得到最终的目标图像序列。
如图8所示，在一个实施例中，模型约束单元333包括：
低秩近似矩阵求解单元3331，用于通过奇异值分解对Hankel矩阵求低秩近似矩阵。
本实施例中，利用公式：[U S V]＝svd(H)对Hankel矩阵做奇异值分解。其中U和V分别是左、右奇异值向量矩阵，对角阵S为奇异值矩阵。令uk和vk分别表示U和V的第k列，σk为第k个奇异值。取前L个最大的奇异值及其对应的奇异向量形成低秩近似矩阵，该低秩近似矩阵可以表示为:
Figure PCTCN2014093168-appb-000012
Hankel矩阵恢复单元3333，用于沿低秩近似矩阵的反对角线方向作平均值运算增强信号，恢复低秩近似矩阵的Hankel结构。
信号提取单元3335，用于从低秩近似Hankel矩阵中提取目标图像序列每一空间位置对应的时间信号。
本实施例中，对低秩近似矩阵的反对角线方向作平均值运算增强信号，恢复低秩近似矩阵的Hankel结构，并从低秩近似Hankel矩阵中提取目标图像序列每一空间位置对应的时间信号。
上述方法通过约束图像序列所形成的Hankel矩阵的低秩特性，将图像序列的指数衰减信号作为先验知识引入图像序列重构过程中，修正了压缩感知图像重建中所产生的误差，提高通过参数估计模型生成参数图像的准确性。以上实施例仅表达了本发明的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对本发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，这些都属于本发明的保护范围。因此，本发明专利的保护范围应以所附权利要求为准。

Claims (8)

  1. 一种磁共振快速参数成像的方法，所述方法包括：
    获取目标的欠采样磁共振信号；
    获取参数模型的先验信息；
    根据所述欠采样磁共振信号和所述先验信息进行目标图像序列重构，得到目标图像序列；
    将所述目标图像序列代入参数估计模型，得到目标参数，并生成参数图像。
  2. 根据权利要求1所述的方法，其特征在于，所述根据所述欠采样磁共振信号和所述先验信息进行目标图像序列重构，得到目标图像序列的步骤包括：
    通过补零傅立叶变换的方法初始化所述目标图像序列；
    通过基于非线性滤波稀疏重建的方法对所述初始化后的目标图像序列进行重构，进而求解得到目标图像序列。
  3. 根据权利要求2所述的方法，其特征在于，所述通过基于非线性滤波稀疏重建的方法对所述初始化后的目标图像序列进行重构，进而求解得到目标图像序列的步骤包括：
    基于稀疏采样理论，对所述目标图像序列在某一稀疏变化域内的系数做软阈值处理；
    根据所述参数估计模型，将目标图像序列每一空间位置对应的时间信号构造成Hankel矩阵；
    对所述Hankel矩阵求低秩近似，以使得所述目标图像序列每一空间位置所对应的时间信号呈指数衰减；
    将所述欠采样磁共振信号代入所述目标图像序列的K空间，以使得目标图像序列与所述欠采样磁共振信号一致；
    重复上述步骤直至算法收敛，得到目标图像序列。
  4. 根据权利要求3所述的方法，其特征在于，所述对所述Hankel矩阵求低秩近似，以使得所述目标图像序列每一空间位置所对应的时间信号呈指数衰减的步骤包括：
    通过奇异值分解对所述Hankel矩阵求低秩近似矩阵；
    沿所述低秩近似矩阵的反对角线方向作平均值运算增强信号，恢复所述低秩近似矩阵的Hankel结构；
    从所述低秩近似Hankel矩阵中提取所述目标图像序列每一空间位置对应的时间信号。
  5. 一种磁共振快速参数成像系统，所述系统包括：
    信号获取模块，用于对目标进行扫描得到欠采样磁共振信号；
    先验信息获取模块，用于获取参数模型的先验信息；
    图像重建模块，用于根据所述欠采样磁共振信号和所述先验信息进行目标图像序列重构，得到目标图像序列；
    参数估计模块，用于将所述目标图像序列代入参数估计模型，得到目标参数，并生成参数图像。
  6. 根据权利要求5所述的系统，其特征在于，所述图像重建模块包括：
    初始化单元，用于通过补零博立叶变换方法初始化所述目标图像序列；
    运算单元，用于通过基于非线性滤波的稀疏重建方法对所述初始化后的目标图像序列进行重构，进而求解得到所述目标图像序列。
  7. 根据权利要求6所述的系统，其特征在于，所述运算单元包括：
    稀疏约束单元，用于基于稀疏采样理论，对所述目标图像序列在某一稀疏变化域内的系数做软阈值处理；
    矩阵构造单元，用于根据所述参数估计模型，将目标图像序列每一空间位置对应的时间信号构造成Hankel矩阵；
    模型约束单元，用于对所述Hankel矩阵求低秩近似，以使得所述目标图像序列每一空间位置所对应的时间信号呈指数衰减；
    信号保真单元，用于将所述欠采样磁共振信号代入所述目标图像序列的 K空间，以使得目标图像序列与所述欠采样磁共振信号一致，所述运算单元通过重复执行所述稀疏约束单元、所述矩阵构造单元、所述模型约束单元和所述信号保真单元的运算，直到所求目标图像序列收敛，得到最终的目标图像序列。
  8. 根据权利要求7所述的系统，其特征在于，所述模型约束单元包括：
    低秩近似矩阵求解单元，用于通过奇异值分解对所述Hankel矩阵求低秩近似矩阵；
    Hankel矩阵恢复单元，用于沿所述低秩近似矩阵的反对角线方向作平均值运算增强信号，恢复所述低秩近似矩阵的Hankel结构；
    信号提取单元，用于从所述低秩近似Hankel矩阵中提取所述目标图像序列每一空间位置对应的时间信号。
PCT/CN2014/093168 2013-12-10 2014-12-05 磁共振快速参数成像方法和系统 WO2015085889A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US15/103,845 US10018698B2 (en) 2013-12-10 2014-12-05 Magnetic resonance rapid parameter imaging method and system

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201310671860.3A CN103654789B (zh) 2013-12-10 2013-12-10 磁共振快速参数成像方法和系统
CN201310671860.3 2013-12-10

Publications (1)

Publication Number Publication Date
WO2015085889A1 true WO2015085889A1 (zh) 2015-06-18

Family

ID=50294424

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2014/093168 WO2015085889A1 (zh) 2013-12-10 2014-12-05 磁共振快速参数成像方法和系统

Country Status (3)

Country Link
US (1) US10018698B2 (zh)
CN (1) CN103654789B (zh)
WO (1) WO2015085889A1 (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105976329A (zh) * 2016-04-27 2016-09-28 厦门大学 一种基于时域信号低秩的频谱恢复方法
CN109633502A (zh) * 2018-12-03 2019-04-16 深圳先进技术研究院 磁共振快速参数成像方法及装置
CN110946578A (zh) * 2018-09-26 2020-04-03 西门子医疗有限公司 自动选择成像参数的方法、相应的成像方法和成像设备
CN115830172A (zh) * 2023-02-20 2023-03-21 泰山学院 迭代细化压缩感知欠采样磁共振图像重构方法、重构系统

Families Citing this family (25)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103654789B (zh) * 2013-12-10 2015-12-30 深圳先进技术研究院 磁共振快速参数成像方法和系统
JP6072723B2 (ja) * 2014-04-21 2017-02-01 株式会社日立製作所 磁気共鳴イメージング装置、及び画像撮像方法
CN104027113B (zh) * 2014-07-01 2016-03-16 中国科学院武汉物理与数学研究所 一种基于先验知识和稀疏采样的肺部快速磁共振成像方法
CN104156994B (zh) * 2014-08-14 2017-03-22 厦门大学 一种压缩感知磁共振成像的重建方法
CN104569880B (zh) * 2014-12-31 2017-04-05 中国科学院深圳先进技术研究院 一种磁共振快速成像方法及系统
KR101667141B1 (ko) * 2015-03-11 2016-10-25 한국과학기술원 소멸필터를 이용한 고속 mr 영상 복원 알고리듬 개발
CN107016656B (zh) * 2017-04-01 2020-06-19 中国科学院光电技术研究所 一种基于压缩感知在图像重建中的小波稀疏基优化方法
CN107423543B (zh) * 2017-04-14 2020-04-07 厦门大学 一种超复数磁共振波谱的快速重建方法
CN107330953B (zh) * 2017-07-06 2020-06-26 桂林电子科技大学 一种基于非凸低秩的动态mri重建方法
US10853977B2 (en) * 2017-08-30 2020-12-01 Korea Advanced Institute Of Science And Technology Apparatus and method for reconstructing image using extended neural network
CN110161442B (zh) * 2018-02-12 2020-04-21 深圳先进技术研究院 磁共振参数成像方法、装置、医学设备及存储介质
CN110146836B (zh) * 2018-02-12 2020-04-21 深圳先进技术研究院 一种磁共振参数成像方法、装置、设备及存储介质
WO2019199960A1 (en) * 2018-04-11 2019-10-17 Cornell University Up-sampling of signals by analytic phase projection
CN108896942B (zh) * 2018-06-27 2019-11-08 中国石油大学(北京) 核磁共振数据压缩方法及装置
CN109165432B (zh) * 2018-08-09 2022-12-13 厦门理工学院 一种基于部分奇异值和的磁共振波谱重建方法
CN109247939B (zh) * 2018-08-29 2020-06-19 中国科学院武汉物理与数学研究所 一种自适应高欠采样超极化气体肺部动态mri重建方法
CN109171727B (zh) * 2018-09-20 2022-03-15 上海东软医疗科技有限公司 一种磁共振成像方法和装置
CN109615675B (zh) * 2018-12-04 2022-05-03 厦门大学 一种多通道磁共振成像的图像重建方法
CN109658468B (zh) * 2018-12-12 2022-12-13 深圳先进技术研究院 磁共振参数成像方法、装置、设备及存储介质
EP3905946A4 (en) * 2019-01-04 2022-08-31 University of Cincinnati SYSTEM AND METHOD FOR RECONSTRUCTING MAGNETIC RESONANCE IMAGES ACQUIRED WITH PARTIAL FOURIER ACQUISITION
CN111856364B (zh) * 2019-04-24 2023-03-28 深圳先进技术研究院 一种磁共振成像方法、装置、系统及存储介质
CN111856365B (zh) * 2019-04-24 2023-03-14 深圳先进技术研究院 磁共振成像方法、装置、系统及存储介质
CN110658484B (zh) * 2019-10-17 2021-07-20 东北大学 一种磁共振波谱重建方法及系统
CN113143243B (zh) * 2021-02-25 2023-08-22 厦门大学 一种基于谱分解的深度学习磁共振波谱重建方法
CN116738764B (zh) * 2023-08-08 2023-10-20 中国海洋大学 一种基于奇异值阈值算法的海洋平台舱室舒适度评估方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101975936A (zh) * 2010-09-03 2011-02-16 杭州电子科技大学 一种基于cs压缩感知技术的快速磁共振成像方法
WO2012127341A1 (en) * 2011-03-24 2012-09-27 Koninklijke Philips Electronics N.V. Compressed sensing mr image reconstruction using constraint from prior acquisition
CN103035017A (zh) * 2011-12-08 2013-04-10 中国科学院深圳先进技术研究院 磁共振参数成像方法和系统
CN103654789A (zh) * 2013-12-10 2014-03-26 深圳先进技术研究院 磁共振快速参数成像方法和系统
CN103705239A (zh) * 2013-12-05 2014-04-09 深圳先进技术研究院 磁共振参数成像方法和系统

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7330027B2 (en) * 2003-07-02 2008-02-12 Universitat Zurich System and method of magnetic resonance imaging for producing successive magnetic resonance images
WO2013067546A1 (en) * 2011-11-06 2013-05-10 Mayo Foundation For Medical Education And Research Method for calibration-free locally low-rank encouraging reconstruction of magnetic resonance images
CN103027682A (zh) * 2011-12-12 2013-04-10 深圳先进技术研究院 动态对比度增强磁共振成像方法和系统
CN103064046B (zh) * 2012-12-25 2015-04-15 深圳先进技术研究院 一种基于稀疏采样的核磁共振成像的图像处理方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101975936A (zh) * 2010-09-03 2011-02-16 杭州电子科技大学 一种基于cs压缩感知技术的快速磁共振成像方法
WO2012127341A1 (en) * 2011-03-24 2012-09-27 Koninklijke Philips Electronics N.V. Compressed sensing mr image reconstruction using constraint from prior acquisition
CN103035017A (zh) * 2011-12-08 2013-04-10 中国科学院深圳先进技术研究院 磁共振参数成像方法和系统
CN103705239A (zh) * 2013-12-05 2014-04-09 深圳先进技术研究院 磁共振参数成像方法和系统
CN103654789A (zh) * 2013-12-10 2014-03-26 深圳先进技术研究院 磁共振快速参数成像方法和系统

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
AN, YISHUO.: "Study of Parameter Mapping in Magnetic Resonance Imaging Based on Enhanced Compressed Sensing", MASTER DEGREE THESIS OF JILIN UNIVERSITY, 15 September 2013 (2013-09-15), pages 39 - 45 *

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105976329A (zh) * 2016-04-27 2016-09-28 厦门大学 一种基于时域信号低秩的频谱恢复方法
CN105976329B (zh) * 2016-04-27 2018-09-25 厦门大学 一种基于时域信号低秩的频谱恢复方法
CN110946578A (zh) * 2018-09-26 2020-04-03 西门子医疗有限公司 自动选择成像参数的方法、相应的成像方法和成像设备
CN110946578B (zh) * 2018-09-26 2023-10-31 西门子医疗有限公司 自动选择成像参数的方法、相应的成像方法和成像设备
CN109633502A (zh) * 2018-12-03 2019-04-16 深圳先进技术研究院 磁共振快速参数成像方法及装置
CN115830172A (zh) * 2023-02-20 2023-03-21 泰山学院 迭代细化压缩感知欠采样磁共振图像重构方法、重构系统

Also Published As

Publication number Publication date
CN103654789A (zh) 2014-03-26
US20160370444A1 (en) 2016-12-22
CN103654789B (zh) 2015-12-30
US10018698B2 (en) 2018-07-10

Similar Documents

Publication Publication Date Title
WO2015085889A1 (zh) 磁共振快速参数成像方法和系统
Zhao et al. Maximum likelihood reconstruction for magnetic resonance fingerprinting
US10671939B2 (en) System, method and computer-accessible medium for learning an optimized variational network for medical image reconstruction
Sumpf et al. Fast T2 mapping with improved accuracy using undersampled spin-echo MRI and model-based reconstructions with a generating function
US8587307B2 (en) Systems and methods for accelerating the acquisition and reconstruction of magnetic resonance images with randomly undersampled and uniformly undersampled data
Hong et al. Compressed sensing MRI with singular value decomposition-based sparsity basis
US8823374B2 (en) System for accelerated MR image reconstruction
CN109615675B (zh) 一种多通道磁共振成像的图像重建方法
US11830111B2 (en) Magnetic resonance imaging method, device, medical device and storage medium
Stolk et al. Understanding the Combined Effect of ${k} $-Space Undersampling and Transient States Excitation in MR Fingerprinting Reconstructions
Majumdar Improved dynamic MRI reconstruction by exploiting sparsity and rank-deficiency
WO2014162300A1 (en) Acceleration of low-rank mri data acquisition
Balachandrasekaran et al. Recovery of damped exponentials using structured low rank matrix completion
Wen et al. A conditional normalizing flow for accelerated multi-coil MR imaging
CN105759233B (zh) 一种快速化学交换饱和转移成像方法和系统
US20100027860A1 (en) Method for determining phase-corrected amplitudes in nmr relaxometric imaging
Zong et al. Fast reconstruction of highly undersampled MR images using one and two dimensional principal component analysis
Huang et al. Accelerating cardiac diffusion tensor imaging combining local low-rank and 3D TV constraint
CN107727678A (zh) 一种核磁共振弛豫高低本征模态耦合方法
Hafftka et al. N-dimensional tensor completion for nuclear magnetic resonance relaxometry
Hu et al. SPICE: Self-supervised learning for MRI with automatic coil sensitivity estimation
Grosser et al. Efficient optimization of mri sampling patterns using the bayesian fisher information matrix
Jelena Comparison of algorithms for compressed sensing of magnetic resonance images
Ye et al. Patch-based directional redundant wavelets in compressed sensing parallel magnetic resonance imaging with radial sampling trajectory
Chaari Bayesian sparse regularized reconstruction in parallel MRI with sensitivity matrix imprecision

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 14868907

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

WWE Wipo information: entry into national phase

Ref document number: 15103845

Country of ref document: US

32PN Ep: public notification in the ep bulletin as address of the adressee cannot be established

Free format text: NOTING OF LOSS OF RIGHTS PURSUANT TO RULE 112(1) EPC (EPO FORM 1205 DATED 29.09.2016)

122 Ep: pct application non-entry in european phase

Ref document number: 14868907

Country of ref document: EP

Kind code of ref document: A1