Add logistic-CDF decay function

Re #52

NAMESPACE

@@ -8,6 +8,7 @@ export(cumulative_interval)
 export(decay_binary)
 export(decay_exponential)
 export(decay_linear)
+export(decay_logistic)
 export(decay_power)
 export(decay_stepped)
 export(fgt_poverty)

R/decay_logistic.R

@@ -0,0 +1,87 @@
+#' Logistic-CDF decay function
+#'
+#' Returns a logist-cdf weighting function parameterized with the
+#' median (inflection point) and standard deviation to be used inside
+#' accessibility calculating functions.
+#'
+#' @template description_generic_cost
+#'
+#' @param cutoff A `numeric` vector. The median or inflection point
+#'   of the logistic CDF in minutes of travel time.
+#' @param sd A `numeric` vector with same length as `cutoff`.
+#'   The standard deviation in minutes of the logistic-CDF decay function
+#'   must be greater than 0 and less than 120.
+#'   Values near 0 result approximate binary decay, values near 120
+#'   approximate linear decay.
+#'
+#' @return A `function` that takes a generic travel cost vector (`numeric`) as
+#'   an input and returns a vector of weights (`numeric`).
+#'
+#' @details When using a function created with `decay_logistic()`, the
+#'   output is named after the combination of cutoff (`"T"`) and sd (`"s"`)
+#'   - e.g. given the cutoff `c(10, 20)` and the sd `c(10, 20)`,
+#'   the outputs will be named `"T10;s10"`, `"T20;s20"`.
+#'
+#' @family decay functions
+#'
+#' @examplesIf identical(tolower(Sys.getenv("NOT_CRAN")), "true")
+#' weighting_function <- decay_logistic(
+#'   cutoff = seq(10, 120, by = 10),
+#'   sd = 10
+#' )
+#'
+#' weighting_function(seq(0, 120, by = 5))
+#'
+#' weighting_function <- decay_logistic(
+#'   c(10, 10, 10, 10, 20, 20, 20, 20),
+#'   c(2, 4, 6, 8, 10, 12, 2, 4, 6, 8, 10, 12)
+#' )
+#'
+#' weighting_function(seq(0, 120, by = 5))
+#'
+#' @export
+decay_logistic <- function(cutoff, sd) {
+  checkmate::assert_numeric(
+    cutoff,
+    lower = 0.001,
+    finite = TRUE,
+    any.missing = FALSE,
+    min.len = 1,
+    unique = TRUE,
+    sorted = TRUE
+  )
+  checkmate::assert_numeric(
+    sd,
+    lower = 0.001,
+    upper = 119.999,
+    any.missing = FALSE,
+    len = length(cutoff)
+  )
+
+  SQRT3 = sqrt(3)
+  g = function(travel_cost, med_m, sd_m) {
+    1 + exp(((travel_cost - med_m) * pi) / (sd_m * SQRT3))
+  }
+
+  weighting_function <- function(travel_cost) {
+    sd_list <- mapply(
+      meds = cutoff,
+      sds = sd,
+      FUN = function(meds, sds) {
+        vapply(
+          travel_cost,
+          function(x) g(0, meds, sds) / g(x, meds, sds),
+          numeric(1)
+        )
+      },
+      SIMPLIFY = FALSE
+    )
+
+    list_names <- sprintf('T%0.0f;s%0.0f', cutoff, sd)
+    names(sd_list) <- list_names
+
+    return(sd_list)
+  }
+
+  return(weighting_function)
+}

tests/testthat/test-decay_logistic.R

@@ -0,0 +1,60 @@
+tester <- function(cutoff = c(20, 30), sd = c(2, 10)) {
+  decay_logistic(cutoff, sd)
+}
+
+test_that("adequately raises errors", {
+  expect_error(tester("a"))
+  expect_error(tester(0))
+  expect_error(tester(c(25, Inf)))
+  expect_error(tester(c(25, NA)))
+  expect_error(tester(integer()))
+  expect_error(tester(c(25, 25)))
+  expect_error(tester(c(25, 20)))
+
+  expect_error(tester(c(25, 50), "a"))
+  expect_error(tester(c(25, 50), 0))
+  expect_error(tester(c(25, 50), 120))
+  expect_error(tester(c(25, 50), c(25, Inf)))
+  expect_error(tester(c(25, 50), c(25, NA)))
+  expect_error(tester(c(25, 50), integer()))
+
+})
+
+test_that("output is a decay function that returns a list of numeric vctrs", {
+  expect_is(tester(), "function")
+
+  output_fn <- tester(c(20, 30), c(2, 10))
+  output_list <- output_fn(0)
+  expect_is(output_list, "list")
+  expect_length(output_list, 2L)
+  expect_named(output_list, c("T20;s2", "T30;s10"))
+
+  expect_equal(
+    output_fn(c(0, 20))[["T20;s2"]],
+    c(1, 0.5)
+  )
+
+  expect_equal(
+    output_fn(c(0, 30))[["T30;s10"]],
+    c(1, 0.502),
+    tolerance = 0.001
+  )
+
+  output_fn <- tester(20, 2)
+  output_list <- output_fn(0)
+  expect_is(output_list, "list")
+  expect_length(output_list, 1L)
+  expect_named(output_list, "T20;s2")
+
+  expect_equal(
+    output_fn(c(0, 20))[["T20;s2"]],
+    c(1, 0.5)
+  )
+})
+
+test_that("output fn returns empty numeric if receives empty numeric/integer", {
+  output_fn <- tester(20, 2)
+
+  expect_identical(output_fn(integer())[["T20;s2"]], numeric())
+  expect_identical(output_fn(numeric())[["T20;s2"]], numeric())
+})

vignettes/decay_functions.Rmd

@@ -3,8 +3,8 @@ title: "Decay functions"
 output: rmarkdown::html_vignette
 bibliography: '`r system.file("REFERENCES.bib", package = "accessibility")`'
 vignette: >
-  %\VignetteIndexEntry{Decay functions} 
-  %\VignetteEngine{knitr::rmarkdown} 
+  %\VignetteIndexEntry{Decay functions}
+  %\VignetteEngine{knitr::rmarkdown}
   %\VignetteEncoding{UTF-8}
 ---
 
@@ -119,6 +119,22 @@ is the total number of steps, $S_{k}$ is the travel cost cutoff that delimits
 the $k^{th}$ step, and $v_{k}$ is the value that the decay function assumes at
 the $k^{th}$ step.
 
+## Logistic-CDF
+
+Weights decay sigmoidally for moderate values of the standard deviation parameter.
+Approximates binary decay at standard deviation values near 0, and linear
+decay at standard deviation values near 120.
+
+$$
+\begin{aligned}
+  g(t_{ij}, m, s) &= 1 + e^\frac{(t_{ij} - m) * \pi}{\sqrt{3}s} \\
+  f(t_{ij}) &= \frac{g(0, m, s)}{g(t_{ij}, m, s)}
+\end{aligned}
+$$
+
+Where $t_{ij}$ is the travel cost between origin *i* and destination *j*, $T$
+is the median or inflection point, and $s$ is the standard deviation of the logistic-CDF.
+
 ### Quick demonstration
 
 All decay functions (`decay_*()`) take decay parameters as input and return a