Topography as a macroscopic index for the dissolved iron productivity of different land cover types in the Amur River Basin

Iron is the limiting nutrient of phytoplankton in the Sea of Okhotsk, and the majority of iron in this system is fed by the Amur River. The recent conversion of wetlands, the main source of iron in the Amur River basin, to agricultural lands will likely impact dissolved iron productivity, which may also influence primary production in the Sea of Okhotsk. Therefore this study was conducted to construct a macroscopic index for use in assessing dissolved iron productivity in the basin. Correlation analysis between climate and topographic parameters and the observed dissolved iron concentration in forests and wetlands revealed that the topographic wetness index (TWI) had a significant correlation with dissolved iron concentration. An exponential curve was found to be the best curve to express this correlation. We assumed that dissolved iron concentration for grasslands and agricultural lands, the other two dominant land cover types, could also be expressed by TWI. Based on this assumption, dissolved iron concentration curves for grasslands and agricultural lands were inversely identified by systematic modification of the curve for forests and wetlands. The results suggest that TWI can describe the average dissolved iron concentration of major land cover types in the basin.


INTRODUCTION
proved that iron was the limiting nutrient of phytoplankton growth in the Northeast Pacific Ocean.The Sea of Okhotsk and the adjacent Oyashio region is also known as a region in which iron is the limiting nutrient of phytoplankton growth (Tsuda et al., 2003).In general, the primary source of oceanic iron was assumed to be aerosols until the 2000s.However, Nakatsuka et al. (2007) proposed an Intermediate-Water Iron Hypothesis based on intensive observations in the Sea of Okhotsk.While there is still some uncertainty associated with this hypothesis, it is highly probable that it explains the major part of iron influx into the Sea of Okhotsk via fresh water from the Amur River (Nishioka, personal communication).Accordingly, recent wetlands conversion to agricultural lands in the basin will likely impact dissolved iron productivity, which may also influence primary production in the Sea of Okhotsk.Thus, it is essential to estimate the amount of dissolved iron produced in the Amur River to assess the impact of the Amur River Basin on primary production in the Sea of Okhotsk.
In general, iron is present in low concentrations within the range of pH and redox conditions of surface water due to the low solubility and thermodynamically stable state of ferric iron.Conversely, the solubility of iron increases under reducing conditions because microbial oxidation consumes electrons from ferric iron.In addition, humic substances strongly interact with both dissolved and particulate iron species (Tipping et al., 1981;Davis, 1982;Warren and Haack, 2001), subsequently forming an iron-humic substance complex, by which dissolved iron is stabilized.Thus, the redox process and interaction with humic substances are two major processes governing iron solubility.
There are two approaches involved in modeling dissolved iron production.One is coupling of a physically based reactive transport model and an iron-binding model to humic substances.Elaborate numerical models that can deal with these processes have been developed during the last several decades (Tipping, 1998;Kinniburgh et al., 1999;Šimunek et al., 2006).However, such models can only deal with a spatial scale of several hundred meters to several kilometers due to spatial variability, computational limitations, and the lack of observed data for validation.Although these models formulate biogeochemical processes explicitly, it is not feasible to apply such models to continental watersheds.
The other possible approach is identifying empirical macroscopic indices that can characterize biogeochemical processes in a simple manner.Topography is one such index because it is an important factor that governs hydrological conditions, which in turn affect various biogeochemical processes.Indeed, studies employing topography as an index have been conducted by several researchers including Vitousek (1977), Ogawa et al. (2006), and Anderson and Nyberg (2009), who all found correlations between topographic parameters and the chemical composition of stream water.Thus, taking the latter approach, this study attempted to identify a macroscopic index including topographic parameters that would represent the dissolved iron concentration of different land cover types in the Amur River Basin.

Study site
The study site was the Amur River Basin (Figure S1), which has a catchment area of 2,050,057 km 2 .The total length of the Amur River is about 4,300 km.The amount of annual fresh water supplied to the Sea of Okhotsk by the river is about 300 km 3 .Average annual precipitation ranges from 300 mm in the west to more than 700 mm in the east.The mean annual temperature also varies from −7°C in the north to 6°C in the south.
One distinguishing spatial pattern of the land use/land cover (LULC) of the basin is the contrast between the Russian side and the Chinese side (Figure S1).Specifically, the majority of the Russian side is forest area, while the land on the Chinese side is primarily used for agriculture.If all forests in the study area are classified as a single forest type, the four most dominant LULC types in the basin are forest (59.5%), agricultural land (dry land and paddy field, 18.3%), grassland (12.2%), and wetland (6.9%).

Data
The observation points for river discharge and dissolved iron concentration are shown in Figure S1.We obtained discharge data from the main course (Stations 6 to 8) and at several large tributaries (Stations 1 to 5).The discharge data was provided by the Federal Service for Hydrometeorology and Environmental Monitoring (ROSHYDROMET) and the Global Runoff Data Center (GRDC) in Koblenz, Germany (http://grdc.bafg.de).The time resolution was daily at Stations 6 to 8, and monthly at Stations 1 to 5.
Dissolved iron concentrations from 1980 to 1995 were also obtained from ROSHYDROMET from a total of 38 sampling points that were sampled about once a month from April to October of each year.The discharge rate was also observed at the same time.Dissolved iron was measured by the colorimetric method with 1,10-phenantroline, which was applied to water filtered through Whatman GF/F filters and acidified to pH < 2 with HCl (Hydrochemical Institute, 2006).
The H08 data set (Hirabayashi et al., 2008) was utilized for climate data such as average, maximum, and minimum air temperature, downward short wave radiation, specific humidity, and precipitation.The spatial and time resolution of the H08 are 0.5° and daily, respectively.SRTM3 data derived from the Shuttle Radar Topography Mission of NASA was used for DEM.A coarser DEM data set with a grid size of 1000 m was produced by averaging SRTM3 for the analysis.

Method
To identify a primal parameter of dissolved iron concentration, correlation analysis was conducted.The average dissolved iron concentration of each watershed was used as the objective variable, while explanatory variables included the spatio-temporal average of climate parameters and spatial average of topographic parameters that might govern the dissolved iron production of each watershed.Specifically, the climate parameters included annual precipitation, summer/winter precipitation, annual average air temperature, and average air temperature during summer/ winter, while the topographic parameters were a/tanβ, slope, and Laplacian.
For these climate parameters, summer was defined as the period from May to August, and winter was defined as September to April of the next year according to the definition by Tachibana et al. (2008).For a/tanβ, a was the watershed area per unit length of the calculation grid, and tanβ was the slope of each grid (Beven and Kirkby, 1979).Since a/tanβ is recognized as a good index of wetness, we hereafter refer to a/tanβ as the topographic wetness index (TWI).Slope was defined as the steepest gradient of each grid, which can be estimated by choosing the steepest gradient among eight surrounding grids for each grid.The Laplacian is an index of land surface roughness, a definition of which is given in document S1.
The calculated period was from 1980 to 1995, during which land cover conditions could be considered the same as for the year 2000 (Chinese Bureau of Statistics, 1980Statistics, -2000)).The watershed area and LULC composition of the measuring points are summarized in Table S1.When calculating the average dissolved iron concentration, the arithmetic average and weighted average were calculated.As a weighting function, measured discharge was used.
Where there was a lack of measured discharge data, monthly precipitation for each watershed was used as an alternative weighting function.The spatio-temporal averages of the climate parameters and spatial averages of the topographic parameters were the arithmetic average of each parameter included in each watershed area.The topographic parameters were calculated by utilizing the DEM.

Correlation analysis
The results shown in Table I clearly indicate that TWI and slope were correlated with dissolved iron concentration, while no other parameters showed a clear correlation with dissolved iron concentration.The watershed area of data used for correlation analysis ranged widely from 100 km 2 to 233,000 km 2 (Table SI).In spite of non-uniformity in spatial scale, topography was found to be a good index of dissolved iron concentration.Because the correlation coefficient of TWI was slightly higher than that of the slope, TWI was adopted to express the dissolved iron concentration.Since no distinct differences between the arithmetic and weighted average were observed, the weighted average was used in the following analysis.

Construction of a concentration curve for forests and wetlands
We attempted to construct a concentration curve for forests and wetlands with respect to TWI.In the following analysis, dissolved iron concentration at a given point was calculated with a spatial resolution of 0.5°.Thus, catchments of which the watershed area was less than 10,000 km 2 were extracted from the original data.In addition, watersheds in which agricultural lands occupied more than 1% of the catchment were excluded.As a result 17 points, of which the dominant LULC types were forest and wetland, were retained.Using this data, mean annual dissolved iron concentrations were plotted against the average TWI.The calculation period was also from 1980 to 1995.Three different types of curve, i.e. linear, power, and exponential curves, were tested as fitting curves.Figure 1 shows the exponential curve that was found to have the highest correlation coefficient among the three different curves.Correlation coefficients of the linear curve and power curve were respectively 0.49 and 0.59 for average, 0.49 and 0.36 for maximum, and 0.53 and 0.57 for minimum.

Uncertainty analysis of concentration curve
The concentration curve inductively generated was only applicable to a watershed primarily covered by forest and wetland.Thus, if it was used to predict the dissolved iron concentration at other points, it is likely that some discrepancies would occur.In addition, as shown in Figure 1, the inter-annual fluctuation range of dissolved iron at each point was large, especially for higher TWI values.Thus, within the range between the fitted curves of the maximum and minimum value shown in Figure 1, Monte Carlo Simulation was implemented to evaluate the uncertainty inherent in the curves.The trial numbers of Monte Carlo Simulation were 1,000 for each curve.The linear congruential method was used to generate random numbers.Results were compared against observed dissolved iron concentrations along Stations a-g.Dissolved iron concentration at a given point can be estimated by summing up the annual discharge weighted average of dissolved iron concentration of grids that are included in the watershed area of the target point.Document S2 provides a calculation procedure used to determine the annual discharge from each grid.
Figure 2 shows a comparison between observed and estimated dissolved iron concentrations during the period from 1980 to 1995.The calculated values were overestimated except for Station c, regardless of the type of fitting curve.Even if overestimation of the discharge of the Songhua River was taken into consideration, discrepancies between the observed and calculated values were significantly large.Because agricultural land and grassland were the two most dominant LULC types following forest   Calculations were conducted using identified concentration curves shown in Figure 1.
and wetland, it is highly probable that the dissolved iron concentrations of grasslands and agricultural lands had lower values than those of forests and wetlands.

Identification of concentration curves for grasslands and agricultural lands
We assumed that the dissolved iron concentration of the grasslands and agricultural lands could also be expressed as a function of TWI.In addition, we presumed that this function could be obtained by multiplying a constant value by the obtained concentration curve.We selected an exponential curve as the concentration curve because the correlation coefficient of the curve was highest.Under these simple assumptions, we introduced two independent constants, b and c, as multipliers for agricultural lands and grasslands, respectively.Moreover, constant a, which is multiplied by the original concentration curve when the LULC is forest, wetland, or other LULC types, was introduced.Multiplying parameters ranging from 0.0 to 2.0 while changing the parameter value at an interval of 0.1 by the concentration curve, a total of 21 × 21 × 21 = 9261 trials were made.The calculation period was 1980 to 1995.The fitness of the calculated result was evaluated by the relative root mean square error (RRMSE).
Figure 3 shows the distribution of the RRMSE in the b-c plane, which represents the intersection of several different values of a.These figures clearly demonstrate that lower RRMSE values were concentrated around the best-fit parameter set, which was as follows: a = 0.8, b = 0.1, c = 0.0, RRMSE = 0.22.These results indicate that the RRMSE surface had only one optimal point in the parameter space and that there was no equifinality problem.The average and variance of parameter sets with good performances ranked in the top 100 were as follows: a = 0.72 ± 0.09, b = 0.17 ± 0.13, c = 0.39 ± 0.28.
By multiplying these average values by the original formula, we obtained new concentration curves for forests/ wetlands, agricultural lands, and grasslands (Figure S3). Figure 4 shows dissolved iron concentration calculated by using the newly developed dissolved iron concentration curves.Most of the discrepancies between the observed and calculated values decreased.In addition, the changing ranges of the observed and calculated values at each site were nearly identical.
Because the dissolved iron data of watersheds dominated by agricultural land or grassland were not obtained, we could not confirm the validity of the curves for agricultural lands and grasslands at this time.However, Yoh et al. (2007) reported that the dissolved iron concentration of dry land and paddy fields was lower than that of wetlands, which indirectly supports the validity of the curve.Annual precipitation on the grasslands for most of the Mongolian high plain is less than 400 mm; thus, it is reasonable to assume that the average dissolved iron concentration of the grasslands in this study area is also low, which supports the validity of the curve for grasslands.

DISCUSSION AND CONCLUSION
This study clarified that TWI can be a good macroscopic index to represent the average dissolved iron concentration of each grid.This means that we can easily assess the dissolved iron concentration of any grid simply by calculating TWI.Since the calculation of TWI only requires DEM data, the calculation procedures are also very simple.Thus, the developed curves will be useful especially for  evaluating the dissolved iron productivity of continentalscale large basins.Moreover, the identified curves can be easily incorporated into a hydrological model as an explicit function.This will open up the possibility of predicting the dissolved iron productivity of any watershed.We therefore consider our results to be a first step toward building a comprehensive terrestrial iron transport model.For the further development of such a model, inter-annual and seasonal changes in dissolved iron should also be formulated.
The key factors governing these temporal changes in dissolved iron are biogeochemical factors such as reducing conditions and the presence of humic substances.
Overall, we obtained a dissolved iron concentration curve for forests and wetlands with respect to TWI.Modifying this function, dissolved iron concentration curves for agricultural lands and grasslands were also identified.The results suggest that TWI can describe the average dissolved iron concentration in areas with different land cover types in the Amur River Basin.Future studies should be conducted to incorporate these concentration curves into a hydrological model to simulate temporal changes in the dissolved iron concentration of the basin.

Figure 1 .
Figure 1.Relationships between mean annual dissolved iron concentrations and average TWI of each basin.Exponentially fitted curves against average, maximum and minimum values are shown.

Figure 2 .
Figure 2. Comparison of the average observed and calculated dissolved iron concentrations during the period between 1980 and 1995 along the main course of the Amur River.Calculations were conducted using identified concentration curves shown in Figure 1.

Figure 3 .
Figure 3. Relative Root Mean Square Error (RRMSE) distributions in the fitting parameter space.Six different planes perpendicular to the a-axis cut at six different a values (a = 0.2, 0.4, 0.6, 0.8, 1.0, 1.5) are shown.

Figure 4 .
Figure 4. Comparison between observed and calculated dissolved iron using identified dissolved iron concentration curves.

Table I .
Pearson's correlation coefficient of the topographic index and the climate index against the average dissolved iron concentration