A bioimage informatics approach to automatically extract complex fungal networks

Bioinformatics, Sep 2012

Motivation: Fungi form extensive interconnected mycelial networks that scavenge efficiently for scarce resources in a heterogeneous environment. The architecture of the network is highly responsive to local nutritional cues, damage or predation, and continuously adapts through growth, branching, fusion or regression. These networks also provide an example of an experimental planar network system that can be subjected to both theoretical analysis and experimental manipulation in multiple replicates. For high-throughput measurements, with hundreds of thousands of branches on each image, manual detection is not a realistic option, especially if extended time series are captured. Furthermore, branches typically show considerable variation in contrast as the individual cords span several orders of magnitude and the compressed soil substrate is not homogeneous in texture making automated segmentation challenging. Results: We have developed and evaluated a high-throughput automated image analysis and processing approach using Phase Congruency Tensors and watershed segmentation to characterize complex fungal networks. The performance of the proposed approach is evaluated using complex images of saprotrophic fungal networks with 105–106 edges. The results obtained demonstrate that this approach provides a fast and robust solution for detection and graph-based representation of complex curvilinear networks. Availability and implementation: The Matlab toolbox is freely available through the Oxford e-Research Centre website: http://www.oerc.ox.ac.uk/research/bioimage/software Contacts: boguslaw.obara{at}oerc.ox.ac.uk

A PDF file should load here. If you do not see its contents the file may be temporarily unavailable at the journal website or you do not have a PDF plug-in installed and enabled in your browser.

Alternatively, you can download the file locally and open with any standalone PDF reader:

https://bioinformatics.oxfordjournals.org/content/28/18/2374.full.pdf

A bioimage informatics approach to automatically extract complex fungal networks

Boguslaw Obara 1 2 Vicente Grau 0 2 Mark D. Fricker 3 Associate Editor: Olga Troyanskaya 0 Institute of Biomedical Engineering 1 Oxford Centre for Integrative Systems Biology 2 Oxford e-Research Centre 3 Department of Plant Sciences, University of Oxford , UK Motivation: Fungi form extensive interconnected mycelial networks that scavenge efficiently for scarce resources in a heterogeneous environment. The architecture of the network is highly responsive to local nutritional cues, damage or predation, and continuously adapts through growth, branching, fusion or regression. These networks also provide an example of an experimental planar network system that can be subjected to both theoretical analysis and experimental manipulation in multiple replicates. For high-throughput measurements, with hundreds of thousands of branches on each image, manual detection is not a realistic option, especially if extended time series are captured. Furthermore, branches typically show considerable variation in contrast as the individual cords span several orders of magnitude and the compressed soil substrate is not homogeneous in texture making automated segmentation challenging. Results: We have developed and evaluated a high-throughput automated image analysis and processing approach using Phase Congruency Tensors and watershed segmentation to characterize complex fungal networks. The performance of the proposed approach is evaluated using complex images of saprotrophic fungal networks with 105-106 edges. The results obtained demonstrate that this approach provides a fast and robust solution for detection and graph-based representation of complex curvilinear networks. Availability and implementation: The Matlab toolbox is freely available through the Oxford e-Research Centre website: http://www.oerc .ox.ac.uk/research/bioimage/software Contacts: The Author 2012. Published by Oxford University Press. All rights reserved. For Permissions, please e-mail: journals.permissions.com. 1 INTRODUCTION Rapid advances in imaging technologies have heralded a new era of quantitative measurements that offer considerable potential to improve our understanding of complex biological systems, particularly those that need to be considered as integrated functional units (Swedlow et al., 2003; Carpenter, 2007; Kvilekval et al., 2010). In particular, enhanced image processing techniques have greatly facilitated characterization of interconnected transport and communication networks in animal systems, including vascular systems, milk ducts, respiratory pathways and neuronal circuitry [e.g. (Chaudhuri et al., 1989; Frangi et al., 1998; Meijering et al., 2004; Rodriguez et al., 2008)]. In contrast, *To whom correspondence should be addressed. there has been remarkably little progress to characterize macroscopic microbial network development for an entire taxonomic group of organisms, namely the fungi. Unlike other biological networks, which exist within an organism, the entire growth form of fungi is as an adaptive network. Foraging saprotrophic fungi play a central role in ecosystem biology as they are the only organisms capable of complete degradation of wood in temperate forests. These fungi form extensive interconnected mycelial networks of multi-hyphal cords that forage for scarce resources that are patchily distributed in time and space (Fricker et al., 2007). The architecture of the network responds to local nutritional cues, damage or predation and the local structure is continuously remodeled through growth, branching, fusion or regression of the hyphal cords. Despite the absence of any centralized control system, these organisms exhibit complex foraging and decision making behavior as a co-ordinated individual. Unlike animal network models, microbial networks also provide an example of an extremely tractable experimental system that can be subjected to both theoretical analysis and experimental manipulation in multiple replicates (Boddy et al., 2009). Early measures of macroscopic mycelial organisation focussed on fractal dimension as a useful tool to capture aspects of the network structure as a metric (Boddy and Donnelly, 2008). However, a single summary statistic, such as the overall fractal dimension provides little opportunity to characterize the network structure and dynamics further or to explore the underlying mechanism leading to network optimisation. There has been some progress at the microscopic level to extract the physical network structure using (semi)-automated image analysis (Barry and Williams, 2011). However, to date, delineation of the network architecture from the larger soil-based systems has only been possibly manually (Fricker et al., 2007), as basic intensity-based segmentation is not effective in this context (M.D.Fricker, L.Boddy and D.P.Bebber, unpublished data). As a result, the total number of macroscopic fungal networks analysed so far is relatively low. Nevertheless, graph-theoretic analysis of the digitised networks has provided evidence that these indeterminate, de-centralized systems yield adaptive networks with both high-transport capacity and robustness to damage, but at a relatively low cost, through a Darwinian process of selective reinforcement of key transport pathways and recycling of redundant routes (Bebber et al., 2007). Furthermore, fungal networks dynamically modify link strengths and local connectivity when subject to experimental attack to readjust the balance between transport capacity, robustness to damage and resource allocation, resulting in increased resilience as the environment becomes more challenging (Boddy et al., 2010). Existing image processing approaches for semi-automated fungal network extraction work best if the mycelium is grown on a suitable substrate, such as cellophane (Ritz et al., 1996) or nitrocellulose (Barry and Williams, 2011) to facilitate acquisition of high-contrast images. The subsequent processing procedure typically involves noise reduction, using median or band-pass filters, intensity-based thresholding, with manual or automated threshold selection, followed by morphological processing to achieve an improved skeletonisation of the colony structure e.g. Barry and Williams (2011). However, these methods require rather specific culture conditions and illumination regimes that cannot readily be adapted for macroscopic networks grown on soil-based media, which is essential to explore the natural behavioural capability of these organisms and the range of their responses. In addition, the branches (termed cords in these systems) are multi-hyphal aggregates and show considerable variation in contrast as the individual cords span several orders of magnitude in diameter from m to mm, within a single colony. Furthermore, the background from the compressed soil substrate is not homogeneous in texture or reflectivity. Under these conditions, automatic characterization of the fungal network is challenging for conventional image processing approaches and more advanced solutions are required. We have developed a high-throughput automated image analysis approach to detect and characterize large complex fungal networks, grown under realistic conditions. We use watershed segmentation to extract the network rapidly, with prior curvilinear feature enhancement to highlight the network structure in noisy or low contrast images. The proposed enhancement method uses a novel brightness- and contrast-invariant approach based on the concept of local phase, called the Phase Congruency Tensor (PCT). The resulting skeleton is then pruned on the basis of local cost functions that incorporate both intensity and tensor direction information. A graph representation of the network is constructed from the pruned skeleton to give the network topology or a weighted network that includes edge lengths. To obtain an estimate of the edge thickness, which is important in predicting the physiological performance of the system, the skeleton is dilated to capture the approximate dimension at each point, and this value is used to set a local sampling region to estimate the reflected cord intensity for comparison with empirically derived calibration curves. This approach ensures that even sub-resolution branches can be analysed correctly. An addition procedure is applied to extract the resource region which is then used to redefine a corresponding part of the graph. Although we apply the method here to fungal networks, it has wide applicability to planar networks in any domain. EXTRACTION OF NETWORK-LIKE STRUCTURES Recent reviews of image processing methods for network-like structure extraction can be found in Dehkordi et al. (2011). According to this review, most image analysis approaches to extract network-like structures are classified into pattern recognition, model-based, tracking, artificial intelligence and neural networks. To set the current work in context, we outline the basic method and the limitations of each approach below. Pattern recognition Many of these methods apply thresholding and object connectivity, followed by a thinning procedure and extraction based on graph description (Kawata et al., 1995). Alternatively, the image intensity is represented as a hypersurface and the curvilinear features are extracted using curvature, with the crest points of the hyper-surface correspond to the center lines of the structure (Prinet et al., 1996). These centerlines can then be used as seeds in a region growing procedure that segments the full curvilinear structure from the images (Higgins et al., 1989). Direct segmentation of curvilinear structures can also be done by multi-resolution segmentation approaches. After segmenting salient structures at low resolution, less prominent network branches can be extracted in the neighborhood of the previously segmented regions at higher resolution (Chwialkowski et al., 1996). Curvilinear structures can also be segmented by mathematical morphology operators using specific structuring elements (Zana and Klein, 1997). Model based Model-based approaches apply explicit curvilinear structure models, or templates, to extract the network. In network extraction applications, the template is usually represented in the form of a series of nodes connected in segments and deformed to fit the image (Summers and Bhalerao, 1995). In curvilinear or tubular object segmentation, objects are commonly described as a set of overlapping template ellipses (2D) or ellipsoids (3D) (Krissian et al., 2000). More general deformable model-based techniques find object contours using parametric curves that deform under the influence of internal and external forces (Klein et al., 1994). Tracking Tracking-based approaches apply local operators on a specific point known to be in an object of interest and track it. Tracking approaches, starting from an initial point, detect network centerlines or boundaries by analyzing the pixels parallel and orthogonal to the tracking direction (Haris et al., 1997; Meijering et al., 2004). Artificial intelligence and neural network-based Artificial Intelligence-based approaches use knowledge to guide the segmentation process and to delineate network structures, e.g. using a general curvilinear model (Smets et al., 1988). A disadvantage of learning methods is that they depend on the quality of the training dataset, which requires significant work and can be susceptible to problems like over-training (Nekovei and Sun, 1995). ENHANCEMENT OF CURVE-LIKE STRUCTURES When the network structure is affected by variations of intensity contrast within the image, image enhancement operators can highlight curvilinear features and reduce background effects. Although other alternatives are available, here we concentrate on methods based on tensors. Intensity-based tensors A tensor representation of an image can give information about how much the image varies along and orthogonal to the dominant orientations within a certain neighborhood (Westin et al., 2001). In particular, the Hessian matrix (a secondorder tensor) has been proposed to describe local shape orientation for elongated structures (Sato et al., 1998). Since curvilinear structures observed in the image can appear in different sizes, scale space representations have been extensively used. This is commonly achieved through the use of a family of Gaussian filters, or their derivatives, with multiple scales achieved by varying the value of the variance. One of the most popular Hessian-based approaches to enhance curvi-linear structures is known as vesselness (Frangi et al., 1998). For a given image I(p) and scale , vesselness is defined as where p x; y T is the spatial location, ; 1; ; 2 are the eigenvalues of H p (the Hessian matrix at scale ) satisfying j ; ij j ; i1j; and c are thresholds that control the sensitivity of the line measurement. Multi-scale vesselness, for a given set of scales f g, can be computed as the maximum of the vesselness values calculated at different scales, and the eigenvectors at the same scale can be used to define local orientation (Obara et al., 2012). An alternative measurement of piecewise linear segments that is also based on the Hessian matrix, called neuriteness, has been proposed in (Meijering et al., 2004). MATERIALS AND METHODS Input images The proposed approach has been tested and validated on images of fungal networks of Phanerochaete velutina or Phallus impudicus that were grown for 3950 days in the dark on compressed soil or sand microcosms and photographed at intervals, as described previously, see figure 1a (Bebber et al., 2007; Boddy et al., 2010). Phase congruency-based tensors Herein, we introduce a novel phase concept for curvilinear feature enhancement called the PCT. In general, local tensor-based representations can be produced by combining the outputs from polar separable quadrature filters, applied on several orientations (Knutsson, 1989). Although tensor representations can be built on purely intensity-based filters (like the Hessian), these have the downside of being sensitive to changes in image contrast. Methods based on local phase have been proposed as a contrast-independent alternative for feature detection. In particular, phase congruency (Kovesi, 2000) is based on the concept that salient features have similar values of local phase when observed at different scales. We exploit the idea that phase congruency values are high in the direction perpendicular to the structure, while they remain close to zero in the direction parallel to the structure. More importantly, the values of phase congruency are minimally affected by contrast changes. Thus, for a given set of scales fs}, a set of orientations fo} and a given set of phase congruency measures PCop (for each orientation o) (Kovesi, 2000), the PCT takes the following form (Obara et al., 2012): where no is the normalized orientation vector in the direction o, 1=m 1, with m being the dimensionality of the image and I is the identity tensor. The local structure can then be described by eigenvectors and eigenvalues of TPC. 4.2.1 PCT vesselness and PCT neuriteness As described in Section 3.2, piecewise curvilinear segments can be detected by analyzing the relations between eigenvalues and eigenvectors of the locally calculated Hessian. In a similar way, the dominant orientation of the surface representing a curvilinear structure is given by the dominant eigenvector of TPC, i.e. the eigenvector corresponding to the eigenvalue of largest magnitude. PCT-based vesselness and PCT-based neuriteness are calculated using the original vesselness and neuriteness equations, with the eigenvalues of TPC substituting those of the Hessian. General approach for network extraction 4.3.1 Curvilinear feature enhancement Let us consider an image of a curvilinear network structure I(p) of a size M N pixels, where p x; y T represents pixel location. When the curvilinear network is not affected by variations of intensity contrast within the image, then the input image I(p) is directly used in the following procedures (see fig. 1a). Otherwise, one of the curvilinear feature enhancement methods, discussed in (Section 4.2), is applied. An enhanced image IFp as well as the corresponding vector field IVFp are calculated. The size of image IVFp is M N 2. 4.3.2 Watershed The watershed transformation was introduced as a tool for segmenting gray-scale images by (Beucher and Lantuejoul, 1979), and is now used as in many segmentation procedures using efficient algorithms (Meyer, 1994). Lets us consider a gray-scale image I(p) (or IFp) as a topographic surface: the gray level of a pixel becomes the elevation of a point, the basins and valleys of the relief correspond to the dark areas, whereas the peaks and crest lines correspond to the bright areas. The watershed line may be intuitively introduced as the set of points where a drop of water may flow down towards several catchment basins of the relief. The watershed algorithm can be used for segmentation by using the calculated catchment basins as segmented regions. In our case, however, we are interested in watershed lines themselves, which delineate ridges in the image and thus correspond to network branches. 4.3.3 Branches and branching points extraction The skeleton defined by watershed lines is subsequently categorised into branches IBp, branching points IBPp and end points IBEp using the topological method proposed by Kong and Rosenfeld (1996). IBp; IBPp and IBEp are images containing unique labels for the pixels that make up a specific object. The performance of this procedure is demonstrated in figure 1b. The majority of the skeleton corresponds visually to the clearest structures of the network. However, the watershed transformation typically produces an over-segmentation of the input image I(p) [or IFp] as extremely subtle variations in image intensity can still yield basins with corresponding watershed boundaries. These extraneous segments can be evaluated using local cost-functions and pruned if they fail to meet the appropriate criteria. Algorithm 1 Pruning procedure. Input: IFp; IVFp; IBp; IBPp; IBEp; , t Output: IBp; IBPp; IBEp NB maxIBp for i 1 ! NB do b findIBp i c CIF; IVF; b; if c5t then IBb 0 IBEbe 0 end if end for IBp > 0 [ IBPp > 0 [ IBEp > 0 BranchesISKp BranchPointsISKp BranchEndsISKp The pruning procedure used here is based on a cost function that includes weighted contributions of both intensity information and the alignment of the local vector field, calculated for every branch of the skeleton. The cost of the branch defined by pixels b is computed using the following formula: where CI is a normalized input or enhanced image intensity-based cost defined by and CV is a vector field-based cost that is calculated as follows: pffi1ffiffiffiffijffijffiffiffiffiffiffiffiffiffiffibffiffiiffiffiffiffi1ffiffi;ffiffibffiffiiffiffiffijffiffij bi; bi1 jIVFbi dbi; bi1j bi=jjbi1 and K is a number of pixels in b. 2 0; 1 determines the relative weight of the CI and CV cost components. |||| is the Euclidean norm. The purpose of the vector field cost CVF is to determine branches aligned with the direction of a vector field given by IVFp. Finally, all branches with a cost value smaller than threshold t and their corresponding branching and end points are removed (see fig. 1c). The procedure is given in Algorithm 1. NB denotes number of labeled branches in image IBp. Functions BranchPoints, BranchEnds and Branches return images containing labeled pixels representing every branching points, ends and branches. Function find returns non-zero pixel positions of an image. 4.3.5 Graph extraction In order to apply network analysis tools, the morphological structure has to be translated into an appropriate Fig. 1. Workflow: (a) Input image as defined by the region of interest outlined in red in Figure 2; (b) watershed-based segmentation: fungal network skeleton (in red), branching points (in green) and endpoints (in blue); (c) pruning (as pointed out by violet arrows); (d) graph built form the extracted structure of the network; (e) pseudo-color coded representation of cord thickness; (f) fungal food source detection (region of interest outlined in green in Figure 2: the boundary is outlined in magenta and (g) graph at the region of the source of food network representation. We assume that the fungal network can be represented as a graph by classifying junctions (branch-points and anastomoses) as nodes and the hyphae or cords between nodes as links (Bebber et al., 2007; Fricker et al., 2007). The graph representation of the network structure is constructed by matching every branch with its two corresponding branching points or with its branching point and endpoint, see figure 1d. The input branches, branching and end points images, IBp, IBPp and IBEp, respectively, are used to construct the graph G. For every branch b in IBp, a set of corresponding branching points fbp} is extracted. When the branch has two branching points then the graphs edge G(i) is defined by fbp(1), bp(2)}, otherwise, G(i) is represented by branching point bp(1) and the end point be. The procedure is given in Algorithm 2. Algorithm 2 Graph extraction procedure. NB maxIBp for i 1 ! NB do b findIBp i bp findIBPp > 0 \ IBp i if dimbp 2 then Gi fbp1; bp2g else if dimbp 1 then be findIBEp > 0 \ IBp i Gi fbp1; beg end if end for 4.3.6 Thickness A very important measure used to characterize biological networks, is the cord thickness, as this has a significant impact on the physiological function of the network. The thickness value for each branch is determined by iterative filtering of the center line image using morphological dilation and reconstruction procedures (Serra, 1982). The input branches and branching points and ends images, IBp; IBPp and IBEp, respectively, are used to construct the skeleton images ISKp. A copy of ISKp, called IWp, is then dilated using a disk structuring element S with 1 pixel radius. The difference between IWp and the result of dilation is stored in the image ISp that contains the surrounding boundaries of the skeleton. Then, for a given threshold value k, thresholding of the image IFp masked by the image ISp is performed and the image IVp is obtained. In order to reconstruct all pixels in the image IVp which are connected to the skeleton in the image IWp, the morphological reconstruction of the image marker IWp under the image mask IMp, a union of ISKp and IVp, is performed producing IRp. The points in IRp are added onto IWp, and the procedure is repeated. The procedure is repeated until no pixels remain to be reconstructed with intensity IFp above k. The distance map IRDp is calculated as the shortest distance from each of the skeleton points ISKp to the background in the reconstructed image IWp. Finally, the thickness map image ITKp is given by the distance map image IRDp weighted by the mean intensity of IFp points defined by the neighborhood of each pixel of the skeleton image ISKp. The output of the procedure is presented in figure 1e and its workflow is shown in Algorithm 3. In this Algorithm, function SIp denotes the morphological dilation of an image I(p) by a structuring element S. SIp; Jp is a function that represents the morphological reconstruction of mask image I(p) from marker image J(p), Jp Ip. Function dist computes the Euclidean distance transform of a binary image. denotes the convolution operation and S(r) is a disk structuring element with r pixel radius. Algorithm 3 Thickness procedure. Application to extraction of fungal networks In all the fungal network extraction experiments, the following parameters have been used: number of scales and orientations 6, weight 0:5, branch cost threshold t 0.2 and thickness threshold k 0.4. In all cases, the PCT-based vesselness approach (Equation 1) is used. An analysis of the sensitivity to different parameters for the PCT calculation can be founded in (Obara et al., 2012). In the particular case of the fungal networks, the graph representation of the network has to be redefined in the region corresponding to the fungal source of food, as the structure of the network within the food source is not visible, but can be approximated as a single node connected to all cords arising from the boundary. 4.4.1 Fungal food source detection In order to extract the fungal food source, the image I(p) is segmented by a global threshold (Otsu, 1979). The segmented image is then filtered using morphological opening and closing using a disk structuring element S, to produce the binary image IFoodp. figure 1f shows the output of this procedure. 4.4.2 Graph for the fungal network The graph G representing the fungal network is finally modified to account for the presence of food resources. The graph is redefined for the vertices that intersect with the segmented resource IFoodp. such that they all connect to a single node in the center of the source of food as the underlying network is not visible within the wood block (fig. 1g). Species such as P. impudicus produce densely cross-linked networks, with relatively well-defined thick cords. To summarise the complete extraction process, images of this species collected at high-resolution and high-contrast (fig. 1a) were rapidly segmented using the watershed method (fig. 1b), with no additional requirement for curvi-linear feature enhancement. Oversegmentation, generally considered as a downside of the watershed transform, allows in this case a rapid delineation of highdensity networks. Pruning based on the average intensity of every branch was then carried out, providing a mechanism to control the overall level of detail of the final network (figure 1c). The graph representation of the network was then calculated by defining the branch or fusion points and end points throughout the network as nodes, and then determining the edges that link them (fig. 1d). The Euclidean length of each edge was also calculated from the segmented skeleton, and the cord thickness estimated from a sequence of morphological dilations of the skeleton until the threshold criterion was reached, combined with filtered pixel intensities to capture the difference in reflectivity for cords of different dimensions (fig. 1e). The final step for the fungal networks was to identify the initial food resource (fig. 1f) and then redefine the graph in this region to connect all cords emanating from the resource with a central node (fig. 1g). The performance of the proposed approach applied to the entire image is presented in figure 2. When no enhancement was applied, the watershed and pruning method was less reliable for low-contrast, low-resolution images, such as those obtained for Phanerochaete velutina, which has thinner, more diffuse cords, particularly at the growing margin. With varying levels of threshold, the network was initially over-segmented, but rapidly became dis-connected as the threshold was raised, making the choice of threshold critical (fig. 3ae). In comparison, even these challenging networks were efficiently extracted following intensity-independent PCT enhancement (see Section 4 for details) to give a very clean network over a broad threshold range (fig. 3fj), that greatly facilitated subsequent pruning operations (fig. 4). To quantify the performance of our approach, five complex regions of interest in the fungal images were selected. Manual tracing of network centerlines was performed exhaustively by an expert, to give a Gold Standard (GS) reference. In total, in all regions of interest, 2256 network branches were annotated. The GS and extracted network were compared using the network network distance measure "d (Gelasca et al., 2009), defined as the average distance between each point on the GS network centerline and the corresponding closest point on the extracted network centerline, and vice versa. The standard deviation "d of the network-network distance measures was also calculated. The distance error evaluation for each ROI, where the network was extracted using the proposed approach based on PCT vesselness, is presented in Table 1. For all annotated regions, the average distance error was "d 0:7393pixel . The recommended setting used to calculate the PCT-based vesselness, can be found in (Obara et al., 2012). Supplementary figure S4 presents normalized histograms of distance error values "d for all analyzed regions of interest in fungal images. Additionally, traced lines were identified as true/false Positives depending on whether a line was found in the gold standard at a distance smaller than two pixels. Precision and recall were calculated for a set of comparable approaches, and are shown in Table 2. Finally, a comparison of the network thickness estimated using optical microscope and our approach, for 10 randomly selected measurement points, is presented in Supplementary figures S5 and S6. Runtime performance The average runtime spent to extract the fungal networks from the grayscale images of 1000 1000 pixels and 60 K branches (links) is 70s (without enhancement) and 200s (with PCT enhancement). This compares with manual segmentation of 3K branches in about 23 days (Fricker et al., 2007). The runtime performance of the implemented methods was tested on a PC with Intel Core 2 Duo (T8300) system with 2GB memory, running Linux and MATLAB R2009b. In this article, an image analysis and processing concept for general curvilinear network detection has been introduced. The core method is based on extraction of the network centerline from a watershed segmentation, followed by a pruning procedure. This approach can be directly applied to high-contrast images, but the combination with a curvilinear feature enhancement based on Phase Congruency Tensors provides a contrast invariant solution for more challenging networks. Using the watershed transform to calculate center lines provides several advantages over other more complex networkanalysis routines. First, it guarantees calculation of connected networks. Second, it requires low computation times in comparison with other segmentation methods. By using a standard morphological watershed transformation on raw or PCTenhanced images, we obtain an oversegmented skeleton image which encompasses the entire biological network, but also a contribution of false edges associated with subtle variations in background intensity that still yield basins in the watershed image. These extraneous features can be removed post-segmentation by judicious selection of cost-function weights or reduced by prior network enhancement and noise-reduction using PCT or other image enhancement techniques described in Section 4. Visual inspection of the results (fig. 2 and Supplementary Figs. S13) confirm the robustness of the proposed approach to extract networks from complex and challenging biological specimens. In particular, the PCT-based enhancement is effective at dealing with highly variable intensity levels of the curvilinear features and is capable of providing high detection responses on low contrast edges against a noisy background (fig. 3). These properties are essential to detect structures in low contrast regions of noisy images that are common in a large number of biomedical images (Meijering et al., 2004). The results demonstrate that the proposed approach provides a fast and robust solution to detect and extract a graph representation of complex curvilinear networks and overcomes a critical bottleneck in biological network analysis. This now permits high-throughput measurements with improved resolution and precision. The approach is generic and can be applied to a wide range of biomedical images. One of the major contributors to errors is optimal threshold selection for the pruning procedure. Therefore, future work will be focused on data-driven optimizing approaches to determine optimal local or global threshold selections. Once the network has been extracted, a wide range of network parameters can be calculated (Fricker et al., 2007). As the data are embedded in Euclidean space, a number of basic morphological measures can also be readily derived. These values either have a straightforward biological meaning in their own right or they provide a comparison with network structures in other domains. In particular, the rich network structures extracted by this approach provide a rapid means to examine network features, 0.3 0.5 0.7 (h) (i) 0.1 Table 1. Distance error evaluation of the proposed method applied to five regions of interest in fungal images No. of branches "d [pixel] Intensity Vesselness Neuriteness PCT vesselness PCT Precision Recall 0.91 0.56 0.95 The assessment of the rate of false positive and false negative segments has been performed within an error diameter of the gold standard segments equal to 2[pixels]. such as transport efficiency, resilience, cost and control complexity. ACKNOWLEDGEMENT The authors thank to L. Boddy, D. ABear, J. Hynes and J. Wood for some of the fungal images analyzed. Funding: BBSRC, EPSRC and NERC, and an Academic Fellowship from Research Councils UK. Conflict of Interest: none declared.


This is a preview of a remote PDF: https://bioinformatics.oxfordjournals.org/content/28/18/2374.full.pdf

Boguslaw Obara, Vicente Grau, Mark D. Fricker. A bioimage informatics approach to automatically extract complex fungal networks, Bioinformatics, 2012, 2374-2381, DOI: 10.1093/bioinformatics/bts364