A data processing method and system applied to magnetic resonance imaging

Abstract

The present invention provides a data processing method and system for magnetic resonance imaging, which relates to the technical field of image processing, and includes obtaining the velocity field time series of each voxel in a three-dimensional space in four-dimensional flow nuclear magnetic resonance imaging data; performing multi-scale frequency domain processing on the velocity field time series of each voxel to obtain a frequency domain feature matrix and a frequency domain energy ratio; based on the frequency domain feature matrix and the frequency domain energy ratio, performing weighted fusion to obtain a comprehensive entropy; combining the comprehensive entropy and a preset regularity index of the frequency domain features to construct an objective optimization function, and solving the objective optimization function to screen out a voxel set that meets the set conditions; performing post-processing on the voxel set that meets the set conditions to determine the image data of the blood vessel region in the four-dimensional flow nuclear magnetic resonance imaging data. The present invention realizes the effective identification of the blood vessel region in 4DFlowMRI data by combining local and global frequency domain features and the regularity index of the frequency domain features.

Description

Technical Field

[0001] The present invention relates to the technical field of image processing, and particularly to a data processing method and system for magnetic resonance imaging data. Background Art

[0002] Magnetic resonance imaging (MRI), as a non-invasive medical imaging technology, plays a crucial role in clinical diagnosis. Among them, four-dimensional flow magnetic resonance imaging (4DFlowMRI) technology, with its ability to capture the dynamic blood flow during the cardiac cycle, provides an unprecedented perspective for the diagnosis and research of cardiovascular diseases. However, accurately extracting the information of the vascular region from a large amount of 4DFlowMRI data has always been a major challenge in this field.

[0003] Traditional 4DFlowMRI data processing methods mainly rely on image segmentation algorithms. Although these algorithms can distinguish vascular and non-vascular regions to a certain extent, they often fail to achieve ideal accuracy when dealing with complex vascular structures or low signal-to-noise ratio data. In addition, due to the high-dimensional characteristics of 4DFlowMRI data (including the time dimension and three-dimensional space dimension), traditional methods also have bottlenecks in processing speed and computational efficiency.

[0004] Therefore, it is necessary to provide a data processing method and system for magnetic resonance imaging to solve the above technical problems. Summary of the Invention

[0005] To solve the above technical problems, the present invention provides a data processing method and system for magnetic resonance imaging, which realizes the effective identification of the vascular region in 4DFlowMRI data by combining local and global frequency domain features and the regularity index of frequency domain features.

[0006] The present invention provides a data processing method for magnetic resonance imaging, and the method includes the following steps:

[0007] Expand the four-dimensional flow magnetic resonance imaging data in the time dimension to obtain the velocity field time series of each voxel;

[0008] Perform multi-scale frequency domain processing on the velocity field time series of each voxel to obtain a frequency domain feature matrix and a frequency domain energy ratio;

[0009] Calculate the local entropy and the global entropy based on the frequency domain feature matrix and calculate the weights of the local entropy and the global entropy respectively in real time based on the frequency domain energy ratio, and perform weighted fusion on the local entropy and the global entropy according to their respective weights to obtain a comprehensive entropy;

[0010] Construct a target optimization function by combining the comprehensive entropy and the regularity index of the preset frequency domain features, and use the gradient descent method to solve the target optimization function to screen out a voxel set that meets the set conditions;

[0011] Perform post-processing on the voxel set that meets the set conditions to determine the image data of the blood vessel region in the four-dimensional flow nuclear magnetic resonance imaging data.

[0012] Preferably, the unfolding of the four-dimensional flow nuclear magnetic resonance imaging data in the time dimension to obtain the velocity field time series of each voxel includes:

[0013] Decompose the four-dimensional flow nuclear magnetic resonance imaging data into three-dimensional velocity field data at multiple moments;

[0014] For each voxel of the three-dimensional velocity field data, extract the velocity values at different time points to form a corresponding velocity field time series, and perform normalization processing on the velocity field time series.

[0015] Preferably, the multi-scale frequency domain processing of the velocity field time series of each voxel to obtain a frequency domain feature matrix and a frequency domain energy ratio includes:

[0016] Perform wavelet transform on the velocity field time series to obtain wavelet coefficients representing the joint time and frequency;

[0017] Extract the energy distribution characteristics of low frequency, medium frequency and high frequency according to the wavelet coefficients to construct a frequency domain feature matrix;

[0018] Define a reference frequency range related to blood vessel flow, and calculate the frequency domain energy ratio of each voxel based on the reference frequency range.

[0019] Preferably, the calculation process of the comprehensive entropy includes:

[0020] Normalize the frequency domain features of each voxel of the frequency domain feature matrix to obtain a probability distribution, and calculate the local entropy based on the probability distribution. The calculation formula of the local entropy is:

[0021] Where, represents the local entropy, represents the voxel at different frequencies on the probability distribution;

[0022] Perform statistical analysis on the entire frequency domain feature matrix to obtain a global probability distribution, and calculate the global entropy based on the global probability distribution. The calculation formula of the global entropy is:

[0023] Where, represents the global entropy; represents the global probability distribution;

[0024] Calculate the weights of the local entropy and the global entropy respectively according to the frequency-domain energy ratio, where the calculation formula of the weight is:

[0025] where, represents the weight of the local entropy, represents the weight of the global entropy, represents the frequency-domain energy ratio, represents a regulation parameter for controlling the sensitivity of the weight to the frequency-domain energy ratio;

[0026] Perform weighted fusion by combining the local entropy, the global entropy and their respective weights to obtain a comprehensive entropy, where the calculation formula of the comprehensive entropy is:

[0027] .

[0028] Preferably, the calculation method of the preset regularity index of the frequency-domain feature includes:

[0029] Extract the energy standard deviation of each frequency sub-band in the frequency-domain feature matrix;

[0030] Define a reference interval according to the physiological frequency range of vascular pulsation, calculate the energy distribution ratio of the energy standard deviation within the reference interval, and use the energy distribution ratio as the regularity index.

[0031] Preferably, the determination process of the voxel set includes:

[0032] Perform linear weighting on the comprehensive entropy and the regularity index to obtain an objective optimization function;

[0033] Use the gradient descent method with an adaptive learning rate to iteratively optimize the objective optimization function, and terminate the optimization when the change amount of the function value is less than a preset threshold;

[0034] Sort all voxels according to the optimized objective function value, and screen out the voxel set with the objective function value higher than the dynamic threshold, where the dynamic threshold is determined according to the percentile of the distribution of the objective function values of all voxels.

[0035] Preferably, the determination method of the dynamic threshold includes the following steps:

[0036] Obtain the optimized objective function values of all voxels and calculate their numerical distribution;

[0037] Based on the preset percentile range, extract the corresponding objective function values from the numerical distribution as candidate dynamic thresholds;

[0038] According to the expected spatial continuity of the vascular region, by statistically analyzing the clustering density of the voxels corresponding to the candidate dynamic thresholds in three-dimensional space, the candidate threshold with the maximum clustering density is selected as the final dynamic threshold, where the calculation method of the clustering density is to statistically analyze the number of connected regions in three-dimensional space with the voxel neighborhood as the unit, and the density is inversely proportional to the number of connected regions.

[0039] Preferably, the post-processing includes performing morphological processing on the set of voxels that meet the set conditions.

[0040] The present invention also provides a data processing system applied to magnetic resonance imaging, which is used to execute a data processing method applied to magnetic resonance imaging. The system includes:

[0041] A sequence acquisition module, which is used to unfold the four-dimensional flow nuclear magnetic resonance imaging data in the time dimension to obtain the time series of the velocity field of each voxel;

[0042] A multi-scale processing module, which is used to perform multi-scale frequency domain processing on the time series of the velocity field of each voxel to obtain a frequency domain feature matrix and a frequency domain energy ratio;

[0043] A comprehensive entropy calculation module, which is used to calculate the local entropy and the global entropy based on the frequency domain feature matrix and calculate the weights of the local entropy and the global entropy in real time based on the frequency domain energy ratio, and perform weighted fusion on the local entropy and the global entropy according to their respective weights to obtain a comprehensive entropy;

[0044] A voxel set screening module, which is used to construct an objective optimization function by combining the comprehensive entropy and a preset regularity index of the frequency domain features, and use the gradient descent method to solve the objective optimization function to screen out a set of voxels that meet the set conditions;

[0045] A post-processing module, which is used to perform post-processing on the set of voxels that meet the set conditions to determine the image data of the vascular region in the four-dimensional flow nuclear magnetic resonance imaging data.

[0046] Compared with the related technology, a data processing method and system applied to magnetic resonance imaging provided by the present invention have the following beneficial effects:

[0047] The present invention first expands 4DFlowMRI data in the time dimension to obtain the time series of the velocity field for each voxel; then performs multi-scale frequency-domain processing on the time series of the velocity field to extract the frequency-domain feature matrix and the frequency-domain energy proportion; then calculates the local entropy and the global entropy based on the frequency-domain feature matrix, and adjusts the weights of the two in real time according to the frequency-domain energy proportion to obtain the comprehensive entropy; then constructs an objective optimization function in combination with the comprehensive entropy and the preset regularity index of the frequency-domain features, and uses the gradient descent method to solve it; finally, post-processes the voxel set that meets the set conditions to determine the image data of the vascular region.

[0048] The proposed invention not only overcomes the limitations of traditional 4DFlowMRI data processing methods in identifying the vascular region, but also improves the processing speed and computational efficiency, providing a more accurate and efficient tool for the diagnosis and research of cardiovascular diseases. Brief Description of the Drawings

[0049] Figure 1 It is a flowchart of a data processing method applied to magnetic resonance imaging provided by the present invention;

[0050] Figure 2 It is a module structure diagram of a data processing system applied to magnetic resonance imaging provided by the present invention. Detailed Embodiments

[0051] The present invention will be further described in detail below with reference to the drawings and embodiments. It can be understood that the specific embodiments described herein are only used to explain the present invention, rather than limiting the present invention. Additionally, it should be noted that for the sake of description, only parts related to the present invention rather than all structures are shown in the drawings. Furthermore, without conflict, the embodiments in the present invention and the features in the embodiments can be combined with each other.

[0052] It should also be noted that for the sake of description, only parts related to the present invention rather than all content are shown in the drawings. Before discussing the exemplary embodiments in more detail, it should be mentioned that some exemplary embodiments are described as processes or methods depicted as flowcharts. Although the flowcharts describe the operations (or steps) as sequential processes, many of the operations can be implemented in parallel, concurrently, or simultaneously. In addition, the order of the operations can be rearranged. The process can be terminated when its operations are completed, but it can also have additional steps not included in the drawings. The process can correspond to a method, function, procedure, subroutine, subprogram, etc.

[0053] Embodiment 1

[0054] Four-dimensional flow magnetic resonance imaging data contains three-dimensional spatial coordinates and blood flow velocity information in the time dimension. Its essence is the quantitative representation of dynamic vascular flow in space-time. Due to the pulsatile characteristics of blood flow in blood vessels being closely related to physiological activities, the velocity changes in the time dimension can reflect hemodynamic characteristics. Therefore, it is necessary to extract the velocity evolution law of each voxel (the smallest unit in three-dimensional space) over time from the original four-dimensional data. In traditional methods, three-dimensional static velocity field analysis is difficult to capture dynamic characteristics, and the direct processing of four-dimensional data has problems such as high computational complexity and sensitivity to noise. It needs to be decomposed into time series through structuring to adapt to subsequent frequency domain analysis and feature extraction. Therefore, the present invention provides a data processing method applied to magnetic resonance images.

[0055] Reference Figure 1 As shown, a data processing method applied to magnetic resonance images provided by the present invention includes the following steps:

[0056] S1: Expand the four-dimensional flow magnetic resonance imaging data according to the time dimension to obtain the velocity field time series of each voxel.

[0057] Among them, step S1 specifically includes the following steps:

[0058] S11: Decompose the four-dimensional flow magnetic resonance imaging data into three-dimensional velocity field data at multiple moments.

[0059] In this embodiment, first read the four-dimensional flow magnetic resonance imaging data file (usually in NIfTI or DICOM format) and parse its data dimension structure. The four-dimensional data consists of X, Y, Z, and T, where X, Y, and Z are three-dimensional spatial coordinates, and T is the number of time points (for example, 20 - 30 phases within a cardiac cycle). After slicing according to the time dimension, each moment corresponds to a three-dimensional velocity field data block, which contains the velocity components of each voxel in the three directions of X, Y, and Z , and . To simplify the calculation, usually the magnitude of the velocity vector is used as the velocity scalar value of the voxel, and finally T independent three-dimensional velocity field data files are generated.

[0060] S12: For each voxel of the three-dimensional velocity field data, extract the velocity values at different time points to form the corresponding velocity field time series, and perform normalization processing on the velocity field time series.

[0061] In this embodiment, step S12 further extracts the velocity values of each voxel at all time points to form a one-dimensional time series with a length of T. To improve the stability of subsequent frequency domain analysis, perform Z-score normalization processing on the time series to eliminate the differences in velocity magnitudes between different voxels.

[0062] Through the above steps, the spatio-temporal characteristics of the four-dimensional data are decoupled into a separate representation of spatial distribution and temporal evolution, which not only retains the dynamic information of vascular pulsation but also reduces the data dimension. The normalization process effectively suppresses the velocity fluctuation deviation caused by magnetic field inhomogeneity or noise, making the frequency-domain characteristics of voxels at different positions comparable, providing a standardized input for subsequent multi-scale frequency-domain analysis (such as wavelet transform), and avoiding the real blood flow signal being overwhelmed by high-frequency noise. This process takes into account both computational efficiency and feature fidelity and is the basis for subsequent entropy value fusion and vascular region screening.

[0063] S2: Perform multi-scale frequency-domain processing on the time series of the velocity field of each voxel to obtain a frequency-domain feature matrix and the proportion of frequency-domain energy.

[0064] The processes of obtaining the frequency-domain feature matrix and the proportion of frequency-domain energy include the following steps:

[0065] S21: Perform wavelet transform on the time series of the velocity field to obtain wavelet coefficients representing the joint time and frequency.

[0066] In this embodiment, the Daubechies wavelet basis is selected to perform discrete wavelet transform on the normalized time series of the velocity field, and the decomposition level is set to 4 - 6 layers to cover the physiologically relevant frequency band of 0.5 - 20 Hz. Multilevel decomposition is performed on the time series of each voxel to obtain approximation coefficients (low-frequency components) and detail coefficients (high-frequency components). For example, the detail coefficients of the first layer correspond to the highest frequency band (10 - 20 Hz, which may contain noise), and the approximation coefficients of the fourth layer correspond to the lowest frequency band (0.5 - 1.25 Hz, reflecting the cardiac cycle fluctuations). By reconstructing the wavelet coefficients of each frequency band, a time-frequency matrix is generated. To enhance the feature stability, the absolute value of the wavelet coefficients of each frequency band is taken and a moving average is performed (window length 3 to 5 time points) to suppress instantaneous noise fluctuations.

[0067] S22: Extract the energy distribution characteristics of low-frequency, medium-frequency, and high-frequency based on the wavelet coefficients to construct a frequency-domain feature matrix.

[0068] In this embodiment, based on the wavelet decomposition level (such as 4 layers), the detail coefficients (D1 - D4) and approximation coefficients (A4) of each layer are divided into three sub-bands of low-frequency, medium-frequency, and high-frequency:

[0069] Low-frequency sub-band (0.5 - 1.25 Hz): Corresponding to the approximation coefficients of the fourth layer (A4), representing the slow changes in blood flow dominated by the cardiac cycle;

[0070] Medium-frequency sub-band (1.25 - 5 Hz): Reconstructed from the detail coefficients of the third and second layers (D3, D2), reflecting local vascular pulsation and turbulence;

[0071] High-frequency sub-band (5 - 20 Hz): Reconstructed from the first-level detail coefficients (D1), may contain noise or rapid blood flow fluctuations.

[0072] For each sub-band, calculate its energy value, which are respectively:

[0073] ;

[0074] ;

[0075] .

[0076] Finally, normalize the three-frequency band energy values of each voxel into a probability distribution, and combine the probability distributions of the three-frequency band energy values of each voxel according to the spatial coordinates into a three-dimensional frequency domain feature matrix, where each channel of the matrix corresponds to the proportion of the frequency band energy.

[0077] S23: Define a reference frequency range related to blood vessel flow, and calculate the proportion of the frequency domain energy of each voxel based on the reference frequency range.

[0078] In this embodiment, define the reference frequency range based on the physiological characteristics of human blood vessel pulsation. The main frequency of arterial blood flow affected by cardiac contraction is 0.8 - 2 Hz (corresponding to a heart rate of 60 - 120 beats per minute), and the main frequency of venous return modulated by respiration is about 0.1 - 0.3 Hz. Considering clinical data comprehensively, the reference frequency range is set to 0.5 - 5 Hz to cover the main blood flow signals. For each voxel, calculate the proportion of the energy of its wavelet coefficients within the reference frequency band:

[0079] Extract all the wavelet coefficients belonging to the reference frequency band (0.5 - 5 Hz) in the time-frequency matrix to form a subset;

[0080] Calculate the total energy and the reference frequency band energy;

[0081] The proportion of the frequency domain energy is defined as the reference frequency band energy divided by the total energy, and smoothing processing is performed to reduce spatial mutations.

[0082] S3: Calculate the local entropy and the global entropy based on the frequency domain feature matrix, and calculate the respective weights of the local entropy and the global entropy in real time based on the proportion of the frequency domain energy, and perform weighted fusion on the local entropy and the global entropy according to their respective weights to obtain the comprehensive entropy.

[0083] In the processing of magnetic resonance imaging data, entropy value analysis of frequency domain features is often used to quantify the complexity and regularity of signals. However, traditional methods often use local entropy (reflecting the uncertainty of the frequency domain distribution of a single voxel) or global entropy (reflecting the statistical characteristics of the overall data set) alone, ignoring their complementarity in blood vessel region detection. Local entropy is sensitive to noise, which may lead to isolated noise voxels being misjudged as blood vessels; while global entropy overly depends on the overall distribution, which may weaken local blood flow dynamic characteristics. In addition, the blood flow frequency domain energy in the blood vessel region is usually concentrated in the frequency bands related to physiological activities (such as the cardiac pulsation frequency), but the existing methods do not dynamically associate the characteristics of the frequency domain energy distribution with the entropy value calculation, resulting in insufficient discrimination between blood vessels and surrounding tissues. Therefore, this application adopts a combined method for processing.

[0084] Among them, the calculation process of the comprehensive entropy includes:

[0085] S31: Normalize the frequency domain features of each voxel in the frequency domain feature matrix to obtain a probability distribution, and calculate the local entropy based on the probability distribution. The calculation formula of the local entropy is:

[0086] Among them, represents the local entropy, represents the voxel at different frequencies on the probability distribution.

[0087] In this embodiment, when the energy proportion of a certain frequency sub-band is significant (such as the main frequency of blood vessel pulsation), the local entropy value is low, indicating that the frequency domain distribution is concentrated; if the energy is dispersed (such as noise interference), the entropy value increases.

[0088] S32: Conduct statistical analysis on the entire frequency domain feature matrix to obtain a global probability distribution, and calculate the global entropy based on the global probability distribution. The calculation formula of the global entropy is:

[0089] Among them, represents the global entropy; represents the global probability distribution.

[0090] In this embodiment, first sum up the energy values of each frequency sub-band of all voxels to generate a global energy distribution, and then calculate the global entropy. The global entropy reflects the frequency domain consistency of the entire data set. For example, if the energy of most voxels is concentrated in the frequency bands related to blood vessels, the global entropy is low.

[0091] S33: Calculate the weights of the local entropy and the global entropy respectively according to the frequency domain energy proportion. The calculation formula of the weight is:

[0092] Among them, represents the weight of the local entropy, represents the weight of the global entropy, represents the proportion of the frequency-domain energy, represents the adjustment parameter, with a range of [0, 1], used to control the sensitivity of the weight to the proportion of the frequency-domain energy.

[0093] In this embodiment, when is larger, the weight is more dependent on the proportion of the frequency-domain energy, that is, it is more inclined to adjust the weight according to the frequency-domain energy distribution of the voxels. When is smaller, the weight is closer to a fixed value, reducing the dependence on the frequency-domain energy distribution.

[0094] S34: Perform weighted fusion by combining the local entropy and the global entropy and their respective weights to obtain the comprehensive entropy. The calculation formula of the comprehensive entropy is:

[0095] .

[0096] In this embodiment, the comprehensive entropy effectively balances the local details and the global statistical characteristics: in the vascular region, a higher proportion of the frequency-domain energy drives the weight to favor the local entropy, highlighting the concentrated frequency-domain distribution of the blood flow pulsation; in the non-vascular region, a lower proportion of the frequency-domain energy makes the weight favor the global entropy, suppressing the local entropy fluctuations caused by random noise.

[0097] S4: Construct an objective optimization function by combining the comprehensive entropy and a preset regularity index of the frequency-domain features, and use the gradient descent method to solve the objective optimization function to screen out the voxel set that meets the set conditions.

[0098] In step S4, the calculation method of the preset regularity index of the frequency-domain features includes:

[0099] First, extract the energy standard deviation of each frequency sub-band in the frequency-domain feature matrix.

[0100] In this embodiment, each voxel in the frequency-domain feature matrix contains the proportion of the energy of three sub-bands: low frequency, medium frequency, and high frequency (as defined in step S22). For each sub-band (low frequency, medium frequency, high frequency), calculate the standard deviation of its energy value over time, which is specifically expressed as follows:

[0101] Low-frequency sub-band (0.5 - 1.25 Hz): energy standard deviation .

[0102] Medium-frequency sub-band (1.25 - 5 Hz): energy standard deviation .

[0103] High-frequency sub-band (5 - 20 Hz): energy standard deviation .

[0104] The standard deviation reflects the intensity of the fluctuation of the frequency energy over time. Due to being driven by the periodic cardiac pulsation, the standard deviation of the main frequency sub-band (such as the medium frequency of 1.25 - 5 Hz) in the vascular region should be relatively low and stable, while the standard deviation of the non-vascular region (such as noise) may show irregular high fluctuations.

[0105] Then, a reference interval is defined according to the physiological frequency range of vascular pulsation, the proportion of the energy distribution of the energy standard deviation within the reference interval is calculated, and the proportion of the energy distribution is used as a regularity index.

[0106] In this embodiment, according to the physiological characteristics of vascular pulsation, the reference frequency range is set to 0.5 - 5 Hz (covering the low-frequency and medium-frequency sub-bands). Calculate the ratio of the sum of the standard deviations of the sub-bands within this range (i.e., the low frequency and the medium frequency) to the sum of the standard deviations of the full frequency band (low + medium + high) , and the specific formula is as follows:

[0107] This ratio The closer it is to 1, the more it indicates that the energy fluctuation is mainly concentrated in the vascular-related frequency bands, and the stronger the regularity.

[0108] Control example Perform three-dimensional Gaussian filtering (kernel size 3×3×3, standard deviation 1.0) to suppress the influence of isolated noise points and enhance the spatial continuity of the vascular region.

[0109] After determining the regularity index of the frequency domain features, all voxels are screened to obtain the required voxel set, and the determination process includes:

[0110] First, linearly weight the comprehensive entropy and the regularity index to obtain the target optimization function.

[0111] In this embodiment, and are the weight coefficients of the comprehensive entropy and the regularity index respectively ( + = 1, by default = 0.6, = 0.4), and are determined by cross-validation optimization.

[0112] Secondly, use the gradient descent method with an adaptive learning rate to iteratively optimize the target optimization function, and terminate the optimization when the change amount of the function value is less than the preset threshold.

[0113] In this embodiment, use the gradient descent method with an adaptive learning rate to optimize the target function: initialize the learning rate to be 0.1, calculate the gradient of the target function for each iteration, if the decreases for 3 consecutive iterations, increase the learning rate ( ← ×1.2), conversely, decrease ( ← ×0.5). The termination condition is the change amount between adjacent iterations variation .

[0114] After the iterative optimization is completed, a dynamic threshold needs to be determined to facilitate the subsequent screening of voxels. The specific steps are as follows:

[0115] Obtain the optimized objective function values of all voxels and calculate their numerical distributions.

[0116] In this embodiment, obtain the optimized objective function values of all voxels (weighted by the comprehensive entropy and the regularity index), and count their numerical distributions. To avoid the interference of extreme values, truncate the objective function values (e.g., remove the outliers in the first 1% and the last 1%) to generate a smooth probability density function.

[0117] Based on the preset percentile range, extract the corresponding objective function values from the numerical distribution as candidate dynamic thresholds.

[0118] In this embodiment, based on the preset percentile range , extract multiple candidate thresholds from the distribution, and each candidate threshold corresponds to a quantile.

[0119] According to the expected spatial continuity of the vascular region, by statistically calculating the clustering density of the voxels corresponding to the candidate dynamic thresholds in the three-dimensional space, select the candidate threshold with the maximum clustering density as the final dynamic threshold. The calculation method of the clustering density is to count the number of connected regions in the three-dimensional space with the voxel neighborhood as the unit, and the density is inversely proportional to the number of connected regions.

[0120] In this embodiment, next, for each candidate threshold, screen out the set of voxels whose objective function values ≥ the candidate threshold, and calculate its clustering density in the three-dimensional space. The definition of the clustering density is: within the 3×3×3 neighborhood centered on each voxel, count the number of connected regions (using 26-neighborhood connectivity judgment). The density is inversely proportional to the number of connected regions. The higher the density, the stronger the spatial continuity of the voxel set, and the more in line with the tubular structure characteristics of blood vessels. The higher the density, the stronger the spatial continuity of the voxel set, and the more in line with the tubular structure characteristics of blood vessels. Finally, select the candidate threshold that maximizes the density as the dynamic threshold.

[0121] Finally, sort all voxels according to the optimized objective function values, and screen out the set of voxels whose objective function values are higher than the dynamic threshold.

[0122] S5: Post-process the set of voxels that meet the set conditions to determine the image data of the vascular region in the four-dimensional flow magnetic resonance imaging data.

[0123] In this embodiment, the post-processing is to perform morphological processing on the voxel set that meets the set conditions. The specific process is as follows:

[0124] Three-dimensional closing operation to fill internal holes:

[0125] Perform a three-dimensional closing operation (first dilate and then erode) on the voxel set using a spherical structuring element (radius 2-3 voxels). The dilation operation marks the area with at least 1 vascular voxel in the voxel neighborhood (3×3×3) as a blood vessel, filling in small breaks caused by noise or threshold deviation; subsequently, the erosion operation removes the edge protrusions introduced by dilation and restores the original size of the blood vessel. For example, for the broken area at the branch of the intracranial artery, the closing operation can connect the discontinuous points with a distance less than or equal to 2 voxels.

[0126] Three-dimensional opening operation to remove isolated noise:

[0127] Perform a three-dimensional opening operation (first erode and then dilate) using the same structuring element. The erosion operation eliminates isolated noise points (single voxels without adjacent vascular voxels around), retaining continuous regions; the dilation operation restores the eroded blood vessel edges. For example, scattered noise voxels generated by respiratory movement in the liver image can be removed while retaining the main structure of the portal vein.

[0128] Multi-scale connected region screening:

[0129] Calculate the volume (number of voxels) of all connected regions and set a dynamic volume threshold:

[0130] Extract the distribution histogram of the volumes of all connected regions and take the 90th percentile as the reference threshold;

[0131] The final volume threshold ensures that small blood vessel branches with at least 10 voxels are retained.

[0132] Delete the connected regions with volumes smaller than the volume threshold to eliminate residual noise. For example, at the carotid artery bifurcation, pseudo-plaque signals with volumes smaller than 15 voxels can be removed.

[0133] Anisotropic diffusion to smooth edges:

[0134] Apply three-dimensional anisotropic diffusion filtering (Perona-Malik model), set the number of iterations to 3-5 times, and the conduction coefficient k = 0.05. This operation smooths the edges along the long axis of the blood vessel while retaining the sharpness at the branches. For example, it can correct the blurring of the blood vessel wall caused by the closing operation and improve the clarity of the coronary artery edges.

[0135] Embodiment 2

[0136] The present invention also provides a data processing system for magnetic resonance imaging, which is used to execute a data processing method for magnetic resonance imaging. Refer to Figure 2 As shown, the system includes:

[0137] A sequence acquisition module 100, configured to unfold four-dimensional flow nuclear magnetic resonance imaging data in the time dimension to obtain a time series of the velocity field of each voxel.

[0138] A multi-scale processing module 200, configured to perform multi-scale frequency domain processing on the time series of the velocity field of each voxel to obtain a frequency domain feature matrix and a frequency domain energy ratio.

[0139] A comprehensive entropy calculation module 300, configured to calculate local entropy and global entropy based on the frequency domain feature matrix, and calculate the weights of local entropy and global entropy respectively in real time based on the frequency domain energy ratio, and perform weighted fusion on local entropy and global entropy according to their respective weights to obtain a comprehensive entropy.

[0140] A voxel set screening module 400, configured to construct an objective optimization function by combining the comprehensive entropy and a preset regularity index of frequency domain features, and use the gradient descent method to solve the objective optimization function to screen out a voxel set that meets the set conditions.

[0141] A post-processing module 500, configured to perform post-processing on the voxel set that meets the set conditions to determine the image data of the blood vessel region in the four-dimensional flow nuclear magnetic resonance imaging data.

[0142] This application is described by referring to the flowcharts and / or block diagrams of methods, devices (systems) and computer program products according to the embodiments of this application. It should be understood that each process and / or block in the flowchart and / or block diagram can be implemented by computer program instructions, and the combination of processes and / or blocks in the flowchart and / or block diagram can also be implemented by computer program instructions. These computer program instructions can be provided to the processor of a general-purpose computer, a special-purpose computer, an embedded processor or other programmable data processing devices to generate a machine, so that the instructions executed by the processor of the computer or other programmable data processing devices generate a device for realizing the functions specified in one Figure 1 one process or multiple processes and / or blocks Figure 1 one block or multiple blocks.

[0143] Those of ordinary skill in the art can understand that all or part of the steps in the various methods of the above embodiments can be completed by instructing relevant hardware through a program. This program can be stored in a computer-readable storage medium, which includes read-only memory (ROM), random access memory (RAM), programmable read-only memory (PROM), erasable programmable read-only memory (EPROM), one-time programmable read-only memory (OTPROM), electrically erasable programmable read-only memory (EEPROM), compact disc read-only memory (CD-ROM), or other optical disc memories, magnetic disk memories, tape memories, or any other medium that can be used to carry or store data and is computer-readable.

[0144] It should also be noted that the term "comprising", "including" or any other variation thereof is intended to cover non-exclusive inclusion, so that a process, method, commodity or device comprising a series of elements not only includes those elements, but also includes other elements not expressly listed, or further includes elements inherent to such process, method, commodity or device. Without further limitation, an element defined by the statement "comprising one..." does not exclude the existence of additional identical elements in the process, method, commodity or device comprising the element.

Claims

1. A data processing method for magnetic resonance imaging, characterized in that: The method comprises the following steps: Expand the 4D flow MRI data according to the time dimension to obtain the velocity field time series of each voxel; Performing multi-scale frequency domain processing on the velocity field time series of each voxel to obtain a frequency domain feature matrix and a frequency domain energy ratio; Calculating local entropy and global entropy based on the frequency domain feature matrix and calculating the weights of the local entropy and the global entropy in real time based on the frequency domain energy proportion, and performing weighted fusion on the local entropy and the global entropy according to their respective weights to obtain a comprehensive entropy; Constructing a target optimization function by combining the comprehensive entropy and a preset regularity index of frequency domain features, and solving the target optimization function by using a gradient descent method to screen out a set of voxels that meet the set conditions; The calculation method of the preset regularity index of the frequency domain feature includes: Extracting the energy standard deviation of each frequency subband in the frequency domain feature matrix; Defining a reference interval according to the physiological frequency range of vascular pulsation, calculating the energy distribution ratio of the energy standard deviation within the reference interval, and using the energy distribution ratio as a regularity indicator; Post-processing is performed on the voxel set that meets the set conditions to determine image data of the blood vessel region in the four-dimensional flow magnetic resonance imaging data.

2. The data processing method for magnetic resonance imaging according to claim 1, characterized in that: The four-dimensional flow magnetic resonance imaging data is expanded according to the time dimension to obtain the velocity field time series of each voxel, including: Decomposing four-dimensional flow MRI data into three-dimensional velocity field data at multiple moments; For each voxel of the three-dimensional velocity field data, velocity values at different time points are extracted to form a corresponding velocity field time series, and the velocity field time series is normalized.

3. The data processing method for magnetic resonance imaging according to claim 2, characterized in that: The multi-scale frequency domain processing is performed on the velocity field time series of each voxel to obtain a frequency domain feature matrix and a frequency domain energy ratio, including: Performing wavelet transform on the velocity field time series to obtain wavelet coefficients jointly representing time and frequency; Extracting low-frequency, medium-frequency and high-frequency energy distribution characteristics according to the wavelet coefficients, and constructing a frequency domain feature matrix; A reference frequency range related to vascular flow is defined, and the frequency domain energy proportion of each voxel is calculated based on the reference frequency range.

4. The data processing method for magnetic resonance imaging according to claim 3, characterized in that: The calculation process of the comprehensive entropy includes: The frequency domain features of each voxel of the frequency domain feature matrix are normalized to obtain a probability distribution, and local entropy is calculated based on the probability distribution, wherein the calculation formula of the local entropy is: ； in, represents the local entropy, Representing voxels At different frequencies The probability distribution on ; Perform statistical analysis on the entire frequency domain feature matrix to obtain a global probability distribution, and calculate the global entropy based on the global probability distribution, wherein the calculation formula of the global entropy is: ； in, represents the global entropy; represents the global probability distribution; The weights of the local entropy and the global entropy are calculated according to the frequency domain energy proportion, wherein the weight calculation formula is: in, represents the weight of local entropy, represents the weight of global entropy, represents the energy ratio in the frequency domain, Represents the adjustment parameter, which is used to control the sensitivity of the weight to the frequency domain energy ratio; The local entropy and the global entropy are combined with their respective weights for weighted fusion to obtain the comprehensive entropy, where the calculation formula of the comprehensive entropy is: 。 5. The data processing method for magnetic resonance imaging according to claim 4, characterized in that: The process of determining the voxel set includes: Linearly weighting the comprehensive entropy and the regularity index to obtain a target optimization function; The target optimization function is iteratively optimized using the gradient descent method with adaptive learning rate, and the optimization is terminated when the function value change is less than the preset threshold; All voxels are sorted according to the optimized objective function values, and a set of voxels whose objective function values are higher than a dynamic threshold is screened out. The dynamic threshold is determined according to the percentile of the objective function value distribution of all voxels.

6. The data processing method for magnetic resonance imaging according to claim 5, characterized in that: The method for determining the dynamic threshold comprises the following steps: Obtain the optimization objective function values of all voxels and calculate their numerical distribution; Based on a preset percentile range, extracting a corresponding objective function value from the numerical distribution as a candidate dynamic threshold; According to the expected spatial continuity of the vascular region, the clustering density of the voxels corresponding to the candidate dynamic thresholds in three-dimensional space is counted, and the candidate threshold with the largest clustering density is selected as the final dynamic threshold. The clustering density is calculated by counting the number of connected areas in three-dimensional space in units of voxel neighborhoods, and the density is inversely proportional to the number of connected areas.

7. The data processing method for magnetic resonance imaging according to claim 1, characterized in that: The post-processing includes performing morphological processing on the voxel set that meets the set conditions.

8. A data processing system for magnetic resonance imaging, configured to execute the data processing method for magnetic resonance imaging according to any one of claims 1 to 7, characterized in that: The system comprises: A sequence acquisition module is used to expand the four-dimensional flow MRI data according to the time dimension to obtain the velocity field time series of each voxel; A multi-scale processing module is used to perform multi-scale frequency domain processing on the velocity field time series of each voxel to obtain a frequency domain feature matrix and a frequency domain energy ratio; A comprehensive entropy calculation module is used to calculate local entropy and global entropy based on the frequency domain feature matrix, and to calculate the weights of the local entropy and the global entropy in real time based on the frequency domain energy proportion, and to perform weighted fusion of the local entropy and the global entropy according to their respective weights to obtain a comprehensive entropy; a voxel set screening module, configured to construct a target optimization function by combining the comprehensive entropy and a preset regularity index of frequency domain characteristics, and solve the target optimization function using a gradient descent method to screen out a voxel set that meets the set conditions; The post-processing module is used to perform post-processing on the voxel set that meets the set conditions to determine the image data of the blood vessel area in the four-dimensional flow magnetic resonance imaging data.

Images