This disclosure is directed to methods for digital subtraction in medical images.
Digital subtraction is a common procedure for many medical imaging applications. For example, before a medical intervention procedure, contrast media is injected into blood vessels to highlight the vessels and then flushed away. A fluoroscopy image with vessels highlighted can be obtained, which is known as a mask image. Then, during the intervention, a guide wire may be inserted into the vessels, and a fluoroscopy image of the vessels is obtained in real-time. To highlight the vessel's location, the mask image is subtracted from a current image that shows the presence of guide-wire. The final result is an image of the guide-wire overlaid with the vessels to help a physician position the guide wire. It can be represented as following equation:
Isub=Icur−Imask
However, since imaging noise on the two images, i.e., the mask image and the current image, is independent of each other, the subtraction result has more noise than either of the two input images:
σsub2=σcur2+σmask2.
Reducing the noise improves image quality. One common way to reduce noise is to apply temporal filtering to average over consecutive frames and hence reduce the noise. An issue with utilizing consecutive frames, however, is preventing ghosting artifacts when there is significant guide wire motion.
Since vessels usually cannot be clearly seen in traditional X-ray fluoroscopy images, contrast media is injected to highlight the vessels. Vessel related intervention procedures are usually divided into several phases. In a first phase, contrast media is injected into the vessels to highlight them. Vessels with contrast media appear darker than the surrounding soft tissues. Then, after a short time, the contrast media is flushed away. During this phase, a mask image can be obtained by taking the minimum intensity of each pixel for all image frames acquired during this phase. A sample mask image is shown in
Exemplary embodiments of the invention as described herein generally include methods and systems for a decomposed temporal filtering for fluoroscopy guided intervention application that utilize the availability of the mask image to decompose the intervention images into moving structures and static structures. Performing temporal filtering only on static structures can prevent ghosting in traditional temporal filtering. A filtering algorithm according to an embodiment of the invention can be performed in real-time during an intervention procedure, and can be extended to other imaging modalities which utilize digital subtraction procedures.
According to an aspect of the invention, there is provided a method for temporally filtering medical images during a fluoroscopy guided intervention procedure, the method including providing a mask image and a fluoroscopy intervention image acquired at a current time during a medical intervention procedure, forming a subtraction image by subtracting the mask image from the intervention image, calculating a motion image of a moving structure in the subtraction image, forming a residual image by subtracting the motion image from the subtraction image, temporally filtering the residual image with a filtered image from a previous time, and adding the motion image to the temporally filtered residual image.
According to a further aspect of the invention, the method includes normalizing contrast of the mask image and the intervention image, before forming the subtraction image.
According to a further aspect of the invention, contrast is normalized by calculating a log of each image.
According to a further aspect of the invention, the moving structure comprises a guide wire being inserted into a blood vessel during the medical intervention procedure.
According to a further aspect of the invention, calculating a motion image comprises calculating an image that minimizes an objective function that includes a mutual information of the subtraction image and the intervention image.
According to a further aspect of the invention, the mutual information is of the form MI(Isub,t(x,y)−Imov,t(x,y), It(x,y)−Imov,t(x,y)), where (x,y) is a pixel location, Imov,t(x,y) is the motion image, Isub,t(x,y) is the subtraction image, and It(x,y) is the intervention image.
According to a further aspect of the invention, the objective function being minimized includes term proportional to a square of the motion image.
According to a further aspect of the invention, the mutual information is of the form
where fk(N1)=N1K, gk(N2)=N2k and N1, N2 represent the first and second arguments of the mutual information, respectively.
According to a further aspect of the invention, the method includes calculating a motion consistent image Ĩmov,t(x,y)=arg(Ĩmov,t)min(∥Ĩmov,t(x,y)−Imov,t(x,y)∥2+λ∥Ĩmov,t(x,y)−Ĩmov,t-1(x,y)∥2) that constrains the motion image to be substantially constant if the moving structure is not moving, where λ is a constant set to 0 when ∥Imov,t(x,y)−Ĩmov,t-1(x,y)∥>ε, where ε is a threshold for motion detection.
According to a further aspect of the invention, the temporal filtering takes the form It(x,y)=wIt-1(x,y)+(1−w)(Isub,t(x,y)−Imov,t(x,y))
where (x,y) is a pixel location, Imov,t(x,y) is the motion image, Isub,t(x,y) is the subtraction image, It-1(x,y) is the filtered image from a previous time, and w is a predetermined weighting factor.
According to another aspect of the invention, there is provided a method of temporally filtering medical images during a fluoroscopy guided intervention procedure, the method including providing a fluoroscopy intervention image of a blood vessel acquired at a current time during a medical intervention procedure and a subtraction image of the blood vessel, calculating a motion image that minimizes an objective function that includes a mutual information of the subtraction image and the intervention image, where the mutual information is of the form MI(Isub,t(x,y)−Imov,t(x,y), It(x,y)−Imov,t(x,y)), where (x,y) is a pixel location, Imov,t(x,y) is the motion image, Isub,t(x,y) is the subtraction image, and It(x,y) is the intervention image, forming a residual image by subtracting the motion image from the subtraction image, and temporally filtering the residual image with a filtered image from a previous time.
According to a further aspect of the invention, the subtraction image comprises providing a mask image and forming the subtraction image by subtracting the mask image from the intervention image, where the mask image is provided before the intervention by injecting a contrast media into the blood vessel and acquiring an image of the highlighted vessels.
According to a further aspect of the invention, the method includes adding the motion image to the temporally filtered residual image.
According to another aspect of the invention, there is provided a program storage device readable by a computer, tangibly embodying a program of instructions executable by the computer to perform the method steps for temporally filtering medical images during a fluoroscopy guided intervention procedure.
a)-(b) compare results of direct subtraction with a subtraction method according to an embodiment of the invention.
a)-(b) compare results of traditional temporal filtering with a decomposed temporary filtering method according to an embodiment of the invention.
Exemplary embodiments of the invention as described herein generally include systems and methods for decomposed temporal filtering for fluoroscopy guided intervention procedures. Accordingly, while the invention is susceptible to various modifications and alternative forms, specific embodiments thereof are shown by way of example in the drawings and will herein be described in detail. It should be understood, however, that there is no intent to limit the invention to the particular forms disclosed, but on the contrary, the invention is to cover all modifications, equivalents, and alternatives falling within the spirit and scope of the invention.
As used herein, the term “image” refers to multi-dimensional data composed of discrete image elements (e.g., pixels for 2-D images and voxels for 3-D images). The image may be, for example, a medical image of a subject collected by computer tomography, magnetic resonance imaging, ultrasound, or any other medical imaging system known to one of skill in the art. The image may also be provided from non-medical contexts, such as, for example, remote sensing systems, electron microscopy, etc. Although an image can be thought of as a function from R3 to R or R7, the methods of the inventions are not limited to such images, and can be applied to images of any dimension, e.g., a 2-D picture or a 3-D volume. For a 2- or 3-dimensional image, the domain of the image is typically a 2- or 3-dimensional rectangular array, wherein each pixel or voxel can be addressed with reference to a set of 2 or 3 mutually orthogonal axes. The terms “digital” and “digitized” as used herein will refer to images or volumes, as appropriate, in a digital or digitized format acquired via a digital acquisition system or via conversion from an analog image.
In the following disclosure, M(x,y) denotes a mask image as a function of pixel location (x,y), It(x,y) refers to the fluoroscopy image acquired during the intervention at time t as a function of pixel location (x,y), and Isub,t(x,y) is the pixel-wise subtraction result of It(x,y) and M(x,y).
Contrast Normalization
Referring now to
D=Aexp(−αs)
When the X-ray finally hits the detector, its energy is measured and represented as an intensity value I. This value reflects the attenuation through the whole path of the X-ray:
I=A·exp(−∫αsds).
The attenuation coefficient αs is different for different tissues. To maintain consistent contrast of the objects (i.e., guide wire or vessels) during the subtraction procedure, the image should be normalized by a log function before subtraction:
Isub,t(x,y)=log(It(x,y))−log(M(x,y)) (1)
This log-normalization makes the image additive, since before normalization, due to the exponential absorption, the intensity is multiplicative. After this normalization step, a guide wire passing through region of different intensity will have the same intensity, different than the neighborhood. According to some embodiments of the invention, the log-conversion may be performed with a look-up-table. According to other embodiments of the invention, the look-up-table is not an exact log function, but is modified to better handle very dark regions.
Decomposed Filtering
Temporal filtering is a common technique for reducing imaging noise by averaging over consecutive image frames. If there is motion during image acquisition, registration algorithms can be used to compensate for the motion before the temporal filtering.
However, motion compensation is challenging in X-ray fluoroscopy imaging. Since X-rays pass though body tissue, objects inside the body, such as organs, guide wires, soft tissue, etc., are projected and overlaid onto one 2D image as semi-transparent layers. Usually the patient stays static, but the guide wire is moving inside a vessel. Registration algorithms cannot typically register images with multiple transparent layers to compensate the motion because most registration algorithms assume one motion vector at each pixel.
According to an embodiment of the invention, a mask image M(x,y) is used at step 42 to decompose the intervention images It(x,y) into two layers, a static body tissue layer and a moving guide wire layer, before the temporal filtering. This scheme can handle the different motion of the body tissue and guide wire, reducing noise while preventing ghosting artifacts caused by the moving guide wire.
Traditionally, direct subtraction, i.e., Isub,t(x,y)=It(x,y)−M(X,y), is used to overlay the moving guide wire with mask image of the vessel tree during the intervention procedure, which is noisy due to the subtraction. Directly applying temporal filtering on It(x,y) before subtraction or on Isub(x,y) is challenging because the images contain both moving guide wire and static structures.
However, according to an embodiment of the invention, utilization of the original mask image and intervention images can separate the moving structures, i.e., Ĩmov,t(x,y), and remove the moving structures from the images before temporal filtering. After the filtering, the moving structures are recombined with the temporal filtered images to form the final output results. A temporal filtering scheme according to an embodiment of the invention can reduce imaging noise while preventing ghosting artifacts generated by the moving structures.
The detailed steps for separating the moving structures and static structures are described in the following sections.
Separate Moving Layer Based on Mutual Information
The direct subtraction image Isub,t(x,y) contains both moving structures, such as the guide wire, and static structures, such as vessels. A moving structure, such as the guide wire, is inserted during intervention, and hence also exists in the intervention image. Therefore, at step 43, it is possible to estimate the moving structure, Imov,t(x,y) based on a mutual information analysis between a subtraction image Isub,t(x,y) and an intervention image It(x,y):
where MI is the mutual information and Imov,t(x,y) is the motion image of the moving structure at time t. Since the motion image Imov,t(x,y) is inserted during intervention, and is the only common signal between the intervention image It(x,y) and the subtraction image Isub,t(x,y), removing the moving structure from the two images should minimize the mutual information between them. The second term is a regularization term that prevents imaging noise from being included in the moving structure.
Mutual information can be computationally intensive to calculate. For more efficient computation, the mutual information term can be approximated according to an embodiment of the invention based on assuming the independence of random variables:
E(f(N1),g(N2))=E(f(N1))·E(g(N2)).
where E is an expectation and N1, N2 are the random arguments. Based on this property, one can measure the mutual information, based on the following objective function:
where
fk(N1)N1K,gk(N2)=N2k.
To find an optimal solution, one can take the derivative of the objective function and set it to zero, which will yield an optimal solution.
Enforce Motion Consistency
Apart from the mutual information constraint on the moving structure, there are also temporal constraints. For example, the intensity of a moving guide wire remains substantially constant if it is not moving. But if the guide wire moves, the intensity of the affected pixels changes. This constraint can be enforced at step 44 by optimizing the Do following cost function:
where the result is a motion consistent image Ĩmov,t(x,y). The first term in the cost function ensures that the motion consistent image looks like the input, Imov,t(x,y). The second term enforces a similar intensity constraint on the motion consistent image. As described above, the similar intensity constraint should be applied only when the guide wire is not moving. So, λ may be to zero if ∥Imov,t(x,y)−Ĩmov,t-1(x,y)∥>ε, where ε is a threshold for motion detection.
Temporal Filtering
An estimated moving structure can be removed from the direct subtraction image so that and temporal filtering may be applied on the residual image Isub,t(x,y)−Ĩmov,t(x,y). This way, the temporal filtering is free of motion artifacts such as “ghosting”. The temporal filter is implemented based on a running average of the input images and the strength of the temporal filter can be controlled by a weighting factor w:
Ĩt(x,y)=wĨt-1(x,y)+(1−w)(Isub,t(x,y)−Ĩmov,t(x,y)), (5)
where Ĩt-1(x,y) is the temporally filtered image from the previous time step, and Isub,t(x,y) is the subtraction image calculated in step 42. The final output can be obtained at step 46 by adding back the motion consistent image from step 44 as follows:
ĨtO(x,y)=Ĩi(x,y)+Ĩmov,t(x,y). (6)
The resulting image ĨtO(x,y) should have significantly reduced image noise while preserving the moving structure without motion artifacts. The temporal filtering of EQ. (5) is exemplary and non-limiting. According to other embodiments of the invention, the temporal filtering involves multiple frames in the sequence. According to other embodiments of the invention, the weighting scheme can be adaptive.
Experimental Results
An algorithm according to an embodiment of the invention was tested on real image sequences. For a first sequence, a direct subtraction result was compared with a subtraction result according to an embodiment of the invention. The results are shown in
Another test compared a filtering method according to an embodiment of the invention with a traditional temporal filtering method. The results are shown in
System Implementations
It is to be understood that embodiments of the present invention can be implemented in various forms of hardware, software, firmware, special purpose processes, or a combination thereof. In one embodiment, the present invention can be implemented in software as an application program tangible embodied on a computer readable program storage device. The application program can be uploaded to, and executed by, a machine comprising any suitable architecture.
The computer system 71 also includes an operating system and micro instruction code. The various processes and functions described herein can either be part of the micro instruction code or part of the application program (or combination thereof) which is executed via the operating system. In addition, various other peripheral devices can be connected to the computer platform such as an additional data storage device and a printing device.
It is to be further understood that, because some of the constituent system components and method steps depicted in the accompanying figures can be implemented in software, the actual connections between the systems components (or the process steps) may differ depending upon the manner in which the present invention is programmed. Given the teachings of the present invention provided herein, one of ordinary skill in the related art will be able to contemplate these and similar implementations or configurations of the present invention.
While the present invention has been described in detail with reference to exemplary embodiments, those skilled in the art will appreciate that various modifications and substitutions can be made thereto without departing from the spirit and scope of the invention as set forth in the appended claims.
This application claims priority from “Decomposed Temporal Filtering Method for Roadmap Application”, U.S. Provisional Application No. 61/104,339 of Chen, et al., filed Oct. 10, 2008, the contents of which are herein incorporated by reference in their entirety.
Number | Name | Date | Kind |
---|---|---|---|
3974386 | Mistretta et al. | Aug 1976 | A |
4482918 | Keyes et al. | Nov 1984 | A |
5018179 | Kaneko | May 1991 | A |
5671263 | Ching-Ming | Sep 1997 | A |
5671743 | Kawai et al. | Sep 1997 | A |
6154518 | Gupta | Nov 2000 | A |
7346198 | Oosawa | Mar 2008 | B2 |
20080025588 | Zhang et al. | Jan 2008 | A1 |
20080037844 | Baumgart | Feb 2008 | A1 |
20080049994 | Rognin et al. | Feb 2008 | A1 |
20080205591 | Ozawa | Aug 2008 | A1 |
20090010512 | Zhu et al. | Jan 2009 | A1 |
Number | Date | Country | |
---|---|---|---|
20100092061 A1 | Apr 2010 | US |
Number | Date | Country | |
---|---|---|---|
61104339 | Oct 2008 | US |