FRACTAL ANALYSIS OF RIVER NEWTWORKS BASED ON REMOTE SENSING DATA: AN EXAMPLE FROM KELANI RIVER BASIN, SRI LANKA

River network patterns are important components of basin geomorphology. The geomorphological complexity of river networks has an impact on various natural processes and hazardous events. In order to find those correlations and establish quantitative relationships, a comprehensive textural examination on river network patterns is crucial. In this study, the fractal geometry approach is applied to analyse river network by the means of their geometric distinctiveness. Since river networks are fractal objects formed by repetition of certain natural processes over a long period of time, fractal geometry provides the most sensitive method to analyse their complex branching structure. Five river networks from Kelani river basin (1. Ambalanpiti oya, 2. Gurugoda oya, 3. Pugoda oya, 4. Pusweli oya, 5. Wakoya) are subjected to fractal analysis. Fixed-size box counting algorithm is applied to obtain fractal measures. The recursive algorithm is applied to the same river network twice; allocating weighted-lines for different orders and then allocating non-weighted equal width line for all tributaries, to test the best suitable option to model their geomorphological complexity. Non-plane-filling behaviour of river networks is confirmed by present values which are greater than 1 and lesser than 2. The largest fractal dimension value is for Gurugodaoya tributary, confirming that it has the most complex geomorphologic pattern. Smallest value is from Pugodaoya concluding the least geomorphological complexity. Multifractal spectra f(α) are constructed for each river network and detailed investigation is required (considering lithological features of the basin) to link f(α) to the physical characteristics of river network.


INTRODUCTION
Computation of spatial patterns and their detailed analysis using nonlinear approach is growing rapidly in many fields including landscape ecology. The multifractal theory was introduced by Mandelbrot (1972Mandelbrot ( , 1985 and developed by Frisch and Parisi (1985) for the study of turbulent dissipative phenomena and by Halsey (1986) for the analysis of strange attractors in the theory of dynamic systems. River networks have long been recognized as scaling, possessing self-similar structures over a considerable range of scales. Horton (1932Horton ( , 1945 provided a set of scaling laws, later refined by Strahler (1952) and Schumm (1956). Mandelbrot (1983) used empirical length-area power law relationship to imply that rivers are fractal objects. Superimposition of different physical and geological processes governing the evolution of rivers, also including random components, produce fractal river networks (Dombradi et al., 2007) as they satisfy the minimum energy requirement. The major objective of this approach is to spot the diversities in natural patterns and study their causes. The linear analyses essentially focus on the secondary parameters (contributing drainage area, channel length, channel slope, elevation, and etc.) and completely overlook the fractal behaviour of drainage network. Typical linear approach uses the slopearea correlation which generates the same results for dissimilar causative effects. In the other hand, a fractal measure is a powerful tool for patterns with dissimilar space filling characteristics. Fractal analysis probes the linearization, heterogeneity and connectivity of the drainage patterns (Mahmood et al., 2011).
Also, river basins and networks were read with a multifractal lens: investigations were effected about the scale properties of some relevant hydrological variables. Particularly, the descriptions of contributing drainage areas of given hierarchical order, slopes, dissipation energy, channel initiation function and width function was realised through the use of the multifractal spectrum formalism proposed by Halsey (1986). The multifractal spectra are index of the river network geometric complexity. The relationships between physical characteristics of the drainage network and multifractal spectra were investigated by Bartolo et al. (2004). The idea of lithological control over the multifractal behaviour of river networks was established by recognizing the influence of source rocks and underlying lithology of river basins (Ijjasz-Vasquez et al., 1992).
Although some studies have been conducted to analyse the river network of Sri Lanka, fractal geometry approach has never been applied to investigate the geomorphological complexity of river network. Therefore, the main objective of this study is to apply the fractal geometry approach to selected local river network and give a quantitative depiction on their geomorphological complexity. This case study will provide a novel method to analyse drainage networks. It is also valuable as a secondary data source in case of investigating the correlations between river network geomorphological complexity and natural processes/hazardous events (flood, erosion, landslides, surface deformation, sediment transportation, etc.). Kelani river basin, one of main four river basins in Sri Lanka and having been suffered from flood recent decades was selected for the case study. As a whole, this study has a great referential value and proposes a new decision supportive tool for natural disaster risk management.

STUDY AREA
Kelani river basin is one of the main river basins of Sri Lanka. Hence the Kelani river network mainly contributes to govern major impacts on local irrigation system. Also severe flood events are being recorded each year from Kelani basin. Kelani River is the fourth longest river (145 km

DATA SOURCES
In order to extract the river network of the selected area, a digital elevation map of Sri Lanka was formed using ASTER Global Digital Elevation Models (ASTER GDEM -30m resolution). ASTER GDEM is a product of METI and NASA.

COMPUTING METHODOLOGY
River network was extracted from ASTER Global Digital Elevation Models (ASTER GDEM -30m resolution). ArcMap 10.4 software package was used to prepare binary images of GDEMs. Arc hydro tools extension of ArcMap 10.4 was used for next steps in image (DEM) processing. All artificial sinks which are not exactly the parts of the image were removed from the prepared binary image of study area. The D8 flow grid algorithm (Jenson and Domingue, 1988;O'Callaghan and Mark, 1984) was used to determine the flow direction. D8 algorithm computes the potential flow guidelines at every pixel towards the neighbouring 8 pixels and the least cost path (whichever the pixel has steepest slope) was fixed as the direction of water flow. The flow accumulation map was constructed according to the fixed flow direction. In the step of defining streams, threshold value was given such that the Strahler order of extracted streams is higher than or equal to 6. The river network of study area was prepared as a binary image where the streams have a pixel value of 1 and rest of the space is given the pixel value zero (Melo et al., 2006). Finally, five selected river networks of Kelani basin ( Figure 2) (1. Ambalanpiti oya, 2. Gurugoda oya, 3. Pugoda oya, 4. Pusweli oya, 5. Wakoya) were prepared as separate networks for the sake of subject to the fractal analysis.

Fig. 3.
Extracted drainage patterns of selected subwatersheds compared to the local elevation.

CONSTRUCTION OF WEIGHTED-LINE NETWORKS AND NON-WEIGHTED-LINE NETWORKS
In fractal analysis, natural fractal object approximates by a fractal structure constructed/ designed in a spatial environment. Normally, river networks extract as poly-line networks for the purpose of analysing. The present work deals with two different fractal structures in case of natural river networks. In first structure, the whole network was given the same line (same pixel width) ( Figure 4). In the other structure, different lines (different pixel widths) were allocated for streams with different Strahler order. In natural river networks, width of the stream gets higher when the order is increasing. Extracted river networks were given line widths considering aforesaid natural behaviour of rivers. Line width was set such that 1 st order streams allocated a line of one-unit width, 2 nd order streams drawn allocated a line of twounit width, 3 rd order streams drawn allocated a line of three-unit width and so on ( Figure  5).

FRACTALS AND THE FRACTALITY OF RIVER NETWORKS
Fractal geometry describes irregular and complex features of natural objects. Fractals are self-similar statistically for a range of scales (Zang et al., 2015). Achieving the planefilling dimension depends on the density of the object. If S is a dense subset of an open region of R n , then its dimension is equal to n (Falconer, 1990), i.e., to 2 in the plane case. But river networks can be assimilated to tree-graphs (connected acyclic graphs), which are not dense when the number of nodes is finite (West, 1996): this prevents the plane filling dimension from being reached. Therefore, the dimension of a river network varies according to its geomorphologic density/ complexity. Traditional morphologic methods for describing the river networks are to consider the pattern characteristics such as, radius of river bend curvature, channel length, and other parameters. Those methods can only investigate the simple geometric configurations and poorly describe natural features of complex river networks (Snow, 1989;Mandelbrot, 1977;Schuller et al., 2001;Hassan and Kurths, 2002;Guillermo et al., 2004;Zang et al., 2015). Many studies have proven the importance of analysing the irregular components of river networks (Leopold and Wolman, 1960;Schumm, 1977;Fredsoe, 1978;Tarboton, 1993Tarboton, , 1995Tokunaga, 1978).
Fractal geometry method includes river shapes characterized by irregularities (Zang et al., 2015). Fractality of a river networks describes its geometrical distribution and complexity (Mahmood et al., 2011).

FRACTAL DIMENSION AND MULTIFRACTAL APPROACH
Fractal dimension (D 0 ) of a river network quantifies the geomorphological complexity of a river and the ability of network to fill a plane. Incorporating the effect of sinuosity, results for a river network fractal dimensions became closer to 2, exhibiting a space-filling behaviour according to the topologically random development of tributaries (Shreve, 1966). But the studies of natural basins were evidenced dimensions being nonspace-filling, indicating that geological constraints prevent channels from developing as purely branching process. In this study, fractal dimension (D 0 ) for each river network was obtained using generalised box counting method and compared.
While a single fractal dimension describes the complexity of a fractal set, multifractal theory (Mandelbrot, 1974;Frisch and Parisi, 1985;Halsey et al., 1986) describes the distribution of a measure over the fractal set and explains dense and sparse areas separately (Ijjasz-Vasquez et al., 1992). i.e., it generates a series of fractals, each with its own scaling exponent (Lipschitz-Hölder/ Singularity exponent/ α), and the corresponding monofractal dimensions (Aharony, 1989). There are two ways to represent multifractals. Either by infinite number of generalised fractal dimensions (D q ) or by singularity/ multifractal spectrum f(α). Lipschitz-Hölder exponents express how the probability of finding neighbouring points belongs to the river network changes with the distance. (Dombradiet al., 2007). In the present work, multifractal spectra were generated for all five river networks (for weighted-line networks and nonweighted line networks) separately. The shape of curves ware compared. Basically, spatial variation of fractal dimension was observed.

GENERALIZED BOX-COUNTING METHOD
The generalized box-counting method is one of fixed-size algorithms. Number of box-counting methods were introduced and developed by many authors. It is applicable in analysing the fractality of river networks (Bartolo et al., 2000(Bartolo et al., , 2004. Box-counting method can be used as a tool to evaluate fractal dimension as well as for the multifractal analysis, with some limitations. The procedure of box-counting method used to evaluate fractal dimensions as follows (Turcotte, 1992): (i) A square box of side length ε was generated. (ii) A grid of those boxes was placed over the binary image of the river network. (iii) The number of boxes, N (ε), which a part of the fractal falls was counted. (iv) Side length ε was decreased and above steps were repeated. (v) Same procedure was repeated for 100 iterations. (Figure 6 and Figure 7) In theory, L = Nε gives the actual length of the river network when ε tends to zero. The area of image for all five river networks was fixed. The relationship between the fractal dimension and the box-counting method is as follows: After performing this method for a series of grids, the relationship between log N (ε) and log (ε) was plotted in a graph. The fractal dimension (D 0 ) of each river network was obtained by the slope of fitted regression line and the corresponding equation can be written as follows:
Equations and relations employed in multifractal analysis are as follows: Where N i is the number of points in the i th cell of given size (ε) and N is total number of points in the fractal set (river network). Thus, μ i gives the probability of finding neighbouring points that belong to the fractal set. The generalized fractal dimensions in box-counting method were computed for q ≠1 as follows:
Summing the q th momentum of the number of points (normalized with the total number of points, μ i = N i / N) through the N(ε) number of cells covering the fractal object at ε resolution.
For the special case of q = 1,

D 1 = (6)
The Lipschitz-Hölder exponent αand the multifractal spectrum can be expressed from generalized fractal dimensions by introducing the mass exponent τ (q):

RESULTS
First, the box-counting method was applied to obtain the fractal dimension, considering each river network as a mono-fractal (Eq.s (1) and (2)). ln N(ε) Vs. ln (ε) graphs were plotted for all fractal sets and the fitted lines were drawn on the same graph (Figure 8. / Figure 9. (A, B, C, D, E)).  In order to compare the geomorphological complexity of selected river networks, D 0 (in monofractal analysis) and D q = 0 (in multifractal analysis) were evaluated. Other generalized fractal dimensions were not investigated in detail. Boxcounting method was supposed to be biased by border effects and inability to construct the right side of multifractal spectra. Even though the right side was constructed, it was not taken part in analysing due to the high probability of error occurrences. All quantitative results were summarized in Table 1. The goodness of fitted regression line was indicated by R 2 value. All R 2 values of river networks drawn with weighted lines are higher than 99%. But some of R 2 of nonweighted line structures were below 99%. The highest value of fractal dimensions was from Gurugodaoya while Pugodaoya gives the lowest values.

DISCUSSION
In the present work, river networks were analysed in the scope of their geomorphological complexity by applying fractal geometry method. Fractal properties of river network textures were revealed and their geometrical complexity distribution was quantified. Results of fractal analysis was corroborated the nonspace filling behaviour of drainage systems as all obtained fractal dimensions were greater than 1 and lesser than 2. The generalized boxcounting method was supposed to suffer from border effects and due to boxes with few points take count, accurate dimensions for negative orders were not able to determine. Therefore, only the left sides of multifractal spectra were considered as the reliable curves. In case of comparing the qualitative characteristics of multifractal spectra of river network geomorphological complexity, considering left side was regarded as appropriate. Area of fractal Box sizes (grid calibres) of 100 and 20 were selected for iterations in boxcounting method for monofractal and multifractal analyses respectively. The set of moment orders (Q set) was set -10 to +10, but in optimized results, moment orders of -1 to +3 were included. Moment orders are set of arbitrary exponents for the calculation of multifractal spectra. Multifractal spectra were calculated using functions of the sums of the mass probability distribution raised to the values of Q and seeing how the function changes with the value of Q. For an instance, higher values of Q were affected the denser parts of a probability distribution for the multifractal spectra. In case of a statistically self-similar fractal set, fractal dimensions for different moment orders (Q) would be the same and consequently, the multifractal spectrum would consist merely a single point (e.g. Ijjasz-Vasquez et al., 1992). Asymmetry of multifractal curves expresses dominancy of low or high fractal exponents in respect to the other (e.g. Telesca et al., 2004). Construction of multifractal spectra, multifractal behaviour and the heterogeneous scaling properties of river networks were observed for the analysed subriver networks of Kelani river basin.
Due to the errors of DEMs, location and distribution of few tributaries of extracted river network was found less obvious (in Figure 1 and Figure 2). But the selected river networks were not affected by those uncertain tributaries and were reliable enough according to the referred maps of Irrigation department and Survey Department, Sri Lanka.
Differences of the multifractal curves were well observed. Singularity/ multifractal spectra f(α) were recognized as index of geomorphological complexity of river networks and the almost all curves were left-skewed (negatively skewed). Necessity of investigating more geophysical characteristics of natural river networks (underlying geological structure, rock type, etc…) was noted as they would helpful to link the behaviour of f(α) and compare the river network evolution and control over those facts.
It is a well-known fact that river networks and drainage patterns are controlled by the geological and structural factors such as relatively soluble or weak layers of the rocks, fractures, faults, and other such structural features. The ultimate objective of this exercise is to demonstrate the strength and the use of fractal approach in modelling the river networks in local context. It was also our intention to experimentally test the most suitable criterion (weighted and non-weighted lines) to model the natural river networks.
A set of quantitative data of selected five river networks (1. Ambalanpiti oya, 2. Gurugoda oya, 3. Pugoda oya, 4. Pusweli oya, 5. Wakoya) from Kelani river basin was generated. This case study was recognized with a high referential value for further river network analyses in aspect of fractal geometry approach. Data points departing from the main system, mainly due to human activities such as change of natural course of the rivers, may indicate the potential points of disaster or vulnerable points. Also, the generated data were established as valuable secondary data source in watershed management, constructing correlations and creating decision supportive tools in disaster risk management as geomorphological complexity of river networks has various direct and indirect impacts on many natural phenomena/ processes/ disasters (flood, landslides, erosion, surface deformation, sediment transportation, etc.).

CONCLUSIONS
Fractal analysis of five sub-river networks (1. Ambalanpiti oya, 2. Gurugoda oya, 3. Pugoda oya, 4. Pusweli oya, 5. Wakoya) from Kelani river basin was performed. Multifractal behaviour of river networks was revealed and the non-plane-filling configuration was established by obtained fractal dimensions. Quantified geometrical complexity was evidenced that the evolution of river networks is not completely affected by a random process and geological constraints including also other geophysical facts govern instead. Results of both weighted-line fractal structures and non-weighted-line fractal structures of river networks were compared and weighted-line structures were pronounced as more objectively real fractal structures as R 2 values of all the weighted-line structured river networks were obtained above 99%. Highest values of fractal dimensions (1.6305 and 1.6026 for weighted and non-weighted lines respectively) were obtained for Gurugodaoya tributary, concluding the highest geomorphological complexity, while the lowest fractal dimensions were obtained for Pugodaoya concluding the lowest geomorphological complexity. All multifractal spectra showed negatively skewed behaviour, which indicated dense areas of the river networks tend to occupy the space more homogeneously while sparse areas are more heterogeneous. This work demonstrate the analysis of river networks using fractal approach can provide essential information on river network complexity, as well as help predicting the vulnerable points and potential for flood in a given river tributary.