Evidence for parameter values

library(spleenrbc, warn.conflicts = FALSE)
suppressPackageStartupMessages(library(ggplot2))
suppressPackageStartupMessages(library(here))

output_dir <- here("outputs")

pf <- baseline_parameters("Pf")
pv <- baseline_parameters("Pv")

ss_pf <- retic_steady_state(pf)
ss_pv <- retic_steady_state(pv)

# Remove non-scalar derived parameters, such as the initial red blood cell
# populations.
remove_nonscalar_derived_parameters <- function(x) {
  x[sapply(x, function(x) length(x) == 1)]
}

# Define how to format parameter values of different magnitudes.
format_value <- function(x, nsmall = 4, digits = 5, min = 1e-4, max = 1e4,
                         delim = TRUE) {
  scientific <- x < min || x > max
  result <- format(x, digits = digits, nsmall = nsmall,
                   scientific = scientific)
  if (delim) {
    start <- "\\("
    end <- "\\)"
  } else {
    start <- ""
    end <- ""
  }
  if (scientific) {
    # Convert from scientific notation to "X x 10^Y".
    parts <- strsplit(result, "e")[[1]]
    digits <- parts[1]
    # Strip the plus sign from exponents.
    power <- sub("^\\+", "", parts[2])
    # Strip leading zeros from positive and negative exponents.
    power <- sub("^(-?)0*", "\\1", power)
    paste0(start, digits, "\\times 10^{", power, "}", end)
  } else {
    # Strip fractional digits for integer values.
    result <- sub("\\.0+$", "", result)
    paste0(start, result, end)
  }
}

# Convenience function for inline R expressions.
pretty <- function(x, digits = 3, nsmall = 0, delim = FALSE) {
  format_value(x, digits = digits, nsmall = nsmall, delim = delim)
}

# Return the baseline parameters for a species as a long data frame.
get_baseline_parameters <- function(species) {
  baseline_parameters(species = species) |>
    remove_nonscalar_derived_parameters() |>
    as.data.frame() |>
    dplyr::mutate(
      species = NULL,
      nospleen = NULL,
      Scenario = NULL
    ) |>
    tidyr::pivot_longer(
      everything(),
      names_to = "Parameter",
    ) |>
    dplyr::mutate(
      value = sapply(value, format_value)
    )
}

# Convert the baseline parameters for both species into data frames.
baseline <- dplyr::inner_join(
  get_baseline_parameters("Pf") |>
    dplyr::rename(Pf = value),
  get_baseline_parameters("Pv") |>
    dplyr::rename(Pv = value),
  by = "Parameter"
) |>
  # Sort by parameter name (case-insensitive).
  dplyr::arrange(tolower(Parameter)) |>
  # Highlight parameters whose values differ between species.
  dplyr::mutate(Parameter = dplyr::case_when(
    .default = Parameter,
    Pf != Pv ~ paste0("**", Parameter, "**")
  ))

# For each parameter, add a link to the relevant section that describes the
# available evidence (if any).
baseline <- baseline |>
  dplyr::mutate(Evidence = dplyr::case_when(
    Parameter == "gamma" ~ "[Normoblast population in the bone marrow]",
    Parameter == "M_t_1" ~ "[Macrophage population in the spleen]",
    Parameter == "MaxFold" ~ "[Maximum fold increase in RBC production]",
    Parameter == "PMF" ~ "[Parasite multiplication]",
    Parameter == "rho0" ~ "[Reticulocyte release from bone marrow]",
    Parameter == "rho_slope" ~ "[Reticulocyte release from bone marrow]",
    Parameter == "rho_inflection" ~ "[Reticulocyte release from bone marrow]",
    Parameter == "kappa" ~ "[Reticulocyte release from bone marrow]",
    Parameter == "T_R" ~ "[Reticulocyte release from bone marrow]",
    Parameter == "T_R_min" ~ "[Reticulocyte release from bone marrow]",
    Parameter == "U_ss" ~ "[Steady-state RBC population]",
    .default = ""
  ))

knitr::kable(
  baseline,
  caption=paste(
    "Evidence for chosen parameter values.",
    "Parameters that differ between species are shown **in bold**."
  )
)