[Dplr-commits] r814 - in pkg/dplR: . R man vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Apr 13 06:11:38 CEST 2014
Author: andybunn
Date: 2014-04-13 06:11:30 +0200 (Sun, 13 Apr 2014)
New Revision: 814
Added:
pkg/dplR/vignettes/xdate-dplR.Rnw
Modified:
pkg/dplR/ChangeLog
pkg/dplR/R/corr.rwl.seg.R
pkg/dplR/R/corr.series.seg.R
pkg/dplR/R/detrend.series.R
pkg/dplR/R/interseries.cor.R
pkg/dplR/R/rwi.stats.running.R
pkg/dplR/TODO
pkg/dplR/man/corr.rwl.seg.Rd
pkg/dplR/man/corr.series.seg.Rd
pkg/dplR/man/detrend.series.Rd
pkg/dplR/man/interseries.cor.Rd
pkg/dplR/man/rwi.stats.running.Rd
pkg/dplR/vignettes/intro-dplR.Rnw
Log:
* Added correlation methods to places that had no specification including rwi.stats.running (not legacy) and the crossdating functions.
* Started a crossdating vignette.
Modified: pkg/dplR/ChangeLog
===================================================================
--- pkg/dplR/ChangeLog 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/ChangeLog 2014-04-13 04:11:30 UTC (rev 814)
@@ -31,7 +31,8 @@
File: detrend.R and detrend.series.R
------------
- Added an Ar detrend method. Revised plotting in detrend.series
-- Added a verbose option to write out parameters used in detrending
+- Added a verbose option to print out parameters used in detrending
+- Added a return.info option to write a list of parameters used in detrending
File: powt.R
------------
@@ -57,11 +58,17 @@
-------------------------
- A vignette to intriduce dplR
+File: corr.series.seg.R
+--------------------
+- Added method argument to specify method for cor.test(). Defauts to
+ "spearman."
+
File: corr.rwl.seg.R
--------------------
- Removed yr.range() function in favor of yr.range() in helpers.R.
They are identical for all practical purposes.
-
+- Added method argument to specify method for cor.test(). Defauts to
+ "spearman."
File: interseries.cor.R
-------------------------
- New function interseries.cor.
@@ -125,6 +132,7 @@
cases in the running window.
- 'c.eff' in the return value is now 0 if no correlations were
computed
+- Added a method argument to change the type of correlation method performed.
- Optimizations
Files: write.compact.R, write.crn.R,
Modified: pkg/dplR/R/corr.rwl.seg.R
===================================================================
--- pkg/dplR/R/corr.rwl.seg.R 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/R/corr.rwl.seg.R 2014-04-13 04:11:30 UTC (rev 814)
@@ -1,5 +1,6 @@
corr.rwl.seg <- function(rwl, seg.length=50, bin.floor=100, n=NULL,
prewhiten = TRUE, pcrit=0.05, biweight=TRUE,
+ method = c("spearman", "pearson", "kendall"),
make.plot = TRUE, label.cex=1,
floor.plus1 = FALSE, master = NULL,
master.yrs = as.numeric(if (is.null(dim(master))) {
@@ -8,7 +9,7 @@
rownames(master)
}),
...) {
-
+ method <- match.arg(method)
## run error checks
qa.xdate(rwl, seg.length, n, bin.floor)
@@ -175,7 +176,7 @@
bin.pval <- NA
} else {
tmp <- cor.test(series[mask], master2[mask],
- method = "spearman", alternative = "greater")
+ method = method, alternative = "greater")
bin.cor <- tmp$estimate
bin.pval <- tmp$p.val
}
@@ -184,7 +185,7 @@
}
## overall correlation
tmp <- cor.test(series, master2,
- method = "spearman", alternative = "greater")
+ method = method, alternative = "greater")
overall.cor[i, 1] <- tmp$estimate
overall.cor[i, 2] <- tmp$p.val
}
Modified: pkg/dplR/R/corr.series.seg.R
===================================================================
--- pkg/dplR/R/corr.series.seg.R 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/R/corr.series.seg.R 2014-04-13 04:11:30 UTC (rev 814)
@@ -1,9 +1,11 @@
corr.series.seg <- function(rwl, series, series.yrs=as.numeric(names(series)),
seg.length=50, bin.floor=100, n=NULL,
prewhiten = TRUE, biweight=TRUE,
+ method = c("spearman", "pearson", "kendall"),
pcrit=0.05, make.plot = TRUE,
floor.plus1 = FALSE, ...) {
+ method <- match.arg(method)
## run error checks
qa.xdate(rwl, seg.length, n, bin.floor)
@@ -92,7 +94,7 @@
bin.cor <- NA
bin.pval <- NA
} else {
- tmp <- cor.test(series2[mask], master[mask], method = "spearman",
+ tmp <- cor.test(series2[mask], master[mask], method = method,
alternative = "greater")
bin.cor <- tmp$estimate
bin.pval <- tmp$p.val
@@ -101,7 +103,7 @@
res.pval[j] <- bin.pval
}
## overall correlation
- tmp <- cor.test(series2, master, method = "spearman",
+ tmp <- cor.test(series2, master, method = method,
alternative = "greater")
overall.cor[1] <- tmp$estimate
overall.cor[2] <- tmp$p.val
Modified: pkg/dplR/R/detrend.series.R
===================================================================
--- pkg/dplR/R/detrend.series.R 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/R/detrend.series.R 2014-04-13 04:11:30 UTC (rev 814)
@@ -381,7 +381,7 @@
xlab=gettext("Age (Yrs)", domain="R-dplR"),
ylab=gettext("RWI", domain="R-dplR"))
abline(h=1)
- mtext(text="Ar residuals are not plotted with raw series",side=3,line=-1)
+ mtext(text="(Not plotted with raw series)",side=3,line=-1,cex=0.75)
}
}
Modified: pkg/dplR/R/interseries.cor.R
===================================================================
--- pkg/dplR/R/interseries.cor.R 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/R/interseries.cor.R 2014-04-13 04:11:30 UTC (rev 814)
@@ -1,8 +1,8 @@
interseries.cor <- function(rwl, n=NULL, prewhiten=TRUE, biweight=TRUE,
- method = c("spearman", "pearson","kendall")) {
+ method = c("spearman", "pearson", "kendall")) {
method <- match.arg(method)
nseries <- length(rwl)
- rho <- numeric(nseries)
+ res.cor <- numeric(nseries)
p.val <- numeric(nseries)
rwl.mat <- as.matrix(rwl)
tmp <- normalize.xdate(rwl=rwl.mat, n=n,
@@ -12,9 +12,11 @@
master <- tmp[["master"]]
for (i in seq_len(nseries)) {
tmp2 <- cor.test(series[, i], master[, i],
- method = method)
- rho[i] <- tmp2[["estimate"]]
+ method = method, alternative = "greater")
+ res.cor[i] <- tmp2[["estimate"]]
p.val[i] <- tmp2[["p.value"]]
}
- data.frame(rho = rho, p.val = p.val, row.names = names(rwl))
+ res <- data.frame(res.cor = res.cor, p.val = p.val, row.names = names(rwl))
+ # change res.cor to r, rho, or tau based on method
+ res
}
Modified: pkg/dplR/R/rwi.stats.running.R
===================================================================
--- pkg/dplR/R/rwi.stats.running.R 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/R/rwi.stats.running.R 2014-04-13 04:11:30 UTC (rev 814)
@@ -2,7 +2,7 @@
### Computes the correlation coefficients between columns of x and y.
### Requires "limit" overlapping values in each pair.
-cor.with.limit <- function(limit, x, y) {
+cor.with.limit <- function(limit, x, y, method) {
n.x <- ncol(x) # caller makes sure that n.x
n.y <- ncol(y) # and n.y >= 1
r.mat <- matrix(NA_real_, n.x, n.y)
@@ -15,7 +15,8 @@
good.both <- which(good.x & good.y)
n.good <- length(good.both)
if (n.good >= limit && n.good > 0) {
- r.mat[i, j] <- cor(this.x[good.both], this.y[good.both])
+ r.mat[i, j] <- cor(this.x[good.both], this.y[good.both],
+ method = method)
}
}
}
@@ -23,7 +24,7 @@
}
### Computes the correlation coefficients between different columns of x.
-cor.with.limit.upper <- function(limit, x) {
+cor.with.limit.upper <- function(limit, x, method) {
n.x <- ncol(x) # caller makes sure that n.x >= 2
r.vec <- rep.int(NA_real_, n.x * (n.x - 1) / 2)
good.x <- !is.na(x)
@@ -34,24 +35,29 @@
k <- k + 1
good.both <- which(good.i & good.x[, j])
if (length(good.both) >= limit) {
- r.vec[k] <- cor(x[good.both, i], x[good.both, j])
+ r.vec[k] <- cor(x[good.both, i], x[good.both, j],
+ method = method)
}
}
}
r.vec
}
-rwi.stats <- function(rwi, ids=NULL, period=c("max", "common"), ...) {
+rwi.stats <- function(rwi, ids=NULL, period=c("max", "common"),
+ method = c("spearman", "pearson", "kendall"),
+ ...) {
args <- list(...)
args[["rwi"]] <- rwi
args[["ids"]] <- ids
args[["period"]] <- period
+ args[["method"]] <- method
args[["running.window"]] <- FALSE
do.call(rwi.stats.running, args)
}
### Main function, exported to user
rwi.stats.running <- function(rwi, ids=NULL, period=c("max", "common"),
+ method = c("spearman", "pearson", "kendall"),
prewhiten=FALSE,n=NULL,
running.window=TRUE,
window.length=min(50, nrow(rwi)),
@@ -61,7 +67,7 @@
round.decimals=3,
zero.is.missing=TRUE) {
period2 <- match.arg(period)
-
+ method <- match.arg(method)
if (running.window) {
if (window.length < 3) {
stop("minimum 'window.length' is 3")
@@ -248,7 +254,8 @@
i.data <- rwi3[year.idx, cores.of.tree[[i]], drop=FALSE]
for (j in (i + 1):n.trees) {
j.data <- rwi3[year.idx, cores.of.tree[[j]], drop=FALSE]
- bt.r.mat <- cor.with.limit(min.corr.overlap, i.data, j.data)
+ bt.r.mat <- cor.with.limit(min.corr.overlap, i.data, j.data,
+ method=method)
bt.r.mat <- bt.r.mat[!is.na(bt.r.mat)]
n.bt.temp <- length(bt.r.mat)
if (n.bt.temp > 0) {
@@ -270,7 +277,8 @@
n.cores.tree[i] <- 1
} else {
these.data <- rwi3[year.idx, these.cores, drop=FALSE]
- wt.r.vec <- cor.with.limit.upper(min.corr.overlap, these.data)
+ wt.r.vec <- cor.with.limit.upper(min.corr.overlap, these.data,
+ method=method)
wt.r.vec <- wt.r.vec[!is.na(wt.r.vec)]
n.wt.temp <- length(wt.r.vec)
if (n.wt.temp > 0) {
Modified: pkg/dplR/TODO
===================================================================
--- pkg/dplR/TODO 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/TODO 2014-04-13 04:11:30 UTC (rev 814)
@@ -4,9 +4,10 @@
- What is the best way to extract the parameters with ModNegExp?
- In detrend() the call to detrend.series doesn't appear to pass the names
in to detrend.series when dopar is invoked. Is that right?
+- Also note that i.detrend will need fixing too for verbose mode.
-* Write more vignettes:
-- Crossdating
+o[andybunn] Write more vignettes:
+- Crossdating (started)
- Spectral and wavelets
- Advanced chronology building (strip.rwl, etc.)
@@ -14,6 +15,9 @@
to specify which correlation method (e.g., spearman). Note that this is
implemented in interseries.cor - but will have to do this throughout.
The default should be spearman.
+- Did this in rwi.stats, corr.rwl.seg, and corr.series.seg but the
+ output names still have "rho" in them. Change?
+- This can't be easily implemented in ccf.series.rwl
* Decide when to use class('rwl') in functions dealing with rwl objects.
Other than the plot and summary S3Method for rwl, are there cases when
Modified: pkg/dplR/man/corr.rwl.seg.Rd
===================================================================
--- pkg/dplR/man/corr.rwl.seg.Rd 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/man/corr.rwl.seg.Rd 2014-04-13 04:11:30 UTC (rev 814)
@@ -7,6 +7,7 @@
\usage{
corr.rwl.seg(rwl, seg.length = 50, bin.floor = 100, n = NULL,
prewhiten = TRUE, pcrit = 0.05, biweight = TRUE,
+ method = c("spearman", "pearson","kendall"),
make.plot = TRUE, label.cex = 1, floor.plus1 = FALSE,
master = NULL,
master.yrs = as.numeric(if (is.null(dim(master))) {
@@ -33,6 +34,9 @@
the correlation test. }
\item{biweight}{ \code{logical} flag. If \code{TRUE} then a robust
mean is calculated using \code{\link{tbrm}}. }
+ \item{method}{Can either "pearson", "kendall", or "spearman" which indicates
+ the correlation coefficient is to be used. Defaults to "spearman." See
+ \code{\link{cor.test}}. }
\item{make.plot}{ \code{logical flag} indicating whether to make a
plot. }
\item{label.cex}{ \code{numeric} scalar for the series labels on the
Modified: pkg/dplR/man/corr.series.seg.Rd
===================================================================
--- pkg/dplR/man/corr.series.seg.Rd 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/man/corr.series.seg.Rd 2014-04-13 04:11:30 UTC (rev 814)
@@ -8,7 +8,9 @@
\usage{
corr.series.seg(rwl, series, series.yrs = as.numeric(names(series)),
seg.length = 50, bin.floor = 100, n = NULL,
- prewhiten = TRUE, biweight = TRUE, pcrit = 0.05,
+ prewhiten = TRUE, biweight = TRUE,
+ method = c("spearman", "pearson","kendall"),
+ pcrit = 0.05,
make.plot = TRUE, floor.plus1 = FALSE, \dots)
}
\arguments{
@@ -30,6 +32,9 @@
whitened using \code{\link{ar}}. }
\item{biweight}{ \code{logical} flag. If \code{TRUE} then a robust
mean is calculated using \code{\link{tbrm}}. }
+ \item{method}{Can either "pearson", "kendall", or "spearman" which indicates
+ the correlation coefficient is to be used. Defaults to "spearman." See
+ \code{\link{cor.test}}. }
\item{pcrit}{ a number between 0 and 1 giving the critical value for
the correlation test. }
\item{make.plot}{ \code{logical} flag indicating whether to make a
Modified: pkg/dplR/man/detrend.series.Rd
===================================================================
--- pkg/dplR/man/detrend.series.Rd 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/man/detrend.series.Rd 2014-04-13 04:11:30 UTC (rev 814)
@@ -193,15 +193,15 @@
\examples{library(stats)
## Using a plausible representation of a tree-ring series
gt <- 0.5 * exp (-0.05 * 1:200) + 0.2
-noise <- c(arima.sim(model = list(ar = 0.7), n = 200, mean = 1, sd = 0.5))
+noise <- c(arima.sim(model = list(ar = 0.7), n = 200, sd = 0.5))+2
series <- gt * noise
-series.rwi <- detrend.series(y = series, y.name = "Foo")
+series.rwi <- detrend.series(y = series, y.name = "Foo", verbose=TRUE)
## Use series CAM011 from the Campito data set
data(ca533)
series <- ca533[, "CAM011"]
names(series) <- rownames(ca533)
# defaults to all four methods
-series.rwi <- detrend.series(y = series, y.name = "CAM011")
+series.rwi <- detrend.series(y = series, y.name = "CAM011", verbose=TRUE)
# see plot with three methods
series.rwi <- detrend.series(y = series, y.name = "CAM011",
method=c("Spline", "ModNegExp","Mean"))
Modified: pkg/dplR/man/interseries.cor.Rd
===================================================================
--- pkg/dplR/man/interseries.cor.Rd 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/man/interseries.cor.Rd 2014-04-13 04:11:30 UTC (rev 814)
@@ -7,7 +7,7 @@
}
\usage{
interseries.cor(rwl,n=NULL,prewhiten=TRUE,biweight=TRUE,
- method = c("spearman", "pearson","kendall"))
+ method = c("spearman", "pearson", "kendall"))
}
\arguments{
\item{rwl}{ a \code{data.frame} with series as columns and years as
@@ -43,13 +43,14 @@
detrending can be seen with \code{\link{series.rwl.plot}}.
This function produces the same output of the "overall" portion of
- \code{\link{corr.rwl.seg}}. The mean rho value given is sometimes referred to as
- the "overall interseries correlation"" or the "COFECHA interseries
- correlation." This output differs from the \code{rbar} statistics given by
- \code{\link{rwi.stats}} in that \code{rbar} is the average pairwise correlation between
- series where this is the correlation between a series and a master chronology.
+ \code{\link{corr.rwl.seg}}. The mean correlation value given is sometimes
+ referred to as the "overall interseries correlation"" or the "COFECHA
+ interseries correlation." This output differs from the \code{rbar}
+ statistics given by \code{\link{rwi.stats}} in that \code{rbar} is
+ the average pairwise correlation between series where this is the
+ correlation between a series and a master chronology.
}
-\value{ a \code{data.frame} with rho values and p-values given from
+\value{ a \code{data.frame} with correlation values and p-values given from
\code{\link{cor.test}}
}
\author{ Andy Bunn, patched and improved by Mikko Korpela }
@@ -58,6 +59,8 @@
foo <- interseries.cor(gp.rwl)
# compare to:
# corr.rwl.seg(rwl=gp.rwl,make.plot=FALSE)$overall
+# using pearson's r
+foo <- interseries.cor(gp.rwl,method="pearson")
# two measures of interseries correlation
# compare interseries.cor to rbar from rwi.stats
Modified: pkg/dplR/man/rwi.stats.running.Rd
===================================================================
--- pkg/dplR/man/rwi.stats.running.Rd 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/man/rwi.stats.running.Rd 2014-04-13 04:11:30 UTC (rev 814)
@@ -12,6 +12,7 @@
}
\usage{
rwi.stats.running(rwi, ids = NULL, period = c("max", "common"),
+ method = c("spearman", "pearson","kendall"),
prewhiten=FALSE,n=NULL,
running.window = TRUE,
window.length = min(50, nrow(rwi)),
@@ -21,7 +22,8 @@
round.decimals = 3,
zero.is.missing = TRUE)
-rwi.stats(rwi, ids=NULL, period=c("max", "common"), \dots)
+rwi.stats(rwi, ids=NULL, period=c("max", "common"),
+ method = c("spearman", "pearson","kendall"), \dots)
rwi.stats.legacy(rwi, ids=NULL, period=c("max", "common"))
}
@@ -42,6 +44,9 @@
complete observations over the period common to all cores (i.e. rows
common to all samples) or the maximum pairwise overlap. Defaults
to \code{"max"}. }
+ \item{method}{Can either "pearson", "kendall", or "spearman" which indicates
+ the correlation coefficient is to be used. Defaults to "spearman." See
+ \code{\link{cor}}. }
\item{n}{ \code{NULL} or an integral value giving the filter length
for the \code{\link{hanning}} filter used for removal of low
frequency variation. }
Modified: pkg/dplR/vignettes/intro-dplR.Rnw
===================================================================
--- pkg/dplR/vignettes/intro-dplR.Rnw 2014-04-12 17:36:14 UTC (rev 813)
+++ pkg/dplR/vignettes/intro-dplR.Rnw 2014-04-13 04:11:30 UTC (rev 814)
@@ -29,34 +29,52 @@
\maketitle
\begin{abstract}
-This document describes basic features of dplR including reading and working
-with ring-width data, detrending and standardization of ring-width data,
-building chronologies, and calculating descriptive statistics. A few
-simple plots are also presented.
+This document describes basic features of dplR by following the inital steps
+that an analyst might follow when working with a new tree-ring data set. The
+vignette starts with reading in ring widths and plotting them. We describe
+a few of the available methods for detrending and then show how to extract
+basic descriptive statistics. We show how to build a and plot simple
+mean-value chronology. We also show how to build a chronology
+using the expressed population signal from the detrended ring widths as an
+example of how more complicated analysis can be done using dplR.
\end{abstract}
\tableofcontents
\newpage
\section{Introduction}
+\subsection{What's 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
used functions in dplR. There is more detailed information
-available in the help files and in the literature including \cite{Bunn2008}
-and \cite{Bunn2010}.
+available in the help files and in the literature \citep{Bunn2008, Bunn2010}.
In this vignette, we will walk through the most basic activities of working
with tree-ring data in roughly the order that a user might follow. E.g.,
reading data, detrending, chronology building, and doing preliminary
exploratory data analysis via descriptive statistics.
+\subsection{Citing dplR and R}
+The creation of dplR is an act of love. We enjoy writing this software and
+helping users. However, neither of us is among the idle rich. Alas. We have
+jobs and occassionally have to answer to our beters. There is a nifty
+\code{citation} function in R that gives you information on how to best
+cite R and, in many cases, its packages. We ask that you please cite dplR
+and R appropraitely in your work. This way when our department chairs and
+deans accuse us of being dilletantes we can point to the use of dplR as a
+partial excuse.
+<<>>=
+citation()
+citation("dplR")
+@
+
\section{Working with Ring-Width Data}
\subsection{Reading Data}
There are, alas, many different ways that tree-ring data are digitally stored.
-These range in sophistication from the simple commonly used
-\href{http://www.ncdc.noaa.gov/paleo/treeinfo.html#formats}{Tucson} (decadal)
-format file of ring widths to the more complex
+These range in sophistication from the simple (and commonly used)
+\href{http://www.ncdc.noaa.gov/paleo/treeinfo.html#formats}{Tucson/decadal}
+format file of ring widths to the more complex (but richer)
\href{http://www.tridas.org/}{TRiDaS format}. We generally refer to these as
\code{rwl} objects for ``ring width length'' but there is no reason these can't be
other types of tree-ring data (e.g., density).
@@ -69,7 +87,7 @@
Throughout this vignette we will use the onboard data set \code{ca533}
which gives the raw ring widths for bristlecone pine \emph{Pinus longaeva} at
-Campito Mountain in California, USA. There are 34 series spanning over 1358
+Campito Mountain in California, USA. There are 34 series spanning 1358
years.
These objects are structured very simply as a \code{data.frame} with the series in
@@ -99,7 +117,7 @@
\begin{figure}[htbp]
\centering
\includegraphics{intro-dplR-a}
-\caption{A spaghetti plot of Campito Mountain ring widths.}
+\caption{A spaghetti plot of the Campito Mountain ring widths.}
\label{fig:rwl.plot}
\end{figure}
@@ -118,7 +136,7 @@
A rwi object has the same basic properties as the \code{rwl} object from which it is
made. I.e., it has the same number of rows and columns, the same names, and so
on. The difference is that each series has been standardized by dividing the
-ringwidths against a growth model (e.g., a stiff spline, a negative
+ring widths against a growth model (e.g., a stiff spline, a negative
exponential, etc.). This gives each series a mean of one (thus referred to
as ``indexed'') and allows a chronology to be built (next section). As
\code{read.rwl} is the primary function for getting data into R,
@@ -129,14 +147,10 @@
\subsection{Common Detrending Methods}
As any dendrochronologist will tell you, detrending is a dark art. In dplR we
have implemented some of the standard tools for detrending but all have
-drawbacks. In all of the detrend methods, the detrending is the estimation and
-removal of the tree's natural biological growth trend. The standardization is
-done by dividing each series by the growth trend to produce units in the
-dimensionless ring-width index (RWI).
-
-We'll discuss detrending via fitting a nonlinear function using
-\code{nls} (the \code{"ModNegExp"} method of \code{detrend}) and detrending
-via cubic smoothing spline (the \code{"Spline"} method of \code{detrend}). Much of the
+drawbacks. In all of the methods, the detrending is the estimation and
+removal of the low frequency variability that is due to biological or stand
+effects. The standardization is done by dividing each series by the growth
+trend to produce units in the dimensionless ring-width index (RWI). Much of the
text that follows is modified from the help page of \code{detrend}.
Probably the most common method for detrending is what is often
@@ -166,26 +180,28 @@
colMeans(ca533.rwi, na.rm=TRUE)
@
-An alternative method in \code{detrend} is to standardize with the \code{"Spline"} approach.
-This method uses a spline as the growth model where the frequency response
-is 0.50 at a wavelength of \(0.67 \times \text{series length}\) (unless specified differently by
-the user). This attempts to remove the low frequency
-variability that is due to biological or stand effects. Rather than detrend the
-entire \code{ca533} \code{rwl} object, we'll illustrate the spline method by detrending a
-single series using the \code{detrend.series} function, which produces a plot by
-default. See Figure~\ref{fig:spline.detrend}.
+When \code{detrend} is run on an \code{rwl} object the function loops through
+each series. It does this by calling a different function
+(\code{detrend.series}) for each column in the \code{rwl} object.
+But, a user can also call \code{detrend.series} and it's useful to do so here
+for educational purposes.
+Let's detrend a single series and apply more than one detrending method when we
+call it. We'll also call \code{detrend.series} using the verbose mode so that
+we can see the parameters applied for each method. The \code{detrend.series}
+function, produces a plot by default. See Figure~\ref{fig:detrend.series}.
+
<<b, fig=TRUE>>=
series <- ca533[, "CAM011"] # extract the series
names(series) <- rownames(ca533) # give it years as rownames
series.rwi <- detrend.series(y = series, y.name = "CAM011",
- method="Spline")
+ verbose=TRUE)
@
\begin{figure}[htbp]
\centering
\includegraphics{intro-dplR-b}
-\caption{Detrending a series via a spline.}
-\label{fig:spline.detrend}
+\caption{Detrending a single series via mutiple methods.}
+\label{fig:detrend.series}
\end{figure}
Often, a user will want to interactively detrend each series and fit a negative
@@ -222,7 +238,7 @@
correlations) as well as the expressed population signal and signal-to-noise
ratio for a data set. These are done in dplR using the \code{rwi.stats}
function so-named because these statistics are typically (but not always)
-carried out on detrended and standardized ring-width indices. If a data set
+carried out on detrended and standardized ring widths (rwi). If a data set
has more than one core taken per tree this information can be used in the
calculations to calculate within vs. between tree correlation. The function
\code{read.ids} is used to identify which trees have multiple cores.
@@ -270,11 +286,12 @@
dim(ca533.crn)
@
-The chronology can be plotted using the \code{crn.plot} function which
-has many arguments for customization. Here we'll just make a simple plot of the
-chronology with a smoothing spline added. See Figure~\ref{fig:crn.plot.spline}.
+An object produced by \code{chron} has a generic S3 moethod for plotting
+which calls the \code{crn.plot} function (which has many arguments for
+customization). Here we'll just make a simple plot of the chronology with
+a smoothing spline added. See Figure~\ref{fig:crn.plot.spline}.
<<c, fig=TRUE>>=
-crn.plot(ca533.crn, add.spline=TRUE, nyrs=20)
+plot(ca533.crn, add.spline=TRUE, nyrs=20)
@
\begin{figure}
\centering
@@ -283,6 +300,7 @@
\label{fig:crn.plot.spline}
\end{figure}
+\section{Prospectus}
In general this vignette aims to give a very cursory overview of basic tasks
that most dendrochronologists will want to be aware of. Know that we are just
scratching the surface of what dplR is capable of. As a small example,
@@ -341,7 +359,6 @@
\label{fig:crn.plot.eps}
\end{figure}
-\section{Prospectus}
We hope that this vignette helps users cover introductory data handling and
processing using dplR and R. As we noted above we are just providing a short
introduction as to what is possible in dplR. There are many other functions in
Added: pkg/dplR/vignettes/xdate-dplR.Rnw
===================================================================
--- pkg/dplR/vignettes/xdate-dplR.Rnw (rev 0)
+++ pkg/dplR/vignettes/xdate-dplR.Rnw 2014-04-13 04:11:30 UTC (rev 814)
@@ -0,0 +1,172 @@
+% -*- mode: noweb; noweb-default-code-mode: R-mode; -*-
+%\VignetteIndexEntry{Crossdating in dplR}
+\documentclass[a4paper,11pt]{article}
+\usepackage{dplR} % dplR settings - needs some work
+\usepackage[utf8]{inputenx} % R CMD build wants this here, not in dplR.sty
+\input{ix-utf8enc.dfu} % more characters supported
+\title{Crossdating in dplR}
+\author{Andy Bunn \and Mikko Korpela}
+\hypersetup{
+ pdfauthor = {Andy Bunn; Mikko Korpela},
+}
+\date{\footnotesize{$ $Processed with dplR
+\Sexpr{packageDescription("dplR", field="Version")}
+in \Sexpr{R.version.string} on \today}}
+
+\begin{document}
+\bibliographystyle{jss}
+
+\setkeys{Gin}{width=1.0\textwidth} % figure width
+\SweaveOpts{concordance=TRUE}
+\SweaveOpts{strip.white=true}
+\SweaveOpts{include=FALSE}
+<<echo=false>>=
+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
+@
+
+\maketitle
+
+\begin{abstract}
+Foo
+\end{abstract}
+\tableofcontents
+
+\newpage
+
+\section{Introduction}
+\subsection{What's Covered}
+The Dendrochronology Program Library in R (dplR) is a package for
+dendrochronologists to handle data processing and analysis. This
+document gives an introduction of some of the crossdating functions in
+dplR. This vignette is essentially a rehashing of \cite{Bunn2010}. Please
+cite that paper if you use dplR for crossdating. There is more detailed
+information on all these functions in the help files.
+
+\subsection{Citing dplR and R}
+The creation of dplR is an act of love. We enjoy writing this software and
+helping users. However, neither of us is among the idle rich. Alas. We have
+jobs and occassionally have to answer to our beters. There is a nifty
+\code{citation} function in R that gives you information on how to best
+cite R and, in many cases, its packages. We ask that you please cite dplR
+and R appropraitely in your work. This way when our department chairs and
+deans accuse us of being dilletantes we can point to the use of dplR as a
+partial excuse.
+
+<<>>=
+citation()
+citation("dplR")
+@
+
+\section{Ruining a Perfectly Good Data Set}
+
+Throughout this vignette we will use the onboard data set \code{co021}
+which gives the raw ring widths for Douglas fir \emph{Pseudotsuga menziesii}
+at Mesa Verde in Colorado, USA. There are 35 series spanning 788 years.
+
+We'll rename the \code{co021} object to \code{dat} because we are going to
+mess around with it and it seems like good practice to rename it. It is a
+beautifully sensitive series with long segment lengths, high standard
+deviation (relative to ring widths), large first-order autocorrelation,
+and a high mean interseries correlation ($\mathrm{r}\approx 0.84$). The data are
+plotted in Figure~\ref{fig:rwl.plot}.
+<<a, fig=TRUE>>=
+library(dplR)
+data(co021)
+dat <- co021
+dat.sum <- summary(dat)
+mean(dat.sum$year)
+mean(dat.sum$stdev)
+mean(dat.sum$median)
+mean(dat.sum$ar1)
+mean(interseries.cor(dat)[,1])
+plot(dat, plot.type="spag")
+@
+\begin{figure}[htbp]
+\centering
+\includegraphics{xdate-dplR-a}
+\caption{A spaghetti plot of the Mesa Verde ring widths.}
+\label{fig:rwl.plot}
+\end{figure}
+
+\textbf{By the way, if this is all new to you - you should stop reading this
+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.}
+
+To demonstrate how crossdating works in dplR, we will take this perfectly
+lovely data set and corrupt the dating of one of the series. By doing so we
+will be able to reenact one of the most common tasks of the dendrochronologist:
+tracking down a misdated core. Here we will take the "641143" core and remove
+one of the years of growth. This simulates a missing ring in the series. We'll
+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
+set.seed(4576)
+i <- sample(x=1:nrow(dat),size=1) # 709
+j <- sample(x=1:ncol(dat),size=1) # 12 "643114"
+tmp <- dat[,j]
+tmp <- c(NA,tmp[-i])
+dat[,j] <- tmp
+@
+We've now deleted the ith observation from the jth core while making sure that
+\code{dat} still has the appropriate numbers of rows. By sticking the NA at the
+start of the series it is as if we missed a ring while measuring.
+
+\section{Crossdating}
+\subsection{Assessing the Dating of a Data Set}
+The primary function for looking the crossdating of a tree-ring data set in
+dplR is \code{corr.rwl.seg}. This function looks at the correlation between
+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
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/dplr -r 814
More information about the Dplr-commits
mailing list