Abstract:
Method and Apparatus for eliminating motion artifacts in magnetic resonance imaging are disclosed according to the present invention. The present invention relates to magnetic resonance imaging field. The method for eliminating motion artifacts in magnetic resonance imaging according to the present invention utilizes the concept of iterative approximation to control the difference between the data lines in the K-space caused by motions and allow the common features between the data lines to be remained, such that the motion artifacts in the reconstructed image are restrained and the motion artifacts caused under various circumstances are well restrained. Accordingly, the quality of the magnetic resonance imaging is improved.

Description:
FIELD OF THE INVENTION  
       [0001]    The present invention relates to Magnetic Resonance Imaging (MRI) field, and more specifically, to method and apparatus for eliminating motion artifact in MRI. 
       BACKGROUND OF THE INVENTION  
       [0002]    In MRI technology, a space stored with magnetic resonance signals is called K-space, which is a Fourier transform space. MRI signals, which are obtained via MRI scan, are filled in the K-space. A magnetic resonance image can be obtained via an inverse Fourier transform of the K-space data. 
         [0003]    Motion interference has always been a headache to the MRI imaging quality. This is even more the case for scans of certain parts of body, such as heart, abdomen. Body breathing or heart beating may result in motion artifacts and thereby affecting the MRI imaging quality. 
         [0004]    Currently, there are some approaches for eliminating motion artifacts in MRI, such as PROPELLER approach, which may to some extent mitigate the motion interference to the imaging and remove the motion artifacts. However, these approaches relates to motion of rigid bodies in a plane. For scans of some parts of the body, motions of non-rigid bodies or out-of-plane motions are the main factors affecting the imaging quality. However, these approaches are unable to solve these problems. 
         [0005]    There are also some other approaches, such as flow compensation technique, which may be used to compensate the artifacts caused by the blood flow. However, these approaches are limited to applications under certain special circumstances (i.e., artifacts caused by the blood flow). 
         [0006]    Overall, the existing approaches for eliminating the motion artifacts in magnetic resonance imaging are limited. 
         [0007]    Accordingly, there is a need for a new method and apparatus for eliminating motion artifacts in magnetic resonance imaging so as to restrain the motion artifacts caused under various circumstances and thus improving the quality of magnetic resonance imaging. 
       SUMMARY OF THE INVENTION  
       [0008]    The problem to be solved by the present invention is that the existing method for eliminating motion artifacts in magnetic resonance imaging is not able to effectively restrain the motion artifacts caused under various circumstances. 
         [0009]    In order to solve the foregoing problem, there is provided a method for eliminating motion artifacts in magnetic resonance imaging according to the present invention. The method includes the following steps: 
         [0010]    1) obtaining fully-sampled initial K-space data corresponding to a plurality of channels; 
         [0011]    2) calculating combination coefficients based on the fully-sampled initial K-space data and determining K 1 -space data based on the combination coefficients; 
         [0012]    3) obtaining K 2 -space data based on the initial K-space data and the K 1 -space data, obtaining a S 2  image by transforming the K 2 -space data to an image domain, obtaining a S 2′  image by extracting part of the errors caused during the determination of K 1 -space data; 
         [0013]    4) obtaining K 3 -space data by transforming the S 2′  image to K-space domain, obtaining K 4 -space data by subtracting the K 3 -space data from the K 1 -space data; 
         [0014]    5) determining whether to perform a next iteration operation; making the K 4 -space data new initial K-space data if a next iteration operation is determined to be performed, and returning to step 2); 
         [0015]    6) obtaining a final image by transforming the K 4 -space data to the image domain. 
         [0016]    Optionally, the step of obtaining K 2 -space data based on the initial K-space data and the K 1 -space data comprises obtaining K 2 -space data by subtracting the initial K-space data from the K 1 -space data. 
         [0017]    Optionally, the step of obtaining K 4 -space data based on the K 3 -space data and the K 1 -space data comprises obtaining K 4 -space data by subtracting the K 3 -space data from the K 1 -space data. 
         [0018]    Optionally, the step of extracting part of the errors caused during the calculation of the K 1 -space includes replacing the values of pixels in the S 2  image with zero so as to obtain the S 2′  image wherein the values of pixels in the S 2  image to be replaced are greater than a first threshold. 
         [0019]    Optionally, the step of extracting part of the errors caused during the step of obtaining the K 1 -space data includes replacing the values of pixels in the S 2  image with zero, and performing mean filtering so as to obtain the S 2 ′ image; wherein the values of pixels in the S 2  image to be replaced are greater than a first threshold. 
         [0020]    Optionally, step 1) can be implemented by directly full-sampling the initial K-spaces for a plurality of channels; or under-sampling K-spaces for a plurality of channels and obtaining the fully-sampled initial K-spaces using a corresponding method. 
         [0021]    Optionally, the corresponding method may be a parallel acquisition reconstruction method or semi-Fourier reconstruction method. 
         [0022]    Optionally, step 2) includes the following sub-steps: 
         [0023]    21) utilizing the initial K-space data to obtain the combination coefficients; 
         [0024]    22) for each channel, using the same rule to divide the initial K-space data into a plurality of sub-K-space data, wherein the sub-K-space data are matrix structures which have the same size as the initial K-space data, and each sub-K-space data comprises parts of the sampled data derived from the initial K-space data and unsampled data to be filled; 
         [0025]    23) utilizing the combination coefficients and sampled data in the sub-K-space data to fill the unsampled datain each sub-K-space; 
         [0026]    24) combining the data in each filled sub-K-space and obtaining the K 1 -space data. 
         [0027]    Optionally, in step 21), the combination coefficients refer to the combination coefficients used in the parallel acquisition reconstruction method. 
         [0028]    Optionally, in step 23), the data to be filled in each sub-K-space is filled using the parallel acquisition reconstruction method. 
         [0029]    Optionally, the parallel acquisition reconstruction method may be GRAPPA method or SPIRiT method. 
         [0030]    Optionally, the step 22) comprises dividing the initial K-space data for each channel into two sub-K-space data, wherein the sub-K-space data are matrix structures which have the same size as the initial K-space. The odd data lines in one of the two sub-K-space data correspond to the odd data lines in the initial K-space data, and the even data lines are the unsampled data to be filled out. The even data lines in the other sub-K-space data correspond to the even data lines in the initial K-space data, and the odd data lines are the unsampled data to be filled out. 
         [0031]    Optionally, step 24) comprises calculating an average of each filled sub-K-space data, and using the average as the K 1 -space data. 
         [0032]    Optionally, the first threshold is the average of the values of the pixels in the S 2  image obtained via the current iteration. 
         [0033]    Optionally, determining whether to perform a next iteration operation in step 5) comprises: pre-setting the number of iteration cycle, proceeding to perform iteration calculation and making the K 4 -space data new initial K-space data and returning to step 2) if the preset number is not reached; stopping the iteration calculation and performing step 6) if the pre-set number is reached. 
         [0034]    Optionally, determining whether to perform a next iteration operation in step 5) comprises setting a second threshold; proceeding the iteration calculation and making the K 4 -space data new initial K-space data and returning to step 2) if a standard deviation of the values of the pixels in the S 2  image obtained via the current iteration is no less than the second threshold; stopping the iteration calculation and performing step 6) if the standard deviation of the values of the pixels in the S 2  image obtained via the current iteration is less than the second threshold. 
         [0035]    Optionally, the first threshold is 30% of the standard deviation of the values of the pixels in the S 2  image obtained via the first iteration. 
         [0036]    Optionally, the K-space data may be transformed to the image domain via inverse Fourier Transformation and the image may be transformed to the K-space domain via the Fourier Transformation. 
         [0037]    There is further provided an apparatus for eliminating motion artifacts in magnetic resonance imaging according to the present invention. The apparatus includes a sampling unit, a K 1 -space obtaining unit, a S 2′  image obtaining unit, a K 4  space obtaining unit, an iteration determination unit and an image obtaining unit. 
         [0038]    The sampling unit is configured to obtain fully-sampled initial K-space data corresponding to a plurality of channels. 
         [0039]    The K 1 -space obtaining unit coupled to the sampling unit is configured to calculate combination coefficients based on the fully-sampled initial K-space data and obtain K 1 -space data based on the combination coefficients. 
         [0040]    The S 2′  image obtaining unit coupled to the K 1 -space obtaining unit is configured to obtain a K 2 -space based on the initial K-space data and the K 1 -space data, obtain a S 2  image by transforming the K 2 -space data to an image domain, and obtain a S 2′  image by extracting part of the errors caused during the step of obtaining of the K 1 -space data. 
         [0041]    The K 4  space obtaining unit coupled to the S 2′  image obtaining unit is configured to obtain K 3 -space data by transforming the S 2′  image to K-space domain and obtain K 4 -space data based on the K 3 -space data and the K 1 -space data. 
         [0042]    The iteration determination unit coupled to the K 1 -space obtaining unit and the K 4 -space obtaining unit is configured to determine whether to perform a next iteration operation. 
         [0043]    The image obtaining unit coupled to the iteration determination unit is configured to transform the K 4 -space data to the image domain to obtain a final image. 
         [0044]    The present invention has the following advantages over the prior art. That is, the present invention utilizes the concept of the iterative approximation to control or even eliminate the difference between the data lines in the K-space caused by the motion and allow the common features between the data lines to be remained, such that the motion artifacts in the reconstructed image are restrained and the motion artifacts caused under various circumstances are well restrained. Accordingly, the quality of the magnetic resonance imaging is improved. 
     
    
     
       BRIEF DESCRIPTION OF THE DRAWINGS  
         [0045]      FIG. 1  illustrates a flowchart of a method for eliminating motion artifacts in magnetic resonance imaging according to the present invention; 
           [0046]      FIG. 2  illustrates a schematic of GRAPPA approach; 
           [0047]      FIG. 3  illustrates a schematic of a detailed process of step S 02  in the method for eliminating motion artifacts in magnetic resonance imaging according to the present invention; 
           [0048]      FIG. 4  illustrates a contrast between the image reconstructed without using the method of the present application and the image reconstructed using the method according to the present invention; and 
           [0049]      FIG. 5  illustrates a schematic of an apparatus for eliminating motion artifacts in magnetic resonance imaging according to the present invention. 
       
    
    
     DETAILED DESCRIPTION OF THE INVENTION  
       [0050]    The purposes, features and advantages concerning the present invention will become more readily appreciated by reference to the following detailed description of the embodiments, when taken in conjunction with the accompanying drawings. 
         [0051]    As shown in  FIG. 1 , the method for eliminating motion artifacts in magnetic resonance imaging includes the following steps. 
         [0052]    Step S 01  is performed by obtaining fully-sampled initial K-spaces for a plurality of channels. The fully-sampled initial K-space can be obtained by the following ways: directly fully-sampling the initial K-spaces for the plurality of channels; or under-sampling the K-spaces for the plurality of channels and obtaining the fully-sampled initial K-spaces using a corresponding method. 
         [0053]    The corresponding method can be an existing method for filling the under-sampled K-space data and can also be other methods. 
         [0054]    Specifically, the corresponding method may be, but not limited to, a parallel acquisition reconstruction method or semi-Fourier reconstruction method. 
         [0055]    In the present embodiment, the fully-sampled initial K-space data is denoted as “S 0 ” for brevity hereinafter. 
         [0056]    Step S 02  is performed by determining a combination coefficient C 0  based on the fully-sampled initial K-space S 0  and determining a K 1 -space based on the combination coefficient. 
         [0057]    Specifically, a parallel acquisition reconstruction method can be used to compute the combination coefficient and obtaining the K 1 -space data. The parallel acquisition reconstruction method may be, but not limited to, GRAPPA (Generalized Auto-calibrating Partially Parallel Acquisitions) method or SPIRiT method. 
         [0058]    The parallel acquisition reconstruction method involves filling the K-space data in the existing process of eliminating the motion artifacts in magnetic resonance imaging, which is known by a skilled person. GRAPPA method is taken below as an example for illustration of the parallel acquisition reconstruction method. 
         [0059]    K-space can be divided into a phase encoding (PE) direction and a read-out (RO) direction. As shown in  FIG. 2 , the horizontal direction represents the phase encoding direction and the vertical direction represents the number of channels. The direction perpendicular to the paper is the read-out direction (a dot in  FIG. 2  can be regarded as a data line in the K-space). The black solid dots  101  denote the actually sampled K-space data. The white hollow dots  102  denote un-sampled data which needs to be filled. The grey solid dots  103  denote data which are to be fitted. In GRAPPA method, any white hollow dot  102  in the figure can be represented by linear addition of the adjacent black solid dots  101 . This equals to the combination of data from a plurality of coils. As shown in  FIG. 2 , the combination coefficient n ij  which corresponds to coil i, position j can be determined by the black solid dots  101  best fitting the grey dots  103 . In GRAPPA method, the data in the fully-sampled area in the K-space can be used as the calibration data (i.e., the block solid dots  101  and the grey solid dots  103  in  FIG. 2  can be used as the calibration data). By calculating the fitting between the data lines in the calibration data, all the combination coefficients can be obtained. After the combination coefficients are determined, coils can be combined based on the combination coefficients and other sampled data, and the un-sampled data (to be filled) in the K-space, as denoted by white hollow dots  102 , are filled. 
         [0060]    Specifically, step S 02  may include the following sub-steps. 
         [0061]    Step S 021  is performed by utilizing the initial K-space data as the calibration data to determine the combination coefficients. 
         [0062]    Step S 022  is performed by dividing the initial K-space into a plurality of sub-K-spaces (also referred to as divided K-spaces) for each channel with a same rule and obtaining a plurality of sub-K-spaces. Each sub-K-space is a matrix structure which has the same size as the initial K-space. Each sub-K-space includes parts of the sampled data points and the unsampled data to be filled. 
         [0063]    Preferably, the initial K-space S 0  for each channels are pre-divided in the following ways. 
         [0064]    As shown in  FIG. 3 , the two-dimensional K-space may be divided into a phase encoding (PE) direction and a read-out (RO) direction. A plurality of data lines arranged sequentially along the PE direction forms a complete K-space data. In the present embodiment, the initial K-space S 0  for each channel is divided into two sub-K-spaces S 01  and S 02 . The sub-K-spaces S 01  and S 02  are matrixes which has the same size as the initial K-space S 0 . The odd data lines in one of the sub-K-spaces (S 01 ) correspond to the odd data lines in the initial K-space S 0 , and the even data lines are the unsampled data to be filled. The even data lines in the other sub-K-space (S 02 ) correspond to the even data lines in the initial K-space S 0 , and the odd data lines are the unsampled data to be filled. 
         [0065]    Step S 023  is performed by utilizing the combination coefficients and sampled data in the sub-K-spaces to fill the unsampled data to be filled in each sub-K-space. 
         [0066]    Preferably, as shown in  FIG. 3 , the even data lines in the sub-K-space S 01  is filled out using the parallel acquisition reconstruction method. The combination coefficient and the even data lines in the sub-K-space S 02  are used to fill out the odd data line in the K-space S 02  in accordance with the parallel acquisition reconstruction method. As such, two fully-sampled sub-K-spaces S 01     —     full  and S 01     —     full  are obtained, respectively. 
         [0067]    Step S 03  is performed by combining each filled sub-K-space data and obtaining the K 1 -space data. 
         [0068]    Preferably, as shown in  FIG. 3 , an average of the two filled sub-K-space data S 01     —     full  and S 01      —     full  is calculated. The average serves as K 1 -space data S 1 , i.e., S 1 =(S 01     —     full +S 02     —     full )/2. 
         [0069]    Step S 03  is performed by obtaining a K 2 -space data S 2  by subtracting the initial K-space data S 0  from the K 1 -space data S 1  (i.e., S 2 =S 1 −S 0 ), obtaining a S 2  image ImS 2  by transforming the K 2 -space data S 2  to an image domain, obtaining a S 2′  image ImS 2′  by extracting part of the errors introduced during the calculation of K 1 -space. 
         [0070]    The method of extracting part of the errors introduced by algorithm during the calculation of the K 1 -space may include, but not limited to, the following two approaches: replacing the values of pixels in the S 2  image, which are greater than a first threshold, with zero so as to obtain the S 2′  image; or replacing the values of pixels in the S 2  image, which are greater than a first threshold, with zero, and performing mean filtering so as to obtain the S 2′  image. 
         [0071]    The method of transforming the K 2 -space data to the image domain may be, but not limited to, an inverse Fourier transformation. 
         [0072]    Preferably, the first threshold is the average of the values of the pixels in the S 2  image. 
         [0073]    Step S 04  is performed by transforming the S 2′  image ImS 2′  to K-space domain to obtain K 3 -space data S 3  and subtracting the K 3 -space data S 3  from the K 1 -space data S 1  to obtain K 4 -space data, i.e., S 4 =S 1 −S 3 . 
         [0074]    The method of transforming the S 2′  image ImS 2′  to K-space domain may be, but not limited to, a Fourier transformation. 
         [0075]    Step S 05  is performed by determining whether to perform iteration. If it is determined to perform the iteration, the K 4 -space data is treated as the new initial K-space data, and step S 02  is returned to. If it is determined not to perform the iteration, step S 06  is performed. 
         [0076]    The determination of whether to perform the iteration can be implemented with one of the following two methods. 
         [0077]    1) pre-sett the number of iteration cycle. If the pre-set number is reached, stop the iteration and perform step S 06 ; 
         [0078]    2) set a second threshold. If a standard deviation of the values of the pixels in the S 2  image obtained via the current iteration is less than the second threshold, stop the iteration and perform step S 06 . 
         [0079]    The second threshold may be a specific value or may be determined by the standard deviation of the values of the pixels in the S 2  image obtained via the first iteration. 
         [0080]    Preferably, the second threshold is 30% of the standard deviation of the values of the pixels in the S 2  image obtained via the first iteration. That is, the iteration is stopped when the standard deviation of the values of the pixels in the S 2  image obtained via the current iteration is less than 30% of the standard deviation of the values of the pixels in the S2 image obtained via the first iteration. 
         [0081]    Step S 06  is performed by transforming the K4-space data S4 to the image domain to obtain a resulting magnetic resonance image. 
         [0082]    Specifically, the method of transforming the K4-space data S4 to the image domain may be, but not limited to, an inverse Fourier transformation. 
         [0083]    Advantages of the present invention and why such advantages can be achieved are detailed below. 
         [0084]    In step S 01 , the obtained initial K-space data can be expressed as S 0 (k)=A(k)+B(k)+n, where K denotes the K th  data line in the K-space; A(k) denotes the signal when the object stands still, which is the desired signal in magnetic resonance image; B(k) denotes the signal generated due to motion interference; n denotes a noise signal. Since the influence caused by motion is changing every moment, the A(k) of the sampled data line of corresponding positions at different channels are consistent with each other and B(K) are different from each other. 
         [0085]    Thus, there is need to a method for obtaining a new K-space where the signals generated due to motion interference can be reduced and thereby motion interference can be restrained. 
         [0086]    Step S 02  is performed so as to obtain the K1-space data S 1 (k)=A(k)+B′(k)+n+C(k). Since in step S 02  a data line in the K1-space is determined by best fitting its neighboring data line, thus B′(k)=Σ i=1 ″ α i B(k+Δ i ). That is, B′(k) represents the sum of convolutions of the motion of the neighboring data line. In magnetic resonance imaging, it can be approximately regarded that the difference made due to the motion of the neighboring data line can be counteracted by each other. Therefore, the motion difference reflected by B′(k) is less than B(k). C(k) denotes that the error introduced during the step of obtaining the K1-space. 
         [0087]    In step S 03 , K 2 -space data is obtained by subtracting the initial K-space data from the K 1 -space data. That is, the K 2 -space data is expressed as S 2 (k)=B′(k)+C(k)−B(k). S 2  image is then obtained by transforming the K 2 -space data to the image domain. By replacing the values of pixels in the S 2  image which are greater than the first threshold with zero (and performing mean filtering), the error C(k) introduced during the step of obtaining the K 1 -space can be extracted. 
         [0088]    Therefore, in step S 04 , part of data in the K4-space is obtained S 4 (k)=A(k)+B′(k)+n. The other part of the data is still S 4 ′(k)=A(k)+B′(k)+n+C(k). Multiple subsequent iterations can be performed to eliminate the error C(k) from S 4 (k). 
         [0089]    With the concept of iterative approximation, steps S 02 ˜S 04  are iterated such that the resulting K 4 -space data S 4 (k) will best approximate A(k). Consequently, the K-space data with motion restrained can be obtained. K-space data which are obtained via iteration are then used to reconstruct an image where motion artifacts are removed. Thus, the quality of the magnetic resonance imaging is improved. 
         [0090]      FIG. 4  ( a ) is an image obtained by a conventional image reconstruction method. That is, the image is obtained by transforming the initial K-space data obtained via sampling method directly to the image domain without using any method for eliminating motion artifacts.  FIG. 4(   b ) is an image obtained by the image reconstruction method for eliminating motion artifacts in accordance with the present invention. As shown in  FIG. 4(   a ) and  FIG. 4(   b ), the image reconstruction method of the present invention has great advantages in eliminating artifacts caused by the motion indicated by the white arrow. 
         [0091]      FIG. 4(   c ) is an image obtained by a conventional image reconstruction method in a heart beating scan.  FIG. 4(   d ) is an image reconstructed by the image reconstruction method for eliminating motion artifacts in this scan image in accordance with the present invention. As shown in  FIG. 4(   c ) and  FIG. 4(   d ), the present invention has remarkable effect on eliminating artifacts caused by heart beating. 
         [0092]      FIG. 4  indicates that the image reconstruction method according to the present invention has sound effect in restraining the motion artifacts under various circumstances. 
         [0093]    There is further provided an apparatus for eliminating motion artifacts in magnetic resonance imaging according to the present invention. As shown in  FIG. 5 , the apparatus  400  for eliminating motion artifacts in magnetic resonance imaging according to the present invention includes a sampling unit  401 , a K 1 -space obtaining unit  402 , the S 2′  image obtaining unit  403 , the K 4 -space obtaining unit  404 , the iteration determination unit  405 , and an image obtaining unit  406 . 
         [0094]    The sampling unit  401  is configured to obtain fully-sampled initial K-spaces corresponding to a plurality of channels. 
         [0095]    The K 1 -space obtaining unit  402  coupled to the sampling unit  401  is configured to calculate combination coefficients based on the fully-sampled initial K-space and obtain a K 1 -space based on the combination coefficients. 
         [0096]    The S 2′  image obtaining unit  403  coupled to the K 1 -space obtaining unit  402  is configured to obtain a K 2 -space by subtracting the initial K-space data from the K 1 -space data, obtain a S 2  image by transforming the K 2 -space to an image domain, and obtain a S 2′  image by extracting part of the errors introduced during the step of obtaining of the K 1 -space. 
         [0097]    The K 4 -space obtaining unit  404  coupled to the S 2′  image obtaining unit  403  is configured to obtain K 3 -space data by transforming the S 2′  image to K-space domain, and obtain K 4 -space data by subtracting the K 3 -space data from the K 1 -space data. 
         [0098]    The iteration determination unit  405  coupled to the K 1 -space obtaining unit  402  and the K 4 -space obtaining unit  404  is configured to determine whether to perform a next iteration. 
         [0099]    The image obtaining unit  406  coupled to the iteration determination unit  405  is configured to transform the K 4 -space data to the image domain to obtain a resulting image. 
         [0100]    Although the present invention has disclosed the above content, the present invention is not limited by such disclosure. Any skilled person in the art is able to make various changes and modifications without departing from the spirit of the present invention. Therefore, any content without departing from the spirit of the present invention and any simple variation or equivalents or modification made to the foregoing embodiments in accordance with the technical substance of the present invention shall all fall within the scope of the present invention.