Skip to content
Chris Guiterman edited this page Jul 21, 2017 · 20 revisions

burnr Cookbook

The official guide to cooking with burnr.

Table of Contents

Make a data.frame with annual rows and columns for each series

Let's say that you want to use an fhx object to create a data.frame that has one row for each year and columns for each series. This is fairly simple with the reshape2 package.

library('reshape2')

data(lgr2)

d <- dcast(lgr2, 'year ~ series', value.var = 'rec_type')

This loads the reshape2 package and then loads the lgr2 fhx object from burnr, which we can use as an example dataset. The last line "casts" the ring data in lgr2 using the year column as rows, the series column as columns and then fills the body of the data.frame with values from rec_type.

Note that when you're casting the data, you may be missing several years from your rows. You can see this in our this example by running head(d) - our first year skips from 1366 to 1411:

year LGR54 LGR44 LGR47 ...
1366 inner_year NA NA ...
1411 NA inner_year NA ...
1423 unknown_fi unknown_fs NA ...
1424 recorder_year recorder_year NA ...
1425 recorder_year recorder_year NA ...
... ... ... ... ...

This is because the series have no recorded observations between 1366 and 1411. In FHX2, and FHX2 files, this is recorded as "null years" but we cut these null years out to save space and computation time. This usually doesn't affect people's code, but it's still something to be aware of.

One dirty way to get around this is

library('plyr')

year_range <- seq(min(lgr2$year), max(lgr2$year))
filler <- data.frame(year = year_range,
                     series = rep('hackishSolution', length(year_range)),
                     rec_type = rep(NA, length(year_range)))

d <- dcast(rbind(lgr2, filler), 
           'year ~ series', value.var = 'rec_type', 
           subset = .(series != 'hackishSolution'))

which should then fill in that null space.

year LGR54 LGR44 LGR47 ...
1366 inner_year NA NA ...
1367 NA NA NA ...
1368 NA NA NA ...
1369 NA NA NA ...
1370 NA NA NA ...
... ... ... ... ...

Annotate a facet plot

This is a quick recipe describing how do add annotations to a facet plot. This example is especially helpful if you want different annotations in each facet.

First, we load burnr, ggplot2, and some example fire history data.

library(ggplot2)

library(burnr)


data(lgr2)

data(lgr2_meta)

Now we use that data to create a facetted plot. Here we are doing one facet for each species in the data.

p <- plot_demograph(lgr2, facet_group = lgr2_meta$SpeciesID, facet_id = lgr2_meta$TreeID)

This is nothing special. Now we need to make a new data frame that holds all of our annotation data.

annotation_df <- data.frame(year = rep(1400, 4),
                            facet_group = c('PIPO', 'PIST', 'POTR', 'PSME'), 
                            label = c(rep('BACON!', 3), 'SPAM!'))

This is the information we need to plot an annotation over the year 1400, for each of the species facets. The label column is a little bit of extra information we need to plot text within these facets. You can see here that we are plotting 'BACON!' in the PIPO, PIST, and POTR facets. PSME, is going to get the text 'SPAM!'.

Now we plug this information into a geom_text and add this layer to our existing plot.

p + geom_text(data = annotation_df, aes(label = label), y = Inf, vjust = 1.5)

We give the function our annotation data and point its label aesthetic to our data frame's label column. Note that ggplot already knows that "year" in our annotation_df is associated with the x-axis. y = Inf is here to tell ggplot2 that we want the text to plot in the upper edge of the plot. We use vjust as a slight vertical adjustment to keep the text from getting cut on the edge of the facet. Generally, with annotations we would use ggplot2's annotate(), but here we need to use the more primitive geom_text() for finer controls. I suspect this might change with future versions of ggplot2, so be sure to keep an eye out for that.

Add site composites to facet plots

A limitation of ggplot graphics is an inability to use composite_rug within a facet generated by plot_demograph(). In this instance, one can make composites and add them to demograph plots within the facets, or stack the composites in a separate facet. We can use the color options of plot_demograph() to help distinguish between trees and site-level composites. This section will demonstrate how.

First, read in three related sites from the Jemez Mountains, New Mexico. These are available on the International Multiproxy Paleofire Database (IMPD).

library(burnr)
library(ggplot2)

url <- "https://www1.ncdc.noaa.gov/pub/data/paleo/firehistory/firescar/northamerica/"
pmr <- read_fhx(paste0(url, 'uspmr001.fhx'))
pme <- read_fhx(paste0(url, 'uspme001.fhx'))
pmw <- read_fhx(paste0(url, 'uspmw001.fhx'))

To designate the different sites and data type (tree vs. composite in this case) we create metadata tables for each site from the fhx objects. If we had species information we color by species as well, and would add "composite" as a sort of species.

pmr.meta <- data.frame(series = series_names(pmr), site = 'PMR', type = 'Tree')
pme.meta <- data.frame(series = series_names(pme), site = 'PME', type = 'Tree')
pmw.meta <- data.frame(series = series_names(pmw), site = 'PMW', type = 'Tree')

Make the composite fire histories for each site.

pmr.comp <- composite(pmr, comp_name = 'PMR.comp')
pme.comp <- composite(pme, comp_name = 'PME.comp')
pmw.comp <- composite(pmw, comp_name = 'PMW.comp')

Make a dataframe for the composites, and combine it with the tree-level data.

comp.meta <- data.frame(series = c('PMR.comp', 'PME.comp', 'PMW.comp'), site = c('PMR', 'PME', 'PMW'), type = 'Composite')

Combine all of the fhx objects and the metadata tables. We add here the sort function to order the series in the fhx object. In this case, the preferred final look of the plot will be trees sorted by their inner-ring date. In the next step we'll flip the series order to place the composites on the bottom of each facet, so here we sort the series in the opposite manner.

all.fhx <- sort(pmr, sort_by = "first_year", decreasing = TRUE) + 
              sort(pme, sort_by = "first_year", decreasing = TRUE) + 
              sort(pmw, sort_by = "first_year", decreasing = TRUE) + 
              pmr.comp + pme.comp + pmw.comp

all.meta <- rbind(pmr.meta, pme.meta, pmw.meta, comp.meta)

To ensure the composite stays at the bottom of the facet, the factor order needs to be reversed.

all.fhx$series <- factor(all.fhx$series, levels = rev(levels(all.fhx$series)))

Make the fire-demography plot.

plot_demograph(all.fhx, facet_group = all.meta$site, facet_id = all.meta$series,
                        color_group = all.meta$type, color_id = all.meta$series,
                        ylabels = FALSE, event_size = c(2.5, 1, 1),
                        plot_legend = TRUE, yearlims = c(1600, 1995)) +
  scale_color_manual(values=c('black', 'red')) +
  theme(legend.position = 'top', 
        legend.direction="horizontal", legend.background=element_rect(fill='white'), 
        legend.box="horizontal")

To place the composites in their own facet, change the site ID to "composite". Make sure the order in which the composites are added to the full fhx object matches the order of the tree-level site facets.

comp.meta <- data.frame(series = c('PMR.comp', 'PME.comp', 'PMW.comp'), site = 'Composite', type = 'Composite')

all.meta <- rbind(pmr.meta, pme.meta, pmw.meta, comp.meta)

plot_demograph(all.fhx, facet_group = all.meta$site, facet_id = all.meta$series,
               ylabels = FALSE, event_size = c(2.5, 1, 1),
               plot_legend = TRUE, yearlims = c(1600, 1995)) +
  theme(legend.position = 'top', 
        legend.direction="horizontal", legend.background=element_rect(fill='white'), 
        legend.box="horizontal")

Simple Recipe Template

_Your name ()_

Short description of what your recipe will do. Keep the language simple and approachable. Don't be pedantic.

All recipes need two parts, regardless of whether you use full section-headings or not: setup and the guts.

Setup

Describe all of the libraries you need to import, data you need to load, and functions you need to define as setup for the actual example in the recipe. Keep it short and digestible.

In general, all recipe code should be explicit rather than implicit. People should be able to reproduce your result based on what you've written in the recipe. Also, be sure your code is readable. This is paramount.

Present code blocks like this:

library('burnr')
    
thisIsAFunction()
look <- 'a_variable'

You can show inline code print('like this').

The Guts (or whatever you want to name this section)

Now that everything has been setup. This is where you present your cook trick or method. Show your code and explain what you're doing and why you're doing it. Assume that the reader has a moderate understanding of R.

Remember that readers should be able to reproduce your results!