Basics of Spatial Data Analysis

class: center, middle, inverse, title-slide

# Basics of Spatial Data Analysis
### Jakub Nowosad
### Workshop, 2019-09-27

---

class: inverse, left, bottom
# Ecosystem(s)

---
# Ecosystem(s)

<div class="figure" style="text-align: center">
<img src="figs/rb.png" alt="From Roger Bivand's slides 'A practical history of R'" width="80%" />
<p class="caption">From Roger Bivand's slides 'A practical history of R'</p>
</div>

---
# Ecosystem(s)

<div class="figure" style="text-align: center">
<img src="figs/rb2.png" alt="From Roger Bivand's slides 'A practical history of R-sig-geo'" width="80%" />
<p class="caption">From Roger Bivand's slides 'A practical history of R-sig-geo'</p>
</div>

---
# Ecosystem(s)

.pull-left[
R packages:

- **sf**, **raster** - spatial classes
- **dplyr**, **rmapshaper** - processing of atttribute tables/geometries 
- **rnaturalearth**, **osmdata**, **getlandsat** - spatial data download
- **rgrass7**, **RQGIS**, **RSAGA**, **link2GI** - connecting with GIS software
- **gstat**, **mlr**, **CAST** - spatial data modeling
- **rasterVis**, **tmap**, **ggplot** - static visualizations
- **leaflet**, **mapview**, **mapdeck** - interactive visualizations
- many more...
]

.pull-right[

]

---
# Ecosystem(s)

.pull-left[
Spatial classes:

- **sf**
- **raster**
]

.pull-right[
Visualizations:

- **tmap**
]

---
# Ecosystem(s)

.pull-left[
Spatial classes:

- **sf**
- **raster**
- **sp**
- **stars**
- **terra**
- **spatial**
- ...
]

.pull-right[
Visualizations:

- **tmap**
- **ggplot** 
- **rasterVis**
- **leaflet**
- **mapview**
- **mapdeck**
- ...
]

---
# Ecosystem(s)

<div class="figure" style="text-align: center">
<img src="figs/sf-infrastructure.png" alt="Pebesma, Edzer. Simple features for R: standardized support for spatial vector data. The R Journal 10.1 (2018): 439-446." width="70%" />
<p class="caption">Pebesma, Edzer. Simple features for R: standardized support for spatial vector data. The R Journal 10.1 (2018): 439-446.</p>
</div>

---
# Prerequisites

.pull-left[
We can install the packages used in the workshop as follows:

```r
pkgs = c(
  "sf",                
  "raster",            
  "dplyr",            
  "tmap",
  "spData",
  "usethis"
)
to_install = !pkgs %in% installed.packages()
if(any(to_install)) {
  install.packages(pkgs[to_install])
}
```
]

.pull-right[
The code below can be used to download all the scripts and data for the workshop:

```r
usethis::use_course("http://bit.ly/whyr19w")
```
]

---
# Exercises

1. Think about how R can help you? 
Why do you want to learn how to use it for spatial data analysis?
2. Have a look at [CRAN Task View: Analysis of Spatial Data](https://cran.r-project.org/web/views/Spatial.html) and check if there are any packages that may be useful to you in your daily work?

---
class: inverse, left, bottom
# Spatial data representation

---
# Spatial data representation

**sf**:
- The **sf** package is the successor of the **sp** package based on the OGC standard Simple Features.
- Combines the functionality of three previous packages: **sp**, **rgeos** and **rgdal**.
- Most of the functions in this package start with a prefix `st_`.
- This package handles additional vector data types (e.g. polygon and multipolygon are two separate classes), allows for easier data processing, and support sfor spatial databases such as PostGIS.
- https://r-spatial.github.io/sf/ and https://github.com/rstudio/cheatsheets/blob/master/sf.pdf

**raster**:
- The **raster** package contains classes and methods representing raster objects and operations.
- It allows raster data to be loaded and saved.
- It allows raster algebra and raster processing.
- It includes a number of additional functions, e.g. for analysis of terrain characteristics.
- It allows you to work on large sets of data.
- ?`raster-package` and http://www.rpubs.com/etiennebr/visualraster

---
# Spatial data reading - vector

.pull-left[

```r
library(sf)
poland = st_read("data/poland.gpkg")
```

```
## Reading layer `polska' from data source `/home/jn/Science/conferences/whyr_19/workshop/data/poland.gpkg' using driver `GPKG'
## Simple feature collection with 1 feature and 3 fields
## geometry type:  POLYGON
## dimension:      XY
## bbox:           xmin: 171842.5 ymin: 132323.8 xmax: 861798.4 ymax: 775305.9
## epsg (SRID):    NA
## proj4string:    +proj=tmerc +lat_0=0 +lon_0=19 +k=0.9993 +x_0=500000 +y_0=-5300000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
```
- The `st_write()` function can be used for spatial vector data writing.
]

.pull-right[

```r
plot(poland)
```

]

---
# Spatial data structure - vector

```r
poland
```

```
## Simple feature collection with 1 feature and 3 fields
## geometry type:  POLYGON
## dimension:      XY
## bbox:           xmin: 171842.5 ymin: 132323.8 xmax: 861798.4 ymax: 775305.9
## epsg (SRID):    NA
## proj4string:    +proj=tmerc +lat_0=0 +lon_0=19 +k=0.9993 +x_0=500000 +y_0=-5300000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
##     name  pop_est gdp_md_est                           geom
## 1 Poland 38476269    1052000 POLYGON ((487928.3 182549, ...
```

---
# Spatial data reading - raster

.pull-left[

```r
library(raster)
dem = raster("data/poland_srtm.tif")
```

- The `writeRaster()` function can be used for spatial raster data writing.
]

.pull-right[

```r
plot(dem)
```

<img src="index_files/figure-html/unnamed-chunk-11-1.png" style="display: block; margin: auto;" />
]

---
# Spatial data structure - raster

```r
dem
```

```
## class      : RasterLayer 
## dimensions : 702, 1202, 843804  (nrow, ncol, ncell)
## resolution : 0.008333333, 0.008333333  (x, y)
## extent     : 14.125, 24.14167, 48.99167, 54.84167  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 
## source     : /home/jn/Science/conferences/whyr_19/workshop/data/polska_srtm.tif 
## names      : polska_srtm 
## values     : -9, 2220  (min, max)
```

---
# Exercises

1. Read the `data/wlkp_pn.gpkg` file containing the borders of the Wielkopolski National Park into R.
Display this object and view its structure.
What can you say about the contents of this file? 
What type of data does it store? 
What is the coordinate system used?
How many attributes does it contain?
What is its geometry?
2. Read the `data/wlkp_pn_dem.tif` file containing the digital elevation model for the Wielkopoland National Park area into R.
Display this object and view its structure.
What can you say about the contents of this file? 
What type of data does it store? 
What is the coordinate system used?
How many attributes does it contain?
How many dimensions does it have?
What is the data resolution?

---
class: inverse, left, bottom
# Spatial visualization

---
# Spatial visualization

```r
library(tmap)
```

- Thematic mapping package
- It allows you to create static maps, animated maps and interactive maps.
- It works by adding subsequent layers to visualize and then modifying them.
- https://github.com/mtennekes/tmap - list of additional materials about the **tmap** package

---
# Spatial visualization

Open the `code/spatial_vis.R` file. 
Run the code and:

1. Change the map title from "My map" to "New Zealand".
2. Update the map credits with your own name and today's date.
3. Change the color palette to "BuGn". 
(You can also try other palettes from http://colorbrewer2.org/)
4. Put the north arrow in the top right corner of the map.
5. Improve the legend title by adding the legend units.
6. Increase the number of breaks in the scale bar.
7. Change the borders' color of the New Zealand's regions to black. 
Decrease the line width.
8. Change the background color to any color of your choice.

---
# Spatial visualization

.left-code[

```r
tm_shape(poland) + 
  tm_polygons()
```
]

.right-plot[
![](index_files/figure-html/tm1-1.png)
]

---
# Spatial visualization

.left-code[

```r
tm_shape(poland) + 
  tm_polygons() + 
  tm_scale_bar(position = c("left", 
                            "bottom")) 
```
]

.right-plot[
![](index_files/figure-html/tm2-1.png)
]

---
# Spatial visualization

.left-code[

```r
tm_shape(poland) + 
  tm_polygons() + 
  tm_scale_bar(position = c("left", 
                            "bottom")) + 
  tm_compass()
```
]

.right-plot[
![](index_files/figure-html/tm3-1.png)
]

---
# Spatial visualization

.left-code[

```r
tm_shape(poland) + 
  tm_polygons() + 
  tm_scale_bar(position = c("left", 
                            "bottom")) + 
  tm_compass() +
  tm_layout(title = "Poland")
```
]

.right-plot[
![](index_files/figure-html/tm4-1.png)
]

---
# Spatial visualization

```r
tmap_mode("view")
```

```r
tm_shape(poland) + 
  tm_polygons() + 
  tm_layout(title = "Poland")
```

```r
tmap_mode("plot")
```

---
# Spatial visualization

.left-code[

```r
tm_shape(dem) + 
  tm_raster()
```
]

.right-plot[
![](index_files/figure-html/tm5-1.png)
]

---
# Spatial visualization

.left-code[

```r
tm_shape(dem) + 
  tm_raster(style = "cont",
            midpoint = NA)
```
]

.right-plot[
![](index_files/figure-html/tm6-1.png)
]

---
# Spatial visualization

.left-code[

```r
tm_shape(dem) + 
  tm_raster(style = "cont", 
            midpoint = NA,
            palette = "-RdYlGn")
```

```r
tmaptools::palette_explorer()
```
]

.right-plot[
![](index_files/figure-html/tm7-1.png)
]

---
# Spatial visualization

.left-code[

```r
tm_shape(dem) + 
  tm_raster(style = "fixed",
            breaks = c(-99, 0, 300,
                       600, 9999), 
            midpoint = NA,
            palette = "-RdYlGn",
            title = "") + 
  tm_layout(legend.position = c("LEFT", 
                                "BOTTOM"))
```
]

.right-plot[
![](index_files/figure-html/tm8-1.png)
]

---
# Spatial visualization

.left-code[

```r
tm_shape(dem) + 
  tm_raster(style = "fixed",
            breaks = c(-99, 0, 300,
                       600, 9999), 
            labels = c("Depressions", 
                       "Plains", 
                       "Hills",
                       "Mountains"),
            midpoint = NA,
            palette = "-RdYlGn",
            title = "") + 
  tm_layout(legend.position = c("LEFT", 
                                "BOTTOM"))
```
]

.right-plot[
![](index_files/figure-html/tm9-1.png)
]

---
# Spatial visualization

.left-code[

```r
tm_shape(dem) + 
  tm_raster(style = "fixed",
            breaks = c(-99, 0, 300,
                       600, 9999), 
            labels = c("Depressions", 
                       "Plains", 
                       "Hills",
                       "Mountains"),
            midpoint = NA,
            palette = c("#5E8B73",
                        "#DAE97A", 
                        "#EADC70", 
                        "#AF8D5C"),
            title = "") + 
  tm_layout(legend.position = c("LEFT",
                                "BOTTOM"))
```
]

.right-plot[
![](index_files/figure-html/tm10-1.png)
]

---
# Spatial visualization

```r
map1 = tm_shape(dem) + 
  tm_raster(style = "fixed",
            breaks = c(-99, 0, 300, 600, 9999),
            labels = c("Depressions", "Plains", "Hills", "Mountains"),
            midpoint = NA, 
            palette = c("#5E8B73", "#DAE97A", "#EADC70", "#AF8D5C"),
            title = "") +
  tm_layout(legend.position = c("LEFT", "BOTTOM"))
```

```r
tmap_save(map1, "my_first_map.png")
```

---
# Exercises

1. Improve the `map1` object by e.g. adding a scale, north arrow, or title.
Also try to add the Polish border to this map.
2. Read the files `data/wlkp_pn.gpkg` and `data/wlkp_png_dem.tif` into R.
Create a map showing the terrain model for the Wielkopolski National Park area.
Save the obtained map to the file "WLKP_YOURNAME.png".
3. Save the creates map to the file "WLKP_YOURNAME.html". 
Open create file. 
How is it different from the ".png" file?

---
class: inverse, left, bottom
# Data manipulation

---
# Data manipulation

```r
library(dplyr)
```

- The **dplyr** package allows for easy processing of tabular data.
- Its capabilities include selecting specific columns or rows, adding new columns, or performing summaries.
- More information can be found at https://dplyr.tidyverse.org/ and https://r4ds.had.co.nz/transform.html

---
# Attributes manipulation

```r
meteo_data = read.csv("data/poland_meteo_2017.csv", encoding = "UTF-8")
head(meteo_data)
```

<table>
 <thead>
  <tr>
   <th style="text-align:right;"> Kod.stacji </th>
   <th style="text-align:left;"> Nazwa.stacji </th>
   <th style="text-align:right;"> Rok </th>
   <th style="text-align:right;"> Miesiac </th>
   <th style="text-align:right;"> Dzien </th>
   <th style="text-align:right;"> tavg </th>
   <th style="text-align:right;"> precip </th>
   <th style="text-align:right;"> vapor_pres </th>
   <th style="text-align:right;"> humidity </th>
   <th style="text-align:right;"> pressure </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:right;"> 354150100 </td>
   <td style="text-align:left;"> KOŁOBRZEG </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 3.6 </td>
   <td style="text-align:right;"> 4.6 </td>
   <td style="text-align:right;"> 7.4 </td>
   <td style="text-align:right;"> 92.5 </td>
   <td style="text-align:right;"> 1013.5 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354150100 </td>
   <td style="text-align:left;"> KOŁOBRZEG </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 3.5 </td>
   <td style="text-align:right;"> 4.7 </td>
   <td style="text-align:right;"> 6.0 </td>
   <td style="text-align:right;"> 76.9 </td>
   <td style="text-align:right;"> 1011.1 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354150100 </td>
   <td style="text-align:left;"> KOŁOBRZEG </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 3 </td>
   <td style="text-align:right;"> 1.9 </td>
   <td style="text-align:right;"> 10.0 </td>
   <td style="text-align:right;"> 6.7 </td>
   <td style="text-align:right;"> 92.3 </td>
   <td style="text-align:right;"> 1007.1 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354150100 </td>
   <td style="text-align:left;"> KOŁOBRZEG </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 4 </td>
   <td style="text-align:right;"> 2.7 </td>
   <td style="text-align:right;"> 2.3 </td>
   <td style="text-align:right;"> 5.8 </td>
   <td style="text-align:right;"> 76.0 </td>
   <td style="text-align:right;"> 993.9 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354150100 </td>
   <td style="text-align:left;"> KOŁOBRZEG </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 5 </td>
   <td style="text-align:right;"> -1.7 </td>
   <td style="text-align:right;"> 6.2 </td>
   <td style="text-align:right;"> 3.9 </td>
   <td style="text-align:right;"> 74.5 </td>
   <td style="text-align:right;"> 1024.3 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354150100 </td>
   <td style="text-align:left;"> KOŁOBRZEG </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 6 </td>
   <td style="text-align:right;"> -4.5 </td>
   <td style="text-align:right;"> 4.0 </td>
   <td style="text-align:right;"> 3.5 </td>
   <td style="text-align:right;"> 77.3 </td>
   <td style="text-align:right;"> 1038.7 </td>
  </tr>
</tbody>
</table>

---
# Attributes manipulation

```r
meteo_stations = st_read("data/poland_stacje.gpkg")
```

```
## Reading layer `polska_stacje' from data source `/home/jn/Science/conferences/whyr_19/workshop/data/poland_stacje.gpkg' using driver `GPKG'
## Simple feature collection with 1998 features and 2 fields
## geometry type:  POINT
## dimension:      XY
## bbox:           xmin: 172691.6 ymin: 139488.3 xmax: 854361.7 ymax: 846760.3
## epsg (SRID):    NA
## proj4string:    +proj=tmerc +lat_0=0 +lon_0=19 +k=0.9993 +x_0=500000 +y_0=-5300000 +ellps=GRS80 +units=m +no_defs
```

```r
plot(st_geometry(meteo_stations))
```

---
# Attributes manipulation

```r
meteo_data_sel = meteo_data %>% 
  filter(Nazwa.stacji %in% unique(meteo_stations$NAZWA_ST))
```

<table>
 <thead>
  <tr>
   <th style="text-align:right;"> Kod.stacji </th>
   <th style="text-align:left;"> Nazwa.stacji </th>
   <th style="text-align:right;"> Rok </th>
   <th style="text-align:right;"> Miesiac </th>
   <th style="text-align:right;"> Dzien </th>
   <th style="text-align:right;"> tavg </th>
   <th style="text-align:right;"> precip </th>
   <th style="text-align:right;"> vapor_pres </th>
   <th style="text-align:right;"> humidity </th>
   <th style="text-align:right;"> pressure </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 3.4 </td>
   <td style="text-align:right;"> 7.0 </td>
   <td style="text-align:right;"> 7.5 </td>
   <td style="text-align:right;"> 95.4 </td>
   <td style="text-align:right;"> 1013.5 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 2.3 </td>
   <td style="text-align:right;"> 0.9 </td>
   <td style="text-align:right;"> 5.8 </td>
   <td style="text-align:right;"> 80.6 </td>
   <td style="text-align:right;"> 1010.7 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 3 </td>
   <td style="text-align:right;"> 1.4 </td>
   <td style="text-align:right;"> 8.0 </td>
   <td style="text-align:right;"> 6.4 </td>
   <td style="text-align:right;"> 93.0 </td>
   <td style="text-align:right;"> 1007.0 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 4 </td>
   <td style="text-align:right;"> 2.3 </td>
   <td style="text-align:right;"> 6.9 </td>
   <td style="text-align:right;"> 5.6 </td>
   <td style="text-align:right;"> 77.0 </td>
   <td style="text-align:right;"> 992.7 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 5 </td>
   <td style="text-align:right;"> -2.9 </td>
   <td style="text-align:right;"> 1.5 </td>
   <td style="text-align:right;"> 2.8 </td>
   <td style="text-align:right;"> 55.6 </td>
   <td style="text-align:right;"> 1023.6 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 6 </td>
   <td style="text-align:right;"> -5.7 </td>
   <td style="text-align:right;"> 0.9 </td>
   <td style="text-align:right;"> 3.2 </td>
   <td style="text-align:right;"> 79.5 </td>
   <td style="text-align:right;"> 1038.8 </td>
  </tr>
</tbody>
</table>

---
# Attributes manipulation

```r
meteo_data_sel = meteo_data %>% 
  filter(Nazwa.stacji %in% unique(meteo_stations$NAZWA_ST)) %>% 
  mutate(data = as.Date(paste0(Rok, "-", Miesiac, "-", Dzien)))
```

<table>
 <thead>
  <tr>
   <th style="text-align:right;"> Kod.stacji </th>
   <th style="text-align:left;"> Nazwa.stacji </th>
   <th style="text-align:right;"> Rok </th>
   <th style="text-align:right;"> Miesiac </th>
   <th style="text-align:right;"> Dzien </th>
   <th style="text-align:right;"> tavg </th>
   <th style="text-align:right;"> precip </th>
   <th style="text-align:right;"> vapor_pres </th>
   <th style="text-align:right;"> humidity </th>
   <th style="text-align:right;"> pressure </th>
   <th style="text-align:left;"> data </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 3.4 </td>
   <td style="text-align:right;"> 7.0 </td>
   <td style="text-align:right;"> 7.5 </td>
   <td style="text-align:right;"> 95.4 </td>
   <td style="text-align:right;"> 1013.5 </td>
   <td style="text-align:left;"> 2017-01-01 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 2.3 </td>
   <td style="text-align:right;"> 0.9 </td>
   <td style="text-align:right;"> 5.8 </td>
   <td style="text-align:right;"> 80.6 </td>
   <td style="text-align:right;"> 1010.7 </td>
   <td style="text-align:left;"> 2017-01-02 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 3 </td>
   <td style="text-align:right;"> 1.4 </td>
   <td style="text-align:right;"> 8.0 </td>
   <td style="text-align:right;"> 6.4 </td>
   <td style="text-align:right;"> 93.0 </td>
   <td style="text-align:right;"> 1007.0 </td>
   <td style="text-align:left;"> 2017-01-03 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 4 </td>
   <td style="text-align:right;"> 2.3 </td>
   <td style="text-align:right;"> 6.9 </td>
   <td style="text-align:right;"> 5.6 </td>
   <td style="text-align:right;"> 77.0 </td>
   <td style="text-align:right;"> 992.7 </td>
   <td style="text-align:left;"> 2017-01-04 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 5 </td>
   <td style="text-align:right;"> -2.9 </td>
   <td style="text-align:right;"> 1.5 </td>
   <td style="text-align:right;"> 2.8 </td>
   <td style="text-align:right;"> 55.6 </td>
   <td style="text-align:right;"> 1023.6 </td>
   <td style="text-align:left;"> 2017-01-05 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2017 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 6 </td>
   <td style="text-align:right;"> -5.7 </td>
   <td style="text-align:right;"> 0.9 </td>
   <td style="text-align:right;"> 3.2 </td>
   <td style="text-align:right;"> 79.5 </td>
   <td style="text-align:right;"> 1038.8 </td>
   <td style="text-align:left;"> 2017-01-06 </td>
  </tr>
</tbody>
</table>

---
# Attributes manipulation

```r
meteo_data_sel = meteo_data %>% 
  filter(Nazwa.stacji %in% unique(meteo_stations$NAZWA_ST)) %>% 
  mutate(data = as.Date(paste0(Rok, "-", Miesiac, "-", Dzien))) %>% 
  dplyr::select(-Rok, -Miesiac, -Dzien)
```

<table>
 <thead>
  <tr>
   <th style="text-align:right;"> Kod.stacji </th>
   <th style="text-align:left;"> Nazwa.stacji </th>
   <th style="text-align:right;"> tavg </th>
   <th style="text-align:right;"> precip </th>
   <th style="text-align:right;"> vapor_pres </th>
   <th style="text-align:right;"> humidity </th>
   <th style="text-align:right;"> pressure </th>
   <th style="text-align:left;"> data </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 3.4 </td>
   <td style="text-align:right;"> 7.0 </td>
   <td style="text-align:right;"> 7.5 </td>
   <td style="text-align:right;"> 95.4 </td>
   <td style="text-align:right;"> 1013.5 </td>
   <td style="text-align:left;"> 2017-01-01 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2.3 </td>
   <td style="text-align:right;"> 0.9 </td>
   <td style="text-align:right;"> 5.8 </td>
   <td style="text-align:right;"> 80.6 </td>
   <td style="text-align:right;"> 1010.7 </td>
   <td style="text-align:left;"> 2017-01-02 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 1.4 </td>
   <td style="text-align:right;"> 8.0 </td>
   <td style="text-align:right;"> 6.4 </td>
   <td style="text-align:right;"> 93.0 </td>
   <td style="text-align:right;"> 1007.0 </td>
   <td style="text-align:left;"> 2017-01-03 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> 2.3 </td>
   <td style="text-align:right;"> 6.9 </td>
   <td style="text-align:right;"> 5.6 </td>
   <td style="text-align:right;"> 77.0 </td>
   <td style="text-align:right;"> 992.7 </td>
   <td style="text-align:left;"> 2017-01-04 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> -2.9 </td>
   <td style="text-align:right;"> 1.5 </td>
   <td style="text-align:right;"> 2.8 </td>
   <td style="text-align:right;"> 55.6 </td>
   <td style="text-align:right;"> 1023.6 </td>
   <td style="text-align:left;"> 2017-01-05 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 354160105 </td>
   <td style="text-align:left;"> KOSZALIN </td>
   <td style="text-align:right;"> -5.7 </td>
   <td style="text-align:right;"> 0.9 </td>
   <td style="text-align:right;"> 3.2 </td>
   <td style="text-align:right;"> 79.5 </td>
   <td style="text-align:right;"> 1038.8 </td>
   <td style="text-align:left;"> 2017-01-06 </td>
  </tr>
</tbody>
</table>

---
# Attributes manipulation

<table>
 <thead>
  <tr>
   <th style="text-align:right;"> tavg </th>
   <th style="text-align:right;"> pressure </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:right;"> 8.685298 </td>
   <td style="text-align:right;"> 960.5246 </td>
  </tr>
</tbody>
</table>

---
# Attributes manipulation

<table>
 <thead>
  <tr>
   <th style="text-align:right;"> Kod.stacji </th>
   <th style="text-align:right;"> tavg </th>
   <th style="text-align:right;"> pressure </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:right;"> 349190600 </td>
   <td style="text-align:right;"> 9.2780822 </td>
   <td style="text-align:right;"> 1017.513 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 349190625 </td>
   <td style="text-align:right;"> 6.5290411 </td>
   <td style="text-align:right;"> 0.000 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 349190650 </td>
   <td style="text-align:right;"> 0.1978082 </td>
   <td style="text-align:right;"> 0.000 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 349200660 </td>
   <td style="text-align:right;"> 9.2180822 </td>
   <td style="text-align:right;"> 1014.784 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 349210670 </td>
   <td style="text-align:right;"> 8.7728767 </td>
   <td style="text-align:right;"> 1017.553 </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 349220690 </td>
   <td style="text-align:right;"> 8.3227397 </td>
   <td style="text-align:right;"> 1017.737 </td>
  </tr>
</tbody>
</table>

---
# Attribute joining

```r
meteo = meteo_stations %>% 
  left_join(meteo_data_sel, by = c("KOD_SZS" = "Kod.stacji"))
```

<table>
 <thead>
  <tr>
   <th style="text-align:right;"> KOD_SZS </th>
   <th style="text-align:left;"> NAZWA_ST </th>
   <th style="text-align:right;"> tavg </th>
   <th style="text-align:right;"> pressure </th>
   <th style="text-align:left;"> geom </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:right;"> 249190020 </td>
   <td style="text-align:left;"> BRZEŹNICA </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:left;"> c(545976.260510018, 233149.457458902) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 249199999 </td>
   <td style="text-align:left;"> OSIEK </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:left;"> c(519109.405047427, 231411.908795745) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 249190040 </td>
   <td style="text-align:left;"> GIERAŁTOWICE </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:left;"> c(527197.290270776, 230151.641154815) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 249190050 </td>
   <td style="text-align:left;"> RADZISZÓW </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:left;"> c(558306.218144104, 230885.909302064) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 249190060 </td>
   <td style="text-align:left;"> SIEPRAW </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:left;"> c(569259.6149073, 226821.435334128) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 249190070 </td>
   <td style="text-align:left;"> KĘTY </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:right;"> NA </td>
   <td style="text-align:left;"> c(516182.656254972, 223533.556193128) </td>
  </tr>
</tbody>
</table>

---
# Attribute joining

```r
meteo = meteo_stations %>% 
  inner_join(meteo_data_sel, by = c("KOD_SZS" = "Kod.stacji"))
```

<table>
 <thead>
  <tr>
   <th style="text-align:right;"> KOD_SZS </th>
   <th style="text-align:left;"> NAZWA_ST </th>
   <th style="text-align:right;"> tavg </th>
   <th style="text-align:right;"> pressure </th>
   <th style="text-align:left;"> geom </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:right;"> 349190625 </td>
   <td style="text-align:left;"> ZAKOPANE </td>
   <td style="text-align:right;"> 6.5290411 </td>
   <td style="text-align:right;"> 0.000 </td>
   <td style="text-align:left;"> c(569804.515915351, 158924.036586795) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 349190650 </td>
   <td style="text-align:left;"> KASPROWY WIERCH </td>
   <td style="text-align:right;"> 0.1978082 </td>
   <td style="text-align:right;"> 0.000 </td>
   <td style="text-align:left;"> c(571455.865783155, 152132.170449405) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 349200660 </td>
   <td style="text-align:left;"> NOWY SĄCZ </td>
   <td style="text-align:right;"> 9.2180822 </td>
   <td style="text-align:right;"> 1014.784 </td>
   <td style="text-align:left;"> c(621910.953224735, 196894.713701549) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 349210670 </td>
   <td style="text-align:left;"> KROSNO </td>
   <td style="text-align:right;"> 8.7728767 </td>
   <td style="text-align:right;"> 1017.553 </td>
   <td style="text-align:left;"> c(699590.639267207, 208053.055658539) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 349220690 </td>
   <td style="text-align:left;"> LESKO </td>
   <td style="text-align:right;"> 8.3227397 </td>
   <td style="text-align:right;"> 1017.737 </td>
   <td style="text-align:left;"> c(742035.646054704, 183036.071495296) </td>
  </tr>
  <tr>
   <td style="text-align:right;"> 350200570 </td>
   <td style="text-align:left;"> KIELCE-SUKÓW </td>
   <td style="text-align:right;"> 8.4849315 </td>
   <td style="text-align:right;"> 1016.895 </td>
   <td style="text-align:left;"> c(619185.547042676, 328424.723335625) </td>
  </tr>
</tbody>
</table>

---
# Reprojections

.pull-left[

```r
plot(dem)
```

<img src="index_files/figure-html/unnamed-chunk-38-1.png" style="display: block; margin: auto;" />
]

.pull-right[

```r
plot(st_geometry(meteo), axes = TRUE)
```

<img src="index_files/figure-html/unnamed-chunk-39-1.png" style="display: block; margin: auto;" />
]

---
# Reprojections

```r
meteo = meteo_stations %>% 
  inner_join(meteo_data_sel, by = c("KOD_SZS" = "Kod.stacji")) %>% 
  st_transform(4326)
```

---
# Reprojections

```r
tm_shape(dem) +
  tm_raster() +
  tm_shape(meteo) +
  tm_symbols(col = "tavg", palette = "RdBu")
```

---
# Vector-raster interactions

```r
meteo$elev = extract(dem, meteo)
head(meteo)
```

```
## Simple feature collection with 6 features and 5 fields
## geometry type:  POINT
## dimension:      XY
## bbox:           xmin: 19.96032 ymin: 49.23253 xmax: 22.3417 ymax: 50.81048
## epsg (SRID):    4326
## proj4string:    +proj=longlat +datum=WGS84 +no_defs
##     KOD_SZS        NAZWA_ST      tavg pressure                      geom
## 1 349190625        ZAKOPANE 6.5290411    0.000 POINT (19.96032 49.29382)
## 2 349190650 KASPROWY WIERCH 0.1978082    0.000 POINT (19.98182 49.23253)
## 3 349200660       NOWY SĄCZ 9.2180822 1014.784  POINT (20.6886 49.62714)
## 4 349210670          KROSNO 8.7728767 1017.553 POINT (21.76917 49.70674)
## 5 349220690           LESKO 8.3227397 1017.737  POINT (22.3417 49.46647)
## 6 350200570    KIELCE-SUKÓW 8.4849315 1016.895 POINT (20.69221 50.81048)
##   elev
## 1  844
## 2 1796
## 3  275
## 4  323
## 5  384
## 6  247
```

---
# Vector-raster interactions

.pull-left[

```r
lc = raster("data/poland_lc.tif")
plot(lc)
```

<img src="index_files/figure-html/unnamed-chunk-43-1.png" style="display: block; margin: auto;" />
]

.pull-right[

```r
wlkp_pn = st_read("data/wlkp_pn.gpkg", quiet = TRUE)
plot(st_geometry(wlkp_pn), axes = TRUE)
```

<img src="index_files/figure-html/unnamed-chunk-44-1.png" style="display: block; margin: auto;" />
]

---
# Vector-raster interactions

```r
wlkp_pn_lc = crop(lc, wlkp_pn)
plot(wlkp_pn_lc)
```

---
# Vector-raster interactions

```r
wlkp_pn_lc2 = mask(wlkp_pn_lc, wlkp_pn)
plot(wlkp_pn_lc2)
```

---
# Exercises

1. Read the 2017 meteorological data and the location of the meteorological stations.
Create a spatial object containing all measurements for three stations: `BIAŁYSTOK`, `KRAKÓW-BALICE`, and `POZNAŃ`. (Tip: use the `right_join()` function.
2. Calculate the average temperature in January and July for these three stations.
3. Determine the land cover of the three stations.
4. Additional: improve the average temperature map presented in the "Reprojections" section.

---
class: inverse, left, bottom
# Resources

---
# Resources

.pull-left[
- http://geocompr.github.io/ - Lovelace R., Nowosad J., Muenchow J. 2019, **Geocomputation with R**. CRC Press
- **sf** - https://r-spatial.github.io/sf/
- **raster** - https://rspatial.org/spatial/index.html
- **tmap** - https://github.com/mtennekes/tmap
- **dplyr** - http://r4ds.had.co.nz/ - Wickham H., Goremund G. 2016. **R for Data Science**. O'Reilly Media
]

.pull-right[
- Blogs (e.g., https://www.r-bloggers.com/, https://rweekly.org/)
- Websites (e.g., https://gis.stackexchange.com/, https://stat.ethz.ch/mailman/listinfo/r-sig-geo)
- Twitter `#rspatial`, `#geocompr`
- Meetups
]

---
class: center, middle, clear

.pull-right[
## About me:

Twitter: <i class="ion  ion-social-twitter "></i> jakub_nowosad

Email: nowosad.jakub@gmail.com

https://nowosad.github.io

https://geocompr.github.io/
]