---
title: "Tutorial (vignette) for the eurostat R package"
date: "`r Sys.Date()`"
output: 
  rmarkdown::html_vignette:
    toc: true
vignette: >
  %\VignetteIndexEntry{eurostat tutorial}
  %\VignetteDepends{Cairo}
  %\VignetteEncoding{UTF-8}
  \usepackage[utf8]{inputenc}
  %\VignetteEngine{knitr::rmarkdown}
editor_options: 
  chunk_output_type: console
---

# R Tools for Eurostat Open Data

This [rOpenGov](http://ropengov.github.io) R package provides tools to access [Eurostat database](http://ec.europa.eu/eurostat/data/database), which you can also browse on-line for the data sets and documentation. For contact information and source code, see the [package website](http://ropengov.github.io/eurostat/).


```{r setup, include=FALSE}
# Global options
library(knitr)
opts_chunk$set(fig.path="fig/")
```

# Installation

Release version [(CRAN)](https://cran.r-project.org/web/packages/eurostat/index.html):

```{r install, eval=FALSE}
install.packages("eurostat")
```

Development version [(Github)](https://github.com/rOpenGov/eurostat):

```{r install2, eval=FALSE}
library(devtools)
install_github("ropengov/eurostat")
```

```{r, echo=FALSE}
library(eurostat)
```


Overall, the eurostat package includes the following functions:

```{r, echo=FALSE,comment=NA}
cat(paste0(library(help = "eurostat")$info[[2]], collapse = "\n"))
```


# Finding data

Function `get_eurostat_toc()` downloads a table of contents of eurostat datasets. The values in column 'code' should be used to download a selected dataset.

```{r get_eurostat_toc, warning=FALSE, message=FALSE}
# Load the package
library(eurostat)
library(rvest)

# Get Eurostat data listing
toc <- get_eurostat_toc()

# Check the first items
library(knitr)
kable(head(toc))
```

With `search_eurostat()` you can search the table of contents for particular patterns, e.g. all datasets related to *passenger transport*. The kable function to produces nice markdown output. Note that with the `type` argument of this function you could restrict the search to for instance datasets or tables.

```{r search_eurostat, warning=FALSE, message=FALSE}
# info about passengers
kable(head(search_eurostat("passenger transport")))
```

Codes for the dataset can be searched also from the [Eurostat
database](http://ec.europa.eu/eurostat/data/database). The Eurostat
database gives codes in the Data Navigation Tree after every dataset
in parenthesis.

# Downloading data 

The package supports two of the Eurostats download methods: the bulk download 
facility and the Web Services' JSON API. The bulk download facility is the 
fastest method to download whole datasets. It is also often the only way as 
the JSON API has limitation of maximum 50 sub-indicators at a time and 
whole datasets usually exceeds that. To download only a small section of the 
dataset the JSON API is faster, as it allows to make a data selection 
before downloading.

A user does not usually have to bother with methods, as both are used via main
function `get_eurostat()`. If only the table id is given, the whole table is 
downloaded from the bulk download facility. If also filters are defined 
the JSON API is used.

Here an example of indicator 'Modal split of passenger transport'. This is the percentage share of each mode of transport in total inland transport, expressed in passenger-kilometres (pkm) based on transport by passenger cars, buses and coaches, and trains. All data should be based on movements on national territory, regardless of the nationality of the vehicle. However, the data collection is not harmonized at the EU level. 

Pick and print the id of the data set to download: 
```{r get_id, warning=FALSE, message=FALSE, results='asis'}
# For the original data, see
# http://ec.europa.eu/eurostat/tgm/table.do?tab=table&init=1&plugin=1&language=en&pcode=tsdtr210
id <- search_eurostat("Modal split of passenger transport", 
        	             type = "table")$code[1]
print(id)
```

Get the whole corresponding table. As the table is annual data, it is more
convient to use a numeric time variable than use the default date format:

```{r get_eurostat, warning=FALSE, message=FALSE, results='asis'}
dat <- get_eurostat(id, time_format = "num")
```

Investigate the structure of the downloaded data set:
```{r str_dat, warning=FALSE, message=FALSE}
str(dat)
```

```{r head_dat, warning=FALSE, message=FALSE, results='asis'}
kable(head(dat))
```

Or you can get only a part of the dataset by defining `filters` argument. It
should be named list, where names corresponds to variable names (lower case) and
values are vectors of codes corresponding desidered series (upper case). For
time variable, in addition to a `time`, also a `sinceTimePeriod` 
and a `lastTimePeriod` can be used.

```{r get_eurostat_json, warning=FALSE, message=FALSE, results='asis', eval=FALSE}
dat2 <- get_eurostat(id, filters = list(geo = c("EU28", "FI"), lastTimePeriod=1), time_format = "num")
kable(dat2)
```

## Replacing codes with labels

By default variables are returned as Eurostat codes, but to get human-readable 
labels instead, use a `type = "label"` argument.  

```{r json_labels, warning=FALSE, message=FALSE, results='asis', eval=FALSE}
datl2 <- get_eurostat(id, filters = list(geo = c("EU28", "FI"), 
                                         lastTimePeriod = 1), 
                      type = "label", time_format = "num")
kable(head(datl2))
```

Eurostat codes in the downloaded data set can be replaced with
human-readable labels from the Eurostat dictionaries with the
`label_eurostat()` function.

```{r labels, warning=FALSE, message=FALSE, results='asis'}
datl <- label_eurostat(dat)
kable(head(datl))
```

The `label_eurostat()` allows conversion of individual variable
vectors or variable names as well.

```{r name_labels, eval = FALSE}
label_eurostat_vars(names(datl))
```


Vehicle information has 3 levels. You can check them now with:

```{r vehicle_levels, eval = FALSE}
levels(datl$vehicle)
```



# Selecting and modifying data

## EFTA, Eurozone, EU and EU candidate countries

To facilitate smooth visualization of standard European geographic areas, the package provides ready-made lists of the country codes used in the eurostat database for EFTA (efta\_countries), Euro area (ea\_countries), EU (eu\_countries) and EU candidate countries (eu\_candidate\_countries). These can be used to select specific groups of countries for closer investigation. For conversions with other standard country coding systems, see the [countrycode](...) R package. To retrieve the country code list for EFTA, for instance, use:

```{r, echo=TRUE, eval=TRUE}
data(efta_countries)
kable(efta_countries)
```


## EU data from 2012 in all vehicles:

```{r eu_12}
dat_eu12 <- subset(datl, geo == "European Union (current composition)" & time == 2012)
kable(dat_eu12, row.names = FALSE)
```

## EU data from 2000 - 2012 with vehicle types as variables:

Reshaping the data is best done with `spread()` in `tidyr`.
```{r eu_vehicles_table}
library("tidyr")
dat_eu_0012 <- subset(dat, geo == "EU28" & time %in% 2000:2012)
dat_eu_0012_wide <- spread(dat_eu_0012, vehicle, values)
kable(subset(dat_eu_0012_wide, select = -geo), row.names = FALSE)
```

## Train passengers for selected EU countries in 2000 - 2012

```{r trains_table}
dat_trains <- subset(datl, geo %in% c("Austria", "Belgium", "Finland", "Sweden")
                     & time %in% 2000:2012 
                     & vehicle == "Trains")

dat_trains_wide <- spread(dat_trains, geo, values) 
kable(subset(dat_trains_wide, select = -vehicle), row.names = FALSE)
```



# Visualization

Visualizing train passenger data with `ggplot2`:

```{r trains_plot, fig.width=6, fig.height=3}
library(ggplot2)
p <- ggplot(dat_trains, aes(x = time, y = values, colour = geo)) 
p <- p + geom_line()
print(p)
```

<a name="triangle"></a>**Triangle plot**

Triangle plot is handy for visualizing data sets with three variables. 

```{r plotGallery, warning=FALSE, message=FALSE, fig.width=6, fig.height=6}
library(tidyr)
library(plotrix)
library(eurostat)
library(dplyr)
library(tidyr)

# Some cleaning of the data is required
 energy3 <- get_eurostat("nrg_114a") %>%
 label_eurostat(dat) %>%
 filter(time == "2013-01-01") %>%
 mutate(geo = gsub(geo, pattern=" \\(.*", replacement="")) %>%
 select(product, geo, values) %>%
 group_by(product, geo) %>%
 summarise(svalue = sum(values)) %>%
 group_by(geo) %>%
 mutate(tvalue = sum(svalue),
 svalue = svalue/sum(svalue)) %>%
 filter(tvalue > 1000) %>% # only large countries
 spread(product, svalue)
 
# Triangle plot
 par(cex=0.75, mar=c(0,0,0,0))
 positions <- plotrix::triax.plot(as.matrix(energy3[, c(3,5,4)]),
                     show.grid = TRUE,
                     label.points= FALSE, point.labels = energy3$geo,
                     col.axis="gray50", col.grid="gray90",
                     pch = 19, cex.axis=0.8, cex.ticks=0.7, col="grey50")

 # Larger labels
 ind <- which(energy3$geo %in%  c("Norway", "Iceland","Denmark","Estonia", "Turkey", "Italy", "Finland"))
 df <- data.frame(positions$xypos, geo = energy3$geo)
 points(df$x[ind], df$y[ind], cex=2, col="red", pch=19)
 text(df$x[ind], df$y[ind], df$geo[ind], adj = c(0.5,-1), cex=1.5)
```



## Maps 

###  Disposable income of private households by NUTS 2 regions at 1:60mln resolution using tmap

The mapping examples below use [`tmap`](https://github.com/mtennekes/tmap) package.


```{r maps1-1, eval=TRUE, fig.width=8, fig.height=8}
library(dplyr)
library(eurostat)
library(sf)
library(tmap)

# Download attribute data from Eurostat
sp_data <- eurostat::get_eurostat("tgs00026",
		                  time_format = "raw",
                                  stringsAsFactors = FALSE) %>% 
  # subset to have only a single row per geo
  dplyr::filter(time == 2010, nchar(geo) == 4) %>% 
  # categorise
  dplyr::mutate(income = cut_to_classes(values, n = 5))

# Download geospatial data from GISCO
geodata <- get_eurostat_geospatial(output_class = "sf",
                                   resolution = "60",
				   nuts_level = 2,
				   year = 2013)

# merge with attribute data with geodata
map_data <- inner_join(geodata, sp_data)
```


Construct the map

```{r map1ex, eval=TRUE, warning=FALSE}
map1 <- tmap::tm_shape(geodata) +
  tmap::tm_fill("lightgrey") +
  tmap::tm_shape(map_data) +
  tmap::tm_grid() +
  tmap::tm_polygons("income", title = "Disposable household\nincomes in 2010",  
                    palette = "Oranges")
print(map1)  
```

Interactive maps can be generated as well

```{r maps1-2, eval=FALSE, fig.width=8, fig.height=8}
# Interactive
tmap_mode("view")
map1

# Set the mode back to normal plotting
tmap_mode("plot")
print(map1)
```

### Disposable income of private households by NUTS 2 regions in Poland with labels at 1:1mln resolution using tmap

```{r maps2, fig.width=8, fig.height=8, warning=FALSE}
library(eurostat)
library(dplyr)
library(sf)
library(RColorBrewer)

# Downloading and manipulating the tabular data
print("Let us focus on year 2014 and NUTS-3 level")
euro_sf2 <- get_eurostat("tgs00026", time_format = "raw",
                         stringsAsFactors = FALSE,
			 filter = list(time = "2014")) %>% 
 
  # Subset to NUTS-3 level
  dplyr::filter(grepl("PL",geo)) %>% 
  # label the single geo column
  mutate(label = paste0(label_eurostat(.)[["geo"]], "\n", values, "€"),
         income = cut_to_classes(values))

print("Download geospatial data from GISCO")
geodata <- get_eurostat_geospatial(output_class = "sf", resolution = "60", nuts_level = 2, year = 2013)

# Merge with attribute data with geodata
map_data <- inner_join(geodata, euro_sf2)

# plot map
map2 <- tm_shape(geodata) +
  tm_fill("lightgrey") +
  tm_shape(map_data, is.master = TRUE) +
  tm_polygons("income", title = "Disposable household incomes in 2014",
              palette = "Oranges", border.col = "white") + 
  tm_text("NUTS_NAME", just = "center") + 
  tm_scale_bar()
map2
```

### Disposable income of private households by NUTS 2 regions at 1:10mln resolution using spplot

```{r maps3, fig.width=8, fig.height=8, dev='CairoPNG'}
library(sp)
library(eurostat)
library(dplyr)
library(RColorBrewer)
dat <- get_eurostat("tgs00026", time_format = "raw", stringsAsFactors = FALSE) %>% 
  # subsetting to year 2014 and NUTS-2 level
  dplyr::filter(time == 2014, nchar(geo) == 4) %>% 
  # classifying the values the variable
  dplyr::mutate(cat = cut_to_classes(values))

# Download geospatial data from GISCO
geodata <- get_eurostat_geospatial(output_class = "spdf", resolution = "10", nuts_level = 2, year = 2013)

# merge with attribute data with geodata
geodata@data <- left_join(geodata@data, dat)

# plot map
sp::spplot(obj = geodata, "cat", main = "Disposable household income",
	   xlim = c(-22,34), ylim = c(35,70), 
           col.regions = c("dim grey", brewer.pal(n = 5, name = "Oranges")),
	   col = "white", usePolypath = FALSE)
```

### Disposable income of private households by NUTS 2 regions at 1:60mln resolution using ggplot2

Meanwhile the CRAN version of `ggplot2` is lacking support for simple features, you can plot maps with `ggplot2` by downloading geospatial data as `data.frame` with `output_class` argument set as `df`.

```{r maps4, fig.width=8, fig.height=8, dev='CairoPNG'}
library(eurostat)
library(dplyr)
library(ggplot2)
dat <- get_eurostat("tgs00026", time_format = "raw", stringsAsFactors = FALSE) %>% 
  # subsetting to year 2014 and NUTS-2 level
  dplyr::filter(time == 2014, nchar(geo) == 4) %>% 
  # classifying the values the variable
  dplyr::mutate(cat = cut_to_classes(values))

# Download geospatial data from GISCO
geodata <- get_eurostat_geospatial(resolution = "60", nuts_level = "2", year = 2013)

# merge with attribute data with geodata
map_data <- inner_join(geodata, dat)

ggplot(data=map_data) + geom_sf(aes(fill=cat),color="dim grey", size=.1) + 
    scale_fill_brewer(palette = "Oranges") +
  guides(fill = guide_legend(reverse=T, title = "euro")) +
  labs(title="Disposable household income in 2014",
       caption="(C) EuroGeographics for the administrative boundaries 
                Map produced in R with a help from Eurostat-package <github.com/ropengov/eurostat/>") +
  theme_light() + theme(legend.position=c(.8,.8)) +
  coord_sf(xlim=c(-12,44), ylim=c(35,70))

```


## SDMX

Eurostat data is available also in the SDMX format. The eurostat R package does not provide custom tools for this but the generic [rsdmx](https://cran.r-project.org/web/packages/rsdmx/index.html) and [rjsdmx](https://github.com/amattioc/SDMX/wiki) R packages can be used to access data in that format when necessary:

```{r rsdmx, fig.width=8, fig.height=8, dev='CairoPNG', eval=FALSE}
library(rsdmx)

# Data set URL
url <- "http://ec.europa.eu/eurostat/SDMX/diss-web/rest/data/cdh_e_fos/..PC.FOS1.BE/?startperiod=2005&endPeriod=2011"

# Read the data from eurostat
d <- readSDMX(url)

# Convert to data frame and show the first entries
df <- as.data.frame(d)

kable(head(df))
```




# Further examples

For further examples, see the [package homepage](http://ropengov.github.io/eurostat/articles/index.html).


# Citations and related work

### Citing the data sources

Eurostat data: cite [Eurostat](http://ec.europa.eu/eurostat/).

Administrative boundaries: cite EuroGeographics


### Citing the eurostat R package

For main developers and contributors, see the [package homepage](http://ropengov.github.io/eurostat).

This work can be freely used, modified and distributed under the
BSD-2-clause (modified FreeBSD) license:

```{r citation, message=FALSE, eval=TRUE, echo=TRUE}
citation("eurostat")
```


### Contact

For contact information, see the [package homepage](http://ropengov.github.io/eurostat).


# Version info

This tutorial was created with

```{r sessioninfo, message=FALSE, warning=FALSE}
sessionInfo()
```
