## E

Considering the image to be formed by a noise-free ideal image /(i, j) and a noise process n(i, j) such that g(i, j) =/(i, j)+ n(i, j), it is straightforward to show that the local variance of /(i, j) is given by oy 2 (i, j) = ffg2(', j) — ff„2(i, jO,

where c„2(i, j) is the nonstationary noise variance assuming that /(i, j) and n(i, j) are independent. The function ff„2(i, j) is assumed to be known from an a priori measurement on the imaging system. An example of local statistics computation of an imaging system is the typical degradation model of Poisson noise in mammographic screening systems:

where Poisson^{ • } is a Poisson random number generator, and X is a proportionality factor such that

The conditional ensemble mean and variance of g*(i, j) given / j) are where gl is the vector whose nonzero element is the output gAp (i, j), and the elements of c are set adaptively according to the weighting factor b(i, j) such that

Var[g *(i, j)|f (i, j )]= Af (i, j). The normalized observation is defined as

Therefore, the noise part has the form

«(i,j) = g(i> j) -/(i>j) =-^^ - /(i>j) (13)

and its variance can be shown to be ff„2 (i, i) = E (f (i, j))/l (14)

From these equations, the ensemble variance of f (i, j) can be obtained as ff2(i,j) = ^(i,j) - (f(i,j)/A) = ^(i,j) - ((i,j)/A). (15)

2.2 The Multiresolution/Multiorientation Wavelet Transform

### Multiresolution Wavelet Transform

The WT utilizes two functions, the mother wavelet n(x) that spans the subspace Wi, and a scaling function 0m„(x) that spans the subspace Vi. The function f is subjected to the functional operations of shifts and dyadic dilation, and the WT may be implemented by using filter banks that have good reconstruction properties and high computational efficiency. Two-, three- and four-channel wavelet transforms were initially used for preliminary studies of medical image enhancement and segmentation [5,9,10,12]. An M-channel WT decomposes the input image into M2 subimages, and the 4-channel WT that we illustrate has 16 decomposed subimages. It is important to select the subimages with significant diagnostic features for reconstruction.

The dilations and translations of the scaling function induce a multiresolution analysis of L2 (.R), in a nested chain of closed subspaces of L2(P)(... c V— c V0 c V c V2 ...) such that

V; and W are related by

where W4 is the subspace spanned by mother wavelet function f.

In order to apply wavelet decompositions to images, two-dimensional extensions of wavelets are required. An efficient way to construct this is to use "separable wavelets" obtained from products of one-dimensional wavelets and scaling functions. The two-dimensional functions can be constructed as

The functions of (17) correspond to the separable two-dimensional filter banks. The basis functions for a 4-channel wavelet transform, where 0(x, y) e V4 and Wj(j = 4, 3,2,1), can be ordered by

03(x) 03 (y) > 03(x) f 3 (y) > 03(x) f 2 (y) > 03(x) f 1 (y) f 3(x) 03 (y) >f 3(x) f 3 (y) >f 3(x) f 2 (y) >f 3(x) f 1 (y)

f 2(x) 03 (y) >f 2(x) f 3 (y) >f 2(x) f 2 (y) >f 2(x) f 1 (y) f 1(x) 03 (y) >f 1(x) f 3 (y) >f 1(x) f 2 (y) >f 1(x) f 1 (y) •

In terms of two dimensional basis functions, the set in (18) can also be represented by the block matrix a3x0(x> y)a3y0(x> y)> a3x0(x> y)b3y0(x> 7)> a3x0(x> 7)52y0(x> 7)>

b3x0 (x> y) ^3y 0 (x> y) > b3x0 (x> y) b3y 0 (x> /) > b3x0 (x> y) b2y 0 (x> /) >

b2x0 (x> y) ^3y 0 (x> y) > b2x0 (x> y) b3y 0 (x> y) > b2x0 (x> y) 52y 0 (x y) >

51x 0(x> y)^37 0(x> y )> 51x 0(x> y )537 0(x> y )> 51x 0(x> y)527 0(x> y )> 51x0(x, y 0(x, y)

where A and B represent linear operators. The wavelet representation of the operator just described has interesting interpretations. First, the representation means that the 2D image is decomposed into 16 subimages, which allows the multichannel image processing techniques to be used. Second, in the foregoing wavelet-based image compression, A3x0(x, y)A3y0(x, y) is the "low-frequency" portion of the operator, while other blocks in Eq. (19) are "high-frequency" components of the operator. This motivates a general and consistent scheme for "scaling down" from a fine to a coarse grid. Compared to a 2-channel wavelet transform, a 4-channel, multiresolution wavelet transform provides more accurate representation, selective reconstruction of the higher order M2 subimages allows better preservation of the enhanced structures, and better detection sensitivity can be obtained on enhanced images. Adaptive linear scaling should be used to adjust image gray scale intensity and allow operator-independent evaluation of the enhanced image.

### The Multiorientation Wavelet Transform

The directional WT (DWT) is a wavelet transform for multiorientation signal decomposition using polar coordinates implemented on a pixel-by-pixel basis [16]. This approach allows enhancement of structures that have a specific orientation in the image. The selection of the wavelet basis functions is important to obtain a high orientation selectivity. This is achieved by introducing a directional sensitivity constraint on the wavelet function [17]. The input image is decomposed by the DWT, yielding two output images. One is a directional texture image, used for directional feature enhancement. The second is a smoothed version of the original image, with directional information removed, and later used for image background correction. This module plays two important roles in mammogram analysis: (1) isolation of the central mass from surrounding parenchymal tissue in the case of stellate tumors as required later for enhancement and segmentation, and (2) direct detection and analysis of spiculations and their differentiation from other directional features in the mammogram using a ray tracing algorithm.

Directional sensitivity can be obtained by retaining only the components that lie in the desired orientation in the WT domain. Selection of the orientation can be achieved with a fan having a desired angular width and oriented in the desired orientation in the WT domain. All orientations can be analyzed by considering multiple such fans, positioned adjacently. For example, four fans with 45° angular widths can cover all orientations. By taking the inverse WT of only those components in a given fan, we can produce an image where structures in the selected orientation have been emphasized. The four fans with 45° width give rise to four images and the dominant orientation corresponds to the image with highest gradient. Since the dominant orientation can change across the image, it is wise to make this decision individually for each pixel. Therefore, for each pixel, we can compute the local gradient magnitude G on each of the four images and determine the highest. The pixel value in the image with the highest local gradient is taken as the orientation selective value of that pixel for the 45° fan.

The appropriate fan width can also change from image to image, and narrower fan widths may be better in many images. Therefore, the DWT is implemented by using multiorientation filter banks with nine fan widths, 45° (4 orientations), 30° (6 orientations), 22.5° (8 orientations), 15° (12 orientations), 11 25° (16 orientations), 9° (20 orientations), 7.2° (25 orientations), 6° (30 orientations), and 5.63° (32 orientations). In each orientation, the maximal gradient magnitude is defined as Gn, i = 1,2, ...9 The technique of the adaptive directional filter bank (ADFB) is shown in Fig. 5 where the output depends on the maximal gradient among all orientations, defined as

Gmax = maxK |Gi|, 0C2 I G2 I, «3 I G3 I, — , «9 I G9 l> (20)

where «(i = 1, 2, 3, 9) are normalizing factors to make all «iGi have the value given by a unit arc area and P(i = 1,2, 3,... 9) are adaptive control parameters If Gmax = «G, then p = 1 and Pj = 0 for j = i. For example, when Gmax = «9G9, the output of the 32-channel directional filter is used as the

 4->CrunnQi DFB ft ,
0 0