1.IntroductionOptical mammography is a technique that employs near-infrared light in the wavelength range of 630 to 1000 nm to image the human breast. In this wavelength range, the main source of image contrast is the absorption of hemoglobin. Angiogenesis, or increased blood vessel formation, and other hemodynamic and oxygenation changes are often associated with the presence of breast tumors.1 These changes typically produce a local increase in optical absorption that allows the use of near-infrared methods to detect the presence of breast tumors. Our goal is to find a way of enhancing the detection of those areas in the breast image that exhibit a local maximum in absorption and might therefore correspond to the location of tumors. We show that a robust way of accomplishing this is by employing the spatial second derivative, an operator commonly used in edge-detection algorithms.2 2.BackgroundWe have previously reported “edge-corrected” optical mammograms, which we call N -images, that enhance the image contrast and detectability of breast tumors in transillumination images.3 An N -image is constructed from raw phase and amplitude data acquired by a frequency-domain instrument that operates in a planar projection geometry. Examples of frequency-domain instruments operating in this geometry are the research prototypes previously developed by Zeiss4 5 and Siemens.6 Figure 1 shows a schematic diagram of the Siemens prototype, the technical specifications of which have been reported elsewhere.7 The need for an edge correction is due primarily to the variable thickness of the breast between the two glass plates that slightly compress it. The N -parameter is defined as where r0 is the separation of the plates, ac0 is the amplitude of the intensity-modulated signal at a pixel where the breast thickness is r0, ac(x,y) is the amplitude measured at pixel (x,y), and r(x,y) is the breast thickness at that pixel derived from the phase information available in the frequency domain. We observe that in this approach, the phase information is solely used to measure the thickness of breast tissue at each pixel, and not to separately measure the absorption and scattering properties of breast tissue. A representative N -image of the left breast, craniocaudal projection, for a 58-year-old patient affected by invasive ductal carcinoma is shown in Fig. 2(a). A line graph of N at y=3.2 cm is shown in Fig. 2(b).Figure 1Schematic diagram of the frequency-domain research prototype developed by Siemens Medical Engineering, Erlangen, Germany, for optical mammography. PMT, photomultiplier tube.  Figure 2(a) Edge-corrected N -image at 690 nm of the left breast, craniocaudal view (lcc), of 58-year-old patient 310. (b) Graph of the N data versus x along the line y=3.2 cm shown in (a). (c) Second-derivative image at 690 nm of the same breast shown in (a). (d) Graph of the negative second derivative of N (−N″) versus x along the line y=3.2 cm shown in (c). The arrows in (a) and (c) indicate the location of a 3-cm invasive ductal carcinoma.  Although N -images represent a significant improvement in quality over the images based on raw intensity data, it is nevertheless often difficult to discern much structure in them. For example, the secondary peak around x=8.8 cm along the line y=3.2 cm is visible in Fig. 2(b), but it is almost completely hidden in the N -image of Fig. 2(a) by the large dynamic range that is due to the main peak at x=12 cm. However, we have found that by calculating the spatial second derivative (N″) and considering only those pixels where the second derivative is negative, we can isolate areas in the N -image characterized by a local maximum in absorption. The resulting N″ -image is presented in Fig. 2(c). In Fig. 2(d) we have plotted the negative of the second derivative along the line y=3.2 cm and set a threshold at N″=0. Notice how the feature around x=8.8 cm, which is most likely due to a blood vessel, is now much more visible, both in the line graph and in the N″ -image. It is also evident from Fig. 2(d) that there is additional structure around x=6 cm, and this, too, is visible in the second-derivative image, even if with relatively low contrast. Overall, the N″ -image is a higher-contrast image that makes it easier to spot those areas in the breast where there is a local maximum in absorption, and therefore potentially smaller tumors that might have been eclipsed by other features. 3.Second-Derivative MethodOur method of generating second-derivative images from edge-corrected optical mammograms (N -images) is summarized in Fig. 3. We discuss each step in detail in the following sections. All of the processing was carried out in Interactive Data Language ((IDL), a higher-level programming language developed by Research Systems Inc.) on a Pentium IV personal computer. Figure 3Conceptual representation of the steps involved in the second-derivative image-processing scheme that converts the input N -image into the second-derivative N″ -image.  3.1.Input N-ImageWe have used the second-derivative method to analyze images collected with a frequency-domain instrument designed and built by Siemens Medical Engineering6 7 (Fig. 1). The instrument consists of four laser diodes (690, 750, 788, and 856 nm) used as optical sources modulated at a frequency of ∼70 MHz. The source and detector fibers are located on opposite sides of the breast, which is slightly compressed between two parallel glass plates. The fibers are scanned in tandem across the breast to yield two-dimensional projection images of the phase and amplitude of the intensity-modulated light. Two projections of each breast were typically acquired—craniocaudal (the geometry illustrated in Fig. 1) and oblique—and obtained by rotating the glass plates by 45 deg. The raw data at each wavelength were converted into an N -image according to Eq. (1). It is necessary to remove edge effects prior to processing with the second derivative in order to avoid the introduction of artifacts near the breast edge. 3.2.Convolve with Smoothing FunctionBefore taking the second derivative, we smoothed the N -image with the 5×5 weighted average matrix shown in Fig. 3 to reduce the effects of noise. Any significant noise spikes or anomalies in the data should be removed prior to smoothing, e.g., by employing traditional image-processing techniques such as median filtering.2 Our 5×5 weighted average function is the result of cascading two 3×3 uniformly weighted average functions. We found that applying a 3×3 uniformly weighted average function once was not sufficiently effective in avoiding noise spikes in the second-derivative images. The amount of smoothing necessary will, of course, depend on the characteristics of the input data. In our case, the pixel size in the input N -image was 2×2 mm, and smoothing over a 1-cm 2 area worked well. 3.3.Calculate Second Derivative in Four DirectionsWe calculated the second derivative of N at each pixel by considering the N values at its nearest neighbors in four directions: horizontal, vertical, and along the two diagonals. Including other intermediate directions in the calculation did not significantly alter the results. We have used a standard forward-difference discrete approximation of the second derivative, but for completeness we define our second-derivative operator in the horizontal direction as follows. Consider three pixels in an N -image that lie in the same row and denote them left (l), center (c), and right (r). The value of N at each pixel is Nl, Nc, and Nr, respectively. Let Δw represent the center-to-center pixel separation in the horizontal direction. Then the first derivatives (N′) to the left and right of the center pixel, respectively, are given by: and From Eqs. (2) and (3), it follows that the second derivative at the center pixel is The second derivative in the vertical and diagonal directions is defined analogously. In the case of the diagonal directions, the center-to-center pixel separation is (2Δw)1/2, so the horizontal and vertical second-derivative coefficients differ from the diagonal ones by a factor of 2. In the matrices of Fig. 3, Δw has been normalized to 1, so only the factor of 2 remains.When implementing the second-derivative calculation as a convolution with the functions of Fig. 3, one must take into account the transition between the image background (area outside the breast) and the breast. In our N -image, the background pixels contain zeros while the breast pixels contain positive N -values. So simply filtering an N -image (or a smoothed version of it) with the second-derivative functions of Fig. 3 yields artifacts along the breast edge whenever background pixels are included in the calculation. One way of correcting this problem is to use an edge-detection algorithm to flag pixels near the breast edge and write a convolution routine that excludes background pixels in this region. A simpler, albeit cruder, solution is to segment the image into inside-breast/outside-breast regions and set the outside-breast region to a value higher than the average inside-breast N -value during the second-derivative calculation. We employed the latter approach in our analysis. To a lesser degree, this same issue arises when convolving the N -image with the smoothing function of Sec. 3.2. To avoid artificially lowering the N -values along the breast edge, we set the outside-breast region to the average inside-breast N -value during this step. 3.4.Take MinimumWe took the minimum of the second derivative along the four directions and assigned that value to the corresponding pixel in the N″ -image. We used the minimum, as opposed to the average along the four directions, in order to avoid missing areas that have a negative second derivative in one direction but are relatively flat or concave in the other directions. This procedure enhances the visualization of directional structures such as blood vessels. 3.5.DisplayThe second-derivative images are normalized by introducing a uniform multiplicative factor that sets to −1 the pixel value corresponding to the most negative second derivative. The normalized second-derivative values are indicated with N norm ″. We set a threshold at zero so that pixels with N norm ″⩾0 are displayed in white. Pixels with a negative second derivative (−1⩽N norm ″<0) are displayed in gray scale, with darker areas corresponding to more negative. We enhanced the visibility of features of interest in the N norm ″ images by performing a linear contrast stretch2 on the gray levels of the image pixels with a negative second derivative. Automatically optimizing the contrast of an image to be interpreted by a human analyst is a topic with a vast literature and a research undertaking outside the scope of the present paper. We did not attempt to generate new algorithms to address this problem; rather, we sought a robust method that could be readily applied to our images to automatically tweak their contrast. Histogram equalization is a popular image-processing technique for contrast manipulation, but its results are often too severe,8 and we found this to be the case with our images as well. So we settled on a linear contrast stretch, which preserves the shape of the image histogram but stretches it between two predetermined points on the display scale. These points are typically selected by the image analyst interactively. We automated the process by basing our choice on the characteristics of the input images. The contrast stretch procedure we employed consists of linearly scaling the gray levels for (100−x) of the pixels with a negative second derivative, where x of the pixels, which are those with the most negative second derivative, are set to the darkest gray level. The percentage x of the pixels excluded from the contrast stretch is given by: where σ2 is the variance of N norm ″ for the pixels with a negative N norm ″. For example, consider an N norm ″ -image with a variance such that x=7, which we wish to display in 256 gray levels (0 to 255). Then the 7 of the negative N norm ″ pixels with the most negative values of N norm ″ would be set to 0 in the contrast-stretched image, and the gray levels of the remaining 93 of the negative N norm ″ pixels would be linearly scaled into the range 0 (black) to 255 (white). The lower and upper limits of 5 and 10, respectively, for x, along with the coefficients of 5747 and 161 in Eq. (5) were chosen empirically from a subset of images and then successfully applied to the rest of the images without the need for further adjustments. There are a number of ways of selecting these parameters for a given database of images. A relatively straightforward approach is to calculate σ2 for all of the images in the database and determine its distribution. Based on this knowledge, one could choose representative images with low and high values of σ2 and decide what level of contrast stretching (i.e., what value of x) works best in each case. These values would then become the end points of a line segment in a two-dimensional plot of σ2 versus x. The remaining images would have a contrast-stretch parameter x, determined by σ2, that falls somewhere along this line. We found that modeling x as a linear function of σ2 worked well for our purposes, but other functions might yield even better results.The contrast stretch we used is asymmetric and therefore somewhat unusual. However, in our case this asymmetry is justified since we lose little by making the darkest features in our images black, while we would be compromising the detectability of the dimmest structures if we made them lighter. With this process, we are in effect using the variance of N norm ″ for those pixels with a negative second derivative to set a new threshold; the pixels whose N norm ″ values are below this threshold (most negative) are set to 0 (black), and those whose N norm ″ values are above it (less negative) are linearly scaled into the range 0 to 255. The motivation for using the variance of N norm ″ at the negative N norm ″ pixels to determine the amount of contrast stretching is simple. Typically, images with the largest variance are those that possess a dark dominant feature and many light areas. Setting a higher threshold in this situation makes the lighter features that much darker, thereby potentially revealing other structures in addition to the most salient one. 4.ResultsWe present four cases that illustrate the advantages of second-derivative images over N -images. On the left side of each figure is the N -image, with the craniocaudal projection on the top and the oblique projection on the bottom. On the right side, the figures show the corresponding craniocaudal and oblique projections of the N″ -image. Arrows indicate the location of the lesion, which we know from X-ray mammography, and whether it is cancerous or a benign mastopathy, which we know from the pathology report. 4.1.Case 197, 2.5-cm CancerPatient 197 is a 72-year-old woman with a 2.5-cm tumor (invasive ductal carcinoma) in the left breast. In this case, the tumor is clearly visible in both projections of the N -image, but notice in Fig. 4 the additional structural information, most likely due to the vasculature, that emerges in the second-derivative image. Some of the blood vessels were faintly visible in the N -image, but they appear as high-contrast features in the N″ -image. We also see that the areas near the edge of the breast in the N″ -image are mostly white, i.e., they correspond to a positive second derivative, so it appears that taking the second derivative does not result in the introduction of artifacts in this region. We have found this to be the case with all N″ -images that we have examined. 4.2.Case 267, 2-cm CancerThe case of patient 267, a 76-year-old woman with an invasive ductal carcinoma in the right breast, is an example of a low-contrast N -image without much visible structure. The tumor is barely discernible in the craniocaudal projection of the N -image [Fig. 5(a)] and cannot be easily differentiated from other dark areas in the image. In the oblique projection [Fig. 5(b)], the tumor is altogether eclipsed by another feature (circle). The N″ -image, however, clearly reveals the location of the tumor in the craniocaudal projection [Fig. 5(c)], while the other dark areas in the N -image now appear to be due to blood vessels. Although the oblique projection [Fig. 5(d)] continues to be dominated by the circled structure, the area corresponding to the location of the tumor is not white (i.e., it has a negative second derivative) and can therefore be distinguished from the background. 4.3.Case 165, MastopathyPatient 165, a 52-year-old woman, has a benign lesion in the right breast. This feature dominates the N -image and as a result not much else is visible. But in the N″ -image, we see a finer structure emerge that is most likely due to blood vessels. In particular, in the oblique projection, we can see that what appeared to be a rather murky and undifferentiated feature in the N -image shows up as a combination of two resolved structures in the N″ -image. 4.4.Case 215, 3-cm CancerThe case of patient 215, a 53-year-old woman with an invasive ductal carcinoma in the left breast, is particularly interesting. In the craniocaudal projection of the N -image [Fig. 7(a)], the cancer is readily apparent, as is another structure: a blood vessel. But even in this situation there is a substantial improvement in the contrast and level of detail visible in the N″ -image [Fig. 7(c)]. For example, in the N″ -image we can now make out what seems to be a network of blood vessels at the tumor site. The N -image oblique projection [Fig. 7(b)] has a lower contrast, and the improvement in the N″ -image [Fig. 7(d)] is therefore more pronounced. The visibility of the tumor in this projection is much greater in the N″ -image than in its counterpart. 5.DiscussionFrom the examples presented, it is clear that applying the second-derivative method enhances the structural information content of N -images. The N″ -images reveal attributes that were just barely visible in N -images or, in some cases, not really visible at all. In those images that are dominated by one feature (Fig. 4 and Fig. 6), we are now able to see additional structure away from the dominant feature and to better resolve finer details of the dominant feature. In cases where the N -image has a low contrast [Figs. 5(a) and 7(b)], or when the lesion is not the dominant structure [Fig. 5(b)], applying this method can lead to the detection of a previously missed tumor. Moreover, since we normalize the N″ -images and select the linear contrast stretch parameters once for all images, the N″ -images can be generated without any user input. The detectability of a particular lesion does not depend on a fortunate choice of threshold or contrast level. Figure 6N -images at 690 nm of the right breast of 52-year-old patient 165 in (a) craniocaudal and (b) oblique projection, and corresponding second-derivative images in (c) and (d). The arrows indicate the location of a benign mastopathy.  Figure 7N -images at 690 nm of the left breast of 53-year-old patient 215 in (a) craniocaudal and (b) oblique projection, and corresponding second-derivative images in (c) and (d). The arrows indicate the location of a 3.0-cm invasive ductal carcinoma.  As Fig. 2 illustrates, the second derivative is more sensitive to smaller departures from the local background absorbance. Furthermore, it is insensitive to the breast edge. The net result is a process capable of bringing out “buried” features [like the peak at x=8.8 cm in Fig. 2(a)] without introducing artifacts close to the breast edge. The greater detail visible in N″ -images improves the sensitivity of tumor detection, but this may come at the price of a reduced specificity, and Fig. 7(d) brings home this point nicely. Were it not for the arrow indicating the location of the tumor in this image, one would be hard pressed to argue that this structure is a lesion and not just part of a blood vessel. Structural information alone will not allow us to properly classify the various features visible in an N″ -image. For this we must make use of data at multiple wavelengths and exploit the potential of optical mammography to provide functional information (in the form of a measurement of oxygenation) from spectral information. There are several groups currently working to develop an accurate way of measuring the oxygenation level of breast lesions in vivo by combining data from two or more wavelengths.9 10 11 12 13 14 Our second-derivative method can play an important role in the interpretation of optical mammograms by robustly selecting areas of interest to be analyzed spectrally in order to discriminate benign from malignant lesions. The second-derivative operator is an ideal choice for this prefiltering step. By setting a threshold at zero and considering only pixels with a negative second derivative, one is, by definition, picking out any area that represents a local absorption maximum, irrespective of its magnitude, size, or shape. It is for this reason that we have elected to use the minimum value of N″ instead of the average along the four directions (Sec. 3.4). Although this tends to emphasize narrow, elongated features (e.g., blood vessels), we felt it best to err on the side of being more inclusive and not to rely on the assumption that breast lesions would be largely symmetrical. The second derivative is analytically straightforward and easy to implement; it is the driving mechanism behind many edge-detection algorithms in use today. Finally, we should comment on the issues of depth discrimination and resolution. From a planar projection, it is not possible to determine the depth of structures visible in N″ -images. But because these images possess a relatively high contrast, they may be particularly well suited for multisource and/or multidetector approaches to depth discrimination.15 16 As for resolution, we should point out that the second-derivative method does not improve the quality of the data available. It does, however, improve the display of the data, thereby allowing the discrimination of nearby structures that are not visibly resolved in the original images. 6.ConclusionWe have developed an image-processing method based on a spatial second-derivative operator to enhance the detection of breast tumors in optical mammograms. Applying this method has revealed the presence of breast tumors and of fine structures, most likely blood vessels, that were not visible in the original images. As a result, this second-derivative method enhances the structural imaging capability of optical mammography, thus leading to a potentially higher sensitivity in the detection of breast cancer. In fact, in addition to the direct visualization of breast tumors, a modified spatial pattern of the vascularization may also provide an indirect indication of the presence of cancer. Furthermore, the identification of suspicious areas of interest provided by second-derivative images may guide the selection of image pixels for a spectral analysis aimed at the assessment of the local oxygen saturation of hemoglobin. This second step, which provides functional information, may improve the specificity of optical mammography by possibly allowing the discrimination of benign and malignant breast lesions. AcknowledgmentsThis work was sponsored by the Department of the Army award No. DAMD17-99-1-9218 and by National Science Foundation award No. BES-93840. REFERENCES
P. Vaupel
,
F. Kallinioski
, and
P. Okunieff
,
“Blood flow, oxygen and nutrient supply, and metabolic microenvironment of human tumors: a review,”
Cancer Res. , 49 6449
–6465
(1989). Google Scholar
S. Fantini
,
M. A. Franceschini
,
G. Gaida
,
E. Gratton
,
H. Jess
,
W. W. Mantulin
,
K. T. Moesta
,
P. M. Schlag
, and
M. Kaschke
,
“Frequency-domain optical mammography: edge-effect corrections,”
Med. Phys. , 23 149
–157
(1996). Google Scholar
M. Kaschke
,
H. Jess
,
G. Gaida
,
J. M. Kaltenbach
, and
W. Wrobel
,
“Transillumination imaging of tissue by phase modulation techniques,” in Advances in Optical Imaging and Photon Migration, R. R. Alfano, Ed.,”
Proc. OSA , 21 88
–92
(1994). Google Scholar
M. A. Franceschini
,
K. T. Moesta
,
S. Fantini
,
G. Gaida
,
E. Gratton
,
H. Jess
,
W. W. Mantulin
,
M. Seeber
,
P. M. Schlag
, and
M. Kaschke
,
“Frequency-domain techniques enhance optical mammography: initial clinical results,”
Proc. Natl. Acad. Sci. U.S.A. , 94 6468
–6473
(1997). Google Scholar
L. Go¨tz
,
S. H. Heywang-Ko¨brunner
,
O. Schu¨tz
, and
H. Siebold
,
“Optische mammographie an pra¨operativen patientinnen,”
Akt. Radiol. , 8 31
–33
(1998). Google Scholar
S. Fantini, E. L. Heffer, M. A. Franceschini, L. Go¨tz, A. Heinig, S. Heywang-Ko¨brunner, Oliver Schu¨tz, and Horst Siebold, “Optical mammography with intensity-modulated light,” Proceedings of Inter-Institute Workshop on In Vivo Optical Imaging at the NIH, A. H. Gandjbakhche, Ed., pp. 111–117, Optical Society of America, Washington, DC (2000).
J. A. Stark
,
“Adaptive image contrast enhancement using generalizations of histogram equalization,”
IEEE Trans. Image Process. , 9
(5), 889
–896
(2000). Google Scholar
S. Zhou
,
C. Xie
,
S. Nioka
,
H. Liu
,
Y. Zhang
, and
B. Chance
,
“Phased array instrumentation appropriate to high precision detection and localization of breast tumor,” in Optical Tomography and Spectroscopy of Tissue: Theory, Instrumentation, Model, and Human Studies II, B. Chance and R. R. Alfano, Eds.,”
Proc. SPIE , 2979 98
–106
(1997). Google Scholar
B. W. Pogue
,
M. Testorf
,
T. McBride
,
U. Osterberg
, and
K. Paulsen
,
“Instrumentation and design of a frequency-domain diffuse optical tomography imager for breast cancer detection,”
Opt. Express , 1 391
–403
(1997). Google Scholar
D. Grosenick
,
H. Wabnitz
,
H. H. Rinneberg
,
K. T. Moesta
, and
P. M. Schlag
,
“Development of a time-domain optical mammography and first in vivo applications,”
Appl. Opt. , 38 2927
–2943
(1999). Google Scholar
V. Ntziachristos
,
A. G. Yodh
,
M. Schnall
, and
B. Chance
,
“Concurrent MRI and diffuse optical tomography of breast after indocyanine green enhancement,”
Proc. Natl. Acad. Sci. U.S.A. , 97 2767
–2772
(2000). Google Scholar
N. Shah
,
A. Cerussi
,
C. Eker
,
J. Espinoza
,
J. Butler
,
J. Fishkin
,
R. Hornung
, and
B. Tromberg
,
“Noninvasive functional optical spectroscopy of human breast tissue,”
Proc. Natl. Acad. Sci. U.S.A. , 98 4420
–4425
(2001). Google Scholar
E. L. Heffer
and
S. Fantini
,
“Quantitative oximetry of breast tumors: A near-infrared method that identifies two optimal wavelengths for each tumor,”
Appl. Opt. , 41 3827
–3839
(2002). Google Scholar
X. D. Li
,
T. Durduran
,
A. G. Yodh
,
B. Chance
, and
D. N. Pattanayak
,
“Diffraction tomography for biochemical imaging with diffuse-photon density waves,”
Opt. Lett. , 22 573
–575
(1997). Google Scholar
|
CITATIONS
Cited by 51 scholarly publications.
Breast
Tumors
Mammography
Image enhancement
Blood vessels
Image processing
Absorption
