BIOLOGICAL INVASION IN SOIL: COMPLEX NETWORK ANALYSIS F Perez-Reche", S. N Taraskini", FM Neri$, C.A. Gilligan", L. da F Costa', M P Viana§, W Otten', D. Grinev' £Department of Chemistry, University of Cambridge, Cambridge, UK tSt. Catharine's College, University of Cambridge, Cambridge, UK $Department of Plant Sciences, University of Cambridge, Cambridge, UK § Instituto de Fisica de Sao Carlos, Universidade de Sao Paulo, Sao Carlos, SP, Brazil +SIMBIOS, University of Abertay Dundee, Dundee UK ABSTRACT A network model for soil pore space is developed and applied to the analysis of biological invasion of microorganisms in soil. The model was parameterized for two soil samples with different compaction (loosely and densely packed) from images derived from an X-ray micro-tomography system. The data were then processed using 3-D imaging techniques, to construct the networks ofpore structures with in the soil samples. The network structure is characterized by the measurement of features that are relevant for biological colonization through soil. These include the distribution of channel lengths, node coordination numbers, location and size of channel bottlenecks, and the topology of the largest connected cluster. The pore-space networks are then used to investigate the spread of a microorganism through soil, in which the transmissibility between pores is defined as a function of the channel characteristics. The same spreading process is investigated in artificially constructed homogeneous networks with the same average properties as the original ones. The comparison shows that the extent of invasion is lower in the original networks than in the homogeneous ones: this proves that inherent heterogeneity and correlations contribute to the resilience of the system to biological invasion. Index Terms- Soil image analysis, complex network, biological invasion, percolation. 1. INTRODUCTION The structure of soil and its transport properties are of significant interest for various physical, geological, biological and agricultural reasons [1, 2, 3, 4, 5]. Usually soil is treated as a medium of pores of various sizes and shapes. It is a challenging experimental task to obtain reliable information about FPR, SNT, FMN and CAG thank BBSRC for funding (Grant No. RG46853), L. da F. Costa thanks FAPESP (05/00587-5) and CNPq (301303/06- 1) for sponsorship and M. P. Viana thanks FAPESP (07/508829) for financial support. WO and DG thank the Scottish Alliance for Geoscience, Environment and Society (SAGES) for support. 978-1-4244-3298-1/09/$25.00 ©2009 IEEE soil structure and, in particular, about the spatial arrangements of pores. Several techniques have been used for soil structure analysis including serial sectioning [6], laser scanning confocal microscopy [7], X-ray tomography [8] and nuclear magnetic resonance imaging [9]. The structure of soil is obtained in the form of an image, which can then be processed by image analysis techniques and used for construction of a network model representing the soil structure [10, 11, 12]. The network model can then be used for studying transport properties of soil, e.g. water flow through the soil (see e.g. [13, 2]). Nowadays commercially available X-ray computerized tomography systems (also known as CT scanners)are capable of resolving micron-sized pores in undisturbed soil samples in three dimensions. Here we introduce the results from a X-ray computerized micro-tomography system that achieves a resolution of 74JLm for contrasting samples of approximately 2 x 2 x 4 cm representing loosely and densely packed field soils. The resulting network structures are deconvoluted into a series probability distributions for channel lengths, node coordination numbers, location and size of channel bottlenecks, and the topology ofthe largest connected cluster defined. The invasion of a microorganism through the largest connected cluster is then modelled using a set of rules in which the transmissibility between each pair of connected nodes in the network is a function of channel properties linking the nodes. Under these rules, the transmission process is equivalent to an epidemiological process described by an SIR (susceptible-infected-removed) model. This, in tum, can be mapped onto isotropic bond percolation. According to this mapping, the system (i.e. the microorganism invading soil) exhibits a second-order phase transition from a non-invasive to an invasive regime. The model enables us to compare artificially constructed homogeneous networks with the same average transmissibilities as the original ones. We show that there are shifts in the invasion curves between the homogeneous and heterogeneous networks (as a function of transmissibility. Those shifts are due to the existence of local correlations in channel properties that reduce the extent of invasion in the DSP 2009 original networks , with greater effect of heterogeneity in the densely packed soil. The paper is organized as follows. In Sec. 2, the experimental technique used for soil sample analysis is described. Sec. 3 describes the details of the network models associated with the soil samples. Sec. 4 deals with biological invasion in the in the network models. Concluding remarks are presented in Sec. 5. 2. EXPERIMENTAL DETAILS Soil aggregates (1 - 2 mm) of an arable sandy loam were packed to attain bulk densities of 1.2 or 1.4 Mg/m', representing two characteristic examples of a loosely and densely packed field soil, respectively for this soil type. Full details can be found in Ref. [14]. We refer to these two soil samples as loosely (1.2 Mg/rn") and densely (1.4 Mg/rn") packed soils hereafter. We scanned these samples with a Metris X-TEK Benchtop micro-tomography system [4, 15], using a molybdenum target, X-ray source settings of 155 kV and 25 ILA, and an aluminum filter (0.25 mm) to reduce beam-hardening artefacts. 2-D radiographs were collected at 1169 angular positions and then reconstructed using a filtered back projection algorithm with a resolution of 74 ILm and isotropic voxel size). Each radiograph was averaged over 32 frames to improve signal-to-noise ratio and ring artefacts were minimized during data acquisition. Numerical corrections were also applied during the reconstruction to minimize remaining beamhardening artefacts. Both 3-D volumes were then imported into VGStudioMax v.I.2 .1 [16] and converted into 260 x 525 8-bit TIFF image stacks with voxel-thick slices. Binary data sets were created by thresholding the grey-scale image stacks in Image] [17]. The choice of the threshold parameter was based on the 3-D statistical analysis of the histogram region corresponding to the pore-solid interface. In order to obtain an optimal threshold value this analysis took into account variation of the grey-scale values in pores of different shapes and sizes. 3. NETWORK ANALYSIS OF THE SOIL STRUCTURE Examples of visualization of loosely and densely packed soil samples are provided in Figs. Ia and Ib, respectively. The filled space represents the pores in the soil. In order to construct the network representing the pore volume, all the pores were processed by a thinning algorithm [19,20,21], required in order to reduce each object (pore) to a respective l-voxel skeleton. The skeleton is a thin structure located at the most central parts of the respective original shape. The skeleton retains all the topological features of the original shape (e.g. branching structure and cycles). The skeletonized pores were then mapped onto a network as follows. The skeleton is a set of intersecting curves with some dead ends (where the skeleton terminates). Each intersection point and all dead ends were associated with the nodes of the network [21]. The pore space around the skeleton between two nodes is called a channel or link (edge) between two nodes. The axis of the channel thus coincides with the skeleton. Figs. 2 (a) and (b) show the pore networks obtained for the loosely and densely packed samples presented in Fig. I. r-\ \ \ I. Fig. 1. Images ofloosely (a) and densely (b) packed soil. The gray colour represents the pore space in soil. The insets show the middle slice (along the vertical axis) of the volumes with the pores shown in yellow. I I I <, ｟ ＮＭ ｾ I 3.1. Network construction After conversion into binary format, pore networks for both soil samples were reconstructed in 3-D by using the Visualization Toolkit (VTK) [18]. Isolated pores with volume smaller than 104 voxels were discarded to aid data visualization and analysis (below we use I pixel = 74JLm as a unit of length). This has negligible impact on the topology of the resulting network on which we perform our analysis . Fig. 2. The network models for the loosely (a) and densely (b) packed soil shown in Fig. I. The complete network (a) consists of Nnode = 10316 nodes and Nlink = 11502 links, while Nnode = 3612 and Nlink = 3846 for the complete network (b). Different colors correspond to distinct isolated clusters with the largest connected component in blue. 3.2. Network characteristics The networks constructed according to the algorithm described in Sec. 3.1 represent the pore (or equivalently soil) topology of two soil samples. The networks are embedded in 3D-space and each node is characterized by a position vector R i . The nodes are connected to each other in a complicated way and, in general, the network consists of isolated (not connected to each other) clusters of different sizes. The largest connected component (LCC) (the blue clusters in Fig. 2) which can be identified for each network is of particular interest for biological invasions as described below. In the case of loosely packed soil, the LCC percolates through the sample and contain ｎ［ｯｾ = 10183 nodes connected by Nkn1c = 11369 channels. The LCC for densely packed soil does not percolate through the sample volume and contains LCC LCC N node = 2613 nodes and N hnk = 2823 channels. One of the standard characteristics of network connectivity is the node degree (Z) distribution, p(Z), where Z is the number of links attached to a node (coordination number). These distributions are shown in Fig.3 (left lower inset) both for complete networks and LCe. The value of Z = 1 corresponds to the number of dead ends in the networks. By construction, there are no nodes with coordination number Z = 2. As seen from Fig. 3, the probability density functions reach a maximum value at Z = 3 and then quickly decay with increasing Z. The complete networks are found rather sparsely connected with mean degrees ((Z ) = 2Nlink I Nnodc), (Z l) ':::' 2.229 and (Z2) ':::' 2.130, for loosely and densely packed soil samples, respectively. The LCC exhibit similar values, (Z l) ':::' 2.233 and (Z2) ':::' 2.161. Therefore, the densely packed soil have smaller degree numbers (pores collapse under compaction). The other important structural characteristics of the network are the arc-length, L ij, of the channel between nodes i and j , and bottleneck diameter of the channel, ¢ j . The arclength is defined as the length of the axis of the channel. The bottleneck diameter, cPij, is defined as the minimal diameter of the maximally inscribed circle for crossection of the channel. The values of cPij were obtained by calculating the maximum distance transform value along the skeletons [22]. Both these characteristics play an important role for biological invasions (see Sec. 4). In the case of limited image resolution, the value of cP* can also be associated with the resolution, i.e. all the pores of sizes less than the resolution are not represented by the network model. The distribution of the arc-lengths for both samples are presented in Fig. 3. It can be seen from Fig. 3 that the numbers of short and very long channels are reduced in the densely packed sample as compared with the loosely packed one. The shape of the channels can be characterized by their relative arc-length as compared with the Euclidean distance between nodes, i.e. by the parameter be = 1- IRj -Ri IIL ij. The value of be describes the deviations of the channel axis 10 5 10-3 • 10"2 10"3 10·.1 40 10- I 2 J z 4 ..-... .. . ODJJ[J] •• [] [] «J : 5 L [] D mo .[] • 6 • •• ••••• 40 20 L 60 Fig. 3. Probability density function of the arc-lengths of channels (in pixels), p(L), in the LCC of the networks representing loosely (solid circles) and densely (open squares) packed soil. The upper right inset magnifies the region of small arclengths. The lower left inset shows the probability density function of the node degrees, p(Z) for networks corresponding to loosely (circles) and densely (squares) packed soil samples. Open and solid symbols are used to denote the complete network and the largest connected component in the network, respectively. Open and solid circles are not distinguishable at this scale. from the straight line connecting two nodes, so that be = 0 in the case of zero deviations. The distribution of be for the LCC in both networks is shown in Fig. 4(a). It is evident from Fig. 4(a) that the majority of channels (which are relatively short in length) do not deviate strongly from the shortest path between two nodes. The deviations (the value of be) increase with the arc-length of the channels . This can be seen in Fig. 4(b) which demonstrates correlations between values of be and the arc-length . The bottleneck diameter distribution is shown in Fig. 5 (a). It follows from this figure that the number of channels with large bottleneck diameters is reduced in the densely packed soil in comparison with the loosely packed one. In order to characterize the position of the bottleneck in a channel {ij} it is convenient to introduce a parameter, bij = IL ib - Lbj II L ij, where Lib and Lbj are the arc-lengths from node i to the bottleneck and from the bottleneck to the node j , respectively. This parameter reflects the relative distance from the bottleneck to the middle of the channel. If the bottleneck is in the middle of the channel then bij = 0 and if it is at one of the ends of the channel then bij = 1. The distribution of the values of parameter bij for both networks is shown in Fig. 5(b). The dip around b = 1 indicates that the bottlenecks mainly occur inside the channels but not close to the nodes. This indicates that the crossing points of several channels (nodes) are characterized by relatively large pore 10 ..."\ . 0.3 ｾ \.. ｉＩｑ (..Q <::» Q. 0.1 0.01 -_. 0. 1 se c c " C ｾ ｾ ｄ c 0. 1 ti•l' •. , C if ", • •'Is,!, CC - aD ｾＱ C- Ｐ ＭＲ " ｾ Ｎ , 40 I I .". amII "",' 10- • . • D D .' 3 ce I- ｾ ｾ ｾ 08 .ｾ · .. Ｎ • IJ Ｎｾ ｾ lJ • ｾ ｣ ｾ 20 . co n Ｇ ｾ ＮＺＭＢＬ ,-..., C ｾ cae •• . 110 _ ••• . . .. • ｜ ＬＮ c ,p ftf: 1Ib 10-1ｬｯＧ ｾ ｣ . ., ' (..Q Q) 0.2 (b) c C \....c c.lf, 2.0 100 0.4 (a) ... •• e c l} ....f"Ｎ ［ ｾ c Ｎ ｩＧ ｾ . ••• 60 L Fig. 4. (a) Probability density function, p(Oe ), of the relafor the LCC of the networks representing tive arc-length, loosely (solid circles) and densely (open squares) packed soil. (b) Dependence of the parameter on the arc-length of the channels, (in pixels). The symbols have the same mean ing as in (a). oe, oe space. The size of the bottleneck is determin ed by the mechanical prop erties of the soil and the topology of the chann el. Fig. 6 demon strates appar ent correlations between the value of the bottleneck diameter and arc-length of the channel , suggesting that cPij ex: Li/' with 0: 1 ｾ 0.41 and 0:2 ｾ 0.46 for loosely and densely packed soil, respe ctively. 4. BIOLOGICAL INVASION IN SOIL The following transport processes can occur in the soil: (i) the spread of liquids such as water; (ii)the spread of gas (e.g. oxygen) and (iii) dispersal of microorganisms. The first two phenom ena have many important applications and have attracted hug e attention from the scientific community (see e.g. [11] and references therein) The last phenomenon, of great importan ce for biological processes including the spread of pathog en, has attracted much less attention [23] and still exhibits several unsolv ed problems. In this section, we address the probl em of invasion ofmicroorg anisms through the system of soil pores represented by a network (see Sec. 3). A question to answer is the following. Given a network of pores and dynamical rules for dispersal, can soil be invaded (occupied) by microorganisms and, if so, what are the conditions for such invasion? We approach this problem using standard methods of statistical physic s for description of non-equilibrium phenomena in complex systems. 4.1. Theoretical background We consider invasion of a microorganism in structured soil, typified by hyphal (filamentous) growth through soil pores of 10- 0 10 20 30 40 50 60 70 <P 0 . "ｾ '" ｾ 'i C [J c DC 0"*. ｉｊ ｾ ... 't 0.5 'i 114 C IJ 0.0 0 0.2 0.4 0.6 0.8 I 8 Fig. 5. Probability density function of (a) bottleneck diameters (in pixels), p(cPij ), and (b) relative distances to the bottleneck, p(0) in the LCC of the networks representing loosely (sol id circles) and densely (open squares) packed soil. a fungal colony expanding from an initial site of introduction with a fixed source of nutrients at the site of introduction. The dynamical rules are as follows: (i) We place an initial source of microorganisms on an arbitrary node i (as defined in Sec. 3) of the LCC of the network representing the pore space in the soil. (ii) The organi sm can spread stochastically along any of Z, channels linked to nod e i and thus , with given probability, Pij , can reach node ( = 1, ... , Zi ), so that node becom es coloni zed (occupied) by the microorganism. (iii) From node the process of invasion continues in a similar mann er with the only exception that the organisms cannot move back to previously colonized nodes. The motion of the organism s along a channe l is assumed to be a Poisson process, so that the probability, dPij (t), to reach node for organisms moving from node i to node along the channel of length L ij in an infinitesimall y small time interval betwe en t and t + dt is [24], dPij (t ) = (1 Pij(t ))f3ijdt , with straightforward solution , Pij(t) = 1 e -{3i t where the coefficient f3ij is the invasion rate along the channel {ij} . The invasion rate can be estimated as f3i j = Vij / ij , with Vij being the typical velocity of motion of organisms through the channel {ij} and ij being the lengths (arc-length as defined in Sec. 3) of the channel betwe en nodes i and Accordin g to the equation above the probability of invasion of node i . dPij (t ), approaches unity when t ----; 00 . This leads to a scenario in which the whole LCC network is invaded in the long-tim e limit, something that is unlikely to happen in practice when nutrient is limited for the invading microorganism. Accordingly, we introduce a time Tij , available for exploration of each {ij} channel , i.e. the microorganism can only move along the chann el for t :::; Tij ' The finite exploration tim e may be a consequence of a limited amount of nutrients inside the channe l or a limited life-tim e of the microorgan isms. Rapid invasion is often of great importance into the chann el. If the bottleneck diameter is greater than the bound value, ¢*, then the organisms can fit into the channel, so that invasion can occur with certain given probability along this channel. Therefore, the transmissibility of channel {ij} is defined as if if • Fig. 6. Dependence of the bottleneck diameter (in pixels) vs arc-length of the channel (in pixels) for loosely (solid circles and solid line) and densely (open squares and dashed curve) packed soil, respectively. The straight lines represent the linear fit by regression. as new niches are colonized predominantly by rapid invaders. Here we assume a local clock for colonisation of each successive pore, which amounts to an assumption of a limited and equal amount of nutrient available for microorganism growth within each pore. Under these assumptions, the probability of invasion of node j from node i in the long-time limit (called transmissibility) is T 'J. . -- 1 _ e - Vi j Ti j / ij . (1) All the parameters affecting the transmissibility T i j depend on chann el characteristics and vary from channel to channel due to their dependence on e.g. channel topology. The parameters Vij and T i j can also depend on the node characteristics (such as amount of nutrient available at the node and the volume of the pore associated with the node). This can result in mutual dependence of transmissibilities of the channels originating from the same node. For simplicity, we assume that the values of Vij and T ij are the same for all channels with the same local invasion scale, k = Vij Ti j being a control parameter in further analysis. The local invasion scale has the meaning of a typical channel length such that the channels with L i j > k are more likely to be closed for invasion while the channels with L i j < k are more likely to be opened. The bottleneck size plays the role of a geometrical cutoff for the channel, so that all the channels with bottleneck sizes less than a critical one, ¢*, are deterministically closed for biological invasion. The value of ¢* can depend on the type of microorganisms and soil properties. It is reasonable to assume, that if the minimal cross-section size along the channel (the diameter of the bottlen eck or throat in the channel), ¢'I/' , is small enough, ¢ijin < ¢ *, then the organisms do not fit ¢ijin > ¢* ¢ijin < ¢* (2) Under these rules, the invasion process is identical to an epidemic process described by an SIR (susceptible-infectedremoved) model which can be mapped onto isotropic percolation [25] by identifying the channel transmissibility with the bond probability. In this mapping, the system (organisms invading soil) exhibits a second-order phase transition from a non-invasive to an invasive regime with an increase of the control parameter, k. The order parameter, Any, invasion probability in infinite system (probability that the initially infected site belongs to the spanning cluster), changes from zero value Any = 0 for k ｾ k c to finite values for k > kc, where k c is the critical value of the local invasion scale. In bond-percolation, the bond probability plays the role of the control parameter [26]. For the sake of comparison with bond percolation, it is convenient to introduce an alternative control parameter for invasion of soil, namely the mean transmissibility, (T) = - e -k / i j ) p ( )d where i j) is an experimentally known probability density function of the channel lengths and averaging is taken only through the channels with ¢ijin > ¢* which are accessible for invasion. Invasion of soil is inherently heterogeneous because the channel transmissibilities and coordination numbers vary through the system. However, if the channel transmissibilities are independent random variables for each channel, then spread on a heterogeneous network characterized by the distribution of T i j is equivalent to spread on the same network with a homogeneous transmissibility, Ti j = (T) [27,28,29], meaning that the function Any ((T) ) is the same in the systems with heterogeneous and homogeneous (meanfield) transmissibilities. One of the aims of this study is to verify this property for the spreading process in soil. The other aim of our analysis is to find how the order parameter depends on both control parameters (k and (T) ) and thus identify the conditions for invasion of soil by microorganisms. 4.2. Results We have investigated spread on two networks describing invasion of loosely and densely packed soil samples . In both cases, the analysis has been performed on the LCC of the network (see Sec. 3). The topology of the largest connected component obviously depends on the value of ¢* , because all the channel s with the bottleneck size less than ¢* are closed to spread. First, we present the data obtained for the networks with all the channels included, i.e. with ¢* = O. Fig. 7(a) shows the dependence of the invasion probability versus the local invasion scale for loosely (solid circles) and densely (solid squares) packed soil samples. The value of Anv increases with increasing k in such a way that Anv ':::' 0 for k ;S ke and Anv > 0 for k ｾ k e , where ke is th e critical value of the local invasion scale (the absence of a sharp transition at k = k e is due to finite-size effects). For loosely packed soil , k; ':::' 15, while k e ':::' 35 for densely packed soil reflecting the shift of the whole curve for Anv(k) to the higher values of This means that the densely packed soil is more resilient with respect to invasion due to the fact that the pore sizes under compression decrease. Therefore, the local invasion scale should be larger in order to activate more channels for the same level of invasion as in the loosely packed soil. The open symbols in Fig. 7(a) represent the invasion probability for topologically the same networks but with homogeneous transmissibility, T = (T) It is clear from Fig . 7(a) that the homogeneous (mean-field in transmissibility) networks are less resilient to invasion (the open symbols are above the solid ones for the same value of This effect is much more pronounced for the densely packed network. Therefore, heterogeneity in transmissibilities makes the networks less acceptable for invasion. A similar effect has been observed in lattice models for spread of epidemics [30]. 0 <1>* = > 0.6 <: ｾＭ 0.4 Fig. 7(b) shows how the invasion probability Anv changes with the mean transmissibility, ( which is analogous to bond probability for a percolation problem. For small transmissibilities (bond probabilities), the invasion probability (the mass of the spanning cluster) is close to zero and it increases with (T) The point where it starts significantly to depart from zero is called the invasion threshold or critical transmissibility, T e . For the percolation problem this value is known as a percolation threshold. It should be noticed that the valu es of the critical transmissibility (percolation threshold) are relatively high (T e ':::' 0.81 for loosely packed soil and T; ':::' 0.96 for the densely packed one) , as compared with T e ':::' 0.25 for bond percolation in simple cubic lattice and T; = 0.5 for bond percolation in square lattice [31], and approach T; = 1 for bond percolation in a linear chain. This is a consequence of the low node coordination in both networks with the mean coordination number ( ':::' 2.23 for loosely and ( ':::' 2.16 for densely packed soil (see Fig. 3). Fig. 7(b) also presents the data (open symbols) for Anv(k) obtained for the effective (mean-field) networks with homogeneous transmissibility (T) If the transmissibilities T are independent random variables then the invasion probability should be the same both in heterogeneous and homogeneous mean-field networks characterized by the mean value of the transmissibility [27, 28, 32, 33] , i.e. the open and solid symbol curves in Fig. 7(b) should coincide. However, as seen from Fig. 7(b), this is not the case for both networks and thus the correlations in transmissibilities playa significant role for both networks. They reduce the transmissibility for a given value of (T) and thus make the real networks more resilient to invasion as compared to the homogeneous ones (the opensymbol curves are above the solid-symbol ones in Fig. 7(b)). 0.2 O. . .Ｎ Ｌ ｾ Ｚ［Ｍ Ｇ Ｍ ｦＬ［ ＺＭＧ ｾ ＭＧ Ｍ［ ｫ 100 ＭＧ ｾ 120 ＭＧ Ｌ ｾ 140 0.8 0.8 f ｾＭ f (b) ｾＭ 0.4 .......... ......... 4 6 8 4-<l 10 0.4 0.8 0.6 0.4 0.2 0.2 0.6 60 80 100 120 140 0 0 k 0.2 Ｐ ｴＮ e-c-e $,=0 06 Ｔｉｩ ｾ ＭＧｉ ｴＭＮ ｾ､ ｾ 0.8 <T> Fig. 7. Invasion probability Anv versus (a) local invasion scale k (in pixels) and (b) mean transmissibility (T) in two networks representing loosely (circles) and densely (squares) soil samples. The solid symbols are used for heterogeneous networks with transmissibilities calculated according to Eq. (2), while the open symbols refer to the effective networks with homogeneous transmissibility equal to the mean value, ( The networks with all the channels (i.e, ¢* = 0) were used for the analysis. ｾ ｾＭ 0.8 06 0.6 0.4 0.4 0.2 0.2 (d) 0.9 T Fig. 8. Invasion probability Anv versus local invasion scal e k (in pixels) (a and b) and transmissibility T (c and d) for loosely (a and c) and densely (b and d) packed soil for different values of ¢* (in pixels) as indicated in the figures. The invasion probability, Finv, has been calculated numerically for several networks corresponding to different values of ¢* as a function of the local invasion scale, k, and mean transmissibility, (T). The results of the analysis are shown in Fig. 8 for several values of ¢*. The LCC in the network corresponding to the loosely packed soil becomes significantly sparser with ｎ［ｯｾ == 6903 and (Z) ｾ 2.148 == 2526 and (Z) ｾ 2.089 for for o; == 4 pixels and ｎ［ｯｾ ¢* == 10 pixels. The LCC in the network for the densely == 2613 packed soil becomes even more sparse with ｎ［ｯｾ and (Z) ｾ 2.161 for o; == 1 pixel and ｎ［ｯｾ == 1983 and (Z) ｾ 2.134 for o; == 2 pixels so that strong finite-size effects influence the invasion probability for greater values of ¢* (not shown). It follows from Fig. 8(a,b) that the removal of narrow channels from soil reduces significantly the invasion probability and thus results in the shift of the critical value of k; to greater values (cf. curves in (a) and (b) for different values of ¢*) and makes the soil more resilient to biological invasion. The reduction in the channel number for the LCC with finite value of ¢* results in a decrease of the mean coordination number and thus in an increase of the critical transmissibility for networks without narrow channels. This can be clearly seen in Fig. 8(c,d) where the probability of invasion curves are shifted to larger values of transmissibilities with increasing value of ¢*. In the analysis above, several assumptions for the invasion process have been made. In particular, we assumed that the invasion rate is inversely proportional to the arc-length of the channels and depends on the channel area in a step-like fashion. Both these assumptions are rather simplistic and proper analysis based on the experimental data should be undertaken (which is the aim of our future work). It would also be desirable to support our findings by analysis of several (rather than two) soil samples using different image resolutions. 5. CONCLUSIONS To conclude, we have presented a network model based on real data for soil structure and used this model to analyze a biological invasion in soil. The main idea of the network model is to reduce the pore space in the soil to a skeleton (network) which is topologically equivalent to the original pore space. The resulting networks corresponding to loosely and densely packed soil exhibit the following topological properties: (i) the mean node coordination number is relatively small and close to two; (ii) the network channels (links between nodes) are distributed in length and deviate from the straight line connecting two nodes with the deviation being proportional to the length of the channel; (iii) the channel bottlenecks are mainly located inside the interior part of the channels rather than close to the nodes; (iv) the size of the channel bottleneck decreases with increasing length of the channel. An invasion of microorganisms is then defined on the largest connected component of the network as a stochastic process on the network with transmission probabilities depending on the properties of the network. Within this approach, biological invasion is a critical phenomenon which either can or cannot occur in the system depending on the values of the control parameter (local invasion length scale). We have found that the networks associated with the pore space in soil are significantly heterogeneous both in topology (e.g. node connectivity and their spatial arrangement) and transmissibility (invasion probability through the channel). The heterogeneity in transmissibilities makes the network more resilient to invasion in comparison with the same network all the channels of which are characterized by the same (mean) transmissibility. The topological heterogeneity is shown to result in correlation effects between transmissibilities of different channels and thus leads to further increase in resilience. It has also been demonstrated that the compaction of the soil changes the soil structure in such a way that the densely packed soil is more resilient to biological invasion as compared to the loosely packed soil. 6. REFERENCES [1] W. Otten and C.A. Gilligan, "Soil structure and soil-borne diseases: using epidemiological concepts to scale from fungal spread to plant epidemics," European Journal of Soil Science, vol. 42, pp. 131-134, 2006. [2] X. Zhang, J.W. Crawford, and I.M. Young, "Does pore water velocity affect the reaction rate of adsorptive solute transport in soils? demonstration with pore-scale modelling.," Adv. Water resour, vol. 31, pp. 425--437, 2007. [3] I.M. Young, J.W. Crawford, N. Nunan, W. Otten, and A. Spiers, "Microbial distribution in soils: physics and scaling.," Advances in Agronomy, vol. 100, pp. 81121,2008. [4] S.N. Johnson, J.W. Crawford, P.J. Gregory, D.V. Grinev, R.W. Mankin, G.J. Masters, P.J. Murray, and X. Zhang, "Non-invasive techniques for investigating and modelling root-feeding insects in managed and natural systems," Agricultural and Forest Entomology, vol. 9,pp.39--46,2007. [5] W. Otten, K. Binnie, I.M. Young, J. Villot, and D. Grinev, "Interconnectivity of habitats in soil: combining x-r ay micro tomography and thin sectioning to reveal fungal-soil structure interactions.," Multitrophic Interactions in Soil IOBC WPRC, vol. 42, pp. 131-134, 2009. [6] MJ. Kwiecien, I.F. Macdonald, and F.A.L. Dullien, "Threedimensional reconstruction of porous media from serial section data," Journal ofMicroscopy, vol. 159,pp.343-359, 1990. [7] J. T. Fredrich, "3d imaging of porous media using laser scanning confocal microscopy with application to microscale transport processes," Phys. Chem. Earth, vol. 24,pp. 551-561,1999. [8] P. Spanne, J. F. Thovert, C. J. Jacquin, W. B. Lindquist, K. W. Jones, and P. M. Adler, "Synchrotron computed microtomography of porous media: Topology and transports," Phys. Rev. Lett., vol. 73, pp. 20012004,1994. [9] L. F. Gladden, M. P. Hollewand, and P. Alexander, "Characterization of structural inhomogeneities in porous media," American Institute ofChemical Engineers Journal, vol. 41, pp. 894-906, 1995. [10] L.A. Ferrand, M.A. Celia, and W.A. SolI, "Percolationbased models for pore-to-lab scale calculations in multifluid porous media," in Dynamics ofFluids in Hierarchical Porous Media, J.H. Cushman, Ed., pp. 463--483. Academic Press, San Diego, California, 1990. [11] B. Berkovitz and R. P. Ewing, "Percolation theory and network modeling applications in soil physics," Surveys in Geophysics, vol. 19, pp. 23-72, 1998. [12] A. Santiago, J. P. Cardenas, J. C. Losada, R. M. Benitol, A. M. Tarquis, and F. Borondo, "Multiscaling of porous soils as heterogeneous complex networks," Nonlin. Processes Geophys., vol. 15, pp. 893 - 902, 2008. [13] A.G. Hunt, "Percolative transport in fractal porous media," Chaos, Solitons and Fractals, vol. 19, pp. 309325,2004. [14] K. Harris, I.M. Young, C.A. Gilligan, W. Otten, and K. Ritz, "Effect of bulk density on the spatial organisation of the fungus rhizoctonia solani in soil.," FEMS Microbio!. Eco!., vol. 44, pp. 45-56, 2003. [15] http://www.simbios.ac.uk. [16] http://www.volumegraphics.com. [17] http://rsbweb.nih.gov/ij/. [18] W. Schroeder, K. Martin, and B. Loresen, The Visualization Toolkit, Prentice Hall, Upper Saddle River, NJ, 1997. [19] R. C. Gonzales and R. E. Woods, Digital Image processing, Addison-Wesley, 1992. [20] L. F. Costa and R. M. Cesar Jr, Shape Analysis and Classification: Theory and Practice, CRC Press Boca Raton, 2001. [21] M. P. Viana, E. Tanck, M. E. Beletti, and L. F. Costa, "Modularity and robustness of bone networks," Molecular Biosystems, vol. submitted, 2009. [22] L. F. Costa, G. Mutinari, and D. Schubert, "Characterizing width uniformity by wave propagation," Physical Review E, vol. 68. [23] C. A. Glasbey, G. W. Horgan, and J.F. Darbyshire, "Image analysis and threedimensional modelling of pores in soil aggregates," Journal ofSoil Science, vol. 42, pp. 479, 1991. [24] D.T. Gillespie, Markov Processes, Academic Press, San Diego, 1992. [25] P. Grassberger, "On the critical behavior of the general epidemic process and dynamical percolation," Math. Biosc., vol. 63, pp. 157-172, 1983. [26] D. Stauffer and A. Aharony, Introduction to Percolation Theory, Taylor and Francis, London, 2nd edition, 1992. [27] L. Sander, C. P. Warren, I. M. Sokolov, C. Simon, and J. Koopman, "Percolation on heterogeneous networks as a model for epidemics," Math. Biosc., vol. 180, pp. 293-305, 2002. [28] L. M. Sander, C. P. Warren, and I.M. Sokolov, "Epidemics, disorder, and percolation.," Physica A, vol. 325, pp. 1-8,2003. [29] M.E.J. Newman, "The spread of epidemic disease on networks," Phys. Rev. E, vol. 66, pp. 016128,2002. [30] F. M. Neri, F. J. Perez-Reche, S. Taraskin, and C. A. Gilligan, "Heterogeneity in sir epidemics on lattices," 2009, To be submitted. [31] M. B. Isichenko, "Percolation, statistical topography, and transport in random media," Rev. Mod. Phys., vol. 64,pp.961-1043, 1992. [32] E. Kenah and J. M. Robins, "Second look at the spread of epidemics on networks," Phys. Rev. E, vol. 76, pp. 036113, 2007. [33] 1. C. Miller, "Epidemic size and probability in populations with heterogeneous infectivity and susceptibility," Phys. Rev. E, vol. 76, pp. 010101(R), 2007.