- Email: [email protected]
Hopfield neural network for the multichannel segmentation of magnetic resonance cerebral images
Pattern Recognition, Vol. 30, No. 6, pp. 921-927, 1997 ~ 1997 Pattern Recognition Society. Published by Elsevier Science Ltd Printed in Great Britain. All rights reserved 0031-3203/97 $17.00+.00
Pergamon
PII: S0031-3203(96)00095-7
HOPFIELD NEURAL NETWORK FOR THE MULTICHANNEL SEGMENTATION OF MAGNETIC RESONANCE CEREBRAL IMAGES RACHID S A M M O U D A t'*, NOBORU NIKI t and HIROMU NISHITANF tDepartment of Optical Science and Technology, University of Tokushima, Minami-josanjima-cho 2-1, Tokushima, Japan ¢Medical School, University of Tokushima, Tokushima, Japan
(Received 8 August 1995; in revised form 15 May 1996: reeeivedJor publication 10 July 1996)
Abstract--In this paper, we present an approach for the segmentation of magnetic resonance images of the brain, based on Hopfield neural network. We formulate the segmentation problem as a minimization of an energy function constructed with two terms, the cost-term, that is a sum of errors' squares, and the second term is a temporary noise added to the cost-term as an excitation to the network to escape from certain local minima and be closer to the global minimum. Also, to ensure the convergence of the network and its utility in the clinic with useful results, the minimization is achieved in a way that after a prespecified period of time the energy function can reach a local minimum close to the global minimum and remain there ever after. We present here, segmentation results of two patients data diagnosed with a metastatic tumor and multiples sclerosis in the brain. ~L' 1997 Pattern Recognition Society. Published by Elsevier Science Ltd. Magnetic resonance images Segmentation Artificial neural networks Optimization Local minimum Global minimum
1. INTRODUCTION Magnetic resonance imaging (MRI) is a multidimensional technique as it provides information about three tissue-dependent parameters: T1 (the spin-lattice relaxation time), T2 (the spin-spin relaxation time) and the proton density. Automating the analysis process of the images of brain obtained using the MRI technique, saves human time and effort and can provide accurate and reproducible results. However, the segmentation of the MR images which is the first step toward a quantitative analysis, is still far from being totally automated because of several artifacts such as: 1. intensity variations introduced by the ratio frequency (RF) coil, 2. imager imposed inter-patient variations, 3. inter-slice intensity variations. To overcome these difficulties of the segmentation process, a more realistic approach is to be provided the user with an integrated computer environment that offers rich display, processing, and classification capabilities to assist the human in his/her task analysis. These environments should be interactive to let the expert controls the results of the classification, but they should reduce as much as possible the amount of interaction required to reach the desired results. A number of studies have explored the potential of various supervised or unsuper-
* Author to whom correspondence should be addressed. Tel.: (81) 886-56-7491; fax: (81) 886-56-7492; e-mail: sammouda @is.tokushima-u.ac.jp. 921
vised pattern recognition techniques for the segmentation of MRI data (Vannier(1), Herskovits (m, Kohn!3)). Typically, supervised classifiers are based on the estimation of a probability density function and a common assumption is that the data can be modeled by a multivariate normal distribution. In recent years, artificial neural networks (ANNs) have been proposed as an attractive alternative solution to a number of pattem recognition problems. Their major advantage is that they do not rely on any assumption about the underlying probability density functions, thus possibly improving the results when the data significantly deviate from normal. Hopfield neural network (HNN) is one of the ANNs which has been used by Armatur C4)for the classification of MR images, based on energy minimization. Generally the results have been found to be acceptable. However, it has been reported ~4) that minimizing the constructed energy function can lead to get trapped in an early local minima in the energy landscape, and the segmentation's resultant images need a post classification filtering which can cost the loss of some small structures of the images. Furthermore, the convergence of the network is data dependent as it needs a preset threshold to detect the stability of the network's inputs. Such limitations decrease the experts confidence in taking the results as a basis for their analysis and diagnosis, and thus makes the segmentation process unpractical in the clinic considering the huge amount of data and time restrictions. To be more practical and closer to real-life time, we present in this paper a contribution of two main wrinkles to the previous work. (4)
R. SAMMOUDA et al.
922
1. Ensure the convergence of the network in prespecifled time. 2. A procedure for performing the hill-climbing required to avoid being trapped in certain local minima and be closer to the global minimum. We show the effectiveness of our contribution by comparing our results to those of the previous work.
2. R E V I E W O F P R E V I O U S W O R K
The segmentation problem is formulated by Armatur (4) as a partition of N pixels of P features among M classes such that the assignment of the pixels minimizes the criterion function: N
1
M
E : ~RklV;i,____ k:l
2
(1)
/:1
where Rkt is considered as the Mahalanobis distance measure between the kth pixel and the centroid of class l, and Vkt is the output of the klth neuron. The minimization is achieved by using HNN and by solving a set of motion equations satisfying:
OUi OE Ot - OVi'
(2)
where Ui and Vi are respectively the input and the output of the ith neuron. Using the maximum neuron rule as input-output function, and if the system satisfies equation (2), the energy function continuously reduces as a function of time and the system converges to a local minimum. However, by implementing the algorithm proposed by Armatur (4), we found that the total energy of the network is not monotically decreasing, and that is because the variation of the input of each neuron not only depends on the variation of the energy with respect to the corresponding neuron's output, but also depends on all the other neurons states. Due to this, the network revisits some previous states and the convergence time becomes longer.
3. O U R C O N T R I B U T I O N S
Our first contribution concems the formulation of the energy function of the segmentation problem. Here, we consider the classification as a partition of N pixels of P features into the best M classes, and we propose the following equation as a cost-term of the energy function: 1
Ec
2=
N
M
/=1
In the case n = l , the cost-term is reduced to equation (1) which is the energy function described by mrmatur (4). In case n=2, the cost-term is the sum of the squares errors (~°), which leads to a better solution as we will show in Section 4. As a second contribution and to increase the convergence speed of the network, the minimization is
achieved by:
OUi Ot --
OE~. #(t) o ~ i ,
(4)
where #(t) is as defined by Jacobs (s), a scalar positive function of time which determines the length of the step to be taken in the direction of the vector d = -~TE,.(v). Furthermore, it is possible to find a convenient step #(t) in order to ensure that the network can reach a local minimum after a specified period of time, say Ts, and remains there ever after. (6'7) The appropriate choice or selection of the step #(t) is some thing of an art. Experimentation (or trial and error) and a familiarity with a given class of optimization problems are often required to find the best function. In our study we found that the function #(t) - t* (Ts t) is convenient for the segmentation of MR data using HNN, where T~ is the prespecified time. The last contribution is how to avoid the network being trapped in an early local minimum. To escape from the local minima, the most widely used methods for global optimization are of a stochastic nature. In these methods, random fluctuations or noise are introduced into the system in order to avoid being trapped in local minima. However, noise in ANNs often seems ubiquitous and apparently useless. On the other hand, by investigating neural cells of living creatures, it has been found that they are essentially stochastic (random) in nature, in that responses of individual neural cells (in the isolated cortex) due to cyclically repeated identical stimuli will never result in identical responses. This implies that noise plays an important role in neural networks. Similarly, in ANNs it has been found that additive noise introduced into the network could be very useful since it increases the probability of convergence to the global minimum. However, in order to get the exact solution this additive noise should be introduced appropriately and then subsequently removed from the network as the global minimum is approached, as the noise destroys valuable information. (s'9) In other words, during the optimization process the additive noise (perturbations) should approach zero in time so that the network itself will become deterministic prior to reaching the final solution. Adding a noise term to the cost-term, the energy function of our system will be described as follows: N
E,. + c(t) Z
M
Z NktV~,,
(5)
k=l /=1
where ~? is called the perturbed form of the energy function, and Nkt is a N x M vector of independent high-frequency white noise sources, c(t) is the parameter controlling the magnitude of noise which must be selected in such away that it provides zero as time t tends to "infinity". The parameter c(t) is usually monotonically decreasing in time. However, in some cases a more complex time schedule can be incorporated in which c(t) can occasionally be increased if the network is trapped in a local minimum. A convenient choice in
Hopfield neural network for the multichannel segmentation our case for the controlling of the parameter c(t) is an exponential schedule which can be expressed mathematically as
c(t) = t~e -c~t,
(6)
where the coefficient/3>0 controls the initial amplitude of the noise at (t=0) and c~>0 determines the rate of noise damping. For a large coefficient c~, the system evolves very rapidly to become deterministic while a low c~ allows the system to revisit chaotic states as needed usually during the optimization process. On the other hand, for a very low value of c~ the convergence process will be very slow. Then we should balance the desire for a possibly high speed of convergence (i.e. large c0 with the need to ensure a solution of good quality (i.e. near to the global minimum).
4. IMPLEMENTATION AND RESULTS
In order to show the effects of our contributions, we have applied the algorithm described by Armatur ¢4)
(a)
923
to the equations (3-5). The segmentation results are given below in comparison with the results of previous work. Figure 1 (a)--(c) shows respectively the TI-, T2-, and the proton density-weighted images of a patient diagnosed with a metastatic tumor in the brain. The T2weighted image shows an abnormal bright region in the white matter. We use these images as an input to the classifier network with a prespecified convergence time T~, for two and three channel segmentation cases. The time Ts defined in the step function #(t), is chosen experimentally as that corresponding to the optimal energy value. Figure 2 shows three logarithmic curves of the energy function of the two channel segmentation using the T1- and the T2-weighted images based on our approach with three different T~ values (150, 200, 230). In this case Ts=200 iterations corresponds to the best solution found. Figure 3 shows the logarithmic curves of the energy function, of the three channel segmentation (using the TI-, T2-, and proton density-weighted images) obtained based on our approach with three T~ values
(b)
(c)
Fig. 1. Magnetic resonance head images of a patient diagnosed with a metastatic tumor in the brain (a) T1weighted (b) T2-weighted and (c) the proton density-weighted.
E(O
E(0 3500
7000 A
6000
B
C
........ .........
A
3000
5OOO
m
13
........
C
.........
2500
aOOO 2000 3000 t500 2000
~~'~:::~.::;............
I000
1000
"::::::I~L~I~ .............
, .
500 i
0 0
50
i
I
100 150 Iteration (t)
i
200
250
Fig. 2. Three logarithmic curves of the energy function of the two channel segmentation using the T1- and the T2-weighted images based on our approach with three T~ values: A(150), B(200) and C(230).
C
A 0
~
1
I
50
i
Ito 150 Iteration (0
i
200
250
Fig. 3. Three logarithmic curves of the energy function of the three channel segmentation based on our approach with three T~ values: A(70), B(120) and C(200).
R. SAMMOUDA et al.
924
Fig. 4. Two channel segmentation results with seven classes using the T1- and the T2-weighted images. A, B, C, and D are obtained using the different approaches described in the text.
(70, 120, 200). Ts=120 corresponds to the best solution found. For the two and three channel segmentation, the parameter c(t) is defined as/3 = ko/(1 + ,4) and a=0.5, with k0 a positive constant and A the difference between the two last energy values of iteration. Figure 4 shows the two channel segmentation results with seven classes using the Tl-weighted and T2weighted images, with the color maps showing the distribution of the pixels in the segmented images along their values in the TI- and T2-weighted images. The result A is obtained using the energy function described by Armatur (4~, which corresponds to ( n = l ) in our proposed cost-term equation. B is obtained based on the energy function as a sum of the cost-term with ( n = l ) and the additional noise term. C is obtained when the energy function is equal to the cost-term with (n=2), which corresponds to the sum of the squares errors. D is obtained using equation (5) including the cost term E,, with (n=2) as an energy function of the classification problem, and using equation (4) for the minimization process. Figure 5 shows the logarithmic curves of the energy function corresponding respectively to the results presented in Fig. 4. Figure 6 shows the three channel segmentation results with seven classes. The results A, B, C, and D are respectively obtained with the approaches described above for Figs 4 and 7 show the logarithmic curves of the energy function corresponding to the results of Fig. 6.
5. DISCUSSION
Comparing the result A obtained by the previous work, and the result B in Fig. 4, we can say that adding the noise term to the system can help the network to give
E(O 7000 ...... "'¢.,.,, \
6000
A 13 .......
5000 4000
': ": "
"'..;".,
A
3000
"""
2000
.
B
- • 2:,-
~::~:":::::~:-ii;);i./~-i .
• ............
l ooo
"
C
:ii: iii.ii ............. D
0
i
0
50
i
i
100 150 Iteration (t)
200
250
Fig. 5. The logarithmic curves of the energy function of the two channel segmentation corresponding respectively to the results presented in Fig. 4.
clearer results by going more deeply into the energy landscape, but it can not ensure the optimal result. However, by considering the energy function of the segmentation problem as the sum of the squares errors and using the algorithm (4), we get the result C which presents more clarity and smoothness than the A and B. The difference between the results C and D in Fig. 4 is that the tumor region colored red in D is not detected in C and is colored yellow, the same as the cerebrospinal fluid (csf). However, comparing the segmentation results to the raw images, and considering the energy value of the network at equilibrium, the result D is considered as the best.
Hopfield neural network for the multichannel segmentation
925
Fig. 6. The three channel segmentation results with seven classes+ A, B, C, and D are obtained respectively with the same approaches described for Fig. 4.
ECt) 3500
i I¸I:
3000
2500
i
2000
'~ ~'
A
1500 •
~
B
1000
::::--
5oo
...................
C
:~::::211~I:IILY21~L:IY::I:I::::2.:~::::: D
0
50
100 150 Iteration (t)
200
Fig. 7. The logarithmic curves of the energy function corresponding respectively to the results of Fig. 6.
In the case of the three channel segmentation, the same remarks can be said with the addition that the energy function is more optimal and the convergence speed is faster than in the case of two channel segmentation. Such behavior is similar to that of the human brain, if we have more information about a problem then we can understand it faster and more clearly. Also, a very important step for the automation process is to fix the convenient convergence time of the network to ensure the optimal result. However, here we fix the class number to seven for the presented case based on some medical information and also on the smoothness of the regions that present the results. In fact the energy value at equilibrium decreases if we increase the class number, but there are not any smooth regions in the segmented image. In fact, this is due to the intensity variation of the data which may have to be smoothed carefully without blurring the image before being used for the segmentation process.
I
(a)
(b)
(c)
Fig. 8, Magnetic resonance head images of a patient diagnosed with a multiple sclerosis in the brain (a) T1weighted (b) T2-weighted and (c) the proton density-weighted.
926
R. SAMMOUDA et al.
Fig. 9. Two channel segmentation results with eight classes using the Tl-weighed and the T2-weighted images obtained with our approach.
Fig. 10. Three channel segmentation results with eight classes obtained with our approach.
To verify the generality of our approach and its parameters on the data used, we operated the program with different subject data. We found that the parameters values of #(t), T~, c~ and/3 used before are also convenient for our database which contains images of several normal and abnormal cases. Figure 8 (a)-(c) shows respectively the Tl-weighted, T2-weighted, and proton densityweighted images of a patient diagnosed with a multiple sclerosis in the white matter. Also, as in the first case, the T2-weighted image shows an abnormal bright region in the white matter. We have used our approach with the same parameters described above, for two and three channel segmentation of these images. Figure 9 shows the two channel segmentation result with eight classes, using the Tl-weighted and the T2-weighted images. The red regions in the segmented image represent the abnormal regions. The color map shows the distribution of the pixels in the segmented image along their values in the T1- and T2- weighted images. Figure 10 shows the three channel segmentation result with eight classes. In this case, we have used eight classes because of the high intensity variations present in the raw images.
6. CONCLUSION In this study, we presented some contributions to a previous work in the use of Hopfield network for the segmentation of MR images. We showed that the performance of the network is mainly dictated by the choice of the energy function which can lead to a good segmentation when it is based on the sum of the squares
errors as a cost-term, and the addition of a noise term as an excitation to the network to detect and escape from certain local minima and be closer to the global minimum. Also, the use of the network with a prespecified time convergence promises a full automatic segmentation in the future if we can determine automatically the class number of the data during the segmentation process. In our current work, we focus on a preclassification process to uniform carefully the intensity values in the raw images without blurring the edges between the different regions such that we can determine the class number automatically during the segmentation process.
REFERENCES
I. M. W. Vannier, R. L. Butterfield, D. L. Rickman, D. ZM. Jordan, W. A. Murphy, R. G. Levit and G. Mohktar, Multispectral analysis of magnetic resonance images, Radiology 154(1), 221-224 (1995). 2. M. I. Kohn, N. K. Tanna and G. T. Herman, Analysis brain and cerebrospinal fluid volumes with MR imaging: Part 1 Methods, reliability and validation, Radiology 178, 115-122 (1991). 3. E. Herkovits, A hybrid classifier for automated radiological diagnosis: Preliminary results and clinical applications, Comput. Methods and Programs in Biomedicine 32(1), 45-52 (1990). 4. S. C. Armatur, D. Piraino and Y. Takefuji, Optimization neural networks for the segmentation of magnetic resonance images, IEEE Trans. Medical Imaging 11(2), 215-220 (1992). 5. R.A. Jacobs, Increased rates of convergence through learning rate adoption, Neural Nem,orks 1, 295-307 (1988). 6. Z. R. Novakovic, Solving systems of non-linear equations using the Lyapunor direct method, Computer Math. Apllic. 20(12), 19-23 (1990). 7. S. Schaffier and H. Warsitz, A trajectory-following method for unconstrained optimization, J. Optimization Theory and Applications 67(1), 133-138 (1990). 8. E Aluffi-pentini, V. Parsi and F. Zirilli, Global optimization and stochastic differential equations, J. Optimization Theory and Applications 47, 1-16 (1985). 9. S. Geman and C.R. Hwang, Diffusions for global optimization, SIAM J. Control and optimization 24(5), 1031 1043 (1986). 10. R. Sammouda, N. Niki and H. Nishitani, Multichannel segmentation of magnetic resonance cerebral images based on neural networks, IEEE ICIP-95, Washington, pp. 484487, October (1995) (1986).
Hopfield neural network for the multichannel segmentation About the Author - - RACHID SAMMOUDA received his B.S. degree in Computer Science from the University of Tunisia in 1991 and the M.S. degree in Information Science and Intelligent Systems from the University of Tokushima in Japan in 1994. He has been a Ph.D. student at the University of Tokushima until March 1996, when he started to teach as Assistant Professor at the Department of Optical Science there. His current interest are medical images analysis and artificial intelligence.
About the A u t h o r - - N O B O R U NIKI received his B.S. and M.S. in Electronic Engineering from Tokushima University in Japan in 1975 and 1977, respectively. Also, he received a Ph.D. degree in Electrical Engineering from Kyoto University in Japan in 1987. He was a Research Assistant in Information Science Department at Tokushima University from 1977 to 1989. In 1989, he became an Associate Professor in same department. Since April 1996 he has been a Professor at the Department of Optical Science and Technology. His current interests are image reconstruction algorithms, medical graphics, pattern recognition, and parallel image processing.
About the A u t h o r - - HIROMU NISHITANI received his M.D. degree from the University of Kyushu in Japan
in 1970, and Ph.D. Medical Science from the University of Kyushu 1984. He has been Professor and chairman of the Department of Radiology, School of Medicine, the University of Tokushima in Japan since 1988. His main interests are diagnostic radiology and medical image analysis.
927
Pergamon
PII: S0031-3203(96)00095-7
HOPFIELD NEURAL NETWORK FOR THE MULTICHANNEL SEGMENTATION OF MAGNETIC RESONANCE CEREBRAL IMAGES RACHID S A M M O U D A t'*, NOBORU NIKI t and HIROMU NISHITANF tDepartment of Optical Science and Technology, University of Tokushima, Minami-josanjima-cho 2-1, Tokushima, Japan ¢Medical School, University of Tokushima, Tokushima, Japan
(Received 8 August 1995; in revised form 15 May 1996: reeeivedJor publication 10 July 1996)
Abstract--In this paper, we present an approach for the segmentation of magnetic resonance images of the brain, based on Hopfield neural network. We formulate the segmentation problem as a minimization of an energy function constructed with two terms, the cost-term, that is a sum of errors' squares, and the second term is a temporary noise added to the cost-term as an excitation to the network to escape from certain local minima and be closer to the global minimum. Also, to ensure the convergence of the network and its utility in the clinic with useful results, the minimization is achieved in a way that after a prespecified period of time the energy function can reach a local minimum close to the global minimum and remain there ever after. We present here, segmentation results of two patients data diagnosed with a metastatic tumor and multiples sclerosis in the brain. ~L' 1997 Pattern Recognition Society. Published by Elsevier Science Ltd. Magnetic resonance images Segmentation Artificial neural networks Optimization Local minimum Global minimum
1. INTRODUCTION Magnetic resonance imaging (MRI) is a multidimensional technique as it provides information about three tissue-dependent parameters: T1 (the spin-lattice relaxation time), T2 (the spin-spin relaxation time) and the proton density. Automating the analysis process of the images of brain obtained using the MRI technique, saves human time and effort and can provide accurate and reproducible results. However, the segmentation of the MR images which is the first step toward a quantitative analysis, is still far from being totally automated because of several artifacts such as: 1. intensity variations introduced by the ratio frequency (RF) coil, 2. imager imposed inter-patient variations, 3. inter-slice intensity variations. To overcome these difficulties of the segmentation process, a more realistic approach is to be provided the user with an integrated computer environment that offers rich display, processing, and classification capabilities to assist the human in his/her task analysis. These environments should be interactive to let the expert controls the results of the classification, but they should reduce as much as possible the amount of interaction required to reach the desired results. A number of studies have explored the potential of various supervised or unsuper-
* Author to whom correspondence should be addressed. Tel.: (81) 886-56-7491; fax: (81) 886-56-7492; e-mail: sammouda @is.tokushima-u.ac.jp. 921
vised pattern recognition techniques for the segmentation of MRI data (Vannier(1), Herskovits (m, Kohn!3)). Typically, supervised classifiers are based on the estimation of a probability density function and a common assumption is that the data can be modeled by a multivariate normal distribution. In recent years, artificial neural networks (ANNs) have been proposed as an attractive alternative solution to a number of pattem recognition problems. Their major advantage is that they do not rely on any assumption about the underlying probability density functions, thus possibly improving the results when the data significantly deviate from normal. Hopfield neural network (HNN) is one of the ANNs which has been used by Armatur C4)for the classification of MR images, based on energy minimization. Generally the results have been found to be acceptable. However, it has been reported ~4) that minimizing the constructed energy function can lead to get trapped in an early local minima in the energy landscape, and the segmentation's resultant images need a post classification filtering which can cost the loss of some small structures of the images. Furthermore, the convergence of the network is data dependent as it needs a preset threshold to detect the stability of the network's inputs. Such limitations decrease the experts confidence in taking the results as a basis for their analysis and diagnosis, and thus makes the segmentation process unpractical in the clinic considering the huge amount of data and time restrictions. To be more practical and closer to real-life time, we present in this paper a contribution of two main wrinkles to the previous work. (4)
R. SAMMOUDA et al.
922
1. Ensure the convergence of the network in prespecifled time. 2. A procedure for performing the hill-climbing required to avoid being trapped in certain local minima and be closer to the global minimum. We show the effectiveness of our contribution by comparing our results to those of the previous work.
2. R E V I E W O F P R E V I O U S W O R K
The segmentation problem is formulated by Armatur (4) as a partition of N pixels of P features among M classes such that the assignment of the pixels minimizes the criterion function: N
1
M
E : ~RklV;i,____ k:l
2
(1)
/:1
where Rkt is considered as the Mahalanobis distance measure between the kth pixel and the centroid of class l, and Vkt is the output of the klth neuron. The minimization is achieved by using HNN and by solving a set of motion equations satisfying:
OUi OE Ot - OVi'
(2)
where Ui and Vi are respectively the input and the output of the ith neuron. Using the maximum neuron rule as input-output function, and if the system satisfies equation (2), the energy function continuously reduces as a function of time and the system converges to a local minimum. However, by implementing the algorithm proposed by Armatur (4), we found that the total energy of the network is not monotically decreasing, and that is because the variation of the input of each neuron not only depends on the variation of the energy with respect to the corresponding neuron's output, but also depends on all the other neurons states. Due to this, the network revisits some previous states and the convergence time becomes longer.
3. O U R C O N T R I B U T I O N S
Our first contribution concems the formulation of the energy function of the segmentation problem. Here, we consider the classification as a partition of N pixels of P features into the best M classes, and we propose the following equation as a cost-term of the energy function: 1
Ec
2=
N
M
/=1
In the case n = l , the cost-term is reduced to equation (1) which is the energy function described by mrmatur (4). In case n=2, the cost-term is the sum of the squares errors (~°), which leads to a better solution as we will show in Section 4. As a second contribution and to increase the convergence speed of the network, the minimization is
achieved by:
OUi Ot --
OE~. #(t) o ~ i ,
(4)
where #(t) is as defined by Jacobs (s), a scalar positive function of time which determines the length of the step to be taken in the direction of the vector d = -~TE,.(v). Furthermore, it is possible to find a convenient step #(t) in order to ensure that the network can reach a local minimum after a specified period of time, say Ts, and remains there ever after. (6'7) The appropriate choice or selection of the step #(t) is some thing of an art. Experimentation (or trial and error) and a familiarity with a given class of optimization problems are often required to find the best function. In our study we found that the function #(t) - t* (Ts t) is convenient for the segmentation of MR data using HNN, where T~ is the prespecified time. The last contribution is how to avoid the network being trapped in an early local minimum. To escape from the local minima, the most widely used methods for global optimization are of a stochastic nature. In these methods, random fluctuations or noise are introduced into the system in order to avoid being trapped in local minima. However, noise in ANNs often seems ubiquitous and apparently useless. On the other hand, by investigating neural cells of living creatures, it has been found that they are essentially stochastic (random) in nature, in that responses of individual neural cells (in the isolated cortex) due to cyclically repeated identical stimuli will never result in identical responses. This implies that noise plays an important role in neural networks. Similarly, in ANNs it has been found that additive noise introduced into the network could be very useful since it increases the probability of convergence to the global minimum. However, in order to get the exact solution this additive noise should be introduced appropriately and then subsequently removed from the network as the global minimum is approached, as the noise destroys valuable information. (s'9) In other words, during the optimization process the additive noise (perturbations) should approach zero in time so that the network itself will become deterministic prior to reaching the final solution. Adding a noise term to the cost-term, the energy function of our system will be described as follows: N
E,. + c(t) Z
M
Z NktV~,,
(5)
k=l /=1
where ~? is called the perturbed form of the energy function, and Nkt is a N x M vector of independent high-frequency white noise sources, c(t) is the parameter controlling the magnitude of noise which must be selected in such away that it provides zero as time t tends to "infinity". The parameter c(t) is usually monotonically decreasing in time. However, in some cases a more complex time schedule can be incorporated in which c(t) can occasionally be increased if the network is trapped in a local minimum. A convenient choice in
Hopfield neural network for the multichannel segmentation our case for the controlling of the parameter c(t) is an exponential schedule which can be expressed mathematically as
c(t) = t~e -c~t,
(6)
where the coefficient/3>0 controls the initial amplitude of the noise at (t=0) and c~>0 determines the rate of noise damping. For a large coefficient c~, the system evolves very rapidly to become deterministic while a low c~ allows the system to revisit chaotic states as needed usually during the optimization process. On the other hand, for a very low value of c~ the convergence process will be very slow. Then we should balance the desire for a possibly high speed of convergence (i.e. large c0 with the need to ensure a solution of good quality (i.e. near to the global minimum).
4. IMPLEMENTATION AND RESULTS
In order to show the effects of our contributions, we have applied the algorithm described by Armatur ¢4)
(a)
923
to the equations (3-5). The segmentation results are given below in comparison with the results of previous work. Figure 1 (a)--(c) shows respectively the TI-, T2-, and the proton density-weighted images of a patient diagnosed with a metastatic tumor in the brain. The T2weighted image shows an abnormal bright region in the white matter. We use these images as an input to the classifier network with a prespecified convergence time T~, for two and three channel segmentation cases. The time Ts defined in the step function #(t), is chosen experimentally as that corresponding to the optimal energy value. Figure 2 shows three logarithmic curves of the energy function of the two channel segmentation using the T1- and the T2-weighted images based on our approach with three different T~ values (150, 200, 230). In this case Ts=200 iterations corresponds to the best solution found. Figure 3 shows the logarithmic curves of the energy function, of the three channel segmentation (using the TI-, T2-, and proton density-weighted images) obtained based on our approach with three T~ values
(b)
(c)
Fig. 1. Magnetic resonance head images of a patient diagnosed with a metastatic tumor in the brain (a) T1weighted (b) T2-weighted and (c) the proton density-weighted.
E(O
E(0 3500
7000 A
6000
B
C
........ .........
A
3000
5OOO
m
13
........
C
.........
2500
aOOO 2000 3000 t500 2000
~~'~:::~.::;............
I000
1000
"::::::I~L~I~ .............
, .
500 i
0 0
50
i
I
100 150 Iteration (t)
i
200
250
Fig. 2. Three logarithmic curves of the energy function of the two channel segmentation using the T1- and the T2-weighted images based on our approach with three T~ values: A(150), B(200) and C(230).
C
A 0
~
1
I
50
i
Ito 150 Iteration (0
i
200
250
Fig. 3. Three logarithmic curves of the energy function of the three channel segmentation based on our approach with three T~ values: A(70), B(120) and C(200).
R. SAMMOUDA et al.
924
Fig. 4. Two channel segmentation results with seven classes using the T1- and the T2-weighted images. A, B, C, and D are obtained using the different approaches described in the text.
(70, 120, 200). Ts=120 corresponds to the best solution found. For the two and three channel segmentation, the parameter c(t) is defined as/3 = ko/(1 + ,4) and a=0.5, with k0 a positive constant and A the difference between the two last energy values of iteration. Figure 4 shows the two channel segmentation results with seven classes using the Tl-weighted and T2weighted images, with the color maps showing the distribution of the pixels in the segmented images along their values in the TI- and T2-weighted images. The result A is obtained using the energy function described by Armatur (4~, which corresponds to ( n = l ) in our proposed cost-term equation. B is obtained based on the energy function as a sum of the cost-term with ( n = l ) and the additional noise term. C is obtained when the energy function is equal to the cost-term with (n=2), which corresponds to the sum of the squares errors. D is obtained using equation (5) including the cost term E,, with (n=2) as an energy function of the classification problem, and using equation (4) for the minimization process. Figure 5 shows the logarithmic curves of the energy function corresponding respectively to the results presented in Fig. 4. Figure 6 shows the three channel segmentation results with seven classes. The results A, B, C, and D are respectively obtained with the approaches described above for Figs 4 and 7 show the logarithmic curves of the energy function corresponding to the results of Fig. 6.
5. DISCUSSION
Comparing the result A obtained by the previous work, and the result B in Fig. 4, we can say that adding the noise term to the system can help the network to give
E(O 7000 ...... "'¢.,.,, \
6000
A 13 .......
5000 4000
': ": "
"'..;".,
A
3000
"""
2000
.
B
- • 2:,-
~::~:":::::~:-ii;);i./~-i .
• ............
l ooo
"
C
:ii: iii.ii ............. D
0
i
0
50
i
i
100 150 Iteration (t)
200
250
Fig. 5. The logarithmic curves of the energy function of the two channel segmentation corresponding respectively to the results presented in Fig. 4.
clearer results by going more deeply into the energy landscape, but it can not ensure the optimal result. However, by considering the energy function of the segmentation problem as the sum of the squares errors and using the algorithm (4), we get the result C which presents more clarity and smoothness than the A and B. The difference between the results C and D in Fig. 4 is that the tumor region colored red in D is not detected in C and is colored yellow, the same as the cerebrospinal fluid (csf). However, comparing the segmentation results to the raw images, and considering the energy value of the network at equilibrium, the result D is considered as the best.
Hopfield neural network for the multichannel segmentation
925
Fig. 6. The three channel segmentation results with seven classes+ A, B, C, and D are obtained respectively with the same approaches described for Fig. 4.
ECt) 3500
i I¸I:
3000
2500
i
2000
'~ ~'
A
1500 •
~
B
1000
::::--
5oo
...................
C
:~::::211~I:IILY21~L:IY::I:I::::2.:~::::: D
0
50
100 150 Iteration (t)
200
Fig. 7. The logarithmic curves of the energy function corresponding respectively to the results of Fig. 6.
In the case of the three channel segmentation, the same remarks can be said with the addition that the energy function is more optimal and the convergence speed is faster than in the case of two channel segmentation. Such behavior is similar to that of the human brain, if we have more information about a problem then we can understand it faster and more clearly. Also, a very important step for the automation process is to fix the convenient convergence time of the network to ensure the optimal result. However, here we fix the class number to seven for the presented case based on some medical information and also on the smoothness of the regions that present the results. In fact the energy value at equilibrium decreases if we increase the class number, but there are not any smooth regions in the segmented image. In fact, this is due to the intensity variation of the data which may have to be smoothed carefully without blurring the image before being used for the segmentation process.
I
(a)
(b)
(c)
Fig. 8, Magnetic resonance head images of a patient diagnosed with a multiple sclerosis in the brain (a) T1weighted (b) T2-weighted and (c) the proton density-weighted.
926
R. SAMMOUDA et al.
Fig. 9. Two channel segmentation results with eight classes using the Tl-weighed and the T2-weighted images obtained with our approach.
Fig. 10. Three channel segmentation results with eight classes obtained with our approach.
To verify the generality of our approach and its parameters on the data used, we operated the program with different subject data. We found that the parameters values of #(t), T~, c~ and/3 used before are also convenient for our database which contains images of several normal and abnormal cases. Figure 8 (a)-(c) shows respectively the Tl-weighted, T2-weighted, and proton densityweighted images of a patient diagnosed with a multiple sclerosis in the white matter. Also, as in the first case, the T2-weighted image shows an abnormal bright region in the white matter. We have used our approach with the same parameters described above, for two and three channel segmentation of these images. Figure 9 shows the two channel segmentation result with eight classes, using the Tl-weighted and the T2-weighted images. The red regions in the segmented image represent the abnormal regions. The color map shows the distribution of the pixels in the segmented image along their values in the T1- and T2- weighted images. Figure 10 shows the three channel segmentation result with eight classes. In this case, we have used eight classes because of the high intensity variations present in the raw images.
6. CONCLUSION In this study, we presented some contributions to a previous work in the use of Hopfield network for the segmentation of MR images. We showed that the performance of the network is mainly dictated by the choice of the energy function which can lead to a good segmentation when it is based on the sum of the squares
errors as a cost-term, and the addition of a noise term as an excitation to the network to detect and escape from certain local minima and be closer to the global minimum. Also, the use of the network with a prespecified time convergence promises a full automatic segmentation in the future if we can determine automatically the class number of the data during the segmentation process. In our current work, we focus on a preclassification process to uniform carefully the intensity values in the raw images without blurring the edges between the different regions such that we can determine the class number automatically during the segmentation process.
REFERENCES
I. M. W. Vannier, R. L. Butterfield, D. L. Rickman, D. ZM. Jordan, W. A. Murphy, R. G. Levit and G. Mohktar, Multispectral analysis of magnetic resonance images, Radiology 154(1), 221-224 (1995). 2. M. I. Kohn, N. K. Tanna and G. T. Herman, Analysis brain and cerebrospinal fluid volumes with MR imaging: Part 1 Methods, reliability and validation, Radiology 178, 115-122 (1991). 3. E. Herkovits, A hybrid classifier for automated radiological diagnosis: Preliminary results and clinical applications, Comput. Methods and Programs in Biomedicine 32(1), 45-52 (1990). 4. S. C. Armatur, D. Piraino and Y. Takefuji, Optimization neural networks for the segmentation of magnetic resonance images, IEEE Trans. Medical Imaging 11(2), 215-220 (1992). 5. R.A. Jacobs, Increased rates of convergence through learning rate adoption, Neural Nem,orks 1, 295-307 (1988). 6. Z. R. Novakovic, Solving systems of non-linear equations using the Lyapunor direct method, Computer Math. Apllic. 20(12), 19-23 (1990). 7. S. Schaffier and H. Warsitz, A trajectory-following method for unconstrained optimization, J. Optimization Theory and Applications 67(1), 133-138 (1990). 8. E Aluffi-pentini, V. Parsi and F. Zirilli, Global optimization and stochastic differential equations, J. Optimization Theory and Applications 47, 1-16 (1985). 9. S. Geman and C.R. Hwang, Diffusions for global optimization, SIAM J. Control and optimization 24(5), 1031 1043 (1986). 10. R. Sammouda, N. Niki and H. Nishitani, Multichannel segmentation of magnetic resonance cerebral images based on neural networks, IEEE ICIP-95, Washington, pp. 484487, October (1995) (1986).
Hopfield neural network for the multichannel segmentation About the Author - - RACHID SAMMOUDA received his B.S. degree in Computer Science from the University of Tunisia in 1991 and the M.S. degree in Information Science and Intelligent Systems from the University of Tokushima in Japan in 1994. He has been a Ph.D. student at the University of Tokushima until March 1996, when he started to teach as Assistant Professor at the Department of Optical Science there. His current interest are medical images analysis and artificial intelligence.
About the A u t h o r - - N O B O R U NIKI received his B.S. and M.S. in Electronic Engineering from Tokushima University in Japan in 1975 and 1977, respectively. Also, he received a Ph.D. degree in Electrical Engineering from Kyoto University in Japan in 1987. He was a Research Assistant in Information Science Department at Tokushima University from 1977 to 1989. In 1989, he became an Associate Professor in same department. Since April 1996 he has been a Professor at the Department of Optical Science and Technology. His current interests are image reconstruction algorithms, medical graphics, pattern recognition, and parallel image processing.
About the A u t h o r - - HIROMU NISHITANI received his M.D. degree from the University of Kyushu in Japan
in 1970, and Ph.D. Medical Science from the University of Kyushu 1984. He has been Professor and chairman of the Department of Radiology, School of Medicine, the University of Tokushima in Japan since 1988. His main interests are diagnostic radiology and medical image analysis.
927