| Title: | Time Series Methods Based on Growth Curves |
|---|---|
| Description: | Provides tools for modelling and forecasting epidemic trajectories using a dynamic Gompertz model within a state space framework, with the Kalman filter for robust estimation of non-linear growth. Includes a reinitialization feature to adapt to new waves, and a leading-indicator extension that uses a related series moving ahead of the variable of interest (e.g. cases ahead of hospitalisations) to improve short-horizon forecasts, with model and lag selection via rolling-origin cross-validation. Applicable to data at daily, monthly, quarterly, or annual frequency, and to non-epidemic trajectories with similar dynamics, such as innovation diffusion and product adoption. Includes functions for data preprocessing, model fitting, forecast visualization, and accuracy evaluation using standard error measures. Methods are described in Harvey and Kattuman (2020) <doi:10.1162/99608f92.828f40de>, Harvey and Kattuman (2021) <doi:10.1098/rsif.2021.0179>, and Ashby, Harvey, Kattuman, Tang, and Thamotheram (2024) <https://www.jbs.cam.ac.uk/wp-content/uploads/2024/03/cchle-tsgc-paper-2024.pdf>. |
| Authors: | Michael Ashby [aut, cre], Paul Kattuman [aut], Andrew Harvey [aut], Edwin Tang [aut], Craig Thamotheram [aut], Guglielmo Secchi [aut], Cambridge Centre for Health Leadership & Enterprise, Cambridge Judge Business School, University of Cambridge [fnd], UK Health Security Agency [fnd], Keynes Fund, Faculty of Economics, University of Cambridge [fnd] (Supported A. Harvey), University of Cambridge Social Science Impact Fund [fnd] (Supported P. Kattuman), Cambridge Mathematics Placements Programme, University of Cambridge [fnd] (Supported E. Tang), Caldecott Bursary Fund, Magdalene College, University of Cambridge [fnd] (Supported E. Tang), Chancellor's Scholarship Scheme, University of Warwick [fnd] (Supported E. Tang), Statistics Centre for Doctoral Training, University of Warwick [fnd] (Supported E. Tang) |
| Maintainer: | Michael Ashby <[email protected]> |
| License: | GPL (>= 3) |
| Version: | 2.0.0 |
| Built: | 2026-10-10 16:13:32 UTC |
| Source: | https://github.com/edwintang903/tsgc |
idx_series by positionSubsets an idx_series using integer positions (in
the same units as x$start), not row numbers within data.
For example, if x$start == 10, then x[10:12] returns the
first three observations of x.
## S3 method for class 'idx_series' x[i, j, ...]## S3 method for class 'idx_series' x[i, j, ...]
x |
An |
i |
Integer vector of positions to select. |
j |
Optional column selector (for matrix-valued series). |
... |
Unused. |
Accessor method to access alphahat, the smoothed state
estimates, from a fitted KFS object.
alphahat(object)alphahat(object)
object |
A |
The smoothed state estimates alphahat.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::alphahat(output(res))library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::alphahat(output(res))
Similar to Python's argmax function.
argmax(x, decreasing = TRUE)argmax(x, decreasing = TRUE)
x |
Object to have its maximum found; either an |
decreasing |
Logical value indicating whether |
If x is an idx_series, a length-1 idx_series
at the position of the maximum. Otherwise the maximum value.
x <- idx_series(cumsum(rpois(30, 5)) + 1, start = 1) argmax(x)x <- idx_series(cumsum(rpois(30, 5)) + 1, start = 1) argmax(x)
idx_series
Coerces plain vectors/matrices to idx_series, or
passes through an idx_series unchanged. Convenience helper so
that internal functions can accept either raw data or an already
constructed idx_series.
as_idx_series(x, start = 1L)as_idx_series(x, start = 1L)
x |
A numeric vector, matrix or |
start |
Integer start position, only used if |
An idx_series object.
idx_series to a matrixCoerce an idx_series to a matrix
## S3 method for class 'idx_series' as.matrix(x, ...)## S3 method for class 'idx_series' as.matrix(x, ...)
x |
An |
... |
Unused. |
idx_series to a plain numeric vector or matrixas.numeric() is a primitive generic in base R and
does not always reliably dispatch S3 methods registered only under the
as.numeric.* name. as.double.idx_series is also provided
for the same reason.
## S3 method for class 'idx_series' as.numeric(x, ...) ## S3 method for class 'idx_series' as.double(x, ...)## S3 method for class 'idx_series' as.numeric(x, ...) ## S3 method for class 'idx_series' as.double(x, ...)
x |
An |
... |
Unused. |
Accessor method to access the filtered state estimates in a
KFS object.
att(object)att(object)
object |
A |
The filtered state estimates.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::att(output(res))library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::att(output(res))
This function performs a walk-forward validation to compare
forecasting performance across different models specified by the user.
It returns a data frame of a user-specified error metric (e.g. MAPE, MAE)
for forecasts n.ahead positions ahead, using the given models with
varying end positions.
cross_val( Y, model_list, est.end, n.ahead = 7, n.estimate = 1, gap = 1, xpred_targ.full = NULL, xpred_lead.full = NULL, LeadIndCol = 1, criterion = "mape" )cross_val( Y, model_list, est.end, n.ahead = 7, n.estimate = 1, gap = 1, xpred_targ.full = NULL, xpred_lead.full = NULL, LeadIndCol = 1, criterion = "mape" )
Y |
An |
model_list |
A list containing |
est.end |
The initial estimation end position for model fitting.
Starting from this position, the function re-estimates the model and
evaluates the performance for each lag in |
n.ahead |
Integer specifying the number of positions to forecast ahead for MAPE evaluation. |
n.estimate |
Integer indicating the total number of walk-forward validation steps to report. |
gap |
Integer specifying the position gap between two successive validations, where the model is re-estimated and evaluated during the walk-forward validation. |
xpred_targ.full |
An |
xpred_lead.full |
(Only required for leading indicator models) An
|
LeadIndCol |
(Only required for leading indicator models) Integer
representing the column number in |
criterion |
A string object indicating how to compare between different models. Available choices are "mape" (by default), "smape", "mae" and "rmse". |
A table summarizing the chosen error metric for each model in
model_list across the specified positions.
Helper method to compute the log growth rates of cumulated variables.
df2ldl(dt)df2ldl(dt)
dt |
Cumulated data series, as an |
An idx_series of log growth rates of the cumulated
variable inputted via the parameter dt.
x <- idx_series(cumsum(rpois(30, 5)) + 1) df2ldl(x)x <- idx_series(cumsum(rpois(30, 5)) + 1) df2ldl(x)
Helper method to compute the successive increments and log growth rates of cumulated variables. It will compute the successive increments and log cumulative growth rate for each column in the 2-column series, which will then be used to predict or estimate with the leading indicator model.
df2ldl_lead(data, LeadIndCol = 1)df2ldl_lead(data, LeadIndCol = 1)
data |
Cumulated data series as an |
LeadIndCol |
Column number of |
A list of idx_series with the original cumulative
variables, successive increments and log growth rates: cLead,
cTarg, newLead, newTarg, LDLlead,
LDLtarg.
Daily Covid-19 metrics for England.
data(england)data(england)
An object of class "xts" with columns:
Cumulative cases of Covid-19
Cumulative hospital admissions with Covid-19 since the start of the pandemic.
Cumulative deaths within 28 days of a positive test for Covid-19 by death date since the start of the pandemic.
Number of patients in hospital with confirmed Covid-19 each day. Note this is NOT cumulative.
Downloaded from https://ukhsa-dashboard.data.gov.uk/topics/covid-19. Now at https://ukhsa-dashboard.data.gov.uk/covid-19-archive-data-download. Full documentation from UKHSA included in downloads from the website.
data(england) # plot daily cases plot(diff(england))data(england) # plot daily cases plot(diff(england))
Daily Weather data of England in 2021.
data(england_weather_2021)data(england_weather_2021)
An object of class "xts" with columns:
Daily mean air temperature in degree Celsius (C) at 2 metres
Daily maximum temperature in degree Celsius (C) at 2 metres
Relative humidity at 2 metres (%)
Wind speed in metres per second (m/s) at 10 metres
Precipitation (millimetres per day; mm/day)
Downloaded from https://cds.climate.copernicus.eu/datasets/reanalysis-era5-single-levels?tab=overview and aggregated across the latitude and longitude corners for England by the authors.
data(england_weather_2021)data(england_weather_2021)
Accessor method that calls the estimate method of a
SSModelDynamicGompertz or SSModelLeadingIndicator object.
estimate(model)estimate(model)
model |
A |
A FilterResults or FilterResultsLI object containing the
estimation results.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) estimate(model)library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) estimate(model)
Computes the (instantaneous) reproduction number
, where is the
filtered or smoothed growth rate of the incidence variable
returned by res$get_gy_ci(), and gen_int is the mean
generation interval of the disease/process being modelled (in the same
integer-position units as res). This estimate, and its confidence
interval, is a deterministic transformation of get_gy_ci()'s
output; res must already have been produced by estimate().
estimate_r0( res, gen_int, n.ahead = 7, smoothed = FALSE, confidence.level = 0.68 )estimate_r0( res, gen_int, n.ahead = 7, smoothed = FALSE, confidence.level = 0.68 )
res |
A |
gen_int |
The mean generation interval, in the same integer-position
units as |
n.ahead |
Number of most recent integer positions to return, taken
from the end of |
smoothed |
Logical value indicating whether to use the smoothed
( |
confidence.level |
Confidence level for the confidence interval
around |
A data frame containing Position, fit,
lower and upper. If res$calendar is available, a
Date column is included before Position. Each row gives the
estimated reproduction number and its confidence interval over one of the
last n.ahead integer positions.
library(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) cal <- idx_calendar(anchor = as.Date("2021-01-01"), anchor_pos = 1L, amount = 1, unit = "days") model <- SSModelDynamicGompertz$new(Y = Y, q = NULL, end = 100, calendar = cal) res <- estimate(model) estimate_r0(res, gen_int = 5, n.ahead = 7)library(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) cal <- idx_calendar(anchor = as.Date("2021-01-01"), anchor_pos = 1L, amount = 1, unit = "days") model <- SSModelDynamicGompertz$new(Y = Y, q = NULL, end = 100, calendar = cal) res <- estimate(model) estimate_r0(res, gen_int = 5, n.ahead = 7)
Monthly downloads of selected leading apps that allow for online share trading in France from January 2017 to July 2021, by app.
data(etrading_apps)data(etrading_apps)
An object of class "xts" with cumulative downloads of each trading app
Taken from https://www.statista.com/statistics/1260008/etrading-app-downloads-france/
data(etrading_apps)data(etrading_apps)
Class for estimated Dynamic Gompertz Curve model and contains
methods to extract smoothed/filtered estimates of the states, the level of
the incidence variable , and forecasts of . The output from the estimate method
of the SSModelDynamicGompertz class is of the class FilterResults.
dataAn idx_series object containing the non-reinitialized
cumulated variable.
xpred_logicalLogical value indicating whether exogenous predictors were used to estimate the FilterResults object.
indexThe integer positions of the observations used in
estimation (of data).
reinit.idxThe reinitialisation position (a single integer) of the
estimated SSModelDynamicGompertz model (if applicable).
ar1Logical value indicating whether an ar1 component should be included in the model.
outputA KFS results object obtained after fitting a
SSModelDynamicGompertz model.
xpred.newAn idx_series object containing exogenous
predictors to be used in prediction. Defaults to NULL, and should
be provided if xpred is used for model estimation.
sea.periodThe period of seasonality, inherited from the estimated
SSModelDynamicGompertz model. For a day-of-the-week
effect with daily data, this would be 7.
calendarAn optional idx_calendar object, inherited
unchanged from the SSModelDynamicGompertz model this was
estimated from. Purely cosmetic: describes how the integer positions in
data map back to calendar time, for use by plotting only.
Defaults to NULL.
get_growth_y(smoothed = FALSE, return.components = FALSE)Returns the growth rate of the incidence () of the cumulated
variable (). Computed as
smoothed Logical value indicating whether to use the
smoothed estimates of and to compute the
growth rate (TRUE), or the contemporaneous filtered estimates
(FALSE). Default is FALSE.
return.components Logical value indicating whether to
return the estimates of and as well as
the estimates of the growth rate, or just the growth rate. Default is
FALSE.
idx_series objects containing
smoothed/filtered growth rates and components ( and
), where applicable.
get_gy_ci(smoothed = FALSE, confidence.level = 0.68)Returns the growth rate of the incidence () of the cumulated
variable (). Computed as
smoothed Logical value indicating whether to use the
smoothed estimates of and to compute the
growth rate (TRUE), or the contemporaneous filtered estimates
(FALSE). Default is FALSE.
confidence.level Confidence level for the confidence
interval. Default is , which is one standard deviation for
a normally distributed random variable.
idx_series object containing smoothed/filtered
growth rates and upper and lower bounds for the confidence intervals.
initialize(
data,
xpred_logical,
index,
reinit.idx,
ar1,
output,
sea.period,
xpred.new = NULL,
calendar = NULL
)Create an instance of the FilterResults class with fields defined
earlier in the fields section.
mapes(n.ahead, Y)Computes five metrics, including Mean Absolute Percentage Error (MAPE), for forecasts against a holdout sample. For more details, please refer to mapes.
sMAPE is computed as mean(100 * abs(Actual - Forecast) / (Actual + Forecast)) i.e. this uses (Actual + Forecast), not the more common (|Actual| + |Forecast|) / 2, in the denominator. As a result this statistic ranges over [0, 100], NOT the [0, 200] range associated with the textbook sMAPE definition that divides by the mean of the absolute values. Do not compare this sMAPE value directly against sMAPE figures computed with the conventional 0-200 formula.
predict_all(
n.ahead,
sea.on = TRUE,
return.all = FALSE,
confidence.level = 0.68
)Returns forecasts of the incidence variable , the state variables
and the conditional covariance matrix
for the states.
n.ahead The number of forecasts you wish to create from
the end of your sample period.
sea.on Logical value indicating whether seasonal
components should be included in the state-space model or not. Default is TRUE.
return.all Logical value indicating whether to return
all filtered estimates and forecasts.
(TRUE) or only the forecasts (FALSE). Default is
FALSE.
confidence.level The confidence level for the log growth
rate that should be used to compute. Confidence intervals only reported
for the incidence variable .
idx_series objects containing the forecast
(and filtered, where applicable) level
of (y.hat), (level.t.t),
(slope.t.t), vector of states including the
seasonals where applicable (a.t.t) and covariance matrix of all
states including seasonals where applicable (P.t.t).
predict_level(
n.ahead,
confidence.level = 0.68,
sea.on = TRUE,
return.diff = TRUE
)Forecast the cumulated variable or the incidence of it. This function returns
the forecast of the cumulated variable , or the forecast of the incidence of the cumulated variable, . For
example, in the case of an epidemic, might be daily new cases of
the disease and
the cumulative number of recorded infections.
n.ahead The number of periods ahead you wish to forecast from
the end of the estimation window.
confidence.level The confidence level for the log growth
rate that should be used to compute
the forecast intervals of .
sea.on Logical value indicating whether to return the prediction
of just the trend or prediction incorporating seasonality. Deafults to TRUE.
return.diff Logical value indicating whether to return the cumulated variable,
, or the incidence of it,
(i.e., the first difference of the cumulated variable). Default is
TRUE.
idx_series object containing the point
forecasts and upper and lower bounds of
the forecast interval.
print()Provides a quick glimpse of model states and standard errors.
print_estimation_results()Prints a table of estimated parameters in a format ready to paste into LaTeX.
summary()Supplies details of the FilterResults object, such as estimated parameter values, start and end positions of estimation.
Harvey, A. C. and Kattuman, P. (2021). A Farewell to R: Time Series Models for Tracking and Forecasting Epidemics, Journal of the Royal Society Interface, vol 18(182): 20210179
library(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) # Specify a model model <- SSModelDynamicGompertz$new(Y = Y, q = 0.005, end = 100) # Estimate a specified model res <- estimate(model) # Show summary of object summary(res) # Print a short description of the object print(res) # Print estimation results res$print_estimation_results() # Forecast 7 periods ahead from the end of the estimation window res$predict_level(n.ahead = 7, confidence.level = 0.68, sea.on=TRUE) # Forecast 7 periods ahead from the model and return filtered states res$predict_all(n.ahead = 7, return.all = TRUE) # Return the filtered growth rate and its components res$get_growth_y(return.components = TRUE) # Return smoothed growth rate of incidence variable and its confidence # interval res$get_gy_ci(smoothed = TRUE, confidence.level = 0.68) # Return MAPE of forecast res$mapes(n.ahead=7,Y)library(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) # Specify a model model <- SSModelDynamicGompertz$new(Y = Y, q = 0.005, end = 100) # Estimate a specified model res <- estimate(model) # Show summary of object summary(res) # Print a short description of the object print(res) # Print estimation results res$print_estimation_results() # Forecast 7 periods ahead from the end of the estimation window res$predict_level(n.ahead = 7, confidence.level = 0.68, sea.on=TRUE) # Forecast 7 periods ahead from the model and return filtered states res$predict_all(n.ahead = 7, return.all = TRUE) # Return the filtered growth rate and its components res$get_growth_y(return.components = TRUE) # Return smoothed growth rate of incidence variable and its confidence # interval res$get_gy_ci(smoothed = TRUE, confidence.level = 0.68) # Return MAPE of forecast res$mapes(n.ahead=7,Y)
A class that holds the information of an estimated Leading Indicator model.
Contains methods to extract smoothed/filtered estimates of the states, the
level of the incidence variable , and forecasts of . The output from the estimate method
of the SSModelLeadingIndicator class is of the class FilterResultsLI.
dataAn idx_series object with cumulated variables: lagged leading
indicator and target variable (plus their increments and log-growth rates - see
df2ldl_lead).
outputA KFS results object obtained after fitting a
SSModelLeadingIndicator model.
n.lagNumber of integer positions the leading indicator is lagged by, inherited
from the estimated SSModelLeadingIndicator model.
sea.periodThe period of seasonality, inherited from the estimated
SSModelLeadingIndicator model. For a day-of-the-week
effect with daily data, this would be 7.
LeadIndColThe column in data that contains the leading
indicator, inherited from the estimated SSModelLeadingIndicator model.
xpred_logicalVector of length 2 with logical values, indicating whether there are exogenous predictors for leading series and target series respectively.
xpred_lead.newAn idx_series object containing the values of exogenous
variables for the leading indicator over the prediction time frame.
xpred_targ.newAn idx_series object containing the values of exogenous
variables for the target variable over the prediction time frame.
startInteger position marking the start of the estimation period.
endInteger position marking the end of the estimation period.
calendarAn optional idx_calendar object, inherited
unchanged from the SSModelLeadingIndicator model this was
estimated from. Purely cosmetic: describes how the integer positions in
data map back to calendar time, for use by plotting only.
Defaults to NULL.
get_growth_y(smoothed = FALSE, return.components = FALSE)Returns the growth rate of the incidence () of the cumulated
variable (). Computed as
smoothed Logical value indicating whether to use the
smoothed estimates of and to compute the
growth rate (TRUE), or the contemporaneous filtered estimates
(FALSE). Default is FALSE.
return.components Logical value indicating whether to
return the estimates of and as well as
the estimates of the growth rate, or just the growth rate. Default is
FALSE.
idx_series objects containing
smoothed/filtered growth rates and components ( and
), where applicable.
get_gy_ci(smoothed = FALSE, confidence.level = 0.68)Returns the growth rate of the incidence () of the cumulated
variable (). Computed as
smoothed Logical value indicating whether to use the
smoothed estimates of and to compute the
growth rate (TRUE), or the contemporaneous filtered estimates
(FALSE). Default is FALSE.
confidence.level Confidence level for the confidence
interval. Default is , which is one standard deviation for
a normally distributed random variable.
idx_series object containing smoothed/filtered
growth rates and upper and lower bounds for the confidence intervals.
initialize(
data,
output,
n.lag,
sea.period,
LeadIndCol,
xpred_logical,
start,
end,
xpred_lead.new = NULL,
xpred_targ.new = NULL,
calendar = NULL
)Create an instance of the FilterResultsLI class with fields defined
earlier in the fields section.
mapes(n.ahead, Y)Computes five metrics, including Mean Absolute Percentage Error (MAPE), for forecasts against a holdout sample. For more details, please refer to mapes.
sMAPE is computed as mean(100 * abs(Actual - Forecast) / (Actual + Forecast)) i.e. this uses (Actual + Forecast), not the more common (|Actual| + |Forecast|) / 2, in the denominator. As a result this statistic ranges over [0, 100], NOT the [0, 200] range associated with the textbook sMAPE definition that divides by the mean of the absolute values. Do not compare this sMAPE value directly against sMAPE figures computed with the conventional 0-200 formula.
predict_all(
n.ahead,
sea.on = TRUE,
return.all = FALSE,
confidence.level = 0.68
)Returns forecasts of the incidence variable , the state variables
and the conditional covariance matrix for the states.
n.ahead The number of forecasts you wish to create from
the end of your sample period.
sea.on Logical value indicating whether to return the prediction
of just the trend or prediction incorporating seasonality. Default is TRUE.
return.all Logical value indicating whether to return
all filtered estimates and forecasts
(TRUE) or only the forecasts (FALSE). Default is
FALSE.
idx_series objects containing the forecast
(and filtered, where applicable) level
of (y.hat), (level.t.t),
(slope.t.t), vector of states including the
seasonals where applicable (a.t.t) and covariance matrix of all
states including seasonals where applicable (P.t.t).
predict_level(n.ahead = n.lag, confidence.level = 0.68, sea.on = TRUE)Forecast the cumulated variable or the incidence of it. This function returns
the forecast of the cumulated variable , or the forecast of the incidence of the cumulated variable, . For
example, in the case of an epidemic, might be daily new cases of
the disease and
n.ahead The number of periods ahead you wish to forecast from
the end of the estimation window. Default is n.lag.
sea.on Logical value indicating whether to return the prediction
of just the trend or prediction incorporating seasonality. Deafults to TRUE.
confidence.level The confidence level for the log growth
rate that should be used to compute the forecast intervals of .
An idx_series object containing the point
forecasts and upper and lower bounds of the forecast interval.
print()Provides a quick glimpse of model states and standard errors.
print_estimation_results()Prints a table of estimated parameters in a format ready to paste into LaTeX.
summary()Supplies details of the FilterResultsLI object, such as estimated parameter values, start and end positions of estimation.
Harvey, A. (2021). TIME SERIES MODELLING OF EPIDEMICS: LEADING INDICATORS, CONTROL GROUPS AND POLICY ASSESSMENT. National Institute Economic Review, 257, 83-100. doi:10.1017/nie.2021.21
library(tsgc) set.seed(1) lead <- cumsum(rpois(150, 6)) + 1 targ <- cumsum(rpois(150, 8)) + 1 Y <- idx_series(cbind(lead, targ), start = 1) # Define and estimate the model model <- SSModelLeadingIndicator(Y = Y, n.lag = 5, q = NULL, LeadIndCol = 1, sea.period = 0, start = 1, end = 100) res <- estimate(model) # Print estimation results res$print_estimation_results() # Forecast 7 positions ahead from the end of the estimation window res$predict_level(n.ahead = 7, confidence.level = 0.68) # Forecast 7 positions ahead from the model and return filtered states res$predict_all(n.ahead = 7, return.all = TRUE) # Return the filtered growth rate and its components res$get_growth_y(return.components = TRUE) # Return smoothed growth rate of incidence variable and its confidence # interval res$get_gy_ci(smoothed = TRUE, confidence.level = 0.68) # Return MAPE of forecast res$mapes(n.ahead = 7, Y = Y)library(tsgc) set.seed(1) lead <- cumsum(rpois(150, 6)) + 1 targ <- cumsum(rpois(150, 8)) + 1 Y <- idx_series(cbind(lead, targ), start = 1) # Define and estimate the model model <- SSModelLeadingIndicator(Y = Y, n.lag = 5, q = NULL, LeadIndCol = 1, sea.period = 0, start = 1, end = 100) res <- estimate(model) # Print estimation results res$print_estimation_results() # Forecast 7 positions ahead from the end of the estimation window res$predict_level(n.ahead = 7, confidence.level = 0.68) # Forecast 7 positions ahead from the model and return filtered states res$predict_all(n.ahead = 7, return.all = TRUE) # Return the filtered growth rate and its components res$get_growth_y(return.components = TRUE) # Return smoothed growth rate of incidence variable and its confidence # interval res$get_gy_ci(smoothed = TRUE, confidence.level = 0.68) # Return MAPE of forecast res$mapes(n.ahead = 7, Y = Y)
Cumulative daily cases of Covid-19 in the South African province of Gauteng.
data(gauteng)data(gauteng)
An object of class "xts" with columns:
Cumulative cases of Covid-19 from 10th March 2020
Downloaded from https://sacoronavirus.co.za/
data(gauteng) # plot daily cases plot(diff(gauteng))data(gauteng) # plot daily cases plot(diff(gauteng))
Daily weather data of the South African province of Gauteng in 2021
data(gauteng_weather_2021)data(gauteng_weather_2021)
An object of class "xts" with columns:
Wind speed in metres per second (m/s) at 10 metres
Relative humidity at 2 metres (%)
Daily mean air temperature in degree Celsius (C) at 2 metres
Daily maximum temperature in degree Celsius (C) at 2 metres
Precipitation (millimetres per day; mm/day)
Downloaded from https://cds.climate.copernicus.eu/datasets/reanalysis-era5-single-levels?tab=overview and aggregated across the latitude and longitude corners for Gauteng province, South Africa by the authors.
data(gauteng_weather_2021)data(gauteng_weather_2021)
idx_series objects given start and end positionsHelper method to subset an idx_series for a
specified range of integer positions.
get_timeframe(df, start, end = NULL)get_timeframe(df, start, end = NULL)
df |
An |
start |
Start position (integer) of the range. |
end |
End position (integer) of the range. Defaults to the last
position in |
An idx_series containing the selected observations, or
NULL if df is NULL.
x <- idx_series(cumsum(rpois(30, 5)) + 1, start = 1) get_timeframe(x, 5, 10) get_timeframe(x, 5)x <- idx_series(cumsum(rpois(30, 5)) + 1, start = 1) get_timeframe(x, 5, 10) get_timeframe(x, 5)
Accessor method to access the non-diffuse part of the error
covariance matrix of the smoothed states in a KFS object.
get_V(object)get_V(object)
object |
A |
The error covariance matrix of the smoothed states.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) get_V(output(res))library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) get_V(output(res))
Accessor method to access the time series y in an
SSModel object.
gety(object)gety(object)
object |
An |
The time series y.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::gety(modelKFS(output(res)))library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::gety(modelKFS(output(res)))
Accessor method to access the prediction y.hat in the
output of predict_all.
gety.hat(object)gety.hat(object)
object |
An object returned by |
The predicted values y.hat.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) all_predictions <- res$predict_all(n.ahead = 7) tsgc:::gety.hat(all_predictions)library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) all_predictions <- res$predict_all(n.ahead = 7) tsgc:::gety.hat(all_predictions)
idx_series
Returns the first n observations of x, in
position order. Unlike indexing x directly with [, n
is a count of observations (as for utils::head), not a set of
absolute positions, so head(x, 6) always returns the first six
observations regardless of x$start.
## S3 method for class 'idx_series' head(x, n = 6L, ...)## S3 method for class 'idx_series' head(x, n = 6L, ...)
x |
An |
n |
Number of observations to return from the start of |
... |
Unused. |
An idx_series of length min(n, length(x)).
idx_series plots is
displayedBundles the small set of cosmetic choices that decide how
plotting functions in this file turn integer idx_series positions
(and an optional idx_calendar) into an x-axis, and whether to
attach a small info box summarising the calendar. Pass an
idx_axis_opts object via the axis argument of the various
plot.*/plot_* functions; NULL (the default
everywhere) reproduces the previous behaviour exactly.
idx_axis_opts( mode = c("auto", "position", "steps", "time_since", "date"), info_box = FALSE, pattern_n = NULL )idx_axis_opts( mode = c("auto", "position", "steps", "time_since", "date"), info_box = FALSE, pattern_n = NULL )
mode |
One of:
|
info_box |
A single logical. If |
pattern_n |
|
An object of class idx_axis_opts.
idx_axis_opts(mode = "steps", info_box = TRUE) idx_axis_opts(mode = "time_since", info_box = TRUE, pattern_n = 5)idx_axis_opts(mode = "steps", info_box = TRUE) idx_axis_opts(mode = "time_since", info_box = TRUE, pattern_n = 5)
idx_calendar
idx_calendar describes how to translate the integer
positions used by idx_series back into calendar time. It
holds no data of its own and is not used by the statistical engine of
this package (estimation, filtering, forecasting); it exists for use at
the edges of the package - plotting, printing, reporting - where a human
needs to see a date rather than an integer position.
An idx_calendar is defined by:
anchor: a reference point in calendar time (e.g. a
Date, or a plain number if the series has no true calendar
meaning - e.g. picoseconds since some experiment start).
anchor_pos: the integer position (in the same units as
idx_series$start) that lines up with anchor.
amount/unit: the size of a single step between
consecutive integer positions, e.g. amount = 2.5, unit = "hours"
or amount = 1, unit = "days". amount may be non-integer.
pattern: a repeating vector of step multipliers, used
when consecutive integer positions do not correspond to evenly spaced
calendar steps. For example, business days (weekends skipped) would
be pattern = c(1, 1, 1, 1, 3): four single-day steps, Mon-Thu,
then a 3-day step from Fri to the following Mon. Defaults to
1, i.e. every step is a single, uniform amount.
pattern_start: the index into pattern that
anchor_pos sits at, i.e. which position in the repeating cycle
the anchor falls on. Defaults to 1.
This object deliberately does no calendar arithmetic itself beyond plain
numeric offsets (number of steps x amount, in units). Any
interpretation of that numeric offset as a true calendar date/time
(including POSIXct-specific concerns such as DST or variable month
lengths) is left to the plotting/translation layer, which is expected to
use idx_calendar_offset() together with the posixct flag on
this object to decide how to turn the offset into an actual date.
idx_calendar( anchor, anchor_pos = 1L, amount = 1, unit = "days", pattern = 1, pattern_start = 1L, posixct = FALSE, anchor_name = NULL )idx_calendar( anchor, anchor_pos = 1L, amount = 1, unit = "days", pattern = 1, pattern_start = 1L, posixct = FALSE, anchor_name = NULL )
anchor |
A single reference point in calendar time. Typically a
|
anchor_pos |
A single integer: the |
amount |
A single positive number: the size of one step between
consecutive integer positions, in |
unit |
A single string naming the unit that |
pattern |
A numeric vector of step multipliers, recycled
cyclically across consecutive integer positions. Defaults to |
pattern_start |
A single positive integer: the index into
|
posixct |
A single logical. If |
anchor_name |
An optional single string giving a human-readable
name for |
An object of class idx_calendar.
# Daily series, anchored so that position 1 is 2024-01-01. cal <- idx_calendar(anchor = as.Date("2024-01-01"), anchor_pos = 1L, amount = 1, unit = "days", posixct = TRUE) idx_to_date(cal, 1:5) # Business-day series (weekends skipped), anchor at a Monday. cal_bd <- idx_calendar(anchor = as.Date("2024-01-01"), anchor_pos = 1L, amount = 1, unit = "days", pattern = c(1, 1, 1, 1, 3), pattern_start = 1L, posixct = TRUE) idx_to_date(cal_bd, 1:10) # Non-calendar series: picosecond timesteps, arbitrary numeric anchor. cal_ps <- idx_calendar(anchor = 0, anchor_pos = 1L, amount = 2.5, unit = "picoseconds") idx_to_date(cal_ps, 1:5)# Daily series, anchored so that position 1 is 2024-01-01. cal <- idx_calendar(anchor = as.Date("2024-01-01"), anchor_pos = 1L, amount = 1, unit = "days", posixct = TRUE) idx_to_date(cal, 1:5) # Business-day series (weekends skipped), anchor at a Monday. cal_bd <- idx_calendar(anchor = as.Date("2024-01-01"), anchor_pos = 1L, amount = 1, unit = "days", pattern = c(1, 1, 1, 1, 3), pattern_start = 1L, posixct = TRUE) idx_to_date(cal_bd, 1:10) # Non-calendar series: picosecond timesteps, arbitrary numeric anchor. cal_ps <- idx_calendar(anchor = 0, anchor_pos = 1L, amount = 2.5, unit = "picoseconds") idx_to_date(cal_ps, 1:5)
idx_calendar from a multi_step_pattern
Third, most general constructor for
idx_calendar. Use this when consecutive positions cycle
through genuinely different kinds of steps - not just different
multiples of one step (that case is covered by the numeric
pattern argument of idx_calendar/
idx_calendar_step). For example, "3 days, 3 days, then 1
month, repeating" needs a multi_step_pattern since two
slots are fixed-duration and one is calendar-relative.
Unlike the numeric pattern mechanism (which pre-computes a
closed-form pattern-weighted step count, see
idx_calendar_offset), a multi_step_pattern is
applied by walking one slot at a time from anchor_pos to the
target position, via idx_step_add. This is required
because heterogeneous steps (e.g. days mixed with months) cannot be
collapsed into a single scalar multiplier the way a uniform-unit
pattern can - so this constructor's calendars are O(distance from
anchor) to convert per position, rather than O(1).
idx_calendar_multi_step( anchor, multi_step, anchor_pos = 1L, pattern_start = 1L, posixct = FALSE, anchor_name = NULL )idx_calendar_multi_step( anchor, multi_step, anchor_pos = 1L, pattern_start = 1L, posixct = FALSE, anchor_name = NULL )
anchor |
A single reference point in calendar time (see
|
multi_step |
A |
anchor_pos, posixct, anchor_name
|
As in |
pattern_start |
A single positive integer: the index into
|
An object of class idx_calendar.
# 3 days, 3 days, then 1 month, repeating. cal <- idx_calendar_multi_step( anchor = as.Date("2024-01-01"), multi_step = multi_step_pattern(idx_step(days = 3), idx_step(days = 3), idx_step(months = 1)), posixct = TRUE ) idx_to_date(cal, 1:5)# 3 days, 3 days, then 1 month, repeating. cal <- idx_calendar_multi_step( anchor = as.Date("2024-01-01"), multi_step = multi_step_pattern(idx_step(days = 3), idx_step(days = 3), idx_step(months = 1)), posixct = TRUE ) idx_to_date(cal, 1:5)
Computes, for one or more integer idx_series
positions, the number of base amount-sized steps that separate
each position from cal$anchor_pos, accounting for the repeating
pattern of step multipliers. This is a pure integer/numeric
computation - no calendar arithmetic and no units are applied here; the
result is in multiples of cal$amount, positive for positions
after the anchor and negative for positions before it. It is intended to
be consumed by idx_to_date (or directly by the plotting
layer's own POSIXct-aware conversion).
idx_calendar_offset(cal, pos)idx_calendar_offset(cal, pos)
cal |
An |
pos |
An integer vector of |
A numeric vector, the same length as pos, giving the
pattern-weighted step offset (in units of cal$amount) of each
position relative to cal$anchor_pos.
idx_calendar from a compound idx_step
Secondary constructor for idx_calendar,
taking a compound idx_step (a combination of years,
quarters, months, weeks, days, hours, minutes, seconds) as the step
size instead of a single amount/unit pair. Use this when
a single unit cannot express the step - e.g. "3 seconds after each
quarter" or "a month and 13 seconds". For a plain single-unit step, the
primary idx_calendar constructor is simpler.
idx_calendar_step( anchor, step, anchor_pos = 1L, pattern = 1, pattern_start = 1L, posixct = FALSE, anchor_name = NULL )idx_calendar_step( anchor, step, anchor_pos = 1L, pattern = 1, pattern_start = 1L, posixct = FALSE, anchor_name = NULL )
anchor |
A single reference point in calendar time (see
|
step |
An |
anchor_pos, pattern, pattern_start, posixct, anchor_name
|
As in
|
An object of class idx_calendar.
# Data recorded 3 seconds after each quarter boundary. cal <- idx_calendar_step( anchor = as.POSIXct("2024-01-01 00:00:03", tz = "UTC"), step = idx_step(quarters = 1, seconds = 3), posixct = TRUE ) idx_to_date(cal, 1:4)# Data recorded 3 seconds after each quarter boundary. cal <- idx_calendar_step( anchor = as.POSIXct("2024-01-01 00:00:03", tz = "UTC"), step = idx_step(quarters = 1, seconds = 3), posixct = TRUE ) idx_to_date(cal, 1:4)
idx_series objectsColumn-binds one or more idx_series objects that
share the same start and length, analogous to cbind for
matrices/xts objects.
idx_cbind(...)idx_cbind(...)
... |
|
A single idx_series with combined columns.
idx_series
Index-based analogue of diff() for calendar series:
returns for each position , dropping the
first lag positions (which have no valid lag within the series).
idx_diff(x, lag = 1L)idx_diff(x, lag = 1L)
x |
An |
lag |
Integer number of positions to lag by. Default |
An idx_series, shorter than x by lag.
idx_series
Index-based analogue of stats::lag(): shifts the
series forward by k positions. The returned series still covers
length(x) observations, but its start is shifted by
k, i.e. position of the result holds the value that
was at position in x.
idx_lag(x, k = 1L)idx_lag(x, k = 1L)
x |
An |
k |
Integer number of positions to lag by. Default |
An idx_series with the same length as x.
idx_series
base::NCOL() is not an S3 generic, so it will not
dispatch on idx_series objects. Use this function instead when
you need the number of columns/variables in an idx_series.
idx_ncol(x)idx_ncol(x)
x |
An |
idx_series integer
positionInverse of idx_to_date's plain-arithmetic
branch (posixct = FALSE, or a non-Date/POSIXct cal$anchor):
converts a value already expressed in cal$anchor's own units/scale
(e.g. a picosecond count) into the integer idx_series position
that corresponds to it.
Only relevant when you actually have an idx_calendar with a
non-Date/POSIXct anchor (e.g. a numeric anchor for sub-second or
otherwise non-calendar timesteps) and want to look up a position from a
value in that anchor's scale. If your series has no calendar meaning at
all - no dates, no anchor, nothing to translate - you don't need an
idx_calendar or this function: idx_series positions
(start/end arguments throughout this package) are already
plain integers, so just use them directly (see idx_range).
idx_offset_to_pos(cal, value)idx_offset_to_pos(cal, value)
cal |
An |
value |
A single number, in the same units as |
An integer (vector, matching value): the idx_series
position(s) corresponding to value.
cal <- idx_calendar(anchor = 0, anchor_pos = 1L, amount = 2.5, unit = "picoseconds") idx_offset_to_pos(cal, 12.5)cal <- idx_calendar(anchor = 0, anchor_pos = 1L, amount = 2.5, unit = "picoseconds") idx_offset_to_pos(cal, 12.5)
idx_series
Returns the full vector of integer positions associated
with each observation in x. This is the index-based analogue of
zoo::index().
idx_positions(x)idx_positions(x)
x |
An |
An integer vector.
idx_series
Returns the first and last integer position covered by
x, i.e. c(x$start, x$start + length(x) - 1).
idx_range(x)idx_range(x)
x |
An |
An integer vector of length 2: c(first, last).
idx_series sequentiallyConcatenates x then y in position order.
y must start exactly one position after x ends.
idx_rbind(x, y)idx_rbind(x, y)
x, y
|
|
A single, combined idx_series.
idx_series is the core data structure used throughout
the analysis parts of this package in place of calendar-indexed objects
(such as xts). A series is stored as a plain numeric vector or
matrix of observations together with a single integer start,
which is the integer position (relative to some external anchor point,
e.g. the first observation ever collected) of the first row of
data. All internal computation - estimation, filtering,
forecasting, reinitialisation - is done purely in terms of these integer
positions. This keeps the statistical engine of the package completely
free of calendar/time-of-day concerns.
Calendar time is reintroduced only at the very edges of the package (e.g. plotting, printing to the user) by a separate, cosmetic translation layer that maps integer positions back to dates given an anchor date. That translation layer is out of scope here.
idx_series(data, start = 1L)idx_series(data, start = 1L)
data |
A numeric vector or matrix of observations. For a matrix, rows are observations (in position order) and columns are variables. |
start |
A single positive integer giving the position of the first
observation in |
An object of class idx_series.
x <- idx_series(cumsum(rpois(50, 5)) + 1, start = 1) length(x) idx_range(x) x[10:20]x <- idx_series(cumsum(rpois(50, 5)) + 1, start = 1) length(x) idx_range(x) x[10:20]
idx_step describes the size of a single step between
consecutive idx_series integer positions as a combination of
calendar-relative and fixed-duration components: years, quarters,
months, weeks, days, hours, minutes, and seconds. Each component is
independent - they are not normalized or collapsed into one another
(e.g. idx_step(months = 1, days = 31) is not simplified to two
months; it is applied as "one month, then 31 days"). This lets a step
express things a single amount/unit pair cannot, such as
"3 seconds after each quarter" (idx_step(quarters = 1, seconds = 3))
or "a month and 13 seconds" (idx_step(months = 1, seconds = 13)).
When applied (see idx_step_add), components are added to
an anchor date in a fixed order, from largest calendar-relative unit to
smallest fixed-duration unit: years, then quarters, then months (each
calendar-aware, i.e. respecting variable month/quarter length), then
weeks, days, hours, minutes, seconds (each fixed-duration). This order
is applied once per full step multiple - i.e. for n steps from
the anchor, every component is scaled by n before being added,
not added once and then repeated.
A step with any of hours, minutes, or seconds
non-zero requires a POSIXct anchor (not a plain Date or
yearqtr/yearmon), since those units need a time-of-day
component to be meaningful.
idx_step( years = 0, quarters = 0, months = 0, weeks = 0, days = 0, hours = 0, minutes = 0, seconds = 0 )idx_step( years = 0, quarters = 0, months = 0, weeks = 0, days = 0, hours = 0, minutes = 0, seconds = 0 )
years, quarters, months, weeks, days, hours, minutes, seconds
|
Single
non-negative, finite numbers giving the count of each component in one
step. All default to |
An object of class idx_step.
# A single quarter, nudged 3 seconds later - e.g. data recorded 3 # seconds after each quarter boundary. idx_step(quarters = 1, seconds = 3) # A month and 13 seconds - the extra seconds accumulate with each step, # same as the calendar-relative part. idx_step(months = 1, seconds = 13) # Equivalent to amount = 1, unit = "days". idx_step(days = 1)# A single quarter, nudged 3 seconds later - e.g. data recorded 3 # seconds after each quarter boundary. idx_step(quarters = 1, seconds = 3) # A month and 13 seconds - the extra seconds accumulate with each step, # same as the calendar-relative part. idx_step(months = 1, seconds = 13) # Equivalent to amount = 1, unit = "days". idx_step(days = 1)
idx_step to an anchor date, scaled by a step countAdds n whole steps of step to anchor,
one component at a time, from largest calendar-relative unit to
smallest fixed-duration unit (years, quarters, months, weeks, days,
hours, minutes, seconds), each scaled by n before being added.
Calendar-relative components (years/quarters/months) use calendar-aware
stepping (seq.Date/seq.POSIXt, or the
equivalent fractional-year arithmetic for yearqtr/yearmon
anchors); fixed-duration components (weeks/days/hours/minutes/seconds)
use plain duration addition.
idx_step_add(anchor, step, n)idx_step_add(anchor, step, n)
anchor |
A |
step |
An |
n |
An integer (or vector of integers): the number of whole steps to apply. May be negative. |
A vector the same length as n, the same class as
anchor (or POSIXct, if a sub-day component forces
promotion from Date).
idx_series positions to calendar timeConverts one or more integer idx_series positions
into calendar time using an idx_calendar.
When cal$posixct is TRUE and cal$anchor is a
Date/POSIXct, and cal$unit is one of the standard
calendar units - "seconds", "minutes", "hours",
"days", "weeks", "months", "quarters",
"years" (singular or plural) - calendar-aware stepping is used via
seq.Date/seq.POSIXt's by argument.
This is exact for every one of these units: seconds through weeks are
fixed-duration and so plain offset arithmetic would already agree with
it, but months/quarters/years are calendar-relative (variable length -
e.g. not every month has the same number of days) and genuinely require
this calendar-aware stepping to land on the correct date (respecting
month-end/DST rules) rather than a fixed number of days.
cal$anchor may also be a zoo yearqtr or
yearmon, for series that are naturally indexed at quarterly or
monthly resolution rather than daily.
Internally, every idx_calendar (whether built via the primary
amount/unit constructor or via
idx_calendar_step) carries a compound cal$step
(see idx_step), and calendar-aware stepping is delegated
uniformly to idx_step_add. This is what allows steps that
mix calendar-relative and fixed-duration components, e.g. "3 seconds
after each quarter" or "a month and 13 seconds" - a single
amount/unit pair cannot express these, but a compound
idx_step can.
Otherwise (posixct = FALSE, or a non-Date/POSIXct/yearqtr/yearmon
anchor) the conversion falls back to plain numeric offset arithmetic: it
multiplies the pattern-weighted step offset (see
idx_calendar_offset) by cal$amount and adds the
result to cal$anchor. This supports arbitrary non-calendar units
(e.g. picoseconds) uniformly, but is not available for calendars built
via idx_calendar_step with posixct = FALSE, since
those have no single amount to multiply by (see
idx_calendar_step); such calendars require
posixct = TRUE.
idx_to_date(cal, pos)idx_to_date(cal, pos)
cal |
An |
pos |
An integer vector of |
A vector, the same length as pos, giving the calendar
time (or plain number, for non-calendar cal$anchor) for each
position.
idx_series integer positionInverse of idx_to_date: converts a calendar
date into the integer position (relative to cal) that corresponds
to it. Intended for locating estimation start/end dates, or for aligning
a second series to an already-converted idx_series (see
xts_to_idx). Requires a calendar-anchored cal
(cal$posixct = TRUE and a Date/POSIXct/
yearqtr/yearmon cal$anchor); use
idx_offset_to_pos instead for non-calendar (e.g.
sub-second/arbitrary numeric) calendars.
cal may be built via any of the three idx_calendar
constructors:
The primary idx_calendar constructor (a single
amount/unit pair). For the calendar-relative units
"months", "quarters", "years", the position is
found by counting whole calendar steps between cal$anchor and
date (the inverse of the calendar-aware stepping used by
idx_to_date), rather than via a fixed-length difference -
a quarter is not a fixed number of days. For the fixed-duration units
"seconds", "minutes", "hours", "days",
"weeks", positions are computed as a plain time difference
scaled by cal$amount; this is exact for these units. A repeating
numeric pattern is respected. Because patterned offsets are
piecewise rather than a single amount that can be divided in closed
form, the inverse is found by an exact binary search over integer
positions. Dates that do not land on a whole step boundary are rejected
rather than rounded.
idx_calendar_step (a compound idx_step) and
idx_calendar_multi_step (a multi_step_pattern)
calendars have no single amount/unit pair to invert in
closed form in general, so date is instead located by search:
idx_calendar_step calendars repeat one compound idx_step
every position, and idx_step_add is strictly monotonic
in the step count, so this uses an O(log distance) binary search;
idx_calendar_multi_step calendars cycle through heterogeneous
steps that cannot be searched this way, so this walks one step at a
time from cal$anchor_pos instead, an O(distance) linear walk.
Both are exact, and both ignore any numeric pattern (a
multi_step_pattern calendar has no numeric pattern to
ignore in the first place - it uses pattern_start only as the
multi_step slot index).
A yearqtr/yearmon cal$anchor is only valid together
with cal$unit of "quarters"/"months" respectively
(matching idx_step_add's restriction that these anchors
have no sub-quarter/sub-month resolution); date is coerced to the
same class as cal$anchor and whole-step boundaries are found via
fractional-year arithmetic, consistent with idx_to_date.
idx_to_pos(cal, date)idx_to_pos(cal, date)
cal |
An |
date |
A single date/time. For |
A single integer: the idx_series position corresponding
to date.
cal <- idx_calendar(anchor = as.Date("2024-01-01"), anchor_pos = 1L, amount = 1, unit = "days", posixct = TRUE) idx_to_pos(cal, "2024-01-10") # A multi_step_pattern calendar (walked, not closed-form). cal_ms <- idx_calendar_multi_step( anchor = as.Date("2024-01-01"), multi_step = multi_step_pattern(idx_step(days = 3), idx_step(days = 3), idx_step(months = 1)), posixct = TRUE ) idx_to_pos(cal_ms, "2024-02-07")cal <- idx_calendar(anchor = as.Date("2024-01-01"), anchor_pos = 1L, amount = 1, unit = "days", posixct = TRUE) idx_to_pos(cal, "2024-01-10") # A multi_step_pattern calendar (walked, not closed-form). cal_ms <- idx_calendar_multi_step( anchor = as.Date("2024-01-01"), multi_step = multi_step_pattern(idx_step(days = 3), idx_step(days = 3), idx_step(months = 1)), posixct = TRUE ) idx_to_pos(cal_ms, "2024-02-07")
idx_series
Portable, always-dispatches accessor for the raw numeric
data underlying an idx_series. Prefer this over as.numeric()
in package-internal code, since it does not depend on S3 dispatch of a
primitive generic.
idx_values(x)idx_values(x)
x |
An |
A plain numeric vector (or matrix, if x wraps a matrix).
idx_calendar
Test whether an object is an idx_calendar
is_idx_calendar(x)is_idx_calendar(x)
x |
Object to test. |
Logical.
idx_series
Test whether an object is an idx_series
is_idx_series(x)is_idx_series(x)
x |
Object to test. |
Logical.
idx_step
Test whether an object is an idx_step
is_idx_step(x)is_idx_step(x)
x |
Object to test. |
Logical.
multi_step_pattern
Test whether an object is a multi_step_pattern
is_multi_step_pattern(x)is_multi_step_pattern(x)
x |
Object to test. |
Logical.
idx_series
Number of observations in an idx_series
## S3 method for class 'idx_series' length(x, ...)## S3 method for class 'idx_series' length(x, ...)
x |
An |
... |
Unused. |
This is a helper function that calculates five error metrics of a forecast generated by time series growth curve (tsgc) models. It compares the forecasted values to a holdout sample, providing a measure of forecast accuracy.
mapes(res, n.ahead, Y)mapes(res, n.ahead, Y)
res |
A |
n.ahead |
Integer specifying the number of periods to forecast ahead. |
Y |
An |
A list containing five error metrics for the forecast, with element names
mape: mean absolute percentage error
smape: symmetric mean absolute percentage error (between 0 to 100)
mae: mean absolute error
rmse: root mean squared error
coverage: Percentage of holdout sample data points that lie inside the confidence interval for predictions
Accessor method to access the matrices used in the
observation, state, and disturbance equations of a KFS object.
matrixKFS(object, matrix)matrixKFS(object, matrix)
object |
A |
matrix |
Character string naming a matrix component of |
The requested matrix.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::matrixKFS(output(res), "Z")library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::matrixKFS(output(res), "Z")
Accessor method to access the fitted SSModel.
modelKFS(object)modelKFS(object)
object |
A |
The SSModel object underlying the KFS fit.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) modelKFS(output(res))library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) modelKFS(output(res))
idx_stepsmulti_step_pattern describes a repeating cycle
where each slot can be a completely different idx_step -
not just a different multiple of one fixed step. This is for cases the
numeric pattern argument of idx_calendar/
idx_calendar_step cannot express, because that mechanism
only ever scales a single amount/unit (or a single
compound idx_step) by different multipliers around the cycle -
it cannot mix step kinds. For example, "3 days, 3 days, then 1 month,
repeating" needs two of the three slots to be fixed-duration and one to
be calendar-relative; multi_step_pattern lets each slot carry
its own idx_step:
multi_step_pattern(idx_step(days = 3), idx_step(days = 3), idx_step(months = 1)).
Use this via idx_calendar_multi_step (a third,
most-general constructor for idx_calendar, alongside the
primary amount/unit constructor and the
idx_calendar_step compound-step constructor). The
existing numeric pattern mechanism is unchanged and remains the
simpler choice whenever every slot in the cycle is just a different
multiple of the same single step.
multi_step_pattern(...)multi_step_pattern(...)
... |
One or more |
An object of class multi_step_pattern: a list of
idx_step objects.
# 3 days, 3 days, then 1 month - repeating. multi_step_pattern(idx_step(days = 3), idx_step(days = 3), idx_step(months = 1))# 3 days, 3 days, then 1 month - repeating. multi_step_pattern(idx_step(days = 3), idx_step(days = 3), idx_step(months = 1))
idx_series (NCOL method)Provided for convenience when calling NCOL()
directly on an idx_series; note that since base::NCOL is
not an S3 generic, this method is only reached if NCOL is called
as NCOL.idx_series() directly, or after NCOL has been
made generic (e.g. via Matrix or similar). Prefer
idx_ncol for portable code within this package.
NCOL.idx_series(x)NCOL.idx_series(x)
x |
An |
Quarterly Sales of Nintendo consoles from 2004 Q4 to 2024 Q1
data(nintendo_sales)data(nintendo_sales)
An object of class "xts" with columns:
Cumulative Wii sales
Cumulative Wii U sales
Cumulative Switch sales (original, lite and OLED combined)
Cumulative original Switch sales
Cumulative Switch Lite sales
Cumulative switch OLED sales
Cumulative 3DS sales
Cumulative NDS sales
Taken from https://vgsales.fandom.com/wiki/Wii
data(nintendo_sales)data(nintendo_sales)
Accessor method to access the fitted KFS model from
FilterResults or FilterResultsLI.
output(object)output(object)
object |
A |
The fitted KFS model.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) output(res)library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) output(res)
Plots forecasts from a list of fitted models on the same axes for visual comparison.
plot_compare_forecast( results, n.ahead = 14, sea.on = TRUE, actual = NULL, title = "Comparison of forecasts", labels = NULL, axis = NULL )plot_compare_forecast( results, n.ahead = 14, sea.on = TRUE, actual = NULL, title = "Comparison of forecasts", labels = NULL, axis = NULL )
results |
A list of |
n.ahead |
Duration of the holdout sample. Default is 14. |
sea.on |
Logical value indicating whether to plot the
seasonality-adjusted forecasts. Defaults to |
actual |
An |
title |
Title for the plot. Defaults to |
labels |
Character vector of legend labels, one per element of
|
axis |
An |
A ggplot2 plot.
Plots actual values of the difference in the cumulated
variable, the forecast of the cumulated variable (including seasonal
components, where specified) and forecast intervals around the forecast.
The forecast intervals are based on the prediction intervals for
.
plot_forecast( res, n.ahead = 14, confidence.level = 0.68, title = NULL, plt.start = NULL, series.name = "target variable", axis = NULL )plot_forecast( res, n.ahead = 14, confidence.level = 0.68, title = NULL, plt.start = NULL, series.name = "target variable", axis = NULL )
res |
A |
n.ahead |
Number of forecasts (i.e. number of positions ahead to forecast from the end of the estimation window). Default is 14. |
confidence.level |
Width of the prediction interval for
|
title |
Title for the plot. |
plt.start |
First integer position of actual data (from the
estimation sample) to plot. |
series.name |
Name of the series the growth rate is being computed
for, e.g. |
axis |
An |
A ggplot2 plot.
library(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) model <- SSModelDynamicGompertz$new(Y = Y, q = 0.005, end = 100) res <- estimate(model) plot_forecast(res, n.ahead = 7, series.name = "cases")library(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) model <- SSModelDynamicGompertz$new(Y = Y, q = 0.005, end = 100) res <- estimate(model) plot_forecast(res, n.ahead = 7, series.name = "cases")
Plots the smoothed/filtered growth rate of the difference
in the cumulated variable () and the associated confidence
intervals.
plot_gy_ci( res, plt.start = NULL, smoothed = FALSE, title = NULL, series.name = NULL, pad.right = NULL, axis = NULL )plot_gy_ci( res, plt.start = NULL, smoothed = FALSE, title = NULL, series.name = NULL, pad.right = NULL, axis = NULL )
res |
A |
plt.start |
Plot start position. |
smoothed |
Logical value indicating whether to use the smoothed
estimates. Default |
title |
Title for the plot. |
series.name |
Name of the series the growth rate is being computed
for, e.g. |
pad.right |
Numerical value for the number of positions of blank space to leave on the right of the graph. |
axis |
An |
A ggplot2 plot.
Plots the smoothed/filtered growth rate of the difference
in the cumulated variable (), the smoothed/filtered growth rate
of the cumulated variable (), and the smoothed/filtered slope of
, . Following Harvey and Kattuman (2021), we
compute as
plot_gy_components( res, plt.start = NULL, smoothed = FALSE, title = NULL, axis = NULL )plot_gy_components( res, plt.start = NULL, smoothed = FALSE, title = NULL, axis = NULL )
res |
A |
plt.start |
Plot start position. |
smoothed |
Logical value indicating whether to use the smoothed
estimates of |
title |
Title for the plot. |
axis |
An |
A ggplot2 plot.
Plots actual values of the difference in the cumulated variable, the forecast of the cumulated variable (including seasonal components, where specified) and forecast intervals around the forecast, plus the actual outcomes from the holdout sample. Also reports the mean absolute percentage error over the holdout sample.
plot_holdout( res, Y, n.ahead = 14, confidence.level = 0.68, series.name = "target variable", title = NULL, caption = NULL, axis = NULL )plot_holdout( res, Y, n.ahead = 14, confidence.level = 0.68, series.name = "target variable", title = NULL, caption = NULL, axis = NULL )
res |
A |
Y |
Values of the cumulated variable, including the holdout sample. |
n.ahead |
Duration of the holdout sample. Default is 14. |
confidence.level |
Width of the prediction interval for
|
series.name |
Name of the variable being forecast, for the y-axis label. |
title |
Title for the plot. |
caption |
Caption for the plot. |
axis |
An |
A ggplot2 plot.
Plots actual and filtered values of the log cumulative
growth rate () in the estimation sample and the forecast
and realised log cumulative growth rate out of the estimation sample.
plot_log_forecast( res, Y, n.ahead = 14, plt.start = NULL, title = "", caption = "", axis = NULL )plot_log_forecast( res, Y, n.ahead = 14, plt.start = NULL, title = "", caption = "", axis = NULL )
res |
A |
Y |
Cumulated dataset ( |
n.ahead |
Number of positions ahead from the end of the sample to forecast. Default is 14. |
plt.start |
Plot start position. |
title |
Plot title. |
caption |
Plot caption. |
axis |
An |
A ggplot2 plot.
Plots , as computed by estimate_r0,
together with its confidence interval, over the most recent
n.ahead integer positions of res. A horizontal reference
line is drawn at , the threshold separating a growing from
a shrinking epidemic/process.
plot_r0( res, gen_int, n.ahead = 7, smoothed = FALSE, confidence.level = 0.68, title = NULL, axis = NULL )plot_r0( res, gen_int, n.ahead = 7, smoothed = FALSE, confidence.level = 0.68, title = NULL, axis = NULL )
res |
A |
gen_int |
The mean generation interval, in the same integer-position
units as |
n.ahead |
Number of most recent integer positions to plot, taken
from the end of |
smoothed |
Logical value indicating whether to use the smoothed
( |
confidence.level |
Confidence level for the confidence interval
around |
title |
Title for the plot. |
axis |
An |
A ggplot2 plot.
library(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) model <- SSModelDynamicGompertz$new(Y = Y, q = NULL, end = 100) res <- estimate(model) plot_r0(res, gen_int = 5, n.ahead = 7)library(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) model <- SSModelDynamicGompertz$new(Y = Y, q = NULL, end = 100) res <- estimate(model) plot_r0(res, gen_int = 5, n.ahead = 7)
Plots the lagged differences of the cumulated dataset
Y stored in a SSModelDynamicGompertz object, optionally
overlaid with a centred moving average.
## S3 method for class 'SSModelDynamicGompertz' plot( x, title = NULL, series.name = "target variable", MA_period = 7, axis = NULL, ... )## S3 method for class 'SSModelDynamicGompertz' plot( x, title = NULL, series.name = "target variable", MA_period = 7, axis = NULL, ... )
x |
A |
title |
Title for the plot. |
series.name |
Name of the series being plotted, used in the y-axis
label. Default is |
MA_period |
Number of positions in the centred moving average to
overlay. |
axis |
An |
... |
Unused. |
A ggplot2 plot.
Plots the daily incidence (new cases) of the leading
indicator and target variable stored in a SSModelLeadingIndicator
object, optionally on a log scale.
## S3 method for class 'SSModelLeadingIndicator' plot( x, title = NULL, series.name.lead = "Leading Indicator", series.name.target = "Target Variable", take.log = TRUE, axis = NULL, ... )## S3 method for class 'SSModelLeadingIndicator' plot( x, title = NULL, series.name.lead = "Leading Indicator", series.name.target = "Target Variable", take.log = TRUE, axis = NULL, ... )
x |
A |
title |
Title for the plot. |
series.name.lead |
Name of the leading indicator series, used in
the legend. Default |
series.name.target |
Name of the target variable series, used in
the legend. Default |
take.log |
Logical value indicating whether to plot the log of the
daily incidence. Default |
axis |
An |
... |
Unused. |
A ggplot2 plot.
summary() on FilterResults/
FilterResultsLI reports estimated variance parameters and
model states, but does not report the model's log-likelihood or
residual diagnostics. print_model_diagnostics() reports these
two additional, complementary diagnostics; it does not replace or
duplicate what summary() already prints.
print_model_diagnostics(res)print_model_diagnostics(res)
res |
A |
NULL, invisibly. Called for its printed output.
Accessor method to print a short description for the objects of
FilterResults class
## S3 method for class 'FilterResults' print(x, ...)## S3 method for class 'FilterResults' print(x, ...)
x |
A |
... |
Additional arguments. |
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) # Specify a model model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) # Estimate a specified model res <- estimate(model) # Return short description of fitted model print(res)library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) # Specify a model model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) # Estimate a specified model res <- estimate(model) # Return short description of fitted model print(res)
Accessor method to print a short description for the objects of
FilterResultsLI class
## S3 method for class 'FilterResultsLI' print(x, ...)## S3 method for class 'FilterResultsLI' print(x, ...)
x |
A |
... |
Additional arguments. |
library(tsgc) data(england, package = "tsgc") conv <- xts_to_idx(england[, 1:2]) out_eng <- SSModelLeadingIndicator( Y = conv$series, n.lag = 4, sea.period = 7, LeadIndCol = 1, calendar = conv$calendar, start = idx_to_pos(conv$calendar, "2021-04-30"), end = idx_to_pos(conv$calendar, "2021-07-24")) res_eng<-estimate(out_eng) print(res_eng)library(tsgc) data(england, package = "tsgc") conv <- xts_to_idx(england[, 1:2]) out_eng <- SSModelLeadingIndicator( Y = conv$series, n.lag = 4, sea.period = 7, LeadIndCol = 1, calendar = conv$calendar, start = idx_to_pos(conv$calendar, "2021-04-30"), end = idx_to_pos(conv$calendar, "2021-07-24")) res_eng<-estimate(out_eng) print(res_eng)
idx_calendar
Print an idx_calendar
## S3 method for class 'idx_calendar' print(x, ...)## S3 method for class 'idx_calendar' print(x, ...)
x |
An |
... |
Unused. |
idx_series
Print an idx_series
## S3 method for class 'idx_series' print(x, ...)## S3 method for class 'idx_series' print(x, ...)
x |
An |
... |
Unused. |
idx_step
Print an idx_step
## S3 method for class 'idx_step' print(x, ...)## S3 method for class 'idx_step' print(x, ...)
x |
An |
... |
Unused. |
multi_step_pattern
Print a multi_step_pattern
## S3 method for class 'multi_step_pattern' print(x, ...)## S3 method for class 'multi_step_pattern' print(x, ...)
x |
A |
... |
Unused. |
Accessor method to print a short description for the objects of
SSModelDynamicGompertz class
## S3 method for class 'SSModelDynamicGompertz' print(x, ...)## S3 method for class 'SSModelDynamicGompertz' print(x, ...)
x |
A |
... |
Additional arguments. |
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) # Specify a model model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-06")) # Print a short description of the model object print(model)library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) # Specify a model model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-06")) # Print a short description of the model object print(model)
Accessor method to print a short description for the objects of
SSModelLeadingIndicator class
## S3 method for class 'SSModelLeadingIndicator' print(x, ...)## S3 method for class 'SSModelLeadingIndicator' print(x, ...)
x |
A |
... |
Additional arguments. |
library(tsgc) data(england, package = "tsgc") conv <- xts_to_idx(england[, 1:2]) # Specify a model out_eng <- SSModelLeadingIndicator( Y = conv$series, n.lag = 4, sea.period = 7, LeadIndCol = 1, calendar = conv$calendar, start = idx_to_pos(conv$calendar, "2021-04-30"), end = idx_to_pos(conv$calendar, "2021-07-24")) # Print a short description of the model object print(out_eng)library(tsgc) data(england, package = "tsgc") conv <- xts_to_idx(england[, 1:2]) # Specify a model out_eng <- SSModelLeadingIndicator( Y = conv$series, n.lag = 4, sea.period = 7, LeadIndCol = 1, calendar = conv$calendar, start = idx_to_pos(conv$calendar, "2021-04-30"), end = idx_to_pos(conv$calendar, "2021-07-24")) # Print a short description of the model object print(out_eng)
Accessor method to access the non-diffuse part of the error
covariance matrix of the filtered states in a KFS object.
Ptt(object)Ptt(object)
object |
A |
The error covariance matrix of the filtered states.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::Ptt(output(res))library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) tsgc:::Ptt(output(res))
Takes a cumulated series and re-bases it so that it starts
from zero at reinit.idx - 1.
reinitialise_dataframe(dt, reinit.idx)reinitialise_dataframe(dt, reinit.idx)
dt |
Cumulated data series, as an |
reinit.idx |
Integer position at which reinitialisation should
occur (i.e. |
The reinitialised series, as an idx_series starting at
reinit.idx.
x <- idx_series(cumsum(rpois(30, 5)) + 1, start = 1) reinitialise_dataframe(x, 10)x <- idx_series(cumsum(rpois(30, 5)) + 1, start = 1) reinitialise_dataframe(x, 10)
Accessor method to access the number of seasonal components
used in a KFS object.
seasonalComp(object)seasonalComp(object)
object |
A |
The seasonal component specification stored on the model terms.
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) seasonalComp(output(res))library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) res <- estimate(model) seasonalComp(output(res))
Class for Dynamic Gompertz Curve State-Space Model Object, which encapsulates model settings and provides methods to obtain a FilterResults object.
The dynamic Gompertz model with an integrated random walk (IRW) trend is defined as:
where is the cumulative variable, , and
The trend component follows:
Here, the observation disturbances and slope disturbances are independent and normally distributed. The signal-to-noise ratio,
,
determines how rapidly the slope adjusts to new observations—higher values lead to faster changes, while lower values induce smoothness.
For models without seasonal terms (sea.period = 0), the priors are given by:
The diffuse prior is defined as with , implemented via the KFAS package (Helske, 2017). For models with a seasonal component (sea.period>1), the prior mean vector and prior covariance matrix are extended accordingly.
See the vignette for details on the state disturbance variance matrix and the observation noise variance .
This class also supports the implementation of the reinitialisation
procedure, described in the vignette and also summarised below.
Let denote the re-initialization position and denote the
position at which the cumulative series is set to 0. As the growth rate of
cumulative cases is defined as , we have:
where is the cumulative cases after re-initialization. We
choose to set the cumulative cases to zero at ,
such that the growth rate of cumulative cases is available from
onwards.
We reinitialise the model by specifying the prior distribution for the
initial states appropriately. See the vignette for details.
YThe cumulated variable, as an idx_series. Must be
strictly increasing.
qThe signal-to-noise ratio (ratio of slope to irregular
variance). Defaults to 'NULL', in which case no
signal-to-noise ratio will be imposed. Instead, it will be estimated.
sea.periodA positive integer specifying the period of seasonality used in the
trigonometric seasonal component of the model. For example, use 7 for daily
data to model day-of-the-week effects. A value of 0 disables the seasonal
component entirely. The default is 7, which is suitable for capturing
weekly seasonality in daily time series.
reinit.idx(Only needed for reinitialization.) The
reinitialisation position , as a single integer. Defaults to
NULL, which represents the non-reinitialized version.
original.results(Only needed for reinitialization.) Rather than re-estimating the model up
to the reinit.idx, a FilterResults class object can be
specified here and the parameters for the reinitialisation will be taken
from this object. Default is NULL. This parameter is optional.
use.presample.info(Only needed for reinitialization.) Logical value denoting whether or
not to use information from before the reinitialisation position in the
reinitialisation procedure. Default is TRUE. If FALSE, the
model is estimated from scratch from the reinitialisation position and no
attempt to use information from before the reinitialisation position is made.
xpredAn idx_series object containing the dataset of exogenous
variables to include in the model. Defaults to NULL.
ar1Logical value indicating whether an ar1 component should be
included in the model. Default is FALSE.
startInteger position marking the start of the estimation period.
endInteger position marking the end of the estimation period.
calendarAn optional idx_calendar object describing how the
integer positions used by Y (and the rest of this model's
idx_series fields) map back to calendar time. Purely cosmetic:
never consulted by estimation/filtering, only carried along so that
plotting can later translate positions to dates. Defaults to
NULL, in which case plots will be labelled with raw integer
positions instead of dates.
estimate()Estimates the dynamic Gompertz curve model when applied to an object of
class SSModelDynamicGompertz.
An object of class FilterResults
containing the result output for the estimated dynamic Gompertz curve
model.
initialize(
Y,
q = NULL,
sea.period = 7,
reinit.idx = NULL,
original.results = NULL,
use.presample.info = TRUE,
xpred = NULL,
ar1 = FALSE,
start = idx_range(Y)[1],
end = idx_range(Y)[2],
calendar = NULL
)Create an instance of the SSModelDynamicGompertz class. Parameters
are defined in 'fields' section.
SSModelDynamicGompertz$new(Y = y, q = 0.005,
reinit.idx = 45)
print()Provides a quick description of the SSModelDynamicGompertz object, providing model states and standard errors.
summary()Supplies details of the SSModelDynamicGompertz object, such as estimated parameter values, start and end positions of the estimation window.
library(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) # Specify a model model <- SSModelDynamicGompertz$new(Y = Y, q = 0.005, end = 100) # Show summary of the model object summary(model) # Print a short description of the model object print(model) # Estimate a specified model res <- estimate(model) reslibrary(tsgc) set.seed(1) Y <- idx_series(cumsum(rpois(120, 8)) + 1, start = 1) # Specify a model model <- SSModelDynamicGompertz$new(Y = Y, q = 0.005, end = 100) # Show summary of the model object summary(model) # Print a short description of the model object print(model) # Estimate a specified model res <- estimate(model) res
A class for specifying the parameters of a leading indicator model. The model settings are stored in the fields of this object, and the class contains methods to obtain a FilterResultsLI object for further analysis.
YA cumulated idx_series with 2 columns: a leading indicator and a target
variable. Both the target variable and the lagged leading indicator must be strictly
increasing within the estimation window.
qThe signal-to-noise ratio (ratio of slope error variance to target variable observation error variance).
Defaults to 'NULL', in which case no
signal-to-noise ratio will be imposed. Instead, it will be estimated.
sea.periodA positive integer specifying the period of seasonality used in the
trigonometric seasonal component of the model. For example, use 7 for daily
data to model day-of-the-week effects. A value of 0 disables the seasonal
component entirely. The default is 7, which is suitable for capturing
weekly seasonality in daily time series.
n.lagNumber of integer positions to lag the leading indicator by.
xpred_leadAn idx_series object containing the values of exogenous
variables for the leading indicator. Dataset must contain values for all positions in the
estimation time frame. Defaults to NULL, indicating no exogenous variables are needed
for the leading indicator.
xpred_targAn idx_series object containing the values of exogenous
variables for the target variable. Dataset must contain values for all positions in the
estimation time frame. Defaults to NULL, indicating no exogenous variables are
needed for the target variable.
LeadIndColThe column in Y that contains the leading indicator.
Defaults to 1.
startInteger position marking the start of the estimation period.
endInteger position marking the end of the estimation period.
calendarAn optional idx_calendar object describing how the
integer positions used by Y map back to calendar time. Since the
leading indicator and target series in Y are already combined
onto one shared position axis (see Y), a single calendar
correctly applies to both. Purely cosmetic: never consulted by
estimation/filtering, only carried along so that plotting can later
translate positions to dates. Defaults to NULL, in which case
plots will be labelled with raw integer positions instead of dates. It
is the user's responsibility to ensure the leading indicator and target
series were aligned onto a shared, consistent calendar before being
combined into Y; this class does not attempt to reconcile two
different calendars.
estimate()Estimates the Leading Indicator model when applied to an object of
class SSModelLeadingIndicator.
An object of class FilterResultsLI
containing the result output for the estimated Leading Indicator
model.
initialize(
Y,
n.lag,
sea.period = 7,
q = NULL,
LeadIndCol = 1,
xpred_lead = NULL,
xpred_targ = NULL,
start = idx_range(Y)[1],
end = idx_range(Y)[2],
calendar = NULL
)Create an instance of the SSModelLeadingIndicator class with the
fields laid out at the beginning of the documentation.
print()Provides a quick description of the SSModelLeadingIndicator object, providing model states and standard errors.
summary()Supplies details of the SSModelLeadingIndicator object, such as estimated parameter values, start and end positions of estimation.
library(tsgc) set.seed(1) lead <- cumsum(rpois(120, 6)) + 1 targ <- cumsum(rpois(120, 8)) + 1 Y <- idx_series(cbind(lead, targ), start = 1) # Specify a model with the estimation timeframe model <- SSModelLeadingIndicator(Y = Y, n.lag = 5, sea.period = 0, LeadIndCol = 1, start = 1, end = 100) # Show summary of the model object summary(model) # Print a short description of the model object print(model) # Estimate a specified model res <- estimate(model) reslibrary(tsgc) set.seed(1) lead <- cumsum(rpois(120, 6)) + 1 targ <- cumsum(rpois(120, 8)) + 1 Y <- idx_series(cbind(lead, targ), start = 1) # Specify a model with the estimation timeframe model <- SSModelLeadingIndicator(Y = Y, n.lag = 5, sea.period = 0, LeadIndCol = 1, start = 1, end = 100) # Show summary of the model object summary(model) # Print a short description of the model object print(model) # Estimate a specified model res <- estimate(model) res
Accessor method to show a summary for the objects of
FilterResults class
## S3 method for class 'FilterResults' summary(object, ...)## S3 method for class 'FilterResults' summary(object, ...)
object |
A |
... |
Additional arguments. |
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) # Specify a model model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) # Estimate a specified model res <- estimate(model) # Return KFS object in output of res summary(res)library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) # Specify a model model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-20")) # Estimate a specified model res <- estimate(model) # Return KFS object in output of res summary(res)
Accessor method to show a summary for the objects of
FilterResultsLI class
## S3 method for class 'FilterResultsLI' summary(object, ...)## S3 method for class 'FilterResultsLI' summary(object, ...)
object |
A |
... |
Additional arguments. |
library(tsgc) data(england, package = "tsgc") conv <- xts_to_idx(england[, 1:2]) out_eng <- SSModelLeadingIndicator( Y = conv$series, n.lag = 4, sea.period = 7, LeadIndCol = 1, calendar = conv$calendar, start = idx_to_pos(conv$calendar, "2021-04-30"), end = idx_to_pos(conv$calendar, "2021-07-24")) res_eng<-estimate(out_eng) summary(res_eng)library(tsgc) data(england, package = "tsgc") conv <- xts_to_idx(england[, 1:2]) out_eng <- SSModelLeadingIndicator( Y = conv$series, n.lag = 4, sea.period = 7, LeadIndCol = 1, calendar = conv$calendar, start = idx_to_pos(conv$calendar, "2021-04-30"), end = idx_to_pos(conv$calendar, "2021-07-24")) res_eng<-estimate(out_eng) summary(res_eng)
Accessor method to show a summary for the objects of
SSModelDynamicGompertz class
## S3 method for class 'SSModelDynamicGompertz' summary(object, ...)## S3 method for class 'SSModelDynamicGompertz' summary(object, ...)
object |
A |
... |
Additional arguments. |
library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) # Specify a model model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-06")) # Show summary of the model object summary(model)library(tsgc) data(gauteng, package = "tsgc") conv <- xts_to_idx(gauteng) # Specify a model model <- SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-06")) # Show summary of the model object summary(model)
Accessor method to show a summary for the objects of
SSModelLeadingIndicator class
## S3 method for class 'SSModelLeadingIndicator' summary(object, ...)## S3 method for class 'SSModelLeadingIndicator' summary(object, ...)
object |
A |
... |
Additional arguments. |
library(tsgc) data(england, package = "tsgc") conv <- xts_to_idx(england[, 1:2]) # Specify a model out_eng <- SSModelLeadingIndicator( Y = conv$series, n.lag = 4, sea.period = 7, LeadIndCol = 1, calendar = conv$calendar, start = idx_to_pos(conv$calendar, "2021-04-30"), end = idx_to_pos(conv$calendar, "2021-07-24")) summary(out_eng)library(tsgc) data(england, package = "tsgc") conv <- xts_to_idx(england[, 1:2]) # Specify a model out_eng <- SSModelLeadingIndicator( Y = conv$series, n.lag = 4, sea.period = 7, LeadIndCol = 1, calendar = conv$calendar, start = idx_to_pos(conv$calendar, "2021-04-30"), end = idx_to_pos(conv$calendar, "2021-07-24")) summary(out_eng)
idx_series
Returns the last n observations of x, in
position order. As with head.idx_series, n is a
count of observations, not a set of absolute positions.
## S3 method for class 'idx_series' tail(x, n = 6L, ...)## S3 method for class 'idx_series' tail(x, n = 6L, ...)
x |
An |
n |
Number of observations to return from the end of |
... |
Unused. |
An idx_series of length min(n, length(x)).
Cumulative daily cases of Covid-19 in Italy and UK, before 14 Dec 2020.
data(ukitaly)data(ukitaly)
An object of class "xts" with columns:
Cumulative cases of Covid-19 in Italy
Cumulative cases of Covid-19 in the UK
Downloaded from https://www.ecdc.europa.eu/en/publications-data/download-todays-data-geographic-distribution-covid-19-cases-worldwide
data(ukitaly) # plot daily cases plot(diff(ukitaly))data(ukitaly) # plot daily cases plot(diff(ukitaly))
Function writes the following results to csv files which get
saved in the location specified in res.dir: forecast new cases or
incidence variable, ; the filtered level and slope of ,
and ; filtered estimates of and the
confidence intervals for these estimates.
write_results(res, res.dir, n.ahead, prefix = "", confidence.level = 0.68)write_results(res, res.dir, n.ahead, prefix = "", confidence.level = 0.68)
res |
Results object of class |
res.dir |
File path to save the results to. A character string. |
n.ahead |
Number of periods ahead to forecast. A positive integer. |
prefix |
The prefix to be added to the file names generated. A character string. |
confidence.level |
Confidence level to use for the confidence interval
on the forecasts |
A number of csv files saved in the directory specified in
res.dir.
# Not run as do not wish to save to local disk when compiling documentation. # Below will run if copied and pasted into console. library(tsgc) library(here) res.dir <- tempdir() data(gauteng,package="tsgc") conv <- xts_to_idx(gauteng) res <- estimate(SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-06"))) tsgc::write_results( res=res, res.dir = res.dir, prefix="dyn_gompertz",n.ahead = 14, confidence.level = 0.68)# Not run as do not wish to save to local disk when compiling documentation. # Below will run if copied and pasted into console. library(tsgc) library(here) res.dir <- tempdir() data(gauteng,package="tsgc") conv <- xts_to_idx(gauteng) res <- estimate(SSModelDynamicGompertz$new(Y = conv$series, q = 0.005, calendar = conv$calendar, end = idx_to_pos(conv$calendar, "2020-07-06"))) tsgc::write_results( res=res, res.dir = res.dir, prefix="dyn_gompertz",n.ahead = 14, confidence.level = 0.68)
xts/zoo object to an idx_series
Converts a calendar-indexed xts or zoo object
into an idx_series, together with an
idx_calendar describing how to translate back to calendar
time. This is the standard entry point for bringing calendar-indexed
data (the user's own xts/zoo/data frame data, or one of
the package's bundled datasets) into the idx_series-based
analysis functions.
By default (detect = TRUE), the step size and any repeating gap
pattern are auto-detected from x's index via
idx_detect_calendar_pattern (searching the full unit
ladder years/quarters/months/weeks/days/hours/minutes/seconds), so
callers do not need to work out amount/unit/pattern
themselves for common cases: a plain daily/weekly/etc. index, a
business-day index (weekends skipped), a monthly/quarterly/yearly
Date index, or a yearqtr/yearmon index. This only
recognises a single amount/unit pair (optionally repeated
via an integer vector pattern); it does not attempt to detect
compound (idx_step) or heterogeneous (multi_step_pattern)
calendars - build those manually via idx_calendar_step/
idx_calendar_multi_step if needed. If detection fails
(the index has no such regular pattern, e.g. genuinely irregular gaps),
xts_to_idx raises an error explaining why, rather than silently
falling back to an incorrect calendar; pass detect = FALSE (and
set amount/unit yourself) for indices detection can't
handle.
xts_to_idx(x, start.pos = 1L, detect = TRUE, amount = 1, unit = "days")xts_to_idx(x, start.pos = 1L, detect = TRUE, amount = 1, unit = "days")
x |
An |
start.pos |
Integer position that the first row of |
detect |
A single logical. If |
amount, unit
|
Used only when |
A list with two elements: series, an idx_series
holding x's values, and calendar, an idx_calendar
anchoring series's positions to x's original dates.
# Daily data - detected automatically. x <- xts::xts(cumsum(rpois(30, 5)) + 1, order.by = Sys.Date() - 29:0) conv <- xts_to_idx(x) conv$series conv$calendar # Business-day data (weekends skipped) - pattern detected automatically. bdays <- seq(as.Date("2024-01-01"), by = "day", length.out = 40) bdays <- bdays[!weekdays(bdays) %in% c("Saturday", "Sunday")] xb <- xts::xts(cumsum(rpois(length(bdays), 5)) + 1, order.by = bdays) conv_b <- xts_to_idx(xb) conv_b$calendar # Monthly Date index - detected as amount = 1, unit = "months" (not days). months <- seq(as.Date("2024-01-01"), by = "month", length.out = 24) xm <- xts::xts(cumsum(rpois(24, 5)) + 1, order.by = months) conv_m <- xts_to_idx(xm) conv_m$calendar# Daily data - detected automatically. x <- xts::xts(cumsum(rpois(30, 5)) + 1, order.by = Sys.Date() - 29:0) conv <- xts_to_idx(x) conv$series conv$calendar # Business-day data (weekends skipped) - pattern detected automatically. bdays <- seq(as.Date("2024-01-01"), by = "day", length.out = 40) bdays <- bdays[!weekdays(bdays) %in% c("Saturday", "Sunday")] xb <- xts::xts(cumsum(rpois(length(bdays), 5)) + 1, order.by = bdays) conv_b <- xts_to_idx(xb) conv_b$calendar # Monthly Date index - detected as amount = 1, unit = "months" (not days). months <- seq(as.Date("2024-01-01"), by = "month", length.out = 24) xm <- xts::xts(cumsum(rpois(24, 5)) + 1, order.by = months) conv_m <- xts_to_idx(xm) conv_m$calendar