R/Ypresent.R

In yymmhaha/PackPaper1: Package related to my first paper

#' Generate linear time lagged response with time series covariates
#'
#' @description The response (univariate) is generated by the formula
#' \deqn{y^{(t)} = X^{(t)}\beta^{(0)} + X^{(t-1)}\beta^{(1)} + ... + \epsilon_t}
#' where each \eqn{X} is a vector of covariates at a certain time point, with corresponding regression coefficients.
#' Also automatically prepares the necessary lagged response matrix for second step (residual fitting)
#'
#' @details The most recent time stamp for covariates is \eqn{t}, the same as the response.
#'
#' The input matrix X is automatically shifted given the number of lags. Due to the shifting scheme for the time lags,
#' the length of the output response and extended matrix (nrow) are (t - lagx).
#'
#'
#' @param beta a list. Default number of list is 5, corresponding to the furthest lag \eqn{x^{t-4}}. Each list is a vector of regression
#' coefficients, corresponding to all covariates at that time stamp.
#' @param sigmay a number. Standard deviation of noise for response.
#' @param X a matrix. Covariate matrix, each column being one \eqn{x}. This corresponds to the output X from \code{XgenSimple, XgenCorr}s.
#' @param lagx a number. Number of lags for covarites, default is 5, meaning that the furthest covariate is \eqn{x^{(t-4)}}.
#' @param lagy a number. Number of lags for response for the preparation of second step regression,
#' default is 5, meaning that the furthest response is \eqn{y^{(t-5)}}.
#'
#' @return a list of components
#'
#'
#' \item{X}{a matrix. The original input matrix X}
#' \item{Y}{a vector. The generated linear time lagged response}
#' \item{ExtendX}{a matrix. Dimension is (t - lagx) by lagx*length(beta), 495 by 50 in our examples}
#' \item{ExtendY}{a matrix. Dimension is (t - lagy) by lagy}
#'
#'
#' @export
#'
#' @examples
#' # see
#'
#'



y.present <- function(beta, sigmay, X, lagx, lagy){

  t <- nrow(X)
  noisey <- rnorm((t - lagx), 0, sigmay)

  # ----------------------------- #
  # 1. generate extended X matrix and compute their product with beta
  # replace extendX1 with extendX2
  extendX2 <- list()
  for (i in 1:(lagx)){   # 6 blocks
    extendX2[[i]] <- X[((lagx+1)-i):(t-i), ]
  }
  ExtendX2 <- do.call(cbind, extendX2)
  betavec <- unlist(beta)
  X.inner <- ExtendX2 %*% betavec



  # ----------------------------- #
  # 2. generate y without AR y
  # length should be 997

  y <- X.inner + noisey

  # ----------------------------- #
  # 3. create lagged y matrix
  # for simplification, do not remove first dummy y0. keep it the same length as ExtendX
  # it's ok just to use 0 in the first several rows

  ExtendY <- matrix(rep(0, (t-lagx)*lagy), (t-lagx), lagy)
  for (j in 1:lagy){
    ExtendY[, j] <- c(rep(0, j), y[1:(length(y) - j)]  )
  }


  return(list('ExtendX' = ExtendX2, 'ExtendY' = ExtendY,
              'Y' = y, 'X' = X ))
}