Segmentation and analysis of colour Doppler images of tumour vasculature

Ultrasound

Pergamon

in Med. & Biol., Vol. 21, No. 5. pp. 635-647, 1995 Copyright 0 1995 Elsevier Science Ltd Printed in the USA. All rights reserved 0301.5629/95 $9.50 + .oO

0301-5629(95)00006-2

aOriginal

Contribution

SEGMENTATION

AND ANALYSIS OF COLOUR OF TUMOUR VASCULATURE

DOPPLER

IMAGES

DAVID S. BELL,JEFFXEY C. BAMBER and ROBERT J. ECISERSLEY Joint Department

of Physics, Institute of Cancer Research and Royal Marsden Hospital, Sutton, Surrey, UK (Received

2 August

1994; in fmal form

22 November

1994)

Abstract-A technique has been developed to segment(separate), from a digitized colour Doppler video image, the colour and greyscaleinformation and then to estimate from the colour information the original mean Doppler frequency shift data from which the image was created. The remapped velocity image is then analysed to extract numerical features of the tumour vasculature. The present version of the software is set-up to work for an Acuson 128 colour Doppler system using the V4 colour scale, although it should work well with any systemwhich modulatesonly two colours for each flow direction and displays a colour calibration scaleat the side of the image. Accuracy of classificationof greyscale, colour and flow direction was estimated as being in the region of 95% for typical breast tumour images.The degree of agreement between the remapped colour velocity values and those stated by the scanner at the sameimage locations was evaluated in terms of the linearity of the relationship ( >99%), precision (better than 5%) and accuracy (better than 7.6% ). We investigated the value, for diagnosisand assessment of responseof, a variety of characteristics of the displayed vascularity. At present, the software calculates the following vascular image features within any region of interest defined by the operator: mean displayed velocity, maximum displayed velocity, standard deviation of displayed velocity, total area occupiedby &our signal, percentagearea occupied by colour signal, area integral of displayed velocity and the total displayedvelocity per unit area. Key Words: Colour Doppler, Image analysis, Image segmentation, Three-dimensional imaging, Tumour blood flow, Tumour vasculature, Ultrasound.

which prompted the development of this flow segmentation and colour map inversion algorithm. The principal reason for using the video signal rather than the raw RF signal is because the program can then be made independent of the scanner. In other words, both the signal connections and acquisition software can be standardised. Also, an auxiliary colour video signal is invariably supplied with a scanner. The same cannot be said for RF signals, where scanner modifications would often be required. The analysis system can therefore be connected to any scanner.

INTRODUCTION Colour Doppler imaging provides real-time two-dimensional visualisation of blood flow in conjunction with the conventional pulse-echo information. An estimate of the instantaneous mean flow velocity, within some limits of spatial and temporal precision, is presented as a colour overlay and is generally assumed to correspond to hypoechoic regions. The resulting composite image contains both anatomical and blood flow information. This article describes a technique to segment (separate), from the digitized colour video image, the colour and greyscale information, and then to estimate from the colour information the original mean Doppler frequency shift or velocity data from which the image was created. The whole procedure was implemented on a commercial image processor so that it could be automated. There were two applications

The automation of a procedure to generate statistical descriptors of tumour blood Jlow to assist in the clinical evaluation of benign and malignant masses using colour Doppler $0~ imaging. Cosgrove et al.

( 1990, 1993) have demonstrated the value of colour Doppler examination for assisting in the differential diagnosis of breast tumours and have described a method for quantifying the percentage area of the screen occupied by displayed vascularity, by manually counting the coincidence of coloured image pixels with

Address correspondence to: Mr. D. S. Bell, Physics Department. Royal Marsden Hospital. Downs Road, Sutton SM2 5PT. UK.

635

636

Ultrasound

in Medicine

and

Biology

randomly positioned dots on a transparent overlay placed on the screen. This procedure, although used later by Shmulewitz et al. ( 1993 ), was time consuming and, for routine diagnostic work, a qualitative visual judgement of the displayed vascularity is still employed. There remains a clinical need for an automated method of quantifying the relative vascular density, which has become particularly pressing since the demonstration that relative vascular density is a useful measure of the degree of response of breast tumours to primary medical therapy ( Kedar et al. 1994). In addition, there remains a desire, once the colour Doppler image data is available in digital form, to investigate the value for diagnosis and assessment of response of additional characteristics of the displayed vascularity. For example, features are described later in this article which might act as relative measures of both overall tissue perfusion and complexity of the vascular pattern. Segmentation of colour and greyscale data as a preprocessing stage for three-dimensional display of tumour vasculature. Current colour Doppler systems display only a two-dimensional representation of anatomy and vasculature, which is inadequate for assessing the complex and disordered pattern of blood vessels associated with tumours. Preliminary work has shown that three-dimensional display of colour Doppler image data can assist visual interpretation and understanding of the pattern of tumour vascularity (Bamber et al. 1992; Eckersley et al. in press), but data for such studies must at present be obtained by scanning the tissue volumes with existing colour Doppler equipment. For effective visual display the set of images covering the volume must be segmented into blood flow and anatomical components, and then rendered two dimensional by a method which enhances the visual appreciation of the three-dimensional information. The method described herein may also be of interest in other areas. For example, there is increasing interest in quality assurance testing of colour Doppler imaging equipment (Evans et al. 1989; Rickey et al. 1992), which may benefit from an ability to extract quantitative information from the video image. The procedure to be described is equivalent to reversing the process of colour coding the displayed image from the estimated Doppler shift frequencies or flow velocities, known as the inverse transformation. In principle this procedure should not be necessary, since the required information is contained within the scanner. However, in general, this information has not been made available to the user by the manufacturers of scanners. Often the only available data are the red, green and blue analogue video signals corresponding to the displayed (composite) image. In fact, even if some manufacturers

Volume

?I,

Number

5, lOY5

were to make the required data available there remains theoretical merit in developing a system which. by working with the video output, would not be restricted to one scanner. There is no implication that identical performance with respect to tumour identification would be obtained, but merely that the analysis system can be connected to any scanner. To appreciate the steps required for segmentation and colour decoding, some aspects of the colour Doppler imaging process need to be understood. From elementary Doppler theory, the measured frequency shift is directly proportional to the velocity component along the acoustic beam axis. The displayed “velocity” therefore relates both to the flow direction and the velocity magnitude. These two terms cannot easily be separated. The sign of the Doppler shift indicates whether flow is towards or away from the transducer. This is the only directional information available. Each of the two flow directions is colour coded separately, the key to which normally takes the form of a pair of colour bars at the side of the displayed image. It is these colour bars which provide the reference for the colourto-velocity inversion mapping. The distance along each colour bar relative to the minimum value is taken to be proportional to the Doppler frequency shift. The validity of this assumption is confirmed later in this study by the results of tests to measure the accuracy of velocity remapping. The remainder of this article will describe the flow segmentation algorithm, the system used for its implementation, the method used to calibrate the algorithm for a particular scanner/probe combination and tests which were carried out to assess the errors in the velocity estimates and the effectiveness of pixel segmentation. Statistical measures intended for future assessment of tumour blood flow characteristics will also be described. APPARATUS Ultrasound scanner Analogue video signals corresponding to the red, green and blue components of a colour Doppler image from an Acuson 128 scanner (Mountain View, CA) were sent to the signal conditioning and digitization circuit of an Imaging Technology series 151 image processor (Boston, MA). These were processed according to the algorithm described later and the resulting velocity image was displayed on a colour monitor. The system has been installed in the clinic, for which a simple software user interface was developed. Both 3.5-MHz (electronic sector S328) and 7.5MHz (linear array L7384) probes were used. A variety of

Colour Doppler

in tumour

Fig. 1. An example of a colour Doppler “capture” image of a breast carcinoma, demonstrating the complexity of tumour vascularity, the wide range of flow velocities often observed (including aliasing of some of the high velocities), the use of the V4 colour scale and the colour scale reference bar on the left of the image. The image was obtained by using the scanner in “colour capture” mode to accumulate the peak values of mean velocity which occurred during a 4-s sweep across the region of highest vascularity.

vasculature

0

D. S.

BELL

637

et ul.

colour-coding schemes were available to display velocity information, not all of which were suitable for use in the flow decoding algorithm. The essential requirement for a colour code to be useful is that no more than two of the colours are needed to represent uniquely each flow direction. The reasons for this will become clear when the algorithm is discussed. The “V4” colour map of the scanner was chosen for use in this study. Figure 1 shows an example colour Doppler image of a breast carcinoma, displayed using the V4 colour map. The reference bar for the velocity-to-colour coding used by the scanner is seen on the left of the image. The velocity scale to which this colour bar refers can be selected interactively to one of about eight different values. The value selected refers to the maximum displayed velocity. The only additional colour Doppler control on this instrument is described as “CD level.” According to our own observations, this appears to function as if it were a gain control for a Doppler signal path used by a comparator to set the binary value of a signal which then controls whether or not to permit the display of the current estimate of mean Doppler frequency in the colour image. Image processing

The image processor is a VME-based modular image processor (Imaging Technology Inc., Series

a Fig. 2. Example (a) B-mode and (b) colour Doppler images of the rotating sponge phantom, showing in which the lower cut-off frequency is visible.

the manner

63X

Ultrasound

in Medicine

and Biology

151), designed to perform simple image processing operations at normal video rates, that is, at 25 frames per second. It is controlled by an IBM-PC fitted with a VME adaptor card. The system is used to digitize the analogue video signals, segment the colour and greyscale data and then to convert colour information to a frequency shift or velocity estimate. The principal features of the system are as follows. For image capture and display an g-bit video analogue-to-digital converter and signal conditioner with four multiplexed analogue input channels is used. Three digital output channels, each with display look-up tables and an 8bit digital-to-analogue converter drive an RGB colour monitor. The system possesses 2 megabytes of image storage memory, organized as eight 5 12 x 5 12 g-bit frame stores, with the ability for hardware pan, scroll and zoom. For image processing there is a general purpose integer processing unit, which can perform vector and scalar operations (some to g-bit and others to 16-bit precision) including add, subtract, multiply, divide, barrel shift, comparisons and look-up table transforms. The system is configured as a pipeline processor and therefore has some parallelism built into it. A basic operation consists of sending the contents of frame memory through a processing pipeline, the output of which returns the data to another frame memory. This takes one video frame time (40 ms), independent of how many distinct modules there are in the processing pipeline. Any algorithm implemented on this system is a set of such operations. There are a number of other commercial image processing systems which operate in a broadly similar manner, which from a programming point of view is not unlike the operation of a SIMD (single instruction multiple data) parallel computer. The algorithm to be described is thus not to be regarded as specific to this hardware. Indeed, the description of the algorithm is given in terms which should permit its implementation on almost any computer with adequate colour frame grabbing and display facilities. Experiments errors

to derive calibration

factors

and assess

A number of experiments are described later in which: (a) measurements were made of factors required to calibrate the colour decoding algorithm for a particular scanner/probe combination; and (b) the errors associated with the segmentation/decoding process were evaluated. These experiments were conducted with apparatus designed to permit repeatable cross-sectional imaging of a rotating speckle phantom. A degassed reticulated foam sponge of the type used by Bamber and Nassiri ( 1987)) with a typical cell size

Volume

2 1. Number

5. 1995

of 0.5 mm, was cut into a cylinder of height 45 mm and diameter 50 mm, immersed in degassed water and mounted on a stainless-steel disc which was driven by a stepping motor geared to provide 1000 steps per revolution. The geometry was arranged so that the clamped ultrasound array produced a cross-sectional image of the sponge. A typical colour Doppler image of this phantom is shown in Fig. 2. METHODS Flow segmentation algorithm

There are two distinct parts to the algorithm: classification and transformation. Classification consists of determining for each pixel whether the Doppler shift frequency is positive, zero (i.e., greyscale only) or negative. Transformation consists of using the classification data to control the conversion of the raw image data into a greyscale velocity magnitude image. Following classification, greyscale pixels are assigned a velocity value of zero, whereas positive and negative flow pixel values are estimated by passing the colour pixel data through look-up tables which are based on analyses of the displayed colour reference bars for the corresponding directions. To be able to create the look-up tables for colour-to-velocity transformation using currently available amounts of memory in small computers, the velocity in each flow direction must be a function of no more than two colours. This condition is satisfied by the so-called V4 colour scale on the Acuson, as shown in Fig. 1. Each of the two phases of the algorithm, classification and transformation, will now be described in more detail. Figure 3 shows the steps taken to achieve the separation between greyscale and both flow directions in the classification phase of the algorithm. An initial (preprocessing) step was to smooth the input images in the horizontal direction using a uniform four-element kernel. This was found to be necessary to reduce ringing artefacts introduced as a result of nonoptimal design of the analogue-to-digital converter input amplifier. Since the colour Doppler image pixel size was measured to be approximately eight greyscale pixels horizontally and four vertically, this horizontal smoothing did not significantly degrade the colour data. The ringing artefacts occurred at boundaries between colour and greyscale pixels and at boundaries between pixels of opposite flow directions. Since the ringing occurs to a different degree for each of the red, green and blue components, failure to apply the smoothing resulted in the existence of colours for which there were no corresponding red, green and blue combinations in the colour reference bars. In fact the smoothing only partially solves this problem and, as

Colour Doppler in tumour vasculature

0

D. S.

BELL

et al.

639

phase of the algorithm, that of decoding the colour data into velocity values. The transformation phase of the algorithm (Fig. 3) consisted of: (a) passing the smoothed red and green images through a specially constructed positive flow look-up table to produce a positive flow image; (b) passing the smoothed green and blue images through a negative flow look-up table to produce a negative flow image; (c) creating an image full of zeroes to represent zero flow; and (d) using the flag image from the classification phase to select which pixel (forward, negative or zero flow) is the correct one. The principal challenge was to generate the required look-up tables for the two colour-to-velocity transforms in stages (a) and (b).

0 N

Look-up table generation

set value to zero @pc&e

aebF’e

+ve flow OP

Fig, 3. Flow chart for the colour segmentation algorithm.

-ve flow O/p and

decoding

will be described later, the residual false colours contributed a source of error in pixel classification. The second stage of classification was to generate magnitude-of-difference images for red minus green and green minus blue. This was to identify forward and reverse flow pixels, the assumption being that red, green and blue were equal for nonflow pixels. Each of these difference images was then compared to a noise threshold, that is, the statistical variation in interchannel difference due to electronic noise in the original digitized and smoothed red, green and blue images. For each pixel in the difference images, a value greater than the threshold was used to classify the pixel as coloured (i.e., nonzero flow). At this stage, only the flow and nonfiow regions had been separated. Figure 4 shows the pixel values for line profiles along and through the centre of the colour reference bars. For the V4 colour scale, the positive flow coding scheme contains no blue component and the negative flow coding scheme contains no red component. Based on this information, it was therefore possible to separate the forward and reverse flows by comparing the red and blue images. This comparison produced a flag image, the pixels of which could take one of three different values, corresponding to zero, positive and negative flow. The flag image was used to control the second

As previously stated, the transformation look-up tables were generated from an analysis of the colour bars displayed by the scanner. The method relied on an assumption that the velocity was proportional to the distance along each colour bar from the centre. To set up a relationship between this distance and velocity required a knowledge of the minimum and maximum displayed velocities, for which separate measurements were made. The experiment conducted to make these measurements is described in the next section. Here we describe the method used for generating the look-up tables, assuming that the relationship between distance along the colour reference bars and velocity has been established. After digitization of a sample image, ten vertical profiles of the colour bars in Fig. 1 were measured at different horizontal offsets. These were then averaged to reduce electronic noise. The profiles (both original and averaged) have a staircase appearance as shown in Fig. 4. Inverse transformation of a particular colour combination to a velocity value involves reading the velocity value corresponding to the distances along the colour bar profiles in a two-dimensional look-up table. The only data that are available to create the two twodimensional matrices of velocity values (one for each direction) are those from the profiles shown in Fig. 4, which can be used to fill elements of the velocity matrices by using them in the forward direction to look-up the colour coordinates for particular velocity values. By displaying the two-dimensional velocity matrices as images in which grey level represents the desired velocity value, and the x and y coordinates are the two colour image pixel values, it may be seen, in Fig. 5, that using only the data from the displayed colour bars to fill the tables results in a very sparse matrix of values. In practise, due to factors such as electronic noise and the ringing artefact described earlier, colour

Ultrasound

in Medicine

and

Biology

-

Volume

I. Numher

5. 1995

180

0

200

Positidf? (pixels)

b

a

Fig. 4. Measured pixel data for colour components of line profiles through the colour reference scales for (a) positive flow and (b) negative flow. Position is measured in terms of distance in pixels from the low velocity end of the colour scale bar.

combinations may be observed in real digitized images for which there will be no entry in the velocity matrices. Several problems can be identified, which prevent such a table from being used to transform a measured colour combination back to a velocity value. These will now be described along with their solution in this implementation of the algorithm.

Green component

a

The first problem is that the staircase nature of the colour bar profiles produces gaps in the range of possible colour combinations. To avoid this phenomenon, the averaged profiles were smoothed before being used. This was done by drawing smooth curves through the averaged profiles by hand, as shown in Fig. 6. Since the look-up tables do not need to be generated

-b

Blue component

*

b

Fig. 5. Velocity look-up tables displayed as colour-space greyscale images for (a) positive and (b ) negative flow, when filled only with data available from the profiles of Fig. 4.

Colour Doppler in tumour vasculature

D.

0

S. BELL

et al.

641

r-240

c

blue

-190 -220 F 8 g2x -5 s -E -140

2 a

- 120

looT .x a

I

0

200

PositioY

I

I

0

(pixels)

a

I

I

I

I

I-

70

200

PositioY(pixels)

b Fig. 6. Colour bar profiles of Fig. 4, after smoothing

very often, this manual method was felt to be adequate. Once the smoothed curves had been generated, they were digitized with a Summasketch II graphics tablet. The second potential problem concerned gaps which could remain in the locus of the path of the reference colours through colour space (i.e., the velocity matrices) even after the velocity matrices were filled using the smoothed profile data. However since the digital dynamic range was maximised for the colour data, the number of such gaps was very small. Measurements made on sample images indicated about 1% missing points. This is small enough to be neglected. A third problem arose because electronic and digitization noise may produce colour combinations which do not fall on the locus of the reference colours. This was solved by filling remaining elements of the velocity matrices, to a distance from the reference locus judged to contain practically all image pixels. For Gaussian noise, this is approximately equal to two standard deviations of the colour noise. The calculations were done using magnitude-of-difference images of red-green and green-blue image pairs. Values were computed from velocity values on the locus by nearest neighbour extrapolation. The basis of this is an assumption that noise will displace a colour from the main reference locus with equal probability for all three colours. Thus, it seemed reasonable that on measuring a colour combination for which there were no reference values, the velocity value for the nearest point on the reference locus should be chosen. Outside

by eye.

of this broadened reference locus the velocity values in the matrices were allowed to remain zero. Finally, the horizontal smoothing used to reduce ringing artefacts (described earlier) itself produced a problem in which artefactual colours were generated at boundaries between flow pixels and greyscale pixels and boundaries between positive and negative flow pixels. These colours placed flow pixels outside the region of even the broadened and interpolated reference locus, so that they would have been given a velocity value of zero. One approach to a solution might have been to extend the nearest neighbour extrapolation of the reference locus to cover the whole of the velocity matrix, rather than just the “noise band.” This would have had the equally unfortunate result of giving extremely high velocity values to the artefactual colours. A compromise was found in which the colour spaces of both velocity matrices were divided in half using the green component. At flow-greyscale boundaries the velocity is low, and therefore the green component is low. At flow-flow boundaries the velocity is high, and therefore the green component is high. This is the basis of using the green component to deal with artefactual colours. A green threshold was set to perform the discrimination, the precise value determined from measurements made on a few sample images. When the green component is high the peak velocity is used. When the green component is low zero velocity is used. One important consequence of this is that, on all images which display a flow reversal across any boundary between two adjacent pixels, the maximum velocity will equal the peak value.

Ultrasound in Medicine and Biology

Green component

Volume 2 I, Number 5. 199s

-b

Blue component

+

E is) 8 E 8

a

b

Fig. 7. (a) The velocity matricesof Fig. 5, after filling, using all stepsdescribedin the text; that is, use of the smoothedprofiles shownin Fig. 6, nearest-neighbour interpolationof the main referencelocus,nearest-neighbour extrapolation to the noiseboundary and binary division of the remainderusingthe greencomponent.(b ) These featureshave beenlabelled.

Figure 7 showsthe final two velocity matrices, generated using all of the steps described above. These are

completely filled and may be used to look-up a velocity value given any measured combination of red and green colours for forward flow and any measured combination of blue and green for reverse flow. Measurement of displayed velocity limits

As mentioned earlier, the generation of the transformation look-up tables required knowledge of the relationship between velocity and distance along the colour reference bars. Assuming a linear relationship, this in turn required knowledge of the maximum and minimum displayed velocities. The velocity scale value displayed on the Acuson was assumed to be equal to the maximum displayed velocity. This was verified by using the rotating sponge and noting that all velocities were less than or equal to the scale value. The scanner does not give any information regarding the minimum displayed velocity. Therefore, it must be measured, since it has a profound effect on the transformation look-up tables. The rotating sponge described in the apparatus section of this study was used for this purpose. For all ultrasound probe/velocity scale configurations of the scanner the following procedure was adopted. The gain and TGC of the scanner were adjusted to produce a uniform B-mode cross-sectional image of the sponge cylinder. With the scanner switched to “RES” mode, which causes this scanner to operate at its maxi-

mum spatial resolution over a selected region of interest, the RES box was adjusted to give maximum resolution over the circular image of the sponge. The sponge was removed from the field of view and the colour Doppler was switched on. The “CD level” control was then adjusted by gradually reducing it from maximum until no colour noise was produced, after which the sponge was then moved back into the field of view. Then, while imaging the sponge in colour Doppler mode, the speed of rotation of the sponge was adjusted to give the full range of colours without significant aliasing. Having set up the scanner, the speed of rotation of the sponge was set so that no colour was displayed within half a radius from the centre of the sponge. This ensured that the velocity increase from the centre was low enough to cause the minimum velocity to be displayed. The Acuson colour cursor was then used to measure this minimum displayed velocity. This value was measured at about five distinct points to provide a cross-check. The measurement was made for both positive and negative velocities. The minimum displayed velocity was found to be directly proportional to the velocity scale and the constant of proportionality independent of the probe type. Statistical features of tumour blood jfow

In the current implementation, after the magnitude-of-velocity images have been derived, a set of local statistical colour Doppler features are calculated

Colour

Doppler

in tumour

within a user-defined rectangular region-of-interest. There are seven features: mean colour velocity (MCV ) , peak colour velocity (PCV) , standard deviation of colour velocity (SDCV), total colour Doppler area (CDA), percentage colour Doppler area ( %CDA), integrated colour velocity (ICV) and integrated colour velocity per unit area (%ICV). Equations for each of these are as follows: MCV = c V(x, Y) NF PCV = max(V(x, SDCV =

y))

c (V(x9 Y))’ _ (MCV)2 NF

CDA = NF %CDA = 100x $ ICV = c V(X, y) %1(-V = c V(x- Y) N where V ( x, y) is the array of colour velocity values, N is the total number of pixels in the analysis box, NF is the number of flow pixels in the analysis box. Summations were over the entire analysis box. Measurement

of errors

Assessment of classi$cation accuracy. There is a facility on the Acuson for turning the colour Doppler data on and off for a frozen image. This is referred to as colour toggle and provides a standard against which the algorithm can be compared. Ten images were acquired comprising two clinical, two colour noise and six sponge images. Each of these was used to produce a “truth table” relating the number of pixels in each of the three classes for reference (colour toggle) and test (classification algorithm) images. For each image, four separate images were identified corresponding to the red, green, blue and greyscale components. The greyscale image was simply the red image with the colour turned off. These four images were used to generate two flag images, the pixels of which contained information relating to velocity direction. Each flag image comprised three values, corresponding to the three classes (zero, positive and negative flow). There were two flag images, corresponding to the colour toggle and the classification algorithm. The colour toggle flag image was generated by subtracting the greyscale image from the

vasculature

643

0 D. S. BELL ef al.

red and blue images and thresholding each of these with a value of ten to eliminate electronic noise. The resulting images were used to define the flag image. The classification algorithm flag image was generated as a byproduct of the algorithm, shown in Fig. 3. Having generated the flag images, tables were constructed to compare corresponding pixels. The tables comprise three rows and three columns. Each column corresponds to a class in the colour toggle flag and each row to a class in the classification algorithm flag. Every pixel pair will map to a specific position in the table, where the value is incremented each time a match is found. The ideal situation is for all nonzero elements to be on the main diagonal, indicating all pixels correctly classified. The sum of all nine table entries must be equal to the number of pixels analysed. Direct velocity comparison. The Acuson has a facility for making point velocity measurements with a crosswire. By displaying this crosswire on the image processor, the digitised image pixel lying closest to a given Acuson colour Doppler pixel could be located and the corresponding calculated velocity measured. The procedure was to set up the rotating sponge to give the maximum displayed velocity at the edge of the sponge, using the L7384 probe. After positioning the Acuson crosswire to read a velocity, a crosswire was displayed by the image processor and overlaid on the Acuson crosswire as displayed on the image processor screen. The velocity magnitude image was then displayed and the value at this point was noted. The procedure was repeated for a number of velocity values covering the range up to the maximum displayed velocity. Results were generated for both positive and negative how. The velocity comparison data were analysed to calculate values for linearity, precision and accuracy. For the linearity measure, a quadratic curve was fitted to each of the datasets. Then the linear coefficient was divided by the sum of the linear and quadratic coefficients to generate the linearity measure. This was expressed as a percentage. The precision was estimated as the root mean square error in the velocity value from the algorithm, estimated from a linear regression assuming no errors in the Acuson velocity value. This was expressed as a percentage of the Acuson velocity scale value. The accuracy was estimated as the slope of the above regression line as compared to the theoretical value of one. This again was expressed as a percentage. Equations describing these calculations are as follows: LINEARITY

= 100

COEFF, COEFF, + COEFF?

llltrasound in Medicine and Biology

Volume

31, Number

5. 19%

Table 2. Worst case classification accuracy (obtained an image of the rotating sponge phantom)

from

Colour toggle (reference clnsses,I Classitication algorithm

Zero flow

Positive

Negative

Zero flow Positive Negative

6327 5686 4181

47 12,610 403

192 227 18,321

Total correctly classified = 77.6%.

Fig. 8. An example velocity magnitude image, obtained by processing the breast carcinoma in Fig. 1.

PRECISION = x

100 VELSCALE

X ( V,p - (GRADIENT.

VA + INTERCEPT)) * N

ACCURACY

= ; c

(

100

where GRADIENT and INTERCEPT refer to the iinear regression line. V, and V,p refer to Acuson and image processorvelocities, respectively. COEFF refers to the first- and second-order quadratic fit coefficients. VELSCALE is the Acuson velocity scale. N is the number of points used for the comparison.

RESULTS Segmentation accuracy Figure 8 illustrates a velocity magnitude image produced by the algorithm and obtained by processing Table

1. Best case classification accuracy (obtained an image of a breast tumour)

from

Colour toggle (reference classes) Classification algorithm

Zero flow

Positive

Negative

Zero flow Positive Negative

15,816 467 854

5 2321 53

0 8 14,396

Total correctly classified = 95.9%.

the image of the breast carcinoma in Fig. 1. Black corresponds to zero flow and white corresponds to maximum flow; that is, equal to the velocity scale value. There is no directional information. It can be seen that all important vascular structures appear to have been preserved and coded in the final flow image, the remaining B-mode information having been successfully excluded. Tables 1 and 2 show the results of the classilication accuracy assessmentfor the best and worst cases, respectively. Each table’s entry is the number of pixels in the specified class. It happened that the best case was for a colour Doppler image of a breast carcinoma and that the overall classification performance for this image was about 96%. The worst case was for an image of the rotating sponge phantom, with an overall performance of 78% correct. Information on the total number of pixels correctly classified for all ten images studied is given in Table 3. Velocity values Figure 9 showsthe observed relationship between velocity calculated by the image processor and that stated by the Acuson scanner. Table 4 summarisesthe statistical properties of these graphs in terms of the computed linearity, precision and accuracy values for the velocity data. The Acuson velocity scale was 6 cm/s. The ratios of minimum displayed velocity to

Table 3. Total numbers of pixels correctly segmentedfor each of the ten test images Object

Classification accuracy

Breast cancer Breast cancer Colour noise Colour noise Rotating sponge Rotating sponge Rotating sponge Rotating sponge Rotating sponge Rotating sponge

95.9% 93.7% 88.6% 84.7% 88.6% 94.3% 88.6% 88.3% 88.2% 77.6%

Colour

Doppler

in tumour

vasculature

0 D. S. BELL

“A

60-

*

-60 I

p20-

-50 1

635

et ul.

-40 I

(mms-‘) -30 -20 I I

-10

L

0

0

--40

5-a

- -50

10l *

0 0

I

I

10

20

I 30

1

I

I

40

50

60

- -60

(mms-‘)

VA

a

b

Fig. 9. Measuredvelocity (V,p), calculatedby the imageprocessor,plotted against“true” velocity (VA) obtained from the Acuson colour cursor.

velocity scaleswere found to be 0.0941 -+ 0.0001 and 0.1252 + 0.0002 for positive and negative velocities, respectively.

DISCUSSION There are undoubtedly many alternative ways in which the problem of segmenting and decoding colour Doppler images could be approached. We make no claim that the method described here is necessarily the best. We have, however, demonstrated that it works and have provided an assessmentof the accuracy and precision of important facets of its operation. We have also described specific problems (and methods for their solution), associated with artefacts arising from the digitising of colour video images, which are likely to be present whatever algorithm is constructed for the segmentation and decoding process. In fact it was partly these artefacts which limited and led us to abandon, in favour of the present method, the development of an earlier attempt by us to accomplish the sametask by first transforming the digitised red, green and blue

Table 4. Resultsof a comparisonbetweenthe velocity as statedby the Acuson scannerand the velocity calculated by the transformationalgorithm Positive Correlation Linearity Precision Accuracy

coeff.

flow

0.981 99.9% 24.1% +1.6%

Negative

flow

0.991 99.1% ?3.2% -4.5%

images into a hue, saturation and intensity (HSI ) colour representation. The intention had been to find a colour mapping schemeused by the scannerfor which the velocity could be uniquely identified by hue value only. A simple method that can be usedfor work, which requires only a segmentation of the colour and greyscale images, is to digitise an appropriate colour component (e.g., green) with and without the scanner switched to colour Doppler mode. Simple subtraction of the images can then be used to isolate the flow component (Humphreys et al. 1992). This method however gives no quantitative flow information and requires manual intervention for each frame to be processed.It was therefore not suitable for our intended purposes. Tamura et al. ( 1991) have analysed colour Doppler images from an Acuson 128 scanner by using the velocity cursor facility to calculate velocity profiles. These measurements were combined with a pulse Doppler spectral analysis. This is not an automated analysis technique. Vattyam et al. ( 1992) have used a PC video frame grabber attached to an Acuson 128 scanner to perform quantitative analysis from video tape. No technical details were given on the method used for velocity segmentation. With regard to the methods chosen for estimation of the errors in the final velocity values, it would have been possible to use a secondrotating sponge, knowing the position of the centre of rotation on the image, to plot measuredvelocity against known velocity. However, such an experiment would be a combined assess-

646

Ultrasound

in Medicine

and

Biology

ment of our velocity decoding algorithm and the accuracy of the scanner’s Doppler system, mean frequency estimator and display processing. It was felt that these latter factors fell beyond the scope of the present study, but that such an experiment might be of considerable interest to those interested in quality assurance testing of colour Doppler flow imaging equipment. Figure 9a shows one point on the horizontal axis which departs significantly from the straight-line fit. It suggests that a calculated velocity of zero results from a nonzero Acuson velocity. The data used to generate Fig. 9 comprised only about 30 spot values of Acuson velocity. It seemed unreasonable to discard this data point on the basis of this limited data set. The intercept on Fig. 9a is clearly not zero. This again may be due to a limited data set. The measurements were extremely time-consuming and as such were not repeated. The speed of operation of the algorithm on the specialised image processing system described here is I .5 s, plus 6 s to calculate the statistical blood flow features for an analysis box equal to the entire image region for the L7384 probe (approximately 143,000 image processor pixels). This is sufficiently fast for routine clinical use. The limiting factor for speed if the algorithm were to be implemented on a general purpose computer using only the main central processing unit to perform all operations would probably be the data transfer rate to and from conventional memory which initiates and terminates the processing chain. The look-up table transformations do not in themselves present a bottleneck once the data are in conventional memory. These transformations in the present system are done in a single video frame by specialised hardware. For implementation on a 50-MHz Intel 80486 processor, the expected execution times are 2 s, plus 0.6 s for the statistics. In this article, we have made no attempt to evaluate the seven statistical features calculated by the analysis software but it may be helpful to describe briefly what we expect each of the features to be related to in terms of soft tissue vascularity and blood flow. Total colour may be thought of as being related to the total amount of tissue vascularization (i.e., vascular volume if the results were to be integrated over many slices through the tissue), and percentage colour the vascularity per unit volume of tissue. These relationships are only true if the flow is detected and subsequently displayed. It has been realized that a large amount of vascularity characterized by very low amplitude and high velocity flow (e.g., from tiny AV shunts) or very slow flow may not be shown in the colour display of many scanners, even though it may contribute significantly to the total tissue vascularization.

Volume

2 I. Number

5. 1995

However, the assumed (approximate ) relationships constitute a useful working model. These two quantities may be expected to depend on the CD gain and velocity scale settings of the scanner to the extent that these settings affect the thresholds of signal strength and llow velocity for display 01‘ a colour signal. They should, however, be independent of whether aliasing occurs. It is important to bear in mind when interpreting these features that the area of colour displayed for a single, small diameter blood vessel is likely to be influenced by the length of vessel lying within the scan slice width, the diameter of the vessel insofar as it affects the signal strength, other factors (such as depth and the presence of overlying attenuating structures) which may affect the signal strength and whether the proportion of flow for which the Doppler component of blood velocity is above the threshold for display. On the other hand, the diameter of the vessel is only likely to affect the colour area in a direct manner when the diameter approaches or exceeds the colour Doppler resolution. Thus, for regions of soft tissue perfused by small vessels (and when signal strengths and flow velocities are broadly similar) comparisons of colour area measures will amount to comparisons of numbers and lengths of blood vessels. Mean velocity, peak velocity and standard deviation of velocity will be directly affected by the presence of aliasing. If these features are to be at all meaningful it will be necessary for the operator to ensure that the Doppler pulse repetition frequency (controlled by the colour Doppler velocity scale) is set so that aliasing does not occur at any point in the image. This would most easily be done by using the most sensitive velocity scale for which the peak colour velocity in the image does not reach the highest value on the colour map of the scanner. These velocity-based features are of interest to us because of their likely diagnostic significance in distinguishing between benign and malignant tumours, particularly in picking up: (a) the existence of shunt vessels in tumours from high peak velocities; (b) simultaneous very slow flow and fast flow in tumours from large standard deviation of velocity; and (c) widely varying vessel orientation from large standard deviation of velocity. With regard to the peak velocity. if there is any aliasing, then the calculated value will always be equal to the velocity scale value; that is, the maximum possible velocity. This is likely to occur for nearly all images. To justify this approach it can be seen that if aliasing has occurred the velocity has exceeded the maximum. Therefore, the peak velocity should be set to this maximum value. The quantities ICV and ICVA, since both combine

Colour

Doppler

in tumour

velocity estimates with vascularization indices, may be thought of as crude relative indices of volume flow and volume flow per unit volume of tissue (i.e., perfusion), respectively. Both, however, will be subject to aliasing artefact. Other difficulties to be wary of in using an automated analysis system such as this are related to the standardisation of the colour Doppler gain setting and the exclusion of unwanted colour signals which may arise from flash artefact or colour noise. Given the complex process by which colour Doppler images are produced and the doubt about the precise relationship between the features calculated and the true tissue vascular characteristics, there is an immediate requirement to test the clinical validity of this colour Doppler image analysis system. This is, in fact, the subject of two independent clinical trials aimed at evaluating measures of tissue perfusion (Hirsch et al. in press) and differential diagnostic criteria (Kedar et al. in press). CONCLUSION We have described a method for separating, from digitised colour Doppler video images, the colour and greyscale components of the information and then to estimate from the colour information the original mean Doppler frequency shifts. This is accomplished without manual intervention, at reasonable speed, and with errors believed to be sufficiently small for routine clinical use. We have also provided examples of tumour blood flow characteristics which may be obtained from an analysis of the separated flow images. Some of these features have been designed with the application of differential diagnosis of benign/malignant tumours in mind while it is hoped that others may provide relative estimates of soft tissue perfusion. The hardware and software system described is in use for collection of clinical data on the colour Doppler characteristics of breast tumour blood flow ( Kedar et al. in press) and has been employed to provide the basic segmentation information for a volumerendered three-dimensional display of tumour vasculature (Bamber et al. 1992; Eckersley et al. in press). Hirsch et al. (in press) have shown that some of the

vasculature

0 D. S. BELL et al.

641

features calculated do indeed provide a good relative measure of tissue perfusion, at least in human calf muscle. The rotating sponge phantom used for calibration and error assessment may be of interest in the area of quality assurance testing of colour Doppler equipment. REFERENCES Bamber, J. C.; Nassiri, D. K. Spatial resolution and information content in echographic texture analysis. In: McAvoy, B. R.. ed. IEEE 1986 Ultrasonics Symposium, IEEE Cat. No. 86CH23754. Piscataway, NJ: IEEE; 1987:937-940. Bamber, J. C.; Eckersley, R. J.; Hubregtse, P.: Bush, N. L.; Bell. D. S.; Crawford, D. C. Data processing for 3-D ultrasound visualisation of tumour anatomy and blood flow. In: SPIE, Vol. 1808. Bellingham, WA: Society of Photo-Optical Instrumentation Engineers; 1992:65 l-663. Cosgrove, D. 0.; Bambcr, J. C.; Davey, J. B.; McKinna, J. A.; Sinnett, H. D. Colour Doppler signals from breast turnours. Radiology 176:175-180; 1990. Cosgrove, D. 0.; Kedar, R. P.; Bamber, J. C., et al. Breast disease: Colour Doppler ultrasound in differential diagnosis. Radiology 189:99-104; 1993. Eckersley, R. J.; Bamber, J. C.; Bell, D. S.; Bush. N. L.: Cosgrove, D. 0. Potential and limitations of colour Doppler ultrasound for three dimensional imaging of tumour vasculature. Eur. J. Ultrasound [in press]. Evans, J. A.; Price, k.; Luhana, F. A novel testing device for Doppler ultrasound eauioment. Phvs. Med. Biol. 34:17011707; 1989. Hirsch, W.; Bamber, J. C.; McCready, V. R.; Bell, D. S. Colour Doppler image analysis for tissue vascularity and perfusion: a preliminary clinical evaluation. Ultrasound Med. Biol. [in press] Humphreys, K.; Reynolds, R. A.; Wilson, K.; O’Keefe, D.; Pitt, T. The development of 3-D ultrasound venographic techniques for the lower limb venous system. BMUS bulletin, November 92. London: The British Medical Ultrasound Society; 1992:24-25. Kedar, R. P.; Cosgrove, D. 0.; Smith, 1. E.; Mansi, J. L.; Bamber, J. C. Breast carcinoma: Measurement of tumor response to primary medical therapy with colour Doppler flow imaging. Radiology 190:825-830; 1994. Kedar, R. P.; Cosgrove. D. 0.; Bamber, J. C.; Bell. D. S. Automated quantification and analysis of colour Doppler signals: A preliminary study in breast tumours. Radiology [in press]. Rickey, D. W.; Rankin, R.; Fenster, A. A velocity evaluation phantom for colour and pulsed Doppler instruments. Ultrasound Med. Biol. 18:479-494; 1992. Shmulewitz, A.; Teefey, S. A.; Codwell, D.; Starr, L. Temperature dependent ultrasound color flow imaging in the study of a VX2 tumor in rabbits: preliminary findings. Ultrasound Med. Biol. 19:221-229; 1993. Tamura, T.; Cobbold, R. S. C.; Johnston, K. W. Quantitative study of steady flow using colour Doppler ultrasound. Ultrasound Med. Biol. l7:595-605; 1991. Vattyam, H. M.; Shu, M. C. S.; Rittgers, S. E. Quantitication of Doppler colour flow images from a stenosed carotid artery model. Ultrasound Med. Biol. 18:195-203: 1992.