diff --git a/NAMESPACE b/NAMESPACE
index 04ac3a520..8c3885e01 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -18,7 +18,6 @@ export(LogNormal)
 export(NonParametric)
 export(Normal)
 export(R_to_growth)
-export(add_day_of_week)
 export(adjust_infection_to_report)
 export(apply_tolerance)
 export(backcalc_opts)
@@ -29,20 +28,9 @@ export(calc_summary_measures)
 export(calc_summary_stats)
 export(clean_nowcasts)
 export(collapse)
-export(construct_output)
 export(convert_to_logmean)
 export(convert_to_logsd)
 export(convolve_and_scale)
-export(copy_results_to_latest)
-export(create_backcalc_data)
-export(create_clean_reported_cases)
-export(create_gp_data)
-export(create_initial_conditions)
-export(create_obs_model)
-export(create_rt_data)
-export(create_shifted_cases)
-export(create_stan_args)
-export(create_stan_data)
 export(delay_opts)
 export(discretise)
 export(discretize)
@@ -54,7 +42,6 @@ export(estimate_delay)
 export(estimate_infections)
 export(estimate_secondary)
 export(estimate_truncation)
-export(estimates_by_report_date)
 export(expose_stan_fns)
 export(extract_CrIs)
 export(extract_inits)
@@ -73,8 +60,6 @@ export(get_parameters)
 export(get_pmf)
 export(get_raw_result)
 export(get_regional_results)
-export(get_regions)
-export(get_regions_with_most_reports)
 export(gp_opts)
 export(growth_to_R)
 export(lognorm_dist_def)
@@ -96,15 +81,11 @@ export(rstan_sampling_opts)
 export(rstan_vb_opts)
 export(rt_opts)
 export(sample_approx_dist)
-export(save_estimate_infections)
-export(save_input)
 export(secondary_opts)
 export(set_dt_single_thread)
 export(setup_default_logging)
-export(setup_dt)
 export(setup_future)
 export(setup_logging)
-export(setup_target_folder)
 export(simulate_infections)
 export(simulate_secondary)
 export(stan_laplace_opts)
@@ -112,11 +93,7 @@ export(stan_opts)
 export(stan_pathfinder_opts)
 export(stan_sampling_opts)
 export(stan_vb_opts)
-export(summarise_key_measures)
-export(summarise_results)
 export(trunc_opts)
-export(update_horizon)
-export(update_list)
 export(update_secondary_args)
 import(Rcpp)
 import(methods)
diff --git a/NEWS.md b/NEWS.md
index df8c3fa15..89a472e5c 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -3,12 +3,20 @@
 ## Major changes
 
 * Delay discretisation is now based on a two-day censoring window (with uniform probability in between), based on recommendations in [Park et al, medRxiv, 2024](https://doi.org/10.1101/2024.01.12.24301247). By @sbfnk in #518 and reviewed by @jamesmbaazam.
+* The interface to generating delay distributions has been completely overhauled. Instead of calling `dist_spec()` users now specify distributions using functions that represent the available distributions, i.e. `LogNormal()`, `Gamma()` and `Fixed()`. See `?EpiNow2::Distributions`. Uncertainty is specified using calls of the same nature, to `Normal()`. More information on the underlying design can be found in `inst/dev/design_dist.md` By @sbfnk in #504 and reviewed by @seabbs.
 
-## OTher breaking changes
+## Deprecations
 
-* The functions `get_dist`, `get_generation_time`, `get_incubation_period` have been deprecated and replaced with examples. By @sbfnk in #481 and reviewed by @seabbs.
+* The functions `sample_approx_dist()`, `report_cases()`, and `adjust_infection_reports()` have been deprecated as the functionality they provide can now be achieved with `simulate_secondary()`. See #597 by @jamesmbaazam and reviewed by @sbfnk.
 * The utility function `update_list()` has been deprecated in favour of `utils::modifyList()` because it comes with an installation of R. By @jamesmbaazam in #491 and reviewed by @seabbs.
 * The `fixed` argument to `dist_spec` has been deprecated and replaced by a `fix_dist()` function. By @sbfnk in #503 and reviewed by @seabbs.
+* The functions `get_dist()`, `get_generation_time()`, `get_incubation_period()` have been deprecated and replaced with examples. By @sbfnk in #481 and reviewed by @seabbs.
+* The function `init_cumulative_fit()` has been deprecated. By @jamesmbaazam in #541 and reviewed by @sbfnk.
+* The model-specific `weigh_delay_priors` argument has been deprecated in favour of delay-specific prior weighting using `weight_priors`. See `generation_time_opts()`, `delay_opts()`, and `trunc_opts()`. By @sbfnk in #630 and reviewed by @jamesmbaazam.
+* All functions now use a `data` argument to pass data. The existing `reported_cases`, `reports`, and `obs` arguments are deprecated and will be removed in v2.0.0. By @jamesmbaazam in #638 and reviewed by @sbfnk.
+
+## Other breaking changes
+
 * Updated `estimate_infections()` so that rather than imputing missing data, it now skips these data points in the likelihood. This is a breaking change as it alters the behaviour of the model when dates are missing from a time series but are known to be zero. We recommend that users check their results when updating to this version but expect this to in most cases improve performance. By @seabbs in #528 and reviewed by @sbfnk.
 * `simulate_infections` has been renamed to `forecast_infections` in line with `simulate_secondary` and `forecast_secondary`. The terminology is: a forecast is done from a fit to existing data, a simulation from first principles. By @sbfnk in #544 and reviewed by @seabbs.
 * A new `simulate_infections` function has been added that can be used to simulate from the model from given initial conditions and parameters. By @sbfnk in #557 and reviewed by @jamesmbaazam.
@@ -23,9 +31,10 @@
 * Removed references to the no longer existing `forecast_infections` function. By @sbfnk in #460 and reviewed by @seabbs.
 * The contribution guide has been improved to include more detail on ways to contribute new features/enhancements, report bugs, and improve or suggest vignettes. By @jamesmbaazam in #464 and reviewed by @seabbs.
 * Updated the code in `inst/CITATION` and added a GitHub Actions workflow to auto-generate `citation.cff` so that the two citation files are always in sync with `DESCRIPTION`. By @jamesmbazam in #467, with contributions from @Bisaloo, and reviewed by @seabbs and @sbfnk.
-* Updated the documentation of the `reported_cases` argument in `estimate_infections()` and `confirm` column in the `obs` argument of `estimate_truncation()` to allow `numeric` types, not just `integer`. See #594, by @jamesmbaazam, and reviewed by @sbfnk.
+* Updated the documentation of the `data` argument in `estimate_infections()` and `confirm` column in the `obs` argument of `estimate_truncation()` to allow `numeric` types, not just `integer`. See #594, by @jamesmbaazam, and reviewed by @sbfnk.
 * Removed the reporting templates that were previously provided. See #604 by @jamesmbaazam, and reviewed by @sbfnk.
 * Clarified how estimated or specified uncertainty around data truncation can be passed to `epinow()`, `regional_epinow()`, and `estimate_infections()` using the `truncation` argument. By @jamesmbaazam in #644 and reviewed by @sbnfk.
+* Internal functions have been removed from the pkgdown index. By @sbfnk in #735.
 
 ## Package
 
@@ -43,9 +52,7 @@
 * Replaced descriptions and plot labels to be more general and clearer. By @sbfnk in #621 and reviewed by @jamesmbaazam.
 * Argument choices have been moved into default arguments. By @sbfnk in #622 and reviewed by @seabbs.
 * `simulate_infections()` gained the argument `seeding_time` to change the seeding time. Additionally, the documentation was improved. By @sbfnk in #627 and reviewed by @jamesmbaazam.
-* The model-specific `weigh_delay_priors` argument has been deprecated in favour of delay-specific prior weighting using `weight_priors`. See `generation_time_opts()`, `delay_opts()`, and `trunc_opts()`. By @sbfnk in #630 and reviewed by @jamesmbaazam.
 * A `cmdstanr` backend has been added. By @sbfnk in #537 and #642 and reviewed by @seabbs.
-* All functions now use a `data` argument to pass data. The existing `reported_cases`, `reports`, and `obs` arguments are deprecated and will be removed in v2.0.0. By @jamesmbaazam in #638 and reviewed by @sbfnk.
 
 ## Model changes
 
diff --git a/R/EpiNow2-package.R b/R/EpiNow2-package.R
index 39781ef5b..63503a710 100644
--- a/R/EpiNow2-package.R
+++ b/R/EpiNow2-package.R
@@ -10,5 +10,7 @@
 # roxygen namespace tags. Modify with care!
 ## usethis namespace: start
 #' @importFrom lifecycle deprecate_soft
+#' @importFrom lifecycle deprecate_warn
+#' @importFrom lifecycle deprecate_stop
 ## usethis namespace: end
 NULL
diff --git a/R/adjust.R b/R/adjust.R
deleted file mode 100644
index 6f9c1dbe1..000000000
--- a/R/adjust.R
+++ /dev/null
@@ -1,372 +0,0 @@
-#' Adjust from Case Counts by Infection Date to Date of Report
-#'
-#' @description  `r lifecycle::badge("deprecated")`
-#' Maps from cases by date of infection to date of report via date of
-#' onset.
-#' @param infections `<data.table>` containing a `date` variable and a numeric
-#' `cases` variable.
-#'
-#' @param delay_defs A list of single row data.tables that each  defines a
-#' delay distribution (model, parameters and maximum delay for each model).
-#' See [lognorm_dist_def()] for an example of the structure.
-#'
-#' @param reporting_effect A numeric vector of length 7 that allows the scaling
-#' of reported cases by the day on which they report (1 = Monday, 7 = Sunday).
-#' By default no scaling occurs.
-#'
-#' @param reporting_model A function that takes a single numeric vector as an
-#' argument and returns a single numeric vector. Can be used to apply stochastic
-#' reporting effects. See the examples for details.
-#'
-#' @return A `data.table` containing a `date` variable (date of report) and a
-#' `cases` variable. If `return_onset = TRUE` there will be a third variable
-#' `reference` which indicates what the date variable refers to.
-#' @export
-#' @inheritParams sample_approx_dist
-#' @importFrom data.table setorder data.table data.table
-#' @importFrom lubridate wday
-#' @importFrom lifecycle deprecate_warn
-#' @examples
-#' \donttest{
-#' # This function is deprecated and its functionality can now be accessed
-#' # from [simulate_secondary()].
-#' # Here are some examples of how to use [simulate_secondary()] to replace
-#' # adjust_infection_to_report().
-#'
-#' # Old (using adjust_infection_to_report()):
-#' # Define example case data
-#' cases <- data.table::copy(example_confirmed)
-#' cases <- cases[, cases := as.integer(confirm)]
-#' # Define a simple reporting delay distribution
-#' delay_def <- lognorm_dist_def(
-#'  mean = 5, mean_sd = 1, sd = 3, sd_sd = 1,
-#'  max_value = 30, samples = 1, to_log = TRUE
-#' )
-#' report <- adjust_infection_to_report(
-#'  cases,
-#'  delay_defs = list(delay_def),
-#'  reporting_model = function(n) rpois(length(n), n)
-#' )
-#' print(report)
-#'
-#' # New (using simulate_secondary()):
-#' cases <- data.table::copy(example_confirmed)
-#' cases <- cases[, primary := as.integer(confirm)]
-#' uncertain_delay <- LogNormal(
-#'  mean = Normal(5, 1), sd = Normal(3, 1),
-#'  max = 30
-#' )
-#' delay <- fix_dist(uncertain_delay, strategy = "sample")
-#' report <- simulate_secondary(
-#'  cases,
-#'  delays = delay_opts(delay),
-#'  obs = obs_opts(family = "poisson")
-#' )
-#' print(report)
-#' }
-adjust_infection_to_report <- function(infections, delay_defs,
-                                       reporting_model, reporting_effect,
-                                       type = "sample",
-                                       truncate_future = TRUE) {
-  deprecate_warn(
-    when = "1.5.0",
-    what = "adjust_infection_to_report()",
-    with = "simulate_secondary()",
-    details = c(
-      "See equivalent examples using `simulate_secondary()`",
-      "in ?adjust_infection_to_report.",
-      "This function will be removed completely in version 2.0.0."
-    )
-  )
-  # Reset DT Defaults on Exit
-  set_dt_single_thread()
-
-  ## deprecated
-  sample_single_dist <- function(input, delay_def) {
-    ## Define sample delay fn
-    sample_delay_fn <- function(n, ...) {
-      EpiNow2::dist_skel(
-        n = n,
-        model = delay_def$model[[1]],
-        params = delay_def$params[[1]],
-        max_value = delay_def$max_value[[1]],
-        ...
-      )
-    }
-
-
-    ## Infection to onset
-    out <- suppressWarnings(EpiNow2::sample_approx_dist(
-      cases = input,
-      dist_fn = sample_delay_fn,
-      max_value = delay_def$max_value,
-      direction = "forwards",
-      type = type,
-      truncate_future = FALSE
-    ))
-
-    return(out)
-  }
-
-  sample_dist_spec <- function(input, delay_def) {
-    ## Define sample delay fn
-    sample_delay_fn <- function(n, dist, cum, ...) {
-      fixed_dist <- discretise(fix_dist(delay_def, strategy = "sample"))
-      if (dist) {
-        fixed_dist[[1]]$pmf[n + 1]
-      } else {
-        sample(seq_along(fixed_dist[[1]]$pmf) - 1, size = n, replace = TRUE)
-      }
-    }
-
-    ## Infection to onset
-    out <- suppressWarnings(EpiNow2::sample_approx_dist(
-      cases = input,
-      dist_fn = sample_delay_fn,
-      max_value = max(delay_def),
-      direction = "forwards",
-      type = type,
-      truncate_future = FALSE
-    ))
-
-    return(out)
-  }
-
-  if (is(delay_defs, "dist_spec")) {
-    report <- sample_dist_spec(infections, extract_single_dist(delay_defs, 1))
-    if (length(delay_defs) > 1) {
-      for (def in seq(2, length(delay_defs))) {
-        report <- sample_dist_spec(report, extract_single_dist(delay_defs, def))
-      }
-    }
-  } else {
-    deprecate_warn(
-      "1.5.0",
-      "adjust_infection_to_report(delay_defs = 'should be a dist_spec')",
-      details = "Specifying this as a list of data tables is deprecated."
-    )
-    report <- sample_single_dist(infections, delay_defs[[1]])
-    if (length(delay_defs) > 1) {
-      for (def in 2:length(delay_defs)) {
-        report <- sample_single_dist(report, delay_defs[[def]])
-      }
-    }
-  }
-  ## Add a weekly reporting effect if present
-  if (!missing(reporting_effect)) {
-    reporting_effect <- data.table::data.table(
-      day = 1:7,
-      effect = reporting_effect
-    )
-
-    report <- report[, day := lubridate::wday(date, week_start = 1)]
-    report <- report[reporting_effect, on = "day"]
-    report <- report[, cases := as.integer(cases * effect)][
-      ,
-      `:=`(effect = NULL, day = NULL)
-    ]
-
-    report <- data.table::setorder(report, date)
-  }
-
-  if (!missing(reporting_model)) {
-    report <- report[, cases := reporting_model(cases)]
-  }
-
-  ## Truncate reported cases by maximum infection date
-  if (type == "sample" && truncate_future) {
-    report <- report[date <= max(infections$date)]
-  }
-  return(report)
-}
-
-#' Approximate Sampling a Distribution using Counts
-#'
-#' @description `r lifecycle::badge("deprecated")`
-#' Convolves cases by a PMF function. This function will soon be removed or
-#' replaced with a more robust stan implementation.
-#'
-#' @param cases A `<data.frame>` of cases (in date order) with the following
-#' variables: `date` and `cases`.
-#'
-#' @param max_value Numeric, maximum value to allow. Defaults to 120 days
-#'
-#' @param direction Character string, defato "backwards". Direction in which to
-#' map cases. Supports either "backwards" or "forwards".
-#'
-#' @param dist_fn Function that takes two arguments with the first being
-#' numeric and the second being logical (and defined as `dist`). Should return
-#' the probability density or a sample from the defined distribution. See
-#' the examples for more.
-#'
-#' @param earliest_allowed_mapped A character string representing a date
-#' ("2020-01-01"). Indicates the earliest allowed mapped value.
-#'
-#' @param type Character string indicating the method to use to transform
-#' counts. Supports either "sample"  which approximates sampling or "median"
-#' would shift by the median of the distribution.
-#'
-#' @param truncate_future Logical, should cases be truncated if they occur
-#' after the first date reported in the data. Defaults to `TRUE`.
-#'
-#' @return A `<data.table>` of cases by date of onset
-#' @export
-#' @importFrom purrr map_dfc
-#' @importFrom data.table data.table setorder
-#' @importFrom lubridate days
-#' @importFrom lifecycle deprecate_warn
-#' @examples
-#' \donttest{
-#' cases <- example_confirmed
-#' cases <- cases[, cases := as.integer(confirm)]
-#' print(cases)
-#'
-#' # total cases
-#' sum(cases$cases)
-#'
-#' delay_fn <- function(n, dist, cum) {
-#'   if (dist) {
-#'     pgamma(n + 0.9999, 2, 1) - pgamma(n - 1e-5, 2, 1)
-#'   } else {
-#'     as.integer(rgamma(n, 2, 1))
-#'   }
-#' }
-#'
-#' onsets <- sample_approx_dist(
-#'   cases = cases,
-#'   dist_fn = delay_fn
-#' )
-#'
-#' # estimated onset distribution
-#' print(onsets)
-#'
-#' # check that sum is equal to reported cases
-#' total_onsets <- median(
-#'   purrr::map_dbl(
-#'     1:10,
-#'     ~ sum(sample_approx_dist(
-#'       cases = cases,
-#'       dist_fn = delay_fn
-#'     )$cases)
-#'   )
-#' )
-#' total_onsets
-#'
-#'
-#' # map from onset cases to reported
-#' reports <- sample_approx_dist(
-#'   cases = cases,
-#'   dist_fn = delay_fn,
-#'   direction = "forwards"
-#' )
-#'
-#'
-#' # map from onset cases to reported using a mean shift
-#' reports <- sample_approx_dist(
-#'   cases = cases,
-#'   dist_fn = delay_fn,
-#'   direction = "forwards",
-#'   type = "median"
-#' )
-#' }
-sample_approx_dist <- function(cases = NULL,
-                               dist_fn = NULL,
-                               max_value = 120,
-                               earliest_allowed_mapped = NULL,
-                               direction = "backwards",
-                               type = "sample",
-                               truncate_future = TRUE) {
-  deprecate_warn(
-    "1.5.0",
-    "sample_approx_dist()",
-    details = "The function will be removed completely in version 2.0.0."
-  )
-  if (type == "sample") {
-    if (direction == "backwards") {
-      direction_fn <- rev
-    } else if (direction == "forwards") {
-      direction_fn <- function(x) {
-        x
-      }
-    }
-    # reverse cases so starts with current first
-    reversed_cases <- direction_fn(cases$cases)
-    reversed_cases[is.na(reversed_cases)] <- 0
-    # draw from the density fn of the dist
-    draw <- dist_fn(0:max_value, dist = TRUE, cum = FALSE)
-
-    # approximate cases
-    mapped_cases <- do.call(cbind, purrr::map(
-      seq_along(reversed_cases),
-      ~ c(
-        rep(0, . - 1),
-        stats::rbinom(
-          length(draw),
-          rep(reversed_cases[.], length(draw)),
-          draw
-        ),
-        rep(0, length(reversed_cases) - .)
-      )
-    ))
-
-
-    # set dates order based on direction mapping
-    if (direction == "backwards") {
-      dates <- seq(min(cases$date) - lubridate::days(length(draw) - 1),
-        max(cases$date),
-        by = "days"
-      )
-    } else if (direction == "forwards") {
-      dates <- seq(min(cases$date),
-        max(cases$date) + lubridate::days(length(draw) - 1),
-        by = "days"
-      )
-    }
-
-    # summarises movements and sample for placement of non-integer cases
-    case_sum <- direction_fn(rowSums(mapped_cases))
-    floor_case_sum <- floor(case_sum)
-    sample_cases <- floor_case_sum +
-      as.numeric((runif(seq_along(case_sum)) < (case_sum - floor_case_sum)))
-
-    # summarise imputed onsets and build output data.table
-    mapped_cases <- data.table::data.table(
-      date = dates,
-      cases = sample_cases
-    )
-
-    # filter out all zero cases until first recorded case
-    mapped_cases <- data.table::setorder(mapped_cases, date)
-    mapped_cases <- mapped_cases[
-      ,
-      cum_cases := cumsum(cases)
-    ][cum_cases != 0][, cum_cases := NULL]
-  } else if (type == "median") {
-    shift <- as.integer(
-      median(as.integer(dist_fn(1000, dist = FALSE)), na.rm = TRUE)
-    )
-
-    if (direction == "backwards") {
-      mapped_cases <- data.table::copy(cases)[
-        ,
-        date := date - lubridate::days(shift)
-      ]
-    } else if (direction == "forwards") {
-      mapped_cases <- data.table::copy(cases)[
-        ,
-        date := date + lubridate::days(shift)
-      ]
-    }
-  }
-
-  if (!is.null(earliest_allowed_mapped)) {
-    mapped_cases <- mapped_cases[date >= as.Date(earliest_allowed_mapped)]
-  }
-
-  # filter out future cases
-  if (direction == "forwards" && truncate_future) {
-    max_date <- max(cases$date)
-    mapped_cases <- mapped_cases[date <= max_date]
-  }
-  return(mapped_cases)
-}
diff --git a/R/checks.R b/R/checks.R
index 68f556453..01dfb7db8 100644
--- a/R/checks.R
+++ b/R/checks.R
@@ -1,6 +1,6 @@
 #' Validate data input
 #'
-#' @description
+#' @description `r lifecycle::badge("stable")`
 #' `check_reports_valid()` checks that the supplied data is a `<data.frame>`,
 #' and that it has the right column names and types. In particular, it checks
 #' that the date column is in date format and does not contain NA's, and that
@@ -54,7 +54,7 @@ check_reports_valid <- function(data,
 
 #' Validate probability distribution for passing to stan
 #'
-#' @description
+#' @description `r lifecycle::badge("stable")`
 #' `check_stan_delay()` checks that the supplied data is a `<dist_spec>`,
 #' that it is a supported distribution, and that is has a finite maximum.
 #'
diff --git a/R/create.R b/R/create.R
index d3e146334..f02bfdc59 100644
--- a/R/create.R
+++ b/R/create.R
@@ -23,9 +23,11 @@
 #' @inheritParams estimate_infections
 #' @importFrom data.table copy merge.data.table setorder setDT frollsum
 #' @return A cleaned data frame of reported cases
-#' @export
+#' @keywords internal
 #' @examples
+#' \dontrun{
 #' create_clean_reported_cases(example_confirmed, 7)
+#' }
 create_clean_reported_cases <- function(data, horizon = 0,
                                         filter_leading_zeros = TRUE,
                                         zero_threshold = Inf,
@@ -126,8 +128,9 @@ create_complete_cases <- function(cases) {
 #' @importFrom stats lm
 #' @importFrom runner mean_run
 #' @return A `<data.frame>` for shifted reported cases
-#' @export
+#' @keywords internal
 #' @examples
+#' \dontrun{
 #' shift <- 7
 #' horizon <- 7
 #' smoothing_window <- 14
@@ -146,6 +149,7 @@ create_complete_cases <- function(cases) {
 #'    cases
 #'  ))
 #' create_shifted_cases(cases, shift, smoothing_window, horizon)
+#' }
 create_shifted_cases <- function(data, shift,
                                  smoothing_window, horizon) {
   shifted_reported_cases <- data.table::copy(data)[
@@ -211,6 +215,7 @@ create_shifted_cases <- function(data, shift,
 #'
 #' @param delay Numeric mean delay
 #' @importFrom rlang arg_match
+#' @keywords internal
 #' @return A list containing a logical called fixed and an integer called from
 create_future_rt <- function(future = c("latest", "project", "estimate"),
                              delay = 0) {
@@ -245,8 +250,9 @@ create_future_rt <- function(future = c("latest", "project", "estimate"),
 #' @seealso rt_settings
 #' @return A list of settings defining the time-varying reproduction number
 #' @inheritParams create_future_rt
-#' @export
+#' @keywords internal
 #' @examples
+#' \dontrun{
 #' # default Rt data
 #' create_rt_data()
 #'
@@ -255,6 +261,7 @@ create_future_rt <- function(future = c("latest", "project", "estimate"),
 #'
 #' # using breakpoints
 #' create_rt_data(rt_opts(use_breakpoints = TRUE), breakpoints = rep(1, 10))
+#' }
 create_rt_data <- function(rt = rt_opts(), breakpoints = NULL,
                            delay = 0, horizon = 0) {
   # Define if GP is on or off
@@ -311,9 +318,7 @@ create_rt_data <- function(rt = rt_opts(), breakpoints = NULL,
 #' @seealso backcalc_opts
 #' @importFrom data.table fcase
 #' @return A list of settings defining the Gaussian process
-#' @export
-#' @examples
-#' create_backcalc_data(backcalc = backcalc_opts())
+#' @keywords internal
 create_backcalc_data <- function(backcalc = backcalc_opts()) {
   data <- list(
    rt_half_window = as.integer((backcalc$rt_window - 1) / 2),
@@ -339,8 +344,9 @@ create_backcalc_data <- function(backcalc = backcalc_opts()) {
 #' @importFrom data.table fcase
 #' @seealso [gp_opts()]
 #' @return A list of settings defining the Gaussian process
-#' @export
+#' @keywords internal
 #' @examples
+#' \dontrun{
 #' # define input data required
 #' data <- list(
 #'   t = 30,
@@ -356,6 +362,7 @@ create_backcalc_data <- function(backcalc = backcalc_opts()) {
 #'
 #' # custom lengthscale
 #' create_gp_data(gp_opts(ls_mean = 14), data)
+#' }
 create_gp_data <- function(gp = gp_opts(), data) {
   # Define if GP is on or off
   if (is.null(gp)) {
@@ -408,8 +415,9 @@ create_gp_data <- function(gp = gp_opts(), data) {
 #' @seealso [obs_opts()]
 #' @return A list of settings ready to be passed to stan defining
 #' the Observation Model
-#' @export
+#' @keywords internal
 #' @examples
+#' \dontrun{
 #' dates <- seq(as.Date("2020-03-15"), by = "days", length.out = 15)
 #' # default observation model data
 #' create_obs_model(dates = dates)
@@ -424,6 +432,7 @@ create_gp_data <- function(gp = gp_opts(), data) {
 #'
 #' # Apply a custom week week length
 #' create_obs_model(obs_opts(week_length = 3), dates = dates)
+#' }
 create_obs_model <- function(obs = obs_opts(), dates) {
   data <- list(
     model_type = as.numeric(obs$family == "negbin"),
@@ -464,12 +473,14 @@ create_obs_model <- function(obs = obs_opts(), dates) {
 #' @importFrom stats lm
 #' @importFrom purrr safely
 #' @return A list of stan data
-#' @export
+#' @keywords internal
 #' @examples
+#' \dontrun{
 #' create_stan_data(
 #'  example_confirmed, 7, rt_opts(), gp_opts(), obs_opts(), 7,
 #'  backcalc_opts(), create_shifted_cases(example_confirmed, 7, 14, 7)
 #' )
+#' }
 create_stan_data <- function(data, seeding_time,
                              rt, gp, obs, horizon,
                              backcalc, shifted_cases) {
@@ -573,7 +584,7 @@ create_delay_inits <- function(data) {
 #' @importFrom purrr map2_dbl
 #' @importFrom truncnorm rtruncnorm
 #' @importFrom data.table fcase
-#' @export
+#' @keywords internal
 create_initial_conditions <- function(data) {
   init_fun <- function() {
     out <- create_delay_inits(data)
@@ -673,13 +684,15 @@ create_initial_conditions <- function(data) {
 #' @importFrom utils modifyList
 #'
 #' @return A list of stan arguments
-#' @export
+#' @keywords internal
 #' @examples
+#' \dontrun{
 #' # default settings
 #' create_stan_args()
 #'
 #' # increasing warmup
 #' create_stan_args(stan = stan_opts(warmup = 1000))
+#' }
 create_stan_args <- function(stan = stan_opts(),
                              data = NULL,
                              init = "random",
@@ -723,6 +736,7 @@ create_stan_args <- function(stan = stan_opts(),
 ##'   determines weight associated with weighted delay priors; default: 1
 ##' @return A list of variables as expected by the stan model
 ##' @importFrom purrr transpose map flatten
+##' @keywords internal
 create_stan_delays <- function(..., time_points = 1L) {
   delays <- list(...)
   ## discretise
diff --git a/R/deprecated.R b/R/deprecated.R
index 73dec0a3d..06f47257e 100644
--- a/R/deprecated.R
+++ b/R/deprecated.R
@@ -1,3 +1,311 @@
+#' Adjust from Case Counts by Infection Date to Date of Report
+#'
+#' @description  `r lifecycle::badge("deprecated")`
+#' Maps from cases by date of infection to date of report via date of
+#' onset.
+#' @param infections `<data.table>` containing a `date` variable and a numeric
+#' `cases` variable.
+#'
+#' @param delay_defs A list of single row data.tables that each  defines a
+#' delay distribution (model, parameters and maximum delay for each model).
+#' See [lognorm_dist_def()] for an example of the structure.
+#'
+#' @param reporting_effect A numeric vector of length 7 that allows the scaling
+#' of reported cases by the day on which they report (1 = Monday, 7 = Sunday).
+#' By default no scaling occurs.
+#'
+#' @param reporting_model A function that takes a single numeric vector as an
+#' argument and returns a single numeric vector. Can be used to apply stochastic
+#' reporting effects. See the examples for details.
+#'
+#' @return A `data.table` containing a `date` variable (date of report) and a
+#' `cases` variable. If `return_onset = TRUE` there will be a third variable
+#' `reference` which indicates what the date variable refers to.
+#' @keywords internal
+#' @export
+#' @inheritParams sample_approx_dist
+#' @importFrom data.table setorder data.table data.table
+#' @importFrom lubridate wday
+#' @examples
+#' \donttest{
+#' # This function is deprecated and its functionality can now be accessed
+#' # from [simulate_secondary()].
+#' # Here are some examples of how to use [simulate_secondary()] to replace
+#' # adjust_infection_to_report().
+#'
+#' # Old (using adjust_infection_to_report()):
+#' # Define example case data
+#' cases <- data.table::copy(example_confirmed)
+#' cases <- cases[, cases := as.integer(confirm)]
+#' # Define a simple reporting delay distribution
+#' delay_def <- lognorm_dist_def(
+#'  mean = 5, mean_sd = 1, sd = 3, sd_sd = 1,
+#'  max_value = 30, samples = 1, to_log = TRUE
+#' )
+#' report <- adjust_infection_to_report(
+#'  cases,
+#'  delay_defs = list(delay_def),
+#'  reporting_model = function(n) rpois(length(n), n)
+#' )
+#' print(report)
+#'
+#' # New (using simulate_secondary()):
+#' cases <- data.table::copy(example_confirmed)
+#' cases <- cases[, primary := as.integer(confirm)]
+#' uncertain_delay <- LogNormal(
+#'  mean = Normal(5, 1), sd = Normal(3, 1),
+#'  max = 30
+#' )
+#' delay <- fix_dist(uncertain_delay, strategy = "sample")
+#' report <- simulate_secondary(
+#'  cases,
+#'  delays = delay_opts(delay),
+#'  obs = obs_opts(family = "poisson")
+#' )
+#' print(report)
+#' }
+adjust_infection_to_report <- function(infections, delay_defs,
+                                       reporting_model, reporting_effect,
+                                       type = "sample",
+                                       truncate_future = TRUE) {
+  deprecate_warn(
+    when = "1.5.0",
+    what = "adjust_infection_to_report()",
+    with = "simulate_secondary()",
+    details = c(
+      "See equivalent examples using `simulate_secondary()`",
+      "in ?adjust_infection_to_report.",
+      "This function will be removed completely in the next version."
+    )
+  )
+  # Reset DT Defaults on Exit
+  set_dt_single_thread()
+
+  ## deprecated
+  sample_single_dist <- function(input, delay_def) {
+    ## Define sample delay fn
+    sample_delay_fn <- function(n, ...) {
+      EpiNow2::dist_skel(
+        n = n,
+        model = delay_def$model[[1]],
+        params = delay_def$params[[1]],
+        max_value = delay_def$max_value[[1]],
+        ...
+      )
+    }
+
+
+    ## Infection to onset
+    out <- suppressWarnings(EpiNow2::sample_approx_dist(
+      cases = input,
+      dist_fn = sample_delay_fn,
+      max_value = delay_def$max_value,
+      direction = "forwards",
+      type = type,
+      truncate_future = FALSE
+    ))
+
+    return(out)
+  }
+
+  sample_dist_spec <- function(input, delay_def) {
+    ## Define sample delay fn
+    sample_delay_fn <- function(n, dist, cum, ...) {
+      fixed_dist <- discretise(fix_dist(delay_def, strategy = "sample"))
+      if (dist) {
+        fixed_dist[[1]]$pmf[n + 1]
+      } else {
+        sample(seq_along(fixed_dist[[1]]$pmf) - 1, size = n, replace = TRUE)
+      }
+    }
+
+    ## Infection to onset
+    out <- suppressWarnings(EpiNow2::sample_approx_dist(
+      cases = input,
+      dist_fn = sample_delay_fn,
+      max_value = max(delay_def),
+      direction = "forwards",
+      type = type,
+      truncate_future = FALSE
+    ))
+
+    return(out)
+  }
+
+  if (is(delay_defs, "dist_spec")) {
+    report <- sample_dist_spec(infections, extract_single_dist(delay_defs, 1))
+    if (length(delay_defs) > 1) {
+      for (def in seq(2, length(delay_defs))) {
+        report <- sample_dist_spec(report, extract_single_dist(delay_defs, def))
+      }
+    }
+  } else {
+    deprecate_warn(
+      "1.5.0",
+      "adjust_infection_to_report(delay_defs = 'should be a dist_spec')",
+      details = "Specifying this as a list of data tables is deprecated."
+    )
+    report <- sample_single_dist(infections, delay_defs[[1]])
+    if (length(delay_defs) > 1) {
+      for (def in 2:length(delay_defs)) {
+        report <- sample_single_dist(report, delay_defs[[def]])
+      }
+    }
+  }
+  ## Add a weekly reporting effect if present
+  if (!missing(reporting_effect)) {
+    reporting_effect <- data.table::data.table(
+      day = 1:7,
+      effect = reporting_effect
+    )
+
+    report <- report[, day := lubridate::wday(date, week_start = 1)]
+    report <- report[reporting_effect, on = "day"]
+    report <- report[, cases := as.integer(cases * effect)][
+      ,
+      `:=`(effect = NULL, day = NULL)
+    ]
+
+    report <- data.table::setorder(report, date)
+  }
+
+  if (!missing(reporting_model)) {
+    report <- report[, cases := reporting_model(cases)]
+  }
+
+  ## Truncate reported cases by maximum infection date
+  if (type == "sample" && truncate_future) {
+    report <- report[date <= max(infections$date)]
+  }
+  return(report)
+}
+
+#' Specify a distribution.
+#'
+#' @description `r lifecycle::badge("deprecated")`
+#' This function is deprecated as a user-facing function (while its
+#' functionality is still used internally). Construct distributions using
+#' the corresponding distribution function such as [Gamma()], [LogNormal()],
+#' [Normal()] or [Fixed()] instead.
+#'
+#' @param distribution Character, defaults to "lognormal". The (discretised)
+#' distribution to be used. Can be "lognormal", "gamma", "normal" or "fixed".
+#' The corresponding parameters (defined in [natural_params()]) are passed
+#' as `params_mean`,  and their uncertainty as `params_sd`.
+#'
+#' @param params_mean Numeric. Central values of the parameters of the
+#' distribution as defined in [natural_params().
+#'
+#' @param params_sd Numeric. Standard deviations of the parameters of the
+#' distribution as defined in [natural_params().
+#'
+#' @param max Numeric, maximum value of the distribution. The distribution will
+#' be truncated at this value. Default: `Inf`, i.e. no maximum.
+#'
+#' @param mean Deprecated; use `params_mean` instead.
+#'
+#' @param sd Deprecated; use `params_mean` instead.
+#'
+#' @param mean_sd Deprecated; use `params_sd` instead.
+#'
+#' @param sd_sd Deprecated; use `params_sd` instead.
+#'
+#' @param pmf Numeric, a vector of values that represent the (nonparametric)
+#' probability mass function of the delay (starting with 0); defaults to an
+#' empty vector corresponding to a parametric specification of the distribution
+#' (using \code{params_mean}, and \code{params_sd}.
+#' @param fixed Deprecated, use [fix_dist()] instead.
+#' @return A list of distribution options.
+#' @importFrom rlang warn arg_match
+#' @export
+#' @keywords internal
+dist_spec <- function(distribution = c(
+                        "lognormal", "normal", "gamma", "fixed", "empty"
+                      ),
+                      params_mean = numeric(0), params_sd = numeric(0),
+                      mean, sd = 0, mean_sd = 0, sd_sd = 0,
+                      max = Inf, pmf = numeric(0), fixed = FALSE) {
+
+  lifecycle::deprecate_warn(
+    "1.5.0",
+    "dist_spec()",
+    details = c(
+      paste0(
+        "Please use distribution functions such as `Gamma()` or `Lognormal()` ",
+        "instead."
+      ),
+      "The function will become internal only in the next version."
+    )
+  )
+  ## check for deprecated parameters
+  if (!missing(fixed)) {
+    lifecycle::deprecate_warn(
+      "1.5.0",
+      "dist_spec(fixed)",
+      "fix_dist()"
+    )
+    params_sd <- NULL
+  }
+  ## check for deprecated parameters
+  if (!all(missing(mean), missing(sd), missing(mean_sd), missing(sd_sd)) &&
+      (length(params_mean) > 0 || length(params_sd) > 0)) {
+    stop("Distributional parameters should not be given as `mean`, `sd`, etc. ",
+         "in addition to `params_mean` or `params_sd`")
+  }
+  distribution <- match.arg(distribution)
+  ## check if distribution is given as empty and warn about deprecation if so
+  if (distribution == "empty") {
+    deprecate_warn(
+      "1.5.0",
+      "dist_spec(distribution = 'must not be \"empty\"')",
+      details = "Please use `Fixed(0)` instead."
+    )
+  }
+
+  if (!all(missing(mean), missing(sd), missing(mean_sd), missing(sd_sd))) {
+    if (sd == 0 && mean_sd == 0 && sd_sd == 0) {
+      distribution <- "fixed"
+    }
+    ## deprecated arguments given
+    if (distribution == "lognormal") {
+      params_mean <- c(meanlog = mean, sdlog = sd)
+      params_sd <- c(meanlog = mean_sd, sdlog = sd_sd)
+    } else if (distribution == "gamma") {
+      temp_dist <- Gamma(
+        mean = Normal(mean, mean_sd),
+        sd = Normal(sd, sd_sd)
+      )
+      params_mean <- vapply(get_parameters(temp_dist), mean, numeric(1))
+      params_sd <- vapply(get_parameters(temp_dist), sd_dist, numeric(1))
+    } else if (distribution == "normal") {
+      params_mean <- c(mean = mean, sd = sd)
+      params_sd <- c(mean = mean_sd, sd = sd_sd)
+    } else if (distribution == "fixed") {
+      params_mean <- mean
+    }
+  }
+  if (length(pmf) > 0) {
+    if (!all(
+      missing(mean), missing(sd), missing(mean_sd), missing(sd_sd),
+      missing(params_mean), missing(params_sd)
+    )) {
+      stop("Distributional parameters should not be given in addition to `pmf`")
+    }
+    distribution <- "nonparametric"
+    parameters <- list(pmf = pmf)
+  } else {
+    if (length(params_sd) == 0) {
+      params_sd <- rep(0, length(params_mean))
+    }
+    parameters <- lapply(seq_along(params_mean), function(id) {
+      Normal(params_mean[id], params_sd[id])
+    })
+    names(parameters) <- natural_params(distribution)
+    parameters$max <- max
+  }
+  return(new_dist_spec(parameters, distribution))
+}
+
 #' Generate a Gamma Distribution Definition Based on Parameter Estimates
 #'
 #' @description `r lifecycle::badge("deprecated")`
@@ -15,6 +323,7 @@
 #'
 #' @importFrom truncnorm rtruncnorm
 #' @return A `<data.table>` defining the distribution as used by [dist_skel()]
+#' @keywords internal
 #' @export
 #' @inheritParams dist_skel
 #' @inheritParams lognorm_dist_def
@@ -43,7 +352,7 @@ gamma_dist_def <- function(shape, shape_sd,
                            max_value, samples) {
   lifecycle::deprecate_warn(
     "1.5.0", "gamma_dist_def()", "Gamma()",
-    "The function will be removed completely in version 2.0.0."
+    "The function will be removed completely in the next version."
   )
 
   if (missing(shape) && missing(scale) && !missing(mean) && !missing(sd)) {
@@ -85,6 +394,85 @@ gamma_dist_def <- function(shape, shape_sd,
   return(dist)
 }
 
+#' Generate initial conditions by fitting to cumulative cases
+#'
+#' @description `r lifecycle::badge("deprecated")`
+#'
+#' Fits a model to cumulative cases. This may be a useful approach to
+#' initialising a full model fit for certain data sets where the sampler gets
+#' stuck or cannot easily be initialised as fitting to cumulative cases changes
+#' the shape of the posterior distribution. In `estimate_infections()`,
+#' `epinow()` and `regional_epinow()` this option can be engaged by setting
+#' `stan_opts(init_fit = "cumulative")`.
+#'
+#' This implementation is based on the approach taken in
+#' [epidemia](https://github.com/ImperialCollegeLondon/epidemia/) authored by
+#' James Scott.
+#'
+#' @param samples Numeric, defaults to 50. Number of posterior samples.
+#'
+#' @param warmup Numeric, defaults to 50. Number of warmup samples.
+#'
+#' @param verbose Logical, should fitting progress be returned. Defaults to
+#' `FALSE`.
+#'
+#' @inheritParams create_initial_conditions
+#' @importFrom rstan sampling
+#' @importFrom futile.logger flog.debug
+#' @importFrom utils capture.output
+#' @inheritParams fit_model_with_nuts
+#' @inheritParams stan_opts
+#' @return A stanfit object
+#' @keywords internal
+init_cumulative_fit <- function(args, samples = 50, warmup = 50,
+                                id = "init", verbose = FALSE,
+                                backend = "rstan") {
+  deprecate_warn(
+    when = "1.5.0",
+    what = "init_cumulative_fit()"
+  )
+  futile.logger::flog.debug(
+    "%s: Fitting to cumulative data to initialise chains", id,
+    name = "EpiNow2.epinow.estimate_infections.fit"
+  )
+  # copy main run settings and override to use only 100 iterations and a single
+  # chain
+  initial_args <- create_stan_args(
+    stan = stan_opts(
+      args$object,
+      samples = samples,
+      warmup = warmup,
+      control = list(adapt_delta = 0.9, max_treedepth = 13),
+      chains = 1,
+      cores = 2,
+      backend = backend,
+      open_progress = FALSE,
+      show_messages = FALSE
+    ),
+    data = args$data, init = args$init
+  )
+  # change observations to be cumulative in order to protect against noise and
+  # give an approximate fit (though for Rt constrained to be > 1)
+  initial_args$data$cases <- cumsum(initial_args$data$cases)
+  initial_args$data$shifted_cases <- cumsum(initial_args$data$shifted_cases)
+
+  # initial fit
+  if (verbose) {
+    fit <- fit_model(initial_args, id = "init_cumulative")
+  } else {
+    out <- tempfile(tmpdir = tempdir(check = TRUE))
+    capture.output(
+      {
+        fit <- fit_model(initial_args, id = "init_cumulative")
+      },
+      type = c("output", "message"),
+      split = FALSE,
+      file = out
+    )
+  }
+  return(fit)
+}
+
 #' Generate a Log Normal Distribution Definition Based on Parameter Estimates
 #'
 #' @description `r lifecycle::badge("deprecated")`
@@ -107,6 +495,7 @@ gamma_dist_def <- function(shape, shape_sd,
 #' @return A `<data.table>` defining the distribution as used by [dist_skel()]
 #' @importFrom truncnorm rtruncnorm
 #' @export
+#' @keywords internal
 #' @inheritParams dist_skel
 #' @examples
 #' def <- lognorm_dist_def(
@@ -131,7 +520,7 @@ lognorm_dist_def <- function(mean, mean_sd,
                              to_log = FALSE) {
   lifecycle::deprecate_warn(
     "1.5.0", "lognorm_dist_def()", "LogNormal()",
-    "The function will be removed completely in version 2.0.0."
+    "The function will be removed completely in the next version."
   )
 
   transform_mean <- function(mu, sig) {
@@ -179,127 +568,322 @@ lognorm_dist_def <- function(mean, mean_sd,
   return(dist)
 }
 
-#' Specify a distribution.
+#' Report case counts by date of report
 #'
 #' @description `r lifecycle::badge("deprecated")`
-#' This function is deprecated as a user-facing function (while its
-#' functionality is still used internally). Construct distributions using
-#' the corresponding distribution function such as [Gamma()], [LogNormal()],
-#' [Normal()] or [Fixed()] instead.
+#' Convolves latent infections to reported cases via an observation model.
+#' Likely to be removed/replaced in later releases by functionality drawing on
+#' the `stan` implementation.
 #'
-#' @param distribution Character, defaults to "lognormal". The (discretised)
-#' distribution to be used. Can be "lognormal", "gamma", "normal" or "fixed".
-#' The corresponding parameters (defined in [natural_params()]) are passed
-#' as `params_mean`,  and their uncertainty as `params_sd`.
+#' @param case_estimates A data.table of case estimates with the following
+#' variables: date, sample, cases
 #'
-#' @param params_mean Numeric. Central values of the parameters of the
-#' distribution as defined in [natural_params().
+#' @param case_forecast A data.table of case forecasts with the following
+#' variables: date, sample, cases. If not supplied the default is not to
+#' incorporate forecasts.
 #'
-#' @param params_sd Numeric. Standard deviations of the parameters of the
-#' distribution as defined in [natural_params().
+#' @param reporting_effect A `data.table` giving the weekly reporting effect
+#'  with the following variables: `sample` (must be the same as in `nowcast`),
+#'  `effect` (numeric scaling factor for each weekday),`day` (numeric 1 - 7
+#'  (1 = Monday and 7 = Sunday)). If not supplied then no weekly reporting
+#'  effect is assumed.
 #'
-#' @param max Numeric, maximum value of the distribution. The distribution will
-#' be truncated at this value. Default: `Inf`, i.e. no maximum.
+#' @return A list of `data.table`s. The first entry contains the following
+#' variables `sample`, `date` and `cases` with the second being summarised
+#' across samples.
 #'
-#' @param mean Deprecated; use `params_mean` instead.
-#'
-#' @param sd Deprecated; use `params_mean` instead.
-#'
-#' @param mean_sd Deprecated; use `params_sd` instead.
-#'
-#' @param sd_sd Deprecated; use `params_sd` instead.
-#'
-#' @param pmf Numeric, a vector of values that represent the (nonparametric)
-#' probability mass function of the delay (starting with 0); defaults to an
-#' empty vector corresponding to a parametric specification of the distribution
-#' (using \code{params_mean}, and \code{params_sd}.
-#' @param fixed Deprecated, use [fix_dist()] instead.
-#' @return A list of distribution options.
-#' @importFrom rlang warn arg_match
+#' @keywords internal
 #' @export
-dist_spec <- function(distribution = c(
-                        "lognormal", "normal", "gamma", "fixed", "empty"
-                      ),
-                      params_mean = numeric(0), params_sd = numeric(0),
-                      mean, sd = 0, mean_sd = 0, sd_sd = 0,
-                      max = Inf, pmf = numeric(0), fixed = FALSE) {
-
-  lifecycle::deprecate_warn(
-    "1.5.0",
-    "dist_spec()",
+#' @inheritParams estimate_infections
+#' @inheritParams adjust_infection_to_report
+#' @importFrom data.table data.table rbindlist
+#' @importFrom future.apply future_lapply
+#' @examples
+#' \donttest{
+#' # This function is deprecated and its functionality can now be accessed
+#' # from [simulate_secondary()].
+#' # Here are some examples of how to use [simulate_secondary()] to replace
+#' # report_cases().
+#' # Old (using report_cases()):
+#' # Define case data
+#' cases <- example_confirmed[1:40]
+#' cases <- cases[, cases := as.integer(confirm)]
+#' cases <- cases[, confirm := NULL][, sample := 1]
+#' reported_cases <- report_cases(
+#'  case_estimates = cases,
+#'  delays = delay_opts(example_incubation_period + example_reporting_delay),
+#'  type = "sample"
+#' )
+#' print(reported_cases$samples)
+#'
+#' # New (using simulate_secondary()):
+#' cases <- example_confirmed[1:40]
+#' cases <- cases[, primary := as.integer(confirm)]
+#' report <- simulate_secondary(
+#'  cases,
+#'  delays = delay_opts(
+#'    fix_dist(example_incubation_period + example_reporting_delay)
+#'  ),
+#'  obs = obs_opts(family = "poisson")
+#' )
+#' print(report)
+#' }
+report_cases <- function(case_estimates,
+                         case_forecast = NULL,
+                         delays,
+                         type = "sample",
+                         reporting_effect,
+                         CrIs = c(0.2, 0.5, 0.9)) {
+  deprecate_warn(
+    when = "1.5.0",
+    what = "report_cases()",
+    with = "simulate_secondary()",
     details = c(
-      paste0(
-        "Please use distribution functions such as `Gamma()` or `Lognormal()` ",
-        "instead."
-      ),
-      "The function will become internal only in version 2.0.0."
+      "See equivalent examples using `simulate_secondary()`",
+      "in ?report_cases.",
+      "This function will be removed completely in the next version."
     )
   )
-  ## check for deprecated parameters
-  if (!missing(fixed)) {
-    lifecycle::deprecate_warn(
-      "1.5.0",
-      "dist_spec(fixed)",
-      "fix_dist()"
+  samples <- length(unique(case_estimates$sample))
+
+  # add a null reporting effect if missing
+  if (missing(reporting_effect)) {
+    reporting_effect <- data.table::data.table(
+      sample = list(1:samples),
+      effect = rep(1, 7),
+      day = 1:7
     )
-    params_sd <- NULL
+    reporting_effect <- reporting_effect[,
+      .(sample = unlist(sample)), by = .(effect, day)
+    ]
   }
-  ## check for deprecated parameters
-  if (!all(missing(mean), missing(sd), missing(mean_sd), missing(sd_sd)) &&
-      (length(params_mean) > 0 || length(params_sd) > 0)) {
-    stop("Distributional parameters should not be given as `mean`, `sd`, etc. ",
-         "in addition to `params_mean` or `params_sd`")
-  }
-  distribution <- match.arg(distribution)
-  ## check if distribution is given as empty and warn about deprecation if so
-  if (distribution == "empty") {
-    deprecate_warn(
-      "1.5.0",
-      "dist_spec(distribution = 'must not be \"empty\"')",
-      details = "Please use `Fixed(0)` instead."
-    )
+  # filter and sum nowcast to use only upscaled cases by date of infection
+  infections <- data.table::copy(case_estimates)
+
+  # add in case forecast if present
+  if (!is.null(case_forecast)) {
+    infections <- data.table::rbindlist(list(
+      infections,
+      case_forecast[, .(date, sample, cases = as.integer(cases))]
+    ), use.names = TRUE)
   }
 
-  if (!all(missing(mean), missing(sd), missing(mean_sd), missing(sd_sd))) {
-    if (sd == 0 && mean_sd == 0 && sd_sd == 0) {
-      distribution <- "fixed"
-    }
-    ## deprecated arguments given
-    if (distribution == "lognormal") {
-      params_mean <- c(meanlog = mean, sdlog = sd)
-      params_sd <- c(meanlog = mean_sd, sdlog = sd_sd)
-    } else if (distribution == "gamma") {
-      temp_dist <- Gamma(
-        mean = Normal(mean, mean_sd),
-        sd = Normal(sd, sd_sd)
+  ## For each sample map to report date
+  report <- future.apply::future_lapply(1:max(infections$sample),
+    function(id) {
+      suppressWarnings(
+        EpiNow2::adjust_infection_to_report(infections[sample == id],
+          delay_defs = delays,
+          type = type,
+          reporting_effect = reporting_effect[sample == id, ]$effect
+        )
       )
-      params_mean <- vapply(get_parameters(temp_dist), mean, numeric(1))
-      params_sd <- vapply(get_parameters(temp_dist), sd_dist, numeric(1))
-    } else if (distribution == "normal") {
-      params_mean <- c(mean = mean, sd = sd)
-      params_sd <- c(mean = mean_sd, sd = sd_sd)
-    } else if (distribution == "fixed") {
-      params_mean <- mean
+    },
+    future.seed = TRUE
+  )
+
+  report <- data.table::rbindlist(report, idcol = "sample")
+
+  out <- list()
+  # bind all samples together
+  out$samples <- report
+  # summarise samples
+  out$summarised <- calc_summary_measures(
+    report[, value := cases][, cases := NULL],
+    summarise_by = "date",
+    order_by = "date",
+    CrIs = CrIs
+  )
+  return(out)
+}
+
+#' Approximate Sampling a Distribution using Counts
+#'
+#' @description `r lifecycle::badge("deprecated")`
+#' Convolves cases by a PMF function. This function will soon be removed or
+#' replaced with a more robust stan implementation.
+#'
+#' @param cases A `<data.frame>` of cases (in date order) with the following
+#' variables: `date` and `cases`.
+#'
+#' @param max_value Numeric, maximum value to allow. Defaults to 120 days
+#'
+#' @param direction Character string, defato "backwards". Direction in which to
+#' map cases. Supports either "backwards" or "forwards".
+#'
+#' @param dist_fn Function that takes two arguments with the first being
+#' numeric and the second being logical (and defined as `dist`). Should return
+#' the probability density or a sample from the defined distribution. See
+#' the examples for more.
+#'
+#' @param earliest_allowed_mapped A character string representing a date
+#' ("2020-01-01"). Indicates the earliest allowed mapped value.
+#'
+#' @param type Character string indicating the method to use to transform
+#' counts. Supports either "sample"  which approximates sampling or "median"
+#' would shift by the median of the distribution.
+#'
+#' @param truncate_future Logical, should cases be truncated if they occur
+#' after the first date reported in the data. Defaults to `TRUE`.
+#'
+#' @return A `<data.table>` of cases by date of onset
+#' @keywords internal
+#' @export
+#' @importFrom purrr map_dfc
+#' @importFrom data.table data.table setorder
+#' @importFrom lubridate days
+#' @examples
+#' \donttest{
+#' cases <- example_confirmed
+#' cases <- cases[, cases := as.integer(confirm)]
+#' print(cases)
+#'
+#' # total cases
+#' sum(cases$cases)
+#'
+#' delay_fn <- function(n, dist, cum) {
+#'   if (dist) {
+#'     pgamma(n + 0.9999, 2, 1) - pgamma(n - 1e-5, 2, 1)
+#'   } else {
+#'     as.integer(rgamma(n, 2, 1))
+#'   }
+#' }
+#'
+#' onsets <- sample_approx_dist(
+#'   cases = cases,
+#'   dist_fn = delay_fn
+#' )
+#'
+#' # estimated onset distribution
+#' print(onsets)
+#'
+#' # check that sum is equal to reported cases
+#' total_onsets <- median(
+#'   purrr::map_dbl(
+#'     1:10,
+#'     ~ sum(sample_approx_dist(
+#'       cases = cases,
+#'       dist_fn = delay_fn
+#'     )$cases)
+#'   )
+#' )
+#' total_onsets
+#'
+#'
+#' # map from onset cases to reported
+#' reports <- sample_approx_dist(
+#'   cases = cases,
+#'   dist_fn = delay_fn,
+#'   direction = "forwards"
+#' )
+#'
+#'
+#' # map from onset cases to reported using a mean shift
+#' reports <- sample_approx_dist(
+#'   cases = cases,
+#'   dist_fn = delay_fn,
+#'   direction = "forwards",
+#'   type = "median"
+#' )
+#' }
+sample_approx_dist <- function(cases = NULL,
+                               dist_fn = NULL,
+                               max_value = 120,
+                               earliest_allowed_mapped = NULL,
+                               direction = "backwards",
+                               type = "sample",
+                               truncate_future = TRUE) {
+  deprecate_warn(
+    "1.5.0",
+    "sample_approx_dist()",
+    details = "The function will be removed completely in the next version."
+  )
+  if (type == "sample") {
+    if (direction == "backwards") {
+      direction_fn <- rev
+    } else if (direction == "forwards") {
+      direction_fn <- function(x) {
+        x
+      }
     }
-  }
-  if (length(pmf) > 0) {
-    if (!all(
-      missing(mean), missing(sd), missing(mean_sd), missing(sd_sd),
-      missing(params_mean), missing(params_sd)
-    )) {
-      stop("Distributional parameters should not be given in addition to `pmf`")
+    # reverse cases so starts with current first
+    reversed_cases <- direction_fn(cases$cases)
+    reversed_cases[is.na(reversed_cases)] <- 0
+    # draw from the density fn of the dist
+    draw <- dist_fn(0:max_value, dist = TRUE, cum = FALSE)
+
+    # approximate cases
+    mapped_cases <- do.call(cbind, purrr::map(
+      seq_along(reversed_cases),
+      ~ c(
+        rep(0, . - 1),
+        stats::rbinom(
+          length(draw),
+          rep(reversed_cases[.], length(draw)),
+          draw
+        ),
+        rep(0, length(reversed_cases) - .)
+      )
+    ))
+
+
+    # set dates order based on direction mapping
+    if (direction == "backwards") {
+      dates <- seq(min(cases$date) - lubridate::days(length(draw) - 1),
+        max(cases$date),
+        by = "days"
+      )
+    } else if (direction == "forwards") {
+      dates <- seq(min(cases$date),
+        max(cases$date) + lubridate::days(length(draw) - 1),
+        by = "days"
+      )
     }
-    distribution <- "nonparametric"
-    parameters <- list(pmf = pmf)
-  } else {
-    if (length(params_sd) == 0) {
-      params_sd <- rep(0, length(params_mean))
+
+    # summarises movements and sample for placement of non-integer cases
+    case_sum <- direction_fn(rowSums(mapped_cases))
+    floor_case_sum <- floor(case_sum)
+    sample_cases <- floor_case_sum +
+      as.numeric((runif(seq_along(case_sum)) < (case_sum - floor_case_sum)))
+
+    # summarise imputed onsets and build output data.table
+    mapped_cases <- data.table::data.table(
+      date = dates,
+      cases = sample_cases
+    )
+
+    # filter out all zero cases until first recorded case
+    mapped_cases <- data.table::setorder(mapped_cases, date)
+    mapped_cases <- mapped_cases[
+      ,
+      cum_cases := cumsum(cases)
+    ][cum_cases != 0][, cum_cases := NULL]
+  } else if (type == "median") {
+    shift <- as.integer(
+      median(as.integer(dist_fn(1000, dist = FALSE)), na.rm = TRUE)
+    )
+
+    if (direction == "backwards") {
+      mapped_cases <- data.table::copy(cases)[
+        ,
+        date := date - lubridate::days(shift)
+      ]
+    } else if (direction == "forwards") {
+      mapped_cases <- data.table::copy(cases)[
+        ,
+        date := date + lubridate::days(shift)
+      ]
     }
-    parameters <- lapply(seq_along(params_mean), function(id) {
-      Normal(params_mean[id], params_sd[id])
-    })
-    names(parameters) <- natural_params(distribution)
-    parameters$max <- max
   }
-  return(new_dist_spec(parameters, distribution))
+
+  if (!is.null(earliest_allowed_mapped)) {
+    mapped_cases <- mapped_cases[date >= as.Date(earliest_allowed_mapped)]
+  }
+
+  # filter out future cases
+  if (direction == "forwards" && truncate_future) {
+    max_date <- max(cases$date)
+    mapped_cases <- mapped_cases[date <= max_date]
+  }
+  return(mapped_cases)
 }
diff --git a/R/dist_spec.R b/R/dist_spec.R
index 6efd2f782..08513d18f 100644
--- a/R/dist_spec.R
+++ b/R/dist_spec.R
@@ -183,8 +183,8 @@ dist_skel <- function(n, dist = FALSE, cum = TRUE, model,
 
 #' Creates a delay distribution as the sum of two other delay distributions.
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' @return A delay distribution representing the sum of the two delays
-
 #' @param e1 The first delay distribution (of type [dist_spec()]) to
 #' combine.
 #'
@@ -210,6 +210,7 @@ dist_skel <- function(n, dist = FALSE, cum = TRUE, model,
 
 #' Combines multiple delay distributions for further processing
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This combines the parameters so that they can be fed as multiple delay
 #' distributions to [epinow()] or [estimate_infections()].
 #'
@@ -245,6 +246,7 @@ c.dist_spec <- function(...) {
 
 #' Returns the mean of one or more delay distribution
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This works out the mean of all the (parametric / nonparametric) delay
 #' distributions combined in the passed [dist_spec()] (ignoring any uncertainty
 #' in parameters)
@@ -308,6 +310,7 @@ mean.dist_spec <- function(x, ..., ignore_uncertainty = FALSE) {
 
 #' Returns the standard deviation of one or more delay distribution
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This works out the standard deviation of all the (parametric /
 #' nonparametric) delay distributions combined in the passed [dist_spec()].
 #'
@@ -365,6 +368,7 @@ sd_dist <- function(x) {
 
 #' Returns the maximum of one or more delay distribution
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This works out the maximum of all the (parametric / nonparametric) delay
 #' distributions combined in the passed [dist_spec()] (ignoring any uncertainty
 #' in parameters)
@@ -402,6 +406,7 @@ max.dist_spec <- function(x, ...) {
 
 #' Discretise a <dist_spec>
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' By default it will discretise all the distributions it can discretise
 #' (i.e. those with finite support and constant parameters).
 #' @title Discretise a <dist_spec>
@@ -465,6 +470,7 @@ discretize <- discretise
 
 #' Collapse nonparametric distributions in a <dist_spec>
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This convolves any consecutive nonparametric distributions contained
 #' in the <dist_spec>.
 #' @param x A `<dist_spec>`
@@ -513,6 +519,7 @@ collapse <- function(x) {
 
 #' Applies a threshold to all nonparametric distributions in a <dist_spec>
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This removes any part of the tail of the nonparametric distributions in the
 #' <dist_spec> where the probability mass is below the threshold level.
 #' @param x A `<dist_spec>`
@@ -545,6 +552,7 @@ apply_tolerance <- function(x, tolerance) {
 
 #' Prints the parameters of one or more delay distributions
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This displays the parameters of the uncertain and probability mass
 #' functions of fixed delay distributions combined in the passed [dist_spec()].
 #' @param x The `<dist_spec>` to use
@@ -565,7 +573,7 @@ apply_tolerance <- function(x, tolerance) {
 print.dist_spec <- function(x, ...) {
   .print.dist_spec(x, indent = 0, ...)
 }
-
+#' @keywords internal
 .print.dist_spec <- function(x, indent, ...) {
   indent_str <- strrep(" ", indent)
   if (length(x) > 1) {
@@ -615,6 +623,7 @@ print.dist_spec <- function(x, ...) {
 
 #' Plot PMF and CDF for a dist_spec object
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This function takes a `<dist_spec>` object and plots its probability mass
 #' function (PMF) and cumulative distribution function (CDF) using `{ggplot2}`.
 #' Note that currently uncertainty in distributions is not plot.
@@ -694,6 +703,7 @@ plot.dist_spec <- function(x, ...) {
 
 #' Extract a single element of a composite `<dist_spec>`
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' @param x A composite `dist_spec` object
 #' @param i The index to extract
 #' @return A single `dist_spec` object
@@ -722,6 +732,7 @@ extract_single_dist <- function(x, i) {
 
 #' Fix the parameters of a `<dist_spec>`
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' If the given `<dist_spec>` has any uncertainty, it is removed and the
 #' corresponding distribution converted into a fixed one.
 #' @return A `<dist_spec>` object without uncertainty
@@ -877,6 +888,7 @@ NonParametric <- function(pmf) {
 
 #' Get the names of the natural parameters of a distribution
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' These are the parameters used in the stan models. All other parameter
 #' representations are converted to these using [convert_to_natural()] before
 #' being passed to the stan models.
@@ -902,6 +914,7 @@ natural_params <- function(distribution) {
 
 #' Get the lower bounds of the parameters of a distribution
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This is used to avoid sampling parameter values that have no support.
 #' @return A numeric vector, the lower bounds.
 #' @inheritParams natural_params
@@ -923,6 +936,8 @@ lower_bounds <- function(distribution) {
   return(ret)
 }
 
+#' Extract parameter names
+#' @description `r lifecycle::badge("experimental")`
 #' Internal function for extracting given parameter names of a distribution
 #' from the environment. Called by `new_dist_spec`
 #'
@@ -945,6 +960,7 @@ extract_params <- function(params, distribution) {
 #' Internal function for generating a `dist_spec` given parameters and a
 #' distribution.
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This will convert all parameters to natural parameters before generating
 #' a `dist_spec`. If they have uncertainty this will be done using sampling.
 #' @param params Parameters of the distribution (including `max`)
@@ -1039,6 +1055,7 @@ new_dist_spec <- function(params, distribution) {
 
 #' Internal function for converting parameters to natural parameters.
 #'
+#' @description `r lifecycle::badge("experimental")`
 #' This is used for preprocessing before generating a `dist_spec` object
 #' from a given set of parameters and distribution
 #' @param params A numerical named parameter vector
diff --git a/R/epinow-internal.R b/R/epinow-internal.R
index a2bee0427..ad25210d7 100644
--- a/R/epinow-internal.R
+++ b/R/epinow-internal.R
@@ -6,7 +6,7 @@
 #' @inheritParams setup_target_folder
 #' @inheritParams estimate_infections
 #' @return Numeric forecast horizon adjusted for the users intention
-#' @export
+#' @keywords internal
 update_horizon <- function(horizon, target_date, data) {
   if (horizon != 0) {
     horizon <- horizon + as.numeric(
@@ -24,7 +24,7 @@ update_horizon <- function(horizon, target_date, data) {
 #' @inheritParams setup_target_folder
 #' @inheritParams epinow
 #' @return No return value, called for side effects
-#' @export
+#' @keywords internal
 save_input <- function(data, target_folder) {
   if (!is.null(target_folder)) {
     latest_date <- data[confirm > 0][date == max(date)]$date
@@ -35,8 +35,6 @@ save_input <- function(data, target_folder) {
   return(invisible(NULL))
 }
 
-
-
 #' Save Estimated Infections
 #'
 #' @description `r lifecycle::badge("stable")`
@@ -52,7 +50,7 @@ save_input <- function(data, target_folder) {
 #' @inheritParams setup_target_folder
 #' @inheritParams  estimate_infections
 #' @return No return value, called for side effects
-#' @export
+#' @keywords internal
 save_estimate_infections <- function(estimates, target_folder = NULL,
                                      samples = TRUE, return_fit = TRUE) {
   if (!is.null(target_folder)) {
@@ -85,7 +83,7 @@ save_estimate_infections <- function(estimates, target_folder = NULL,
 #'
 #' @return A list of samples and summarised estimates of estimated cases by
 #' date of report.
-#' @export
+#' @keywords internal
 #' @importFrom data.table := rbindlist
 estimates_by_report_date <- function(estimates, CrIs = c(0.2, 0.5, 0.9),
                                      target_folder = NULL, samples = TRUE) {
@@ -130,7 +128,7 @@ estimates_by_report_date <- function(estimates, CrIs = c(0.2, 0.5, 0.9),
 #' @inheritParams setup_target_folder
 #'
 #' @return No return value, called for side effects
-#' @export
+#' @keywords internal
 copy_results_to_latest <- function(target_folder = NULL, latest_folder = NULL) {
   if (!is.null(target_folder)) {
     ## Save all results to a latest folder as well
@@ -170,7 +168,7 @@ copy_results_to_latest <- function(target_folder = NULL, latest_folder = NULL) {
 #' @inheritParams save_estimate_infections
 #'
 #' @return A list of output as returned by `epinow`
-#' @export
+#' @keywords internal
 construct_output <- function(estimates,
                              estimated_reported_cases,
                              plots = NULL,
diff --git a/R/epinow.R b/R/epinow.R
index 075f7638d..906c1f803 100644
--- a/R/epinow.R
+++ b/R/epinow.R
@@ -1,6 +1,6 @@
 #' Real-time Rt Estimation, Forecasting and Reporting
 #'
-#' @description `r lifecycle::badge("maturing")`
+#' @description `r lifecycle::badge("stable")`
 #' This function wraps the functionality of [estimate_infections()] in order
 #' to estimate Rt and cases by date of infection and forecast these infections
 #' into the future. In addition to the functionality of
@@ -109,7 +109,7 @@ epinow <- function(data,
       "1.5.0",
       "epinow(reported_cases)",
       "epinow(data)",
-      "The argument will be removed completely in version 2.0.0."
+      "The argument will be removed completely in the next version."
     )
     data <- reported_cases
   }
@@ -313,4 +313,4 @@ epinow <- function(data,
     return(invisible(NULL))
   }
 }
-# nolint end: cyclocomp_linter
\ No newline at end of file
+# nolint end: cyclocomp_linter
diff --git a/R/estimate_infections.R b/R/estimate_infections.R
index ecef78dab..b5e9f6f6d 100644
--- a/R/estimate_infections.R
+++ b/R/estimate_infections.R
@@ -138,7 +138,7 @@ estimate_infections <- function(data,
       "1.5.0",
       "estimate_infections(reported_cases)",
       "estimate_infections(data)",
-      "The argument will be removed completely in version 2.0.0."
+      "The argument will be removed completely in the next version."
     )
     data <- reported_cases
   }
@@ -290,328 +290,6 @@ estimate_infections <- function(data,
   return(format_out)
 }
 
-
-#' Generate initial conditions by fitting to cumulative cases
-#'
-#' @description `r lifecycle::badge("deprecated")`
-#'
-#' Fits a model to cumulative cases. This may be a useful approach to
-#' initialising a full model fit for certain data sets where the sampler gets
-#' stuck or cannot easily be initialised as fitting to cumulative cases changes
-#' the shape of the posterior distribution. In `estimate_infections()`,
-#' `epinow()` and `regional_epinow()` this option can be engaged by setting
-#' `stan_opts(init_fit = "cumulative")`.
-#'
-#' This implementation is based on the approach taken in
-#' [epidemia](https://github.com/ImperialCollegeLondon/epidemia/) authored by
-#' James Scott.
-#'
-#' @param samples Numeric, defaults to 50. Number of posterior samples.
-#'
-#' @param warmup Numeric, defaults to 50. Number of warmup samples.
-#'
-#' @param verbose Logical, should fitting progress be returned. Defaults to
-#' `FALSE`.
-#'
-#' @inheritParams create_initial_conditions
-#' @importFrom rstan sampling
-#' @importFrom futile.logger flog.debug
-#' @importFrom utils capture.output
-#' @inheritParams fit_model_with_nuts
-#' @inheritParams stan_opts
-#' @return A stanfit object
-init_cumulative_fit <- function(args, samples = 50, warmup = 50,
-                                id = "init", verbose = FALSE,
-                                backend = "rstan") {
-  deprecate_warn(
-    when = "1.5.0",
-    what = "init_cumulative_fit()"
-  )
-  futile.logger::flog.debug(
-    "%s: Fitting to cumulative data to initialise chains", id,
-    name = "EpiNow2.epinow.estimate_infections.fit"
-  )
-  # copy main run settings and override to use only 100 iterations and a single
-  # chain
-  initial_args <- create_stan_args(
-    stan = stan_opts(
-      args$object,
-      samples = samples,
-      warmup = warmup,
-      control = list(adapt_delta = 0.9, max_treedepth = 13),
-      chains = 1,
-      cores = 2,
-      backend = backend,
-      open_progress = FALSE,
-      show_messages = FALSE
-    ),
-    data = args$data, init = args$init
-  )
-  # change observations to be cumulative in order to protect against noise and
-  # give an approximate fit (though for Rt constrained to be > 1)
-  initial_args$data$cases <- cumsum(initial_args$data$cases)
-  initial_args$data$shifted_cases <- cumsum(initial_args$data$shifted_cases)
-
-  # initial fit
-  if (verbose) {
-    fit <- fit_model(initial_args, id = "init_cumulative")
-  } else {
-    out <- tempfile(tmpdir = tempdir(check = TRUE))
-    capture.output(
-      {
-        fit <- fit_model(initial_args, id = "init_cumulative")
-      },
-      type = c("output", "message"),
-      split = FALSE,
-      file = out
-    )
-  }
-  return(fit)
-}
-
-#' Fit a Stan Model using the NUTs sampler
-#'
-#' @description `r lifecycle::badge("maturing")`
-#' Fits a stan model using [rstan::sampling()]. Provides the optional ability to
-#' run chains using `future` with error catching, timeouts and merging of
-#' completed chains.
-#'
-#' @param args List of stan arguments.
-#'
-#' @param future Logical, defaults to `FALSE`. Should `future` be used to run
-#' stan chains in parallel.
-#'
-#' @param max_execution_time Numeric, defaults to Inf. What is the maximum
-#' execution time per chain in seconds. Results will still be returned as long
-#' as at least 2 chains complete successfully within the timelimit.
-#'
-#' @param id A character string used to assign logging information on error.
-#' Used by [regional_epinow()] to assign errors to regions. Alter the default to
-#' run with error catching.
-#'
-#' @importFrom futile.logger flog.debug flog.info flog.error
-#' @importFrom R.utils withTimeout
-#' @importFrom future.apply future_lapply
-#' @importFrom purrr compact
-#' @importFrom rstan sflist2stanfit sampling
-#' @importFrom rlang abort cnd_muffle
-#' @return A stan model object
-fit_model_with_nuts <- function(args, future = FALSE, max_execution_time = Inf,
-                                id = "stan") {
-  args$method <- NULL
-  args$max_execution_time <- NULL
-  args$future <- NULL
-
-  futile.logger::flog.debug(
-    paste0(
-      "%s: Running in exact mode for ",
-      ceiling(args$iter - args$warmup) * args$chains,
-      " samples (across ", args$chains,
-      " chains each with a warm up of ", args$warmup, " iterations each) and ",
-      args$data$t, " time steps of which ", args$data$horizon, " are a forecast"
-    ),
-    id,
-    name = "EpiNow2.epinow.estimate_infections.fit"
-  )
-
-  if (exists("stuck_chains", args)) {
-    stuck_chains <- args$stuck_chains
-    args$stuck_chains <- NULL
-  } else {
-    stuck_chains <- 0
-  }
-
-  fit_chain <- function(chain, stan_args, max_time, catch = FALSE) {
-    stan_args$chain_id <- chain
-    if (inherits(stan_args$object, "stanmodel")) {
-      sample_func <- rstan::sampling
-    } else if (inherits(stan_args$object, "CmdStanModel")) {
-      sample_func <- stan_args$object$sample
-      stan_args$object <- NULL
-    }
-    if (catch) {
-      fit <- tryCatch(
-        withCallingHandlers(
-          R.utils::withTimeout(do.call(sample_func, stan_args),
-            timeout = max_time,
-            onTimeout = "silent"
-          ),
-          warning = function(w) {
-            futile.logger::flog.warn(
-              "%s (chain: %s): %s - %s", id, chain, w$message, toString(w$call),
-              name = "EpiNow2.epinow.estimate_infections.fit"
-            )
-            rlang::cnd_muffle(w)
-          }
-        ),
-        error = function(e) {
-          error_text <- sprintf(
-            "%s (chain: %s): %s - %s", id, chain, e$message, toString(e$call)
-          )
-          futile.logger::flog.error(error_text,
-            name = "EpiNow2.epinow.estimate_infections.fit"
-          )
-          return(NULL)
-        }
-      )
-    } else {
-      fit <- R.utils::withTimeout(do.call(sample_func, stan_args),
-        timeout = max_time,
-        onTimeout = "silent"
-      )
-    }
-
-    if ((inherits(fit, "stanfit") && fit@mode != 2L) ||
-        inherits(fit, "CmdStanMCMC")) {
-      return(fit)
-    } else {
-      return(NULL)
-    }
-  }
-
-  if (future) {
-    chains <- args$chains
-    args$chains <- 1
-    args$cores <- 1
-    fits <- future.apply::future_lapply(1:chains,
-      fit_chain,
-      stan_args = args,
-      max_time = max_execution_time,
-      catch = TRUE,
-      future.seed = TRUE
-    )
-    if (stuck_chains > 0) {
-      fits[1:stuck_chains] <- NULL
-    }
-    fit <- purrr::compact(fits)
-    if (length(fit) == 0) {
-      fit <- NULL
-      if (is.null(fit)) {
-        rlang::abort(
-          "all chains failed - try inspecting the output for errors or",
-          " increasing the max_execution_time"
-        )
-      }
-    } else {
-      failed_chains <- chains - length(fit)
-      if (failed_chains > 0) {
-        futile.logger::flog.warn(
-          "%s: %s chains failed or were timed out.", id, failed_chains,
-          name = "EpiNow2.epinow.estimate_infections.fit"
-        )
-        if ((chains - failed_chains) < 2) {
-          rlang::abort(
-            "model fitting failed as too few chains were returned to assess",
-            " convergence (2 or more required)"
-          )
-        }
-      }
-      fit <- rstan::sflist2stanfit(fit)
-    }
-  } else {
-    fit <- fit_chain(seq_len(args$chains),
-        stan_args = args, max_time = max_execution_time,
-        catch = !id %in% c("estimate_infections", "epinow")
-    )
-    if (stuck_chains > 0) {
-      fit <- NULL
-    }
-    if (is.null(fit)) {
-      rlang::abort("model fitting was timed out or failed")
-    }
-  }
-  return(fit)
-}
-
-#' Fit a Stan Model using an approximate method
-#'
-#' @description `r lifecycle::badge("maturing")`
-#' Fits a stan model using variational inference.
-#'
-#' @inheritParams fit_model_with_nuts
-#' @importFrom futile.logger flog.debug flog.info flog.error
-#' @importFrom purrr safely
-#' @importFrom rstan vb
-#' @importFrom rlang abort
-#' @return A stan model object
-fit_model_approximate <- function(args, future = FALSE, id = "stan") {
-  method <- args$method
-  args$method <- NULL
-  futile.logger::flog.debug(
-    paste0(
-      "%s: Running in approximate mode for ", args$iter,
-      " iterations (with ", args$trials, " attempts). Extracting ",
-      args$output_samples, " approximate posterior samples for ",
-      args$data$t, " time steps of which ",
-      args$data$horizon, " are a forecast"
-    ),
-    id,
-    name = "EpiNow2.epinow.estimate_infections.fit"
-  )
-
-  if (exists("trials", args)) {
-    trials <- args$trials
-    args$trials <- NULL
-  } else {
-    trials <- 1
-  }
-
-  fit_approximate <- function(stan_args) {
-    if (inherits(stan_args$object, "stanmodel")) {
-      if (method == "vb") {
-        sample_func <- rstan::vb
-      } else {
-        stop("Laplace approximation only available in the cmdstanr backend")
-      }
-    } else if (inherits(stan_args$object, "CmdStanModel")) {
-      if (method == "vb") {
-        sample_func <- stan_args$object$variational
-      } else if (method == "laplace") {
-        sample_func <- stan_args$object$laplace
-      } else {
-        sample_func <- stan_args$object$pathfinder
-      }
-      stan_args$object <- NULL
-    }
-    fit <- do.call(sample_func, stan_args)
-
-    if (length(names(fit)) == 0) {
-      return(NULL)
-    } else {
-      return(fit)
-    }
-    return(fit)
-  }
-  safe_fit <- purrr::safely(fit_approximate) # nolint
-  fit <- NULL
-  current_trials <- 0
-
-  while (current_trials <= trials && is.null(fit)) {
-    fit <- safe_fit(args)
-
-    error <- fit[[2]]
-    fit <- fit[[1]]
-    if (is(fit, "CmdStanFit") && fit$return_codes() > 0) {
-      error <- tail(capture.output(fit$output()), 1)
-      fit <- NULL
-    }
-    current_trials <- current_trials + 1
-  }
-
-  if (is.null(fit)) {
-    futile.logger::flog.error(
-      paste(
-        "%s: Fitting failed - try increasing stan_args$trials or inspecting",
-        "the model input"
-      ),
-      id,
-      name = "EpiNow2.epinow.estimate_infections.fit"
-    )
-    rlang::abort(paste("Approximate inference failed due to:", error))
-  }
-  return(fit)
-}
-
 #' Format Posterior Samples
 #'
 #' @description `r lifecycle::badge("stable")`
@@ -633,6 +311,7 @@ fit_model_approximate <- function(args, future = FALSE, id = "stan") {
 #' @importFrom lubridate days
 #' @importFrom futile.logger flog.info
 #' @return A list of samples and summarised posterior parameter estimates.
+#' @keywords internal
 format_fit <- function(posterior_samples, horizon, shift, burn_in, start_date,
                        CrIs) {
   format_out <- list()
@@ -657,7 +336,6 @@ format_fit <- function(posterior_samples, horizon, shift, burn_in, start_date,
       )
   ]
 
-
   # remove burn in period if specified
   if (burn_in > 0) {
     futile.logger::flog.info(
diff --git a/R/estimate_secondary.R b/R/estimate_secondary.R
index d215515f8..e3eeb1e5b 100644
--- a/R/estimate_secondary.R
+++ b/R/estimate_secondary.R
@@ -166,7 +166,7 @@ estimate_secondary <- function(data,
       "1.5.0",
       "estimate_secondary(reports)",
       "estimate_secondary(data)",
-      "The argument will be removed completely in version 2.0.0."
+      "The argument will be removed completely in the next version."
     )
     data <- reports
   }
diff --git a/R/estimate_truncation.R b/R/estimate_truncation.R
index f09ff5383..d628de996 100644
--- a/R/estimate_truncation.R
+++ b/R/estimate_truncation.R
@@ -136,7 +136,7 @@ estimate_truncation <- function(data, max_truncation, trunc_max = 10,
       "1.5.0",
       "estimate_truncation(obs)",
       "estimate_truncation(data)",
-      "The argument will be removed completely in version 2.0.0."
+      "The argument will be removed completely in the next version."
     )
     data <- obs
   }
@@ -152,7 +152,7 @@ estimate_truncation <- function(data, max_truncation, trunc_max = 10,
       "1.5.0",
       "estimate_truncation(weigh_delay_priors)",
       "trunc_opts(weight_prior)",
-      detail = "This argument will be removed completely in version 2.0.0"
+      detail = "This argument will be removed completely in the next version"
     )
   }
    # Validate inputs
@@ -165,7 +165,7 @@ estimate_truncation <- function(data, max_truncation, trunc_max = 10,
   assert_logical(weigh_delay_priors)
   assert_logical(verbose)
 
-  ## code block to remove in EpiNow2 2.0.0
+  ## code block to remove in next EpiNow2 version
   construct_trunc <- FALSE
   if (!missing(trunc_max)) {
     if (!missing(truncation)) {
@@ -178,7 +178,7 @@ estimate_truncation <- function(data, max_truncation, trunc_max = 10,
         "`max_truncation` and `trunc_max` arguments are both given. ",
         "Use only `truncation` instead.")
     }
-    deprecate_warn(
+    deprecate_stop(
       "1.4.0",
       "estimate_truncation(trunc_max)",
       "estimate_truncation(truncation)"
@@ -191,7 +191,7 @@ estimate_truncation <- function(data, max_truncation, trunc_max = 10,
         "`max_truncation` and `truncation` arguments are both given. ",
         "Use only `truncation` instead.")
     }
-    deprecate_warn(
+    deprecate_stop(
       "1.4.0",
       "estimate_truncation(max_truncation)",
       "estimate_truncation(truncation)"
@@ -206,7 +206,7 @@ estimate_truncation <- function(data, max_truncation, trunc_max = 10,
         "`trunc_dist` and `truncation` arguments are both given. ",
         "Use only `truncation` instead.")
     }
-    deprecate_warn(
+    deprecate_stop(
       "1.4.0",
       "estimate_truncation(trunc_dist)",
       "estimate_truncation(truncation)"
@@ -250,7 +250,7 @@ estimate_truncation <- function(data, max_truncation, trunc_max = 10,
     confirm := NULL
   ])
   obs <- purrr::reduce(obs, merge, all = TRUE)
-  obs_start <- max(nrow(obs) - trunc_max - sum(is.na(obs$`1`)) + 1, 1)
+  obs_start <- max(nrow(obs) - max(truncation) - sum(is.na(obs$`1`)) + 1, 1)
   obs_dist <- purrr::map_dbl(2:(ncol(obs)), ~ sum(is.na(obs[[.]])))
   obs_data <- obs[, -1][, purrr::map(.SD, ~ ifelse(is.na(.), 0, .))]
   obs_data <- as.matrix(obs_data[obs_start:.N])
@@ -318,7 +318,7 @@ estimate_truncation <- function(data, max_truncation, trunc_max = 10,
     ]
   link_obs <- function(index) {
     target_obs <- dirty_obs[[index]][, index := .N - 0:(.N - 1)]
-    target_obs <- target_obs[index < trunc_max]
+    target_obs <- target_obs[index < max(truncation)]
     estimates <- recon_obs[dataset == index][, c("id", "dataset") := NULL]
     estimates <- estimates[, lapply(.SD, as.integer)]
     estimates <- estimates[, index := .N - 0:(.N - 1)]
diff --git a/R/extract.R b/R/extract.R
index e3ad05c9c..dc6e0d44d 100644
--- a/R/extract.R
+++ b/R/extract.R
@@ -14,6 +14,7 @@
 #' @return A `<data.frame>` containing the parameter name, date, sample id and
 #' sample value.
 #' @importFrom data.table melt as.data.table
+#' @keywords internal
 extract_parameter <- function(param, samples, dates) {
   param_df <- data.table::as.data.table(
     t(
@@ -43,6 +44,7 @@ extract_parameter <- function(param, samples, dates) {
 #' @inheritParams extract_parameter
 #' @return A `<data.frame>` containing the parameter name, sample id and sample
 #' value
+#' @keywords internal
 extract_static_parameter <- function(param, samples) {
   data.table::data.table(
     parameter = param,
@@ -145,6 +147,7 @@ extract_samples <- function(stan_fit, pars = NULL, include = TRUE) {
 #' parameter
 #' @importFrom rstan extract
 #' @importFrom data.table data.table
+#' @keywords internal
 extract_parameter_samples <- function(stan_fit, data, reported_dates,
                                       reported_inf_dates,
                                       drop_length_1 = FALSE, merge = FALSE) {
diff --git a/R/fit.R b/R/fit.R
new file mode 100644
index 000000000..6f2e53257
--- /dev/null
+++ b/R/fit.R
@@ -0,0 +1,244 @@
+#' Fit a Stan Model using the NUTs sampler
+#'
+#' @description `r lifecycle::badge("maturing")`
+#' Fits a stan model using [rstan::sampling()]. Provides the optional ability to
+#' run chains using `future` with error catching, timeouts and merging of
+#' completed chains.
+#'
+#' @param args List of stan arguments.
+#'
+#' @param future Logical, defaults to `FALSE`. Should `future` be used to run
+#' stan chains in parallel.
+#'
+#' @param max_execution_time Numeric, defaults to Inf. What is the maximum
+#' execution time per chain in seconds. Results will still be returned as long
+#' as at least 2 chains complete successfully within the timelimit.
+#'
+#' @param id A character string used to assign logging information on error.
+#' Used by [regional_epinow()] to assign errors to regions. Alter the default to
+#' run with error catching.
+#'
+#' @importFrom futile.logger flog.debug flog.info flog.error
+#' @importFrom R.utils withTimeout
+#' @importFrom future.apply future_lapply
+#' @importFrom purrr compact
+#' @importFrom rstan sflist2stanfit sampling
+#' @importFrom rlang abort cnd_muffle
+#' @return A stan model object
+#' @keywords internal
+fit_model_with_nuts <- function(args, future = FALSE, max_execution_time = Inf,
+                                id = "stan") {
+  args$method <- NULL
+  args$max_execution_time <- NULL
+  args$future <- NULL
+
+  futile.logger::flog.debug(
+    paste0(
+      "%s: Running in exact mode for ",
+      ceiling(args$iter - args$warmup) * args$chains,
+      " samples (across ", args$chains,
+      " chains each with a warm up of ", args$warmup, " iterations each) and ",
+      args$data$t, " time steps of which ", args$data$horizon, " are a forecast"
+    ),
+    id,
+    name = "EpiNow2.epinow.estimate_infections.fit"
+  )
+
+  if (exists("stuck_chains", args)) {
+    stuck_chains <- args$stuck_chains
+    args$stuck_chains <- NULL
+  } else {
+    stuck_chains <- 0
+  }
+
+  fit_chain <- function(chain, stan_args, max_time, catch = FALSE) {
+    stan_args$chain_id <- chain
+    if (inherits(stan_args$object, "stanmodel")) {
+      sample_func <- rstan::sampling
+    } else if (inherits(stan_args$object, "CmdStanModel")) {
+      sample_func <- stan_args$object$sample
+      stan_args$object <- NULL
+    }
+    if (catch) {
+      fit <- tryCatch(
+        withCallingHandlers(
+          R.utils::withTimeout(do.call(sample_func, stan_args),
+            timeout = max_time,
+            onTimeout = "silent"
+          ),
+          warning = function(w) {
+            futile.logger::flog.warn(
+              "%s (chain: %s): %s - %s", id, chain, w$message, toString(w$call),
+              name = "EpiNow2.epinow.estimate_infections.fit"
+            )
+            rlang::cnd_muffle(w)
+          }
+        ),
+        error = function(e) {
+          error_text <- sprintf(
+            "%s (chain: %s): %s - %s", id, chain, e$message, toString(e$call)
+          )
+          futile.logger::flog.error(error_text,
+            name = "EpiNow2.epinow.estimate_infections.fit"
+          )
+          return(NULL)
+        }
+      )
+    } else {
+      fit <- R.utils::withTimeout(do.call(sample_func, stan_args),
+        timeout = max_time,
+        onTimeout = "silent"
+      )
+    }
+
+    if ((inherits(fit, "stanfit") && fit@mode != 2L) ||
+        inherits(fit, "CmdStanMCMC")) {
+      return(fit)
+    } else {
+      return(NULL)
+    }
+  }
+
+  if (future) {
+    chains <- args$chains
+    args$chains <- 1
+    args$cores <- 1
+    fits <- future.apply::future_lapply(1:chains,
+      fit_chain,
+      stan_args = args,
+      max_time = max_execution_time,
+      catch = TRUE,
+      future.seed = TRUE
+    )
+    if (stuck_chains > 0) {
+      fits[1:stuck_chains] <- NULL
+    }
+    fit <- purrr::compact(fits)
+    if (length(fit) == 0) {
+      fit <- NULL
+      if (is.null(fit)) {
+        rlang::abort(
+          "all chains failed - try inspecting the output for errors or",
+          " increasing the max_execution_time"
+        )
+      }
+    } else {
+      failed_chains <- chains - length(fit)
+      if (failed_chains > 0) {
+        futile.logger::flog.warn(
+          "%s: %s chains failed or were timed out.", id, failed_chains,
+          name = "EpiNow2.fit"
+        )
+        if ((chains - failed_chains) < 2) {
+          rlang::abort(
+            "model fitting failed as too few chains were returned to assess",
+            " convergence (2 or more required)"
+          )
+        }
+      }
+      fit <- rstan::sflist2stanfit(fit)
+    }
+  } else {
+    fit <- fit_chain(seq_len(args$chains),
+        stan_args = args, max_time = max_execution_time,
+        catch = !id %in% c("estimate_infections", "epinow")
+    )
+    if (stuck_chains > 0) {
+      fit <- NULL
+    }
+    if (is.null(fit)) {
+      rlang::abort("model fitting was timed out or failed")
+    }
+  }
+  return(fit)
+}
+
+#' Fit a Stan Model using an approximate method
+#'
+#' @description `r lifecycle::badge("maturing")`
+#' Fits a stan model using variational inference.
+#'
+#' @inheritParams fit_model_with_nuts
+#' @importFrom futile.logger flog.debug flog.info flog.error
+#' @importFrom purrr safely
+#' @importFrom rstan vb
+#' @importFrom rlang abort
+#' @return A stan model object
+#' @keywords internal
+fit_model_approximate <- function(args, future = FALSE, id = "stan") {
+  method <- args$method
+  args$method <- NULL
+  futile.logger::flog.debug(
+    paste0(
+      "%s: Running in approximate mode for ", args$iter,
+      " iterations (with ", args$trials, " attempts). Extracting ",
+      args$output_samples, " approximate posterior samples for ",
+      args$data$t, " time steps of which ",
+      args$data$horizon, " are a forecast"
+    ),
+    id,
+    name = "EpiNow2.epinow.estimate_infections.fit"
+  )
+
+  if (exists("trials", args)) {
+    trials <- args$trials
+    args$trials <- NULL
+  } else {
+    trials <- 1
+  }
+
+  fit_approximate <- function(stan_args) {
+    if (inherits(stan_args$object, "stanmodel")) {
+      if (method == "vb") {
+        sample_func <- rstan::vb
+      } else {
+        stop("Laplace approximation only available in the cmdstanr backend")
+      }
+    } else if (inherits(stan_args$object, "CmdStanModel")) {
+      if (method == "vb") {
+        sample_func <- stan_args$object$variational
+      } else if (method == "laplace") {
+        sample_func <- stan_args$object$laplace
+      } else {
+        sample_func <- stan_args$object$pathfinder
+      }
+      stan_args$object <- NULL
+    }
+    fit <- do.call(sample_func, stan_args)
+
+    if (length(names(fit)) == 0) {
+      return(NULL)
+    } else {
+      return(fit)
+    }
+    return(fit)
+  }
+  safe_fit <- purrr::safely(fit_approximate) # nolint
+  fit <- NULL
+  current_trials <- 0
+
+  while (current_trials <= trials && is.null(fit)) {
+    fit <- safe_fit(args)
+
+    error <- fit[[2]]
+    fit <- fit[[1]]
+    if (is(fit, "CmdStanFit") && fit$return_codes() > 0) {
+      error <- tail(capture.output(fit$output()), 1)
+      fit <- NULL
+    }
+    current_trials <- current_trials + 1
+  }
+
+  if (is.null(fit)) {
+    futile.logger::flog.error(
+      paste(
+        "%s: Fitting failed - try increasing stan_args$trials or inspecting",
+        "the model input"
+      ),
+      id,
+      name = "EpiNow2.fit"
+    )
+    rlang::abort(paste("Approximate inference failed due to:", error))
+  }
+  return(fit)
+}
diff --git a/R/get.R b/R/get.R
index 891501fd5..913a1bfea 100644
--- a/R/get.R
+++ b/R/get.R
@@ -6,7 +6,7 @@
 #'  are stored (as produced by [regional_epinow()]).
 #'
 #' @return A named character vector containing the results to plot.
-#' @export
+#' @keywords internal
 get_regions <- function(results_dir) {
   # regions to include - based on folder names
   regions <- list.dirs(results_dir,
@@ -20,6 +20,7 @@ get_regions <- function(results_dir) {
   names(regions) <- regions
   return(regions)
 }
+
 #' Get a Single Raw Result
 #'
 #' @description `r lifecycle::badge("stable")`
@@ -34,7 +35,7 @@ get_regions <- function(results_dir) {
 #' directory.
 #'
 #' @return An R object read in from the targeted `.rds` file
-#' @export
+#' @keywords internal
 get_raw_result <- function(file, region, date,
                            result_dir) {
   file_path <- file.path(result_dir, region, date, file)
@@ -90,7 +91,7 @@ get_regional_results <- function(regional_output,
     # find all regions
     regions <- get_regions(results_dir)
 
-    load_data <- purrr::safely(EpiNow2::get_raw_result) # nolint
+    load_data <- purrr::safely(get_raw_result) # nolint
 
     # get estimates
     get_estimates_file <- function(samples_path, summarised_path) {
@@ -152,6 +153,7 @@ get_regional_results <- function(regional_output,
   }
   return(out)
 }
+
 #' Get a Literature Distribution
 #'
 #'
@@ -175,6 +177,7 @@ get_regional_results <- function(regional_output,
 #'
 #' @seealso [dist_spec()]
 #' @export
+#' @keywords internal
 get_dist <- function(data, disease, source, max_value = 14, fixed = FALSE) {
   lifecycle::deprecate_warn(
     "1.5.0", "get_dist()",
@@ -183,7 +186,7 @@ get_dist <- function(data, disease, source, max_value = 14, fixed = FALSE) {
         "Please use distribution functions such as `Gamma()` or `Lognormal()`",
         "instead."
       ),
-      "The function will be removed completely in version 2.0.0."
+      "The function will be removed completely in the next version."
     )
   )
   target_disease <- disease
@@ -220,6 +223,7 @@ get_dist <- function(data, disease, source, max_value = 14, fixed = FALSE) {
 #' @inherit get_dist
 #' @export
 #' @seealso [dist_spec()]
+#' @keywords internal
 get_generation_time <- function(disease, source, max_value = 14,
                                 fixed = FALSE) {
   lifecycle::deprecate_warn(
@@ -229,7 +233,7 @@ get_generation_time <- function(disease, source, max_value = 14,
         "Please use distribution functions such as `Gamma()` or `Lognormal()`",
         "instead."
       ),
-      "The function will be removed completely in version 2.0.0.",
+      "The function will be removed completely in the next version.",
       paste(
         "To obtain the previous estimate by Ganyani et al. (2020) use",
         "`example_generation_time`."
@@ -254,6 +258,7 @@ get_generation_time <- function(disease, source, max_value = 14,
 #' @inheritParams get_dist
 #' @inherit get_dist
 #' @export
+#' @keywords internal
 get_incubation_period <- function(disease, source, max_value = 14,
                                   fixed = FALSE) {
   lifecycle::deprecate_warn(
@@ -263,7 +268,7 @@ get_incubation_period <- function(disease, source, max_value = 14,
         "Please use distribution functions such as `Gamma()` or `Lognormal()`",
         "instead."
       ),
-      "The function will be removed completely in version 2.0.0.",
+      "The function will be removed completely in the next version.",
      paste(
       "To obtain the previous estimate by Ganyani et al. (2020) use",
       "`example_incubation_period`."
@@ -294,7 +299,7 @@ get_incubation_period <- function(disease, source, max_value = 14,
 #'
 #' @importFrom data.table copy setorderv
 #' @importFrom lubridate days
-#' @export
+#' @keywords internal
 get_regions_with_most_reports <- function(data,
                                           time_window = 7,
                                           no_regions = 6) {
@@ -321,6 +326,7 @@ get_regions_with_most_reports <- function(data,
 ##' to be at least the maximum generation time
 ##' @inheritParams estimate_infections
 ##' @return An integer seeding time
+##' @keywords internal
 get_seeding_time <- function(delays, generation_time, rt = rt_opts()) {
   # Estimate the mean delay -----------------------------------------------
   seeding_time <- sum(mean(delays, ignore_uncertainty = TRUE))
diff --git a/R/opts.R b/R/opts.R
index 3ea60fed0..3d212614d 100644
--- a/R/opts.R
+++ b/R/opts.R
@@ -11,8 +11,6 @@
 #' @param source deprecated; use `dist` instead
 #' @param max deprecated; use `dist` instead
 #' @param fixed deprecated; use `dist` instead
-#' @param prior_weight deprecated; prior weights are now specified as a
-#' model option. Use the `weight_prior` argument instead
 #' @param weight_prior Logical; if TRUE (default), any priors given in `dist`
 #'   will be weighted by the number of observation data points, in doing so
 #'   approximately placing an independent prior at each time step and usually
@@ -45,8 +43,7 @@
 #' generation_time_opts(example_generation_time)
 generation_time_opts <- function(dist = Fixed(1), ...,
                                  disease, source, max = 14, fixed = FALSE,
-                                 prior_weight, tolerance = 0.001,
-                                 weight_prior = TRUE) {
+                                 tolerance = 0.001, weight_prior = TRUE) {
   deprecated_options_given <- FALSE
   dot_options <- list(...)
 
@@ -85,13 +82,6 @@ generation_time_opts <- function(dist = Fixed(1), ...,
     }
     deprecated_options_given <- TRUE
   }
-  if (!missing(prior_weight)) {
-    deprecate_warn(
-      "1.4.0", "generation_time_opts(prior_weight)",
-      "generation_time_opts(weight_prior)",
-      "This argument will be removed in version 2.0.0."
-    )
-  }
   if (deprecated_options_given) {
     warning(
       "The generation time distribution should be given to ",
@@ -99,7 +89,7 @@ generation_time_opts <- function(dist = Fixed(1), ...,
       "This behaviour has changed from previous versions of `EpiNow2` and ",
       "any code using it may need to be updated as any other ways of ",
       "specifying the generation time are deprecated and will be removed in ",
-      "version 2.0.0. For examples and more ",
+      "the next version. For examples and more ",
       "information, see the relevant documentation pages using ",
       "`?generation_time_opts`")
   }
@@ -239,7 +229,7 @@ delay_opts <- function(dist = Fixed(0), ..., fixed = FALSE, tolerance = 0.001,
       "This behaviour has changed from previous versions of `EpiNow2` and ",
       "any code using it may need to be updated as any other ways of ",
       "specifying delays are deprecated and will be removed in ",
-      "version 2.0.0. For examples and more ",
+      "the next version. For examples and more ",
       "information, see the relevant documentation pages using ",
       "`?delay_opts`."
     )
@@ -298,7 +288,7 @@ trunc_opts <- function(dist = Fixed(0), tolerance = 0.001,
       "This behaviour has changed from previous versions of `EpiNow2` and ",
       "any code using it may need to be updated as any other ways of ",
       "specifying delays are deprecated and will be removed in ",
-      "version 2.0.0. For examples and more ",
+      "the next version. For examples and more ",
       "information, see the relevant documentation pages using ",
       "`?trunc_opts`"
     )
@@ -935,7 +925,7 @@ rstan_opts <- function(object = NULL,
 #' James Scott.
 #'
 #' This argument is deprecated and the default (NULL) will be used from
-#' version 2.0.0.
+#' the next version.
 #'
 #' @param return_fit Logical, defaults to TRUE. Should the fit stan model be
 #' returned.
@@ -1015,7 +1005,7 @@ stan_opts <- function(object = NULL,
       when = "1.5.0",
       what = "stan_opts(init_fit)",
       details = paste("This argument is deprecated and the default (NULL)",
-                      "will be used from version 2.0.0.")
+                      "will be used from the next version.")
     )
     if (is.character(init_fit)) {
       init_fit <- arg_match(init_fit, values = "cumulative")
diff --git a/R/regional_epinow.R b/R/regional_epinow.R
index ad8642b50..070e75c32 100644
--- a/R/regional_epinow.R
+++ b/R/regional_epinow.R
@@ -124,7 +124,7 @@ regional_epinow <- function(data,
       "1.5.0",
       "regional_epinow(reported_cases)",
       "regional_epinow(data)",
-      "The argument will be removed completely in version 2.0.0."
+      "The argument will be removed completely in the next version."
     )
     data <- reported_cases
   }
@@ -418,6 +418,7 @@ run_region <- function(target_region,
 #' @seealso [regional_epinow()]
 #' @importFrom futile.logger flog.info
 #' @return A list of processed output
+#' @keywords internal
 process_region <- function(out, target_region, timing,
                            return_output = TRUE, return_timing = TRUE,
                            complete_logger = "EpiNow2.epinow") {
@@ -453,6 +454,7 @@ process_region <- function(out, target_region, timing,
 #' @importFrom futile.logger flog.trace flog.info
 #' @seealso [regional_epinow()] [epinow()]
 #' @return A list of all regional estimates and successful regional estimates
+#' @keywords internal
 process_regions <- function(regional_out, regions) {
   # names on regional_out
   names(regional_out) <- regions
diff --git a/R/report.R b/R/report.R
index d3b3ce1bb..797e227b9 100644
--- a/R/report.R
+++ b/R/report.R
@@ -1,132 +1,3 @@
-#' Report case counts by date of report
-#'
-#' @description `r lifecycle::badge("deprecated")`
-#' Convolves latent infections to reported cases via an observation model.
-#' Likely to be removed/replaced in later releases by functionality drawing on
-#' the `stan` implementation.
-#'
-#' @param case_estimates A data.table of case estimates with the following
-#' variables: date, sample, cases
-#'
-#' @param case_forecast A data.table of case forecasts with the following
-#' variables: date, sample, cases. If not supplied the default is not to
-#' incorporate forecasts.
-#'
-#' @param reporting_effect A `data.table` giving the weekly reporting effect
-#'  with the following variables: `sample` (must be the same as in `nowcast`),
-#'  `effect` (numeric scaling factor for each weekday),`day` (numeric 1 - 7
-#'  (1 = Monday and 7 = Sunday)). If not supplied then no weekly reporting
-#'  effect is assumed.
-#'
-#' @return A list of `data.table`s. The first entry contains the following
-#' variables `sample`, `date` and `cases` with the second being summarised
-#' across samples.
-#'
-#' @export
-#' @inheritParams estimate_infections
-#' @inheritParams adjust_infection_to_report
-#' @importFrom data.table data.table rbindlist
-#' @importFrom future.apply future_lapply
-#' @importFrom lifecycle deprecate_warn
-#' @examples
-#' \donttest{
-#' # This function is deprecated and its functionality can now be accessed
-#' # from [simulate_secondary()].
-#' # Here are some examples of how to use [simulate_secondary()] to replace
-#' # report_cases().
-#' # Old (using report_cases()):
-#' # Define case data
-#' cases <- example_confirmed[1:40]
-#' cases <- cases[, cases := as.integer(confirm)]
-#' cases <- cases[, confirm := NULL][, sample := 1]
-#' reported_cases <- report_cases(
-#'  case_estimates = cases,
-#'  delays = delay_opts(example_incubation_period + example_reporting_delay),
-#'  type = "sample"
-#' )
-#' print(reported_cases$samples)
-#'
-#' # New (using simulate_secondary()):
-#' cases <- example_confirmed[1:40]
-#' cases <- cases[, primary := as.integer(confirm)]
-#' report <- simulate_secondary(
-#'  cases,
-#'  delays = delay_opts(
-#'    fix_dist(example_incubation_period + example_reporting_delay)
-#'  ),
-#'  obs = obs_opts(family = "poisson")
-#' )
-#' print(report)
-#' }
-report_cases <- function(case_estimates,
-                         case_forecast = NULL,
-                         delays,
-                         type = "sample",
-                         reporting_effect,
-                         CrIs = c(0.2, 0.5, 0.9)) {
-  deprecate_warn(
-    when = "1.5.0",
-    what = "report_cases()",
-    with = "simulate_secondary()",
-    details = c(
-      "See equivalent examples using `simulate_secondary()`",
-      "in ?report_cases.",
-      "This function will be removed completely in version 2.0.0."
-    )
-  )
-  samples <- length(unique(case_estimates$sample))
-
-  # add a null reporting effect if missing
-  if (missing(reporting_effect)) {
-    reporting_effect <- data.table::data.table(
-      sample = list(1:samples),
-      effect = rep(1, 7),
-      day = 1:7
-    )
-    reporting_effect <- reporting_effect[,
-      .(sample = unlist(sample)), by = .(effect, day)
-    ]
-  }
-  # filter and sum nowcast to use only upscaled cases by date of infection
-  infections <- data.table::copy(case_estimates)
-
-  # add in case forecast if present
-  if (!is.null(case_forecast)) {
-    infections <- data.table::rbindlist(list(
-      infections,
-      case_forecast[, .(date, sample, cases = as.integer(cases))]
-    ), use.names = TRUE)
-  }
-
-  ## For each sample map to report date
-  report <- future.apply::future_lapply(1:max(infections$sample),
-    function(id) {
-      suppressWarnings(
-        EpiNow2::adjust_infection_to_report(infections[sample == id],
-          delay_defs = delays,
-          type = type,
-          reporting_effect = reporting_effect[sample == id, ]$effect
-        )
-      )
-    },
-    future.seed = TRUE
-  )
-
-  report <- data.table::rbindlist(report, idcol = "sample")
-
-  out <- list()
-  # bind all samples together
-  out$samples <- report
-  # summarise samples
-  out$summarised <- calc_summary_measures(
-    report[, value := cases][, cases := NULL],
-    summarise_by = "date",
-    order_by = "date",
-    CrIs = CrIs
-  )
-  return(out)
-}
-
 #' Provide Summary Statistics for Estimated Infections and Rt
 #' @description `r lifecycle::badge("questioning")`
 #' Creates a snapshot summary of estimates. May be removed in later releases as
diff --git a/R/setup.R b/R/setup.R
index 39453cf5e..107a1ece2 100644
--- a/R/setup.R
+++ b/R/setup.R
@@ -187,7 +187,7 @@ setup_future <- function(data,
 #' maps input to be a `{data.table}`
 #' @inheritParams estimate_infections
 #' @return A data table
-#' @export
+#' @keywords internal
 setup_dt <- function(data) {
   suppressMessages(data.table::setDTthreads(threads = 1))
   data <- data.table::setDT(data)
@@ -205,7 +205,7 @@ setup_dt <- function(data) {
 #' create if not present).
 #'
 #' @return A list containing the path to the dated folder and the latest folder
-#' @export
+#' @keywords internal
 setup_target_folder <- function(target_folder = NULL, target_date) {
   if (!is.null(target_folder)) {
     latest_folder <- file.path(target_folder, "latest")
diff --git a/R/simulate_infections.R b/R/simulate_infections.R
index 64562b048..0bd4631bd 100644
--- a/R/simulate_infections.R
+++ b/R/simulate_infections.R
@@ -11,7 +11,7 @@
 #' A previous function called [simulate_infections()] that simulates from a
 #' given model fit has been renamed [forecast_infections()]. Using
 #' [simulate_infections()] with existing estimates is now deprecated. This
-#' option will be removed in version 2.0.0.
+#' option will be removed in the next version.
 #' @param R a data frame of reproduction numbers (column `R`) by date (column
 #'   `date`). Column `R` must be numeric and `date` must be in date format. If
 #'   not all days between the first and last day in the `date` are present,
@@ -41,7 +41,6 @@
 #' @inheritParams estimate_infections
 #' @inheritParams rt_opts
 #' @inheritParams stan_opts
-#' @importFrom lifecycle deprecate_warn
 #' @importFrom checkmate assert_data_frame assert_date assert_numeric
 #'   assert_subset assert_integer
 #' @importFrom data.table data.table merge.data.table nafill rbindlist
@@ -82,7 +81,7 @@ simulate_infections <- function(estimates, R, initial_infections,
       "forecast_infections()",
       details = paste0(
         "The `estimates` option will be removed from [simulate_infections()] ",
-        "in version 2.0.0."
+        "in the next version."
       )
     )
     return(forecast_infections(estimates = estimates, ...))
diff --git a/R/summarise.R b/R/summarise.R
index 948421176..c52f43869 100644
--- a/R/summarise.R
+++ b/R/summarise.R
@@ -22,7 +22,7 @@
 #' @importFrom purrr safely map_chr map_dbl map_chr
 #' @importFrom data.table setorderv melt merge.data.table dcast
 #' @return A list of summary data
-#' @export
+#' @keywords internal
 summarise_results <- function(regions,
                               summaries = NULL,
                               results_dir = NULL,
@@ -202,7 +202,7 @@ regional_summary <- function(regional_output = NULL,
   if (is.null(regional_output)) {
     if (!is.null(results_dir)) {
       futile.logger::flog.info("Extracting results from: %s", results_dir)
-      regions <- EpiNow2::get_regions(results_dir)
+      regions <- get_regions(results_dir)
       if (is.null(target_date)) {
         target_date <- "latest"
       }
@@ -429,7 +429,7 @@ regional_summary <- function(regional_output = NULL,
 #' @seealso regional_summary
 #' @return A list of summarised Rt, cases by date of infection and cases by
 #' date of report
-#' @export
+#' @keywords internal
 #' @importFrom data.table setnames fwrite setorderv
 summarise_key_measures <- function(regional_results = NULL,
                                    results_dir = NULL, summary_dir = NULL,
@@ -491,6 +491,7 @@ summarise_key_measures <- function(regional_results = NULL,
   save_variable(out$cases_by_report, "cases_by_report")
   return(out)
 }
+
 #' Summarise Regional Runtimes
 #'
 #' @description `r lifecycle::badge("maturing")`
@@ -503,6 +504,7 @@ summarise_key_measures <- function(regional_results = NULL,
 #' @export
 #' @importFrom data.table data.table fwrite
 #' @importFrom purrr map safely
+#' @keywords internal
 #' @examples
 #' regional_out <- readRDS(system.file(
 #'     package = "EpiNow2", "extdata", "example_regional_epinow.rds"
diff --git a/R/utilities.R b/R/utilities.R
index c3e3be1a4..54177c107 100644
--- a/R/utilities.R
+++ b/R/utilities.R
@@ -150,6 +150,7 @@ R_to_growth <- function(R, gamma_mean, gamma_sd) {
 #' @param delay_var List of numeric delays
 #' @param no_delays Numeric, number of delays
 #' @return A numeric array
+#' @keywords internal
 allocate_delays <- function(delay_var, no_delays) {
   if (no_delays > 0) {
     out <- unlist(delay_var)
@@ -169,6 +170,7 @@ allocate_delays <- function(delay_var, no_delays) {
 #' empty if missing.
 #' @param n Numeric, number of samples to assign an empty array
 #' @return A list of parameters some allocated to be empty
+#' @keywords internal
 allocate_empty <- function(data, params, n = 0) {
   for (param in params) {
     if (!exists(param, data)) {
@@ -196,6 +198,7 @@ allocate_empty <- function(data, params, n = 0) {
 #'
 #' @return A logical vector of named output arguments
 #' @importFrom  futile.logger flog.info flog.debug
+#' @keywords internal
 match_output_arguments <- function(input_args = NULL,
                                    supported_args = NULL,
                                    logger = NULL,
@@ -340,41 +343,23 @@ discretised_gamma_pmf <- function(mean, sd, max_d, zero_pad = 0,
   return(pmf)
 }
 
-#' Update a List
-#'
-#' @description `r lifecycle::badge("deprecated")`
-#' Used to handle updating settings in a list. For example when making
-#' changes to [opts_list()] output.
-#' @param defaults A list of default settings
-#' @param optional A list of optional settings to override defaults
-#' @importFrom lifecycle deprecate_stop
-#' @return A list
-#' @export
-update_list <- function(defaults = list(), optional = list()) {
-  deprecate_stop(
-    when = "1.4.1",
-    what = "update_list()",
-    with = "utils::modifyList()",
-    details = "Please file an issue if deprecating this \
-      function has caused any issues."
-    )
-}
-
 #' Adds a day of the week vector
 #'
 #' @param dates Vector of dates
 #' @param week_effect Numeric from 1 to 7 defaults to 7
 #'
 #' @return A numeric vector containing the period day of the week index
-#' @export
+#' @keywords internal
 #' @importFrom lubridate wday
 #' @examples
+#' \dontrun{
 #' dates <- seq(as.Date("2020-03-15"), by = "days", length.out = 15)
 #' # Add date based day of week
 #' add_day_of_week(dates, 7)
 #'
 #' # Add shorter week
 #' add_day_of_week(dates, 4)
+#' }
 add_day_of_week <- function(dates, week_effect = 7) {
   if (week_effect == 7) {
     day_of_week <- lubridate::wday(dates, week_start = 1)
@@ -407,7 +392,6 @@ add_day_of_week <- function(dates, week_effect = 7) {
 #' data.table::getDTthreads()
 #' }
 #' @export
-
 set_dt_single_thread <- function() {
   a <- list2env(x = list(dt_previous_threads = data.table::getDTthreads()))
 
@@ -422,7 +406,6 @@ set_dt_single_thread <- function() {
 }
 
 #' @importFrom stats glm median na.omit pexp pgamma plnorm quasipoisson rexp
-#' @importFrom lifecycle deprecate_warn
 #' @importFrom stats rlnorm rnorm rpois runif sd var rgamma pnorm
 globalVariables(
   c(
diff --git a/_pkgdown.yml b/_pkgdown.yml
index 2abaf1338..2f74b0679 100644
--- a/_pkgdown.yml
+++ b/_pkgdown.yml
@@ -73,40 +73,22 @@ reference:
       - forecast_secondary
       - estimate_delay
       - estimate_truncation
-  - title: Internal epinow Functions
-    desc: Functions used internally by epinow
-    contents:
-      - update_horizon
-      - contains("save_")
-      - estimates_by_report_date
-      - copy_results_to_latest
-      - construct_output
-  - title: Internal estimate_infections Functions
+  - title: Specify an Argument
     desc: Functions used by estimate_infections
     contents:
       - contains("_opts")
       - opts_list
-      - contains("create_")
-      - allocate_delays
-      - allocate_empty
-  - title: Model Fitting
-    desc: Functions for model fitting
-    contents:
-      - contains("fit_")
-      - init_cumulative_fit
   - title: Summarise Across Regions
     desc: Functions used for summarising across regions (designed for use with regional_epinow)
     contents:
       - regional_summary
       - regional_runtimes
+      - get_regional_results
   - title: Summarise results
-    desc: Functions from summarising results
+    desc: Functions for summarising results
     contents:
       - contains("summary.")
-      - contains("summarise")
       - contains("calc_")
-      - contains("process_")
-      - format_fit
       - make_conf
       - map_prob_change
   - title: Plot Results
@@ -116,7 +98,6 @@ reference:
   - title: Report results
     desc: Functions to report results
     contents:
-      - report_cases
       - report_plots
       - report_summary  
   - title: Define, Fit and Parameterise Distributions
@@ -137,7 +118,6 @@ reference:
       - simulate_infections
       - simulate_secondary
       - convolve_and_scale
-      - adjust_infection_to_report
   - title: Data
     desc: Package datasets that may be used to parameterise other functions or in examples
     contents:
@@ -164,13 +144,13 @@ reference:
   - title: Utilities
     desc: Utility functions
     contents:
-      - match_output_arguments
       - run_region
       - expose_stan_fns
       - convert_to_logmean
       - convert_to_logsd
       - growth_to_R
       - R_to_growth
-      - add_day_of_week
       - update_secondary_args
-      - update_list
+  - title: internal
+    contents:
+      - has_keyword("internal")
diff --git a/man/add_day_of_week.Rd b/man/add_day_of_week.Rd
index 633dfa0de..41828262c 100644
--- a/man/add_day_of_week.Rd
+++ b/man/add_day_of_week.Rd
@@ -18,6 +18,7 @@ A numeric vector containing the period day of the week index
 Adds a day of the week vector
 }
 \examples{
+\dontrun{
 dates <- seq(as.Date("2020-03-15"), by = "days", length.out = 15)
 # Add date based day of week
 add_day_of_week(dates, 7)
@@ -25,3 +26,5 @@ add_day_of_week(dates, 7)
 # Add shorter week
 add_day_of_week(dates, 4)
 }
+}
+\keyword{internal}
diff --git a/man/adjust_infection_to_report.Rd b/man/adjust_infection_to_report.Rd
index 9c952d6de..decaa3a95 100644
--- a/man/adjust_infection_to_report.Rd
+++ b/man/adjust_infection_to_report.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/adjust.R
+% Please edit documentation in R/deprecated.R
 \name{adjust_infection_to_report}
 \alias{adjust_infection_to_report}
 \title{Adjust from Case Counts by Infection Date to Date of Report}
@@ -85,3 +85,4 @@ report <- simulate_secondary(
 print(report)
 }
 }
+\keyword{internal}
diff --git a/man/allocate_delays.Rd b/man/allocate_delays.Rd
index beacd41a3..62ca0a9cf 100644
--- a/man/allocate_delays.Rd
+++ b/man/allocate_delays.Rd
@@ -18,3 +18,4 @@ A numeric array
 \ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 Allocate delays for stan. Used in \code{\link[=delay_opts]{delay_opts()}}.
 }
+\keyword{internal}
diff --git a/man/allocate_empty.Rd b/man/allocate_empty.Rd
index e269e2aac..8246dfbdd 100644
--- a/man/allocate_empty.Rd
+++ b/man/allocate_empty.Rd
@@ -22,3 +22,4 @@ A list of parameters some allocated to be empty
 Allocate missing parameters to be empty two dimensional arrays. Used
 internally by \code{\link[=forecast_infections]{forecast_infections()}}.
 }
+\keyword{internal}
diff --git a/man/apply_tolerance.Rd b/man/apply_tolerance.Rd
index c74e3623c..b7e2c81b6 100644
--- a/man/apply_tolerance.Rd
+++ b/man/apply_tolerance.Rd
@@ -16,6 +16,7 @@ A \verb{<dist_spec>} where probability masses below the threshold level
 have been removed
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This removes any part of the tail of the nonparametric distributions in the
 <dist_spec> where the probability mass is below the threshold level.
 }
diff --git a/man/c.dist_spec.Rd b/man/c.dist_spec.Rd
index 7a9bfed87..9c8415679 100644
--- a/man/c.dist_spec.Rd
+++ b/man/c.dist_spec.Rd
@@ -13,6 +13,7 @@
 Combined delay distributions (with class \verb{<dist_spec>})
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This combines the parameters so that they can be fed as multiple delay
 distributions to \code{\link[=epinow]{epinow()}} or \code{\link[=estimate_infections]{estimate_infections()}}.
 }
diff --git a/man/check_reports_valid.Rd b/man/check_reports_valid.Rd
index 71798df1a..665019568 100644
--- a/man/check_reports_valid.Rd
+++ b/man/check_reports_valid.Rd
@@ -26,6 +26,7 @@ This is used to determine which checks to perform on the data input.}
 Called for its side effects.
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 \code{check_reports_valid()} checks that the supplied data is a \verb{<data.frame>},
 and that it has the right column names and types. In particular, it checks
 that the date column is in date format and does not contain NA's, and that
diff --git a/man/check_stan_delay.Rd b/man/check_stan_delay.Rd
index 75bf99d90..5839ad74a 100644
--- a/man/check_stan_delay.Rd
+++ b/man/check_stan_delay.Rd
@@ -13,6 +13,7 @@ check_stan_delay(dist)
 Called for its side effects.
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 \code{check_stan_delay()} checks that the supplied data is a \verb{<dist_spec>},
 that it is a supported distribution, and that is has a finite maximum.
 }
diff --git a/man/collapse.Rd b/man/collapse.Rd
index 26424f862..815ee0419 100644
--- a/man/collapse.Rd
+++ b/man/collapse.Rd
@@ -14,6 +14,7 @@ A \verb{<dist_spec>} where consecutive nonparametric distributions
 have been convolved
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This convolves any consecutive nonparametric distributions contained
 in the <dist_spec>.
 }
diff --git a/man/construct_output.Rd b/man/construct_output.Rd
index 2700aa77d..9a87d6de0 100644
--- a/man/construct_output.Rd
+++ b/man/construct_output.Rd
@@ -31,3 +31,4 @@ A list of output as returned by \code{epinow}
 \ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 Combines the output produced internally by \code{epinow} into a single list.
 }
+\keyword{internal}
diff --git a/man/convert_to_natural.Rd b/man/convert_to_natural.Rd
index ec1596efc..6b72ebf03 100644
--- a/man/convert_to_natural.Rd
+++ b/man/convert_to_natural.Rd
@@ -16,6 +16,7 @@ A list with two elements, \code{params_mean} and \code{params_sd}, containing
 mean and sd of natural parameters.
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This is used for preprocessing before generating a \code{dist_spec} object
 from a given set of parameters and distribution
 }
diff --git a/man/copy_results_to_latest.Rd b/man/copy_results_to_latest.Rd
index e9fb5d46b..0f912f880 100644
--- a/man/copy_results_to_latest.Rd
+++ b/man/copy_results_to_latest.Rd
@@ -21,3 +21,4 @@ No return value, called for side effects
 Copies output from the dated folder to a latest folder. May be undergo
 changes in later releases.
 }
+\keyword{internal}
diff --git a/man/create_backcalc_data.Rd b/man/create_backcalc_data.Rd
index 87190e25d..f87e778ca 100644
--- a/man/create_backcalc_data.Rd
+++ b/man/create_backcalc_data.Rd
@@ -18,9 +18,7 @@ A list of settings defining the Gaussian process
 Takes the output of \code{\link[=backcalc_opts]{backcalc_opts()}} and converts it into a list understood
 by stan.
 }
-\examples{
-create_backcalc_data(backcalc = backcalc_opts())
-}
 \seealso{
 backcalc_opts
 }
+\keyword{internal}
diff --git a/man/create_clean_reported_cases.Rd b/man/create_clean_reported_cases.Rd
index 988c9494f..f782e588e 100644
--- a/man/create_clean_reported_cases.Rd
+++ b/man/create_clean_reported_cases.Rd
@@ -47,5 +47,8 @@ which point 0 cases are replaced with a user supplied value (defaults to
 \code{NA}).
 }
 \examples{
+\dontrun{
 create_clean_reported_cases(example_confirmed, 7)
 }
+}
+\keyword{internal}
diff --git a/man/create_future_rt.Rd b/man/create_future_rt.Rd
index 03c51d654..86b476274 100644
--- a/man/create_future_rt.Rd
+++ b/man/create_future_rt.Rd
@@ -25,3 +25,4 @@ A list containing a logical called fixed and an integer called from
 Converts the \code{future} argument from \code{\link[=rt_opts]{rt_opts()}} into arguments that can be
 passed to stan.
 }
+\keyword{internal}
diff --git a/man/create_gp_data.Rd b/man/create_gp_data.Rd
index 427fe693e..94b666e01 100644
--- a/man/create_gp_data.Rd
+++ b/man/create_gp_data.Rd
@@ -23,6 +23,7 @@ Takes the output of \code{\link[=gp_opts]{gp_opts()}} and converts it into a lis
 stan.
 }
 \examples{
+\dontrun{
 # define input data required
 data <- list(
   t = 30,
@@ -39,6 +40,8 @@ create_gp_data(NULL, data)
 # custom lengthscale
 create_gp_data(gp_opts(ls_mean = 14), data)
 }
+}
 \seealso{
 \code{\link[=gp_opts]{gp_opts()}}
 }
+\keyword{internal}
diff --git a/man/create_initial_conditions.Rd b/man/create_initial_conditions.Rd
index 140d8cf5d..308a1b8ba 100644
--- a/man/create_initial_conditions.Rd
+++ b/man/create_initial_conditions.Rd
@@ -19,3 +19,4 @@ used to sample from the prior distributions (or as close as possible) for
 parameters. Used in order to initialise each stan chain within a range of
 plausible values.
 }
+\keyword{internal}
diff --git a/man/create_obs_model.Rd b/man/create_obs_model.Rd
index 6295cbcf4..736385fcf 100644
--- a/man/create_obs_model.Rd
+++ b/man/create_obs_model.Rd
@@ -22,6 +22,7 @@ Takes the output of \code{\link[=obs_opts]{obs_opts()}} and converts it into a l
 by stan.
 }
 \examples{
+\dontrun{
 dates <- seq(as.Date("2020-03-15"), by = "days", length.out = 15)
 # default observation model data
 create_obs_model(dates = dates)
@@ -37,6 +38,8 @@ create_obs_model(
 # Apply a custom week week length
 create_obs_model(obs_opts(week_length = 3), dates = dates)
 }
+}
 \seealso{
 \code{\link[=obs_opts]{obs_opts()}}
 }
+\keyword{internal}
diff --git a/man/create_rt_data.Rd b/man/create_rt_data.Rd
index 573b0aba6..4a02c1664 100644
--- a/man/create_rt_data.Rd
+++ b/man/create_rt_data.Rd
@@ -27,6 +27,7 @@ Takes the output from \code{\link[=rt_opts]{rt_opts()}} and converts it into a l
 stan.
 }
 \examples{
+\dontrun{
 # default Rt data
 create_rt_data()
 
@@ -36,6 +37,8 @@ create_rt_data(rt = NULL)
 # using breakpoints
 create_rt_data(rt_opts(use_breakpoints = TRUE), breakpoints = rep(1, 10))
 }
+}
 \seealso{
 rt_settings
 }
+\keyword{internal}
diff --git a/man/create_shifted_cases.Rd b/man/create_shifted_cases.Rd
index 80484fed4..1a1dff5a6 100644
--- a/man/create_shifted_cases.Rd
+++ b/man/create_shifted_cases.Rd
@@ -42,6 +42,7 @@ projected to the end of the forecast horizon. The initial part of the data
 any non-integer resulting values rounded up.
 }
 \examples{
+\dontrun{
 shift <- 7
 horizon <- 7
 smoothing_window <- 14
@@ -61,3 +62,5 @@ cases <- data.table::rbindlist(list(
  ))
 create_shifted_cases(cases, shift, smoothing_window, horizon)
 }
+}
+\keyword{internal}
diff --git a/man/create_stan_args.Rd b/man/create_stan_args.Rd
index 1c25217f9..065737d1e 100644
--- a/man/create_stan_args.Rd
+++ b/man/create_stan_args.Rd
@@ -44,9 +44,12 @@ Initialisation defaults to random but it is expected that
 \code{\link[=create_initial_conditions]{create_initial_conditions()}} will be used.
 }
 \examples{
+\dontrun{
 # default settings
 create_stan_args()
 
 # increasing warmup
 create_stan_args(stan = stan_opts(warmup = 1000))
 }
+}
+\keyword{internal}
diff --git a/man/create_stan_data.Rd b/man/create_stan_data.Rd
index 3e63c7db7..e755fca8d 100644
--- a/man/create_stan_data.Rd
+++ b/man/create_stan_data.Rd
@@ -52,8 +52,11 @@ construct a single list for input into stan with all data required
 present.
 }
 \examples{
+\dontrun{
 create_stan_data(
  example_confirmed, 7, rt_opts(), gp_opts(), obs_opts(), 7,
  backcalc_opts(), create_shifted_cases(example_confirmed, 7, 14, 7)
 )
 }
+}
+\keyword{internal}
diff --git a/man/create_stan_delays.Rd b/man/create_stan_delays.Rd
index 35e017c08..ea10c9118 100644
--- a/man/create_stan_delays.Rd
+++ b/man/create_stan_delays.Rd
@@ -18,3 +18,4 @@ A list of variables as expected by the stan model
 \description{
 Create delay variables for stan
 }
+\keyword{internal}
diff --git a/man/discretise.Rd b/man/discretise.Rd
index a24848ebd..cbf2cbc45 100644
--- a/man/discretise.Rd
+++ b/man/discretise.Rd
@@ -22,12 +22,13 @@ A \verb{<dist_spec>} where all distributions with constant parameters are
 nonparametric.
 }
 \description{
-Discretise a <dist_spec>
-}
-\details{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 By default it will discretise all the distributions it can discretise
 (i.e. those with finite support and constant parameters).
 }
+\details{
+Discretise a <dist_spec>
+}
 \examples{
 # A fixed gamma distribution with mean 5 and sd 1.
 dist1 <- Gamma(mean = 5, sd = 1, max = 20)
diff --git a/man/dist_spec.Rd b/man/dist_spec.Rd
index a2086b5f6..0af001f4a 100644
--- a/man/dist_spec.Rd
+++ b/man/dist_spec.Rd
@@ -57,3 +57,4 @@ functionality is still used internally). Construct distributions using
 the corresponding distribution function such as \code{\link[=Gamma]{Gamma()}}, \code{\link[=LogNormal]{LogNormal()}},
 \code{\link[=Normal]{Normal()}} or \code{\link[=Fixed]{Fixed()}} instead.
 }
+\keyword{internal}
diff --git a/man/epinow.Rd b/man/epinow.Rd
index 37742a5f5..b191cc5e1 100644
--- a/man/epinow.Rd
+++ b/man/epinow.Rd
@@ -118,7 +118,7 @@ A list of output from estimate_infections with additional elements
 summarising results and reporting errors if they have occurred.
 }
 \description{
-\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#maturing}{\figure{lifecycle-maturing.svg}{options: alt='[Maturing]'}}}{\strong{[Maturing]}}
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 This function wraps the functionality of \code{\link[=estimate_infections]{estimate_infections()}} in order
 to estimate Rt and cases by date of infection and forecast these infections
 into the future. In addition to the functionality of
diff --git a/man/estimates_by_report_date.Rd b/man/estimates_by_report_date.Rd
index 288c73b91..7d749567b 100644
--- a/man/estimates_by_report_date.Rd
+++ b/man/estimates_by_report_date.Rd
@@ -30,3 +30,4 @@ date of report.
 Either extracts or converts reported cases from an input data table. For
 output from \code{estimate_infections} this is a simple filtering step.
 }
+\keyword{internal}
diff --git a/man/extract_parameter.Rd b/man/extract_parameter.Rd
index e17227f9d..1ddf27206 100644
--- a/man/extract_parameter.Rd
+++ b/man/extract_parameter.Rd
@@ -23,3 +23,4 @@ sample value.
 Extracts a single from a list of stan output and returns it as a
 \verb{<data.table>}.
 }
+\keyword{internal}
diff --git a/man/extract_parameter_samples.Rd b/man/extract_parameter_samples.Rd
index e657732c4..85d59c3e9 100644
--- a/man/extract_parameter_samples.Rd
+++ b/man/extract_parameter_samples.Rd
@@ -39,3 +39,4 @@ parameter
 Extracts a custom set of parameters from a stan object and adds
 stratification and dates where appropriate.
 }
+\keyword{internal}
diff --git a/man/extract_params.Rd b/man/extract_params.Rd
index 34387b5f2..69e6c65b5 100644
--- a/man/extract_params.Rd
+++ b/man/extract_params.Rd
@@ -2,8 +2,7 @@
 % Please edit documentation in R/dist_spec.R
 \name{extract_params}
 \alias{extract_params}
-\title{Internal function for extracting given parameter names of a distribution
-from the environment. Called by \code{new_dist_spec}}
+\title{Extract parameter names}
 \usage{
 extract_params(params, distribution)
 }
@@ -16,6 +15,7 @@ extract_params(params, distribution)
 A character vector of parameters and their values.
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 Internal function for extracting given parameter names of a distribution
 from the environment. Called by \code{new_dist_spec}
 }
diff --git a/man/extract_single_dist.Rd b/man/extract_single_dist.Rd
index 9e79c5361..e3407cf6b 100644
--- a/man/extract_single_dist.Rd
+++ b/man/extract_single_dist.Rd
@@ -15,7 +15,7 @@ extract_single_dist(x, i)
 A single \code{dist_spec} object
 }
 \description{
-Extract a single element of a composite \verb{<dist_spec>}
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 }
 \examples{
 dist1 <- LogNormal(mean = 1.6, sd = 0.5, max = 20)
diff --git a/man/extract_static_parameter.Rd b/man/extract_static_parameter.Rd
index b99b6a934..5a39a6030 100644
--- a/man/extract_static_parameter.Rd
+++ b/man/extract_static_parameter.Rd
@@ -18,3 +18,4 @@ value
 \description{
 Extract Samples from a Parameter with a Single Dimension
 }
+\keyword{internal}
diff --git a/man/fit_model_approximate.Rd b/man/fit_model_approximate.Rd
index e5a2e9fc9..efa1b80b0 100644
--- a/man/fit_model_approximate.Rd
+++ b/man/fit_model_approximate.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/estimate_infections.R
+% Please edit documentation in R/fit.R
 \name{fit_model_approximate}
 \alias{fit_model_approximate}
 \title{Fit a Stan Model using an approximate method}
@@ -23,3 +23,4 @@ A stan model object
 \ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#maturing}{\figure{lifecycle-maturing.svg}{options: alt='[Maturing]'}}}{\strong{[Maturing]}}
 Fits a stan model using variational inference.
 }
+\keyword{internal}
diff --git a/man/fit_model_with_nuts.Rd b/man/fit_model_with_nuts.Rd
index e57b492f9..104c28064 100644
--- a/man/fit_model_with_nuts.Rd
+++ b/man/fit_model_with_nuts.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/estimate_infections.R
+% Please edit documentation in R/fit.R
 \name{fit_model_with_nuts}
 \alias{fit_model_with_nuts}
 \title{Fit a Stan Model using the NUTs sampler}
@@ -34,3 +34,4 @@ Fits a stan model using \code{\link[rstan:stanmodel-method-sampling]{rstan::samp
 run chains using \code{future} with error catching, timeouts and merging of
 completed chains.
 }
+\keyword{internal}
diff --git a/man/fix_dist.Rd b/man/fix_dist.Rd
index 4a78dd6a9..a4c400b6e 100644
--- a/man/fix_dist.Rd
+++ b/man/fix_dist.Rd
@@ -17,6 +17,7 @@ standard deviation from uncertainty given in the \verb{<dist_spec>}}
 A \verb{<dist_spec>} object without uncertainty
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 If the given \verb{<dist_spec>} has any uncertainty, it is removed and the
 corresponding distribution converted into a fixed one.
 }
diff --git a/man/format_fit.Rd b/man/format_fit.Rd
index cb222e6c6..b22935b2c 100644
--- a/man/format_fit.Rd
+++ b/man/format_fit.Rd
@@ -27,3 +27,4 @@ A list of samples and summarised posterior parameter estimates.
 \ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 Summaries posterior samples and adds additional custom variables.
 }
+\keyword{internal}
diff --git a/man/gamma_dist_def.Rd b/man/gamma_dist_def.Rd
index 39680f9d3..760e60bb4 100644
--- a/man/gamma_dist_def.Rd
+++ b/man/gamma_dist_def.Rd
@@ -65,3 +65,4 @@ def <- gamma_dist_def(
 print(def)
 def$params[[1]]
 }
+\keyword{internal}
diff --git a/man/generation_time_opts.Rd b/man/generation_time_opts.Rd
index a687024d6..844e8157c 100644
--- a/man/generation_time_opts.Rd
+++ b/man/generation_time_opts.Rd
@@ -11,7 +11,6 @@ generation_time_opts(
   source,
   max = 14,
   fixed = FALSE,
-  prior_weight,
   tolerance = 0.001,
   weight_prior = TRUE
 )
@@ -30,9 +29,6 @@ distribution is given a fixed generation time of 1 will be assumed.}
 
 \item{fixed}{deprecated; use \code{dist} instead}
 
-\item{prior_weight}{deprecated; prior weights are now specified as a
-model option. Use the \code{weight_prior} argument instead}
-
 \item{tolerance}{Numeric; the desired tolerance level.}
 
 \item{weight_prior}{Logical; if TRUE (default), any priors given in \code{dist}
diff --git a/man/get_dist.Rd b/man/get_dist.Rd
index 46ce221e1..36a20357c 100644
--- a/man/get_dist.Rd
+++ b/man/get_dist.Rd
@@ -30,3 +30,4 @@ using functions such as \code{\link[=Gamma]{Gamma()}} or \code{\link[=LogNormal]
 \seealso{
 \code{\link[=dist_spec]{dist_spec()}}
 }
+\keyword{internal}
diff --git a/man/get_generation_time.Rd b/man/get_generation_time.Rd
index 7b9e4f84a..6dc2db796 100644
--- a/man/get_generation_time.Rd
+++ b/man/get_generation_time.Rd
@@ -29,3 +29,4 @@ using functions such as \code{\link[=Gamma]{Gamma()}} or \code{\link[=LogNormal]
 \seealso{
 \code{\link[=dist_spec]{dist_spec()}}
 }
+\keyword{internal}
diff --git a/man/get_incubation_period.Rd b/man/get_incubation_period.Rd
index 98f97f023..accfaaa1d 100644
--- a/man/get_incubation_period.Rd
+++ b/man/get_incubation_period.Rd
@@ -29,3 +29,4 @@ using functions such as \code{\link[=Gamma]{Gamma()}} or \code{\link[=LogNormal]
 \seealso{
 \code{\link[=dist_spec]{dist_spec()}}
 }
+\keyword{internal}
diff --git a/man/get_raw_result.Rd b/man/get_raw_result.Rd
index e14fdf660..1344afd6f 100644
--- a/man/get_raw_result.Rd
+++ b/man/get_raw_result.Rd
@@ -22,3 +22,4 @@ An R object read in from the targeted \code{.rds} file
 \description{
 \ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 }
+\keyword{internal}
diff --git a/man/get_regions.Rd b/man/get_regions.Rd
index abb8d811f..65e3a938a 100644
--- a/man/get_regions.Rd
+++ b/man/get_regions.Rd
@@ -16,3 +16,4 @@ A named character vector containing the results to plot.
 \description{
 \ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 }
+\keyword{internal}
diff --git a/man/get_regions_with_most_reports.Rd b/man/get_regions_with_most_reports.Rd
index 0fd24c028..4c8ed438f 100644
--- a/man/get_regions_with_most_reports.Rd
+++ b/man/get_regions_with_most_reports.Rd
@@ -23,3 +23,4 @@ A character vector of regions with the highest reported cases
 Extract a vector of regions with the most reported cases in a set time
 window.
 }
+\keyword{internal}
diff --git a/man/get_seeding_time.Rd b/man/get_seeding_time.Rd
index 98df076c1..5bf04b32c 100644
--- a/man/get_seeding_time.Rd
+++ b/man/get_seeding_time.Rd
@@ -26,3 +26,4 @@ An integer seeding time
 The seeding time is set to the mean of the specified delays, constrained
 to be at least the maximum generation time
 }
+\keyword{internal}
diff --git a/man/init_cumulative_fit.Rd b/man/init_cumulative_fit.Rd
index b9ececde4..5e3c60e73 100644
--- a/man/init_cumulative_fit.Rd
+++ b/man/init_cumulative_fit.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/estimate_infections.R
+% Please edit documentation in R/deprecated.R
 \name{init_cumulative_fit}
 \alias{init_cumulative_fit}
 \title{Generate initial conditions by fitting to cumulative cases}
@@ -47,3 +47,4 @@ This implementation is based on the approach taken in
 \href{https://github.com/ImperialCollegeLondon/epidemia/}{epidemia} authored by
 James Scott.
 }
+\keyword{internal}
diff --git a/man/lognorm_dist_def.Rd b/man/lognorm_dist_def.Rd
index b36209cac..5eb9b031a 100644
--- a/man/lognorm_dist_def.Rd
+++ b/man/lognorm_dist_def.Rd
@@ -49,3 +49,4 @@ def <- lognorm_dist_def(
 print(def)
 def$params[[1]]
 }
+\keyword{internal}
diff --git a/man/lower_bounds.Rd b/man/lower_bounds.Rd
index 3b54b3e6b..28e62c5f1 100644
--- a/man/lower_bounds.Rd
+++ b/man/lower_bounds.Rd
@@ -13,6 +13,7 @@ lower_bounds(distribution)
 A numeric vector, the lower bounds.
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This is used to avoid sampling parameter values that have no support.
 }
 \examples{
diff --git a/man/match_output_arguments.Rd b/man/match_output_arguments.Rd
index 47123aa69..826aac4bd 100644
--- a/man/match_output_arguments.Rd
+++ b/man/match_output_arguments.Rd
@@ -31,3 +31,4 @@ A logical vector of named output arguments
 Match user supplied arguments with supported options and return a logical
 list for internal usage.
 }
+\keyword{internal}
diff --git a/man/max.dist_spec.Rd b/man/max.dist_spec.Rd
index 9db3a2d40..5e4aed9e9 100644
--- a/man/max.dist_spec.Rd
+++ b/man/max.dist_spec.Rd
@@ -15,6 +15,7 @@
 A vector of means.
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This works out the maximum of all the (parametric / nonparametric) delay
 distributions combined in the passed \code{\link[=dist_spec]{dist_spec()}} (ignoring any uncertainty
 in parameters)
diff --git a/man/mean.dist_spec.Rd b/man/mean.dist_spec.Rd
index b508ce7c8..fe38f65bf 100644
--- a/man/mean.dist_spec.Rd
+++ b/man/mean.dist_spec.Rd
@@ -16,6 +16,7 @@ parameters. If set to FALSE (the default) then the mean of any uncertain
 parameters will be returned as NA.}
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This works out the mean of all the (parametric / nonparametric) delay
 distributions combined in the passed \code{\link[=dist_spec]{dist_spec()}} (ignoring any uncertainty
 in parameters)
diff --git a/man/natural_params.Rd b/man/natural_params.Rd
index 6e7812b9f..1fc7bf594 100644
--- a/man/natural_params.Rd
+++ b/man/natural_params.Rd
@@ -13,6 +13,7 @@ natural_params(distribution)
 A character vector, the natural parameters.
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 These are the parameters used in the stan models. All other parameter
 representations are converted to these using \code{\link[=convert_to_natural]{convert_to_natural()}} before
 being passed to the stan models.
diff --git a/man/new_dist_spec.Rd b/man/new_dist_spec.Rd
index 1c649873f..0a97ff863 100644
--- a/man/new_dist_spec.Rd
+++ b/man/new_dist_spec.Rd
@@ -16,6 +16,7 @@ new_dist_spec(params, distribution)
 A \code{dist_spec} of the given specification.
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This will convert all parameters to natural parameters before generating
 a \code{dist_spec}. If they have uncertainty this will be done using sampling.
 }
diff --git a/man/plot.dist_spec.Rd b/man/plot.dist_spec.Rd
index c6b80f5fb..05cef9dce 100644
--- a/man/plot.dist_spec.Rd
+++ b/man/plot.dist_spec.Rd
@@ -12,6 +12,7 @@
 \item{...}{Additional arguments to pass to \code{{ggplot}}.}
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This function takes a \verb{<dist_spec>} object and plots its probability mass
 function (PMF) and cumulative distribution function (CDF) using \code{{ggplot2}}.
 Note that currently uncertainty in distributions is not plot.
diff --git a/man/plus-.dist_spec.Rd b/man/plus-.dist_spec.Rd
index 6a010bab4..7692bafbe 100644
--- a/man/plus-.dist_spec.Rd
+++ b/man/plus-.dist_spec.Rd
@@ -17,7 +17,7 @@ combine.}
 A delay distribution representing the sum of the two delays
 }
 \description{
-Creates a delay distribution as the sum of two other delay distributions.
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 }
 \examples{
 # A fixed lognormal distribution with mean 5 and sd 1.
diff --git a/man/print.dist_spec.Rd b/man/print.dist_spec.Rd
index 337b4781e..21dbba464 100644
--- a/man/print.dist_spec.Rd
+++ b/man/print.dist_spec.Rd
@@ -15,6 +15,7 @@
 invisible
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This displays the parameters of the uncertain and probability mass
 functions of fixed delay distributions combined in the passed \code{\link[=dist_spec]{dist_spec()}}.
 }
diff --git a/man/process_region.Rd b/man/process_region.Rd
index 590b1610f..0ec9888e3 100644
--- a/man/process_region.Rd
+++ b/man/process_region.Rd
@@ -39,3 +39,4 @@ logging information.
 \seealso{
 \code{\link[=regional_epinow]{regional_epinow()}}
 }
+\keyword{internal}
diff --git a/man/process_regions.Rd b/man/process_regions.Rd
index 47161d458..29ec0bd41 100644
--- a/man/process_regions.Rd
+++ b/man/process_regions.Rd
@@ -23,3 +23,4 @@ adds summary logging information.
 \seealso{
 \code{\link[=regional_epinow]{regional_epinow()}} \code{\link[=epinow]{epinow()}}
 }
+\keyword{internal}
diff --git a/man/regional_runtimes.Rd b/man/regional_runtimes.Rd
index 45d683b8b..945a0435a 100644
--- a/man/regional_runtimes.Rd
+++ b/man/regional_runtimes.Rd
@@ -41,3 +41,4 @@ regional_runtimes(regional_output = regional_out$regional)
 \seealso{
 regional_summary regional_epinow
 }
+\keyword{internal}
diff --git a/man/report_cases.Rd b/man/report_cases.Rd
index 64b2b921e..a729a7c6b 100644
--- a/man/report_cases.Rd
+++ b/man/report_cases.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/report.R
+% Please edit documentation in R/deprecated.R
 \name{report_cases}
 \alias{report_cases}
 \title{Report case counts by date of report}
@@ -79,3 +79,4 @@ report <- simulate_secondary(
 print(report)
 }
 }
+\keyword{internal}
diff --git a/man/sample_approx_dist.Rd b/man/sample_approx_dist.Rd
index b579cf9af..6ac8adfdb 100644
--- a/man/sample_approx_dist.Rd
+++ b/man/sample_approx_dist.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/adjust.R
+% Please edit documentation in R/deprecated.R
 \name{sample_approx_dist}
 \alias{sample_approx_dist}
 \title{Approximate Sampling a Distribution using Counts}
@@ -101,3 +101,4 @@ reports <- sample_approx_dist(
 )
 }
 }
+\keyword{internal}
diff --git a/man/save_estimate_infections.Rd b/man/save_estimate_infections.Rd
index 91cc68bc5..16e4b02e4 100644
--- a/man/save_estimate_infections.Rd
+++ b/man/save_estimate_infections.Rd
@@ -32,3 +32,4 @@ Saves output from \code{estimate_infections} to a target directory.
 \seealso{
 estimate_infections
 }
+\keyword{internal}
diff --git a/man/save_input.Rd b/man/save_input.Rd
index fcf533520..2f1ae2747 100644
--- a/man/save_input.Rd
+++ b/man/save_input.Rd
@@ -20,3 +20,4 @@ No return value, called for side effects
 \ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 Saves observed data to a target location if given.
 }
+\keyword{internal}
diff --git a/man/sd_dist.Rd b/man/sd_dist.Rd
index 22f5ca1f1..be8ef17f6 100644
--- a/man/sd_dist.Rd
+++ b/man/sd_dist.Rd
@@ -13,6 +13,7 @@ sd_dist(x)
 A vector of standard deviations.
 }
 \description{
+\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#experimental}{\figure{lifecycle-experimental.svg}{options: alt='[Experimental]'}}}{\strong{[Experimental]}}
 This works out the standard deviation of all the (parametric /
 nonparametric) delay distributions combined in the passed \code{\link[=dist_spec]{dist_spec()}}.
 }
diff --git a/man/setup_dt.Rd b/man/setup_dt.Rd
index a930cdf50..e3811555c 100644
--- a/man/setup_dt.Rd
+++ b/man/setup_dt.Rd
@@ -18,3 +18,4 @@ A data table
 Convenience function that sets the number of \code{{data.table}} cores to 1 and
 maps input to be a \code{{data.table}}
 }
+\keyword{internal}
diff --git a/man/setup_target_folder.Rd b/man/setup_target_folder.Rd
index 0872e260d..e24440214 100644
--- a/man/setup_target_folder.Rd
+++ b/man/setup_target_folder.Rd
@@ -20,3 +20,4 @@ A list containing the path to the dated folder and the latest folder
 \ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#stable}{\figure{lifecycle-stable.svg}{options: alt='[Stable]'}}}{\strong{[Stable]}}
 Sets up a folders for saving results
 }
+\keyword{internal}
diff --git a/man/simulate_infections.Rd b/man/simulate_infections.Rd
index 17941eba1..b3ac1147a 100644
--- a/man/simulate_infections.Rd
+++ b/man/simulate_infections.Rd
@@ -98,7 +98,7 @@ Uncertain parameters are not allowed.
 A previous function called \code{\link[=simulate_infections]{simulate_infections()}} that simulates from a
 given model fit has been renamed \code{\link[=forecast_infections]{forecast_infections()}}. Using
 \code{\link[=simulate_infections]{simulate_infections()}} with existing estimates is now deprecated. This
-option will be removed in version 2.0.0.
+option will be removed in the next version.
 }
 \examples{
 \donttest{
diff --git a/man/stan_opts.Rd b/man/stan_opts.Rd
index bcb7474a3..828f1940b 100644
--- a/man/stan_opts.Rd
+++ b/man/stan_opts.Rd
@@ -45,7 +45,7 @@ This implementation is based on the approach taken in
 James Scott.
 
 This argument is deprecated and the default (NULL) will be used from
-version 2.0.0.}
+the next version.}
 
 \item{return_fit}{Logical, defaults to TRUE. Should the fit stan model be
 returned.}
diff --git a/man/summarise_key_measures.Rd b/man/summarise_key_measures.Rd
index 913911229..1b1001980 100644
--- a/man/summarise_key_measures.Rd
+++ b/man/summarise_key_measures.Rd
@@ -41,3 +41,4 @@ Used internally by \code{regional_summary}.
 \seealso{
 regional_summary
 }
+\keyword{internal}
diff --git a/man/summarise_results.Rd b/man/summarise_results.Rd
index a317277eb..d621752e9 100644
--- a/man/summarise_results.Rd
+++ b/man/summarise_results.Rd
@@ -36,3 +36,4 @@ A list of summary data
 Used internally by \code{regional_summary} to produce a summary table of results.
 May be streamlined in later releases.
 }
+\keyword{internal}
diff --git a/man/update_horizon.Rd b/man/update_horizon.Rd
index 9b68027ad..387646e0b 100644
--- a/man/update_horizon.Rd
+++ b/man/update_horizon.Rd
@@ -24,3 +24,4 @@ Numeric forecast horizon adjusted for the users intention
 Makes sure that a forecast is returned for the user specified time period
 beyond the target date.
 }
+\keyword{internal}
diff --git a/man/update_list.Rd b/man/update_list.Rd
deleted file mode 100644
index 4ef385d1b..000000000
--- a/man/update_list.Rd
+++ /dev/null
@@ -1,21 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/utilities.R
-\name{update_list}
-\alias{update_list}
-\title{Update a List}
-\usage{
-update_list(defaults = list(), optional = list())
-}
-\arguments{
-\item{defaults}{A list of default settings}
-
-\item{optional}{A list of optional settings to override defaults}
-}
-\value{
-A list
-}
-\description{
-\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#deprecated}{\figure{lifecycle-deprecated.svg}{options: alt='[Deprecated]'}}}{\strong{[Deprecated]}}
-Used to handle updating settings in a list. For example when making
-changes to \code{\link[=opts_list]{opts_list()}} output.
-}
diff --git a/tests/testthat/test-estimate_truncation.R b/tests/testthat/test-estimate_truncation.R
index 9c2e7abe5..ccd7382ac 100644
--- a/tests/testthat/test-estimate_truncation.R
+++ b/tests/testthat/test-estimate_truncation.R
@@ -85,7 +85,6 @@ test_that("estimate_truncation works with zero_threshold set", {
 })
 
 test_that("deprecated arguments are recognised", {
-  options(warn = 2)
   expect_error(estimate_truncation(example_truncated,
     verbose = FALSE, trunc_max = 10
   ), "deprecated")