Skip to content

Commit

Permalink
version 0.2.8
Browse files Browse the repository at this point in the history
  • Loading branch information
Inder Tecuapetla-Gómez authored and cran-robot committed Jun 29, 2023
0 parents commit a54f301
Show file tree
Hide file tree
Showing 19 changed files with 1,149 additions and 0 deletions.
26 changes: 26 additions & 0 deletions DESCRIPTION
@@ -0,0 +1,26 @@
Package: dbacf
Version: 0.2.8
Date: 2023-06-27
Title: Autocovariance Estimation via Difference-Based Methods
Authors@R: person("Inder", "Tecuapetla-Gómez",
email = "itecuapetla@conabio.gob.mx",
role = c("aut", "cre"))
Author: Inder Tecuapetla-Gómez [aut, cre]
Maintainer: Inder Tecuapetla-Gómez
<itecuapetla@conabio.gob.mx>
Description: Provides methods for (auto)covariance/correlation function estimation
in change point regression with stationary errors circumventing the pre-estimation
of the underlying signal of the observations. Generic, first-order, (m+1)-gapped,
difference-based autocovariance function estimator is based on M. Levine and I. Tecuapetla-Gómez (2023) <doi:10.48550/arXiv.1905.04578>. Bias-reducing, second-order, (m+1)-gapped,
difference-based estimator is based on I. Tecuapetla-Gómez and A. Munk (2017)
<doi:10.1111/sjos.12256>. Robust autocovariance estimator for change point regression with autoregressive errors is based on S. Chakar et al. (2017) <doi:10.3150/15-BEJ782>.
It also includes a general projection-based method for covariance matrix estimation.
License: GPL (>= 2)
Encoding: UTF-8
Imports: Matrix
Depends: R (>= 2.15.3)
NeedsCompilation: no
Packaged: 2023-06-27 18:18:26 UTC; itecuapetla
RoxygenNote: 7.2.0
Repository: CRAN
Date/Publication: 2023-06-29 14:30:16 UTC
18 changes: 18 additions & 0 deletions MD5
@@ -0,0 +1,18 @@
203b21e30195ed237523821dbd8e77da *DESCRIPTION
fc11d1f112f8f16b55466037def1a9f5 *NAMESPACE
7f51d8ae5b38847b82d76f8b052fbdcf *R/auxFunctions.R
3f8a578850de6857e49781b9bb9b95ce *R/dbacf-package.R
b6f72fd6450caf35746ca35759e546bc *R/dbacf.R
f2d696c41ef4f5c08249eb90588e42ce *R/dbacf_AR1.R
ec1232f0b34c6267d5186ed07d6d0b1a *R/nearPDToeplitz.R
a41ebc5ef22af87a4ca8a15eec44b34f *R/projectToeplitz.R
58c5a81da30de3a1a717645a3556f0e3 *R/symBandedToeplitz.R
f5b9c678eb53ee6a35808013c10ca91e *man/dbacf-package.Rd
42730fe590423206643591ed3b953e85 *man/dbacf.Rd
87aa6dc768e9e93f5158a79d24d62797 *man/dbacf_AR1.Rd
19f04c8173e68b22c572f32142f2b097 *man/nearPDToeplitz.Rd
ea587fc714ae32951b084a6deca1420d *man/plot.dbacf.Rd
1a68b4bd231611b9f609dbee2fb59f6a *man/projectToeplitz.Rd
798653164b5d501e5ba3bc08a5d9d7d4 *man/symBandedToeplitz.Rd
a2c960adc8a8118340cfd9970f525eab *tests/test.dbacf.R
c726251372f60e62e4d74c012a7ff182 *tests/test.projections.R
10 changes: 10 additions & 0 deletions NAMESPACE
@@ -0,0 +1,10 @@
# Generated by roxygen2: do not edit by hand

S3method(plot,dbacf)
export(dbacf)
export(dbacf_AR1)
export(nearPDToeplitz)
export(projectToeplitz)
export(symBandedToeplitz)
importFrom(graphics,abline)
importFrom(stats,median)
100 changes: 100 additions & 0 deletions R/auxFunctions.R
@@ -0,0 +1,100 @@
#-----------------------------------------------------------------------
# implements Eq. (1.5) of Munk and Tecuapetla-Gómez (2015)
dbe <- function(h, data) {
mean(diff(data, h)^2) / 2
}
#-----------------------------------------------------------------------
# implements Eq. (1.6) of Munk and Tecuapetla-Gómez (2015)
gammaZeroEst <- function(m, d, data) {
n <- length(data)
sum((data[1:(n - 2 * (m + 1))] - (1 + d) * data[(m + 2):(n - m - 1)] + d * data[(2 * m + 3):n])^2) /
(2 * (1 + d + d^2) * (n - 2 * (m + 1)))
}
#
#-----------------------------------------------------------------------
# implements Eq. (1.7) of Munk and Tecuapetla-Gómez (2015)
gh <- function(m, h, d, data) {
ifelse(h <= m, gammaZeroEst(m=m, d=d, data=data) - dbe(h=h, data=data), 0)
}
#

gh_first <- function(m, h, data) {
ifelse(h <= m, dbe(h=m+1, data=data) - dbe(h=h, data=data), 0)
}


#-----------------------------------------------------------------------
# computes first row of projection onto space of Toeplitz matrices
# see section 2.2.1 of Munk and Tecuapetla-Gómez (2015)
computeFirstRow <- function(matrix) {
n <- ncol(matrix)
firstRow <- numeric(n)

for (i in 1:n) {
for (j in 1:(n + 1 - i)) {
firstRow[i] <- firstRow[i] + matrix[j, j - 1 + i]
}
firstRow[i] <- firstRow[i] / (n + 1 - i)
}

firstRow
}
#
#-----------------------------------------------------------------------
# tests a numeric for an integer value
is.wholenumber <- function(x, tol = .Machine$double.eps^0.5) abs(x - round(x)) < tol

# -----------------------------------------------------------------------------
# --- Added on June 23, 2023
dbacfSecondOrder <- function(data, m, d){

if (missing(d)) {
d <- numeric(m + 1L)
d[1] <- 1
if (m > 0) {
for (i in 1:m) {
aux <- i^2 - 4 * (m + 1 - i)^2
if (aux >= 0) {
d[i + 1L] <- (i + sqrt(aux)) / (2 * (m + 1 - i))
} else {
d[i + 1L] <- 1
}
}
}
} else {
if (length(d) != 1L && length(d) != m + 1L) {
stop("length of d must be 1 or m + 1")
}

if ( !all(is.finite(d))) {
stop("d must be a finite numeric")
}

if (length(d) == 1L) {
d <- rep(d, m + 1L)
}
}

acf <- numeric(m)
if (m > 0) {
for (h in 1:m)
acf[h] <- gh(m=m, h=h, d=d[h + 1L], data=data)
}
acf <- c(gammaZeroEst(m=m, d=d[1], data=data), acf)

list(acf=acf, d=d)
}

# ----------------------------------------------------------------------------
# --- Added on June 23, 2023
dbacfFirstOrder <- function(data, m){
acf <- numeric(m)
if (m > 0) {
for (h in 1:m)
acf[h] <- gh_first(m=m+1, h=h, data=data) #dbe(k=h,data=data)
}
acf <- c(dbe(h=m+1, data=data), acf)

list(acf=acf, d=NULL)
}

55 changes: 55 additions & 0 deletions R/dbacf-package.R
@@ -0,0 +1,55 @@
#' Autocovariance function estimation via difference-based methods
#'
#' Difference-based (auto)covariance/correlation estimation in change point
#' regression with stationary errors.
#'
#' Provides bias-reducing methods for (auto)covariance-correlation
#' estimation in change point regression with stationary \eqn{m}-dependent errors
#' without having to pre-estimate the underlying signal of the observations.
#' In the same spirit, provides a robust estimator of the autorregressive coefficient
#' in change point regression with stationary, \eqn{AR(1)} errors.
#' It also includes a general projection-based method for covariance matrix estimation.
#'
#' @name dbacf-package
#' @author Tecuapetla-Gómez, I. \email{itecuap@@conabio.gob.mx}
#'
#' @section Autocovariance Estimation:
#' \code{\link{dbacf}} returns \emph{and plots by default} (auto)covariance/correlation
#' estimates without pre-estimating the underlying \emph{not necessarily smooth}
#' signal of observations with \emph{stationary} \eqn{m}-dependent errors. The corresponding
#' plot method \code{\link[dbacf]{plot.dbacf}} allows for adjusting graphical
#' parameters to users' liking. This method is based on \code{\link[stats]{plot.acf}}.
#'
#' \code{\link{dbacf_AR1}} returns (auto)covariance/correlation estimates while
#' circumventing the difficult estimation of the underlying change point regression
#' function from observations with stationary \eqn{AR(1)} errors.
#'
#' @section Covariance Matrix Estimation:
#' Given a matrix estimate, \emph{not necessarily positive definite}, of
#' the covariance matrix of a stationary process,
#' \code{\link{nearPDToeplitz}} returns the nearest, \emph{in the Frobenius norm},
#' covariance matrix to the initial estimate. See \code{\link{projectToeplitz}}
#' for the projection of a given symmetric matrix onto the space of Toeplitz matrices.
#' See also \code{\link{symBandedToeplitz}} for creating a (stationary process'
#' large covariance) matrix by specifying its dimension and values of its
#' autocovariance function.
#'
#' @references Grigoriadis, K.M., Frazho, A., Skelton, R. (1994).
#' \emph{Application of alternating convex projection methods for computation
#' of positive Toeplitz matrices}, IEEE Transactions on signal processing
#' \bold{42(7)}, 1873--1875.
#'
#' @references N. Higham (2002). \emph{Computing the nearest correlation matrix - a
#' problem from finance}, Journal of Numerical Analysis \bold{22}, 329--343.
#'
#' @references Tecuapetla-Gómez, I and Munk, A. (2017). \emph{Autocovariance
#' estimation in regression with a discontinuous signal and \eqn{m}-dependent errors: A
#' difference-based approach}. Scandinavian Journal of Statistics, \bold{44(2)}, 346--368.
#'
#' @references Levine, M. and Tecuapetla-Gómez, I. (2023). \emph{Autocovariance
#' function estimation via difference schemes for a semiparametric change point model
#' with \eqn{m}-dependent errors}. Submitted.
#'
#' @keywords package
#' @docType package
NULL
156 changes: 156 additions & 0 deletions R/dbacf.R
@@ -0,0 +1,156 @@
#' Difference-based (auto)covariance/correlation function estimation
#'
#' Computes \emph{and by default plots} the (auto)covariance/correlation function
#' estimate without pre-estimating the underlying \emph{piecewise constant signal}
#' of the observations. To that end, a class of second-order
#' \emph{difference-based estimators} is implemented according to Eqs.(2.5)-(2.6)
#' of \cite{Tecuapetla-Gómez and Munk (2017)}. By default, this function computes
#' a subclass of estimates with minimal bias according to Eqs.(2.12)-(2.14) of the
#' aforementioned paper.
#'
#' @param data numeric vector or a univariate object of class
#' \code{\link[stats]{ts}} of length at least \code{2(m + 1)}.
#' @param m integer scalar giving the underlying level of dependency.
#' @param d numeric vector giving the weights used in difference-based
#' estimation method. Only pertinent when \code{order=second}.
#' If missing, the weights \code{d} are calculated according
#' to Eqs.(2.12)-(2.14) of \cite{Tecuapetla-Gómez and Munk (2017)}.
#' When a single value \eqn{d^\ast}{d*} is specified,
#' \code{d = rep(}\eqn{d^\ast}{d*}\code{, m + 1)}.
#' @param type character string specifying whether covariance (default)
#' or correlation must be computed.
#' @param order character specifying whether a \code{first} (default)
#' or a \code{second} difference-based estimate should be employed.
#' @param plot logical. If \code{TRUE} (default) the acf is plotted.
#' @param \dots further arguments passed to \code{\link{plot.dbacf}}.
#'
#' @note Although the theoretical properties of the methods implemented
#' in this function were derived for change point regression with stationary
#' \emph{Gaussian} \eqn{m}-dependent errors, these methods have proven robust against
#' non-normality of the errors and as efficient as other methods in which
#' pre-estimation of an underlying smooth signal is required. For further
#' details see Section 6 of \cite{Tecuapetla-Gómez and Munk (2017)}.
#'
#' The first-order difference-based estimator was implemented following Eqs.(4)-(5)
#' of \cite{Levine and Tecuapetla-Gómez (2023)}. For the robustness of this estimator
#' see Section 4 of the just mentioned paper.
#'
#' @export
#'
#' @examples
#' ma2 <- arima.sim(n = 50, model = list(ma = c(0.4, -0.4), order = c(0, 0, 2)),
#' sd = 0.25)
#' dbacf(data=ma2, m = 2)
#' dbacf(data=ma2, m = 2, order="first")
#'
#' @seealso \code{\link[stats]{acf}}, \code{\link[dbacf]{plot.dbacf}}
#
#' @return An object of class "dbacf" containing:
#' \item{acf}{numeric vector of length \code{m + 1} giving estimated
#' (auto)covariance-correlation.}
#' \item{m}{integer giving underlying level of dependency.}
#' \item{d}{numeric vector containing the weights used to estimate acf.}
#' \item{acfType}{string indicating whether \code{covariance} or
#' \code{correlation} has been computed.}
#' \item{n}{integer giving \code{length(data)}.}
#' \item{series}{string with name of variable \code{data}.}
#'
#' @references Tecuapetla-Gómez, I and Munk, A. (2017). \emph{Autocovariance
#' estimation in regression with a discontinuous signal and \eqn{m}-dependent errors: A
#' difference-based approach}. Scandinavian Journal of Statistics, \bold{44(2)}, 346--368.
#'
#' @references Levine, M. and Tecuapetla-Gómez, I. (2023). \emph{Autocovariance
#' function estimation via difference schemes for a semiparametric change point model
#' with \eqn{m}-dependent errors}. Submitted.
#'
dbacf <- function(data, m, d, type = c("covariance", "correlation"),
order = c("second", "first"), plot = TRUE,
...) {

if (!is.numeric(data))
stop("'data' must be numeric")

if (length(data) > 2^31 - 1) {
stop("length of 'data' must be smaller than 2^31 - 1")
}
sampleT <- length(data)

series <- deparse(substitute(data))

if (missing(m)) {
stop("'m' must be specified")
} else {
if ( !is.wholenumber(m) || m < 0)
stop("m must be a positive integer")
}

m <- as.integer(m + 2 * .Machine$double.eps ^ 0.5)

if (length(data) <= 2L * (m + 1L))
stop("length of data must be greater than 2 * (m + 1)")

type <- match.arg(type)

order <- match.arg(order)

if(order == "first"){
output <- dbacfFirstOrder(data=data, m=m)
} else {
output <- dbacfSecondOrder(data=data, m=m, d=d)
}

if (type == "correlation"){
acf <- output$acf / output$acf[1]
} else {
acf <- output$acf
}

dbacf.out <- structure(list(acf = acf, m = m,
d = output$d, acfType = type,
n = sampleT, series = series),
class = "dbacf")

if (plot) {
plot.dbacf(dbacf.out, ...)
invisible(dbacf.out)
}
else dbacf.out
}
#
#' Plot autocovariance and autocorrelation functions
#'
#' This function returns the plot method for objects of class "dbacf".
#'
#' @param x an object of class "dbacf".
#' @param type what type of plot should be drawn. For possible types see
#' \code{\link[graphics]{plot}}.
#' @param xlab the x label of the plot.
#' @param ylab the y label of the plot.
#' @param xlim numeric vector of length 2 giving the \code{x} coordinates
#' range.
#' @param main an overall title for the plot.
#' @param ltyZeroLine type of line used to draw horizontal line passing at 0.
#' @param colZeroLine string indicating color of horizontal line passing at 0.
#' @param \dots extra arguments to be passed to plot.
#'
#' @rdname plot.dbacf
#' @method plot dbacf
#' @export
#'
#' @note \code{\link[dbacf]{dbacf}} documents the structure of objects of class "dbacf".
#'
#' @return No return value
#'
#' @importFrom graphics abline
#'
#' @seealso \code{\link[stats]{acf}}, \code{\link[dbacf]{dbacf}}.
#'
plot.dbacf <- function(x, type = "h", xlab = "Lag",
ylab = paste("ACF", ifelse(x$acfType == "covariance",
"(cov)", " ")),
xlim = c(0, x$m + 1), main = paste("Series", x$series),
ltyZeroLine = 3, colZeroLine = "blue", ...) {
plot(0:x$m, x$acf, type = type, xlab = xlab,
ylab = ylab, xlim = xlim, main = main, ...)
abline(h = 0, lty = ltyZeroLine, col = colZeroLine)
}

0 comments on commit a54f301

Please sign in to comment.