MATERIALS & METHODS

The present study was carried out in the Chugoku district of Japan that is located in the West of Honshu island within (130^o 55' 16'’ and 133 ^o 12' 11'’) longitude and (33 ^o 57' 40'’ and 35 ^o 23' 34'’) latitude. It includes five Prefectures (Hiroshima, Yamaguchi, Tottori, Shimane and Okayama), and covers 32,000 km2. The underlying geology is largely volcanic (andesite and rhyolite) in central and the northern parts, Mesozoic sedimentary formations (sandstone/ shale/ pudding stone) in the west and quaternary sedimentary formations (gravel and clay) in the low lands of the study area. Dystric regosols, fluvic gleysols, gleysols, humic cambisols, ochric cambisols and rhodic arcsols are the dominant soil types that can be observed in each catchment. The spatial analysis of land cover map indicated that a high percentage (79.34%) of the study area is covered by forest. Other land cover classes including urban, agriculture, grassland and water body make up 5.15, 8.33, 6.63 and 0.46% of the study area, respectively (Table 1). There are 7,732,499 inhabitants in the study area. There are 54 major rivers in the Chugoku district, out of which 21 were selected as our study area.

Total Nitrogen (TN) and Total Phosphorus (TP) were selected as the factors representing the water quality of the rivers (Table 2). The water quality data were obtained from the five Prefecture offices (Hiroshima, Yamaguchi, Shimane, Tottori and Okayama) in the study area that could be defined as secondary database (Sliva and William, 2001). Water sampling of the stream and analysis are carried out monthly based on Japanese Industrial Standard (TN/TP: JIS K0102 45.2, 45.3). Details of sampling method and analysis procedures would be found in JSA JIS K 0102 (Japanese Standard Association, 1998). Annual mean of water quality data in 2001 were calculated using monthly data without data transformations for this study. Population data was obtained from the population census in 2000 that carried out by Japanese Statistic Bureau of Ministry of Internal Affairs and Communication. The present study was conducted using wholecatchments land use approach (Sliva and Williams, 2001). This approach uses the results of spatial analysis for the land use or land covers in whole catchments as variable in the multivariate analysis. In order to facilitate the spatial analysis and determination of morphological attributes, Geographical Information System was established using theArcView 3.2 (ESRI, 1999).For allwater sampling stations, upstream catchmentsboundaries were hand digitized by using JGSI topographic quadrangle maps (1:200 000). All databases were transformed into a common digital format, projected onto a common coordinate system (UTM, zones 52 and 53). Human Population density map was firstly generated by linking county-scale population database with the digital map of counties. It was then overlaid by the catchments map and aggregated in order to find the catchments-scale population density map in the study area and was transformed into spatial data. For generating land cover map in the study area, two scenes of NASA Landsat-5 TM (2000) were used. The satellite images were georeferenced by the affine procedure. The supervised classification method was applied to classify land uses, which included forest, agriculture, grassland, urban and water body (wetland, natural and artificial lake). Finally, the generated land cover map was verified using JGSI maps. Preparation, interpretation and analysis of the satellite data were carried out using Integrated Land and Water Information System (ILWIS, 2004). The generated land cover map was then superimposed by catchments map for calculating the real extend of each land cover for each catchment, which subsequently divided by the catchments area to drive the percentage of the catchments covered by each type.

Water quality data including TN and TP and area (%) of land cover types were tested for normality using the Sharpio-Wilk test with a pvalue of less than 0.05 (Table 3). For determination of linkage between land cover and stream water qualities, Multiple Linear Regression (MLR) modelling was applied using stepwise approach. Inter-variable co-linearity of the models was investigated referring to Variance Inflation Factor (VIF). It measures the impact of co-linearity among the X’s in a regression model on the precision of estimation. VIF expresses the degree to which co-linearity among the predictors degrades the precision of an estimate. (Chatterjee et al., 2000). Normality of residuals of the models were then examined by Anderson-Darling test with a p-value<0.05. For a given water quality variable, the appropriate model was selected based on regression statistics (r2, p-value) and considering the significance of the coefficients of the model provided that the residuals of the model was normally distributed. Finally, the goodness-of-fit of the statistically significant models was evaluated by scatter plot and simple linear regression of observed versus equivalent model prediction (Ahearn et al., 2005). All descriptive statistics and other statistical analysis were completed using Excel Add-ins (XLSTATTM 2006) and SPSS for Windows Release10. In environmental modelling, it is very important to know to what extend the independent variables contri1bute in explaining the variations of the dependent variable (Lek et al., 1999). Sensitivity analysis was considered by some researchers (e.g. see, Maier and Dandy 2001; Lek et al.,1996; Lek et al.,1999 ;Spitz and Lek, 1999) whom took advantage of some techniques in environmental informatics such as artificial neural networks for specify how much changes in the input variables would contribute in explaining the variations in the model output.

As a part of the sensitivity analysis, each of the inputs is increased by a certain percentage (5%) in turn, and the change in the output caused by the change in the input is calculated (Maier and Dandy, 2001). In applying this approach, which is classified as mathematical sensitivity analysis, a critical point should be kept in mind that all the variables of the model should be independent. In our case, mathematical sensitivity analysis was not seemed to be possible since structure of land use types within a given river basin is inter-related so that making change in the percentage of a land use type (e.g. increase of the urban area) should be in expense of other one (e.g. decrease in agricultural area) in a given river basin. Therefore, a stochastic sensitivity analysis was paid attention for this study. Considering the number of rivers that their water quality data were available (21 basins), it seemed that number of rivers that should be kept for testing the models would be very few (3 basins). In this case, reliable information will not be provided for judging on the calibrated models. Application of testing data subset was avoided for testing the calibrated models, which was real but few in numbers. For dealing with the shortage of the testing data set, it was decided to use Monte Carlo Method for sensitivity analysis of the calibrated models. According to the Monte Carlo Method, 5000 random numbers were generated using random number generator command in Microsoft Excel for six variables consisting of area (%) urban, forest, agriculture, grassland and water body, and human population density. Table 4 summarized the statistics of this randomly generated data set. The calibrated models were tested using the generated random number data set based on the Monte Carlo procedure. In generating random number, three conditions were considered as follows:

1) random numbers should vary between 0-100; 2) summation of five variables for land use type to be equal to 100%, and 3) random numbers for each variable should follow a known distribution.

Initially, six variables were defined in the Microsoft Excel. 5000 random numbers were generated for each variable using the related command (rand between (bottom)(top)) in it for meeting the first condition. Each row stands for a hypothetical river basin, which consisted of five variables for land use types (urban, forest, agriculture, grassland and water body areas in percent) and the last variable for human population density. For meeting the second condition, summation of five land use variables was calculated and each land use variable divided by the summation since in a real situation for a given river basin, land use types as land use variable can not be independent and any change in one land use type should be in expense of other ones. Based on the Monte Carlo Procedure, statistical distribution of variable(s) should be known before they are used for further analysis. Distribution fitting approach was applied in order to specify what statistical distribution every randomly generated variable follows. It was carried out by computing the Maximum Likelihood test with maximum 100 iterations and accuracy of 10-4 . For goodness of fit, the Kolmogrov-Smirnov test was carried out. The results of the distribution fitting of the input variables were summarized in Table 6. Following on from that, this data set was entered into LR models that were developed for the TN and TP separately. The output of each model was then computed. For the output of MLR models, the distribution fitting was carried out one more time in order to determine what statistical distribution the output of the models follow? Moreover, another purpose for distribution fitting was to specify the probability of the event such as Pr (output< 0), since only positive values make sense in relation with the water quality variables in general and for TN and TP in particular. For materializing aforementioned purpose, the estimated cumulative distribution curve was plotted for the output of TN and TP models, which were developed.

TN = 0.219+1.375*10^-3P+5.151*10-²U (1)

where TN is concentration of total Nitrogen in mg/L in the river, and P stands for human population density (person/km2) in the catchment, and U is urban area (%) in the catchment. Eq. 4 indicates that increasing the proportion of urban area (%) would increase total nitrogen level in the rivers. It is confirmed by the findings of Tong & Chen (2002) on the strong relationship between the total nitrogen levels and residential (urbanized area). The TP model (Eq. 2), (r²=0.47, p<0.01), only compositional attribute (%) of urban area was introduced into the regression model as significant variable. Others were excluded during the stepwise regression modelling process.

TP=1.362*10^-2+6.495*10^-3U. (2)

where TP is concentration of total phosphorus in mg/L in the river, and U is urban area (%) in the catchment. Eq. 2 suggests a positive relationship between the proportion of urban area (%) and concentration of total phosphorus in the rivers. It is confirmed by the findings of Tong & Chen (2002) and Ngoye & Machiwa (2004) who reported strong relationship between TP levels and urban area. The statistics value (r² ) for the TP regression model is not relatively high. It might be related to effects of other factors, which were not considered in the present study. These factors could be background sources of chemical loads from either catchment geology and soils or catchment hydrology, which plays important role in the transport of pollutants from the catchment surface to the rivers (Jarvie et al, 2003). Intervariable co-linearity of the models was investigated by referring to VIF (Table 5). While a VIF>10 could be considered as severe co-linearity within variables in the models (Neter et al., 1996; Chatterjee et al., 2000), all models have revealed no co-linearity (VIF<2). Normality of residuals of the models were tested using the Anderson-Darling test with a p<0.05 whether they follow a normal distribution. The results of the test were summarized in Table 5. It suggests the residuals of all models were normally distributed at significance level of p<0.05. For validating of the goodness-of-fit of the statistically significant models, simple linear regression analysis of observed value versus the values predicted by the relevant models was carried out. Sensitivity analysis of TN model suggested that stream water TN concentration might be predicted between 0.4 to 3.60 mg/L using the developed model. On the other hands, upper limit, which were specified by sensitivity analysis, might not be responded accurately by the model in real situation. Besides, it seems to be no probability for Pr (TN<0) since the cumulative distribution curve achieved a turning point at 0.4 mg/L. It might be mentioned that the TN model would not be able to predict the targeted water quality variable in a level more than 3.6 mg/L since the cumulative distribution curve have reached a turning point at this concentration and a steady state was revealed afterward. Concerning the TP regression model, the sensitivity analysis revealed that the model might be able to predict the stream water TP concentration within a range of 0.04 to 0.4 mg/L. This model could not be used for prediction of the TP concentration beyond the upper and lower limits. Moreover, no probability for occurring the Pr (TP<0) was observed for the output of the TP model. To know the probability for Pr(x<0), which x stands for any water quality variable, would be very important since the generation of negative values by water quality predictive model would not make any sense.

Application of the sensitivity analysis, which used the randomly generated data, could provide reliable information with us for judging on the prediction abilities of the developed models stochastically. On the other hands, using the sensitivity analysis, which was based on the Monte Carlo Method could dealt with the issue of lack of access to the sufficient database in numbers. The results of this study indicated that the Monte Carlo method-based sensitivity analysis would be helpful for determination of the upper and lower limits of the prediction abilities of the models, which are developed, in case of facing with insufficient database in number. It might be noted that the prediction of the behaviour of the TN and TP regression models, which were conducted in present study, would be valid provided that input data that will use for the prediction of stream water TN and TP concentrations to be the range shown in Table 4. It will be obvious that the behaviour of the model will be in expose of uncertainty if the input data to be out of the range. Results of this study indicated that the sensitivity analysis might be considered as one of inevitable steps in model development. Its importance should not be underestimated once any effort is taking into action for investigating the linkage between land cover and stream water quality. Information, which is provided with us by the sensitivity analysis about output space of the model that was developed would be more reliable than those which is brought up by judging from initial conditions based on which the model was built.