[Dplr-commits] r830 - in pkg/dplR: . R man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Mon Apr 21 22:20:48 CEST 2014


Author: andybunn
Date: 2014-04-21 22:20:48 +0200 (Mon, 21 Apr 2014)
New Revision: 830

Added:
   pkg/dplR/R/xskel.plot.R
   pkg/dplR/man/xskel.plot.Rd
Modified:
   pkg/dplR/ChangeLog
   pkg/dplR/NAMESPACE
   pkg/dplR/R/ccf.series.rwl.R
   pkg/dplR/R/xskel.ccf.plot.R
   pkg/dplR/man/xskel.ccf.plot.Rd
Log:
* made an xskel.plot function to go with xskel.ccf.plot
* fussed with the colors a bit in ccf.series.rwl

Modified: pkg/dplR/ChangeLog
===================================================================
--- pkg/dplR/ChangeLog	2014-04-21 19:34:45 UTC (rev 829)
+++ pkg/dplR/ChangeLog	2014-04-21 20:20:48 UTC (rev 830)
@@ -13,10 +13,10 @@
 - Added summary.rwl as an S3Method.
 - Added insert and delete.ring functions.
 
-File: skel.ccf.plot.R
+File: xskel.ccf.plot.R and xskel.plot.R and 
 ---------------
 
-- New and amitious plotting function to help cross date with skeleton plot
+- New plotting functions to help crossdate with skeleton plot
   and cross correlation plots.
 
 File: ccf.series.rwl.R
@@ -25,6 +25,7 @@
 - Switched the order of x and y in the call to ccf(). This makes a great deal
   more logical sense now as a missing ring shows up with a positive lag rather
   than a negative lag.
+- Changed color scheme a bit to look less harsh
 
 Files: ccf.series.rwl.R, corr.series.seg.R, series.rwl.plot.R
 -------------------------------------------------------------

Modified: pkg/dplR/NAMESPACE
===================================================================
--- pkg/dplR/NAMESPACE	2014-04-21 19:34:45 UTC (rev 829)
+++ pkg/dplR/NAMESPACE	2014-04-21 20:20:48 UTC (rev 830)
@@ -38,7 +38,7 @@
        tridas.vocabulary, uuid.gen, wavelet.plot, wc.to.po,
        write.compact, write.crn, write.rwl, write.tridas,
        write.tucson, plot.rwl, interseries.cor, summary.rwl,
-       plot.crn, insert.ring, delete.ring, xskel.ccf.plot)
+       plot.crn, insert.ring, delete.ring, xskel.ccf.plot, xskel.plot)
 
 S3method(print, redfit)
 S3method(plot, rwl)

Modified: pkg/dplR/R/ccf.series.rwl.R
===================================================================
--- pkg/dplR/R/ccf.series.rwl.R	2014-04-21 19:34:45 UTC (rev 829)
+++ pkg/dplR/R/ccf.series.rwl.R	2014-04-21 20:20:48 UTC (rev 830)
@@ -117,10 +117,16 @@
                        panel.abline(v=lag.vec, lty="solid", col="gray")
                        panel.abline(h=0, v=0, lwd=2)
                        panel.abline(h=sig, lwd=2, lty="dashed")
-                       col <- ifelse(y > 0, "#E41A1C", "#377EB8")
+                       #col <- ifelse(y > 0, "#E41A1C", "#377EB8")
+                       col <- ifelse(y > 0, "darkred", "darkblue")
+                       bg <- ifelse(y > 0, "lightsalmon", "lightblue")
                        ## segments, dots for all r
-                       panel.segments(x1=x, y1=0, x2=x, y2=y, col=col, lwd=2)
-                       panel.dotplot(x, y, col = col, ...)
+                       #panel.segments(x1=x, y1=0, x2=x, y2=y, col=col, lwd=2)
+                       #panel.dotplot(x, y, col = col, ...)
+                       panel.segments(x1=x, y1=0, x2=x, y2=y, 
+                                      col=col, lwd=2)
+                       panel.dotplot(x, y, col = col, fill=bg,
+                                     pch=21,...)
                    }, ...)
         trellis.par.set(strip.background = list(col = "transparent"),
                         warn = FALSE)

Modified: pkg/dplR/R/xskel.ccf.plot.R
===================================================================
--- pkg/dplR/R/xskel.ccf.plot.R	2014-04-21 19:34:45 UTC (rev 829)
+++ pkg/dplR/R/xskel.ccf.plot.R	2014-04-21 20:20:48 UTC (rev 830)
@@ -192,16 +192,16 @@
   
   # text  
   grid.text(master.yrs.sig, x=unit(master.yrs.sig,"native"), 
-            y = unit(-12, "native"), rot = 90,
+            y = unit(0, "npc"), rot = 90,just="right",
             gp=gpar(fontsize=12))
   grid.text(series.yrs.sig, x=unit(series.yrs.sig,"native"), 
-            y = unit(12, "native"), rot = 90,
+            y = unit(1, "npc"), rot = 90,just="left",
             gp= gpar(fontsize = 12))
-  grid.text("Master",x=unit(min(yrs)-1,"native"),
-            y=unit(-10,"native"),just = "left",rot=90,
+  grid.text("Master",x=unit(0,"npc"),
+            y=unit(0,"npc"),hjust = 0,vjust = 0,rot=90,
             gp= gpar(fontsize = 12))
-  grid.text("Series",x=unit(min(yrs)-1,"native"),
-            y=unit(10,"native"),just = "right",rot=90,
+  grid.text("Series",x=unit(0,"npc"),
+            y=unit(1,"npc"),hjust=1,vjust=0,rot=90,
             gp= gpar(fontsize = 12))
   
   upViewport(3) # back to bnd

Added: pkg/dplR/R/xskel.plot.R
===================================================================
--- pkg/dplR/R/xskel.plot.R	                        (rev 0)
+++ pkg/dplR/R/xskel.plot.R	2014-04-21 20:20:48 UTC (rev 830)
@@ -0,0 +1,152 @@
+xskel.plot <- function(rwl,series,series.yrs = as.numeric(names(series)),
+         win.start, win.end=win.start+100, n = NULL, prewhiten = TRUE, 
+         biweight = TRUE) {      
+  
+  ## Handle different types of 'series'
+  tmp <- pick.rwl.series(rwl, series, series.yrs)
+  rwl <- tmp[[1]]
+  series <- tmp[[2]]
+  
+  master.yrs <- as.numeric(rownames(rwl))
+  series.yrs <- as.numeric(names(series))
+  yrs <- seq(from=win.start,to=win.end)
+  nyrs <- length(yrs)
+  
+  if(nyrs > 101){
+    warning("These plots get crowded with windows longer than 100 years.")
+  }
+  # check window overlap with master and series yrs
+  if (!all(yrs %in% series.yrs)) {
+    cat("Window Years: ", min(yrs), "-", max(yrs)," & ",
+        "Series Years: ", min(series.yrs), "-", max(series.yrs),
+        "\n", sep="")
+    stop("Fix window overlap")
+  }
+  if (!all(yrs %in% master.yrs)) {
+    cat("Window Years: ", min(yrs), "-", max(yrs)," & ",
+        "Master Years: ", min(master.yrs), "-", max(master.yrs),
+        "\n", sep="")
+    stop("Fix window overlap")
+  }
+  
+  # normalize.
+  names(series) <- series.yrs
+  tmp <- normalize.xdate(rwl, series, n, prewhiten, biweight)
+  
+  # master
+  master <- tmp$master
+  master.yrs <- as.numeric(names(master))
+  master <- master[master.yrs%in%yrs]
+  master.yrs <- as.numeric(names(master))
+  # series
+  series <- tmp$series
+  series.yrs <- as.numeric(names(series))
+  series <- series[series.yrs%in%yrs]
+  series.yrs <- as.numeric(names(series))
+  
+  
+  # skeleton
+  master.skel <- cbind(master.yrs,xskel.calc(master))
+  master.skel <- master.skel[master.skel[,1]%in%yrs,]
+  master.yrs.sig <- master.skel[!is.na(master.skel[,2]),1]
+  series.skel <- cbind(series.yrs,xskel.calc(series))
+  series.skel <- series.skel[series.skel[,1]%in%yrs,]
+  series.yrs.sig <- series.skel[!is.na(series.skel[,2]),1]
+  
+  # cor and skel agreement
+  overall.r <- round(cor(series,master),3)    
+  overall.agree <- sum(series.yrs.sig%in%master.yrs.sig)/length(master.yrs.sig)
+  overall.agree <- round(overall.agree*100,1)
+  
+  # build page for plotting
+  grid.newpage()
+  # 1.0 a bounding box for margins  
+  bnd.vp <- plotViewport(margins=rep(0.5,4),name = "bnd.vp") # 1/2 line margin
+  # go from bottom up.
+  
+  # 4.1 bounding box for skeleton plot. 55% of area  
+  skel.bnd.vp <- viewport(x = 0, y = 0, width = 1, height = 0.95,
+                          just = c("left", "bottom"), name = "skel.bnd.vp")
+  # 4.2 plotting region for skeleton plot. 2 lines left and right. 
+  # 3 lines on top and bottom
+  skel.region.vp <- plotViewport(margins=c(3,2,3,2),
+                                 xscale=c(min(yrs)-0.5,max(yrs)+0.5), 
+                                 yscale=c(-10,10),
+                                 name = "skel.region.vp")
+  # 5.0 a box for overall text. 5%
+  overall.txt.vp <- viewport(x = 0, y = 0.95, width = 1, height = 0.05,
+                             just = c("left", "bottom"), 
+                             name = "overall.txt.vp")
+  
+  # actual plotting
+  pushViewport(bnd.vp) # inside margins
+  pushViewport(skel.bnd.vp) # inside skel
+  pushViewport(skel.region.vp) # inside margins
+  grid.rect(gp = gpar(col="lightgreen", lwd=1))
+  grid.grill(h = unit(seq(-10, 10, by=1), "native"),
+             v = unit(yrs-0.5, "native"),
+             gp = gpar(col="lightgreen", lineend = "square", 
+                       linejoin = "round"))
+  # rw plot
+  master.tmp <- master*-2
+  for(i in 1:length(yrs)){
+    xx <- c(yrs[i]+0.5,yrs[i]-0.5,yrs[i]-0.5,yrs[i]+0.5)
+    yy <- c(0,0,master.tmp[i],master.tmp[i])
+    grid.polygon(xx,yy,default.units="native", 
+                 gp=gpar(fill='lightgreen',col='darkgreen')) 
+  }
+  series.tmp <- series*2
+  for(i in 1:length(yrs)){
+    xx <- c(yrs[i]+0.5,yrs[i]-0.5,yrs[i]-0.5,yrs[i]+0.5)
+    yy <- c(0,0,series.tmp[i],series.tmp[i])
+    grid.polygon(xx,yy,default.units="native", 
+                 gp=gpar(fill='lightgreen',col='darkgreen')) 
+  }
+  
+  #master
+  grid.segments(x0=master.yrs.sig,y0=0,
+                x1=master.yrs.sig,y1=-10,
+                default.units="native",
+                gp=gpar(lwd=1,col='black',lineend="butt"))
+  grid.segments(x0=master.skel[,1],y0=0,
+                x1=master.skel[,1],y1=master.skel[,2]*-1,
+                default.units="native",
+                gp=gpar(lwd=5,col='black',lineend="butt"))
+  #series
+  grid.segments(x0=series.yrs.sig,y0=0,
+                x1=series.yrs.sig,y1=10,
+                default.units="native",
+                gp=gpar(lwd=1,col='black',lineend="butt"))
+  grid.segments(x0=series.skel[,1],y0=0,
+                x1=series.skel[,1],y1=series.skel[,2],
+                default.units="native",
+                gp=gpar(lwd=5,col='black',lineend="butt"))
+  
+  # text  
+  grid.text(master.yrs.sig, x=unit(master.yrs.sig,"native"), 
+            y = unit(0, "npc"), rot = 90,just="right",
+            gp=gpar(fontsize=12))
+  grid.text(series.yrs.sig, x=unit(series.yrs.sig,"native"), 
+            y = unit(1, "npc"), rot = 90,just="left",
+            gp= gpar(fontsize = 12))
+  grid.text("Master",x=unit(0,"npc"),
+            y=unit(0,"npc"),hjust = 0,vjust = 0,rot=90,
+            gp= gpar(fontsize = 12))
+  grid.text("Series",x=unit(0,"npc"),
+            y=unit(1,"npc"),hjust=1,vjust=0,rot=90,
+            gp= gpar(fontsize = 12))
+  
+  upViewport(3) # back to bnd
+  pushViewport(overall.txt.vp) # description
+  tmp.txt <- paste("Period: ",min(yrs),"-",max(yrs), 
+                   ", r(lag0)= ", overall.r, 
+                   ". Skeleton Agreement ", overall.agree, "%",sep="")
+  tmp.txt <- bquote("Period:" ~ .(min(yrs)) * "-" * 
+                      .(max(yrs)) * ","~r[lag0] * "=" * .(overall.r)*
+                      ","~"Skeleton Agreement"~.(overall.agree)*"%")
+  grid.rect(gp=gpar(col=NA,fill="white"))
+  grid.text(tmp.txt,y=unit(0.5,"npc"),x=unit(0.5,"npc"),
+            just = "center",
+            gp= gpar(fontsize = 12))
+  
+}
\ No newline at end of file

Modified: pkg/dplR/man/xskel.ccf.plot.Rd
===================================================================
--- pkg/dplR/man/xskel.ccf.plot.Rd	2014-04-21 19:34:45 UTC (rev 829)
+++ pkg/dplR/man/xskel.ccf.plot.Rd	2014-04-21 20:20:48 UTC (rev 830)
@@ -1,6 +1,6 @@
 \name{xskel.ccf.plot}
 \alias{xskel.ccf.plot}
-\title{ Skeleton Plot with Cross Correlation }
+\title{ Skeleton Plot for Series and Master with Cross Correlation }
 \description{
   ...
 }

Added: pkg/dplR/man/xskel.plot.Rd
===================================================================
--- pkg/dplR/man/xskel.plot.Rd	                        (rev 0)
+++ pkg/dplR/man/xskel.plot.Rd	2014-04-21 20:20:48 UTC (rev 830)
@@ -0,0 +1,80 @@
+\name{xskel.plot}
+\alias{xskel.plot}
+\title{ Skeleton Plot for Series and Master }
+\description{
+  ...
+}
+\usage{
+xskel.plot(rwl,series,series.yrs = as.numeric(names(series)),
+                       win.start, win.end=win.start+100, n = NULL, 
+                       prewhiten = TRUE, biweight = TRUE)
+
+}
+\arguments{
+  \item{rwl}{ a \code{data.frame} with series as columns and years as rows
+    such as that produced by \code{\link{read.rwl}}. }
+  \item{series}{ a \code{numeric} or \code{character} vector.  Usually a
+    tree-ring series.  If the length of the value is 1, the
+    corresponding column of \code{\var{rwl}} is selected (by name or
+    position) as the series and ignored when building the master
+    chronology.  Otherwise, the value must be \code{numeric}. }
+  \item{series.yrs}{ a \code{numeric} vector giving the years of
+    \code{\var{series}}.  Defaults to
+    \code{as.numeric(names(\var{series}))}. }
+  \item{win.start}{ year to start window }
+  \item{win.end}{ year to end window } 
+  \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. }
+  \item{prewhiten}{ \code{logical} flag. If \code{TRUE} each series is
+    whitened using \code{\link{ar}}. }
+  \item{biweight}{ \code{logical} flag. If \code{TRUE} then a robust
+    mean is calculated using \code{\link{tbrm}}. }
+}
+\details{
+This function produces a plot that is a mix of a skeleton plot and a 
+cross-correlation plot. It's used in crossdating.
+
+The top panel shows the normalized values for the master chronology
+(bottom half) and the series (top half) in green. The values are the 
+detrended and standardized data (e.g., RWI).
+
+Similarly, the black lines are a skeleton plot for the master and 
+series with the marker years annotated for the master on the bottom axis and 
+series 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 years of skeleton bars for the series and master. 
+I.e., if all the black lines occur in the same years the percentage would be
+100\%. 
+
+The bottom panels show cross correlations for the first half (left) and second 
+half of the time series using function \code{\link{ccf}} as
+\code{ccf(x=series,y=master,lag.max=5}. 
+
+The plot is built using the \code{\link[grid]{Grid}} package which allows for
+great flexibility in building complicated plots. However, these plots look best
+when they don't cover too wide a range of years (unless the plotting device
+is wider than is typical). For that reason the user 
+will get a warning if \code{win.width} is greater than 100 years.
+
+}
+\value{
+  None. Invoked for side effect (plot).
+}
+\author{ Andy Bunn. Patched and improved by Mikko Korpela. }
+\seealso{ \code{\link{ccf}}
+}
+\examples{data(co021)
+dat <- co021
+#corrupt a series
+bad.series <- dat$"641143"
+names(bad.series) <- rownames(dat)
+bad.series <- delete.ring(bad.series,year=1825)
+# good match
+xskel.plot(rwl=dat,series=bad.series,win.start=1850)
+# overlap missing ring
+xskel.plot(rwl=dat,series=bad.series,win.start=1800)
+}
+
+\keyword{ hplot }
+



More information about the Dplr-commits mailing list