Apparatus and method for estimating velocity field

Abstract

The invention relates to a method of determining a two-dimensional characteristic of a propagating displacement wave in a patient, the method comprising: estimating a plurality of velocity diagrams of a target region over time, filtering (30) the propagating displacement wave within the plurality of velocity diagrams to obtain a filtered plurality of velocity diagrams, and determining the two-dimensional characteristic of the propagating displacement wave within the filtered plurality of velocity diagrams. And determining two-dimensional characteristics of a transmissible displacement wave of each point of interest A in the plurality of velocity diagrams, the two-dimensional characteristics comprising a propagation velocity and a propagation direction of the transmissible displacement wave to obtain at least one two-dimensional characteristic diagram, the point of interest A having the same coordinates in the plurality of velocity diagrams.

Description

Technical Field

[0001] The present invention relates to a method for determining two-dimensional (2D) properties of propagating displacement waves in a patient.

[0002] More particularly, the present invention relates to methods for measuring properties of biological tissue of interest.

[0003] It is applicable to, but not limited to, measuring viscoelastic parameters of the liver or spleen of humans or animals, which measurements are correlated with the degree of fibrosis present in the liver. Background Art

[0004] Shear wave elastography is a well-known technique for measuring the viscoelastic properties (eg, stiffness) of a medium (eg, tissue, organ, or other sample).

[0005] This technique involves measuring the propagation velocity of shear waves in a medium, which is directly related to the viscoelastic properties of the medium being analyzed.

[0006] Shear wave elastography uses ultrasound radiation force (acoustic radiation force pulse technique) or mechanical vibrators (vibration controlled transient elastography VCTE) to TM technology, or the technology described in document WO2022084502) causes medium displacement.

[0007] In both cases, the propagation direction of the shear wave front is known a priori, and the induced shear wave propagation is analyzed in an analysis direction chosen a priori to be perpendicular to the assumed shear wave front.

[0008] Shear waves are also generated naturally by living structures such as the heart, arteries, veins, and muscles.

[0009] 1. Generation of shear waves by acoustic radiation force

[0010] In the case of shear wave generation by acoustic radiation force, shear waves are generated using acoustic radiation force impulses (ARFI) transmitted as push pulses: ultrasonic energy is transmitted to a focal area of ​​the medium to generate shear waves, causing the medium to displace around the focal area.

[0011] An ultrasound scanning step is performed to track the displacement of the medium over time within a region of interest (ROI) contained within the medium. Specifically, an ultrasonic compression wave is emitted at a faster rate after generating a shear wave (see document WO0055616). This ultrasound scanning step enables a series of "echometry maps" of the medium to be obtained during the propagation of the shear wave. Specifically, successive images of incremental displacement are determined by correlating successive echometry maps.

[0012] A processing step is then performed on these incremental displacement images to determine the properties of the shear wave (e.g., propagation velocity along the analysis direction). These properties represent the specific viscoelastic properties of the medium.

[0013] However, the shear wave characteristics inferred from the displacements may be biased if the analyzed direction of shear wave propagation differs from the actual direction of shear wave propagation. In fact, shear waves may be refracted or reflected due to, for example, inhomogeneities in the medium.

[0014] This means that the processing steps will not be performed in the correct direction, leading to an overestimation of the shear wave propagation velocity.

[0015] 2. Generation of shear waves by mechanical excitation

[0016] Vibration-controlled transient elastography (VCTE) TM ) or the technology described in WO 2022 / 084502, the shear wave source is mechanically coupled to the patient's skin. When the shear wave source is driven, the vibrations generated by the shear wave source on the patient's skin at a low frequency (typically between 30 Hz and 200 Hz) can excite shear waves that propagate within the medium. The shear wave propagation direction is assumed to be perpendicular to the interface between the skin and the shear wave source.

[0017] In the subsequent ultrasound scanning step, ultrasound pulses are emitted to track the displacement of the medium caused by this low-frequency periodic vibration along one or more scan lines that should be parallel to the shear wave propagation direction (assuming it is known a priori). A one-dimensional displacement along the scan line is obtained. This displacement is the projection of the displacement caused by the propagation of the shear wave at different points within the ROI at the same time along the scan line.

[0018] A processing step of these scan lines is performed to determine the shear wave propagation velocity, given the given assumption that the propagation direction of the shear wave is parallel to the scan lines. The shear wave propagation velocity is used to derive the specific viscoelasticity of the medium.

[0019] However, there are cases where the assumption about the shear wave propagation direction is incorrect, causing the determined shear wave propagation velocity to be overestimated.

[0020] For example, when imaging a patient's organ (e.g., the liver), the shear wave source must be applied to the patient's skin at a location directly above the intercostal space between the patient's ribs. This can cause the undesirable effect of changing the propagation direction of the shear wave (mechanical coupling effects and shear wave shielding effects, discussed in more detail below). If this change in propagation direction is not taken into account in the processing steps, the shear wave propagation velocity may be overestimated.

[0021] 3. Natural generation of shear waves

[0022] Finally, a technique called "passive elastography" relies on shear waves generated naturally by the body, for example, by the vibrations caused by a beating heart.

[0023] In passive elastography, the propagation direction of the shear wave may not be known a priori.

[0024] Therefore, any solution such as those described above for ARFI and transient elastography, which makes an assumption about the shear wave propagation direction before performing the processing steps, will lead to biased estimates.

[0025] 4. Object of the invention

[0026] The object of the present invention is to propose a method for determining 2D characteristics of a propagating displacement wave that overcomes at least one of the above-mentioned drawbacks.

[0027] More specifically, the object of the present invention is to provide a method for determining the 2D characteristics of a propagating displacement wave without requiring prior knowledge of the propagation direction of said propagating displacement wave. Summary of the Invention

[0028] To this end, the present invention proposes a method for determining 2D characteristics of propagating displacement waves in a patient, the method comprising:

[0029] - Estimating multiple velocity images / maps of the target area over time, the estimation stage includes the following steps:

[0030] o Construct multiple echogenicity maps by repeatedly performing the following sub-steps for each echogenicity map:

[0031] ■ using a transducer element array comprising a group of transducer elements to transmit an ultrasound signal,

[0032] ■ The backscattered echo signal is received by the transducer element array. Each transducer element can collect each time domain signal in a set of time domain signals. The time domain signal corresponds to the amplitude of the received backscattered echo signal on the transducer element array.

[0033] ■Build an echo intensity map based on a set of time domain signals,

[0034] ○ Estimation of multiple velocity images / maps by correlating and processing multiple echo intensity maps,

[0035] - filtering the propagating displacement waves within the plurality of velocity images / maps to obtain a plurality of filtered velocity images / maps,

[0036] - determining 2D characteristics of a propagating displacement wave for each point of interest A within a plurality of velocity images / maps, the two-dimensional characteristics comprising a propagation velocity and a propagation direction of the propagating displacement wave, to obtain at least one 2D characteristic map, wherein the point of interest A has the same coordinates within the plurality of velocity images / maps:

[0037] o associating at least two points B, C within the plurality of velocity images / maps such that points A, B, C are not collinear, the at least two points B, C having the same coordinates within the plurality of velocity images / maps,

[0038] ○ For each of three points A, B, C (point of interest A and at least two points B, C), extract the time domain propagation displacement signal from the filtered multiple velocity images / maps,

[0039] o Estimate the first and second apparent flight times (specifically, the “longitudinal” apparent flight time between point A and point B, and the “transverse” apparent flight time between point A and point C) of the propagating displacement wave between the three points A, B, and C (point A of interest and at least two points B, C) using the time-domain propagating displacement signals associated with the three points A, B, and C (point A of interest and at least two points B, C) in the following manner:

[0040] Cross-correlating the time-domain propagation displacement signals associated with points A and B (point A of interest and one (B) of at least two points B and C) to estimate a first apparent flight time (longitudinal apparent flight time) of the propagation displacement wave from point A to point B (from point A of interest to one (B) of at least two points B and C),

[0041] Cross-correlating the time-domain propagative displacement signals associated with points A and C (point A of interest and the other (C) of the at least two points B and C) to estimate a second apparent flight time (lateral apparent flight time) of the propagative displacement wave from point A to point C (from point A of interest to the other (C) of the at least two points B and C),

[0042] o Determine the propagation direction and propagation velocity of the propagating displacement wave based on the first and second apparent flight times.

[0043] In the context of the present invention, the expression "propagating displacement waves" refers to shear waves generated naturally (ie, passive elastography) or artificially (ie, transient elastography).

[0044] Preferred but non-limiting aspects of the method according to the invention are as follows:

[0045] - the method may further comprise the step of estimating the stiffness of the tissue comprised within the target area using the 2D characteristic map;

[0046] - The steps for estimating the stiffness of the tissue may include:

[0047] o extracting the velocity of the propagating displacement wave from the 2D characteristic map for a subset of points representing tissue included in the target region,

[0048] ○ Calculate the median or average value of the velocity of the propagating displacement wave for the extracted subset of points,

[0049] ○ Estimate the hardness based on the calculated median or average value;

[0050] - The steps for estimating the stiffness of the tissue may include:

[0051] ○Calculate the hardness map of the region of interest based on the 2D characteristic map,

[0052] o extracting the hardness of a subset of points from the hardness map, said subset of points representing tissue included in the target region, and

[0053] ○ Calculate the median or average hardness of the extracted point subset;

[0054] - The method may further comprise the step of determining a quality factor associated with the mean propagative displacement direction, said step of determining the quality factor comprising:

[0055] ○ Extract the direction of the propagating displacement wave of a subset of points from the 2D characteristic map,

[0056] ○ Calculate the average propagation direction of the propagation displacement wave based on the direction of the propagation displacement wave of the extracted point subset,

[0057] ○ comparing the average propagation direction of the propagating displacement wave with a predefined range of directions to obtain comparative data, and

[0058] ○Derive quality factor based on comparative data;

[0059] - The method may further comprise the step of determining a coefficient representative of the pressure applied to the patient by the transducer array, said step of determining a coefficient comprising the following sub-steps:

[0060] ○Calculate the curvature of the propagating displacement wave front based on the 2D characteristic diagram,

[0061] If the calculated curvature is a concave shape corresponding to converging propagating displacement waves, the coefficient is assigned a value representing light pressure,

[0062] ○ If the calculated curvature is a convex shape corresponding to a diverging propagating displacement wave,

[0063] Then the coefficient is assigned a value representing a high pressure greater than a light pressure;

[0064] The method may further comprise the step of determining a parameter representative of the quality of the propagated 2D displacement estimate, said step of determining the parameter comprising the following sub-steps:

[0065] Calculate the amplitude of the change in the apparent transverse flight time (the amplitude of the change in the flight time between points A and C), which represents the spatial variation in the direction of the propagating displacement wave.

[0066] o comparing the amount with a predefined threshold value and determining a parameter based on the result of the comparison sub-step;

[0067] - The step of estimating a plurality of velocity images / maps by correlating and processing a plurality of complex echo intensity maps comprises:

[0068] ○ Cross-correlate multiple complex echo intensity maps with a certain time lag to obtain multiple phase shift maps,

[0069] ○Derivation of multiple velocity maps based on the associated multiple phase shift maps;

[0070] - The sub-step of constructing an echogenicity map from a set of time domain signals comprises solving an inverse problem to generate said echogenicity map.

[0071] The present invention also describes a method for estimating shear wave velocity V in biological tissue, the method comprising:

[0072] ● Detecting shear waves excited in biological tissue by a shear wave source, where the local characteristics of the shear waves are manifested as shear wave fronts propagating along the propagation direction,

[0073] ● Calculate the time domain displacement signal representing the propagation of the shear wave front within the region of interest of the biological tissue,

[0074] Determine the maximum cross-correlation time of the shear wave time-domain displacement between at least three non-collinear points located at different positions within the region of interest:

[0075] o a first point and a second point extending along a first line (assumed to be roughly parallel to the direction of the shear wave front) within the region of interest,

[0076] o a first point and a third point extending along a second line within the region of interest (assumed to be approximately perpendicular to the direction of propagation of the shear wave), and

[0077] • Determine the local shear wave velocity V based on the time.

[0078] Preferred but non-limiting aspects of the method according to the invention are as follows:

[0079] The step of determining the shear wave velocity comprises the following sub-steps:

[0080] o Calculate the first apparent shear wave front velocity V1 along a first line assumed to be roughly parallel to the propagation direction of the shear wave,

[0081] o Calculate a second apparent shear wave front velocity V2 along a second line assumed to be approximately perpendicular to the direction of propagation of the shear wave,

[0082] The step of determining the local shear wave velocity V comprises the sub-step of estimating said local shear wave velocity V based on the first shear wave front velocity V1 and the second shear wave front velocity V2;

[0083] Advantageously:

[0084] To calculate the first apparent shear wave front velocity V1, the method includes:

[0085] ■ selecting a first point A and a second point B located at different positions along the first line, the first point A being closer to the transducer array than the second point B, and deriving a first distance between the first point A and the second point B,

[0086] ■ Determine the first apparent propagation time of the shear wave front between the first point A and the second point B using the time domain propagation displacement signal,

[0087] ■ deriving the first apparent shear wave front velocity V1 by determining the ratio of the first distance divided by the first propagation time,

[0088] ○ To calculate the second apparent shear wave front velocity:

[0089] ■ selecting a third point C located along the second line and deriving a second distance between the first point A and the third point C,

[0090] ■ Determine the second apparent propagation time of the shear wave front between the first point A and the third point C using the time domain propagation displacement signal,

[0091] ■ The second apparent shear wave front velocity V2 is derived by determining the ratio of the second distance divided by the second propagation time.

[0092] The present invention also relates to a process for assisting a user in positioning a probe and controlling pressure on a surface of an anatomical structure to be imaged, the probe comprising a transducer array (T1-T2) for imaging the anatomical structure. n ) and at least one vibrator for generating shear waves through anatomical structures, the process comprising implementing the above method. BRIEF DESCRIPTION OF THE DRAWINGS

[0093] The present invention may be more fully understood upon consideration of the following detailed description of various embodiments in conjunction with the accompanying drawings, in which:

[0094] · Figure 1is a schematic diagram of the stages implemented in a method for determining 2D properties (propagation velocity and propagation direction) of a propagating displacement wave in a patient,

[0095] · Figure 2 Is used to implement Figure 1 A schematic diagram of the processing components of the illustrated method,

[0096] · Figure 3 yes Figure 1 A schematic diagram of the steps implemented within the tracking phase is shown,

[0097] · Figure 4 yes Figure 1 A schematic diagram of the steps performed within the processing stage is shown,

[0098] · Figure 5 is a schematic diagram showing the radiation pattern of shear waves generated by a point source,

[0099] · Figure 6 It is shown in Schematic diagram of the ultrasonic pulse elastography probe implemented in the system,

[0100] · Figure 7 is Schematic diagram of the radiation pattern of the shear waves generated by the system,

[0101] · Figure 8 is a schematic diagram showing an ultrasonic pulse elastography probe implemented in a Hepatoscope,

[0102] · Figure 9 is a schematic diagram showing the shear wave generation area on the ultrasound probe of the Hepatoscope,

[0103] · Figure 10a and Figure 10b When used separately Schematic diagram of the shear wave shielding effect in the plane perpendicular to the patient's ribs when using a hepatoscope.

[0104] · Figure 11a and Figure 11b When used separately Schematic diagram of the shear wave shielding effect in a plane parallel to the patient's ribs when using a Hepatoscope.

[0105] · Figure 12a and Figure 12b When used separately Schematic diagram of the mechanical coupling effect in the plane perpendicular to the patient's ribs when using a Hepatoscope.

[0106] · Figure 13a andFigure 13b When used separately Schematic diagram of the mechanical coupling effect in a plane parallel to the patient's ribs when using a Hepatoscope.

[0107] · Figure 14 is a geometric diagram showing shear wave front analysis,

[0108] · Figure 15 is a schematic diagram illustrating the extraction of a time-propagating displacement signal from a filtered velocity map. DETAILED DESCRIPTION

[0109] Different examples of the method and apparatus according to the present invention will now be described with reference to the accompanying drawings. In these different drawings, equivalent elements are represented by the same reference numerals.

[0110] 1. Method for determining 2D characteristics of propagative displacement waves

[0111] 1.1 SUMMARY

[0112] refer to Figure 1 , showing steps of a method configured to determine 2D characteristics of a propagating displacement wave propagating within a region of interest (ROI) of a patient.

[0113] As will be described in more detail below, the 2D characteristics of the propagating displacement wave may be:

[0114] The propagation velocity of the propagating displacement wave,

[0115] The direction of propagation of the propagating displacement wave, or

[0116] • More generally, any 2D representation (Cartesian, polar, ...) of the displacement wave velocity vector.

[0117] In the context of the present invention, a propagating displacement wave may be a shear wave generated by:

[0118] - naturally (ie, passive elastography), by propagating displacement waves generated by an organ such as the heart during the cardiac cycle, or

[0119] - Artificially (ie transient elastography), the propagating displacement waves are generated by an external source, for example an inertial vibration exciter of the type described in WO 2022 / 084502.

[0120] The method consists of the following stages:

[0121] a) an optional excitation phase 10, in which a propagating displacement wave is generated,

[0122] b) a tracking phase 20, in which the displacement of the medium within the ROI, caused in part by the propagating displacement wave, is tracked; this tracking phase enables the estimation of a plurality of velocity maps of the ROI over time,

[0123] c) a filtering phase 30 for separating the propagating displacement wave from the plurality of velocity maps,

[0124] d) a processing phase 40, during which the 2D characteristics of the propagating displacement wave are determined from the plurality of velocity maps.

[0125] In what follows, the analysis method will be described with reference to the processing of data acquired with an ultrasound probe, which is able to:

[0126] - generate at least one low-frequency elastic wave, called "propagating displacement wave", in the medium, and

[0127] - while generating the low-frequency elastic wave:

[0128] o emit high-frequency ultrasound waves, and

[0129] o receive acoustic echoes resulting from the reflection of the high-frequency ultrasound waves in the medium,

[0130] so as to observe the propagation of the low-frequency elastic wave in the medium.

[0131] However, the person skilled in the art will understand that this analysis method can be implemented in the context of passive elastography techniques in which the shear wave is naturally generated by the human body. In this case, the ultrasound probe comprises an array of transducer elements for emitting high-frequency ultrasound waves (1-20 MHz) and receiving acoustic echoes, so as to observe the propagation of the low-frequency elastic wave(s) naturally generated by the human body.

[0132] 1.2 Excitation phase (optional)

[0133] With reference to Figure 2 , a processing assembly for implementing the method according to Figure 1 is shown.

[0134] This processing assembly comprises:

[0135] • a signal acquisition probe S, and

[0136] • a driving and processing unit Uc for:

[0137] • controlling the probe S, and

[0138] • processing the signals acquired by the probe S.

[0139] The probe S comprises an array of transducer elements for emitting ultrasonic waves and receiving acoustic echoes, and an actuator for generating shear waves. The actuator can be an inertial vibration actuator of the type described in WO 2022 / 084502. The fixed portion of the actuator is mechanically attached to the transducer array to transmit vibrations to the probe, thereby mechanically generating propagating displacement waves.

[0140] The driving and processing unit Uc is connected to the probe S via wired or wireless communication means. The driving and processing unit Uc is configured to control the transducer elements of the probe S and process the data collected by the transducers of the probe S. The driving and processing unit Uc is further configured to activate the inertial vibration exciter to excite one (or more) propagating displacement waves in the medium. More specifically, the driving and processing unit Uc is configured to control the inertial vibration exciter to generate vibrations capable of generating propagating displacement waves in the medium, command the transducer elements to emit high-frequency ultrasonic waves in the medium, command the transducer unit to receive echoes reflected by the medium, convert the echoes into electrical reception signals, and process the reception signals.

[0141] The transmit and receive electronics are located between the drive and processing unit Uc and the probe S and provide the necessary signal generation, amplification, digitization, and transmission between the probe and Uc. This hardware can be physically located close to Uc and connected to the probe using a micro coaxial cable that transmits analog signals, or it can be located in the probe handle, in which case communication between the probe and Uc utilizes a purely digital link.

[0142] The driving and processing unit Uc may consist of one or more different physical entities that may be remote from the probe S. For example, the driving and processing unit Uc may be composed of one or more different physical entities that may be remote from the probe S. c This can include:

[0143] - one (or more) controllers 11, such as a smartphone and / or electronic tablet (e.g. ) and / or a personal digital assistant (PDA), or a chip on an electronic board (e.g., a microcontroller, an FPGA), or may be of any other type known to those skilled in the art, and

[0144] - One (or more) computers 12 (beamformer, processor, etc.), such as personal computers and / or workstations, etc.

[0145] One (or more) storage units 13 comprising at least one memory, such as random access memory (RAM) and / or read only memory (ROM) and / or a USB key, cloud storage, etc. The storage unit may be part of the controller 11 or the computer 12 .

[0146] In addition to storing acquired / processed data, the memory unit 13 can store programming code instructions intended to execute Figure 1 Stages of the analytical method are shown.

[0147] Figure 2 The working principle of the processing components shown is as follows.

[0148] To perform the excitation phase 10, the user contacts the probe S with the patient's skin. For example, if the user wishes to analyze the characteristics of the patient's liver, the transducer element array contacts the patient's skin in the right rib area so that the transducer element array extends substantially between two adjacent ribs.

[0149] Drive and processing unit U C The actuator is driven to generate vibrations, which excite propagating displacement waves in the patient's body.

[0150] 1.3 Tracking phase

[0151] During the tracking phase 20, the displacement of the medium caused by the propagation of the propagating displacement wave is measured within the ROI by pulse echo ultrasound (e.g., tissue Doppler). Specifically, multiple velocity maps of the ROI are estimated over time during the tracking phase. These velocity maps correspond to "displacement images" and include information representing the propagation of the propagating displacement wave in the ROI.

[0152] To implement the tracking phase, the driving and processing unit Uc controls the array of transducer elements of the probe S to emit a series of ultrasonic waves (for example, each having a central frequency between 1 MHz and 20 MHz).

[0153] These ultrasonic waves penetrate the medium and are reflected by scattering particles (eg, collagen particles) contained in the medium, thereby being able to track the displacement of these scattering particles in the medium.

[0154] Specifically, refer to Figure 3 , the tracking phase 20 includes the following sub-steps:

[0155] - the driving and processing unit Uc controls the transducer element array to transmit a series of ultrasonic waves in the medium at a rate of 10 to 10,000 transmissions per second (step 201),

[0156] - the driving and processing unit Uc controls the reception of the array of transducer elements and records (in real time) in the storage unit 13 the time-dependent acoustic signals received by the transducer array, said time-dependent acoustic signals representing the echoes generated by the ultrasound waves interacting with the scattering particles of the ROI (step 202),

[0157] -Drive and processing unit Uc:

[0158] o constructing a plurality of complex echogenicity maps of the ROI using the time-dependent acoustic signal (step 203), and

[0159] o Estimate a plurality of velocity maps of the ROI over time based on the plurality of complex echogenicity maps (step 204).

[0160] More precisely, the sub-step of constructing a plurality of echogenicity maps comprises determining the images of the ROIs included in the tracking field (i.e., the area penetrated by the ultrasound waves) by solving an inverse problem that can be simplified to a conventional "beamforming" problem. Each echogenicity map is generated by a composite processing of one or more consecutive transmissions.

[0161] In the subsequent sub-step of estimating multiple velocity maps, these consecutive ROI images are processed by filtering and correlation processing (in particular by cross-correlation processing):

[0162] - Delay between two frames (time domain) to perform inter-frame correlation processing,

[0163] - Or perform inter-frame correlation processing with a reference image.

[0164] 1.4 Filtering phase

[0165] As described above, the velocity map represents the propagation of the propagating displacement wave in the ROI. However, the velocity may also include components that may seriously interfere with the propagating displacement wave.

[0166] The filtering stage 30 separates the portion of the displacement due to the propagation of the propagating displacement wave from other artifact components (eg static deformations, compression waves, out-of-plane shear waves or reflection waves).

[0167] The filtering stage can be based on different techniques known in the art. For example, it can be implemented using the solution disclosed in document FR2210011 entitled “PROCEDE D'ANALYSE D'UN MILIEU PERMETTANT DE REDUIRE LESEFFETSD'ARTEFACTSDUS A DES DEFORMATION STATIQUES DANS LE MILIEU”.

[0168] The filtering stage can also be based on temporal filtering (widely used for harmonic elastography, as it enables band-pass filtering of the velocity map in the frequency band of interest), or on spatial filtering, or on spatio-temporal joint filtering. An example of such a technique is disclosed in the article by Manduca et al. (A. Manduca, D. S. Lake, S. A. Kruse, R. L. Ehman, “Spatio-temporal directional filtering for improved inversion of MR elastography images”, Med. Image Anal., Vol. 7, No. 4, pp. 465-473, December 2003).

[0169] In any case, the filtering stage enables obtaining a plurality of filtered velocity maps, from which some unwanted components have been removed.

[0170] 1.5 Processing phase

[0171] Once filtered, the plurality of filtered velocity maps is processed to determine, for each point of interest within the velocity map of the plurality of filtered ROIs, the 2D characteristics of the propagating displacement wave. The 2D characteristics of the propagating displacement wave can include:

[0172] - the propagation velocity of the propagating displacement wave, and

[0173] - the propagation direction of the propagating displacement wave.

[0174] These 2D characteristics of the propagating displacement wave are estimated for each point of interest within the plurality of velocity maps. This enables obtaining a map of the 2D characteristics of the propagating displacement wave.

[0175] In the following, the processing stage will be described for a single point of interest A of a ROI, it being understood that the different steps disclosed below can be repeated for any point of interest (A', A", etc.) of the ROI.

[0176] 1.5.1. Correlation of neighboring points

[0177] The processing stage comprises associating two (or more than two) neighboring points B, C (step 402) with the point of interest A.

[0178] The two neighboring points B, C are chosen such that the three points A, B, C are not collinear.

[0179] The skilled person will understand that the point A, the point B, the point C have the same coordinates within all the filtered velocity maps, as the plurality of filtered velocity maps show the same ROI.

[0180] 1.5.2. Extraction of time domain displacement signals

[0181] The processing stage includes extracting the time domain propagation displacement signal S of point A, point B and point C from the filtered multiple velocity images / maps A 、S B 、S C step (step 403).

[0182] It should be noted that each filtered velocity map:

[0183] - represents the velocity of the propagating displacement wave at each point of the ROI, and

[0184] - corresponds to each time t during the propagation of said propagating displacement wave.

[0185] like Figure 15 As shown at the midpoint A, the time domain propagation displacement signal S is determined A Including: the velocity V at point A estimated at different times t1, t2, t3, t4 in the velocity diagrams VM1, VM2, VM3, VM4 A To express.

[0186] By extracting the speed V included in the plurality of speed maps VM1, VM2, VM3, VM4 B 、V C , the time domain propagation displacement signal S is obtained in the same way B 、S C , each speed diagram corresponds to each time t1, t2, t3, t4.

[0187] 1.5.3. Estimation of apparent time of flight

[0188] After obtaining the time domain propagation displacement signal S A 、S B 、S C Afterwards, the processing stage includes the time domain propagation displacement information S A 、S B 、S C The step of estimating the apparent flight time (step 404) determines the propagation time of the propagating displacement wave between point A and point B on the one hand, and between point A and point C on the other hand.

[0189] Each time domain propagation displacement signal S A 、S B 、S C Represents the displacement of the propagating displacement wave at each point A, B, and C in the ROI.

[0190] In other words, the time domain propagation displacement signal S A 、SB 、S C It represents the change in the relative position of scattering particles in the medium caused by the propagation of the propagating displacement wave in the medium.

[0191] Since points A, B, and C are separated from each other in space, and since it is assumed that the waveform of the propagating displacement wave remains unchanged during its propagation in the medium, the time-domain propagating displacement signal S A 、S B 、S C The waveforms are the same, but there is a time shift (for example, the propagating displacement wave reaches point A before reaching point B or point C).

[0192] In order to determine the apparent flight time of the propagating wave between point A and point B (longitudinal apparent flight time), the time domain propagating displacement signal S A and S B Perform cross-correlation processing to determine the value of the A ) relative to another (e.g., S B ) to measure the signal S A and S B This allows us to evaluate the similarity of the signal S A With S B The time delay between points A and B corresponds to the apparent longitudinal flight time of the propagating displacement wave between points A and B.

[0193] In order to determine the apparent flight time (transverse apparent flight time) of the propagating wave between point A and point C, the time domain propagating displacement signal S A and S C Perform cross-correlation processing to determine the value of the A ) relative to another (e.g., S C ) to measure the time shift of the signal S A and S C This allows us to evaluate the similarity of the signal S A With S C The time delay between points A and C corresponds to the apparent lateral flight time of the propagating displacement wave between points A and C.

[0194] 1.5.4. Determination of 2D characteristics of propagative displacement waves

[0195] As described above, select adjacent points B and C so that:

[0196] - Assume that the line segment [A, B] is approximately parallel to the direction of propagation of the propagating displacement wave,

[0197] - Assume that the line segment [A, C] is approximately perpendicular to the direction of propagation of the propagating displacement wave,

[0198] - Assume that the angle θ between B, A and C (or ) is approximately equal to π / 2 (i.e., line segment [A, B] is perpendicular to line segment [B, C]).

[0199] Based on the longitudinal and transverse apparent flight times, determination of the 2D characteristics of the propagating displacement wave (eg, propagation speed and propagation direction) may be achieved (step 405 ).

[0200] 1.5.4.1. Calculation of longitudinal and transversal apparent velocities

[0201] The longitudinal apparent shear wave front velocity V1 between points A and B can be determined by calculating the ratio between the following two terms:

[0202] -The distance separating point A and point B (hereinafter referred to as the "first distance") divided by

[0203] - The apparent longitudinal flight time between points A and B.

[0204] In fact, since the positions of point A and point B are known (from the association of adjacent points B in the association step), the first distance is known and the longitudinal apparent flight time is calculated in step 1.5.3.

[0205] In a similar manner, the transverse apparent shear wave front velocity V2 between points A and C can be derived by calculating the ratio between the following two terms:

[0206] - the second distance separating point A from point C (known a priori from the association step where point A is associated with point C) divided by

[0207] - The apparent lateral flight time between points A and C.

[0208] 1.5.4.2. Calculation of real velocity and direction of propagative displacement waves

[0209] Using the calculated longitudinal apparent shear wave front velocity V1 and transverse apparent shear wave front velocity V2, the propagating displacement wave can be derived:

[0210] - actual shear wave front velocity V, and

[0211] - Actual direction of propagation.

[0212] See Section 3.1.8 for details.

[0213] 1.5.4.3. SUMMARY

[0214] By selecting neighboring points B and C of the point of interest A (said points B and C being selected so that (AB) and (AC) are approximately perpendicular), by estimating the longitudinal and transversal apparent flight times of the propagating wave between points A, B and C (using the time-domain propagating displacement signal S A B C ), by deriving the longitudinal apparent velocity VI and the transversal apparent velocity V2 of the propagating wave between points A, B and C, and by inferring the actual velocity and direction of the propagating displacement wave using the longitudinal apparent velocity VI and the transversal apparent velocity V2, the above method is able to determine:

[0215] - the actual propagation velocity of the propagating displacement wave, and

[0216] - the actual propagation direction of the propagating displacement wave.

[0217] The processing stages disclosed above are implemented for a plurality of points of interest (A, A', A", etc.). This allows obtaining a 2D profile of the propagating displacement wave within the ROI.

[0218] The 2D profile of the propagating displacement wave can then be used for a variety of applications.

[0219] 2. Applications

[0220] 2.1. Hardness estimation

[0221] For example, the 2D profile of the propagating displacement wave can be used to determine an estimate of the stiffness of the medium within the ROI. Indeed, the propagation rate (velocity) of the propagating displacement wave varies according to the stiffness of the target region. The velocity of the propagating displacement wave in a fibrotic liver (typically up to 5 m / s) is significantly higher than in a healthy liver (typically 1 m / s).

[0222] To determine the stiffness of the tissue located within the region of interest, the method can comprise a step of estimating the stiffness of the target tissue using the 2D profile (more specifically, the velocity information comprised in said profile).

[0223] The estimate of the stiffness of a given subset of points (corresponding to the target tissue) within the region of interest can comprise the following sub-steps:

[0224] - computing a median or mean of the propagation velocities of said subset of points, and

[0225] - estimating the stiffness from said median or mean of the propagation velocities, for example using a simple relationship between stiffness and velocity (for example, the formula E = 3 x c 2 , E being the stiffness and c being the velocity) applicable to quasi-incompressible media such as biological tissues.

[0226] ​​Alternatively, estimating the hardness of a given subset of points within the region of interest (corresponding to the target tissue) may include the following sub-steps:

[0227] - Calculating a stiffness map from the 2D characteristic map (more specifically the velocity information included in said map), for example using a simple relationship between stiffness and velocity (e.g. the formula E=3×c for quasi-incompressible media such as biological tissue) 2 , E is hardness, c is speed), then

[0228] - Calculate the median or mean hardness value of a subset of points based on the hardness map.

[0229] 2.2. Determination of quality factors related to orientation and position of the probe

[0230] The 2D characteristic diagram of the propagating displacement wave can also be used to prompt the user whether the direction or position of the probe on the patient's body surface is correct.

[0231] To this end, the method may include the step of determining a quality factor representing the orientation or position of the probe on the patient's body surface. The step of determining the quality factor may include the following sub-steps: determining an average shear wave propagation direction (or angle) for a subset of points of the 2D characteristic map, comparing the average shear wave propagation angle with a preferred angle range (which is predefined), and providing a quality factor for the probe orientation and probe position based on the result of the comparison:

[0232] - If the average shear wave propagation angle lies within a predefined range, the probe's orientation or position is considered correct and the quality factor indicates to the user that the probe's orientation or position is correct,

[0233] - If the average shear wave propagation angle is outside a predefined range, the probe is considered to be oriented or positioned incorrectly, and the quality factor indicates to the user that the probe is oriented or positioned incorrectly; the user can change the probe's orientation or position.

[0234] 2.3. Determination of a coefficient representative of the pressure exerted by the probe on the patient

[0235] Furthermore, in the context of transient elastography, a 2D characteristic map of the propagating displacement wave can be utilized to provide the user with information related to the pressure exerted by the probe on the patient's skin.

[0236] To this end, the method may comprise the step of determining a coefficient representative of the pressure exerted by the probe on the patient, said determination step comprising the following sub-steps:

[0237] - Calculate the curvature of the propagating displacement wave front based on the 2D characteristic map,

[0238] - Compare the calculated curvature to the baseline information to determine whether the curvature is concave or convex,

[0239] - If the calculated curvature is concave, assign the coefficient a value representing light pressure,

[0240] - If the calculated curvature is convex, assign the coefficient a value representing high pressure.

[0241] 2.4. Determination of parameters related to the quality of contact between the probe and the patient

[0242] In addition, the 2D characteristic diagram of the propagating displacement wave can be used to prompt the user whether the contact quality between the probe and the patient's skin is good.

[0243] To this end, the method may include the step of determining a parameter representing the contact quality between the probe and the patient's skin, the determining step including the following sub-steps: estimating a variance in the direction of the propagating displacement wave based on the 2D characteristic map, and determining a value for the parameter based on the determined variance. An exemplary embodiment may be comparing the determined variance with a predefined threshold. If the variance exceeds the threshold, the contact quality between the probe and the patient is considered low.

[0244] In another embodiment, we can calculate the amplitude of the variation in the transverse time of flight (corresponding to the time of flight between points A and C), which is a quantity that represents the amplitude of the spatial variation in the direction of the propagating displacement wave.

[0245] 3. Theory related to the invention

[0246] 3.1 Detailed theoretical description of the invention

[0247] 3.1.1 General problem and proposed solution

[0248] In ultrasound elastography, tissue displacement can be induced by:

[0249] ○Mechanical vibrator (Fibroscan VCTE technology),

[0250] ○Ultrasound radiation force (ARFI, ultrasonic shear wave elastography),

[0251] o Natural vibrations of the body (e.g., heartbeat).

[0252] In both cases, the main assumption is that the displacement direction is known and the analysis of the induced shear wave propagation is performed in a direction chosen a priori to be perpendicular to the shear wave front.

[0253] There are cases where this prior assumption is wrong:

[0254] o In transient elastography (i.e. when the shear wave propagation is induced by an external vibration), if there is a mechanical coupling between the elastography probe (shear wave source) and the rib, the shear wave can originate from the rib edge. In this case, the propagation direction can be at an angle with respect to the probe symmetry axis, leading to an overestimation of the liver stiffness,

[0255] o In ARFI or ultrasound shear wave imaging, if the medium is not homogeneous, shear waves can be refracted and reflected. Thus, if they are analyzed along their propagation direction, it will again lead to an overestimation of the stiffness,

[0256] o In passive elastography (where shear waves are generated with the natural vibrations of the body), the propagation direction of the shear waves is not known beforehand, and again analyzing the shear waves with a pre-selected propagation direction will lead to biased results.

[0257] We propose to solve the above problems by performing a 2D analysis of the propagation of the generated shear waves, enabling us to:

[0258] o estimate the shear wave propagation direction in the 2D imaging plane,

[0259] o use the real shear wave propagation direction to estimate the correct, not overestimated stiffness value,

[0260] o in the case of liver elastography, the shear wave propagation direction measured in the 2D imaging plane helps the operator to position and orient the probe and to apply the right amount of pressure to the patient.

[0261] 3.1.2 Doppler processing of velocity map estimation

[0262] Doppler processing aims at extracting the displacement of structures such as blood and tissue from continuously transmitted ultrasound (US) beams. In US imaging, displacement extraction is performed by analyzing the variations in a number of echo intensity maps (which can be in-phase quadrature (IQ) or radio-frequency (RF) beamformed maps (i.e. echo intensity maps) along the slow time scale). Therefore, all techniques start first with the transmission / reception and beamforming of a set of IQ / RF images to obtain echo intensity maps represented as slow time samples. Such images are usually acquired at a constant rate represented as the pulse repetition frequency (PRF).

[0263] To extract the Doppler frequency shifts from the RF or IQ image slow time samples, the academic community has developed a number of techniques. However, different techniques are based on the same common principle: tracking the displacement of the recorded backscattered signals between consecutive pulses.

[0264] The difference between the different methods lies in the way this motion is tracked, which mainly depends on the nature of the received backscattered echoes (RF or IQ, narrowband or broadband).

[0265] 3.1.2.1 Phase-based estimator

[0266] Well-known techniques include calculating the phase shift using a lagged 1-order temporal autocorrelation (often denoted as Kasai autocorrelation (or one-dimensional (1D) autocorrelation)), as described by Angelsen (BA Angelsen, "Instantaneous frequency, mean frequency, and variance of mean frequency estimators for ultrasonic blood velocity Doppler signals," IEEE Trans. Biomed. Eng., Vol. 28, No. 11, November 1981, pp. 733-741) and Kasai and Nagekawa (C. Kasai and K. Namekawa, "Real-time two-dimensional blood flow imaging using an autocorrelation technique," IEEE 1985 Ultrasonics Symposium, San Francisco, CA, USA, 1985. doi: 10.1109 / ultsym.1985.198654). The Kasai 1D autocorrelation method is briefly reviewed below.

[0267] Formally, considering that we can use N EL (set length) consecutive slow time samples, each slow time sample is composed of Nr×N θ The beamformed image (i.e., echo intensity map) is formed on a grid of pixels. Slow time samples are collected at time intervals TP (pulse repetition period), T P =1 / f P , where f P is the pulse repetition frequency (PRF). Doppler phase shift It can be estimated by the following equation:

[0268]

[0269] Among them, γ k is the slow time sample, w k is the window coefficient.

[0270] Note that the set length can be reduced to 2 consecutive frames. The estimated phase shift and Doppler shift are related by the following equations:

[0271] 2πδf p ×T P =ΔΦ

[0272] Furthermore, we can relate velocity to Doppler shift as

[0273]

[0274] where c corresponds to the average speed of sound and f0 is the demodulation frequency.

[0275] Therefore, we can intuitively see that there is a maximum speed associated with a given PRF, corresponding to the maximum phase shift between consecutive signals, which can be derived from the above equation as:

[0276]

[0277] The main difference between the Loupas estimator and the Kasai 1D correlator is that the mean frequency used to extract the velocity from the Doppler shift is calculated locally using fast time samples, thus accounting for local variations in frequency caused by random fluctuations in the speckle.

[0278] One can generalize the Kasai estimator defined above to arbitrary lags, which may be useful when the slow-time samples include IQ data with beamforming at different angles. In this case, the Kasai estimator can be directly expressed as follows:

[0279]

[0280] Here, m represents the lag order of the correlation.

[0281] The parameters for the velocity extraction step are as follows:

[0282] ο set length,

[0283] o range gate length (used by the Loupas estimator to calculate the average frequency),

[0284] o transverse average length (used to filter the correlation along the azimuth dimension for noise reduction purposes),

[0285] ο time window coefficient,

[0286] ο radial window coefficient,

[0287] ο horizontal window coefficient,

[0288] ο boundary conditions in different dimensions,

[0289] οConvolution modes of different dimensions,

[0290] ο Displacement extraction method: Kasai or Loupas,

[0291] o Additional lag order (in case of lagged m-order correlation).

[0292] 3.1.2.2. Time-based estimator

[0293] Time-based estimators (C. Kasai and K. Namekawa, “Real-time two-dimensional blood flow imaging using an autocorrelation technique,” ​​IEEE 1985 Ultrasonics Symposium, San Francisco, CA, USA, 1985. doi: 10.1109 / ultsym.1985.198654) represent an alternative to phase-based estimators that are suitable for broadband US pulses and large displacements where aliasing may occur with phase-based estimators.

[0294] Consider two slow time samples γ for the same radial line (or equivalently axial line) k (r) and γ k+1 (r), where the azimuthal (or equivalently transverse) coordinates are neglected for clarity. Assume that the displacement of the medium is radial and the velocity v r Constant. It can be expressed as:

[0295] γ k+1 (r) = γ k (r+v r T P ),

[0296] Among them, T P Indicates the pulse repetition period.

[0297] Therefore, the slow-time samples are only shifted in time with respect to each other, and the time shift depends on the radial velocity. To estimate the time shift, a windowed cross-correlation is used such that:

[0298]

[0299] When Δr=v r T P· Therefore, maximizing the cross-correlation gives an estimate of the radial velocity.

[0300] 3.1.3. Filtering techniques of velocity maps

[0301] This section aims to provide scientific details related to the filtering stage of the analysis method according to the present invention.

[0302] Filtering techniques aim to separate the portion of the measured displacement that comes from shear wave propagation from other components (e.g., static deformations, compressional waves, out-of-plane shear waves, or reflected waves). They can rely on signal processing techniques (e.g., filtering) and clever acquisition or shear wave generation procedures. A few examples are given below.

[0303] In transient elastography, when shear wave generation relies on the mechanical vibration of the probe, the propagating displacement wave is contaminated by static deformations, compression waves, and potential interferences.

[0304] Different methods are used to reduce static deformations (see document entitled “PROCEDE d'ANALYSE d'UN MILIEUPERMETTANT DE REDUIRE LES EFFETS d'ARTEFACTS DUS A DES DEFORMATION STATIQUESDANS LE MILIEU”, application number FR2210011).

[0305] Symmetrical sampling (DC Melema et al., “Probe Oscillation Shear Elastography (PROSE)”, IEEE Trans. Med. Imaging, Vol. 35, No. 9, September 2016, pp. 2098-2106) is based on acquiring the displacement field at exactly the same phase instant as the static deformation.

[0306] Impulse vibration generation (S. Catheline, F. Wu and M. Fink, “A solution to diffraction biases in sonoelasticity: the acoustic impulse technique”, J. Acoustic Soc. Am., Vol. 105, No. 5, May 1999, pp. 2941-2950) is based on generating shorter mechanical vibrations to separate static deformation from shear wave propagation.

[0307] Low-pass IIR filtering (D.C. Melema et al., “Probe Oscillation Shear Elastography (PROSE): A High Frame-Rate Method for Two-Dimensional Ultrasound Shear Wave Elastography”, IEEE Trans. Med. Imaging, Vol. 35, No. 9, pp. 2098-2106, September 2016) and empirical mode decomposition (D.C. Melema et al., “Probe Oscillation Shear Wave Elastography: Initial In Vivo Results in Liver”, IEEE Trans. Med. Imaging, Vol. 37, No. 5, pp. 1214-1223, May 2018) are filtering techniques that act directly on the displacement field.

[0308] To further filter propagating displacement waves from more common artifacts (e.g. compression waves, interference, in-plane shear wave reflections), time, spatial and spatio-temporal (i.e. directional) filtering techniques were developed and are commonly used in almost all types of dynamic elastography: magnetic resonance elastography (A. Manduca, D. S. Lake, S. A. Kruse, and R. L. Ehman, “Spatio-temporal directional filtering for improved inversion of MR elastography images,” Med. Image Anal., vol. 7, no. 4, pp. 465-473, Dec. 2003), harmonic elastography (H. Zhao et al., “External vibration multi-directional ultrasound shear wave elastography (EVMUSE): application in liver fibrosis staging,” IEEE Trans. Med. Imaging, vol. 33, no. 11, pp. 2140-2148, Nov. 2014), shear wave elastography (T. Defieux, J.-L. Gennisson, J. Bercoff, and M. Tanter, “On the effects of reflected waves in transient shear wave elastography,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 58, no. 10, pp. 2032-2035, Oct. 2011). These methods are based on filtering the displacement field directly.

[0309] Time domain filtering is widely used in harmonic elastography as it enables band-pass filtering of the displacement field around the frequency of interest. Spatial and spatio-temporal techniques are used in almost all elastography methods.

[0310] 3.1.4 Time shift estimation of shear wave propagation based on tracking of velocity maps

[0311] Once the displacement (i.e. multiple velocity maps) is computed and filtered, the shear wave velocity can be computed and thus the stiffness can be inferred. The computation of the shear wave velocity (i.e. group velocity or shear wave front velocity) can be done by considering two points in the medium with a known distance between them and performing a time shift estimation of the velocity signals of these two points along the slow time scale.

[0312] Estimating time shift from displacement is based on well-known techniques in signal processing. The principle is to compare slow-time velocity signals corresponding to two points oriented parallel to the shear wave propagation direction. These velocity signals are therefore phase-shifted or time-shifted replicas of each other. For example, the maximum value of the generalized cross-correlation function can be used to recover the time / phase shift. Given the distance and time / phase shift, the shear wave velocity can be easily inferred.

[0313] 3.1.5. Instantaneous elastography for the generation of mechanical shear waves

[0314] 3.1.5.1. Theory and far field beam pattern related to point source excitation

[0315] Figure 5 The radiation amplitude 1 of a shear wave generated by a point source 2 on a patient's skin 3 is shown. For example, the point source is a vibrating needle, which can excite shear waves in the patient's body when the vibrating needle acts on the patient's skin.

[0316] The shear wave front is approximately hemispherical centered at the excitation point (corresponding to the contact point between the point source and the patient's skin 1).

[0317] As can be seen by those skilled in the art, the direction of the maximum amplitude of the shear wave is at an oblique angle relative to the normal line of the skin surface. Unfortunately, the amplitude is close to zero in the normal direction.

[0318] 3.1.5.2 Green function for point source excitation

[0319] In an infinite medium, the Green's function associated with a point source excitation can be decomposed into four terms (see S. Catheline and N. Benech, "Longitudinal shear wave and transverse dilatational wave in solids," J. Acoust. Soc. Am., Vol. 137, No. 2, February 2015, pp. EL200-5):

[0320]

[0321] Where m and n represent the direction of the point source force and the direction of observation, and r is the distance from the point source.

[0322] and are Green's function terms associated with the well-known far-field compressional and shear waves. These waves are longitudinally and transversely polarized, respectively.

[0323] The other two Green's function terms and denote the near-field shear wave term and the near-field compression wave term with longitudinal and transverse polarization, respectively. Note that the polarization of the near-field term is opposite to that of the far-field term. These terms are near-field terms, which are expressed in terms of r -2 Decrease.

[0324] The dominant terms are far-field compression (P) and shear (S) waves. These terms are observed to vanish in the directions perpendicular and parallel to the point source force, respectively. A weak near-field term can be observed in the region where the far-field terms vanish.

[0325] In soft media, compression waves can be neglected, simplifying the Green's function to the sum of two shear wave terms.

[0326] Now, assuming that the point source force is applied in the vertical direction (n=3), and we observe the field in the direction parallel to the point source force (m=3), we have

[0327] and

[0328] Where ρ is the density, k represents the shear wave number, and β is the shear wave velocity.

[0329] Therefore, in the direction parallel to the point source force direction, the only remaining propagation term is the near-field term, which has the same characteristics as the far-field term in terms of phase velocity and frequency, but with a longitudinal polarization.

[0330] In a semi-infinite medium (i.e., when a point source is located at the surface of the medium), the expression for the Green's function in the medium is similar to that for an infinite medium, except that a coefficient of 2 is introduced because the interface reflects the originally backward-propagating wave like a mirror (see L. Sandrin, D. Cassereau, and M. Fink, "The role of the coupling termin transient elastography," J. Acoust. Soc. Am., Vol. 115, No. 1, January 2004, pp. 73-83).

[0331] Therefore, in soft media, the only remaining term is the very weak shear wave near-field term in the direction parallel to the direction of the point source force.

[0332] 3.1.6. Example of a clinical system for liver stiffness measurement relying on instantaneous elastography

[0333] 3.1.6.1.

[0334] The system relies on the damped vibrations of a single-element ultrasound probe. Depending on the patient's size (children, adults with different BMIs), different probes with different diameters (7mm, 9mm, 12mm) can be used.

[0335] refer to Figure 6 , the probe P1 includes:

[0336] ○ Ultrasonic transducer 14, which is used to transmit ultrasonic waves and collect echoes, and

[0337] o A vibrator 15 forming a non-point source of shear waves, said vibrator 15 comprising a cylindrical rod, the free end of which (a disc of a few millimeters in diameter) is adapted to come into contact with the patient's skin.

[0338] The transducer 14 is attached to the end of a vibrator 15. The vibrator 15 vibrates the transducer 14 to generate shear waves. The operating principle of a medical pulse elastography device is as follows. The vibrator 15 is activated, causing the transducer 14 to move, generating low-frequency shear waves in the tissue being analyzed. During the propagation of the low-frequency shear waves, the transducer 14 transmits and receives high-frequency ultrasound waves, enabling the study of the propagation of the low-frequency shear waves.

[0339] Figure 7 1 shows the shear wave amplitude of the shear wave generated by the vibrator 15. Figure 7 As shown in Figure 1, the location of the shear wave source corresponds to the periphery of the free end of the rod. This is the convolution of the point source Green's function described above with the shape of the ultrasonic transducer tip.

[0340] The shear wave on the axis is no longer near zero due to constructive interference of shear waves from the rod's periphery (which have non-zero angles). This phenomenon is due to the non-zero diameter of the rod; otherwise, only a very weak near-field term would be seen. The larger the rod diameter, the larger the angle at a given depth.

[0341] If the angle is θ, then the apparent velocity of the shear wave measured on the axis will be greater than the true shear wave velocity by a factor For example, if the rod diameter is 10 mm and the analysis depth is 40 mm, the shear wave velocity will be overestimated by 0.8% and the hardness will be overestimated by 1.6%. Such fixed overestimations related to geometry are small and can be corrected.

[0342] 3.1.6.2. Hepatoscope

[0343] Other types of ultrasonic pulse elastography probes have been proposed in the art. For example, document WO 2022 / 084502 describes a probe called "Hepatoscope" that can measure the elasticity of a medium. Figure 8 As shown, the probe P2 includes an inertial vibration exciter 31, a transducer array 32 and other components (housing, electronic card 33, etc.).

[0344] Apart from the moving parts of the inertial vibration exciter 31 (not shown), all components are mechanically integrated so that when the inertial vibration exciter 31 is activated, the entire probe vibrates.

[0345] This ultrasonic pulse elastography probe has a more complex shear wave generation mode. Figure 5 It is obtained by convolving the Dirac point source excitation pattern (Dirac Stylus diagram) with the shape of the vibrating probe.

[0346] The intensity of the shear waves generated from each probe surface location will increase as the local curvature of the interface between the probe and the patient's skin increases.

[0347] Figure 9 Shown is a bottom view of the probe P2 and its lens and the shear wave generation area:

[0348] o The first region 21 has a medium curvature, corresponding to the lens in contact with the patient's skin; the first region 21 has a medium shear wave generating capability,

[0349] ○ The second region 22 is the longitudinal edge of the contact surface between the probe and the patient's skin; the second region has a high shear wave generating capability,

[0350] o A third region 23 corresponds to the lateral edges of the contact surface between the probe and the patient's skin; this third region has normal shear wave generating capabilities.

[0351] 3.1.7. Bad effects related to instantaneous elastography for liver stiffness measurement

[0352] As mentioned above, a probe of the type previously described ( A hepatoscope must be placed in contact with the patient's skin to generate shear waves. When such a probe is used to measure the viscoelasticity of an organ (e.g., a patient's liver), the probe is in contact with the skin area above the patient's ribs.

[0353] In this case, different adverse effects may occur:

[0354] ○ On the one hand, there is the shear wave shielding effect, and

[0355] ○On the other hand, there is the mechanical coupling effect.

[0356] 3.1.7.1. Shear wave shadowing effect

[0357] When the user applies "light" pressure (e.g., less than 5 Newtons) to bring the UPE probe into contact with the patient's skin, a shear wave shielding effect is triggered. This shear wave shielding effect is caused by the shear waves generated by the UPE probe being reflected off the ribs.

[0358] Figure 10a and Figure 10b Shown in a plane perpendicular to the ribs, when using ( Figure 6 ) and when using a Hepatoscope( Figure 8 ) occurs when the shear wave shielding effect occurs.

[0359] like Figure 10a As shown, in In this case, only a portion of the free end of the rod is able to generate shear waves that propagate in the medium between two consecutive ribs. In fact, the shear waves generated by the portion of the free end of the rod located above the ribs are shielded and reflected onto the ribs.

[0360] like Figure 10b As shown, in the case of the hepatoscope, the only shear waves that are obscured and reflected onto the ribs are the shear waves generated by the second region (corresponding to the longitudinal edge of the contact surface between the probe and the patient's skin). The shear waves generated by the first and third regions propagate between the ribs and within the medium.

[0361] Figure 11a and Figure 11b Shown in a plane parallel to the ribs, when using ( Figure 6 ) and when using Hepatoscope( Figure 8 ) propagation of shear waves in a medium. A person skilled in the art will understand that:

[0362] ○In In the case of Figure 11a This is why Fibroscan offers different probe sizes, in order to minimize shadowing effects.

[0363] o In the case of the Hepatoscope, the shear waves generated by the first and third regions travel through the medium beneath the ribs; the shear wave front generated by their interference has a concave shape (see Figure 11b ).

[0364] 3.1.7.2. Mechanical coupling effect

[0365] When the user applies a "high" pressure (e.g., a force greater than or equal to 5 Newtons) to press the UPE probe against the rib, a mechanical coupling effect is induced. In this case, the displacements of the probe and rib are identical, the probe and rib move as a single mechanical assembly, and the shear generation region is limited to the rib edge and the lateral edge of the probe between the ribs.

[0366] Figure 12a 、 Figure 13a and Figure 12b 、 Figure 13b Shows when using ( Figure 12a 、 Figure 13a ) and when using a Hepatoscope( Figure 12b 、 Figure 13b ) occurs when the mechanical coupling effect Figure 12a and Figure 12b In the plane perpendicular to the ribs, Figure 13a and Figure 13b In a plane parallel to the ribs. In both cases, the main shear waves originate from the outer edges of the ribs, the region where the maximum possible shear forces are generated between the ribs and the tissue during displacement. In the region where the probe moves with the tissue, no shear waves are generated. The probe generates weak shear waves only in the region where it is displaced relative to the tissue: the probe edge between the ribs (the third region in the case of a hepatoscope).

[0367] The shear waves excited from the ribs are efficient, but the excited shear waves excited by this vibration have such a propagation direction:

[0368] ○For outside the axis, and

[0369] ○ For Hepatoscope, outside the imaging plane.

[0370] Since the geometry (rib dimensions and positions) is unknown, the projection angle onto the axis or imaging plane is incorrect: the shear wave velocity estimate will always overestimate the unknown amount.

[0371] In the case of a hepatoscope, the angle of the shear wave front projection on the imaging plane can be determined by processing multiple velocity maps. When the probe is pressed symmetrically against two ribs, the excited wave front is convex, which makes it possible to distinguish between high-pressure situations (convex wave front) and low-pressure situations (concave wave front).

[0372] However, the user may press the probe against only one rib, and the analysis shows that the induced shear wave front has a large angle relative to the normal (greater than 10 degrees or less than -10 degrees).

[0373] 3.1.7.3. CONCLUSION

[0374] When mechanical coupling effects occur, no quantitative value for the shear wave velocity in the liver can be derived. It is therefore preferred to avoid (or at least limit) such mechanical coupling effects.

[0375] exist In the case of , due to its low shear wave excitation efficiency, it is necessary to apply a relatively high pressure between the probe and the patient's skin. In order to avoid mechanical coupling effects, it is necessary to use a rod with a free end diameter smaller than the rib spacing. This is The reason for providing the user with a set of three probes, each comprising a free end of a rod of a different diameter, is that this allows for optimization of the shear wave intensity on the one hand and the shear wave penetration depth on the other hand.

[0376] For hepatoscopes, where shear wave excitation is more efficient (due to a larger mechanical aperture), a “light” pressure can be applied to bring the ultrasound pulse elastography probe into contact with the patient’s skin. This limits mechanical coupling effects while ensuring that shear waves with sufficient amplitude are excited to propagate in the medium.

[0377] Analyzing shear wave propagation angles in the 2D imaging plane enables:

[0378] ○ Correction of bias in hardness estimation caused by unknown shear wave propagation direction in the imaging plane,

[0379] ○ Distinguish between geometric features specific to light pressure and those specific to rib coupling; this will allow the rejection of measurements of shear waves excited by out-of-plane ribs and provide a visual aid to the user to avoid this.

[0380] 3.1.8. 2D shear wave front angle analysis performed from velocity signals of points A, B, C

[0381] Figure 14 The shear wave front is shown, with a propagation velocity of V, and three analysis points in the image: point A, point B, and point C.

[0382] The principle of the present invention is to calculate:

[0383] o The time of flight of the shear wave between points A and B, which are on a line chosen to be approximately parallel to the direction of propagation; knowing this time of flight enables calculation of the projected velocity V1 of the shear wave on this line,

[0384] o The time of flight of the shear wave between points A and C, which are on a line chosen to be approximately perpendicular to the direction of propagation; knowing this time of flight enables calculation of the projected velocity V2 of the shear wave on this line.

[0385] Then, from the projected velocities V1 and V2, the corrected actual propagation velocity V of the shear wave at point A can be estimated based on the following derivation.

[0386] Since: θ=θ1+θ2, and V=V1cosθ1=V2cosθ2,

[0387] We have

[0388]

[0389] Here, if the propagation displacement wave propagates from A to B, V2 is considered to be positive, and if the propagation displacement wave propagates from B to A, V2 is considered to be negative.

[0390] If θ is chosen to be close to pi / 2, the expression becomes:

[0391]

[0392] Then

[0393]

[0394] if Then it becomes:

[0395]

[0396] In a practical setup, to obtain the local shear wave velocity, points A and B can be selected to be adjacent on a radial line originating from the center of the convex probe. To obtain a more robust estimate of the velocity V2 for velocity correction, points A and C can be further apart.

[0397] Therefore, the proposed shear wave front angle analysis allows us to:

[0398] ○ Eliminate deviations in the local value of shear wave velocity at each point in the medium

[0399] o Use the angle θ1 and the corresponding azimuthal coordinates of the scan line between A and B to obtain the local value of the shear wave front angle.

[0400] 3.1.9. Use of shear wave front angle as a quality factor and as an indicator of the pressure of the probe on the tissue

[0401] From the above, it is clear that in order to avoid the excitation of out-of-plane shear waves from the rib edges (due to mechanical coupling effects), a “light” probe pressure on the tissue is required.

[0402] One possible way to ensure this condition (without excluding other possibilities) is as follows:

[0403] To differentiate between “light” pressure and “high” pressure (causing rib coupling), we can use the following facts:

[0404] In some cases, the curvature (convexity) of the shear wave front is an indicator of pressure, especially when equal pressure is applied to the upper and lower ribs.

[0405] If the absolute value of the average shear wave front angle within the ROI is large (e.g., greater than 10 degrees), the pressure on one rib may be greater (asymmetric pressure between the upper and lower ribs),

[0406] ○ We can then perform analysis using the ROI located to the right or left of the probe center and determine, for example:

[0407] If the average shear wave angle within the ROI is within a certain range (e.g., between 0 and -10 degrees, the sign of the angle depends on the coordinate system used), then the pressure is appropriate ("light" pressure, no rib coupling),

[0408] o If the average angle is outside this range, the pressure may be too "high".

[0409] Those skilled in the art will note that this indicator can be combined with other quality factors, such as:

[0410] Probe deflection, which can be monitored based on displacement before shear wave extraction, and / or

[0411] ●Shear wave amplitude.

[0412] Those skilled in the art will understand that:

[0413] ○Quality factor related to the orientation and position of the probe,

[0414] ○ represents the coefficient of pressure applied by the probe to the patient, and

[0415] ○Parameters related to the contact quality between the probe and the patient,

[0416] The above-mentioned items described herein may be combined in any manner to obtain composite figures of merit and / or indices.

[0417] Those skilled in the art will appreciate that many modifications may be made to the invention described above without materially departing from the novel teachings and advantages described herein. Accordingly, all modifications of this type are intended to be incorporated within the scope of the appended claims.

Claims

1. A method for determining two-dimensional characteristics of a propagating displacement wave in a patient, the method comprising: ○ Estimating multiple velocity maps of the target area over time, the estimation stage includes the following steps: ○ Construct multiple echogenicity maps by repeatedly performing the following substeps for each echogenicity map: o transmitting an ultrasound signal (201) using a transducer element array comprising a set of transducer elements, o receiving a backscattered echo signal (202) by an array of transducer elements, each transducer element being capable of collecting a respective time domain signal from a set of time domain signals, the time domain signal corresponding to the amplitude of the received backscattered echo signal on the array of transducer elements, ○ constructing an echo intensity map (203) based on a set of time domain signals, o estimating multiple velocity maps (204) by correlating and processing the multiple echogenicity maps, ○ filtering (30) the propagating displacement waves within the plurality of velocity maps to obtain a plurality of filtered velocity maps, Determining a two-dimensional characteristic of a propagating displacement wave for each point of interest A within a plurality of velocity maps, the two-dimensional characteristic comprising a propagation velocity and a propagation direction of the propagating displacement wave, to obtain at least one two-dimensional characteristic map, wherein the point of interest A has the same coordinates within the plurality of velocity maps, by: o associating at least two points B, C within the plurality of velocity images / maps such that points A, B, C are not collinear (402), said points B, C having the same coordinates within the plurality of velocity images / maps, ○ For each of the three points A, B, and C, extract the time domain propagation displacement signal (403) from the filtered velocity maps, o Estimate the longitudinal apparent flight time and transverse apparent flight time (404) of the propagating displacement wave between points A, B, and C using the time domain propagating displacement signals associated with points A, B, and C by: ■ Cross-correlate the time-domain propagation displacement signals associated with points A and B to estimate the longitudinal apparent flight time of the propagation displacement wave from point A to point B, ■ Cross-correlate the time-domain propagation displacement signals associated with points A and C to estimate the apparent lateral flight time of the propagation displacement wave from point A to point C, ○ Determine the propagation direction and propagation speed of the propagating displacement wave based on the transverse apparent flight time and the longitudinal apparent flight time (405).

2. The method according to claim 1, further comprising the step of estimating the hardness of tissue included in the target area using the two-dimensional characteristic map.

3. The method according to claim 2, wherein: The steps to estimate tissue stiffness include: o extracting the velocity of the propagating displacement wave for a subset of points representing tissue included in the target region from the two-dimensional characteristic map, o Calculate the median or average of the velocity of the propagating displacement wave for the extracted subset of points, and o Estimate the hardness based on the calculated median or mean.

4. The method according to claim 2, wherein: The steps to estimate tissue stiffness include: ○Calculate the hardness map of the region of interest based on the two-dimensional characteristic map, o extracting the hardness of a subset of points from the hardness map, said subset of points representing tissue included in the target region, and ○ Calculate the median or average hardness of the extracted point subset.

5. The method according to any one of claims 1 to 4, further comprising the step of determining a quality factor associated with the mean propagation displacement direction, the step of determining the quality factor comprising: ○ Extract the direction of the propagating displacement wave of a subset of points from the two-dimensional characteristic map, ○ Calculate the average propagation direction of the propagation displacement wave based on the direction of the propagation displacement wave of the extracted point subset, ○ comparing the average propagation direction of the propagating displacement wave with a predefined range of directions to obtain comparative data, and ○ Derive quality factors based on comparative data.

6. The method according to any one of claims 1 to 5, further comprising determining a value representing the energy distribution of the transducer array (T1-T n ) the step of determining a coefficient of the pressure applied to the patient, said step of determining the coefficient comprising the following sub-steps: - Calculate the curvature of the wavefront of the propagating displacement wave based on the two-dimensional characteristic diagram, - if the calculated curvature is a concave shape corresponding to converging propagating displacement waves, the coefficient is assigned a value representing a light pressure, If the calculated curvature is a convex shape corresponding to a diverging propagating displacement wave, the coefficient is assigned a value representing a high pressure greater than a light pressure.

7. The method according to any one of claims 1 to 6, further comprising the step of determining a parameter representative of the quality of the propagated two-dimensional displacement estimate, said step of determining the parameter comprising the following sub-steps: - Calculate the amplitude of the change in the apparent transverse flight time, which represents the spatial variation in the direction of the propagating displacement wave, - comparing said amount with a predefined threshold value and determining a parameter as a function of the result of said comparison sub-step.

8. The method according to any one of claims 1 to 7, wherein The steps of estimating multiple velocity images / maps by correlating and processing multiple complex echo intensity maps include: - cross-correlating multiple complex echo intensity maps with a certain time lag to obtain multiple phase shift maps, -Derivation of multiple velocity maps from the correlated multiple phase shift maps.

9. The method according to any one of claims 1 to 8, wherein The sub-step of constructing an echogenicity map (a two-dimensional map) from a set of time-domain signals comprises solving an inverse problem to generate said echogenicity map.

Images