Nothing
#' @title Check the stationarity and identification conditions of specified GMAR, StMAR, or G-StMAR model.
#'
#' @description \code{is_stationary_int} checks the stationarity condition and \code{is_identifiable} checks the identification condition
#' of the specified GMAR, StMAR, or G-StMAR model.
#'
#' @inheritParams loglikelihood_int
#' @details \code{is_stationary_int} does not support models imposing linear constraints. In order to use it for a model imposing linear
#' constraints, one needs to expand the constraints first to obtain an unconstrained parameter vector.
#'
#' Note that \code{is_stationary_int} returns \code{FALSE} for stationary parameter vectors if they are extremely close to the boundary
#' of the stationarity region.
#'
#' \code{is_identifiable} checks that the regimes are sorted according to the mixing weight parameters and that there are no duplicate
#' regimes.
#' @return Returns \code{TRUE} or \code{FALSE} accordingly.
#' @section Warning:
#' These functions don't have any argument checks!
#' @references
#' \itemize{
#' \item Kalliovirta L., Meitz M. and Saikkonen P. 2015. Gaussian Mixture Autoregressive model for univariate time series.
#' \emph{Journal of Time Series Analysis}, \strong{36}(2), 247-266.
#' \item Meitz M., Preve D., Saikkonen P. 2023. A mixture autoregressive model based on Student's t-distribution.
#' \emph{Communications in Statistics - Theory and Methods}, \strong{52}(2), 499-515.
#' \item Virolainen S. 2022. A mixture autoregressive model based on Gaussian and Student's t-distributions.
#' Studies in Nonlinear Dynamics & Econometrics, \strong{26}(4) 559-580.
#' }
#' @keywords internal
is_stationary_int <- function(p, M, params, restricted=FALSE) {
M <- sum(M)
if(restricted == FALSE) {
pars <- matrix(params[1:(M*(p + 2))], ncol=M)
for(i1 in 1:M) { # Go through regimes
if(any(abs(polyroot(c(1, -pars[2:(p + 1), i1]))) <= 1 + 1e-8)) { # The stationarity condition with num tol 1e-8
return(FALSE)
}
}
} else { # restricted -> common AR coefficients in all regimes
absroots <- abs(polyroot(c(1, -params[(M + 1):(M + p)])))
if(any(absroots <= 1 + 1e-8)) {
return(FALSE)
}
}
TRUE
}
#' @rdname is_stationary_int
is_identifiable <- function(p, M, params, model=c("GMAR", "StMAR", "G-StMAR"), restricted=FALSE, constraints=NULL) {
# Pick parameters etc
model <- match.arg(model)
params <- remove_all_constraints(p=p, M=M, params=params, model=model, restricted=restricted, constraints=constraints)
M_orig <- M
if(model == "G-StMAR") {
M1 <- M[1]
M2 <- M[2]
M <- sum(M)
if(M1 == 1 & M2 == 1) {
return(TRUE)
}
}
if(M == 1) {
return(TRUE)
}
pars <- pick_pars(p=p, M=M_orig, params=params, model=model, restricted=FALSE, constraints=NULL)
alphas <- pick_alphas(p=p, M=M_orig, params=params, model=model, restricted=FALSE, constraints=NULL)
dfs <- pick_dfs(p=p, M=M_orig, params=params, model=model)
if(model == "StMAR") {
pars <- rbind(pars, dfs)
} else if(model == "G-StMAR") {
pars1 <- as.matrix(pars[,1:M1])
pars2 <- as.matrix(pars[,(M1 + 1):M])
pars2 <- rbind(pars2, dfs)
alphas1 <- alphas[1:M1]
alphas2 <- alphas[(M1 + 1):M]
}
# Functions for determining are mixing weight parameters in the decreasing order
# or are there regimes with identical parameter values
alphas_unsorted <- function(alps) is.unsorted(rev(alps), strictly=TRUE)
pars_dublicates <- function(prs0) anyDuplicated(t(prs0)) != 0
if(model == "G-StMAR") { # Check GMAR parameters and StMAR parameters separately
if(M1 > 1) {
if(alphas_unsorted(alphas1)) {
return(FALSE)
} else if(pars_dublicates(pars1)) {
return(FALSE)
}
}
if(M2 > 1) {
if(alphas_unsorted(alphas2)) {
return(FALSE)
} else if(pars_dublicates(pars2)) {
return(FALSE)
}
}
} else { # if model != "G-StMAR
if(alphas_unsorted(alphas)) {
return(FALSE)
} else if(pars_dublicates(pars)) {
return(FALSE)
}
}
TRUE
}
#' @title Check the stationary condition of specified GMAR, StMAR, or G-StMAR model.
#'
#' @description \code{is_stationary} checks the stationarity condition of the specified GMAR, StMAR, or G-StMAR model.
#'
#' @inheritParams loglikelihood
#' @details This function falsely returns \code{FALSE} for stationary models when the parameter is extremely close
#' to the boundary of the stationarity region.
#' @inherit is_stationary_int return references
#' @examples
#' # GMAR model
#' params22 <- c(0.4, 0.39, 0.6, 0.3, 0.4, 0.1, 0.6, 0.3, 0.8)
#' is_stationary(p=2, M=2, params=params22)
#'
#' # StMAR model
#' params12t <- c(-0.3, 1, 0.9, 0.1, 0.8, 0.6, 0.7, 10, 12)
#' is_stationary(p=1, M=2, params=params12t, model="StMAR")
#'
#' # G-StMAR model
#' params12gs <- c(1, 0.1, 1, 2, 0.2, 2, 0.8, 20)
#' is_stationary(p=1, M=c(1, 1), params=params12gs, model="G-StMAR")
#'
#' # Restricted GMAR model
#' params13r <- c(0.1, 0.2, 0.3, -0.99, 0.1, 0.2, 0.3, 0.5, 0.3)
#' is_stationary(p=1, M=3, params=params13r, restricted=TRUE)
#'
#' # Such StMAR(3, 2) that the AR coefficients are restricted to be the
#' # same for both regimes and that the second AR coefficients are
#' # constrained to zero.
#' params32trc <- c(1, 2, 0.8, -0.3, 1, 2, 0.7, 11, 12)
#' is_stationary(p=3, M=2, params=params32trc, model="StMAR", restricted=TRUE,
#' constraints=matrix(c(1, 0, 0, 0, 0, 1), ncol=2))
#' @export
is_stationary <- function(p, M, params, model=c("GMAR", "StMAR", "G-StMAR"), restricted=FALSE, constraints=NULL) {
# Pick and check parameter values
model <- match.arg(model)
check_pM(p=p, M=M, model=model)
check_params_length(p=p, M=M, params=params, model=model, restricted=restricted, constraints=constraints)
if(!is.null(constraints)) {
check_constraint_mat(p=p, M=M, restricted=restricted, constraints=constraints)
params <- reform_constrained_pars(p=p, M=M, params=params, model=model, restricted=restricted, constraints=constraints)
}
is_stationary_int(p=p, M=M, params=params, restricted=restricted) # Call the internal function to determine whether the process is stationary
}
#' @title Check that the data is set correctly and correct if not
#'
#' @description \code{check_and_correct_data} checks that the data is set correctly and
#' throws an error if there is something wrong with the data.
#'
#' @inheritParams loglikelihood_int
#' @return Returns the data as a class 'ts' object.
#' @keywords internal
check_and_correct_data <- function(data, p) {
if(anyNA(data)) { # Check if there are problems with the data
stop("The data contains NA values")
} else if(length(data) < p + 2) {
stop("The data must contain at least p+2 values")
} else if(!is.numeric(data)) {
stop("The data should be numeric")
}
if(is.matrix(data)) {
if(ncol(data) > 1) {
stop("Only univariate time series are supported! For multivariate analysis, try the package 'gmvarkit'.")
}
}
if(!is.ts(data)) { # Make data a class ts object it is not already
data <- as.ts(data)
}
data
}
#' @title Check the parameter vector is specified correctly
#'
#' @description \code{parameter_checks} checks dimension, restrictions, and stationarity of the given parameter
#' of a GMAR, StMAR, or G-StMAR model. Throws an error if any check fails. Does NOT consider the identifiability
#' condition!
#'
#' @inheritParams loglikelihood_int
#' @return Throws an informative error if any check fails. Does not return anything.
#' @keywords internal
parameter_checks <- function(p, M, params, model=c("GMAR", "StMAR", "G-StMAR"), restricted=FALSE, constraints=NULL) {
# Pick the parameter values
model <- match.arg(model)
check_params_length(p=p, M=M, params=params, model=model, restricted=restricted, constraints=constraints)
params <- remove_all_constraints(p=p, M=M, params=params, model=model, restricted=restricted, constraints=constraints)
pars <- pick_pars(p=p, M=M, params=params, model=model, restricted=FALSE, constraints=NULL)
alphas <- pick_alphas(p=p, M=M, params=params, model=model, restricted=FALSE, constraints=NULL)
dfs <- pick_dfs(p=p, M=M, params=params, model=model)
# Check degrees of freedom parameters for StMAR and G-StMAR models
if(model == "StMAR" | model == "G-StMAR") {
if(any(dfs <= 2)) {
stop("The degrees of freedom parameters have to be larger than 2")
} else if(any(dfs > 1e+5)) {
stop("We have set an upper bound of 1e+5 for the degrees of freedom parameters
in order to avoid numerical problems. This is not, however, restrictive
since t-distribution with such a large degrees of freedom parameter value
already strongly resembles the Gaussian distribution by its shape. Moreover,
the conditional variance of the corresponding regime is close to the constant
parameter 'sigma^2' if the degrees of freedom parameter is very large.")
}
}
# Check assumptions regarding the rest of the parameters
if(sum(M) >= 2 & sum(alphas[-sum(M)]) >= 1) {
stop("The mixing weight parameters don't sum to one")
} else if(any(alphas <= 0)) {
stop("The mixing weight parameters have to be strictly larger than zero")
} else if(any(pars[p + 2,] <= 0)) {
stop("The variance parameters have to be strictly larger than zero")
} else if(!is_stationary_int(p=p, M=M, params=params, restricted=FALSE)) {
stop("The model doesn't satisfy the stationarity condition")
}
}
#' @title Check the constraint matrices
#'
#' @description \code{check_constraint_mat} checks for some parts that the constraint matrices are correctly set.
#' @inheritParams loglikelihood_int
#' @return Doesn't return anything but throws an informative error if finds out that something is wrong.
#' @keywords internal
check_constraint_mat <- function(p, M, restricted=FALSE, constraints=NULL) {
if(!is.null(constraints)) { # Check only if constraints are imposed
M <- sum(M)
if(restricted == TRUE) { # 'constraints' is a single matrix
if(!is.matrix(constraints)) {
stop("The constraint matrix has to be a matrix (not in a list)")
} else if(nrow(as.matrix(constraints)) != p) {
stop("The constraint matrix has wrong dimension")
} else if(ncol(as.matrix(constraints)) > p) {
stop("The constraint matrix has more columns than rows!? Please make sure your constraints make sense!")
} else if(qr(as.matrix(constraints))$rank != ncol(as.matrix(constraints))) {
stop("The constraint matrix does not have full column rank")
}
} else { # restricted == FALSE, so 'constraints' is a list of matrices
if(!is.list(constraints) | length(constraints) != M) {
stop("The argument constraints should be a list of M constraint matrices - one for each mixture component")
}
for(i1 in 1:M) { # Go through the regimes and check constraints for each of them separately
C0 <- as.matrix(constraints[[i1]])
q <- ncol(C0)
if(nrow(C0) != p) {
stop(paste("The constraint matrix for regime", i1 ,"has wrong dimension"))
} else if(q > p) {
stop(paste("The constraint matrix for regime", i1, "has more columns than rows!? Please make sure your constraints make sense!"))
} else if(qr(C0)$rank != ncol(C0)) {
stop(paste("The constraint matrix for regime", i1 ,"does not have full column rank"))
}
}
}
}
}
#' @title Check that p and M are correctly set
#'
#' @description \code{check_pM} checks that the arguments p and M are correctly set.
#' @inheritParams loglikelihood_int
#' @return Doesn't return anything but throws an informative error if something is wrong.
#' @keywords internal
check_pM <- function(p, M, model=c("GMAR", "StMAR", "G-StMAR")) {
model <- match.arg(model)
if(model == "G-StMAR") { # Check that M is correctly set for G-StMAR model
if(length(M) != 2 | !all_pos_ints(M)) {
stop("For a G-StMAR model the argument M should be a length 2 positive integer vector")
}
} else if(!all_pos_ints(M) | length(M) != 1) { # Check that M is correctly set for GMAR and StMAR models
stop("Argument M has to be positive integer")
}
if(!all_pos_ints(p) | length(p) != 1) { # Check that p is a positive integer
stop("Argument p has to be positive integer")
}
}
#' @title Calculate the number of parameters
#'
#' @description \code{n_params} calculates the number of parameters that should be in the parameter vector.
#' @inheritParams loglikelihood_int
#' @return Returns the number of parameters.
#' @keywords internal
n_params <- function(p, M, model=c("GMAR", "StMAR", "G-StMAR"), restricted=FALSE, constraints=NULL) {
model <- match.arg(model)
# We go through case by case as the number of parameters is different with different types of constraints
if(restricted == FALSE) {
# Function to calculate the number of AR parameters in each regime when constraints are employed
n_const_pars <- function(M, constraints) sum(vapply(1:M, function(i1) ncol(as.matrix(constraints[[i1]])), numeric(1)))
if(model == "StMAR") {
if(is.null(constraints)) {
d <- M*(p + 4) - 1
} else {
d <- 4*M - 1 + n_const_pars(M=M, constraints=constraints)
}
} else if(model == "G-StMAR") {
M1 <- M[1]
M2 <- M[2]
M <- sum(M)
if(is.null(constraints)) {
d <- M*(p + 3) + M2 - 1
} else {
d <- 3*M + M2 - 1 + n_const_pars(M=M, constraints=constraints)
}
} else { # If model == "GMAR"
if(is.null(constraints)) {
d <- M*(p + 3) - 1
} else {
d <- 3*M - 1 + n_const_pars(M=M, constraints=constraints)
}
}
} else { # if restricted == TRUE
if(model == "StMAR") {
if(is.null(constraints)) {
d <- 4*M + p - 1
} else {
d <- 4*M + ncol(as.matrix(constraints)) - 1
}
} else if(model == "G-StMAR") {
M1 <- M[1]
M2 <- M[2]
M <- sum(M)
if(is.null(constraints)) {
d <- 3*M + M2 + p - 1
} else {
d <- 3*M + M2 + ncol(as.matrix(constraints)) - 1
}
} else { # if model == "GMAR"
if(is.null(constraints)) {
d <- 3*M + p - 1
} else {
d <- 3*M + ncol(as.matrix(constraints)) - 1
}
}
}
d
}
#' @title Check whether all arguments are strictly positive natural numbers
#'
#' @description \code{all_pos_ints} tells whether all the elements in a vector
#' are strictly positive natural numbers.
#'
#' @param x a vector containing the elements to be tested.
#' @return Returns \code{TRUE} or \code{FALSE} accordingly.
#' @keywords internal
all_pos_ints <- function(x) {
all(vapply(x, function(x) x %% 1 == 0 & length(x) == 1 & x >= 1, logical(1)))
}
#' @title Check that the argument 'model' is correctly specified.
#'
#' @description \code{check_model} checks that the argument 'model' is correctly specified.
#'
#' @inheritParams loglikelihood_int
#' @return Doesn't return anything but throws and error if something is wrong.
#' @keywords internal
check_model <- function(model) {
if(!model %in% c("GMAR", "StMAR", "G-StMAR")) {
stop("The argument 'model' has to be 'GMAR', 'StMAR', or 'G-StMAR'")
}
}
#' @title Check that the parameter vector has the correct dimension
#'
#' @description \code{check_model} checks that the parameter vector has the correct dimension.
#' @inheritParams loglikelihood_int
#' @inherit check_model return
#' @keywords internal
check_params_length <- function(p, M, params, model=c("GMAR", "StMAR", "G-StMAR"), restricted=FALSE, constraints=NULL) {
if(length(params) != n_params(p=p, M=M, model=model, restricted=restricted, constraints=constraints)) {
stop("The parameter vector has wrong dimension")
}
}
#' @title Check that given object has class attribute 'gsmar'
#'
#' @description \code{check_gsmar} checks that the given object has class attribute 'gsmar'.
#'
#' @param object an object to be tested
#' @param object_name the name of the tested object
#' @inherit check_model return
#' @keywords internal
check_gsmar <- function(object, object_name="gsmar") {
if(!any(class(object) == "gsmar")) stop(paste("The argument", object_name, "has to be an object of class 'gsmar', often created with function 'fitGSMAR' or 'GSMAR'."))
}
#' @title Check that given object contains data
#'
#' @description \code{check_data} checks that a given object contains data.
#'
#' @inheritParams check_gsmar
#' @inherit check_gsmar return
#' @keywords internal
check_data <- function(object) {
if(is.null(object$data)) stop("The model has to contain data! Data can be added with the function add_data.")
}
#' @title Warn about large degrees of freedom parameter values
#'
#' @description \code{warn_dfs} warns if the model contains large degrees of freedom parameter values.
#'
#' @inheritParams check_gsmar
#' @inheritParams loglikelihood_int
#' @param warn_about warn about inaccurate derivatives or standard errors?
#' @details Either provide a class 'gsmar' object or specify the model by hand.
#' @return Doesn't return anything but throws a warning if any degrees of freedom parameters have value
#' larger than 100.
#' @keywords internal
warn_dfs <- function(object, p, M, params, model=c("GMAR", "StMAR", "G-StMAR"), warn_about=c("derivs", "errors")) {
if(!missing(object)) { # If the model is given as class gsmar object, extract the required information from it.
M <- object$model$M
params <- object$params
model <- object$model$model
}
if(model %in% c("StMAR", "G-StMAR")) { # Check whether there are large degrees of freedom parameters in some regime
M2 <- ifelse(model == "G-StMAR", M[2], M)
d <- length(params)
dfs <- params[(d - M2 + 1):d]
if(any(dfs > 100)) {
warning("The model contains overly large degrees of freedom parameter values. Consider switching to a G-StMAR model by setting the corresponding regimes to be GMAR type with the function 'stmar_to_gstmar'.")
}
}
}
#' @title Warn about near-unit-roots in some regimes
#'
#' @description \code{warn_ar_roots} warns if the model contains near-unit-roots in some regimes
#'
#' @inheritParams simulate.gsmar
#' @param tol if some root is smaller that \code{1 + tol}, a warning is thrown
#' @details Warns if some moduli of the autoregressive polynomial's roots are close to one.
#' @return Doesn't return anything.
#' @keywords internal
warn_ar_roots <- function(gsmar, tol=0.005) {
ar_roots <- get_ar_roots(gsmar) # Roots of the AR-polynomials
near_nonstat <- vapply(1:sum(gsmar$model$M), function(i1) any(abs(ar_roots[[i1]]) < 1 + tol), logical(1)) # Which are near the boundary of the stationarity region
if(any(near_nonstat)) { # If some are near the boundary, throw a warning
my_string <- ifelse(sum(near_nonstat) == 1,
paste("Regime", which(near_nonstat),"has near-unit-roots!"),
paste("Regimes", paste(which(near_nonstat), collapse=" and ") ,"have near-unit-roots!"))
warning(paste(my_string, "Consider building a model from the next-largest local maximum with the function 'alt_gsmar' by adjusting its argument 'which_largest'."))
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.