introduction-to-r-part-3.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### R in a nutshell\n",
    "\n",
    "- Statistical programming environments\n",
    "- Originally designed and implemented by statisticians\n",
    "- Widely popular due to its extensive collection of community-contributed packages\n",
    "- Quickly gaining places among traditional proprietary tools such as SAS and STATA for data analytics\n",
    "\n",
    "#### Learning Objectives\n",
    "\n",
    "- Understand basic programming concepts: variables, assignment, functions, loops, conditions\n",
    "- Understand core R concepts: data loading, data types, data access, libraries\n",
    "- Understand advanced R concepts: data manipulation, visualization\n",
    "- Understand HPC concepts: running R codes on the Palmetto supercomputer via batch submission scripts\n",
    "\n",
    "#### Materials on this notebook is based on two lessons by Software Carpentry and Data Carpentry:\n",
    "\n",
    "- Introduction to Programming using R\n",
    "- Data Analysis and Visualization in R for Ecology"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Where am I?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "getwd()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "surveys_complete <- read.csv(\"output_data/surveys_complete.csv\")\n",
    "head(surveys_complete)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data visualization"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "** Installation: **"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "setupLibrary <- function(libraryName){\n",
    "  if (!require(libraryName, character.only = TRUE)){\n",
    "    install.packages(libraryName, dep = TRUE)\n",
    "    if (!require(libraryName, character.only = TRUE)){\n",
    "      print('Package not found')\n",
    "    }\n",
    "  } else {\n",
    "    print('Package is loaded')\n",
    "  }\n",
    "}\n",
    "\n",
    "setupLibrary('dplyr')\n",
    "setupLibrary('ggplot2')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*ggplot:* graphical presentations are described as a combination of elements and built by adding new elements"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "** Steps to build a ggplot **"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Bind the plot to a specific data frame using the `data` argument*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data=surveys_complete)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*define aesthetics (**`aes`**), by selecting the variables to be plotted and the variables to define presentation such as plotting size, shape, color, etc ...*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = surveys_complete, aes(x = weight, y = hindfoot_length))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*add **`geoms`** – graphical representation of the data in the plot (points, lines, bars). To add a geom to the plot use **`+`** operator:*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = surveys_complete, aes(x = weight, y = hindfoot_length)) +\n",
    "    geom_point()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "The **`+`** enables using plot templates to explore various plot designs:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "surveys_plot <- ggplot(data = surveys_complete, aes(x = weight, y = hindfoot_length))\n",
    "surveys_plot + geom_point()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Notes: **\n",
    "\n",
    "- Any configurations defined inside `ggplot()` are also visible by all `geom_` layers. This includes x and y axis set up in `aes()`. \n",
    "- Individual `aes()` for a given `geom_` can be set independely of the global `aes()` in `ggplot()`.\n",
    "- In the case of multi-line presentation, the `+` sign must be placed at the end of each line containing layer(s)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Challenge: **\n",
    "\n",
    "Scatter plots can be useful exploratory tools for small datasets. For data sets with large numbers of observations, such as the surveys_complete data set, overplotting of points can be a limitation of scatter plots. One strategy for handling such settings is to use hexagonal binning of observations. The plot space is tesselated into hexagons. Each hexagon is assigned a color based on the number of observations that fall within its boundaries. To use hexagonal binning with ggplot2, first install the R package hexbin from CRAN:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "install.packages(\"hexbin\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Then use `geom_hex()` function from the ggplot2 package to visualize data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- What are the relative strengths and weaknesses of a hexagonal bin plot compared to a scatter plot? \n",
    "- Examine the above scatter plot and compare it with the hexagonal bin plot that you created."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Building plots iteratively: **"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*Start with a template:*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "surveys_plot <- ggplot(data = surveys_complete, aes(x = weight, y = hindfoot_length))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*Choose a geom:*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "surveys_plot + geom_point()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*Customize geom:*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# add transparency\n",
    "surveys_plot + geom_point(alpha = 0.1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# add colors for all the points:\n",
    "surveys_plot + geom_point(alpha = 0.1, color = \"blue\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# color each species in the plot differently:\n",
    "surveys_plot + geom_point(alpha = 0.1, aes(color=species_id))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*Change to a different geom, `geom_boxplot()`*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "surveys_plot <- ggplot(data = surveys_complete, aes(x = species_id, y = hindfoot_length))\n",
    "surveys_plot + geom_boxplot()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# adding points to boxplot to understand number of measurements and distribution:\n",
    "surveys_plot + geom_boxplot(alpha = 0) +\n",
    "    geom_jitter(alpha = 0.3, color = \"tomato\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Challenge: **\n",
    "\n",
    "How can we show the box plot?     "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Challenge: **\n",
    "\n",
    "Boxplots are useful summaries, but hide the shape of the distribution. For example, if there is a bimodal distribution, it would not be observed with a boxplot. An alternative to the boxplot is the violin plot (sometimes known as a beanplot), where the shape (of the density of points) is drawn.\n",
    "\n",
    "- Replace the box plot with a violin plot; see `geom_violin()`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "In many types of data, it is important to consider the scale of the observations. For example, it may be worth changing the scale of the axis to better distribute the observations in the space of the plot. Changing the scale of the axes is done similarly to adding/modifying other components (i.e., by incrementally adding commands). Try making these modifications:\n",
    "\n",
    "- Represent weight on the log10 scale; see `scale_y_log10()`\n",
    "- Create boxplot for the scaled weights.\n",
    "- Add color to the datapoints on your boxplot according to the species from which the sample was taken (species_id)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Plotting time-series data: **"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Number of counts per year for each species:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "yearly_counts <- surveys_complete %>%\n",
    "                 group_by(year, species_id) %>%\n",
    "                 tally\n",
    "head(yearly_counts)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Plot everything!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = yearly_counts, aes(x = year, y = n)) +\n",
    "     geom_line()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*We can improve the clarity of the graphs by separating them into individual lines for individual species*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = yearly_counts, aes(x = year, y = n, group = species_id)) +\n",
    "    geom_line()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*We can make this even better!*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = yearly_counts, aes(x = year, y = n, group = species_id, colour = species_id)) +\n",
    "    geom_line()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Faceting: ** \n",
    "\n",
    "It is possible to create multiple plots within a larger plot frame based on a factor variable within the data set. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "str(yearly_counts)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = yearly_counts, aes(x = year, y = n, group = species_id, colour = species_id)) +\n",
    "    geom_line() +\n",
    "    facet_wrap(~ species_id)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Further customization on faceted plots can be done. For example, we can split data within each individual plots into lines presenting male and female:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "yearly_sex_counts <- surveys_complete %>%\n",
    "                      group_by(year, species_id, sex) %>%\n",
    "                      tally"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "head(yearly_sex_counts)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    " ggplot(data = yearly_sex_counts, aes(x = year, y = n, color = species_id, group = sex)) +\n",
    "     geom_line() +\n",
    "     facet_wrap(~ species_id)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "To improve this presentation, we can customize the theme layers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    " ggplot(data = yearly_sex_counts, aes(x = year, y = n, color = species_id, group = sex)) +\n",
    "     geom_line() +\n",
    "     facet_wrap(~ species_id) +\n",
    "     theme_bw() +\n",
    "     theme(panel.grid.major.x = element_blank(),\n",
    "       panel.grid.minor.x = element_blank(),\n",
    "       panel.grid.major.y = element_blank(),\n",
    "       panel.grid.minor.y = element_blank())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Since we already separate species into individual plots, we do not need to color the plot, but we will need to differentiate between the male and femal lines within each plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = yearly_sex_counts, aes(x = year, y = n, color = sex, group = sex)) +\n",
    "    geom_line() +\n",
    "    facet_wrap(~ species_id) +\n",
    "    theme_bw()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Challenge: **\n",
    "\n",
    "Create a plot that shows how the average weight of each species changes over the years. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- The `facet_wrap` geometry extracts plots into an arbitrary number of dimensions to allow them to cleanly fit on one page. \n",
    "- The `facet_grid` geometry allows you to explicitly specify how you want your plots to be arranged via formula notation (rows ~ columns; a . can be used as a placeholder that indicates only one row or column)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# One column, facet by rows\n",
    "yearly_sex_weight <- surveys_complete %>%\n",
    "    group_by(year, sex, species_id) %>%\n",
    "    summarize(avg_weight = mean(weight))\n",
    "ggplot(data = yearly_sex_weight, aes(x=year, y=avg_weight, color = species_id, group = species_id)) +\n",
    "    geom_line() +\n",
    "    facet_grid(sex ~ .)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# One row, facet by column\n",
    "ggplot(data = yearly_sex_weight, aes(x=year, y=avg_weight, color = species_id, group = species_id)) +\n",
    "    geom_line() +\n",
    "    facet_grid(. ~ sex)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Customization: **\n",
    "\n",
    "Cheatsheet: https://www.rstudio.com/wp-content/uploads/2016/11/ggplot2-cheatsheet-2.1.pdf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*Customize title and axis' titles*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = yearly_sex_counts, aes(x = year, y = n, color = sex, group = sex)) +\n",
    "    geom_line() +\n",
    "    facet_wrap(~ species_id) +\n",
    "    labs(title = 'Observed species in time',\n",
    "         x = 'Year of observation',\n",
    "         y = 'Number of species') +\n",
    "    theme_bw()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*Change font size*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = yearly_sex_counts, aes(x = year, y = n, color = sex, group = sex)) +\n",
    "    geom_line() +\n",
    "    facet_wrap(~ species_id) +\n",
    "    labs(title = 'Observed species in time',\n",
    "        x = 'Year of observation',\n",
    "        y = 'Number of species') +\n",
    "    theme_bw() +\n",
    "    theme(text=element_text(size=16, family=\"Arial\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*Change label orientation:*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data = yearly_sex_counts, aes(x = year, y = n, color = sex, group = sex)) +\n",
    "    geom_line() +\n",
    "    facet_wrap(~ species_id) +\n",
    "    labs(title = 'Observed species in time',\n",
    "        x = 'Year of observation',\n",
    "        y = 'Number of species') +\n",
    "    theme_bw() +\n",
    "    theme(axis.text.x = element_text(colour=\"grey20\", size=12, angle=90, hjust=.5, vjust=.5),\n",
    "                        axis.text.y = element_text(colour=\"grey20\", size=12),\n",
    "          text=element_text(size=16, family=\"Arial\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "** Save plots to files: **"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "my_plot <- ggplot(data = yearly_sex_counts, aes(x = year, y = n, color = sex, group = sex)) +\n",
    "    geom_line() +\n",
    "    facet_wrap(~ species_id) +\n",
    "    labs(title = 'Observed species in time',\n",
    "        x = 'Year of observation',\n",
    "        y = 'Number of species') +\n",
    "    theme_bw() +\n",
    "    theme(axis.text.x = element_text(colour=\"grey20\", size=12, angle=90, hjust=.5, vjust=.5),\n",
    "                        axis.text.y = element_text(colour=\"grey20\", size=12),\n",
    "          text=element_text(size=16, family=\"Arial\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "current_dir <- getwd()\n",
    "output_graph_dir <- 'output_graph'\n",
    "\n",
    "if (!file.exists(output_graph_dir)){\n",
    "    dir.create(file.path(current_dir, output_graph_dir))\n",
    "} else {\n",
    "    print (\"Directory already exists\")\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ggsave(file.path(current_dir, output_graph_dir,\"yearly_sex_counts.png\"), \n",
    "       my_plot, width=15, height=10)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R 3.4.2 (Juniper)",
   "language": "R",
   "name": "juniper_r3.4.2"
  },
  "language_info": {
   "codemirror_mode": "r",
   "file_extension": ".R",
   "mimetype": "text/x-r-source",
   "name": "R",
   "pygments_lexer": "r",
   "version": "3.4.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}