Accounting for DEM Error in Sea Level Rise Assessment within Riverine Regions; Case Study from the Shatt Al-Arab River Region

Abstract

Global elevation datasets such as the Shuttle Radar Topography Mission (SRTM) digital elevation model (DEM) are the best available terrain data in many parts of the world. Consequently, SRTM is widely used for understanding the risk of coastal inundation due to climate change-induced sea level rise. However, SRTM elevations are prone to error, giving rise to uncertainty in the quality of the inundation projections. This study investigated the error propagation model for the Shatt al-Arab River region (SARR) to understand the impact of DEM error on an inundation model in this sensitive, low-lying coastal region. The analysis involved three stages. First, a multiple regression model, parameterized from the Mississippi River delta region, was used to generate an expected DEM error surface for the SARR. This surface was subtracted from the SRTM DEM for the SARR to adjust it. Second, residuals from this model were simulated for the SARR. Modelled residuals were subtracted from the adjusted SRTM to produce 50 DEM realizations capturing potential elevation variation. Third, the DEM realizations were each used in a geospatial “bathtub” inundation model to estimate flooding area in the region given 1 m of sea level rise. Across all realizations, the area predicted to flood covered about 50% of the entire region, while predicted flooding using the raw SRTM covered only about 28%, indicating substantial underprediction of the affected area when error was not accounted for. This study can be an applicable approach within such environments worldwide.

1. Introduction

Sea level rise (SLR) varies regionally around the world: in the Persian (Arabian) Gulf, the rate is similar to that in south Asia, which is 2.2 ± 0.5 mm/year [1,2]. As the Persian Gulf experiences high temperatures in the summer, the mean sea level is lower in winter and higher in the summer with an inter-annual range of 26 cm [1,3]. As a consequence, tidal wetlands and nearby low elevation areas are extremely sensitive to even slight changes in elevation [4] or sea level. Hydrology plays a key role in understanding the ecology of these environments: predicting and managing changes in coastal environments requires an understanding of the linkages between estuarine wetland ecology and hydrology. The long-term effects of climate change, especially sea level rise and human hydraulic modification of water flow, are also critical for these systems [5]. In the Shatt al-Arab River region (SARR), low-lying areas support important human activities such as large urban systems and agriculture, and the influence of SLR may have a substantial human impact in this region. The SARR has supported the largest cultivated date palm forest in the world, which stretches from the riverbanks towards the desert for distances varying from a few hundred meters to almost six kilometers. A fifth of the world’s palm trees were cultivated in this region in the mid-1970s, amounting to 17–18 million date palms. By 2002, salinization wiped out more than 80% of the SARR palm forest, a total of 14 million palms [6,7]. Environmental problems have been a major focus of several studies in this region, particularly regarding land change due to increased salinity in both the river water and the soil [8,9,10].

Decreasing levels of freshwater in the Shatt Al-Arab River is a primary driver for many problems in the region [11]. There are two main reasons for the decrease of fresh water. The first involves the concentration of Shatt Al-Arab salinity, which increases towards the southern source of salinity in the Persian Gulf. The northern part of the Shatt Al-Arab River is less affected by salinity since it is more distant from the Gulf. A symptom of the second problem involves the vegetated area around the river, which is notably disappearing from the south while some of the riparian vegetation still remain in the northern part of the river [12]. Recently, SLR has played a major role in environmental and vegetation changes, as well as the associated land area changes in this region, and there is a need to assess its impact on the future of the SARR.

In this study, we used a digital elevation model (DEM) product, the Shuttle Radar Topography Mission 1 arc-second global product (SRTM). This DEM was produced from data collected during a US Space Shuttle mission in February 2000 that covered 80% of the global land surface [13]. It is the first global DEM with 30 m resolution, but these DEM data have large vertical errors in densely vegetated and urban areas [14,15,16]. To apply a geospatial inundation analysis, which is widely used to model SLR e.g., [17,18,19,20], while accounting for these errors, error propagation modelling must be applied. We used geostatistical simulations and the Monte Carlo method to understand the effect of SRTM error on inundation likelihood. More specifically, kriging was used to develop linear estimates of parameters for the SRTM error distribution for each cell in the DEM, then that distribution was sampled to generate error, and hence true elevation realizations. We did not have reference elevations for this region, so we used parameters from an error model constructed on the SRTM data in the Mississippi River Delta region (MRDR). The study created correlated, statistical DEM error realizations for the study region to achieve the main aims of this research, which were to:

• Develop an error propagation model to enhance the SRTM 1 Arc-Second Global elevation model for the region.
• Transfer parameters such as vegetation cover and slope for error modelling of the SRTM.
• Apply this model to identify areas that are potentially vulnerable to inundation in the SARR.
• Understand the impact of error on inundation models in low-lying coastal regions.

2. Data and Methods

2.1. Study Area

The SARR is in southeastern Iraq and southwestern Iran and includes land on both sides of the Shatt al-Arab River from north of Basra City to the Persian Gulf, as depicted in Figure 1.

Its total area is 2189 km² [21]. Land elevation along the river descends gradually from five meters above sea level near Al Qurnah in the north to zero in the south. The banks of the river generally have higher elevations than adjacent interior areas; elevation typically ranges from 0 to 5 m in the lowlands away from the river. The Shatt Al-Arab River flows through Iraq for about 80 km from just above the confluence with the Karun River (from the east) to the Persian Gulf; the Shatt Al-Arab forms the border between Iran and Iraq. Cities and towns, which together are home to more than three million people as shown in Figure 1, are interspersed with agricultural land and marshes along its banks. The date-palm groves along the river are naturally irrigated by tidal action [22]. Various studies found that Shatt Al-Arab is distinctive in that it has two tide sequences during the day so that the water level is directly affected by the water in the Persian Gulf [11]. The tidal movement regulates the hydrology of the river, especially since the reduction in freshwater approaching from the Euphrates, Tigris, and Karun Rivers [11] during the phase evaluated in this study (1975 to 2017). Tides come into the river from the Persian Gulf and their effects spread at least as far as Basra city [23]. The region includes a portion of the Governate of Basra in southeastern Iraq and the Abadan County, located in western Iran as shown in Figure 1a. This region is rich in economically important agricultural resources, including crops, vegetables, and palm trees [8,9].

2.2. SRTM 1 Arc-Second Global Coverage (30 m)

In February 2000, interferometric synthetic aperture radar (IFSAR) data were collected by the Shuttle Radar Topography Mission (SRTM) covering 80% of the Earth’s land surface [15]. The data collected brought substantial significance with it, especially in regions with little or no free terrain data at medium to high resolution. SRTM DEMs are not without their issues, of which accuracy is particularly important here [14]. SRTM 1 arc second data were downloaded from the USGS Earth Explorer, then clipped and mosaiced in ArcGIS to acquire a DEM for the SARR. This DEM was then projected to UTM zone 38 and 39 N referenced to WGS 1984 datum; final cell size was 30 m as shown in Figure 1b.

2.3. Vegetation Data

In regional and global studies, vegetation cover fraction (VCF) has been used as a variable in global land process models for earth surface change and climate change assessment. This study also required the use of VCF as a variable since SRTM error is associated with vegetation cover [14,24,25,26,27]. VCF Landsat 5 TM Collection 1 products are available for download from the USGS Earth Explorer; four images, which were acquired on 19 February 2002, were used in this study (Path 165, 166, and Row 39). Since the SRTM resolution was 30 m, the downloaded images were also projected and resampled at 30 m resolution. All of the images also had the same time observation of the Earth as shown in Figure 1d. VCF is defined as the percentage of green vegetation area relative to the total pixel area [24,27,28]. As part of this study, the normalized difference vegetation index (NDVI) image was developed first, and it was determined to have two values, which are bare soil and NDVI value full canopy. The percentage VCF was estimated by the following equation:

$$\mathrm{VCF}=\frac{\left(NDVI-NDV{I}_s\right)}{\left(NDV{I}_c-NDV{I}_b\right)}$$

(1)

where NDVI is the normalized difference vegetation index, NDVI_c is the NDVI for a full canopy pixel, NDVI_b is the vegetation pixel value, and NDVI_S is the bare soil NDVI pixel value [29].

2.4. Methods

2.4.1. Monte Carlo Error Propagation

We used Monte Carlo simulation, which is one of the main approaches to error propagation modelling, as a technique to determine how spatially autocorrelated data errors affect applications using spatial data. The Monte Carlo simulation operation was repeated many times to generate a sample for the input probability distribution [30,31]. The SRTM error for the SARR was expected to have a mean component that was associated with available covariates (SRTM elevation, VCF and SRTM slope) and a spatially autocorrelated, mean-zero residual component. We developed the error model (raw error ~SRTM elevation + VCF + SRTM slope) in a separate study for the Mississippi River delta region (MRDR). We then conducted validation on a separate MRDR dataset to contrast the model-fitted error estimates versus actual SRTM error; finally, we constructed a variogram of the model residuals. We used regression kriging to fit the model and handle its spatially autocorrelated residuals (after [32]):

$$\widehat{z}RK\left(s_0\right)={\displaystyle \sum }_{k=0}^p{\widehat{\beta }}_k\times q_k\left(s_0\right)+{\displaystyle \sum }_{i=1}^n {\omega }_i\left(s_0\right)\times e\left(s_i\right)$$

(2)

where $\widehat{z}$RK(s₀) is the predicted SRTM error, q_k(s₀) is the fitted linear regression of SRTM error using VCF, SRTM elevation, and slope as independent variables, and e(s₀) is the simple kriging interpolated residual at the prediction location s₀. The k represents the regression coefficients, which are estimated using generalized least squares [31], and ${\widehat{\beta }}_k$ represents explanatory variables. The ω_i are kriging weights for the residuals, where e(s_i) is the residual at the measurement location s_i, that is, the difference between the elevation error and the fitted trend at s_i.

To better understand the effect of DEM data uncertainty on inundation modelling, DEM realizations were constructed from the SRTM data and the error model as shown below:

$$DEM real=\mathrm{SRTM}-\mathrm{Regression}\mathrm{predicted}\mathrm{error}-\mathrm{Residual}\mathrm{realization}$$

(3)

The regression kriging predicted SRTM error for SARR was subtracted from the original SRTM data. Unconditional error model residual realizations (their derivation is explained in the next section) were generated and subtracted as well. QQ plots were used to check the normality assumptions of the error model. Predictions for the residual error surface in the SARR were necessarily unconditional, meaning that residual realizations matched spatial properties of the MRDR model but did not rely on ground truth, since we did not have reference elevation data for the SARR. The error realizations were used to enhance the SARR SRTM 30 m DEM, and we applied a SLR inundation model to detect and characterize the flooding impact of SLR for the study region, as shown in Figure 2.

2.4.2. Unconditional Simulation with Convolution Filtering

This study had a large amount of data that was too computationally intensive for standard methods. We employed process convolution [33] to develop the residual simulations in an efficient manner. The convolution filter reproduced a spatial Gaussian covariance function:

$$\rho \left(h\right)=exp\left[-3\left(\frac{h^2}{r^2}\right)\right]$$

(4)

where h is the correlation coefficient and r is the specified range of the DEM error. The convolution filter method worked as follows. We simulated standardized Gaussian noise for pixels throughout a raster. We then calculated a moving window filter using the convolution algorithm proposed by Oksanen and Sarjakoski [33] to reproduce the desired Gaussian covariance function; an example is shown in Figure 3.

Finally, this filter was passed across the random noise raster to generate a realization with the desired spatial structure. In this study, the geostatistical modelling of the error elevation was accomplished by applying the R statistical package; the inundation model used the following libraries: raster, grdevices, and gstat [34,35,36].

2.4.3. Inundation Modelling

Many previous studies on geospatial inundation modelling (bathtub modelling) used either high resolution LIDAR or SRTM with error modelling [18,31,37,38]. In SARR, we did not have reference elevation data, so this study used the error propagation model from the MRDR and applied it to the SARR. We employed the sea level rise of 1 m and applied the inundation model to the area. The inundation model function is:

$$Floode{d}_j=\{\begin{array}{c}2\\ 1\\ 0\end{array} \begin{array}{c}if DE{M}_j>f\\ if f\ge DE{M}_j\ge 0\\ if DE{M}_j\le 0\end{array}$$

(5)

where Flooded_j represents cell j in the DEM, the cell elevation value for j is represented by DEM_j, and f is the inundation elevation. The inundation model was applied to every DEM realization, resulting in dozens of output maps. The distribution of inundated area across these maps was calculated and portrayed in a histogram. Finally, the probability of inundation for each cell was computed and a map of inundation probability for the SARR was produced.

3. Results

3.1. Regression Kriging-Based Error Modelling

The error model was developed and validated using data from the MRDR and transferred to the SARR. Validation in the MRDR was positive: the regression model showed high significance and explanatory power, and the residuals were spatially autocorrelated.

The variogram model of SRTM errors was a Gaussian model with a range = 800, partial sill = 0.1, and nugget = 0.04, as shown in Figure 4a and

Table 1. Regression coefficients for regression model variables.

Parameter	Estimate	Std. Error	t Value	Pr (>|t|)
Intercept	−0.1382	0.01884	−7.337	2.36 × 10−13
SRTM Elevation	0.8008	0.0050	159.976	<2 × 10−16
VCF	0.0054	0.0007	7.991	1.49 × 10−15
SRTM Slope	−0.0803	0.0104	−7.723	1.24 × 10−14

.

We transferred the error model to the SARR to model SRTM error for this region and produce realizations of the potential DEM surface, accounting for uncertainty in SRTM associated with vegetation, slope, elevation values, and spatially autocorrelated residuals. Figure 4b illustrates the predicted SRTM error surface from the regression model for MRDR, while Figure 4d shows the predicted error surface for the SARR. Histograms of both are shown in Figure 4c,e, illustrating that these are centered just above zero with a bell-shaped distribution with long tails.

3.2. Realizations

A sample realization generated following the approach from Section 2.4.1 is depicted in Figure 5.

The QQ plot for an initial Gaussian random noise realization as shown in Figure 5a,b was straight and, therefore, matched expectations, while that for the filtered realization was not as normal, but its shorter tails were acceptable. We generated 50 realizations of error model residual maps and used the equation above to produce DEM realizations for the SARR; one of these is depicted in Figure 5c.

3.3. Inundation Modelling

Consistent to the literature and according to the 1 m scenario in this study, the spatial extent of inundation for a 1 m sea level rise in the SARR was calculated using the original SRTM DEM and is mapped in Figure 6c. Compound model results, developed by running the inundation model on all fifty realizations, are mapped in Figure 6e. In this comparison, each map shows three classes: water, where areas are at or below 0 m of elevation; land, where areas are above 1 m; and flooded, where areas are between 0 and 1 m due to being inundated by this amount of flooding. While only 28% of the land surface was indicated as flooded using the raw SRTM, the DEM realizations showed a substantial difference with nearly 50% of the land area being flooded across all realizations. Figure 6d,f depict histograms of the proportion of land area that remained dry for these realizations. Figure 6f also shows the raw SRTM data depicted by the vertical dotted line positioned well to the right of any bars in the histogram.

4. Discussion

In this study, we used the inaccurate but globally available SRTM as the base elevation data for SARR and transferred a SRTM regression-kriging error model from the MRDR. This model could account for SRTM error caused by vegetation and topographic factors, while also characterizing the spatial autocorrelation of the regression residuals. We applied a standard geospatial bathtub model to multiple realizations of this elevation model to assess and visualize the impact of one meter of SLR in the SARR. The key principle was to employ regression-kriging modelling with unconditional simulation of the residual process via convolution filtering to create and propagate the impact of SRTM error to the geospatial inundation model in the SARR. To ensure robust results, 50 realizations of the elevation surface were generated. This study is the first one to transfer an error propagation model for bathtub modelling between different places in the world and to evaluate the effects of SRTM error on SLR inundation of near-coastal riverine regions. We consider that the use of the MRDR error model in SARR is justified for several reasons:

• SRTM 30 m is a global DEM product that covers both regions. We assume that, because the sensors and algorithms that produced it are the same, the error properties in both areas should be similar.
• The topography of both regions is similar. Both regions have a low elevation of between 0 to 35 m. Elevation is relatively high on the levees along the river banks and lower away from the rivers. The northwestern parts of both regions are higher than the southeastern, and both end in marine gulfs.
• We can acquire the parameters to explain SRTM error in both regions such as vegetation cover and the slope from SRTM itself and from Landsat imagery.
• We have surface reference elevation data for the MRDR, but we do not have any reference elevation data for the SARR. Geostatistical simulation is computationally expensive, especially over large regions. In this study, a Gaussian convolution filter was applied to the random noise raster, which rapidly generated spatially autocorrelated DEM error realizations.

The results showed that the methodology improved simulation efficiency while guaranteeing a relative normality of residuals in each simulation. The comparison of geospatial inundation modelling for raw SRTM 30 m and adjusted SRTM 30 m, as shown in Figure 6f for a 1 m SLR, and the impact of SRTM error on flooding estimation is stark. The raw SRTM flooded around 28% of the adjacent land, whereas the adjusted SRTM indicated that more than 45% would be flooded. This finding has potential implications for other studies that use SRTM to estimate SLR. Just one previous study by [18] considered the effect of SRTM error on inundation modelling. That study used SRTM to assess the performance of NOAA’s Global Land 1-km Base Elevation Project (GLOBE) in the US but failed to reach any definitive conclusion because of the spatial heterogeneity in SRTM 30 m modelling quality [18].

Limitations

A major limitation for conducting SLR modelling in the SARR is the lack of local data. In comparison, the MRDR has high-resolution, high-accuracy 10 m DEMs derived from LIDAR, and tide station gauges in the MRDR provide regular sea-level elevations in this region at high temporal frequency over many years [39]. The National Elevation Data (NED) is known to have random errors when it comes to absolute vertical accuracy statistics. Varying the consequences for vertical error in the elevation data causes a difference of 1.64 m within the near-shore topography and bathymetry, particularly in marshes and wetlands [40,41,42]. However, this study focused on the actual land area around the SARR, not the shorelines or the transition point between land and water. Thus, this study did not conduct its own error analysis, rather limiting the analysis to the generation of best estimates. In contrast, in the SARR several studies have been concerned with the intrusion of seawater in the region and how the river became salty because of a reduced fresh water input from upstream, local demand, and seawater intrusion. However, there is no reliable, long-term gauge or station to provide the data we can use to measure SLR, or a high-quality, freely available DEM. This study was the first to consider SLR impact in the SARR, where the best freely available DEM is SRTM, a global data product with substantial error. After developing the regression-kriging error model for the MRDR, we transferred this model to the SARR. The residual error surfaces for each realization were unconditional, and residual realizations matched properties of the model but did not have local reference data for validation. We must assume that error properties truly are similar between the two regions; if this is not the case, then the effects of error on the inundation model are mis-specified. Second, the resulting surface in the regression kriging, such as any interpolation procedure, does not represent slope and aspect. In both regions, the surface is smoother than depicted, as slopes and aspect are poorly represented.

Another limitation of the study is that the methods described in this study can only be applied in coastal riverine regions. Therefore, any differences between enforcement methods in the propagation regression-kriging modelling can only be said to be transferable to other areas with similar geographical characteristics, especially topography and vegetation cover. Finally, there are multiple SLR models available for the data-rich MRDR, including more formal hydrological modelling [43,44,45], but this study is limited to bathtub inundation modelling in order provide a fair comparison between the MRDR and the SARR.

5. Conclusions

The results presented in this paper furnished valuable information on the error propagation modelling of SRTM 30 m data to characterize sea level rise in the SARR. It is the first to apply error propagation modelling in this context. Model output included a probability of inundation for each location in the SARR, and provided a comparison of total inundated area. The approach was effective for determining and mapping the land that would be inundated in the region from a 1 m rise in sea level based on bathtub modelling. Most critically, in the SARR, error propagation modelling found much larger areas were at risk of flooding than would be estimated with raw SRTM. SRTM, with its global coverage, can be used to model SLR if error propagation modelling is applied properly. The approach will identify reasonable flooding risk for equivalent regions facing SLR. This research provides an appropriate example of the successful application of transferring a geostatistical error model for error propagation modelling to assess SLR inundation on near-coastal riverine regions.