2.1 Shoreline displacement curve for Ørland

Geologically derived displacement curves are in Norway often reported with a uniform probability between the lower and upper limit of the possible elevation of the sea-level over time. This is also the case for the displacement curve published for Ørland. At present, the required format for displacement curves to be used with shoredate is therefore a data frame with the columns bce, upperelev, lowerelev and name. Here, upperelev and lowerelev denote the highest and lowest possible elevation for the sea-level for each year BCE. Years BCE, defined by the column bce, should be given at an interval with equal or higher resolution than that which is desired for the resulting shoreline date. The probability between the upper and lower elevation limits are here assumed to be uniformly distributed.

The displacement curve for Ørland is provided with shoredate and can be loaded with the following code:

# The packages ggplot2 and sf are explictly loaded here, instead of using the
# package::function() syntax with each function call below.
library(ggplot2)
library(sf)

# Load shoredate
library(shoredate)

# Load the displacement curve for Ørland
orland_disp <- get(load(system.file("extdata/orland_displacement_curve.rda",
                  package = "shoredate")))

## Warning in readChar(con, 5L, useBytes = TRUE): cannot open compressed file '', probable reason 'No such file or directory'

## Error in readChar(con, 5L, useBytes = TRUE): cannot open the connection

# Print the last few rows of the Ørland displacement curve
tail(orland_disp)

## Error in eval(expr, envir, enclos): object 'orland_disp' not found

This can then be plotted with displacement_plot(). The name is here adjusted with the argument target_name to something more informative, and the geologically derived displacement curves from south-eastern Norway are excluded by specifying the displacement_curves argument, as these are not of relevance here.

displacement_plot(target_curve = orland_disp,
                  target_name = "Displacement curve for Ørland, Central Norway",
                  displacement_curves = NA)

## Error in displacement_plot(target_curve = orland_disp, target_name = "Displacement curve for Ørland, Central Norway", : unused arguments (target_curve = orland_disp, target_name = "Displacement curve for Ørland, Central Norway")

A note should be made that the Ørland curve was originally published with reference to the highest astronomical tide in the region (see Romundset and Lakeman 2019:66). To adjust this to mean sea-level, the difference in elevation between the highest astronomical tide and the mean sea-level has been been subtracted from the elevation values.

Furthermore, while variable uplift rates are also relevant for this area, this is not corrected for in this example. The example therefore focuses on the small area for which the curve was developed, where it can be assumed to be directly applicable. Functionality to adjust for variable displacement rates with shoredate is only accommodated for the Skagerrak area through the function interpolate_curve(). Thus, if similar adjustments are to made in other areas, this will at present therefore have to be done independently of shoredate, to which the adjusted curves can then be passed.

2.2 Creating maps of Ørland

Having loaded the displacement curve, this can then be directly used to shoreline date sites in locations where it applies, provided the elevation of the sites above present sea-level is known. First we will create a couple of fictitious site examples, each represented by a point, and maps of their location to demonstrate the extendibility of the target_plot() function.

# Create example sites
target_points <-  st_sfc(st_point(c(532719, 7065723)),
                         st_point(c(532896, 7066260)))

# Set CRS
target_points <- st_as_sf(target_points, crs = 32632)

# Add site names
target_points$name <- c("Example 1", "Example 2")

To create a map of where Ørland and these target points are located, one can adapt the target_plot() function. First, set up the geometries to be plotted:

# Load in the limit of the spatial coverage in south-eastern Norway,
# which is provided with the package
senorway <- st_read(system.file("extdata/spatial_limit.gpkg",
                                    package = "shoredate"), quiet = TRUE)

# Assign a name to this for the map legend
senorway$name <- "Skagerrak limit"

# Retrieve the first of the example points, to represent Ørland
orland <- target_points[1,]
orland$location = "Ørland"

Once this has been done we can use target_plot() to set up a plot. Setting naturalearth_basemap to TRUE downloads a world map from https://www.naturalearthdata.com/ using the rnaturalearth package. This is stored in a temporary folder and is deleted when the current R session is ended. The argument naturalearth_zoom specifies the amount of cropping that is done on this world map, with the provided targets as the focal point. The argument crs_epsg is here the same as the default, but is explicitly called to highlight that different coordinate reference systems can be used. Setting the argument isobases to NA means that the default isobases pertaining to south-eastern Norway are excluded from the plot. Finally, setting target_labels to FALSE excludes the labelling of the target points in the plot, which will instead be handled with a legend in the code to follow below.

overview_map <- target_plot(targets = orland,
                      naturalearth_basemap = TRUE,
                      naturalearth_zoom = c(1000000, 1000000),
                      crs_epsg = 32632,
                      base_col = "black",
                      base_fill = "grey85",
                      isobases = NA,
                      target_labels =  FALSE)

## Error in target_plot(targets = orland, naturalearth_basemap = TRUE, naturalearth_zoom = c(1e+06, : unused arguments (naturalearth_basemap = TRUE, naturalearth_zoom = c(1e+06, 1e+06), crs_epsg = 32632, base_col = "black", base_fill = "grey85", target_labels = FALSE)

Having created a base plot, this can now be manipulated using other functions from the package ggplot2.

overview_map <- overview_map +
        geom_sf(data = senorway, aes(col = name), fill = NA,
                   lwd = 0.5, show.legend = "polygon") +
        # Replotting the point for Ørland to add it to the legend
        geom_sf(data = orland, aes(fill = location),
                   col = "black", shape = 21,
                   size = 3, show.legend = "point") +
        scale_fill_manual(values = c("Ørland" = "red")) +
        scale_colour_manual(values = c("Skagerrak limit" = "red"),
                               guide = guide_legend(
                                          override.aes = list(shape = NA))) +
        theme(legend.position = "bottom",
              legend.title = element_blank()) +
        ggtitle(paste("Location of Ørland relative to the",
                       "spatial limit in south-eastern Norway"))

## Error in eval(expr, envir, enclos): object 'overview_map' not found

When the map displaying the location of Ørland is created, we can create a second map that shows the location of the example points, zoomed in at Ørland, which is located at the tip of the Fosen peninsula. While the map is fairly simple, it can nonetheless be useful to perform this exercise to make sure everything looks as it should.

# Create basemap
examples_map <- target_plot(targets = target_points,
                      naturalearth_basemap = TRUE,
                      naturalearth_zoom = c(15000, 10000),
                      base_col = "black",
                      base_fill = "grey85",
                      isobases = NA)

## Error in target_plot(targets = target_points, naturalearth_basemap = TRUE, : unused arguments (naturalearth_basemap = TRUE, naturalearth_zoom = c(15000, 10000), base_col = "black", base_fill = "grey85")

# Add axis labels and ticks, which are not returned with target_plot()
examples_map <- examples_map +
      theme(
        axis.text.y = element_text(),
        axis.text.x = element_text(),
        axis.ticks = element_line()) +
      coord_sf(datum = st_crs(target_points), expand = FALSE) +
      ggtitle("Location of example points in Ørland")

## Error in eval(expr, envir, enclos): object 'examples_map' not found

# Call overview map
overview_map

## Error in eval(expr, envir, enclos): object 'overview_map' not found

# Call map displaying the location of the example points
examples_map

## Error in eval(expr, envir, enclos): object 'examples_map' not found

2.3 Dating example sites in Ørland

For this example, a raster retrieved from Amazon Web Service Terrain Tiles is used for finding the elevation of the example sites. This follows the procedure that is outlined in the main vignette.

# Retrieve raster
target_wgs84 <- st_transform(target_points, crs = 4326)
elev_raster <- elevatr::get_elev_raster(target_wgs84, z = 14,
                                        verbose = FALSE, src = "aws")
elev_raster <- terra::project(terra::rast(elev_raster), "epsg:32632")

# Plot the raster and sites for inspection
terra::plot(elev_raster)
plot(target_points, col = "black", add = TRUE)
text(st_coordinates(target_points) - 50, labels = target_points$name)

plot of chunk examples-plot

After the raster has been loaded we can find the elevation of the sites:

terra::extract(elev_raster, terra::vect(target_points))[, -1]

## [1] 17.74340 19.53946

This elevation is retrieved by shoreline_date() if the raster is passed to its elevation argument and the sites to be dated are provided as spatial geometries. However, to illustrate the point that all that is needed to use shoreline_date() is a displacement curve and knowledge of these elevations, the function is here called by only providing the displacement curve and providing a character vector with the name of the sites and a numerical vector with their elevations:

target_dates <- shoreline_date(sites = c("Example 1", "Example 2"),
                               target_curve = orland_disp,
                               elevation = c(17.7, 19.5))

## Error in shoreline_date(sites = c("Example 1", "Example 2"), target_curve = orland_disp, : unused argument (target_curve = orland_disp)

# Plot the results
shoredate_plot(target_dates, multiplot = TRUE)

plot of chunk orland-dateplot

2.4 Finding the earliest possible date

As mentioned, the above implementation assumes that the relationship between the sites and the shoreline is characterised by the same gamma function as that identified for sites in south-eastern Norway. At present, the only adjustment that is possible to do to this assumed relationship is changing the parameters for the gamma function when calling shoreline_date(), using the argument model_parameters.

However, as an alternative it is also possible to perform the dating procedure without accounting for any distance between the site and the contemporaneous shoreline. This is done by setting the parameter model to “none” when calling shoreline_date(). This thus effectively provides a termnius post quem date, under the assumption that the earliest possible date for when the site was in use is when the location of the site emerged from the sea.

While a terminus post quem date limits the further inferential steps that can be taken, it might be more appropriate to apply shoreline dating in this manner in regions where the relationship between sites and the shoreline is unknown. This can potentially also be extended to the dating of other phenomena such as rock art, or other cases where this relationship is less certain (see e.g. Sognnes 2003).

Furthermore, it could also be possible to reverse this logic in regions that have instead been subject to relative sea-level rise, where the date for when a location was inundated can provide a terminus ante quem date – the latest possible date for the use of a site. However, implementation of shoredate to regions which have experienced continuous or disjoint phases of relative sea-level rise remains to be developed.

# Finding the earliest possible date for the first example point
earliest_date <- shoreline_date(target_points[1,],
                                target_curve = orland_disp,
                                elevation = elev_raster,
                                model = "none",
                                hdr_prob = 1)

## Error in shoreline_date(target_points[1, ], target_curve = orland_disp, : unused argument (target_curve = orland_disp)

# Call to plot
shoredate_plot(earliest_date, site_name = TRUE)

## Error in eval(expr, envir, enclos): object 'earliest_date' not found