Assignment Instructions 📖

Create a map of rainfall in Ireland for the 25 weather stations, colour coding the symbol for each station according to its median rainfall level in January.


1. Introduction 🖥️

This blog will look at how to create a rainfall map in Ireland using the weather_stations.geojson data set which covers the period 1850-2014. The process entails the development of an interactive map using 25 weather stations around Ireland and colour coding them based on their median rainfall during the month of January.

Data sets used

Alongside the weather_stations.geojson data set, the counties.geojson data set is also used to help build the interactive map and its subsequent layers.


2. Methodology 📊

How the Map was made 🗺️

We will start by loading in the required packages used in the script below.

library(tidyverse) # sf package uses this
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.1     ✔ stringr   1.5.2
## ✔ ggplot2   4.0.0     ✔ tibble    3.3.0
## ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
## ✔ purrr     1.1.0     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(sf) # The 'sf' package handles geographical data
## Linking to GEOS 3.13.0, GDAL 3.8.5, PROJ 9.5.1; sf_use_s2() is TRUE
library(tmap) 
library(leaflet)
library(leafpop)

Loading data sets

#read in rainfall 
rainfall <- st_read('weather_stations.geojson')  %>% st_transform(4326)
## Reading layer `weather_stations' from data source 
##   `/Users/willryanuni/Library/CloudStorage/OneDrive-MaynoothUniversity/GY672[A] — Analysing Spatial and Temporal Data Using R/Class Materials-20250924/weather_stations.geojson' 
##   using driver `GeoJSON'
## Simple feature collection with 49500 features and 11 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: 46633.44 ymin: 59756.08 xmax: 330201.3 ymax: 457227.9
## Projected CRS: TM65 / Irish Grid
rainfall # looking at our rainfall data using geojson
## Simple feature collection with 49500 features and 11 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: -10.23045 ymin: 51.79042 xmax: -5.991011 ymax: 55.35998
## Geodetic CRS:  WGS 84
## First 10 features:
##    Year Month Rainfall Station Elevation  Easting Northing  County Abbreviation
## 1  1850   Jan    169.0  Ardara        15 180787.7   394679 Donegal           AR
## 2  1851   Jan    236.4  Ardara        15 180787.7   394679 Donegal           AR
## 3  1852   Jan    249.7  Ardara        15 180787.7   394679 Donegal           AR
## 4  1853   Jan    209.1  Ardara        15 180787.7   394679 Donegal           AR
## 5  1854   Jan    188.5  Ardara        15 180787.7   394679 Donegal           AR
## 6  1855   Jan     32.3  Ardara        15 180787.7   394679 Donegal           AR
## 7  1856   Jan    151.6  Ardara        15 180787.7   394679 Donegal           AR
## 8  1857   Jan    179.1  Ardara        15 180787.7   394679 Donegal           AR
## 9  1858   Jan    110.1  Ardara        15 180787.7   394679 Donegal           AR
## 10 1859   Jan    157.8  Ardara        15 180787.7   394679 Donegal           AR
##    Source coast_dist                   geometry
## 1  Briffa   17.36709 POINT (-8.290712 54.79005)
## 2  Briffa   17.36709 POINT (-8.290712 54.79005)
## 3  Briffa   17.36709 POINT (-8.290712 54.79005)
## 4  Briffa   17.36709 POINT (-8.290712 54.79005)
## 5  Briffa   17.36709 POINT (-8.290712 54.79005)
## 6  Briffa   17.36709 POINT (-8.290712 54.79005)
## 7  Briffa   17.36709 POINT (-8.290712 54.79005)
## 8  Briffa   17.36709 POINT (-8.290712 54.79005)
## 9  Briffa   17.36709 POINT (-8.290712 54.79005)
## 10 Briffa   17.36709 POINT (-8.290712 54.79005)
#read in counties 
counties <- st_read('counties.geojson') %>% st_transform(4326)
## Reading layer `counties' from data source 
##   `/Users/willryanuni/Library/CloudStorage/OneDrive-MaynoothUniversity/GY672[A] — Analysing Spatial and Temporal Data Using R/Class Materials-20250924/counties.geojson' 
##   using driver `GeoJSON'
## Simple feature collection with 40 features and 2 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: 17528.59 ymin: 19537.25 xmax: 366403.6 ymax: 466923.2
## Projected CRS: TM65 / Irish Grid
counties # looking at our rainfall data using geojson
## Simple feature collection with 40 features and 2 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -10.66221 ymin: 51.41943 xmax: -5.432885 ymax: 55.44661
## Geodetic CRS:  WGS 84
## First 10 features:
##            County mean_rain                       geometry
## 1          Antrim  77.71545 MULTIPOLYGON (((-6.054397 5...
## 2          Armagh  77.71545 MULTIPOLYGON (((-6.353978 5...
## 3   Carlow County  92.14104 MULTIPOLYGON (((-6.708988 5...
## 4    Cavan County  72.24133 MULTIPOLYGON (((-7.923589 5...
## 5    Clare County 100.73245 MULTIPOLYGON (((-9.090629 5...
## 6       Cork City 100.00000 MULTIPOLYGON (((-8.476869 5...
## 7     Cork County 107.23233 MULTIPOLYGON (((-10.2472 51...
## 8  Donegal County 104.23685 MULTIPOLYGON (((-8.832665 5...
## 9            Down  77.71545 MULTIPOLYGON (((-6.079602 5...
## 10    Dublin City  62.23451 MULTIPOLYGON (((-6.224726 5...

The st_transform function allows us to convert the geometry points in the data sets to latitude and longitude co-coordinates. Which will help us manipulate and map the data later.

Filtering and calculating the data

Median rainfall for January in each station was filtered and calculated from the data in the period 1850-2014.

rainfall %>% filter(Month == 'Jan') %>% group_by(Station, County) %>%
  summarise(median_rain_jan = median(Rainfall, na.rm = TRUE)) -> median_jan_rainfall
## `summarise()` has grouped output by 'Station'. You can override using the
## `.groups` argument.

After this process of data manipulation, visually checking the new median_jan_rainfall data set is a good way to identify any irregularities or errors that may have occurred in the previous step.

head(median_jan_rainfall)
## Simple feature collection with 6 features and 3 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: -8.290712 ymin: 52.19037 xmax: -5.991011 ymax: 54.79005
## Geodetic CRS:  WGS 84
## # A tibble: 6 × 4
## # Groups:   Station [6]
##   Station    County    median_rain_jan             geometry
##   <chr>      <chr>               <dbl>          <POINT [°]>
## 1 Ardara     Donegal             172.  (-8.290712 54.79005)
## 2 Armagh     Armagh               75   (-6.640924 54.35011)
## 3 Athboy     Meath                87.1  (-6.930877 53.6002)
## 4 Belfast    Antrim              102.   (-5.991011 54.5001)
## 5 Birr       Offaly               77.5  (-7.88075 53.08026)
## 6 Cappoquinn Waterford           147.  (-7.800753 52.19037)
tail(median_jan_rainfall)
## Simple feature collection with 6 features and 3 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: -10.23045 ymin: 51.79042 xmax: -7.090842 ymax: 53.75018
## Geodetic CRS:  WGS 84
## # A tibble: 6 × 4
## # Groups:   Station [6]
##   Station                   County    median_rain_jan             geometry
##   <chr>                     <chr>               <dbl>          <POINT [°]>
## 1 Roches Point              Cork                111   (-8.240696 51.79042)
## 2 Shannon Airport           Clare                92.9 (-8.920616 52.68031)
## 3 Strokestown               Roscommon           102.  (-8.100728 53.75018)
## 4 University College Galway Galway              124.  (-9.060602 53.27024)
## 5 Valentia                  Kerry               166    (-10.23045 51.9304)
## 6 Waterford                 Waterford           103.  (-7.090842 52.29037)

Create a filename for each station

Following this, a file name was created for each station. This will later be added to median_jan_rainfall as this will be used for the embedded popups as part of the finalised interactive map object. To prevent errors in joining the data sets, we will ensure that the sf object geometry is dropped, and only include Station and Filename variables in the file object. This prevents issues that arise when joining two data sets which have sf objects in them, and also where variable names can become duplicated.

# make a filename for each station
median_jan_rainfall %>% st_drop_geometry() %>% mutate(Filename=file.path(getwd(),paste0(Station,'.png'))) %>% select(Station, Filename) -> files
files %>% head
## # A tibble: 6 × 2
## # Groups:   Station [6]
##   Station    Filename                                                           
##   <chr>      <chr>                                                              
## 1 Ardara     /Users/willryanuni/Library/CloudStorage/OneDrive-MaynoothUniversit…
## 2 Armagh     /Users/willryanuni/Library/CloudStorage/OneDrive-MaynoothUniversit…
## 3 Athboy     /Users/willryanuni/Library/CloudStorage/OneDrive-MaynoothUniversit…
## 4 Belfast    /Users/willryanuni/Library/CloudStorage/OneDrive-MaynoothUniversit…
## 5 Birr       /Users/willryanuni/Library/CloudStorage/OneDrive-MaynoothUniversit…
## 6 Cappoquinn /Users/willryanuni/Library/CloudStorage/OneDrive-MaynoothUniversit…

Construct the interactive map base

tmap_mode('view')
## ℹ tmap modes "plot" - "view"
## ℹ toggle with `tmap::ttm()`
tm_start_map <- tm_basemap(leaflet::providers$CartoDB.Positron) + tm_shape(counties) + tm_borders(lwd = 0.5)
tm_start_map

Produce the plots for each station

In the loop below, bar chart plots for each station were created so that they can be assigned to the Filenames that were created above (In the Section: Create a filename for each station). It is designed so that as the loop is going through each row in the files object, it will identify the current station and the median rainfall value associated with it.

This allows us to assign a colour that highlights the specific station in the graph for all 25 stations. Enabling us to visually see in the bar chart how the particular station compares to others. As we are using a median rainfall value for each station, and not values across a time series, it is important to find clear way which show the single variable, and how it changes spatially.

# Loop through each station to create a bar chart highlighting the specific weather station

for (i in 1:nrow(files)) {
  
  current_station <- files$Station[i]
  
  # Get rainfall value
  current_rainfall <- median_jan_rainfall %>%
    filter(Station == current_station) %>%
    pull(median_rain_jan)
  
  # Assign fill colour for highlighting
  median_jan_rainfall$fill_colour <- ifelse(
    median_jan_rainfall$Station == current_station,
    "#006778",
    "lightgrey"
  )
  
  # Open PNG device
  png(files$Filename[i], width = 500, height = 300)
  
  # Plot
  station_plot <- ggplot(median_jan_rainfall, aes(x = Station, y = median_rain_jan, fill = fill_colour)) +
    geom_bar(stat = "identity") +
    scale_fill_identity() +
    theme_minimal() +
    coord_cartesian(ylim = c(0, 200)) +
    labs(
      title = paste("Median January Rainfall for", current_station,
                    paste0("(", round(current_rainfall, 2), " mm)")),
      x = "Station",
      y = "Median Rainfall (mm)"
    ) +
    theme(
      plot.title = element_text(size = 13, face = "bold"), 
      axis.title.x = element_text(size = 13, face = "bold"),
      axis.title.y = element_text(size = 13, face = "bold"),
      axis.text.x = element_text(angle = 45,size = 10, face = "bold", hjust = 1),
      axis.text.y = element_text(size = 10, face = "bold")
      )
  
  # Print and save
  print(station_plot)
  dev.off()
  
}

Join Filenames with the January median rainfall variables in median_jan_rainfall

median_jan_rainfall %>% left_join(files, by = "Station") %>%
  select(Station, median_rain_jan, Filename) -> Station_jan_med_rainfall

Use leaflet to convert tmap objects and add Popups

tmap_leaflet(tm_start_map) %>% addCircleMarkers(data = Station_jan_med_rainfall, 
                                                popup=popupImage(Station_jan_med_rainfall$Filename,embed=TRUE)) 

3. Bringing it all together 📈

Now that the Filename, interactive map object and Popups have been created, we can bring it all together in the final interactive map of this blog. Note that tmap has been overhauled in its newest version (v4), and uses updated syntax compared to v3. Two useful resources that were invaluable in overcoming this learning curve were the tmap website on github and the book Spatial Data Visualization with tmap. Both are excellent resources, and highly recommended for working with spatial data in R. I will rely on them in the future for sure!

tm_start2 <- tm_shape(Station_jan_med_rainfall) + 
               tm_dots(
                 fill = "median_rain_jan",
                 popup.vars = FALSE,
                 fill.scale = tm_scale(values = "brewer.blues"), 
                 fill.legend = tm_legend(title = "Median January Rainfall: 1850-2014 (mm) ", reverse = TRUE),
                 size = 0.8
                 ) +
               tm_shape(counties) +
                tm_borders(lwd = 0.5)
               

tmap_leaflet(tm_start2) %>% 
  addCircleMarkers(data=Station_jan_med_rainfall,stroke = NA,
                   popup=popupImage(Station_jan_med_rainfall$Filename,embed=TRUE)) 

4. Spatial patterns Identified

We can see from the interactive map (Section: [3. Bringing it all together 📈]) that coastal weather stations experienced higher median January rainfall levels compared to inland weather stations. The three stations that recorded the highest median levels for this period were Killarney (177.7 mm), Ardara (171.6 mm) and Valentia (166 mm); all of which are located on the west coast. The image below shows the annual rainfall (mm) in Ireland from Ireland’s Climate Averages 1991-2020 report. It reveals that patterns identified of higher precipitation values in coastal regions, and lower values inland are reflected in empirical data. The report also notes that in addition to the coastal and western location, elevation also plays a role in increasing rainfall.


Conclusion

This blog shows the process of how to create an interactive map of median rainfall in Ireland using 25 weather stations. It goes through the steps on how to manipulate data sets to filter and select particular variables, how to create file names, Popup images and embed them successfully into an interactive map object using the updated tmap package and leaflet. Further work to improve the map could be to include how the Elevation and coastal distance (coast_dist) variables in the weather_stations.geojson data set can further explain the spatial patterns identified from this blog! Thanks for reading the first entry for January 2026, maybe the next blog can show how these observations compare to the rainfall records collected during this month?

Happy mapping!