Pages (Desktop)

Pages (Mobile)

Pollen diagrams in R using rioja - Part 3

In the third part of this guide series, which looks at plotting pollen diagrams using the "rioja" package, I will show you how you can combine different plot styles into a single pollen diagram figure. Now, your pollen diagram (or other stratigraphic diagram) could consist of bar plots, line plots and/or silhouette plots, rather than just a single plot type...

! Never used R before? Check out my getting started with R guide first.

Guide Information

Title Pollen diagrams in R using rioja - Part 3
Author Benjamin Bell
Published February 16, 2018
Last updated
R version 3.4.2
Packages base; rioja; vegan

This is a 3 part guide:

Part 1: Pollen diagrams using rioja How to plot pollen diagrams using rioja in R.
Part 2: Pollen diagrams using rioja "Hacks" to further customise the plot.
Part 3: Pollen diagrams using rioja How to use different plot types on the same plot.

Pollen diagrams using rioja - Part 3

This is part 3 of a guide series looking at the plotting pollen diagrams using the rioja package in R. If you missed part 1 or part 2, you can find them here. Ideally, you need to have read part 1 to follow this guide, or you need to be aware of how to use rioja.

In part 1, I showed you how to create different style pollen diagrams - where the entire pollen diagram would be a single style (e.g. bar plot, line plot, silhouette plot etc.). But, what if you wanted to create a pollen diagram that used multiple styles on a single figure?

The default options for the strat.plot() function in the rioja package only allow you to define the plot style for the entire diagram using logical (TRUE or FALSE) arguments. (Correct as at 15-Feb-2018, and rioja v0.9-15.1).

So, to create a bar plot pollen diagram you would use plot.line=FALSE, plot.poly=FALSE, plot.bar=TRUE, or to create a silhouette pollen diagram you would use plot.line=FALSE, plot.poly=TRUE, plot.bar=FALSE etc.

You might also remember that the option to add exaggerated curves allows you to create a vector of logical arguments to control which pollen types (or species) have an exaggerated curve. If you tried to do this with the plot type arguments, for example plot.bar=c(TRUE, TRUE, TRUE, FALSE, FALSE) etc. it would result in the following warning message, and would only plot the diagram using the first argument:

There were 28 warnings (use warnings() to see them)
> warnings()
1: In if (plot.bar) { ... :
  the condition has length > 1 and only the first element will be used

To use different plot styles on the same pollen diagram, there is a workaround! This involves using the add=TRUE argument to plot another diagram on top of the previous one, and defining different colour schemes for each diagram.

© Benjamin Bell. All Rights Reserved. http://www.benjaminbell.co.uk

Multiple plot styles on a single pollen diagram

For this example, we'll use the same data from part 1 of this guide, which can be downloaded from Mendeley Data - have a look at part 1 for detailed instructions on downloading the data (if you don't already have it).

Once you have the data, we'll import it into R, and then subset the data to remove unnecessary columns, although this time we'll keep the total arboreal (tree) pollen column and move it next to the other tree pollen columns.

# Import data
ma.pollen.raw <- read.csv("Morocco Pollen Surface Sample Data.csv", header=TRUE, row.names=1, sep=",", check.names=FALSE)
### Subset imported data
# Remove the first four columns
ma.pollen <- ma.pollen.raw[-1:-4]
# Remove total NAP and pollen concentration columns (but leave AP column)
ma.pollen <- ma.pollen[-50:-51]
# Remove columns where total abundanace is less than 10%
ma.sum <- colSums(ma.pollen) # Calculate column sums
ma.pollen1 <- ma.pollen[, which(ma.sum > 10)] # Subset data
# Re-order data frame to move AP column next to tree pollen columns
ma.pollen1 <- ma.pollen1[c(1:7, 28, 8:27)]

For an explanation of this code, please see part 1.


For the pollen diagram, we will use three plot styles: bar plots for the first 7 columns (the tree pollen), a line plot for the 8th column (total tree pollen), and silhouette plots for the remaining 20 columns (the non-tree pollen).

For the workaround to work, we'll plot this pollen diagram 3 times, and define the colour scheme for each plot type to "hide" the columns we do not want to plot.

So, we first need to create separate colour schemes for each plot type. You'll need to create a vector which defines a colour value for every column (pollen type or species) in the pollen diagram. In order to "hide" a particular plot style from appearing on a particular column, you should replace the colour value with NA. Use the following code to define different colour schemes for the different plot types:

# Bar plot colour scheme
c.bar <- c(rep("forestgreen", times=7), rep(NA, times=21))
# Line plot colour scheme
c.line <- c(rep(NA, times=7), "green2", rep(NA, times=20))
# Silhouette plot colour scheme
c.poly <- c(rep(NA, times=8), rep("gold2", times=20))
c.polyline <- c(rep(NA, times=8), rep("black", times=20))

Here, we have created several vectors which contain a string of colour values. The rep() command saves you from having to type the same thing multiple times! Type "c.bar" into the R console to see:

> c.bar
 [1] "forestgreen" "forestgreen" "forestgreen" "forestgreen" "forestgreen"
 [6] "forestgreen" "forestgreen" NA            NA            NA           
[11] NA            NA            NA            NA            NA           
[16] NA            NA            NA            NA            NA           
[21] NA            NA            NA            NA            NA           
[26] NA            NA            NA    

The important thing to remember, is that the position in the colour vector relates to the order of the columns in the pollen diagram. So, column 1 of the pollen diagram would use the first colour value in the vector, column 2 would use the second colour value and so on. No colour will be used where there is a NA value, which effectively hides the column on the diagram.

To plot the "first" pollen diagram, use the following code:

# Bar plot
pol.plot <- strat.plot(ma.pollen1, yvar=y.scale, y.tks=y.scale, y.rev=TRUE, plot.line=FALSE, plot.poly=FALSE, plot.bar=TRUE, col.bar=c.bar, lwd.bar=10, scale.percent=TRUE, xSpace=0.01, x.pc.lab=TRUE, x.pc.omit0=TRUE, las=2, srt.xlabel=45, y.axis=FALSE)

For each pollen diagram to line up correctly, you need to specify the same y axis and x axis arguments for each plot. However, you also need to use y.axis=FALSE to stop R from repeatedly drawing the y axis. We'll add this manually at the end.

This will give us the "first" part of the pollen diagram:

Next, we are going to plot the "second" pollen diagram:

# Line plot
pol.plot <- strat.plot(ma.pollen1, yvar=y.scale, y.tks=y.scale, y.rev=TRUE, plot.line=TRUE, plot.poly=FALSE, plot.bar=FALSE, col.line=c.line, lwd.line=3, scale.percent=TRUE, xSpace=0.01, x.pc.lab=TRUE, x.pc.omit0=TRUE, srt.xlabel=45, las=2, x.axis=FALSE, x.names=NA, y.axis=FALSE, add=TRUE)

In this code, we have changed the colour scheme to "c.line". We have also used additional arguments to stop R from re-drawing the x axes: x.axis=FALSE stops the values along the bottom of the diagram from being drawn, and x.names=NA stops the pollen names at the top from being drawn.

The add=TRUE argument means that R will draw this pollen diagram "on top" of the first diagram. So the pollen diagram should now look like this:

Next, we'll add the "third" pollen diagram:

# Silhouette plot with black lines
pol.plot <- strat.plot(ma.pollen1, yvar=y.scale, y.tks=y.scale, y.rev=TRUE, plot.line=FALSE, plot.poly=TRUE, plot.bar=TRUE, col.poly=c.poly, col.poly.line=c.polyline, col.bar=c.polyline, scale.percent=TRUE, xSpace=0.01, x.pc.lab=TRUE, x.pc.omit0=TRUE, srt.xlabel=45, las=2, x.axis=FALSE, x.names=NA, y.axis=FALSE, add=TRUE)

Finally, we'll add the y axis to the pollen diagram:

# Add y axis 
par(fig=c(0.01, 1, 0.07, 0.8))
axis(2, at=y.scale, labels=row.names(ma.pollen1), las=2)

par(fig) is used to change the position of the plot on the figure, using relative values between 0 and 1, for xLeft, xRight, yBottom, yTop. In this case, we need to change the position and size of the y axis to match the position of the pollen diagram.

The final "complete" pollen diagram should look like this:

Thanks for reading this guide! Please leave any comments or questions below.

This is a 3 part guide:

Part 1: Pollen diagrams using rioja How to plot pollen diagrams using rioja in R.
Part 2: Pollen diagrams using rioja "Hacks" to further customise the plot.
Part 3: Pollen diagrams using rioja How to use different plot types on the same plot.
© Benjamin Bell. All Rights Reserved. http://www.benjaminbell.co.uk

Further reading

Pollen diagrams using rioja - Part 1 - In the first part of this guide, I show you how to plot pollen diagrams using rioja.

Pollen diagrams using rioja - Part 2 - In the second part of this guide, I show you some additional options available for plotting pollen diagrams.

25 comments

  1. Is there any way to make a stacked line graph plot? For example the classic stacked plot showing trees, shrubs and herbs all on the same section at the end of pollen diagrams?

    ReplyDelete
    Replies
    1. Hi,

      As far as i know, there is no option to do this using rioja - unless you were to use the "fun" argument to add your own custom curves to a plot, although this is not something i've investigated yet..

      My suggestion would be to make the stacked bar graph using the base plot() function in R, and use an image editor to add it to your pollen diagram.

      Best wishes,
      Ben

      Delete
  2. Hi Ben,

    Is there any way to graph diatoms density plots and PCA plot? I have tried to do it following your procedure, but I do not know how to change the x-axis scale.
    Thank you.

    ReplyDelete
    Replies
    1. Hi,

      Apologies for delay in the reply. I unfortunately do not know enough about diatom density plots to be able to advise whether this is possible using rioja. And, as far as i am aware, you cannot do PCA with rioja either - although it is of course possible to use R for this, and i have a guide to doing PCA on this blog.

      The x-axis options are rather limited, either plotting by percentage or absolute values - you can adjust the labels, but rioja tends to control the scaling.

      Best wishes,
      Ben

      Delete
  3. Hey there, is it possible to include both negative and positive values on the plot?

    ReplyDelete
    Replies
    1. Hi,

      I am not sure actually - give it a try and see if it works :)

      Best wishes,
      Ben

      Delete
  4. Is there an option to plot more than 200 taxa in the same plot in a nice way?

    ReplyDelete
    Replies
    1. Hi,

      I am not sure as I have never plotted so many. You might need to split your taxa up, and plot several separate diagrams, and then stitch them together using a graphics editor.

      Best wishes,
      Ben

      Delete
  5. Excellent, Thansk so much, now I know another way to plot interactions diagrams between plants and bees

    ReplyDelete
  6. Hi Ben,

    Do you know how one would code to have the y-axis as 'Age' and scale the axis relative to an age-depth model, such as one outputted by rBacon? This can be done in Tilia, so hoping it can be done in r, too.

    Cheers,
    Jonny

    ReplyDelete
    Replies
    1. Hi Jonny,

      Just use the age values for each sample as the y-axis argument i.e. replace depths with ages, and it should plot fine. Best using a line graph for this.

      Best wishes

      Delete
  7. Hi Ben,
    thank you very much for this tutorial. I have two quiestions.
    1) The pollen is merge to bigger groups like Trees, Woody Schrubs and Herb with clamps in your article. It is there some code in R to do it, or you need to use some othre graphical editor.
    2) How can I cite this tutorial in my article?
    Best regards and stay heathy
    Lukas

    ReplyDelete
    Replies
    1. Hi Lukas,

      1) I assume by "clamps" you are referring to the black lines used to group the taxa above their names. Yes, this was done in a graphics editor, I do not believe rioja can draw these for you.

      2) No need to cite the tutorial, but it is good practice to cite the package used. A good way to get the information you need is use the citation() command in the R terminal. See below:

      > library(rioja)
      This is rioja 0.9-26
      > citation("rioja")

      To cite package 'rioja' in publications use:

      Juggins, S. (2020) rioja: Analysis of Quaternary Science Data, R package version
      (0.9-26). (https://cran.r-project.org/package=rioja).

      A BibTeX entry for LaTeX users is

      @Manual{,
      title = {rioja: Analysis of Quaternary Science Data},
      author = {Steve Juggins},
      year = {2020},
      note = {R package version 0.9-26},
      url = {https://cran.r-project.org/package=rioja},
      }


      Best wishes,
      Ben

      Delete
  8. Is it possible to put the lithology and the names of the zones? Also add carbon dates

    ReplyDelete
    Replies
    1. Hi Jorge, apologies for my slow response - i think this is something you would have to do manually using a graphics editor, as i dont think you can add zone names using rioja.

      Delete
  9. Hi Ben!

    Fantastic series of guides, really helped me understand how to plot with this package.

    I have some data that occasionally dips into negative values. I cannot figure out how to get these display.

    Any ideas how I might achieve this, or is it simply not possible with the strat.plots?

    Like your ABBOREAL data, it is only a single diagram that would need to show negatives.

    Thank you!

    ReplyDelete
    Replies
    1. Hi Virgel,

      You can plot negative values if you change the scale.percent argument to FALSE (scale.percent=FALSE), and set the plot style to a line plot (plot.line=TRUE, plot.poly=FALSE, plot.bar=FALSE).

      Although, not sure i fully understand how you have negative values :D

      Best wishes,
      Ben

      Delete
    2. Thank you Ben!

      That certainly allows the negative values to display. I need only one of the diagrams to have negative numbers.
      What I understand from your guides makes me suspect I cannot do both.

      My own plots would have a line plot like your ABBOREAL data, but this could dip into negative numbers at times.

      I have my count data as you do, but display some water-table data along with it, which occasionally dips into negatives to indicate the depth is above 0cm, that is to say the surface.

      Maybe you know otherwise, but I could not seem to find a way to do both in one single plot.

      Best,
      V.

      Delete
  10. Hi Ben,
    Would this multi-plot method work for plotting different types of data from different sample depths? For example, I have biogenic silica data and macrofossil data from the same core that were analyzed at different sample depths/ages that I would like to plot all on one stratigraphic diagram. Also wondering if I could add data from a different core and chronology from the same lake all on the same age y axis? Thank you! I am excited to try this. I am an old Tilia user and always ran into some limitations.
    Laurie

    ReplyDelete
  11. Hi! i would like to know how you did the groupings as shown in the PDF version of your data. How did you group them as trees, shrubs and woody shrubs? is there a code for it? if there is do let me know.

    ReplyDelete
    Replies
    1. Hi Gaya,

      I think you are referring to the diagram in the published paper - and i think you are talking about the group names that appear above the taxa names with lines which group them?

      If so - that was done manually in a graphics editor - i don't think its possible with rioja.

      Best wishes,
      Ben

      Delete
    2. Thank you for your quick response.

      Delete
  12. Hi! Ben, do you have a r code for Canonical correspondence analysis (CCA) or Correspondence analysis ?

    ReplyDelete
    Replies
    1. Hi Gaya,

      Take a look at the vegan package which is able to do CCA and CA: https://cran.r-project.org/web/packages/vegan/index.html

      Unfortunately, i do not have a tutorial for this yet.

      Best wishes,
      Ben

      Delete

Comments are moderated. There may be a delay until your comment appears.