The indicative meaning calculator – quantification of paleo sea-level relationships by using global wave and tide datasets

Abstract

The study of past sea levels relies largely on the interpretation and quantification of sea-level indicators. These are fossil coastal landforms, bioconstructions or deposits that have a quantifiable relation to paleo sea level (called the indicative meaning) and can be assigned an age of formation. The calculation of the indicative meaning should always rely on the quantitative comparison between the paleo sea-level indicator and the elevational range of the same feature in the modern environment with respect to modern sea level. This may prove difficult to quantify when no site-specific data is available, for example when compiling databases with large geographic scope. In this paper, we provide a method to quantify the indicative meaning of sea-level indicators using simple hydro- and morphodynamic equations with inputs from global wave and tide datasets. We apply and compare our approach exemplary to data published in a relative sea-level database on the Marine Isotope Stage 5e highstand and show how this approach can help in estimating the indicative meaning if no site-specific data on the modern analog is available. For the easier usability of the described approach in this paper, we also present a stand-alone Java-app (IMCalc), which can be used as a calculator of the indicative meaning.

Introduction

Past interglacials are one of the main targets of paleoclimate science, as the climate patterns that characterized them can serve as partial analogs for future climate conditions. In particular, the generally higher-than-today global and polar temperatures that characterized these periods induced smaller ice sheets and subsequently higher sea levels [1, 3, 7, 8]. Constraining global sea-level patterns during such warming climates can help improving models predicting the extent of future ice-sheets and consequently modern sea-level rise [1].

Every feature that has a quantifiable relation to the paleo sea level during the time of their formation (e.g. fossil coastal landforms, bioconstructions or deposits), and for which an age can be established, can be used as relative sea-level (RSL) indicator [18, 23]. RSL indicators are the only direct proxies that can be used to assess paleo sea levels. A typical example of RSL indicators are coral reef terraces that formed in direct connection with sea level in the past, and are preserved today as fossil reef terraces (e.g., [20]). It is possible to measure the modern elevation of a fossil reef terrace and establish its age through U-series dating of the fossil corals collected from it [6]. In order to quantify the paleo RSL associated with the reef terrace, it is additionally necessary to quantify the water depth at which it was located during the time of formation. This relation is called the indicative meaning (IM) [18, 23] and has to be established for each site reporting a sea-level indicator. The IM is an important measure to describe the vertical offset to a common tide level, such as the mean tide level (MTL), as only very few RSL indicators form exactly at this point.

The most rigorous approach to establish the indicative meaning of a RSL indicator is to measure a nearby modern analog of the indicator, and transfer the elevation range at which the modern analog is occurring today to the paleo environment (see [15], for an illustrated example applied to the Last Interglacial). Another approach relies on the interpretation of facies or biological elements (e.g., coral species and their living depth range, [6]) within the RSL indicator to reconstruct the paleo environment, and hence the location of paleo RSL. In a focused study on Last Interglacial (Marine Isotope Stage (MIS) 5e, ~ 130–116 ka) beach deposits in the island of Mallorca, Lorscheid et al. [11] showed that the indicative meaning of RSL indicators can be calculated using hydrodynamic models and data on local wave conditions. It is important to note that the concept of the IM is based on the assumption that the same hydro- and morphodynamic processes that dictated their formation in the past are the same or at least very similar to the processes shaping their modern analogs. It is likely that the uncertainties increase for older time scales, such as the MIS 5e highstand, nevertheless, the concept of the IM can still be applied, but for certain limits, e.g. tidal levels, more specific investigations on the changes over long time spans have to be undertaken.

In this study, we apply the same concept adopted by Lorscheid et al. [11] to the most common types of RSL indicators used to study the sea level during past interglacials, as listed by Rovere et al. [15]. We use global wave and tide datasets as input to a series of simple hydro- and morphodynamic equations in order to calculate the indicative meaning for several types of RSL indicators at any point along the global shorelines. We make our proposed method available for other researchers by presenting an operating-system independent software tool for calculating these values. We propose that this methodology and software can be used in absence of site-specific data on modern analogs to gain a first-order estimate on the local indicative meaning.

Methods

Rationale and coastal points database

The main aim of this study was to develop a methodology to calculate, at any point along the global coastlines, the indicative meaning of the seven most common RSL indicator types as described by Rovere et al. [15]. For digitizing the global shorelines, we used the Level-1-shapefile of the ‘Global Self-consistent, Hierarchical, High-resolution Geography Database’ in its full resolution (GSHHG_f_L1, [24]). These polygons of land masses were converted into vector lines and those again sampled in a regular distance of 5 km by using the ‘Generate Points Along Lines’ tool (ESRI ArcMap), resulting in a point-shapefile representing the global shorelines. This sampling distance should ensure a sufficient data distribution, while keeping a relatively small file size. Therefore, coastal features below 5 km in size might not be fully captured. In case such a feature is present in the study area, a more detailed study on the effects on the IM should be undertaken. Only coastline-points between 60°N and 60°S were retained in order to exclude Polar Regions, where RSL indicators are rarely preserved [12] and wave data is often not available (Fig. 1). To each point, we assigned values on tides and waves by using data as described below. This collection of points is called hereafter ‘Coastal Points Database’ (CPD).

Fig. 1

Global map of the Coastal Points Database (CPD). Blue coastlines represent locations included in the CPD, where tide and wave conditions were extracted and where the local indicative meaning can be calculated. Points where no data could be extracted is also shown here (red coastlines), but are not included in the CPD. Points above 60°N or below 60°S were removed from the CPD, in order to exclude polar regions. Additionally, locations of the K09 database ([8], yellow pentagons) are displayed

Tide and wave datasets

In this study, tidal levels were quantified using a three-stepped approach. First, water levels were modelled on a global 0.5 arc-degree grid (~ 55 km latitudinal distance) with the OSU Tidal Prediction Software (OTPS, [4, 5]) and the Tidal Model Driver toolbox for a period of 19-years (the time covering an entire lunar cycle) from 1998 to 2017. Second, we used common mathematical definitions [13, 18] to calculate seven tidal datums from these time series: the Mean Higher High Water (MHHW, Fig. 2c) and the Mean Lower Low Water (MLLW, Fig. 2d), the Highest (HAT) and Lowest (LAT) Astronomical Tide, the Mean High (MHW) and Low (MLW) Water, as well as the Mean Tide Level (MTL). Finally, we used a bilinear interpolation (ESRI ArcMap tool ‘Add Surface Information’) to assign the values of all seven datums to the CPD.

Fig. 2

Global grids of the wave conditions and tidal datums. a Mean significant wave height (Hs); b Mean wave period (Tp); c Mean Higher High Water (MHHW); d Mean Lower Low Water (MLLW)

Wave characteristics were downloaded from the CAWCR (Collaboration for Australian Weather and Climate Research) wave hindcast [2]. This dataset is based on the NOAA WaveWatch III wave model [21, 22] and the NCEP CFSR surface winds and sea ice data [16]. Compared to the NOAA WaveWatch III dataset [21, 22], the CAWCR hindcast has a finer grid size (0.4°/~ 44 km), a longer record of data (January 1979 to June 2016) and includes also the Mediterranean Sea and the Persian Gulf. The downloaded data spans the entire available period (except for January 1979, which represents a model spin-up month) on a global grid. From the provided hindcast data, the hourly significant wave height (Hs) and wave period (Tp) were extracted (Fig. 2a and b, respectively). For both attributes, we calculated four raster datasets corresponding to the arithmetic mean and the standard deviation, as well as the maximum and minimum value for each grid point. From the eight resulting grids, we extracted values for each CPD-point using the same method as described for the tidal values.

Calculation of limits determining the indicative range

The classic definition of the indicative meaning [18, 23] foresees that it is determined by two values: the indicative range (IR) and the reference water level (RWL). The IR measures the uncertainty associated with the modern analog and is defined as the range between upper and lower limit of formation, whereas the RWL represents the mathematically averaged position of sea level, i.e. the midpoint of the IR [17, 18, 23]:

$$ IR= upper\ limit- lower\ limit $$

(1)

$$ RWL=\frac{upper\ limit+ lower\ limit}{2} $$

(2)

Together with the modern measured elevation and measurement error of the indicator, these two values are used to calculate the paleo RSL and the associated uncertainty [15, 17, 18, 23]:

$$ Paleo\ RSL= elevation\ measurment- RWL $$

(3)

$$ Paleo\ RSL\ error=\sqrt{{\left(\frac{IR}{2}\right)}^2+{\left(\frac{measurement\ precision}{2}\right)}^2} $$

(4)

Please note that for a full estimation of the paleo RSL error also other uncertainties have to be included. Especially changes in the environmental conditions, e.g. tidal range or wave energy, or post-depositional effects like sediment compaction, need to be included to reduce the associated uncertainties. Here we only calculate the paleo RSL error using the IR and measurement error as these are in most cases the main influences on the uncertainty.

Following these formulas the values for the upper and lower limits of formation of the RSL indicator need to be quantified in order to determine the IM and further the paleo RSL. Note that it is important that both the RWL and the altitude measurement need to be referenced to the same datum in order to ensure the correct IM and paleo sea-level indication. Rovere et al. [15] already defined qualitatively the geomorphological limits of ten sea-level indicators, from which we picked the seven most common (Table 1) for being used in the IMCalc software.

Table 1 List of sea-level indicators and their descriptive limits (modified from [15]) as used in IMCalc

In order to quantify these values, we use the following formulas and relationships. The breaking depth of waves db, i.e. the lowest point of interaction between sediment and water, was calculated with the formula from Rattanapitikon et al. [14], as reported in Lee and Mizutani [9]:

$$ \frac{db}{L_0}=\left(3.86{s}^2-1.98s+0.88\right){\left(\frac{Hs}{L_0}\right)}^{0.84}\kern0.5em $$

(5)

In this equation, s is the beach slope, Hs the mean deepwater significant wave height (extracted from the CAWCR dataset) and L₀ is the mean deepwater wave length. For this study, the slope of beaches s was calculated averaging the slope of beaches reported in a global compilation by Liu et al. [10], resulting in an average slope of 0.08. In the software, this value is only used if no site-specific value for the beach slope is entered. The deepwater wave length L₀ was calculated using the formula proposed by Stockdon et al. [19]:

$$ {L}_0=\frac{\mathrm{g}\ast {\mathrm{Tp}}^2}{2\ast \uppi}\kern0.75em $$

(6)

with the gravity constant g (9.81 m/s²) and the mean deepwater wave period Tp (also extracted from the CAWCR dataset).

The ordinary berm Ob can be approximated by formulas calculating the run-up of waves. We use the Stockdon et al. [19] formulation for the highest 2% of wave run-up (R₂), adding it to the value of MHHW to account for the run-up occurring during regular high tides:

$$ Ob=1.1\ast \left(0.35\ast s\ast {\left( Hs\ast {L}_0\right)}^{0.5}+\frac{{\left[{H}_s{L}_0\left(0.563\ast {s}^2+0.004\right)\right]}^{0.5}}{2}\right)+ MHHW $$

(7)

In this equation, the Mean Higher High Water MHHW was calculated by the OTPS model as described above, s is again the average beach slope (0.08), Hs is the mean deepwater significant wave height (extracted from the CAWCR dataset) and L₀ the mean deepwater wave length (also calculated with Eq. 6).

To calculate the storm wave swash height SWSH, i.e. the highest point wave splash is able to reach under storm conditions, we also used Eq. 7, but instead of using only the mean value for Hs and Tp, we used the mean value elevated by two standard deviations (2σ) extracted from the CAWCR wave hindcast.

For the deepest point of lagoons ld a uniform value of − 2 m was used. This value is derived by averaging 40 maximum depths of modern lagoons distributed around the world reported in Rovere et al. [15]. Further details on the location and depths can be found in Fig. 11 and supplementary material in Rovere et al. [15].

The IMCalc app

In order to make the data and the approach presented here easily accesible, we created a stand-alone software, which computes the indicative meaning and paleo RSL at any given coordinate, if that coordinate falls within the radius of 50 km distance to the closest data point of our Coastal Point Database (CPD).

The global datasets of wave and tidal conditions were, as described above, sampled in a regular distance, which resulted in the CPD and used as the main data input to our software. The fundamental necessary user input is limited to the coordinates and the RSL indicator. From this input (red lines in Fig. 3), the sofware uses the relations shown in the previous section to calculate the indicative meaning. If the measured elevation and elevation error are given, the software also adds to the outputs the elevation of paleo RSL and the associated uncertainty using Eqs. 1, 2, 3 and 4. All elevations entered and returned by the software, as well as all elevation values mentioned in this section, are related to mean sea level (MSL).

Fig. 3

Flow diagram showing the relations between the input and output values. Red lines represent the calculation using the Coastal Point Database (CPD, see text for details), blue dotted lines the calculation with an own modern analog

As an example, we will show here the calculation of the IM and paleo RSL for a beach deposit reported by Rovere et al. [15] at Cala Millor, Mallorca (3.385249 E; 39.590490 N). The elevation was measured with differential GPS at 1.22 ± 0.3 m. Now, these values on location, indicator and elevation can be entered into the IMCalc software. From the given coordinates the nearest values for the average significant wave height (Hs; 0.87 m) and wave period (Tp; 4.88 s) as well as the Mean Higher High Water (MHHW; 0.07 m) are extracted from the CPD. Since the RSL indicator is reported as a beach deposit, the limits can be calculated from Eq. 5 for the breaking depth (db; − 1.17 m) and from Eq. 7 for the ordinary berm (Ob; 0.54 m). Following Eqs. 1 and 2, the program will then calculate the IR (1.71 m) and the RWL (− 0.31 m), which together represent the indicative meaning of this location for a beach deposit. Finally, taking into account the measured elevation, Eqs. 3 and 4 are used to calculate the paleo RSL (1.53 m) as well as the associated uncertainty (0.91 m).

This simple example should only show the intended basal functionality of the software. The software also allows the user to provide data from local wave or tide models (if they are available) or for an own modern analog (if this was measured on-site). If this data is inserted, the respective CPD information is skipped and the indicative meaning and paleo RSL are calculated accordingly (blue, dotted lines in Fig. 3). Optionally, the software can show the wave and tide conditions, as well as a graph showing a sketch for the selected RSL indicator. In order to calculate larger datasets, a batch-processing option was also implemented. More details on the functionality and usage of the software are provided in an accompanying manual.

Results and discussion

General remarks

The CPD counts 414,296 points, where wave and tidal data conditions could be attributed (Fig. 1), making it possible to calculate the indicative meaning for any of the seven mentioned paleo RSL indicators. Besides the excluded Polar Regions (see above) and many smaller embayments, for which points near the open coastline may still be usable, data gaps are mainly restricted to the Black Sea and the Red Sea. In these areas, we recommend to input the local wave and/or tide conditions from available sources.

The IMCalc app can be used to quantify the indicative meaning of RSL indicators using a globally applicable, systematic approach that employs known relationships of wave and tide conditions with morpho- and hydrodynamic boundaries of landforms and deposits. We highlight that the results can only be used as first-order estimate in sea-level reconstructions, and should wherever possible be refined by local studies on the indicative meaning, as already highlighted by several authors.

We therefore stress that, compared to the approach used by our software, more detailed studies on the limits describing the indicators as well as on the equations for quantifying the limits can lead to more precise values for the indicative meaning. This is especially applicable when the limits can be detailed using local characteristics and through on-site measurements. Our approach only approximates the indicative meaning, but can be used as systematic and standardized method, if no data on specific modern analogs is available, as it is often the case for regional or global compilations of sea-level indicators. In the next section, we discuss the performance of our proposed approach towards a widely used dataset on Last Interglacial proxies, and we describe the caveats that characterize our work.

Comparison with a global database

The approach used by our IMCalc software can be assessed by comparing the calculated IR and RWL values to those obtained from independent literature. Hereafter, we compare the data extracted from IMCalc using the coordinates of 30 non-polar sites (for locations see Fig. 1) contained in the global MIS 5e RSL database of Kopp et al. [8]. These authors compiled information on the IM of MIS 5e RSL indicators at each site based on a qualitative and quantitative expert assessment of the original descriptions contained in the source papers. We only included those locations that could be correlated to one of the seven indicators described here. The difference between the RWL calculated with IMCalc to that reported in K09 is mostly in the range of ±2 m, the IR instead has a larger spread, and results in a difference in the range of ±4 m (Fig. 4). Despite this relatively high uncertainties, we want to highlight, that this approach can help in getting reliable estimations of the IM, although they might be high. Nevertheless, for the analysis of global databases of the paleo sea level reliable uncertainties are important, and too low uncertainties would lead to more problematic results than too high estimates (i.e. the statistical interpretation will be more affected by results with very narrow or no uncertainty, than by results with high uncertainty).

Fig. 4

Comparison between the indicative meanings as reported in the database of Kopp et al. [8] and as calculated by IMCalc. For all sites reported by K09 (location indicated in Fig. 1), the values have been calculated by IMCalc and the difference per site is shown here. Histograms (number of occurrence on the y-axis) show differences (x-axis, i.e. negative numbers show smaller values from IMCalc, positive numbers smaller values in K09) in the RWL (e-h) and IR (a-d) for coral reef terraces (a, e), beach deposits (b, f), marine terraces (c, g) and lagoonal deposits (d, h)

Comparing the values based on the type of indicator, the RWL of coral reef terraces has in average a quite small difference, although there is some diversion towards higher values calculated by IMCalc. The IR instead has also a large spread, but mainly gives smaller values than reported by K09. Regarding beach deposits, the RWL shows nearly a Gaussian distribution around 0, indicating a good correlation between the calculated and reported values. For these indicators the IR is quite irregular and shows mainly higher values when calculated by IMCalc. A similar pattern with fairly good constrained RWL and higher calculated IR can be seen when regarding marine terraces. Lagoonal deposits instead show a contrary pattern, where the IR is calculated similar to the reported values, but the RWL is always much smaller than the reported values by K09.

In general, this comparison shows that the approach presented here can calculate the values for the RWL and IR fairly good when compared to values reported in the literature. Nevertheless, it is important to mention that the number of locations and indicators compared above is very small. Especially the comparison of lagoonal deposits (n = 5) and marine terraces (n = 6) has to be taken with caution, whereas the beach deposits (n = 17) and coral reef terraces (n = 33) reach at least a notable number.

Limitations of the IMCalc app

The equations used to quantify the limits for the different RSL indicators are general applied equations for the quantification of the different hydrodynamic limits. A more detailed investigation on their indication of the specific limit as well as more comparisons with site-specific morphodynamic models might help in the further understanding of the limitations of this approach in calculating the indicative meaning and further refining of these equations.

In order to obtain a globally applicable approach, some aspects needed simplification. First, the beach slope used in IMCalc (Eqs. 5 and 7) was set to a globally uniform value (0.08), representing the average value of a global compilation of beach slopes reported by Liu et al. [10]. This value would need to be significantly higher (> 0.5) in order to have a major effect on the resulting limit elevations. Nevertheless, detailed values on the local beach slope will also result in more precise estimations of the local limits. For lagoonal deposits, we also use a globally uniform lower limit (− 2 m), which was calculated as average of worldwide lagoonal depth reported by Rovere et al. [15]. Therefore, studies using lagoonal deposits as sea-level indicator, need to account for this, or better should use a local modern analogous lagoonal depth as lower limit of the indicative range.

The comparison between the indicative meaning values of Kopp et al. [8] and those obtained with IMCalc show that our software can be used for a first estimate on the indicative meaning. In general, this comparison shows that we constrain well the RWL but overestimate the IR, which in turn leads to larger uncertainties in the final paleo sea-level estimate. Another aspect to consider is that IMCalc uses a general geomorphological description of the limits describing the IR, which was adopted from Rovere et al. [15]. In order to limit these effects, the paleo RSL information may be detailed, at many sites, through the combined application of different indicators (e.g. investigating two different indicators occurring at the same site) or through more precise information on the indicator itself (e.g. identifying sedimentological elements in a beach deposit that characterize it as formed in an intertidal environment). Therefore, a site-specific investigation of the modern analog should, if possible, always be undertaken and reported alongside the original study on the RSL indicators.

Conclusions

In this paper, we present a standardized and globally applicable method together with an accompanying software to calculate the indicative meaning of MIS 5e RSL indicators. This method is especially helpful when no data on a site-specific modern analog can be obtained. Also for the (re-)analysis of regional and global databases of paleo RSL indicators this method can be very helpful in order to quantify the IR on a local scale.

Nevertheless, there are two aspects that need to be regarded if this approach is used in sea-level reconstructions. First, the underlying datasets have a global scale and are relatively coarse. As the availability of detailed local data is often restricted or not available, this approach allows the user to get a first good estimate on the indicative meaning for every location around the world. Depending on the specific study site there might be more detailed data available, in this case more detailed data should be used as it would also improve the resulting values of the indicative meaning. Secondly, the formulas used in this study derive from coastal engineering studies usually not focusing on sea-level reconstructions. We showed above that these formulas can still be applied for this kind of study, but nevertheless, a more specialized study on these formulas in respect to sea-level indicators can also lead to a more detailed calculation of the respective limits and therefore a lower indicative range.

In the past, studies dealing with paleo sea-level reconstructions, often lacked to report the indicative meaning of the investigated sea-level indicators on a local scale. This often happened due to the lack of easy accessible site-specific data. In order, to overcome this problem, we developed the software presented here. Despite the proposed approach might still be improved by further research on this topic and needs to be validated with local studies, it should help in easily getting a first-order quantitative estimate on the IM for a specific study site.