Level Set Segmentation of Retinal OCT Images
Bashir Isa Dodo
1
, Yongmin Li
1
, XiaoHui Liu
1
and Muhammad Isa Dodo
2
1
Department of Computer Science, Brunel University, U.K.
2
Katsina State Institute of Technology and Management, Katsina, Nigeria
Keywords:
Image Segmentation, Level Set, Evolution Constrained Optimisation, Optical Coherence Tomography.
Abstract:
Optical coherence tomography (OCT) yields high-resolution images of the retina. Reliable identification of
the retinal layers is necessary for the extraction of clinically useful information used for tracking the progress
of medication and diagnosis of various ocular diseases. Many automatic methods have been proposed to aid
with the analysis of retinal layers, mainly, due to the complexity of retinal structures, the cumbersomeness of
manual segmentation and variation from one specialist to the other. However, a common drawback suffered by
existing methods is the challenge of dealing with image artefacts and inhomogeneity in pathological structures.
In this paper, we embed prior knowledge of the retinal architecture derived from the gradient information, into
the level set method to segment seven (7) layers of the retina. Mainly, we start by establishing the region
of interest (ROI).The gradient edges obtained from the ROI are used to initialise curves for the layers, and
the layer topology is used in constraining the evolution process towards the actual layer boundaries based on
image forces. Experimental results show our method obtains curves that are close to the manual layers labelled
by experts.
1 INTRODUCTION
Optical coherence tomography first introduced by
(Huang et al., 1991) is a noninvasive imaging tech-
nique that provides cross-sectional images of the
retina with an acquisition speed of approximately
25,000 A-scans per second, and an axial image reso-
lution of approximately 5-7µm (Raftopoulos and Trip,
2012; Jaffe, 2012; Adhi and Duker, 2013). In cur-
rent ophthalmology, identifying various layers of the
retina on OCT has become a vital tool for diagnos-
ing and tracking the progress of medication of various
visual impairments (Adhi and Duker, 2013). Since
the manual segmentation is not only tedious but also
subjected to intra- and inter-grader variance (Yazdan-
panah et al., 2011), many automatic methods have
been proposed to assist with the segmentation pro-
cess(Sun et al., 2016). Current research in retinal
OCT segmentation focuses on improving various as-
pects ranging from computational time, the number of
layers segmented, use of prior knowledge and compu-
tational complexity to mention a few. In general, seg-
mentation is partitioning based on some image char-
acteristics. How these characteristics are defined de-
termines the computation burden of the algorithm.
In some cases, this computational burden is reduced
through dynamic programming (Chiu et al., 2010) or
topology modification (Liu et al., 2018). The number
of questions (conditions) an algorithm has to check or
satisfy is usually the most crucial factor.
Particularly, Markov Boundary Model
(Koozekanani et al., 2001), later extended by
(Boyer et al., 2006), geodesic distance (Duan et al.,
2017), level sets (Wang et al., 2015; Novosel et al.,
2013), graph-based methods (Chiu et al., 2010;
Garvin, 2008; Dodo et al., 2017; Dodo et al., 2018),
and machine learning (Lang et al., 2013) have been
used in the segmentation of retinal OCT images.
Although the level sets method has automatic topo-
logical handling, the steps can be computationally
expensive (Shi and Karl, 2005), while adding
complex constraints in the segmentation method
usually increases the complexity of an algorithm.
In this work, we incorporate simple yet efficient
topological constraints to the evolution process
of the level set method, to improve accuracy and
reduce computational complexity for OCT image
segmentation.
The method we propose is based on the fol-
lowing considerations: 1. The image gradients are
used to initialise curves in order to handle under-
segmentation and over-segmentation of the image; 2.
Dodo, B., Li, Y., Liu, X. and Dodo, M.
Level Set Segmentation of Retinal OCT Images.
DOI: 10.5220/0007577600490056
In Proceedings of the 12th International Joint Conference on Biomedical Engineering Systems and Technologies (BIOSTEC 2019), pages 49-56
ISBN: 978-989-758-353-7
Copyright
c
2019 by SCITEPRESS Science and Technology Publications, Lda. All rights reserved
49
The evolution of a curve is explicitly based on layer
arrangements and implicitly based on OCT topology.
This means for each image the initial contours are
specific to the image under investigation, while the
forces in the normal direction and the topology con-
strains guide the contours towards layer boundaries.
Our method segments an OCT image into 7 seg-
ments, relating to: Nerve Fibre Layer (NFL); Gan-
glion Cell Layer + Inner Plexiform Layer + In-
ner Nuclear Layer (GCL+IPL+INL); Outer Plexiform
Layer (OPL); Outer Nuclear Layer to Inner Segment
(ONL+IS), Outer Segment (OS) and Retinal Pigment
Epithelium (RPE). Locations of these layers on an
OCT image are shown in Figure 1. The rest of the
paper is organised as follows. Section 2 discusses the
proposed method in details. Experimental results and
discussions are treated in Section 3. Finally, conclu-
sions are drawn in Section 4.
2 PROPOSED METHOD
This section describes our approach of segmenting
retinal OCT images. A schematic representation of
the method is illustrated in Figure 2, and details of
each step are elaborated in the subsequent subsec-
tions.
2.1 Pre-processing
The pre-processing steps are illustrated in Figure 3.
Each OCT B-scan image I is first enhanced with a
Gaussian filter to reduce the image noise. The lay-
ers targeted in our study lies within the total retinal
thickness (TRT), which starts from the Internal Limit-
ing Membrane(ILM) to the posterior boundary of the
Retinal pigment epithelium (RPE), i.e. the boundary
between the retinal nerve fibre layer and the vitreous,
and the boundary between the RPE and the choroid
regions respectively. It is commonly accepted that the
NFL, IS-OS and RPE exhibits high reflectivity in an
OCT image (Chiu et al., 2010; Lu et al., 2011; Tian
et al., 2015), based on experiments the ILM and RPE
exibits the highest transitions from dark-bright and
bright-dark, respectively (Dodo et al., 2018). Based
on this understanding of the retinal structure, the ILM
and RPE are identified using the shortest path (Dijk-
stra, 1959), by searching for the highest transitions on
two separate adjacency matrices(Chiu et al., 2010).
Using the ILM and RPE points the image is
cropped to I
cropped
, such that only the Region of in-
terest (ROI) with useful layer information is remain-
ing. This helps in dealing with layer-like structures
outside the ROI and the computational cost associ-
ated with handling image background in segmenta-
tion. The pre-processing is vital in our segmentation
process because only the actual layer properties im-
pact the evolution of the curve. The next operation on
the image is to generate a mask I
mask
of the cropped
image, and then multiply it by the original image I.
The examples of resulting images from this step are
shown in figure 3, column 4, and expressed by the
equation below:
I
processed
= I
mask
I (1)
The process in this subsection is essential because
only the layer structures are obtained when the gradi-
ent of the image is acquired. One of the significant
roles of the pre-processing is to eliminate the need for
handling background as depicted in figure 4 and fur-
ther discussed in the next few subsections 2.2 and 2.3.
From Figure 4 (a) the background and image noise
will affect the segmentation processes (because the
area highlighted in red will be initialised), and for this
reason we establish the ROI only based on part of the
image that contains useful information Figure 4 (b).
The estabilishment of the ROI also complements the
thresholding and refinement processes in the layer ini-
tialisation stage. We use the size of the cropped image
to reduce computational time further and to eliminate
the need for storing idle points.
2.2 Boundary Initialisation
To initialise contours we obtain the vertical gradient
I
processed
of the processed image and threshold it by
a constant T , in our case T = 0.0018. The value of
T should ideally be low, to avoid getting more com-
ponents in the GCL - IPL regions which will nega-
tively impact the segmentation results. We obtain the
edges of the thresholded gradient T G image 5(a), and
then refine it in two simple steps: first, by area open-
ing ((Vincent, 1994), where any component less than
P pixels (P = 30pixels) is deleted to remove small
objects from the image (most especially the GCL to
IPL regions); second, we extrapolate incomplete layer
lines to span the image horizontally. This is carried
out by linking broken lines to the closet neighbouring
points, where each line starts from the first column
and ends at the last column of the image, such that
only complete layers are initialised, figure 5(b). The
initialisation and evolution of complete layer lines ex-
clusively is an essential factor, which further ensures
accurate segmentation, i.e. no merging or splitting of
boundaries is allowed, which prevents under- or over-
segmentation. Without losing context, the refinement
step can be ignored, if the layers from the GCL to INL
BIOIMAGING 2019 - 6th International Conference on Bioimaging
50
Figure 1: Location of Nerve Fibre Layer (NFL); Ganglion Cell Layer + Inner Plexiform Layer + Inner Nuclear Layer
(GCL+IPL+INL); Outer Plexiform Layer (OPL); Outer Nuclear Layer to Inner Segment (ONL+IS), Outer Segment (OS),
Retinal Pigment Epithelium (RPE) and the total retinal thickness, on an OCT image.
Figure 2: Schematic representation of the proposed level set approach.
are targets of the method. However, this will require
a condition for handling the splitting and merging of
boundaries or an alternate measure to correctly iden-
tify which layers are segmented. Moving further, the
edges of the refined image serve as initial curves such
that the number of identified regions in the final out-
put cannot exceed the number of the initial curves.
An important point to note is that, if T is reduced,
then P is to be increased, mainly because the size of
the small objects in 5 (b) will increase with a smaller
Level Set Segmentation of Retinal OCT Images
51
Figure 3: Preprocessing steps showing: Column 1 - Enhanced images; Column 2 - identified ILM (red) and RPE (Green);
Column 3 - image masks I
mask
; and Column 4 - Cropped images I
cropped
. Row 1 - Nasal region; Row 2 - Foveal Region; and
Row 3- Temporal Region.
Figure 4: Gradient of full image ( Figure - 3 Column 1, Row 1) with background noise and layer-like structures in red (a), and
thresholded gradient of preprocessed (Figure 3 - Column 4, Row 1) image TG with ROI only (b).
T . Each boundary curve C
b
is therefore represented
by a collection of C
b
(x,y) on the image.
2.3 Segmentation of OCT
Next, each initial boundary curve C
b
is evolved de-
pending on a speed field F based on the following
differential equation (Shi and Karl, 2005) :
dC
b
dt
= F
~
N (2)
Where
~
N is the normal of the curve pointing out-
ward. The speed field F is made of an external speed
derived from the image data and a characteristic speed
BIOIMAGING 2019 - 6th International Conference on Bioimaging
52
Figure 5: Edges before refinement (a), and refined edges used for contour initialisation (b).
based on C
b
. We associate F and the ensuing evolu-
tion with a gradient descent solution to solve the min-
imisation problem based on a Mumford-shah model
evolution perspective (Tsai et al., 2001). This means
the curve C
b
will evolve until it gets to a local minima
C
bmin
of the energy, i.e. static point of the dynamic
equation (2). Adapting from (Shi and Karl, 2005), we
represent a layer boundary uniquely through two lists
of inside L
in
and outside L
out
points of C
b
. Which are
defined as:
L
out
= {x|φ(x) > 0 and y N(x) : φ(y) < 0}
L
in
= {x|φ(x) < 0 and y N(x): φ(y) > 0}
where N(x) is a distinct neighbourhood of x, in the
level set function φ at pixel x. Based on this definition,
a positive force switches a point from L
in
to L
out
and a
negative force switches a point from L
out
to L
in
. Each
point (x, y) in level set function is defined in relation
to the curve C
b
as follows (Shi and Karl, 2005):
φ =
3, if x,y is outside C
b
and x,y 6∈ L
out
;
1, if x,y L
out
;
3, if x,y is inside C
b
and x,y 6∈ L
in
;
1, if x,y L
in
.
(3)
Based on the definition of φ in equation (3), it
is undemanding to recognise the location of a point
(x,y) on the image in relation to C
b
. In our formu-
lation, we use a 2D list to represent initial bound-
ary points, and to save the positions of final bound-
ary points for straightforward mapping in generating
the final image output. Alternatively, a 1-D list can
be used to save the boundaries point at position φ, as
suggested in (Liu et al., 2018). A boundary position
(x,y) of φ, can either expand or shrink based on:
Expand(x,y) : C
b
(x,y) := C
b
(x,y) + 1
Shrink(x,y) : C
b
(x,y) := C
b
(x,y) 1
(4)
The evolution of each point is influenced by the
image forces computed by a fast gradient vector
field(Wang and Boyer, 2012) and the topology con-
straints to be described in the next subsection 2.4.
2.4 Topology Constrains
As highlighted earlier in subsection 2.1 the ordering
of the layers must be preserved. Taking into account
the architecture of the OCT image, boundary C
b2
is
always below C
b1
for any given boundary points C
b1
and C
b2
, i.e. a point (x,y) on the curve will neither
Shrink nor Expand if it makes C
b1
(x,y) C
b2
(x,y).
Hence, with our appropriate initialisation, we enforce
the topology requirement by carrying out this simple
topology validation before either shrinking or expand-
ing a boundary. Finally, we employ an intuitive ap-
proach to ensure this topology is preserved, by addi-
tionally refining the topology constraint in the vertical
direction:
1. Because each layer boundary spans the image hor-
izontally (one boundary point per column) we add
a condition for evolving a boundary point C
b
(x,y)
to a new boundary point C
bmin
(x,y). We restrict
Expand(x, y) if its neighbour points are u con-
secutive points above it; do not Shrink(x, y) if its
neighbour points are u consecutive points below
it;
Level Set Segmentation of Retinal OCT Images
53
2. Looking at the sample of initial layer boundaries
in figure 5, a boundary point C
b
(x,y) is limited
to a maximum of v operations (either Expand or
Shrink) consecutively in the vertical direction).
The parameters u and v are two prior constants, in
our experiments u and v are set to 3 and 20 respec-
tively. The parameter u aids with boundary smooth-
ness and avoiding peaks for Expand(x, y) or valleys
for Shrink(x, y) on the boundaries, while v further en-
sures the layered architecture is preserved. Addition-
ally, this is why our layer initialisation is ideal be-
cause the starting points are based on the individual
image. The topology constraints facilitate the evolu-
tion because the validation is performed before ex-
panding or shrinking a boundary. Perhaps, this might
not be ideal for abnormal structures. However, con-
sidering the ordering of the layers where C
b
2 will al-
ways be below C
b
1 the layers will move together even
in the case of abnormal retinal structure. The pseudo-
code of our algorithm is illustrated in Algorithm 1.
3 RESULTS AND DISCUSSIONS
We applied our method to 200 macula OCT images.
The original size of each image is 512 by 992 pixels,
with a resolution of 16 bits per pixel. We crop the
image in the pre-processing stage to improve the re-
sults of our method. Comparison to the ground truth
labelling is carried out on full image size. In our ex-
periments N(x) = 8 neighbourhood, mainly, because
only the layers are remaining in the cropped image
and the effect of inhomogeneity is reduced. Experi-
mental results show that our method successfully seg-
ments seven (7) layers of the retina. Samples of the
method output are shown in figure 6.
Table 1 shows the mean and standard deviation
for the performance of the method compared to the
labelling of manual graders. The values show the
promising performance of our method in converging
at curves C
bmin
very close to the actual layer bound-
aries. Notably, the RNL thickness is used for diag-
nosing major eye diseases such as glaucoma, and the
mean (0.951) and standard deviation (±0.022) of dice
coefficient for this layer is reassuring.
Moreover, it can be deduced that our method is
consistent in identifying the layer boundaries from
the distribution of the values in figure 7. Consider-
ing the second quarterlies of the NFL, IS and RPE
begin at 0.900 further attests to the optimum per-
formance of our method, except for few instances in
the GCL+IPL+INL and OPL layers, where the dice
score is below 0.800. However, in few instances the
method could not properly identify the GCL-INL and
Figure 6: From top to Bottom: Sample results from Nasal,
Foveal and Temporal regions respectively.
Table 1: Performance evaluation: Mean and Standard Devi-
ation (STDEV) of Dice Coefficient on 200 B-Scan images
(Units in pixels).
RetinalLayer Mean ST DEV
NFL 0.951 ± 0.022
GCL+IPL+INL 0.879 ± 0.031
OPL 0.892 ± 0.032
ONL 0.907 ± 0.030
IS 0.932 ± 0.017
OS 0.920 ± 0.028
RPE 0.934 ± 0.021
the OPL due to some of the small components not
been removed. The method avoids over and under
segmentation, due to our layer initialisation and topo-
logical constraint, which prevents merging or splitting
of boundaries. On the other hand, the method will
come in short when compared to studies targeting the
choroidal region, which is believed to provide details
on some of the visual impairments (Sun et al., 2016).
BIOIMAGING 2019 - 6th International Conference on Bioimaging
54
Algorithm 1: Boundary Evolution.
1: Initialise Boundaries
2: loop:
3: if evolution will not make C
b1
(x,y) C
b2
(x,y) then
4: %% Shrink Boundary
5: if neighbours of C
b
(x,y) 6= consecutive u points below it then
6: if C
b
(x,y) has not moved v consecutive points in the vertical direction then
7: if force at point is negative then
8: Shrink(x, y)
9: %% Expand Boundaries.
10: if neighbours of C
b
(x,y) 6= consecutive u points above it then
11: if C
b
(x,y) has not moved v consecutive points in the vertical direction then
12: if force at point is positive then
13: Expand(x, y)
14: %% C
b
(x,y) are at relative local minima.
15: if no changes made to all C
b
(x,y) then
16: break
Figure 7: Box plot of mean Dice Coefficient distribution for the seven (7) layers.
4 CONCLUSIONS
We have presented an automatic level set method
for retinal OCT segmentation. The proposed work
separates retinal OCT images into seven (7) non-
overlapping layers. Our approach has explored image
segmentation using level set from the point of initial-
isation, and the constraining of curve evolution based
on retinal layer topology explicitly. Refined edges
of gradients images are used to initialise the curves.
Image forces constrained by the topological architec-
ture of the OCT are used to guide the evolution of
the curve to its minimum layer boundaries. These
two components ensure the boundaries obtained by
the method are close to the actual features of interest.
The proposed method takes advantage of handling the
obstruction of image background noise in the pre-
processing stage, consequently making the segmen-
tation process to suffer less from the image artefacts.
Additionally, refinement of the gradient edge infor-
mation ensures only the targeted layers are initialised
at the beginning of the evolution process. Experimen-
tal results show that the proposed approach success-
fully segmented the target layers from OCT images,
and the segmentation results are close to the manu-
ally labelled ground-truth.Future work will seek to in-
clude the GCL to IPL and choroid regions in the ROI,
which has to do with implicitly assigning the param-
eters T and P based on the image as opposed to con-
stants used in our approach.
REFERENCES
Adhi, M. and Duker, J. S. (2013). Optical coherence
tomography–current and future applications. Current
Level Set Segmentation of Retinal OCT Images
55
opinion in ophthalmology, 24(3):213.
Boyer, K. L., Herzog, A., and Roberts, C. (2006). Auto-
matic recovery of the optic nervehead geometry in op-
tical coherence tomography. IEEE Transactions on
Medical Imaging, 25(5):553–570.
Chiu, S. J., Li, X. T., Nicholas, P., Toth, C. a., Izatt, J. a., and
Farsiu, S. (2010). Automatic segmentation of seven
retinal layers in SDOCT images congruent with expert
manual segmentation. Optics express, 18(18):19413–
19428.
Dijkstra, E. W. (1959). A note on two problems in connex-
ion with graphs. Numerische mathematik, 1(1):269–
271.
Dodo, B. I., Li, Y., Eltayef, K., and Liu, X. (2018). Graph-
cut segmentation of retinal layers from oct images. In
Proceedings of the 11th International Joint Confer-
ence on Biomedical Engineering Systems and Tech-
nologies - Volume 2: BIOIMAGING,, pages 35–42.
INSTICC, SciTePress.
Dodo, B. I., Li, Y., and Liu, X. (2017). Retinal oct im-
age segmentation using fuzzy histogram hyperboliza-
tion and continuous max-flow. In 2017 IEEE 30th In-
ternational Symposium on Computer-Based Medical
Systems (CBMS), pages 745–750. IEEE.
Duan, J., Tench, C., Gottlob, I., Proudlock, F., and Bai,
L. (2017). Automated segmentation of retinal lay-
ers from optical coherence tomography images using
geodesic distance. Pattern Recognition, 72:158 – 175.
Garvin, M. K. (2008). Automated 3-D segmentation and
analysis of retinal optical coherence tomography im-
ages. PhD thesis - The University of Iowa.
Huang, D., Swanson, E. A., Lin, C. P., Schuman, J. S.,
Stinson, W. G., Chang, W., Hee, M. R., Flotte, T.,
Gregory, K., and Puliafito, C. A. (1991). Optical
coherence tomography. Science (New York, N.Y.),
254(5035):1178–81.
Jaffe, G. J. (2012). OCT of the Macula: An expert provides
a primer on useful scans, identifying artifacts and time
domain vs. spectral domain technology. Reinal Physi-
cian, pages 10–12.
Koozekanani, D., Boyer, K., and Roberts, C. (2001). Reti-
nal thickness measurements from optical coherence
tomography using a Markov boundary model. IEEE
Transactions on Medical Imaging, 20(9):900–916.
Lang, A., Carass, A., Hauser, M., Sotirchos, E. S., Cal-
abresi, P. a., Ying, H. S., and Prince, J. L. (2013).
Retinal layer segmentation of macular OCT images
using boundary classification. Biomedical optics ex-
press, 4(7):1133–52.
Liu, Y., Carass, A., Solomon, S. D., Saidha, S., Calabresi,
P. A., and Prince, J. L. (2018). Multi-layer fast level
set segmentation for macular oct. In 2018 IEEE
15th International Symposium on Biomedical Imaging
(ISBI 2018), pages 1445–1448.
Lu, S., Yim-liu, C., Lim, J. H., Leung, C. K.-s., and Wong,
T. Y. (2011). Automated layer segmentation of op-
tical coherence tomography images. Proceedings -
2011 4th International Conference on Biomedical En-
gineering and Informatics, BMEI 2011, 1(10):142–
146.
Novosel, J., Vermeer, K. A., Thepass, G., Lemij, H. G., and
Vliet, L. J. V. (2013). Loosely Coupled Level Sets
For Retinal Layer Segmentation In Optical Coherence
Tomography. IEEE 10th International Symposium on
Biomedical Imaging, pages 998–1001.
Raftopoulos, R. and Trip, A. (2012). The Application of
Optical Coherence Tomography ( OCT ) in Neurolog-
ical Disease. Advances In Clinincal Neuroscience and
Rehabilitation, 12(2):30–33.
Shi, Y. and Karl, W. C. (2005). A fast level set method
without solving pdes [image segmentation applica-
tions]. In Proceedings. (ICASSP ’05). IEEE Interna-
tional Conference on Acoustics, Speech, and Signal
Processing, 2005., volume 2, pages ii/97–ii100 Vol. 2.
Sun, Y., Zhang, T., Zhao, Y., and He, Y. (2016). 3d auto-
matic segmentation method for retinal optical coher-
ence tomography volume data using boundary surface
enhancement. Journal of Innovative Optical Health
Sciences, 9(02):1650008.
Tian, J., Varga, B., Somfai, G. M., Lee, W. H., Smiddy,
W. E., and DeBuc, D. C. (2015). Real-time automatic
segmentation of optical coherence tomography vol-
ume data of the macular region. PLoS ONE, 10(8):1–
20.
Tsai, A., Yezzi, A., and Willsky, A. S. (2001). Curve evolu-
tion implementation of the mumford-shah functional
for image segmentation, denoising, interpolation, and
magnification. IEEE Transactions on Image Process-
ing, 10(8):1169–1186.
Vincent, L. (1994). Morphological area openings and clos-
ings for grey-scale images. In Shape in Picture, pages
197–208. Springer.
Wang, C., Wang, Y., Kaba, D., Wang, Z., Liu, X., and Li., Y.
(2015). Automated layer segmentation of 3d macular
images using hybrid methods. In Proc. International
Conference on Image and Graphics. Tianjing, China.,
volume 9217, pages 614–628.
Wang, Q. and Boyer, K. L. (2012). The active geometric
shape model: A new robust deformable shape model
and its applications. Computer Vision and Image Un-
derstanding, 116(12):1178 – 1194.
Yazdanpanah, A., Hamarneh, G., Smith, B. R., and Sarunic,
M. V. (2011). Segmentation of intra-retinal layers
from optical coherence tomography images using an
active contour approach. IEEE Transactions on Med-
ical Imaging, 30:484–496.
BIOIMAGING 2019 - 6th International Conference on Bioimaging
56