Copyright © IFAC Modelling and Control in Biomedical Systems, Melbourne, Australia, 2003
PUBLICATIONS www.elsevier.com/locale/i
NOISY HEART SEQUENCES FILTERING BY OPTIMAL MEDIAN BASED OPERATORS Rasto Lukae!, Bogdan Smolka 2, Pavol Zavarsky3, Ko~tas N. Plataniotis and Aoastasios N. Venetsanopoulos
4 ,
I Slovak Image Processing Center. Jarkova 343. 04925 Dobsina. Slovak Republic.
[email protected] 2 Department 0/ Automatic Control. Silesian University o/Technology. Akademicka 16 Str.. 44-101 Gliwice. Poland.
[email protected] 3 Nagaoka University of Technology. 1603-1 Kamitomioka. 940-2 I 88 Nagaoka. Japan.
[email protected] 4 The Edward S. Rogers Sr. Department 0/ Electrical and Computer Engineering. University o/Toronto. 10 King's College Road. Toronto. Canada. {kostas. anv}@dsp. toronto. edu
Abstract: Bit errors and impulses introduced into cardiographic image sequences prohibit to process and evaluate the heart dynamics correctly. Thus, biomedical imaging such as vascular imaging and quantification of heart dynamics is closely related to digital filtering. We provide three-dimensional impulse detection based adaptive median filters designed to perform the optimal filtering situation, where the robust median filter affects only noisy samples, whereas the desired image features remains unchanged. The proposed method holds excellent preserving characteristics and takes the advantage of the superposition property despite the nonlinearity of the median filtering. Copyright © 2003 IFAC Keywords: Detection, Impulses, Median Filters, Optimal Filtering, Sequences, Smoothing
This paper focuses on three-dimensional (3-D) adaptive median filters based on the impulse detection approach designed to effectively remove the impulsive component from noisy cardiographic image sequences. Impulse noise affects the useful information in the form of bit errors and it introduces to the image high frequency changes that prohibit to process and to evaluate the heart dynamics correctly. Therefore biomedical imaging such as vascular imaging and quantification of heart dynamics is closely related to digital filtering. In order to suppress impulse noise effectively, well-known nonlinear filters based on the robust order-statistic theory provide interesting results. Although median filters hold good impulse noise attenuation characteristics. their performance is often accompanied with undesired processing of noise-free samples, which results in edge and texture blurring.
I. INTRODUCTION In general, cardiographic image sequences (Bosnjak et aI., 2001; Plebe and Gallo.2001; Lukac,2002b; Lukac et aI., 2002) are often affecated by noise preventing an immediate and accurate tracking of the heart motion. In order to restore the heart dynamics, which is the most important a-priori knowledge (Bosnjak, 200 I) available for cardiographic image sequence. we provide a set of adaptive median filters. Since image sequences, i.e. the time-varying images, can be considered as a time sequence of twodimensional (2-0) images, the spatiotemporal filters (Arce, 200 I; Kleihorst et aI., 1995) utilising both temporal correlation of frames and inherent spatial correlation present in frames of image sequences, can produce the best estimate.
311
/
T7/
r7 L7
LJLJ r7 L7
rIU
LJ TI
CITY
o/J7J
/J~I r1 2T ICl 7 err]
I
/1//
1~7J
,/1//
I
m'I·,er U 4:10
k
la)
/7/7
image sequences, these filters process each frame independently (Fig. Ib), so that the excluded information about the motion trajectory can result in a large motion blurring, (Arce, 2001).
Ifl 1 17/7
nz r
rT7
.,
ITI 1 ITI I I 1 I ] 1 r7 I
FT} J
/777
It can be easily seen, that 3-D filters (Kleihorst et al. 1995, Kokaram,1998; Nikolaidis and Pitas, 2001) also called spatiotemporal filters (Fig. Ic) utilize both temporal correlation and spatial correlation present in image sequences. For that reason, spatiotemporal filters constitute a natural filter class for noisy image sequences and it is possible to observe the significant reduction of the spatial and motion blurring in the dependence on the chosen filtering algorithm.
)
(0)
(c)
Fig. I. Image sequence filtering divided according to dimensionality of the filter window: temporal (a), spatial (b) and spatiotemporal filters (c). The reason is that the median filters do not satisfy the superposition property and thus, the optimal filtering situation can never be fully obtained.
Since the degree of noise corruption significantly influences the motion estimation precision, the image sequence filtering algorithms are usually referred with no motion estimation. Thus, a 3 x 3 x 3 cube spatiotemporal filter window represents the most frequently used window shape, because it respects all kinds of the correlation present in image sequences.
In the case of the impulse noise corruption, the aim of the optimal filtering (Astola and Kuosmanen, 1997) is to design the filtering algorithm that will affect only corrupted samples, whereas the desired (noise-free) samples will be invariant to the filtering operation. This is realised by adaptive impulsedetection based median filters, which replace the noisy samples by the median of the input set and perform the identity operation related to noise-free of the output samples. The measure distortiondepends on the capability to detect atypical image samples - impulses and outliers, that can be very similar to samples belonging to an edge.
3. IMPULSE DETECTORS In the case of heavy tailed noise, e.g. impulse noise with uniform distribution of random values, the popular class of linear filters exhibits worse noise attenuation characteristics. In order to suppress the impulse noise effectively, a wide class of nonlinear order-statistic filters (Astola and Kuosmanen, 1997), based on the ordering (Pitas and Venetsanopoulos, 1992) of the input samples spanned by a filter window is preferred. The ordering operation removes atypical image samples, often noise, to the borders of the ordered set and the mid-positioned samples, e.g. the median value, of the ordered sequence represent the robust estimates.
In this paper, the performance of the proposed approaches is tested for the heart image sequence of 38 frames with a wide range of the noise corruption. The results are evaluated in terms of mean absolute error, mean square error and cross correlation. 2. IMAGE SEQUENCE FILTERING In terms of signal dimensionality and correlation present in the input set, the filtering techniques for noisy image sequences or generally image sequences can be divided into three classes (Arce, 1991; Kokaram, 1998; Lukac and Marchevsky, 200 Ia; Nikolaidis and Pitas.2001): temporal filters, spatial filters and spatiotemporal filters. This division is based on the fact that image sequences represent time sequences of 2-D images, i.e. the spatiotemporal data.
Probably the most popular nonlinear filter is the wellknown median (Astola and Kuosmanen, 1997; Pitas and Venetsanopou!os, 1992; Lukac and Marchevsky, 200 Ib) characterised by excellent robustness against impulse noise. However, it often introduces to an image too much smoothing resulting in a blurring that can be more perceivable and disturbing than the original noise. Since the nonlinear filters (Peltonen et aI., 2001, Lukac 2002a) do not satisfy the superposition property, the optimal filtering situation, where only noisy samples are filtered whereas the desired image features will be invariant to a filtering operation, can never be fully obtained. For that reason there were developed some adaptive impulse detection-based median filters (Fig.2) taking the advantage of the optimal filtering situation.
The class of temporal filters (Fig. I a) is referred to the temporal correlation of frames. These onedimensional (I-D) filters remove noise without impairing the spatial resolution in stationary areas. In order to improve a filter performance, the temporal filtering is usually connected with a motion compensation method (Kleihorst et aI., 1995; Klima et at.. 1994; Kokaram, 1998; Zahradnik, 1994) so as to filter objects along their motion trajectory. However, this way is very computationally complex and it often does not work well because of the spatial warping and scene changes. In general, 2-D filters (Astola and Kuosmanen, 1997) are probably the most frequently used image filtering techniques used in image processing. In the case of
Fig. 2. Impulse detector based filtering.
312
The sigma filtering concept (Lee, 1983), i.e. the filtering based on the simple statistical measures such as the mean value and standard deviation is given by the operation value defined as
In general, the detection-filtering algorithm can be stated (Smolka et ai, 2002) as follows:
IF
Val ~ Tol
TIIEN ELSE
y =
JIIF
(I)
y=X,V+1I12 Val=lfi-x(N+I)/21
where Val is the detector operation value that is compared with the threshold value Tol . If the condition Val ~ Tol is satisfied, the central sample X(N+I)/2 of the filter window represents the impulse and it is estimated by median filter, i.e. the filter output y = y MF . Otherwise, the adaptive median filter performs no smoothing or identity operation, where the central sample X(N+I)/2 remains unchanged, i.e. Y = X(N+I)/2'
where
I Tol =
Let x"x2 "",XN be a discrete-time continuos-valued input set determined by a filter window and x(I)'x(2)"",X(N) an ordered set so that X(I) ::; x(2) ::; ... ::; X(N)
fori=I,2, ...,N
N
-
3.4 Local contrast probability based detector (Lep)
The idea of the LCP detector (Beghdadi and Khellaf, 1997) is based on fixed threshold Tol = 1/ N computed from the window size N and the operation value given by
(2)
(3)
The median filter output is given by YMF = X«(N+I)/2)
Val =
(4)
where X(N+II/2) is the middle positioned orderstatistic, i.e. the median.
x(N+I)'21
The original heart image sequence of 38 frames of the size 256 x 256 image samples and 8 bit per sample representation was used for the simulations, (Fig.3 a,b). The original signal was cOffilpted by 5% and 10% (Fig.3c) random-valued impulse noise (Astola and Kuosmanen, 1997; Boncellet, 2000) defined as
Now, let us consider the set of n mid-positioned ordered samples x((n+I)/21' x(I+(n+II/2)' ..• , x(,V-(n+I)/21 and the simple trimmed mean defined as
L:
X'I
(6)
Since, the extreme order-statistics (usually outliers) are excluded and the trimmed mean fin is determined by probably no cOffilpted samples, it will provide the precise value necessary to be compared with the central sample X(N+II,2 . Thus, the detector value (Park and Kurz, 1996) can be stated as follows -XIN + 1)121
={
o',) u
with probability I - P..
.
with probability
(11)
Pv
where X"j is the noisy image sample, 0,.) describes the sample original image, i,j are indices of the sample location, u is the random value from <0,255> and Po is the impulse probability. Note that impulse noise frequently occurs as bit errors (Astola and Kuosmanen, 1997):
1=(n ... I).l2
Val = lfin
(10)
4. EXPERIMENTAL RESULTS
(5)
3.2 Order-statistic detector (OD)
n
fil
where fi is the sample mean of observed data x, ,x2 '''''X N ' If the equation (10) results in the value greater than or equal to the threshold value Tol = 1/ N, the central sample is considered as the impulse and the output value is determined by the median.
i.e. the design parameter Val is equal to the absolute difference between the median output YMF and the central sample X(N+I)/2 . Our experiments indicate, that for the majority of noisy scenarios the threshold value Tol is optimally set as Tol = 20.
I V-t n .1)12 fin =XI')
IX(N+I)/2 N
L:lxi -fil
The advantage of the median operator (pitas and Venetsanopulos, 1992) lies in the robust estimate in the environments cOffilpted by impulse noise, because the median value has the highest probability to be the noise-free sample. For that reason, the median value can be used successfully in the equation (I) so that the parameter Val is given by Val = IYMF -
(9)
L:(xi - fi)2
N i=1 If the detection value is greater than or equal to the standard deviation or the adaptive threshold Tol, the central sample is probably distorted, because it is more different from other input samples.
3.1 Median detector (MD)
X(i)E{X"X2 , ...,XN }
(8)
is the sample mean of observed data x" x 2' ..., x" and the adaptive threshold value Tol is given by fi
'kM '.)
=
r { 1- k,~
with probability 1- Pv
'I
with probability
Pv
(12)
where Po is the bit change probability, k'~1 and 'ki~) are binary values {O, I} of B -bit original sample 0,) and noisy sample x,,/ given by o = k'I.J 2 8 - 1 + k' 2 8 - ' + ... + k 8 -'2 + k 8 (13)
(7)
In the case of cube filter window, the OD optimally works for the threshold Tol = 15 and the number of considered mid-positioned samples n = 13 .
I.}
I.)
t,)
I.}
(14) The objective measure of the filters performance were evaluated according to three error criteria
3,3 Standard deviation based detector (SD)
313
(a)
(b)
(d)
(e)
(c)
(f)
Fig. 3. Achieved results. (a) Original 5th frame, (b) Original 35th frame, (c) Noisy (10% impulse noise) 5th frame, (d-f)image Fig. 3c filtered by: (d) 2-D median, (e) 3-D median, (f) 3-D OD + 2-D median Spatial Median Median MD based median SD based median OD based median LCP based median
(Lukac and Marchevsky, 200Ia): mean absolute error (MAE), mean square error (MSE) and cross correlation M , computed for whole sequences. Mathematically, the MAE and MSE are defined as I
N
M
K
MAE=-""""""Io. NMK f::~f:: l,;.! -xi,j,/ 1
N
M
K
MSE=--""""""(o. NMK f::~f:: I.J.t
I
NM L.L. '.; ,:,
-x) t.}.t
Method / Criteria MAE MSE Identity 5.105 572.3 Spatial Median 1.219 5.9 Median 1.473 7.3 MD based median 0.311 4.1 SD based median 0.386 6.6 OD based median 0.379 2.7 LCP based median 0.418 3.7 Table 3 Achieved results using 10% impulse
(16)
XI.'-EIEI·'I
Method / Criteria Identity Spatial Median Median MD based median SD based median OD based median LCP based median
'.1
(17)
;=1
at
0"+1
where E' is the mean value and a ' is the standard deviation of the tth frame and E'·' is the mean value and a'" is the standard deviation of the (t + I) th frame. In general, MAE is a mirror of the signal-details preservation, MSE well evaluates the noise suppression and M expresses the preservation of the motion trajectory in an image sequence. Thus, using the three measures, the quality of the filtered cardiographic image sequences is quantified with a high accuracy.
MAE
MSE
MAE 9.940 1.342 1.531 0.375 0.757 0.489 0.564
MSE 1506.8 8.2 7.7 4.7 25.9 4.1 12.0
~R
0.345 0.001 0.003 0.002 0.003 0.001 0.001 noise ~R
0.504 0.001 0.003 0.001 0.015 0.001 0.006
The achieved results are depicted in Tables 1-3 and shown in Fig.3d-f and Fig.4. It can be seen that the median (Fig.3d,e and FigAc) allows to suppress noise (Fig.4b) effectively, however it results in significant blurring of fine edges. The adaptive impulse detector-based median filters such as MD and OD can achieve excellent signal-detail preservation with the simultaneous noise suppression (Fig.3f and Fig.4e), however other approaches
Table I Achieved results related to original sequence Method / Criteria Identity
0.001 0.003 0.001 0.001 0.001 0.001
Table 2 Achieved results using 5% impulse noise 2
filtered (noisy) image, i,j,t are indices of the sample position and N, M, K characterise the sequence size. The cross correlation M is given by
l _l_~~XI
5.3 6.9 0.9 3.0 2.4 3.5
(15)
where {o.,j} is the original image, {X,./} is the
R' =
1.136 1.432 0.023 0.520 0.268 0.652
~R
314
(a)
(b)
(c)
(d)
(e)
(f)
Fig. 4. Detailed view (a) Original, (b) 10% Impulse noise, (c) 2-D median, (d) SD, (e) OD, (f) LCP
Bosnjak, A., Torrealba, V., Montilla, G., Villegas, H., Burdin, V., Solaiman, B., Roux, C. (2001). Segmentation, Reconstruction, Modeling and 3D Visualization of the Ventricles in Echocardiographics Images. Proc. ISPA'OI, 260-265. Boncelet, C. (2000), Image Noise Models, in Handbook of Image & Video Processing (ed. A. Bovik), Academic Press, 325-336. Kleihorst, R.P., Lagendijk, R.L., and Biemond, 1. (1995). Noise Reduction of Image Sequences Using Motion Compensation and Signal Decomposition. IEEE Trans. on lP, 4,274-284. Klima, M.. Dvorak, P., Zahradnik, P., Kolar, 1., and Kott, P. (1994). Motion Detection and Target Tracking in a TV Image for Security Purposes, Proc. IEEE Carnahan Conference on Security Technology 1994 in Albuquerque USA, 43-44. Kokaram, A. (1998). Motion Picture Restoration. Springer-Verlag, London.
(Fig.4d,f) provide insufficient noise attenuation capability because these approaches are more appropriate for standard video applications. 5. CONCLUSION This paper focused on the reconstruction of the heart image sequences corrupted by bit errors or impulse noise. Since standard robust nonlinear filtering approaches, e.g. well-known median filters, blur important edges of the heart images and affect the details of heart sequences, in order to allow immediate and accurate tracking of the heart motion, we have tested four impulse detector based adaptive median filters that provides much better signal detail preservation than the standard median filtering. The most appropriate adaptive approaches are represented by OD and MD medians. In addition, these filters can be effectively implemented through binary operations (Lukac, 2002a).
Lee, l.S., (1983). Digital Image Smoothing and the Sigma Filter. Computer Vision, Graphics. and Image Processing, 24, 255-269. Lukac, R., and Marchevsky, S. (2001a). LUM Smoother with Smooth Control for Noisy Image Sequences. EURASIP J. ofASP, 2001, 110-120. Lukac, R. and Marchevsky, S. (200Ib). Boolean Expression of LUM Smoothers. IEEE Signal Processing Letters, 8, 292-294. Lukac, R. (2002a). Binary LUM Smoothing. IEEE Signal Processing Letters, 9, 400-403. Lukac, R. (2002b). Reconstruction of Noisy Hearth Dynamics by Optimised Permutation Filters. Proc. BIOSIGNAL '02. 289-291.
REFERENCES Arce, G.R., (1991). Mulstistage Order Statistic Filters for Image Sequence Processing. IEEE Trans. on Signal Processing, 39,1146-1163. Astola, 1., and Kuosmanen, P. (1997) Fundamentals ofNonlinear Digital Filtering. CRC Press. Beghdadi. A., and Khellaf, A. (1997). A NoiseFiltering Method Using a Local Information Measure, IEEE Trans. on lP, 6, 879-882.
315
Statistics in Digital Image Processing. Proceedings ofthe IEEE, 80,1892-1919. Plebe, A., and Gallo, G. (2001). Filtering Echocardiographic Image Sequences in Frequency Domain. Proc.ISPA'OI, 238-243. Smolka, 8., Studer, M., Stomrnel M., and Wojciechowski, K. (2002). Vector Median Based Color Gamma Filter. Proc. ICCVG 2002, 677-684. Zahradnik, P., Klirna, M., and Cuda, J. (1995). Real Time Video Motion Detection based on the TMS 320C50 Signal Processor, Proceedings of The Texas Instruments Educators Conf, 51-53.
Lukac, R., Smolka, B., and Galajda, P. (2002). Hearth Dynamics Filtered by Optimised Adaptive Median. Proc. EMBEC'02, 3, 954-955. Nikolaidis, N., and Pitas, 1. (2000) 3-D Image Processing Algorithms. Wiley. Park, 1., and Kurz, L. (1996). Image Enhancement Using the Modified ICM Method. IEEE Transactions on Image Processing, 5, 765-771. Peltonen, S., Gabbouj, M., and Astola, J. (2001). Nonlinear Filter Design: Methodologies and Challenges. Proc.ISPA'OI, 102-106. Pitas, 1. and Venetsanopoulos, A.N. (1992). Order
316