Abstract

Hypsometric curve and hypsometric integral are important watershed health indicators that can used be during the assessment of erosion status of a watershed. The recent advancements in Geographic Information Systems and the availability of LIDAR data facilitating the production of accurate Digital Elevation Models have made the calculation of morphometric indices a tractable process and this allows analysts to readily integrate hypsometric analysis in the landscape planning and analysis studies.

Despite this easy calculation of morphometric indices, the determination of the hypsometric integral for a watershed remains a relatively more protracted process due to its iterative nature and this hampers its use in watershed analysis. Many software tools have been developed to help in this process of calculating such hypsometric indices for a single catchment, but all fall short of facilitating such calculation for sub-catchments within a watershed and this is a required workflow during catchment prioritization studies.

hypsoLoop has therefore been developed as tool to automate the drawing of hypsographic curves and the calculation of hypsometric curves for a watershed that contains several sub-catchments. It has been tested on Yanze watershed, a small basin in central Rwanda, and was found to provide an accurate estimation of hypsometric integrals with a simple and accessible interface that can also be easily integrated into larger watershed analysis workflows.

Installation

You can install the development version of hypsoLoop from GitHub with:

# devtools::install_github("fgashakamba/hypsoLoop")

and you can install the published beta version at CRAN with:

# Still in beta stage (version 0.1.0)
# install.packages("hypsoLoop")

1.Introduction

Morphometry - defined as a formal or mathematical analysis of the configuration of the earth’s surface, shape, and dimension of its landforms (Altaf et al., 2013) – is powerful tool used by geomorphologists, hydrologists, and environmentalists in general to characterize a landform and predict its response to external factors such as precipitation, land cover change or farming.

Hypsometry is the most important aspect of morphometry and its analysis constitutes the starting point of any study aimed at characterizing a watershed, especially the studies about the effects of erosion in a landscape. A watershed is the most appropriate unit to study for several processes of land surface and the definition of morphometric indices at the watershed level, rather than at individual streams level is more natural and useful in revealing relevant patterns.

Hypsometric analysis of a watershed (area-elevation analysis) has generally been used to reveal the stages of geomorphic development (stabilized, mature and young) and it describes the elevation distribution across an area of land surface. The two most considered metrics during hypsometric analysis are the “Hypsometric Curve” (also called the hypsographic curve) and the “Hypsometric Integral”. The former is a graph that shows the proportion of surface area of earth at various elevations by plotting relative area against relative elevation. The non-dimensional hypsometric curve is the standardized form usually used by hydrologists and it is drawn by scaling elevation and area by the maximum values. The Hypsometric Integral is the area beneath the hypsometric curve whose value varies from 0 to 1 (S. K. Sharma et al., 2018). In terms of its interpretation, given a hypothetical Hypsometric Integral of 0.35 for example, that would mean that only 35 per cent of the land masses in that basin remain to be eroded.

In addition, the valuable suggestion was made that, under certain conditions, this It has been argued that the hypsometric integral can provide a quantitative expression of “stage” of basin denudation (Chorley & Morley, 2016) and therefore is important for estimation of erosion status of watershed. Thus, it is used during the prioritization of watersheds for proposing soil and water conservation activities (S. Sharma & Mahajan, 2020). the hypsometric integral also provides a simple morphological index that can be used in surface runoff and sediment yield prediction from watersheds (S. K. Sharma et al., 2018; S. Sharma & Mahajan, 2020).

Given the dimension-less nature of the above-mentioned hypsometric parameters, their use in the characterization of a watershed has an added advantage of allowing comparison between the different sub-watersheds of the basin irrespective of their true absolute scale (Singh et al., 2008).

The measurement techniques used to calculate the morphometric parameters of landforms have evolved over the years but it’s the recent advancements in and widespread use of Geographic Information Systems (GIS) and the availability of LIDAR data facilitating the production of accurate Digital Elevation Models (DEMs) that made the calculation of morphometric indices a tractable process and helped establish geomorphometry as the science of quantitative land-surface analysis.

Notwithstanding this easy calculation of morphometric indices, the determination of the hypsometric integral for a watershed remains a relatively more protracted process due to its iterative nature and this hampers its use in watershed analysis. There are many tools that have been developed to automate this process and make it amenable to its integration into general watershed analysis workflows. However, most of these tools, among other things, fall short of facilitating such calculation for sub-catchments within a watershed which is the required workflow during catchment prioritization studies.

This project therefore aimed at developing a tool to be used to automate the drawing of hypsographic curves and the calculation of hypsometric curves for a watershed that contains several sub-catchments. The tool was extensively tested on Yanze watershed, a small watershed located in the north-west of Kigali city, Rwanda. For comparison purposes, the tool was also applied to the whole country of Rwanda to demonstrate the assumption that the dimensionless nature of hypsometric integral allows for comparison between the different sub-watersheds of the basin irrespective of their true absolute scale.

2. Problem statement

There exist several software tools that can be used during the process of calculating hypsometric curves and integrals. However, for studies related to watershed analysis as part of watershed management planning, these tools have gaps in terms of (1) not being simple enough for use by non-specialist practitioners, (2) not being open-source for the purpose of transparency and ease of modification, (3) being provided as part of larger programs (sometimes expensive packages) which becomes a bloated solution for a simple problem, and (4) not being automated enough for easy integration into watershed analysis workflows conducted in an accessible environment such as R.

Before planning to develop hypsoLoop, 5 tools were tested and found to be lacking in one or more of the above-mentioned areas. For instance, calHypso is a relatively easy-to-use tool developed by J.V. Pérez-Pena et.al (Pérez-Peña et al., 2009) as part of the MAT project (Gudowicz & Paluszkiewicz, 2021). It is offered as an extension of ESRI’s ArcGIS software. However we found its usefulness limited in the way that its reliance on the ArcMap environment, a commercial, non-open source package, makes it not accessible to everybody.

At the other end of the spectrum, the hydroTSM tool (Mauricio Zambrano-Bigiarini, n.d.) also can be used to draw hypsometric curves and, as an R package, is accessible to everybody and transparent in terms of the algorithms it uses. However, being developed as a general tool for hydrologists thus incorporating many other functions for time series management and hydrological modelling for instance, its hypsometric analysis function is rough and lacks many useful ingredients such as good visual presentation of the curve plot and proper labeling. In addition, its general purpose nature makes it bloated and clutters the working environment of the analyst with functions that are not needed.

In the middle ground, we find tools such as the Hypsometry Toolbox (a Python script) (Davis, n.d.) and the Hypsometric Integral Toolbox for ArcGIS (Matos, A.; Dilts, 2019) that are fine but lack the required iterative function to process all the sub-catchments at once.

Finally, some other analysts have tried to leverage the functionality of some tools that were not specifically designed for hypsometric analysis but can generate the watershed elevation tables that can be fed into a statistical package to generate the curve and the integral. This approach is also so cumbersome and not scalable. hypsoLoop was therefore developed to close that gap by providing a way to draw hypsographic curves and calculate hypsometric integrals of a watershed with multiple sub-catchments in an interactive and customizable manner. The tool is used in R, an open-source data analysis environment of choice for researchers which ensures it can easily be integrated into further watershed analysis workflows.

3. Materials and Methods

3.1. Input data considerations

hypsoLoop relies heavily on the raster package for generating elevation tables required for drawing hypsographic curves and calculation of hypsometric integrals.

The two datasets required as inputs are:

A Digital Elevation Model (DEM) covering the study area. The DEM can be supplied as a raster object from the raster, terra, or stars packages. A note should be taken of the coordinate reference system used by the DEM because it will affect the accuracy of area calculation and therefore influence the results of the estimation of hypsometric integrals. In general, planar coordinate reference systems are preferred over angular ones (lon/lat). Notwithstanding this recommendation, appropriate steps have been taken to minimize issues during area calculation regardless of which type of reference coordinate system is used. Furthermore, hypsoLoop assumes the values in the input DEM are in meters (m) and the output elevation tables give contour areas in hectares. In this regard, the analyst should ensure appropriate unit transformations on the DEM are done beforehand. Finally, the preprocessing of the DEM dataset should include filling all the sinks and leveling spurious peaks therein; otherwise, the interpretation of the resulting elevation tables will not be easily interpretable.
Delineation of the sub-catchments in the area of interest. The boundaries of sub-catchments in study area are given as a spatial vector object of class sf dataframe or SpatialPolygonsDataframe (from sp package). This object can be created from an existing shapefile or GeoJSON data using st_read() function of sf package or readOGR() function of rgdal package. When preparing the layer, a column called Name should be included to hold the name or code of each sub-catchment. The user should also note that it is recommended that the DEM and the sub-catchments object be in the same coordinate reference system. If not, the sub-catchments object will be re-projected to the coordinate reference system used by the DEM.

3.2. Description of the algorithm

hypsoLoop works by first producing elevation tables of each sub-catchment in the input sub-catchments boundaries input data. These tables are generated by extracting the DEM of each sub-catchment and reclassifying it into 30 equal elevation range classes and then calculating the area covered by each range class.

The hypsometric integrals are then calculated by fitting a 3^rd degree polynomial equation to the normalized contour area-elevation values and integrating it within the limits of 1 and 0 (S. K. Sharma et al., 2018). The choice of 3^rd degree as the polynomial fit was based on the prevalent S-shape profile of hypsometric curves of watersheds at inequilibrium (young) and equilibrium (mature) stages. This approximation doesn’t necessarily work for watersheds at old (monadnock) stage (Vanderwaal & Ssegane, 2013). Therefore, the coefficients of the fitted equation should closely be examined to ensure the right equation is being fitted.

3.3. Structure of the tool

Two functions do the heavy-lifting work of hypsoLoop. First, the generateHypsoTables() function does the extraction of DEM for each sub-catchment and generates an elevation table for each sub-catchment while also compiling the minimum and maximum elevation values for each sub-catchment in another separate table.

For input data sanitization, a utility function called check_arguments() is used and it checks whether the provided data are of the right types. If necessary, and the DEM is converted into a raster object (raster package) and the vector layer is converted into an sf dataframe (sf package).

The function then checks whether the two datasets use the same coordinate reference system and re-projects the vector layer if it’s not the case. For the calculation of areas, another utility function called calc_area() is called and it uses the appropriate area calculation method depending on whether the coordinate reference system used by the provided DEM is planar or angular.

After this data sanitization step, the process to draw the hypsometric curves and calculate the hypsometric integrals for each sub-catchment is done through the following steps:

Normalization of the area-elevation tables: For elevation, the minimum elevation is subtracted from each contour elevation and then divided by the range. For area, the cumulative area at each contour starting from the top contour is divided by the total area of the sub-catchment.

$$normalizedElevation=\frac{countourElevation-minElevation}{maxElevation-minElevation}$$ $$normalizedArea=\frac{cumulativeSumArea}{totalArea}$$

Plotting of the hypsometric curve: Using ggplot2 package, the normalized data is plotted on a 2D chart showing the relative area on the X-axis and the relative elevation on the Y-axis.
Calculation of the hypsometric integral: Using PolynomF package, a 3rd degree polynomial equation is fitted to the data and is integrated within the limits of 0 and 1.
Outputting the summary table and charts: Each of the generated charts along with a summary table are saved in a folder named “HYPSO_OUTPUT” located on the working directory

Pseudo-code flowchart of hypsoLoop

3.4. Case Study: Yanze watershed

Yanze is a relatively small watershed located in north-west of the City of Kigali, Rwanda. It is part of a bigger Nyabugogo watershed (Figure 1) which was identified by the Rwanda Water Resources Board as one of the five most important watersheds in the country (Rob et al., 2018).

The watershed is composed of 3 main sub-catchments namely Cyonyonyo (3667Ha), Mulindi (2572Ha), and Yanze Downstream (3442Ha); and together they total an area of 9,681Ha (Figure 2.a). The elevation in Yanze watershed ranges between 1370m and 2225m above sea level and it forms a depression between two high elevation zones namely Mount Jali and Shyorongi highland on the eastern and western sides of the watershed respectively. The watershed is characterized by steep slopes on hillsides separated by V-shaped valleys with an average slope of 38.7% (Figure 2.b). The drainage pattern in Yanze watershed is characterized by dendritic streams with a three-order permanent stream flow system. The main tributaries, Cyonyonyo and Mulindi, feed into Yanze stream and they are themselves fed by streams that take source in the surrounding hillsides. These third-order streams include Nyakabingo, Ruhonwa, Ntakaro, Munyarwanda, and Kinywamagana.

Map of Yanze watersheds and its location in Rwanda

4. Results and discussion

The drawHypsoCurves function was run with the DEM of Yanze watershed and the level-4 yanze catchments vector dataset as the input arguments.

The following is the summary table representing the results of the hypsometric analysis of the watershed:

Sub Catchment Name	Minimum Elevation (m)	Maximum Elevation (m)	Total Area (Ha)	Hypsometric Integral
Mulindi	1,616	2,225	2,571.53	0.736
Cyonyonyo	1,614	2,217	3,666.87	0.791
Yanze	1,368	2,071	3,431.98	0.857

This means that all of Yanze sub-catchments are still at the young stage of their geomorphological development and therefore have a potential to undergo intensive erosion if not enough soil conservation measures are undertaken.

Hypsographic curves of Yanze’s level-4 sub-catchments

As it can be observed in Figure 3, the hypsometric curves of both Mulindi and Cyonyonyo catchments are clearly convex upward which normally indicates the youth stage and hints to the high susceptibility to erosion and land sliding. On the contrary, the downstream catchment boasts a characteristically S-shaped curve that is concave upward at high elevations and convex downward at low elevations. This normally indicates a watershed that is at a mature stage of its geomorphologic development which is understandable given the position of much of that sub-catchment at the peneplain of the watershed.

Regarding the suitability of using the 3^rd degree polynomial integration approach when calculating the hypsometric integral, hypsoLoop gives estimates of the significance of the fitted coefficients. If the p-value of any of the coefficients is too big, then the user is advised not to use this this tool for that particular catchment. For Yanze, the p-values of estimates of the coefficients fitted for all the sub-catchments were within the acceptable limits with the coefficient of determination (R-squared) values of 0.998 for mulindi and cyonyonyo sub-catchments and 0.991 for the Yanze downstream sub-catchment.

Coefficients fitted by hypsoLoop during the calculation of hypsometric integrals of Yanze sub-catchments

5. Conclusion

Hypsometric analysis of watersheds expresses the complexity of denudational processes and the rate of morphological changes. Therefore, it is useful to comprehend the erosion status of watersheds and prioritize them for undertaking soil and water conservation measures.

Among the various methods used to determine the hypsometric integral, the integration of the polynomial equation fitted to the hypsographic curve of the watershed is considered to be the best given its closeness to the very definition of hypsometric integral. However, this method is time consuming and necessitates mathematical integration procedures and subsequent calculation within the desired limits of hypsographic curve.

hypsoLoop has been developed to automate this process and to allow full hypsometric analysis of a watershed with multiple sub-catchments (hundreds if needed) in one go. The application of this tool to Yanze watershed revealed some interesting insights where it was shown that the two upstream sub-catchments have convex upward curves typical to catchments that are at a young stage, whereas the curve of the sub-catchment located in the downstream part of the watershed had the S-shaped profile typical to catchments at mature or equilibrium stage.

6. References

[1] Altaf, F., Meraj, G. and Romshoo, S.A. 2013. Morphometric Analysis to Infer Hydrological Behaviour of Lidder Watershed, Western Himalaya, India. Geography Journal. 2013, May 2014 (2013), 1–14. https://doi.org/10.1155/2013/178021.

[2] Chorley, R.J. and Morley, L.S.D. 2016. A Simplified Approximation for the Hypsometric Integral. The University of Chicago Press. 67, 5 (2016), 566–571. https://doi.org/10.1086/626608.

[3] Davis, J. Hypsometry Toolbox.

[4] Gudowicz, J. and Paluszkiewicz, R. 2021. MAT: GIS-based morphometry assessment tools for concave landforms. Remote Sensing (MDPI). 13, 14 (2021). https://doi.org/10.3390/rs13142810.

[5] Hypsometric Integral Toolbox for ArcGIS: 2019. https://www.arcgis.com/home/item.html?id=23a2dd9d127f41c195628457187d4a54. Accessed: 2021-09-23.

[6] Mauricio Zambrano-Bigiarini hydroTSM: Time Series Management, Analysis and Interpolation for Hydrological Modelling.

[7] Pérez-Peña, J. V., Azañón, J.M. and Azor, A. 2009. CalHypso: An ArcGIS extension to calculate hypsometric curves and their statistical moments. Applications to drainage basin analysis in SE Spain. Computers and Geosciences. 35, 6 (2009), 1214–1223. https://doi.org/10.1016/j.cageo.2008.06.006.

[8] Rob, N., Benon, Z. and Ebel, S. 2018. IWRM Programme Rwanda, Catchment Plan of Nyabugogo Watershed 2017-2023 – Version 1.0.

[9] Sharma, S. and Mahajan, A.K. 2020. GIS-based sub-watershed prioritization through morphometric analysis in the outer Himalayan region of India. Applied Water Science. 10, 7 (2020), 1–11. https://doi.org/10.1007/s13201-020-01243-x.

[10] Sharma, S.K., Gajbhiye, S., Tignath, S. and Patil, R.J. 2018. Hypsometric Analysis for Assessing Erosion Status of Watershed Using Geographical Information System. (2018), 263–276. https://doi.org/10.1007/978-981-10-5801-1_19.

[11] Singh, O., Sarangi, A. and Sharma, M.C. 2008. Hypsometric integral estimation methods and its relevance on erosion status of North-Western Lesser Himalayan watersheds. Water Resources Management. 22, 11 (2008), 1545–1560. https://doi.org/10.1007/s11269-008-9242-z.

[12] Vanderwaal, J.A. and Ssegane, H. 2013. Do Polynomials Adequately Describe the Hypsometry of Monadnock Phase Watersheds? Journal of the American Water Resources Association. 49, 6 (2013), 1485–1495. https://doi.org/10.1111/jawr.12089.

hypsoLoop: A tool used to conduct hypsometric analysis of a watershed