WO2024235017A1 - 变加速因子k空间轨迹的联合设计方法 - Google Patents

变加速因子k空间轨迹的联合设计方法 Download PDF

Info

Publication number
WO2024235017A1
WO2024235017A1 PCT/CN2024/090880 CN2024090880W WO2024235017A1 WO 2024235017 A1 WO2024235017 A1 WO 2024235017A1 CN 2024090880 W CN2024090880 W CN 2024090880W WO 2024235017 A1 WO2024235017 A1 WO 2024235017A1
Authority
WO
WIPO (PCT)
Prior art keywords
echo
trajectory
image frame
acceleration factor
under
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Ceased
Application number
PCT/CN2024/090880
Other languages
English (en)
French (fr)
Inventor
张祎
祖涛
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University
Original Assignee
Zhejiang University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang University filed Critical Zhejiang University
Priority to JP2025567713A priority Critical patent/JP2026517238A/ja
Publication of WO2024235017A1 publication Critical patent/WO2024235017A1/zh
Anticipated expiration legal-status Critical
Ceased legal-status Critical Current

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/05—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves
    • A61B5/055—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves involving electronic [EMR] or nuclear [NMR] magnetic resonance, e.g. magnetic resonance imaging
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • G01R33/56—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • G01R33/56—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution
    • G01R33/565—Correction of image distortions, e.g. due to magnetic field inhomogeneities
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/58—Calibration of imaging systems, e.g. using test probes, Phantoms; Calibration objects or fiducial markers such as active or passive RF coils surrounding an MR active material
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Definitions

  • the invention belongs to the field of magnetic resonance imaging and can design a k-space trajectory for a magnetic resonance parameter measurement sequence that requires acquisition with a variable acceleration factor.
  • Magnetic resonance imaging (MRI) technology can measure tissue parameters including nuclear spin density, longitudinal relaxation time (T1), transverse relaxation time (T2), apparent diffusion coefficient (ADC), cerebral blood flow rate, magnetization transfer rate, etc. These parameters can come from the nuclei of endogenous compounds present in the body (such as 1 H in water) or from exogenous substances or tracers. Directly quantified MRI parameters can provide more valuable clinical information about various pathologies in neurological, musculoskeletal, liver and myocardial imaging. This includes relaxation time imaging, diffusion imaging, perfusion imaging, functional magnetic resonance imaging (fMRI) and chemical exchange saturation transfer (CEST) imaging, which generally require the acquisition of multi-frame image data sets with pulse sequence parameter modulation. Target parameter estimation can be achieved by fitting the signal evolution in the multi-frame image data set with the corresponding mathematical model.
  • the Keyhole method combines the acquired low-resolution image and the high-resolution reference image to reduce the data required for each acquisition, but it will sacrifice the high-frequency information of the image.
  • the SLAM method can directly generate regional CEST measurements from the region of interest of arbitrary shape, with a much higher signal-to-noise efficiency than the traditional single-voxel and multi-voxel methods, but it will lose information about the internal tissue heterogeneity of the region.
  • Model-based methods such as k-wROSA, can directly estimate the asymmetric Z spectrum of interest from complete or incomplete measurements, but its acceleration capability is limited.
  • the CS-CEST method based on compressed sensing can use the sparsity of the image in the transform domain to speed up, but complete random undersampling has not yet been achieved in practical applications. These methods can improve the efficiency of data acquisition, reduce scanning time, and also help improve image quality.
  • Multi-frame imaging usually acquires multiple frames at the same location, ranging from a few to dozens. Full k-space sampling of all frames would be too time-consuming, and using the same high acceleration factor for each frame would result in strong artifacts and SNR loss. Therefore, it is a reasonable approach to use a variable acceleration factor based on the importance of different frames. This can maximize the balance between signal-to-noise ratio and scanning speed.
  • the variable acceleration sensitivity encoding method vSENSE
  • the vSENSE method adopts a variable acceleration factor strategy to achieve higher acceleration and higher image quality.
  • the KIPI method implements automatic calibration, does not require the acquisition of additional sensitivity maps, does not require full sampling frames, and allows a more flexible variable acceleration strategy.
  • variable acceleration factor strategy on multi-frame imaging sequences with long echo chain readout, and high-quality actual acceleration has not yet been achieved. Therefore, in the field of multi-frame imaging, it is of great significance to propose a method for k-space trajectory design for variable acceleration factor acquisition to achieve a consistent T2 filtering effect for all image frames.
  • the purpose of the present invention is to solve the problems of additional contrast difference caused by variable acceleration factor sampling in multi-frame imaging technology, difficulty in implementing the variable acceleration factor strategy, and to provide a joint design method of k-space trajectory with variable acceleration factor (English name: Joint design of trajectory with variable acceleration factor, hereinafter referred to as JOINT method), which can automatically design the sampling trajectory of under-sampled image frames with variable acceleration factor and eliminate the contrast difference between different frames caused by sampling.
  • JOINT method Joint design of trajectory with variable acceleration factor
  • a joint design method of a variable acceleration factor k-space trajectory is used to design a k-space trajectory of a multi-frame imaging sequence when acquiring undersampled image frames, wherein the undersampled image frames include a first undersampled image frame and a second undersampled image frame with an acceleration factor not less than 2, and the acceleration factor of the first undersampled image frame is not higher than the acceleration factor of the second undersampled image frame, and the readout module of the multi-frame imaging sequence is an echo chain readout;
  • S2 determining the number of first echo chains to be collected for the first under-sampled image frame according to the number of acquisition points in the sampling mask and the length of the echo chain, wherein the number of first echo chains is the rounded-up value of the ratio of the number of acquisition points in the sampling mask to the length of the echo chain; then sequentially designing the acquisition trajectory of each echo chain, and forming a first echo chain trajectory queue by the trajectories of all echo chains of the first under-sampled image frame;
  • echo chain trajectories satisfying the second echo chain number are selected in sequence from the first echo chain trajectory queue, and a second echo chain trajectory queue of a second under-sampled image frame is formed, thereby completing the echo chain trajectory design.
  • the first under-sampled image frame and the second under-sampled image frame are both two-dimensional images; the acceleration factor of the first under-sampled image frame is preferably 2, and the acceleration factor of the second under-sampled image frame is preferably 2 to 4.
  • both the first under-sampled image frame and the second under-sampled image frame are three-dimensional images; the acceleration factor of the first under-sampled image frame is preferably 4, and the acceleration factor of the second under-sampled image frame is preferably 4 to 12.
  • the sampling mask is a uniform downsampling mask or a variable density Poisson downsampling mask.
  • a full sampling area is reserved in the center of the sampling mask.
  • the density of the variable density Poisson downsampling mask decreases from the center to the edge.
  • the acquisition trajectory of each echo chain is designed to be read outward from the center of the k-space.
  • the acquisition trajectory of each designed echo chain should be dispersed in the entire k-space.
  • both the first under-sampled image frame and the second under-sampled image frame are three-dimensional images, and the acquisition trajectory design method of each echo chain of the first under-sampled image frame is as follows:
  • the j-th echo is sequentially assigned to the acquisition point p in the j-th group, 1 ⁇ j ⁇ L, where the acquisition point p in the j-th group is the one with the polar angle ⁇ closest to all the unassigned acquisition points in the j-th group.
  • the L acquisition points thus obtained by the allocation constitute the trajectory of the t-th echo chain;
  • the jth echo is assigned to the only remaining unassigned acquisition point in the jth group, 1 ⁇ j ⁇ L; if M>0, the jth echo is assigned to the only remaining unassigned acquisition point in the jth group, 1 ⁇ j ⁇ M.
  • the second echo chain trajectory queue is composed of the front echo chain trajectory in the first echo chain trajectory queue that meets the second echo chain number.
  • the second under-sampled image frame has multiple frames, different second under-sampled image frames are allowed to have different acceleration factors, and each second under-sampled image frame needs to design a second echo chain trajectory queue according to its own acceleration factor.
  • the present invention has the following beneficial effects:
  • the prior art designs acquisition trajectories for undersampled image frames of different acceleration factors separately.
  • the difference in trajectories leads to inconsistent T2 filtering between different frames, which introduces additional T2 contrast between different frames and interferes with the measurement of magnetic resonance parameters.
  • the present invention allows the echo chain trajectory of high acceleration factor frames to be completely derived from the echo chain trajectory queue of low acceleration factor frames, so that the T2 filtering effect is consistent between all frames and there is no additional T2 contrast, thereby achieving optimal acquisition.
  • the present invention allows arbitrary changes in the acceleration factor while jointly designing the acquisition trajectory of the image frame with the variable acceleration factor, thereby improving the quantitative effect of magnetic resonance parameters.
  • the invention enables the variable acceleration factor undersampling strategy to be implemented to provide shorter scanning time and higher image quality.
  • the method can be applied to any acquisition technique in magnetic resonance imaging that requires echo train readout, including but not limited to TSE (FSE), EPI, and variable flip angle TSE sequence (also known as SPACE sequence).
  • Multi-frame imaging is any technique in magnetic resonance imaging that requires multiple acquisitions of the same imaging object, including but not limited to relaxation time (T 1 , T 2 ) imaging, functional imaging (fMRI), magnetic resonance spectroscopy (MRS) imaging, chemical exchange saturation transfer (CEST) imaging, diffusion imaging, and perfusion imaging.
  • FIG1 is a part of the echo chain acquisition trajectory of the first undersampled image frame in the embodiment.
  • the first row is the result of the traditional method
  • the second row is the result of the JOINT method; from left to right, each column is the trajectory of the 1st, 3rd, 5th, and 7th echo chains respectively.
  • the color depth represents the acquisition order in the echo chain
  • the lightest color represents the k-space line acquired by the first echo in the echo chain
  • the darkest color represents the k-space line acquired by the last echo in the echo chain.
  • FIG2 shows the acquisition trajectory of all echo chains designed for the first undersampled image frame and the second undersampled image frame by the conventional method and the JOINT method in the embodiment.
  • (C), (F) The difference in the order of echoes in the echo chain at the same position in the trajectory of the second undersampled image frame and the trajectory of the first undersampled image frame, where (C) is the conventional method and (F) is the JOINT method.
  • FIG3 is a layer of the 3.5ppm and -3.5ppm 3D-CEST source images of healthy volunteers in the embodiment.
  • the first column is the 3.5ppm source image
  • the second column is the -3.5ppm source image
  • the third column is the difference between the 3.5ppm source image and the -3.5ppm source image.
  • Both the traditional method and the JOINT method use the KIPI method for reconstruction.
  • Figure 4 is an APTw parameter image of a healthy volunteer in the embodiment.
  • Both the traditional method and the JOINT method use KIPI for reconstruction.
  • FIG5 is a schematic flow chart of the JOINT trajectory design method of the present invention.
  • variable acceleration factor strategy uses a low acceleration factor for image frames with low signal-to-noise ratio and a high acceleration factor for image frames with high signal-to-noise ratio, which can achieve the best balance between signal-to-noise ratio and acquisition speed.
  • acquisition sequences with long echo chain readout modules are widely used in the field of multi-frame imaging, such as EPI, TSE (FSE) and variable flip angle TSE sequences.
  • variable acceleration factor strategy When the variable acceleration factor strategy is used in the above sequences, the T2 attenuation effect between different frames varies greatly due to the difference in echo chain trajectories of frames with different acceleration factors, which will introduce additional T2 - related contrast and interfere with the quantification of magnetic resonance parameters.
  • the present invention realizes that the echo chain trajectory of the high acceleration factor frame is completely derived from the echo chain trajectory of the low acceleration factor frame by jointly designing the echo chain trajectory of different acceleration factor frames, thereby ensuring that the T2 filtering effects of different acceleration factor frames are the same.
  • the process of the joint design method of variable acceleration factor k-space trajectory provided by the present invention is shown in FIG5 , and the specific implementation process thereof is described in detail below.
  • the design method performs joint design of k-space acquisition trajectory for undersampled image frames with variable acceleration factors on magnetic resonance machines.
  • the undersampled image frames are divided into first undersampled image frames and second undersampled image frames.
  • the first undersampled image frames are undersampled with a lower acceleration factor, while the second undersampled image frames are undersampled with a higher acceleration factor (AF ⁇ 2).
  • the acceleration factor of the second undersampled image frame is not lower than that of the first undersampled image frame.
  • the multi-frame image data collected by the magnetic resonance machine can be 2D data or 3D data, and the variable acceleration factor is changed according to the actual situation. If the first under-sampled image frame and the second under-sampled image frame are 2D images, the acceleration factor of the first under-sampled image frame is preferably 2, and the acceleration factor of the second under-sampled image frame is preferably 2 to 4; if the first under-sampled image frame and the second under-sampled image frame are 3D images, the acceleration factor of the first under-sampled image frame is preferably 4, and the acceleration factor of the second under-sampled image frame is preferably 4 to 12.
  • the present invention does not involve the design of acquisition in the readout direction. If it is a 2D image, the present invention designs the acquisition trajectory in the phase encoding direction; if it is a 3D image, the present invention designs the acquisition trajectory in the phase encoding and layer selection encoding directions.
  • the joint design method of the variable acceleration factor k-space trajectory in the present invention comprises the following specific steps:
  • a sampling mask is designed according to the acceleration factor and echo chain length set in advance.
  • the sampling mask design for the first under-sampling frame in the present invention can adopt a uniform downsampling mask or a variable density Poisson downsampling mask.
  • the present invention recommends that a smaller full sampling area be reserved at the center of the mask. Further, if a variable density Poisson downsampling mask is used, the present invention recommends that the density of the variable density Poisson downsampling decreases from the center to the edge.
  • the first echo chain number is the rounded-up value of the ratio of the number of acquisition points in the sampling mask to the echo chain length, i.e.
  • N represents the number of acquisition points in the sampling mask
  • L represents the echo chain length
  • the acquisition trajectory of each echo chain designed in the present invention is preferably read outward from the center of the k-space, and the acquisition trajectory of each designed echo chain should be dispersed in the entire k-space as much as possible.
  • the acquisition trajectory of each echo chain is read outward from the center of the k-space, other methods such as random arrangement can also be adopted according to actual conditions.
  • a design method for the echo chain acquisition trajectory of the first under-sampled image frame based on the golden angle is designed for the phase encoding ky and layer selection encoding kz direction trajectories of three-dimensional imaging.
  • the coordinate positions (ky, kz) of all k-space points to be acquired i.e., all N acquisition points in the sampling mask of the first under-sampled image frame, are converted into polar coordinates (r, ⁇ ) and arranged in ascending order according to the polar radius r.
  • the sorted acquisition point sequence is grouped based on the echo chain length L and the number of the first echo chain ShotNum 1.
  • the grouping rules are divided into two forms based on the remainder M of N divided by L:
  • the N acquisition points in the acquisition point sequence are directly grouped into one group with one acquisition point per ShotNum, which is exactly L groups.
  • N is not divisible by L (M>0): first divide all N collection points into two queues, A and B. The first ShotNum 1 *M collection points are divided into queue A, and the rest of the collection points are divided into queue B. In queue A, each ShotNum 1 collection point is divided into a group, which is exactly M groups; in queue B, each ShotNum 1 -1 collection point is divided into a group, which is exactly LM groups. Therefore, this situation is also divided into L groups in total, of which the first M groups each have ShotNum 1 collection point, and the last LM groups each have ShotNum 1 -1 collection points.
  • the trajectories of the ShotNum 1 echo chains are designed in sequence, where:
  • First determine the trajectory of the first echo chain assign the jth echo (1 ⁇ j ⁇ L) to the acquisition point p in the jth group, where the acquisition point p in the jth group is the acquisition point with the polar angle ⁇ closest to (j-1) ⁇ among all the unassigned acquisition points in the jth group, so that the L allocated acquisition points constitute the trajectory of the first echo chain.
  • After determining the trajectory of the first echo chain mark the allocated acquisition points in each group, and the number of unassigned acquisition points in each group will be reduced by 1.
  • the allocated acquisition points in each group are marked, and the unallocated acquisition points in each group will be reduced by 1.
  • the second echo chain number is the product of the ratio of the first under-sampled image frame acceleration factor to the second under-sampled image frame acceleration factor multiplied by the first echo chain number rounded up, that is, AF1 and AF2 represent a first under-sampled image frame acceleration factor and a second under-sampled image frame acceleration factor, respectively.
  • echo chain lengths of the first and second under-sampled image frames of the present invention are the same.
  • the second echo chain trajectory queue ⁇ 2 when constructing the second echo chain trajectory queue ⁇ 2, it is necessary to select ShotNum 2 echo chain trajectories in sequence from the first echo chain trajectory queue ⁇ 1 of the first under-sampled image frame echo chain, wherein the selected position in the queue can be adjusted according to actual needs, and preferably, the front ShotNum 2 echo chain trajectories in the ⁇ 1 queue are selected as the ⁇ 2 queue, but other methods can also be selected according to actual conditions, for example, ShotNum 2 echo chain trajectories can be randomly selected from the ⁇ 1 queue as the ⁇ 2 queue.
  • the second under-sampled image frame in the present invention has multiple frames, and different second under-sampled image frames may have different acceleration factors.
  • Each second under-sampled image frame needs to determine the second echo chain number ShotNum 2 and design the second echo chain trajectory queue according to its own acceleration factor.
  • the acceleration factor of a second under-sampled image frame is X
  • the number of echo chains to be collected is ShotNum 2 _X
  • the acceleration factor of another second under-sampled image frame is Y
  • the number of echo chains to be collected is ShotNum 2 _Y
  • the above S1 to S4 constitute the joint design method of the variable acceleration factor k-space trajectory of the present invention.
  • variable acceleration factor k-space trajectory of the present invention will further illustrate the aforementioned variable acceleration factor k-space trajectory of the present invention through an embodiment.
  • the technical effects that can be achieved by the combined design method are presented so that those skilled in the art can better understand the essence of the present invention.
  • This embodiment applies the joint design method of k-space trajectory with variable acceleration factor shown in S1 to S4 (i.e., JOINT method) to a specific example.
  • JOINT method joint design method of k-space trajectory with variable acceleration factor shown in S1 to S4
  • the downsampling mask designed in this way has 1276 sampling points.
  • each echo chain is designed to be read from the inside to the outside of the k-space, and each echo chain is evenly distributed in the entire k-space.
  • the specific method is as follows:
  • First determine the trajectory of the first echo chain fill the jth (1 ⁇ j ⁇ 142) echo into the position in the jth group whose angle ⁇ is closest to ((j-1)* ⁇ ). After determining the trajectory of the first echo chain, remove the determined positions, and reduce the remaining positions in each group by 1.
  • the jth echo is filled into the only remaining position in the jth group, where 1 ⁇ j ⁇ 140.
  • the human experiments were performed on a 3T Siemens MRI scanner (MAGNETOM Prisma) using a 64-channel head receiving coil. The human experiments were approved by the local institutional review board. The sequence used was the CEST imaging sequence, and the MRI parameter measured was the amide proton transfer effect amplitude.
  • a total of 7 CEST saturation offset frames were collected for amide proton transfer weighted (APTw) imaging, including (S0), ⁇ 3, ⁇ 3.5, and ⁇ 4-ppm.
  • the SPACE sequence is a variable flip angle TSE sequence, and the trajectory designed by the above JOINT method was used for acquisition, with an acquisition time of 2.5 minutes.
  • the acquisition data of the trajectory designed by the JOINT method of the present invention and the acquisition data of the traditional method are both reconstructed using the KIPI method: the first under-sampled image frame is first reconstructed using the GRAPPA method, and then the coil sensitivity map and the correction factor map of the second under-sampled image frame are derived therefrom; the second under-sampled image frame is reconstructed using CG-SENSE, and the image generated by CG-SENSE is multiplied point by point with the correction factor map to obtain the final reconstructed image.
  • the APTw parametric image was calculated as follows. First, the source image was registered to the first undersampled frame image 3.5ppm. Second, the phase difference of the GRE images acquired at two different TEs was calculated as a B0 map. Third, based on the calculated B0 map, the corrected +3.5-ppm and -3.5-ppm signal values were generated for each voxel. Finally, the corrected 3.5-ppm and +3.5-ppm images were subtracted to obtain the APTw parametric image.
  • FIG1 shows the trajectory design for the first undersampled image frame.
  • the trajectory of each echo chain in the traditional method is relatively concentrated (A-D), while the method used in this example ensures that the trajectory of each echo chain is relatively uniform and dispersed (E-H).
  • FIG2 shows the trajectory design results of the conventional method and the JOINT method of the present invention for the first and second undersampled image frames.
  • the conventional method designs the trajectories of the two undersampled image frames separately, resulting in different echo ordinal numbers at the same acquisition position, which can differ by up to 20 (C);
  • the JOINT method of the present invention uses a joint design, and the echo chain trajectory of the second undersampled image frame is completely derived from the first undersampled image frame, so the echo ordinal numbers at the same acquisition position of the first and second undersampled image frames are exactly the same (F), with no difference.
  • Figure 3 shows a layer from the 3.5ppm and -3.5ppm 3D-CEST source images of a healthy volunteer.
  • the second row is the data collected by the trajectory designed by the traditional method (DF)
  • the third row is the data collected by the trajectory designed by the JOINT method of the present invention (GI)
  • the trajectory is designed jointly.
  • the collected data obtained by the traditional method and the JOINT method are reconstructed using KIPI. It can be seen that due to the different T2 filtering effects of the first and second undersampled frames of the traditional method, the difference (F) of the two frames shows obvious T2- related signals, which is very different from the reference (C); while the JOINT method of the present invention makes the T2 filtering effect consistent in each frame through joint trajectory design, and the difference (I) of the two frames is consistent with the reference method (C), without additional artifacts.
  • Figure 4 shows the whole-brain APTw parameter image of a healthy volunteer.
  • the conventional method results in obvious artifacts (DF) due to the failure to maintain a consistent T2 filtering effect between all frames;
  • the JOINT method of the present invention makes the T2 filtering effect of all frames consistent, and the APTw image (GI) shows a high degree of consistency with the reference result (AC);

Landscapes

  • Physics & Mathematics (AREA)
  • Health & Medical Sciences (AREA)
  • High Energy & Nuclear Physics (AREA)
  • General Physics & Mathematics (AREA)
  • Condensed Matter Physics & Semiconductors (AREA)
  • Engineering & Computer Science (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • General Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Signal Processing (AREA)
  • Radiology & Medical Imaging (AREA)
  • Animal Behavior & Ethology (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Pathology (AREA)
  • Biophysics (AREA)
  • Public Health (AREA)
  • Veterinary Medicine (AREA)
  • Medical Informatics (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Biomedical Technology (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

本发明公开了一种变加速因子k空间轨迹的联合设计方法，属于磁共振成像领域。该方法针对具有回波链读出的多帧成像序列进行k空间轨迹设计：首先确定低加速因子帧的采样掩膜，然后根据此掩膜设计每个回波链的轨迹，再确定高加速因子帧需要采集的回波链个数，然后取低加速因子帧的相应个数的回波链轨迹作为高加速因子帧的回波链轨迹。本发明能够实现多帧磁共振成像中变加速因子采集的轨迹设计，保证不同的加速因子图像帧之间没有额外的T2对比度差异，从而使得重建后的图像具有一致性。该方法作为一种自动的轨迹设计方案，尤其适合应用于三维多帧磁共振成像。

Description

变加速因子k空间轨迹的联合设计方法 技术领域
本发明属于磁共振成像领域，可以针对需要变加速因子采集的磁共振参数测量序列进行k空间轨迹设计。
背景技术
磁共振成像(MRI)技术可以测量包括核自旋密度、纵向弛豫时间(T1)、横向弛豫时间(T2)、表观扩散系数(ADC)、脑血流速率、磁化转移速率等组织参数。这些参数可以来自于生物体内存在的内源性化合物的核(如水中的1H)，也可以来自于外源性物质或示踪剂。直接量化的MRI参数可以提供有关神经、肌肉骨骼、肝脏和心肌成像中各种病理的更有价值的临床信息。这包括弛豫时间成像、扩散成像、灌注成像、功能磁共振(fMRI)和化学交换饱和转移(CEST)成像等技术，它们一般需要采集具有脉冲序列参数调制的多帧图像数据集。目标参数估计则可以通过将多帧图像数据集中的信号演化与相应的数学模型进行拟合来实现。
为了加速多帧成像中的参数测量，学者们已经提出了很多方法。例如，Keyhole方法将获取的低分辨率图像以及高分辨率的参考图像相结合以减少每次所需采集的数据，但是会牺牲图像的高频信息。SLAM方法可以直接从任意形状的感兴趣区域产生区域CEST测量值，具有比传统单体素和多体素方法高得多的信噪比效率，但是会丢失区域内部组织异质性的信息。基于模型的方法，如k-wROSA可以直接从完整或不完整的测量中估计感兴趣的不对称Z谱，但是其加速能力有限。基于压缩感知的CS-CEST方法可以利用图像在变换域的稀疏性进行提速，但在实际应用中尚未实现完全的随机欠采样。这些方法可以提高数据采集的效率，减少扫描时间，也有利于提高图像质量。
多帧成像通常在同一位置采集多个帧，范围从几个到几十个不等。所有帧的全k空间采样将过于耗时，而对每一帧使用相同的高加速因子会导致强烈的伪影和SNR损失。因此，根据不同帧的重要性采用变加速因子是一种合理的方法， 这可以最大化平衡信噪比和扫描速度。例如，最近变加速敏感度编码方法(vSENSE)被提出以加速参数定量。vSENSE方法采用变加速因子策略，可以实现较高倍数的加速和较高的图像质量。作为对vSENSE方法的改进，KIPI方法实现了自动校准，不需要采集额外的灵敏度图，也无需全采样帧，且允许更为灵活的变加速策略。
另一方面，由于长回波链读出的序列在采集效率和信噪比方面的优势，其被广泛应用于多帧成像领域，例如TSE、EPI等。然而由于信号在时间维度上的T2或者T2 *衰减，k空间在沿回波编码的方向被施加了指数衰减滤波。对于从中心向外读出的回波轨迹，该效应等效于k空间中的低通滤波。当使用回波链读出序列进行多帧成像时，变加速因子策略会导致不同加速因子图像帧的降采样掩膜不同。传统方法对每个加速因子的降采样掩膜分别计算回波链轨迹，这导致不同帧的T2滤波效应不同，引入了额外的T2对比度，干扰了磁共振参数测量。
因此在具有长回波链读出的多帧成像序列上使用变加速因子策略仍有一定困难，高质量的实际加速尚未实现。故而在多帧成像领域，提出一种针对变加速因子采集进行k空间轨迹设计的方法以实现所有图像帧一致的T2滤波效应，具有重要意义。
发明内容
本发明的目的在于解决多帧成像技术中变加速因子采样导致的额外对比度差异、变加速因子策略落地应用困难等问题，并提供一种变加速因子k空间轨迹的联合设计方法(英文名字为Joint design of trajectory with variable acceleration factor，后续可简称为JOINT方法)，可实现自动设计变加速因子欠采样图像帧的采样轨迹，消除不同帧之间由于采样引起的对比度差异。
本发明所采用的具体技术方案如下：
一种变加速因子k空间轨迹的联合设计方法，用于设计多帧成像序列在采集欠采样图像帧时的k空间轨迹，所述欠采样图像帧包含第一欠采样图像帧和加速因子不低于2的第二欠采样图像帧，且第一欠采样图像帧的加速因子不高于第二欠采样图像帧的加速因子，多帧成像序列的读出模块为回波链读出；
所述设计方法步骤如下：
S1：针对第一欠采样图像帧，根据提前设定的加速因子和回波链长度设计采 样掩膜；
S2：根据采样掩膜中采集点的个数和回波链长度，确定第一欠采样图像帧需要采集的第一回波链个数，其中第一回波链个数为采样掩膜中采集点的个数与回波链长度之比的向上取整值；然后依次设计每个回波链的采集轨迹，并由第一欠采样图像帧所有回波链的轨迹构成第一回波链轨迹队列；
S3：针对第二欠采样图像帧，确定要采集的第二回波链个数，其中第二回波链个数为第一欠采样图像帧加速因子与第二欠采样图像帧加速因子之比乘积第一回波链个数后的向上取整值；
S4：从第一回波链轨迹队列中依次选取满足第二回波链个数的回波链轨迹，并构成第二欠采样图像帧的第二回波链轨迹队列，从而完成回波链轨迹设计。
作为优选，所述第一欠采样图像帧和第二欠采样图像帧均为二维图像；所述第一欠采样图像帧的加速因子优选为2，所述第二欠采样图像帧的加速因子优选为2到4。
作为优选，所述第一欠采样图像帧和第二欠采样图像帧均为三维图像；所述第一欠采样图像帧的加速因子优选为4，所述第二欠采样图像帧的加速因子优选为4到12。
作为优选，所述S1中，采样掩膜为均匀降采样掩膜或者变密度泊松降采样掩膜。
作为优选，所述S1中，采样掩膜的中心保留有一块全采样区域。
作为优选，所述变密度泊松降采样掩膜的密度为从中心向边缘递减。
作为优选，所述S2中，设计的每个回波链的采集轨迹为从k空间中心向外读出。
作为优选，所述S2中，设计的每个回波链的采集轨迹应分散在整个k空间中。
作为优选，所述S2中，所述第一欠采样图像帧和第二欠采样图像帧均为三维图像，第一欠采样图像帧的每个回波链的采集轨迹设计方法如下：
S21、将第一欠采样图像帧的采样掩膜中所有N个采集点的坐标转换为极坐标(r,θ)形式，然后按照极径r对所有采集点进行升序排列，然后将排序后的采集点序列基于回波链长度L和第一回波链个数ShotNum1进行分组，且分组规则基 于N除以L的余数M分为两种形式：
若M＝0，则将采集点序列中的N个采集点直接以每ShotNum1个采集点为一组按序分为L组；
若M>0，则先将采集点序列中的前ShotNum1*M个采集点以每ShotNum1个采集点为一组按序分M组，然后再将剩余的采集点以每ShotNum1-1个采集点为一组按序分L-M组，一共分为L组；
S22、基于分组结果以及黄金角φ对ShotNum1个回波链的轨迹依次进行设计，其中：
对于1≤t≤(ShotNum1-1)，设计第t个回波链的轨迹时，依次将第j个回波分配到第j组里面的采集点p，1≤j≤L，其中第j组里面的采集点p是第j组的所有未分配采集点中极角θ最接近于的采集点，从而由分配得到的L个采集点构成第t个回波链的轨迹；
设计最后的第t＝ShotNum1个回波链的轨迹时，若M＝0，则将第j个回波分配到第j组里面的仅剩的1个未分配采集点，1≤j≤L；若M>0，则将第j个回波分配到第j组里面的仅剩的1个未分配采集点，1≤j≤M。
作为优选，所述S4中，所述第二回波链轨迹队列由第一回波链轨迹队列中最靠前的满足第二回波链个数的回波链轨迹构成。
作为优选，所述第二欠采样图像帧具有多帧，不同的第二欠采样图像帧允许具有不同的加速因子，每帧第二欠采样图像帧需要根据其自身的加速因子设计第二回波链轨迹队列。
本发明相对于现有技术而言，具有以下有益效果：
现有技术对不同的加速因子的欠采样图像帧单独设计采集轨迹，轨迹的差异导致不同帧之间T2滤波不一致，这在不同的帧之间引入了额外的T2对比度、干扰了磁共振参数测量。本发明通过让高加速因子帧的回波链轨迹完全来自于低加速因子帧的回波链轨迹队列，使得T2滤波效应在所有帧之间具有一致性且没有额外的T2对比度，可实现最优采集。
相比于现有技术，本发明在允许任意变化加速因子的同时，对变加速因子图像帧的采集轨迹进行联合设计，提升了磁共振参数定量效果。本发明的提出使得变加速因子欠采样策略可以落地应用，以提供较短的扫描时间和较高的图像质量。
该方法可适用于磁共振成像中任何需要回波链读出的采集技术，包括但不限于例如TSE(FSE)、EPI、变翻转角TSE序列(又称SPACE序列)。多帧成像为磁共振成像中任何需要对同一成像对象进行多次采集的技术，包括但不限于弛豫时间(T1，T2)成像，功能成像(fMRI)，磁共振波谱(MRS)成像，化学交换饱和转移(CEST)成像，扩散成像和灌注成像。
附图说明
图1为实施例中第一欠采样图像帧的一部分回波链采集轨迹。第一行是传统方法的结果，第二行是JOINT方法的结果；从左到右每一列分别为第1、3、5、7个回波链的轨迹。颜色深浅代表其在回波链中的采集顺序，颜色最浅代表回波链中第一个回波采集的k空间线，颜色最深代表回波链中最后一个回波采集的k空间线。
图2为实施例中传统方法和JOINT方法对第一欠采样图像帧和第二欠采样图像帧设计的所有回波链构成的采集轨迹。对于第一欠采样图像帧，传统方法和JOINT方法都是加速因子(AF)＝2×2(分别在相位编码和层选编码方向)，中心有15×15全采样区域(A，D)；对于第二欠采样图像帧，传统方法使用AF＝4×2且回波链数量为4(B)，JOINT方法的AF＝9且回波链数量为4(E)。(C)，(F)：第二欠采样图像帧轨迹和第一欠采样图像帧轨迹中相同位置在回波链中的回波次序之差，其中(C)为传统方法，(F)为JOINT方法。
图3为实施例中健康志愿者的3.5ppm和-3.5ppm 3D-CEST源图像中的一层。第一列为3.5ppm源图像，第二列为-3.5ppm源图像，第三列为3.5ppm源图像与-3.5ppm源图像的差。第一行是参考方法，AF＝2×2的常规GRAPPA；第二行是传统方法轨迹采集的数据，3.5ppm帧为第一欠采样图像帧AF＝2×2，-3.5ppm帧为第二欠采样图像帧AF＝4×2，并单独设计轨迹；第三行是JOINT方法轨迹采集的数据，3.5ppm帧为第一欠采样图像帧AF＝2×2，-3.5ppm帧为第二欠采样图像帧AF＝9，并联合设计轨迹。传统方法和JOINT方法都使用KIPI方法进行重建。
图4为实施例中健康志愿者的APTw参数图像。第一行是参考方法(A-C)，使用AF＝2×2的常规GRAPPA；第二行是传统方法设计的轨迹采集的数据(D-F)，3.5ppm帧为第一欠采样图像帧AF＝2×2，其他6帧为第二欠采样图像帧AF＝4 ×2，并单独设计轨迹；第三行是JOINT方法轨迹采集的数据，3.5ppm帧为第一欠采样图像帧AF＝2×2，其他6帧为第二欠采样图像帧AF＝9，并联合设计轨迹。传统方法和JOINT方法都使用KIPI进行重建。
图5为本发明JOINT轨迹设计方法的流程示意图。
具体实施方式
下面结合附图和具体实施方式对本发明做进一步阐述和说明。本发明中各个实施方式的技术特征在没有相互冲突的前提下，均可进行相应组合。
多帧成像由于需要对同一成像物体采集多次，因此采集时间较长。变加速因子策略对信噪比低的图像帧使用低加速因子，对信噪比高的图像帧使用高加速因子，能够实现信噪比和采集速度的最佳平衡。另一方面，多帧成像领域广泛使用具有长回波链读出模块的采集序列，例如EPI、TSE(FSE)和变化翻转角TSE序列。当在上述序列中使用变加速因子策略时，由于不同加速因子帧的回波链轨迹差异导致不同帧之间的T2衰减效应有较大差异，这会引入额外的T2相关的对比度，给磁共振参数定量带来了干扰。
本发明通过对不同加速因子帧联合设计回波链轨迹，实现了高加速因子帧的回波链轨迹完全来自于低加速因子帧回波链轨迹，保证了不同加速因子帧的T2滤波效应是相同的。本发明提供的变加速因子k空间轨迹的联合设计方法流程如图5所示，下面对其具体实现过程进行展开描述。
该设计方法针对磁共振机器上变加速因子欠采样图像帧进行k空间采集轨迹的联合设计。欠采样图像帧分为第一欠采样图像帧和第二欠采样图像帧，第一欠采样图像帧以较低的加速因子欠采样，而第二欠采样图像帧以较高的加速因子(AF≥2)欠采样。第二欠采样图像帧的加速因子不低于第一欠采样图像帧。第一欠采样图像帧只有一帧，其余皆为第二欠采样图像帧。
需要说明的是，磁共振机器采集的多帧图像数据可以是2D数据也可以是3D数据，变加速因子根据实际情况进行改变。若第一欠采样图像帧和第二欠采样图像帧为2D图像，则第一欠采样图像帧的加速因子优选为2，而第二欠采样图像帧的加速因子优选为2到4；若第一欠采样图像帧和第二欠采样图像帧为3D图像，则第一欠采样图像帧的加速因子优选为4，而第二欠采样图像帧的加速因子优选为4到12。
需要说明的是，本发明不涉及对读出方向采集的设计，若为2D图像，本发明设计的是相位编码方向上的采集轨迹；若为3D图像，本发明设计的是相位编码和层选编码方向上的采集轨迹。
本发明中该变加速因子k空间轨迹的联合设计方法，具体步骤如下：
S1：针对第一欠采样图像帧，根据提前设定的加速因子和回波链长度设计采样掩膜。
本发明中对第一欠采样帧的采样掩膜设计可以采用均匀降采样掩膜或者变密度泊松降采样掩膜。本发明推荐在掩膜中心保留有较小的全采样区域。进一步，若使用变密度泊松降采样掩膜，本发明推荐变密度泊松降采的密度为从中心向边缘递减。
但需要说明的是，在本步骤中，虽然推荐使用均匀降采样掩膜或者变密度泊松降采样掩膜，但是根据情况也可以采取其他降采样掩膜，例如带有部分傅里叶的降采样掩膜或者随机降采样掩膜。
S2：根据采样掩膜中采集点的个数和回波链长度，确定第一欠采样图像帧需要采集的回波链个数(记为第一回波链个数ShotNum1)，第一回波链个数为采样掩膜中采集点的个数与回波链长度之比的向上取整值即其中N表示采样掩膜中采集点的个数，L表示回波链长度，表示向上取整。然后依次设计每个回波链的采集轨迹trai，第一欠采样图像帧所有ShotNum1个回波链的轨迹构成的队列记为第一回波链轨迹队列Ω1＝{trai},i＝1...ShotNum1。
本发明中设计的每个回波链的采集轨迹优选为从k空间中心向外读出，设计的每个回波链的采集轨迹应尽量分散在整个k空间。但需要说明的是，在本步骤中，虽然推荐每个回波链的采集轨迹为从k空间中心向外读出，但也可以根据实际情况采取其他方式，例如随机排列。
另外，本发明的实施例中，为了保证设计的每个回波链的采集轨迹应尽量分散在整个k空间中，还针对三维成像的相位编码ky和层选编码kz方向轨迹设计了一种基于黄金角的第一欠采样图像帧的回波链采集轨迹的设计方法。
为了便于描述，定义M为N除以L的余数，黄金角或基于该定义，当第一欠采样图像帧和第二欠采样图像帧均为三 维图像时，第一欠采样图像帧的每个回波链的采集轨迹设计方法如下：
首先，对于三维成像，将所有需要采集的k空间点即第一欠采样图像帧的采样掩膜中所有N个采集点的坐标位置(ky,kz)转换为极坐标形式(r,θ)，并按照极径r进行升序排列，然后将排序后的采集点序列基于回波链长度L和第一回波链个数ShotNum1进行分组。且分组规则基于N整除L的余数M分为两种形式：
对于N可以被L整除的情况(M＝0)，将采集点序列中的N个采集点直接每ShotNum1个采集点分为一组，恰好分为L组。
对于N不可以被L整除的情况(M>0)：先将所有N个采集点分为A，B两个队列，前面的ShotNum1*M个采集点分到A队列，其余的采集点分为B队列；在A队列中以每ShotNum1个采集点分为一组，恰好分为M组；在B队列中以每ShotNum1-1个采集点分为一组，恰好分为L-M组。因此这种情况一共也分为L组，其中前M组每组有ShotNum1个采集点，后L-M组每组有ShotNum1-1个采集点。
然后，基于分组结果以及黄金角φ对ShotNum1个回波链的轨迹依次进行设计，其中：
先确定第1个回波链的轨迹：将第j个回波(1≤j≤L)分配到上述第j组里面的采集点p，其中第j组里面的采集点p是第j组的所有未分配采集点中极角θ最接近于(j-1)φ的采集点，从而由分配得到的L个采集点构成第1个回波链的轨迹。确定好第1个回波链轨迹后，将各组中已分配的采集点进行标记，每组的未分配采集点将减少1个。
依次类推，确定第t个回波链的轨迹(2≤t≤(ShotNum1-1))：将第j个回波(1≤j≤L)分配到上述第j组里面的采集点p，其中第j组里面的采集点p是第j组的所有未分配采集点中极角θ最接近于的采集点，确定好第t个回波链轨迹后，将各组中已分配的采集点进行标记，每组的未分配采集点将减少1个。
再确定最后一个回波链即第t＝ShotNum1个回波链的轨迹，也分为两种情况：
对于N可以被L整除(M＝0)的情况，此时L组中每组剩余1个采集点，因此将第j个回波分配到第j组里面的仅剩的1个未分配采集点，其中1≤j≤L；
对于N不可以被L整除(M>0)的情况，此时前M组中每组剩余1个未分配采集点，后L-M组中每组剩余0个采集点，因此将第j个回波分配到第j组里面的仅剩的一个采集点，其中1≤j≤M。
S3：针对第二欠采样图像帧，确定要采集的回波链个数(记为第二回波链个数ShotNum2)，第二回波链个数为第一欠采样图像帧加速因子与第二欠采样图像帧加速因子之比乘积第一回波链个数后的向上取整值即 AF1和AF2分别表示第一欠采样图像帧加速因子和第二欠采样图像帧加速因子。
需要说明的是，本发明的第一、二欠采样图像帧的回波链长度相同。
S4：从第一欠采样图像帧的第一回波链轨迹队列中依次选取ShotNum2个回波链轨迹构成第二欠采样图像帧的第二回波链轨迹队列Ω2，至此所有欠采样图像帧的回波链轨迹设计完成。
本发明中，构建第二回波链轨迹队列Ω2时，需要从第一欠采样图像帧回波链的第一回波链轨迹队列Ω1中依次选取ShotNum2个回波链轨迹，其中在队列中选取的位置可以根据实际需要调整，优选为选取Ω1队列中最靠前的ShotNum2个回波链轨迹作为Ω2队列，但也可根据实际情况选择其他方式，例如也可从Ω1队列中随机选取ShotNum2个回波链轨迹作为Ω2队列。
需要注意的是，本发明中的第二欠采样图像帧具有多帧，并且不同的第二欠采样图像帧可以具有不同的加速因子，每帧第二欠采样图像帧需要根据其自身的加速因子确定第二回波链个数ShotNum2并设计第二回波链轨迹队列。具体而言，如果一帧第二欠采样图像帧的加速因子为X，需要采集的回波链个数为ShotNum2_X时，那么需要从第一欠采样图像帧回波链的轨迹队列中选取ShotNum2_X个回波链轨迹作为此帧的Ω2队列；当另外一帧第二欠采样图像帧的加速因子为Y时，需要采集的回波链个数为ShotNum2_Y时，那么需要从第一欠采样图像帧回波链的轨迹队列中选取ShotNum2_Y个回波链轨迹作为此帧的Ω2队列。
由此，上述S1～S4即构成了本发明的变加速因子k空间轨迹的联合设计方法。
下面将通过一个实施例来进一步展示本发明前述的变加速因子k空间轨迹 的联合设计方法所能够实现的技术效果，以便本领域技术人员更好地理解本发明的实质。
实施例
本实施例将上述S1～S4所示的变加速因子k空间轨迹的联合设计方法(即JOINT方法)应用于具体实例中，下面详细描述本方法的具体数据来源以及各步骤的具体实现方式如下：
1、轨迹设计
所有的处理和分析都是使用在PC电脑(3.2GHz)上编写的MATLAB(MathWorks,Natick,MA)软件离线进行的。
预先设定相位编码和层选编码方向的采集矩阵大小＝80×72，回波链长度＝142，第一欠采样图像帧，AF＝4；第二欠采样图像帧只有一个加速因子，AF＝9。
本实施例的轨迹设计分为以下几步：
1、设计第一欠采样图像帧的采样掩膜，采用AF＝2×2的均匀降采样掩膜(分别在相位编码和层选编码方向)，并保留中心15×15的全采样区域，并以80、72作为椭圆的长轴和短轴，椭圆之外不采集数据。据此设计的降采样掩膜有1276个采样点。
2、根据降采样掩膜计算得需要9个回波链。注意由于1276除以142等于8余140，因此一共需要9个回波链，但是第9个回波链实际采集的回波个数是140。设计每个回波链的轨迹为从k空间由内向外读出，并且每个回波链较均匀的分布在整个k空间。具体方法如下：
将所有需要采集的k空间位置(ky，kz)转换为极坐标形式(r，θ)；按照r进行升序排列，然后按序每9个分为一组，先得到140组，最后16个位置按序每8个分为一组，有2组；一共142组。这样第一组的9个位置就对应所有9个回波链的第一个回波，第X组的9个位置就对应所有9个回波链的第X个回波(1≤X≤140)。
下面描述每个回波链上的所有回波如何分配到相应的位置以保证每个回波链的轨迹较为均匀：定义黄金角
先确定第1个回波链的轨迹：将第j个(1≤j≤142)回波填入第j组里面角度θ离((j-1)*φ)最近的位置，确定好第一个回波链轨迹后移除已经确定的位置，每组剩余减少1个位置；
对于第2个回波链的轨迹：将第j(1≤j≤142)个回波填入第j组里面角度θ离((j-1)*φ+φ)最近的位置，确定好第二个回波链轨迹后移除已经确定的位置，每组剩余继续减少1个位置；
依次类推，对于第m(3≤m≤8)个回波链的轨迹：在未分配的位置中进行选择，将第j(1≤j≤142)个回波填入第j组里面角度θ离((j-1)*φ+m*(m-1)/2*φ)最近的位置，然后移除已经确定的位置。
依次确定好8个回波链的轨迹之后，第1到140组每组剩余一个位置，141、142组里没有剩余位置。因此对于第9个回波链的轨迹，将第j个回波填入第j组里仅剩的一个位置，其中1≤j≤140。
至此，第一欠采样图像帧的所有回波链轨迹设计完成。
3、第二欠采样图像帧的回波链数计算为4/9*9＝4个。因此选择第一欠采样图像帧回波链轨迹队列中的第1到4个回波链轨迹作为第二欠采样图像帧的采集轨迹。
由此，第二欠采样图像帧的回波链轨迹设计完成。
2、数据采集
人体实验在一台3T西门子磁共振扫描仪(MAGNETOM Prisma)上进行，使用64通道头部接收线圈。人体实验得到了当地机构审查委员会的批准。所使用的序列为CEST成像序列，测量的MRI参数为酰胺质子转移效应幅度。
对于3D人体脑实验，使用SPACE-CEST序列进行采集数据，运行参数如下:TE＝17ms,TR＝3s,FOV＝232×232×209mm3，分辨率＝2.9×2.9×2.9mm3，采集矩阵大小＝80×80×72(分别为相位编码，读出编码，层选编码)，回波链长度＝142,NSA＝1.3。共采集7个CEST饱和偏移帧用于酰胺质子转移加权(APTw)成像，包括(S0),±3,±3.5,and±4-ppm。选择+3.5ppm帧作为第一欠采样图像帧(AF＝4)，其余6帧为第二欠采样图像帧(AF＝9)。SPACE序列即是变翻转角TSE序列，使用上述JOINT方法设计的轨迹进行采集，采集时间为2.5分钟。
为了比较，进行了传统方法采集变加速数据的实验，其中+3.5ppm帧AF＝2 ×2，其他帧AF＝4×2，并且+3.5ppm帧中心保留15×15的ACS区域。同时将常规的GRAPPA＝2×2作为参考，所有7帧：AF＝2×2且中心有24×24的ACS区域，采集时间为5min。
为了对酰胺质子转移效应幅度进行评估，需要采集B0校正数据，采用与CEST序列相同视场、方位和分辨率的GRE序列，TR为30ms。GRE序列在双回波模式下运行，TE分别为4.92ms和9.84ms。
3、图像重建
使用本发明JOINT方法设计轨迹的采集数据和传统方法的采集数据都使用KIPI方法进行重建：先使用GRAPPA方法重建第一欠采样图像帧，然后从中得出线圈灵敏度图和第二欠采样图像帧的校正因子图；对于第二欠采样图像帧使用CG-SENSE进行重建，并且将CG-SENSE生成的图像与校正因子图逐点相乘以得到最终的重建图像。
APTw参数图像计算如下。首先，将源图像配准到第一欠采样帧图像3.5ppm。其次，将两种不同TE获取的GRE图像的相位差计算为B0地图。第三，根据计算出的B0地图，为每个体素生成校正后的+3.5-ppm和-3.5-ppm信号值。最后，将校正后的3.5-ppm和+3.5-ppm图像相减得到APTw参数图像。
4、结果分析
图1显示了对于第一欠采样图像帧的轨迹设计，传统方法每个回波链的轨迹较为集中(A-D)，而本实例中采用的方法保证了每个回波链的轨迹较为均匀和分散(E-H)。
图2展示了传统方法和本发明的JOINT方法对于第一、二欠采样图像帧的轨迹设计结果。对于第一欠采样图像帧，传统方法和提出方法都是AF＝2×2，中心有15×15全采样区域(A，D)；对于第二欠采样图像帧，传统方法使用AF＝4×2且回波链数量为4(B)，JOINT方法的AF＝9且回波链数量为4(E)。传统方法分别设计两个欠采样图像帧的轨迹，导致相同采集位置的回波序数不同，相差可达到20(C)；本发明的JOINT方法使用联合设计，第二欠采样图像帧的回波链轨迹完全来自于第一欠采样图像帧，因此第一、二欠采样图像帧的相同采集位置的回波序数完全相同(F)，没有差别。
图3展示了健康志愿者的3.5ppm和-3.5ppm 3D-CEST源图像中的一层。第 一行是参考方法，AF＝2×2的常规GRAPPA(A-C)；第二行是传统方法设计的轨迹采集的数据(D-F)，3.5ppm帧为第一欠采样图像帧AF＝2×2，-3.5ppm帧为第二欠采样图像帧AF＝4×2，并单独设计轨迹；第三行是本发明的JOINT方法设计的轨迹采集的数据(G-I)，3.5ppm帧为第一欠采样图像帧AF＝2×2，-3.5ppm帧为第二欠采样图像帧AF＝9，并联合设计轨迹。传统方法和JOINT方法得到的采集数据都使用KIPI进行重建。可以看出，由于传统方法第一、二欠采样帧的T2滤波效应不同，两帧的差(F)显示出明显的T2相关的信号，与参考(C)相差极大；而本发明的JOINT方法通过联合轨迹设计使得T2滤波效应在各帧一致，两帧的差(I)与参考方法(C)一致，没有额外的伪影。
图4显示了健康志愿者的全脑APTw参数图像。由于未能在所有帧之间保持一致的T2滤波效应，传统方法导致明显的伪影(D-F)；3.5ppm帧为第一欠采样图像帧AF＝2×2，其他6帧为第二欠采样图像帧AF＝2×4，并单独设计轨迹。而本发明的JOINT方法使得所有帧的T2滤波效应是一致的，APTw图像(G-I)与参考结果(A-C)显示出高度一致性；3.5ppm帧为第一欠采样图像帧AF＝2×2，其他6帧为第二欠采样图像帧AF＝9，并联合设计轨迹。
以上所述的实施例只是本发明的一种较佳的方案，然其并非用以限制本发明。有关技术领域的普通技术人员，在不脱离本发明的精神和范围的情况下，还可以做出各种变化和变型。因此凡采取等同替换或等效变换的方式所获得的技术方案，均落在本发明的保护范围内。

Claims (10)

  1. 一种变加速因子k空间轨迹的联合设计方法，用于设计多帧成像序列在采集欠采样图像帧时的k空间轨迹，其特征在于，所述欠采样图像帧包含第一欠采样图像帧和加速因子不低于2的第二欠采样图像帧，且第一欠采样图像帧的加速因子不高于第二欠采样图像帧的加速因子，多帧成像序列的读出模块为回波链读出；
    所述设计方法步骤如下：
    S1：针对第一欠采样图像帧，根据提前设定的加速因子和回波链长度设计采样掩膜；
    S2：根据采样掩膜中采集点的个数和回波链长度，确定第一欠采样图像帧需要采集的第一回波链个数，其中第一回波链个数为采样掩膜中采集点的个数与回波链长度之比的向上取整值；然后依次设计每个回波链的采集轨迹，并由第一欠采样图像帧所有回波链的轨迹构成第一回波链轨迹队列；
    S3：针对第二欠采样图像帧，确定要采集的第二回波链个数，其中第二回波链个数为第一欠采样图像帧加速因子与第二欠采样图像帧加速因子之比乘积第一回波链个数后的向上取整值；
    S4：从第一回波链轨迹队列中依次选取满足第二回波链个数的回波链轨迹，并构成第二欠采样图像帧的第二回波链轨迹队列，从而完成回波链轨迹设计。
  2. 如权利要求1所述的变加速因子k空间轨迹的联合设计方法，其特征在于，所述第一欠采样图像帧和第二欠采样图像帧均为二维图像；所述第一欠采样图像帧的加速因子优选为2，所述第二欠采样图像帧的加速因子优选为2到4。
  3. 如权利要求1所述的变加速因子k空间轨迹的联合设计方法，其特征在于，所述第一欠采样图像帧和第二欠采样图像帧均为三维图像；所述第一欠采样图像帧的加速因子优选为4，所述第二欠采样图像帧的加速因子优选为4到12。
  4. 如权利要求1所述的变加速因子k空间轨迹的联合设计方法，其特征在于，所述S1中，采样掩膜为均匀降采样掩膜或者变密度泊松降采样掩膜。
  5. 如权利要求1所述的变加速因子k空间轨迹的联合设计方法，其特征在于，所述S1中，所述采样掩膜的中心保留有一块全采样区域。
  6. 如权利要求4所述的变加速因子k空间轨迹的联合设计方法，其特征在 于，所述变密度泊松降采样掩膜的密度为从中心向边缘递减。
  7. 如权利要求1所述的变加速因子k空间轨迹的联合设计方法，其特征在于，所述S2中，设计的每个回波链的采集轨迹为从k空间中心向外读出，且设计的每个回波链的采集轨迹应分散在整个k空间中。
  8. 如权利要求1所述的变加速因子k空间轨迹的联合设计方法，其特征在于，所述S2中，所述第一欠采样图像帧和第二欠采样图像帧均为三维图像，第一欠采样图像帧的每个回波链的采集轨迹设计方法如下：
    S21、将第一欠采样图像帧的采样掩膜中所有N个采集点的坐标转换为极坐标(r,θ)形式，然后按照极径r对所有采集点进行升序排列，然后将排序后的采集点序列基于回波链长度L和第一回波链个数ShotNum1进行分组，且分组规则基于N除以L的余数M分为两种形式：
    若M＝0，则将采集点序列中的N个采集点直接以每ShotNum1个采集点为一组按序分为L组；
    若M>0，则先将采集点序列中的前ShotNum1*M个采集点以每ShotNum1个采集点为一组按序分M组，然后再将剩余的采集点以每ShotNum1-1个采集点为一组按序分L-M组，一共分为L组；
    S22、基于分组结果以及黄金角φ对ShotNum1个回波链的轨迹依次进行设计，其中：
    对于1≤t≤(ShotNum1-1)，设计第t个回波链的轨迹时，依次将第j个回波分配到第j组里面的采集点p，1≤j≤L，其中第j组里面的采集点p是第j组的所有未分配采集点中极角θ最接近于的采集点，从而由分配得到的L个采集点构成第t个回波链的轨迹；
    设计最后的第t＝ShotNum1个回波链的轨迹时，若M＝0，则将第j个回波分配到第j组里面的仅剩的1个未分配采集点，1≤j≤L；若M>0，则将第j个回波分配到第j组里面的仅剩的1个未分配采集点，1≤j≤M。
  9. 如权利要求1所述的变加速因子k空间轨迹的联合设计方法，其特征在于，所述S4中，所述第二回波链轨迹队列由第一回波链轨迹队列中最靠前的满足第二回波链个数的回波链轨迹构成。
  10. 如权利要求1所述的变加速因子k空间轨迹的联合设计方法，其特征在 于，所述第二欠采样图像帧具有多帧，不同的第二欠采样图像帧允许具有不同的加速因子，每帧第二欠采样图像帧需要根据其自身的加速因子设计第二回波链轨迹队列。
PCT/CN2024/090880 2023-05-16 2024-04-30 变加速因子k空间轨迹的联合设计方法 Ceased WO2024235017A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2025567713A JP2026517238A (ja) 2023-05-16 2024-04-30 可変加速係数ｋ空間軌跡の共同設計方法

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN202310549861.4 2023-05-16
CN202310549861.4A CN116643224B (zh) 2023-05-16 2023-05-16 变加速因子k空间轨迹的联合设计方法

Publications (1)

Publication Number Publication Date
WO2024235017A1 true WO2024235017A1 (zh) 2024-11-21

Family

ID=87614611

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2024/090880 Ceased WO2024235017A1 (zh) 2023-05-16 2024-04-30 变加速因子k空间轨迹的联合设计方法

Country Status (3)

Country Link
JP (1) JP2026517238A (zh)
CN (1) CN116643224B (zh)
WO (1) WO2024235017A1 (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116643224B (zh) * 2023-05-16 2025-12-30 浙江大学 变加速因子k空间轨迹的联合设计方法

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20150022207A1 (en) * 2013-07-16 2015-01-22 Siemens Aktiengesellschaft Method and magnetic resonance apparatus to acquire mr data and to reconstruct an mr image
CN106680748A (zh) * 2016-12-30 2017-05-17 上海联影医疗科技有限公司 磁共振图像获取方法及磁共振系统
CN107209239A (zh) * 2015-01-30 2017-09-26 通用电气公司 磁共振中的回波链长度的规范
CN109839607A (zh) * 2019-01-10 2019-06-04 浙江大学 一种基于变加速敏感度编码的cest图像重建方法和装置
CN113009398A (zh) * 2021-04-08 2021-06-22 浙江大学 结合k空间和图像空间重建的成像方法和装置
CN113066141A (zh) * 2020-01-02 2021-07-02 通用电气公司 Mr图像数据的重建
CN116643224A (zh) * 2023-05-16 2023-08-25 浙江大学 变加速因子k空间轨迹的联合设计方法

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2018184056A1 (en) * 2017-04-05 2018-10-11 The University Of Queensland Magnetic resonance imaging method and apparatus
CN107219481B (zh) * 2017-05-26 2019-09-27 深圳先进技术研究院 一种通用的三维欠采样轨迹设计方法
US10823806B2 (en) * 2017-11-22 2020-11-03 Siemens Healthcare Gmbh Magnetic resonance imaging using dataset undersampling
US11143730B2 (en) * 2019-04-05 2021-10-12 University Of Cincinnati System and method for parallel magnetic resonance imaging
CN113384260B (zh) * 2020-03-13 2023-04-25 上海联影医疗科技股份有限公司 加速因子调节方法、磁共振成像扫描方法、装置和设备

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20150022207A1 (en) * 2013-07-16 2015-01-22 Siemens Aktiengesellschaft Method and magnetic resonance apparatus to acquire mr data and to reconstruct an mr image
CN107209239A (zh) * 2015-01-30 2017-09-26 通用电气公司 磁共振中的回波链长度的规范
CN106680748A (zh) * 2016-12-30 2017-05-17 上海联影医疗科技有限公司 磁共振图像获取方法及磁共振系统
CN109839607A (zh) * 2019-01-10 2019-06-04 浙江大学 一种基于变加速敏感度编码的cest图像重建方法和装置
CN113066141A (zh) * 2020-01-02 2021-07-02 通用电气公司 Mr图像数据的重建
CN113009398A (zh) * 2021-04-08 2021-06-22 浙江大学 结合k空间和图像空间重建的成像方法和装置
CN116643224A (zh) * 2023-05-16 2023-08-25 浙江大学 变加速因子k空间轨迹的联合设计方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
ZU TAO, SUN YI, WU DAN, ZHANG YI: "Joint K-space and Image-space Parallel Imaging (KIPI) for accelerated chemical exchange saturation transfer acquisition", MAGNETIC RESONANCE IN MEDICINE, WILEY-LISS, US, vol. 89, no. 3, 1 March 2023 (2023-03-01), US , pages 922 - 936, XP093240704, ISSN: 0740-3194, DOI: 10.1002/mrm.29480 *

Also Published As

Publication number Publication date
JP2026517238A (ja) 2026-05-28
CN116643224B (zh) 2025-12-30
CN116643224A (zh) 2023-08-25

Similar Documents

Publication Publication Date Title
Chen et al. High-resolution 3D MR fingerprinting using parallel imaging and deep learning
Zhao et al. Applications of a deep learning method for anti-aliasing and super-resolution in MRI
Deshmane et al. Parallel MR imaging
US7576539B2 (en) Method and apparatus for reconstructing parallel-acquired image
US12523726B2 (en) Joint k-space and image-space reconstruction imaging method and device
CN100581456C (zh) 磁共振成像装置
US11047942B2 (en) Magnetic resonance imaging method and system with optimal variable flip angles
US10338178B2 (en) System and method for high-resolution spectroscopic imaging
US7224163B2 (en) Magnetic resonance imaging device
US9488710B2 (en) MRI pulse sequence based on Q-space trajectory technique
US20080119721A1 (en) Magnetic resonance imaging apparatus
Kasten et al. Data-driven MRSI spectral localization via low-rank component analysis
Solano-Castiella et al. Parcellation of human amygdala in vivo using ultra high field structural MRI
US9538936B2 (en) MRI apparatus acquires first and second MR data and generates therefrom third image data having higher contrast between blood and background tissues
Van Steenkiste et al. Super‐resolution T1 estimation: quantitative high resolution T1 mapping from a set of low resolution T1‐weighted images with different slice orientations
Hutter et al. Slice-level diffusion encoding for motion and distortion correction
CN110604571A (zh) 一种分段编码的双核同步磁共振成像方法
Lee et al. Rapid dual‐RF, dual‐echo, 3D ultrashort echo time craniofacial imaging: a feasibility study
Krämer et al. Segmentation and visualization of the human cranial bone by T2* approximation using ultra-short echo time (UTE) magnetic resonance imaging
Weigel et al. Spinal cord imaging using averaged magnetization inversion recovery acquisitions
EP2645119A1 (de) Verfahren für schnelle Spin-Echo MRT Bildgebung
CN116643224B (zh) 变加速因子k空间轨迹的联合设计方法
KR102170977B1 (ko) 자기공명분광 기반 뇌 대사물질에 대한 시변함수를 이용한 뇌 대사물질 네트워크 생성 시스템 및 방법
US12529746B2 (en) Method and system for improved magnetic resonance imaging
Yang et al. Evaluation of the inner ear by 3D fast asymmetric spin echo (FASE) MR imaging: phantom and volunteer studies

Legal Events

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

Ref document number: 24806373

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2025567713

Country of ref document: JP

Kind code of ref document: A

WWE Wipo information: entry into national phase

Ref document number: 2025567713

Country of ref document: JP

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 24806373

Country of ref document: EP

Kind code of ref document: A1