- Email: [email protected]
Spectrogram analysis of arterial Doppler signals for off-line automated HITS detection
Ultrasound in Med. & Biol., Vol. 25, No. 3, pp. 349 –359, 1998 Copyright © 1999 World Federation for Ultrasound in Medicine & Biology Printed in the USA. All rights reserved 0301-5629/99/$–see front matter
PII S0301-5629(98)00173-2
● Original Contribution SPECTROGRAM ANALYSIS OF ARTERIAL DOPPLER SIGNALS FOR OFF-LINE AUTOMATED HITS DETECTION EMMANUEL ROY,*† SILVIO MONTRE´ SOR,* PIERRE ABRAHAM† and JEAN-LOUIS SAUMET† *Laboratoire d’Informatique de l’Universite´ du Maine, rue Lae¨nec, BP 535, 72085 Le Mans Cedex 9 France; and † Laboratoire de Physiologie et d’Explorations Fonctionnelles Vasculaires, Centre Hospitalier Universitaire, 49033 Angers Cedex 01 France (Received 7 July 1998; in final form 27 October 1998)
Abstract—Recently, a time processing of arterial Doppler signals was proposed to detect automatically highintensity transient signals (HITS). This technique provided satisfactory detection results, but was not always constantly accurate, particularly with high-resistance blood velocity profiles. A time-frequency processing, based on the spectrogram, is presented to detect the presence of emboli in the arterial Doppler signals. The method uses the narrow-band hypothesis and extracts the detection criterion from the time-frequency representation (TFR). A first database of 560 peripheral arterial Doppler HITS was created to study microemboli and to define the normal limits to be used in our method. A threshold was experimentally defined using this database, and then applied to 38 recordings from 12 patients. Using another database, 6 human expert Doppler users reported 140, 176, 155, 161, 161 and 146 HITS, corresponding to a total of 197 different observed HITS. When an event was detected by at least 6, 5, 4, 3, 2 and 1 of the observers, sensitivity of the automatic detection was 93.9, 91.7, 89.6, 88.7, 84.7 and 73.1%, respectively. The sensitivity of our automatic detection is, thus, highly associated with the number of observers in agreement. A preliminary experiment has been performed to test the method in the case of long recording duration. In 15 patients, 6 h 24 min of recordings have been analyzed. The proposed automated processing provided an overall sensibility of 91.5%. The present work shows that this detection scheme preserves good sensibility and improves the positive predictive value compared with the time-processing recently proposed. © 1999 World Federation for Ultrasound in Medicine & Biology. Key Words: Spectrogram, Doppler ultrasound, Embolism, Detection.
depending on their size, characteristics, embolus-to– blood ratio and the ultrasonic carrier frequency. Within the appropriate dynamic range of bidirectional Doppler equipment, an embolic signal is unidirectional within the Doppler velocity spectrum. HITS are always audible, and their sound depends on their velocity. Arterial embolisms can cause serious damage in the embolized organ in both cerebral (Spencer 1997; Tegeler 1994; Siebler et al. 1992) and peripheral (Abraham et al. 1997) circulation. Arterial emboli may involve the extremities, especially the legs and feet. Less common sites include the kidneys, intestines and other areas. Thus, the continuous monitoring of patients with significant risk of thromboembolic complications is of major interest. Several commercial software systems have been developed to detect automatically cerebral embolic signals with pulsed Doppler devices (Siebler et al. 1994; Van Zuilen et al. 1996). Automated techniques are needed for detecting peripheral emboli with continuous Doppler systems. Recently, Roy et al. (1998) proposed a
INTRODUCTION The possibility of detecting microemboli using Doppler ultrasound was reported in the 1960s by Spencer and co-workers during hyperbaric decompression experiments (Spencer and Campbell 1968; Spencer et al. 1969). The passage of a solid embolus through the beam of a Doppler ultrasound instrument appears in the Doppler spectrum as a high intensity transient signal (HITS). Embolic signals, when recorded optimally, have been defined based on four main features (Consensus Committee NICHS 1995; Nicholls et al. 1996; Spencer 1996). They are transient (usually , 300 ms) and their duration depends on their time of passage through the Doppler sample volume. Their amplitude is usually at least 3 dB higher than that of the background blood-flow signal, Address correspondence to: Professor J. L. Saumet, Centre Hospitalier Universitaire, Laboratoire de Physiologie et d’Explorations Fonctionnelles Vasculaires, 49033 Angers Cedex 01, France. E-mail: [email protected] 349
350
Ultrasound in Medicine and Biology
time-processing approach to perform the detection. They assumed that the Doppler signal could be considered as a narrow-band signal when an embolus occurs. The detection procedure consisted in retaining only instants where both the variation of the signal envelope was greater and the variation of the instantaneous frequency estimation was lower than respective variations of the background Doppler signal. This preliminary technique provided satisfactory detection results compared with detections made by six human experts. Nevertheless, the velocity profiles in the peripheral circulation vary widely compared with the those in the cerebral circulation, depending on the artery recorded, from one subject to another, or in the same subject over time. As mentioned by the authors (Roy et al. 1998), the originally reported method was not always stable, particularly with high-resistance blood-velocity profiles, which entailed false-positives. In this work, time-frequency processing is presented to detect the presence of emboli in arterial Doppler signals. The method uses the narrow-band hypothesis previously reported, but extracts the detection criterion from the representation. The same experimental scheme as in the previous work was followed to show the improvement of the detection results. A training experiment was first realized to set the optimal detection parameters. Thereafter, we compared the results of the detection obtained through our automatic method to those of six human experts. A complementary experiment was performed to evaluate the behavior of the new proposed method with long recording durations.
Volume 25, Number 3, 1999
the Fourier transform. For a discrete time signal, the subsampled form of the spectrogram can be written as (Allen and Rabiner 1977; Portnoff 1981):
SP~n, f ! 5
U
O
m5N21
m52N11
Time-frequency representation (TFR) of the Doppler signal Spectrogram representations are often used for nonstationary signal analysis, such as speech processing, acoustic analysis or biomedical signal processing (Hlawatsch and Boudreaux-Bartels 1992). The classic definition of the spectrogram corresponds to the application of a specific short time-window on the signal around each time sample before computing the square module of
U
2
, (1)
where n represents the discrete time axis, f represents the frequency axis, Dn is the analysis stride (the time distance between successive applications of the window to the data) and h(m) is an analysis window of 2N 2 1 points. This expression can be implemented efficiently using a fast Fourier transform (FFT) algorithm. The FFT size (Nfft) can be chosen equal to the next power of two of the analysis window lengths using zero-padding technique (Oppenheim and Schafer 1989). Changing the duration of the analysis window can modify the time and frequency resolution of the spectrogram. As a result, Hanning windows with different durations have been used to determine the optimal window length that provides the best results. Features extraction From the TFR previously computed, two characteristic functions have been extracted: the root-mean of the local power wpectrum (RMPS) and the modal frequency (fmode). Because the signals were real, only the positive Nfft frequencies of the TFR were investigated ( 1 1 2 frequencies):
DETECTION METHOD As described by Roy et al. (1998), the passage of a particle under the ultrasonic probe provided a shrinkage of the Doppler spectrum within a dominant frequency, accompanied by an increase of the Doppler power spectrum. Thus, the information carried by the passage of the embolus was concentrated around this dominant frequency. The different quantities described below will define regions of interest in the time-frequency representation of the Doppler signal, as well as detection criterion.
h~m! z s~nDn 1 m! z e
2j2pfm
RMPS~n! 5
Î
2 z Nfft 1 2
O SP~n,f !
f5fs/ 2
f50
f mode~n! 5 argmax0#f#fs/ 2~SP~n,f !! , (2) where fs is the sampling frequency. Supplementary conditions were imposed to avoid errors in the modal frequency estimation, particularly when the signal-to–noise ratio (SNR) was low. Furthermore, the ambient noise was not a white noise, which gave rise to an important positive bias. To reduce the effects of the colored noise and the low SNR, the modal frequency was revalued ~f˜mode ! in two steps. First, the global statistics (histogram) of RMPS were calculated. The histogram of RMPS typically has a bimodal distribution (Fig. 1a, b). This distribution can be approximated with two Gaussian densities: N1(m1,s1) and N2(m2,s2). To estimate the distribution parameters, an expectation-maximize procedure (EM algorithm) was
Automated detection of HITS signals ● E. ROY et al.
351
This threshold corresponded to the smallest modal frequency value when the estimation was unbiased by the noise. Thus, the second stage consisted in revaluing fmode at each time nI where RMPS was inferior to TRMPS. To attenuate the noise effect in the modal frequency estimation, ˜fmode was sought in the interval [0, Tfmode]: ˜f mode~n j! 5 argmax0#f#Tf ~SP~n j, f !!. mode
(6)
Figure 2 illustrates the different steps of the procedure. Assuming HITS are narrow-band signals, the Doppler spectrum must be concentrated around ˜fmode with a smaller frequency spreading when an embolus occurs than for the remainder of the Doppler signal (Fig. 3). An edge line flim(n) was defined, above which no significant important embolus information was present in the timefrequency plane: f lim~n! 5 min~2 z ˜f mode~n!, fs/ 2!.
(7)
(The minimum value between fs/2 was chosen to avoid flim results superior to the Nyquist frequency). To characterize frequency spreading, the first and second time moments of the TFR (Cohen 1995) were computed with the frequency summation limit defined by flim(n): Fig. 1. ( a) Root-mean of the local power spectrum; (b) corresponding histogram. The results of the EM algorithm are represented by N1(m1,s1) and N2(m2,s2).
used (Redner and Walker 1984). The Gaussian distribution N1(m1,s1) can be viewed as the Doppler-signal-in– noise part. Then, a threshold was set to twice the RMPS mean value m1 (; 6 dB above the mean background signal level):
O f z SP~n, f !
f5flim~n!
f m~n! 5
f50
(8)
O SP~n, f !
f5flim~n!
f50
O ~ f 2 f ~n!! z SP~n, f !
f5flim~n!
2
T RMPS 5 2 z m1,
m
(3)
We have made the assumption that the frequency modal estimation was biased when the RMPS was lower to this threshold, and unbiased otherwise. Then at each instant ni where RMPS was superior to this threshold, the frequency modal estimation was considered unbiased by the noise. Also the revalued modal frequency ˜fmode~ni !. remained unchanged: ˜f mode~n i! 5 f mode~n i!,
(4)
and a second threshold was fixed as follows: T fmode 5 min~f˜mode~n i!!.
(5)
B 2~n! 5
f50
O SP~n, f !
f5flim~n!
,
(9)
f50
where fm(n) and B(n) provided information on the average frequency and the bandwidth under flim(n), respectively. To remove unwanted fluctuations, the secondorder moment was subjected to a sixth-order Butterworth low-pass filter (fc 5 15 Hz). Detection criterion The HITS detection can be formulated as the following binary hypothesis problem (Scharf 1991; Helstrom 1995):
352
Ultrasound in Medicine and Biology
Volume 25, Number 3, 1999
Fig. 2. Illustration of the corrective method for the determination of the modal frequency.
Automated detection of HITS signals ● E. ROY et al.
353
Then, the detection test consists in comparing RMPS(n) to the reference: H1
RMPS~n! D~n! 5 RMPS ref
. ,
h,
(14)
H0
where h is a threshold that will be experimentally determined. To avoid multiple detections for the same embolus, all consecutive detections occurring within 10 ms from the previous detection were considered to be related to a single episode. If two consecutive detections were separated by more than this particular interval, the second event was considered as a new one. EXPERIMENTS Fig. 3. Root-nean of local power spectrum (upper) and smooth bandwidth (bottom) of typical arterial Doppler signal. This figure illustrates the narrow-band hypothesis when a HITS occurs.
H 0 : s(n) 5 b(n)
(10)
H 1 : s(n) 5 e(n) 1 b(n),
(11)
where s(n) is the observed signal, b(n) is the additive noise, and e(n) is the signal to be detected. The hypothesis H0 represents the cases where no embolic signal is present, b(n) is the Doppler signal without HITS and is considered as a nonstationary noise. H1 is the situation in which an embolus occurs within the Doppler signal. As described previously, when an embolus is present in the insonated blood flow, RMPS(n) increases and B(n) decreases compared to the Doppler signal without HITS. Generally, for the Doppler signal without embolus, these two quantities (RMPS(n), B(n)) simultaneously reach their greatest values around the peak systole (Fig. 3). Thus, we take as reference the instant where B(n) is the highest and RMPS(n) is a local maximum: n ref 5 argmaxn~B~n!!.
(12)
Similarly, the RMPS reference is computed as follows: RMPS ref 5 RMPS~n ref!.
(13)
Doppler signals have been collected on the posterotibial artery using a continuous wave 4-MHz ultrasound Doppler device (Angiodop 481, DMS, France). Because the dynamic range of the continuous-wave ultrasound system could not be changed, particular attention was paid to avoid saturation of the Doppler signal by using a medium recording power. Signals were digitized at 16 bits with a sampling frequency of 11,025 Hz with a compatible sound blaster card (ESS1488 Audiodrive, ESS Technology, Inc.), using a portable computer (Pentium 75, Notebook Computer, Taiwan). Recordings were performed during clinical procedures. The results of this study were reported elsewhere (Abraham et al. 1997). The nature of the embolic signals are explained in Roy et al. (1998). From Doppler recordings, three databases were analyzed as follows. The analysis technique and the detection procedure were implemented in MATLABt (version 4.2c, The MathWorks Inc., Natick, MA) using the Signal Processing Toolboxt. To analyze the performance of the proposed method, several criteria were used. Positive predictive value (good detection) is the ratio of detected events to the number of events detected by human analyses. Falsepositive (false detection) is the ratio of the number of detected events that were not found by human experts to the total number of events found by the automated method. False-negative (no detection) is the ratio of human detected events that were not found by the automated method to the total number of human-detected events. Thereby, the sum of good detection, false-positive detection and no detection can be greater than 100%. First experiment The first database consisted of 560 HITS distributed within 190 files of 6 s each. These HITS were detected
354
Ultrasound in Medicine and Biology
Volume 25, Number 3, 1999
by three human experts and have been chosen to offer a wide variety of embolic appearances (e.g., sound, duration, isolated or as showers, Doppler spectrum). Each event noted by at least one of the three experts was considered as a HITS. The detection procedure was tested using this database with different observation window duration and with different values of the threshold h. Thus, the spectrogram was computed with distinct time-window duration of 10, 12, 14, 16, 18, 20, 22, 24 and 26 ms, with a 50% window overlap. The FFT size (Nfft) was matched to the next power of two of the analysis window lengths. For each time window duration, the detection procedure was applied with various threshold values: h ranging from 1 to 3.5 with an increment of 0.1. hopt was defined as the h value that was the closest to the threshold corresponding to an equal error rate (false-positive rate 5 false-negative rate). In the second experiment, the optimal duration window and the hopt threshold were chosen to compare results between the automated detection and human experts. Second experiment The second database came from 12 other patients. Because very long interprocedure periods occurred in the original recordings, an independent expert has focused the analysis periods in the areas of interest of these patients, resulting in 38 files of 15–30 s each. The experts were asked independently to analyze the 38 recordings and to record the time of each embolic signal they heard. The experts were blinded to the results of the others. The recordings were subdivided into periods of 0.5 s to localize the moments when HITS occurred. The study of the agreement between the six experts has been reported elsewhere (Roy et al. 1998). The detection ability of the automated technique was compared with the six other detections obtained by the human experts. For different human “gold standards” (number of experts in agreement), Kappa values (Roy et al. 1998) obtained by the proposed method (M2) were compared with those obtained by the previously reported method (M1). Third experiment A preliminary experiment was performed to test the method in the case of long recording duration (clinical condition). Fifteen patients (20 –30 min each) led to a 6-h 24-min recording that has been analyzed. Signals have been subdivided into periods of 1 s to localize the moments at which embolic signals occurred (as suggested by Van Zuilen et al. 1996). Each event noted by all of the three experts was considered as a HITS. The automated method was then applied to these recordings with the
Fig. 4. ROC curves from the first experiment.
window duration and the h value determined as in the first experiment. RESULTS First experiment The results of the first part of this experiment are shown in Fig. 4a– c. These figures represent the evolution of the good detection as a function of the false-positive detection for different values of the threshold h. Time
Automated detection of HITS signals ● E. ROY et al.
355
the automated method (Fig. 6) and each observer yielded a percentage of true-positives ranging from 80.1 to 90.9%, and a percentage of false-positives ranging from 3.4 to 15.1%. On average, 86.2% true-positives were obtained (a percentage of false-negatives equal to 13.8%), corresponding to 7.9% of false-positives. Table 2 reported the concordance between the methods (M1:M2) and the different human “gold standard” (HGS). For M1, the Kappa values (and the standard error) ranged from 0.79 (0.03) to 0.88 (0.02). For the proposed method, the Kappa values ranged from 0.82 (0.02) to 0.92 (0.02). As presented in Fig. 7, 93.9% of 115 HITS heard by all observers were detected by the proposed automated method (M2). This proportion decreases as the number of observers in agreement decreases. If the cut-off level is set on significant HITS heard by at least one, two, three, four or five experts, M2 obtains 73.1, 84.7, 88.7, 89.6 and 91.7% of good detections with 1.4, 1.4, 2.7, 6.2 and 9.6 % false-negative detections, respectively. No experts heard two events of the 146 detected automatically (1.37%).
Fig. 5. Graphic determination of the optimal time window duration.
window durations of 22 ms, 20 ms and 18 ms give similar ROC curves, and short time durations obtain more distinct results. Figure 5 represents the evolution of the normalized area of these nine ROC curves as a function of the time window duration. The window size that maximizes this area is that size chosen for the next experiments. For each time window duration, we have reported the threshold that is the closest to the threshold that realizes an equal error rate. Percentages of false and good detections that correspond to this threshold have been summarized in Table 1. Estimated values of this threshold are close to each other whatever the window length, but they give equal error rates that increase when the window size decreases. Thus, a 22-ms Hanning window duration and a threshold of 2.1 were used as standard setting for the two other experiments.
Third experiment As reported in Table 3, for these 6-h 24-min recordings, the automated method found 54 HITS (59 detected by the human experts) and 31 signals in excess. Among these 31 extra detections, 4 are “true” false-positives, and 19 events correspond to artefact. And, after a second human analysis of these 31 events, 8 signals were reported as possible HITS. An overall sensibility of 91.5% was achieved. DISCUSSION
Second experiment On the recordings of the second database, the six experts reported a total of 197 different events. Observers Number 1 to 6 reported 140, 176, 155, 161, 161 and 146 events, respectively; 115 events were found by all observers, 29 by five observers, nine by four observers, seven by three observers, ten by two observers and 27 events by one observer. The automated method detected 146 events. Comparison between detections realized by
Exploiting the narrow-band hypothesis (Roy et al. 1998), the proposed method is based on the time-frequency processing of the arterial Doppler signal, and the previously reported method used temporal processing to perform the automated detection. The chosen time-frequency representation is the spectrogram that is used for classical Doppler signal analysis (Jensen 1996). Recently, several studies (Wang
Table 1. Estimated threshold values that are the closest to the threshold that realizes an equal error rate for each time window duration Time window duration
Threshold % False-positive % False-negative
26 ms
24 ms
22 ms
20 ms
18 ms
16 ms
14 ms
12 ms
10 ms
1.7 11.5 13.7
1.7 15.2 14.1
2.1 7.9 6.6
2.1 8 7.5
2.1 8.5 9.6
2.0 10.8 11.6
2.0 18.8 16.6
2.2 19.8 20.4
2.9 22.3 20.4
356
Ultrasound in Medicine and Biology
Volume 25, Number 3, 1999
Fig. 6. Comparisons of automated detection with human detection for each observer.
and Fish 1996; Cardoso et al. 1996; Guo et al. 1994) have been performed to determine what was the best time-frequency distribution for the analysis of the Doppler blood flow signal. Authors concluded that Bessel distribution (Guo et al. 1994) or Choi–Williams distribution (Wang and Fish 1996; Cardoso et al. 1996) provided higher performances than the spectrogram for the analysis of simulated Doppler ultrasound signals. How-
ever, the authors specified that their Doppler signals were less complex than those found in clinical practice and that their results were strictly valid only for the simulated signals used. Although the main shortcoming of the spectrogram is its tradeoff between time and frequency resolution, it has the advantage of providing nonnegative values, fast computation and easy interpretation. Furthermore, the presented detection scheme can be realized
Table 2. Results from the third experiment Automated detection False-positive
Tape
Duration
Number of human detection
True detection
Possible HITS
Artefact
“True” falsepositive
Sensibility
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 All
0:27:50 0:21:52 0:23:23 0:21:14 0:23:15 0:29:02 0:22:10 0:22:43 0:27:51 0:30:00 0:32:40 0:21:43 0:27:21 0:24:29 0:28:17 6:23:50
4 8 4 14 4 – – – 18 – – – – 7 – 59
4 8 4 12 4 – – – 16 – – – – 6 – 54
5 2 – – – – – – 1 – – – – – – 8
5 – 3 1 2 – – – – 3 2 – – 2 1 19
– – – – – – 3 – 1 – – – – – – 4
100 100 100 85.7 100 – – – 88.9 – – – – 85.7 – 91.5
Automated false-positives were subdivided into three categories.
Automated detection of HITS signals ● E. ROY et al.
Fig. 7. Percentage of automated false-positive detection and percentage of automated true detection as a function of human “gold standard” (M1 5 previously reported method; M2 5 proposed method).
with any TFR of the Cohen class (Cohen 1995). The only condition is to take the positive or absolute values of the TFR. Therefore, further studies should be realized to compare the ability of each Cohen class member to detect HITS. However, Choi–Williams, Bessel or smoothed-pseudoWigner–Ville distributions can be powerful tools for a fine characterization of already detected HITS. Their high time and frequency resolutions could allow a more precise measurement (duration, intensity, chirp-rate,. . .) but would require more timeconsuming computation for the detection than the spectrogram (Smith et al. 1994). Because it is assumed that emboli are high-intensity narrow-band signals, the root-mean of the local power spectrum and the spectral bandwidth allow for the extraction of a high-intensity wide band signal from the nonstationary background blood flow. The main test used for the automatic detection is the comparison between the signal spectral contents to this reference. That is close
357
to the general method, which used embolic-to– blood ratio (EBR) to size and characterize embolic events (Moehring and Klepper 1994; Moehring and Ritcey 1996; Moehring et al. 1996). However, our approach endeavors to look for embolic-free regions of the signal with no human intervention and with no a priori knowledge of the Doppler signal. Nevertheless, the peripheral blood flow is quite complex: distribution of blood velocities in human vessels depends on different conditions, such as arterial and venous sites, vessel dimensions or concomitant pathology, thermal or exercise stress. Furthermore, the unknown noise level and the colored nature of the noise can engender large bias in the spectral bandwidth estimation. Because flim defines the summation limit in bandwidth computation, the determination of ˜fmode and flim is an important step in the detection scheme presented in this study. This bandwidth estimation allows the characterization of the energy concentration around a dominant frequency, rather than for the exact Doppler bandwidth. A more sophisticated method could use maximum Doppler frequency estimators (the percentile method, the threshold crossing methods, the hybrid method or the geometrical method) (Mo et al. 1988; Marasek and Nowicki 1994) to define the summation limit. But, it was shown that their behavior strongly depends on SNR, noise level estimation and the shape of spectra (Marasek and Nowicki 1994). Nevertheless, the summation limit defined by flim does not truncate the information carried by embolic signals. A particularly interesting part of the first experiment is that our result differs from the window length selection usually reported. Multiple reports have proved that the spectrogram shows an optimum window length at 8 ms (Guo et al. 1994) or 10 ms (Wang and Fish 1996; Cardoso et al 1996) that is in agreement with conventional use. In our study, 18, 20 or 22 ms provide better results than smaller window durations, although ana-
Table 3. Percentages of true-positives (tp), percentages of false-positives (fp) and Kappa values (K) obtained by the previously reported method (M1) and the presented method (M2) for different human “gold standard” Methods M1
M2
Human “gold standard”
% tp
% fp
K
% tp
% fp
K
At least 1 expert agreed At least 2 experts agreed At least 3 experts agreed At least 4 experts agreed At least 5 experts agreed Exactly 6 experts agreed Mean
72.08 83.3 86.88 88.89 90.97 94.78 86.15
6.58 6.58 8.55 10.53 13.82 28.29 12.39
0.79 (0.03) 0.86 (0.02) 0.88 (0.02) 0.88 (0.02) 0.87 (0.02) 0.80 (0.03) 0.85 (0.02)
73.1 84.71 88.75 89.56 91.67 93.91 86.95
1.37 1.37 2.74 6.16 9.59 26.03 7.88
0.82 (0.02) 0.90 (0.02) 0.92 (0.02) 0.91 (0.02) 0.90 (0.02) 0.81 (0.03) 0.88 (0.02)
Values in parentheses are the standard error, for each Kappa value p , 0.00001.
358
Ultrasound in Medicine and Biology
lyzed Doppler signals can no longer be presumed stationary (5–20 ms; Jensen 1996). This may, in part, be ascribed to the summation limit described above. Furthermore, consistent with the narrow band hypothesis (Roy et al. 1998), the proposed method needs to provide a high-frequency resolution to differentiate HITS bandwidth from the background blood flow bandwidth. Thus, the best frequency resolution is obtained with a long duration time window. The choice of the method for the optimal threshold determination is arbitrary, but corresponds to a minimax rule ( Scharf 1991), where the cost assigned to the probability of a false alarm matches the cost assigned to the probability of a miss. For our study, this choice is suggested as being the best compromise for clinical use. Although 18-, 20- and 22-ms time windows show close results, the 22-ms duration has been chosen as the optimal value. In this case, the chosen threshold provides 6.6% false-negatives for 7.9% false-positives when the previously reported method allowed an equal error rate about 14.5%. In the second experiment, detected events by the automated technique have been compared to the events detected by the six human experts. The previous work (Roy et al. 1998) reported that the agreement between the six human experts was lower than the one mentioned by Van Zuilen et al. (1996). A possible explanation is that our second database shows fewer correct negative counts than those used by Van Zuilen and co-workers. This could have some influence on the agreement measurement. Another explanation is that our database contains quite different HITS intensity, duration and Doppler spectrum. Nevertheless, our proposed method obtains an overall mean rate of true-positives (tp) of 86.2%, for an overall mean rate of false-positives (fp) of 7.9%. Kappa values between the proposed method and the different human gold standard (overall mean: 0.88 [0.02]) are greater than Kappa values between the previously reported method and HGS (overall mean: 0.85 [0.02]). These results confirm the improvement for the agreement between human experts and the proposed method. Although the percentage of true-positives provided by M2 is superior to the percentage of true-positives provided by M1, the improvement is not statistically significant. Whereas the x2 test showed a significant difference in false-positives between M1 and M2 (p , 0.05) when the HGS have at least 1, 2 or 3 experts in agreement. The previously reported method provided an high percentage of true-positives; thererfore, it is difficult to improve this proportion significantly. As we specified in the introduction, our goal is rather to reduce significantly the number of false-positive detections from that of the previously reported method. The overall sensitivity (91.5%) achieved in the third experiment is consistent with that provided in the first
Volume 25, Number 3, 1999
experiment by the experimental determination of the optimal threshold (tp 5 92.1%, fp 5 6.6%). We believe the automatic detection of embolic signals should be performed as in 24-h electrocardiogram Holters. This requires a high sensitivity, and postprocessing analysis of the detected events should be performed by human experts to differentiate HITS from artefacts. Because an inverse relationship exists between sensitivity and specificity, and before a human postprocessing analysis is performed, automated detection should provide a high sensitivity, rather than high specificity. The fact that, in the postprocessing human analysis, some of the automatically detected events were probably to be considered as true HITS not detected by human experts in the 6-h 25-min recording is consistent with this approach. Furthermore, the human postanalysis is a time-consuming process, and the fact that the proposed method decreases the number of false detections is an interesting result. Thereby, extra detections have been divided into three groups: artefact, events not heard by the three experts (true false-positives) and missing HITS. In the present study, artefacts were not classified as “true” false-positives because our method does not include a specific methodology for artefact discrimination from HITS. Generally, the unidirectional property of embolic signals is considered as differentiation criterion (Consensus Committee NICHS 1995; Markus et al. 1993). But, some large gaseous emboli often overload the Doppler frontend, which leads to apparently bidirectional signals like artefacts. A more recent technique has been developed using dual-gated Doppler ultrasound (Smith et al. 1996) to differentiate emboli from artefact, but it requires the use of a pulsed Doppler device that cannot easily be adapted to the peripheral limb circulation. Thus, the artefact differentiation problem, with the continuous Doppler system, remains a new direction of research. CONCLUSION We have presented a new scheme for automated HITS detection. The narrow band hypothesis was exploited and characterized by time-frequency processing. The present work shows that this scheme improves the positive predictive value compared with the time processing recently proposed. Moreover, the detection scheme perdmits the use of the conventional spectrogram as for any Cohen class member. For real-time clinical routine use, the off-line property of this method could appear as the major drawback. Nevertheless, some modifications could be conceivable to perform online detection. When starting the clinical exam, the detection scheme could be used to extract initialization setups. As a result, the online detection could consist in comparing instantaneous computed in-
Automated detection of HITS signals ● E. ROY et al.
dex with the initial setup. Keeping in mind that clinical recording conditions may change during the investigation, initial settings could be reassessed after some min and/or after each automated detection. Acknowledgements—The authors gratefully acknowledge B. Vielle for his help in the statistical computation.
REFERENCES Abraham P, L’Hoste P, Pasco-Papon A, Roy E, Leftheriotis G, Saumet JL. Asymptomatic peripheral embolism is almost constant during lower limb diagnostic arteriography. J Vasc Invest 1997;3(4):193–196. Allen JB, Rabiner LR. A unified approach to short-time Fourier analysis and synthesis. Proc IEEE 1977;65(11):1558 –1564. Cardoso JC, Ruano MG, Fish PJ. Nonstationary broadening reduction in pulsed Doppler spectrum measurements using time-frequency estimators. IEEE Trans Biomed Eng 1996;43(12):1176 –1186. Cohen L. Time-frequency analysis. Englewood Cliffs, NJ: PrenticeHall, 1995. Consensus Committee, NICHS (Ninth International Cerebral Hemodynamic Symposium). Basic identification criteria of Doppler microembolic signals. Stroke 1995;26(6):1123–. Guo Z, Durand LG, Lee HC. Comparison of time-frequency distributions techniques for analysis of simulated Doppler ultrasound signals of the femoral artery. IEEE Trans Biomed Eng 1994;41(4): 332–342. Helstrom C. Elements of signal detection and estimation. Englewood Cliffs, NJ: Prentice-Hall, 1995. Hlawatsch F, Boudreaux-Bartels GF. Linear and quadratic time-frequency signal representations. IEEE Signal Processing Mag 1992; 9:21– 67. Jensen JA. Estimation of blood velocities using ultrasound: a signal processing approach. Cambridge, UK: Cambridge University Press, 1996. Marasek K, Nowicky A. Comparison of the performance of three maximum Doppler frequency estimators coupled with different spectral estimation methods. Ultrasound Med Biol 1994;20(7): 629 – 638. Markus H, Loh A, Brown MM. Computerized detection of cerebral emboli and discrimination from artifact using Doppler ultrasound. Stroke 1993;24:1667–1672. Mo L, Yun L, Cobbold R. Comparison of four digital maximum frequency estimators for Doppler ultrasound. Ultrasound Med Biol 1988;14(5):355–363. Moehring M, Klepper J. Pulse Doppler ultrasound detection, characterization and size estimation of emboli in flowing blood. IEEE Trans Biomed Eng 1994;41(1):35– 44. Moehring M, Ritcey J. Sizing emboli in blood using pulse Doppler ultrasound-1: verification of the EBR model. IEEE Trans Biomed Eng 1996;43(6):572–580.
359
Moehring M, Ritcey J, Ishimaru A. Sizing emboli in blood using pulse Doppler ultrasound-2: effect of beam refraction. IEEE Trans Biomed Eng 1996;43(6):581–588. Nicholls SC, O’Brien JK, Sutton MG. Venous thromboembolism: Detection by duplex scanning. J Vasc Surg 1996;:511–516. Oppenheim AV, Schafer RW. Discrete-time signal processing. Englewood Cliffs, NJ: Prentice-Hall, 1989. Portnoff MR. Short-time Fourier analysis of sampled speech. IEEE Trans Acoust Speech Sig Proc 1981;29(3):364 –373. Redner R, Walker H. Mixture densities, maximum likelihood, and the EM algorithm. SIAM Rev 1984;26:195–239. Roy E, Abraham P, Montre´sor S, Baudry M, Saumet JL. The narrow band hypothesis: an interesting approach for high-intensity transient signals (HITS) detection. Ultrasound Med Biol 1998;24(3): 375–382. Scharf L. Statistical signal processing: detection, estimation and time series analysis. Reading, MA: Addison-Wesley, 1991. Siebler M, Rose G, Sitzer M, Bender A, Steinmetz H. Real-time identification of cerebral microemboli with US feature detection by a neural network. Radiology 1994;192(3):739 –742. Siebler M, Sitzer M, Steinmetz H. Detection of intracranial emboli in patients with symptomatic extracranial carotid artery disease. Stroke 1992;23:1652–1654. Smith JL, Evans DH, Fan L, Bell PRF, Naylor AR. Differentiation between emboli and artefacts using dual gated transcranial Doppler ultrasound. Ultrasound Med Biol 1996;22:1031–1036. Smith JL, Evans DH, Fan L, Thrush AJ, Naylor AR. Processing Doppler ultrasound signals from blood-borne emboli. Ultrasound Med Biol 1994;20:455– 462. Spencer MP. Detection of embolism with Doppler ultrasound: a review. Echocardiography 1996;13(5):519 –527. Spencer MP. Transcranial Doppler monitoring and causes of stroke from carotid endarterectomy. Stroke 1997;28:685– 691. Spencer MP, Campbell SD. Development of bubbles in venous and arterial blood during hyperbaric decompression. Bull Mason Clin 1968;22:26 –32. Spencer MP, Campbell SD, Sealey JL, Henry FCH, Lindbergh J. Experiments on decompression bubbles in the circulation using ultrasonic and electromagnetic flowmeters. J Occup Med 1969; 11(5):238 –245. Tegeler C. High-intensity transient signals detected by Doppler ultrasonography: searching for answers. Cerebrovasc Dis 1994;4:379 – 382. Van Zuilen EV, Mess HW, Jansen C, Van Der Tweel I, Van Gijn J, Ackerstaff RGA. Automatic embolus detection compared with human experts. Stroke 1996;27:1840 –1843. Wang Y, Fish PJ. Comparison of Doppler signal analysis techniques for velocity waveform, turbulence and vortex measurement: a simulation study. Ultrasound Med Biol 1996;22(5):635– 649.
PII S0301-5629(98)00173-2
● Original Contribution SPECTROGRAM ANALYSIS OF ARTERIAL DOPPLER SIGNALS FOR OFF-LINE AUTOMATED HITS DETECTION EMMANUEL ROY,*† SILVIO MONTRE´ SOR,* PIERRE ABRAHAM† and JEAN-LOUIS SAUMET† *Laboratoire d’Informatique de l’Universite´ du Maine, rue Lae¨nec, BP 535, 72085 Le Mans Cedex 9 France; and † Laboratoire de Physiologie et d’Explorations Fonctionnelles Vasculaires, Centre Hospitalier Universitaire, 49033 Angers Cedex 01 France (Received 7 July 1998; in final form 27 October 1998)
Abstract—Recently, a time processing of arterial Doppler signals was proposed to detect automatically highintensity transient signals (HITS). This technique provided satisfactory detection results, but was not always constantly accurate, particularly with high-resistance blood velocity profiles. A time-frequency processing, based on the spectrogram, is presented to detect the presence of emboli in the arterial Doppler signals. The method uses the narrow-band hypothesis and extracts the detection criterion from the time-frequency representation (TFR). A first database of 560 peripheral arterial Doppler HITS was created to study microemboli and to define the normal limits to be used in our method. A threshold was experimentally defined using this database, and then applied to 38 recordings from 12 patients. Using another database, 6 human expert Doppler users reported 140, 176, 155, 161, 161 and 146 HITS, corresponding to a total of 197 different observed HITS. When an event was detected by at least 6, 5, 4, 3, 2 and 1 of the observers, sensitivity of the automatic detection was 93.9, 91.7, 89.6, 88.7, 84.7 and 73.1%, respectively. The sensitivity of our automatic detection is, thus, highly associated with the number of observers in agreement. A preliminary experiment has been performed to test the method in the case of long recording duration. In 15 patients, 6 h 24 min of recordings have been analyzed. The proposed automated processing provided an overall sensibility of 91.5%. The present work shows that this detection scheme preserves good sensibility and improves the positive predictive value compared with the time-processing recently proposed. © 1999 World Federation for Ultrasound in Medicine & Biology. Key Words: Spectrogram, Doppler ultrasound, Embolism, Detection.
depending on their size, characteristics, embolus-to– blood ratio and the ultrasonic carrier frequency. Within the appropriate dynamic range of bidirectional Doppler equipment, an embolic signal is unidirectional within the Doppler velocity spectrum. HITS are always audible, and their sound depends on their velocity. Arterial embolisms can cause serious damage in the embolized organ in both cerebral (Spencer 1997; Tegeler 1994; Siebler et al. 1992) and peripheral (Abraham et al. 1997) circulation. Arterial emboli may involve the extremities, especially the legs and feet. Less common sites include the kidneys, intestines and other areas. Thus, the continuous monitoring of patients with significant risk of thromboembolic complications is of major interest. Several commercial software systems have been developed to detect automatically cerebral embolic signals with pulsed Doppler devices (Siebler et al. 1994; Van Zuilen et al. 1996). Automated techniques are needed for detecting peripheral emboli with continuous Doppler systems. Recently, Roy et al. (1998) proposed a
INTRODUCTION The possibility of detecting microemboli using Doppler ultrasound was reported in the 1960s by Spencer and co-workers during hyperbaric decompression experiments (Spencer and Campbell 1968; Spencer et al. 1969). The passage of a solid embolus through the beam of a Doppler ultrasound instrument appears in the Doppler spectrum as a high intensity transient signal (HITS). Embolic signals, when recorded optimally, have been defined based on four main features (Consensus Committee NICHS 1995; Nicholls et al. 1996; Spencer 1996). They are transient (usually , 300 ms) and their duration depends on their time of passage through the Doppler sample volume. Their amplitude is usually at least 3 dB higher than that of the background blood-flow signal, Address correspondence to: Professor J. L. Saumet, Centre Hospitalier Universitaire, Laboratoire de Physiologie et d’Explorations Fonctionnelles Vasculaires, 49033 Angers Cedex 01, France. E-mail: [email protected] 349
350
Ultrasound in Medicine and Biology
time-processing approach to perform the detection. They assumed that the Doppler signal could be considered as a narrow-band signal when an embolus occurs. The detection procedure consisted in retaining only instants where both the variation of the signal envelope was greater and the variation of the instantaneous frequency estimation was lower than respective variations of the background Doppler signal. This preliminary technique provided satisfactory detection results compared with detections made by six human experts. Nevertheless, the velocity profiles in the peripheral circulation vary widely compared with the those in the cerebral circulation, depending on the artery recorded, from one subject to another, or in the same subject over time. As mentioned by the authors (Roy et al. 1998), the originally reported method was not always stable, particularly with high-resistance blood-velocity profiles, which entailed false-positives. In this work, time-frequency processing is presented to detect the presence of emboli in arterial Doppler signals. The method uses the narrow-band hypothesis previously reported, but extracts the detection criterion from the representation. The same experimental scheme as in the previous work was followed to show the improvement of the detection results. A training experiment was first realized to set the optimal detection parameters. Thereafter, we compared the results of the detection obtained through our automatic method to those of six human experts. A complementary experiment was performed to evaluate the behavior of the new proposed method with long recording durations.
Volume 25, Number 3, 1999
the Fourier transform. For a discrete time signal, the subsampled form of the spectrogram can be written as (Allen and Rabiner 1977; Portnoff 1981):
SP~n, f ! 5
U
O
m5N21
m52N11
Time-frequency representation (TFR) of the Doppler signal Spectrogram representations are often used for nonstationary signal analysis, such as speech processing, acoustic analysis or biomedical signal processing (Hlawatsch and Boudreaux-Bartels 1992). The classic definition of the spectrogram corresponds to the application of a specific short time-window on the signal around each time sample before computing the square module of
U
2
, (1)
where n represents the discrete time axis, f represents the frequency axis, Dn is the analysis stride (the time distance between successive applications of the window to the data) and h(m) is an analysis window of 2N 2 1 points. This expression can be implemented efficiently using a fast Fourier transform (FFT) algorithm. The FFT size (Nfft) can be chosen equal to the next power of two of the analysis window lengths using zero-padding technique (Oppenheim and Schafer 1989). Changing the duration of the analysis window can modify the time and frequency resolution of the spectrogram. As a result, Hanning windows with different durations have been used to determine the optimal window length that provides the best results. Features extraction From the TFR previously computed, two characteristic functions have been extracted: the root-mean of the local power wpectrum (RMPS) and the modal frequency (fmode). Because the signals were real, only the positive Nfft frequencies of the TFR were investigated ( 1 1 2 frequencies):
DETECTION METHOD As described by Roy et al. (1998), the passage of a particle under the ultrasonic probe provided a shrinkage of the Doppler spectrum within a dominant frequency, accompanied by an increase of the Doppler power spectrum. Thus, the information carried by the passage of the embolus was concentrated around this dominant frequency. The different quantities described below will define regions of interest in the time-frequency representation of the Doppler signal, as well as detection criterion.
h~m! z s~nDn 1 m! z e
2j2pfm
RMPS~n! 5
Î
2 z Nfft 1 2
O SP~n,f !
f5fs/ 2
f50
f mode~n! 5 argmax0#f#fs/ 2~SP~n,f !! , (2) where fs is the sampling frequency. Supplementary conditions were imposed to avoid errors in the modal frequency estimation, particularly when the signal-to–noise ratio (SNR) was low. Furthermore, the ambient noise was not a white noise, which gave rise to an important positive bias. To reduce the effects of the colored noise and the low SNR, the modal frequency was revalued ~f˜mode ! in two steps. First, the global statistics (histogram) of RMPS were calculated. The histogram of RMPS typically has a bimodal distribution (Fig. 1a, b). This distribution can be approximated with two Gaussian densities: N1(m1,s1) and N2(m2,s2). To estimate the distribution parameters, an expectation-maximize procedure (EM algorithm) was
Automated detection of HITS signals ● E. ROY et al.
351
This threshold corresponded to the smallest modal frequency value when the estimation was unbiased by the noise. Thus, the second stage consisted in revaluing fmode at each time nI where RMPS was inferior to TRMPS. To attenuate the noise effect in the modal frequency estimation, ˜fmode was sought in the interval [0, Tfmode]: ˜f mode~n j! 5 argmax0#f#Tf ~SP~n j, f !!. mode
(6)
Figure 2 illustrates the different steps of the procedure. Assuming HITS are narrow-band signals, the Doppler spectrum must be concentrated around ˜fmode with a smaller frequency spreading when an embolus occurs than for the remainder of the Doppler signal (Fig. 3). An edge line flim(n) was defined, above which no significant important embolus information was present in the timefrequency plane: f lim~n! 5 min~2 z ˜f mode~n!, fs/ 2!.
(7)
(The minimum value between fs/2 was chosen to avoid flim results superior to the Nyquist frequency). To characterize frequency spreading, the first and second time moments of the TFR (Cohen 1995) were computed with the frequency summation limit defined by flim(n): Fig. 1. ( a) Root-mean of the local power spectrum; (b) corresponding histogram. The results of the EM algorithm are represented by N1(m1,s1) and N2(m2,s2).
used (Redner and Walker 1984). The Gaussian distribution N1(m1,s1) can be viewed as the Doppler-signal-in– noise part. Then, a threshold was set to twice the RMPS mean value m1 (; 6 dB above the mean background signal level):
O f z SP~n, f !
f5flim~n!
f m~n! 5
f50
(8)
O SP~n, f !
f5flim~n!
f50
O ~ f 2 f ~n!! z SP~n, f !
f5flim~n!
2
T RMPS 5 2 z m1,
m
(3)
We have made the assumption that the frequency modal estimation was biased when the RMPS was lower to this threshold, and unbiased otherwise. Then at each instant ni where RMPS was superior to this threshold, the frequency modal estimation was considered unbiased by the noise. Also the revalued modal frequency ˜fmode~ni !. remained unchanged: ˜f mode~n i! 5 f mode~n i!,
(4)
and a second threshold was fixed as follows: T fmode 5 min~f˜mode~n i!!.
(5)
B 2~n! 5
f50
O SP~n, f !
f5flim~n!
,
(9)
f50
where fm(n) and B(n) provided information on the average frequency and the bandwidth under flim(n), respectively. To remove unwanted fluctuations, the secondorder moment was subjected to a sixth-order Butterworth low-pass filter (fc 5 15 Hz). Detection criterion The HITS detection can be formulated as the following binary hypothesis problem (Scharf 1991; Helstrom 1995):
352
Ultrasound in Medicine and Biology
Volume 25, Number 3, 1999
Fig. 2. Illustration of the corrective method for the determination of the modal frequency.
Automated detection of HITS signals ● E. ROY et al.
353
Then, the detection test consists in comparing RMPS(n) to the reference: H1
RMPS~n! D~n! 5 RMPS ref
. ,
h,
(14)
H0
where h is a threshold that will be experimentally determined. To avoid multiple detections for the same embolus, all consecutive detections occurring within 10 ms from the previous detection were considered to be related to a single episode. If two consecutive detections were separated by more than this particular interval, the second event was considered as a new one. EXPERIMENTS Fig. 3. Root-nean of local power spectrum (upper) and smooth bandwidth (bottom) of typical arterial Doppler signal. This figure illustrates the narrow-band hypothesis when a HITS occurs.
H 0 : s(n) 5 b(n)
(10)
H 1 : s(n) 5 e(n) 1 b(n),
(11)
where s(n) is the observed signal, b(n) is the additive noise, and e(n) is the signal to be detected. The hypothesis H0 represents the cases where no embolic signal is present, b(n) is the Doppler signal without HITS and is considered as a nonstationary noise. H1 is the situation in which an embolus occurs within the Doppler signal. As described previously, when an embolus is present in the insonated blood flow, RMPS(n) increases and B(n) decreases compared to the Doppler signal without HITS. Generally, for the Doppler signal without embolus, these two quantities (RMPS(n), B(n)) simultaneously reach their greatest values around the peak systole (Fig. 3). Thus, we take as reference the instant where B(n) is the highest and RMPS(n) is a local maximum: n ref 5 argmaxn~B~n!!.
(12)
Similarly, the RMPS reference is computed as follows: RMPS ref 5 RMPS~n ref!.
(13)
Doppler signals have been collected on the posterotibial artery using a continuous wave 4-MHz ultrasound Doppler device (Angiodop 481, DMS, France). Because the dynamic range of the continuous-wave ultrasound system could not be changed, particular attention was paid to avoid saturation of the Doppler signal by using a medium recording power. Signals were digitized at 16 bits with a sampling frequency of 11,025 Hz with a compatible sound blaster card (ESS1488 Audiodrive, ESS Technology, Inc.), using a portable computer (Pentium 75, Notebook Computer, Taiwan). Recordings were performed during clinical procedures. The results of this study were reported elsewhere (Abraham et al. 1997). The nature of the embolic signals are explained in Roy et al. (1998). From Doppler recordings, three databases were analyzed as follows. The analysis technique and the detection procedure were implemented in MATLABt (version 4.2c, The MathWorks Inc., Natick, MA) using the Signal Processing Toolboxt. To analyze the performance of the proposed method, several criteria were used. Positive predictive value (good detection) is the ratio of detected events to the number of events detected by human analyses. Falsepositive (false detection) is the ratio of the number of detected events that were not found by human experts to the total number of events found by the automated method. False-negative (no detection) is the ratio of human detected events that were not found by the automated method to the total number of human-detected events. Thereby, the sum of good detection, false-positive detection and no detection can be greater than 100%. First experiment The first database consisted of 560 HITS distributed within 190 files of 6 s each. These HITS were detected
354
Ultrasound in Medicine and Biology
Volume 25, Number 3, 1999
by three human experts and have been chosen to offer a wide variety of embolic appearances (e.g., sound, duration, isolated or as showers, Doppler spectrum). Each event noted by at least one of the three experts was considered as a HITS. The detection procedure was tested using this database with different observation window duration and with different values of the threshold h. Thus, the spectrogram was computed with distinct time-window duration of 10, 12, 14, 16, 18, 20, 22, 24 and 26 ms, with a 50% window overlap. The FFT size (Nfft) was matched to the next power of two of the analysis window lengths. For each time window duration, the detection procedure was applied with various threshold values: h ranging from 1 to 3.5 with an increment of 0.1. hopt was defined as the h value that was the closest to the threshold corresponding to an equal error rate (false-positive rate 5 false-negative rate). In the second experiment, the optimal duration window and the hopt threshold were chosen to compare results between the automated detection and human experts. Second experiment The second database came from 12 other patients. Because very long interprocedure periods occurred in the original recordings, an independent expert has focused the analysis periods in the areas of interest of these patients, resulting in 38 files of 15–30 s each. The experts were asked independently to analyze the 38 recordings and to record the time of each embolic signal they heard. The experts were blinded to the results of the others. The recordings were subdivided into periods of 0.5 s to localize the moments when HITS occurred. The study of the agreement between the six experts has been reported elsewhere (Roy et al. 1998). The detection ability of the automated technique was compared with the six other detections obtained by the human experts. For different human “gold standards” (number of experts in agreement), Kappa values (Roy et al. 1998) obtained by the proposed method (M2) were compared with those obtained by the previously reported method (M1). Third experiment A preliminary experiment was performed to test the method in the case of long recording duration (clinical condition). Fifteen patients (20 –30 min each) led to a 6-h 24-min recording that has been analyzed. Signals have been subdivided into periods of 1 s to localize the moments at which embolic signals occurred (as suggested by Van Zuilen et al. 1996). Each event noted by all of the three experts was considered as a HITS. The automated method was then applied to these recordings with the
Fig. 4. ROC curves from the first experiment.
window duration and the h value determined as in the first experiment. RESULTS First experiment The results of the first part of this experiment are shown in Fig. 4a– c. These figures represent the evolution of the good detection as a function of the false-positive detection for different values of the threshold h. Time
Automated detection of HITS signals ● E. ROY et al.
355
the automated method (Fig. 6) and each observer yielded a percentage of true-positives ranging from 80.1 to 90.9%, and a percentage of false-positives ranging from 3.4 to 15.1%. On average, 86.2% true-positives were obtained (a percentage of false-negatives equal to 13.8%), corresponding to 7.9% of false-positives. Table 2 reported the concordance between the methods (M1:M2) and the different human “gold standard” (HGS). For M1, the Kappa values (and the standard error) ranged from 0.79 (0.03) to 0.88 (0.02). For the proposed method, the Kappa values ranged from 0.82 (0.02) to 0.92 (0.02). As presented in Fig. 7, 93.9% of 115 HITS heard by all observers were detected by the proposed automated method (M2). This proportion decreases as the number of observers in agreement decreases. If the cut-off level is set on significant HITS heard by at least one, two, three, four or five experts, M2 obtains 73.1, 84.7, 88.7, 89.6 and 91.7% of good detections with 1.4, 1.4, 2.7, 6.2 and 9.6 % false-negative detections, respectively. No experts heard two events of the 146 detected automatically (1.37%).
Fig. 5. Graphic determination of the optimal time window duration.
window durations of 22 ms, 20 ms and 18 ms give similar ROC curves, and short time durations obtain more distinct results. Figure 5 represents the evolution of the normalized area of these nine ROC curves as a function of the time window duration. The window size that maximizes this area is that size chosen for the next experiments. For each time window duration, we have reported the threshold that is the closest to the threshold that realizes an equal error rate. Percentages of false and good detections that correspond to this threshold have been summarized in Table 1. Estimated values of this threshold are close to each other whatever the window length, but they give equal error rates that increase when the window size decreases. Thus, a 22-ms Hanning window duration and a threshold of 2.1 were used as standard setting for the two other experiments.
Third experiment As reported in Table 3, for these 6-h 24-min recordings, the automated method found 54 HITS (59 detected by the human experts) and 31 signals in excess. Among these 31 extra detections, 4 are “true” false-positives, and 19 events correspond to artefact. And, after a second human analysis of these 31 events, 8 signals were reported as possible HITS. An overall sensibility of 91.5% was achieved. DISCUSSION
Second experiment On the recordings of the second database, the six experts reported a total of 197 different events. Observers Number 1 to 6 reported 140, 176, 155, 161, 161 and 146 events, respectively; 115 events were found by all observers, 29 by five observers, nine by four observers, seven by three observers, ten by two observers and 27 events by one observer. The automated method detected 146 events. Comparison between detections realized by
Exploiting the narrow-band hypothesis (Roy et al. 1998), the proposed method is based on the time-frequency processing of the arterial Doppler signal, and the previously reported method used temporal processing to perform the automated detection. The chosen time-frequency representation is the spectrogram that is used for classical Doppler signal analysis (Jensen 1996). Recently, several studies (Wang
Table 1. Estimated threshold values that are the closest to the threshold that realizes an equal error rate for each time window duration Time window duration
Threshold % False-positive % False-negative
26 ms
24 ms
22 ms
20 ms
18 ms
16 ms
14 ms
12 ms
10 ms
1.7 11.5 13.7
1.7 15.2 14.1
2.1 7.9 6.6
2.1 8 7.5
2.1 8.5 9.6
2.0 10.8 11.6
2.0 18.8 16.6
2.2 19.8 20.4
2.9 22.3 20.4
356
Ultrasound in Medicine and Biology
Volume 25, Number 3, 1999
Fig. 6. Comparisons of automated detection with human detection for each observer.
and Fish 1996; Cardoso et al. 1996; Guo et al. 1994) have been performed to determine what was the best time-frequency distribution for the analysis of the Doppler blood flow signal. Authors concluded that Bessel distribution (Guo et al. 1994) or Choi–Williams distribution (Wang and Fish 1996; Cardoso et al. 1996) provided higher performances than the spectrogram for the analysis of simulated Doppler ultrasound signals. How-
ever, the authors specified that their Doppler signals were less complex than those found in clinical practice and that their results were strictly valid only for the simulated signals used. Although the main shortcoming of the spectrogram is its tradeoff between time and frequency resolution, it has the advantage of providing nonnegative values, fast computation and easy interpretation. Furthermore, the presented detection scheme can be realized
Table 2. Results from the third experiment Automated detection False-positive
Tape
Duration
Number of human detection
True detection
Possible HITS
Artefact
“True” falsepositive
Sensibility
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 All
0:27:50 0:21:52 0:23:23 0:21:14 0:23:15 0:29:02 0:22:10 0:22:43 0:27:51 0:30:00 0:32:40 0:21:43 0:27:21 0:24:29 0:28:17 6:23:50
4 8 4 14 4 – – – 18 – – – – 7 – 59
4 8 4 12 4 – – – 16 – – – – 6 – 54
5 2 – – – – – – 1 – – – – – – 8
5 – 3 1 2 – – – – 3 2 – – 2 1 19
– – – – – – 3 – 1 – – – – – – 4
100 100 100 85.7 100 – – – 88.9 – – – – 85.7 – 91.5
Automated false-positives were subdivided into three categories.
Automated detection of HITS signals ● E. ROY et al.
Fig. 7. Percentage of automated false-positive detection and percentage of automated true detection as a function of human “gold standard” (M1 5 previously reported method; M2 5 proposed method).
with any TFR of the Cohen class (Cohen 1995). The only condition is to take the positive or absolute values of the TFR. Therefore, further studies should be realized to compare the ability of each Cohen class member to detect HITS. However, Choi–Williams, Bessel or smoothed-pseudoWigner–Ville distributions can be powerful tools for a fine characterization of already detected HITS. Their high time and frequency resolutions could allow a more precise measurement (duration, intensity, chirp-rate,. . .) but would require more timeconsuming computation for the detection than the spectrogram (Smith et al. 1994). Because it is assumed that emboli are high-intensity narrow-band signals, the root-mean of the local power spectrum and the spectral bandwidth allow for the extraction of a high-intensity wide band signal from the nonstationary background blood flow. The main test used for the automatic detection is the comparison between the signal spectral contents to this reference. That is close
357
to the general method, which used embolic-to– blood ratio (EBR) to size and characterize embolic events (Moehring and Klepper 1994; Moehring and Ritcey 1996; Moehring et al. 1996). However, our approach endeavors to look for embolic-free regions of the signal with no human intervention and with no a priori knowledge of the Doppler signal. Nevertheless, the peripheral blood flow is quite complex: distribution of blood velocities in human vessels depends on different conditions, such as arterial and venous sites, vessel dimensions or concomitant pathology, thermal or exercise stress. Furthermore, the unknown noise level and the colored nature of the noise can engender large bias in the spectral bandwidth estimation. Because flim defines the summation limit in bandwidth computation, the determination of ˜fmode and flim is an important step in the detection scheme presented in this study. This bandwidth estimation allows the characterization of the energy concentration around a dominant frequency, rather than for the exact Doppler bandwidth. A more sophisticated method could use maximum Doppler frequency estimators (the percentile method, the threshold crossing methods, the hybrid method or the geometrical method) (Mo et al. 1988; Marasek and Nowicki 1994) to define the summation limit. But, it was shown that their behavior strongly depends on SNR, noise level estimation and the shape of spectra (Marasek and Nowicki 1994). Nevertheless, the summation limit defined by flim does not truncate the information carried by embolic signals. A particularly interesting part of the first experiment is that our result differs from the window length selection usually reported. Multiple reports have proved that the spectrogram shows an optimum window length at 8 ms (Guo et al. 1994) or 10 ms (Wang and Fish 1996; Cardoso et al 1996) that is in agreement with conventional use. In our study, 18, 20 or 22 ms provide better results than smaller window durations, although ana-
Table 3. Percentages of true-positives (tp), percentages of false-positives (fp) and Kappa values (K) obtained by the previously reported method (M1) and the presented method (M2) for different human “gold standard” Methods M1
M2
Human “gold standard”
% tp
% fp
K
% tp
% fp
K
At least 1 expert agreed At least 2 experts agreed At least 3 experts agreed At least 4 experts agreed At least 5 experts agreed Exactly 6 experts agreed Mean
72.08 83.3 86.88 88.89 90.97 94.78 86.15
6.58 6.58 8.55 10.53 13.82 28.29 12.39
0.79 (0.03) 0.86 (0.02) 0.88 (0.02) 0.88 (0.02) 0.87 (0.02) 0.80 (0.03) 0.85 (0.02)
73.1 84.71 88.75 89.56 91.67 93.91 86.95
1.37 1.37 2.74 6.16 9.59 26.03 7.88
0.82 (0.02) 0.90 (0.02) 0.92 (0.02) 0.91 (0.02) 0.90 (0.02) 0.81 (0.03) 0.88 (0.02)
Values in parentheses are the standard error, for each Kappa value p , 0.00001.
358
Ultrasound in Medicine and Biology
lyzed Doppler signals can no longer be presumed stationary (5–20 ms; Jensen 1996). This may, in part, be ascribed to the summation limit described above. Furthermore, consistent with the narrow band hypothesis (Roy et al. 1998), the proposed method needs to provide a high-frequency resolution to differentiate HITS bandwidth from the background blood flow bandwidth. Thus, the best frequency resolution is obtained with a long duration time window. The choice of the method for the optimal threshold determination is arbitrary, but corresponds to a minimax rule ( Scharf 1991), where the cost assigned to the probability of a false alarm matches the cost assigned to the probability of a miss. For our study, this choice is suggested as being the best compromise for clinical use. Although 18-, 20- and 22-ms time windows show close results, the 22-ms duration has been chosen as the optimal value. In this case, the chosen threshold provides 6.6% false-negatives for 7.9% false-positives when the previously reported method allowed an equal error rate about 14.5%. In the second experiment, detected events by the automated technique have been compared to the events detected by the six human experts. The previous work (Roy et al. 1998) reported that the agreement between the six human experts was lower than the one mentioned by Van Zuilen et al. (1996). A possible explanation is that our second database shows fewer correct negative counts than those used by Van Zuilen and co-workers. This could have some influence on the agreement measurement. Another explanation is that our database contains quite different HITS intensity, duration and Doppler spectrum. Nevertheless, our proposed method obtains an overall mean rate of true-positives (tp) of 86.2%, for an overall mean rate of false-positives (fp) of 7.9%. Kappa values between the proposed method and the different human gold standard (overall mean: 0.88 [0.02]) are greater than Kappa values between the previously reported method and HGS (overall mean: 0.85 [0.02]). These results confirm the improvement for the agreement between human experts and the proposed method. Although the percentage of true-positives provided by M2 is superior to the percentage of true-positives provided by M1, the improvement is not statistically significant. Whereas the x2 test showed a significant difference in false-positives between M1 and M2 (p , 0.05) when the HGS have at least 1, 2 or 3 experts in agreement. The previously reported method provided an high percentage of true-positives; thererfore, it is difficult to improve this proportion significantly. As we specified in the introduction, our goal is rather to reduce significantly the number of false-positive detections from that of the previously reported method. The overall sensitivity (91.5%) achieved in the third experiment is consistent with that provided in the first
Volume 25, Number 3, 1999
experiment by the experimental determination of the optimal threshold (tp 5 92.1%, fp 5 6.6%). We believe the automatic detection of embolic signals should be performed as in 24-h electrocardiogram Holters. This requires a high sensitivity, and postprocessing analysis of the detected events should be performed by human experts to differentiate HITS from artefacts. Because an inverse relationship exists between sensitivity and specificity, and before a human postprocessing analysis is performed, automated detection should provide a high sensitivity, rather than high specificity. The fact that, in the postprocessing human analysis, some of the automatically detected events were probably to be considered as true HITS not detected by human experts in the 6-h 25-min recording is consistent with this approach. Furthermore, the human postanalysis is a time-consuming process, and the fact that the proposed method decreases the number of false detections is an interesting result. Thereby, extra detections have been divided into three groups: artefact, events not heard by the three experts (true false-positives) and missing HITS. In the present study, artefacts were not classified as “true” false-positives because our method does not include a specific methodology for artefact discrimination from HITS. Generally, the unidirectional property of embolic signals is considered as differentiation criterion (Consensus Committee NICHS 1995; Markus et al. 1993). But, some large gaseous emboli often overload the Doppler frontend, which leads to apparently bidirectional signals like artefacts. A more recent technique has been developed using dual-gated Doppler ultrasound (Smith et al. 1996) to differentiate emboli from artefact, but it requires the use of a pulsed Doppler device that cannot easily be adapted to the peripheral limb circulation. Thus, the artefact differentiation problem, with the continuous Doppler system, remains a new direction of research. CONCLUSION We have presented a new scheme for automated HITS detection. The narrow band hypothesis was exploited and characterized by time-frequency processing. The present work shows that this scheme improves the positive predictive value compared with the time processing recently proposed. Moreover, the detection scheme perdmits the use of the conventional spectrogram as for any Cohen class member. For real-time clinical routine use, the off-line property of this method could appear as the major drawback. Nevertheless, some modifications could be conceivable to perform online detection. When starting the clinical exam, the detection scheme could be used to extract initialization setups. As a result, the online detection could consist in comparing instantaneous computed in-
Automated detection of HITS signals ● E. ROY et al.
dex with the initial setup. Keeping in mind that clinical recording conditions may change during the investigation, initial settings could be reassessed after some min and/or after each automated detection. Acknowledgements—The authors gratefully acknowledge B. Vielle for his help in the statistical computation.
REFERENCES Abraham P, L’Hoste P, Pasco-Papon A, Roy E, Leftheriotis G, Saumet JL. Asymptomatic peripheral embolism is almost constant during lower limb diagnostic arteriography. J Vasc Invest 1997;3(4):193–196. Allen JB, Rabiner LR. A unified approach to short-time Fourier analysis and synthesis. Proc IEEE 1977;65(11):1558 –1564. Cardoso JC, Ruano MG, Fish PJ. Nonstationary broadening reduction in pulsed Doppler spectrum measurements using time-frequency estimators. IEEE Trans Biomed Eng 1996;43(12):1176 –1186. Cohen L. Time-frequency analysis. Englewood Cliffs, NJ: PrenticeHall, 1995. Consensus Committee, NICHS (Ninth International Cerebral Hemodynamic Symposium). Basic identification criteria of Doppler microembolic signals. Stroke 1995;26(6):1123–. Guo Z, Durand LG, Lee HC. Comparison of time-frequency distributions techniques for analysis of simulated Doppler ultrasound signals of the femoral artery. IEEE Trans Biomed Eng 1994;41(4): 332–342. Helstrom C. Elements of signal detection and estimation. Englewood Cliffs, NJ: Prentice-Hall, 1995. Hlawatsch F, Boudreaux-Bartels GF. Linear and quadratic time-frequency signal representations. IEEE Signal Processing Mag 1992; 9:21– 67. Jensen JA. Estimation of blood velocities using ultrasound: a signal processing approach. Cambridge, UK: Cambridge University Press, 1996. Marasek K, Nowicky A. Comparison of the performance of three maximum Doppler frequency estimators coupled with different spectral estimation methods. Ultrasound Med Biol 1994;20(7): 629 – 638. Markus H, Loh A, Brown MM. Computerized detection of cerebral emboli and discrimination from artifact using Doppler ultrasound. Stroke 1993;24:1667–1672. Mo L, Yun L, Cobbold R. Comparison of four digital maximum frequency estimators for Doppler ultrasound. Ultrasound Med Biol 1988;14(5):355–363. Moehring M, Klepper J. Pulse Doppler ultrasound detection, characterization and size estimation of emboli in flowing blood. IEEE Trans Biomed Eng 1994;41(1):35– 44. Moehring M, Ritcey J. Sizing emboli in blood using pulse Doppler ultrasound-1: verification of the EBR model. IEEE Trans Biomed Eng 1996;43(6):572–580.
359
Moehring M, Ritcey J, Ishimaru A. Sizing emboli in blood using pulse Doppler ultrasound-2: effect of beam refraction. IEEE Trans Biomed Eng 1996;43(6):581–588. Nicholls SC, O’Brien JK, Sutton MG. Venous thromboembolism: Detection by duplex scanning. J Vasc Surg 1996;:511–516. Oppenheim AV, Schafer RW. Discrete-time signal processing. Englewood Cliffs, NJ: Prentice-Hall, 1989. Portnoff MR. Short-time Fourier analysis of sampled speech. IEEE Trans Acoust Speech Sig Proc 1981;29(3):364 –373. Redner R, Walker H. Mixture densities, maximum likelihood, and the EM algorithm. SIAM Rev 1984;26:195–239. Roy E, Abraham P, Montre´sor S, Baudry M, Saumet JL. The narrow band hypothesis: an interesting approach for high-intensity transient signals (HITS) detection. Ultrasound Med Biol 1998;24(3): 375–382. Scharf L. Statistical signal processing: detection, estimation and time series analysis. Reading, MA: Addison-Wesley, 1991. Siebler M, Rose G, Sitzer M, Bender A, Steinmetz H. Real-time identification of cerebral microemboli with US feature detection by a neural network. Radiology 1994;192(3):739 –742. Siebler M, Sitzer M, Steinmetz H. Detection of intracranial emboli in patients with symptomatic extracranial carotid artery disease. Stroke 1992;23:1652–1654. Smith JL, Evans DH, Fan L, Bell PRF, Naylor AR. Differentiation between emboli and artefacts using dual gated transcranial Doppler ultrasound. Ultrasound Med Biol 1996;22:1031–1036. Smith JL, Evans DH, Fan L, Thrush AJ, Naylor AR. Processing Doppler ultrasound signals from blood-borne emboli. Ultrasound Med Biol 1994;20:455– 462. Spencer MP. Detection of embolism with Doppler ultrasound: a review. Echocardiography 1996;13(5):519 –527. Spencer MP. Transcranial Doppler monitoring and causes of stroke from carotid endarterectomy. Stroke 1997;28:685– 691. Spencer MP, Campbell SD. Development of bubbles in venous and arterial blood during hyperbaric decompression. Bull Mason Clin 1968;22:26 –32. Spencer MP, Campbell SD, Sealey JL, Henry FCH, Lindbergh J. Experiments on decompression bubbles in the circulation using ultrasonic and electromagnetic flowmeters. J Occup Med 1969; 11(5):238 –245. Tegeler C. High-intensity transient signals detected by Doppler ultrasonography: searching for answers. Cerebrovasc Dis 1994;4:379 – 382. Van Zuilen EV, Mess HW, Jansen C, Van Der Tweel I, Van Gijn J, Ackerstaff RGA. Automatic embolus detection compared with human experts. Stroke 1996;27:1840 –1843. Wang Y, Fish PJ. Comparison of Doppler signal analysis techniques for velocity waveform, turbulence and vortex measurement: a simulation study. Ultrasound Med Biol 1996;22(5):635– 649.