"A fully integrated multiplatform image analysis utility"

Search provided by Google

9. Tissue Classification

9.1 Accessing the Segmentation Tool

Figure: The Segmentation Menu. Left: FSL found, Right: FSL not found - see Section A.1 for details on installing and configuring FSL.

Image Segmentation is the process of identifying features in images and marking them as distinct from one another. These features can be things like grey and white matter in an MRI of the brain or individual organs in an MR or CT image. We describe here the segmentation tools in BioImage Suite, some of which leverage the functionality of the FSL library of image analysis tools (See http://www.fmrib.ox.ac.uk/fsl/). See also Appendix B for instructions as to how to install FSL.

The Segmentation tools are accessed under the Segmentation Menu in most viewers - see Figure 9.1. The entries in the menu change depending on whether BioImage Suite successfully finds an installation of FSL on the computer. If FSL is found, two additional options appear for calling the Brain Extraction Tool (Skull Stripping) and the FSL - Gray/White Segmentation methodology respectively.

Many of the graphical user interfaces described below make use of a simple embedded image GUI control shown in figure 9.2. The appearance of this is somewhat variable but the basic concept is that it stores an image (usually an algorithm output image) and enables easy manipulation of this.

Figure 9.2: The image GUI control for manipulating an image. The filename is at the top (empty) with the dimensions (1 1 1) at the top right. The textbox below gives a longer description of the image. Common functionality (Load, Save, Info, Display) is available from the buttonbar below. Sometimes additional buttons appear here depending on the viewer. For example, in the Brain Register applications there are two buttons ``Display Ref'' and ``Display Trn'' instead of the ``Display'' button for displaying the image in either of the two viewers.

9.2 Math Morphology

Figure 9.3: Left: The Segmentation Control with the Morphology Tab selected. Right: Some examples of the results of the operations in the Math Morphology Tab. The original is at the top left. The next image in the cascade is the mask image generated by thresholding. This mask is then eroded, yielding the next image in the series. The threshold result mask was then dilated to generate the fourth image. The threshold result mask was also passed through the BMedian function repeatedly ( 4x) to remove outlying ``islands'' (note the midbrain). These provide rough examples of the results possible using these tools.

This section of the Image Segmentation toolbox provides some functions for simple mathematical morphology. These functions require you to threshold the image first, yielding either a mask in addition to the original image (in the objectmap viewer) or a binary image in the 3-component viewers. They are accessed via the (Segmentation | Math_Morphology) menu selection, as shown in Figure 9.3 above.

The top section of the Math_Morpology section of the dialog box contains two sliders for setting upper and lower thresholds. Once they are set, simply press the Threshold button, and if the Display Mask checkbox is checked, the thresholded image will be displayed in the viewer. If not, simply hit the Display button to send the result of the threshold operation to the viewer. If, instead of thresholding, you want to generate a mask based on the output of the brain extraction tool, which is useful for generating whole-brain masks, hit the large Generate Mask from Brain Extract Tool Output button.

The Simple Segment button is just a macro for segmenting that calls the threshold function followed by a couple other functions described below (The order is: Threshold; Erode; Connectivity; Dilate)..more on this function and its constituents below. See Fig 9.3 for some examples of a rough segmentation of an MRI brain image into grey and white matter regions. The red regions are the mask as shown in the 4-component viewer.

Math Morphology Operations

In the lower half of the Math_Morphology tab, there are a group of buttons that operate on the mask generated either by the thresholding or brain extraction tool import done in the upper half of the tab. If the ``2D'' checkbutton is enabled all operations, apart from Connectivity, are performed in 2D (i.e. for each slice), otherwise the operations (e.g. Dilate/Erode/Median) are performed in 3D which is the default. Their functions are summarized below:

Note: The erode and dilate operations utilize as a parameter a kernel size, which sets how many of each voxels neighbors in each dimension should be examined. This is set with the ``Kernel Size'' menu. If the 2D option is selected, only the within-slice dimensions are considered. Essentially, the kernel size parameter controls how aggressive the erosion or dilation is.

Erode: Erode shrinks the mask by looking at the surrounding neighborhood of each voxel. If any of the neighbors have a value of zero, then the voxel in question is set to zero.

Dilate: Dilate grows the mask in the exact opposite manner as erode shrinks it: if any neighboring voxel has value 1, then the voxel in question is set to 1.

BMedian: BMedian polls the surrounding voxels and if the majority are non-zero, then the voxel in question is set to 1; if most of the neighbors are zero, then the voxel is set to 0.

Connect Foreground: This operation is dependent on the location of the cursor. It sets to zero any voxels in the mask that are not contiguous with the voxel indicated by the seed. The position of the seed can either be obtained using the viewer cursor or specified manually. (This voxel must be equal to one, otherwise the whole image goes to zero!) This is useful for separating out regions of interest.

Connect Background: This operation is dependent on the location of the cursor. It sets to one any voxels in the background that are not contiguous with the voxel indicated by the seed (This voxel must be equal to zero!) This is useful for eliminating small holes inside foreground regions.

9.3 Histogram Segmentation

This tab provides access to two algorithms, Histogram Segmentation and Region-based Segmentation. Both algorithms take an input image and generate an objectmap of values 0 to N-1, where N is the number of classes. Histogram Segmentation operates exclusively on the image histogram which makes it very fast, the downside is that there are no regional constraints on the segmentation. It is essentially a customized implementation of k-means clustering optimized for this application. The region segmentation extends Histogram segmentation (which it uses for initialization) to allow for regional homogeneity constraints using a Markov Random Field (MRF) smoothness model.

Histogram Segmentation This algorithm divides the image into N normally distributed intensity classes using a k-means clustering algorithm. It has three standard parameters. (Additional parameters can be set in the ``experimental'' tab which are not discussed here). The first parameter is ``Classes'' whish sets the number of classes (N) in which the image is to be divided. The second parameter is sigma ratio which sets the maximum value of the ratio of the maximum standard deviation to the minimum standard deviation. Setting this to 1 ensures that all classes have the same standard deviation. A small non-zero value is useful to ensure that small classes do not disappear. The last option is the preprocessing mode which describes the pre-processing applied to the image prior to generating a histogram. It has three settings: None - no preprocessing, Positive - all non-positive voxels are ignored and Autorange which eliminates voxel outside the 1%-99% span of the cumulative histogram.

Figure 9.4: The Histogram Segmentation Tab.

The class means for upto seven classes are also displayed. This can be initialized manually - set the ``Use Initial'' checkbutton to use the initial parameters. On completion of the segmentation - invoked by the ``Segment'' button the values in the mean text boxes reflect the results. The ``Show Params'' button provides additional details.

Region Segmentation This algorithm essentially adds a Markov Random Field smoothness model to the previous algorithm. This has the result that the segmentation can no longer be performed on the intensity histogram alone hence this algorithm is significantly slower, as each voxel needs to be examined individually. It uses all the parameters from the Histogram Segmentation as well as the following additional four parameters (plus additional ones in the experimental tab which are not described here). (a) Iterations - this is the maximum number of iterations, (b) Smoothness - this sets the weight of the spatial smoothness component, setting this to zero reverts to histogram segmentation, (c) Convergence - sets the maximum percentage of voxels that are allowed to change in a given iteration and the algorithm considered to have converged and (d) Noise sigma - which is an estimate of the image noise standard deviation. This algorithm is very similar to the FSL/FAST algorithm with one key difference. FSL/FAST attempts to perform intensity-inhomogeneity correction at the same time as segmentation (which works OK at 1.5T, less at higher field strengths) whereas ``Region Segmentation'' is a pure segmentation method and assumes that the Bias field has been corrected previous to its use.

Segmentation Output This control enables the output of the histogram/region segmentations to be used to generate a mask for either bias field correction (see below) or math morphology operations above. By default the output of the last histogram or region segmentation run is stored in the ``Segmentation Output'' control which has options for Loading/Saving and Displaying this objectmap.

Below this, one can select which classes to use to create a mask. Clicking ``Generate'' takes the ROI defined by these classes and generates a binary mask for use in the Math Morphology Tool. Clicking ``Do Bias Field'' executes the polynomial bias field correction algorithm using the mask defined in this control and the parameters defined in the Bias Field Correction tab.

9.4 FSL Brain Extraction Tool

The brain extraction tool is used to remove the skull from an image, leaving only the region occupied by actual brain tissue. It segments these by using the dark space between the skull and brain, occupied by the Cerebro-spinal Fluid (CSF). This tool comes from the external program FSL's toolkit (see the FSL webpage. Thus, it will only appear if FSL is installed. (See also Section A.1 for instructions as to how to install FSL.)

Figure 9.5: The Brain Extraction Tool Tab - this calls the external program BET from FSL.
Figure 9.6: Brain Skull Overlay.

9.4.1 The BioImage Suite Interface to BET

In the ``Brain Extraction Tool'' tab in the Segmentation Control window, or via the ( Segmentation | Brain Extraction Tool ) menu selection, you have just a couple parameters to potentially set, and a button to execute the tool. The options are fairly straightforward.

To make the algorithm more aggressive, removing a larger portion of the image, set the Fractional Intensity Threshold slider to a higher value. To skew the stripping towards either the top or bottom of the image, use the Gradient Threshold slider. A positive value here yields a brain that has been stripped more aggressively at the top. If you check the ``Use Current Crosshairs'' checkbox, the location of the crosshairs in the viewer will be taken into account, and a sphere centered on that point will be the starting point for the brain extraction. This can help if there is a feature in the image that you want to be sure to remove. Set the crosshairs to a location away from this spot before executing the extraction, and check this box. The ``Initial Radius'' field only operates with the ``Use Current Crosshairs'' option selected, and sets the initial size of the sphere that i used to initialize the stripping. To execute simply press the Execute Brain Extraction Tool button. An output dialog box will pop up and inform you of the progress.

Figure 9.7: Methods used for the First Brain Extraction.

1. Choose (Segmentation | Brain Extraction Tool) - See Figure 9.7.
2. Put crosshairs in center of brain above corpus callosum along the midline
3. Click use current crosshairs (make the box to the left, red in color)
4. Set Fractional Intensity Threshold to 0.6 and the Gradient Threshold to 0.0
5. Click execute brain extraction tool - BioImageSuite will automatically save new stripped brain as [filename]_stripped.hdr
6. Click Display Under Brain/Skull Overlay to check if any tissue is erroneously included or excluded
7. Examine outline of stripping by clicking and holding the left mouse button and dragging it back and forth over one view space (coronal) while watching the other (sagittal).
   • If part of the brain is excluded: Decrease Fractional Intensity Threshold (in increments of 0.1 or 0.5) and repeat steps 2-7.
   • If too much skull/throat/neck is included: Increase Fraction Intensity Threshold (in increments of 0.1 or 0.5) and repeat steps 2-7.

Figure 9.8: Methods used for the Second Brain Extraction.

1. Choose (File | Load), then (see Figure 9.8)
2. Choose [filename]_stripped.hdr and click Open
3. Choose (Segmentation | Brain Extraction Tool). As above in steps 2 and 3, put the crosshairs in the center of the brain and make sure that the box to the left of use current crosshairs is red.
4. Change Fractional Intensity Threshold to 0.4 or 0.3 and the Gradient Threshold to 0.2
5. Click Execute Brain Extraction Tool. BioImageSuite will automatically save second stripped brain as [filename]_stripped_stripped.hdr.
6. Click Display under Brain/Skull Overlay to check if any brain is erroneously included or excluded.
7. Examine outline of stripping as stated above in Step number 7.

When the brain extraction is done, the results will be automatically saved, with ``_stripped'' appended to the original filename. The filename is placed into the ``Stripped Brain'' box at the bottom of the window, where you can save it somewhere new (Save ), load another image (Load ), get a popup with some information about the image dimensions and orientation (Info ), and display it in the viewer (Display ). The image is displayed in the Results view of whichever viewer contained the original. A brain/skull overlay image is also generated (example at right), which shows where the boundary between brain and non-brain has been set. This lets you quickly see if any regions in your image are being erroneously included or excluded from the brain, and change parameters accordingly. Here too, you can save and load the image, display its properties, and send it to the viewer.

Sometimes doing multiple rounds of brain extraction on the same image yields better results than a single extremely aggressive operation. To do this, simply execute the brain extraction tool, and copy the contents of the results view to the image view (see Image Display vs. Results Display), then repeat, with slightly more aggressive parameters.

9.4.2 A Recipe for Brain Extraction - Removing the Skull from a 3D image

The brain extraction tool is used to remove the skull from an image, leaving only the region occupied by actual brain tissue. It segments these by using the dark space between the skull and brain, occupied by the CSF. This tool comes from the external program FSL's toolkit (see the FSL webpage). Thus, it will only appear if FSL is installed. (Such integration is easier on Unix systems at this point).

9.4.2.0.1 First Brain Extraction:

The goal of the first brain extraction process is to remove most of the skull/throat/neck without loosing any brain tissue. This extraction will not be extremely accurate and may include more skull than desired in order to keep all of the brain tissue. The remainder of the skull can be removed in the second brain extraction. This process is described in Figure 9.7.

9.4.2.0.2 Second Brain Extraction:

The goal of the second brain extraction process is to be more accurate around the brain and remove the meninges and any remaining skull. This process is described in Figure 9.8

9.5 Grey/White Segmentation - using FSL Fast

Note: The Grey/White Segmentation control assumes that you have already performed a brain extraction (Described previously), yielding a brain with no skull.

Also note: The algorithm is optimized for a 1.5T scanner. Results may be different for other field strengths.

Upon execution of the Grey/White segmentation tool, by default the image is segmented into three classes: Grey matter, White matter, and Cerebro-spinal Fluid (CSF). You may specify more or less classes with ``Number of Classes'' menu on the right side of the window. When you click the Execute Automated Segmentation Tool button, a console pops up, showing the output of FSL program that performs this segmentation function. When the brain extraction is finished, the console will show a table of volumes for each class in the image.

Figure 9.9: Grey/White Segmentation Control.

The volumes are reported in $mm^3$, except for ``brain percentage''. These values can also be examined in more detail by clicking the Statistics button in the ``Post Process Output'' box. The image output is a classification map (filename: *_segm_restore.hdr), which has n + 1 frames, where n is the number of classes selected. The first frame of the image shows all classes, with a value between 0 and 1 for each class. The remaining frames are binary images that correspond to each class.

This algorithm assumes that the Brain is already skull stripped, use the Brain Extraction Tool to do so. The algorithm is very similar to the ``Region Segmentation'' algorithm described above with the addition of an integrated Bias Field Correction step. This works well for 1.5T brain images, its performance is worse at higher field strengths.

The outputs generated are a classification map which is binary objectmap with labeled tissue classes and a restored brain which is the result of applying the FSL/FAST bias field correction algorithm to the original image.

9.6 Bias Field Correction

Often MRI images are corrupted by intensity inhomogeneity or shading artifacts. This tab, shown in Figure 9.10, provides access to two algorithms, Slice Homogeneity Correction and Volumetric Bias Field Correction. Both aim to remove image intensity inhomogeneities caused by MR Scanner receiver coil sensitivity variations. Both algorithms may use a Mask Image to set the Region of Interest (ROI) in which the optimization is performed - this is stored in Mask Frame (labeled as C in Figure 9.10). The output of both algorithms is a Bias Corrected Image and optionally an estimate of the Bias Field which go to the Outputs Frame (D in Figure 9.10) in the bottom.

9.6.1 Slice Homogeneity Correction

Figure 9.10: The Bias Field Correction Tab.

This algorithm aims to eliminate slice-by-slice intensity variations. This is particularly useful in multi breathhold multislice acquisitions such as in the case of abdominal imaging. This very simple but surprisingly effective algorithm performs a best straight line fit between adjacent slices (using neighboring voxels as data point pairs - but excluding edge points) to estimate the scaling and (optionally if the ``Pure Scaling'' checkbutton is off) the offset between the two slices. If the ``Use Mask'' checkbutton is set then the Mask defined in the Mask Frame (C) is set. If the ``Median'' checkbutton is set then the line fitting is done using a least absolute value criterion (median - which can be more robust) whereas if it is off a standard least squares fit is performed. There are two buttons to invoke the algorithms. Run Slice - performs homogeneity correction only in the Z-direction (slice acquisition) Run All Slices - performs homogeneity correction in all three orthogonal axis directions (X,Y,Z). This is very useful as a pre-processing step for the polynomial (volumetric) bias field correction algorithm, described next. It is especially useful in cases of large inhomogeneities. See Frame A of Figure 9.10.

9.6.2 Polynomial Bias Field Correction

This algorithm is a custom implementation of the PABIC algorithm (Styner et al, TMI 2000 [111]). It aims to improve image inhomogeneity by estimating a polynomial model of the inhomogeneity that will result in the image intensities being more clustered with respect to their class means - the numbers of classes and the means invoke the Histogram Segmentation algorithm and its current parameter settings. The algorithm has 4 key parameters: (a) Resolution: 1 = Native, 2=Half Native etc. Since the bias field is a low frequency effect using a reduced resolution can be just as accurate and results in significant computational savings. (b) Degree: 1-3, determines the maximum degree of the polynomial used 1=linear, 2=quadratic, 3=cubic. (c) ``Use Mask'' - if on restricts the optimization to the ROI defined by the Mask Image (in frame C of Figure 9.10). (d) ``Use Initial Parameters'' - if on the user can specify initial values for the linear and some of the quadratic parameters to manually optimize the field estimate. Trial and error can be used to set this, clicking ``Apply'' simply applies these values with no optimization. Manual initialization may be necessary in the case of large inhomogeneities such as those present at high fields (e.g. 4T as in the example shown in Figure 4). An additional three buttons (Print,Load, Save) enable direct access to the polynomial coefficients underlying the bias field estimate.

The algorithm is invoked with the Optimize button. The Clear button resets the parameters to their default values (all zeros)

The Mask frame is used to manipulate the mask restricting the optimization of the algorithm. The Grab Morph button copies the current mask used in the Math Morphology tab to this control, whereas the Set Morph does the opposite, i.e. it transfers the current mask to the Math Morphology control for further editing.

Figure 9.11: Bias Field Correction Example. Left: The left panel shows an image before bias field correction. Note the inhomogeneity in the intensities of the image (the frontal and temporal lobes are darker than the rest of the brain). Right In the image shown in this panel, this artifact has been corrected. The bias field in this example was particularly difficult to parameterize, and the data was of poor quality, but the result image is significantly more uniform in intensity. Often, when the bias field is linear or parabolic in nature, the result is even better.

9.6.2.0.1 Outputs:

The corrected image from the last slice homogeneity or bias field correction operation is stored in the Bias Corrected Image control in frame D of the tab as shown in Figure 9.10. The estimated bias field is stored in the Bias Field control which has similar functionality.

9.6.3 A Recipe for Bias Field Correction for a Brain T1-Weighted Image

At higher field strengths, sometimes structural images acquire an intensity gradient across the image making some parts of the image brighter than others. This intensity gradient can influence segmentation algorithms erroneously, therefore a method has been developed to remove this intensity gradient from the image, it is known as bias field correction.

Figure 9.12: Methods used to correct the bias field of a structural image.

1. Choose (Segmentation | Brain Extraction Tool). A new Segmentation Control window will appear. Put the crosshairs in the center of the brain above the corpus callosum along the midline.
2. On the Segmentation Control Window, set the Fractional Intensity Threshold to 0.6 and the Gradient Threshold to 0.0
3. On the Segmentation Control Window, click use current crosshairs (make the box to the left, red in color). Click execute brain extraction tool - BioImageSuite will automatically save new stripped brain as [filename]_stripped.hdr. Click Display Under Brain/Skull Overlay to check if any tissue is erroneously included or excluded.
4. Once a satisfactory stripping of the brain has been achieved, click the Math Morphology tab of the Segmentation Control Window.
5. On the (Segmentation Control | Math Morphology) tab click Generate Mask from Brain Extract Tool Output. This will create a white mask of the stripped brain.
6. On the (Segmentation Control | Math Morphology) tab click Dilate to increase the size of the mask slightly.
7. On the Segmentation Control Window click the Bias Field Correction tab.
8. On the (Segmentation Control | Bias Field Correction) tab click Grab Morph to input the previous mask into the bias field correction algorithm.
9. On the Viewer window choose (File | Load)
10. Choose the stripped version of the brain image.
11. On the (Segmentation Control | Bias Field Correction) tab click Optimize. This will generate the bias filed correction and apply it to the image.
12. The results shown in (Figure 10, image II), show that before the correction the histogram did not differentiate between gray and white matter, while after there are two clear peaks in the histogram. Also, a simple segmentation is more precise after the correction rather than before.

9.7 Appendix: A Little bit of Theory

9.7.1 The K-means Clustering Algorithm

The K-means clustering algorithm (see Duda & Hart) can be used to essentially estimate the thresholds above. In our implementation we assume that the intensity in each tissue class can be modeled as having a normal distribution with mean $\mu$ and standard deviation $\sigma$. We will use the notation $p(i\vert l=c) = {\cal N}(\mu_c,\sigma_c)$ to define this, where $l$ is the label and $c$ is the class number. The goal of the version of the k-means algorithm that we will describe is to estimate the class parameters ($\mu,\sigma$) for all classes and then to assign each voxel to the class that maximizes the function:

$${\hat l({\bf x})} = {\arg\max \atop l } p (i({\bf x})\vert l({\bf x})=c)$$

$$= {\arg\max \atop l } \frac{1}{\sqrt{2\pi\sigma_c^2}}e^{\frac{-(i({\bf x})-\mu_c)^2}{2\sigma_c^2}}$$ (9.1)

A simpler form of the algorithm assumes that all $\sigma_c$'s have the same value. This reduces the problem to estimating the means only. The procedure can be described in recipe form as:

1. Decide on the number of classes $M$
2. Assign class centroids $\mu_c$ and optionally standard deviations $\sigma_c$ for each class. The most common way to do this is to equally space the $\mu_i$'s and set all $\sigma_i$'s to some constant value.
3. Estimate the labels $l({\bf x})$ using equation 9.1. This is an exhaustive optimization - compute $p(i\vert l=c)$ for all $l$'s and pick the $l$ that maximizes this probability.
4. Estimate a new set of $\mu_c$'s and $\sigma_c$'s by computing the mean and standard deviation of the intensities of the voxels labeled as having class $c$.
5. Repeat steps 3-4 until the parameters $\mu_c$ and $\sigma_c$ converge.

Note that, since the spatial position of the voxels ${\bf x}$ does not enter into the calculation, a quick way to implement this method is by first forming the image histogram and performing all operations on the histogram itself. This can speed up the iterations by a factor of 15000 or so in the case of an $128\times128\times 128$ image whose histogram can be described by 128 bins.

9.7.2 Imposing Local Smoothness using Markov Random Fields

The key weakness of the previous method is that, as noted, the spatial position of each voxel is not used during the segmentation. Spatial homogeneity is often a powerful constraint that can be used to improve the segmentation in the presence of image noise. This can be loosely thought of as finding the label at each voxel that is an optimal constraint between (i) the intensity at that voxel and (ii) the labels assigned to its neighboring voxels.

9.7.2.0.1 Markov Random Fields:

The probabilistic structure used, most frequently, to capture such homogeneity is to model the label (classification) image as a Markov Random Field. This (and we are skipping a lot of math here) basically reduces to describing the probability of each voxel belonging to class $l$, as having a Gibbs distribution of the form:

$$p(l({\bf x})) = k_1\exp(-W(L({\cal R}_{\bf x}),l))$$ (9.2)

where $k_1$ is a normalization constant, $L({\cal R}_{\bf x})$ is the set of labels of the neighboring voxels and $W$ is a positive semi-definite function. This is a way of converting an ``energy-function'' like smoothness term into a probability distribution for integration into a probabilistic framework. The function $W$ can take many forms, the one we will use here is from Zhang et al TMI 2001 [129] (the method in FSL):

$$W(l({\bf X}_n)= \sum_{{\bf x'} \in {\cal R}_{\bf x}} \delta(l({\bf x'}) - l({\bf x}))$$ (9.3)

This essentially counts the number of voxels in the neighborhood of the voxel at location ${\bf x}$ that have labels different from it.

9.7.2.0.2 Overall Formulation using the Expectation-Minimization Framework:

We first define the vector $\Theta$ as the collection of the means and standard deviations of all $C$ classes, i.e. $\Theta = [ \mu_0,\sigma_0,\ldots,\mu_{c-1},\sigma_{c-1} ]$. The goal of the segmentation is to estimate both the set of labels $L$ and the class parameters $\Theta$, given the image $I$. We can express this mathematically as:

$${\hat L},{\hat \Theta} = {\arg\max \atop L,\Theta } p(L,\Theta \vert I)$$ (9.4)

As is commonly done, this can be solved iteratively in the same spirit as the EM-framework as:

$$\hspace*{-0.1in}\mbox{\textbf{E-Step: }}\mbox{ } \Theta^k = ... ...ep: }}\mbox{ } L^k = {\arg\max \atop L } p(L\vert I,\Theta^k)$$ (9.5)

where $k$ is the iteration number. In the E-Step we estimate a new set of parameters $\Theta^k$ given the current classification $L^{k-1}$. In the M-Step, using the newly estimated $\Theta^k$ we estimate a new classification $L^k$.

E-Step: This is straightforward. For each class $i$ we estimate the mean and standard deviation of the intensities $I$ of all the voxels where $M=i$. This is identical to the procedure used in the k-means algorithm above.

M-Step: This takes the form of a Bayesian a-posterior maximization. First we express

$${\hat l({\bf x})} = {\arg\max \atop l } \log p(l({\bf x})\vert i({\bf x}),\Theta^k, L({\cal R}_{\bf x})$$

$$k_1 + \underbrace{\log p(i({\bf x}),\Theta^k \vert l({\bf x}) )}_... ...} + \underbrace{log p(l\vert L({\cal R}_{\bf x}))}_{\mbox{MRF Smoothness Term}}$$ (9.6)

where $k_1$ is a constant. This equation is easily maximized by a greedy search strategy as $M$ can only take values of $1 \ldots C$. The prior term on the classification, $p(L)$, can be defined by modeling $L$ as a Markov random field (see discussion above and equation 9.3). We express the likelihood (or data-adherence) term for each possible value of $l({\bf x})=c$ as:

$$p(i({\bf x}), \Theta^k \vert l({\bf x})=c ) = p(i({\bf x}) \vert \Theta^k , l({\bf x})=c )$$ (9.7)

which is similar to the model previously used in equation 9.1.

Webpage design by Marcello DiStasio.


This document was generated using the LaTeX2HTML translator Version 2002-2-1 (1.71)