This document collects sources of evidence for baseline 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")
# 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**."
)
)
Parameter | Pf | Pv | Evidence |
---|---|---|---|
a50_beta | 80 | 80 | |
bM | 0.051398 | 0.051398 | |
consUl | 0.3300 | 0.3300 | |
delta_iR | 0.5620 | 0.5620 | |
delta_iR_prime | 0.01686 | 0.01686 | |
delta_iS | 1.1240 | 1.1240 | |
delta_iS_prime | 0.01124 | 0.01124 | |
delta_iS_scale | 2 | 2 | |
deltaIa_c50 | 26 | 26 | |
deltaIa_slope | 10 | 10 | |
deltaU_A | 0.74303 | 0.74303 | |
deltaU_c50 | 2954.3060 | 2954.3060 | |
deltaU_g | 43.7334 | 43.7334 | |
deltaU_kmax | 1.2069 | 1.2069 | |
deltaU_kmin | 2.1641 × 10−5 | 2.1641 × 10−5 | |
gamma | 7.203 × 109 | 7.203 × 109 | Normoblast population in the bone marrow |
gamma_M | 0.2500 | 0.2500 | |
gdUloss | 1 | 1 | |
I_t_1 | 100 | 100 | |
k_iR | 0.0300 | 0.0300 | |
k_iS | 0.0100 | 0.0100 | |
kappa | 1 × 10−9 | 1 × 10−9 | Reticulocyte release from bone marrow |
kM | 0.0100 | 0.0100 | |
knuIRBC | 3 | 3 | |
knuURBC | 1 | 1 | |
lambdaI.sel | 1.5 × 10−11 | 1.5 × 10−11 | |
lambdaU.sel | 5 × 10−7 | 5 × 10−7 | |
M_t_1 | 1.2 × 109 | 1.2 × 109 | Macrophage population in the spleen |
mag_deltaUr | 10 | 10 | |
MaxFold | 10 | 10 | Maximum fold increase in RBC production |
mu_deltaUr | 3.6500 | 3.6500 | |
nu | 0 × 10 | 0 × 10 | |
p_to_spleen | 1 | 1 | |
PMF | 8 | 8 | Parasite multiplication |
Pv_slope | 4.5000 | 4.5000 | |
rho0 | 0.0010 | 0.0010 | Reticulocyte release from bone marrow |
rho_inflection | 0.5000 | 0.5000 | Reticulocyte release from bone marrow |
rho_slope | 10 | 10 | Reticulocyte release from bone marrow |
sigma_deltaUr | 0.0025 | 0.0025 | |
sl_beta | 20 | 20 | |
slope_e | 16 | 16 | |
T_irbc | 48 | 48 | |
T_M | 108 | 108 | |
T_n | 120 | 120 | |
T_R | 84 | 84 | Reticulocyte release from bone marrow |
T_R_min | 24 | 24 | Reticulocyte release from bone marrow |
T_urbc | 2880 | 2880 | |
U_ss | 1.752 × 1013 | 1.752 × 1013 | Steady-state RBC population |
UIdelta50iRBC | 0.0001 | 0.0001 | |
UIdelta50uRBC | 1 × 10−7 | 1 × 10−7 | |
Ul | 5.7815 × 1012 | 5.7815 × 1012 | |
w_rc | 0.1000 | 0.1000 | |
zeta_a50 | 26 | 26 | |
zeta_sl | 10 | 10 |
Steven Kho:
I cannot find a reference showing this exact figure. Our own estimates of total macrophages in the spleen based on counts by histology is 8.7 × 109 cells, though our spleens are infected and enlarged, therefore may not represent baseline steady-state. I have done similar calculations from published histological counts in control spleens (Urban et al., 2005), and the number comes to 1.2 × 109 cells.
As per Watson et al. (2017):
In extreme anaemia [RBC production] can be increased fivefold or more, for example.
Figure 3 of Appendix 1 shows the haematocrit response and reticulocyte count response for max-fold values of 1..10.
Steven Kho:
The reference we have previously used for bone marrow cellularity and content is (Harrison, 1962) which states a mean of 11.1 × 109 nucleated cells/kg body weight, of which 28.4% are normoblasts. Therefore, in an average 60kg Papuan male, this data suggests the bone marrow contains 6.7 × 1011 normoblasts, which is about 2 log-folds higher than the current value in the model. If I look at other studies cited in this paper (Table 5), values do fall in this range when you adjust to total bodyweight.
At steady-state we obtain γ = 7.2 × 109, which results in a normoblast population of 8.64 × 1011. This is 29% higher than 6.7 × 1011.
Presumably we can’t assume that the splenectomised Papuans had an average bodyweight of ≈ 96 kg.
We can adjust the normoblast population at homeostasis by adjusting model parameters. However, adjusting the baseline rate parameter ρ0 and scaling parameter κ for the reticulocyte release rate results in a reduction of < 1%. Here are the results obtained by adjusting other parameters:
But the decrease in δumin also causes the circulating RBC population to remain stable during infections. With the current baseline parameters, the final circulating RBC count is ≈ 73.7% of the steady-state count. Decreasing δumin 100-fold results in a final circulating RBC count that is ≈ 96.6% of the steady-state count.
Note that the median total RBC count in asymptomatic participants in the HHS differs very little between uninfected (1.75 × 1013) and infected (1.73 × 1013 for Pf, 1.66 × 1013 for Pv).
Steven Kho:
We reported in the AJH paper (Kho et al., 2023) that asymptomatic infection resulted in a 8.1% loss of circulating RBCs, which translates to 91.9% of steady-state. Perhaps the rate of removal of RBCs (except young and old) in the spleen could be tweaked to a smaller fold-change to address this?
If we scale δumin by 0.15, the steady-state normoblast population is 8.64 × 1011, and the final circulating RBC count is ≈ 8.3% for Pf and ≈ 9.7% for Pv.
From Simpson et al. (2002):
The mean population estimate of “PMR” was approximately 8, and was highly dependent on the P. falciparum “strain”. PMR also varied significantly between patients with a 90% prediction interval varying from 5.5 to 12.3-fold.
As per the caption of Figure 2 in Koepke and Koepke (1986):
With increasing anaemia (and erythropoietin production) the maturation time of the erythroid marrow normoblasts and reticulocytes progressively shortens from a normal 3.5 days to 1.5 days or less. Conversely the reticulocytes in the peripheral blood persist for a longer time when the patient is anaemic.
The manuscript text also states:
The maturation time of the reticulocyte in the peripheral blood is taken as 1 day when the PCV is 0.45±0.5, 1.5 days when the PCV is 0.35±0.05, 2 days when the PCV is 0.25±0.05, and 3 days when the PCV is 0.15±0.05 (Hillman & Finch 1967).
Accordingly, the residence time in the bone marrow decreases from 3.5 days when the PCV is 0.45±0.5, down to 1.5 day when the PCV is 0.15±0.05.
The following function provides a reasonable characterisation, while maintaining a release age of 3.5 days in the absence of anaemia:
orig_retic_release_age <- function(u_frac) {
T_R <- 3.5 * 24
T_R_min <- 24
pkg_data <- new.env()
utils::data("rbc_steady_state", package = "spleenrbc", envir = pkg_data)
u_ss <- pkg_data$rbc_steady_state
u_rbc <- pmin(u_frac * u_ss, u_ss)
scale <- log(u_ss) - log(u_rbc)
t_release <- T_R_min + (T_R - T_R_min) * exp(-100 * scale)
t_release
}
new_retic_release_age <- function(u_frac, inflection, slope) {
T_R <- 3.5 * 24
T_R_min <- 24
u_frac <- pmin(u_frac, 1)
scale <- T_R - T_R_min
T_R_min + scale / (1 + exp(- slope * (u_frac - inflection)))
}
plot_retic_release_age <- function() {
u_frac <- seq(from = 0, to = 1.0, by = 1e-3)
slope <- 10
inflection <- 0.5
orig_age <- orig_retic_release_age(u_frac)
new_age <- new_retic_release_age(u_frac, inflection, slope)
df <- data.frame(
u_pcnt = rep(100 * u_frac, 2),
age = c(orig_age, new_age),
equation = rep(c("Original", "New"), each = length(u_frac))
)
df_koepke <- data.frame(
hematocrit = c(45, 35, 25, 15),
age = 24 * c(3.5, 3.0, 2.5, 1.5)
) |> dplyr::mutate(
u_pcnt = 100 * hematocrit / 45,
equation = "Koepke & Koepke (1986)"
)
ggplot(df, aes(u_pcnt, age / 24)) +
geom_line(mapping = aes(colour = equation)) +
geom_point(data = df_koepke) +
expand_limits(y = 0) +
scale_x_continuous(
"uRBC Population (% of steady-state)",
breaks = c(0, 33, 66.7, 100),
labels = c("0%", "33%", "67%", "100%")
) +
ylab("Minimum age of release (days)") +
scale_colour_brewer(
"Equation: ",
palette = "Dark2"
) +
theme(legend.position = "top")
}
plot_retic_release_age()
We used the median RBC count from Papuans with no malaria infection, who had not reported fever in the past 24 hours, as reported by Pava et al. (2016).