Exercise10.Rmd

# Advancing Further in R {#ch10}

```{r html-doc, child="html_doc.Rmd", echo = FALSE}
```

```{r setup, include=FALSE}
source("setup.R")
library(tidyverse)
```

Throughout these tutorials, we have introduced you to evolutionary genetic concepts using R. By now, you should have some familiarity with the versatility of R for data analysis and what is possible with it. We have gone from manipulating vectors and dataframes to processing genome-resequencing data and calculating population genomic statistics. With a course as broad as this, we understand it is difficult to feel like you fully understand every aspect of the analyses you are conducting - there is obviously a lot to learn! Ultimately, it is impossible for you to achieve this within a single course, using R properly takes experience and practice. Indeed, we'd argue that there is no real mastery of R, it is a tool with which you are always able to learn new things. We learn all the time in R and we have been using it for quite some time!

Nonetheless, there are obviously basics that you can master and that you can build upon in your own work, research and analysis. You are already familiar with many of the most important of these - data structures, how to use functions, how to manipulate and visualise data. How can we go further, to giving you an introduction to some more advanced techniques in R? In this final session, will take a step back from population and evolutionary genetics to focus once more on how R itself works, except this time we will focus mainly on programming. We will return to some R programming topics you have actually already encountered but in more detail, with more of a focus on explaining them piece-by-piece.

### What to expect {.unnumbered}

In this section we will:

-   learn about some advanced features of RStudio
-   learn about joining data sets
-   learn more about vectorisation, and how to use `lapply()` and `sapply()`

### Getting started {.unnumbered}

As always, we need to set up our R environment. We'll load `tidyverse` as usual and that's the only package we will use today.

```{r, eval = FALSE, echo = TRUE, results = "hide", message = FALSE, warning = FALSE}
# clear the R environment
rm(list = ls())

# load packages
library(tidyverse)
```

## Advanced features of RStudio

```{r html-doc, child="html_doc.Rmd", echo = FALSE}
```

Back in the [introductory session](#ch01), we learned that RStudio is a front-end for R. More specifically, RStudio is an **integrated development environment** that makes it more straightforward to do a lot [more than just interact with the console](https://www.rstudio.com/products/rstudio/features/). Much of these features are specifically designed with programming in mind and so we will give a short introduction to them.

### Projects, projects, projects

One of the most useful features of RStudio is it's ability to separate your work into different projects. Throughout the tutorials, you might have amassed a large number of files (i.e. datasets, scripts and so on) with little obvious organisation. You might have also had to set your working directory repeatedly.

This can become very tedious when you are working on multiple things at once in R. Instead, you can very easily subdivide your work into different RStudio projects, which are essentially separate working directories. You can make brand new ones or associate them with existing work.

You can see here for [a more detailed explanation on how to manage projects in RStudio.](https://support.rstudio.com/hc/en-us/articles/200526207-Using-Projects)

### Everyone has a history

Another useful feature of RStudio (and also the standard R distribution) is the fact that every command, function and input you type to the console is stored as a history. You can access this very easily using the function `history`. For example:

```{r, eval = FALSE, echo = TRUE, results = 'hidden', message = FALSE}
# show recent command history
history()
```

In RStudio, this command will actually open the `history` pane and show you previous commands that you have run. You can even use the buttons at the top of this pane to reread commands into the console (i.e. `To Console`) or copy them into an R script (`To Source`).

There is also an even easier way for you to quickly access the history or previously run commands when you are in the console. You can simply press the `up` key on your keyboard to rerun the previous lines. This can save you a lot of time and retyping.

### Tab complete and other hotkeys

When typing code into a script or the R console, RStudio allows you to do something called **tab completion**. This simply means if you hit the `tab` key while typing a function or object name, it will give you a list of options. For example, try typing `pl` into the R console and then pressing tab. You will see a list of available functions and you can then select the function you want from them.

Tab completion also works within a function too. Try using tab complete to call the function `plot` and then pressing it again inside the brackets. You will see the arguments the function requires and if you hover over them, a small pop-up will give you a brief explanation of what they are.

There are actually a large number of keyboard based shortcuts you can make use of with Rstudio - these include useful things like `Ctrl + 1` or `Ctrl + 2` to switch between the console and script panes. You can see [a whole list of them here](https://support.rstudio.com/hc/en-us/articles/200711853-Keyboard-Shortcuts) or by selecting `Tools > Keyboard Shortcuts Help` from the menu bar at the top of the program.


## More on data handling: categorising and joining

```{r html-doc, child="html_doc.Rmd", echo = FALSE}
```


Handling large amounts of data with few lines of code is one of R's strong points. In this section we will show how you can use `ifelse()` and the `left_join()` function from dplyr to make categories and add information to your data

### `ifelse()` for making categories

Let's first review what the `ifelse()` function does. As you learned back in [chapter 5](#ifelse), `ifelse()` takes a logical statement, and returns something different depending on whether the condition is `TRUE` or `FALSE`:

```{r, eval = FALSE, echo = TRUE, results = 'hide', message = FALSE}
# make a small vector
y <- c(20, 30, 50)
# use ifelse to evaluate it
ifelse(y > 25, "Greater than 25", "Not greater than 25")
```

This function can be extremely useful for creating new variables in datasets. Let's return to the familiar `starwars` data from `dplyr` in order to use the function in this way.

```{r, eval = TRUE, echo = TRUE, results = 'hidden', message = FALSE}
starwars
```

Now we can take a look at the `starwars$species` vector. There are a lot of different species, so what if we wanted to create a vector that simply states whether an individual is a droid or not?

```{r, eval = FALSE, echo = TRUE, results = 'hidden', message = FALSE}
# if else to identify droids and non-droids
ifelse(starwars$species == "Droid", "Droid", "Non-Droid")
```

This can be useful e.g. for comparing droids and non-droids in a plot or table. Say we want to label our species based on whether they are a droid, human or neither of the two. A useful thing with `ifelse()` is that the third argument can be *another* `ifelse()` function! So we can actually chain `ifelse()` commands like this:

```{r, eval = FALSE}
ifelse(starwars$species == "Droid", "Droid", ifelse(starwars$species == "Human", "Human", "Neither human nor droid"))
```

This is useful, but quickly becomes convoluted. Imagine how the code would look if we threw in a third and fourth category there! In cases like this, remember to use linebreaks to make the code more readable. You can have a linebreak anywhere after starting a function, and R will still understand that it's part of the same function. A suggestion for better formatting than above:

```{r, eval = FALSE}
ifelse(
  starwars$species == "Droid", "Droid", 
  ifelse(
    starwars$species == "Human", "Human", 
    "Neither human nor droid"
    )
  )
```

Still, if you have more than, say, four-five categories, this becomes difficult to read and time-consuming. For e.g. adding more information to a data frame, *joining* may be a better alternative, which we will go through next.

### Joining

For this section we will revisit the [copepods.txt](https://BIOS1140.github.io/data/copepods.txt) data that we encountered way back in [week 2](#w02). Start by reading in this data. You should know enough by now to do this by yourself, so we won't show you how.

```{r echo = FALSE}
copepods <- read.table("docs/data/copepods.txt", header = TRUE)
```

```{r}
copepods
```

Next, we use `pivot_longer()` to get all taxa in a single column, i.e., convert to long format. See if you manage to do this yourself before looking at my code below.

::: {.fold .c}
```{r}
copepods_long <- pivot_longer(copepods, 
                              -depth, 
                              names_to = "taxon", 
                              values_to = "count")
```
:::

```{r}
copepods_long
```

Now, say that you have recorded the temperature at each depth, and want to add that information to your copepod data. How would you go about doing that? First, here is the data in a data frame:

```{r}
temps <- data.frame(
  depth = c(0,2,4,6,8,10,12,14,16),
  temp = c(15.5, 15.4, 15.2, 14.7, 11.4, 8.3, 7.6, 7.0, 6.8)
)
temps
```

One way would be using nested `ifelse()` functions, like we learned in the previous section. This is a lot of work and doesn't look good, but it's written it out below just to show you it's possible:

::: {.fold .c}
```{r}
copepods_long$depthtemp <- ifelse(
  copepods_long$depth == 0, 15.5,
  ifelse(
    copepods_long$depth == 2, 15.4,
    ifelse(
      copepods_long$depth == 4, 15.2,
      ifelse(
        copepods_long$depth == 6, 14.7,
        ifelse(
          copepods_long$depth == 8, 11.4,
          ifelse(
            copepods_long$depth == 10, 8.3,
            ifelse(
              copepods_long$depth == 12, 7.6,
              ifelse(
                copepods_long$depth == 14, 7.0,
                ifelse(
                  copepods_long$depth == 16, 6.8,
                  NA
                )
              )
            )
          )
        )
      )
    )
  )
)


```
:::

Instead, you can use the `left_join()` function from `dplyr`. You have to supply it the original data, the data you want to join with, and a vector of column names to join by (here "depth").

```{r}
copepods_temp <- left_join(copepods_long, temps, by = "depth")
copepods_temp
```

`left_join()` matches one or more columns in your two data sets, and add rows from the second data set into the first data set in the correct place. You see that the `temp` column is equal to the `depthtemp` we created earlier, but it's so much easier to work with! Keep in mind that it is this simple in our case because `depth` has the exact same name in both data frames. Remember this when recording data in the future!

::: {.green}
**Important concept:**\
If you have your data spread out over multiple files, remember to name columns appropriately. All columns that contain the same kind of data should have the *exact same name* across all data sets. Similarly, the data should be entered in the same way in both data sets (e.g., don't record depth as "2" in one data set and "2m" in the other). If you do this, you can easily join data sets with the `left_join()` function.
:::

## Lists

In this course, we have stored large amounts of data in two types of objects: vectors and data frames. These have some inherent limitations:

-   A vector can only contain one type of value. This means, for example, that one vector can't contain both numeric values AND character values.
-   In a data frame, all columns need to be of the same length. In addition, a column can only contain a single type of value.

A list in R, on the other hand, has none of these limitations. It can contain values of multiple types and different lengths. This makes it ideal for storing more complex objects. In fact, much of the complex data you have encountered in these tutorials, like the phylogenetic trees last week, are based on lists.

You can make a list with the `list()` function:

```{r}
my_list <- list(names = c("Ask", "Embla"),
                numbers = c(11, 15, 6, 8),
                condition = TRUE)
my_list
```

You can access list elements with `$`, or with double square brackets `[[]]`:

```{r}
# access list element 1
my_list$names
my_list[[1]]
```

You can access elements within the vector within the list by adding a second set of (single) square brackets.

```{r}
# access element 1 of the vector in list element 1
my_list[[1]][1]
```

Since a list can contain any type of data (even data frames and other lists), there are times when they are very useful. We won't go into much more detail about lists here, but it is important that you know that this structure exists and a bit about how to work with them.

## Vectorisation {#vectorisation2}

```{r html-doc, child="html_doc.Rmd", echo = FALSE}
```

During the tutorials in this course, you have encountered the term "vectorisation" a few times(e.g. [chapter 5](#vectorisation)). In short, this is the concept of using functions and doing operations on a vector of values, in contrast to looping over each value and calculating each time. In many cases, it can make calculations faster and more readable.

Just to remind you, here is the example we used in chapter 5:

```{r}
# define function to calculate volume of a cylinder
cyl_vol <- function(r, h){
  return(pi * r^2 * h)
}

# use the function on single values
r <- 3
h <- 8
cyl_vol(r, h)

# use the function on vectors of values
r <- 10:1
h <- 1:10
cyl_vol(r, h)
```

An equivalent for-loop to the last example would be:

```{r}
# create empty vector
vols <- rep(NA, length(r))
# do calculations in loop
for (i in 1:length(r)){
  vols[i] <- cyl_vol(r[i], h[i])

}
vols
```

Notice how much simpler the vectorised solution looks. For simple operations like these, vectorising instead of looping is a no-brainer. Now, we will look at some more ways to vectorise, with the `sapply()` and `lapply()` functions.

### Using `lapply()` to vectorise

In [chapter 5](#vectorisation), you learned about the `apply()` function, for applying a function to a data frame or matrix. The `lapply()` function work similarly, but on vectors and lists instead.

Say you have a simple function that takes a vector as input, and rescales each element to the proportion of the total by dividing on the sum of the vector. In other words, calculate and return $\frac{x}{sum(x)}$.

**Exercise: Make a function that does the above, and test it with a couple of different vectors.**

::: {.fold .s}
```{r}
# define function
calc_proportion <- function(x){
  return(x/sum(x))
}

# test it
calc_proportion(c(1, 3, 6))
calc_proportion(c(17, 36, 24, 55))
```
:::

Now, what if we want to do this with several vectors? We may for example have a list of vectors, and want to apply this function to all of them:

```{r}
numbers_list <- list(
  c(1, 3, 6),
  c(17, 36, 24, 55),
  c(100, 500, 400, 38, 75)
)
numbers_list
```

If we simply try to use the function on the list directly, we get an error. R has no way of using `sum()` on a list!

```{r, error = TRUE}
calc_proportion(numbers_list)
```

This is where `lapply()` comes in. `lapply()` (or "list apply") takes a list or vector as input, and applies a function to each element of the list/vector. The syntax is `lapply(list, function)`. Like with `apply()`, the function should not have parentheses after its name. Let's use this to apply `calc_proportion()` to `numbers_list`:

```{r}
lapply(numbers_list, calc_proportion)
```

Now we get a list of our rescaled vectors! Note that `lapply()` always returns a list (hence "list apply"). Like earlier examples we've seen, it doesn't matter if `numbers_list` has a single value or a million, the code would still look the same[^exercise10-1].

[^exercise10-1]: Of course, this can also be achieved with a for-loop, which looks like this:

    ```{r}
    # create empty list, confusingly with the vector() function ...
    proportion_list <- vector(mode = "list", length = length(numbers_list))

    # loop and calculate
    for (i in 1:length(numbers_list)){
      proportion_list[[i]] <- calc_proportion(numbers_list[[i]])
    }
    proportion_list
    ```

    \ 

::: {.green}
**Important concept:**\
`lapply()` allows for fast and simple vectorisation when your function cannot be applied to your data directly. It takes a list/vector and a function as arguments, and returns a list of values.
:::

### `sapply()`

`sapply()` (or "simplified apply") is very similar to `lapply()` in that it takes a list or vector and a function as arguments, and applies the function to the list/vector. The key difference is that where `lapply()` always returns a list, `sapply()` returns the simplest possible object. For example, say we want to use `sum()` on each vector in `numbers_list`. With `lapply()` we get:

```{r}
lapply(numbers_list, sum)
```

However, each list element only contains one value. Wouldn't it be easier to have it stored in a vector? This is what we get with `sapply()`:

```{r}
sapply(numbers_list, sum)
```

::: {.green}
**Important concept:**\
`sapply()` can be convenient to simplify the output of a vectorised operation. Use `sapply()` if you want the simplest data structure possible, and `lapply()` when you want to be sure that your operation returns a list.
:::

### Anonymous functions

A concept that is often used together with `lapply()` and `sapply()` is *anonymous functions*. Like the term implies, this is a function without a name. For example, instead of creating the `calc_proportion()` function above, we could have defined the function **inside** our `lapply()` function like this[^exercise10-2]:

[^exercise10-2]: Note that the function here is defined without the use of curly brackets and `return()`. This is not unique to anonymous functions (but most often used there), but an equivalent way to define a function in R. This means that all these four functions are exactly equal:

    ```{r}
    f1 <- function(x) x/sum(x)
    f2 <- function(x) return(x/sum(x))
    f3 <- function(x){
      x/sum(x)
    }
    f4 <- function(x){
      return(x/sum(x))
    }
    ```

    When the function is more than a single line, you should always use the curly brackets. I personally prefer to use an explicit `return()` in more complicated functions, but you will encounter both in code you read.

```{r}
lapply(numbers_list, function(x) x/sum(x))
```

It works in the exact same way, but you don't have to define a function beforehand, and is easy to understand for the reader right away. Only use this for very simple functions that you are only going to use once. Otherwise, defining a named function is better. It is good to know about anonymous functions anyway, as they are frequently used.

## Concluding remarks

We hope the more indepth focus on programming and more advanced use of R in this tutorial has given you more insight into how these aspects of R work. Ultimately, the best way to become familiar and experienced with these methods is to repeatedly use them, often on your own data, to solve your own issues. Programming is a challenge for many people in any language - and it can be especially challenging for biologists. The key thing to keep in mind is that you are now familiar with the basics. You might still feel daunted at the prospect of sitting down and writing a function or set of code completely from scratch. However as you have seen from these tutorials and also those linked at the end of each section, there is a wealth of resources available online that can demystify many of these concepts further. The most important thing for your own future work is that you can grapple with these tutorials and then take things a step further each time. Even for experienced users, it takes time and thought to code in R - we hope that these tutorials have set you on the trajectory to becoming more advanced users!

Good luck!

## Study questions

The study questions for week 9 are found [here](#w10). These are meant only to test yourself. Hand-in is not required and will not count on the required number of approved hands-in. We will, however, provide a hand-in folder for those who would like feedback from one of the group teachers.

## Going further

There are many, many resources out there for learning more about programming in R. We have collated a few good ones here:

-   [A series of tutorials exploring aspects of R Studio useful for coding](https://www.rstudio.com/online-learning/)
-   [Some good examples of control flow in R](http://resbaz.github.io/2014-r-materials/lessons/30-control-flow/)
-   [A short description of control flow, functions and vectorisation in R](https://ramnathv.github.io/pycon2014-r/learn/controls.html)
-   [A more advanced look at why vectorisation is worth pursuing](http://www.noamross.net/blog/2014/4/16/vectorization-in-r--why.html)
-   [A software carpentry course on creating functions in R](https://swcarpentry.github.io/r-novice-inflammation/02-func-R/)