Magnetic Resonance Imaging Apparatus, Image Processing Apparatus, and Image Creation Method

Abstract

The present invention relates to a magnetic resonance imaging apparatus, an image processing apparatus, and an image creation method. The following problem is solved: When image reconstruction is performed on measurement data obtained without full sampling by successive approximation calculation in an MRI apparatus, a large number of repetitions are required, and the reconstruction time becomes long. Density is estimated from intermediate results of the repetitions, and density correction is performed using the estimated density. Density correction is performed not only for the initial value but also in each repetition process, whereby convergence can be achieved with fewer times.

Description

Technical Field

[0001] The present invention relates to a magnetic resonance imaging apparatus (hereinafter referred to as an MRI apparatus), and particularly to a technique for reconstructing an image by successive approximation using insufficient data. Background Art

[0002] There is a technique for reconstructing an image by successive approximation when data is insufficient, which is used in magnetic resonance imaging (MRI) in, for example, compressed sensing (CS) (Non-Patent Document 1, Patent Document 1). By reducing the data acquired during imaging, the imaging time can be shortened. However, in the reconstruction of an image, a large number of repetitive processes involving wavelet transform, curvelet transform, and their inverse transforms are required, and there is a problem that the reconstruction time becomes long.

[0003] There are several methods for speeding up the repetitive process. For example, there is an excellent method called the Alternating Split Bregman method, which enables high-speed processing even when multiple transforms such as wavelet transform and curvelet transform are used in the repetitive process (Non-Patent Document 2). In addition, the following method is also carried out: at the initial value of repetition, density correction is performed to accelerate convergence. Density correction is a process of weighting data corresponding to the density distribution of measurement data.

[0004] Prior Art Documents

[0005] Patent Documents

[0006] Patent Document 1: U.S. Patent Publication No. 2015 / 0126850

[0007] Non-Patent Documents

[0008] Non-Patent Document 1: Michael Lustig, David Donoho, and John M. Pauly. Sparse MRI: The Application of Compressed Sensing for Rapid MR Imaging. Magnetic Resonance in Medicine 58: 1182-1195 (2007).

[0009] Non-Patent Document 2: Gerlind Plonka, Jianwei Ma. Curvelet-Wavelet Regularized Split Bregman Iteration for Compressed Sensing. International Journal of Wavelets, Multiresolution and Information Processing | Vol.09, No.01, pp.79-110 (2011)

[0010] Although the efficiency of the iterative process has been considerably improved by the alternating split Bregman method and pre-density correction described in Non-Patent Document 2, a large number of iterations are still required, leaving the problem of a long reconstruction time. There is a desire to further shorten the reconstruction time. Summary of the Invention

[0011] To solve the above problems, the present invention performs density correction on the data that is the intermediate result of the iteration. Density correction is premised on the density of the data to be corrected, but for the data that is the intermediate result of the iterative calculation, since it is data that has not been obtained, the density changes. In the present invention, the density correction of the intermediate result is performed corresponding to the changing density.

[0012] That is, the MRI apparatus of the present invention includes: an imaging unit that acquires measurement data including nuclear magnetic resonance signals; and an image generation unit that generates an image by using the measurement data obtained by undersampling by the imaging unit and performing successive approximation calculations. The image generation unit includes: a successive approximation calculation unit that updates the measurement data by performing successive approximation calculations on the measurement data; and a density correction unit that performs density correction on the measurement data. The density correction unit performs density correction while changing the density correction amount each time during the repetition of the successive approximation calculation.

[0013] The image processing apparatus of the present invention has the functions of the above-described image generation unit.

[0014] In addition, the image processing method of the present invention is an image processing method for generating an image by successive approximation calculation of measurement data, and includes the following processes: transforming the measurement data into processing space data for successive approximation calculation; transforming the processed space data that has undergone successive approximation calculation into measurement data to update the measurement data; and repeatedly executing the process of performing the transformation and the process of performing the update using the updated measurement data. At this time, density correction is performed on the updated measurement data, and the density correction amount is made different each time during the repetition.

[0015] Advantages of the Invention

[0016] According to the present invention, not only density correction using the initial value of the density of measurement data is performed, but also density correction corresponding to the density of the data to be targeted is performed in each repeated process. As a result, the successive approximation calculation can converge with fewer repetitions.

[0017] In addition, by finally performing density estimation and density correction, the present invention can form a reconstructed image closer to the correct answer. Regarding this part, the number of repetitions can be reduced to achieve high speed. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] Figure 1 It is a diagram showing the overall structure of an example of an MRI apparatus to which the present invention is applied.

[0019] Figure 2 It is a functional block diagram of a signal processing unit of the MRI apparatus according to the first embodiment.

[0020] Figure 3 It is a flowchart showing an outline of the processing of the MRI apparatus according to the embodiment.

[0021] Figure 4 It is a flowchart showing an example of the order of iterative reconstruction according to the first embodiment.

[0022] Figure 5 It is a diagram showing an example of a density estimation method according to the first embodiment.

[0023] Figure 6 It is a diagram for explaining the effect of the first embodiment.

[0024] Figure 7 It is a flowchart showing an example of the order of iterative reconstruction according to Modification 5 of the first embodiment.

[0025] Figure 8 It is a diagram showing an example of a UI screen according to Modification 5.

[0026] Figure 9 It is a flowchart showing the order of iterative reconstruction based on the alternating split Bregman method.

[0027] Figure 10 It is a flowchart showing an example of the order of iterative reconstruction according to the second embodiment.

[0028] Figure 11 It is a flowchart showing another example of the order of iterative reconstruction according to the second embodiment.

[0029] Figure 12 It is a flowchart showing an example of the order of iterative reconstruction according to a modification of the second embodiment.

[0030] Figure 13This is a flowchart showing another example of the order of repetitive reconstruction of a modified example of the second embodiment.

[0031] Figure 14 This is a structural diagram of an image processing apparatus.

[0032] Figure 15 This is a diagram showing an example of the UI screen of the first embodiment.

[0033] Explanation of Reference Numerals

[0034] 2: Static magnetic field generation unit, 3: Gradient magnetic field generation unit, 4: Sequence generator, 5: Transmission unit, 6: Reception unit, 7: Signal processing unit, 8: Central processing unit (CPU), 14a: Transmission coil, 14b: Reception coil, 18: External storage device, 19: Display, 20: Operation unit, 71: Imaging setting unit, 73: Image creation unit, 731: Successive approximation calculation unit, 732: Density estimation unit, 733: Density correction unit, 100: MRI apparatus, 200: Image processing apparatus. Detailed Embodiment

[0035] First, the overall outline of the MRI apparatus to which the invention is applied will be described with reference to the accompanying drawings. Figure 1 This is a block diagram showing the overall structure of an embodiment of the MRI apparatus according to the present invention. This MRI apparatus uses the NMR phenomenon to obtain tomographic images of a subject, and includes a static magnetic field generation unit 2, a gradient magnetic field generation unit 3, a transmission unit 5, a reception unit 6, a signal processing unit 7, a sequence generator 4, and a central processing unit (CPU) 8. In addition, in the following description, the static magnetic field generation unit 2, the gradient magnetic field generation unit 3, the sequence generator 4, the transmission unit 5, and the reception unit 6 will also be collectively referred to as the imaging unit.

[0036] The static magnetic field generation unit 2 generates a uniform static magnetic field in the space around the subject 1, and includes a static magnetic field generation device of a permanent magnet type, a normal conduction type, or a superconducting type. Depending on the direction of the generated magnetic field, there are a vertical magnetic field type, a horizontal magnetic field type, etc., and the present invention can be applied to any of these types.

[0037] The inclined magnetic field generation unit 3 includes an inclined magnetic field coil 9 that applies an inclined magnetic field to the X, Y, and Z axes of the coordinate system (stationary coordinate system) of the MRI apparatus, and an inclined magnetic field power supply 10 that drives each inclined magnetic field coil. By driving the inclined magnetic field power supplies 10 of the respective coils according to commands from the sequence generator 4, an inclined magnetic field Gx, Gy, Gz is applied to the X, Y, and Z axes. By combining the inclined magnetic fields in these three axial directions, an inclined magnetic field in an arbitrary direction can be generated. For example, during imaging, a slice direction inclined magnetic field pulse (Gs) is applied in a direction orthogonal to the slice plane (imaging section) to set the slice plane for the subject 1, and a phase encoding direction inclined magnetic field pulse (Gp) and a frequency encoding direction inclined magnetic field pulse (Gf) are applied in the remaining two mutually orthogonal directions orthogonal to the slice plane, and position information in each direction is encoded for the echo signal.

[0038] The sequence generator 4 is a control unit that repeatedly applies high-frequency magnetic field pulses (hereinafter referred to as "RF pulses") and inclined magnetic field pulses in a given pulse sequence, operates under the control of the CPU 8, and sends various commands required for collecting data of tomographic images of the subject 1 to the transmission unit 5, the inclined magnetic field generation unit 3, and the reception unit 6.

[0039] The transmission unit 5 irradiates the subject 1 with an RF pulse in order to cause nuclear magnetic resonance of the atomic nuclei of the biological tissues constituting the subject 1, and includes a high-frequency oscillator 11, a modulator 12, a high-frequency amplifier 13, and a transmission-side high-frequency coil (transmission coil)14a. The modulator 12 amplitude-modulates the RF pulse output from the high-frequency oscillator 11 at a timing based on an instruction from the sequence generator 4, and after amplifying the amplitude-modulated RF pulse in the high-frequency amplifier 13, supplies it to the transmission coil 14a arranged close to the subject 1, thereby irradiating the subject 1 with the RF pulse.

[0040] The reception unit 6 detects an echo signal (NMR signal) emitted by nuclear magnetic resonance of the atomic nuclei of the biological tissues constituting the subject 1, and includes a reception-side high-frequency coil (reception coil)14b, a signal amplifier 15, a quadrature phase detector 16, and an A / D converter 17. The NMR signal of the response of the subject 1 induced by the electromagnetic wave irradiated from the transmission coil 14a is detected in the reception coil 14b arranged close to the subject 1, amplified in the signal amplifier 15, and then divided into two orthogonal systems of signals at a timing based on an instruction from the sequence generator 4 by the quadrature phase detector 16, and are respectively converted into digital quantities in the A / D converter 17 and sent to the signal processing unit 7.

[0041] The signal processing unit 7 performs various data processing operations and displays and stores the processing results. CPU 8 implements some of these functions. The processing performed by the signal processing unit 7 (CPU 8) includes image reconstruction of the subject 1, calculation of numerical values representing the subject's characteristics, and correction of signals and processing results using device characteristics such as receiving coil sensitivity.

[0042] The signal processing unit 7 includes an external storage device 18 such as an optical disk or magnetic disk; a display 19; and an operating unit 20 comprised of a trackball, a mouse, a keyboard, or the like. When data from the receiving unit 6 is input to the CPU 8, the CPU 8 performs signal processing, image reconstruction, and other processing, displaying the resulting tomographic image of the subject 1 on the display 19 and recording it to the magnetic disk or the like in the external storage device 18. The CPU 8 functions as a calculation unit that performs the aforementioned processing and also as a control unit that controls the sequencer 4 and the entire apparatus.

[0043] The operation unit 20 inputs various control information of the MRI apparatus and control information of processing performed by the signal processing unit 7 and is arranged close to the display 19 . The operator interactively controls various processing of the MRI apparatus through the operation unit 20 while viewing the display 19 .

[0044] In addition, Figure 1 In the embodiment, the high-frequency coil 14a on the transmitting side and the gradient magnetic field coil 9 are arranged within the static magnetic field space of the static magnetic field generating unit 2, into which the subject 1 is inserted. If the vertical magnetic field method is used, they are arranged opposite the subject 1, and if the horizontal magnetic field method is used, they are arranged to surround the subject 1. Furthermore, the high-frequency coil 14b on the receiving side is arranged opposite the subject 1 or to surround it.

[0045] The transmitting coil 14a and receiving coil 14b may be separate high-frequency coils, but a single high-frequency coil may serve as both. Alternatively, a whole-body coil that performs both transmission and reception functions may be used in combination with a separate high-frequency coil, such as a local coil.

[0046] The MRI apparatus of this embodiment includes a function of performing successive approximation operations such as CS as the function of the above-mentioned signal processing unit 7, especially as an image generation function. Figure 2 The following schematically shows the functions of the signal processing unit 7 related to image generation.

[0047] like Figure 2 As shown, the signal processing unit 7 generally comprises an imaging setting unit 71 that sets imaging conditions and creates a pulse sequence, and an image generating unit 73. The image generating unit 73 includes a successive approximation calculation unit (sometimes simply referred to as a calculation unit) 731, a density estimation unit 732, and a density correction unit 733.

[0048] The imaging setting unit 71 calculates a pulse sequence according to the imaging method and imaging conditions set by the user via the operation unit 20 or the like, and sets the sequence generator 4. The pulse sequence includes not only the imaging method of scanning along the axes of the Cartesian coordinate system, but also non-Cartesian coordinate system scans such as radial scans. The interval method of the measurement space, the density of the data points to be measured, etc. are determined according to the scanning method and the doubling rate set as the imaging condition.

[0049] The successive approximation calculation unit 731 updates the measurement data while repeating the inverse Fourier transform, wavelet transform, etc. of the measurement data to transform to a regularized space, minimization in the regularized space, Fourier transform, and inverse transform to the regularized space, etc., to create an image.

[0050] The density estimation unit 732 estimates the density of the data updated during the repeated calculation. The density estimation methods include methods of estimating based on the difference of the data, methods based on a pre-determined density change curve, etc. The details of the estimation method will be described in the embodiments below. The density correction unit 733 corrects the density of the data during the repeated processing using the density estimated by the density estimation unit 732.

[0051] The functions of the above signal processing unit 7 can be realized by executing software such as programs stored in the memory (storage unit) and CPU 8, which are built in the CPU 8. In addition, a part of its functions can also be realized by hardware such as ASIC and FPGA.

[0052] Next, refer to Figure 3 the following flow to illustrate the outline of the imaging and signal processing of the above MRI device.

[0053] First, the subject 1 is placed in the static magnetic field space of the MRI device (S101), and imaging is performed (S102). The imaging unit applies an RF pulse and a gradient magnetic field pulse to the subject according to the pulse sequence calculated by the imaging setting unit 71, and collects the NMR signal generated from the subject, thereby performing imaging. The NMR signal becomes the measurement data (also referred to as k-space data) of the matrix configured in the k-space. In this embodiment, instead of collecting data that fills the entire k-space, undersampled measurement data is obtained. The density of the obtained measurement data is determined according to the imaging parameters such as the encoding step based on the gradient magnetic field pulse, and imaging conditions such as the undersampling rate (or doubling rate).

[0054] Next, the image generation unit 73 generates the image that is ultimately targeted (S103) while updating the unmeasured data through iterative calculations by the successive approximation calculation unit 731 using the measurement data (NMR signal) obtained by imaging. The image targeted includes proton density images, enhanced images such as diffusion-weighted images, and calculated value images derived from these images. In this image generation step S103, the density correction unit 733 performs density correction of the data during processing (hereinafter referred to as density correction processing) in accordance with the density of the data that continuously changes during the iterative calculations. The density of the changing data is the density estimated by the density estimation unit 732. Methods of density estimation include methods that estimate based on the difference of data and methods based on a predetermined density change curve.

[0055] Hereinafter, a specific implementation mode of the image generation method of the present invention will be described. In addition, Figure 1 and Figure 2 the structure of the device shown in Figure 3 and the outline of the operations shown in

[0056] <First Embodiment>

[0057] In the present embodiment, in the case of reconstructing by updating the unacquired data through iterative calculations on the measurement data including the skipped, i.e., unacquired, measurement points, in each iteration, the density estimation unit 732 estimates the density based on the progress of the update of the unacquired data. The density correction unit 733 performs density correction based on the estimated density.

[0058] In Figure 4 an example of the iterative processing in the image generation process (S103) of the present embodiment is shown. The data s0(k) acquired in imaging (S102) is density-corrected (S1031), and an image is formed by inverse Fourier transform (S1032). After performing wavelet transform (S1033) on the image and extracting only the large coefficients in the threshold processing (S1034), the image is restored by inverse wavelet transform (S1035) and restored to the k-space (the space of the acquired data) by Fourier transform (S1036). Density estimation (S1037) is performed based on the updated data. The data acquired by imaging among the updated data is replaced with s0(k) (S1038) to obtain the final updated data s(k). If the updated data does not converge, the process is repeated starting from the process of density-correcting s(k) with the estimated density (S1031). If the updated data converges, the iteration ends, and an image is generated by inverse Fourier transform (S1040), and the process ends.

[0059] Here, the density correction (S1031) performed in the first repetition is based on the density (initial value) of the data s0(k). After the second time, it is based on the density of the updated data s(k), that is, the density estimated in S1037. The density correction is a process of weighting the data corresponding to the density distribution of the data. The higher the density, the smaller the weight, and the lower the density, the larger the weight. The reciprocal of the density is called the density correction amount.

[0060] The processing of each step other than step S1037 is the same as that of the known successive approximation operations such as compressive sensing, and the detailed description thereof is omitted here. In addition, in Figure 4 , a wavelet transform is shown as an example of the transform to the space for regularization, but a basis function other than the wavelet transform, a curvelet transform, etc. can also be used, or multiple basis functions can be combined.

[0061] Hereinafter, the density estimation (S1037) in the repeated processing will be specifically described with reference to Figure 5 . Figure 5 (a) is a diagram showing the existence of measurement data in the k-space 500, and the positions of the acquired data (acquired data) are shown in black. Initially, it is determined according to the sampling pattern at the time of imaging. The density ρ0(k) of this sampling pattern can be directly obtained from the density of the data around the position k. If the probability of acquiring data at each position k is given as Figure 5 p0(k) as shown by the contour line 502 in (b), then this probability p can also be set as the density ρ.

[0062] Regarding the data during the update, the density is estimated corresponding to the progress of the update. Therefore, first, the degree of inadequacy α(k) of the data update is defined as follows.

[0063] [Equation 1]

[0064] α(k) = (b(k) - c(k)) ÷ b(k) (1)

[0065] In Equation (1), b(k) is the average obtained by dividing the sum of the absolute values of the signal values s0(k) of the data acquired during imaging ( Figure 5 the black area in (a), called the data acquisition area) by the number of acquired data points, and c(k) is the average obtained by dividing the sum of the absolute values of the non-acquired data area (the white area, called the data non-acquisition area) in the signal value s(k) of the updated data by the number of data non-acquisition area points. α(k) becomes 1 in the state where the data is not updated, and becomes 0 if the average value of the data non-acquisition area is updated to the size of the average value of the acquired data.

[0066] Next, using the updated inadequacy α(k) obtained by Equation (1) and the initial value ρ0(k) of the density distribution at the time of data acquisition, the density ρ(k) that has changed due to data update is estimated by the following Equation (2).

[0067] [Equation 2]

[0068] ρ(k) = α(k) × (ρ0(k) - 1) + 1 (2)

[0069] As can be seen from Equation (2), when α(k) = 1 (data not updated), the density is estimated to be ρ0(k), and when α(k) = 0 (data is estimated to be completely updated), the density is estimated to be 1.

[0070] In addition, in this example, it is described that b(k) and c(k) are averaged over the entire data and α(k) is a fixed value independent of k, but b(k) and c(k) can also be defined for each small region that divides the k space. That is, as Figure 5 shown by the small square in (a), it is set as the average obtained by dividing the sum of the absolute values of the signal values s0(k) of the acquired data contained in the region 501 around the position k by the number of data points acquired in the region 501, and the definition of c(k) is set as the average obtained by dividing the sum of the absolute values of the signal values s(k) of the updated data in the region 501 by the number of updated data points in the region 501. In this case, α(k) becomes a value that varies according to k.

[0071] In this case, as the region 501 for averaging when obtaining b(k) and c(k), it can also be not the region around k, but set as Figure 5 the vicinity of the side of the square passing through k (the region shown in black) 503 as shown in (c). When the data acquisition probability p0(k) as shown in Figure 5 (b) is given, it can be set as the region 504 where the value of the probability p is close, as shown in Figure 5 (d), and can be freely determined.

[0072] The density correction unit 733 uses the density estimated in the above step S1037 to correct the density of the data in the middle of the iterative calculation (S1031). That is, the density of the measurement data is corrected with a density correction amount corresponding to the density, and each process of the iteration is advanced.

[0073] According to this embodiment, by estimating the density according to the progress of data update, in the middle of the iteration, the electric energy can correct the density with high accuracy according to the density of the confidence at that time.

[0074] In Figure 6Shows the result of comparing the image with density correction and the image without density correction in this embodiment. In Figure 6 , (A) is the image in the case where there is no data shortage (full sampling is performed), (B) is the image reconstructed by performing successive approximation calculation 10 times without density correction in the state of data shortage, and (C) is the same as (B) in the state of data shortage. While performing density correction in the middle of the repetition, the image is reconstructed by successive approximation calculation, and the number of repetitions is 5 times.

[0075] (D) and (E) are the differences between (A) and (B), and between (A) and (C). From the comparison between (B) and (C) or the comparison between (D) and (E), it can be seen that in this embodiment, with half the number of repetitions in the case of not performing density correction, it is close to the image obtained by full sampling. That is, it can be known that through this embodiment, the accuracy of image reconstruction is good, and it approaches the original image in a shorter time than the conventional image reconstruction.

[0076] The application of this embodiment can be preset, but the user can also select whether to perform density correction during repetition via the user interface ( Figure 1 : operation unit 20). In Figure 15 Shows an example of the user interface screen. In this example, regarding density correction during repetition, the user can select whether to implement it in a drop-down selection manner.

[0077] <Modification Example 1 of the First Embodiment>

[0078] In the first embodiment, as the density estimation method, the progress of data update is directly obtained based on the update data of the data non-acquisition area, but instead of using the update data of the data non-acquisition area, the progress of data update can be indirectly obtained based on the update data in the data acquisition area.

[0079] In this modification example, the inadequacy α(k) of data update is defined by the following formula (3), and the density ρ(k) is estimated using formula (2) of the first embodiment.

[0080] [Mathematical Formula 3]

[0081] α(k) = Δb(k) ÷ Δb0(k) (3)

[0082] Here, Δb(k) and Δb0(k) are represented by the following equations (4) and (5) respectively. Δb(k) represents the value obtained by summing up the absolute values of the differences between the signal value s(k) of the updated data at position k in the data acquisition region and the signal value s0(k) of the acquired data at the same position k in the data acquisition region. Δb0(k) represents the value obtained by summing up the absolute values of the differences between the value obtained by density correction by dividing the signal value s0(k) of the acquired data at position k in the data acquisition region by the density ρ0(k) at that position and the signal value s0(k) in the data acquisition region.

[0083] [Equation 4]

[0084] Δb(k) = ∑|s(k) - s0(k)| (4)

[0085] [Equation 5]

[0086] Δb0(k) = ∑|s0(k)÷ρ0(k) - s0(k)| (5)

[0087] The data acquisition region can be the data acquisition region in the entire k-space, or it can be set as a partial data acquisition region including position k.

[0088] In Equation (3), when the update is insufficient, since s(k) deviates from s0(k), α(k) becomes larger and approaches 1. When the data is completely updated, since s(k) does not deviate from s0(k), α(k) becomes smaller and approaches 0, and Equation (2) can be used to estimate the density.

[0089] <Variant Example 2 of the First Embodiment>

[0090] In the first embodiment, the density ρ(k) estimated by Equation (2) has no difference between the data acquisition region and the data non-acquisition region, but the density can also be made different between the data acquisition region and the data non-acquisition region. For example, the density can be set to 1 without change in the data acquisition region, and only in the data non-acquisition region, the degree of update insufficiency α(k) defined by Equation (1) can be used to set it to

[0091] [Equation 6]

[0092] ρ(k) = 1 - α(k) (6).

[0093] Alternatively, the density can be defined as shown in the following Equation (7) to vary the density corresponding to the region.

[0094] [Equation 7]

[0095] ρ(k) = w(k)×α(k)×((1 - w(k))×ρ0(k) - 1) + 1 (7)

[0096] In Equation (7), w(k) is the weight corresponding to the region. For example, if w(k) = 0 is set in the data acquisition region and w(k) = 1 is set in the data acquisition region, then it becomes

[0097] ρ(k) = 1 (data acquisition region)

[0098] ρ(k) = 1 - α(k) (data not acquired region)

[0099] It becomes the case of Equation (6).

[0100] <Modification Example 3 of the First Embodiment>

[0101] It is also possible to further adjust the density estimated in the first embodiment or its modification examples 1 and 2. As an adjustment method, the following methods can be cited. Adjust α(k) to be 0 or more and 1 or less, and adjust ρ(k) to always be in the range from ρ0(k) to 1. Alternatively, it can also be adjusted so that ρ(k) increases monotonically, that is, if the density estimated at the nth repetition is set as ρn, then with respect to the density ρn-1 estimated at the previous repetition in the repetition, ρ(k) is adjusted so that ρn ≥ ρn-1.

[0102] By performing such adjustment, even when the acquired data and the signal values of the updated data contain abnormal values due to artifacts, noise, etc., density estimation and density correction based on it can be performed without being affected by the abnormal values.

[0103] <Modification Example 4 of the First Embodiment>

[0104] In the first embodiment, the progress status α of the update is calculated using the data during the update, and based on this, the density of the data during the update is estimated. However, it is also possible to define the progress status α of the update using a function with the number of repetitions n of the operation as a variable, and use it to estimate the density of the data during the update. As the function, for example, the following Equation (8) can be used.

[0105] [Mathematical Formula 8]

[0106] α(n) = 0.5^n (8)

[0107] In this case, the density ρ(k) can also be estimated using the above Equation (2). In addition, as long as the function has an initial value of 1 and gradually decreases, it is not limited to Equation (8).

[0108] In this modification example, the above modification example 2 can also be applied. For example, the density in the data acquisition region can be fixed at 1, and the density can be changed only in the data not acquired region.

[0109] According to this modification example, since the density can be estimated according to a pre-determined function without using the updated data, the density estimation can be performed simply.

[0110] <Modification Example 5 of the First Embodiment>

[0111] In the first embodiment, it is shown that after the repetition ends, the inverse Fourier transform (S1040) is performed on s(k) to obtain the final image. However, it is also possible to perform density estimation on s(k), perform density correction using the estimated density, and then perform the inverse Fourier transform to obtain the final image.

[0112] In Figure 7 shows the processing flow of this modification example. Up to S1038 is the same as Figure 4 However, in this modification example, density estimation is performed on the measurement data obtained after the repeated operation ends, and density correction based on the estimated density is performed. Any method among the methods of the first embodiment and its modification examples 1 to 4 can be used as the method of density estimation.

[0113] The application of this modification example can be set in advance, but it is also possible for the user to select whether to append density correction after the repetition ends via the user interface ( Figure 1 : operation unit 20). In Figure 8 shows an example of the user interface screen. In this example, as shown in (A) and (B), for the density correction during repetition and the density correction after the repetition ends, the user can select whether to perform them respectively in a pull-down selection manner.

[0114] According to this modification example, since the image quality can be improved by performing density correction even if the convergence is not sufficient during the repeated processing, a good approximate image can be obtained with fewer repetitions. In addition, by enabling the user to make a selection, image reconstruction corresponding to the user's preference such as whether to prioritize image quality or the time until the image is formed can be performed.

[0115] <Second Embodiment>

[0116] In the first embodiment, for the general successive approximation repeated operation, a reconstruction is performed in which density correction is appended to the data s(k) every time the repetition is performed. However, in this embodiment, the successive approximation operation unit 731 uses the fast alternating split Bregman method during the repeated processing. In this embodiment, appending density correction during the repeated processing is the same as in the first embodiment, but it is different in that the entire processing is changed instead of appending the process of density correction to the data s(k) in the process as Figure 4 described.

[0117] First, the successive approximation operation based on the alternating split Bregman method adopted in this embodiment will be described.

[0118] The image reconstruction based on the successive approximation operation is a process of minimizing the L1 norm of measurement data that has been set as image space data through inverse Fourier transform using wavelet transform and curvelet transform to obtain a reconstructed image, and is a process of solving the problem represented by the following equation (9).

[0119]

Mathematical formula 9

[0120]

[0121] Satisfying: There exists a where f = Φu,

[0122] In equation (9), f is the imaging data, u is the reconstructed image, Φ is the transform that performs Fourier transform and extracts only in the data acquisition region, Ψc is the curvelet transform, Ψw is the wavelet transform, θc is the curvelet transform coefficient, θw is the wavelet transform coefficient, Λc is the weight set for the curvelet coefficients, and Λc is the weight set for the wavelet coefficients.

[0123] This problem is replaced by the Bregman method to repeatedly solve the following equation (10). In the following equations, the superscript character "n" represents the number of repetitions (the same hereinafter).

[0124]

Mathematical formula 10

[0125]

[0126] Here u n , is the u of the nth repetition, the solution of, is the Bregman distance, defined by .

[0127] In addition defined by . is the column gradient of E in. In addition, μ is the weight of the term and takes a value greater than 0.

[0128] Equation (10) is further replaced by the split Bregman method as follows.

[0129]

Mathematical formula 11

[0130]

[0131] Here

[0132]

[0133] Equation (11) is further replaced as follows by the alternating split Bregman method.

[0134]

Mathematical formula 12

[0135]

[0136] The optimization problems for u, θc, and θw characterized by equations (12-1) to (12-3) can be solved analytically.

[0137] The image generation process based on this alternating split Bregman method ( Figure 3 : S103) becomes Figure 9 as follows. First, initial values of u, θw, θc, bc, and bw are set (S2031). The initial values can be 0, for example. Next, u is updated by equation (12-1) (S2032). Next, θc is updated by equation (12-2) and θw is updated by equation (12-3) (S2033). Next, bc is updated by equation (12-4) and bw is updated by equation (12-5) (S2034). Next, f is updated by equation (12-6) (S2035). If the solution u does not converge, the process returns to S2032 and is repeated. If the solution u converges, the process ends.

[0138] In this embodiment, density correction is introduced into this process. Therefore, first, equations (10) for obtaining u, θc, and θw are changed as follows to equation (13).

[0139]

Mathematical formula 13

[0140]

[0141] In equation (13), F represents the Fourier transform, Φ is changed to a transform that is extracted only in the data acquisition region without including the Fourier transform, and (Ω n ) -1 represents the density correction process. The method of density correction is the same as the method described in the first embodiment and is based on the density estimated from the data (Fun) obtained in the nth repetition or the estimated density determined by the number of repetitions n. Since Ω n can be arbitrarily determined, it is set to be independent of u. μ is the weight of the term and takes a value greater than 0.

[0142] Here, in the density estimation for determining the density correction amount of (Ω n ) -1 , if it is assumed that after a certain n, the density becomes 1 and the value for which density correction is not performed, it will converge to the same solution as equation (9), and equation (9) can be solved by equation (13). If, regardless of how much n is increased, the density does not become 1, it is also assumed that after a certain n, Ω n becomes Ω∞ will converge to the solution of Equation (14), which is a variation of Equation (9). That is, as long as Ω is made to converge as in Modification 3 of the First Embodiment, n the convergence of the solution of Equation (13) can be ensured.

[0143]

Mathematical Formula 14

[0144]

[0145] satisfies: there exists a where f = Φu,

[0146]

[0147] In this embodiment, by further replacing Equation (13) with the following Equation (15), a combination (sequence) {(u 1 , θ 1 c, θ 1 w)(u 2 , θ 2 c, θ 2 w)(u 3 , θ 3 c, θ 3 w)…} of the same solution as Equation (13) can be obtained. Additionally, u, θc, θw, f, bc, and bw can be set as real numbers and thus vectors, or they can be set as complex numbers and thus represented as an n-row 2-column matrix using the matrix representation of complex numbers.

[0148]

Mathematical Formula 15

[0149]

[0150] In addition, the fixed point based on this repetition must be a solution of Equation (13).

[0151] Furthermore, Equation (15-1) is changed as in Equations (16-1) to (16-3), and optimization is performed in sequence for u, θc, and θw. Equations (16-4) to (16-6) are the same as Equations (15-2) to (15-4).

[0152]

Mathematical Formula 16

[0153]

[0154] Equation (16-1) specifically becomes

[0155]

Mathematical Formula 17

[0156]

[0157] Here, F T , Φ TThe superscript character T of etc. indicates a transposed matrix, and I N represents the identity matrix and becomes

[0158]

[0159]

[0160] Equations (16-2) and (16-3) specifically become

[0161]

Mathematical formula 18

[0162]

[0163]

[0164] Here (Ψ c u n+1 )l, λ c，l 、 (Ψ w u n+1 ) l 、 λ w，l are respectively Ψ c u n+1 、 λ c 、 Ψ w u n+1 、 λ w the l-th element of.

[0165] In Figure 10 shows the flow of the image creation process (S103) of the present embodiment in which the above density correction process (the process based on Equation (16)) is introduced.

[0166] First, initial values of u, θw, θc, bc, and bw are set (S3031). This step is the same as Figure 9 step S2031, and the initial values can be, for example, 0. Next, u is updated with Equation (16-1) including density correction (S3032). Next, θc is updated with Equation (16-2) including density correction, and θw is updated with Equation (16-3) (S3033). Next, bc is updated with Equation (16-4) including density correction, and bw is updated with Equation (16-5) (S3034). Next, f is updated with Equation (16-6) including density correction (S3035). If the solution u does not converge, return to S3032 and repeat the process. If the solution u converges, end.

[0167] The variables (θc, θw, bc, bw) updated in each iteration can also be changed as follows. The left - hand sides of the third and fourth equations are denoted as b(~) in this article n c 、b(~) n w 。

[0168]

Mathematical Formula 19

[0169]

[0170]

[0171]

[0172]

[0173] In this case, equations (16 - 1) to (16 - 6) are as follows

[0174]

Mathematical Formula 20

[0175]

[0176] In addition, equation (17) is as follows

[0177]

Mathematical Formula 21

[0178]

[0179] In addition, equation (18) is as follows

[0180]

Mathematical Formula 22

[0181]

[0182]

[0183]

[0184]

[0185] Here λ c，l 、 λ w，l are respectively λ c 、 λ w 's l - th element

[0186] In Figure 11 shows the flow of the image generation process (S103) in this case. First, set u, uw, uc, b(~)n c and b(~) n w The initial value (S4031). The initial value can be 0, for example. Next, update u (S4032) with Equation (20-1) including density correction. Next, update uc with Equation (20-2) including density correction, and update uw with Equation (20-3) including density correction (S4033). Next, update b(~) n c with Equation (20-5) including density correction n w (S4034). Next, update f with Equation (20-6) including density correction (S4035). If the solution u does not converge, return to S4032 and repeat the process. If the solution u converges, end the process.

[0187] <Modification Example 1 of the Second Embodiment>

[0188] In the second embodiment, as shown in Equation (14), density correction (Ω) is performed before wavelet transform and curvelet transform (i.e., transform to the space that minimizes the L1 norm). -1 However, density correction can also be performed on the measurement data under the condition f = Φu that is consistent with the measurement data. In this case, instead of adding the density correction Ω in Equation (13) of the second embodiment to the terms of curvelet transform and wavelet transform, it is added to the term of ||ΦFu - f|| as follows.

[0189]

Mathematical Formula 23

[0190]

[0191] In this case, the processing flow is also the same as that Figure 10 、 Figure 11 shown.

[0192] <Modification Example 2 of the Second Embodiment>

[0193] In the second embodiment, similar to Modification Example 5 of the first embodiment, density correction can be performed on the measurement data obtained at the end of repetition before obtaining the final image from the measurement data after the end of repetition, and then inverse Fourier transform is performed to obtain the final image.

[0194] In this case, as Figure 12 、 Figure 13 shown, in Figure 10 、 Figure 11After the repetition ends, processing for performing a Fourier transform on u (S5038), density correction (S1039), and inverse Fourier transform (S1040) are appended. As the density estimation method, any one of the methods of the first embodiment and its modification examples 1 to 4 can be used.

[0195] According to this modification example, since the image quality can be improved by performing density correction even if sufficient convergence is not achieved in the repetition process, a good approximate image can be obtained with further fewer repetitions.

[0196] As described above, the embodiment of the image processing of the present invention has been described as being implemented in the signal processing unit 7 (image creation unit) of the MRI apparatus, but as shown Figure 14 also, the function of the image creation unit 73 can be implemented in an image processing apparatus 200 independent of the MRI apparatus 100. The image processing apparatus 200 and the MRI apparatus 100 are connected by wire, wirelessly, or via an Internet line or the like, or receive the measurement data obtained by the imaging unit of the MRI apparatus 100 via a removable medium, and perform image reconstruction calculation. The created image can be displayed on a display attached to the image processing apparatus 200 or stored in a storage device, or can be sent to the MRI apparatus 100.

Claims

1. A magnetic resonance imaging device, characterized in that, Comprising: An imaging unit that acquires measurement data including nuclear magnetic resonance signals; and An image generation unit that generates an image by successive approximation using the measurement data obtained by undersampling performed by the imaging unit, The image generation unit comprises: A successive approximation operation unit that updates the measurement data by successive approximation for the measurement data; and A density correction unit that corrects the density of the measurement data, In the repetition of the successive approximation operation, the density correction unit corrects the density while changing the density correction amount each time, The density correction unit divides the measurement data into a plurality of regions, determines the density correction amount based on the density obtained for each region, and corrects the density.

2. A magnetic resonance imaging apparatus, characterized in that, Comprising: An imaging unit that acquires measurement data including nuclear magnetic resonance signals; and An image generation unit that generates an image by successive approximation using the measurement data obtained by undersampling performed by the imaging unit, The image generation unit comprises: A successive approximation operation unit that updates the measurement data by successive approximation for the measurement data; and A density correction unit that corrects the density of the measurement data, In the repetition of the successive approximation operation, the density correction unit corrects the density while changing the density correction amount each time, The density correction unit corrects the density of the data at the measurement points acquired by the imaging unit and the data at the measurement points estimated by successive approximation among the measurement data to be corrected with different weights.

3. The magnetic resonance imaging apparatus according to claim 1 or 2, wherein The image generation unit further comprises: A density estimation unit that estimates the density of the measurement data each time in the successive approximation operation, The density correction unit corrects the density based on the density estimated by the density estimation unit.

4. The magnetic resonance imaging apparatus according to claim 1 or 2, wherein The density correction unit corrects the density based on a preset change in density.

5. The magnetic resonance imaging apparatus according to claim 1 or 2, wherein The density correction unit adjusts the density correction amount to correct the density under the condition that the degree of insufficiency of the update of the measurement data updated by successive approximation is 0 or more and 1 or less.

6. The magnetic resonance imaging apparatus according to claim 1 or 2, wherein The density correction unit changes the density correction amount to correct the density under the condition that the change in density is monotonically increasing.

7. The magnetic resonance imaging apparatus according to claim 1 or 2, wherein The density correction unit determines the initial value of the density correction amount corresponding to the data collection method in the imaging unit.

8. The magnetic resonance imaging apparatus according to claim 1 or 2, wherein The density correction unit corrects the density after the repetition of the successive approximation operation in the operation unit ends.

9. The magnetic resonance imaging apparatus according to claim 1 or 2, wherein The successive approximation operation unit comprises: A Fourier transform unit that performs the transformation of real space data and its inverse transformation from the measurement data; and A spatial transformation unit that performs transformation from the real space data to the regularized space data and its inverse transformation.

10. The magnetic resonance imaging apparatus according to claim 1 or 2, characterized in that the successive approximation operation unit performs an operation based on the alternating split Bregman method as the successive approximation operation.

11. An image processing apparatus including an operation unit that inputs measurement data and creates an image through successive approximation operations, the image processing apparatus being characterized by including: a density correction unit that corrects the density of the measurement data, wherein the density correction unit corrects the density while varying the density correction amount each time during the repetition of the successive approximation operation, and the density correction unit divides the measurement data into a plurality of regions, determines the density correction amount based on the density obtained for each region, and corrects the density.

12. An image processing apparatus including an operation unit that inputs measurement data and creates an image through successive approximation operations, the image processing apparatus being characterized by including: a density correction unit that corrects the density of the measurement data, wherein the density correction unit corrects the density while varying the density correction amount each time during the repetition of the successive approximation operation, and the density correction unit corrects the density of the data of the measurement points obtained by the imaging unit and the data of the measurement points estimated through successive approximation operations among the measurement data to be corrected with different weights.

13. The image processing apparatus according to claim 11 or 12, characterized in that the measurement data is measurement data including nuclear magnetic resonance signals measured by a magnetic resonance imaging apparatus.

14. An image processing method for creating an image through successive approximation operations on measurement data, the image processing method being characterized by including the following processes: transforming the measurement data into processing space data for successive approximation operations; transforming the processed space data on which the successive approximation operation has been performed into measurement data and updating the measurement data; and repeatedly executing the process of performing the transformation and the process of performing the update using the updated measurement data, and at this time, correcting the density of the updated measurement data, wherein the density correction makes the density correction amount different each time during the repetition, and the density correction divides the measurement data into a plurality of regions, determines the density correction amount based on the density obtained for each region, and corrects the density.

15. An image processing method for creating an image through successive approximation operations on measurement data, the image processing method being characterized by including the following processes: transforming the measurement data into processing space data for successive approximation operations; transforming the processed space data on which the successive approximation operation has been performed into measurement data and updating the measurement data; and repeatedly executing the process of performing the transformation and the process of performing the update using the updated measurement data, and at this time, correcting the density of the updated measurement data, wherein the density correction makes the density correction amount different each time during the repetition, The density correction corrects the density of the measurement data to be corrected, which is the data of the measurement points obtained by the imaging unit and the data of the measurement points estimated by successive approximation calculation, with different weights.

16. The image processing method according to claim 14 or 15, wherein the image processing method further includes the following processing: estimating the density of the updated measurement data, and the density correction is performed based on the estimated density.

17. The image processing method according to claim 14 or 15, wherein the density correction is performed based on a preset change in density.

18. The image processing method according to claim 14 or 15, wherein the density correction adjusts the density correction amount to perform density correction under the condition that the degree of insufficiency of the update of the measurement data updated by successive approximation calculation is 0 or more and 1 or less.

19. The image processing method according to claim 14 or 15, wherein the density correction changes the density correction amount to perform density correction under the condition that the change in density is monotonically increasing.

20. The image processing method according to claim 14 or 15, wherein the density correction determines the initial value of the density correction amount corresponding to the data collection method in the imaging unit.

21. The image processing method according to claim 14 or 15, wherein the density correction is performed after the repetition of the successive approximation calculation ends.

22. The image processing method according to claim 14 or 15, wherein the successive approximation calculation includes the following processing: transforming the real space data from the measurement data and its inverse transformation; and transforming the real space data into the regularized space data and its inverse transformation.

23. The image processing method according to claim 14 or 15, wherein performing an operation based on the alternating split Bregman method as the successive approximation calculation.

Images