diff --git a/.github/workflows/pkgdown.yaml b/.github/workflows/pkgdown.yaml
index 5aee4cb6..3b7e4ae2 100644
--- a/.github/workflows/pkgdown.yaml
+++ b/.github/workflows/pkgdown.yaml
@@ -33,10 +33,10 @@ jobs:
           extra-packages: any::pkgdown, local::.
           needs: website
           
-      - name: Get dev terra for dev pastclim
-        if: ${{ github.base_ref == 'dev' }} || ${{ github.ref_name == 'dev' }}
-        run: install.packages('terra', repos='https://rspatial.r-universe.dev')
-        shell: Rscript {0}
+#      - name: Get dev terra for dev pastclim
+#        if: ${{ github.base_ref == 'dev' }} || ${{ github.ref_name == 'dev' }}
+#        run: install.packages('terra', repos='https://rspatial.r-universe.dev')
+#        shell: Rscript {0}
 
       - name: Build site
         run: pkgdown::build_site_github_pages(new_process = FALSE, install = FALSE)
diff --git a/DESCRIPTION b/DESCRIPTION
index 53e7df5c..c65cbae7 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Package: pastclim
 Type: Package
 Title: Manipulate Time Series of Climate Reconstructions 
-Version: 2.2.0.9000
+Version: 2.2.0.9001
 Authors@R: c(
        person("Michela", "Leonardi", role = "aut"),
        person(c("Emily","Y."), "Hallet", role = "ctb"),
@@ -25,7 +25,7 @@ BugReports: https://github.com/EvolEcolGroup/pastclim/issues
 Encoding: UTF-8
 LazyData: true
 LazyDataCompression: xz
-RoxygenNote: 7.3.2
+RoxygenNote: 7.3.3
 Roxygen: list(markdown = TRUE)
 Depends: 
   R (>= 4.0.0),
diff --git a/R/bathy_to_spatraster.R b/R/bathy_to_spatraster.R
index 9c6910d0..7306323e 100644
--- a/R/bathy_to_spatraster.R
+++ b/R/bathy_to_spatraster.R
@@ -1,9 +1,9 @@
 #' Cast `bathy` to `SpatRaster`
 #'
-#' This function converts a [`marmap::bathy`][`marmap::as.bathy()`] object to
+#' This function converts a `marmap::bathy` object to
 #' a [`terra::SpatRaster`].
 #'
-#' @param bathy a [`marmap::bathy`][`marmap::as.bathy()`] to convert
+#' @param bathy a `marmap::bathy` to convert
 #' @returns a [`terra::SpatRaster`] with the relief for the chosen region
 #'
 #' @keywords internal
diff --git a/R/download_etopo_subset.R b/R/download_etopo_subset.R
index 8da2ff3e..486fb6bb 100644
--- a/R/download_etopo_subset.R
+++ b/R/download_etopo_subset.R
@@ -8,7 +8,7 @@
 #' fly from the NOAA server. If you plan to use the ETOPO2022 dataset
 #' extensively, it is worthwhile downloading it permanently to your computer
 #' with [download_etopo()], but beware that it is a large file (>1Gb). This
-#' function uses [marmap::getNOAA.bathy()] to download the data, and then
+#' function uses `marmap::getNOAA.bathy()` to download the data, and then
 #' converts them into a [`terra::SpatRaster`] formatted to be compatible with
 #' `pastclim`. NOTE: this function does not save the relief, it returns a
 #' [`terra::SpatRaster`]. If you plan to reuse this relief multiple times, it
@@ -17,7 +17,7 @@
 #' @param rast_template a [`terra::SpatRaster`] providing the extent and
 #'   resolution to be downloaded. This raster needs to have identical vertical
 #'   and horizontal resolution, and standard lat/long projection.
-#' @param ... additional parameters to be passed to [marmap::getNOAA.bathy()] to
+#' @param ... additional parameters to be passed to `marmap::getNOAA.bathy()` to
 #'   customise how files are stored. See the manpage for that function for
 #'   details
 #' @returns a [`terra::SpatRaster`] with the relief for the chosen region
diff --git a/R/download_paleoclim.R b/R/download_paleoclim.R
index 0e00fd77..48b69a8e 100644
--- a/R/download_paleoclim.R
+++ b/R/download_paleoclim.R
@@ -47,7 +47,7 @@ download_paleoclim <- function(dataset, bio_var, filename = NULL) {
   paleoclim_path[1] <- file.path(paleoclim_path[1], resolution)
   paleoclim_path[8] <- file.path(paleoclim_path[8], resolution)
   # create a vrt for each variable
-  for (i in seq_len(length(band_vector))) {
+  for (i in seq_along(band_vector)) {
     # build the vsizip paths
     paleoclim_vsizip <- paste0("/vsizip/", file.path(
       paleoclim_path,
diff --git a/R/download_worldclim_future.R b/R/download_worldclim_future.R
index 80f7840a..3c34a4c6 100644
--- a/R/download_worldclim_future.R
+++ b/R/download_worldclim_future.R
@@ -78,7 +78,7 @@ download_worldclim_future <- function(dataset, bio_var, filename = NULL) {
 
 
   # create a vrt for each variable
-  for (i in seq_len(length(band_vector))) {
+  for (i in seq_along(band_vector)) {
     vrt_path <- file.path(get_data_path(), paste0(
       dataset, "_",
       band_vector[i], "_v",
diff --git a/R/download_worldclim_present.R b/R/download_worldclim_present.R
index 1ac47139..aa680490 100644
--- a/R/download_worldclim_present.R
+++ b/R/download_worldclim_present.R
@@ -64,7 +64,7 @@ download_worldclim_present <- function(dataset, bio_var, filename) {
   }
 
   # create a vrt for each variable
-  for (i in seq_len(length(band_vector))) {
+  for (i in seq_along(band_vector)) {
     # build the vsizip paths
     if (!grepl("altitude", bio_var)) {
       worldclim_vsizip <- paste0("/vsizip/", file.path(
diff --git a/R/get_biome_classes.R b/R/get_biome_classes.R
index cabfbc81..696f10e6 100644
--- a/R/get_biome_classes.R
+++ b/R/get_biome_classes.R
@@ -34,7 +34,7 @@ get_biome_classes <- function(dataset) {
     biomes_string <- biomes_string[-length(biomes_string)]
     biomes_string <- substr(biomes_string, 4, nchar(biomes_string))
     biome_categories <- data.frame(
-      id = seq_len(length(biomes_string)),
+      id = seq_along(biomes_string),
       category = biomes_string
     )
   }
diff --git a/R/get_land_mask.R b/R/get_land_mask.R
index 25ada190..9591c99a 100644
--- a/R/get_land_mask.R
+++ b/R/get_land_mask.R
@@ -68,7 +68,6 @@ get_land_mask <- function(time_bp = NULL, time_ce = NULL, dataset) {
   ), terra::nlyr(land_mask))
 
 
-
   if (is.null(time_ce)) {
     names(land_mask) <- paste("land_mask", time_bp(land_mask), sep = "_")
   } else {
diff --git a/R/koeppen_geiger.R b/R/koeppen_geiger.R
index 7b74b132..4a176bf2 100644
--- a/R/koeppen_geiger.R
+++ b/R/koeppen_geiger.R
@@ -251,7 +251,6 @@ methods::setMethod(
 )
 
 
-
 #' @param filename filename to save the raster (optional).
 #' @rdname koeppen_geiger-methods
 #' @export
diff --git a/R/location_series.R b/R/location_series.R
index 5851cbe0..cf31b5ef 100644
--- a/R/location_series.R
+++ b/R/location_series.R
@@ -65,15 +65,21 @@ location_series <-
 
     check_dataset_path(dataset = dataset, path_to_nc = path_to_nc)
 
-    # if we are using standard datasets, check whether a variables exists
-    # and get the times
-    if (dataset != "custom") {
-      check_var_downloaded(bio_variables, dataset)
-      times <- get_time_bp_steps(dataset = dataset, path_to_nc = path_to_nc)
-    } else { # else check that the variables exist in the custom nc
-      check_var_in_nc(bio_variables, path_to_nc)
-      times <- get_time_bp_steps(dataset = "custom", path_to_nc = path_to_nc)
+    # get the region series for this dataset
+    climate_brick <- region_series(
+      bio_variables = bio_variables,
+      dataset = dataset,
+      path_to_nc = path_to_nc
+    )
+
+    # get all available times
+    times <- time_bp(climate_brick)
+
+    # if time_bp is NULL, get all times from the region series
+    if (is.null(time_bp)) {
+      time_bp <- time_bp(climate_brick)
     }
+
     time_bp_i <- time_bp_to_i_series(
       time_bp = time_bp,
       time_steps = times
@@ -105,21 +111,30 @@ location_series <-
 
     # now copy over the times to match the coordinates
     time_bp <- rep(time_bp, each = n_loc)
-    # and now feed the info to location_slice
-    location_ts <- location_slice(
-      x = x, time_bp = time_bp, coords = coords, bio_variables = bio_variables,
-      dataset = dataset, path_to_nc = path_to_nc,
-      nn_interpol = nn_interpol, buffer = buffer,
+
+    # now simply wrap around location_slice_from_region_series
+    location_ts <- location_slice_from_region_series(
+      x = x,
+      time_bp = time_bp,
+      time_ce = NULL,
+      coords = coords,
+      region_series = climate_brick,
+      nn_interpol = nn_interpol,
+      buffer = buffer,
       directions = directions
     )
 
-    # TODO if we had time_ce, we should convert back from time_bp
+    # if we had time_ce, we should convert back from time_bp
+    if (!is.null(time_ce)) {
+      location_ts$time_ce <- location_ts$time_bp + 1950
+      # remove the time_bp column
+      location_ts <- location_ts[, !names(location_ts) %in% "time_bp"]
+    }
 
     return(location_ts[, !names(location_ts) %in% "time_bp_slice"])
   }
 
 
-
 #' Extract a time series of bioclimatic variables for one or more locations.
 #'
 #' Deprecated version of [location_series()]
diff --git a/R/location_slice_from_region_series.R b/R/location_slice_from_region_series_fast.R
similarity index 66%
rename from R/location_slice_from_region_series.R
rename to R/location_slice_from_region_series_fast.R
index 381dced1..7a85cd7a 100644
--- a/R/location_slice_from_region_series.R
+++ b/R/location_slice_from_region_series_fast.R
@@ -55,10 +55,10 @@ location_slice_from_region_series <- # nolint
     bio_variables <- names(region_series)
 
     time_bp <- check_time_vars(time_bp = time_bp, time_ce = time_ce)
-    # boolean whether we will need to readd time_ce instead of time_bp
-    readd_ce <- FALSE
+    # boolean whether we will need to read time_ce instead of time_bp
+    read_ce <- FALSE
     if (any(!is.null(time_ce), "time_ce" %in% names(x))) {
-      readd_ce <- TRUE
+      read_ce <- TRUE
     }
 
     # if we have a data.frame
@@ -123,78 +123,97 @@ location_slice_from_region_series <- # nolint
       time_bp = locations_data$time_bp, time_steps = times
     )
     locations_data$time_bp_slice <- times[time_indeces]
-    unique_times <- unique(locations_data$time_bp_slice)
 
-    for (i_time in unique_times) {
-      this_slice <- slice_region_series(climate_brick,
-        time_bp = i_time
+    # if we are not using a buffer, we try to get all the values from
+    # specific locations directly
+    if (!buffer) {
+      locations_climate <- terra::extract(climate_brick,
+        y = locations_data[coords],
+        layer = time_indeces
       )
+      # for each bio_variable (the list element), extract the value column from
+      # the data.frame and assign to locations_data
+      for (var in bio_variables) {
+        locations_data[[var]] <- locations_climate[[var]]$value
+      }
+    } else {
+      # if we are using a buffer, set to NA as we will compute them later
+      for (var in bio_variables) {
+        locations_data[[var]] <- NA
+      }
+    }
 
-      this_slice_indeces <- which(locations_data$time_bp_slice == i_time)
-      if (!buffer) { # get the specific values for those locations
-        this_climate <- terra::extract(
-          x = this_slice,
-          y = locations_data[locations_data$time_bp_slice == i_time, coords]
+    # if we interpolate or use a buffer, we have to find the incomplete cases
+    if (nn_interpol || buffer) {
+      loc_id_to_move <- which(!stats::complete.cases(locations_data))
+      # only do something if we have some locations for which we have no data
+      if (length(loc_id_to_move) != 0) {
+        # for each id, get the appropriate list of neighbours
+        rast_id_to_move <- terra::cellFromXY(
+          climate_brick[[1]],
+          as.matrix(coords_df[loc_id_to_move, ])
         )
-        # factors don't behave nicely when adding new elements, cast to
-        # character
-        if ("biome" %in% names(this_climate)) {
-          this_climate$biome <- as.character(this_climate$biome)
-        }
-        # sort out the indexing here
-        locations_data[locations_data$time_bp_slice == i_time, bio_variables] <-
-          this_climate[
-            ,
-            bio_variables
-          ]
-      } else { # set to NA as we will compute them with a buffer
-        locations_data[this_slice_indeces, ] <- NA
-      }
+        neighbours_list <- terra::adjacent(
+          climate_brick[[1]],
+          rast_id_to_move,
+          directions = directions,
+          pairs = FALSE
+        )
+        # convert from matrix (one row per location) to data.frame where first
+        # column is focal location, second column is id of each neighbour
+        neighbours_list <- data.frame(
+          focal_id = rep(as.numeric(rownames(neighbours_list)),
+            each = ncol(neighbours_list)
+          ),
+          neighbour_id = as.vector(t(neighbours_list)),
+          layer = rep(time_indeces[loc_id_to_move],
+            each = ncol(neighbours_list)
+          ),
+          unique_id = rep(seq_len(nrow(neighbours_list)),
+            each = ncol(neighbours_list)
+          )
+        )
+        # there is a BUG terra does not seem to cope with using layer when
+        # extracting from sds if y is cellID, y has to be a data.frame
+        neighbours_coords <- as.data.frame(terra::xyFromCell(
+          climate_brick[[1]],
+          neighbours_list$neighbour_id
+        ))
+        # extract the climate for all neighbours
+        neighbours_values <-
+          terra::extract(
+            x = climate_brick,
+            y = neighbours_coords,
+            layer = neighbours_list$layer
+          )
 
-      if (nn_interpol || buffer) {
-        locations_to_move <- this_slice_indeces[
-          this_slice_indeces %in%
-            which(!stats::complete.cases(locations_data))
-        ]
-        if (length(locations_to_move) == 0) {
-          next
-        }
-        for (i in locations_to_move) {
-          if (inherits(x, "data.frame")) {
-            cell_id <-
-              terra::cellFromXY(this_slice, as.matrix(coords_df[
-                i,
-              ]))
+        # for each variable, compute the mean across neighbours
+        # and replace the values in locations_data
+        for (i_var in bio_variables) {
+          if (i_var == "biome") {
+            # for factors, compute the mode across neighbours
+            neighbours_mode <- tapply(
+              neighbours_values[[i_var]]$value,
+              neighbours_list$unique_id,
+              mode
+            )
+            locations_data[loc_id_to_move, i_var] <-
+              neighbours_mode
           } else {
-            cell_id <- coords_df[i]
-          }
-          neighbours_ids <-
-            terra::adjacent(this_slice, cell_id,
-              directions = directions, pairs = FALSE
+            neighbours_mean <- tapply(
+              neighbours_values[[i_var]]$value,
+              neighbours_list$unique_id,
+              mean,
+              na.rm = TRUE
             )
-
-          neighbours_values <-
-            terra::extract(
-              x = this_slice,
-              y = neighbours_ids[1, ]
-            ) # [, bio_variables]
-
-          neighbours_values_mean <- colMeans(
-            neighbours_values[, !names(neighbours_values) %in% "biome",
-              drop = FALSE
-            ],
-            na.rm = TRUE
-          )
-          if ("biome" %in% bio_variables) {
-            neighbours_values_mean["biome"] <-
-              mode(as.character(neighbours_values[, "biome"]))
+            locations_data[loc_id_to_move, i_var] <-
+              neighbours_mean
           }
-          locations_data[i, bio_variables] <-
-            neighbours_values_mean[bio_variables]
         }
       }
     }
-    # is.nan has not method for a data.frame
+
+    # is.nan has no method for a data.frame
     # nolint start
     is.nan.data.frame <- function(x) {
       do.call(cbind, lapply(x, is.nan))
@@ -205,7 +224,7 @@ location_slice_from_region_series <- # nolint
 
     locations_data <- locations_data[order(orig_id), ]
 
-    if (readd_ce) {
+    if (read_ce) {
       locations_data$time_ce <- locations_data$time_bp + 1950
       locations_data$time_ce_slice <- locations_data$time_bp_slice + 1950
       locations_data <- locations_data[
diff --git a/R/region_series.R b/R/region_series.R
index 200139b0..2a92d276 100644
--- a/R/region_series.R
+++ b/R/region_series.R
@@ -104,22 +104,26 @@ region_series <-
         this_var_longname <- NULL
         this_var_units <- NULL
       }
+
+      # retrieve time axis for virtual file
+      var_brick <- pastclim_rast(
+        x = this_file, bio_var_orig = this_var_orig,
+        bio_var_pastclim = this_var, var_longname = this_var_longname,
+        var_units = this_var_units
+      )
+
       # figure out the time indeces the first time we run this
       if (is.null(time_index)) {
-        # as we have the file name, we can us the same code for custom and
-        # standard datasets.
-        times <- get_time_bp_steps(dataset = "custom", path_to_nc = this_file)
+        # get times from the var_brick
+        times <- time_bp(var_brick)
+        # convert it to indeces
         time_index <- time_bp_to_i_series(
           time_bp = time_bp,
           time_steps = times
         )
       }
-      # retrieve time axis for virtual file
-      var_brick <- pastclim_rast(
-        x = this_file, bio_var_orig = this_var_orig,
-        bio_var_pastclim = this_var, var_longname = this_var_longname,
-        var_units = this_var_units
-      )
+
+
       # subset to time steps
       if (!is.null(time_bp)) {
         var_brick <- terra::subset(var_brick, subset = time_index)
diff --git a/R/sample_region_series.R b/R/sample_region_series.R
index 32022045..f01d6f21 100644
--- a/R/sample_region_series.R
+++ b/R/sample_region_series.R
@@ -130,7 +130,7 @@ sample_rs_variable <- function(x, size, method = "random", replace = FALSE,
   # create list to store samples for each time step
   sample_list <- list()
   t_steps <- time_bp(x[1])
-  for (i in seq_len(length(size))) {
+  for (i in seq_along(size)) {
     if (size[i] > 0) {
       x_step <- slice_region_series(x, t_steps[i])
       # spatSample samples additional points to make sure it has enough points
diff --git a/R/time_bp.R b/R/time_bp.R
index 41ec8ae0..df7530a8 100644
--- a/R/time_bp.R
+++ b/R/time_bp.R
@@ -65,7 +65,6 @@ setMethod(
 )
 
 
-
 #' @rdname time_bp
 #' @export
 setGeneric("time_bp<-", function(x, value) {
diff --git a/R/vrt_set_meta.R b/R/vrt_set_meta.R
index 1785ebc0..d923fc2d 100644
--- a/R/vrt_set_meta.R
+++ b/R/vrt_set_meta.R
@@ -37,7 +37,7 @@ vrt_set_meta <- function(vrt_path, description, time_vector, time_bp = TRUE) {
     )
     return(FALSE)
   }
-  for (i_node in seq_len(length(band_nodes))) {
+  for (i_node in seq_along(band_nodes)) {
     # add a unique description label for this band (variable_time combination)
     xml2::xml_add_child(band_nodes[i_node], "Description",
       paste(description, time_vector[i_node], sep = "_"),
diff --git a/data-raw/to_update/repackage_beyer2020 rewrite.R b/data-raw/to_update/repackage_beyer2020 rewrite.R
index f4bd1d3a..caaf5805 100644
--- a/data-raw/to_update/repackage_beyer2020 rewrite.R	
+++ b/data-raw/to_update/repackage_beyer2020 rewrite.R	
@@ -220,7 +220,6 @@ ncdf4::nc_close(nc_in)
 unlink("Beyer2020_monthly0*")
 
 
-
 ################################
 # Now remove ice sheets and internal seas
 
diff --git a/data-raw/to_update/repackage_krapp2021.R b/data-raw/to_update/repackage_krapp2021.R
index 2a6d632a..b8c37a15 100644
--- a/data-raw/to_update/repackage_krapp2021.R
+++ b/data-raw/to_update/repackage_krapp2021.R
@@ -158,8 +158,6 @@ monthly_dirs <- file.path(wkdir, "original_files", monthly_dirs)
 # ncatted -a pastclim_version,global,m,c,'1.4.0' Krapp2021_prec_monthly_v1.4.0.nc
 
 
-
-
 # NOTE this needs to be changed to work for uncut as well.
 for (file_id in list.files("./repackaged/cut", full.names = TRUE)) {
   ncatted(paste0("-a author,global,d,, -a history,global,d,, -a description,global,d,, -a command_line,global,d,, -h ", file_id))
diff --git a/man/bathy_to_spatraster.Rd b/man/bathy_to_spatraster.Rd
index d8cbb2fa..0d8a5f25 100644
--- a/man/bathy_to_spatraster.Rd
+++ b/man/bathy_to_spatraster.Rd
@@ -7,13 +7,13 @@
 bathy_to_spatraster(bathy)
 }
 \arguments{
-\item{bathy}{a \code{\link[marmap:as.bathy]{marmap::bathy}} to convert}
+\item{bathy}{a \code{marmap::bathy} to convert}
 }
 \value{
 a \code{\link[terra:SpatRaster-class]{terra::SpatRaster}} with the relief for the chosen region
 }
 \description{
-This function converts a \code{\link[marmap:as.bathy]{marmap::bathy}} object to
+This function converts a \code{marmap::bathy} object to
 a \code{\link[terra:SpatRaster-class]{terra::SpatRaster}}.
 }
 \keyword{internal}
diff --git a/man/download_etopo_subset.Rd b/man/download_etopo_subset.Rd
index ab5c2898..05eeb123 100644
--- a/man/download_etopo_subset.Rd
+++ b/man/download_etopo_subset.Rd
@@ -11,7 +11,7 @@ download_etopo_subset(rast_template, ...)
 resolution to be downloaded. This raster needs to have identical vertical
 and horizontal resolution, and standard lat/long projection.}
 
-\item{...}{additional parameters to be passed to \code{\link[marmap:getNOAA.bathy]{marmap::getNOAA.bathy()}} to
+\item{...}{additional parameters to be passed to \code{marmap::getNOAA.bathy()} to
 customise how files are stored. See the manpage for that function for
 details}
 }
@@ -28,7 +28,7 @@ relatively low resolution. This function fetches the necessary subset on the
 fly from the NOAA server. If you plan to use the ETOPO2022 dataset
 extensively, it is worthwhile downloading it permanently to your computer
 with \code{\link[=download_etopo]{download_etopo()}}, but beware that it is a large file (>1Gb). This
-function uses \code{\link[marmap:getNOAA.bathy]{marmap::getNOAA.bathy()}} to download the data, and then
+function uses \code{marmap::getNOAA.bathy()} to download the data, and then
 converts them into a \code{\link[terra:SpatRaster-class]{terra::SpatRaster}} formatted to be compatible with
 \code{pastclim}. NOTE: this function does not save the relief, it returns a
 \code{\link[terra:SpatRaster-class]{terra::SpatRaster}}. If you plan to reuse this relief multiple times, it
diff --git a/man/location_slice_from_region_series.Rd b/man/location_slice_from_region_series.Rd
index 2c32d03f..2d7988b4 100644
--- a/man/location_slice_from_region_series.Rd
+++ b/man/location_slice_from_region_series.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/location_slice_from_region_series.R
+% Please edit documentation in R/location_slice_from_region_series_fast.R
 \name{location_slice_from_region_series}
 \alias{location_slice_from_region_series}
 \title{Extract local climate for one or more locations for a given time slice.}
diff --git a/man/mode.Rd b/man/mode.Rd
index d6dad4f1..9fc55b83 100644
--- a/man/mode.Rd
+++ b/man/mode.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/location_slice_from_region_series.R
+% Please edit documentation in R/location_slice_from_region_series_fast.R
 \name{mode}
 \alias{mode}
 \title{Mode}
diff --git a/man/pastclim-package.Rd b/man/pastclim-package.Rd
index 8cea9219..26bb98ee 100644
--- a/man/pastclim-package.Rd
+++ b/man/pastclim-package.Rd
@@ -20,7 +20,7 @@ Useful links:
 
 }
 \author{
-\strong{Maintainer}: Andrea Manica \email{am315@cam.ac.uk}
+\strong{Maintainer}: Andrea Manica \email{am315@cam.ac.uk} [copyright holder]
 
 Authors:
 \itemize{
diff --git a/tests/testthat/test_filenames_chelsa_future.R b/tests/testthat/test_filenames_chelsa_future.R
index c9121239..c361cc61 100644
--- a/tests/testthat/test_filenames_chelsa_future.R
+++ b/tests/testthat/test_filenames_chelsa_future.R
@@ -16,7 +16,7 @@ test_that("filenames_chelsa_future are correct", {
     "0.5m",
     sep = "_"
   )
-  for (i in seq_len(length(target_datasets))) {
+  for (i in seq_along(target_datasets)) {
     bio_files <- filenames_chelsa_future(
       dataset = target_datasets[i],
       bio_var = sample(variables, 1)
@@ -29,7 +29,7 @@ test_that("filenames_chelsa_future are correct", {
     paste0("precipitation_", sprintf("%02d", 1:12))
   )
 
-  for (i in seq_len(length(target_datasets))) {
+  for (i in seq_along(target_datasets)) {
     bio_files <- filenames_chelsa_future(
       dataset = target_datasets[i],
       bio_var = sample(variables, 1)
diff --git a/tests/testthat/test_filenames_worldclim_future.R b/tests/testthat/test_filenames_worldclim_future.R
index 4853454c..af056ac1 100644
--- a/tests/testthat/test_filenames_worldclim_future.R
+++ b/tests/testthat/test_filenames_worldclim_future.R
@@ -28,7 +28,7 @@ test_that("filenames_worldclim_future are correct", {
   ]
   # nolint end
 
-  for (i in seq_len(length(target_datasets))) {
+  for (i in seq_along(target_datasets)) {
     bio_files <- filenames_worldclim_future(
       dataset = target_datasets[i],
       bio_var = sample(variables, 1)
@@ -42,7 +42,7 @@ test_that("filenames_worldclim_future are correct", {
     paste0("precipitation_", sprintf("%02d", 1:12))
   )
 
-  for (i in seq_len(length(target_datasets))) {
+  for (i in seq_along(target_datasets)) {
     bio_files <- filenames_worldclim_future(
       dataset = target_datasets[i],
       bio_var = sample(variables, 1)
diff --git a/tests/testthat/test_location_series.R b/tests/testthat/test_location_series.R
index 6f0b4b10..7d02f39a 100644
--- a/tests/testthat/test_location_series.R
+++ b/tests/testthat/test_location_series.R
@@ -47,6 +47,25 @@ test_that("location_series", {
   )
   expect_true(nrow(locations_ts) == 12)
 
+  # do the same for time-ce
+  locations_ts_ce <- location_series(
+    x = locations[, c("longitude", "latitude")],
+    time_ce = c(-20000, -10000, -5000) + 1950,
+    bio_variables = c("bio01", "bio12"),
+    dataset = "Example"
+  )
+  expect_true(nrow(locations_ts_ce) == 12)
+  # test that time_ce is a column in the output
+  expect_true("time_ce" %in% names(locations_ts_ce))
+  # test that time_bp is not a column in the output
+  expect_false("time_bp" %in% names(locations_ts_ce))
+  # test that the values are the same as when using time_bp
+  expect_equal(
+    locations_ts$bio01,
+    locations_ts_ce$bio01
+  )
+
+
   locations_ts <- location_series(
     x = locations[, c("longitude", "latitude")],
     time_bp = list(min = -10000, max = -5000),