Session1.Rmd

---
title: 'Geospatial Carpentry: Session 1: Intro to R'
author: "Aleksandra Wilczynska"
date: "2022-11-18"
output: 
  html_document:
     toc: true
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(
	echo = TRUE,
	message = FALSE,
	warning = FALSE
)

# Load libraries ----------------------------------------------------------
# Package names
packages <- c("tidyverse", "here")

# Install packages not yet installed
installed_packages <- packages %in% rownames(installed.packages())
if (any(installed_packages == FALSE)) {
  install.packages(packages[!installed_packages])
}

# Packages loading
invisible(lapply(packages, library, character.only = TRUE))

```


# Resources

- [Coordination document](https://docs.google.com/document/d/1x4mRcsSTyUueN_mCiZHE1TOTvF-OZRwq2gV_J9svPRI/edit#heading=h.rv3i711oni3d)
- [Introduction to R for Geospatial Data](https://datacarpentry.org/r-intro-geospatial/)

<!-- need to reduce this part by  15 minutes -> to 80 minutes -->

# [Introduction to  RStudio](https://datacarpentry.org/r-intro-geospatial/01-rstudio-intro/index.html) and [Project management](https://datacarpentry.org/r-intro-geospatial/02-project-intro/index.html) 


## Project management 
RStudio is an integrated development environment (IDE). It provides a (much prettier) interface for the R software. R is integrated into RStudio, so you never actually have to open R.


R Studio gives a functionality of creating projects: self-contained working space (i.e. working directory), to which R will refer to, when looking for and saving files. You can create projects in existing  directories( folders) or create a directory anew. 

We’re going to create a  project in RStudio in the existing directory:

- `File` 
- `New Project` 
- `Existing directory` 
- browse for directory directory you created when downloading the data:  `gdc`
-`Create project`


<!-- ### Tips for project organisation -->

<!-- + **Raw data** is read only -->
<!-- + **Cleaned data** as read only -->
<!-- + **Generated output** is disposable -->
<!-- + **Related data together** (Some GIS file formats are really 3-6 files that need to be kept together and have the same name, e.g. shapefiles) -->
<!-- + Keep consistent **naming schema** -->
<!-- + Name all files to reflect their content or function -->
<!-- + Stage your scripts: Creating separate R scripts or Rmarkdown documents for different stages of a project -->

### Organising project/working directory

This is one suggestion of how your R project can look like. Your data folder is already there. Let's go ahead and create the other folders. 

![](assets/img/rstudio_project_files.jpeg)

### Two main ways to interact with R

- Test and play within the **interactive R console** (chat)
    - **Pros:** immediate results
    - **Cons:** work lost once we close RStudio
    
- Start writing in an **.R file** ( email )
    - File still will be executed in the console
    - **Pros:** complete record of what you did!
    - **Cons:** Can be messy if we're just want to print things out

### Running the code

- In the console: `Enter`
- In the script:
   - `Ctrl` + `Enter`  (for MAC users: `Command` + `Enter`)
   - `Run` button on right left - current line or selection

### Creating a script
We're going to work with a script. Let's create one now and save it in the `scripts` directory.

- `File` 
- `New File` 
- `R Script` 
- A new `Untitled` script will appear in the source pane. Save it using floppy disc icon.
- Name it `intro-to-r.R`

<!-- ### Escaping  -->

<!-- The console shows it's ready to get new commands with `>` sign. It will show `+` sign if it still requires input for the command to be executed.  -->

<!-- Sometimes you don't know what is missing/ you change your mind and want to run something else, or your code is running much too long and you just want it to stop. The way to do it is to hit `Esc`. -->

### Packages 

A great power of R lays in **packages: add-on sets of functions** that are build by the community and once they go through a quality process they are available to download from a repository called CRAN. They need to be explicitly activated. 
Now, we will be using `tidyverse` package, which is actually a collection of useful packages. Another package that will be useful for us is `here`.

If you have have not installed this package earlier, please do so. You can check if you have it installed in the `Packages` pane in the bottom-right window.

```{r install-package, eval=FALSE}
# install.packages('tidyverse')
install.packages('here')
```

You need to install package only once, but you will need to load it each time you want to use its functionalities. To do that you use `library()` command:

```{r load-package}
library(tidyverse)
library(here)

```

### Handling paths 

![](assets/img/relative_root.png)

<font size="3">*Credit*:[kaggle.com](https://www.kaggle.com/code/rtatman/reproducibility-tips-absolute-vs-relative-paths/notebook)</font>


You have created a project which is your working directory, and a number of subfolders, that will help you organise your project better. But now each time you will save or retrieve a file from those folders, you will need to specify the path from the folder you are in (most likely `scripts`). 


That can become complicated and can become a reproducibility problem if the person using your code (e.g. future you) is working in a different subfolder. 


`here()` to the rescue! This package provides absolute paths from the root (main directory) of your project. 

![](assets/img/here.png)

<font size="3">*Credit*:[Allison horst](https://github.com/allisonhorst)</font>


```{r here}
here('data')
```


### Download files  

We still need to download data for the first part of the workshop. You can do with with the function `download.file()`. We will save it in the `data` folder, where the **raw** data should go. 

```{r download-files}
download.file('bit.ly/GeospatialGapminder', here('data','gapminder_data.csv'), mode = 'wb')

```

<!-- #### Part 2:  -->
<!-- + Download zip file by clicking on [this link](https://ndownloader.figshare.com/articles/2009586/versions/10) -->
<!-- + Move the downloaded zip file to the data directory. -->
<!-- + Once the data have been moved, unzip all files.  -->

# [Intro to R](https://datacarpentry.org/r-intro-geospatial/01-rstudio-intro/index.html)

## Use R as a calculator

```{r calculator}
1+100


```

## Variables and assignment

We can store values in variables using the assignment operator `<-`, like this:

```{r asignment-operator}
x <- 1/40
```

Notice that assignment does not print a value. Instead, we stored it for later in something called a variable. `x` now contains the value `0.025`:
```{r asignment-operator2}
x
```

Look for the `Environment` tab in one of the panes of RStudio, and you will see that `x` and its value have appeared. 
Our variable `x` can be used in place of a number in any calculation that expects a number, e.g. when caclulating a square root:
```{r use-variable}
sqrt(x)
```

Variables can be also reassigned:
```{r reassign}
x <- 100
x
```

You can use the 'old' value when reassigning the value
```{r reassign2}
y <- sqrt(x) # you can use value stored in object x to create y
y
```


# [Data Structures](https://datacarpentry.org/r-intro-geospatial/03-data-structures-part1/index.html) 

## Vectors 
So far we've looked on individual values. Now we will move to a data structure called vectors. Vectors are arrays of values of a same data type (will explain in a second :) ) .

You can create a vector with a `c()` function. 

```{r vectors}

numeric_vector <- c(2, 6, 3) # vector of numbers - numeric data type.
numeric_vector

character_vector <- c('banana', 'apple', 'orange') # vector of words - more precisely strings of characters- character data type
character_vector

logical_vector <- c(TRUE, FALSE, TRUE) # vector of logical values (is something true or false?)- logical data type.
logical_vector

```

### Combining vectors 

The combine function, `c()`, will also append things to an existing vector:

```{r combine-vectors}

ab_vector <- c('a', 'b')
ab_vector

abcd_vector <- c(ab_vector, 'c', 'd')
abcd_vector

```

### Missing values  

A common operation you want to perform is to remove all the missing values (in R denoted as `NA`). Let's have a look how to do it: 

```{r remove-na}
with_na <- c(1, 2, 1, 1, NA, 3, NA ) # vector including missing value
```

First, let's try to calculate mean for the values in this vector
```{r remove-na1}
mean(with_na) # mean() function cannot interpret the missing values

mean(with_na, na.rm = T) # You can add the argument na.rm=TRUE to calculate the result while ignoring the missing values.
```

However, sometimes, you would like to have the `NA` completely removed from your vector. for this you need to identify which elements of the vector hold missing values with `is.na()` function. 

```{r remove-na2}
is.na(with_na) #  This will produce a vector of logical values, stating if a statement 'This element of the vector is a missing value' is true or not
!is.na(with_na) # # The ! operator means negation ,i.e. not is.na(with_na)

```

We know which elements in the vectors are `NA`. Now we need to retrieve the subset of the `with_na` vector that is not `NA`.
Any subsetting in `R` is done with square brackets`[ ]`. 


```{r remove-na3}

without_na <- with_na[!is.na(with_na)] # this notation will return only the elements that have TRUE on their respective positions

without_na

```


## Factors (adapted from [Starting with Data](https://datacarpentry.org/r-socialsci/02-starting-with-data/index.html))

Another important data structure is called a factor. Factors look like character data, but are used to represent categorical information.

Factors create a structured relation between the different levels (values) of a categorical variable, such as days of the week or responses to a question in a survey. While factors look (and often behave) like character vectors, they are actually treated as integer vectors by R. So you need to be very careful when treating them as strings.

### Create factors
Once created, factors can only contain a pre-defined set of values, known as levels. 

```{r factor-create}

nordic_str <- c('Norway', 'Sweden', 'Norway', 'Denmark', 'Sweden')
nordic_str # regular character vectors printed out

nordic_cat <- factor(nordic_str) # factor() function converts a vector to factor data type
nordic_cat # With factors, R prints out additional information - 'Levels'

```

### Inspect factors
R will treat each unique value from a factor vector as a **level** and (silently) assign numerical values to it. This come in handy when performing statistical analysis. You can inspect and adapt levels of the factor. 

```{r factor-inspect}
levels(nordic_cat) # returns all levels of a factor vector.  

nlevels(nordic_cat) # returns number of levels in a vector
```

### Reorder levels
Note that R sorts the levels in the alphabetic order, not in the order of occurrence in the vector. R assigns value of 1 to level 'Denmark', 2 to 'Norway' and 3 to 'Sweden'. This is important as it can affect e.g. the order in which categories are displayed in a plot or which category is taken as a baseline in a statistical model.

You can reorder the categories using `factor()` function.

```{r factor-reorder}

nordic_cat <- factor(nordic_cat, levels = c('Norway' , 'Denmark', 'Sweden')) # now Norway should be the first category, Denmark second and Sweden third

nordic_cat

str(nordic_cat) # you can also inspect vectors with str() function. In facto vectors, it shows the underlying values of each category. You can also see the structure in the environment tab of RStudio.
```
There is more than one way to reorder factors. Later in the lesson, we will use `fct_relevel()` function from `forcats` package to do the reordering.


### Note of caution 
Remember that once created, factors can only contain a pre-defined set of values, known as levels. 
It means that whenever you try to add something to the factor vector outside of this set, it will become an unknown/missing value detonated by R as `NA`.

```{r factor-missing-level}
nordic_str
nordic_cat2 <- factor(nordic_str, levels = c('Norway', 'Denmark'))
nordic_cat2 # since we have not included Sweden in the list of factor levels, it has become NA.
```


<!-- need to reduce this part by by 55 minutes -->
# [Exploring Data frames](https://datacarpentry.org/r-intro-geospatial/04-data-structures-part2/index.html)

Now we turn to the bread-and-butter of working with R: working with tabular data. In R data are stored in a data structure called **data frames**.  

A data frame is the representation of data in the format of a **table** where the columns are **vectors** that all have the **same length**. 


Because columns are vectors, each column must contain a **single type of data** (e.g., characters, integers, factors). 
For example, here is a figure depicting a data frame comprising a numeric, a character, and a logical vector.

![](assets/img/data-frame.svg)

## Reading data

`read.csv()` is a function used to read coma separated data files (`.csv` format)). There are other functions for files separated with other delimiters. 
We're gonna read in the gap minder data set with information about countries' size, GDP and average life expectancy in different years.

```{r reading-data}
gapminder <- read.csv(here('data','gapminder_data.csv') )

```

## Exploring dataset
Let’s investigate the gapminder data frame a bit; the first thing we should always do is check out what the data looks like.

It is important to see if all the variables (columns) have the data type that we require. Otherwise we can run into trouble.

```{r inspecting-data-str}
str(gapminder) 

```
We can see that the `gapminder` object is a data.frame with `r nrow(gapminder)` observations/ rows and `r ncol(gapminder)` variables/columns. 
In each line after a `$` sign, we see the name of each column, its type and first few values. 


There are multiple ways to explore a data set. Here are just a few examples

```{r}
head(gapminder) # see first 5  rows of the data set

summary(gapminder) # gives basic statistical information about each column. Information format differes by data type.

```


When you're analyzing a data set, you often need to access its specific columns.

One handy way to access a column is using it's name and a dollar sign `$`: 
```{r subset-dollar-sign}
country_vec <- gapminder$country  # Notation means: From dataset gapminder, give me column country. You can see that the column accessed in this way is just a vector of characters. 

head(country_vec)

```
Note that the calling a column with a `$` sign will return a vector, it's not a data frame anymore.

# [Data frame Manipulation with dplyr](https://datacarpentry.org/r-intro-geospatial/06-dplyr/index.html) 

## Select
Let's start manipulating the data. 

First we will adapt our dataset, by keeping only the columns we're interested in using the `select()` function from `dplyr` package:

```{r dplyr-select}
year_country_gdp <- select(gapminder, year, country, gdpPercap) 

head(year_country_gdp)

```

## Pipe
Now, this is not the most common notation when working with `dplyr` package. `dplyr` offers an operator `%>%` called a pipe, which allows you build up a very complicated commands in a readable way.


In newer installation of `R` you can also find a notation `|>` . This pipe does exactly the same, the only difference is that you don't need to load any pacakges to have it available.


The `select()` statement with pipe would look like that:

```{r dplyr-pipe}

year_country_gdp <- gapminder %>% 
  select(year,country,gdpPercap)

head(year_country_gdp)

```

First we define data set, then with the use of pipe we pass it on to the `select()` function. This way we can chain multiple functions together, which we will be doing now. 

## Filter

We already now how to select only the needed columns. But now, we also want to filter the data set via certain condition with `filter()` function. Instead doing it in separate steps , we can do it all together. 

In the `gapminder` data set, we want to see the results only for Europe for 21st century. 
```{r}
year_country_gdp_euro <- gapminder %>% 
  filter(continent == "Europe" & year> 2000) %>%
  select(year, country, gdpPercap)

head(year_country_gdp_euro)
```
Let's now find all the observations from Eurasia: 

```{r}
year_country_gdp_eurasia <- gapminder %>% 
  filter(continent == "Europe" | continent == "Asia") %>%
  select(year, country, gdpPercap)

head(year_country_gdp_eurasia)
```


<div class="alert alert-info">
<strong>Challenge</strong> 
Write a single command (which can span multiple lines and includes pipes) that will produce a dataframe that has the African values for life expectancy, country and year, but not for other Continents. How many rows does your dataframe have and why?

`countdown::countdown_fullscreen(minutes = 3)`

</div>
```{r ex5, class.source="bg-info", include = F}
year_country_lifeExp_Africa <- gapminder %>%
                           filter(continent=="Africa"  ) %>%
                           select(year,country,lifeExp)

``` 


## Group and summarize
So far, we have created a dataset for one of the continents represented in the `gapminder` dataset. But rather than doing that, we want to know statistics about all of the continents, presented by group.

```{r dplyr-group}
gapminder %>% # select the dataset
  group_by(continent) %>% # group by continent
  summarize(avg_gdpPercap = mean(gdpPercap)) # summarize function creates statistics for the data set 

```

<div class="alert alert-info">
<strong>Challenge</strong> 
Calculate the average life expectancy per country. Which country has the longest average life expectancy and which has the shortest average life expectancy?

<strong>Hint</strong> Use `max()`  and `min()` functions to find minimum and maximum.

`countdown::countdown_fullscreen(minutes = 3)`

</div>
```{r ex6 , class.source="bg-info", include=F}
lifeExp_bycountry <- gapminder %>%
   group_by(country) %>%
   summarize(avg_lifeExp=mean(lifeExp))


lifeExp_bycountry %>%
   filter(avg_lifeExp == min(avg_lifeExp) | avg_lifeExp == max(avg_lifeExp))
```


You can also group by multiple columns:

```{r dplyr-group-multi}

gapminder %>%
  group_by(continent, year) %>%
  summarize(avg_gdpPercap = mean(gdpPercap))

```

On top of this, you can also make multiple summaries of those groups:
```{r dplyr-summ}
gdp_pop_bycontinents_byyear <- gapminder %>%
  group_by(continent,year) %>%
  summarize(
    avg_gdpPercap = mean(gdpPercap),
    sd_gdpPercap = sd(gdpPercap),
    avg_pop = mean(pop),
    sd_pop = sd(pop),
    n_obs = n()
    )

```

## Frequencies

If you need only a number of observations per group, you can use the `count()` function
```{r dplyr-count}

gapminder %>%
    group_by(continent) %>%
    count()
```
 

## Mutate

Frequently you’ll want to create new columns based on the values in existing columns, for example to do unit conversions, or to find the ratio of values in two columns. For this we’ll use `mutate()`.

```{r dplyr-mutate}
gapminder_gdp <- gapminder %>%
  mutate(gdpBillion = gdpPercap*pop/10^9)

head(gapminder_gdp)

```


# [Introduction to Visualisation](https://datacarpentry.org/r-intro-geospatial/07-plot-ggplot2/index.html) 
Package `ggplot2` is a powerful plotting system. I will introduce key features of `ggplot`. Later today/ on Monday you will use this package to visualize geospatial data.
`gg` stands for grammar of graphics, the idea that three components needed to create a graph are:
- data
- aesthetics - coordinate system on which we map the data ( what is represented on x axis, what on y axis)
- geometries - visual representation of the data (points, bars, etc.)

fun part about `ggplot` is that you can then add additional layers to the plot providing more information and make it more beautiful. 

First, lets plot distribution of life expectancy in the `gapminder` dataset. 

```{r ggplot}
  ggplot(data =gapminder,  aes(x = lifeExp)) + # aesthetics layer 
  geom_histogram() # geometry layer

```

You can see that in `ggplot` you use `+` as a pipe, to add layers. Within `ggplot` call, it is the only pipe that will work. 
But, it is possible to chain operations on a dataset with a pipe that we have already learned: `%>%` ( or `|>`) and follow them by ggplot grammar. 

Let's create another plot, this time only on a subset of observations:

```{r ggplot-col}
gapminder %>%  # we select a dataset
  filter(year == 2007, 
         continent == 'Americas') %>% # and filter it to keep only one year and one continent
  ggplot(aes(x = country, y = gdpPercap)) + # we create aesthetics, both x and y axis represent values of  columns
  geom_col() # we select a column graph as a geometry
```

Now, you can iteratively improve how the plot looks. For example, you might want to flip it, to better display the labels.

```{r ggplot-coord-flip}
gapminder %>%  
  filter(year == 2007, 
         continent == 'Americas') %>% 
  ggplot(aes(x = country, y = gdpPercap)) + 
  geom_col()+ 
  coord_flip()
```

One thing you might want to change here is the order in which countries are displayed. It would be easier to compare GDP per capita, if theY were showed in order. 
To do that, we need to reorder factor levels (you remember, we've already done this before).
the order of the levels will depend on another variable - GDP per capita.

```{r ggplot-color}
gapminder %>%  
  filter(year == 2007, 
         continent == 'Americas') %>% 
  mutate(country = fct_reorder(country, gdpPercap )) %>%
  ggplot(aes(x = country , y = gdpPercap)) + 
  geom_col() +
  coord_flip()

```

Let's make things more colorful - let's represent the average life expectancy of a country by color

```{r ggplot-colors}
gapminder %>%  
  filter(year == 2007, 
         continent == 'Americas') %>% 
  mutate(country = fct_reorder(country, gdpPercap )) %>%
  ggplot(aes(x = country, y = gdpPercap, fill = lifeExp   )) + # fill argument for coloring surfaces, color for points and lines
  geom_col()+ 
  coord_flip()


```

We can also adapt the color scale. Common choice that is used for its colorblind-proofness is `viridis` package.
```{r ggplot-colors-adapt}
plot_c <-
  gapminder %>%  
  filter(year == 2007, 
         continent == 'Americas') %>% 
  mutate(country = fct_reorder(country, gdpPercap )) %>%
  ggplot(aes(x = country, y = gdpPercap, fill = lifeExp   )) + 
  geom_col()+ 
  coord_flip()+
  scale_fill_viridis_c() # _c stands for continous scale 

```

Maybe we don't need that much information about the life expectancy. We only want to know if it's below or above average.

```{r ggplot-colors-discrete}
plot_d <-  # this time let's save the plot in the object.
  gapminder %>%  
  filter(year == 2007, 
         continent == 'Americas') %>% 
  mutate(country_reordered = fct_reorder(country, gdpPercap ),
         lifeExpCat = if_else(lifeExp >= mean(lifeExp), 'high', 'low' )
         ) %>%
  ggplot(aes(x = country_reordered, y = gdpPercap, fill = lifeExpCat   )) + 
  geom_col()+ 
  coord_flip()+
  scale_fill_manual(values = c('light blue', 'orange')) 

```

Since we saved a plot as an object, nothing has been printed out. Just like with any other object in `R`, if you want to see it, you need to call it.  

```{r ggplot-call}
plot_d

```

Now we can make use of the saved object and add things to it.

Let's also give it a title and name the axes:
```{r ggplot-titles}
plot_d <- 
  plot_d +
  ggtitle('GDP per capita in Americas', subtitle = 'Year 2007') +
  xlab('Country')+
  ylab('GDP per capita')

plot_d
```

# [Writing data](https://datacarpentry.org/r-intro-geospatial/08-writing-data/index.html)

Once we are happy with our plot we can save it in a format of our choice. Remember to save it in the dedicated folder. 

```{r save plot}
ggsave(plot = plot_d, filename = here('fig_output','plot_americas_2007.pdf') ) # By default, ggsave() saves the last displayed plot, but you can also explicitly name the plot you want to save

```


Another output of your work you want to save is a cleaned dataset. In your analysis, you can then load directly that dataset. Say we want to save the data only for Australia:

```{r writing-data}
gapminder_amr_2007 <-  gapminder %>%  
  filter(year == 2007, 
         continent == 'Americas') %>% 
  mutate(country_reordered = fct_reorder(country, gdpPercap ),
         lifeExpCat = if_else(lifeExp >= mean(lifeExp), 'high', 'low' )) 

write.csv(gapminder_amr_2007, 
          here('data_output', 'gapminder_americas_2007.csv'), 
          row.names=FALSE)


```