Computers in Biology and Medicine 35 (2005) 915 – 931 www.intl.elsevierhealth.com/journals/cobm
Reconstruction and representation of caudal vasculature of zebrafish embryo from confocal scanning laser fluorescence microscopic images Jun Fenga , Shuk Han Chengb , Po K. Chanb , Horace H.S. Ipa,∗ a Image Computing Group, Department of Computer Science, City University of Hong Kong, 83 Tat Chee Avenue, Kowloon,
Hong Kong b Department of Biology and Chemistry, City University of Hong Kong, 83 Tat Chee Avenue, Kowloon, Hong Kong
Received 23 December 2003; accepted 19 May 2004
Abstract Three-dimensional (3D) reconstruction from a series of sections is an important technique in medical imaging, particularly for visualization of blood vessels from angiography. Here, we present a framework for automatic segmentation and registration of different kind of blood vessels from 2-day-old zebrafish embryos. Series of optical sections were acquired from confocal microscopy with the blood vessels labeled by fluorescent microbeads (0.02 m) injected into blood stream of 2-day-old zebrafish embryos. Blood vessels were extracted and their morphological parameters, including length and diameter, were calculated. At the same time, individual blood vessels were registered automatically. Vasculature was represented by attributed vessel represent graph (AVRG), which contained morphological data and connectivity of every blood vessel. Using AVRG to represent a vasculature made the comparison between vasculatures of different embryos more easy. Visualization, as well as quantification, of reconstructed 3D model of AVRG was presented in an interactive interface. The framework was implemented by Visual C++ as Windows-based program. 䉷 2004 Elsevier Ltd. All rights reserved. Keywords: Medical imaging; 3D vascular reconstruction; 3D vascular tree analysis; Volume visualization; 3D medical image processing
∗ Corresponding author. Tel.: +852-27888641; fax: +852-27884916.
E-mail address:
[email protected] (H.H.S. Ip). 0010-4825/$ - see front matter 䉷 2004 Elsevier Ltd. All rights reserved. doi:10.1016/j.compbiomed.2004.05.003
916
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
1. Introduction The vascular system is an organ system for oxygen and nutrient delivery as well as a conduit for communication between distant tissues. The importance of the vascular system in embryos was demonstrated by the lethality of animals with vascular defects [1–4]. The zebrafish is a model organism in the study of vascular development with some unique benefits, such as external fertilization and transparent embryos [5,6]. Several techniques have taken advantage of the optical clarity of zebrafish embryos to observe circulatory patterns and the developing vasculature. Microangiography is one of the techniques, in which fluorescent microbeads with size of 0.02 m are injected into the circulation of a zebrafish embryo [7]. The microangiograph is visualized using a confocal scanning laser microscopy (CSLM) (Fig. 1a). This technique has been used to characterize aberrant vascular patterns that are associated with phenotypes of zebrafish gridlock mutant [7] and provide a detailed description of wild-type circulatory patterns during early development [8]. Recently, we adopted microangiography to study the adverse effects of cadmium on the early development of vascular system in zebrafish embryos [9]. Our approach is based on serial images obtained from microangiography. Although these serial images provides three-dimensional (3D) volumetric information, to our knowledge, no 3D quantification and analysis of such images of the zebrafish embryo has been done. In this paper, we present techniques for automatic analysis and quantification in three-dimensions of the vasculature of the zebrafish embryo, which will be applicable to mass screening of such embryos in toxicology and applied genomic studies. The origin of 3D reconstruction from serial sections dated back to 1880 in the study of human embryos [10]. It was not until early the 1970, when computerized techniques were sufficiently well developed to assist in such work, that 3D modeling of small structures became possible [11,12]. With the advent of CT and MRI imaging, computerized object reconstruction from sectional images became a common practice in medicine, particularly in the area of surgical planning [13]. Nowadays, very realistic model and visualization can be achieved for human anatomy [14]. Recent work on simulated experiments includes simulating blood flow in human vessel by a system called ViVa [16]. Initially, 3D reconstructions were purely done manually, but over the years and with the advent of computer image processing, semi-automatic techniques were introduced to aid image registration and
Fig. 1. (a) Vasculature of a 2-day-old zebrafish embryo. Red box indicates the caudal region of vasculature. (b) Typical images of three parts of caudal vasculature of 2-day-old zebrafish embryo.
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
917
segmentation to make the reconstruction less time consuming [15]. There also exists commercially available software, which can stack the series of cross-sections together mainly for 3D visualization such as Imaris and Volocity. Some recent systems construct the 3D mesh model of vessels for quantification and visualization [16]. However, such systems do not provide a method for automatic morphological analysis and comparison of the vessel connectivity patterns and their junctions. In this paper, we present techniques and a computational framework for reconstruction, representation and visualization the zebrafish vasculature from confocal microangiographic images. Serial silhouette optical sections of embryonic vasculature are prepared with fluorescent microbeads and recorded by CSLM [9]. Based on anatomical atlas of the zebrafish vasculature [8], we have developed automatic techniques to segment, track and register individual blood vessel. Specifically, the branching pattern of the zebrafish vasculature is represented by an attributed vessel representation graph (AVRG), which captures all the connectivity as well as volumetric information of each individual blood vessel within the vascular network. The topological structure of embryonic vasculatures and the associated volumetric data are subsequently stored in a relational database for generally available for further analysis and comparison by researchers such as developmental biologists. Based on a set of serial confocal images, our system is able to automate the entire process of 3D image pre-processing, image analysis, vasculature reconstruction and quantification and the visualization of reconstructed embryonic vasculature in three-dimensions and, thus, reduce significantly the time for extracting and analyzing information on the zebrafish vasculature from raw confocal images. Due to the space constraint in this paper, instead of elaborating on every step of the vasculature reconstruction, we focus on the framework and detail of the reconstruction and morphological analysis as well as the application of this framework in genomic experiments, and present the image pre-processing parts (e.g. image enhancement and binarization) briefly. Further detail can be found in our previous references.
2. Specimen preparation and image processing 2.1. Specimen preparation and data acquisition Microangiography in zebrafish embryos using fluorescent microbeads was first adopted by Weinstein et al. [7]. Procedures of microinjection and microangiography were described in our previous paper [9]. Fluoresceinated carboxylated latex micro-beads of diameter 0.02 m were obtained commercially (Molecular probes, USA). Fifty microliter bead suspension was diluted (1:1) with 2% bovine serum albumin in double distilled water. Diluted fluorescent beads were sonicated for 1 min at maximum power on a sonicator and then centrifuged for 5 min in an Eppendorf microcentrifuge. Dechorionated 2-dayold embryos were an anesthetised by incubating in embryo medium containing 16.8 mg/100 ml Tricane (Sigma) for 3 min at 28.5 ◦ C. Embryos were then transferred to a 25 mm Petri dish for micromanipulation. Holding pipettes and microinjection pipettes were made from borosilicate glass capillaries (Cat# MTW100-4, World Precision Instruments Inc.) using a vertical puller (Narishige, Japan). The tip of each pipette was broken and then beveled at a 30◦ angle with a microgrinder (Narishige, Japan). Sharpened tip diameters were 300 m for holding pipettes and 20 m for microinjection pipettes. Five microliter of bead suspension was taken up by the microinjection pipette and approximately one-tenth volume was injected into the posterior region of sinus venosus. Injected embryos were scanned with a confocal
918
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
microscope equipped with LSM 5 (Zeiss, Germany) within 15 min. To get images at high detail, images were acquired by 10× objective lens, with addition of 2× digital zoom. At such magnification, short portion of dorsal aorta segment with only six or seven intersegmental vessels could be viewed. Therefore, vasculature of caudal region was subdivided into three parts for image acquisition (Fig. 1b). To rebuild the caudal vasculature, there was at least one intersegmental vessel overlapped in each part of images. The resolution was 0.9 m in x, y and z direction. The sequences of confocal images were then exported in TIFF format with 8-bit intensity depth. Image size was 512 pixels by 512 pixels 2.2. Image restoration, enhancement and binarization The resolution of the optical images obtained from confocal microscopy is orientation dependent [24] and for the specimen stained with fluorescent probes, distortion arises from the fact that the incoherent fluorescent light, which is emitted in all directions, is very difficult to collect fully by passing the beam through the finite aperture of the objective lens. To improve the signal to noise ratio (SNR), we increased the microscopy aperture, which allows more out-of-focus light to pass through, but at the same time, the image is degraded as a result of superposition of blurred images from out-of-focus planes. To address this problem, we have developed a fast method to enhance the images from CSLM [17]. Specifically, we split the 3D enhancement process into two separate stages and adopt different methods, respectively, to deconvolve the images in the lateral and optical directions. The resulting algorithm increased both lateral and optical resolution and reduced the effects due to out-of-focus light in the images. The enhanced volumetric image is characterized by an improvement both in contrast and in edge definition. Assuming the imaging system to be linear, we derived in [17] the theoretical expression for 3D point spread function of CSLM in the idealized case as 4 1 2 2 2 2 h(x, y, z) = 2 sin x + y /|z| x + y /|z| . (1) z The resulting 3D image (volumetric data) is a convolution of the object with the 3D impulse response of the CSLM. Next, we approximate the 3D specimen by a stack of N object planes separated at equal intervals z along the z-axis and obtain the finite summation g(x, y, j z) ≈
N
f (x, y, i z)h(x, y, j z − i z)z.
(2)
i=1
Eq. (2) implies that, for CLSM, the current image plane can be enhanced by considering and subtracting only the blurring effects due to the two immediate adjacent out-of-focus image planes. After that, we model the imaging as a Wiener estimator with an additive noise source n(x, y), which can be used to get the approximate immediate adjacent out-of-focus image planes for a in-focus one. The transfer function T (u, v) of the optimal deconvolution filter in the mean-square sense [23] is T (u, v) =
H ∗ (u, v)Pf (u, v) . |H (u, v)|2 Pf (u, v) + Pn (u, v)
(3)
Further detail of application of this filter to enhancing CSLM images by the authors can be found in [17].
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
919
Fig. 2. Examples of images after enhancement and thresholding.
To binarize the potential vessel regions in the enhanced images, we have developed an adaptive thresholding method that makes use of the 3D information available in a stack of serial images. Our approach first estimates the principal orientation of the major vessels within the image volume and analyses the intensity profiles perpendicular to the principal vessel orientation to obtain an initial threshold for each image sub-region adaptively. To exploit the spatial coherence of vascular intensities, the threshold values of those neighboring regions near the current one both along the x–y and the optical directions were referenced to ensure the continuity of the binarized content in the resulting thresholded images. The developed technique is characterized by suppressing both noise and non-interesting structure such as tissues conglutinated to vessels, while improving the vessel junctions and edge definitions. Examples of images resulting from the enhancement and thresholding are shown in Fig. 2. The detail of the binarization technique can be found in [18]. The main artery and vein of the zebrafish vasculature are the caudal artery (CA) and the caudal vein (CV), respectively. To estimate and transform the entire vasculature so that these major vessels lie along the horizontal axis of the image plane, we take the most left-down point (p1 ) and the most right-down point (p2 ) of the caudal vein in the middle image of the stack and rotate the entire image stack by an angle obtained from (4) in which y is the vertical direction and x is the horizontal direction. This alignment process is done automatically. The alignment process allows us to digitally re-slice the image stack from different directions which are orthogonal to the major and minor vessels in our vessel segmentation process. = arctg[(p2 y − p1 y)/(p2 x − p1 x)].
(4)
2.3. Vasculature segmentation and registration This section presents the automated techniques for segmenting, tracking and identifying the individual vessels within the caudal vascular network of embryo. Our approach digitally dissects the voxel data and obtained multi-orientation serial sections along different cutting planes of the zebrafish to simplify the vessel segmentation and tracking processes. The resulting algorithm successfully identified and labeled not only major vessels such as Caudal Vein and Caudal Artery but also minor ones such as the intersegmental vessels (ISVs). Our technique first detects the outline of the caudal artery and caudal vein using edge tracking and curve fitting as shown in Fig. 3(a). These initial longitudinal vessel outlines then serve to provide the expectation windows for segmenting the vessel cross-sections along the transverse planes. These transverse sections are obtained by digitally cutting along the dorsal–ventral axis of the zebrafish. From these transverse sections, we can easily distinguish the cross-sections of CV, CA or plexus based on the expectation regions derived from the longitudinal outlines as shown in Fig. 3(b).
920
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
Fig. 3. (a) Outlines of CA and CV were detected from the maximum projected thresholded image by edge detection algorithm. (b) Identification of CA and CV. Their outlines shown in (a) were then served as expectation windows (red lines) for segmentation CA and CV from transverse cross-sections. (c) Classification of blood vessels with different colors at horizontal planes. ISVs: blue; CA: red; DLAV: green; VTA: yellow.
After that, intersegmental vessels and dorsal longitudinal anastomotic vessels (DLAV) are segmented from the horizontal sections along the ventral to dorsal axis of zebrafish embryo. From these horizontal section images, we can observe that the ISVs connects the VTA (vertebral artery) and CV at roughly perpendicular angles, its shapes and areas as seen on the images are significantly different from the VTA and CV. To identify the different types of vessels based on their shape, size and order of appearance through the image sequence, we develop a classification algorithm based on the following a priori knowledge as we step through the image sequence from the ventral to the dorsal direction: (1) The cross-sections of CV are the ones to be observed first from the stack of horizontal images and their elongation values: elongation = length/width > T1 and areas > T2 . (2) ISV cross-sections can be observed after that of the CV, elongation < T1 and area < T2 . (3) VTA sections emerge half way towards the stack of horizontal images and they connect ISV. (4) DLAV is the last vessel that could be seen in the horizontal image stack. It connects with the dorsal end of ISV. (4) T1 and T2 could be easily obtained from statistic values of several fish images. Fig. 3 (c) shows the result of the classification, where different vessel classes have been labeled with different colors. In this figure, red represents CV, blue means ISV, green stands for DLAV and VTA is represented as yellow. It can be seen that the VTA has been segmented and classified successfully from ISVs even though they are connected together prior to the processing. Detailed information could be found in [19]. After segmentation, the resulting contour of each vessel provides the basis for reconstructing the surface model of the vasculature in the next stage. 3. Vasculature reconstruction 3.1. Construction of 3D surface models from volumetric data In computer graphics, 3D solid objects can be represented by means of their enclosing surface. Such a model is hierarchical in the sense that the surfaces are in turn represented by means of their enclosing curves (loops) and so on. The geometry, that is, the actual description of the object shape, is represented
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
921
Fig. 4. (a) Cross-sections of a segment of blood vessel. (b) Triangle meshes between two consecutive cross-sections. (c) Connections between blood vessels (blue points) in caudal vasculature. Red lines represented the medial axes of blood vessels.
by a collection of points, curves and surfaces, while the topology, that is, the description of how the geometrical elements are joined to form the global model, is constituted by vertices, edges, loops, faces, shells and blocks. The 3D geometrical models or their wire frames consume much less time and space than volumetric-based solid objects in both visualization and quantification. Since the objects of interests in this work are vessels with an inherently tubular structure, we can model these vessels in terms of elliptical tubes winding in a 3D space. The hierarchical geometrical elements in our global model are described as follows: Voxels: Voxels are the spatial points located at the contours of the vessels. Each voxel has a location vector [x, y, z], which represents its local position, with the local reference with respect to the local reference coordinate system. Cross-sections: A cross-section is the resulting curve intersecting a vessel tube with the plane perpendicular to the orientation of the current piece of the vessel. Because of the inherent nature of blood vessels, the cross-sections can be considered as ellipses. In general, an ellipse in 3D space can be defined in terms of a center voxel C, and two axis radii, the major and minor axis NR1 and NR2 . Each cross-section therefore has a normal vector N [nx, ny, nz], which defines the direction of the current piece of the vessel. −−−−−−→ We have N = (C2 − C1 ), where C1 is center point of current piece of vessels while C2 is center points of adjacent one. The spatial relationship between NR1 , NR2 , and N is NR1 ⊥ NR2 ⊥ N. The contours of cross-sections in our system are represented by collections of voxels (see Fig. 4(a)). Triangle meshes: Triangle meshes are built between successive cross-sections. Every triangle is formed by three voxels: one is an element included in the first cross-section CS1 , the other two are located in the adjacent cross-section CS2 . The spatial relationship between these three voxels must be V1 = Min(|V , (V2 − V3 )|)(V is collections of voxels of CS2 ).
922
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
That means V1 is the voxel which has the minimum distance to the line (V2 , V3 ). This voxel selection specification can ensure no interlacing triangle existed in the surface of the vessel. The set of triangle meshes generated this way constitute the whole vascular surface and the subsequent wire frame of the vasculature and rendering of the vasculature are based on these meshes (see Fig. 4(b)). Triangle is preferred instead of rectangle or polygon is because four points in two cross-sections (two in CS1 and other two in CS2 ) might not be in the same plane, since two consecutive cross-sections may have different radius and orientation and a 3D spline surface will increase the computation load. Vessels: Vessels can be defined as the curve surfaces of elliptical tubes constructed from a series of cross-sections. Each vessel is represented and stored independently, and therefore can be visualized and analyzed individually. Junctions: Junctions are voxels which are the intersections of the center axis of two or more vessels. These junctions in a vasculature provide very important information of the blood vessels structures and circulation pattern that can be used to characterize the vasculature. Besides the properties of general voxels which can be represented by a 3D location of (x, y, z), the internal representation of a junction also holds information of the vessels which it connects with. Some main junctions in a caudal vasculature of zebrafish embryo are shown in Fig. 4(c). Based on the above hierarchical elements, the 3D geometric model of the vasculature is automatically built using the following procedures: (1) From the volumetric model obtained from segmenting and vessel labeling the process described in the previous sections, the central point for each cross-section is computed. The contour points of every cross-section are treated as “Voxels” and labeled as belonging to a particular vessel. (2) Assuming each contour of the cross-section to be ellipse, we compute its major and minor axes and the normal direction of the cross-section. The normal vector (N) is estimated by the segment linking the central points of two consecutive cross-sections. The major axis direction is the vector (NR1 ) which has the maximum distance cross the ellipse. Therefore, the direction of minor axis is: NR2 = NR1 × AN. Here, “×” is the cross-product of two vectors. (3) Resample a suitable number of contour points of each cross-section. For every two adjacent points, find the other points in the next section based on the criterion described above and link them to form a triangular tile. It should be noted that these triangular meshes are generated in the real time during visualization, so the number of the contour points as well as the number of cross-sections of the visualization can be modified interactively by users.
3.2. Vasculature visualization A graphical interactive interface for vasculature visualization has been developed for a general PC platform (see Fig. 5). Although during image capturing, the caudal vascular has been recorded separately in three different portions, the visualization subsystem can combine them together to form the whole caudal vasculature. Therefore, the entire vasculature can be viewed from different angles and perspective interactively and rendered in color. Biologists can select different visualization effects such as shaded display, wire frame or just the vessel cross-sections. The number of the cross-sections and voxels per section to be displayed can be adjusted in the real time dependent on how much details or how fast the users expect to view the models. Furthermore, quantification information of the entire vasculature or an individual vessel being selected is also computed and displayed in real time during visualization.
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
923
Fig. 5. Program interface of ZEBVARS. (a) Whole caudal vasculature was reconstructed and visualized. (b) CV was selected and browsed with its morphological data. (c) The interface for vasculatures comparison.
4. Vasculature quantification and representation 4.1. Morphological measurements of the vessels Morphological description of a vascular system is related to its physiological function.Any change in the morphology, e.g. length and diameter, may affect the efficiency of delivery in vascular system. Therefore, it is important to represent the morphological parameter of a vasculature in systematic manner for further data analysis. Here, we enumerate several quantitative morphological and topological measurements of the vessels and group them under different categories as follow: (1) Number and connectivity: Vascular structure and blood vessel amount are the most important anatomical characteristics. Missing one or more vessels will directly influence the perfused area. Disconnection between vessels could cause the alteration in the flow of blood and thus enlargement of surrounding blood vessel to compensate the requirement of blood delivery. To get a better description of the vasculature and to facilitate the automatic analysis of the effects of mutations, we represent a vascular in terms of AVRG, which will be described in the next subsection. (2) Location and separation: Vessel position and its distance to neighboring vessels are also meaningful in the study of circulation. For example, the ISVs sprout many pairs almost perpendicular to the horizontal planes, but their locations in the artery and the interval between consecutive pairs do not necessarily follow a fixed pattern from normal embryos to their mutations. To measure these alterations, we find the corresponding voxels between two adjacent ISVs pairs. Average location and distance are obtained along all the vessel curves. Intersection voxels of ISVs pairs and artery are also identified to the junctions. (3) Size, length and volume: Since we model blood vessels to elliptical tubes, it is natural to measure the length and the average radius of these tubes. Actual length in 3D is computed by accumulating the
924
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
Fig. 6. (a and b) Two mean axis of ISVs (red) and their linkage vectors (white lines). The MLVs of these two ISVs were showed in (c).
Euclidean distance of every two central points of consecutive cross-sections. Furthermore, we derive a parameter: bend = Real Length/mean axis to estimate the bend degree of each vessel while mean axis is the direct length from the starting point of one vessel to the end of the vessel. Volume is calculated by summation of areas of all sections. (4) Twist and orientation: Twist degrees and orientations are also important measurements of blood vessels since different branching patterns would induce various blood flow rates and directions. However, they are difficult to be accurately quantified since the vessels wriggle in the 3D space randomly. To get a rough estimation, straight lines which have the minimum distance from center points to the mean axis are generated and called “Linkage Segment (LS)” (see Fig. 6). These segment directions and lengths can be used to measure the degree of vessel bending and twist orientations. For example, in Fig. 6, two vessels in the left and right twist towards different angles, respectively. We compute a vector called the maximum linkage vector (MLV), which has the longest distance among those LS. We can see the directional angles of two vessels in x and z direction are almost right about. Average length of all the LS could measure the whole variance distance from vessels to their mean axis and the length of MLV shows whether a piece of vessel winds drastically. 4.2. AVRG Attributed hypergraph have been applied to describing 3D objects and 3D scene [20–22]. However, to our knowledge, the application of graph theory as a formal framework for representing and analyzing embryonic vasculature has not been reported. In this work, we use attributed graph as a generic and hierarchical representation for vessels pattern structure. With this representation, analysis and synthesis of vasculature structures and vessel patterns can be carried out at different resolution levels. More importantly, this representation provides a fundamental and mathematical abstract model for vasculature comparison and for describing and quantifying anatomical changes. Definition 1. A graph (G) is an ordered pair G = (V , E), where V = {Vi |1 i n} is the finite set of vertices and E = {Ej |1 j n} is the finite set of edges. Definition 2. A vessel graph (VG) is a graph which defines vessels as vertices and edges as connections. Two vertices (vessels) u and v in a VG are called neighbors, If e = {u, v} is an edge of VG. The edge e means there is a connection between u and v.
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
925
Definition 3. An AVRG is a VG defined as an ordered pair G = (V , E) associated with Av × AE , where V = {Vi |1 i n} is the finite set of vessels and E = {Ej |1 j n} is the finite set of connections. Av or AE are the finite set of vessel and connection attribute values, respectively. Vi in V may be assigned a set of values from Av and Ej in E may be assigned a set of values from AE . Definition 4. An attribute distance(AD) between two vessels V1 and V2 with attribute set AV = {AVi | 1 i n} is defined as: AD(V1 , V2 , AV) = wi × |AVi − AV2i |)1 i n, where wi is the weight value of each attribute in AV . Definition 5. Suppose AVH is a set of physical attribute values (AVH ∈ AV), Two vertices or vessels V1 and V2 are homogeneous if AD(V1 , V2 , AVH) = 0. Homogeneous is a one-to-one and onto function from V1 to V2 , we will call it the H function. Definition 6. Two vertices or vessels V1 and V2 with geometrical attributes AVG in VARG are correspondence if AD(V1 , V2 , AVG) dif, where dif is a value which is very small. Definition 7. TwoVARGs VG1 ={V1 , E1 }, E1 ={Ej |1 j n1 } and VG2 ={V2 , E2 }, E2 ={Ej |1 j n2 } (n1 n2 ) are structure isomorphic if there are two vertices a and b adjacent in VG1 if and only if H (a) andH (b) are adjacent in VG2 , for all a and b in V1 . Definition 8. Two VARGs VG1 = {V1 , E1 }, and VG2 = {V2 , E2 } are attributed isomorphic if all of the homogeneous vertices in VG1 and VG2 are correspondence. The two isomorphic operators defined in Definitions 7 and 8 above are used to measure the anatomical changes among the vasculatures. However, to obtain the set of suitable weight values and the specific physical and geometrical attribute values, a significant numbers of investigations of the normal and abnormal vasculatures must be conducted which is outside the scope of this paper. The main purpose of VARG is to characterize the anatomical changes between two embryo vasculatures. In graph theory, two graphs have exactly the same pattern are called isomorphic [23], in the sense that there is a one-to-one strict correspondence between their vertex sets that preserves edges. However, as natural biological organs, vasculatures of the zebrafish embryo rarely were exactly the same patterns. To differentiate between natural occurrence of individual variation and the anatomical abnormality due to treatment, we categorize vessel features into two classes: physical and geometrical attributes. Physical attributes characterize the rudimental properties of vessels, such as vessel types (which the system detects as a result of vessel segmentation and registration), relevant location (which distinguishes different vessels such as ISVs pairs) and so on. Other features which describe the vessel shape, size or length, etc. are generalized to geometrical attributes. 4.3. Entity-relation of vascular data—database design For a genomic experiment, a large amount of vessel data and anatomical information must be collected, stored and accessed efficiently and conveniently. A hierarchical data structure is built and the entire vasculature of an embryo consists of four layers: cross-sections, blood vessels, vessel plexus and subnetwork.
926
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
Cross-sections are digital sections perpendicular to the orientation of current piece of vessels as described in Section 2, they are the basic elements of blood vessels. The latter clusters to form vessel plexus, which are sets of vessels belong to the same vessel type. For example, there are always two DLAV and 15–20 ISVs in the caudal vasculature of the zebrafish embryo. Some vessel plexus just have a single vessel, such as the main artery and the main vein. A specific number of vessel plexus constitute a sub-network of the vasculature. In our project, the whole vasculature of the zebrafish embryo is made up from three parts: caudal, trunk and head. These four levels of vessels data are stored as the tables in a relational database as well as the connections. Their attributed values are saved in the fields of each table. Since these tables implicitly record hierarchical relationships, a single vessel such as “caudal ISV-2” in the Embryo_10 can be found by tracing the path “Embryo_10–caudal–ISV–ISV-2”. Finally, the 3D structure of zebrafish embryonic vasculature was visualized by rebuilding these data from the database. Based on the database, data of the vasculature between individual embryos were analyzed and compared. Any change in the morphological structure of a vasculature was identified.
5. Experimental results of ZEBVARS The above computational framework for representing and analyzing vasculature was applied to the reconstruction and analysis of a number of zebrafish embryos to establish some base-line parameters for our future work and to validate our automated approach. The framework was implemented on a PC platform using the Visual C++ language. The resulting system, which we called ZEBVARS (ZEBrafish VAscular Reconstruction System), was applied successfully to segment all blood vessels in the caudal part of 48-hpf zebrafish embryos. Four kinds of blood vessels were tracked and identified automatically, including dorsal aorta and caudal vein in the ventral side, ISV running between each pair of somites and DLAV in the dorsal side connecting to ISV. Morphological data, such as volume and length, are extracted and measured in three-dimensions. Totally, 13 embryos have been analyzed. The length of dorsal aorta was 1260.00 ± 76.86 m and the diameter was 20.27 ± 3.18 m. The volume of dorsal aorta was 1, 409, 795.08 ± 608, 400.38 m3 . We compared results from ZEBVARS with the results obtained from the third-party image analysis software which can only conduct manual interactive measurements in two-dimensions. The diameter of dorsal aorta at the 13th pair of ISV was 20.58 ± 2.43 m, which was comparable to the data extracted automatically from ZEBVARS. Since the volume is expected to be correlated with the diameter, we plotted the volume against the diameter measured from the third-party software (Fig. 7(a)) and ZEBVARS (Fig. 7(b)). Diameter and volume of ISV on the left-hand side of embryos, i.e. odd numbered ISV in the database, were measured. The graph showed that diameter measured from the third-party software was not related to the volume measurements (R 2 = 0.2097), while better correlation is obtained for the diameter and volume measurements from ZEBVARS (R 2 = 0.7366). It was because the measurement from the thirdparty software was derived from 2D projection image of 3D confocal images, while the measurement from ZEBVARS was derived from 3D confocal images directly. This highlights the importance and values of obtaining true 3D data for morphological analysis of vasculatures. The growth of ISV in zebrafish is guided by the process of angiogenesis. Normally there are 13 ISVs in the caudal part of zebrafish embryos. One way to describe the pattern of the ISVs in caudal portion of zebrafish embryos is the distance between each ISV. Table 1 summarized the distance between each
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
927
30 25
2
R = 0.2097
Diameter
20 15 10 5 0 0.00E+00 5.00E+05 1.00E+06 1.50E+06 2.00E+06 2.50E+06 3.00E+06
(a)
Volume 30 R2 = 0.7366
Diameter
25 20 15 10 5 0 0.00E+00 5.00E+05 1.00E+06 1.50E+06 2.00E+06 2.50E+06 3.00E+06
(b)
Volume
Fig. 7. Correlation between the volume and average vessel diameter of CA measured by a third-party commercial image analysis software (a) and ZEBVARS (b) Diameter of vessel was measured XY plane of images in the third-party commercial software while diameter was measured in 3D reconstructed model from the database in ZEBVARS. The volume was calculated and measured by ZEBVARS. Table 1 Summary of the length and distance of 13 ISV in caudal vasculature from 10 embryos ISV
1 2 3 4 5 6 7 8 9 10 11 12 13
Length
Distance
Our program
MetaMorph
Our Program
MetaMorph
67.85 ± 14.32 71.46 ± 12.07 79.38 ± 10.15 82.46 ± 13.91 88.00 ± 11.72 91.85 ± 13.65 95.00 ± 15.04 104.23 ± 8.15 109.85 ± 8.77 113.00 ± 10.24 118.46 ± 10.92 116.91 ± 13.02 114.00 ± 11.83
64.99 ± 14.16 71.91 ± 4.79 75.38 ± 4.95 80.25 ± 5.25 83.79 ± 6.09 90.46 ± 7.81 95.48 ± 6.37 99.66 ± 6.54 103.86 ± 6.12 105.75 ± 6.66 110.02 ± 6.62 116.15 ± 6.71 116.88 ± 6.27
58.95 ± 14.32 60.67 ± 18.52 67.85 ± 11.46 71.70 ± 5.18 72.82 ± 5.23 70.73 ± 8.34 71.15 ± 12.81 77.50 ± 5.71 77.74 ± 3.79 80.49 ± 5.65 82.92 ± 4.01 84.30 ± 4.36
62.30 ± 5.53 64.30 ± 6.02 70.33 ± 5.78 67.88 ± 6.39 71.91 ± 6.97 72.95 ± 7.35 73.06 ± 6.45 76.81 ± 7.82 75.20 ± 7.93 81.83 ± 7.12 82.93 ± 7.05 82.01 ± 5.52
Data from both the third-party commercial software and ZEBVARS were listed.
928
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931 100
Length of ISV
80
60 ZEBVARS 3rd Party Software
40
20
0
1
3
5
7
9
11
13
15
17
19
21
23
ISV
(a) 140
Distance to next
120 100 80
ZEBVARS 3rd Party Software
60 40 20 0
(b)
1
3
5
7
9
11
13
15
17
19
21
23
25
ISV
Fig. 8. Morphological analysis of ISV in their length (a) and distance between two consecutive vessels (b). Data from the third-party commercial software and ZEBVARS were compared and both showed similar results. The data were obtained from 13 embryos, in which only the left-hand side of ISV were measured (1, 3, 5, . . . , 25).
ISV. The number of ISV was counted from left-hand side, posterior to anterior, i.e. the ISV was counted from the tail tip. Distance between ISV increased from posterior to anterior. Similarly, the length of ISV increased from posterior to anterior. Since ISV is running between each pair of somites, the distance between ISVs and their length are the size of somites. The increase in the distance between ISV and length from posterior to anterior is due to the size of somites. This is because the wave of development of somites is from the anterior to posterior and somites at the posterior are still developing. Again, we compared the distance between ISV and the length measured by ZEBVARS and that measured with MetaMorph (see Fig. 8). The two sets of results are comparable and similar trend was also observed in the data from the third-party software that distance between ISV increased from posterior to anterior.
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
929
6. Conclusion The use of confocal microscopy provides a means to obtain 3D information from biological samples. However, automatic reconstruction and analysis of biology sample in 3D is not easy. Here, we propose a computational framework to automatically segment, register, quantify and reconstruct the 3D structure of zebrafish embryonic vasculature. Complex networking structure of the caudal vasculature was represented by AVRG, by which makes the comparison between embryos easier. Furthermore, we have implemented the framework in PC-based program, ZEBVARS, and make use of it to measure the caudal vasculature of 13 embryos. Data from ZEBVARS were compared with the third-party image analysis software which performs measurement and analysis in 2D approach. There was difference in the measurement of 3D parameter, volume. The pattern of growth of ISV was assessed by measuring the length and distance between two consecutive ISV. ZEBVARS, as well as the third-party software, showed similar trend from posterior to anterior that both the length and distance between two consecutive ISV were increased. This is constant with the developmental process of somites and ISV from anterior to posterior. Nevertheless, ZEBVARS permits automatic and reconstruction, data representation of 3D confocal images and allows rapid analysis of vasculature in 3D approach in the study of alteration of vasculature due to genetic modification and toxicant exposure.
Acknowledgements The work described in this paper was supported by CityU grant (9010005) and HKSAR Research Grant Council CERG 1210/03E and 1150/01E.
References [1] P. Carmeliet, V. Ferreira, G. Breier, S. Pollefeyt, L. Kieckens, M. Gertsenstein, M. Fahrig, A. Vandenhoeck, K. Harpal, C. Eberhardt, C. Declercq, J. Pawling, L. Moons, D. Collen, W. Risau, A. Nagy, Abnormal blood vessel development and lethality in embryos lacking a single VEGF allele, Nature 380 (1996) 435–439. [2] N. Ferrara, K. Carver-Moore, H. Chen, M. Dowd, L. Lu, K.S. O’Shea, L. Powell-Braxton, K.J. Hillan, M.W. Moore, Heterozygous embryonic lethality induced by targeted inactivation of the VEGF gene, Nature 380 (1996) 439–442. [3] F. Shalaby, J. Rossant, T.P. Yamaguchi, M. Gertsenstein, X.F. Wu, M.L. Breitman, A.C. Schuh, Failure of blood-island formation and vasculogenesis in Flk-1-deficient mice, Nature 376 (1995) 62–66. [4] F. Shalaby, J. Ho, W.L. Stanford, K.D. Fischer, A.C. Schuh, L. Schwartz, A. Bernstein, J. Rossant, A requirement for Flk1 in primitive and definitive hematopoiesis and vasculogenesis, Cell 89 (1997) 981–990. [5] N.D. Lawson, B.M. Weinstein, Arteries and veins: making a difference with zebrafish, Nat. Rev. Genet. 3 (2002) 674–682. [6] B.M. Weinstein, Plumbing the mysteries of vascular development using the zebrafish, Semin. Cell Dev. Biol. 13 (2002) 515–522. [7] B.M. Weinstein, D.L. Stemple, W.D. Driever, M.C. Fishman, Gridlock, a localized heritable vascular patterning defect in the zebrafish, Nat. Med. 11 (1995) 1143–1147. [8] S. Isogai, M. Horiguchi, B.M. Weinstein, The vascular anatomy of the developing zebrafish: an atlas of embryonic and early larval development, Dev. Biol. 230 (2001) 278–301. [9] S.H. Cheng, P.K. Chan, R.S. Wu, The use of microangiography in detecting aberrant vasculature in zebrafish embryos exposed to cadmium, Aquat. Toxicol. 52 (2001) 61–71. [10] W. His, Anatomie Menschlicher, Embryonen, vol.1, F.C.W. Vogel, Leipzip, 1880, pp. 7–12. [11] C. Levinthal, R. Ware, Three-dimensional reconstruction from serial sections, Nature 236 (1972) 207–210.
930
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
[12] D.F. Wann, et al., An on-line digital computer system for the semiautomatic analysis of Golgi-impregnated neurons, IEEE Trans. Biomed. Eng. BME-20 (1973) 233–247. [13] B. Geiger, Three dimensional modeling of human organs and its application to diagnosis and surgical planning, Ph.D. Dissertation, INRIA, 1993. [14] N. Thalmann, F. Cordier, Construction of a human topological model from medical data, IEEE Trans. Inform. Technol. Med. 4 (2) (2000) 137–143. [15] I. Carlbom, D. Terzopoulos, K.M. Harris, Computer-assisted registration, segmentation, and 3D reconstruction from images of neuronal tissue sections, IEEE Trans. Med. Imaging 13 (1994) 351–362. [16] G. Aboulaev, S.D. Cadeddu, G. Delussu, M. Donizelli, et al., ViVa: The Virtual Vascular Project, IEEE Trans. Inform. Technol. Biomed. 2 (4) (1998) 268–273. [17] H. Ip, J.J. Feng, S.H. Cheng, Three-dimensional enhancement of confocal scanning laser fluorescence microscope images for vascular reconstruction, in: 16th International Conference on Pattern Recognition, 2002, pp. 863–866. [18] J. Feng, H. Ip, S.H. Cheng, Adaptive thresholding of image series from fluorescence confocal scanning laser microscope using object intensity profiles, SPIE Medical Imaging, 2004. [19] H. Ip, J.J. Feng, S.H. Cheng, Automatic registration of vasculature from confocal scanning laser fluorescence microscope images sequences using multi-orientation dissection sections, in: First IEEE International Symposium on Biomedical Imaging, Washington DC, USA, July 7–10, 2002, pp. 249–252. [20] A.K.C. Wong, S.W. Lu, M. Rioux, Recognition and shape synthesis of 3-D objects based on attributed hypergraphs, IEEE Trans. Pattern Anal. Mach. Intell. 11 (1989) 279–290. [21] L. Rong, A.K.C. Wong, Attributed hypergraph representation for dynamic vision, IEEE Trans. Pattern Anal. Mach. Intell. 11 (1989) 279–290. [22] A.K.C. Wong, B.A. McArthur, in: H. Nasr (Ed.), Random graph Representation for 3-D Object Models, SPIE Milestone Series, vol. 72, 1991, pp. 229–238. [23] K.H. Rosen, Discrete Mathematics and its Application (English Version 4), China Machine Press, 1999. [24] X. Bertero, P. Boccacci, G.J. Brakenhoef, F. Malfanti, H.T.M. Van Der Vooret, Three-dimensional image restoration and super-resolution in fluorescence confocal microscopy, J. Microsc. 157 (1989) 3–20.
Ms. Jun Feng received her B.Sc. degree in Information Engineering from XiDian University and M.Phil degree in Computer Science from Northwest University respectively, both in Xi’an, China. She was an instructor in Northwest University. From 2000–2001, she worked as a Research Assistant in Image Computing Group, City University of HongKong where she is now a Ph.D. student. Her research interests include image processing, 3D reconstruction and pattern recognition. Dr. Shuk Han Cheng is an Associate Professor in the Department of Biology and Chemistry at the City University of Hong Kong. After graduating from the University of Hong Kong, she was awarded her Ph.D. from the Royal Postgraduate Medical School, University of London. She then received her postdoctoral training at the Ontario Cancer Institute, University of Toronto. Her current research interests in aquatic toxicology include gene expression analysis and developmental toxicity in zebrafish embryos. Dr. Po Kwok, Chan is a Research Fellow in the Department of Biology and Chemistry at the City University of Hong Kong. He was awarded his Ph.D. degree from the City University of Hong Kong on 2002. His current research projects in developmental toxicology of the zebrafish embryo include molecular mechanisms of toxicants as well as application of computational image analysis for quantification assessment of morphology and gene expression. Professor Horace H.S. Ip received his B.Sc. (First Class Honours) degree in Applied Physics and Ph.D. degree in Image Processing from University College London, United Kingdom, in 1980 and 1983, respectively. Presently, he is the Chair Professor of the Computer Science Department and the founding director of the AIMtech Centre (Centre for Innovative Applications of Internet and Multimedia Technologies) at City University of Hong Kong. His research interests include image processing and analysis, pattern recognition, hypermedia computing systems and computer graphics. Prof. Ip is a member of the Editorial Boards of the Pattern Recognition Journal (Elsevier), The Visual Computer (Springers), the International Journal of Multimedia Tools and Applications (Kluwer Academic), the Chinese Journal of CAD and Computer Graphics (The Chinese Academy of Science),
J. Feng et al. / Computers in Biology and Medicine 35 (2005) 915 – 931
931
Journal of Infrared and Millimeter Waves, (Chinese Optical Society) and a guest editor for the international journals Real-Time Imaging (Academic Press) and Real-Time Systems (Kluwer Academic Publishers). Prof. Ip serves on the International Association for Pattern Recognition (IAPR) Governing Board and served as founding cochair of its Technical Committee on Multimedia Systems, he is currently the Vice-Chairman of the China Computer Federation, Technical Committee on CAD & Computer Graphics. He was the Chairman of the IEEE (Hong Kong Section) Computer chapter, Council member of the Hong Kong Computer Society and the Founding President of the Hong Kong Society for Multimedia and Image Computing. He has published over 130 papers in international journals and conference proceedings. Prof. Ip is a Fellow of the Hong Kong Institution of Engineers (HKIE), a Fellow of the Institution of Electrical Engineers (IEE), United Kingdom, and a Fellow of IAPR.