[Dplr-commits] r836 - pkg/dplR/vignettes

[noreply at r-forge.r-project.org]

Author: andybunn
Date: 2014-04-25 00:46:45 +0200 (Fri, 25 Apr 2014)
New Revision: 836

Modified:
   pkg/dplR/vignettes/intro-dplR.Rnw
   pkg/dplR/vignettes/timeseries-dplR.Rnw
Log:
* Vignette edits.

Modified: pkg/dplR/vignettes/intro-dplR.Rnw
===================================================================
--- pkg/dplR/vignettes/intro-dplR.Rnw	2014-04-24 07:25:43 UTC (rev 835)
+++ pkg/dplR/vignettes/intro-dplR.Rnw	2014-04-24 22:46:45 UTC (rev 836)
@@ -43,7 +43,7 @@
 \newpage
 
 \section{Introduction}
-\subsection{What Is Covered}
+\subsection{What is Covered}
 The Dendrochronology Program Library in R (dplR) is a package for 
 dendrochronologists to handle data processing and analysis. This 
 document gives just a brief introduction of some of the most commonly 

Modified: pkg/dplR/vignettes/timeseries-dplR.Rnw
===================================================================
--- pkg/dplR/vignettes/timeseries-dplR.Rnw	2014-04-24 07:25:43 UTC (rev 835)
+++ pkg/dplR/vignettes/timeseries-dplR.Rnw	2014-04-24 22:46:45 UTC (rev 836)
@@ -112,15 +112,20 @@
 \label{fig:crn.plot}
 \end{figure}
 
-The \code{co021.crn} object has two columns, the first giving the chronology
+\section{Characterizing the Data}
+
+Let's start with a quick exploratory data analysis into the time-series
+proThe \code{co021.crn} object has two columns, the first giving the chronology
 and the second the sample depth during that year. We will start our analysis
 on the chronology by looking at its autocorrelation structure using R's 
 \code{acf} and \code{pacf} functions.
 <<c, fig=TRUE>>=
 dat <- co021.crn[, 1]
+op <- par(no.readonly = TRUE) # Save to reset on exit
 par(mfcol=c(1, 2))
 acf(dat)
 pacf(dat)
+par(op)
 @
 \begin{figure}[h]
 \centering
@@ -157,30 +162,93 @@
 }
 @ 
 \ifforecastUsable% Conditional: If "forecast" is available
-<<d, fig=TRUE>>=
+<<>>=
 if (require("forecast", character.only = TRUE)) {
     dat.arima <- auto.arima(dat, ic="bic")
     summary(dat.arima)
     head(residuals(dat.arima))
     coef(dat.arima)
-    acf(residuals(dat.arima))
+    acf(residuals(dat.arima),plot=FALSE)
 }
 @
-\begin{figure}[h]
-  \centering
-  \includegraphics{timeseries-dplR-d}
-  \caption{ACF plot of the ARIMA(1,1) residuals.}
-  \label{fig:acf.resid}
-\end{figure}
 Instead of an AR(4) model, \code{auto.arima} went for an 
 ARMA(1,1) model -- or an ARIMA(1,0,1). The parsimony principle certainly likes 
 a nice simple ARMA(1,1) model! Note that we could look at the residuals
 (just the first few), model coefficients, etc. quite easily. And indeed the 
-residuals are quite clean as we would expect (Figure~\ref{fig:acf.resid}).
+residuals are quite clean as we would expect.
 \else% If "forecast" is not available
 An example was dropped because \code{"forecast"} is not available.
 \fi% End of conditional
+\section{Frequency Domain}
+There is, at times, and almost manic desire to better characterize the 
+spectral aspects of a tree-ring series. In dplR, we've implemented two 
+of the most common ways that dendrochronologists go about this and there 
+are a host of other approaches in R that we won't get to in this vignette.
 
+The redfit function in dplR is a port of Schulz's REDFIT (version 3.8e) 
+program and estimates the red-noise spectra of a time series. 
+<<d, fig=TRUE>>=
+
+redf.dat <- redfit(dat, nsim = 1000)
+
+par(tcl = 0.5, mar = rep(2.2, 4), mgp = c(1.1, 0.1, 0))
+
+plot(redf.dat[["freq"]], redf.dat[["gxxc"]],
+     ylim = range(redf.dat[["ci99"]], redf.dat[["gxxc"]]),
+     type = "n", ylab = "Spectrum (dB)", xlab = "Frequency (1/yr)",
+     axes = FALSE)
+grid()
+lines(redf.dat[["freq"]], redf.dat[["gxxc"]], col = "#1B9E77")
+lines(redf.dat[["freq"]], redf.dat[["ci99"]], col = "#D95F02")
+lines(redf.dat[["freq"]], redf.dat[["ci95"]], col = "#7570B3")
+lines(redf.dat[["freq"]], redf.dat[["ci90"]], col = "#E7298A")
+freqs <- pretty(redf.dat[["freq"]])
+pers <- round(1 / freqs, 2)
+axis(1, at = freqs, labels = TRUE)
+axis(3, at = freqs, labels = pers)
+mtext(text = "Period (yr)", side = 3, line = 1.1)
+axis(2); axis(4)
+legend("topright", c("dat", "CI99", "CI95", "CI90"), lwd = 2,
+       col = c("#1B9E77", "#D95F02", "#7570B3", "#E7298A"),
+       bg = "white")
+box()
+par(op)
+@
+
+
+\begin{figure}[h]
+  \centering
+  \includegraphics{timeseries-dplR-d}
+  \caption{Spectra of Mesa Verde chronology using redfit}
+  \label{fig:redfit}
+\end{figure}
+Using the Mesa Verde chronology we see that there are frequencies in that
+time series that are significantly different from a red-noise assumption
+in the interannual (<3 years) and at low frequencies (multidecadal). These
+are plotted in Figure~\ref{fig:redfit}.
+
+
+Another popular way to visualize a tree-ring chronology in the frequency 
+domain is through a continuous wavelet transform. In dplR, there is are
+functions for calculalting the transform via \code{wavelet} and plotting
+the result via \code{wavelet.plot}.
+
+<<e, fig=TRUE>>=
+yrs <- as.numeric(rownames(co021.crn))
+out.wave <- morlet(y1 = dat, x1 = yrs, p2 = 8, dj = 0.1,
+                   siglvl = 0.99)
+wavelet.plot(out.wave)
+@
+\begin{figure}[h]
+  \centering
+  \includegraphics{timeseries-dplR-e}
+  \caption{Continuous wavelet of the Mesa Verde chronology}
+  \label{fig:wavelet}
+\end{figure}
+The wavelet plot (Figure~\ref{fig:wavelet}) shows a similar story as the
+plot from \code{redfit} in Figure~\ref{fig:redfit} with significant 
+variation at interannual to multidecadal scales.
+
 \bibliography{dplR}
 
 \end{document}