Prostate Segmentation for Medical Interventions

Prostate Segmentation for Medical
Interventions
by
Seyedeh Sara Mahdavi
B.A.Sc., Sharif University of Technology, 2004
M.A.Sc., University of Tehran, 2007
A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF
THE REQUIREMENTS FOR THE DEGREE OF
DOCTOR OF PHILOSOPHY
in
The Faculty of Graduate Studies
(Electrical and Computer Engineering)
THE UNIVERSITY OF BRITISH COLUMBIA
(Vancouver)
November 2012
c Seyedeh Sara Mahdavi 2012
Abstract
Prostate cancer is the most prevalent type of cancer among men. Accurate
delineation and appropriate visualization of the prostatic region can greatly
affect treatment of prostate cancer and has the potential to reduce some of
the treatment side-effects. The main goal of this research is to develop a
prostate segmentation tool which is suitable to replace manual delineation.
Manual segmentation, the current standard in procedures such as low dose
rate prostate brachytherapy, is tedious, time consuming and observer dependent. We propose a 3D semi-automatic segmentation tool to overcome these
limitations. To show the clinical value of this method we perform extensive
dosimetric evaluation on in-vivo ultrasound images. This tool is currently
being clinically used as part of the prostate brachytherapy treatment procedure at the BC Cancer Agency.
Ultrasound is the common modality for imaging the prostate. Although
safe and simple to use, it can not always allow the prostate to be reliably
delineated. Vibro-elastography is a relatively new imaging method which
is used to characterize mechanical properties of tissue. We investigate the
suitability of vibro-elastography for visualizing the prostate. We compare
in-vivo B-mode ultrasound and vibro-elastograpy images with the gold standard MRI in terms of contrast, edge visibility and the shape and size of the
gland as seen in these images. Based on our results we develop a 3D automatic prostate segmentation tool in which, in addition to B-mode, information from vibro-elastography images is used in an iterative model-based
segmentation approach.
We conclude this work by studying the visibility of cancer itself in vibroelastography images. Areas suspected for cancer are manually marked on
the images and the results are compared to the marked cancer in registered
pathology slices. Our preliminary results show that vibro-elastography has
the potential to be used for detecting prostate cancer; however, we suggest a combined use of various modalities or image types to improve cancer
detection.
ii
Preface
This thesis is primarily based on four manuscripts, resulting from the collaboration between multiple researchers. All publications have been modified
to make the thesis coherent. Ethical approval for conducting this research
has been provided by the Clinical Research Ethics Board, certificate numbers: H08-02696, H06-70146.
A version of Chapter 2 has been published in:
• S. Sara Mahdavi, Nick Chng, Ingrid Spadinger, William J. Morris,
and Septimiu E. Salcudean. Semi-automatic segmentation for prostate
interventions, Medical Image Analysis, 15(2011)226-237
The contribution of the author was developing, implementing, and evaluating the method. Dr. Chng developed a user interface in order to simplify
the use of the original algorithm for evaluation and to prepare the method
as a tool to be used at the Vancouver Cancer Center. Dr. Spadinger helped
develop the idea of using nine sectors in evaluation. Dr. Morris and Prof.
Salcudean helped with their valuable suggestions in improving the methodology. All co-authors contributed to the editing of the manuscript.
A version of Chapter 3 is in press for publication in:
• S.S. Mahdavi, I. Spadinger, N. Chng, S.E. Salcudean, W.J. Morris.
Semiautomatic segmentation for prostate brachytherapy: Dosimetric
Evaluation, Brachytherapy, 2011
The contribution of the author was developing and implementing the dosimetric evaluation methods. The co-authors helped with their suggestions
in the evaluation and Dr. Spadinger created the treatment plans used in
the evaluation process. All co-authors contributed to the editing of the
manuscript.
A version of Chapter 4 has been published in:
iii
• S.S. Mahdavi, M. Moradi, X. Wen, W.J. Morris, S.E. Salcudean. Evaluation of visualization of the prostate gland in vibro-elastography images, Medical Image Analysis, 15(2011)589-600
Prof. Salcudean provided the data acquisition system, including a motorized 3D ultrasound imaging system and the associated software. The sections ’Edge strength’ is the contribution of Dr. Moradi in addition to his
help in the preparation and editing of the manuscript. Dr. Morris assissted
with patient recruiting, hospital access and data collection in the operating
room. Dr. Moradi and Dr. Wen contributed in the data collection process
and contributed in the initial vibro-elastography imaging and data processing software.
A version of Chapter 5 has been accepted for publication in:
• S.S. Mahdavi, M. Moradi, W.J. Morris, S.L. Goldenberg, S.E. Salcudean. Fusion of ultrasound B-mode and vibro-elastography images
for automatic 3D segmentation of the prostate, IEEE Trans. Medical
Imaging.
Drs. Morris and Goldenberg contributed to the medical aspects of the work,
such as providing access to the operating room and recruiting patients. Dr.
Moradi contributed to editing the manuscript and data collection.
Portions of the above papers, mainly text and figures from background
and literature review, also appear in Chapter 1.
In addition to the above Prof. Salcudean, as my supervisor, has helped
me with his valuable ideas and suggestions in the course of improving the
methods and algorithms, in addition to editing all the manuscripts.
iv
Table of Contents
Abstract
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
ii
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
iii
Table of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . .
v
Preface
List of Tables
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viii
List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
List of Abbreviations
x
. . . . . . . . . . . . . . . . . . . . . . . . . xviii
Acknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . xix
1 Introduction . . . . . . . . . . . . . . .
1.1 Motivation . . . . . . . . . . . . . .
1.2 Background . . . . . . . . . . . . .
1.2.1 Prostate Segmentation . . .
1.2.2 LDR Prostate Brachytherapy
1.2.3 Ultrasound Elastography . .
1.3 Thesis Objectives . . . . . . . . . .
1.4 Thesis Outline . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
. .
. .
. .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
1
4
7
7
10
14
19
20
2 Semi-automatic Segmentation of the Prostate
2.1 Introduction . . . . . . . . . . . . . . . . . . .
2.2 Methodology . . . . . . . . . . . . . . . . . . .
2.2.1 Algorithm . . . . . . . . . . . . . . . .
2.2.2 Evaluation . . . . . . . . . . . . . . . .
2.3 Results . . . . . . . . . . . . . . . . . . . . . .
2.3.1 Accuracy . . . . . . . . . . . . . . . . .
2.3.2 Repeatability . . . . . . . . . . . . . . .
2.3.3 Performance . . . . . . . . . . . . . . .
2.4 Discussion and Conclusions . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
22
22
23
27
33
38
38
41
43
46
v
3 Segmentation of the Prostate: Dosimetric Evaluation
3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . .
3.2 Patient Characteristics and Treatment Planning . . . .
3.3 Evaluation Measures . . . . . . . . . . . . . . . . . . . .
3.4 Dosimetric Evaluation of the Algorithm . . . . . . . . .
3.4.1 Impact of Planning using TES Contours . . . .
3.4.2 The Effect of Contour Variability on Planning .
3.5 Results . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.6 Discussion . . . . . . . . . . . . . . . . . . . . . . . . .
3.7 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
49
49
50
51
52
52
53
53
63
63
4 Visualization of the Prostate Gland in
4.1 Introduction . . . . . . . . . . . . . .
4.2 Materials and Methods . . . . . . . .
4.2.1 Data Acquisition . . . . . . . .
4.2.2 Evaluation Methods
. . . . .
4.3 Results . . . . . . . . . . . . . . . . .
4.3.1 Image Comparison . . . . . . .
4.3.2 Volume Comparison . . . . . .
4.4 Discussion and Conclusions
. . . . .
VE Images
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
65
65
65
65
68
77
77
80
81
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
87
87
88
88
90
98
100
102
102
102
105
6 Visibility of Cancer in VE Images of the Prostate
6.1 Introduction . . . . . . . . . . . . . . . . . . . . . .
6.2 Methods and Materials . . . . . . . . . . . . . . . .
6.2.1 Data Collection and Processing . . . . . . .
6.2.2 Cancer Visibility in VE Images . . . . . . .
6.3 Discussion and Conclusions . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
109
109
109
109
115
118
5 Automatic Segmentation of the Prostate
5.1 Introduction . . . . . . . . . . . . . . . .
5.2 Methods and Materials . . . . . . . . . .
5.2.1 Data Collection and Processing .
5.2.2 2D Segmentation Algorithm . . .
5.2.3 3D Segmentation . . . . . . . . .
5.2.4 Evaluation . . . . . . . . . . . . .
5.3 Results . . . . . . . . . . . . . . . . . . .
5.3.1 2D Segmentation . . . . . . . . .
5.3.2 3D Segmentation . . . . . . . . .
5.4 Discussion and Conclusions . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
vi
7 Conclusions and Future Research . . . . . . . . . . . . . . . . 123
7.1 Contributions . . . . . . . . . . . . . . . . . . . . . . . . . . 123
7.2 Future work . . . . . . . . . . . . . . . . . . . . . . . . . . . 125
Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128
vii
List of Tables
2.1
2.2
2.3
2.4
2.5
2.6
2.7
4.1
4.2
4.3
4.4
5.1
5.2
5.3
Comparison of prostate sectors created from Raw TES and
RO reviewed TES contours . . . . . . . . . . . . . . . . . . .
Comparison of RO reviewed TES and manual prostate sectors.
Comparison of Raw TES and manual prostate sectors. . . . .
Inter-observer variability in manual and TES contouring - as
characterized by the volume error Verr (Ant.: anterior, Lat.:
lateral, Post.: posterior, B: base, M: mid-gland, A: apex, Tot.:
total gland). . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Inter-observer variability in manual and TES contouring - as
characterized by the volume difference Vdif f (Ant.: anterior,
Lat.: lateral, Post.: posterior, B: base, M: mid-gland, A:
apex, Tot.: total gland). . . . . . . . . . . . . . . . . . . . . .
Intra-observer variability in manual and TES contouring - as
characterized by the volume error Verr (Ant.: anterior, Lat.:
lateral, Post.: posterior, B: base, M: mid-gland, A: apex, Tot.:
total gland). . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Intra-observer variability in manual and TES contouring - as
characterized by the volume difference Vdif f (Ant.: anterior,
Lat.: lateral, Post.: posterior, B: base, M: mid-gland, A:
apex, Tot.: total gland). . . . . . . . . . . . . . . . . . . . . .
Description of the data used in this chapter . . . . . . . . . .
CNR comparison of VE and B-mode images . . . . . . . . . .
The percentage of non-stationary prostate edge profiles in
different areas of the B-mode and VE images. Standard deviations reported characterize inter-patient variations. . . . .
Average edge continuity measure, K, for the nine regions of
the gland. . . . . . . . . . . . . . . . . . . . . . . . . . . . .
39
40
40
42
43
44
45
68
77
79
80
The main steps in the ASM algorithm. . . . . . . . . . . . . . 92
2D segmentation results . . . . . . . . . . . . . . . . . . . . . 102
3D segmentation results . . . . . . . . . . . . . . . . . . . . . 103
viii
6.1
The sensitivity, percentage of undetected tumors (false negatives/number of tumors) and the percentage of incorrect
guesses (false positives/number of guessed tumors) for the
four regions of the matching slices and the whole gland. . . . 117
ix
List of Figures
1.1
1.2
1.3
1.4
An illustration of the prostate and its periphery, in the coronal
(left) and sagittal (right) views, displaying the bladder (BL),
rectum (R), seminal vesicles (SV) and vas deferens (VD). The
three zones are also shown. The central zone (CZ) surrounds
the ejaculatory ducts (ed), the transitional zone (TZ) surrounds the prostatic urethra (u) proximal to the verumontanum (V) and the peripheral zone (PZ) forms a majority of
the gland. The fibromascular stroma (FS) lies in the midanterior part of the prostate and lacks glandular components.
2
MRI images of the prostate in the coronal (left) and transverse (right) views displaying the outline of the prostate (Pdashed line), the bladder (BL) and the penile bulb (PB). In
the transverse image an endorectal coil (EC) is placed inside the rectum for higher image quality, better showing the
central zone (CZ), peripheral zone (PZ) and neurovascular
bundles (NV). . . . . . . . . . . . . . . . . . . . . . . . . . . .
4
Transverse ultrasound B-mode images of the prostate at the
base (left) and mid-gland (right). The seminal vesicles (SV)
and the shadow of the bladder (BL) are visible at the base.
At the mid-gland, the location of the transrectal ultrasound
(TRUS) probe inside the rectum (R) is labeled. The various zones of the prostate (P) are not always distinguishable
in TRUS images. However, in this case the neurovascular
bundles can be faintly seen at locations 1-B and 1-F of the grid. 5
A snapshot of the software used for treatment planning (VariSeedT M ),
showing the B-mode images of the prostate taken at 10 depths
(smaller images at right). In the larger image, the CTV (red)
and PTV (cyan), needle locations (hollow circles), and the
seed locations (cyan circles) are shown. The 200%, 150% and
100% isodose lines are also displayed. . . . . . . . . . . . . . . 13
x
1.5
1.6
1.7
1.8
1.9
2.1
2.2
2.3
2.4
2.5
2.6
2.7
Schematic representation of vibro-elastography imaging. From
left to right, a broadband external mechanical excitation x0 (t)
is applied to the sample; the axial motion xjk (t) of a tissue
block indexed axially by j and laterally by k within the ultrasound image is obtained from RF images at time t and t − ∆t;
transfer functions from a reference to each of the tissue blocks
are estimated, with the illustration showing the reference being the motion xjk (t) of block jk; images are produced by
displaying the L2 -norm of the difference in transfer functions.
System setup for vibro-elastography data acquisition . . . . .
Details of the motorized cradle . . . . . . . . . . . . . . . . .
An illustration of the longitudinal planes at which data is collected, with respect to the prostate and TRUS probe. Here,
the mid-sagittal plane, at 0◦ , and two other planes at angles
of -20◦ and +10◦ are shown. The orientation of a transverse
plane in the mid-gland is also displayed. . . . . . . . . . . . .
A mid-sagittal (left) and transverse (right) VE image of the
prostate showing the urethra (dashed red). The prostate
boundary is marked by the dashed yellow lines. The section from which the transverse image is displayed is marked
by the dashed blue line. . . . . . . . . . . . . . . . . . . . . .
A snapshot of the graphical user interface for automatic prostate
segmentation used at the Vancouver Cancer Centre. . . . . .
The main steps of the segmentation algorithm. . . . . . . . .
a) Initialization points, b) Image un-warping, IMMPDA edge
detection and tapered ellipse fitting (dashed), c) Image untapering, IMMPDA edge detection and ellipse fitting (dashed).
28
The effect of changing the tapering parameter from -1 (thick
dashed line) to 1 (solid line). A tapering value of zero corresponds to an ellipse (dotted line). . . . . . . . . . . . . . . .
a) Final contours on the TRUS B-mode images, b) Final semiautomatic 3D volume (thick lines) compared to manual segmentation (thin lines). . . . . . . . . . . . . . . . . . . . . .
Division of the gland into nine sectors. . . . . . . . . . . . .
Schematic of changes made to each sector of the TES contours
after modification. . . . . . . . . . . . . . . . . . . . . . . . .
15
17
17
18
19
25
26
30
34
36
38
xi
3.1
3.2
Plot illustrates the means and 95% confidence intervals of
the paired difference in the V100 when the treatment plans
generated based on Raw TES PTVs are overlaid on RO reviewed TES PTVs (i.e. V100 of the plan generated on Raw
TES PTVs with respect to the Raw TES PTVs, the observed
clinical baselines, subtracted from V100 of the plan generated
on Raw TES PTVs with respect to the RO reviewed TES
PTVs). A negative result indicates a decrease in V100 . Each
column of data represents a particular region, labeled at the
bottom. Each data point represents the result of each of the
41 cases analyzed. The p-values for a null-hypothesis of no
effect in each region are presented at the top of the figure.
Two data points fall outside the range of this figure, and are
indicated by ’+1’. Their values are -46.5% (anterior base) and
-50.4% (posterior apex). A statistically significant reduction
in V100 is evident in most regions, although the mean magnitude of the reduction is less than 5% (Ant.: anterior, Lat.:
lateral, Post.: posterior, B: base, M: mid-gland, A: apex). . .
Plot illustrates the means and 95% confidence intervals of
the paired difference in the CI100 when the treatment plans
generated for Raw TES PTV’s are overlaid on RO reviewed
TES PTVs (i.e. CI100 of the plans created on Raw TES PTVs
with respect to the Raw TES PTVs, the observed clinical
baselines, subtracted from the CI100 of plans created on Raw
TES PTVs and overlaid on RO reviewed TES PTVs). Each
column of data represents a particular region, labeled at the
bottom. Each data point represents the result of one of the 41
cases analyzed. A negative value result indicates a CI100 less
than that of the observed clinical baseline. The number of
data points falling outside the displayed range in each region,
are indicated. The range of values for these four sectors: Ant.
B, Ant. A, Lat. A, and Post. A, are: [-1.70 0.73], [-2.32 1.36],
[-2.10 0.61], and [-1.24 1.49]. The lateral base and posterior
base exhibit a statistically significant increase in conformity
(Ant.: anterior, Lat.: lateral, Post.: posterior, B: base, M:
mid-gland, A: apex). . . . . . . . . . . . . . . . . . . . . . . .
57
58
xii
3.3
3.3
3.4
3.4
4.1
4.2
4.3
The means and 95% confidence intervals of the V100 when (a)
a reference plan, created for one set of contours, is overlaid on
the manual contours of 10 other observers in 5 patients. The
typical goal V100 goal of 97% is shown as a dashed horizontal
line. For comparison, the results of performing the test using
the TES contours or their derived plans are displayed using
the triangle symbols (Ant.: anterior, Lat.: lateral, Post.: posterior, B: base, M: mid-gland, A: apex). . . . . . . . . . . . .
Continued. The means and 95% confidence intervals of the
V100 when (b) a reference PTV is overlaid on the treatment
plans created for the manual contours of 10 radiation oncologists, for a single patient (second of the five patients in
(a). The typical goal V100 goal of 97% is shown as a dashed
horizontal line. For comparison, the results of performing
the test using the TES contours or their derived plans are
displayed using the triangle symbols (Ant.: anterior, Lat.:
lateral, Post.: posterior, B: base, M: mid-gland, A: apex). . .
The means and 95% confidence intervals of the CI100 when (a)
a reference plan, created for one set of contours, is overlaid on
the manual contours of 10 other observers in 5 patients. For
comparison, the results of performing the test using the TES
contours or their derived plans are displayed using the triangle
symbols (Ant.: anterior, Lat.: lateral, Post.: posterior, B:
base, M: mid-gland, A: apex). . . . . . . . . . . . . . . . . . .
Continued. The means and 95% confidence intervals of the
CI100 when (b) a reference PTV is overlaid on the treatment
plans created for the manual contours of 10 radiation oncologists, for a single patient (second of the five patients in (a)).
For comparison, the results of performing the test using the
TES contours or their derived plans are displayed using the
triangle symbols (Ant.: anterior, Lat.: lateral, Post.: posterior, B: base, M: mid-gland, A: apex). . . . . . . . . . . . . .
59
60
61
62
Transverse B-mode (left), VE (middle) and MRI (right) prostate
images of two patients. The boundary of the prostate is partially segmented in the second set of images. . . . . . . . . . . 69
VE (top) and B-mode (bottom) sagittal images of the prostate
of three different patients. The boundary of the prostate is
outlined in one of the patients. . . . . . . . . . . . . . . . . . 69
Division of the prostate into nine sectors . . . . . . . . . . . . 70
xiii
4.4
4.5
4.6
4.7
4.8
4.9
4.10
4.10
4.11
4.11
5.1
Selection of inside (yellow boxes) and outside (blue boxes)
regions of the prostate for CNR computation in (a) B-mode
and (b) VE images. . . . . . . . . . . . . . . . . . . . . . . . .
(a) The radii used for edge profile extraction in a VE image,
originating from C and with angles θi = i × 30◦ , i = 1, ..., 12.
For one of the radii, the two neighboring radii are also illustrated. (b) A close-up view of one of the rays used for extracting the edge profile Iθi (r). The measures of edge strength are
calculated for the window of r ∈ [rθi −∆d, rθi +∆d], where rθi
is a manually selected edge point along the ray. The neighboring edge profiles Iθi ±δθ , are extracted similarly and used
along with Iθi in the edge continuity measure. . . . . . . . . .
Registration of B-mode/VE prostate surfaces to MRI . . . . .
The values of the normalized gradient-based edge strength
measure (M ) in arbitrary units vs. the distance from the
edge in cm, separated for the nine regions of the prostate gland.
The values of the gradient-based edge strength measure (M )
along the edge profiles for VE and B-mode images. Error bars
represent the inter-patient standard deviation of the M values.
Comparison of VE (magenta) vs. MRI (blue) 3D surfaces, on
the left, and B-mode (green) vs. MRI (blue) 3D surfaces, on
the right, from one of the patients. . . . . . . . . . . . . . . .
A comparison between VE vs. MRI volume error and B-mode
vs. MRI volume error, showing the mean and inter-patient
standard deviation of Verr % for three observers. Sample data
points are also shown as gray dots. Figure continues on the
next page (Ant.: anterior, Lat.: lateral, Post.: posterior, B:
base, M: mid-gland, A: apex). . . . . . . . . . . . . . . . . . .
continued . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
A comparison between VE vs. MRI volume difference and
B-mode vs. MRI volume difference, showing the mean and
inter-patient standard deviation of Vdif f % for three observers.
Sample data points are also shown as gray dots. Figure continues on the next page (Ant.: anterior, Lat.: lateral, Post.:
posterior, B: base, M: mid-gland, A: apex). . . . . . . . . . .
continued . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
84
85
Illustration showing the span of the RF data planes with respect to the TRUS probe. . . . . . . . . . . . . . . . . . . . .
89
71
72
76
78
79
81
82
83
xiv
5.2
5.3
5.4
5.5
5.6
A B-mode image (top left) and the corresponding vibro-elastography
image (top right) of the prostate mid-gland along with the
phase symmetry of the B-mode (middle left) and vibro-elastography
(middle right). The sum of the B-mode image and its phase
symmetry is shown in the lower left. . . . . . . . . . . . . . . 91
The 30 boundary points used to model the prostate boundary (left) highlighting the 10 main points as red circles and
blue squares and the aligned and scaled, manually segmented
prostate boundaries in the training set, showing only the 10
main points (right). We define the posterior points as the
main points shown as blue squares and the four points in
between them. . . . . . . . . . . . . . . . . . . . . . . . . . . 93
Illustration of how for point i on the model edge of training
image j, a normalized edge derivative profile, gij normal to
the boundary and centered at the model point and of length
2np is extracted. The average of gij over the number of training images results in g¯i . . . . . . . . . . . . . . . . . . . . . . 94
Illustration showing the computation of gray level similarity
and edge continuity for two boundary points. In the upper
point the edge derivative profile hi (d) is compared with that
of the model, g¯i (shown in the rectangle) along a line normal to the edge and of length 2l. dgE is the distance from
the current boundary point where the maximum similarity is
obtained. In the lower point the gray level intensity profile,
ei (d), is compared to its two neighbors, e1i (d) and e2i (d), to
obtain correlation functions, C1i (d) and C2i (d). dKE is the
distance from the current boundary point to where a maximum correlation exists between the center profile and its
neighbors. A similar procedure is applied to B-mode images
to obtain dgU S and dKU S . . . . . . . . . . . . . . . . . . . . . 97
The main steps of the 3D semi-automatic prostate segmentation algorithm proposed in [60] including the applied modifications. Manual initialization is replaced by automatic segmentation of the mid-gland and vibro-elastography images
are included in the segmentation process. . . . . . . . . . . . 99
xv
5.7
5.8
6.1
6.2
6.3
2D segmentation results shown on corresponding B-mode (left)
and vibro-elastography (right) mid-gland images. The initial
contour Xinit is shown in dashed yellow and the final contour
in solid red. The final 3D segmented surface (triangulated
blue) is also compared to the manually created surface (red)
in the lower image. . . . . . . . . . . . . . . . . . . . . . . . . 104
The anterior peri-prostatic venous plexus is observed as a
dark region in B-mode and VE images. When its contrast
is similar to the prostate-to-background contrast, it can be
mistaken as the gland itself. In this figure the dashed yellow is
the initial automatic contour, and solid red is the final. Solid
blue indicates the manually segmented prostate boundary. . . 107
Illustration of the transfer function improvement by removing
data with low CC. Approximately 80 frames of RF data are
collected from the prostate, for each angle, using the TRUS
probe (left). For each region in the probe’s field of view (red
square) the tissue displacement (lower right plot) is estimated
through time, as the tissue is vibrated. Along with the displacement estimation, the correlation coefficient (CC) is also
computed (upper right plot). Instead of computing the transfer function (TF) using all the displacement data, a threshold
on the CC (in this case 0.85) is set to remove erroneous displacement data. TFs are computed from portions of correct
data (TF1 and TF2) which are finally averaged to provide
the final TF. The average value of the magnitude of the TF,
over the vibration frequency range, defines the intensity of
the corresponding pixel (red square). . . . . . . . . . . . . . . 111
The effect of removing data with low correlation (CC). The
two figures show the magnitude of the transfer function from
the same prostate in the sagittal view, before (left) and after
(right) removing low CC data. In this case a threshold of
CC=0.85 was used. . . . . . . . . . . . . . . . . . . . . . . . . 112
An example of pathology slides of one case with cancerous
regions and grades defined. Starting from the apex, the number defined as ‘level’ corresponds to the depth at which the
sample was taken. Cancer is marked by dotted lines. The left
and right side of the prostate is marked as ‘L’ and ‘R’, the
anterior and posterior marked as ‘A’ and ‘P’, respectively. . . 113
xvi
6.4
6.5
6.6
6.7
An example of manual correction of the orientation in pathology slides based on the gland internal structure. In this
case, guided by the prostate boundary and internal structures (marked as blue dashed and black solid lines), the slice
on the right will be rotated counter-clockwise to align with
its prior slice on the left. . . . . . . . . . . . . . . . . . . . . .
The set of pathology images are converted to a Stradwin file
format. The prostate and the marked cancer regions are then
segmented in this software. . . . . . . . . . . . . . . . . . . .
Examples of how evaluation of cancer visibility in VE was performed in five pairs of transverse VE (left column) and pathology (right column) images. Corresponding VE and pathology
images were found via registration and guessed tumors were
marked on VE images and compared to those in pathology.
Matching tumors (true positives) are linked with red lines and
the false positive and false negatives are marked by FP and
FN, respectively. . . . . . . . . . . . . . . . . . . . . . . . . .
An example of an anteriorly located tumor detected as a dark
region in the corresponding VE image. . . . . . . . . . . . . .
114
114
116
118
xvii
List of Abbreviations
2D/3D
ASM
BCCA
B-mode
CC
CF
CNR
CT
CTV
DRE
DSC
ICP
IMMPDA
LDR
MRI
OR
PSA
PTV
RF
RO
ROI
TDPE
TES
TF
TRUS
US
VCC
VE
2-/3-Dimensional
Active Shape Model
British Columbia Cancer Agency
Brightness Mode
Cross Correlation
Coherence Function
Contast to Noise Ratio
Computed Tomography
Clinical Target Volume
Digital Rectal Examination
Dice Similarity Coefficient
Iterative Closest Point
Interacting Multiple Model Probabilistic Data Association
Low Dose Rate
Magnetic Resonance Imaging
Operating Room
Prostate Specific Antigen
Planning Target Volume
Radio Frequency
Radiation Oncologist
Region Of Interest
Time Domain cross-correlation with Prior Estimates
Tapered Ellipsoid Segmentation
Transfer Function
Trans-Rectal Ultrasound
Ultrasound
Vancouver Cancer Center
Vibro-elastography
xviii
Acknowledgements
I would like to thank my supervisor Prof. Tim Salcudean for his support
and valuable guidance, not only throughout my research, but in many other
aspects during my PhD years. I owe special thanks to Dr. Jim Morris
and Dr. Ingrid Spadinger for their insightful feedback in the clinical aspects of my thesis, Dr. Nick Chng for his help and creating a user interface
of the semi-automatic prostate segmentation algorithm, and Dr. Mehdi
Moradi, the lead collaborator in a majority of my work. Special thanks to
Dr. Orcun Goksel for sharing his knowledge in computer programming and
Drs. Hani Eskandari and Reza Zahiri-Azar for their help in the theory of
elastography. Many thanks to my colleagues and the faculty and staff of
UBC, VCC and VGH. I cannot name them all, but here are the name of
some people who helped me during my PhD: Drs. Larry Goldenberg, Tom
Pickles, Mira Keyes, Ted Jones, Michael McKenzie, Rowan Casey, Anthony
Koupparis, Louis Olivier, Nicholas Buchan, Piotr Kozlowski, Ali Baghani,
and Xu Wen, Mr Guy Nir and Mrs Farheen Taquee. Special thanks to
the following institutions for their financial support during the course of my
studies: the University of British Columbia, the Prostate Cancer Foundation
of British Columbia (PCFBC), the National Science and Engineering Council of Canada (NSERC), and the Canadian Institutes of Health Research
(CIHR). I offer my gratitude to my parents whose unconditional support
helped me accomplish many milestones in my life including this. Last but
not least I would like to thank my husband, Mahdi, who has been with me
in every step, and the very new addition to our family, Parsa, who has given
a new meaning to my life.
xix
Chapter 1
Introduction
The prostate is an exocrine gland of the male reproductive system. It weighs
approximately 20gr in young adults and has the approximate shape and size
of a walnut [38]. The normal prostate of an adult has average dimensions of
33mm in length (superior-inferior), 24mm in thickness (anterior-posterior)
and 41mm in width (left-right). This may vary in different geographical
locations or due to malignancies such as tumors or benign prostatic hyperplasia. The prostate is enveloped by a thin but firm fibrous capsule. Two
neurovascular bundles lie posterolateral to the gland which are collections of
vessels and neural structures and are responsible for proper erectile function.
Hence, it is desirable to preserve them during surgery. The superior region
of the prostate, or the base, encircles the bladder neck, while the inferior region of the gland, or the apex, is bounded by the muscles of the pelvic floor.
Posterior to the prostate lies the rectum. The prostate is transversally at its
largest in the superior mid-gland region and tapers down toward the apex.
Originating from the bladder and passing through the gland lies the urethra.
At the midportion of the urethra a raised ridge called the verumontanum
exists and at this location the two ejaculatory ducts join the urethra. At
the distal end of the ejaculatory ducts, at the posterior base of the prostate,
lie the seminal vesicles and vas deferens which provide a majority of the
ejaculate volume.
Based on McNeal’s definition [68] the prostate is divided into three zones:
the central zone, the peripheral zone and the transitional zone. The central
zone, which composes 25% of the gland, is a cone shaped region that surrounds the ejaculatory ducts. The peripheral zone forms a majority of the
gland and lies in the posteriolateral region. The transitional zone surrounds
the prostatic urethra proximal to the verumontanum and forms only 5-10%
of the gland in healthy tissue. This region grows with age and is responsible
for benign prostatic hyperplasia (BPH). Finally, the fibromuscular stroma,
which is composed of only fibrous and muscle tissue and lacks glandular
components, lies in the mid-anterior region of the gland and accounts for
approximately 5% of the gland volume. Figure 1.1 illustrates the prostate,
its periphery and the zones.
Chapter 1. Introduction
Figure 1.1: An illustration of the prostate and its periphery, in the coronal (left) and sagittal (right) views, displaying the bladder (BL), rectum
(R), seminal vesicles (SV) and vas deferens (VD). The three zones are also
shown. The central zone (CZ) surrounds the ejaculatory ducts (ed), the
transitional zone (TZ) surrounds the prostatic urethra (u) proximal to the
verumontanum (V) and the peripheral zone (PZ) forms a majority of the
gland. The fibromascular stroma (FS) lies in the mid-anterior part of the
prostate and lacks glandular components.
Out of various prostate disorders, prostatitis (inflammation of the gland),
benign prostatic hyperplasia (or BPH - enlargement of the transitional
zone leading to difficulty in urination) and prostate cancer can be named.
Prostate cancer is the most prevalent type of cancer among men and is projected to affect 26,500 men in Canada [12] and 241,740 in the U.S. [74] in
2012. It is the most numerous cancer diagnosed in European men (382,000
cases in 2008, as reported by Ferlay et al.) [28]. Prostate cancer is multifocal in nature. According to McNeal et al. [64] 68% of prostate cancer cases
occur in the peripheral zone, 8% in the central zone and 24% in the transitional zone. The small vessels and nerves which enter the prostate from
the neurovascular bundles also create a potential pathway for the spread of
prostate cancer to outside the prosatatic capsule.
Based on the malignancy of the cancer, the patient’s condition (e.g.
age, other diseases, etc), and patient preference, various cancer treatment
methods exists which may be performed alone, or in combination. Some of
these methods include:
• Watchful waiting: Closely monitoring the patient’s condition without
2
Chapter 1. Introduction
giving treatment until symptoms change. It is used in cases where
prostate cancer has a slow growth.
• Prostatectomy: Removal of the prostate, either by laparoscopic or
open surgery. The surgical margin within which the prostate is removed depends on the stage of the disease.
• External beam radiation therapy (EBRT): Utilizing external radiation
focused on the prostate to kill cancer cells. It is recommended for
patients with a significant risk of extraprostatic cancer [69].
• Low dose rate brachytherapy: Generally used for early stage, intracapsular prostate cancer and has rapidly gained acceptance due to its
highly successful clinical results [10, 36, 66, 71, 124]. In this treatment,
small radioactive seeds (Iodine-125 or Palladium-103) are inserted
through the perineum and permanently implanted into the prostate
and periprostatic tissue.
• High dose rate (HDR) brachytherapy: Temporary catheters are placed
inside the prostate which allow treatment with a high dose rate radioactive source, usually iridium-192, for 5-15 minutes, delivering the
radiation treatment over a number of fractions, typically over a few
days.
• Hormone therapy: Lowering testosterone levels to stop the growth of
cancer cells.
• Chemotherapy: Given as pills or by injection, it is usually used when
cancer is in its advanced stage.
• High-intensity focused ultrasound (HIFU): The use of ultrasound to
destroy cancer cells.
Currently, common early detection methods are screening the prostate
specific antigen (PSA) levels and digital rectal examination (DRE). PSA is a
protein produced by normal prostate cells. A rise in levels of PSA (measured
in ng/mL) is associated with prostate cancer, benign prostatic hyperplasia
or acute bacterial prostatitis. DRE is a physical examination in which the
prostate is palpated through the rectum for bumps or abnormalities. Both
methods are not accurate and further tests are required. Further confirmation of cancer can be obtained by biopsy (performed either transrectal,
transurethral or transperineal). However, false negative rates are relatively
high in prostate biopsy.
3
Chapter 1. Introduction
Figure 1.2: MRI images of the prostate in the coronal (left) and transverse
(right) views displaying the outline of the prostate (P- dashed line), the
bladder (BL) and the penile bulb (PB). In the transverse image an endorectal
coil (EC) is placed inside the rectum for higher image quality, better showing
the central zone (CZ), peripheral zone (PZ) and neurovascular bundles (NV).
Commonly used modalities for imaging the prostate include transrectal
ultrasound (TRUS), computed tomography (CT), and magnetic resonance
imaging (MRI). Figure 1.2 shows the coronal and transverse view of the
prostate in MRI and Fig. 1.3 displays transverse ultrasound B-mode images
of the prostate at two levels: the base and the mid-gland.
1.1
Motivation
Effective localized treatment of prostate cancer, regardless of the treatment
method used, requires accurate visualization of the gland and its surrounding
region.
Reliably detecting the prostate boundary can aid proper removal of the
prostate in prostatectomy. In prostate radiation therapy the prostate boundary is required both prior to treatment, during planning and delivery and,
in the case of low dose rate brachytherapy, after treatment for dosimetric
assessment. Accurate delineation of the prostate during treatment planning
and delivery has the potential to reduce some of the possible side-effects of
the current treatment methods. These complications include urinary incontinence, impotence, and damage to the rectum and urethra [107].
Ultrasound is the most commonly used modality for imaging the prostate.
This is due to its availability, safety and ease of use. However, ultrasound
4
Chapter 1. Introduction
Figure 1.3: Transverse ultrasound B-mode images of the prostate at the base
(left) and mid-gland (right). The seminal vesicles (SV) and the shadow of
the bladder (BL) are visible at the base. At the mid-gland, the location of
the transrectal ultrasound (TRUS) probe inside the rectum (R) is labeled.
The various zones of the prostate (P) are not always distinguishable in TRUS
images. However, in this case the neurovascular bundles can be faintly seen
at locations 1-B and 1-F of the grid.
B-mode images do not always delineate the prostate reliably. As a result,
manual extraction of the prostate boundary becomes a highly subjective
process, depending on the experience, skill, and technique of the contourer.
This is particularly evident in the base and apical regions of the gland, where
the prostate boundary blends into the surrounding anatomy in TRUS images, becoming more ambiguous and difficult to delineate. Consequently,
there is substantial inter-observer variation in the contours for a given patient in these regions. Smith et al. [98] have evaluated the reproducibility
of prostate contouring, after brachytherapy implants, using 3DTRUS, T2weighted MR and CT imaging. Prostates from 10 patients with early-stage
prostate cancer (T2b or less) were segmented twice by seven observers. Their
results showed high contouring variability of the anterior base and apex in
3DTRUS, whereas the prostate-rectum interface had the smallest variability.
A study by [15, 108] showed that prostate volume measurement by transrectal ultrasound may vary among observers when patients have large prostates
(≥30 cm3 ). The average volume difference between 101 prostates measured
by two experienced observers was reported as 6.00 cm3 for prostates with
a mean measured volume ≥30 cm3 and 1.51 cm3 for prostates with a mean
measured volume ≤30 cm3 . These numbers increased to 6.84 cm3 and
5
Chapter 1. Introduction
3.99 cm3 , respectively, when measurements were performed by one experienced and one less experienced observer (110 prostate volumes measured
in this case). Manual contouring is also time consuming. This is a major
drawback if an intraoperative treatment procedure is to be performed.
The motivation for an automatic segmentation method is to reduce user
dependence of contouring, replace a repetitive and time-consuming task, and
develop a technique which is fast enough to be used intraoperatively. Automated contouring can also be integrated into other intraoperative treatment
methods such as for LDR brachytherapy where imaging the prostate, segmenting the boundary and creating the treatment plan can all be carried
out inside the operating room.
Improving the visibility of the prostate gland can further aid segmentation. To this end, various attempts have been made and reported in the
literature. These vary from processing ultrasound B-mode images [85, 89],
to the use of additional information from other modalities such as MRI [20].
Another recently developed option is the use of ultrasound elastography.
Elastography has the ability to highlight differences in the mechanical properties of tissue, such as stiffness, therefore, it can be used for highlighting
the prostate boundary since the gland and its periphery are known to have
various stiffness values. However, it is important to evaluate the appearance of the gland in such images to determine the usefulness of these images
in diagnosis and treatment of malignancies. For example, the volume of
the prostate is an important measure for post-operative assessment in LDR
brachytherapy. To be able to use elastography images for segmenting the
prostate, it is required that the gland volume extracted from the images is
valid, either being similar to that shown in a standard modality, such as
B-mode, or the difference being known.
Furthermore, the ability to distinguish internal structures of the gland,
especially cancerous tumors, can aid diagnostic and treatment methods.
Guided prostate biopsy, with knowledge of lesion locations, can reduce the
false negative rate in detecting cancer and the need for repeating the procedure. Additionally, localized treatment (e.g. surgical removal of the gland or
LDR brachytherapy) as opposed to general treatment can be pursued when
the location of disease is known. This can potentially reduce treatment
side-effects.
In this thesis, which is based on the importance of developing tools
which aid diagnosis and treatment of prostate cancer, two main topics are
studied: (i) Segmentation of the prostate boundary and (ii) Visibility of
the gland in elastography images. Chapters 2 and 3 present a prostate
segmentation tool and Chapter 4 addresses the visibility of the gland in
6
Chapter 1. Introduction
vibro-elastography images. In Chapter 5 the methods and results of the
aforementioned chapters are utilized to develop an automatic segmentation
tool and finally, in Chapter 6 initial results of the appearance of prostate
cancer in vibro-elastography images is reported and compared to other work.
1.2
Background
A literature review on various prostate segmentation methods is first provided in Section 1.2.1. Our proposed prostate segmentation algorithm is
applied to and evaluated for prostate brachytherapy, hence, an introduction
to low dose rate prostate brachytherapy is deemed necessary and is provided
in Section 1.2.2. As an attempt to improve the visibility of the prostate, ultrasound elastography, more specifically vibro-elastography, which has been
proposed in our group, is used for imaging the prostate. Section 1.2.3 provides an introduction and a brief literature review of ultrasound elastography and describes the data collection apparatus and the procedure for this
method.
1.2.1
Prostate Segmentation
Various prostate segmentation methods have been proposed in recent years [76].
These can be classified into those which solely rely on image data, and
those which incorporate prior information about the expected shape of the
prostate. The advantage of not constraining the solution to certain geometrical classes may result in more robust segmentation of irregular prostate
shapes and some prostate surface abnormal deformations that may be caused
by cancer. Examples of image-based segmentation methods used in the
literature are contrast enhancement, image smoothing and Canny’s edge
detector [85], multi-resolution discrete wavelet pyramids for coarse-to-fine
segmentation [11] and the use of phase symmetry for outlining dominant
edges which are then linked to create a final contour [123].
The disadvantage of methods that rely solely on image information is
that they are more sensitive to factors such as image quality and noise, and
generally require more effort to achieve reasonable results. Because they do
not use prior shape data, a majority of the proposed general segmentation
methods do not work on ultrasound images of the prostate. Deformable
models such as active contour models (ACM) and snake models have been
used widely for medical image segmentation. They are curves or surfaces of
which the deformation (e.g. with the means of Gabor filters) is guided by
internal forces (e.g. the constraint of the curve being smooth) and external
7
Chapter 1. Introduction
forces (e.g. edges in the image). Because they can not tolerate the large
amount of noise in ultrasound images very well, for such methods to work,
additional constraints must be imposed. For example, the deformation must
be limited by an a priori shape, or alternatively, significant user interaction
must be used [76].
Active shape models (ASM) [18] use shape models that deform within
some constraints. These constraints and the initial shape model are derived,
statistically, from a training set. [25, 39, 42, 45, 55, 72, 73, 87, 95, 105, 119,
125, 126]. Active shape models are used in the 2D segmentation method
proposed by Hodge et al. [39]. This method is extended to 3D by using a
rotational-based slicing method. After each 2D image is segmented, manual editing (required in 26.3% of the 2D images in each case) is allowed,
after which the modified points are clamped and the model is re-deformed.
They report an average mean absolute distance and maximum distance of
1.09 ± 0.49 mm and 7.27 ± 2.32 mm between manual (the average of three
repeated manual delineations performed by three trained graduate students)
and automatic contours (on 36 patient data sets). The average run time is
reported to be 6.9 ± 2.1 minutes which includes the time needed for manual
editing of the 2D images.
Yan et al. [117] use adaptive learning during segmentation, in addition
to global population-based shape statistics for 2D segmentation of the entire
gland in TRUS video sequences. Their adaptive learning method is reported
to improve segmentation of the base and apex. They report a mean absolute
distance of 1.65±0.47 mm and a Dice similarity coefficient value of 0.91±0.03
between automatic and manually delineated TRUS video sequences of 19
patients.
The use of ellipses, ellipsoids, superellipses and other similar shapes
has been a relatively attractive approach for prostate segmentation as most
prostates conform well to these representations, and the reduced parameter
space results in fast segmentation algorithms. In the work by Kachouie et al.
[46] the velocity of an evolving ellipse initially placed on the center of the
prostate is guided by the first and second moments of the Gaussian probability density function fitted to the image histogram. However, further work is
said to be needed for robustness to poorer images. Gong et al. [33] have used
deformable superellipses; ellipses that can be deformed by adjusting their
aspect ratio, squareness, tapering and bending. Their extensive comparison
of manual and computer-generated contours on 125 prostate ultrasound images resulted in a mean error of less than 2 mm between computer generated
and manual contours. The segmentation duration was reported to be less
than 5 seconds for each 256 × 256 image on a Pentium 4 PC running at 2
8
Chapter 1. Introduction
GHz. However, this method only generates 2D contours.
In the method proposed by Badiei et al. [6] an ellipsoidal volume is generated. The method is based on fitting ellipses and ellipsoids with the aid of
initial points defined by the user and the Interacting Multiple Model Probabilistic Data Association (IMMPDA) [2] edge detection technique. It has
the benefit of being fast, owing to solving a convex problem, but often results in relatively large false positive regions mainly in the anterolateral and
inferior regions of the gland. This is due to the tapering of most prostates
both in the 2D transverse plane and in the 3D volume along the main axis
from the base to the apex.
The need for using 3D models as opposed to 2D is specifically seen when
the quality drops in some image slices, specially at the base and the apex.
With the use of a 3D model and the higher quality mid-region images, contours can be generated for the lower quality base and apical images based
on the shape of the mid-gland. Additionally, surface smoothness from one
image to the other is more easily maintained. Amongst the many proposed
methods for segmentation of the prostate in ultrasound images, only a small
number offer 3D segmentation of the gland [3, 39, 41, 86, 113, 118, 125]. The
method proposed by Penna et al. [86] uses Fourier ellipsoids as the 3D surface model. This method requires 90 seconds to create the prostate surface
model and generate the solid models necessary for HIFU therapy planning.
Of this duration, only 10 seconds are required to model the surface. Manual
tracing of approximately 5 transverse and 3 sagittal images of the prostate
is needed to initialize this algorithm, which, as reported, requires another
1.5 minutes and introduces operator variability that has not been quantified. In the work of Tutar et al. [113], 3D semi-automatic segmentation is
based on fitting the best surface to a set of images under shape constraints.
These constraints are derived by modeling the shape of the prostate using
spherical harmonics. A measure of percent volume overlap (the intersection
divided by the union of two volumes) between automatic and manual (average of three observers’ manual contours) of 83.5% is reported on a data
set of 30 patients, while the inter-observer variability is 82.8%. However,
this method also requires the user to manually segment the mid-gland axial
and sagittal contours for initialization. After initialization, the 3D prostate
shape is identified in 1-4 minutes, of which the authors report initialization
times of approximately one minute. In [125] Zhan et al. use a deformable
model based on statistical matching of both texture and shape on 3D TRUS
images (256 × 256 × 176 volume size, 0.312 mm/voxel). They report a
mean overlap volume error of 4.16% compared to manual segmentation on
six patients and a segmentation time of 3 minutes. Since their method is
9
Chapter 1. Introduction
applied to 3D TRUS images, how it would fare relative to a conventional
prostate volume acquired with a brachytherapy stepper, where the slices are
at 5 mm intervals (resulting in 8-12 images per case), is not discussed in the
paper.
To conclude this section, it is worth noting that when evaluating a segmentation algorithm, it is important to measure the intra- and inter-observer
variability in addition to accuracy. Unfortunately this is not commonly reported in the literature. These two measures provide information about the
repeatability of a method when various individuals perform segmentation
or an individual produces contours at different times. It also indirectly provides information about the data set used for evaluation. Authors each have
their own set of images for assessment of their algorithm, hence direct comparison of reported accuracy numbers is not reasonable. Obviously a high
quality set of images provides better results compared to a lower quality
set. One set may have more challenging images than another and therefore,
lower accuracy results. If the accuracy of a method, compared to the gold
standard (usually manual contours) falls within the inter- and intra-observer
variability measured from the same data set, this method can be regarded
as a suitable replacement for manual contouring (needless to say, in terms
of accuracy).
In this thesis, the presented 3D semi-automatic segmentation algorithm
is based on fitting warped and tapered ellipsoids as the shape prior. Volumetric and dosimetric evaluation of this method as applied to LDR brachytherapy is then performed. Also, intra- and inter-observer variability of the
method is studied. This method is partially automated with the use active
shape models and stiffness information from vibro-elastography images.
1.2.2
LDR Prostate Brachytherapy
LDR brachytherapy is generally recommended as a treatment option for
early stage and organ confined cancer [69]. Details of the procedure may
differ among various centers, therefore, a more general description will be
provided here and details specific to the method at the Vancouver Cancer
Center (VCC) will be given in Chapter 3. Prior to the procedure, a volume
study is performed in which transverse ultrasound B-mode images of the
prostate are collected, usually in the transverse plane. A grid with marker
spacing of 0.5 cm, a horizontal axis labeled as ’A a B a · · · G’ and a vertical axis labeled ‘1.5 2.0 2.5 · · · 6.0’ is overlaid on the images such that
the prostate is symmetrical around the mid-line D and the most posterior
aspect of the gland lies on or slightly below one of the horizontal grid lines
10
Chapter 1. Introduction
(usually 1.0, 1.5, or 2.0). The boundary of the prostate is then identified
and a treatment plan is generated based on the prostate boundary. In some
centers (such as VCC) a planning target volume (PTV) is created from
the delineated prostate boundary (also known as the clinical target volume
(CTV)), and is used to generate the treatment plan. The treatment plan
specifies the locations of 50-160 radioactive seeds (depending on the PTV
volume) placed selected grid points. Most or all of the seeds are placed at
grid points falling within the boundary of the PTV. Two commonly used
radioisotopes are 125 I [36] and 103 Pd [10]. The goal is to cover the PTV with
the 100% isodose surface while maintaining tolerable dose to sensitive tissues such as the urethra and rectum, thereby avoiding complications such as
urinary incontinence and erectile dysfunction. However, since an ideal dose
distribution is not physically achievable, guidelines for allowable exceptions
and maximum margins outside the PTV are given.
Commonly used parameters for evaluating the dose distribution are the
V100 , V150 , V200 , where for example a PTV V100 of 98% means that the isodose
surface providing 100% or more of the prescribed dose (a prescription being
144 Gy) covers 98% of the PTV. Ideally, the PTV V100 should be 100%
whereas to prevent overdose, the V150 should typically be around 50-60%.
Other parameters are the D80 and D90 which indicate the minimum dose
recieved by 80% or 90% of the PTV. The shape of the isodose surface is
also important. For example, at some centers, the 150% isodose generally
has an upright horseshoe shape with very little margin at the posterior (to
avoid overdosing the rectum) and the hollow region around the urethra.
The seeds are delivered via pre-loaded radioactive needles (typically, for
monotherapy, 144 Gy for 125 I and 115-120 Gy for 103 Pd) through the perineum and guided by the TRUS probe and a template with guide holes that
match the grid that was overlaid on the images during planning, placed in
front of the perineum. Fluoroscopy may be used sparingly throughout the
procedure to evaluate the position of the implanted seeds and evaluate the
implantation process.
If treatment plans are created pre-operatively, the position of the TRUS
probe is adjusted intra-operatively to provide images of the prostate similar to the pre-operative images. This is necessary due to the fact that
the orientation (mainly due to patient positioning), shape and size of the
prostate may change from pre-operative to intra-operative imaging sessions.
Occasionally, despite the adjustments, the treatment plan may need to be
updated by the physician by adding or relocating seeds to adhere to the
changed prostate boundary. The adjustment of the images and the plan
are usually performed manually, which is time consuming and not free from
11
Chapter 1. Introduction
error. These adjustments are of less concern if intra-operative planning,
which has been established and is being used in some centers, is performed
[1, 4, 71].
Many factors cause the implanted seed positions to deviate from the
treatment plan. These include motion of the prostate during implantation
of the seeds and changes in the size, position and shape of the prostate, which
are effected by inflammation and progress of the disease (although this is not
common), position of the patient between the volume study and treatment
sessions, and pressure of the TRUS probe. Factors such as interference of
the pelvic bone may also require the plan to be changed. As a result, postimplant dosimetry is recommended by the American Brachytherapy Society
for all patients undergiong LDR permanent brachytherapy [70]. This is
required to identify regions where less or more than the prescribed dose is
delivered to the tissue. These regions mainly include the prostate boundary
(for coverage) as well as the anterior rectum and prostatic urethra (for rectal
and urinary morbidities), the neurovascular bundle and penile bulb (for
erectile dysfunction) [19]. The standard imaging modality for post-implant
assessment is CT where the seed locations, prostate and peripheral organ
boundaries are defined. Parameters such as the prostate D100 , D90 and V100
are then obtained and compared to the intended dose distribution.
Figure 1.4 shows a snapshot of the software used at VCC for LDR
brachytherapy treatment planning (VariSeedT M , Varian Medical Systems
Inc., Palo Alto, CA, USA) and 10 B-mode images of a sample prostate
taken at various depths. In the selected image, at depth 5, the planning
grid is overlaid on the image. Both the contours and the plan are symmetric
about the D column of the template through which the center of the ultrasound probe passes. Guided by the CTV (red) and PTV (cyan), needles
(hollow circles) and seeds (cyan) are placed at grid locations. Isodose lines
are calculated and displayed based on the prescribed dose and the location
of the seeds, here showing the isodose lines corresponding to 200%, 150%
and 100% of the prescribed dose (in this case 144 Gy). Notice that the seeds
are placed such that large doses are avoided in the center of the prostate,
where the urethra lies, and in the posterior, close to the rectum.
12
13
Chapter 1. Introduction
Figure 1.4: A snapshot of the software used for treatment planning (VariSeedT M ), showing the B-mode images of
the prostate taken at 10 depths (smaller images at right). In the larger image, the CTV (red) and PTV (cyan),
needle locations (hollow circles), and the seed locations (cyan circles) are shown. The 200%, 150% and 100%
isodose lines are also displayed.
Chapter 1. Introduction
1.2.3
Ultrasound Elastography
Elastography [78, 79] is a promising technique for imaging soft tissues, it
relies upon measuring tissue strain in response to excitation, and and is
based on the fact that, when compressed, softer tissue will compress more
than stiffer tissue, and therefore experience larger strain. Ultrasound, as a
means for measuring tissue deformation, is a low-cost, portable, real-time
and safe modality and therefore, widely used. Another popular modality for
measuring tissue motion is MRI [67]. MR elastography (MRE) offers high
sensitivity, a large field of view, and the ability to obtain full 3D displacement information throughout the 3D volume. Its applications are vast and
still growing [44, 97, 120]. However, our main focus in this thesis will be
ultrasound-based elastography.
Tissue motion can be achieved using external excitation, such as the
inward motion of the ultrasound transducer in the image axial direction
(usually strain levels of less than 2%), or internal excitation, such as from
the natural motion of the cardiac muscle in myocardial elastography [52].
The type of excitation can be static/quasi-static, in which the tissue is
slowly compressed and only the static behavior of the tissue, the stiffness,
is measured [37]. In transient (impulsive or short tone burst) elastography,
the propagation of a shear wave is imaged with parallel receive ultrasound
[14, 103]. When the excitation is dynamic, the shear modulus can be estimated as a complex function of frequency, thus providing information on
the viscoelastic properties of tissue [40].
During tissue excitation, ultrasound radio-frequency (RF) signals are collected and tissue motion is estimated. This is done in various ways such as
cross-correlation techniques [79], peak searching methods [26], zero-crossing
[100], and speckle tracking [77]. Alternatively, tissue vibration induced by
the exciter can be measured with Doppler ultrasound [56], with larger vibrations corresponding to softer tissue. Motion estimation can be performed
in 1D, in the direction of beam propagation where ultrasound resolution is
higher, or the lateral and/or elevational directions resulting in 2D and 3D
estimation [51, 53, 58, 103]. For more details related to the basics of elastography, various excitation techniques, and computation methods we refer the
reader to [30, 35, 40, 80, 84, 88]. A concise review on the fundamental principles of elastography, developments in various categories of elastographic
imaging in the last 20 years and recent clinical results, can be read in [83].
Clinical applications of elastography include, but are not limited to,
breast lesions [32, 47, 57, 94, 97], liver fibrosis [13, 43, 120], vascular vulnerable plaque [93], elastic properties of skeletal muscle [23], thyroid gland
14
Chapter 1. Introduction
Figure 1.5: Schematic representation of vibro-elastography imaging. From
left to right, a broadband external mechanical excitation x0 (t) is applied to
the sample; the axial motion xjk (t) of a tissue block indexed axially by j
and laterally by k within the ultrasound image is obtained from RF images
at time t and t − ∆t; transfer functions from a reference to each of the tissue
blocks are estimated, with the illustration showing the reference being the
motion xjk (t) of block jk; images are produced by displaying the L2 -norm
of the difference in transfer functions.
tumors [59], assessment of thermal tissue ablation [116], and detection of
prostate cancer [16, 29, 34, 48, 65, 82, 99, 127]. A recent review on clinical
applications of elastography is presented in [31].
Vibro-Elastography
Vibro-elastography is a dynamic ultrasound elastography method which
models viscoelastic properties of tissue [112]. The approach is illustrated
in Fig. 1.5. The technique relies on the continuous real time acquisition of
unprocessed ultrasound echo data as a time series of ‘radio-frequency’ (RF)
data images, while, simultaneously, tissue is externally vibrated with a mechanical excitation that generates low-frequency broad-band compression
waves. A time-series of tissue displacement or strain images is computed
from consecutive RF data images. The tissue displacement as a function
of time at a given spatial location can be viewed as the output of a linear
system whose input is the motion of the exciter as a function of time. Therejk
(ω), that relates
fore a frequency response or Transfer Function (TF), Href
15
Chapter 1. Introduction
the tissue motion at any spatial location, xjk , with the exciter motion as a
reference, xref , can be computed in the frequency domain. Alternatively,
if the exciter motion is not measured, a tissue region, typically in the focal
area of the ultrasound beam, can be selected as the reference. The transfer function can provide elastic
properties
of tissue such as viscosity and
jk
stiffness. At zero frequency,Href (0) gives the amount of axial tissue deformation at location xjk as a response to unit displacement at location xref .
In order to reduce the displacement estimation errors, usually an average
of the TF at low frequencies is computed instead of at zero frequency. In
this sense, computing the stiffness based on the TF is superior to the static
method. The reliability of the transfer function can be assessed using the
jk
(ω). Its value, being between 0 and
corresponding coherence function, Cref
1, indicates what portion of the input energy at location xref at frequency
ω appears as the output energy at location xjk at frequency ω. A value of
1 indicates that the system is linear and the signal-to-noise ratio is high.
Figure 1.6 shows the transrectal ultrasound (TRUS) VE actuation system developed in our research group for VE imaging of the prostate. A LDR
brachytherapy stepper (EXII, CIVCO Medical Solutions) was modified to
enable acquisition of RF data volumes. The rotation of the cradle was motorized to enable smooth rotation of the TRUS transducer and a shaker was
mounted on the transducer cradle in order to radially vibrate the TRUS
probe. Figure 1.7 displays a close-up view of the motorized cradle. A PCbased control interface allows the user to control the rotation motion range
(-45 to 50 degrees, 0 degrees being the mid-sagittal plane), the amplitude,
and frequency range of the transducer vibration. Synchronized with the
motion of the probe, frames of RF data are collected from a SonixRP ultrasound machine with the sagittal array of a dual-plane linear/microconvex
broadband 5 - 9 MHz endorectal transducer (Ultrasonix Medical Corp.).
The RF sampling frequency is 20 MHz.
For the collection of a majority of our data, the vibrating TRUS probe,
with a vibration range of 2-10 Hz and vibration amplitude of approximately
1 mm, was swept from −45◦ to 45◦ covering a fan of longitudinal planes
passing through the axial transducer axis with a predefined angular separation. At each angle, one B-mode image and frames of RF data were collected
continuously from which vibro-elastography images were created. Figure 1.8
displays the orientation of the fan of longitudinal planes at which data is collected, with respect to the prostate and TRUS probe. These planes, which
we will simply refer to as ‘longitudinal’, pass through the transducer axis
and include the mid-sagittal plane at angle 0◦ .
16
Chapter 1. Introduction
Figure 1.6: System setup for vibro-elastography data acquisition
Figure 1.7: Details of the motorized cradle
Further details of the data collection and processing will be explained in
Chapter 4.
17
Chapter 1. Introduction
Figure 1.8: An illustration of the longitudinal planes at which data is collected, with respect to the prostate and TRUS probe. Here, the mid-sagittal
plane, at 0◦ , and two other planes at angles of -20◦ and +10◦ are shown.
The orientation of a transverse plane in the mid-gland is also displayed.
Appearance of the prostate in vibro-elastography images
Figure 1.9 shows a mid-sagittal and transverse image of the prostate as
seen in VE. In general, the prostate appears as a darker (stiffer) structure
compared to its periphery, with the prostate apex exhibiting more stiffness
contrast with respect to the surrounding tissue. On the other hand, the
base of the prostate has less contrast. This could be due to the anatomy of
this region or the base merging into the bladder. Another possibility is that
since the base is imaged using the distal end of the probe and the probe
in motion acts like a cantilever, the desired excitation frequency may not
properly propagate to the end of the probe. Hence, tissue in this region
is not excited with the desired motion. Similar to B-mode images, it is
sometimes possible to see the path of the urethra. This is seen as a softer
(brighter) region passing from the base to the apex. A contrast between
central and peripheral region of the prostate can be seen in some cases.
18
Chapter 1. Introduction
Figure 1.9: A mid-sagittal (left) and transverse (right) VE image of the
prostate showing the urethra (dashed red). The prostate boundary is marked
by the dashed yellow lines. The section from which the transverse image is
displayed is marked by the dashed blue line.
1.3
Thesis Objectives
The specific objectives of this thesis are as follows:
1. Developing a clinically acceptable prostate segmentation method.
2. Evaluating ultrasound vibro-elastography as an imaging method for
visualization of the prostate boundary.
3. Incorporating ultrasound vibro-elastography in the proposed B-modebased prostate segmentation method to enhance segmentation.
4. Studying the feasibility of using ultrasound vibro-elastography as a
tool for detecting prostate cancer.
In the course of achieving these objectives, the following contributions
were made:
• A semi-automatic prostate segmentation tool was developed and is
currently being used as part of the clinical procedure for prostate
brachytherapy at Vancouver Cancer Center. It has also been introduced to other centers of the BC Cancer Agency.
• The proposed segmentation algorithm has been thoroughly evaluated
both volumetrically and dosimetrically for nine regions of the prostate
and has shown to be suitable for replacing manual delineation, which
is currently the standard.
19
Chapter 1. Introduction
• Ultrasound vibro-elastography (VE) has been evaluated for prostate
imaging. The shape and volume of the prostates shown in VE have
been compared to ultrasound B-mode and with MRI as the goldstandard. Additionally visibility of the prostate boundary and the
contrast of the gland in VE images have been compared to those of
ultrasound B-mode.
• The fusion of information from vibro-elastography and B-mode images
has led to the development of a fully automatic prostate segmentation
algorithm.
• Groundwork for studying the visibility of prostate cancer in vibroelastography images has been laid.
c was the main environment in which algorithms are impleMATLAB
c was also used for processing the RF data.
mented. Visual C++
1.4
Thesis Outline
This chapter introduced the research topic, motivation for performing this
research and the thesis objectives and contributions. Related background
accompanied by an overview of the current literature of prostate segmentation, prostate brachytherapy, ultrasound elastography and prostate cancer
detection methods has also been presented.
In Chapter 2, a semi-automatic segmentation method for delineating
the prostate in ultrasound B-mode images is described. Volumetric accuracy evaluation of this method is performed by comparing semi-automatic
contours with manual contours and physician corrected semi-automatic contours. Repeatability of the method is assessed by computing inter- and
intra-observer variability in semi-automatic segmentation and comparing
the results to inter- and intra-observer variability in manual contouring.
The proposed semi-automatic segmentation algorithm is clinically used
for delineating the prostate in LDR brachytherapy. Since the goal of prostate
brachytherapy is to deliver a sufficient amount of radioactive dose to the
diseased prostate while sparing critical regions, it is important to assess the
dosimetric impact of using the semi-automatic segmentation algorithm, as
opposed to conventional manual contouring. This is the subject of Chapter 3.
In Chapter 4, visualization of the prostate in ultrasound vibro-elastography
images is evaluated. Ultrasound B-mode imaging is well known for its poor
20
Chapter 1. Introduction
quality, especially at the prostate base and apex. Methods that improve
visibility of the gland can greatly aid delineation and cancer detection. In
this chapter delineated prostate surfaces from vibro-elastography images are
compared to those of B-mode with MRI as the gold standard. The quality
of edges and contrast of the images is also compared to those of B-mode.
Based on results obtained from Chapter 4, we concluded that vibroelastography images are superior to B-mode in some aspects. In Chapter 5 the combination of information from ultrasound B-mode and vibroelastography is used for developing a fully automatic 2D prostate segmentation algorithm, which is then used to automate the 3D semi-automatic
segmentation algorithm proposed in Chapter 2.
Finally, in Chapter 6 we provide preliminary results of the visibility
of cancer in vibro-elastography images. We suggest that although, at this
stage, VE alone can not always be used for detecting cancer, it can potentially be combined with other image types for developing an automatic
cancer detection method. This final chapter does not report peer reviewed
contributions but provides tools and the foundation for future work in the
area of cancer research.
21
Chapter 2
Semi-automatic
Segmentation of the
Prostate:
Method Description and
Volumetric evaluation
2.1
Introduction
In this chapter a semi-automatic prostate segmentation method for use in
LDR prostate brachytherapy is proposed.
For an algorithm to improve the efficiency and consistency of the brachytherapy treatment, it needs to satisfy the following important requirements.
First, it should be able to produce contours that are not distinguishable
from those generated by medical practitioners. Second, it should produce
contours that are amenable to the design of uncomplicated treatment plans.
Third, it should not require considerable changes to the conventional clinical
procedure. Moreover, a fast and automatic solution to manual contouring
would greatly facilitate intra-operative planning, where additional and significant gains in treatment quality are likely to be realized.
In the proposed method, the boundary of the prostate in TRUS images
is delineated based on prior knowledge of the shape of the gland, resulting
in smooth, symmetric and less user dependent contours. The 3D geometric
model of the prostate is created based on the assumption that the prostate
has a tapered ellipsoidal shape and is slightly warped posteriorly due to the
presence of the TRUS probe. Using an a-priori shape – in this case, a tapered, warped ellipsoid – aids segmentation in the less visible base and apex
of the gland. Manual initialization of the algorithm makes use of the practitioner’s experience. The simplicity of the algorithm and the formulation
Chapter 2. Semi-automatic Segmentation of the Prostate
of the most intensive part of the computation as a convex problem results
in a fast, close to real-time, and repeatable segmentation method.
Symmetry and smoothness are two other desired features which greatly
aid treatment planning. The use of symmetric contours with respect to
the sagittal plane for treatment planning is an accepted practice at the BC
Cancer Agency and some other centers. Symmetric contours lead to simple treatment plans that are also easier to change to ensure adequate dose
coverage, should the shape, size and position of the prostate change significantly with respect to the volume study. A population study of biochemical
and survival outcomes on a total of 1006 consecutive prostate brachytherapy interventions (July 1998 to October 2003) performed at the BC Cancer
Agency, the institution where we collected the patient data, shows excellent
oncologic outcomes for men with favorable-risk disease, with very low rates
of biochemical or metastatic recurrence [66].
2.2
Methodology
Given a series of 2D trans-rectal transverse B-mode ultrasound images of
the prostate, from the base to the apex, the goal is to generate a 3D volume of the gland in a real (or close to real) time, preferably using the least
user interaction. The images are collected using a B&K Pro-Focus System
B-Series machine (BK Medical, Peabody, MA), with the MFI Biplane Transducer. The image specifications are: image size 640 × 480 pixels, 0.155 mm
× 0.155 mm 2D image pixel size, 5 mm image spacing. It is assumed that
the prostate is positioned in the TRUS images such that symmetry is maintained with respect to the mid-sagittal plane. This is usually met since the
current treatment planning procedure requires such a condition. If such a
condition is not met, an initial rotation and translation can be easily applied
to the images.
TRUS images of the prostate show posterior warping of the gland due to
the presence of the TRUS probe. Additionally, tapering is often seen both
transversally, toward the anterior aspect of the gland (creating a pointed
appearance), and axially, with the gland narrowing toward the apex. Our
method is based on the assumption of separability between the probe-induced
warping, the tapering, and ellipsoidal fit, so that each can be handled independently. In doing so, the fitting and hence, the segmentation problem,
is simplified into the convex problem of fitting an ellipsoid to the preprocessed data. Figure 2.1 shows a snapshot of the graphical user interface used
in the Vancouver Cancer Centre for prostate contouring using the proposed
23
Chapter 2. Semi-automatic Segmentation of the Prostate
method. Figure 2.2 illustrates the main steps of the segmentation algorithm.
Details of each step are presented below. We will hereafter and in Chapter 3
refer to the algorithm and the resulting contours as the ‘Tapered Ellipsoid
Segmentation algorithm’ (TES) and ‘TES contours’.
24
Chapter 2. Semi-automatic Segmentation of the Prostate
25
Figure 2.1: A snapshot of the graphical user interface for automatic prostate segmentation used at the Vancouver
Cancer Centre.
26
Chapter 2. Semi-automatic Segmentation of the Prostate
Figure 2.2: The main steps of the segmentation algorithm.
Chapter 2. Semi-automatic Segmentation of the Prostate
2.2.1
Algorithm
Initialization
A prostate volume study used for brachytherapy treatment planning consists
of a set of transversal prostate images. Selection of the mid-gland, apex, and
base images is the first step in the semi-automatic segmentation algorithm.
The mid-gland image is where the initial 2D segmentation will be carried
out. It contains the largest and most visible section of the gland. The
base and apex are the extreme superior and inferior images of the prostate,
respectively. They are not always visible in TRUS images, and will not be
used in the TES algorithm. Their contours will be defined by projection of
the final 3D shape on the respective planes. However, the depth of the base
and apex images is important in order to extend the segmentation from 2D,
on the mid-gland, to 3D. The term ‘slices’, hereafter, will refer to the images
starting from image base + 1 to image apex − 1.
On the mid-gland image, the user selects six boundary points which,
in addition to p1 - the TRUS probe center, initiate the algorithm. These
boundary points are: p2 -lowest posterior lateral, p3 -extreme right, p4 -midposterior, p5 -mid-anterior, p6 - intersection of the mid-perpendicular line
between p2 and p3 with the boundary, and p7 - intersection of the midperpendicular line between p3 and p5 with the boundary (Fig. 2.3a). Points
p3 , p4 , and p5 are primarily determined by the size of the prostate. Points
p1 , p2 , p4 are used to calculate the amount of warping. Points p6 and p7 have
the main contribution in determining the amount of tapering. The aim is to
extract the most information from the image while keeping the variability
of the point selection low by directing the user to specific regions; either to
extremes (e.g. points p3 , p4 , and p5 ) or by guiding lines (points p6 and p7 ).
These points along with their symmetric reflections across the medial line
will be referred to as the ‘initial points’.
Image un-warping
Based on the initial points, the posterior region of the mid-gland image
is un-warped to reduce the deformation caused by the TRUS probe using
Eq.2.1 below. In this equation r is the current distance of an image pixel
on a radial line starting from the probe center with angle φ ( φ = 90◦ being
the medial line) and rnew is the distance of the re-located pixel. In this
sinusoidal Gaussian function, the maximum deformation is achieved when
27
Chapter 2. Semi-automatic Segmentation of the Prostate
p5
p7
p
3
p4
p2
p6
p
1
(a)
(b)
(c)
Figure 2.3: a) Initialization points, b) Image un-warping, IMMPDA edge detection and tapered ellipse fitting (dashed), c) Image un-tapering, IMMPDA
edge detection and ellipse fitting (dashed).
φ = 90◦ and reduces as r increases.
rnew = r − r sin(φ) exp(−r 2 /2σ 2 )
(2.1)
28
Chapter 2. Semi-automatic Segmentation of the Prostate
σ is a variable which represents the amount of radial stretch and is calculated
by solving Eq. 2.1 for φ = 90o to obtain:
σ=
s
−r 2
,
2 ln(1 − rnew /r)
rnew < r
(2.2)
For this angle, r is set as the distance between p1 and p4 and rnew as the
distance between p1 and the reflection of p4 about a horizontal line passing
through p2 . If rnew > r no un-warping is required.
Eq.2.1 is also used to un-warp the initial points. Assuming that the
presence of the TRUS probe causes uniform deformation along the prostate
images, all the other TRUS images are also un-warped. Therefore, the effect
of the TRUS probe on the gland is largely removed.
Image un-tapering and mid-gland ellipse fitting
The un-warped initial points are used to fit an initial tapered ellipse on the
now un-warped mid-gland slice (Fig. 2.3b). The tapered ellipse parameters
P = (x0 , y0 , ax , ay , t1 ) are found by solving the following problem using the
recursive Levenberg-Marquardt algorithm [7]:
minP {eT e | ei =
√
x
′
y
′
ax ay [( axi )2 + ( ayi )2 − 1]}
(2.3)
xi = (xi − x0 )/( at1y (yi − y0 ) + 1),
′
′
yi = (yi − y0 )
with ax and ay being the radii along the axes, −1 ≤ t1 ≤ 1 the tapering
parameter and [x0 , y0 ] the center of the shape. [xi , yi ], i = 1 . . . N , are the
coordinates of the boundary points with N being the number of initial points
as previously defined. Figure 2.4 shows the effect of changing the tapering
value, t, in a tapered ellipse from -1 (thick dashed line) to +1 (thick line).
Throughout this algorithm, 2D segmentation is carried out with the aid
of the IMMPDA edge detector [2]. In this edge detector, the evolution of
the radius from an arbitrary seed point inside the contour to the contour
edge is modeled as a dynamic system in the radius angle. Multiple models
can be incorporated in the dynamic system to accommodate abrupt changes
in the edge. Two constant velocity models (Eq. 2.4 ) have been used in our
approach.
29
Chapter 2. Semi-automatic Segmentation of the Prostate
Figure 2.4: The effect of changing the tapering parameter from -1 (thick
dashed line) to 1 (solid line). A tapering value of zero corresponds to an
ellipse (dotted line).
Xj (k + 1) =
Zj (k + 1) =
"
1 δθ
0 1
h
1 0
i
#
Xj (k) + Vj (k)
(2.4)
Xj (k) + ωj (k)
where j = 1, 2 is the
h number of itrajectory models used to describe the
boundary, Xj (k) = dj (k) dθj (k)T are the system states at the kth radius
in which dj (k) is the distance of the prostate boundary from the arbitrary
seed point inside the prostate and dθj (k) is its derivative with respect to
angle θ. Vj (k) is the process noise vector with covariance Qj (k). Zj (k) is the
process output (measured boundary location) and ωj (k) is the measurement
noise with covariance Rj (k). The two models used in our implementation
have the following process noise vector and measurement noise covariances:
Q1 (k) =
"
δθ 4 /4 δθ 3 /2
δθ 3 /2 δθ 2
#
101
Q2 (k) =
"
δθ 4 /4 δθ 3 /2
δθ 3 /2 δθ 2
#
105
R1,2 = 20
where the numerical values have been obtained by trial and error in prior
work [5, 6] and have not been adjusted for any of the patient data they
30
Chapter 2. Semi-automatic Segmentation of the Prostate
were used in. Each trajectory model is associated with a Kalman filter.
The output of these filters is combined with a probabilistic data association
filter for more accurate contour extraction. Since no numerical optimization
technique is used in the IMMPDA method, it is fast enough to be used
repetitively within a segmentation algorithm. Meanwhile, the interacting
multiple model (IMM) estimator increases its accuracy and robustness when
noise is present in the images. For further details please refer to [2]
The initial tapered ellipse fitted to the initial points is used to guide
the IMMPDA edge detector by setting limits on how far from this contour
the edge detector can search. These limits prevent the edge detection from
drifting away from the prostate boundary in regions with low image contrast.
The resulting edge points are once again fed to the Levenberg-Marquardt
algorithm to obtain an improved fit of the tapered ellipse. The tapering
value of this contour, t1 , is used to un-taper the ultrasound images. We
assume that the prostate is most tapered at the mid-gland and the tapering
linearly reduces to zero superiorly and inferiorly. Using the negative of the
tapering value for each slice, all images along with the initial points on the
mid-gland slice are un-tapered. The result of this step is transverse images
in which the prostate has an elliptical shape.
An ellipse is fitted to the un-warped and un-tapered initial points. This
can be done in many ways; we use the simple and fast generalized eigenvector
solution of [114]. The ellipse guides the IMMPDA edge detector and to the
extracted edges, a second ellipse is fitted on the mid-gland slice (Fig. 2.3c).
Semi-ellipsoid fitting for slice contour propagation
An IMMPDA edge detection and second ellipse fitting similar to that applied
to the mid-gland slice is carried out on the rest of the slices. However, since
there are no user defined points on these images, the initial ellipses are
created by fitting two semi-ellipsoids; one that extends superiorly toward
the base, and another inferiorly toward the apex (see Fig. 2.2). Each is
fitted to the union of the data points on the mid-gland contour and the
point on the intersection of the axial line passing through the center of the
mid-gland contour with the respective base − 1, or apex + 1 slice. This is
again solved as a generalized eigenvector problem. Two semi-ellipsoids were
found to give a better initial approximation than a single ellipsoid because
the prostate tapers more rapidly toward the apex.
The intersections of these two semi-ellipsoids with each of the slices are
used to guide the subsequent IMMPDA edge detection on each slice. A line
31
Chapter 2. Semi-automatic Segmentation of the Prostate
is fitted to the centers of the resulting 2D contours. This line represents the
main axis of the prostate and will be later used for re-aligning the final 2D
contours.
Tapered ellipsoid fitting
At this stage, 2D contours have been generated from all un-warped and untapered images. Yet there is no guarantee that smoothness and continuity is
maintained from one contour to the next, in the direction of the TRUS probe.
Hence, a tapered ellipsoid with an elliptical cross-section and tapering along
its main axis is fitted to these contours. Similarly to the 2D tapered ellipse
fitting, P = (x0 , y0 , z0 , ax , ay , az , t2 , t3 ) is obtained by solving the following
problem using the Levenberg-Marquardt algorithm:
minP {eT e | ei =
√
ax ay az [f (xi , yi , zi , P ) − 1]}
f (xi , yi , zi , P ) = (
′
xi
ax
)2 + (
′
yi
ay
)2 + (
′
zi
az
(2.5)
)2
xi = (xi − x0 )/( at2z (zi − z0 ) + 1),
′
yi = (yi − y0 )/( at3z (zi − z0 ) + 1),
′
′
zi = (zi − z0 ),
where ax , ay , az are the radii along the axes, [x0 , y0 , z0 ] is the position of
the center of the volume and t2 and t3 are the respective x and y tapering
values in the direction of the TRUS probe. [xi , yi , zi ] , i = 1 . . . M are the
coordinates of the M boundary points generated by segmenting all image
slices.
The fitting of this 3D shape to the boundary points is no longer a convex problem and is the most time consuming part of the algorithm. Suitable
selection of the initial parameters can greatly reduce the search period and
prevent the optimization algorithm from converging to local minima. To aid
the optimization algorithm, the center of the 3D shape, [x0 , y0 , z0 ], and the
axes, ax , ay , az , are determined by first fitting an ellipsoid to the data cloud
consisting of the ellipse contours of all slices. Since this is a generalized
eigenvalue problem, the one and only minimum is found almost instantly.
The six derived parameters are used along with two tapering parameters, t2
and t3 , initially set to zero to define an ellipsoid, as starting values for the
optimization algorithm (Eq.2.5).
32
Chapter 2. Semi-automatic Segmentation of the Prostate
Contour tapering and warping
We have assumed that the prostate is parallel to the TRUS probe and not
rotated about its main axis. This reasonable assumption simplifies the fitting problem since introducing three additional rotation parameters into the
optimization algorithm will increase both the running time and the possibility of the algorithm converging to local minima. To take account of the
possibility that the prostate may be slightly pitched or yawed, the best line
fitted to the centers of the extracted 2D contours of all images before 3D
shape fitting is used as the actual axis of the final 3D shape. After slicing
the 3D tapered ellipsoid at the corresponding image depths, the centers of
the final 2D contours are repositioned to align with this axis.
Finally, the repositioned contours are tapered and warped to match the
original images (Fig. 2.5a, 2.5b). The negative of the same tapering values
initially used to un-taper the images is used to taper the contours. However,
to ensure that the posterior aspects of the contours do not overlap with the
rectum, the warping parameter may need to be modified. When the midposterior point on the final contour is lower than the initial posterior point
selected by the user, the σ is re-calculated using (2.2) but with r and rnew
respectively set to the distance between p1 and p4 and the distance between
p1 and the mid-posterior point on the final mid-gland contour. This change
in σ is only made if r > rnew , otherwise, the previous value of σ is used.
2.2.2
Evaluation
The presented semi-automatic prostate segmentation algorithm is currently
being used by all the radiation oncologists who practice prostate brachytherapy at the Vancouver Cancer Centre, BCCA. Prior to treatment planning,
initial delineations of the prostate, hereinafter called ‘Raw TES contours’,
are approved and modified, if deemed necessary, by the radiation oncologist (RO) in charge of the patient treatment. We will call these approved
contours, whether modified or not, ‘RO reviewed TES contours’. The RO
reviewed TES contours are used by medical physicists to generate a treatment plan, which in turn is again checked by the radiation oncologist before
the actual implant takes place. Thus, while our semi-automated prostate
segmentation is used to provide initial contours, these are not used for treatment without the treating radiation oncologists having the final say.
Modifications applied to the contours are not always due to segmentation errors, but are mainly required for brachytherapy treatment planning
purposes, to ensure that the treatment plans created on these contours are
33
Chapter 2. Semi-automatic Segmentation of the Prostate
(a)
(b)
Figure 2.5: a) Final contours on the TRUS B-mode images, b) Final semiautomatic 3D volume (thick lines) compared to manual segmentation (thin
lines).
robust and implantable. This is also observed in manual segmentation where
the physician may not necessarily follow the prostate boundary, e.g. to avoid
needle interference with the pubic arch or manage dose distribution in a certain region.
To evaluate this algorithm we have carried out two sets of studies: evalu34
Chapter 2. Semi-automatic Segmentation of the Prostate
ation of the accuracy of the algorithm, and evaluation of the repeatability of
the algorithm. We start by comparing the semi-automatic contours before
and after modification to indicate to what extent and in which regions the
performance of the algorithm was not satisfactory for brachytherapy treatment planning. In order to understand how biased the physicians’ modifications are by the initial semi-automatic contours, we compare the RO reviewed TES results with those of manual segmentation on the raw images.
In order to evaluate our method based on a commonly used approach in the
literature, we compare the Raw TES results with that of manual, generally
used as the reference in the literature. In the next step we measure interand intra-observer variability of both manual and semi-automatic segmentation. A comparison of these two provides a judgment on the repeatability of
the algorithm. The acceptable range for the segmentation error is provided
by the intra- and inter-observer variability of manual contouring.
The term ‘reference’, used throughout this work refers to the prostate
geometric shape relative to which the comparison is carried out. Depending on the study, this can be manual segmentation results, RO reviewed
TES contours, etc. The term ‘case’ refers to a set of patient images. An
‘observer’ is an individual carrying out segmentation. All observers who perform contouring in this chapter have adequate knowledge in manual and/or
semi-automatic segmentation and include an expert radiation oncologist, a
radiation therapist and a graduate student with significant training in ultrasound prostate segmentation.
The comparison measures used are:
- Mean Absolute Distance, M AD : the average absolute radial distance
between contours CA and CB , in a slice.
- Maximum Distance, M axD : the maximum absolute radial distance
between contours CA and CB in a slice.
Since the boundary of the prostate in the mid-gland slice is visible enough
for accurate manual segmentation, these two measures are calculated for this
slice only.
- Percent volume difference, Vdif f : the difference between the volumes
of two delineated prostates defined as:
Vdif f = (VA − Vref )/Vref × 100
(2.6)
in which ref in this equation is the reference and V denotes the volume.
- Percent volume error, Verr : the volume of the non overlapping region
between two delineated prostates defined as:
Verr = |(VA + VB − 2(VA∩B ))| /(VA + VB ) × 100
(2.7)
35
Chapter 2. Semi-automatic Segmentation of the Prostate
where Verr is in fact 1 − Dice similarity coef f icient [21]. Vdif f and Verr
provide measures of size similarity and shape similarity, respectively. They
are either calculated for the entire prostate, or for each of the nine sectors
of the gland (Fig. 2.6). The nine sectors are created by first subdividing
the prostate into posterior, anterior, and two lateral sectors, the latter of
which are considered as a single region. The axis of division is the axis
of the reference shape. These regions are then partitioned according to
whether they are in the base, mid-gland or apex (respectively 0.3, 0.4, 0.3
of the length of the base-apex axis), forming a total of nine sectors for the
purposes of analysis.
Figure 2.6: Division of the gland into nine sectors.
Quadrant-based division of the prostate is a common approach in the
literature [96, 106]. However, for our application, this subdivision scheme
provides more detailed reporting of results for different regions of the gland
and is clinically motivated by the different consequences of segmentation
errors with respect to treatment planning in these regions. For example,
because the posterior aspect of the prostate is adjacent to the rectum, overestimating the boundary there can result in high doses to the radiosensitive
rectal wall and subsequently higher rates of rectal morbidity. This analysis
aims to consider contouring performance in the context of treatment.
Accuracy
This consists of a comparison between 3D shapes generated from: (i) Raw
TES contours and RO reviewed TES contours (the reference), (ii) manually
segmented contours (the reference) and RO reviewed TES contours, (iii)
manually segmented contours (the reference) and Raw TES contours.
36
Chapter 2. Semi-automatic Segmentation of the Prostate
Comparison between Raw TES and RO reviewed TES contours can give
a measure of how satisfied the physicians are with the results of the algorithm and which regions of the prostate need the most modifications. It
can indicate the degree to which semi-automatic contouring alone is successful in providing reasonable prostate contours for treatment planning. In
this analysis, a total of 40 cases (randomly selected) was semi-automatically
segmented by various radiation therapists and then modified by radiation
oncologists (a similar study on a larger population of 140 cases is presented
in [63]).
Modifications by the physicians to the TES contours may be biased by
the initially given contours. To measure this bias, we compare manually
created contours for 10 prostate image sets with RO reviewed TES contours,
generated by three observers (one expert and two trained by experts) for the
same image sets.
Finally, to include a commonly used method of evaluation reported in
the literature, a comparison between prostate shapes generated by the TES
algorithm and by manual segmentation, is also carried out. For this purpose, 21 cases were manually segmented by an expert and two individuals
trained by an expert and the average of their contours was compared to TES
segmentations of an observer experienced with the algorithm and blind to
the manual contours.
Repeatability
This analysis is performed to evaluate the consistency of the contours across
different observers and at different times. It consists of a comparison between (i) manual contours generated by different observers vs. the ‘average’
manual contour (the reference) and Raw TES contours generated by different observers vs. the ‘average’ TES contour (the reference) and (ii) initial
Raw TES (the reference) vs. repeated Raw TES and initial manual (the reference) vs. repeated manual contours. For this aim, 10 cases were segmented
by different observers once manually and once using the algorithm. Five of
these cases were segmented again, both manually and using the algorithm,
after approximately 2 weeks. All observers were blind to each others contours, their previous segmentations and patient data. (i) gives a measure of
inter-observer variability and in (ii) intra-observer variability is quantified.
Manual intra-observer and inter-observer variability also provide means of
evaluating the accuracy obtained in the previous analysis.
The ‘average’ manual/TES contours are the average of the manually/TES
delineated gland, by each observer, on each slice.
37
Chapter 2. Semi-automatic Segmentation of the Prostate
Figure 2.7: Schematic of changes made to each sector of the TES contours
after modification as computed in Tables 2.1a and 2.1b (B: base, M: midgland, A: apex).
2.3
Results
The results of each evaluation study are as follows:
2.3.1
Accuracy
The mean and standard deviation of the percent volume error and volume
difference between Raw TES (segmentation algorithm’s results) and RO
reviewed TES contours (expert approved contours used for treatment planning) are calculated for the nine sectors and the total gland and presented
in Tables 2.1a and 2.1b. The average prostate volume created from Raw
TES contours and RO reviewed TES contours is 42.2 ml and 41.5 ml, respectively. The average whole gland Verr of 5.82 ± 4.15% corresponds to a
non-overlapping volume of 4.7 ml.
Whether these errors are clinically considered large or not will be determined in the following section. Based on these two tables, a schematic
of the modifications made in each sector is drawn in Fig. 2.7. From the
sagittal view it appears that on average, the algorithm over-estimates the
mid-anterior, anterior-apex and posterior-base and under-estimates anteriorbase, mid-posterior and posterior-apex. This may be a result of the algorithm not entirely capturing the tilt of the prostate. The coronal view shows
that laterally, the size of the base is increased and the mid and apex are reduced in size after modifications.
The average MAD and MaxD between Raw TES and RO reviewed
TES contours for the 40 cases on the mid-gland slice is 0.71±0.75 mm and
38
Chapter 2. Semi-automatic Segmentation of the Prostate
%
Ant.
Lat.
Post
Total
Base
6.22 ± 8.18
5.39 ± 5.10
7.97 ± 8.07
Mid
Apex
5.14 ± 4.78
9.38 ± 9.40
4.65 ± 4.86 10.18 ± 8.99
4.52 ± 4.44 11.31 ± 14.42
5.82 ± 4.15
(a) Verr
%
Ant.
Lat.
Post
Total
Base
-1.09 ± 19.67
-2.22 ± 11.94
11.13 ± 14.80
Mid
6.75 ± 13.31
4.66 ± 13.74
-4.33 ± 9.35
2.11 ± 9.93
Apex
16.75 ± 30.02
12.91 ± 28.67
-13.05 ± 24.35
(b) Vdif f
Table 2.1: Comparison of prostate sectors created from Raw TES and RO
reviewed TES contours showing the volume error (a) and volume difference
(b).
2.00±1.87 mm with 11 out of the 40 mid-gland contours requiring absolutely
no modifications.
Tables 2.2a and 2.2b show the average and inter-observer standard deviation of the absolute percent volume error and volume difference between
manual and RO reviewed TES contours of 10 cases created by three observers. This comparison measures the amount of bias in the physicians’
contouring when segmentation is done entirely manually as opposed to when
initial TES contours are given to be approved. Most of the difference between the manual and RO reviewed TES contours is seen in the base and
apex where visibility of the gland is low or absent. In these regions, the
observer tends to rely more on the given TES contours. In the mid-gland
region (which consists of over 50% of the prostate volume) the bias is lower
and similar to the other computed errors for this region (Raw TES vs. RO
reviewed TES and Raw TES vs. manual contours). A negative Vdif f in all
regions indicates that the prostate volume created from manual contours
tends to be larger than that of RO reviewed TES. The average manually
created volume is 55.6 ml and the Verr of 7.25 ± 0.39% corresponds to 8.2 ml.
Tables 2.3a and 2.3b show the absolute percent volume error and volume
difference between manual and Raw TES contours created on 21 cases. The
average MAD and MaxD between the manual and TES contours on the
39
Chapter 2. Semi-automatic Segmentation of the Prostate
%
Ant.
Lat.
Post
Total
Base
7.65 ± 0.16
10.52 ± 3.71
13.97 ± 1.28
Mid
5.68 ± 1.50
5.33 ± 0.52
5.27 ± 0.54
7.25 ± 0.39
Apex
8.66 ± 1.46
10.17 ± 2.35
14.05 ± 2.64
Mid
-0.98 ± 7.53
-4.23 ± 3.59
-2.44 ± 4.93
-6.64 ± 2.36
Apex
-6.75 ± 3.42
-9.92 ± 3.99
-13.71 ± 6.18
(a) Verr
%
Ant.
Lat.
Post
Total
Base
-1.15 ± 7.33
-12.15 ± 4.14
-11.56 ± 8.24
(b) Vdif f
Table 2.2: Comparison of RO reviewed TES and manual prostate sectors
showing the volume error (a) and volume difference (b).
mid-gland slice are respectively 1.38±0.61 mm and 3.49±1.10 mm.
%
Ant.
Lat.
Post
Total
Base
7.77 ± 4.96
6.45 ± 2.30
7.71 ± 3.59
Mid
Apex
6.16 ± 4.08 11.27 ± 5.95
5.18 ± 1.40
9.59 ± 2.90
5.49 ± 2.66 11.83 ± 7.87
6.63 ± 0.90
(a) Verr
%
Ant.
Lat.
Post
Total
Base
-0.92 ± 16.84
0.23 ± 10.54
3.63 ± 15.65
Mid
8.31 ± 13.57
2.24 ± 6.98
6.77 ± 9.92
2.43 ± 6.08
Apex
2.79 ± 25.80
1.08 ± 16.85
0.19 ± 29.86
(b) Vdif f
Table 2.3: Comparison of Raw TES and manual prostate sectors showing
the volume error (a) and volume difference (b).
Since the initial boundary points are selected on the mid-gland slice, and
the mid-gland choice is not unique, we measured the sensitivity of the algorithm to the mid-gland slice selection. For this aim, 11 randomly selected
40
Chapter 2. Semi-automatic Segmentation of the Prostate
cases were semi-automatically segmented twice by an experienced individual. In the first round, the mid-gland slice was chosen as the first candidate,
i.e. the largest and most visible image of the gland. In the second round
it was selected as the next best candidate, one slice above or below the selection in the first round. The whole gland volume error between manually
created surfaces and surfaces created using the best mid-gland candidate was
5.56±1.21% and between manually created surfaces and surfaces created using the next best mid-gland candidate was 6.96±1.81%. Finally, the whole
gland volume error between the two TES surfaces with different mid-gland
slices was 6.16± 2.16%.
2.3.2
Repeatability
Tables 2.4-2.7 show the inter- and intra-observer variability for manual and
TES segmentation. In the presented tables, the mean and standard deviation (shown by the bars and error bars respectively) of Verr and Vdif f are
derived from the average performance of each observer over all cases (i.e.
from the mean value of Verr and Vdif f for each observer over all cases). It is
the standard deviation that determines the observer variability in segmenting the prostate.
As shown in Table 2.4c and Table 2.5c, the inter-observer variability of
TES contouring is less than that of manual segmentation in most of the
sectors. The only sector in which manual inter-observer variability is less
(as seen in both the Verr and the Vdif f bar graphs) is the posterior-apex
sector. Additionally, the relatively small manual Vdif f mean values along
with the large standard deviation values (Table 2.5) compared to that of
TES, indicates that most regions of the prostate can be simultaneously overestimated by some observers while under-estimated by the others in manual
segmentation, whereas in TES segmentation, there is more agreement between observers in under-estimating or over-estimating different regions of
the gland.
Figures 2.6c and 2.7c show the intra-observer variability in manual and
TES contouring. As in inter-observer variability, in most sectors, intraobserver variability is less in TES contouring compared to that of manual.
The sectors in which manual intra-observer variability is noticeably better
are the anterior-apex and lateral-apex. However, similar to inter-observer
variability, intra-observer variability is less on the entire gland in TES segmentation.
A comparison of manual Verr values in Table 2.4 and Table 2.6 with Table 2.1a shows that the error between the TES contours and those modified
41
Chapter 2. Semi-automatic Segmentation of the Prostate
%
Ant.
Lat.
Post
Total
Base
9.62 ± 2.98
5.77 ± 1.00
6.73 ± 1.15
Mid
4.40 ± 0.83
3.60 ± 0.69
3.90 ± 1.40
4.65 ± 0.77
Apex
5.47 ± 1.38
5.36 ± 0.97
4.93 ± 0.89
Mid
1.86 ± 0.30
2.45 ± 0.48
3.96 ± 1.02
3.04 ± 0.30
Apex
3.59 ± 0.64
4.64 ± 0.48
6.21 ± 1.38
(a) manual
%
Ant.
Lat.
Post
Total
Base
3.32 ± 0.62
3.26 ± 0.30
4.95 ± 0.47
(b) TES
14
Manual
Semi−auto
12
10
Verr
8
6
4
2
0
B
M
Ant.
A
B
M
A
B
Lat.
M
Post.
A
Tot.
(c) Verr (%), manual vs. TES
Table 2.4: Inter-observer variability in manual and TES contouring - as
characterized by the volume error Verr (Ant.: anterior, Lat.: lateral, Post.:
posterior, B: base, M: mid-gland, A: apex, Tot.: total gland).
for treatment planning is in the order of the difference in the segmentations
done by different observers or by an individual at different times. Therefore, based on the above results, it is reasonable to claim that if the TES
contours were not to be modified by physicians and were to be used directly
for treatment planning, they are just as likely to provide adequate results
as manual or modified contours.
42
Chapter 2. Semi-automatic Segmentation of the Prostate
%
Ant.
Lat.
Post
Total
Base
0.59 ± 16.91
-0.01 ± 6.02
-0.19 ± 6.13
Mid
-0.42 ± 3.89
-0.41 ± 4.32
-0.69 ± 7.80
-0.31 ± 4.10
Apex
-0.67 ± 6.75
-0.24 ± 7.84
-0.39 ± 7.43
Mid
-1.80 ± 1.11
-1.10 ± 2.14
-1.36 ± 4.93
-1.23 ± 2.02
Apex
-2.08 ± 2.87
-0.65 ± 4.75
-1.29 ± 7.86
(a) manual
%
Ant.
Lat.
Post
Total
Base
-1.97 ± 3.12
-1.07 ± 1.14
-1.24 ± 4.59
(b) TES
10
Manual
Semi−auto
8
6
4
Vdiff
2
0
−2
−4
−6
−8
−10
B
M
Ant.
A
B
M
Lat.
A
M
B
Post.
A
Tot.
(c) Vdif f (%), manual vs. TES
Table 2.5: Inter-observer variability in manual and TES contouring - as
characterized by the volume difference Vdif f (Ant.: anterior, Lat.: lateral,
Post.: posterior, B: base, M: mid-gland, A: apex, Tot.: total gland).
2.3.3
Performance
The graphical interface used at the Vancouver Cancer Centre for semiautomatic segmentation gives the user the opportunity to modify the initial
points to best fit a 2D contour to the mid-gland image. After approval
of this contour, the algorithm proceeds to the 3D shape fitting. This can
be repeated until a satisfactory delineation is achieved. The contours are
then exported for further brachytherapy planning to the VariSeed software
43
Chapter 2. Semi-automatic Segmentation of the Prostate
%
Ant.
Lat.
Post
Total
Base
9.77 ± 7.04
6.65 ± 2.45
7.20 ± 0.91
Mid
6.01 ± 3.19
4.55 ± 1.36
3.90 ± 1.43
5.95 ± 1.59
Apex
7.99 ± 1.66
7.73 ± 0.69
9.75 ± 1.91
Mid
2.16 ± 0.40
3.21 ± 0.86
2.46 ± 0.90
3.48 ± 0.95
Apex
5.20 ± 2.41
6.18 ± 2.57
6.04 ± 2.23
(a) manual
%
Ant.
Lat.
Post
Total
Base
4.30 ± 2.01
4.07 ± 1.06
4.23 ± 1.05
(b) TES
18
Manual
Semi−auto
16
14
Verr
12
10
8
6
4
2
0
B
M
Ant.
A
B
M
A
B
Lat.
M
Post.
A
Tot.
(c) Verr (%), manual vs. TES
Table 2.6: Intra-observer variability in manual and TES contouring - as
characterized by the volume error Verr (Ant.: anterior, Lat.: lateral, Post.:
posterior, B: base, M: mid-gland, A: apex, Tot.: total gland).
(Varian Medical Systems, Palo Alto, CA).
The average duration of the TES segmentation per case (calculated for
the 40 cases used in the accuracy study), from after initialization until the
final contours are created is 14.36 ± 1.39 s on a standard PC (Intel Xeon
CPU, Intel; Santa Clara, CA, 2.27 GHz, 3.23 GB RAM). Of this duration,
the most time-consuming sections are the 3D tapered ellipsoid fitting, being an iterative process (2.72 ± 0.27 s), and the image un-warping of all
44
Chapter 2. Semi-automatic Segmentation of the Prostate
%
Ant.
Lat.
Post
Total
Base
11.59 ± 7.51
-4.22 ± 5.98
-4.60 ± 5.15
Mid
8.09 ± 5.17
3.88 ± 4.03
2.57 ± 5.45
2.29 ± 3.74
Apex
8.84 ± 2.97
2.59 ± 8.02
-5.31 ± 6.38
Mid
1.60 ± 1.13
1.42 ± 3.98
-1.11 ± 2.08
0.90 ± 2.83
Apex
0.98 ± 8.07
3.11 ± 9.71
-2.76 ± 2.94
(a) manual
%
Ant.
Lat.
Post
Total
Base
3.55 ± 7.20
1.20 ± 5.48
-0.74 ± 6.85
(b) TES
20
Manual
Semi−auto
15
10
Vdiff
5
0
−5
−10
−15
B
M
Ant.
A
B
M
Lat.
A
M
B
Post.
A
Tot.
(c) Vdif f (%), manual vs. TES
Table 2.7: Intra-observer variability in manual and TES contouring - as
characterized by the volume difference Vdif f (Ant.: anterior, Lat.: lateral,
Post.: posterior, B: base, M: mid-gland, A: apex, Tot.: total gland).
TRUS images (4.23 ± 0.38 s). However, with further code optimization,
these durations can be reduced. The selection of the initial slices and the
7 initial points requires 32.3 ± 14.0 s for someone familiar with TRUS images of the prostate. The physicians using the results of this algorithm have
reported an average modification time of 1-3 minutes. Their modifications
are reported to mainly consist of shifting or changing the overall size of
the contour (especially the base and apex) which is done with ease in the
45
Chapter 2. Semi-automatic Segmentation of the Prostate
VariSeed software. Based on the above, the total segmentation duration, including initialization, is less than 1 minute. Including contour modifications
for prostate brachytherapy purposes, it is expected to fall within the range
of 2-4 minutes.
The time required for manual segmentation varies between users, and
depends on their experience. An experienced radiation oncologist requires
approximately 5-10 minutes per case to perform manual segmentation. Clinical fellows, during training, require up to 30 minutes for manual segmentation, but reach the 5 to 15 minute range after 1-3 months of brachytherapy
training. However, a study to measure manual prostate segmentation time
has not been performed
By visual observation of the pre- and post-modified contours, it was
seen that out of the 369 segmented images of 40 patients, 28% needed either
absolutely no modification (26%) or very little modification (2% - contour
displacement of 1 mm or less).
2.4
Discussion and Conclusions
In this chapter, we presented a semi-automatic prostate segmentation algorithm in TRUS images and performed various clinical studies to evaluate the algorithm. Clinical results show that the inter-observer and intraobserver variability of the semi-automatic contours are less than that of
manual contours in most sectors. The regions in which the variability is
higher in semi-automatic segmentation is mainly the apical region of the
gland. Comparison of Raw TES and RO reviewed TES contours shows a
percent volume error of 5.82±4.15% on the entire gland which is comparable
to both the manual inter-observer (4.65 ± 0.77%) and manual intra-observer
(5.95 ± 1.59%) variability. Comparison of RO reviewed TES contours and
manual contours, a measure of physician bias when modifying the contours,
shows a percent volume error of 7.25±0.39% for the whole gland. The duration of segmentation after initialization has been reduced from a few minutes
(as seen in the literature) to less than 15 seconds (14.36± 1.39 s). By including the initialization and possible modification time, the total segmentation
time is less than 4 minutes. With the above results, we conclude that the
proposed semi-automatic prostate segmentation method is accurate and consistent enough to replace manual segmentation of the gland. By applying
slight modifications, such as removing the need for manual initialization,
this method has the potential to be used as a real-time intra-operative segmentation method. This is the subject of Chapter 5.
46
Chapter 2. Semi-automatic Segmentation of the Prostate
We observed that poor image quality could in some cases lead to unsatisfactory results. However, the algorithm is guided by the manually selected
initial mid-gland boundary points and the positions of the base and apex,
from which initial contours and surfaces are produced. Since the edge detection is performed within a certain limit of these initial contours and surfaces,
artefacts inside the prostate, such as calcifications should not pose a problem
and, as long as the image quality at mid-gland is adequate for the observer
to perform initialization, our method should provide consistent results.
The region-based volume measure of physicians’ modifications applied
to our semi-automatic contours (Fig. 2.7) suggests that our tapered ellipsoid shape assumption is reasonable. Other models based on priors could
also be implemented, and may improve segmentation in terms of accuracy.
For example, a statistically obtained prostate model, which includes possible prostate shape abnormalities (e.g. due to tumors), would be a good
choice. However, in addition to possibly increasing the segmentation time,
such a model will complicate the process of treatment planning and treatment modification, as the plans may be assumed to change continuously
with the parameters describing the prostate shape. A tapered ellipsoid parameterizes the shape with only a few parameters that are intuitive and easy
to understand.
As per the BCCA protocol, the contouring algorithm assumes that a
smooth and symmetric clinical target volume is the aim of the oncologist, who consequently positions the prostate symmetrically across the midsagittal plane during TRUS image acquisition. The use of symmetric contours for treatment planning is widely practiced as part of the popular Seattle preplanning technique [102]. Symmetric contours lead to simple treatment plans that are also simple to modify to ensure adequate dose coverage
should the shape, size, or position of the prostate change significantly with
respect to the volume study. By maintaining symmetry during the preoperative volume study, reproducing the prostate image intra-operatively is
relatively simple since the body’s long axis can be identified easily in the
dorsolithotomy position and does not change over time or in response to
shifting leg positions and tissue relaxation. However, replicating a specific
arrangement of misalignment is not easily accomplished since there are numerous ways to misalign the axes of the TRUS probe and the prostate, each
of which creates a somewhat different visual pattern of asymmetry on the
TRUS images. We emphasize the need to maintain proper body alignment
throughout both the TRUS image acquisition and intra-operatively since, in
most cases, maintaining this is sufficient to achieve symmetry on all slices.
Effective implementation of a symmetric planning approach is demonstrated
47
Chapter 2. Semi-automatic Segmentation of the Prostate
by the results of a population-based analysis that showed only 35 recurrence
events among the first 1006 consecutive BCCA prostate brachytherapy patients who underwent implant between July 1998 to October 2003 [66].
Since the semi-automatic segmentation algorithm is currently being used
in practice (to this date, more than 1000 patients have been treated based
on our segmentation method), it is possible to analyze the growing collection
of segmented prostate images. The existance of any pattern or bias in the
modifications could also be studied. Finally, we suggest the formation of a
complete, randomly selected standard clinical data set. A comparison of the
available segmentation methods on such a data set may be a suitable topic
for further work.
48
Chapter 3
Semi-automatic
Segmentation of the
Prostate:
Dosimetric evaluation
3.1
Introduction
In this chapter we provide a clinical validation of the tapered ellipsoid segmentation (TES) method from chapter 2.Currently, the semi-automatic contour is first approved and modified, if required, prior to treatment planning.
The volumetric results in chapter 2 suggested that such modifications are
so minor that they may not be necessary in many cases. Indeed, a volumetric study showed that the semi-automatic segmentation error is within
the range of inter- and intra-observer variability of manual contours in most
regions of the prostate, which suggests that on average, no greater variation is introduced by utilizing the algorithm than would be expected if a
different oncologist performed the contour. Here, we aim to show that the
segmentation error leads to a dose error that is negligible. The dosimetric
comparisons were designed to investigate what the impact on coverage of the
RO reviewed PTV would have been if planning had been performed directly
on the Raw TES PTV. To do this, treatment plans were originally created
on Raw TES contours, while satisfying the BC Cancer Agency (BCCA)
planning goals, and subsequently superimposed on the corresponding RO
reviewed TES contours. Plans derived from Raw TES PTVs were also compared to the plans created on the manual contours of different observers
on the same image set. Details of each of the evaluation methods are described in the following sections. The terminolgy used below is defined in
Chapter 1.2.2.
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
3.2
Patient Characteristics and Treatment
Planning
LDR brachytherapy is indicated at the BCCA for low and intermediate-risk
prostate cancer (all of: pretreatment prostate-specific antigen [iPSA]≤20,
Gleason score [GS]≤7, clinical stage [CS]≤T2c [International Union Against
Cancer (UICC) 1997]). Three to four weeks before the implant, a radiation oncologist performs a volume study in which 2D ultrasound images
are obtained at 5 mm intervals with the use of a trans-rectal ultrasound
probe (B&K Pro-Focus System B-series ultrasound machine; BK Medical,
Peabody,MA, with the MFI Biplane Transducer, 640×480 pixels image size,
0.155 mm × 0.155 mm pixel size). The patient is in the dorsal lithotomy
position during imaging. For applying the TES algorithm on these clinical
images appropriate institutional and ethics committee approval have been
obtained. The TES algorithm is initiated by a radiation therapist to produce
a clinical target volume (CTV) called the ‘Raw TES CTV’. The contours
are then transferred to the treatment planning system (VariSeed, Varian
Medical Systems, Palo Alto, CA) where they are modified if necessary and
then approved for planning by a radiation oncologist (RO). This volume will
be referred to as the ‘RO reviewed TES CTV’, which is used to produce the
planning target volume (PTV). For the purposes of comparing the dosimetric effect of the RO modifications, a second PTV was also generated directly
from the Raw TES CTV, which will be referred to as the ‘Raw TES PTV’.
The guidelines for the creation of the PTV at this institution recommend
applying 0.3 to 0.5 cm lateral, 0 to 0.3 cm anterior, and 0.5 cm superior
margins to the CTV. No planning margins are added posteriorly or inferiorly to spare the rectum and penile bulb. Although small variations in the
size of the margins were present among clinically generated PTVs, the margins applied to generate the Raw TES PTVs for this study complied with
the guideline recommendations (0.3 cm lateral, 0.2 cm anterior, and 0.5 cm
superior). An additional component of this study involved the use of contours that were generated completely manually (i.e. without the presence of
any preliminary contours on the image sets) by multiple blinded observers
(radiation oncologists, radiation therapists and/or individuals trained by
experts). We will describe these contours and their derivative structures
as ‘manually’ generated to distinguish them from the ‘RO reviewed TES’
contours which are initially produced by the TES algorithm.
Brachytherapy treatment plans were developed for the PTVs by a single medical physicist. These plans adhered to the standard BCCA plan-
50
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
ning algorithm, which can be generally described as following a low-activity
(0.424 U) modified peripheral loading strategy using custom loaded, stranded
seeds (RAPIDStrand, Oncura, Arlington Heights, IL). The goal is to have
the 100% isodose (the surface at which the radiation dose measures 100%
of the prescribed dose) cover the PTV on all slices. Since this is not usually
achieveable, guidelines are provided for maximum and minimum margins as
well as allowable exceptions. These guidelines may vary from one center to
the other. At VCC, each plan is designed to provide ≥97% coverage of the
PTV and ≥99% coverage of the CTV by the 100% (144 Gy) isodose, with
a CTV V150 between 56 and 65% and PTV V150 between 50 and-60%. The
V150 is geometrically biased to the posterolateral aspects of the target. The
volume that does not reach prescription dose in planning is confined to a
small region of the anterior base of the PTV whenever possible.
3.3
Evaluation Measures
The following measures were used for dosimetric evaluation of the TES
method. These measured were calculated for nine sectors of the prostate
as described in chapter 2.
The standard dose parameter, V100 : the volume of the PTV receiving ≥ 100% of the prescribed dose, computed for the nine sectors of the PTV
and the whole PTV. These values were calculated by the VariSeed software.
The External Index 150, EI150 : To characterize extraprostatic dose,
the External Index (EI) [92], defined in Equation 3.1, measures the amount
of tissue external to the PTV that receives doses 150% of the prescribed
dose.
EI150 = (isoV150 − V150 )/V
(3.1)
isoV150 is the total volume of the 150% isodose surface, V150 is the volume of the PTV receiving ≥150% of the prescribed dose (the volume of the
intersection between the isoV150 and PTV surfaces) and V is the volume of
the PTV. Ideally, EI150 is zero.
The Conformity Index 100, CI100 : A 3D extension of the conformity
index (CI) defined by Otto and Clark [81] is used, which measures both
under-coverage of the target as well as over-treatment of the normal tissues.
51
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
CI100 % = 100 ×
Vr − (Vu + Vh )
Vr
(3.2)
where Vr is the volume of the PTV (or one of its nine sectors), Vu is the
volume of the PTV (or one of its nine sectors) that is receiving less than
100% of the prescribed dose and Vh is the volume of the region outside the
PTV (or one of its nine sectors) that is receiving ≥ 100% of the prescribed
dose. A maximum conformity index of 1 shows perfect conformity of the
100% isodose to the region being observed.
We would like to note that although the above mentioned dose parameters are computed to evaluate the TES method, the planning process places
quantitative constraints only on the whole prostate and whole PTV V100 ,
V150 , and V200 (i.e. a constraint of V150 between 50 and-60% is for the whole
PTV, not a region). On the other hand, regional PTV coverage by the 100%
isodose, location and shape of the 150% isodose, and the extent of extraprostatic dose are evaluated and constrained on the basis of visual (qualitative)
inspection of the isodose distributions, rather than by adherence to regional
constraints on V100 or CI100 . Hence, in this study, we refer to the values
of the dosimetric parameters in each region, after planning, as the observed
clinical baselines (which do not necessarily adhere with the whole PTV constraints) and compute the dosimetric paired differences from these observed
clinical baselines in each region rather than reporting the absolute values.
3.4
3.4.1
Dosimetric Evaluation of the Algorithm
Impact of Planning using TES Contours
The aim of the dosimetric evaluation is to examine the clinical impact of
planning using Raw TES contours. This helps to put differences in volumetric coincidence in perspective, because if such differences do not result
in a significant degradation in dosimetry when a TES-derived plan is used
to treat a reference contour, then it is reasonable to suppose that the TES
and reference contour are of equivalent utility for planning purposes.
To investigate this, 41 anonymized consecutive patients (seen between
January-April 2009) had treatment plans generated using their Raw TES
PTVs as described in the ”Patient characteristics and treatment planning”
section. The aforementioned dose parameters for these plans were calculated
for the PTV and the nine sectors and used as the observed clinical baselines.
52
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
These plans were then overlaid on the reference (RO reviewed TES) contours, and the resulting dose parameters calculated for the PTV and the
nine sectors. The distribution of paired differences in the dose parameters
was calculated (i.e. dose parameter of the plan generated using Raw TES
PTVs and overlaid on RO reviewed TES PTVs minus the observed clinical
baseline values).
3.4.2
The Effect of Contour Variability on Planning
Although the impact of TES-based planning is readily calculated, establishing a sensible threshold for the acceptable amount of dosimetric degradation
below which the adoption of TES-based planning is unacceptable is challenging. For example, a plan with a whole PTV V100 below 97% would not be
accepted for implant at our institution, so it may seem natural to set this
as a target for TES-based planning. However, the patient might have been
seen by any of a number of oncologists, none of whose plans are explicitly
required to meet the 97% criterion on the contours of their colleagues. To
avoid a double standard, the evaluation of any automatic contouring algorithm cannot ignore the implicitly accepted differences in dosimetry which
arise from the endemic variability in target definition between observers.
Consequently, TES-based planning in this study is evaluated in this context, by elucidating the range of variability in dosimetric quality that is
entirely attributable to manual inter-observer variability, and then examining how the impact of TES-based planning compares to this range. This
is performed in two tests. First, in a subset of five (of the 41) cases, the
treatment plan produced on the set of contours originally used in the patient’s treatment was overlaid on the contours of all 10 observers with the
exception of the implanting radiation oncologist. In the second test, in one
of the 41 cases, the set of plans produced on the 10 observers’ PTVs were
mapped back on to the original planning PTV. In all of these tests, the observers were radiation oncologists, blinded to their colleagues’ contours. In
this study we argue that if TES-based plans fall within the range of manual
variability, it is reasonable to conclude that planning on the Raw TES CTVs
is as reliable, in a statistical sense, as planning on the contours drawn by a
colleague.
3.5
Results
Figures 3.1 and 3.2 display the paired differences in the V100 and CI100 when
the plans created on the Raw TES PTVs are mapped to the RO reviewed
53
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
TES PTVs. In both figures, each data point represents the difference in the
respective dose parameters for each of the 41 cases, for each sector of the
gland as well as the whole PTV. Overlaid on the raw data are the means
and 95% confidence intervals. Where this interval does not include zero, the
impact is considered to be statistically significant, and the corresponding
p-value is displayed for each region. The greatest impact on V100 was seen
in the anterior base, anterior apex, posterior base and posterior apex. In all
regions except the anterior base and apex, a statistically significant decrease
in V100 was found (p<0.05). For the whole gland, the mean PTV V100 fell
from 98.62±0.12% (observed clinical baseline), to 96.45±0.70% when the
Raw TES derived plans were applied to the RO reviewed TES contours.
With respect to CI100 , variability in the CI100 was most pronounced in
the apex and lowest in the mid-gland. The greatest mean decrease was observed in the anterior apex, which is consistent with the volumetric analysis
establishing a tendency of TES to over-contour this region (see Table 2.1).
However, in neither this nor any other region was there a statistically significant impact on the CI100 (p>0.05). For the PTV as a whole, the mean
CI100 of 0.68±0.02 fell to only 0.66±0.3 when the Raw TES derived plans
were mapped to the RO reviewed contours.
The mean and 95% confidence interval of the PTV 150% isodose external index EI150 (data not shown), increased from 0.065±0.004 (range
0.037-0.109) to 0.072±0.010 (range 0.025-0.160), a statistically insignificant increase in extra-target dose (p=0.22). The most significant increases
(p<0.05) in the EI150 were in the mid-anterior (0.01±0.004 to 0.02±0.01,
p=0.03) and lateral apex (0.21±0.02 to 0.27±0.06, p=0.04). However, significant decrease (p<0.05) in the extra-target dose was observed in the lateral
base (0.18±0.02 to 0.15±0.02, p=0.00) and posterior base (0.10±0.01 to
0.07±0.01, p=0.000). No significant changes were observed in other regions.
The planning goals in our center require a CTV V100 of 99% or greater
and a CTV V150 between 56 and 65%. Out of the 41 cases, 11 (27%) had a
CTV V100 less than 99%, 3 of which were less than 98% (96.0%, 97.8% and
97.3%). In 6 of these 11 cases the CTV V150 was also below 56% (range
50.3-55.9%).
54
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
Substantial variability in dosimetric coverage and conformity arising
from manual variability in target delineation is evident in Figs. 3.3 and 3.4,
which look at the V100 and CI100 parameters respectively. The subfigures
in each case indicate whether the reference plan was (a) mapped to other
observers’ PTVs, or (b) other observers’ plans were mapped back to a reference PTV. The reference PTV and plan were those of the oncologist who
treated the patient. For the test in which there was a reference plan, the
figures show the mean and 95% confidence intervals of the dose parameter
resulting from the application of the reference plan to each of 10 alternate
contours produced by the other observers. This was repeated for five cases,
and a mean and interval is shown for each of these cases, for each region of
the gland. For the test in which there was a reference PTV, only one of these
five cases was analyzed (the second of the five cases shown in Figs. 3.3a and
3.4a), and so only a single interval is shown in each region. The raw data
points associated with the data derived from manual contours are hidden in
order to highlight the relationship between the special case (marked by the
triangle symbol) in which the test involved the Raw TES PTV or its derived
plan, with the distribution of manual variability (i.e., using the Raw TES
PTV or its derived plan instead of the manual PTV or its derived plan, in
each test).
Extensive inter-observer and inter-case variability of the V100 in the anterior base, anterior apex, posterior base, and posterior apex and of the CI100
in the apex is noticeable. It is clear from the figures that in most of the
examined situations, the impact on dosimetric quality resulting from utilizing the TES algorithm is indistinguishable from the mean impact expected
when using another observer’s contours. In many cases where the impact is
not within this range, the degradation is less pronounced when using TES
contours than its manual alternatives.
A one-way analysis of variance test confirmed that in most regions of
the prostate in the test of a reference plan (Figs. 3.3a and 3.4a), the dose
distribution accuracy of the plans created on Raw TES PTVs, in terms of
the V100 parameter, was better than, or indistinguishable from that of the
manual distribution. The exceptions are in the anterior base and anterior
mid regions in three of the five cases (p<0.05). In terms of the CI100 , the
TES results are superior in almost all regions of the prostate for all five
cases (p<0.05) with the exceptions of the anterior base in two cases and
the mid-posterior and posterior apex sectors in another. For the tests in
which there was a reference PTV (Figs. 3.3b and 3.4b), the majority of the
TES results are either superior to or fall within the manual variability of
the manual results. The exceptions are in the anterior base for the V100 and
55
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
in the anterior base and mid posterior sectors for the CI100 .
It is clear from the figures that the dose parameters computed from
overlaying the reference treatment plan on contours from different observers
greatly differs from overlaying their plans on the reference treatment contour
(compare Figs. 3.3b and 3.4b with the second of the five cases in Figs. 3.3a
and 3.4a). For this case, the V100 values in Fig. 3.3b are in general less than
those in Fig. 3.3a. However, the opposite is observed for the CI100 values
in Fig. 3.4. This was expected since the radiation oncologist who created
the reference treatment plan for this case tends to create larger PTV’s.
Thus, the plans created on the other radiation oncologists’ contours, can
not completely cover the reference treatment PTV, resulting in lower V100
coverage. However, the reference plan created on the relatively large PTV,
when overlaid on other radiation oncologists’ manual contours, will result in
overdose. This is confirmed by more regions having negative CI100 values
in Fig. 3.4a. Therefore, it is important to consider the dynamics of both
parameters in evaluating impact, especially if only one of the above two
tests are performed. Looking at V100 in isolation obscures the inherent bias
towards over-treatment, as a plan generated for a high volume target is more
likely encompass the volumes of other observers and result in good coverage.
56
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
Figure 3.1: Plot illustrates the means and 95% confidence intervals of the
paired difference in the V100 when the treatment plans generated based on
Raw TES PTVs are overlaid on RO reviewed TES PTVs (i.e. V100 of the
plan generated on Raw TES PTVs with respect to the Raw TES PTVs, the
observed clinical baselines, subtracted from V100 of the plan generated on
Raw TES PTVs with respect to the RO reviewed TES PTVs). A negative
result indicates a decrease in V100 . Each column of data represents a particular region, labeled at the bottom. Each data point represents the result
of each of the 41 cases analyzed. The p-values for a null-hypothesis of no
effect in each region are presented at the top of the figure. Two data points
fall outside the range of this figure, and are indicated by ’+1’. Their values are -46.5% (anterior base) and -50.4% (posterior apex). A statistically
significant reduction in V100 is evident in most regions, although the mean
magnitude of the reduction is less than 5% (Ant.: anterior, Lat.: lateral,
Post.: posterior, B: base, M: mid-gland, A: apex).
57
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
Figure 3.2: Plot illustrates the means and 95% confidence intervals of the
paired difference in the CI100 when the treatment plans generated for Raw
TES PTV’s are overlaid on RO reviewed TES PTVs (i.e. CI100 of the
plans created on Raw TES PTVs with respect to the Raw TES PTVs, the
observed clinical baselines, subtracted from the CI100 of plans created on
Raw TES PTVs and overlaid on RO reviewed TES PTVs). Each column
of data represents a particular region, labeled at the bottom. Each data
point represents the result of one of the 41 cases analyzed. A negative value
result indicates a CI100 less than that of the observed clinical baseline. The
number of data points falling outside the displayed range in each region,
are indicated. The range of values for these four sectors: Ant. B, Ant. A,
Lat. A, and Post. A, are: [-1.70 0.73], [-2.32 1.36], [-2.10 0.61], and [-1.24
1.49]. The lateral base and posterior base exhibit a statistically significant
increase in conformity (Ant.: anterior, Lat.: lateral, Post.: posterior, B:
base, M: mid-gland, A: apex).
58
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
(a)
Figure 3.3: The means and 95% confidence intervals of the V100 when (a)
a reference plan, created for one set of contours, is overlaid on the manual
contours of 10 other observers in 5 patients. The typical goal V100 goal
of 97% is shown as a dashed horizontal line. For comparison, the results
of performing the test using the TES contours or their derived plans are
displayed using the triangle symbols (Ant.: anterior, Lat.: lateral, Post.:
posterior, B: base, M: mid-gland, A: apex).
59
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
(b)
Figure 3.3: Continued. The means and 95% confidence intervals of the V100
when (b) a reference PTV is overlaid on the treatment plans created for
the manual contours of 10 radiation oncologists, for a single patient (second
of the five patients in (a). The typical goal V100 goal of 97% is shown as
a dashed horizontal line. For comparison, the results of performing the
test using the TES contours or their derived plans are displayed using the
triangle symbols (Ant.: anterior, Lat.: lateral, Post.: posterior, B: base, M:
mid-gland, A: apex).
60
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
(a)
Figure 3.4: The means and 95% confidence intervals of the CI100 when (a)
a reference plan, created for one set of contours, is overlaid on the manual
contours of 10 other observers in 5 patients. For comparison, the results
of performing the test using the TES contours or their derived plans are
displayed using the triangle symbols (Ant.: anterior, Lat.: lateral, Post.:
posterior, B: base, M: mid-gland, A: apex).
61
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
(b)
Figure 3.4: Continued. The means and 95% confidence intervals of the CI100
when (b) a reference PTV is overlaid on the treatment plans created for the
manual contours of 10 radiation oncologists, for a single patient (second
of the five patients in (a)). For comparison, the results of performing the
test using the TES contours or their derived plans are displayed using the
triangle symbols (Ant.: anterior, Lat.: lateral, Post.: posterior, B: base, M:
mid-gland, A: apex).
62
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
3.6
Discussion
In this chapter we presented a dosimetric evaluation of our semi-automatic
prostate segmentation algorithm (TES) for ultrasound images [60]. For 41
cases we measured the difference between the V100 and CI100 dose parameters of treatment plans created for the Raw TES PTV, used as the baseline,
and treatment plans created for the Raw TES PTV’s but overlaid on RO
reviewed TES PTV’s. The mean decrease in V100 and CI100 was less than
5% and 0.2 respectively, in all regions of the gland. The greatest degradation in quality occurred in the posterior base and apex, and anterior base
and apex for the V100 and in the apex for the CI100 . However, this study
has demonstrated, in a subset analysis of 5 cases with 10 blinded observers,
that any differences in the distribution of dose when planning using TES
contours are largely comparable with manual dosimetric variability between
observers. Moreover, this variability only considered a single institution,
and may be even greater between experts at different institutions due to
diversity in training backgrounds and treatment strategies.
Our program regards the reproducibility of the alignment between the
prostate, the probe and the patient’s craniocaudal axes to be important,
as the accurate registration of the pre-planned PTV with the prostate as
visualized on the day of the implant to be a vital component in streamlining
the procedure and reducing setup complications. This is facilitated by ensuring that the prostate is positioned so as to have mid-sagittal symmetry
in the planning images. In addition, a smooth CTV topography is preferred
to reduce the sensitivity of coverage to misalignment. The TES algorithm
achieves these two goals with a minimum of operator assistance.
In our experience, the algorithm greatly reduces the time necessary to
arrive at an acceptable CTV. The initialization of the algorithm and generation of a smooth and symmetric 3D surface, which is tedious to accomplish
by hand, requires less than a minute by a radiation therapist. Once this (the
Raw TES) CTV is complete, only 2-4 minutes of review and modification
are required by the radiation oncologist to arrive at what we have described
as the RO reviewed TES CTV, which is currently used for planning.
3.7
Conclusions
The results of this study suggest that many of the modifications to the Raw
TES PTVs prior to planning are superfluous, in the sense that the impact of
not performing the modifications will result in a planned dose distribution
63
Chapter 3. Segmentation of the Prostate: Dosimetric Evaluation
not dissimilar in quality to that which would have been delivered if the
patient had been treated by a colleague. On the basis of this finding, we
conclude that the proposed TES algorithm is a suitable replacement for
manual prostate segmentation in a pre-planned treatment methodology.
64
Chapter 4
Visualization of the Prostate
Gland in Vibro-elastography
Images
4.1
Introduction
In this chapter, vibro-elastography (VE), an ultrasound-based method that
creates images of tissue viscoelasticity contrast, is evaluated as an imaging
modality to visualize and segment the prostate. We report a clinical study
to characterize the visibility of the prostate in VE images and the ability
to detect the boundary of the gland. Measures for contrast, edge strength
characterized by gradient and statistical intensity change at the edge, and
the continuity of the edges are proposed and computed for VE and B-mode
ultrasound images. Furthermore, using MRI as the gold standard, we compare the error in the computation of the volume of the gland from VE and
B-mode images.
4.2
4.2.1
Materials and Methods
Data Acquisition
The B-mode ultrasound, VE, and MR images used in this study were acquired from patients going through the standard LDR prostate brachytherapy procedure at the Vancouver Cancer Center, BC Cancer Agency, or radical prostatectomy at the Vancouver General Hospital. Institutional ethics
approval and informed patient consent were obtained prior to data collection.
The MR images were collected between one and two weeks prior to treatment. T2-weighted transverse MR images (slice spacing 4 mm, pixel size
0.27 mm ×0.27 mm) were collected at the UBC MRI research centre, Vancouver, with a Philips Medical Systems Achieva 3.0 Tesla MRI scanner. To
Chapter 4. Visualization of the Prostate Gland in VE Images
minimize the deformation of the gland and for patient comfort, imaging was
carried out in a supine position and a cardiac coil was used.
The ultrasound data were collected intra-operatively, just prior to the
actual brachytherapy or prostatectomy intervention, using the hardware described in Section 1.2.3 and the signal processing described in Section 4.2.1.
The vibrating TRUS probe, with a vibration range of 2-10 Hz (for the
first 14 cases) and 2-5 Hz (for the last 6 cases) and vibration amplitude
of approximately 1 mm, was swept from −45◦ to 50◦ , covering a fan of
longitudinal planes passing through the transducer axis with a predefined
angular separation of 2◦ . At each angle, one B-mode image and frames of
RF data were collected continuously for 2 seconds at an approximate frame
rate of 40 fps. Imaging depth was set to 5 cm (or 6 cm for large prostates).
The entire procedure, including setting up the system, adjusting the
TRUS probe in order to see the entire prostate in the ultrasound image
(performed by a physician) and collecting the data requires approximately
20 minutes of the OR time.
VE system signal processing
The RF data were processed to obtain longitudinal VE transfer function
images. Single DOF axial motion estimation was used in processing the
RF data. Each RF data frame has 128 lines of RF data, each having 1296
samples (1424 for larger prostates). The collected RF data lines are split
into blocks of 26 samples, or equivalently 1 mm, with 50% block overlap.
The axial displacement xjk (t) of each block at axial location k and lateral
location (line) j, and subsequently the axial strain, were computed from one
RF data frame to the next by using a correlation-based method, as described
in [122], resulting in axial displacement images of 128 lines by 100 blocks.
An image-based reference was used to compute the transfer function
images. It was computed as the average strain, at half the tissue imaging
depth, of all the 128 lines. This specific depth was chosen for to two reasons. First, the ultrasound imaging focal point is usually set at this depth,
resulting in a more accurate motion estimation. Second, in our images, the
half depth line is more likely to be enclosed entirely within the prostate,
resulting in more uniform mechanical properties along the line.
jk
(j2πf ) from this reference to each of the
The transfer functions (Href
blocks in the strain images were computed using standard signal processing
methods described in [90]. Vibro-elastography images used in this chapter
66
Chapter 4. Visualization of the Prostate Gland in VE Images
were generated by computing:
V Ejk
1
=
f2 − f1
Z
f2
f1
jk
j−1,k
Href (j2πf ) − Href (j2πf ) df
(4.1)
where V Ejk is the pixel intensity value at axial location k and lateral location
(line) j and [f1 , f2 ] describes the frequency range of interest, which in this
work coincides with the range of the broad-band vibration applied to tissue.
If we assume that the Fourier transform of the reference described above
is unity, then the contrast in the transfer function image from one spatial
location (e.g. a reference) to another corresponds to the strain energy difference between these locations in the frequency range f1 to f2 . By computing
the ‘difference’ between transfer functions of consecutive blocks j and j − 1
on line k, with respect to reference ref (equation 4.1), changes in stiffness,
including the prostate boundary, are highlighted. Throughout this chapter,
and this chapter only, we will refer to these vibro-elastography ‘difference’
images as VE images.
From these sets of longitudinal images, 3D VE volumes were generated
by interpolation (slice spacing in transverse direction, 0.43 mm, pixel size
0.5 mm ×0.5 mm). Similarly, the 3D B-mode volume (slice spacing in transverse direction 0.43 mm, pixel size 0.37 mm ×0.37 mm) was created from
the collected longitudinal B-mode images. Various approaches to creating
3D volumes are described in detail in [27]. Our method of constructing a
3D volume by interpolating longitudinal images suffers from a decrease in
image accuracy as the depth increases. However, the accuracy can be increased by using finer angle increments at the expense of scan time. Another
alternative is to use the transverse (convex) array of the TRUS probe to generate a 3D sector view by changing the depth of the transducer with respect
to the prostate. This approach also has the issue of decreasing accuracy
at increased depth because the scan is acquired with a convex array. Furthermore, such a transverse plane sweep changes the TRUS position in the
rectum causing the prostate to move as a function of probe depth, which is
clearly undesirable. The analysis performed in this work uses the transverse
plane of the constructed 3D volumes for evaluation. This choice was made
because of the preference of the transverse view by the clinicians, and also
for agreement with the transverse MRI images.
Figure 4.1 shows transverse VE, B-mode and MRI images of the prostate
mid-gland of two patients. Half of the prostate boundary is delineated in
one of the image sets. In Fig. 4.2 mid-sagittal VE (top row) and B-mode
(bottom row) prostate images of three different patients can be seen. The
prostate is the stiffer region as outlined in the VE image.
67
Chapter 4. Visualization of the Prostate Gland in VE Images
Table 4.1 describes the image types and number of cases used in this
chapter. VE and B-mode images were obtained from 20 patients. Among
these 20 cases, which consist of 178 co-registered VE and B-mode image
pairs, 107 images are in the transverse plane and are used for edge evaluation. The remaining 71 images are in the longitudinal plane at multiples
of 5◦ intervals with the mid-sagittal plane and are used along with the 107
transverse images, for CNR analysis. MRI data were also available for nine
of the patients. These nine data sets are used for volume evaluation. The
initial study started by recruiting patients for VE imaging at the time of
the volume examination. However, we found that the patients were nervous
and very uncomfortable with the additional time required to take the images
with a different ultrasound machine. In some cases we were able to collect
US data, but the same patients were not available for MRI, due to illness, ineligibility for MRI, and difficulty traveling for the additional medical exam.
For the data acquired during the volume study, the acquisition time had to
be minimized. Therefore the number of MRI scans that we have used in
our study is limited. Nevertheless, this number is similar to that used in
[98] that compares B-mode based prostate segmentation to MRI. Furthermore, in [98] and most of the studies cited therein, post-operative MRI is
compared with TRUS. However, tissue will change as a result of treatment.
We believe that comparison based on pre-operative images is necessary for
general assessment of the gland’s visibility and treatment planning.
Table 4.1: Description of the data used in this chapter
Data set size
9 cases
Image type
MRI, VE,
B-mode
Image plane
transverse
CNR analysis
178 images
(20 cases)
VE, B-mode
107 transverse,
71 longitudinal
Edge analysis
107 images
(11 cases)
VE, B-mode
transverse
Volume analysis
4.2.2
Evaluation Methods
In order to assess the visibility of the prostate in VE images, we utilize
several evaluation measures. These include image-based measures, which
consist of contrast-to-noise ratio, edge continuity and edge strength, and
volume-based measures, which consist of volume error and volume difference.
68
Chapter 4. Visualization of the Prostate Gland in VE Images
Figure 4.1: Transverse B-mode (left), VE (middle) and MRI (right) prostate
images of two patients. The boundary of the prostate is partially segmented
in the second set of images.
Figure 4.2: VE (top) and B-mode (bottom) sagittal images of the prostate
of three different patients. The boundary of the prostate is outlined in one
of the patients.
For more accurate region-based evaluation, the results of these measures will
be reported for nine sectors of the gland, similar to those used in the previous
69
Chapter 4. Visualization of the Prostate Gland in VE Images
Figure 4.3: Division of the prostate into nine sectors
chapters and as shown in Fig. 4.3.
Image-based measures
Accurate delineation of the prostate is affected by how the prostate is visualized. A “good” image of the prostate is one in which the prostate has high
contrast relative to the background and has well-defined edges. A good edge
should be both “continuous” and “strong”. In this analysis, VE images of
the prostate are compared to B-mode ultrasound images using the measures
of contrast, edge strength, and edge continuity as described below.
The Contrast-to-noise ratio
To compare the contrast of VE and B-mode images, the contrast-to-noise
ratio (CNR) was calculated using the following equation [9]:
CN R =
2(mt − mb )2
σt2 + σb2
(4.2)
in which m and σ 2 are the mean and variance of pixel intensities of the
target, t, and background, b, in a region of interest (ROI). The ROI of
the target and background were manually selected in four regions of the
prostate (lateral left, lateral right, anterior and posterior). The target is an
area inside the prostate while the background is an area outside the prostate
close to the target ROI (respectively yellow and blue boxes in Fig. 4.4). The
physical size and position of the ROI in VE and B-mode images were similar.
70
Chapter 4. Visualization of the Prostate Gland in VE Images
(a) B-mode ultrasound
(b) VE
Figure 4.4: Selection of inside (yellow boxes) and outside (blue boxes) regions of the prostate for CNR computation in (a) B-mode and (b) VE images.
Additionally, to ensure similarity between the intensity range in both set of
images, histogram stretching [0-255] was initially applied to the images.
The contrast between the prostate and its background is typically low
in B-mode ultrasound images (see Fig. 4.1). However, due to the presence
of a visible edge, delineation of the prostate is still possible. Therefore,
CNR alone can not evaluate the visibility of an object in an image and edge
evaluation is also required.
Edge strength
An edge can be defined as the boundary between two regions in an image
that are different from each other with respect to some local property [50].
While conceptually clear, this definition is difficult to quantify. However,
for the purpose of this study we need quantitative measures. Two different measures of edge strength are presented here. One is a conventional
gradient-based edge filter and the second one is a test that identifies the
statistical changes in image intensity. The edge strength measures are computed for radial edge profiles formed as follows: On each transverse image,
we extended radii at θi = i × 30◦ , i = 1, ..., 12 in polar coordinates originating from a manually selected center point, C, inside the prostate. The
intersections rθi , i = 1, ..., 12, of these radii with the prostate boundary were
manually identified for each angle θi . For each edge point, a radial edge
intensity profile Iθi (r) was extracted. The measures of edge strength were
calculated for the window of r ∈ [rθi − ∆d, rθi + ∆d], where ∆d is half the
length of the edge intensity profile (see Fig. 4.5).
71
Chapter 4. Visualization of the Prostate Gland in VE Images
(a)
(b)
Figure 4.5: (a) The radii used for edge profile extraction in a VE image,
originating from C and with angles θi = i × 30◦ , i = 1, ..., 12. For one of
the radii, the two neighboring radii are also illustrated. (b) A close-up view
of one of the rays used for extracting the edge profile Iθi (r). The measures
of edge strength are calculated for the window of r ∈ [rθi − ∆d, rθi + ∆d],
where rθi is a manually selected edge point along the ray. The neighboring
edge profiles Iθi ±δθ , are extracted similarly and used along with Iθi in the
edge continuity measure.
Gradient-based measure of edge strength: The following gradient formulation, also used in [2], was used as the edge filter acting on a radial edge
intensity I(r), where we remove, for convenience, the θi index of I:
fedge (r) = 1/3×[I(r+2∆r)+I(r+∆r)+I(r)−I(r−∆r)−I(r−2∆r)−I(r−3∆r)]2
(4.3)
where ∆r is the physical size of the image pixel.
M , our measure of edge strength, is the sum of the filter outputs on windows
of size n pixels (n an odd integer) normalized to the sum of edge filter values
on the entire edge profile:
j=rθi +∆d
j=r+(n−1)∆r/2
M (r) =
X
j=r−(n−1)∆r/2
fedge (j) /
X
fedge (j)
(4.4)
j=rθi −∆d
The normalization is performed so that we can compare results for dif72
Chapter 4. Visualization of the Prostate Gland in VE Images
ferent areas. M is computed along r ∈ [rθi − ∆d, rθi + ∆d], and is expected
to have a strong peak at the edge. In our implementation, parameter values
were set to ∆d = 0.5 cm and n = 5 pixels.
DF statistical test: The performance of the gradient-based edge detector can be plagued by local minima in US images. Therefore, we present
a new approach that models the difference of the radial edge intensity profile as an autoregressive process. The edge strength is characterized based
on the degree of stationarity of this process. The stationarity of the edge
profile is tested using the statistical test proposed by [22]. Each radial intensity profile, was considered as a time series Iθi (k) := Iθi (k∆r), where the
discretized radial distance k, such that k∆r ∈ [rθi − ∆d, rθi + ∆d], replaces
the usual time index.
The edge profiles, Iθi (k), can be modeled as a first order autoregressive
AR(1) processes as follows:
Iθi (k) = ρIθi (k − 1) + er
(4.5)
where ρ is a real number and er is a sequence of independent normal variables with mean 0 and variance σ 2 . In order to show that an AR(1) model
is sufficient for modeling the edge profiles, we computed the partial autocorrelation function (PACF) of Iθi (k). The PACF of an AR(1) process has
significant values only at lag=1. In 78% of the edge profiles extracted from
both B-mode and VE images, at the significance level of 0.05, the PACF
function only has significant values at lag=1.
Iθi (k) is stationary if |ρ| < 1. If a unit root exists (|ρ| = 1), then
the variance of Iθi (k) is rσ 2 , and therefore Iθi (k) is non-stationary. In many
economics applications, the existence of the unit root, which is an indication
of a “trend” or a “shock”, is important for modeling and forecasting the
future observations of a time series. Dickey et al. in [22] provided a statistical
method to test an AR model for the null hypothesis of the existence of a
unit root. If we re-write (4.5) as follows
∆Iθi (k) = Iθi (k) − Iθi (k − 1)
= (ρ − 1)Iθi (k − 1) + er
= γIθi (k − 1) + er
then the DF test is formulated as follows:
H0 : ρ = 1 ↔ H0 : γ = 0
(4.6)
73
Chapter 4. Visualization of the Prostate Gland in VE Images
H1 : ρ < 1 ↔ H0 : γ < 0
(4.7)
Note that the test is performed on the residuals and not the time series
samples. Therefore, the standard t-distribution cannot be used. Dickey and
Fuller provide a non-standard statistic τ , which depends on the number of
observations, and provide tables of critical values for it. In other words,
based on the calculated value of τ , they provide the significance level at
which the null hypothesis above can be rejected. We used the implementation reported in [49] to perform the DF test.
The statistical properties of the intensity profile are altered at the edge.
The existence of a strong edge in the time series is an obvious trend. A signal
with a trend cannot be stationary, since its statistical moments depend on
time, or in our model, on distance. In other words, if one traces the image
in the radial direction, the intensity profile tends to become non-stationary
upon passing through an edge. To evaluate the edge quality, we compute
and report the percentage of edge profiles for which, according to the DF
test, the stationarity hypothesis is strongly rejected (p < 0.05) for the edge
profile. For an image with strong edges, one expects to see a high percentage
of edge profiles for which the unit root exists and the stationarity hypothesis
is strongly rejected.
Edge continuity
Detecting an edge point on a single edge profile can not guarantee the presence of an edge. Continuity of the presence of such a point within a neighborhood is an important factor for the visibility of the edge. We measure
this continuity by measuring the similarity of the edge region within a neighborhood.
For a continuous edge, two neighboring edge profiles are expected to
be similar, although slight differences may be present due to image noise
or in the case of ultrasound, speckle. The normalized cross-correlations
Rθi ,θi ±δθ (r) of the two neighboring edge intensity profiles Iθi (r) and Iθi ±δθ (r),
r ∈ [rθi − ∆d, rθi + ∆d] and the average c(θi )(r) = 21 [Rθi ,θi +δθ + Rθi ,θi −δθ ]
were calculated for each point, rθi , at which the radius rθi intersects the
edge (see Fig. 4.5). The parameters we used for this implementation were
θi = i × 30o , i = 1, ..., 12, ∆d = 0.2 cm, and δθ such that the arc length
between the two adjacent edge profiles is 0.2 cm.
For a large similarity between adjacent edge profiles, cθi should have a
shape similar to the shape of a Gaussian function with large peak at the edge
74
Chapter 4. Visualization of the Prostate Gland in VE Images
point and a small standard deviation. Thus, we propose that the following
edge continuity measure, K(θi ), be calculated:
Kθi =
P 2 (θi )
σ(θi )
(4.8)
in which P (θi ) and σ(θi ) are the peak and standard deviation of a Gaussian
function fitted to cθi .
Volume-based measures
The volume of the prostate is an important parameter in planning the dose
distribution in LDR brachytherapy. A well defined 3D shape can also aid
prostatectomy. To evaluate the shape and size of the prostate created from
VE images, we compare them with the shape and size extracted from the
commonly used B-mode ultrasound images. Since MR images of the prostate
provide more anatomical details including visualizing the boundaries, the 3D
surface extracted from MR images is used as the gold standard. Contouring
in all three image types was performed manually and by three observers: one
radiation oncologist and two trained by experts. The ‘volume difference’ and
‘volume error’ between surfaces created from MR and VE/B-mode images
(i.e. MR vs. VE or MR vs. B-mode) provide shape and size similarity errors
that can be used to compare to the gold standard provided by MRI.
The percent volume difference provides a measure of the difference in
the size of the total gland for each of the nine sectors in VE/B-mode images
compared to the gold standard MRI. It is defined as:
Vdif f % = 100 ×
VV E/B−mode − VM RI
VM RI
(4.9)
The percent volume error is the volume of the non-overlapping region
between the surface of the prostate generated from VE/B-mode images and
that of MR, divided by the sum of the two volumes. This is defined as:
Verr % = 100 ×
(VV E/B−mode + VM RI − 2(VV E/B−mode T M RI ))
VV E/B−mode + VM RI
(4.10)
In other words, it provides a measure of the difference in the shape of
the gland compared to the gold standard. This definition is equivalent to
1 − DSC, in which DSC is the Dice Similarity Coefficient [21].
The 3D shape of the prostate, in each modality, was extracted by manually segmenting the 2D images with the use of the Stradwin [109] software.
75
Chapter 4. Visualization of the Prostate Gland in VE Images
Figure 4.6: Registration of B-mode/VE prostate surfaces to MRI
This software was then used to transform the 2D contours into a 3D surface
c.
which could be opened in MATLAB
To register the MRI and ultrasound images we opted for a rigid registration approach, as opposed to a deformable registration one. Indeed, these
images were not acquired at the same time and the patients were in the
dorsal lithotomy position during B-mode and VE acquisition and supine position during MR acquisition. Furthermore, the patients were anaesthetized
and relaxed prior to the procedure, while they were awake and possibly
tense during the MRI exam. A non-rigid (deformable) registration is not a
suitable registration option since the goal is to understand the differences
between the manually segmented surfaces. A deformable registration which
maps one surface to the other will conceal these differences.
We assume that the main cause of mis-registration of the surfaces is the
angle of the prostates’ base-apex axis with the TRUS probe. This is mainly
a result of the patients’ orientation during imaging. We also assume that
the size of the prostate does not change between data collection sessions and
the prostate does not rotate around its base-apex axis by a large angle. For
each surface, this axis was obtained by fitting a line to the centers of the
prostate contours in each image slice.
With these assumptions, to register the prostate surfaces PB (the VE or
B-mode surface) to PA (the MRI surface), a translation ~tAB is first applied
in order to match their geometric centers CA and CB . A rotation Rθ around
the matched geometric center is then applied to rotate the base-apex axis
of the translated prostate PB to the main axis of PA . θ is the angle between
the base-apex axis of the two surfaces, obtained as described above. Finally
a six degree-of-freedom Iterative Closest Point (ICP) method [8] translates
76
Chapter 4. Visualization of the Prostate Gland in VE Images
(~tICP ) and rotates (RICP ) the resulting surface to fine-tune the registration,
resulting in PB−registered (Fig. 4.6), described as:
PB−registered = RICP [Rθ (PB + ~tAB ) + ~tICP ]
(4.11)
The volume error and volume difference is calculated for the rigidly registered surfaces.
4.3
Results
4.3.1
Image Comparison
The Contrast
The CNR of VE and B-mode images, averaged over the three sections
of the prostate, is shown in Table 4.2. Data from 178 images (from 20 patients) were used in this analysis. In all three regions the CNR of VE is
significantly higher than that of B-mode (p < 0.05 for base, and p < 0.001
for the mid-gland and apex regions).
Table 4.2: CNR comparison of VE and B-mode images
CNR VE
CNR B-mode
Base
10.25 ± 12.83
2.07 ± 1.06
Mid
13.73 ± 5.85
1.43 ± 0.75
Apex
20.51 ± 23.13
1.56 ± 1.16
Edge Strength: Gradient-based measure
The gradient-based measure of edge strength, M , was computed for edge
profiles in the nine regions described in Section 4.2.2. The VE and B-mode
images have different resolutions. This is due to the fact that the window
size used for displacement estimation in VE images is larger than the Bmode pixel. Therefore, the choice of n (the window size), affects the value
computed for M . We examined values of window sizes in the range of n = 2
to n = 10 and in all cases, the results of edge analysis were consistent. The
reported results were acquired for n = 5. The result, extracted from 107
images of 11 patient datasets, is illustrated in Fig. 4.7 for the nine regions
and in Fig. 4.8 for all regions combined. As illustrated in Fig. 4.7, the VE
edges are notably stronger than B-mode edges with the exception of the
anterior base and the posterior apex. The other observation of note is that
the values of M show an evident peak in all areas of VE images, with the
77
Chapter 4. Visualization of the Prostate Gland in VE Images
exception of the anterior base. In the case of B-mode images, in some areas
such as the lateral and anterior apex and the lateral mid-gland, the edges
are very weak and are not represented by a dominant peak in the gradient
values. In the posterior region of the prostate (the last row in Fig. 4.7) the
B-mode edge strength appears to be similar or slightly higher than the VE
edge strength.
It is evident from Fig. 4.8 that, overall, edges in VE images are stronger.
At the edge point, the normalized M value is 2.2 times larger in VE images
compared to B-mode images.
M, (AU)
base−anterior
mid−anterior
1
1
0.5
0.5
−0.5
0
0.5
−0.5
base−lateral
M, (AU)
M, (AU)
0.5
0
0.5
0
0
−0.5
0.5
−0.5
0
0.5
apex−lateral
1
0.5
0.5
base−posterior
1
0.5
1
0.5
0
−0.5
1
mid−lateral
1
B−mode
VE
apex−anterior
0
0.5
0
−0.5
mid−posterior
0
0.5
apex−posterior
1
1
0.5
0.5
0
0
−0.5
0
0.5
−0.5
0
0.5
−0.5
0
0.5
distance from the edge point (cm) distance from the edge point (cm) distance from the edge point (cm)
Figure 4.7: The values of the normalized gradient-based edge strength measure (M ) in arbitrary units vs. the distance from the edge in cm, separated
for the nine regions of the prostate gland.
Edge Strength: DF test of stationarity
Table 4.3 shows the percentage of edge profiles for which the DF test
strongly shows a unit root (p < 0.05), and therefore, non-stationarity. 107
images from 11 patient datasets and 12 edge profiles within each image,
were used in this analysis. This result is reported separately for the nine
regions. The edges appear stronger in VE images in all areas, with the exception of the anterior base (matching the result from the gradient-based
edge strength). In general, the edges in the base region are relatively weak
in both modalities.
78
Chapter 4. Visualization of the Prostate Gland in VE Images
edge filter values in B−mode vs. VE
normalized edge strength M, (AU)
1.5
VE
B−mode
1
0.5
0
−0.5
−0.4
−0.3
−0.2
−0.1
0
0.1
0.2
distance from the edge point (cm)
0.3
0.4
0.5
Figure 4.8: The values of the gradient-based edge strength measure (M )
along the edge profiles for VE and B-mode images. Error bars represent the
inter-patient standard deviation of the M values.
Table 4.3: The percentage of non-stationary prostate edge profiles in different areas of the B-mode and VE images. Standard deviations reported
characterize inter-patient variations.
%
Ant.
Lat.
Post.
Base
67.9±4.6
79.5±3.1
78.6±3.5
VE
Mid
99.2±1.1
98.8±0.9
99.2±1.9
Apex
99.9±0.2
98.1±1.7
99.8±0.2
B-mode ultrasound
Base
Mid
Apex
71.4±2.8 80.8±4.1 76.6±7.5
66.1±3.5 81.3±2.5 88.9±5.1
67.9±3.2 80.8±4.1 84.1±8.1
79
Chapter 4. Visualization of the Prostate Gland in VE Images
Edge Continuity
Table 4.4 shows the edge continuity value, K, computed for VE and
B-mode images in nine regions of the gland. 107 images from 11 patients
were used in this analysis. The edge continuity of VE in all regions except
for the anterior base and posterior apex is shown to be superior to that of
B-mode images. This is also in agreement with the results of the gradient
based measure of edge strength (Fig. 4.7).
Table 4.4: Average edge continuity measure, K, for the nine regions of the
gland.
Ant.
Lat.
Post.
4.3.2
Base
0.64±0.37
0.84±0.39
0.99±0.44
K(θi ) VE
Mid
0.87±0.33
1.15±0.36
0.82±0.35
Apex
0.84±0.37
0.85±0.27
0.60±0.55
K(θi )
Base
0.71±0.52
0.63±0.40
0.52±0.38
B-mode ultrasound
Mid
Apex
0.69±0.41 0.50±0.32
0.48±0.23 0.39±0.29
0.51±0.32 0.75±0.55
Volume Comparison
Figure 4.9 shows an example of VE (magenta) and B-mode (green) 3D
surfaces compared to that of MRI (blue) for one of the patients. Figures 4.10
and 4.11 compare the percent volume error and volume difference between
B-mode and MRI and between VE and MRI prostate surfaces manually
created by one radiation oncologist (Figs. 4.10b, 4.11b) and two individuals
trained by experts (Figs. 4.10a, 4.11a and Figs. 4.10c, 4.11c). Images from
nine patient data sets were used. From Fig. 4.10 it can be seen that in most
regions of the prostate, the mean volume error between VE and MRI is less
than that between B-mode and MRI.
Over the total gland, a volume error of 8.8 ± 2.5% for VE vs. MRI and
10.3 ± 4.6% for B-mode vs. MRI, and a volume difference of −4.6 ± 11.1%
for VE vs. MRI and −4.1 ± 17.1% for B-mode vs. MRI, averaged over nine
patients and three observers, were obtained. However, a one sided analysis of variance does not show any statistically significant difference, which
may also be due to the limited number of patients. In general because the
prostate is not discernible at the base and apex, a larger error in shape
(volume error) and size (volume difference) between VE or B-mode surface
and the gold standard in these two regions compared to the mid-gland is
expected. However, in general, and over the three observers, VE errors are
smaller than those from B-mode.
80
Chapter 4. Visualization of the Prostate Gland in VE Images
Figure 4.9: Comparison of VE (magenta) vs. MRI (blue) 3D surfaces, on
the left, and B-mode (green) vs. MRI (blue) 3D surfaces, on the right, from
one of the patients.
4.4
Discussion and Conclusions
In this chapter ultrasound vibro-elastography was evaluated as an imaging
modality for the visualization of the prostate. We have quantitatively shown
that VE is a promising imaging modality for delineation of the prostatic
region and the use of such data along with B-mode ultrasound can improve
the visualization of the gland.
VE transfer function (TF) images were qualitatively and quantitatively
compared with the commonly used B-mode ultrasound. The evaluation measures used were both image based (CNR, edge continuity and edge strength)
and volume based (volume error and volume difference). A nine sector analysis was used for more detailed characterization. The results on 178 images
suggest that the VE images are significantly superior to B-mode images
in terms of contrast of the gland. The gradient-based and the DF-based
measures of edge strength and the edge continuity measure, on 107 images,
all show that on average VE images provide stronger edges as well. It is
important to emphasize that the statistically significant outcomes of our
analysis of edge continuity and gradient-based edge strength are consistent:
the edges in base-lateral, mid-lateral, mid-anterior, and lateral apex are significantly stronger and more continuous (p < 0.05 for the one-sided analysis
of variance) in VE images compared to B-mode images. The K and M values, in addition to the visual inspection of the images, show that on average
the anterior base and posterior apex of the prostate has weak edges on VE
81
Chapter 4. Visualization of the Prostate Gland in VE Images
70
VE vs. MRI
B−mode vs MRI
60
Verr− Obs 1
50
40
30
20
10
0
Ant.B
Ant.M
Ant.A
Lat.B
Lat.M
Lat.A
Post.B Post.M Post.A
Total
(a) Observer 1: trained by expert
70
VE vs. MRI
B−mode vs MRI
60
Verr− Obs 2
50
40
30
20
10
0
Ant.B
Ant.M
Ant.A
Lat.B
Lat.M
Lat.A
Post.B Post.M Post.A
Total
(b) Observer 2: Radiation Oncologist
Figure 4.10: A comparison between VE vs. MRI volume error and B-mode
vs. MRI volume error, showing the mean and inter-patient standard deviation of Verr % for three observers. Sample data points are also shown as
gray dots. Figure continues on the next page (Ant.: anterior, Lat.: lateral,
Post.: posterior, B: base, M: mid-gland, A: apex).
82
Chapter 4. Visualization of the Prostate Gland in VE Images
70
VE vs. MRI
B−mode vs MRI
60
Verr− Obs 3
50
40
30
20
10
0
Ant.B
Ant.M
Ant.A
Lat.B
Lat.M
Lat.A
Post.B Post.M Post.A
Total
(c) Observer 3: trained by expert
Figure 4.10: continued
images.
The B-mode data used in this work was collected simultaneously with the
RF data used for creating the VE images. In other words, B-mode and VE
images were acquired under exactly similar conditions and are co-registered.
The comparison of volumes calculated from VE and B-mode prostate
images show that VE volumes are closer to the MRI gold standard in most
regions of the prostate, both in terms of shape and size. This confirms that
the outlined region used in the image-based evaluation of VE, is indeed the
prostate. In this comparison, prostate images from 9 patients were manually
delineated by three observers (one expert and two trained by experts). The
total gland volume error for VE vs. MRI was 8.4 ± 2.9%, 8.3 ± 1.8%, and
9.8 ± 2.8% for the three observers. For B-mode vs. MRI, these values were
8.6 ± 4.8%, 11.3 ± 4.7%, and 11.0 ± 4.4%. The total gland volume difference
for VE vs. MRI was −5.3 ± 11.7%, −8.1 ± 8.8%, and −0.3 ± 12.1% for
the three observers. For B-mode vs. MRI, these values were −8.7 ± 13.8%,
−1.4 ± 20.2%, and −2.1±17.7%.
It is also worth noting that Vdif f (see Equation 4.9) measures the prostate
volume ratios for B-mode/MRI and VE/MRI (Vdif f = VB−mode /VM RI − 1).
83
Chapter 4. Visualization of the Prostate Gland in VE Images
80
VE vs. MRI
B−mode vs MRI
60
Vdiff− Obs 1
40
20
0
−20
−40
−60
−80
Ant.B Ant.M Ant.A Lat.B Lat.M Lat.A Post.B Post.M Post.A Total
(a) Observer 1: trained by expert
80
VE vs. MRI
B−mode vs MRI
60
Vdiff− Obs 2
40
20
0
−20
−40
−60
−80
Ant.B Ant.M Ant.A Lat.B Lat.M Lat.A Post.B Post.M Post.A Total
(b) Observer 2: Radiation Oncologist
Figure 4.11: A comparison between VE vs. MRI volume difference and
B-mode vs. MRI volume difference, showing the mean and inter-patient
standard deviation of Vdif f % for three observers. Sample data points are
also shown as gray dots. Figure continues on the next page (Ant.: anterior,
Lat.: lateral, Post.: posterior, B: base, M: mid-gland, A: apex).
84
Chapter 4. Visualization of the Prostate Gland in VE Images
80
VE vs. MRI
B−mode vs MRI
60
Vdiff− Obs 3
40
20
0
−20
−40
−60
−80
Ant.B Ant.M Ant.A Lat.B Lat.M Lat.A Post.B Post.M Post.A Total
(c) Observer 3: trained by expert
Figure 4.11: continued
[98] report a MR/US prostate volume ratio of 1.11 ± 0.10 averaged over
prostate surfaces from 10 patients outlined by 7 observers, repeated twice.
Similar results from other work have been listed in a table therein. The
MR/US prostate volume ratios for the 9 patients in our database are 1.12 ±
0.18, 1.05 ± 0.21, and 1.05 ± 0.19, for the three observers, which support
previous work. For MR/VE, these values are 1.07 ± 0.13, 1.10 ± 0.10, and
1.02± 0.13. However, the volume ratio or Vdif f only report differences in the
size of the prostate, whereas two prostate surfaces can be greatly different
in shape but have equal volumes. Verr , which we have included here, can
provide useful characterization of the shape difference.
Currently, one dimensional axial motion estimation is used for creating
the VE images. Any out of plane motion, including the lateral motions
caused by the slippage between the protective sheath on the probe and the
surface of the rectum during VE data collection, can affect the VE images.
This can be improved upon by using 2D motion estimation techniques to
account for lateral tissue motion. We have not used 2D motion tracking
because most of the displacement estimates between consecutive frames of
RF data are sub-sample. It has been shown that the gain based from 2D
motion estimation in such cases is not significant enough to warrant its
use [121] unless beam-steering is employed, such that the lateral motion
is estimated from two axial measurements at different beam angles, e.g.
85
Chapter 4. Visualization of the Prostate Gland in VE Images
±10 degrees. The increase in the number of axial measurements lowers the
sampling rate and we opted not to do that in this set of experiments.
One should also keep in mind that the presence of the TRUS probe
results in deformation of the posterior region of the gland in both VE and
B-mode surfaces. For more accurate registration of these surfaces with those
from MRI, this deformation could be accounted for. However, since the
B-mode and VE images that we use are exactly registered to each other,
this deformation will only result in a bias in the volume error and volume
difference, especially in the posterior region, with MRI surfaces. The result
of such a bias can be seen to be a generally larger mean and standard
deviation in the volume error and volume difference of the posterior region
in all three observers, compared to that of the lateral and anterior regions.
As a final note, the best validation of the prostate gland segmentation
could be a comparison to pathology, the undisputed gold standard. Such a
study requires access to the whole mount pathology, and viable solutions to
the open problem of how to register the shrunk, misaligned pathology slices
to pre-operative images.
86
Chapter 5
Automatic prostate
segmentation using fused
ultrasound B-mode and
vibro-elastography images
5.1
Introduction
In this chapter an automatic method for segmenting the prostate is proposed. In this method, which is an extension of the semi-automatic method
of Chapter 2, a combination of vibro-elastography and B-mode images is
used to eliminate the need for manual initialization.
The elimination of user interaction in a prostate segmentation algorithm
is desirable for various reasons. Manual interaction introduces user variability due to different user preferences and level of experience. An automatic
algorithm can also be more time-efficient. Furthermore, it can be integrated
into an automatic treatment method and be used intra-operatively.
In LDR prostate brachytherapy if treatment plans are created pre-operatively,
the placement of the TRUS probe is initially adjusted intra-operatively to
provide images of the prostate similar to the pre-operative volume study
images. This is necessary due to the fact that the orientation (mainly due
to patient positioning), shape and size of the prostate may change from
the pre-operative volume study to the intra-operative imaging sessions. Occasionally, despite the adjustments, the treatment plan may need to be
updated by the physician by adding or relocating seeds to adhere to the
changed prostate boundary. The adjustment of the images and the plan is
performed manually, which is time consuming and not free from uncertainty.
If intra-operative planning, an established procedure used in some centers,
is performed [1, 4, 71], the need for such adjustments is of less concern.
Nevertheless, delineation of the prostate in a time-efficient way is needed for
Chapter 5. Automatic Segmentation of the Prostate
intra-operative planning due to time limitations inside the operating room.
Its application is not only limited to LDR brachytherapy, but can also be
useful in high dose-rate brachytherapy or external beam radiation treatment
and biopsy planning.
In Chapter 4 we showed that images produced by the vibro-elastography
method developed in our group [90, 112] have superior prostate-background
contrast compared to regular B-mode images, a result of the prostate tissue
being generally stiffer than the surrounding tissue. However in some regions,
and especially in the posterior region, the boundary of the prostate is better
defined in B-mode images. In this chapter, we make use of the advantages
of each image type and show that with the combination of elastography
images and B-mode images we can achieve automatic 3D segmentation of
the gland. The approach we use replaces the semi-automatic segmentation
of the mid-gland slice in the method described in Chapter 2 with a fully
automatic approach using an Active Shape Model (ASM).
To the best of our knowledge, this is the first report on fusion of vibroelastography and B-mode images for prostate segmentation.
5.2
5.2.1
Methods and Materials
Data Collection and Processing
The B-mode and vibro-elastography images used in this chapter were acquired from patients going through the standard LDR prostate brachytherapy procedure at the Vancouver Cancer Center, British Columbia Cancer
Agency or the radical prostatectomy procedure at the Vancouver General
Hospital. Institutional ethics approval and informed patient consent were
obtained prior to data collection.
Intra-operatively, prior to the procedure, RF data and B-mode images
were simultaneously collected using the system described in [62]. RF data
were collected using a Sonix RP ultrasound machine (Ultrasonix Medical
Corp., Richmond, BC, Canada) at approximately 40 frames per second with
the sagittal array of a vibrating transrectal ultrasound (TRUS) probe (dualplane linear/microconvex broadband, 5-9MHz). The TRUS probe was vibrated with an amplitude of 0.5-2 mm peak-to-peak and frequency range of
2-10 Hz while being rotated along its longitudinal axis from −45◦ to 45◦ to
collect a series of RF data planes spanning a 3D volume in the shape of a
cylinder with a sector base as shown in Fig. 5.1. The angle step size was 2◦ ,
and data were collected for 2 seconds per angle, resulting in approximately
80 frames at each angular location. The total duration for collecting the 3D
88
Chapter 5. Automatic Segmentation of the Prostate
Figure 5.1: Illustration showing the span of the RF data planes with respect
to the TRUS probe.
data was approximately 1.5 minutes.
Each RF data frame has 128 axial lines corresponding to the 128 probe elements. These frames were processed using the method described in [90], resulting in a time series of strain images of size 128 lines by 100 blocks for each
angle for the period of vibration. For every angle one vibro-elastography image was then obtained as follows:
V Ejk
1
=
f2 − f1
Z
f2
f1
jk
Href (j2πf ) df
(5.1)
jk
(j2πf ) is the transfer function (TF) of the strain computed using
where Href
standard signal processing methods from a reference point in the image
(computed as the spatially averaged strain, at half the tissue imaging depth)
to a block at axial location k and lateral location j in the strain image.
[f1 , f2 ] is the range of the broad-band vibration applied to the tissue. Further
details are provided in [90] and [61].
By interpolation, 3D matrices with voxel intensities representing tissue
stiffness were achieved from the longitudinal TF images.
The B-mode and vibro-elastography 3D volumes are matrices of size
162×275×128 pixels and 100×169×128 pixels respectively. The 128 transverse images, corresponding to the 128 transducer elements, are spaced at
0.43 mm and the pixel sizes in B-mode and vibro-elastography are 0.31 mm
and 0.5 mm, respectively. For 3D segmentation, to increase the speed of the
algorithm, we only used 32 out of the 128 transverse images, each spaced at
89
Chapter 5. Automatic Segmentation of the Prostate
1.73 mm.
Throughout the algorithm we make use of B-mode and vibro-elastography
phase symmetry images. The code and description of the calculation of
phase symmetry images are available in [54]). The B-mode and VE phase
symmetry images are used in the initialization of the algorithm while the
weighted sum of pixel intensities in B-mode phase symmetry and B-mode
image is used in the deformation of the posterior part of the prostate contour. We observed that this combination highlights the boundary of the
prostate particularly well in the posterior region. Figure 5.2 shows a Bmode image, its phase symmetry image and the weighted sum of the two.
A vibro-elastography image and its phase symmetry image are also shown.
5.2.2
2D Segmentation Algorithm
We use an Active Shape Model (ASM) approach [18] for 2D segmentation
of the prostate. To guide the deformation of the model at each iteration
we use a combination of information from vibro-elastography and B-mode
ultrasound images. The algorithm steps are summarized in Table 5.1 and
each step is further discussed below.
The ASM algorithm consists of an offline phase in which the training set
is created and the statistics of the model are extracted and an online segmentation phase where the model is evolved based on image characteristics.
In the offline training phase:
1. Each contour in the training set is created from 30 manually selected
specific points on the prostate boundary, X = [x1 · · · x30 ; y1 · · · y30 ]T . The
training set consists of N = 25 contours from the mid-gland images of 7
patients. The images were chosen to cover a broad range of prostate shapes
and none of the training set images were used in the segmentation step.
Figure 5.3 shows the location of the 30 points on the prostate boundary. 10
of these points are selected from specific anatomical locations (mid-posterior,
mid-anterior, lateral left and right, lowest points and lower left (right) points
of the lobes and two points between the mid-anterior and lateral points).
The remaining 20 points are selected to be at equal distances between every
adjacent pair of the ten points specified above.
2. The manually segmented contours are then aligned and scaled using
a least-squares-based iterative scaling, rotation and translation method and
the mean shape, Xmean , is the mean position of all corresponding boundary
points. The distribution of the prostate contours and the 10 main points in
the training set, after alignment, are shown in Fig. 5.3.
90
Chapter 5. Automatic Segmentation of the Prostate
Figure 5.2: A B-mode image (top left) and the corresponding vibroelastography image (top right) of the prostate mid-gland along with the
phase symmetry of the B-mode (middle left) and vibro-elastography (middle right). The sum of the B-mode image and its phase symmetry is shown
in the lower left.
91
Chapter 5. Automatic Segmentation of the Prostate
Table 5.1: The main steps in the ASM algorithm.
Offline training
1. Create a training set of manually marked gland boundary points.
2. Align and scale the prostate contours and compute the mean shape Xmean
3. Extract the shape statistics from the training set by computing the modes
of variation of the prostate shape.
Online segmentation of the target image
4. Obtain the initial shape Xinit using Xmean and target image information.
Set Xcurrent = Xmean
5. Find the desired shape deformation:
dXdesired = Xdesired − Xcurrent
5.(i) For non-posterior boundary points use a combination of gray
level similarity and edge continuity measures from B-mode and vibroelastography images.
5.(ii) For posterior boundary points use intensity gradient information
from B-mode and its phase symmetry image.
6. Calculate the required adjustments in the pose parameters which best
map Xcurrent to Xdesired .
For that, find the rotation, scale and translation adjustments dθ, ds, and dt
which minimize the weighted sum
E = (Xcurrent − Tds,dθ,dt (Xdesired ))T (Xcurrent − Tds,dθ,dt (Xdesired ))
using a least-squares approach, where if, for example applied to a point
[x, y]T
#
#"
# "
"
x
s cos(θ) −s sin(θ)
tx
T
+
Ts,θ,t ([x, y] ) =
y
s sin(θ) s cos(θ)
ty
The residual shape adjustments, dx, can then be made by deforming the
shape model.
7. Calculate db = P T dx to obtain the required shape adjustments in the
model parameter space.
8. Update the shape parameters: b → b + db and the pose parameters:
s → s(1 + ds), θ → θ + dθ, t → t + dt .
9. Update the shape: Xnext = Ts,θ,t(Xinit ) + P b, which consists of applying
the updated pose parameters (Tθ,s,t (.)) and the shape parameters (b) to the
initial shape.
10. Set Xcurrent =Xnext and repeat steps 5-10 until convergence.
92
Chapter 5. Automatic Segmentation of the Prostate
Figure 5.3: The 30 boundary points used to model the prostate boundary
(left) highlighting the 10 main points as red circles and blue squares and the
aligned and scaled, manually segmented prostate boundaries in the training
set, showing only the 10 main points (right). We define the posterior points
as the main points shown as blue squares and the four points in between
them.
3. To capture the statistics of the training set, the covariance matrix is
computed according to (5.2):
S=
N
1 X
dxi dxTi
N i=1
(5.2)
where dxi , in the training phase, is the vector between each point on the
manual contour and the corresponding point on Xmean . The modes of variation of the shape are described by the eigenvalues λi , and eigenvectors pi ,
of S, from which the nλ largest eigenvalues are selected as the most significant modes of variation. We selected nλ such that for i = 1, · · · , nλ ,
P
λi / N
i=1 λi > 5%. For our training set nλ = 3 satisfied this condition.
A shape instance consistent with the training set can thus be created
from an initial contour Xinit using (5.3):
X = Xinit + P b
(5.3)
where P = (p1 · · · pnλ ) is the matrix of the first nλ eigenvectors of S and
b = (b1 · · · bnλ )T is a vector of weights, also called the shape parameters.
For each point i of each image j of the training set, a normalized edge
derivative profile gij , normal to the boundary, centered at the model point,
and of length 2np is extracted (Fig. 5.4). This is the derivative of the edge
93
Chapter 5. Automatic Segmentation of the Prostate
intensity profile divided by the sum of its absolute values over the length of
the profile. From gij , the mean edge derivative profile g¯i , (the average over
all corresponding boundary points in the training set) is computed and the
difference between the derivative profiles and g¯i is used to compute the covariance matrix Sgi [17]. These provide a statistical description of the gray
level appearance for every location of the prostate boundary, as illustrated
in Fig. 5.4, and are used in computing the edge gray level similarity measure,
dg , described further in the chapter in step 5.(i).
Figure 5.4: Illustration of how for point i on the model edge of training
image j, a normalized edge derivative profile, gij normal to the boundary
and centered at the model point and of length 2np is extracted. The average
of gij over the number of training images results in g¯i .
In the online segmentation phase:
The goal is to find the pose (scaling, translation and rotation) and shape
parameters (deformations within the model constraints) which cause the
model to deform into the structure of interest. That is, starting from the
initial shape Xinit , at each iteration, we aim to deform the current shape
points Xcurrent into a desired shape Xdesired obtained from image features.
4. To obtain the initial shape Xinit , the mean prostate shape, Xmean ,
from the training set is scaled and translated to provide a closer shape
match to the shape instance of the prostate to be segmented and hence,
faster algorithm convergence. This is done by first fitting ellipses to the
B-mode and vibro-elastography phase symmetry images (Fig. 5.2, second
row). Using the size and center of the smaller of the two fitted ellipses we
scale and center Xmean to obtain Xinit . We use the scale and center position
of the the smaller ellipse since it was seen that, in some cases, the presence
of anatomical features in the prostate anterior may falsely be detected as
94
Chapter 5. Automatic Segmentation of the Prostate
part of the prostate during the ellipse fitting process resulting in a larger
than expected ellipse.
5. The computation of dXdesired , the desired change in the boundary
point locations to obtain a shape which best fits the prostate, is described
next.
6. Having computed dXdesired , first the pose parameters are adjusted
using a least-squares approach described in [18] to find the rotation, scaling
and translation which best maps Xcurrent to Xdesired , the difference being
the residual dx.
7. Then, the shape parameters, db, are adjusted given dx and using (5.4)
to obtain the shape closest to Xdesired which adheres to the model:
db = P T dx
(5.4)
Based on our image characteristics, we compute dXdesired for the posterior and non-posterior points using two different methods. A combination
of candidate motions from both B-mode and vibro-elastography images is
used for the non-posterior points. In the posterior portion of the prostate
boundary vibro-elastography images show low contrast but this boundary
is more visible in B-mode images. We attribute this low contrast mainly to
the slippage between the probe, its protective sheath and the soft surface
of the rectum and the use of one dimensional axial tissue strain estimation.
Therefore in this region we use a gradient-based motion on the weighted
sum of the pixel intensities in B-mode and their phase symmetry images.
Posterior points include the three posterior main points and the four points
in between them as shown in Fig. 5.3 and the non-posterior points include
the remaining 23 boundary points.
5.(i) To calculate the movement, dXdesired , for the non-posterior model
points, we use a combination of four candidate motions, as defined in (5.5),
based on the Mahalanobis distance [18] and an edge continuity measure that
we derived in [62]. These two measures are computed from both B-mode
and vibro-elastography images for each non-posterior boundary point i:
dXdesired = df un
(5.5)
df = αgU S dgU S + αgE dgE + αKU S dKU S + αKE dKE
where un is the unit vector normal to the boundary at the boundary point,
dgU S and dgE are the motions of the boundary point along un derived
from the gray level similarity in B-mode and vibro-elastography images,
95
Chapter 5. Automatic Segmentation of the Prostate
and dKU S and dKE are the motions of the boundary point along un computed from the edge continuity measures in B-mode and vibro-elastography
images, respectively (index i dropped for clarity). A set of weights α =
[αgU S αgE αKU S αKE ]= [0.15 0.25 0.25 0.35] that works well was found
by trial and error.
Details of how each of the two measures dg and dK are computed are as
follows:
The edge gray level similarity measure, dg
During the shape fitting phase, at each iteration and for every point,
sample edge derivative profiles hi (d) (d = −l · · · l) of length 2np are extracted
at the current boundary point similar to the model edge derivative profile
(g¯i ) which was computed in the offline phase. The square of the Mahalanobis
distance between these profiles from the model profile, provides a measure
of edge gray level similarity, and is maximized:
n
o
(hi (d) − g¯i )
dg = arg max = (hi (d) − g¯i )T Sg−1
i
d
(5.6)
for each point, to compute dgU S (or dgE ) (see Fig. 5.5).
The physical values of l and np are set to 8 and 2.5 mm, respectively.
These particular values were used to ensure both adequate coverage of the
edges and reasonable profile lengths in pixels (l=16 pixels and np =5 pixels
for a vibro-elastography image).
The edge continuity measure, dK
We observed that the gray level edge similarity measure alone gives many
false positives due to ultrasound speckle or sharp edge-like structures since
only 1D information (normal to the edge) is being analyzed. Therefore, in
our approach, we incorporate the edge continuity measure proposed in [62],
where an edge must not only be detected along the direction of the edge
profile, but must also be present in its neighborhood.
To compute the edge continuity at a distance d (d = −l · · · l) from each
boundary point i and along un , the edge intensity profile ei (d) of length 2np
is extracted. The two neighboring edge intensity profiles, e1i (d) and e2i (d),
are of the same length 2np , parallel to and spaced at 1.5 mm on either sides
of ei (d), as shown in Fig. 5.5. We compute the average normalized crosscorrelations, C1i (d) and C2i (d), between the edge intensity profile ei (d) and
its two neighboring edge intensity profiles.
96
Chapter 5. Automatic Segmentation of the Prostate
Figure 5.5: Illustration showing the computation of gray level similarity
and edge continuity for two boundary points. In the upper point the edge
derivative profile hi (d) is compared with that of the model, g¯i (shown in
the rectangle) along a line normal to the edge and of length 2l. dgE is the
distance from the current boundary point where the maximum similarity is
obtained. In the lower point the gray level intensity profile, ei (d), is compared to its two neighbors, e1i (d) and e2i (d), to obtain correlation functions,
C1i (d) and C2i (d). dKE is the distance from the current boundary point to
where a maximum correlation exists between the center profile and its neighbors. A similar procedure is applied to B-mode images to obtain dgU S and
dKU S .
For a continuous edge (i.e. large similarity between ei (d) and its neighboring profiles), C1i (d) and C2i (d) should each have a shape similar to a
Gaussian function with a large peak and a small standard deviation. We
define the edge continuity measure as follows:
h
i
Ki (d) = PC21i + PC22i /2σi (d)
(5.7)
in which PC1i and PC2i characterize the peak and σi is the standard deviation
of a Gaussian function fitted to C1i (d) and C2i (d). For each point, we
compute:
97
Chapter 5. Automatic Segmentation of the Prostate
dK = arg max {K(d)}
d
(5.8)
Figure 5.5 displays how the grey level similarity and edge continuity measures are computed for two boundary points on the prostate contour, respectively. This is shown for a vibro-elastography image but is computed in the
same manner for B-mode images. Due to the difference in image resolution
we ensured that corresponding parameters, such as profile length, match in
terms of physical dimensions in both B-mode and vibro-elastography images.
5.(ii) To compute dXdesired for the posterior points we use the weighted
sum of the pixel intensities of the B-mode and their phase symmetry images:
IBmode (m, n) + 0.5Iphase (m, n) where IBmode (m, n) and Iphase (m, n) are the
pixel intensities of the B-mode image and its phase symmetry at location
(m, n) (assuming an identical intensity range in both images). (see Fig. 5.2).
In the posterior region the edge intensity profile, ei (d), is extracted as discussed before and the desired location of the boundary point in the next
iteration is set to the point along the profile normal to the edge where the
gradient of ei (d) is maximized:
dXdesired = dgradU S un
dgradU S = arg maxd {grad(ei (d))}
(5.9)
where grad(ei ) is the gradient of the edge intensity profile.
The algorithm stops when either the maximum number of allowed iterations has been reached (set to 50) or when at least 96% (29 out of 30) of the
boundary points for 2D evaluation, and 94% (28 out of 30) of the boundary
points for 3D evaluation, have a small dXdesired (we chose a value less than
2.5 mm).
5.2.3
3D Segmentation
For 3D segmentation we use a modification of the method in [60], summarized here for completeness (see Fig. 5.6):
98
99
Chapter 5. Automatic Segmentation of the Prostate
Figure 5.6: The main steps of the 3D semi-automatic prostate segmentation algorithm proposed in [60] including
the applied modifications. Manual initialization is replaced by automatic segmentation of the mid-gland and
vibro-elastography images are included in the segmentation process.
Chapter 5. Automatic Segmentation of the Prostate
The prostate is modeled as a warped and tapered ellipsoid. Warping
is present in the prostate posterior due to the pressure from the TRUS
probe. Tapering is present in the anatomical shape of the gland, both in
the transverse plane and from the base to the apex. The aim of the 3D
segmentation method is to fit such a shape to a set of 8-12 transverse Bmode images that are collected from the base to the prostate apex prior to
the procedure.
The user manually selects six specific boundary points on the mid-gland
image. With the aid of these points, the transverse images are first unwarped and then untapered, resulting in prostate cross section images that
are elliptical in shape, and the mid-gland image is segmented. The prostate
is now modeled as an egg shape with an elliptical cross section. The remaining images are segmented by first fitting two semi-ellipsoids to the mid-gland
contour and the base and apex positions (usually defined during imaging),
slicing them at image depths and using the obtained ellipses as guides for the
IMMPDA edge detector [2]. Finally an ellipsoid tapered along the base-apex
axis is fitted to all the detected edges, sliced at image depths and tapered
and warped (using the inverse of the same untapering and unwarping values
initially used), to match the original images. Figure 5.6 displays the main
steps of the 3D original algorithm along with the modifications applied in
this chapter.
In order to automate the aforementioned algorithm, we replace the manual initialization step with the 2D automatic segmentation method proposed
in Section 5.2.2 by replacing the six manually selected boundary points with
the 30 boundary points of the ASM-based contour. Additionally, we use
vibro-elastography images, in addition to B-mode, for the ASM-based 2D
segmentation and when segmenting the remaining slices. That is, after slicing the two fitted semi-ellipsoids to the mid-gland contour and the base and
apex locations, the resulting contours are used to guide the IMMPDA edge
detector in the vibro-elastography images first and then, the results are used
to guide the IMMPDA in the B-mode images.
5.2.4
Evaluation
We compare the automatic segmentation results with manually segmented
contours from ultrasound B-mode images. Manual segmentations are performed by two trained individuals (author and colleague) , blind to the
automatic results. A total of 61 co-registered pairs of B-mode and vibroelastography images from 17 cases are used for 2D evaluation. Since the
ultimate goal is to automate the 3D semi-automatic segmentation method
100
Chapter 5. Automatic Segmentation of the Prostate
which requires initial segmentation of the mid-gland slice, the 2D evaluation is performed on mid-gland images. Out of the 17 cases available, the
whole prostate is reasonably visible in 11 cases. We use these 11 cases for
evaluation of the 3D algorithm. In the remaining 6 cases the prostate was
not completely visible either due to artifacts (e.g. air bubbles covering the
TRUS probe or insufficient contact of the probe with the prostate) or the
prostate being out of the field of view (when the prostate is large or the
TRUS probe is not appropriately positioned causing the full length of the
gland not being covered).
For 2D evaluation, the measures used for comparing the manually created contour Cman , with that of the automatic contours Cauto , are as follows:
- MaxD: The maximum absolute radial distance between Cman and Cauto ,
in a slice.
- MAD : The average absolute radial distance between contours Cman and
Cauto , in a slice.
- DSC : The Dice Similarity Coefficient [21] equal to
T
2 (Aman Aauto ) / (Aman + Aauto ), where A denotes area of contours.
For 3D evaluation, the measures used for comparing the manually created volume Vman , with that of the automatic Vauto , are as follows:
- Verr : Volume error. The absolute volume of the non overlapping region
between two volumes defined as:
Verr = |(Vauto + Vman − 2(Vauto∩man ))| /(Vauto + Vman )
(5.10)
- Vdif f : Volume difference. The difference between the volumes of two
delineated prostates defined as:
Vdif f = (Vauto − Vman )/Vman
(5.11)
These measures are computed for both the whole gland and the base,
mid-gland, and apex regions (being respectively, 0.3, 0.4, 0.3 of the length
of the gland).
101
Chapter 5. Automatic Segmentation of the Prostate
Table 5.2: 2D segmentation results
5.3
5.3.1
MaxD (mm)
MAD (mm)
DSC
Obs. 1
6.49±3.48
2.55±1.41
0.87±0.07
Obs. 2
7.04±3.36
2.66±1.50
0.87±0.08
Results
2D Segmentation
Table 5.2 shows the results of comparing the 2D automatic and manual
contours (two observers) for the 61 mid-gland images. The duration of the
algorithm was 13.3±12.3 s and the percentage of cases which converged
within the maximum number of iterations (50) was 80.3%.
To show the benefits of using vibro-elastography images, we repeated
the evaluation using the weights α =[0.5 0 0.5 0], that is using gray level
similarity and edge continuity in only B-mode images for the non-posterior
points (motion of the posterior points stays the same). In only 19.7% of
the images the algorithm converged when vibro-elastography images were
not used. The MaxD and MAD values of the 19.7% converged cases were
11.41±4.22 mm and 4.59±1.43 mm for observer 1 and 12.09±4.07 mm and
4.93±1.56 mm for observer 2.
5.3.2
3D Segmentation
Table 5.3a and 5.3b display the volume error and volume difference between
the 3D automatic segmentation surfaces and manually created surfaces by
both observers for three regions of the prostate (base, mid-gland and apex)
and the whole gland. In Table5.3c, the volume of the manually and automatically created surfaces, and the volume of the non-overlapping region
between them are shown in cm3 . The Stradwin software [109] was used to
create the surfaces from the segmented images. The 3D segmentation duration, that is the time required for segmenting the 3D volume of the prostate
starting from the 2D segmented mid-gland image, was 31.7±3.3 s.
Figure 5.7 shows an example of 2D automatic segmentation of the midgland. The initial contour Xinit which is the average of the training set
102
Chapter 5. Automatic Segmentation of the Prostate
Table 5.3: 3D segmentation results
(a) Verr (%)
Base
Mid-gland
Apex
Whole gland
Obs. 1
12.7±3.0
7.7±2.7
13.1±3.8
10.2±2.2
Obs. 2
13.4±3.2
10.2±3.6
21.4±9.6
13.5±4.1
(b) Vdif f (%)
Base
Mid-gland
Apex
Whole gland
Obs. 1
-2.5±12.1
-9.2±5.5
-7.4±16.8
-7.2±9.1
Obs. 2
-7.3±11.1
-14.0±9.5
-14.0±30.2
-13.3±12.6
(c) Volume error, and volume of automatic and manual segmentation (cm3 )
Obs. 1
Base
Mid-gland
Apex
Whole gland
Vol. err.
Vol. man.
3.0±1.0
12.4±5.2
3.2±1.6
22.3±8.3
2.3±1.3
8.7±3.8
8.5±3.7
43.5±17.1
Obs. 2
Base
Mid-gland
Apex
Whole gland
Vol. err.
Vol. man.
3.3±1.3
13.1±5.5
4.3±1.7
23.3±7.8
3.2±2.2
9.9±4.3
11.4±4.4
46.3±16.9
Vol. auto
11.8±4.2
20.1±7.1
7.8±2.7
39.7±13.9
contours scaled and centered according to the phase symmetry B-mode and
vibro-elastography images, is shown in dashed yellow. The final contour (in
this case, converged after 6 iterations) is shown in solid red. As shown in
the figure, the deformed contour can extend past the field of view. This was
103
Chapter 5. Automatic Segmentation of the Prostate
allowed so that the missing lateral regions of the prostate have less affect on
segmentation. For this purpose, in every iteration, any boundary point that
falls outside the field of view (i.e. falls in the black region) is given a lower
weight when obtaining the deformed shape. The 3D segmentation result
initialized by the 2D automatic contour is also shown (blue triangulated
surface) and compared to the manually segmented surface (red).
Figure 5.7: 2D segmentation results shown on corresponding B-mode (left)
and vibro-elastography (right) mid-gland images. The initial contour Xinit
is shown in dashed yellow and the final contour in solid red. The final
3D segmented surface (triangulated blue) is also compared to the manually
created surface (red) in the lower image.
104
Chapter 5. Automatic Segmentation of the Prostate
5.4
Discussion and Conclusions
In this chapter we presented an automatic 2D and 3D prostate segmentation algorithm using a combination of both ultrasound B-mode and vibroelastography images. The high contrast of vibro-elastography images aids
segmentation especially since B-mode image quality is known to be relatively poor. To the best of our knowledge the combination of these two
image types and an edge continuity measure have not been used for prostate
segmentation before.
As we will discuss below, with similar volumetric errors and algorithm
duration, the advantage of the current method over the previously proposed
semi-automatic method is the elimination of manual initialization by replacing the manually selected initial boundary points with an automatically
created mid-gland contour. The use of vibro-elastography images in the
segmentation procedure was found to be necessary as we showed in our 2D
evaluation that when only B-mode images are used, algorithm convergence
is achieved in only 19.7% of the cases.
Various methods are proposed for computing dXdesired . In [39] the goal
is to minimize a cost function based on a grayscale profile normal to the
edge. This method, although simple and fast, can fail where the boundary
is absent or weaker than a non-boundary grayscale transition. In [17] the
square of the Mahalanobis distance is used, a measure of gray level similarity
at the boundary point and normal to the boundary.
Even though due to the difference in data sets and evaluation methods,
direct comparison is not reasonable, we will provide some reported results
from the literature. The 2D ASM-based method proposed in [39] results in
an average mean absolute distance and maximum distance of 1.09 ± 0.49 mm
and 7.27 ± 2.32 mm between manual (the average of three repeated manual
delineations of three trained graduate students) and automatic contours (on
36 patient data sets) in their ASM-based method. In the 3D extension of
this method, after each 2D image is segmented, manual editing (required in
26.3% of the 2D images in each case) is allowed, after which the modified
points are clamped and the model is re-deformed. In the 3D semi-automatic
method proposed by [113] a 83.5% percent volume overlap (the intersection
divided by the union of two volumes) between automatic and manual (average of three observers’ manual contours) is reported on a data set of 30
patients. In [125] a mean overlap volume error of 4.16% compared to manual
segmentation on six patients is reported. [117] report a mean absolute distance of 1.65±0.47 mm and a Dice similarity coefficient value of 0.91±0.03
between automatic and manually delineated TRUS video sequences of 19
105
Chapter 5. Automatic Segmentation of the Prostate
patients.
Our results are comparable to the above, however, with standard quality TRUS images we expect lower errors. In the 3D semi-automatic method
that we proposed in [60] higher quality B-mode images, equivalent to those
displayed on the ultrasound machines’ screen, were used. But the B-mode
images that we were able to save along with the RF data to use in this work
were of lower quality. That is because we did not have access to the additional post processing steps applied to the images before they are displayed
on the machine. Nevertheless, the setup used for collecting RF data is a
modified standard brachytherapy setup and the vibro-elastography images
are created from RF data collected simultaneously with the B-mode images.
Therefore, the B-mode and vibro-elastography images are co-registered.
For a better comparison and to see the effect of image quality on the
results, we ran the semi-automatic segmentation algorithm on the same 11
cases used in this chapter and performed volumetric comparison with the
manual surfaces created by the two observers. The whole gland volume error between the semi-automatic surfaces and manual surfaces created by the
two observers were 10.0±2.6% and 13.6±3.9% and the volume differences
were -2.2±12.3% and -8.4±16.9% for 11 cases. By comparing these results
with those reported in [60] we observe higher error values (volume error of
6.63±0.90% and volume difference of 2.43±6.08%). We attribute this to the
lower quality of the images compared to those used for evaluation in our
previous work. On the other hand, these errors are similar to those that we
report in Table 5.3 (volume errors of 10.2±2.2% and 13.5±4.1% for the two
observers and volume differences of -7.2±9.1% and -13.3±12.6%). In [63]
our dosimetric analysis on treatment plans created based on semi-automatic
volumes showed that the semi-automatic method is a suitable replacement
for manual segmentation in prostate brachytherapy. The similar errors suggest that if higher quality B-mode images, similar to those used in [60], were
used, the surfaces created from the automatic method can also be a suitable
replacement.
One of the sources of error, seen in our 2D results, is the presence of the
anterior peri-prostatic venous plexus. This can be a fairly large structure,
which collapses when circulation is stopped. It is seen in our in-vivo data
as a dark region in VE and B-mode images and therefore, can sometimes be
mistaken as the prostate anterior during automatic segmentation, especially
when it has a contrast similar to that of the gland. For example, see Fig. 5.8.
In this case, the anatomical structure, forming a shadow anterior to the
prostate, has resulted in a large segmentation error. This structure is also
present in Fig. 5.7, but has not caused a problem. A second source of error is
106
Chapter 5. Automatic Segmentation of the Prostate
the limited field of view (-45 to 45◦ in the transverse view, due to mechanical
limitations in our setup) which in some cases does not entirely cover the
lateral portions of the prostate, especially in the mid-region (for example,
see Fig. 5.7). As a result, in some cases we observe automatic volumes
smaller than that from manual segmentation. Even though this missing
information has been accounted for to some extent in the 2D segmentation
method (seen as the contours extending beyond the field of view, as Fig. 5.7
shows) it results in smaller 3D surfaces. By increasing the field of view, this
source of error can be eliminated.
Figure 5.8: The anterior peri-prostatic venous plexus is observed as a dark
region in B-mode and VE images. When its contrast is similar to the
prostate-to-background contrast, it can be mistaken as the gland itself. In
this figure the dashed yellow is the initial automatic contour, and solid red
is the final. Solid blue indicates the manually segmented prostate boundary.
In [60] a 3D segmentation time of 17.17±2.71 s for segmenting the
11 cases using the semi-automatic method was reported (code written in
c ) and an additional 32±14 s manual initialization time was reMATLAB
quired resulting in a total of approximately 49 s. A similar value is seen in
the current automatic method; once the training set is created (which can
be done offline), approximately 45 s is required for segmenting the gland.
However, optimization of the algorithm can further reduce the segmentation
time.
Based on the above results and discussions we conclude that the proposed
3D automatic segmentation algorithm has the potential to be a suitable replacement for manual segmentation of the prostate. We showed the importance of adding vibro-elastography images to the conventional B-mode, for
improving prostate segmentation. The proposed algorithm can be optimized
107
Chapter 5. Automatic Segmentation of the Prostate
to be used as a real-time segmentation method for intra-operative prostate
interventions. -
108
Chapter 6
Visibility of Cancer in
Vibro-elastography Images of
the Prostate
6.1
Introduction
In this chapter, we provide initial results from an ongoing study on cancer
visibility in vibro-elastography images. Vibro-elastography and pathology
data are collected, processed and registered, and manual inspection with
the aim of finding cancer is performed by two observers on the VE images,
the results of which are compared to pathology. In the following sections
the methodology is explained and results are reported. In the Discussions
and Conclusions section results from similar work are reported, and the
limitations of the current work, suggestions for further improvements, and
future pathways are discussed.
6.2
6.2.1
Methods and Materials
Data Collection and Processing
RF longitudinal data frames were obtained from 5 cases using the setup
described in Chapter 4. Pathology results were also available. Computation
of the vibro-elastography images was performed similarly to that of previous
chapters, however, additional changes were made to improve image quality
as follows:
Computation of the vibro-elastography images
The Correlation Coeffiecient (CC) [115] is a measure computed for each
block of tissue, as a function of time, during tissue displacement estimation.
The CC is a measure of how accurate the displacement estimation is, ideally
being 1, and is computed as below:
Chapter 6. Visibility of Cancer in VE Images of the Prostate
CCF1,2
CC1,2 = p
ACF1,1 ACF2,2
(6.1)
where CCF1,2 is the maximum value of the cross-correlation function between signals 1 and 2 and ACF1,1 and ACF2,2 are the maximum values of
the auto-correlation functions. To improve the quality of the images, we
monitor the CC. If the CC falls below a certain threshold set by the user,
the corresponding displacement estimations for that portion of time are not
used for computing the transfer function. If more than a specified number
of frames of data have a CC less than the defined threshold, a ‘gap’ exists
and data from the frames before and after this gap, assuming that they have
an acceptable length, are separately used to compute two transfer functions.
Accordingly, if more gaps exist, multiple transfer functions are computed,
which are ultimately averaged. We set the acceptable length of data to 8
frames and define a gap as 2 or more consecutive unacceptable frames. For
only one unacceptable frame, the corrupt data is replaced by interpolating
between the frames before and after. Currently the CC threshold is set by
trial and error for each case in order to obtain visually improved images (a
value in the range 0.7-0.9). Fig. 6.1 illustrates the concept. This method
has shown to improve the resulting VE images as seen in Fig. 6.2.
The pixel intensity value, V Ejk , at axial location k and lateral location
(line) j of a longitudinal vibro-elastography image used in this chapter is
computed as the weighted average of the magnitude of the transfer function:
V Ejk =
N
1 X
N
n=1
R f2 f1
jk
jk
(j2πf ) Href
CFref
(j2πf ) df
(f2 − f1 )
R f2
f1
jk
(j2πf )df
CFref
(6.2)
jk
where (Href
(j2πf ) is the transfer function computed with respect to a refjk
(j2πf ) is the coherence function, used as the
erence location ref , CFref
weight, and [f1 , f2 ] describes the frequency range of interest, which in this
work coincides with the range of the broad-band vibration applied to tissue,
and N is the number of transfer functions computed from portions of correct
data, after removing low CC frames.
If at location jk no acceptable length of data exists for which the corresponding CC is above the defined threshold, a value of zero is returned and
the stiffness value of that block of tissue is considered unreliable. In order
to distinguish between this unreliable value and a very stiff region (i.e. a TF
value of zero), for each TF image a corresponding binary image of the same
size is produced. In these images pixel values of zero indicate unreliable
110
Chapter 6. Visibility of Cancer in VE Images of the Prostate
Figure 6.1: Illustration of the transfer function improvement by removing
data with low CC. Approximately 80 frames of RF data are collected from
the prostate, for each angle, using the TRUS probe (left). For each region
in the probe’s field of view (red square) the tissue displacement (lower right
plot) is estimated through time, as the tissue is vibrated. Along with the
displacement estimation, the correlation coefficient (CC) is also computed
(upper right plot). Instead of computing the transfer function (TF) using
all the displacement data, a threshold on the CC (in this case 0.85) is set
to remove erroneous displacement data. TFs are computed from portions
of correct data (TF1 and TF2) which are finally averaged to provide the
final TF. The average value of the magnitude of the TF, over the vibration
frequency range, defines the intensity of the corresponding pixel (red square).
data.
Finally, from the set of fan-shaped longitudinal images, 3D VE matrices
or ‘volumes’ were constructed via interpolation, as described in Chapter 4,
Section 4.2.1.
111
Chapter 6. Visibility of Cancer in VE Images of the Prostate
Figure 6.2: The effect of removing data with low correlation (CC). The two
figures show the magnitude of the transfer function from the same prostate
in the sagittal view, before (left) and after (right) removing low CC data.
In this case a threshold of CC=0.85 was used.
Pathology data
In addition to RF data, whole-mount pathology analysis has also been performed on the extracted prostates. The whole-mount sections were produced
using the multi-bladed cutting device described in [24]. Each section is 4 mm
thick from which, in a routine manner, a fine slice was cut and mounted on
a glass slide for hematoxylin and eosin staining. From each prostate, 9-13
slides, depending on the size of the gland, were produced on which cancer
regions were detected and marked by a certified pathologist (see Fig. 6.3).
The slides were scanned and aligned to correct for the in-plane translation
and rotation of each slide with respect to its adjacent slides. To correct
the in-plane translation, bounding boxes were created for each prostate slice
and the centers of these bounding boxes were aligned. Correction of the
in-plane rotation of the slides was performed manually by inspecting the
internal structures of the prostate, including the urethra and the texture of
the central zone. Figure 6.4 shows an example of manual correction of the
orientation of two consecutive pathology slices based on the gland boundary
and internal structures. The aligned pathology images were converted into
a Stradwin file format [109]. Using the Stradwin software we observed and
manually segmented the prostate boundary and cancerous regions to create
prostate and cancer surfaces, an example shown in Fig. 6.5.
112
Chapter 6. Visibility of Cancer in VE Images of the Prostate
Figure 6.3: An example of pathology slides of one case with cancerous regions and grades defined. Starting from the apex, the number defined as
‘level’ corresponds to the depth at which the sample was taken. Cancer is
marked by dotted lines. The left and right side of the prostate is marked as
‘L’ and ‘R’, the anterior and posterior marked as ‘A’ and ‘P’, respectively.
113
Chapter 6. Visibility of Cancer in VE Images of the Prostate
Figure 6.4: An example of manual correction of the orientation in pathology
slides based on the gland internal structure. In this case, guided by the
prostate boundary and internal structures (marked as blue dashed and black
solid lines), the slice on the right will be rotated counter-clockwise to align
with its prior slice on the left.
Figure 6.5: The set of pathology images are converted to a Stradwin file
format. The prostate and the marked cancer regions are then segmented in
this software.
114
Chapter 6. Visibility of Cancer in VE Images of the Prostate
6.2.2
Cancer Visibility in VE Images
To evaluate the visibility of cancer in VE images the VE volumes were observed by two individuals who were blind to the pathology results. The first
observer (the author of this thesis) was more familiar with the appearance
of cancer in VE images, while the second had less experience. Transverse
slices in the VE volume corresponding to the pathology slices were identified
by the registration procedure proposed by Taquee et al. in [104]. In these
slices, all regions suspected for cancer were segmented by the two individuals. We will refer to these identified tumors as the ‘guessed tumors’, with the
marked tumors in pathology being the ‘gold standard’. Suspicious regions
were chosen as regions with low intensity (or dark areas, corresponding to
relatively stiffer tissue) which were continuously seen in multiple transverse
VE images. We assume that dark regions seen in only one or two images
are artefacts.
These suspected tumors were then compared with pathology results in a
region-based manner. For this, the outlined tumors in each of the pathology
slices were labeled as being in one of the four regions: left and right lateral,
anterior and posterior. A guessed tumor in VE which also exists in the same
region and has approximately the same size, is marked as a true positive.
Figure 6.6 shows a few examples of how the comparison was performed.
In this figure, matching tumors are linked with red lines and the false positive
and false negatives are marked by FP and FN, respectively. We considered
only tumors with an area larger than 7 mm2 (3 mm diameter) in both
pathology and VE since the coarse resolution of the VE images and the
presence of noise complicates the detection of small tumors.
In the five studied cases, 25% of the tumors were located in the anterior,
25% in the left lateral, 19% in the right lateral and 31% in the posterior. For
the first observer, the sensitivity of detecting cancer in the correct region
of the VE image corresponding to the pathology slice was 72.5%. 24.4%
of tumors could not be correctly detected in VE (false negatives). 37.3%
of the guessed tumors in VE were in fact not cancerous (false positives).
The sensitivity of VE for detecting cancer with a Gleason score of 7 and
above was 85%. In terms of size, the undetected tumors (considering only
those >7 mm2 in size) had an average area of 24±20 mm2 ( 5.4±5 mm
diameter). It was observed that many of these false negatives were in the
mid-posterior region of the prostate, above the TRUS probe. In this area
we have observed that the VE image quality (i.e. visibility of the prostate
boundary) degrades. Our results showed a higher the detection rate in VE
for higher Gleason scores. Gleason scores for 51 tumors were available. 20
115
Chapter 6. Visibility of Cancer in VE Images of the Prostate
Figure 6.6: Examples of how evaluation of cancer visibility in VE was performed in five pairs of transverse VE (left column) and pathology (right
column) images. Corresponding VE and pathology images were found via
registration and guessed tumors were marked on VE images and compared
to those in pathology. Matching tumors (true positives) are linked with red
lines and the false positive and false negatives are marked by FP and FN,
respectively.
116
Chapter 6. Visibility of Cancer in VE Images of the Prostate
out of the 31 tumors with Gleason scores of 3+3 (64.5%), 13 of the 16 tumors
with Gleason scores of 3+4 (81.25%), both tumors with Gleason scores of
4+3+5 and both tumors with Gleason scores of 4+5 were detected.
However, observer variability in visual detection of cancer is high and the
level of experience of the observer greatly affects cancer detection results.
For the second observer (with less experience with cancer in VE images)
the sensitivity of detecting cancer in the correct region drops to 43.1%. The
percentage of undetected tumors was 56.9% and the percentage of incorrect
guesses was 47.6%. The average area of undetected tumors (again, considering only those >7 mm2 in size) was 63±50 mm2 ( 9±8 mm diameter),
larger than that of the first observer.
Table 6.1 shows the sensitivity, percentage of undetected tumors (false
negatives/number of tumors) and the percentage of incorrect guesses (false
positives/number of guessed tumors) for the four prostate regions of the
registered slices and the whole gland.
Ant.
LL
RL.
Post.
Whole
Sens.%
Obs. 1 Obs. 2
78.3
28.6
77.8
31.3
82.3
68.8
64.3
26.1
75.8
43.1
% undetected
Obs. 1 Obs. 2
21.7
71.4
22.3
68.8
17.6
31.3
35.7
73.9
24.2
56.9
% incorrect guess
Obs. 1 Obs. 2
30.8
16.7
44.0
61.5
50.0
50.0
41.9
33.4
37.3
47.6
Table 6.1: The sensitivity, percentage of undetected tumors (false negatives/number of tumors) and the percentage of incorrect guesses (false positives/number of guessed tumors) for the four regions of the matching slices
and the whole gland.
We noticed that in some cases cancerous tumors which were located in
pathology on the anterior of the prostate were also visible in the VE images.
Such tumors may be missed by TRUS based biopsy. The reason for that is
the depth (and therefore the anterior-posterior position of the biopsy core)
is determined by the length of the core which rarely exceeds 10 mm. Thus
in most cases biopsy only samples the posterior 1/4 to 1/2 of the gland
depending on the anterior-posterior diameter. This is usually adequate,
since anteriorly located tumors are generally restricted to the apex where
the anterior-posterior diameter is small. However, when cancer is located
out of the reach of biopsy it is most likely to be missed and therefore, VE
has the potential to be used in image guided prostate biopsy. Figure 6.7
117
Chapter 6. Visibility of Cancer in VE Images of the Prostate
Figure 6.7: An example of an anteriorly located tumor detected as a dark
region in the corresponding VE image.
shows an example of an anteriorly located tumor detected as a dark region
in the corresponding VE image.
6.3
Discussion and Conclusions
In recent years, prostate cancer detection with the use of elastography has
been of interest to many groups, but is still in an early stage. A general
agreement exists that in the prostate, cancer is stiffer than healthy tissue
(with a stiffness ratio of 1.5-3.5:1). This knowledge has been used as the
basis for detecting cancer. Several groups have performed similar studies on
the visibility of cancer in elastography images.
Salomon et al. [91] investigate the sensitivity and specificity of cancer
detection in the prostate using ultrasound elastography in patients going
through radical prostatectomy. The elastography images they used were
strain images created by manually applying compression and decompression
of the prostate using the transrectal ultrasound probe. An EUB-6500HV
Hitachi ultrasound system with a V53W 7.5 MHz transrectal end-fire probe
was used (Hitachi Medical, Kashiwa, Japan). In a single institution, single
observer study on 109 patients with positive biopsies, B-mode and elastography images were recorded in the transverse plane and areas suspicious for
cancer were marked on the elastography images. Histopathological evaluation was performed on slices 3 mm apart. Sensitivity and specificity was
reported for the 6 regions of the gland (right and left side of the apex, midgland and base) with an average of 75.4% (range 68.1-84.0%) and 76.7%
(range 67.0-84.0%) respectively. They report a positive predictive and neg118
Chapter 6. Visibility of Cancer in VE Images of the Prostate
ative predictive value of 87.8% and 59% and an accuracy of 76%. The
detection rate of cancer foci was best in the apex where approximately 90%
of the tumors were detected and worst at the base where only 75% of the
tumors were detected.
In a paper by Pallwein et al. [82] the prostate of 15 patients with organ
confined cancer confirmed by biopsy and scheduled for radical prostatectomy
were studied. Real-time elastography was performed by two radiologists,
blind to the biopsy results, and the images were interpreted by consensus.
A Voluson 730 (GE Healthcare, Chalfont St Giles, UK) with an endorectal
side fire 7.5MHz transducer was used. Displacement estimates and strain
images were computed from radio frequency data produced from manual
compression and decompression of the transducer. Similar to [91] six outer
gland regions (left and right of base, mid, and apex) were evaluated for cancer. Only hard lesions with a diameter of ≥ 5 mm were considered malignant
in elastography images. 35 tumors of this size were found on histopathology
slides, 32 of which were in the outer gland and 3 in the inner gland. Elastography detected 28 of the 35 tumors and of the 7 undetected foci, 3 were in
the inner gland, and 3 at the base. The best sensitivity for elastography was
at the apex (100%) and the worst at the base (detection of only two of the
five foci), with the mid-gland sensitivity being 94%. Four sites with false
positive findings on elastography were reported. The overall outer gland
sensitivity, specificity and accuracy values of 0.88, 0.93 and 0.91 were obtained. The positive and negative predictive values were, respectively, 0.88
and 0.93.
In a work by Sumura et al. [101] real-time elastography findings were
compared to digital rectal examination, color Doppler ultrasonography, grayscale transrectal ultrasonography and T2-weighted and dynamic contrast
enhanced MRI. Histopathology was used as the gold standard (microslides
spaced at 5 mm intervals) on which 26 tumors with volumes larger than
0.1mL were marked in two regions: anterior (including the transition zone
and/or anterior side of the peripheral zone) and posterior. Strain images
were obtained using an EUB8500 ultrasound system (Hitachi, Japan) with
an end-firing 7.5MHZ transrectal probe by manual compression and decompression of the TRUS probe. They report that the cancer detection rate is
superior to all other methods listed above, with an overall detection rate
of 75% in the anterior and 73.7% in the posterior of the prostate. As expected, the larger the tumor and the higher the Gleason score, the higher
the detection rate.
Tsutsumi et al. [111] collected strain images in the transverse plane
from the base to the apex from 51 patients going through radical prostate119
Chapter 6. Visibility of Cancer in VE Images of the Prostate
ctomy and compared to the corresponding pathological transverse sections
of the prostate. Strain images were obtained by manual compression of a
7.5MHz biplanar probe and an EUB6500 Hitachi ultrasound system (Hitachi, Japan). The detection rate of tumors and relationship of tumor locations in elastography, and relationship between Gleason score and cancer
stage and elastographic findings was observed. For the anterior, middle and
posterior part of the prostate, they report a sensitivity of 94%, 76% and
57% and a positive predictive value of 83%, 72%, and 70%, respectively.
The minimum detectable tumor diameter in their images was 4 mm. Surprisingly, the detection of low grade tumors was more accurate than high
grade tumors (100% for Gleason score 6, 85% for Gleason score 7, 85% for
Gleason score 8 and 63% for Gleason score 9 or 10) and the detection of
T1c tumors was more accurate than T2 or T3. They attribute these results
to the predominance of the peripheral location of high-grade tumors, where
their sensitivity results were low, compared to the anterior location of the
low-grade tumors, where they report better sensitivity results.
In a similar study [65] but on a larger population (311 patients), Miyagawa et al. show that the sensitivity of elastography and elastography+TRUS
imaging in detecting cancer (confirmed by biopsy) is higher than that of
DRE or TRUS only. They report sensitivity of cancer detection in DRE
(37.9%), TRUS (59.0%), elastography (72.6%) and TRUS+elastography
(89.5%). Higher prostate-specific antigen (PSA) level and smaller prostate
volume are reported to increase the sensitivity of elastography and elastography+TRUS. However, the high frequency of false-positive elastography
results and difficulty in the detection of cancer in the peripheral zone are
two main problems reported in their work. They report no significant difference in the sensitivity related to Gleason score. Biopsy findings were used
as the gold standard and regional evaluation was not performed.
In a recent study, Tsutsumi et al. [110] introduce a new excitation
technique called real-time balloon inflation elastography (RBIE) in which
manual inflation and deflation of a balloon covering the TRUS probe by
a piston-type injector is used instead of manual compression of the tissue.
This results in more stable images and reduces operator dependence. They
presented the accuracy and feasibility of this technique in the detection of
prostate cancer on 55 cases in which pathologic specimens were compared
to RBIE images. In this method only 1% of the images were affected by
artifacts due to slippage in the compression plane whereas in manual compression this number was 32%. They report a sensitivity of 84% in the
anterior, 85% in the middle and 60% in the posterior of the prostate and
positive predictive values of 75%, 79%, and 93% for the three regions. A
120
Chapter 6. Visibility of Cancer in VE Images of the Prostate
specificity of 80%, 91% and 96% for the three regions was obtained; an improvment in the posterior region compared to manual compression. They
report that with the RBIE method the detection rate of Gleason scores 7 or
more is larger than cases with a Gleason score of 6 or less.
We showed that it is possible to visually detect cancer in VE images
and our sensitivity values are comparable to some of the above work, but
our results are yet to be improved. In the above studies (except for that of
Tsutsumi et al. [110]) elastography images were obtained by manually exciting the tissue. This increases user dependability, whereas, in our method
excitation is performed automatically. However, in all the aforementioned
work, including ours, detection of cancer is manually performed where an
individual observes the images and decides whether it is cancer or not. Results from this form of detection greatly depend on the experience of the
observer, as we showed. Unfortunately, in the above studies, results from
multiple observers were not reported.
There are many interesting pathways which can be followed from here.
Next, we enlist some current limitations of our work and potential research
paths.
Current limitations and potential research areas
An exact shape and/or area-based comparison between pathology and
suspected cancer regions was not performed due to the various sources of
error that currently exist, some of which are as follows:
Regarding the registration of the VE volume to the pathology volume,
in an ideal situation it is expected that the cancer contours extracted from
pathology images, when overlaid on the corresponding registered VE images,
exactly delineate the cancer as seen in VE. However this is not the case for
various reasons.
Firstly, the registration method is surface-based. Therefore, any prostate
segmentation error will directly affect the registration outcomes. Although
purely images based registration errors are not unaffected by image artifacts
and noise, a combined image- and surface-based method may improve the
results.
The second issue is that there are no obvious and consistent landmarks
seen inside both the VE and pathology volume. Therefore identifying, for
example, which slice in the VE volume corresponds to the first slice of the
pathology is not an easy task. We assumed that the volume of the prostate
or the area of the prostate in each slice does not change significantly. We
segment a portion of the prostate in the VE volume that has the same
superior-inferior length of the prostate in pathology and that the base and
apex areas approximately match.
121
Chapter 6. Visibility of Cancer in VE Images of the Prostate
Thirdly, the registration algorithm used assumes that the TRUS probe
is parallel to the superior-inferior axis of the prostate or, in other words,
the registered transverse pathology slices can not have an orientation where
a pitch (rotation around the left-right axis of the prostate) exists. The
solution provided by Nir et al.[75], considers this and is therefore, deemed
more appropriate for use in the future.
Some error is introduced during the preparation of the pathology slides.
The exact location from which pathology specimens are created is not known.
The ends of the prostate are first removed before creating the whole-mounts
but it is not clear of what size. The specimens are then cut from 4 mm
thick whole-mount slices and placed on slides. Their exact location in the
whole-mount is not known. Placing the specimens of the slides also results
in a difference in the orientation of each specimen with respect to the other.
As a result of these errors, the location of cancer with respect to the prostate
surface from pathology is not exact but only approximate.
Additional study of the false positives is also required. The pathology
results confirmed only tumors with Gleason scores of 6 and above. There
could be a possibility that the stiff regions falsely detected as cancer in VE
could eventually become cancer or were seen as stiff regions due to other
changes in the prostate tissue, such as prostatic intraepithelial neoplasia
or PIN, which is an abnormality of the prostatic gland and believed to
precede the development of prostate adenocarcinoma. Alternatively, the
higher stiffness of the tissue could be due to reasons other than cancer, such
as the presence of calcifications. A pathology analysis more detailed than a
clinical report is therefore essential.
In summary, in this chapter we provided groundwork and preliminary
tools for pursuing cancer detection in VE images. We listed some of the
current limitations of our work and proposed solutions. Our current results
show that cancer can sometimes be visually detected in VE images, however,
as other work also suggest, VE alone can not guarantee a high sensitivity
and specificity. We suggest that the combined use of different ultrasound
image types or modalities can possibly lead to better detection results and
furthermore, lead to the development of an automated method for in-vivo
cancer detection.
122
Chapter 7
Conclusions and Future
Research
In Chapters 2, 3, 4, and 5 novel algorithms and methods were proposed with
the ultimate goal of segmenting the prostate boundary. In Chapter 2 a 3D
semi-automatic prostate segmentation method was presented. The method
was clinically tested with low dose rate brachytherapy as its application in
Chapter 3. Vibro-elastography was shown to improve the visibility of the
prostate in Chapter 4, the results of which were used to devise an automatic
segmentation method in Chapter 5. As a further step, visibility of cancer
was studied in Chapter 6. All proposed methods were tested on in-vivo
patient data.
7.1
Contributions
The contributions of this thesis are summarized as follows:
• Developing a clinically acceptable prostate segmentation tool. This
tool is currently being clinically used as part of the low dose rate
prostate brachytherapy procedure in BC Cancer Agency, Vancouver
Cancer Center. This tool simplifies delineation, reduces user variability and produces smooth and symmetrical prostate contours which
adhere to the brachytherapy guidelines used in this center. By using
this method, individuals other than physicians can create the contours
and physicians are only required to provide corrections.
• Modeling the prostate boundary. A warped and tapered ellispoid was
found to be a suitable shape for modeling the prostate. The benefits
of using a model is less sensitivity of the segmentation algorithm to
image noise. On the other hand, irregular shapes, such as prostates
with extracapsular cancer, can not be properly segmented by such
a pre-defined model. However, for low dose rate brachytherapy, the
Chapter 7. Conclusions and Future Research
application in which our method was used, cancer is usually intracapsular, and smooth and symmetric contours are required, therefore, this
will not cause an issue for the majority of cases.
• Extensive evaluation of the semi-automatic segmentation method. We
would like to emphasize the need for proper evaluation of proposed
segmentation methods. In most segmentation methods seen in the literature, evaluation is limited to computing the accuracy of the method
compared to a gold standard, usually manual contours. Since there
is no standard set of segmented ultrasound images of the prostate,
there is no ground truth to compare various segmentation methods
against. Obviously, a set of high quality prostate images will lead to
good automatic segmentation results. A reported error is meaningful
only when a measure of the quality of the images used for evaluation
or the accepted range for that specific application is also provided.
We provide this acceptability range by measuring the inter- and intraobserver variability in manual segmentation.
• Region-based evaluation. Since segmentation errors in various regions have different consequences with respect to treatment planning,
a region-based evaluation provides more information about the suitability of a segmentation method. Our region-based evaluation, where
error values are reported for nine sectors of the prostate, provides this.
• The feasibility of using vibro-elastography images for visualizing the
prostate boundary was quantitatively studied. Vibro-elastography images were compared to those of B-mode, the commonly used image
type, with MRI as the gold standard. The assessment was performed
both on the boundary visibility and contrast of the gland with respect
to its surrounding tissue. Additionally, a volume-based comparison
between surfaces created from B-mode and VE images and those from
MRI images was performed.
• Developing a new edge measure: the edge continuity. This measure
has been used for evaluating visibility of the prostate boundary in
B-mode and VE images. It has also been used for segmenting the
prostate. We showed that the addition of this measure to our active
shape model-based segmentation algorithm improves edge detection.
124
Chapter 7. Conclusions and Future Research
• Developing an automatic 3D algorithm for segmentation of the prostate
using a combination of B-mode and vibro-elastography images.
7.2
Future work
Novel methods have been presented in this thesis for segmenting and assessment of visualization of the prostate gland. A number of interesting areas
of research can be suggested as follows:
• Extending the semi-automatic segmentation algorithm to be used in
prostate interventions other than LDR brachytherapy. Whether this
algorithm is suitable for use with other prostate treatment methods or
on other modalities is yet to be tested. For example, our preliminary
tests showed that MRI images are too detailed for the IMMPDA edge
detector to detect the correct prostate edge.
• Integrating the prostate segmentation algorithm into an intra-operative
treatment method. Currently in LDR brachytherapy at the Vancouver Cancer Center the prostate is imaged a few days or weeks prior
to the treatment, this being called the volume study. The volume
is segmented, manually or automatically, and based on the obtained
planning target volume, a treatment plan is created. In the operating
room the B-mode images of the gland are manually registered to the
images taken at the volume study and treatment is carried out. The
segmentation algorithm, used intra-operatively, can aid registration.
Once the two image sets are segmented, the results of registration of
the contours can be used to guide the motion of the TRUS probe to
achieve registered pre- and intra-operative images. Furthermore, if
used with an intra-operative treatment planning algorithm the additional pre-operative imaging can be bypassed, avoiding the necessity
to register the images, needless to say the extra costs.
• The automatic method can be used in segmentation of other organs,
e.g. the kidney, since any model can be created and used from a
training set.
• Active shape models only use shape constraints with some information
about the image collected near landmarks (e.g. edge points). Methods
such as active appearance models on both image types can improve
segmentation results as they incorporate other information such as
image texture.
125
Chapter 7. Conclusions and Future Research
• The quality of vibro-elastography images still has room for improvement. Areas where improvements could be made include the data
collection procedure and image processing. The data collection setup
used in this thesis is favorable in the sense that very little modifications to a standard brachytherapy setup are needed. However, a few
flaws exist:
The roll of the rotation motor covers a range of approximately -45
to 45 degrees. Large prostates fall outside this field of view laterally
especially in the mid-region.
The TRUS probe is not a completely rigid object and acts as a cantileverlike structure. Therefore, the desired motion cannot be exactly produced throughout the length of the probe. Direct measurement of the
probe motions via sensors embedded in the probe, or modeling the
motion of the probe during data processing can increase the accuracy
of our results.
Proper placement of the TRUS probe is important for obtaining high
quality images. Applying too much pressure on the prostate leads
to excessive deformation and compression of the gland while too little
pressure will not guarantee adequate propagation of vibration throughout the prostate. Currently the position of the TRUS probe, prior to
collecting RF data, is adjusted by observing the quality of the B-mode
image and is based on the experience of the individual positioning the
probe. Providing visual guides for probe placement, such as with the
help of force sensors placed along the transducer, is suggested.
The resolution of the images that we use is coarse compared to the
resolution of standard medical images such as ultrasound B-mode and
MRI. With the current code, increasing the resolution is simple. However, processing time will increase which is not in favor of our goal of
obtaining a real-time imaging method.
In computing the VE images, 1D motion estimation is performed in
the axial direction. It is clear that lateral motion is also present during
vibration of the tissue. We observed that the quality of VE images
is relatively poor in the posterior region. This is possibly due to the
presence of soft tissue and sliding of the probe on the rectal wall which,
we believe, results in excessive lateral motion in this region.
In the VE images used in this work stiffness is computed ‘relative’ to a
reference location. Based on our qualitative assessment of VE images
obtained using various reference points we chose the average motion
126
of the tissue along a line passing through the focal point, being set
to half the tissue depth. Further investigation on the selection of the
reference may improve image quality.
• The area of cancer detection in elastography images is a relatively
new area of research. Currently a suspicious prostate is screened for
cancer using biopsy. Biopsy is an invasive procedure with a high false
negative rate. Repeats are commonly required. An elastography image
which reliably displays the prostate in a timely manner can potentially
eliminate the need for biopsy. Or, at the least, be used for guided
biopsy to improve detection rates. A detailed list of potential research
projects in this area is presented in Chapter 6.
Bibliography
[1] S. H. Aaltomaa, V. V. Kataja, T. Lahtinen, J. Palmgren, and
T. Forsell. Eight years experience of local prostate cancer treatment
with permanent I125 seed brachytherapy–morbidity and outcome results. Radiotherapy and Oncology, 91(2):213–216, May 2009.
[2] P. Abolmaesumi and M. R. Sirouspour. An interacting multiple model
probabilistic data association filter for cavity boundary extraction
from ultrasound images. IEEE Transactions on Medical Imaging,
23(6):772–784, 2004.
[3] H. Akbari, X. Yang, L. V. Halig, and B. Fei. 3D segmentation of
prostate ultrasound images using wavelet transform. Proceedings of
SPIE, 7962:79622K, 2011.
[4] D. Ash, A. Flynn, J. Battermann, T. de Reijke, P. Lavagnini, L. Blank,
E.S.T.R.A/E.A.U. Urological Brachytherapy Group, and E.O.R.T.C.
Radiotherapy Group. Estro/eau/eortc recommendations on permanent seed implantation for localized prostate cancer. Radiotherapy
and Oncology, 57(3):315–321, Dec 2000.
[5] S. Badiei. Prostate segmentation in ultrasound images using image
warping and ellipse fitting. Master’s thesis, University of British
Columbia, 2007.
[6] S. Badiei, S. E. Salcudean, J. Varah, and W. J. Morris. Prostate
segmentation in 2D ultrasound images using image warping and ellipse
fitting. International Conference on Medical Image Computing and
Computer Assisted Interventions, LNCS 4191:17–24, 2006.
[7] M.S. Bazaraa, H.D. Sherali, and C.M. Shetty. Nonlinear programming:
theory and algorithms. Wiley-Interscience, 2006.
[8] P.J. Besl and H.D. McKay. A method for registration of 3D shapes.
IEEE Transactions on Pattern Analysis and Machine Intelligence,
14(2):239–256, 1992.
[9] M. Bilgen and M. F. Insana. Predicting target detectability on acoustic
elastography. In IEEE Ultrasonics Symposium, pages 1427–1430, 1997.
[10] J. C. Blasko, P. D. Grimm, J. E. Sylvester, K. R. Badiozamani,
D. Hoak, and W. Cavanagh. Palladium-103 brachytherapy for prostate
carcinoma. International Journal of Radiation Oncology Biology
Physics, 46(4):839–850, 2000.
[11] D. Boukerroui, A. Baskurt, J.A. Noble, and O. Basset. Segmentation
of ultrasound images–multiresolution 2D and 3D algorithm based on
global and local statistics. Pattern Recognition Letters, 24(4-5):779–
790, 2003.
[12] Canadian Cancer Society. Canadian cancer statistics. Technical report, 2012.
[13] L. Castéra, J. Vergniol, J. Foucher, B. Le Bail, E. Chanteloup,
M. Haaser, M. Darriet, P. Couzigou, and V. de Lédinghen. Prospective
comparison of transient elastography, fibrotest, apri, and liver biopsy
for the assessment of fibrosis in chronic hepatitis C. Gastroenterology,
128(2):343–350, 2005.
[14] S. Catheline, J. L. Gennisson, and M. Fink. Measurement of elastic
nonlinearity of soft solid with transient elastography. Journal of the
Acoustical Society of America, 114(6 Pt 1):3087–3091, Dec 2003.
[15] Y. J. Choi, J. K. Kim, H. J. Kim, and K. S. Cho. Interobserver
variability of transrectal ultrasound for prostate volume measurement
according to volume and observer experience. American Journal of
Roentgenology, 192(2):444–449, 2009.
[16] D. L. Cochlin, R. H. Ganatra, and D. F. R. Griffiths. Elastography
in the detection of prostatic cancer. Clinical Radiology, 57(11):1014–
1020, 2002.
[17] T. F. Cootes, A. Hill, C. J. Taylor, and J. Haslam. The use of active
shape models for locating structures in medical images. Image and
Vision Computing, 12(6):355–366, 1994.
[18] T. F. Cootes, C. J. Taylor, D. H. Cooper, and J. Graham. Active
shape models-their training and application. Computer Vision and
Image Understanding, 61:38–59, 1995.
[19] J. M. Crook, L. Potters, R. G. Stock, and M. J. Zelefsky. Critical
organ dosimetry in permanent seed prostate brachytherapy: defining
the organs at risk. Brachytherapy, 4(3):186–194, 2005.
[20] V. Daanen, J. Gastaldo, J. Y. Giraud, P. Fourneret, J. L. Descotes,
M. Bolla, D. Collomb, and J. Troccaz. MRI/TRUS data fusion for
brachytherapy. International Journal of Medical Robotics, 2(3):256–
261, 2006.
[21] L. R. Dice. Measures of the amount of ecologic acssociation between
species. Ecology, 26(3):297–302, 1945.
[22] D. A. Dickey and W. A. Fuller. Distribution of the estimators for
autoregressive time series with a unit root. Journal of the American
Statistical Association, 74(366):427–431, 1979.
[23] M. A. Dresner, G. H. Rose, P. J. Rossman, R. Muthupillai, A. Manduca, and R. L. Ehman. Magnetic resonance elastography of skeletal
muscle. Journal of Magnetic Resonance Imaging, 13(2):269–276, Feb
2001.
[24] B. Drew, E. C. Jones, S. Reinsberg, A. C. Yung, S. L. Goldenberg,
and P. Kozlowski. Device for sectioning prostatectomy specimens to
facilitate comparison between histology and in vivo MRI. Journal of
Magnetic Resonance Imaging, 32(4):992–996, Oct 2010.
[25] O. Ecabert, J. Peters, H. Schramm, C. Lorenz, J. von Berg, M. J.
Walker, M. Vembar, M. E. Olszewski, K. Subramanyan, G. Lavi, and
J. Weese. Automatic model-based segmentation of the heart in CT images. IEEE Transactions on Medical Imaging, 27(9):1189–1201, 2008.
[26] H. Eskandari, S. E. Salcudean, and R. Rohling. Tissue strain imaging using a wavelet transform-based peak search algorithm. IEEE
Transactions on Ultrasonics, Ferroelectrics and Frequency Control,
54(6):1118–1130, 2007.
[27] A. Fenster, D. B. Downey, and H. N. Cardinal. Three-dimensional
ultrasound imaging. Physics in Medicine and Biology, 46(5):R67–R99,
May 2001.
[28] J. Ferlay, D. M. Parkin, and E. Steliarova-Foucher. Estimates of cancer incidence and mortality in Europe in 2008. European Journal of
Cancer, 46(4):765–781, 2010.
[29] I. N. Fleming, H. Rivaz, K. Macura, L. M. Su, U. Hamper, G. A.
Lagoda, A. L. Burnett II, T. Lotan, R. H. Taylor, G. D. Hager, and
E. M. Boctor. Ultrasound elastography - enabling technology for image guided laparoscopic prostatectomy. In Proceedings of the SPIE
Medical Imaging 2009: Visualization, Image-Guided Procedures, and
Modeling., volume 7261, page 72612I, 2009.
[30] L. Gao, K. J. Parker, R. M. Lerner, and S. F. Levinson. Imaging of
the elastic properties of tissue–a review. Ultrasound in Medicine and
Biology, 22(8):959–977, 1996.
[31] B. S Garra. Elastography: current status, future prospects, and making it work for you. Ultrasound Q, 27(3):177–186, Sep 2011.
[32] B. S. Garra, E. I. Cespedes, J. Ophir, S. R. Spratt, R. A. Zuurbier,
C. M. Magnant, and M. F. Pennanen. Elastography of breast lesions:
initial clinical results. Radiology, 202(1):79–86, Jan 1997.
[33] L. Gong, S. D. Pathak, D. R. Haynor, P. S. Cho, and Y. Kim. Parametric shape modeling using deformable superellipses for prostate segmentation. IEEE Transactions on Medical Imaging, 23(3):340–349,
2004.
[34] S. Gravas, C. Mamoulakis, J. Rioja, V. Tzortzis, T. de Reijke, H. Wijkstra, and J. de la Rosette. Advances in ultrasound technology in
oncologic urology. Urologic Clinics of North America, 36(2):133–45,
vii, May 2009.
[35] J. F. Greenleaf, M. Fatemi, and M. Insana. Selected methods for imaging elastic properties of biological tissues. Annual Review of Biomedical Engineering, 5:57–78, 2003.
[36] P. D. Grimm, J. C. Blasko, J. E. Sylvester, R. M. Meier, and W. Cavanagh. 10-year biochemical (prostate-specific antigen) control of
prostate cancer with (125)I brachytherapy. International Journal of
Radiation Oncology Biology Physics, 51(1):31–40, 2001.
[37] T. J. Hall, Y. Zhu, and C. S. Spalding. In vivo real-time freehand
palpation imaging. Ultrasound in Medicine and Biology, 29(3):427–
435, Mar 2003.
[38] E. J. Halpern, D. L. Cochlin, and B. B. Goldberg. Imaging of the
prostate. Martin Dunitz, 2002.
[39] A. C. Hodge, A. Fenster, D. B. Downey, and H. M. Ladak. Prostate
boundary segmentation from ultrasound images using 2D active shape
models: optimisation and extension to 3D. Computer Methods and
Programs in Biomedicine, 84(2-3):99–113, 2006.
[40] K. Hoyt, K. J. Parker, and D. J. Rubens. Real-time shear velocity
imaging using sonoelastographic techniques. Ultrasound in Medicine
and Biology, 33(7):1086–1097, Jul 2007.
[41] N. Hu, D. B. Downey, A. Fenster, and H. M. Ladak. Prostate
boundary segmentation from 3D ultrasound images. Medical Physics,
30(7):1648–1659, Jul 2003.
[42] X. Huang and D. N. Metaxas. Metamorphs: deformable shape and
appearance models. IEEE Transactions on Pattern Analysis and Machine Intelligence, 30(8):1444–1459, 2008.
[43] L. Huwart, C. Sempoux, E. Vicaut, N. Salameh, L. Annet, E. Danse,
F. Peeters, L. C. ter Beek, J. Rahier, R. Sinkus, Y. Horsmans, and
B. E. Van Beers. Magnetic resonance elastography for the noninvasive
staging of liver fibrosis. Gastroenterology, 135(1):32–40, 2008.
[44] A. Di Ieva, F. Grizzi, E. Rognone, Z. Tsz Ho Tse, T. Parittotokkaporn,
F. Rodriguez Y Baena, M. Tschabitscher, C. Matula, S. Trattnig, and
R. Rodriguez Y Baena. Magnetic resonance elastography: a general
overview of its current and future applications in brain imaging. Neurosurgical Review, 33(2):137–45; discussion 145, Apr 2010.
[45] A. Jendoubi, J. Zeng, and M. F. Chouikha. Top-down approach to
segmentation of prostate boundaries in ultrasound images. In Applied Imagery Pattern Recognition Workshop, 2004. Proceedings. 33rd,
pages 145–149, 2004.
[46] N. N. Kachouie, P. Fieguth, and S. Rahnamayan. An elliptical level
set method for automatic TRUS prostate image segmentation. In Proceedings of the IEEE International Symposium on Signal Processing
and Information Technology, pages 191–196, Aug. 2006.
[47] M. J. Kadour and J. A. Noble. Assisted-freehand ultrasound elasticity imaging. IEEE Transactions on Ultrasonics, Ferroelectrics and
Frequency Control, 56(1):36–43, 2009.
[48] K. Kamoi, K. Okihara, A. Ochiai, O. Ukimura, Y. Mizutani,
A. Kawauchi, and T. Miki. The utility of transrectal real-time elastography in the diagnosis of prostate cancer. Ultrasound in Medicine
and Biology, 34(7):1025–1032, Jul 2008.
[49] L. Kanzler. A study of the efficiency of the foreign exchange market through analysis of ultra-high frequency data. PhD thesis, Oxford
University, Oxford, UK, 1998.
[50] L. Kitchen and A. Rosenfeld. Edge evaluation using local edge
coherence. IEEE Transactions on Systems, Man and Cybernetics,
11(9):597–605, 1981.
[51] E. Konofagou and J. Ophir. A new elastographic method for estimation and imaging of lateral displacements, lateral strains, corrected
axial strains and Poisson’s ratios in tissues. Ultrasound in Medicine
and Biology, 24(8):1183–1199, Oct 1998.
[52] E. E. Konofagou, J. D’hooge, and J. Ophir. Myocardial elastography–
a feasibility study in vivo. Ultrasound in Medicine and Biology,
28(4):475–482, Apr 2002.
[53] E. E. Konofagou and J. Ophir. Precision estimation and imaging of
normal and shear components of the 3D strain tensor in elastography.
Physics in Medicine and Biology, 45(6):1553–1563, Jun 2000.
[54] P. D. Kovesi. MATLAB and Octave functions for computer vision and
image processing. Centre for Exploration Targeting, School of Earth
and Environment, The University of Western Australia. Available
from: <http://www.csse.uwa.edu.au/∼pk/Research/research.html>.
[55] H. M. Ladak, F. Mao, Y. Wang, D. B. Downey, D. A. Steinman,
and A. Fenster. Prostate boundary segmentation from 2D ultrasound
images. Medical Physics, 27(8):1777–1788, Aug 2000.
[56] R.M. Lerner, K.J. Parker, J. Holen, R. Gramiak, and R.C. Waag.
Sono-elasticity: Medical elasticity images derived from ultrasound signals in mechanically vibrated targets. Acoustic Imaging, 16:317–327,
1988.
[57] J. Li, Y. Cui, R. E. English, and J. A. Noble. Ultrasound estimation
of breast tissue biomechanical properties using a similarity-based nonlinear optimization approach. Journal of Strain Analysis, 44(5):363–
374, 2009.
[58] R. G. P. Lopata, M. M. Nillesen, H. H. G. Hansen, I. H. Gerrits, J. M.
Thijssen, and C. L. de Korte. Performance evaluation of methods for
two-dimensional displacement and strain estimation using ultrasound
radio frequency data. Ultrasound in Medicine and Biology, 35(5):796–
812, May 2009.
[59] A. Lyshchik, T. Higashi, R. Asato, S. Tanaka, J. Ito, J. J. Mai,
C. Pellot-Barakat, M. F. Insana, A. B. Brill, T. Saga, M. Hiraoka,
and K. Togashi. Thyroid gland tumor diagnosis at us elastography.
Radiology, 237(1):202–211, Oct 2005.
[60] S. S. Mahdavi, N. Chng, I. Spadinger, W. J. Morris, and S. E. Salcudean. Semi-automatic segmentation for prostate interventions. Medical
Image Analysis, 15(2):226–237, Apr 2011.
[61] S. S. Mahdavi, M. Moradi, X. Wen, W. J. Morris, and S. E. Salcudean.
Vibro-elastography for visualization of the prostate region. In Medical
Image Computing and Computer Assisted Intervention, volume 5762
of LNCS, pages 339–347, 2009.
[62] S. S. Mahdavi, M. Moradi, X. Wen, W. J. Morris, and S. E. Salcudean.
Evaluation of visualization of the prostate gland in vibro-elastography
images. Medical Image Analysis, 15(4):589–600, Aug 2011.
[63] S. S. Mahdavi, I. Spadinger, N. Chng, S. E. Salcudean, and W. J. Morris. Semiautomatic segmentation for prostate brachytherapy: Dosimetric evaluation. Brachytherapy, 2011.
[64] J. E. McNeal, E. A. Redwine, F. S. Freiha, and T. A. Stamey. Zonal
distribution of prostatic adenocarcinoma. Correlation with histologic
pattern and direction of spread. American Journal of Surgical Pathology, 12(12):897–906, 1988.
[65] T. Miyagawa, M. Tsutsumi, T. Matsumura, N. Kawazoe, S. Ishikawa,
T. Shimokama, N. Miyanaga, and H. Akaza. Real-time elastography
for the diagnosis of prostate cancer: evaluation of elastographic moving
images. Japanese Journal of Clinical Oncology, 39(6):394–398, Jun
2009.
[66] W. J. Morris, M. Keyes, D. Palma, I. Spadinger, M. R. McKenzie,
A. Agranovich, T. Pickles, M. Liu, W. Kwan, J. Wu, E. Berthelet, and
H. Pai. Population-based study of biochemical and survival outcomes
after permanent 125I brachytherapy for low- and intermediate-risk
prostate cancer. Urology, 73(4):860–5; discussion 865–7, Apr 2009.
[67] R. Muthupillai, D. J. Lomas, P. J. Rossman, J. F. Greenleaf, A. Manduca, and R. L. Ehman. Magnetic resonance elastography by direct visualization of propagating acoustic strain waves. Science,
269(5232):1854–1857, Sep 1995.
[68] R. P. Myers. Structure of the adult prostate from a clinician’s standpoint. Clinical Anatomy, 13(3):214–215, 2000.
[69] S. Nag, D. Beyer, J. Friedland, P. Grimm, and R. Nath. American
brachytherapy society (ABS) recommendations for transperineal permanent brachytherapy of prostate cancer. International Journal of
Radiation Oncology Biology Physics, 44(4):789–799, 1999.
[70] S. Nag, W. Bice, K. DeWyngaert, B. Prestidge, R. Stock, and Y. Yu.
The american brachytherapy society recommendations for permanent
prostate brachytherapy postimplant dosimetric analysis. International
Journal of Radiation Oncology Biology Physics, 46(1):221–230, 2000.
[71] S. Nag, J. P. Ciezki, R. Cormack, S. Doggett, K. DeWyngaert, G. K.
Edmundson, R. G. Stock, N. N. Stone, Y. Yu, M. J. Zelefsky, and
American Brachytherapy Society Clinical Research Committe. Intraoperative planning and evaluation of permanent prostate brachytherapy: report of the american brachytherapy society. International Journal of Radiation Oncology Biology Physics, 51(5):1422–1430, 2001.
[72] N. D. Nanayakkara, J. Samarabandu, and A. Fenster. Prostate
segmentation by feature enhancement using domain knowledge and
adaptive region based operations. Physics in Medicine and Biology,
51(7):1831–1848, Apr 2006.
[73] J. C. Nascimento and J. S. Marques. Robust shape tracking with
multiple models in ultrasound images. IEEE Transactions on Image
Processing, 17(3):392–406, 2008.
[74] National Cancer Institute, 2012.
[75] G. Nir and S. E. Salcudean. Registration of whole-mount histology
and tomography of the prostate using particle filtering. In Proceedings
of SPIE Medical Imaging, volume 8676, 2013.
[76] J. A. Noble and d. Boukerroui. Ultrasound image segmentation: a
survey. IEEE Transactions on Medical Imaging, 25(8):987–1010, .
2006.
[77] M. O’Donnell, A.R. Skovoroda, B.M. Shapo, and S.Y. Emelianov. Internal displacement and strain imaging using ultrasonic speckle tracking. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency
Control, 41(3):314–325, 1994.
[78] J. Ophir, I. Céspedes, B. Garra, H. Ponnekanti, Y. Huang, and
N. Maklad. Elastography: ultrasonic imaging of tissue strain and
elastic modulus in vivo. European Journal of Ultrasound, 3(1):49–70,
1996.
[79] J. Ophir, I. Céspedes, H. Ponnekanti, Y. Yazdi, and X. Li. Elastography: a quantitative method for imaging the elasticity of biological
tissues. Ultrasonic Imaging, 13(2):111–134, Apr 1991.
[80] J. Ophir, B. Garra, F. Kallel, E. Konofagou, T. Krouskop, R. Righetti,
and T. Varghese. Elastographic imaging. Ultrasound in Medicine and
Biology, 26 Suppl 1:S23–S29, May 2000.
[81] K. Otto, B. G. Clark, and C. Huntzinger. Exploring the limits of spatial resolution in radiation dose delivery. Medical Physics, 29(8):1823–
1831, Aug 2002.
[82] L. Pallwein, M. Mitterberger, P. Struve, G. Pinggera, W. Horninger,
G. Bartsch, F. Aigner, A. Lorenz, F. Pedross, and F. Frauscher. Realtime elastography for detecting prostate cancer: preliminary experience. British Journal of Urology International, 100(1):42–46, 2007.
[83] K. J. Parker, M. M. Doyley, and D. J. Rubens. Imaging the elastic
properties of tissue: the 20 year perspective. Physics in Medicine and
Biology, 56(1):R1–R29, Jan 2011.
[84] K. J. Parker, L. S. Taylor, S. Gracewski, and D. J. Rubens. A unified view of imaging the elastic properties of tissue. Journal of the
Acoustical Society of America, 117(5):2705–2712, May 2005.
[85] S. D. Pathak, V. Chalana, D. R. Haynor, and Y. Kim. Edge-guided
boundary delineation in prostate ultrasound images. IEEE Transactions on Medical Imaging, 19(12):1211–1219, 2000.
[86] M. A. Penna, K. A. Dines, R. Seip, R. F. Carlson, and N. T. Sanghvi.
Modeling prostate anatomy from multiple view TRUS images for
image-guided HIFU therapy. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control, 54(1):52–69, 2007.
[87] C. Pluempitiwiriyawej, J. M. F. Moura, Y. J. L. Wu, and C. Ho.
STACS: new active contour scheme for cardiac MR image segmentation. IEEE Transactions on Medical Imaging, 24(5):593–603, 2005.
[88] H. Rivaz, E. Boctor, P. Foroughi, R. Zellars, G. Fichtinger, and
G. Hager. Ultrasound elastography: a dynamic programming approach. IEEE Transactions on Medical Imaging, 27(10):1373–1377,
2008.
[89] F. Sahba, H. R. Tizhoosh, and M. M. Salama. A coarse-to-fine
approach to prostate boundary segmentation in ultrasound images.
Biomedical Engineering Online, 4(58):58, 2005.
[90] S. E. Salcudean, D. French, S. Bachmann, R. Zahiri-Azar, X. Wen,
and W. J. Morris. Viscoelasticity modeling of the prostate region using
vibro-elastography. International Conference on Medical Image Computing and Computer Assisted Interventions, 9(Pt 1):389–396, 2006.
[91] G. Salomon, J. Köllerman, I. Thederan, F. K. H. Chun, L. Budus,
T. Schlomm, H. Isbarn, H. Heinzer, H. Huland, and M. Graefen.
Evaluation of prostate cancer detection with ultrasound real-time elastography: a comparison with step section pathological analysis after
radical prostatectomy. European Urology, 54(6):1354–1362, 2008.
[92] C. B. Saw and N. Suntharalingam. Quantitative assessment of interstitial implants. International Journal of Radiation Oncology Biology
Physics, 20(1):135–139, Jan 1991.
[93] J. A. Schaar, C. L. de Korte, F. Mastik, C. Strijder, G. Pasterkamp,
E. Boersma, P. W. Serruys, and A. F. W. van Der Steen. Characterizing vulnerable plaque features with intravascular elastography.
Circulation, 108(21):2636–2641, 2003.
[94] F. K. W. Schaefer, I. Heer, P. J. Schaefer, C. Mundhenke, S. Osterholz, B. M. Order, N. Hofheinz, J. Hedderich, M. Heller, W. Jonat,
and I. Schreer. Breast ultrasound elastography–results of 193 breast lesions in a prospective study with histopathologic correlation. Eurpean
Journal of Radiology, 77(3):450–456, 2011.
[95] D. Shen, Y. Zhan, and C. Davatzikos. Segmentation of prostate boundaries from ultrasound images using statistical shape model. IEEE
Transactions on Medical Imaging, 22(4):539–551, 2003.
[96] S. Sidhu, W. J. Morris, I. Spadinger, M. Keyes, M. Hilts, R. Harrison,
K. Otto, M. McKenzie, and A. Agranovich. Prostate brachytherapy
postimplant dosimetry: a comparison of prostate quadrants. International Journal of Radiation Oncology Biology Physics, 52(2):544–552,
Feb 2002.
[97] R. Sinkus, J. Lorenzen, D. Schrader, M. Lorenzen, M. Dargatz, and
D. Holz. High-resolution tensor MR elastography for breast tumour
detection. Physics in Medicine and Biology, 45(6):1649–1664, Jun
2000.
[98] W. L. Smith, C. Lewis, G. Bauman, G. Rodrigues, D. D’Souza, R. Ash,
D. Ho, V. Venkatesan, D. Downey, and A. Fenster. Prostate volume
contouring: a 3D analysis of segmentation using 3DTRUS, CT, and
MR. International Journal of Radiation Oncology Biology Physics,
67(4):1238–1247, Mar 2007.
[99] R. Souchon, V. Hervieu, A. Gelet J. Ophir, and J. Y. Chapelon. Human prostate elastography: in vitro study. In IEEE Symposium on
Ultrasonics, volume 2, pages 1251–1253, Oct. 2003.
[100] S. Srinivasan and J. Ophir. A zero-crossing strain estimator for elastography. Ultrasound in Medicine and Biology, 29(2):227–238, Feb
2003.
[101] M. Sumura, K. Shigeno, T. Hyuga, T. Yoneda, H. Shiina, and
M. Igawa. Initial evaluation of prostate cancer with real-time elastography based on step-section pathologic analysis after radical prostatectomy: a preliminary study. International Journal of Urology,
14(9):811–816, Sep 2007.
[102] J. E. Sylvester, P. D. Grimm, S. M. Eulau, R. K. Takamiya, and
D. Naidoo. Permanent prostate brachytherapy preplanned technique:
The modern Seattle method step-by-step and dosimetric outcomes.
Brachytherapy, 8(2):197–206, 2009.
[103] M. Tanter, J. Bercoff, L. Sandrin, and M. Fink. Ultrafast compound
imaging for 2-D motion vector estimation: application to transient
elastography. IEEE Transactions on Ultrasonics, Ferroelectrics and
Frequency Control, 49(10):1363–1374, 2002.
[104] F. Taquee, O. Goksel, S. S. Mahdavi, M. Keyes, W. J. Morris,
I. Spadinger, and S. Salcudean. Deformable prostate registration from
MR and TRUS images using surface error driven FEM models. In
Proceedings of the SPIE 8316, 831612, 2012.
[105] P. Thevenaz and M. Unser. Snakuscules. IEEE Transactions on Image
Processing, 17(4):585–593, 2008.
[106] C. W. Thomas, A. Kruk, C. E. McGahan, I. Spadinger, and W. J.
Morris. Prostate brachytherapy post-implant dosimetry: a comparison
between higher and lower source density. Radiotherapy and Oncology,
83(1):18–24, Apr 2007.
[107] I. Thompson, J. B. Thrasher, G. Aus, A. L. Burnett, E. D. CanbyHagino, M. S. Cookson, A. V D’Amico, R. R. Dmochowski, D. T. Eton,
J. D. Forman, S. L. Goldenberg, J. Hernandez, C. S. Higano, S. R.
Kraus, J. W. Moul, C. M. Tangen, and A. U. A. Prostate Cancer
Clinical Guideline Update Panel. Guideline for the management of
clinically localized prostate cancer: 2007 update. Journal of Urology,
177(6):2106–2131, Jun 2007.
[108] S. Tong, H. N. Cardinal, R. F. McLoughlin, D. B. Downey, and A. Fenster. Intra- and inter-observer variability and reliability of prostate
volume measurement via two-dimensional and three-dimensional ultrasound imaging. Ultrasound in Medicine and Biology, 24(5):673–681,
Jun 1998.
[109] G. M. Treece, R. W. Prager, A. H. Gee, and L. Berman. Surface
interpolation from sparse cross sections using region correspondence.
IEEE Transactions on Medical Imaging, 19(11):1106–1114, 2000.
[110] M. Tsutsumi, T. Miyagawa, T. Matsumura, T. Endo, S. Kandori,
T. Shimokama, and S. Ishikawa. Real-time balloon inflation elastography for prostate cancer detection and initial evaluation of clinicopathologic analysis. American Journal of Roentgenology, 194(6):W471–
W476, 2010.
[111] M. Tsutsumi, T. Miyagawa, T. Matsumura, N. Kawazoe, S. Ishikawa,
T. Shimokama, T. Shiina, N. Miyanaga, and H. Akaza. The impact
of real-time tissue elasticity imaging (elastography) on the detection
of prostate cancer: clinicopathological analysis. International Journal
of Clinical Oncology, 12(4):250–255, 2007.
[112] E. Turgay, S. Salcudean, and R. Rohling. Identifying the mechanical properties of tissue by ultrasound strain imaging. Ultrasound in
Medicine and Biology, 32(2):221–235, Feb 2006.
[113] I. B. Tutar, S. D. Pathak, L. Gong, P. S. Cho, K. Wallner, and
Y. Kim. Semiautomatic 3-D prostate segmentation from TRUS images
using spherical harmonics. IEEE Transactions on Medical Imaging,
25(12):1645–1654, 2006.
[114] J. M. Varah. Least squares data fitting with implicit functions. BIT,
36(4):842–854, 1996.
[115] T. Varghese and J. Ophir. Estimating tissue strain from signal decorrelation using the correlation coefficient. Ultrasound in Medicine and
Biology, 22(9):1249–1254, 1996.
[116] T. Wu, J. P. Felmlee, J. F. Greenleaf, S. J. Riederer, and R. L. Ehman.
Assessment of thermal tissue ablation with MR elastography. Magnetic
Resonance in Medicine, 45(1):80–87, Jan 2001.
[117] P. Yan, S. Xu, B. Turkbey, and J. Kruecker. Adaptively learning
local shape statistics for prostate segmentation in ultrasound. IEEE
Transactions on Biomedical Engineering, 58(3):633–641, 2011.
[118] X. Yang and B. Fei. 3D prostate segmentation of ultrasound images
combining longitudinal image registration and machine learning. In
SPIE, volume 8316, page 83162O, 2012.
[119] A. Yezzi, S. Kichenassamy, A. Kumar, P. Olver, and A. Tannenbaum.
A geometric snake model for segmentation of medical imagery. IEEE
Transactions on Medical Imaging, 16(2):199–209, 1997.
[120] M. Yin, J. A. Talwalkar, K. J. Glaser, A. Manduca, R. C. Grimm, P. J.
Rossman, J. L. Fidler, and R. L. Ehman. Assessment of hepatic fibrosis
with magnetic resonance elastography. Clinical Gastroenterology and
Hepatology, 5(10):1207–1213.e2, 2007.
[121] R. Zahiri-Azar, O. Goksel, and S. E. Salcudean. Sub-sample displacement estimation from digitized ultrasound RF signals using multidimensional polynomial fitting of the cross-correlation function. IEEE
Transactions on Ultrasonics, Ferroelectrics and Frequency Control,
57(11):2403–2420, 2010.
[122] R. Zahiri-Azar and S. E. Salcudean. Motion estimation in ultrasound
images using time domain cross correlation with prior estimates. IEEE
Transactions on Biomedical Engineering, 53(10):1990–2000, 2006.
[123] A. Zaim. An edge-based approach for segmentation of prostate ultrasound images using phase symmetry. In ISCCSP 2008, pages 10–13,
2008.
[124] M. J. Zelefsky, D. A. Kuban, L. B. Levy, L. Potters, D. C. Beyer,
J. C. Blasko, B. J. Moran, J. P. Ciezki, A. L. Zietman, T. M. Pisansky, M. Elshaikh, and E. M. Horwitz. Multi-institutional analysis
of long-term outcome for stages T1-T2 prostate cancer treated with
permanent seed implantation. International Journal of Radiation Oncology Biology Physics, 67(2):327–333, Feb 2007.
[125] Y. Zhan and D. Shen. Deformable segmentation of 3-D ultrasound
prostate images using statistical texture matching method. IEEE
Transactions on Medical Imaging, 25(3):256–272, 2006.
[126] Y. Zhan, D. Shen, J. Zeng, L. Sun, G. Fichtinger, J. Moul, and C. Davatzikos. Targeted prostate biopsy using statistical image analysis.
IEEE Transactions on Medical Imaging, 26(6):779–788, 2007.
[127] M. Zhang, P. Nigwekar, B. Castaneda, K. Hoyt, J. V. Joseph,
A. di Sant’Agnese, E. M. Messing, J. G. Strang, D. J. Rubens, and
K. J. Parker. Quantitative characterization of viscoelastic properties
of human prostate correlated with histology. Ultrasound in Medicine
and Biology, 34(7):1033–1042, Jul 2008.