Microvascular image evaluation method, device and equipment and storage medium

Abstract

The invention discloses a capillary image evaluation method, device and equipment and a medium, and the method comprises the steps: obtaining IQ data of a to-be-detected region, and carrying out the beam forming of the IQ data; performing clutter filtering processing on the IQ data of each continuous data frame after beam forming, and extracting a blood flow signal; removing background noise from the blood flow signal to obtain a blood vessel structure image, and segmenting the blood vessel structure image to obtain a distribution position of a blood vessel region; screening a peak region based on the distribution position of the blood vessel region, and determining a peak center position; trajectory tracking is carried out based on the peak value center position, and all trajectories are superposed and displayed to generate a capillary image; and analyzing the generated capillary image, and extracting related parameters of the capillary. Therefore, the method is not limited by image resolution, and can be applied to microvascular images with limited resolution and super-resolution microvascular images, and vascular morphology and hemodynamic parameters can be extracted from the microvascular images with limited resolution and super-resolution.

Description

Technical Field

[0001] The present invention relates to a microvascular image evaluation method, device, equipment and storage medium, and relates to the technical field of ultrasonic imaging and image analysis, and in particular to ultrasonic imaging technology and quantitative evaluation of microvessels. Background Art

[0002] Microvascular imaging and analysis technology plays a vital role in evaluating physiological states and pathological diseases, such as tumors, chronic liver diseases, and cerebrovascular diseases. Ultrasound imaging is a safe, convenient, low-cost and widely used clinical method that can be used to detect and evaluate blood flow in vivo. In recent years, ultrasound microvascular imaging, as an emerging blood flow imaging technology, can achieve sub-millimeter-level high-resolution dynamic blood flow imaging with a wider imaging field of view (tens of millimeters) and a deeper imaging depth (tens of millimeters). Moreover, after intravenous injection of an appropriate amount of ultrasound microbubble contrast agent, super-resolution localization microscopy technology can be used to further obtain vascular structure and blood flow velocity images at the order of tens of microns. Therefore, compared with traditional imaging modes such as color Doppler imaging and CT angiography, ultrasound microvascular imaging can achieve dynamic monitoring of blood flow in small and microvessels with better blood flow sensitivity, higher imaging resolution and faster imaging speed. Further analysis of ultrasound microvascular images can more intuitively and clearly display information such as small vessel structure and microblood flow characteristics, which is convenient for effective evaluation of the physiological and pathological conditions of the imaging site.

[0003] The prior art discloses a method for evaluating the efficacy of neoadjuvant chemotherapy for breast cancer using ultrasound positioning microscopy, which extracts vascular morphological parameters (including vascular diameter, vascular spacing, vascular density, vascular density ratio and vascular tortuosity) and hemodynamic parameters (including average velocity, blood perfusion index, perfusion ratio and directional entropy) from microvascular images of breast cancer masses to evaluate the dynamic changes of microvascular reconstruction during neoadjuvant chemotherapy for breast cancer. However, this evaluation method relies on the continuous circulation perfusion of microbubble contrast agents in the blood vessels in the imaging target area, and there are problems such as insufficient microbubble perfusion, long imaging time and contraindications to contrast. The scope of application of the method is limited, and the quality of the microvascular image obtained after being restricted by external factors is significantly reduced, which is not conducive to the subsequent extraction of quantitative parameters to accurately evaluate the microvascular state.

[0004] There is a need for a microvascular imaging technique that can be widely applicable to different organs and tissues and can provide a comprehensive and systematic assessment of microvascular morphology and hemodynamics. Summary of the invention

[0005] The present invention aims to solve at least one of the technical problems existing in the prior art. To this end, in view of the above problems, the purpose of the present invention is to provide a microvascular image evaluation method, device, equipment and storage medium that can accurately extract vascular morphology and hemodynamic parameters.

[0006] In order to achieve the above-mentioned invention object, the technical solution adopted by the present invention is:

[0007] In a first aspect, the present invention provides a microvascular image assessment method, the method comprising:

[0008] Obtaining IQ data of the area to be measured and performing beam synthesis on the IQ data;

[0009] Perform clutter filtering on the IQ data of each continuous data frame after beamforming to extract blood flow signals;

[0010] The background noise of the blood flow signal is removed to obtain a vascular structure image, and the vascular structure image is segmented to obtain the distribution position of the vascular area;

[0011] The peak area is screened based on the distribution position of the blood vessel area, and the peak center position is determined;

[0012] Track the trajectory based on the peak center position, and overlay all trajectories to generate a microvascular image;

[0013] The generated microvascular images were analyzed and microvascular related parameters were extracted.

[0014] Some possible implementations further include a step of performing motion estimation on the IQ data after beamforming to obtain a plurality of continuous data frames whose relative motion meets the requirements, including:

[0015] Calculate the mutual correlation coefficient between two frames of data of all set interval frames, and select any one of the two frames of data with the largest mutual correlation coefficient as the reference stable frame;

[0016] The structural similarity between all data frames and the reference stable frame is calculated in sequence. If the structural similarity is lower than the threshold, it means that the frame data has a large amplitude and difficult to compensate motion, so it is eliminated; if the structural similarity is higher than the threshold, the continuous multi-frame data is retained as a data block;

[0017] For the reserved data blocks, global rigid motion estimation and local non-rigid motion estimation are used to obtain the motion of each frame relative to the reference stable frame.

[0018] In some possible implementations, performing clutter filtering on the IQ data of each continuous data frame after beamforming to extract the blood flow signal includes:

[0019] The IQ data is divided into multiple blocks of data with a set size and a preset overlap between them;

[0020] The singular value decomposition algorithm is used to decompose each block of data into a spatial singular matrix U, a singular value matrix S and a temporal singular matrix V;

[0021] Based on the singular value matrix S, a singular value variation curve on the diagonal of the singular value matrix is ​​obtained, and a low-order singular value threshold and a high-order singular value threshold used for filtering are determined according to the singular value variation curve, wherein the low-order singular value threshold is used to separate background tissue clutter and blood flow signals, and the high-order singular value threshold is used to separate blood flow signals and noise signals, and the singular values ​​between the low-order singular value threshold and the high-order singular value threshold are considered to be blood flow signals;

[0022] Set all singular values ​​outside the singular value threshold range in the singular value matrix S to zero, and obtain the data after clutter filtering through U*S*V;

[0023] The filtered data blocks are superimposed according to their original positions, and the signal at each spatial position is the average value of the signals of all the overlapping data blocks at that position.

[0024] In some possible implementations, background noise is removed from the blood flow signal to obtain a vascular structure image, and the vascular structure image is segmented to obtain the distribution position of the vascular region, including:

[0025] The power of the blood flow signal is accumulated along the time direction to obtain a power Doppler image, and the power Doppler image is normalized;

[0026] Processing the normalized image to obtain an estimated value of the background noise of the image;

[0027] The background noise is eliminated by subtracting the background noise estimation map from the normalized image to obtain a vascular structure image;

[0028] The local pixel threshold is determined according to the average value of a pixel neighborhood of a local set size in the vascular structure image, and the vascular structure image is segmented using the local adaptive threshold to obtain the position of the vascular area.

[0029] In some possible implementations, trajectories are tracked based on the peak center position, and all trajectories are superimposed and displayed to generate a microvascular image, including:

[0030] Pairing the center positions of multiple peak signals obtained in each frame of the image to generate a trajectory;

[0031] Unreasonable or invalid trajectories are removed based on the set maximum peak distance and shortest trajectory length.

[0032] In some possible implementations, all trajectories are superimposed and displayed to generate a microvascular image, including:

[0033] All trajectories are superimposed and displayed on a super-resolution grid, with the value of each grid node equal to the number of blood flow signal trajectories passing through the node, generating a vascular density map reflecting the microvascular structure;

[0034] The velocity trajectory is superimposed and displayed on the super-resolution grid. The value of each grid node is equal to the average velocity of all velocity trajectory segments passing through the node. Color coding is performed according to the velocity magnitude, and a linear mapping relationship is established between the velocity value and the RGB color value to obtain a blood flow velocity map reflecting the microblood flow velocity.

[0035] In some possible implementations, extracting microvascular related parameters includes:

[0036] Extracting vascular morphological parameters based on microvascular density maps, including vessel diameter, vessel density index, vessel spacing, vessel tortuosity, fractal dimension, average bifurcation number and / or average branch length;

[0037] Hemodynamic parameters including blood flow velocity, perfusion index and / or local flow uniformity are extracted based on the micro-blood flow velocity map.

[0038] In a second aspect, the present invention further includes a microvascular image evaluation device, comprising:

[0039] A beamforming unit is configured to obtain IQ data of the area to be measured and perform beamforming on the IQ data;

[0040] The blood flow signal extraction unit is configured to perform clutter filtering on the IQ data of each continuous data frame after beamforming to extract the blood flow signal;

[0041] The blood vessel region segmentation unit is configured to remove background noise from the blood flow signal to obtain a blood vessel structure image, and to segment the blood vessel structure image to obtain the distribution position of the blood vessel region;

[0042] A peak extraction unit is configured to screen the peak area based on the distribution position of the blood vessel area and determine the peak center position;

[0043] A microvascular image generation unit is configured to perform trajectory tracking based on the peak center position and to superimpose all trajectories to generate a microvascular image;

[0044] The image analysis unit is configured to analyze the generated microvessel image to obtain microvessel related parameters.

[0045] In a third aspect, the present invention further includes an electronic device, comprising:

[0046] at least one processor; and

[0047] A memory communicatively connected to the processor; wherein,

[0048] The memory stores instructions executable by the processor, and the instructions are executed by the processor to enable the processor to perform any one of the methods described above.

[0049] In a fourth aspect, the present invention also includes a computer-readable storage medium storing one or more programs, wherein the one or more programs include computer instructions, and the computer instructions are used to enable a computer to execute any one of the methods described.

[0050] The present invention adopts the above technical solution, and has the following characteristics:

[0051] 1. The present invention is not limited to microbubble contrast agents. In the absence of contrast, microvascular images can be obtained by locating red blood cells in the blood. In the presence of contrast, microvascular images with higher signal-to-noise ratio, clarity and resolution can be obtained by locating microbubble contrast agents in blood vessels.

[0052] 2. The clutter filtering of the present invention adopts a block singular value decomposition algorithm, and uses different singular value thresholds for signals in different areas of the image to adaptively extract blood flow signals, thereby avoiding the selection of singular value thresholds during global filtering being affected by the background signals of local highlighted tissues in the image and being unable to effectively remove background clutter.

[0053] 3. The motion estimation of the present invention adopts a two-step estimation method of global rigid motion and local non-rigid motion, thereby reducing the search window range during non-rigid motion estimation and improving estimation accuracy, shortening the total time of motion estimation and improving robustness.

[0054] 4. The peak location of the present invention uses the blood vessel area segmented from the power Doppler image to constrain the peak position finally located, thereby reducing blood vessel artifacts caused by erroneous positioning and improving the contrast of microvascular images.

[0055] 5. The present invention uses Kalman filtering and adaptive trajectory interpolation methods to make the trajectory smoother and more continuous during trajectory tracking, thereby enhancing the visualization of blood vessels in the microvascular density map and reducing the disorder of blood flow velocity in the microblood flow velocity map.

[0056] In summary, the present invention is not limited by image resolution, and can be applied to both resolution-limited microvascular images and super-resolution microvascular images to extract vascular morphology and hemodynamic parameters, and can be widely used in the analysis and processing of microvascular images. BRIEF DESCRIPTION OF THE DRAWINGS

[0057] Various other advantages and benefits will become apparent to those of ordinary skill in the art by reading the detailed description of the preferred embodiments below. The accompanying drawings are only for the purpose of illustrating the preferred embodiments and are not to be considered as limiting the present invention. Throughout the accompanying drawings, the same reference numerals are used to represent the same components. In the accompanying drawings:

[0058] Figure 1 Flow chart of the ultrasonic microvascular assessment method according to an embodiment of the present invention.

[0059] Figure 2 1 is a blood vessel density map before and after threshold segmentation according to an embodiment of the present invention, wherein Figure (a) is the original image and Figure (b) is the blood vessel image after threshold segmentation.

[0060] Figure 3 1 and 2 are images of the blood vessel skeleton extraction process and result according to an embodiment of the present invention, wherein (a) is a distance transformation image and (b) is a blood vessel central axis.

[0061] Figure 4 This is a diagram showing the branching points and end points of the main branch of a blood vessel and its branches according to an embodiment of the present invention.

[0062] Figure 5 1 and 2 are blood flow velocity and blood flow images according to an embodiment of the present invention, wherein FIG. (a) is a blood flow velocity diagram and FIG. (b) is a blood flow perfusion distribution diagram.

[0063] Figure 6 FIG. 2 is a diagram of an electronic device according to an embodiment of the present invention. DETAILED DESCRIPTION

[0064] It should be understood that the terms used herein are only for the purpose of describing specific example embodiments and are not intended to be limiting. Unless the context clearly indicates otherwise, the singular forms "one", "an" and "said" as used herein may also be meant to include plural forms. The terms "include", "comprise", "contain", and "have" are inclusive, and therefore specify the existence of stated features, steps, operations, elements and / or parts, but do not exclude the existence or addition of one or more other features, steps, operations, elements, parts, and / or combinations thereof. The method steps, processes, and operations described herein are not interpreted as necessarily requiring them to be performed in the specific order described or illustrated, unless the execution order is clearly indicated. It should also be understood that additional or alternative steps may be used.

[0065] Although the terms first, second, third, etc. can be used in the text to describe multiple elements, components, regions, layers and / or sections, these elements, components, regions, layers and / or sections should not be limited by these terms. These terms can only be used to distinguish an element, component, region, layer or section from another region, layer or section. Unless the context clearly indicates, terms such as "first", "second" and other numerical terms do not imply order or sequence when used in the text. Therefore, the first element, component, region, layer or section discussed below can be referred to as the second element, component, region, layer or section without departing from the teaching of the example embodiments.

[0066] For ease of description, spatially relative terms may be used herein to describe the relationship of one element or feature relative to another element or feature as shown in the figures, such as "inside", "outside", "inner side", "outer side", "below", "above", etc. Such spatially relative terms are intended to include different orientations of the device in use or operation in addition to the orientation depicted in the figures.

[0067] At present, there is a need for a microvascular imaging technology that can be widely applied to different organs and tissues, and can comprehensively and systematically evaluate the morphology and hemodynamics of microvessels. The microvascular image evaluation method, device, equipment and medium provided by the present invention include: obtaining IQ data of the area to be measured, and beamforming the IQ data; performing clutter filtering on the IQ data of each continuous data frame after beamforming to extract blood flow signals; removing background noise from the blood flow signal to obtain a vascular structure image, and segmenting the vascular structure image to obtain the distribution position of the vascular area; screening the peak area based on the distribution position of the vascular area to determine the peak center position; tracking the trajectory based on the peak center position, and superimposing all trajectories to generate a microvascular image; analyzing the generated microvascular image to extract microvascular related parameters. Therefore, the present invention is not limited to image resolution, and can be applied to both resolution-limited microvascular images and super-resolution microvascular images to extract vascular morphology and hemodynamic parameters.

[0068] The exemplary embodiments of the present invention will be described in more detail below with reference to the accompanying drawings. Although the exemplary embodiments of the present invention are shown in the accompanying drawings, it should be understood that the present invention can be implemented in various forms and should not be limited by the embodiments described herein. On the contrary, these embodiments are provided in order to enable a more thorough understanding of the present invention and to fully convey the scope of the present invention to those skilled in the art.

[0069] Embodiment 1: Figure 1 As shown, the microvascular image assessment method provided in this embodiment includes:

[0070] S1. Obtain ultrasonic radio frequency data, and obtain in-phase / quadrature (IQ) data after orthogonal demodulation of the ultrasonic radio frequency data.

[0071] In this embodiment, ultrafast plane wave imaging technology is used to collect raw channel data of the imaging target area to obtain ultrasonic radio frequency (RF) data, and the RF data is orthogonally demodulated to obtain in-phase / quadrature (IQ) data.

[0072] Furthermore, ultrasound data may be acquired using non-contrast imaging or contrast imaging, which is not limited herein.

[0073] S2. Beamforming.

[0074] In this embodiment, beam synthesis based on a delay-and-addition algorithm is performed on the IQ data offline to obtain beam-synthesized IQ data. The delay-and-addition algorithm can be implemented on a CPU or a GPU.

[0075] S3. Motion estimation.

[0076] In this embodiment, the relative motion of each frame of image data is estimated from the data before clutter filtering, and the motion estimation is divided into two steps: global rigid motion estimation and local non-rigid motion estimation. The global rigid motion estimation uses a Fourier domain phase correlation algorithm to quickly estimate the overall large-scale motion, thereby reducing the search range and processing time of subsequent local motion estimation; the local non-rigid motion estimation uses a two-dimensional normalized cross-correlation algorithm based on a local sliding window to accurately estimate the local small-scale motion.

[0077] In this embodiment, the specific implementation process of motion estimation is as follows:

[0078] A two-dimensional normalized cross-correlation algorithm is used to calculate the cross-correlation coefficient between all two frames of data with an interval of, for example, 10 frames, and any frame of the two frames of data with the largest cross-correlation coefficient is selected as a reference stable frame. The structural similarity between all data frames and the reference stable frame is calculated in turn, for example, 0.7 to 0.8 is selected as the screening threshold. If the structural similarity is lower than the threshold, it means that the frame data has a large amplitude and difficult to compensate motion, and needs to be eliminated. Continuous multiple frames of data with structural similarity higher than the threshold are retained as a data block with small motion. For these data blocks, a global rigid motion algorithm based on phase correlation is first used to estimate the motion between each frame of the image and the reference stable frame; then a sliding image window of, for example, 1mm×1mm is selected to perform sub-pixel non-rigid motion estimation based on two-dimensional cross-correlation. Each group of image windows is first corrected for rigid motion and then cross-correlated with the corresponding reference window in the reference frame. Finally, the motion of each frame of the image relative to the reference stable frame is obtained by superimposing the global rigid motion and the local non-rigid motion.

[0079] Furthermore, the calculation method of the two-dimensional normalized mutual correlation coefficient between images A and B of the same size is given by the following formula:

[0080]

[0081] Among them, A i,j and B i,j Represent the pixel values ​​of the i-th row and j-th column in image A and image B respectively, and Represent the average value of all pixels in image A and image B respectively.

[0082] S4, noise filtering.

[0083] In this embodiment, a block singular value decomposition algorithm is used to adaptively remove background clutter and noise signals in ultrasound data and extract blood flow signals. The specific implementation process is as follows:

[0084] The image data is first divided into multiple blocks of data with a set size and a preset overlap, and clutter filtering is performed separately. The singular value decomposition algorithm decomposes the original data into a spatial singular matrix U, a singular value matrix S, and a temporal singular matrix V. The singular values ​​on the diagonal of the singular value matrix S are calculated to their common logarithm with a base of 10, and the curve of the logarithm value as the singular value order changes is the singular value curve. The low-order singular value threshold and high-order singular value threshold used for filtering are determined according to the singular value change curve on the diagonal of the singular value matrix. Among them, the low-order singular value threshold is used to separate background tissue clutter and blood flow signals, and the high-order singular value threshold is used to separate blood flow signals and noise signals. The singular values ​​between the low-order singular value threshold and the high-order singular value threshold are considered to be mainly blood flow signal components. The low-order singular value threshold is selected by finding the descending inflection point of the logarithmic singular value curve, and is generally taken at about one-tenth of the number of data frames. The high-order singular value threshold is determined by finding the divergence point between the fitting straight line at the end of the logarithmic singular value curve and the singular value curve, and is generally taken to be about half of the number of data frames. All singular values ​​outside the singular value threshold range in the singular value matrix S are set to zero, and then the data after clutter filtering is obtained through U*S*V. Finally, the filtered data blocks are superimposed according to their original positions, and the signal at each spatial position is the average value of the signals of all overlapping data blocks at that location.

[0085] S5. Background noise equalization.

[0086] In this embodiment, the background noise intensity distribution is estimated from the data before clutter filtering, and then the intensity of the blood flow signal is adjusted based on the noise intensity level. After performing singular value decomposition on the entire block of data before clutter filtering, a spatial singular matrix U, a singular value matrix S, and a temporal singular matrix V are obtained, and all elements except the last element on the diagonal of the singular value matrix S are set to zero, and then an estimated value of the background noise is obtained through U*S*V.

[0087] S6. Blood vessel region segmentation.

[0088] In this embodiment, the power of the blood flow signal is accumulated along the time direction, and then the ultrafast power Doppler image is obtained after logarithmic compression. The distribution position of the blood vessel area is obtained using the adaptive threshold segmentation method. The specific implementation process is as follows:

[0089] The power of the obtained blood flow signal is accumulated along the time direction to obtain a power Doppler image, which is normalized with the maximum value of the image and logarithmically compressed. Then, the first set value is taken as the lower limit of the image intensity, and the second set value is taken as the upper limit of the image intensity. The image value is normalized again to between 0 and 1.

[0090] A top-hat transform is performed on the normalized image to obtain an estimate of the background noise of the image.

[0091] The background noise estimation map was subtracted from the normalized image to eliminate the background noise and obtain the vascular structure image.

[0092] The local pixel threshold is determined according to the average value of the pixel neighborhood of a local set size in the vascular structure image, and the local adaptive threshold is used to segment the vascular structure image to obtain the location of the vascular area for subsequent screening of the peak signal in the blood vessel.

[0093] S7, peak location.

[0094] In this embodiment, several point spread function areas with intensities close to Gaussian distribution that are far from each other are manually selected to estimate the point spread function template of Gaussian intensity distribution, and then two-dimensional normalized cross-correlation is performed within the entire image range for each frame of blood flow intensity data to screen out the microbubble or red blood cell scattering peak area. The peak area with a cross-correlation coefficient higher than a specific value is the target area, and then the local image peak detection algorithm is used to determine the center position of the peak signal, and then the area around the peak signal is taken to determine the exact position of the peak (pixel position with decimals) by weighted averaging. The peak position is corrected according to the movement of each frame of the image relative to the reference frame (the peak position before correction plus the relative motion displacement value at the position is the peak position after correction), and then the peak position that is not in the blood vessel is removed according to the blood vessel area position obtained in S6.

[0095] S8. Trajectory tracking.

[0096] In this embodiment, the specific process of trajectory tracking is:

[0097] For the multiple peak sub-pixel positions of each frame of the image obtained in the previous step, the Kuhn-Munkres algorithm is used to pair and generate trajectories: first, pair the peak positions with the closest distance in two adjacent frames, and then use the paired peak distance as the weight, based on the idea of ​​bipartite graph matching to minimize the total weight of the pairing in adjacent image frames.

[0098] Based on the set maximum peak distance and shortest trajectory length, unreasonable or invalid trajectory connections are removed. Then, the Kalman filter algorithm is used to smooth all valid trajectories, and an appropriate number of trajectory points are inserted between the peak points to make the final generated trajectory smoother and more continuous.

[0099] S9. Microvascular image generation.

[0100] In this embodiment, the specific process of generating the microvessel image is as follows:

[0101] All trajectories are superimposed and displayed on the super-resolution grid. The value of each grid node is equal to the number of blood flow signal trajectories passing through the node (among the multiple trajectories obtained after the previous step of trajectory tracking, each time a trajectory passes through the grid node, the number of trajectories counts +1), thereby generating a vascular density map reflecting the microvascular structure;

[0102] The speed at different time points on each trajectory is calculated according to the trajectory and image frame rate (the spatial distance between two adjacent trajectory points on the trajectory divided by the frame interval of the trajectory points and the image frame rate). The speed trajectory is superimposed and displayed on the super-resolution grid. The value of each grid node is equal to the average speed of all speed trajectory segments passing through the node, and color coding is performed according to the speed size. A linear mapping relationship is established between the speed values ​​within a certain range and multiple groups of RGB color values ​​to obtain a blood flow velocity map reflecting the microblood flow velocity.

[0103] S10. Analyze microvascular images to obtain relevant parameters.

[0104] In this embodiment, relevant parameters are obtained, including:

[0105] The vascular morphological parameters were extracted based on the microvessel density map, including vessel diameter, vessel density index, vessel spacing, vessel tortuosity, fractal dimension, average number of bifurcations and / or average branch length.

[0106] Hemodynamic parameters including blood flow velocity, perfusion index and / or local flow uniformity are extracted based on the micro-blood flow velocity map.

[0107] In this embodiment, a central axis extraction algorithm based on distance transformation is used to calculate the shortest distance from the pixel to the point where the surrounding pixel value is zero, and the pixel value of the point is replaced by the distance. The central axis position of the blood vessel is determined by calculating the local extreme value of the distance. The distance value of the pixel point in the blood vessel reflects the shortest distance from the point to the boundary of the blood vessel, while the pixel point outside the blood vessel is still zero after the distance transformation. The image after distance transformation is converted into a quantitative image that characterizes the relative position of the pixel point in the blood vessel, wherein the distance value of the point on the central axis represents the radius of the cross section of the blood vessel.

[0108] In this embodiment, the extracted vascular skeleton is labeled and analyzed to determine the branches and end points in the structure. This process uses a node and path labeling algorithm based on an undirected graph. After performing morphological operations to determine the branch points and end points in the skeleton, the branch points are removed and the branches of the blood vessels are marked separately. The mapping relationship between the endpoints and the skeleton nodes is detected branch by branch, and the path length of the branch is calculated. Finally, an undirected graph is constructed based on the node coordinates and the distance between nodes, thereby obtaining various quantitative information of the vascular skeleton.

[0109] The application of the microvascular image evaluation method of the present invention is described in detail below through specific embodiments.

[0110] The microvascular image evaluation method proposed in the present invention is used to continuously monitor the microvascular structure and hemodynamic changes of rat liver, including:

[0111] 1. First, a 9-week-old healthy SD rat was anesthetized by intraperitoneal injection of 0.1 mL / 100 g body weight of Shutai 50 combined with xylazine, and then respiratory anesthesia was performed with an isoflurane concentration of 0.5% and an air flow rate of 500 mL / min to maintain the depth of anesthesia.

[0112] 2. The ultrasonic data of rat liver were collected through the abdominal wall using the Vantage 256 research ultrasound data acquisition system and L20-10 linear array transducer from Verasonics, USA. The transmission center frequency was 15 MHz, the sampling frequency was 62.5 MHz, the imaging depth was 15 mm, the deflection angles of the five plane waves emitted were -4°, -2°, 0°, 2° and 4°, the pulse repetition frequency (PRF) of plane wave imaging was 5 kHz, and the effective frame rate after angle compounding was 1000 Hz.

[0113] Furthermore, for contrast-free imaging, the imaging voltage is 20 V, and 1000 frames of raw ultrasound data are continuously acquired in a single section.

[0114] Furthermore, for contrast imaging, the imaging voltage was 10 V, 0.2 mL of SuperVue microbubble solution diluted 10 times with normal saline was injected through the femoral vein for angiography, and then 6000 frames of raw ultrasound data were continuously acquired in a single section.

[0115] 3. The linear array probe emits plane waves at different angles in sequence, and the ultrasound system collects the corresponding echo signals to obtain the original channel RF data. After orthogonal demodulation of the RF data, the delay and superposition algorithm is used for beam synthesis, and then the plane wave imaging data at different deflection angles are composited frame by frame to obtain the in-phase / quadrature (IQ) data after beam synthesis.

[0116] 4. Use the global cross-correlation algorithm to detect stable frames in the continuous multi-frame data, remove the continuous data frames with out-of-plane motion in the middle part, and leave multiple continuous data frames with relatively small motion, including:

[0117] Firstly, a global rigid motion algorithm based on phase correlation is used to estimate the motion between each frame and the reference stable frame.

[0118] Then, a 1mm×1mm image window is selected by sliding to perform sub-pixel non-rigid motion estimation based on cross-correlation. Each group of image windows is first corrected for rigid motion and then cross-correlated with the corresponding reference window in the reference frame.

[0119] Finally, the motion of each frame relative to the reference stable frame is obtained by superimposing the global rigid motion on the local non-rigid motion.

[0120] 5. Use the time-space domain block singular value decomposition algorithm to perform clutter filtering on each continuous data frame. The specific process is as follows:

[0121] The image data is first divided into multiple blocks of data with a spatial size of 3mm×3mm and an overlap of 80% and clutter filtering is performed separately.

[0122] The singular value decomposition algorithm decomposes the original data into a spatial singular matrix U, a singular value matrix S, and a temporal singular matrix V. The low-order singular value threshold and the high-order singular value threshold used for filtering are determined according to the singular value change curve on the diagonal of the singular value matrix. Among them, the low-order singular value threshold is used to separate background tissue clutter and blood flow signals, and the high-order singular value threshold is used to separate blood flow signals and noise signals. The singular values ​​between the low-order singular value threshold and the high-order singular value threshold are considered to be mainly blood flow signal components. The low-order singular value threshold is selected by finding the descending inflection point of the logarithmic singular value curve, and the value is generally about one-tenth of the number of data frames. The high-order singular value threshold is determined by finding the divergence point between the fitting straight line at the end of the logarithmic singular value curve and the singular value curve, and the value is generally about one-half of the number of data frames. Finally, the filtered data blocks are superimposed according to their original positions, and the signal at each spatial position is the average value of all overlapping data block signals at that location.

[0123] 6. Accumulate the power of the obtained blood flow signal along the time direction to obtain a power Doppler image, normalize it with the maximum value of the image and perform logarithmic compression, then take -30dB as the lower limit of the image intensity and 0dB as the upper limit of the image intensity, and normalize the image value to between 0 and 1 again. Perform a top-hat transformation on the normalized image to obtain an estimated value of the background noise of the image. The structural element used in the top-hat transformation is a circle with a diameter of 3mm. This is taken as an example and is not limited to this. Use the normalized image to subtract the background noise estimation map to eliminate the background noise and obtain a vascular structure image. Determine the local pixel threshold based on the average value of the local 2mm×2mm pixel neighborhood in the vascular structure image, use the local adaptive threshold to segment the vascular structure image, and obtain the location of the vascular area for subsequent screening of the peak signal in the blood vessel.

[0124] 7. Convert blood flow IQ data into signal intensity data through envelope extraction to determine the position with the maximum intensity of the microbubble signal center. Manually select several typical point spread function areas in the image to estimate the point spread function template of Gaussian intensity distribution, and then perform cross-correlation in the entire image range for each frame of blood flow intensity data to calculate the cross-correlation coefficient of each pixel position, for example, use 0.6 as the cross-correlation coefficient threshold to screen the signal peak area. Then use the local area extreme value search algorithm to determine the center position of the peak signal, and sort the top 200 peak areas with strong signals in each frame of the image according to the maximum intensity of the peak signal, and then select the signal area of ​​0.2mm×0.2mm around the peak to determine the sub-pixel position of the peak using the weighted average method. Correct the peak position based on the movement of each frame of the image relative to the reference frame, and then remove the peak positioning points outside the vascular area based on the position of the segmented vascular area.

[0125] 8. For multiple peak positioning points in each frame, the Kuhn-Munkres algorithm is used to pair the closest peak positions in two adjacent frames, and then the paired peak distance is used as the weight. Based on the idea of ​​bipartite graph matching, the total weight of the paired points in adjacent image frames is minimized to generate trajectories. The maximum peak distance of 0.1mm and the shortest trajectory length of 5 are used to remove unqualified trajectories. Then the Kalman filter algorithm is used to smooth all valid trajectories with a smoothing window size of 5, and an appropriate number of trajectory points are inserted between the peak points to make the final generated trajectory smoother and continuous. Trajectory interpolation will ensure that the distance between adjacent trajectory points is less than 0.01mm.

[0126] 9. All trajectories are superimposed and displayed on a super-resolution grid with a pixel size of 5μm×5μm. The value of each pixel is equal to the number of blood flow signal trajectories passing through the node, thereby generating a vascular density map reflecting the microvascular structure; the speed at different time points on each trajectory is calculated according to the trajectory and image frame rate, and the speed trajectory is superimposed and displayed on a super-resolution grid with a pixel size of 5μm×5μm, and color-coded according to the speed size to obtain a blood flow velocity map reflecting the microblood flow velocity.

[0127] 10. For the microvessel density map, global threshold segmentation was performed with an image intensity threshold of 2. Figure 2 The vascular density map before and after threshold segmentation is shown. For the segmented vascular structure, the central axis of the vessel needs to be further extracted to represent the morphological direction of the vessel.

[0128] In this embodiment, a central axis extraction algorithm based on distance transformation is used. The algorithm calculates the shortest distance from the pixel to the surrounding zero pixel points pixel by pixel, and replaces the pixel value of the point with the distance. The central axis position of the blood vessel is determined by calculating the local extreme value of the distance. The distance value of the pixel point in the blood vessel reflects the shortest distance from the point to the boundary of the blood vessel, while the pixel point outside the blood vessel is still zero after the distance transformation. The image after distance transformation is converted into a quantitative image that characterizes the relative position of the pixel point in the blood vessel, wherein the distance value of the point on the central axis represents the radius of the cross section of the blood vessel. Figure 3 The vascular structure image and the extracted skeleton after distance transformation are shown.

[0129] In this embodiment, the extracted vascular skeleton is labeled and analyzed to determine the branches and end points in the structure. This process uses a node and path labeling algorithm based on an undirected graph. After performing morphological operations to determine the branch points and end points in the skeleton, the branch points are removed and the branches of the blood vessels are marked separately. The mapping relationship between the endpoints and the skeleton nodes is detected branch by branch, and the path length of the branch is calculated. Finally, an undirected graph is constructed based on the node coordinates and the distance between nodes, thereby obtaining various quantitative information of the vascular skeleton. Figure 4 The main vascular skeleton and its nodes used for labeling analysis are shown.

[0130] Furthermore, after the above analysis and processing, the following morphological parameters can be obtained from the processing results:

[0131] (1) Vascular density index: The vascular density index can be obtained by calculating the ratio of the number of pixels in the vascular area after threshold segmentation to the total number of pixels in the image. It is used to reflect the distribution number and density of blood vessels in the imaging field of view.

[0132] (2) Blood vessel diameter: Twice the distance transformation value of a point on the central axis of a blood vessel is the diameter of the cross section of the blood vessel at that location, which is used to distinguish between large and small blood vessels and facilitate subsequent grading analysis of blood vessels.

[0133] (3) Vascular tortuosity: By selecting any two nodes on an undirected graph, the total length of the shortest path between the nodes can be calculated by finding the shortest path between the nodes, and then the Euclidean distance between the two nodes is calculated. The ratio of the path length to the Euclidean distance is the tortuosity of the blood vessel between the two nodes, which reflects the course and structural characteristics of the blood vessel.

[0134] (4) Fractal dimension: Using the box dimension method, the image is placed on a uniform grid and the number of grids required to cover the vascular structure is calculated. As the grid becomes denser, the size of a single grid decreases, and the ratio of the number of grids required to cover the vascular structure to the total number of grids gradually increases. The ratio of the logarithm of this ratio to the logarithm of the grid size will approach the fractal dimension of the blood vessel, which is used to reflect the complexity of vascular branching.

[0135] (5) Vascular spacing: For each point on the vascular skeleton, the shortest Euclidean distance from the point to other branches of the vessel is calculated, which is the shortest distance between the vessel at that point and the surrounding vessels. Then, the average vascular spacing of the imaging site is obtained by taking the statistical average of all vascular skeleton points, which is used to reflect the density of blood vessels in the imaging field of view.

[0136] (6) Average branch length: After the primary blood vessels are determined in the vascular image by means of vascular diameter or manual identification, the vascular branches at each level can be marked according to the step-by-step bifurcation of the blood vessels. The length of the vascular path of each level of vascular branch is counted and averaged to obtain the average length of the blood vessels at that level, which reflects the structural characteristics of the blood vessels.

[0137] (7) Average number of bifurcations: After determining the primary blood vessel, mark the bifurcation points along the primary blood vessel one by one. Count the number of vascular branches branching out of each bifurcation point and average them to get the average number of bifurcations at that level, which reflects the branching status of each branch of the blood vessel.

[0138] Furthermore, the following hemodynamic parameters can be extracted based on the microblood flow velocity map and the above morphological analysis:

[0139] (1) Blood flow velocity: The blood flow velocity of a certain blood vessel is directly obtained from the image data. After the region of interest is circled, the mean, median, standard deviation and other statistical parameters of the blood flow velocity in the region can be calculated.

[0140] (2) Local flow uniformity: The blood flow direction of a certain blood vessel is directly obtained from the image data. After the region of interest is circled, the entropy of the direction distribution is calculated. This parameter indicates the degree of turbulence in the velocity direction. The smaller the entropy value, the more consistent the flow velocity uniformity.

[0141] (3) Perfusion index: The perfusion index of each part of the blood vessel is calculated by multiplying the blood vessel diameter with the average blood flow velocity in the cross section along the vascular skeleton, such as Figure 5 This parameter reflects the blood transport capacity of the blood vessels and the perfusion of the tissues supplied by the blood vessels.

[0142] Embodiment 2: Embodiment 1 above provides a microvascular image evaluation method, and correspondingly, this embodiment provides a microvascular image evaluation device. The device provided in this embodiment can implement the microvascular image evaluation method of embodiment 1, and the device can be implemented by software, hardware, or a combination of software and hardware. For the convenience of description, this embodiment is described by dividing the functions into various units and describing them separately. Of course, the functions of each unit can be implemented in the same or more software and / or hardware during implementation. For example, the device may include integrated or separate functional modules or functional units to execute the corresponding steps in each method of embodiment 1. Since the device of this embodiment is basically similar to the method embodiment, the description process of this embodiment is relatively simple, and the relevant parts can refer to the partial description of embodiment 1. The embodiment of the microvascular image evaluation device provided by the present invention is only schematic.

[0143] Specifically, the present invention also includes a microvascular image evaluation device, comprising:

[0144] A beamforming unit is configured to obtain IQ data of the area to be measured and perform beamforming on the IQ data;

[0145] The blood flow signal extraction unit is configured to perform clutter filtering on the IQ data of each continuous data frame after beamforming to extract the blood flow signal;

[0146] The blood vessel region segmentation unit is configured to remove background noise from the blood flow signal to obtain a blood vessel structure image, and to segment the blood vessel structure image to obtain the distribution position of the blood vessel region;

[0147] A peak extraction unit is configured to screen the peak area based on the distribution position of the blood vessel area and determine the peak center position;

[0148] A microvascular image generation unit is configured to perform trajectory tracking based on the peak center position and to superimpose all trajectories to generate a microvascular image;

[0149] The image analysis unit is configured to analyze the generated microvessel image to obtain microvessel related parameters.

[0150] Embodiment 3: This embodiment provides an electronic device corresponding to the microvascular image evaluation method provided in Embodiment 1. The electronic device may be an electronic device for a client, such as a mobile phone, a laptop computer, a tablet computer, a desktop computer, etc., to execute the method of Embodiment 1.

[0151] like Figure 6As shown, the electronic device includes a processor, a memory, a communication interface and a bus, and the processor, the memory and the communication interface are connected through the bus to complete mutual communication. The bus can be an Industry Standard Architecture (ISA) bus, a Peripheral Component Interconnect (PCI) bus or an Extended Industry Standard Component (EISA) bus, etc. The memory stores a computer program that can be run on the processor, and the processor executes the method of embodiment 1 when running the computer program. The implementation principle and technical effect are similar to those of embodiment 1, and will not be repeated here. Those skilled in the art can understand that Figure 6 The structure shown in the figure is only a block diagram of a part of the structure related to the scheme of the present application, and does not constitute a limitation on the computing device to which the scheme of the present application is applied. The specific computing device may include more or fewer components than shown in the figure, or combine certain components, or have a different arrangement of components.

[0152] In a preferred embodiment, the logic instructions in the above-mentioned memory can be implemented in the form of a software functional unit and can be stored in a computer-readable storage medium when it is sold or used as an independent product. Based on this understanding, the technical solution of the present application can be essentially or partly embodied in the form of a software product that contributes to the prior art. The computer software product is stored in a storage medium, including several instructions to enable a computer device (which can be a personal computer, a server, or a network device, etc.) to perform all or part of the steps of the method described in each embodiment of the present application. The aforementioned storage medium includes: U disk, mobile hard disk, read-only memory (ROM, Read-Only Memory), random access memory (RAM, Random Access Memory), optical disk and other media that can store program codes.

[0153] In a preferred embodiment, the processor may be a central processing unit (CPU), a digital signal processor (DSP) or other general-purpose processors of various types, which are not limited here.

[0154] Embodiment 4: This embodiment provides a computer-readable storage medium storing one or more programs, wherein the one or more programs include computer instructions. When the computer instructions are executed by a computer, the computer executes the method provided in the above-mentioned embodiment 1.

[0155] Embodiment 5: This embodiment provides a computer program product. The computer program product may include a computer program stored on a computer-readable storage medium. The computer program includes program instructions. When the program instructions are executed by a computer, the computer can execute the method provided in the above-mentioned embodiment 1. Its implementation principle and technical effects are similar to those of embodiment 1 and will not be repeated here.

[0156] In a preferred embodiment, the computer-readable storage medium may be a tangible device that holds and stores instructions used by the instruction execution device, such as but not limited to an electrical storage device, a magnetic storage device, an optical storage device, an electromagnetic storage device, a semiconductor storage device, or any combination thereof. The computer-readable storage medium stores computer program instructions that cause the computer to execute the method provided in the first embodiment.

[0157] The present application is described with reference to the flowcharts and / or block diagrams of the methods, devices (apparatus), and computer program products according to the embodiments of the present application. It should be understood that each process and / or box in the flowchart and / or block diagram, as well as the combination of the processes and / or boxes in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to generate a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowchart and / or block diagram. Figure 1 A process or multiple processes and / or boxes Figure 1 A device that provides the functions specified in a block or multiple blocks.

[0158] These computer program instructions may also be stored in a computer readable memory capable of directing a computer or other programmable data processing device to operate in a specific manner, so that the instructions stored in the computer readable memory produce an article of manufacture including an instruction device, which implements the process Figure 1 A process or multiple processes and / or boxes Figure 1 A function specified in one or more boxes.

[0159] These computer program instructions can also be loaded onto a computer or other programmable data processing device so that a series of operating steps are executed on the computer or other programmable device to produce a computer-implemented process, thereby providing instructions for implementing the process. Figure 1 A process or multiple processes and / or boxes Figure 1 The steps for the functions specified in one or more boxes.

[0160] Each embodiment in this specification is described in a progressive manner, and the same and similar parts between the embodiments can be referred to each other, and each embodiment focuses on the differences from other embodiments. In the description of this specification, the description of reference terms "a preferred embodiment", "further", "specifically", "in the present embodiment", etc. means that the specific features, structures, materials or characteristics described in conjunction with the embodiment or example are included in at least one embodiment or example of the embodiment of this specification. In this specification, the schematic representation of the above terms does not necessarily refer to the same embodiment or example. Moreover, the specific features, structures, materials or characteristics described can be combined in any one or more embodiments or examples in a suitable manner. In addition, those skilled in the art can combine and combine the different embodiments or examples described in this specification and the features of the different embodiments or examples without contradiction.

[0161] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or make equivalent replacements for some of the technical features therein. However, these modifications or replacements do not deviate the essence of the corresponding technical solutions from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A microvascular image assessment method, characterized in that: The method includes: Obtaining IQ data of the area to be measured and performing beam synthesis on the IQ data; Perform clutter filtering on the IQ data of each continuous data frame after beamforming to extract blood flow signals; The background noise of the blood flow signal is removed to obtain a vascular structure image, and the vascular structure image is segmented to obtain the distribution position of the vascular area; The peak area is screened based on the distribution position of the blood vessel area, and the peak center position is determined; Track the trajectory based on the peak center position, and overlay all trajectories to generate a microvascular image; The generated microvascular images were analyzed and microvascular related parameters were extracted.

2. The microvascular image evaluation method according to claim 1, characterized in that: The method also includes the steps of performing motion estimation on the IQ data after beamforming to obtain a plurality of continuous data frames whose relative motion meets the requirements, including: Calculate the mutual correlation coefficient between two frames of data of all set interval frames, and select any one of the two frames of data with the largest mutual correlation coefficient as the reference stable frame; The structural similarity between all data frames and the reference stable frame is calculated in sequence. If the structural similarity is lower than the threshold, it means that the frame data has a large amplitude and difficult to compensate motion, so it is eliminated; if the structural similarity is higher than the threshold, the continuous multi-frame data is retained as a data block; For the reserved data blocks, global rigid motion estimation and local non-rigid motion estimation are used to obtain the motion of each frame relative to the reference stable frame.

3. The microvascular image evaluation method according to claim 1, characterized in that: Perform clutter filtering on the IQ data of each continuous data frame after beamforming to extract the blood flow signal, including: The IQ data is divided into multiple blocks of data with a set size and a preset overlap between them; The singular value decomposition algorithm is used to decompose each block of data into a spatial singular matrix U, a singular value matrix S and a temporal singular matrix V; Based on the singular value matrix S, a singular value variation curve on the diagonal of the singular value matrix is ​​obtained, and a low-order singular value threshold and a high-order singular value threshold used for filtering are determined according to the singular value variation curve, wherein the low-order singular value threshold is used to separate background tissue clutter and blood flow signals, and the high-order singular value threshold is used to separate blood flow signals and noise signals, and the singular values ​​between the low-order singular value threshold and the high-order singular value threshold are considered to be blood flow signals; Set all singular values ​​outside the singular value threshold range in the singular value matrix S to zero, and obtain the data after clutter filtering through U*S*V; The filtered data blocks are superimposed according to their original positions, and the signal at each spatial position is the average value of the signals of all the overlapping data blocks at that position.

4. The microvascular image evaluation method according to claim 1, characterized in that: The background noise of the blood flow signal is removed to obtain the vascular structure image, and the vascular structure image is segmented to obtain the distribution position of the vascular area, including: The power of the blood flow signal is accumulated along the time direction to obtain a power Doppler image, and the power Doppler image is normalized; Processing the normalized image to obtain an estimated value of the background noise of the image; The background noise is eliminated by subtracting the background noise estimation map from the normalized image to obtain a vascular structure image; The local pixel threshold is determined according to the average value of a pixel neighborhood of a local set size in the vascular structure image, and the vascular structure image is segmented using the local adaptive threshold to obtain the position of the vascular area.

5. The microvascular image evaluation method according to claim 1, characterized in that: Track the trajectory based on the peak center position, and overlay all trajectories to generate microvascular images, including: Pairing the center positions of multiple peak signals obtained in each frame of the image to generate a trajectory; Unreasonable or invalid trajectories are removed based on the set maximum peak distance and shortest trajectory length.

6. The microvascular image evaluation method according to claim 1, characterized in that: All trajectories are superimposed to generate microvascular images, including: All trajectories are superimposed and displayed on a super-resolution grid, with the value of each grid node equal to the number of blood flow signal trajectories passing through the node, generating a vascular density map reflecting the microvascular structure; The velocity trajectory is superimposed and displayed on the super-resolution grid. The value of each grid node is equal to the average velocity of all velocity trajectory segments passing through the node. Color coding is performed according to the velocity magnitude, and a linear mapping relationship is established between the velocity value and the RGB color value to obtain a blood flow velocity map reflecting the microblood flow velocity.

7. The microvascular image evaluation method according to claim 6, characterized in that: Extract microvascular related parameters, including: Extracting vascular morphological parameters based on microvascular density maps, including vessel diameter, vessel density index, vessel spacing, vessel tortuosity, fractal dimension, average bifurcation number and / or average branch length; Hemodynamic parameters including blood flow velocity, perfusion index and / or local flow uniformity are extracted based on the micro-blood flow velocity map.

8. A microvascular image evaluation device, characterized in that: The device includes: A beamforming unit is configured to obtain IQ data of the area to be measured and perform beamforming on the IQ data; The blood flow signal extraction unit is configured to perform clutter filtering on the IQ data of each continuous data frame after beamforming to extract the blood flow signal; The blood vessel region segmentation unit is configured to remove background noise from the blood flow signal to obtain a blood vessel structure image, and to segment the blood vessel structure image to obtain the distribution position of the blood vessel region; A peak extraction unit is configured to screen the peak area based on the distribution position of the blood vessel area and determine the peak center position; A microvascular image generation unit is configured to perform trajectory tracking based on the peak center position and to superimpose all trajectories to generate a microvascular image; The image analysis unit is configured to analyze the generated microvessel image to obtain microvessel related parameters.

9. An electronic device, characterized in that: include: at least one processor; as well as A memory communicatively connected to the processor; wherein, The memory stores instructions executable by the processor, and the instructions are executed by the processor to enable the processor to perform the method according to any one of claims 1 to 7.

10. A computer-readable storage medium storing one or more programs, characterized in that: The one or more programs include computer instructions for causing a computer to execute the method according to any one of claims 1-7.

Images