Skip to content

Commit

Permalink
Updated wtc function to handle user-provided lag1 coefficients
Browse files Browse the repository at this point in the history
  • Loading branch information
tgouhier committed Jul 12, 2016
1 parent bea1161 commit fc06ac9
Show file tree
Hide file tree
Showing 5 changed files with 48 additions and 36 deletions.
4 changes: 2 additions & 2 deletions DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,8 +1,8 @@
Package: biwavelet
Type: Package
Title: Conduct Univariate and Bivariate Wavelet Analyses
Version: 0.20.8
Date: 2016-06-25
Version: 0.20.9
Date: 2016-07-12
Author: Tarik C. Gouhier, Aslak Grinsted, Viliam Simko
Maintainer: Tarik C. Gouhier <[email protected]>
Description: This is a port of the WTC MATLAB package written by Aslak Grinsted
Expand Down
53 changes: 28 additions & 25 deletions R/wtc.R
Original file line number Diff line number Diff line change
@@ -1,11 +1,11 @@
#' Compute wavelet coherence
#'
#'
#' @author Tarik C. Gouhier (tarik.gouhier@@gmail.com)
#'
#'
#' Code based on WTC MATLAB package written by Aslak Grinsted.
#'
#' @param d1 time series 1 in matrix format (\code{n} rows x 2 columns). The
#' first column should contain the time steps and the second column should
#'
#' @param d1 time series 1 in matrix format (\code{n} rows x 2 columns). The
#' first column should contain the time steps and the second column should
#' contain the values.
#' @param d2 time series 2 in matrix format (\code{n} rows x 2 columns). The
#' first column should contain the time steps and the second column should
Expand All @@ -27,14 +27,14 @@
#' 2, then do a scale-average test.
#' @param nrands number of Monte Carlo randomizations. Default is 300.
#' @param quiet Do not display progress bar. Default is \code{FALSE}
#'
#'
#' @return Return a \code{biwavelet} object containing:
#' \item{coi}{matrix containg cone of influence}
#' \item{wave}{matrix containing the cross-wavelet transform}
#' \item{wave.corr}{matrix containing the bias-corrected cross-wavelet transform
#' using the method described by \code{Veleda et al. (2012)}}
#' \item{power}{matrix of power}
#' \item{power.corr}{matrix of bias-corrected cross-wavelet power using the method described
#' \item{power.corr}{matrix of bias-corrected cross-wavelet power using the method described
#' by \code{Veleda et al. (2012)}}
#' \item{rsq}{matrix of wavelet coherence}
#' \item{phase}{matrix of phases}
Expand All @@ -54,40 +54,40 @@
#'
#' @references
#' Cazelles, B., M. Chavez, D. Berteaux, F. Menard, J. O. Vik, S. Jenouvrier,
#' and N. C. Stenseth. 2008. Wavelet analysis of ecological time series.
#' and N. C. Stenseth. 2008. Wavelet analysis of ecological time series.
#' \emph{Oecologia} 156:287-304.
#'
#' Grinsted, A., J. C. Moore, and S. Jevrejeva. 2004. Application of the cross
#' wavelet transform and wavelet coherence to geophysical time series.
#'
#' Grinsted, A., J. C. Moore, and S. Jevrejeva. 2004. Application of the cross
#' wavelet transform and wavelet coherence to geophysical time series.
#' \emph{Nonlinear Processes in Geophysics} 11:561-566.
#'
#' Torrence, C., and G. P. Compo. 1998. A Practical Guide to Wavelet Analysis.
#'
#' Torrence, C., and G. P. Compo. 1998. A Practical Guide to Wavelet Analysis.
#' \emph{Bulletin of the American Meteorological Society} 79:61-78.
#'
#'
#' Torrence, C., and P. J. Webster. 1998. The annual cycle of persistence in the
#' El Nino/Southern Oscillation. \emph{Quarterly Journal of the Royal
#' Meteorological Society} 124:1985-2004.
#'
#'
#' Veleda, D., R. Montagne, and M. Araujo. 2012. Cross-Wavelet Bias Corrected by
#' Normalizing Scales. \emph{Journal of Atmospheric and Oceanic Technology}
#' 29:1401-1408.
#'
#'
#' @note The Monte Carlo randomizations can be extremely slow for large
#' datasets. For instance, 1000 randomizations of a dataset consisting of 1000
#' samples will take ~30 minutes on a 2.66 GHz dual-core Xeon processor.
#'
#' @examples
#' t1 <- cbind(1:100, rnorm(100))
#' t2 <- cbind(1:100, rnorm(100))
#'
#'
#' ## Wavelet coherence
#' wtc.t1t2 <- wtc(t1, t2, nrands = 10)
#'
#'
#' ## Plot wavelet coherence and phase difference (arrows)
#' ## Make room to the right for the color bar
#' par(oma = c(0, 0, 0, 1), mar = c(5, 4, 4, 5) + 0.1)
#' plot(wtc.t1t2, plot.cb = TRUE, plot.phase = TRUE)
#'
#'
#' @export
wtc <- function(d1, d2, pad = TRUE, dj = 1 / 12, s0 = 2 * dt,
J1 = NULL, max.scale = NULL, mother = "morlet",
Expand All @@ -111,18 +111,21 @@ wtc <- function(d1, d2, pad = TRUE, dj = 1 / 12, s0 = 2 * dt,
J1 <- round(log2(max.scale / s0) / dj)
}

# Get AR(1) coefficients for each time series
d1.ar1 <- arima(d1[,2], order = c(1, 0, 0))$coef[1]
d2.ar1 <- arima(d2[,2], order = c(1, 0, 0))$coef[1]
if (is.null(lag1)) {
# Get AR(1) coefficients for each time series
d1.ar1 <- arima(d1[,2], order = c(1, 0, 0))$coef[1]
d2.ar1 <- arima(d2[,2], order = c(1, 0, 0))$coef[1]
lag1 <- c(d1.ar1, d2.ar1)
}

# Get CWT of each time series
wt1 <- wt(d = d1, pad = pad, dj = dj, s0 = s0, J1 = J1,
max.scale = max.scale, mother = mother, param = param,
sig.level = sig.level, sig.test = sig.test, lag1 = lag1)
sig.level = sig.level, sig.test = sig.test, lag1 = lag1[1])

wt2 <- wt(d = d2, pad = pad, dj = dj, s0 = s0, J1 = J1,
max.scale = max.scale, mother = mother, param = param,
sig.level = sig.level, sig.test = sig.test, lag1 = lag1)
sig.level = sig.level, sig.test = sig.test, lag1 = lag1[2])

# Standard deviation for each time series
d1.sigma <- sd(d1[,2], na.rm = T)
Expand Down Expand Up @@ -156,7 +159,7 @@ wtc <- function(d1, d2, pad = TRUE, dj = 1 / 12, s0 = 2 * dt,
# Phase difference
phase <- atan2(Im(CW), Re(CW))
if (nrands > 0) {
signif <- wtc.sig(nrands = nrands, lag1 = c(d1.ar1, d2.ar1),
signif <- wtc.sig(nrands = nrands, lag1 = lag1,
dt = dt, n, pad = pad, dj = dj, J1 = J1, s0 = s0,
max.scale = max.scale, mother = mother,
sig.level = sig.level, quiet = quiet)
Expand Down
4 changes: 2 additions & 2 deletions inst/CITATION
Original file line number Diff line number Diff line change
Expand Up @@ -2,12 +2,12 @@ bibentry("Manual",
title="biwavelet: Conduct univariate and bivariate wavelet analyses",
author="Tarik C. Gouhier, Aslak Grinsted, Viliam Simko",
year="2016",
note="(Version 0.20.8)",
note="(Version 0.20.9)",
url="http://github.com/tgouhier/biwavelet",

textVersion =
paste("T.C. Gouhier, A. Grinstead and V. Simko (2016). ",
"biwavelet: Conduct univariate and bivariate wavelet analyses (Version 0.20.7). ",
"biwavelet: Conduct univariate and bivariate wavelet analyses (Version 0.20.9). ",
"Available from http://github.com/tgouhier/biwavelet",
sep=""),

Expand Down
9 changes: 9 additions & 0 deletions inst/NEWS.Rd
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,15 @@
\title{News for Package 'biwavelet'}
\encoding{UTF-8}

\section{Changes in biwavelet version 0.20.9 (2016-07-12)}{
\subsection{fixed}{
\itemize{
\item Fixed handling of \code{lag1} coefficients in \code{wtc}.
}
}
}


\section{Changes in biwavelet version 0.20.8 (2016-06-25)}{
\subsection{fixed}{
\itemize{
Expand Down
14 changes: 7 additions & 7 deletions man/wtc.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

0 comments on commit fc06ac9

Please sign in to comment.