This function estimates parameters for a dual-logistic model applied to RR interval (RRi) signals. The model is designed to capture both the rapid drop and subsequent recovery of RRi values during exercise and recovery periods, as described in Castillo-Aguilar et al. (2025). A robust Huber loss function (with a default \(\delta\) of 50 ms) is used to downweight the influence of outliers, ensuring that the optimization process is robust even in the presence of noisy or ectopic beats.
estimate_RRi_curve(
time,
RRi,
start_params = c(alpha = 800, beta = -380, c = 0.85, lambda = -3, phi = -2, tau = 6,
delta = 3),
lower_lim = c(alpha = 300, beta = -750, c = 0.1, lambda = -10, phi = -10, tau =
min(time), delta = min(time)),
upper_lim = c(alpha = 2000, beta = -10, c = 2, lambda = -0.1, phi = -0.1, tau =
max(time), delta = max(time)),
method = "L-BFGS-B"
)A list containing:
A data frame with columns for time, the original RRi values, and the fitted values obtained from the dual-logistic model.
The optimization method used.
The estimated parameters from the model.
The final value of the objective (Huber loss) function.
An integer code indicating convergence (0 indicates success).
A numeric vector of time points.
A numeric vector of RR interval values.
A named numeric vector or list of initial parameter estimates. Default is
c(alpha = 800, beta = -380, c = 0.85, lambda = -3, phi = -2, tau = 6, delta = 3).
A named numeric vector specifying the lower bound for each parameter.
Default is c(alpha = 300, beta = -750, c = 0.1, lambda = -10, phi = -10, tau = min(time), delta = min(time)).
A named numeric vector specifying the upper bound for each parameter.
Default is c(alpha = 2000, beta = -10, c = 2.0, lambda = -0.1, phi = -0.1, tau = max(time), delta = max(time)).
A character string specifying the optimization method to use with optim(). The default is
"L-BFGS-B", which allows for parameter constraints.
The dual-logistic model, as described in Castillo-Aguilar et al. (2025), is defined as:
$$ RRi(t) = \alpha + \frac{\beta}{1 + e^{\lambda (t - \tau)}} + \frac{-c \cdot \beta}{1 + e^{\phi (t - \tau - \delta)}} $$
where:
is the baseline RRi level.
controls the amplitude of the drop.
modulates the steepness of the drop phase.
represents the time at which the drop is centered.
scales the amplitude of the recovery relative to \(\beta\).
controls the steepness of the recovery phase.
shifts the recovery phase in time relative to \(\tau\).
The function first removes any missing cases from the input data and then defines the dual-logistic model, which represents the dynamic behavior of RR intervals during exercise and recovery. The objective function is based on the Huber loss (with a default threshold of 50 ms), which provides robustness against outliers by penalizing large deviations less harshly than the standard squared error. This objective function quantifies the discrepancy between the observed RRi values and those predicted by the model.
Parameter optimization is performed using optim() with box constraints when the "L-BFGS-B" method is
used. These constraints ensure that the parameters remain within physiologically plausible ranges. For other optimization
methods, the bounds are ignored by setting the lower limit to -Inf and the upper limit to Inf.
It is important to note that the default starting parameters and bounds provided in the function are general
guidelines and may not be optimal for every dataset or experimental scenario. Users are encouraged to customize the
starting parameters (start_params) and, if necessary, the lower and upper bounds (lower_lim and upper_lim)
based on the specific characteristics of their RRi signal. This customization is crucial for achieving robust convergence
and accurate parameter estimates in diverse applications.
Castillo-Aguilar, et al. (2025). Enhancing Cardiovascular Monitoring: A Non-linear Model for Characterizing RR Interval Fluctuations in Exercise and Recovery. Scientific Reports, 15(1), 8628.
true_params <- c(alpha = 800, beta = -300, c = 0.80,
lambda = -3, phi = -1, tau = 6, delta = 3)
time_vec <- seq(0, 20, by = 0.01)
set.seed(1234)
# Simulate an example RRi signal:
RRi_simulated <- dual_logistic(time_vec, true_params) +
rnorm(length(time_vec), sd = 30)
# Estimate the model parameters:
fit <- estimate_RRi_curve(time = time_vec, RRi = RRi_simulated)
# Print method
print(fit)
# Summary method
summary(fit)
# Plot method
plot(fit)
Run the code above in your browser using DataLab