[Dplr-commits] r833 - in pkg/dplR: . vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Wed Apr 23 11:34:09 CEST 2014
Author: mvkorpel
Date: 2014-04-23 11:34:08 +0200 (Wed, 23 Apr 2014)
New Revision: 833
Modified:
pkg/dplR/DESCRIPTION
pkg/dplR/vignettes/intro-dplR.Rnw
pkg/dplR/vignettes/timeseries-dplR.Rnw
pkg/dplR/vignettes/xdate-dplR.Rnw
Log:
In DESCRIPTION: Added "forecast" to Suggests because the package is
needed for a complete build of the timeseries vignette.
In vignettes:
* Fixed some typos
* Did some formatting
- added spaces
- inserted \code where missing
- wrapped some long code lines
- added second # to comment lines (Emacs likes that)
* The timeseries vignette can now be built without the "forecast"
package. If the package is missing, the following note is shown in
the document instead of the auto.arima() example: An example was
dropped because "forecast" is not available.
* Moved around various sections of the timeseries vignette around the
auto.arima() example. This was done to facilitate the conditional
inclusion of the example: it is easier to exclude one big chunk of
an example than have the conditional text and figure intertwined
with contents that are always included.
Modified: pkg/dplR/DESCRIPTION
===================================================================
--- pkg/dplR/DESCRIPTION 2014-04-22 23:48:13 UTC (rev 832)
+++ pkg/dplR/DESCRIPTION 2014-04-23 09:34:08 UTC (rev 833)
@@ -3,7 +3,7 @@
Type: Package
Title: Dendrochronology Program Library in R
Version: 1.6.0
-Date: 2014-04-21
+Date: 2014-04-23
Authors at R: c(person("Andy", "Bunn", role = c("aut", "cph",
"cre", "trl"), email = "andy.bunn at wwu.edu"), person("Mikko",
"Korpela", role = c("aut", "trl")), person("Franco", "Biondi",
@@ -21,7 +21,7 @@
Imports: gmp (>= 0.5-2), graphics, grDevices, grid, stats, utils,
digest (>= 0.2.3), lattice (>= 0.13-6), stringr (>= 0.4), XML
(>= 2.1-0)
-Suggests: foreach, iterators, RUnit (>= 0.4.25)
+Suggests: foreach, forecast, iterators, RUnit (>= 0.4.25)
Description: This package contains functions for performing tree-ring
analyses, IO, and graphics.
LazyData: no
Modified: pkg/dplR/vignettes/intro-dplR.Rnw
===================================================================
--- pkg/dplR/vignettes/intro-dplR.Rnw 2014-04-22 23:48:13 UTC (rev 832)
+++ pkg/dplR/vignettes/intro-dplR.Rnw 2014-04-23 09:34:08 UTC (rev 833)
@@ -20,7 +20,7 @@
\SweaveOpts{concordance=TRUE}
\SweaveOpts{strip.white=true}
\SweaveOpts{include=FALSE}
-<<echo=false>>=
+<<echo=FALSE,results=hide>>=
options(width=62) # width of paper (number of characters)
options(useFancyQuotes=FALSE) # fancy quotes not included in fixed-width font?
Sys.setenv(LANGUAGE="en") # no translations to languages other than English
@@ -128,7 +128,7 @@
negative exponential curves, and so on. There are also dplR functions for
less commonly used detrending methods like regional curve standardization.
-\textbf{By the way, if this is all new to you - you should stop reading this
+\textbf{By the way, if this is all new to you -- you should stop reading this
vignette and proceed immediately to a good primer on dendrochronology like
\cite{Fritts2001}. This vignette is not intended to teach you about how to do
tree-ring analysis. It's intended to teach you how to use the package.}
@@ -263,7 +263,7 @@
<<>>=
ca533.rho <- interseries.cor(ca533.rwi, prewhiten=TRUE,
method="spearman")
-ca533.rho[1:5,]
+ca533.rho[1:5, ]
mean(ca533.rho[, 1])
@
Modified: pkg/dplR/vignettes/timeseries-dplR.Rnw
===================================================================
--- pkg/dplR/vignettes/timeseries-dplR.Rnw 2014-04-22 23:48:13 UTC (rev 832)
+++ pkg/dplR/vignettes/timeseries-dplR.Rnw 2014-04-23 09:34:08 UTC (rev 833)
@@ -20,7 +20,7 @@
\SweaveOpts{concordance=TRUE}
\SweaveOpts{strip.white=true}
\SweaveOpts{include=FALSE}
-<<echo=false>>=
+<<echo=FALSE,results=hide>>=
options(width=62) # width of paper (number of characters)
options(useFancyQuotes=FALSE) # fancy quotes not included in fixed-width font?
Sys.setenv(LANGUAGE="en") # no translations to languages other than English
@@ -82,7 +82,7 @@
mean(co021.sum$stdev)
mean(co021.sum$median)
mean(co021.sum$ar1)
-mean(interseries.cor(co021)[,1])
+mean(interseries.cor(co021)[, 1])
plot(co021, plot.type="spag")
@
\begin{figure}[h]
@@ -92,7 +92,7 @@
\label{fig:rwl}
\end{figure}
-\textbf{By the way, if this is all new to you - you should
+\textbf{By the way, if this is all new to you -- you should
proceed immediately to a good primer on dendrochronology like
\cite{Fritts2001}. This vignette is not intended to teach you about how to do
tree-ring analysis. It's intended to teach you how to use the package.}
@@ -101,9 +101,9 @@
each series with a frequency response of 50\% at a wavelength of 2/3 of
each series's length. The chronology is plotted in Figure~\ref{fig:crn.plot}.
<<b, fig=TRUE>>=
-co021.rwi <- detrend(co021,method="Spline")
-co021.crn <- chron(co021.rwi,prefix="MES")
-plot(co021.crn,add.spline=TRUE,nyrs=64)
+co021.rwi <- detrend(co021, method="Spline")
+co021.crn <- chron(co021.rwi, prefix="MES")
+plot(co021.crn, add.spline=TRUE, nyrs=64)
@
\begin{figure}[h]
\centering
@@ -112,13 +112,13 @@
\label{fig:crn.plot}
\end{figure}
-The \code{co021.crn} obejct has two columns, the first giving the chronology
+The \code{co021.crn} object has two columns, the first giving the chronology
and the second the sample depth during that year. We'll 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]
-par(mfcol=c(1,2))
+dat <- co021.crn[, 1]
+par(mfcol=c(1, 2))
acf(dat)
pacf(dat)
@
@@ -130,42 +130,56 @@
\end{figure}
The ACF function indicates significant autocorrelation out to a lag of
about 10 years (which is not uncommon in tree-ring data) while the PACF plot
-suggests that the persitence after lag 4 is due to the propagation of the
+suggests that the persistence after lag 4 is due to the propagation of the
autocorrelation at earlier lags (Figure~\ref{fig:acf.plot}). And one could
very well argue that the best model here is an AR(2) model given the
marginal significance of the PACF value at lags 3 and 4. After all, you can get
three opinions by asking one statistician to look a time series. But we digress.
We now have the first bit of solid information about the time-series
-properties of these data, it looks like they fits and AR(4) model.
+properties of these data, it looks like they fit an AR(4) model.
But, R being R, there are many other ways to check this.
-The easiest way to to use the \code{ar} function which fits an autoregressive
-and selects the order by AIC. We can do the same sort of automatically fitting
-an ARMA model using the \code{auto.arima} function in the library
-\code{forecast}.
+The easiest way is to use the \code{ar} function which fits an autoregressive
+model and selects the order by AIC.
<<>>=
dat.ar <- ar(dat)
dat.ar
-require(forecast)
-dat.arima <- auto.arima(dat,ic="bic")
-@
-Indeed, \code{ar} produces an AR(4) model. Yet, auto.arima went for an
-ARMA(1,1) model - or an ARIMA(1,0,1). The parsimony priciple certainyl 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
-residials are quite clean as we'd expect Figure~\ref{fig:acf.resid}.
+@
+Indeed, \code{ar} produces an AR(4) model.
+We can do the same sort of analysis by automatically fitting
+an ARMA model using the \code{auto.arima} function in the package
+\code{"forecast"}.
+\newif\ifforecastUsable% Define boolean variable
+<<echo=FALSE,results=tex>>=
+## Test if forecast can be loaded
+if (require("forecast", character.only = TRUE)) {
+ cat("\\forecastUsabletrue\n\n")# output to LaTeX
+}
+@
+\ifforecastUsable% Conditional: If "forecast" is available
<<d, fig=TRUE>>=
-summary(dat.arima)
-head(residuals(dat.arima))
-coef(dat.arima)
-acf(residuals(dat.arima))
+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))
+}
@
\begin{figure}[h]
-\centering
-\includegraphics{timeseries-dplR-d}
-\caption{ACF plot of the ARIMA(1,1) residuals.}
-\label{fig:acf.resid}
+ \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'd expect (Figure~\ref{fig:acf.resid}).
+\else% If "forecast" is not available
+An example was dropped because \code{"forecast"} is not available.
+\fi% End of conditional
\bibliography{dplR}
Modified: pkg/dplR/vignettes/xdate-dplR.Rnw
===================================================================
--- pkg/dplR/vignettes/xdate-dplR.Rnw 2014-04-22 23:48:13 UTC (rev 832)
+++ pkg/dplR/vignettes/xdate-dplR.Rnw 2014-04-23 09:34:08 UTC (rev 833)
@@ -20,7 +20,7 @@
\SweaveOpts{concordance=TRUE}
\SweaveOpts{strip.white=true}
\SweaveOpts{include=FALSE}
-<<echo=false>>=
+<<echo=FALSE,results=hide>>=
options(width=62) # width of paper (number of characters)
options(useFancyQuotes=FALSE) # fancy quotes not included in fixed-width font?
Sys.setenv(LANGUAGE="en") # no translations to languages other than English
@@ -84,7 +84,7 @@
mean(dat.sum$stdev)
mean(dat.sum$median)
mean(dat.sum$ar1)
-mean(interseries.cor(dat)[,1])
+mean(interseries.cor(dat)[, 1])
plot(dat, plot.type="spag")
@
\begin{figure}[h]
@@ -94,7 +94,7 @@
\label{fig:rwl.plot}
\end{figure}
-\textbf{By the way, if this is all new to you - you should
+\textbf{By the way, if this is all new to you -- you should
proceed immediately to a good primer on dendrochronology like
\cite{Fritts2001}. This vignette is not intended to teach you about how to do
tree-ring analysis. It's intended to teach you how to use the package.}
@@ -107,14 +107,15 @@
pick a random year in the core to give us a bit of a challenge in finding it.
<<>>=
-# create a missing ring by deleting a random year of
-# growth in a random series
+## create a missing ring by deleting a random year of
+## growth in a random series
+RNGversion("2.15.0")
set.seed(4576)
-i <- sample(x=1:nrow(dat),size=1)
-j <- sample(x=1:ncol(dat),size=1)
-tmp <- dat[,j]
-tmp <- c(NA,tmp[-i])
-dat[,j] <- tmp
+i <- sample(x=nrow(dat), size=1)
+j <- sample(x=ncol(dat), size=1)
+tmp <- dat[, j]
+tmp <- c(NA, tmp[-i])
+dat[, j] <- tmp
@
We've now deleted the $i^{th}$ observation from the $j^{th}$ core while making sure that
\code{dat} still has the appropriate numbers of rows. By sticking the NA at the
@@ -126,7 +127,7 @@
each tree-ring series and a master chronology built from all the other series
in the rwl object (leave-one-out principle). These correlations are calculated
on overlapping segments (e.g., 50-year segments would be overlapped by
-25-years). By default, each of the series is filtered to remove low-frequency
+25 years). By default, each of the series is filtered to remove low-frequency
variation prior to the correlation analysis. The help file has abundant
details. Here will will look at overlapping 60 year segments. A plot is
produced by default with \code{corr.rwl.seg} (Figure~\ref{fig:corr.rwl.plot}).
@@ -142,7 +143,7 @@
the time period and thus have no correlations calculated. Our modified data
set indicates one series with dating problems.
<<b, fig=TRUE>>=
-rwl.60 <- corr.rwl.seg(dat,seg.length=60,pcrit=0.01)
+rwl.60 <- corr.rwl.seg(dat, seg.length=60, pcrit=0.01)
@
\begin{figure}[h]
\centering
@@ -153,23 +154,24 @@
\label{fig:corr.rwl.plot}
\end{figure}
-The low correlation between series ``643114'' and the master indicates a dating
+The low correlation between series \code{643114} and the master indicates a dating
problem (Figure~\ref{fig:corr.rwl.plot}). Now that we suspect a dating problem,
let's take closer look at this problem child. Figure~\ref{fig:corr.seg.plot}
-shows that series 643114 begins to lose correlation with the master at the
+shows that series \code{643114} begins to lose correlation with the master at the
end of the 19th century.
<<c, fig=TRUE>>=
-# look at this series with a running correlation
-seg.60 <- corr.series.seg(rwl=dat,series="643114",seg.length=60)
+## look at this series with a running correlation
+seg.60 <- corr.series.seg(rwl=dat, series="643114",
+ seg.length=60)
@
\begin{figure}[h]
\centering
\includegraphics{xdate-dplR-c}
-\caption{Correlations between series 643114 and the master chronology are
+\caption{Correlations between series \code{643114} and the master chronology are
shown with horizontal lines according (60-year segments lagged by 30 years).
A centered running correlation with a length of 60 years complements the
-segment correlations. The critical level is show with a dashed line.}
+segment correlations. The critical level is shown with a dashed line.}
\label{fig:corr.seg.plot}
\end{figure}
@@ -182,9 +184,9 @@
<<d, fig=TRUE>>=
win <- 1800:1960
dat.yrs <- as.numeric(rownames(dat))
-dat.trunc <- dat[dat.yrs%in%win,]
-ccf.30 <- ccf.series.rwl(rwl=dat.trunc,series="643114",
- seg.length=30,bin.floor=50)
+dat.trunc <- dat[dat.yrs %in% win, ]
+ccf.30 <- ccf.series.rwl(rwl=dat.trunc, series="643114",
+ seg.length=30, bin.floor=50)
@
\begin{figure}[h]
@@ -192,7 +194,7 @@
\includegraphics{xdate-dplR-d}
\caption{Cross-correlations between series \code{643114} and the master
chronology are shown for each segment (30-year segments lagged by 15 years).
-The series correlates well at lag 0 until the 1865-1894 bin and then at
+The series correlates well at lag 0 until the 1865--1894 bin and then at
lag +1 prior to 1865.}
\label{fig:ccf.plot}
\end{figure}
@@ -202,7 +204,7 @@
\code{643114} is missing a ring as it better correlates to the master
chronology with a one-year offset. \footnote{As of dplR version 1.60,
the cross correlations in \code{ccf.series.rwl} are calculated calling
-as \code{ccf(x=series, y=master, lag.max=lag.max, plot=FALSE)}. Note that
+\code{ccf(x=series, y=master, lag.max=lag.max, plot=FALSE)}. Note that
prior to dplR version 1.60, the \code{master} was set as \code{x} and
the \code{series} as \code{y}. This was changed to be more in line with
user expectations so that a missing ring in a series produces a positive
@@ -219,25 +221,25 @@
<<e, fig=TRUE>>=
win <- 1850:1900
-dat.trunc <- dat[dat.yrs%in%win,]
-ccf.20 <- ccf.series.rwl(rwl=dat.trunc,series="643114",
- seg.length=20,bin.floor=0)
+dat.trunc <- dat[dat.yrs %in% win, ]
+ccf.20 <- ccf.series.rwl(rwl=dat.trunc, series="643114",
+ seg.length=20, bin.floor=0)
@
\begin{figure}[h]
\centering
\includegraphics{xdate-dplR-e}
\caption{Cross-correlations between series \code{643114} and the master
-chronology at 20-year segments lagged by 10 years over 1850-1900.}
+chronology at 20-year segments lagged by 10 years over 1850--1900.}
\label{fig:ccf.plot2}
\end{figure}
By 1879 the correlation between series \code{643114} and the master is solidly at
lag +1 (Figure~\ref{fig:ccf.plot2}). The 1870 to 1889 correlation is marginal
-while the dating at 1880-1899 seems accurate (lag 0). This suggests that
+while the dating at 1880--1899 seems accurate (lag 0). This suggests that
the dating error is between 1879 and 1889.
-We have strong inference now that series 643114 is misdated somewhere in a ten
+We have strong inference now that series \code{643114} is misdated somewhere in a ten
year period around 1885. On final tool that dplR has is the ability to combine
the visual style of crossdating using skeleton plots with the statistical
approach of cross-correlation analysis. The \code{skel.ccf.plot} function does
@@ -245,8 +247,8 @@
suspected dating error (1885):
<<f, fig=TRUE>>=
-xskel.ccf.plot(rwl=dat,series="643114",win.start=1865, win.width=40)
-
+xskel.ccf.plot(rwl=dat, series="643114",
+ win.start=1865, win.width=40)
@
\begin{figure}[h]
@@ -258,9 +260,9 @@
series with the marker years annotated for the master on the bottom axis and
series \code{643114} on the top. The text at the top of the figure gives the
correlation between the series and master (green bars) as well as the percentage
-of agreement betwen the skeleton bars for the series and master. The bottom
+of agreement between the skeleton bars for the series and master. The bottom
panels show cross correlations for the first half (left) and second half of
-the time series. In this case, the early period (1865-1884) shows a mismatch
+the time series. In this case, the early period (1865--1884) shows a mismatch
of the skeleton plot by one year coupled with a strong lag +1 correlation. }
\label{fig:ccf.skel.plot}
\end{figure}
@@ -271,9 +273,9 @@
one could go to to figure out exactly where the dating problem might be
but nothing ever takes the place of looking at the sample!
-We have strong inference now that series 643114 is misdated somewhere in a ten
+We have strong inference now that series \code{643114} is misdated somewhere in a ten
year period around 1885. We have still not revealed whether this is correct
-or not. Let's look at the values for i and j and see how we did:
+or not. Let's look at the values for \code{i} and \code{j} and see how we did:
<<>>=
j
colnames(co021)[j]
@@ -281,8 +283,8 @@
rownames(co021)[i]
@
-Our sluething indicated that our dating error was around the year 1885. In
-fact, i was the year 1884. This is pretty spectacular!
+Our sleuthing indicated that our dating error was around the year 1885. In
+fact, \code{i} was the year 1884. This is pretty spectacular!
\bibliography{dplR}
More information about the Dplr-commits
mailing list