Magnetic resonance imaging device and sensitivity distribution calculation program
Patent Information
- Application Number
- JP2021085645
- Authority / Receiving Office
- JP · JP
- Patent Type
- Patents
- Current Assignee / Owner
- Filing Date
- 2021-05-20
- Publication Date
- 2025-01-06
- Estimated Expiration
- 2041-05-20
AI Technical Summary
【0011】 本発明によれば、k空間で感度分布を求める際に、感度分布をより安定に求めることができる。また、安定に感度分布を求められることを担保しつつ、カーネルサイズを大きくでき、感度分布の高周波成分も表現できるようになる。
Smart Images

Figure 0007608268000006 
Figure 0007608268000007 
Figure 0007608268000008
Abstract
Description
[Technical field]
[0001] The present invention relates to a magnetic resonance imaging apparatus (hereinafter referred to as an MRI apparatus), and in particular to a technique for calculating sensitivity distribution in k-space. [Background technology]
[0002] There is a technique in which data is thinned out and received by a receiver coil with multiple channels, and an image is reconstructed using information on the sensitivity distribution of the coil (Non-Patent Document 1, Non-Patent Document 2). To perform this image reconstruction, it is necessary to accurately calculate the sensitivity distribution of the coil, but there is a problem that it is difficult to accurately calculate the boundary part of the subject. This is because, for stable reconstruction, it is necessary to calculate the sensitivity distribution by slightly extrapolating to the area where the subject does not exist, but this extrapolation is difficult and because the received data that is the basis for sensitivity calculation becomes small at the boundary of the subject, resulting in a low SNR and large errors. To solve this problem, there is a technique in which the sensitivity distribution is obtained in the k-space (frequency domain) (Patent Document 1, Patent Document 2).
[0003] In the calculation to obtain the sensitivity distribution in k-space, the image of each channel and an image (reference image) that does not depend on the sensitivity distribution are used, and the frequency components of the sensitivity distribution convolved in k-space with the reference image are deconvolved in k-space to obtain the sensitivity. In such a calculation in k-space, the sensitivity distribution in the boundary region where the signal is small is obtained using the surrounding regions where the signal is large, which prevents the sensitivity distribution from becoming extremely inaccurate in the boundary region. In addition, based on the assumption that the sensitivity distribution does not change significantly, the size of the kernel obtained by deconvolution can be made small, which limits large changes and makes it possible to stably obtain a natural sensitivity distribution. [Prior art documents] [Patent documents]
[0004] [Patent Document 1] Patent No. 6820876 specification [Patent Document 2] U.S. Pat. No. 10,353,023 [Non-patent literature]
[0005] [Non-Patent Document 1] Klaas P.Pruessmann, et al. Magnetic Resonance in Medicine 42:952-962(1999) [Non-Patent Document 2] Mark A. Griswold, et al. Magnetic Resonance in Medicine 47:1202-1210(2002) Summary of the Invention [Problem to be solved by the invention]
[0006] In the technology disclosed in Patent Documents 1 and 2 for calculating sensitivity distribution in k-space, it is necessary to take measures against spectral leakage. If the measures against spectral leakage are insufficient, the sensitivity distribution becomes unstable, for example, vibrating in areas with small signals. In order to reduce the effect of spectral leakage, there is a method of zero-filling the data for calculating the sensitivity distribution in the image space and expanding it to a specified FOV, but when the FOV is expanded by zero-filling, the degree of freedom of the solution increases, which causes a problem of large vibration in non-signal areas.
[0007] By reducing the kernel size during convolution, vibrations can be suppressed; however, reducing the kernel size narrows the range of frequencies that can be handled, which creates the problem that high frequencies cannot be expressed. [Means for solving the problem]
[0008] In order to solve the above problem, the present invention calculates a sensitivity distribution while suppressing spectral leakage by expanding data, which is the basis for calculating the sensitivity distribution, like a mirror image.
[0009] That is, the MRI apparatus of the present invention includes a receiver coil having multiple channels, a measurement unit that measures the nuclear magnetic resonance signal of the subject for each channel of the receiver coil, and an image calculation unit that creates an image of the subject using the sensitivity distribution for each channel of the receiver coil and a channel image generated from the nuclear magnetic resonance signal for each channel measured by the measurement unit. The image calculation unit includes a sensitivity distribution calculation unit that calculates the sensitivity distribution for each channel of the receiver coil in k-space using a sensitivity distribution image obtained for the sensitivity distribution, and the sensitivity distribution calculation unit has an unfolded image creation unit that performs one or more inversion processes on at least a part of the sensitivity distribution image to generate an unfolded image combining the generated inverted image and the original sensitivity distribution image, and calculates the sensitivity distribution using the unfolded image.
[0010] In the present invention, the "sensitivity distribution" includes not only the sensitivity distribution in image space and k space (frequency space), but also information equivalent to the sensitivity distribution such as the weighting coefficient in k space used in the image reconstruction method (GRAPPA method) in Non-Patent Document 2, and the present invention can be applied to obtaining such sensitivity distribution information. Effect of the Invention
[0011] According to the present invention, when a sensitivity distribution is obtained in k-space, the sensitivity distribution can be obtained more stably. In addition, while ensuring that the sensitivity distribution can be obtained stably, the kernel size can be increased, and high frequency components of the sensitivity distribution can also be expressed. [Brief description of the drawings]
[0012] [Figure 1] Functional block diagram of an MRI device. [Diagram 2] 1 is a diagram showing the overall configuration of an MRI apparatus according to an embodiment. [Diagram 3] 4A to 4C are diagrams showing a photographing procedure according to the first embodiment; [Figure 4] FIG. 4 is a diagram showing an example of a sensitivity distribution calculation pulse sequence in the first embodiment. [Diagram 5] FIG. 4 is a diagram showing a procedure for calculating sensitivity distribution. [Figure 6] FIG. 4 is a functional block diagram of a sensitivity distribution calculation unit. [Figure 7] FIG. 6 is a diagram showing details of the slice extraction process in the flow of FIG. 5. [Figure 8] (A) A diagram showing details of the process of unfolding into a mirror image, and (B) a diagram showing details of the process of cutting out from the mirror image. [Figure 9] FIG. 4 is a diagram showing a process according to the first embodiment. [Figure 10] 5A to 5C are diagrams showing the effects of the first embodiment. [Figure 11] FIG. 11 is a diagram for explaining a process of Modification 1 of the first embodiment. [Figure 12] 13A to 13C are diagrams illustrating generation of an expanded image according to the third modification of the first embodiment. [Figure 13] 13A to 13C are diagrams illustrating generation of an exfoliated image according to the fourth modification of the first embodiment. [Figure 14] 13A to 13C are diagrams illustrating generation of an expanded image according to the second embodiment. [Figure 15] FIG. 13 is a diagram showing a process according to a third embodiment. [Figure 16] FIG. 13 is a diagram showing an example of a user interface. DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS
[0013] First, an embodiment of an MRI apparatus common to each embodiment of the present invention described later will be described.
[0014] 1, the MRI apparatus of this embodiment includes a measurement unit 10 that acquires nuclear magnetic resonance signals, a calculation unit (image calculation unit) 20 that performs calculations such as image reconstruction using the nuclear magnetic resonance signals (hereinafter referred to as MR signals) acquired by the measurement unit, and a control unit 30 that controls the operation of each unit constituting the MRI apparatus, such as the measurement unit 10 and the calculation unit 20. The control unit 30 may be a part of the function of the calculation unit 20. Furthermore, the MRI apparatus may include a user interface unit 40 that includes a device for a user to input commands such as starting / ending imaging and setting imaging conditions to the measurement unit 10 and the control unit 30, and to present the processing results of the calculation unit 20 to the user.
[0015] The measurement unit 10 has a configuration similar to that of a general MRI device, and includes a static magnetic field generating unit 11 that generates a uniform magnetic field in the space in which the subject is placed, a gradient magnetic field generating unit 12 that applies a magnetic field gradient to the static magnetic field, a transmitting unit 13 that irradiates radio frequency magnetic field pulses (hereinafter referred to as RF pulses) that excite the nuclei of atoms that constitute the subject's tissue, typically protons, and a receiving unit 14 that receives nuclear magnetic resonance signals generated from the subject.
[0016] Specifically, as shown in FIG. 2, the static magnetic field generating unit 11 includes a static magnetic field coil 102, a shim coil 104 for correcting the magnetic field uniformity, and a power supply unit 113 for the same; the gradient magnetic field generating unit 12 includes a gradient magnetic field coil 103 and a power supply unit 113 for the same; the transmitting unit 13 includes a transmitting coil 105 and a transmitter 107 for sending an RF signal to the transmitting coil 105; and the receiving unit includes a receiving coil 106 and a receiver 108 for detecting the MR signal received by the receiving coil 106.
[0017] The static magnetic field generating unit 11 can use a static magnetic field coil 102 of a normal conducting magnet type or a superconducting magnet type, or a permanent magnet. Depending on the direction of the magnetic field, there are a vertical magnetic field type, a horizontal magnetic field type, and a type in which the magnetic field direction is inclined relative to the horizontal direction, and any of these may be used.
[0018] 2, the gradient magnetic field coil 103 is composed of gradient magnetic field coils in three mutually orthogonal axial directions, each of which is connected to a respective power supply unit 112. By applying a gradient magnetic field in each axis, it is possible to select an excitation region of the subject and to add position information to the MR signal.
[0019] The transmission coil 105 is generally fixed in the static magnetic field space together with the shim coil 104 and the gradient magnetic field coil 103 .
[0020] Although the receiving coil 106 is shown as one block in Fig. 2, it is composed of receiving coils with multiple channels and has a function of receiving signals for each channel. As a receiving coil with multiple channels, a known multiple channel receiving coil can be used, such as a receiving coil combining sub-coils each corresponding to a channel, or a receiving coil in which multiple loop coils are arranged in a predetermined arrangement, and here, each of the sub-coils or loop coils corresponds to an individual receiving coil. Also, a body coil that also serves as a transmitting coil may be used as the receiving coil 106, and reception may be performed by switching between them.
[0021] The measurement unit 10, under the control of the control unit 30, operates the transmission coil 105 and the gradient magnetic field coil 103 according to the imaging conditions and a predetermined pulse sequence designated by the user to irradiate RF pulses to the subject 101 in a static magnetic field, and the receiving coil 106 receives MR signals generated from the subject. Of the controls performed by the control unit 30, the control of the measurement unit 10 by the pulse sequence is realized by a sequence control device 114. Unless otherwise specified, the imaging method performed by the measurement unit 10 can adopt known MRI technology.
[0022] The calculation unit 20 includes a sensitivity distribution calculation unit 21 that calculates the sensitivity distribution of the receiving coil using the nuclear magnetic resonance signal measured by the measurement unit 10, and an image generation unit 23 that generates an image using the sensitivity distribution calculated by the sensitivity distribution calculation unit 21 and k-space data, which is measurement data. In this embodiment, the sensitivity distribution calculation unit 21 uses the sensitivity distribution calculation data to calculate the frequency components of the sensitivity distribution for each channel of the receiving coil by calculation in the k-space, and calculates the coil sensitivity distribution as real space data or weights in the k-space equivalent to the sensitivity distribution (collectively referred to as sensitivity distribution) using the frequency components. At this time, an inversion process and an expanded image creation process, which will be described later in detail, are performed on the image used for the sensitivity distribution, and the sensitivity distribution is calculated using the expanded image. The image generation unit 23 performs image reconstruction based on the parallel imaging method using the calculated sensitivity distribution or weights in the k-space.
[0023] The functions of the calculation unit 20 can be realized, for example, by installing and executing a program that defines the procedure for calculating sensitivity distribution and generating images on a computer equipped with a CPU or GPU and memory. However, it is also possible to configure a part or all of the calculations with hardware such as an ASIC (Application Specific Integrated Circuit) or an FPGA (Field Programmable Gate Array), or to realize the calculations with a computer other than the MRI apparatus or a calculation means on the cloud.
[0024] The control unit 30 includes a sequence control device 114 that controls the measurement unit 10 according to a predetermined pulse sequence, and controls the operation of the measurement unit 10 and the calculation unit 20. When the calculation unit 20 is realized by a CPU of a computer, some or all of its functions may be realized by loading a control program into the same CPU, or some or all of it may be realized by a separate CPU or hardware.
[0025] In the configuration shown in FIG. 2, the functions of the calculation unit 20 and the control unit 30, which are realized by software, are incorporated into a computer 109, and the computer 109 is connected to a display 110 which serves as a user interface device 40, an input device 115, an external storage device 111, and the like.
[0026] The following describes an outline of the imaging procedure when performing parallel imaging using the MRI apparatus of this embodiment. As shown in Fig. 3, the parallel imaging includes the steps of measurement (S1001), sensitivity distribution calculation (S1002), signal separation processing (S1003), and image display (S1004).
[0027] In the measurement step S1001, a multi-channel receiver coil is used to measure NMR signals with each receiver coil to obtain image data for each channel. NMR signals become k-space data of a given matrix size determined by phase encoding and the number of samples, but in parallel imaging, signals are thinned out (undersampling) to perform high-speed measurement.
[0028] In the sensitivity distribution calculation step S1002, the sensitivity distribution of each channel is calculated. The sensitivity distribution is calculated using image data of each channel and image data of the reference from which the influence of the sensitivity distribution is eliminated. These images for calculating the sensitivity distribution may be obtained by performing a measurement separately from the measurement (main measurement) in step S1001, or the sensitivity distribution data can be obtained during the main measurement. When a wide-area receiving coil such as a body coil is provided, it is also possible to obtain a reference image using the wide-area receiving coil. In the MRI apparatus of this embodiment, the sensitivity distribution of each channel or the weighting coefficient in the k-space is calculated by using these images for the sensitivity distribution through a calculation in the k-space. At this time, an inversion process is repeated on an image (original image) generated from the data for the sensitivity distribution to generate a developed image, and the sensitivity distribution is calculated in the frequency domain as the data for the sensitivity distribution. Details of the inversion process and the calculation of the sensitivity distribution will be described later.
[0029] In the signal separation processing step S1003, signals that would otherwise overlap in the image space due to undersampling are separated (expanded) using sensitivity distribution to generate an image without aliasing. The signal separation processing can be performed using a known parallel imaging operation such as the SENSE method, the SMASH method, or the GRAPPA method.
[0030] In an image display step S1004, the image of the subject generated by the signal separation process is displayed on the display 110 or the like.
[0031] According to the MRI apparatus of this embodiment, the problem of spectral leakage is not solved by widening the FOV of the sensitivity distribution image by zero-filling as in the conventional method, but by using an expanded image (mirror image) of the sensitivity distribution data, the problem of spectral leakage is solved and the problem of the sensitivity distribution oscillating in the non-signal region is solved, and a stable sensitivity distribution can be obtained. Furthermore, according to this embodiment, the kernel size can be relatively large when performing deconvolution calculation in the frequency domain, and a sensitivity distribution that represents a wide frequency range can be obtained. The reasons for this are briefly explained below.
[0032] When imaging an imaging region with a multi-channel coil, the sensitivity is higher the closer to the coil, so if the sensitivity distribution is high at one end of the FOV, it is low at the opposite end, and is discontinuous at the end of the FOV. When a Fourier transform (conversion of image data to k-space data) is performed in this state, spectral leakage occurs, and the sensitivity distribution has a wide range of frequencies. In the conventional method, this problem is addressed by expanding the FOV of the sensitivity distribution image by a factor of two or more with zero filling, allowing any sensitivity distribution in the zero-filled region and making it possible to connect the ends of the FOV. In that case, there is an arbitrary sensitivity distribution in the zero-filled region, which increases the degree of freedom of the solution and results in large values throughout the kernel, causing the sensitivity distribution in the zero-filled region and background region to oscillate significantly.
[0033] In this embodiment, instead of expanding the sensitivity distribution image with zero fill, a mirror image is created and used as an expanded image, which makes it possible to eliminate discontinuities at the ends of the FOV while suppressing the degree of freedom of the solution, and to suppress vibrations in the sensitivity distribution in the background region.
[0034] In the conventional method, the vibration of the sensitivity distribution in the background region can be suppressed by reducing the kernel size, but reducing the kernel size means that the system cannot handle high frequencies, which results in errors. In this embodiment, the vibration can be suppressed, so the kernel size can be made relatively large, and a sensitivity distribution that also represents high frequencies can be obtained.
[0035] Next, an embodiment in which the function of the sensitivity distribution calculation unit is further specified will be described. In the following embodiment, an embodiment of the sensitivity distribution calculation method will be described with a focus on the function of the sensitivity distribution calculation unit 21, but the functions of the measurement unit 10 and the image generation unit 23 will also be mentioned as necessary. In addition, among the configurations shown in Fig. 1 and Fig. 3, for configurations common to each embodiment, refer to these drawings as appropriate.
[0036] First Embodiment In this embodiment, a case where a two-dimensional (2D) sensitivity distribution is calculated will be described. In order to calculate the 2D sensitivity distribution, a process of inverting the sensitivity distribution data (original sensitivity distribution image) obtained in the sensitivity distribution sequence is repeated to generate a 2D developed image, which is used to calculate the sensitivity distribution.
[0037] The processing of the MRI apparatus of this embodiment will be described again with reference to Fig. 3. Here, as an example, a case will be described in which a two-dimensional (2D) sensitivity distribution is calculated using data (sensitivity distribution data) obtained from a three-dimensional (3D) sensitivity distribution sequence.
[0038] First, in a measurement step (FIG. 3, S1001), the measurement unit 10 performs measurement to obtain data for sensitivity distribution. In this embodiment, the data for sensitivity distribution is obtained by executing a pulse sequence for sensitivity distribution separately from a pulse sequence for actual imaging (actual imaging sequence) for obtaining an image of the subject. The data for sensitivity distribution is obtained by both a receiver coil having multiple channels (hereinafter referred to as an array coil) and a body coil.
[0039] The type of pulse sequence is not particularly limited, but FIG. 4 shows a typical pulse sequence 400 of the 3D-GrE (gradient echo) method as an example of a sequence for sensitivity distribution. In the figure, from the top to bottom, RF indicates the irradiation timing of a radio frequency magnetic field pulse 402 and the reception timing of a nuclear magnetic resonance signal 407, and Gs, Gp, and Gr indicate the polarity / intensity and application timing of gradient magnetic field pulses in three mutually orthogonal axial directions. The same three-digit number indicates the same type of pulse, and the same number followed by a hyphen after the three-digit number indicates that the same number is a pulse applied within the same repetition time. In addition, Gs is a slice selection gradient magnetic field, Gp is a phase encoding gradient magnetic field, and Gr is a readout gradient magnetic field, and the vertical arrow next to the pulse indicates that the gradient magnetic field intensity changes with each repetition of the pulse sequence.
[0040] In the imaging sequence of FIG. 4, first, an RF pulse 402 and a slice selection gradient magnetic field 401 are applied to excite a predetermined region, and then a rephase pulse 403 of a slice gradient magnetic field Gs and a phase encoding pulse 404 are applied. At this time, a dephase pulse 405 of a readout gradient magnetic field Gr is applied, and a rephase pulse with the polarity reversed is applied while sampling a nuclear magnetic resonance signal 407. Finally, a rephase pulse 409 is applied in the phase encoding direction. The application of the RF pulse 402 to the rephase pulse 409 is repeated at a predetermined repetition time TR to collect nuclear magnetic resonance signals required for image reconstruction. Digital data of the sampled nuclear magnetic resonance signals becomes k-space data. In this sequence, 3D-data is obtained. The 3D-data is obtained for each channel of the array coil and for each of the body coils.
[0041] The actual imaging sequence and the sensitivity distribution sequence may be the same or different in type, but in the case of the sensitivity distribution sequence, the number of encoding steps is small for both slice encoding and phase encoding, and a low-resolution sensitivity distribution image is acquired. This makes it possible to reduce the measurement time required for the prescan when executing the sensitivity distribution sequence as a prescan separate from the actual imaging sequence. However, since the resolution and field of view (FOV) of the low-resolution sensitivity distribution do not match those of the high-resolution image of the actual imaging, in order to apply the coil sensitivity distribution to image reconstruction, a process is performed to match it with the image of the actual imaging in the sensitivity distribution calculation process described later.
[0042] Next, in a sensitivity distribution calculation step (FIG. 3: S1002), the sensitivity distribution calculation unit 21 calculates a coil sensitivity distribution to be applied to the image of the actual imaging (target image) using data for sensitivity distribution (hereinafter referred to as Smap data) collected by the measurement unit 10. Since parallel imaging is performed in this embodiment, the coil sensitivity distribution is calculated for each of a plurality of receiving coils (channels).
[0043] Details of the sensitivity distribution calculation step S1002 are shown in Fig. 5. The sensitivity distribution calculation is roughly divided into a slice extraction process S151, an unfolded image creation process S152, a conversion of the Smap data that became the unfolded image to k-space data S153, a calculation of frequency components of each receiving coil using the k-space data S154, a conversion of the frequency components (k-space) of the sensitivity distribution to real space data S155, a process S156 of extracting the original region from the sensitivity distribution that became the unfolded image, and a process S157 of multiplying a mask image as necessary.
[0044] To achieve these processes, the sensitivity distribution calculation unit 21 includes, for example, a slice cutout unit 211, an expanded image creation unit 212, a Fourier transform unit 213, a sensitivity distribution frequency component calculation unit 214, a mirror image cutout unit 215, and a mask image multiplication unit 216, as shown in the functional block diagram of Fig. 6. Note that Fig. 6 is merely one example of the sensitivity distribution calculation unit 21, and it is possible to omit some of the functional units as appropriate, or to add other functional units. Below, details of the processes performed by the sensitivity distribution calculation unit 21 will be described as the operation of each unit.
[0045] [ slice Cutout process S151] The slice extraction process is performed on both the Smap data received by the array coil (the same number of Smap data as the number of channels) and the Smap data received by the body coil. Note that the Smap data obtained by the body coil is used as a reference in the subsequent sensitivity distribution calculation, so in the following explanation, it is referred to as Ref data (Ref image data in the image space).
[0046] Details of the slice extraction process S151 are shown in FIG. In the slice cutout process S151, the Fourier transform unit 213 first performs Fourier transform on the 3D-Smap data and the 3D-Ref data, respectively, and converts the k-space data into image data (real space data) (S1511). Then, the slice cutout unit 211 cuts out 2D image data (slice) from the 3D image data, which is real space data (S1512). At this time, the 2D image data of the slice to be cut out is cut out so that its FOV is the same as that of the target image (image obtained in this imaging sequence). In general, in sensitivity distribution measurement, a sensitivity distribution image is acquired with a FOV larger than that of the target image, so that the FOV of the sensitivity distribution image is made to match the FOV of the target image. However, the image matrix size of the sensitivity distribution image is a low-resolution image smaller than that of the target image.
[0047] [Expanded image creation process S152] The exfoliated image generating unit 212 repeats the inversion process on the 2D-Smap image data and the 2D-Ref image data obtained by the slice extracting unit 211 in the slice extracting process S151 to generate an exfoliated image. The details of the inversion process are shown in FIG. 8(A). In FIG. 8(A), the image 1A shown in the upper left is the 2D-Smap image data 801 obtained in the slice extracting process S151, and is set as the original image. First, the original image 1A is copied downward so as to be line-symmetrical to generate a first mirror image 1B. Next, the image (1A+1B) obtained by combining the original image 1A and the mirror image 1B is copied leftward so as to be line-symmetrical to generate a second mirror image (1C+1D), and the image (1A+1B) and the second mirror image (1C+1D) are combined to obtain an exfoliated image 802. The same process is performed on the 2D-Ref image data.
[0048] [Fourier Transform S153] The Fourier transform unit 213 performs a Fourier transform on the exfoliated image that has been mirror-image expanded by the exfoliated image creation unit 212 in the exfoliated image creation process S152, and converts it into data in frequency space, that is, k-space.
[0049] [Frequency component calculation S154] The sensitivity distribution frequency component calculation unit 214 calculates the frequency components of the sensitivity distribution from the k-space data of the Smap image. The calculation of the frequency components of the sensitivity distribution can be performed by applying the method described in Patent Document 1. Hereinafter, the calculation method of the frequency components of the sensitivity distribution will be described with reference to FIG.
[0050] In the unfolded image creation process S152, Smap image data created by mirror image unfolding from the data received by the i-th channel of the array coil is designated as Si, and Ref image data created by mirror image unfolding from the data received by the body coil is designated as Sb.
[0051] In real space, if the vector of each position in real space is r, the Smap image Si(r) of the i-th channel of the array coil is expressed as the product of the sensitivity distribution Ci(r) of each coil and the Ref image Sb(r) as shown in equation (1) (Figure 9: 901). (Number 1) S i (r) = C i (r) × S b (r) (1)
[0052] When both sides of equation (1) are converted to k-space data, the relationship in equation (1) is expressed as the convolution integral of equation (2) (Figure 9: 902).
[0053] (Number 2) Si'(k)=Ci'(k)*Sb'(k) (2) In the formula, k is the k-space position vector, and Si'(k), Ci'(k), and Sb'(k) are the real space images of Si(r), Ci(r), and Sb(r) converted into k-space data, respectively (same below).
[0054] As shown in equation (2), Si'(k) is a convolution of Sb'(k) and the frequency component Ci'(k) of the sensitivity distribution.
[0055] The above-mentioned equation (2) can be expressed as equation (3) by using a weighting coefficient wi of a predetermined kernel size by ignoring the high frequency regions of Si'(k) and Ci'(k).
[0056]
number
[0057] If the elements of Si' and wi in equation (3) are arranged vertically as vectors, and Sb' is re-expressed as a matrix in which the elements correspond to each other, equation (3) can be rewritten as equation (4), and the weight vector wi (i.e., the weight coefficient) can be calculated from the pseudo-inverse matrix shown in equation (5).
[0058]
number
[0059]
number
[0060] This weighting coefficient wi corresponds to the frequency components of the sensitivity distribution of the receiving coil (channel), but since the image used for the calculation is an unfolded image, it is not the frequency components of the sensitivity distribution itself (Figure 9: 903). [Conversion to real space S155]
[0061] The Fourier transform unit 213 adjusts the weighting coefficient wi thus obtained to the resolution of the actual imaging by zero-filling (zero-filling in the k-space) or the like, and then converts it into image data by inverse Fourier transform. The resulting image is an image 904 obtained by expanding the sensitivity distribution image, reflecting S1 and Sb used to calculate the frequency components of the sensitivity distribution, as shown in Fig. 9.
[0062] [Mirror image cutout process S156] As shown in FIG. 8(B), the mirror image clipping unit 215 clips out the original region (the region corresponding to the Smap image 801 in FIG. 8) before the original Smap image is expanded into its mirror image from the inverse Fourier transformed image 904, and obtains a sensitivity distribution (sensitivity map) 905, which is real space data.
[0063] [Mask Multiplication S157] Thereafter, the mask image multiplication unit 216 multiplies the sensitivity distribution image 905 after the Fourier transform by a predetermined mask image (not shown) in order to remove noise around the sensitivity distribution. For this reason, the mask image multiplication unit 216 creates a mask image from the low-resolution image acquired as the sensitivity distribution image. The mask image can be created at any time after a slice is cut out from the 3D sensitivity image (FIG. 7: 1512).
[0064] The above-described steps S151 to S157 complete the sensitivity distribution calculation process S1002 in Fig. 3. These processes are performed for each channel of the receiving coil and for each slice position extracted in S151, and a sensitivity distribution for each slice is obtained for each channel.
[0065] Thereafter, the process of generating an image by parallel imaging calculation using the sensitivity distribution calculated by the above process (FIG. 3: signal separation step S1003, image display step S1004) is as described above. However, in this embodiment, since the sensitivity distribution is obtained for each slice, if the image obtained by this imaging sequence is a 3D image, signal separation processing is performed for each slice on the image at the corresponding slice position.
[0066] If the image obtained in the main imaging is a multi-slice image, in slice extraction (S1512), the same slice positions as those of the main imaging are extracted and sensitivity distribution calculation is performed, and signal separation processing is performed for each slice of the corresponding multi-slice image using the sensitivity distribution of the corresponding slice position.
[0067] According to this embodiment, when calculating the sensitivity distribution in the frequency domain, an image obtained by expanding the original image is used as the sensitivity distribution data, so that the vibration of the non-signal region occurring in the sensitivity distribution obtained when the original image is used as it is can be suppressed, and a stable sensitivity distribution can be obtained. In addition, since the vibration suppression can be effectively performed, the kernel size can be relatively large, which makes it possible to express the high-frequency components of the sensitivity distribution, and improves the accuracy of the signal separation processing in parallel imaging.
[0068] In order to confirm the effect of this embodiment, the results of calculating a sensitivity distribution image using an actual image of a human head are shown in Fig. 10. In Fig. 10, (a) is an Smap image received on a certain channel, (b) is a sensitivity distribution image calculated with a kernel size of 61 × 61 by the sensitivity distribution calculation method of this embodiment, (c) is an image Smac(r) obtained by synthesizing each channel image of the Smap image using the sensitivity distribution in the following formula (6), and (d) is the absolute value of the difference between the synthesized image in (c) multiplied by the sensitivity distribution in (b) and the image of each channel in (a), which is obtained by weighting the error in the sensitivity distribution of each channel by the magnitude of the pixel value, and serves as a measure of the magnitude of the error occurring in image reconstruction.
[0069]
number
[0070] On the other hand, (e), (f), and (g) are the results calculated by zero-filling the FOV of the original image by approximately twice its original size using a conventional method, rather than the method of expanding to a mirror image of this embodiment, and correspond to (b), (c), and (d), respectively. The kernel size is the same, 61 x 61. (h), (i), and (j) show the results when the kernel size is further reduced to 41 x 41 using the conventional method.
[0071] It can be seen that the sensitivity distribution according to this embodiment (b) has less vibration in the background region than the sensitivity distribution (e) according to the conventional method. Also, in the object boundary region, the error that occurs in the region indicated by the arrow in (g) does not occur in (d), so it can be seen that vibration in the object boundary region is also suppressed. Comparing (d) and (g), it can be seen that the errors in the region inside the object are equivalent. In the conventional method, if the kernel size is reduced from 61 × 61 to 41 × 41, as shown in (h), the oscillations in the sensitivity distribution can be suppressed to the same level as in (b), but the errors in (j) are generally larger than in (d).
[0072] Although the first embodiment of the MRI apparatus of the present invention has been described above, the method of acquiring the sensitivity distribution image data and the like are not limited to the method of the first embodiment, and various modifications are possible.
[0073] <Modification 1, Modification 2> For example, in the first embodiment, data received by the body coil is used as Ref image data, but as shown in FIG. 11, data synthesized using Smap data of an array coil using a method such as that described in Patent Document 1 (a synthesis method using a virtual coil sensitivity map) without receiving the data by the body coil may be used instead (variant example 1).
[0074] In the first embodiment, the sensitivity distribution sequence is performed separately from the actual imaging sequence. However, the sensitivity distribution calculation can also be performed (self-calibration) using a portion of the measurement data obtained in the actual imaging sequence (variation example 2). Furthermore, the exfoliated image creation process S152 can be modified as follows.
[0075] <Variation 3> In the first embodiment, in the exfoliated image creation process S152, the images after cutting out slices from the 3D-Smap image data are directly mirrored to obtain exfoliated images, but in this modification, for each of the Smap image and the Ref image, the periphery of the original image is zero-filled, and then the same inversion process as in the first embodiment is repeated to generate exfoliated images in which the original images after zero-filling are exfoliated. An overview of the process in this modification is shown in FIG.
[0076] In this modification, as shown in Fig. 12, the periphery of an original image 801 (1A) is zero-filled to widen the FOV to produce an image 1A0. In conventional zero-filling, such a zero-filled image is used as is for calculating the sensitivity distribution in the frequency domain, but in this modification, the image 1A0 is copied twice to be line-symmetric, as in the process (S152) of the first embodiment, to generate an exfoliated image (1A0+1B0+1C0+1D0) 803. The sensitivity distribution frequency components can be calculated using this exfoliated image in the same manner as in the first embodiment.
[0077] According to this modified example, by combining images with slight zero-filling around them, it is possible to smoothly combine mirror images, and while maintaining the effect of zero-filling in the conventional method (suppression of spectral leakage), it is possible to suppress the occurrence of vibrations in non-signal areas caused by the increased freedom in calculating the sensitivity distribution.
[0078] <Variation 4> In the first embodiment, the entire original image 1A is inverted to generate the mirror image 1B, but it is also possible to invert only a part of the image and combine them. The processing of this modification will be explained with reference to Fig. 13. For ease of understanding, Fig. 13 shows an example of a 4 x 4 image. The numbers 1 to 16 are numbers for distinguishing pixels.
[0079] In the first inversion process, only pixels 5 to 12 in the central region of original image 1A consisting of 16 pixels are copied (1B') and combined so as to be linearly symmetrical with respect to pixels 1 to 4 of original image 1A. Next, only pixels in the central two rows of combined image (1A+1B') are copied and combined so as to be linearly symmetrical with respect to pixels in the leftmost row of image (1A+1B') to generate exfoliated image 804. In this manner, exfoliated image 804 formed in this modification is an image in which the upper left region (1A) of exfoliated image 804 corresponds to original image 801, and the matrix size is smaller than that of the exfoliated image of the first embodiment.
[0080] Using this unfolded image 804 (Smap image and Ref image), the frequency components of the sensitivity distribution are calculated as coefficient w1 in k-space. The calculation method is the same as in the first embodiment, but in this modification, the calculated coefficients in k-space are mirror images with no overlapping at the edge regions as in the image 804 in Fig. 13, so it is only necessary to calculate the coefficients in the range corresponding to region 1A of the image 804. That is, the range of bx and by of the variable wi(bx, by) calculated by equation (3) is only the range of 0 or more as in equation (7), and the scale of calculation can be reduced.
[0081]
number
[0082] As described above, according to this modification, an expanded image is generated by combining images in which only a part of the sensitivity distribution data is inverted, and the expanded image is used to calculate the sensitivity distribution in the frequency domain, thereby making it possible to obtain an effect of reducing the amount of calculation in addition to the same effect as in the first embodiment.
[0083] Second Embodiment In the first embodiment, the case has been described in which the main imaging is multi-slice imaging and 2D sensitivity distribution for each slice is obtained. However, in the present embodiment, in 3D imaging, 3D data is acquired as sensitivity distribution data, and the 3D sensitivity distribution is calculated.
[0084] In this embodiment, in the measurement step (FIG. 3: S1001), 3D imaging is performed in both the main imaging sequence and the sequence for sensitivity distribution, and 3D data is obtained as data for sensitivity distribution.
[0085] In the sensitivity distribution calculation step (FIG. 3: S1002), the sensitivity distribution data (3D-Smap data) received by each channel of the receiving coil and the sensitivity distribution data (3D-Ref data) received by the body coil are used. At this time, the 3D-image data after being converted into the image space is expanded on each axis as shown in FIG. 14 to obtain an expanded image. That is, the original image (2A) 1401 is inverted in the x-axis direction to obtain a mirror image (2B) that is symmetrical with respect to the yz plane. Next, the image 1402 obtained by combining the mirror image 2B and the original image 2A is inverted in the y-axis direction to obtain a mirror image (2C+2D) that is symmetrical with respect to the xz plane. The image 1403 obtained by combining the mirror image (2C+2D) and the image 1402 is inverted in the z-axis direction to obtain a mirror image 1404 that is symmetrical with respect to the xy plane. The mirror image 1404 is used as the expanded image, and the subsequent sensitivity distribution calculation process (S153 to S157 in FIG. 5) is performed.
[0086] The sensitivity distribution can be calculated by using an equation in which each element of equation (3) is expanded three-dimensionally, and the frequency components of the sensitivity distribution are calculated in the same manner as in the first embodiment.
[0087] In this embodiment, as in variants 3 and 4 of the first embodiment, it is possible to employ variants for generating an expanded image, such as zero-filling the periphery of the original image 1401 or inverting and combining only the central region of the original image 1401, and to obtain the same effects as variants 3 and 4.
[0088] <Third embodiment> In the first and second embodiments, the sensitivity distribution for each channel is calculated on the premise that signal separation processing (S1003 in FIG. 3) is performed for each channel. In the present embodiment, however, information corresponding to the sensitivity distribution of each channel is calculated on the premise that parallel imaging of the GRAPPA method (method described in Non-Patent Document 2) is performed as the image generation method.
[0089] The processing of this embodiment will be described by taking the case of 2D imaging as an example. Fig. 15 shows the processing flow of this embodiment. In Fig. 15, the same processes as those in Fig. 5 are indicated with the same reference numerals. As shown in the figure, this embodiment is also similar to the first embodiment in that slices are extracted from 3D data obtained by a 3D sensitivity distribution sequence (S151), the exfoliated images shown in the first embodiment and its modified examples 3 and 4 are used as sensitivity distribution data to be used for calculation (S152), and the exfoliated images are converted into k-space data (S153).
[0090] The GRAPPA method is a method for reconstructing an image by interpolating non-measurement lines in the k-space for undersampled k-space data, and when interpolating non-measurement lines, data of multiple measurement lines including other coils are interpolated by weighted addition. The weights of the multiple measurement lines correspond to the weights of each channel, and in this embodiment, the weights of the data of each channel are calculated as information corresponding to the sensitivity distribution (S158).
[0091] Specifically, in calculating the frequency components of the sensitivity distribution, the frequency components of the first embodiment are Calculation S154 Instead of the equation (2) used in (1), the weight c of the i-th channel for the j-th channel is calculated based on the following equation (8): i 'Find (k, j).
[0092] (Number 8) S i '(k)=Σjc i '(k, j)*S j '(k) (8) This weight c i '(k, j) corresponds to the frequency component (weighting wi in k space) of the sensitivity distribution of each channel obtained in the first embodiment. In the first embodiment, the sensitivity distribution in the image space was obtained by performing an inverse Fourier transform on this weighting wi. In this embodiment, i Using (k, j) as a weight for the measurement data in the k-space, non-measurement lines are interpolated and an image is reconstructed.
[0093] The above describes the embodiments and modified examples of the MRI apparatus of the present invention, but these are merely examples, and it is possible to add elements or omit some elements from the structures and processes shown in the figures.
[0094] In addition, although the above-described embodiment does not specifically mention user involvement, the user may specify via a user interface such as that shown in FIG. 16 whether to expand to a mirror image (the method described in the first embodiment) or whether to perform zero filling (variant example 3 or variant example 4).
[0095] Although each embodiment of the present invention has been described above, the function (sensitivity distribution calculation program) of the sensitivity distribution calculation unit 21 in these embodiments can be realized not only by an MRI apparatus but also by any apparatus (image processing apparatus) having a means for enabling sensitivity distribution calculation in frequency space, and the function realized by such an image processing apparatus is also included in the present invention.
[0096] The present invention can be applied not only to parallel imaging that thins out k-space data, but also to multiple slice selective excitation (SMS) parallel imaging and multi-channel image synthesis (MAC synthesis). In general MAC synthesis, the sum of squares of sensitivity distributions is used, so phase information disappears, but by using the sensitivity distribution obtained by the method of this embodiment, an image with phase information (complex image) can be obtained. [Explanation of symbols]
[0097] 10 ··· measurement unit, 20 ··· calculation unit, 21 ··· sensitivity distribution calculation unit, 23 ··· image generation unit, 30 ··· control unit, 40 ··· user interface, 211 ··· slice extraction unit, 212 ··· developed image creation unit, 213 ··· Fourier transform unit, 214 ··· sensitivity distribution frequency component calculation unit, 215 ··· mirror image extraction unit, 216 ··· mask image multiplication unit.
Claims
1. a measurement unit including a receiving coil having a plurality of channels, the measurement unit measuring a nuclear magnetic resonance signal of a subject for each channel of the receiving coil; an image calculation unit that creates an image of the subject using a sensitivity distribution for each channel of the receiving coil and a channel image generated from the nuclear magnetic resonance signal for each channel measured by the measurement unit, the image calculation unit includes a sensitivity distribution calculation unit that calculates a sensitivity distribution for each channel of the receiving coil in k-space by using a sensitivity distribution image obtained for a sensitivity distribution; The sensitivity distribution calculation unit includes an exfoliated image creation unit that performs one or more inversion processes to invert at least a part of the sensitivity distribution image to form a mirror image, and generates an exfoliated image by combining the sensitivity distribution image and an inverted image, and calculates the sensitivity distribution using the exfoliated image.
2. 2. The magnetic resonance imaging apparatus according to claim 1, The magnetic resonance imaging apparatus according to claim 1, wherein the developed image generating unit generates, as the inverted image, a mirror image obtained by inverting an entire image.
3. 2. The magnetic resonance imaging apparatus according to claim 1, The magnetic resonance imaging apparatus is characterized in that the unfolded image creation unit forms, as the inverted image, a partial mirror image by inverting an area excluding an end area of the image, and generates the unfolded image by combining the partial mirror image with an end area of the original sensitivity distribution image.
4. 2. The magnetic resonance imaging apparatus according to claim 1, The magnetic resonance imaging apparatus according to claim 1, wherein the expanded image generating unit performs one or more inversion processes after expanding the periphery of the sensitivity distribution image with zero-fill.
5. 2. The magnetic resonance imaging apparatus according to claim 1, The magnetic resonance imaging apparatus, wherein the sensitivity distribution calculation unit calculates the sensitivity distribution by using a region in k-space corresponding to a region of the original sensitivity distribution image in the unfolded image.
6. 2. The magnetic resonance imaging apparatus according to claim 1, The sensitivity distribution image is a three-dimensional image, The magnetic resonance imaging apparatus, wherein the sensitivity distribution calculation unit further includes a slice extraction unit that extracts a two-dimensional sensitivity distribution image from the three-dimensional sensitivity distribution image, and calculates the sensitivity distribution for each slice.
7. 7. The magnetic resonance imaging apparatus according to claim 6, The magnetic resonance imaging apparatus according to claim 1, wherein the slice cutout unit cuts out the two-dimensional sensitivity distribution image so that the two-dimensional sensitivity distribution image has the same FOV as an FOV of the image of the subject created by the image calculation unit.
8. 2. The magnetic resonance imaging apparatus according to claim 1, The sensitivity distribution image is a three-dimensional image, The developed image generating unit performs an inversion process on the three-dimensional image, A magnetic resonance imaging apparatus for generating a three-dimensional expanded image by expanding the three-dimensional image in three-dimensional directions.
9. 2. The magnetic resonance imaging apparatus according to claim 1, The sensitivity distribution image includes a reference image and a sensitivity distribution channel image generated from a nuclear magnetic resonance signal for each channel, The magnetic resonance imaging apparatus according to claim 1, wherein the measurement unit further comprises a wide-area receiving coil as a receiving coil for measuring a nuclear magnetic resonance signal for a reference image.
10. 2. The magnetic resonance imaging apparatus according to claim 1, The sensitivity distribution image includes a reference image and a sensitivity distribution channel image generated from a nuclear magnetic resonance signal for each channel, The magnetic resonance imaging apparatus according to claim 1, wherein the image calculation unit further comprises an image synthesis unit that synthesizes the channel images to generate the reference image.
11. 2. The magnetic resonance imaging apparatus according to claim 1, the sensitivity distribution calculation unit calculates, as the sensitivity distribution, a weight for all channels of the sensitivity distribution data for each channel; The magnetic resonance imaging apparatus according to claim 1, wherein the image calculation unit performs image reconstruction based on the GRAPPA method using a sensitivity distribution including a weight in the k-space of each channel for the nuclear magnetic resonance signals measured for each channel.
12. A sensitivity distribution calculation program for calculating a sensitivity distribution for each channel by processing a nuclear magnetic resonance signal measured in each channel of a magnetic resonance imaging apparatus having a receiving coil consisting of multiple channels, the program being executed on a computer. A step of reading a sensitivity distribution image including a reference image and a channel image, or generating the sensitivity distribution image from a nuclear magnetic resonance signal acquired by the magnetic resonance imaging apparatus; performing one or more inversion processes to invert at least a portion of the sensitivity distribution image to form a mirror image, and combining the generated inverted image with the original sensitivity distribution image to generate an expanded image for each of the reference image and the channel image; A step of performing an inverse Fourier transform on the developed image to obtain k-space data; a sensitivity distribution calculation program for executing a step of calculating a sensitivity distribution for each channel from the unfolded image in k-space.
13. The sensitivity distribution calculation program according to claim 12, The sensitivity distribution image is three-dimensional data, The sensitivity distribution calculation program further comprising a step of extracting a two-dimensional sensitivity distribution image from the three-dimensional data.
14. The sensitivity distribution calculation program according to claim 12, 11. A sensitivity distribution calculation program, wherein the step of calculating a sensitivity distribution further comprises a step of removing peripheral noise from the sensitivity distribution calculated for each channel.
Citation Information
Patent Citations
Magnetic resonance imaging apparatus
JP2002315731A
Magnetic resonance imaging apparatus
JP2008259749A
Magnetic resonance imaging apparatus
JP2014236987A
Calculation of mri rf coil sensitivity using interpolation to the extended field of view
JP2017502727A
Magnetic resonance imaging apparatus and body motion correction method
JP2018196661A
