R/TwoBlockUnitRoots.R

In YueyunZhu/UnitRootsofPIAR: Eliminate unit roots for periodically integrated autoregressive models

#' @title Function to compute the coefficients of periodic integrated filter
#' @description Function to compute the coefficients of periodic integrated filter.
#'              Assume there are two block unit roots in the process.
#' @param Xsim similarity matrix for Jordan decomposition of the multi-companion matrix F
#' @param d period
#'
#' @return coefficients of $(1-\beta_sL)(1-\alpha_sL)$
#' @export
#'
#' @examples
#' library(mcompanion)
#' Fd <- sim_mc(dim = 4, mo = 4, mo.col = 2, eigval = c(1), len.block = c(2))
#' icoef <- piar2unit.block(Xsim = Fd$eigvec, d = 4)
piar2unit.block <- function(Xsim, d = 4) {
  vecX1 <- Xsim[, 1]
  vecX2 <- Xsim[, 2]
  alpha <- rep(NA, d)
  beta <- rep(NA, d)
  Delta <- matrix(NA, d, d)
  for (i in 1:d) {
    {
      for (j in 1:d) {
        Delta[i, j] <- vecX1[i] * vecX2[j] - vecX1[j] * vecX2[i]
      }
    }
  }
  alpha[1] <- vecX1[4]/vecX1[1]
  for (i in 2:d) {
    alpha[i] <- vecX1[d - i + 1]/vecX1[d - i + 2]
  }
  beta[1] <- -(vecX1[2] * (vecX1[1] * (vecX1[4] + vecX2[4]) - vecX1[4] * vecX2[1]))/
    (vecX1[1] * Delta[1, 2])
  beta[2] <- (vecX1[1] * Delta[4, 3])/
    (vecX1[4] * (vecX1[1] * (vecX1[4] + vecX2[4]) - vecX1[4] * vecX2[1]))
  beta[3] <- (vecX1[4] * Delta[2, 3])/(vecX1[3] * Delta[3, 4])
  beta[4] <- (vecX1[3] * Delta[1, 2])/(vecX1[2] * Delta[2, 3])
  return(par = list(alpha = alpha, beta = beta))
}