In this post, I will explain how to create interactive HTML reports including widgets such as maps and interactive graphs using R. I assume that the reader is familiar with basic R. I will use two public APIs from the Environment Agency to download data, the Water Quality Archive for water quality data (phosphate in this example) and the Rainfall API to look for rainfall monitoring locations.
First, let’s load the packages we need and get the current date.
library(readr)
library(sf)
library(dplyr)
library(leaflet)
library(rgdal)
library(tidyr)
library(ggplot2)
library(plotly)
library(DT)
library(htmltools)
###Get Current Date and Time
now <- format(as.Date(Sys.Date()), '%d/%m/%Y')
Then, we can retrieve phosphate information for two sites (Powick and Leintwardine on the River Teme, West Midlands). The information to retrieve determinand and site codes are in the documentation linked above. Once we have retrieved the data, we can rename the columns to something a bit more sensible than the original output and create a new concentration column that converts concentrations below the limit of detection to half of the limit of detection. We also add latitude and longitude columns converting from UK-specific coordinates (eastings and northings).
# Set URL parameters
base_url = "https://environment.data.gov.uk/water-quality/data/measurement.csv?_limit=10000"
sampling_points = "&samplingPoint=MD-13598380&samplingPoint=MD-13624310"
start_date = "&startDate=2019-01-01"
detcode = "&determinand=0118"
# retrieve data from API and rename columns
all_records <- read_csv(paste0(base_url, sampling_points, detcode, start_date))
colnames(all_records) <- c("id", "smpt_url", "smpt_code", "smpt_name", "sample_datetime", "determinand", "determinand_def", "detcode", "less_than", "concentration", "result_interpretation",
"unit", "smpt_type", "is_compliance", "sample_purpose", "easting", "northing")
# Add concentration and latitude/longitude columns
all_records$less_than_bool <- ifelse(is.na(all_records$less_than), FALSE, TRUE)
all_records$concentration_adj <- as.numeric(ifelse(all_records$less_than_bool == "TRUE", all_records$concentration/2, all_records$concentration))
latlon <- all_records %>%
st_as_sf(coords = c("easting", "northing"), crs = 27700) %>%
st_transform(4326) %>%
st_coordinates() %>%
as_tibble()
all_records$lat <- latlon$Y
all_records$lon <- latlon$X
We then create an HTML table with site names, codes and coordinates.
sites <- all_records %>% select(smpt_code, smpt_name, easting,
northing, lat, lon) %>%
distinct(.keep_all = TRUE)
# Sites table
table_df <- sites %>% select(smpt_code, smpt_name, easting, northing) %>%
rename("Site code" = smpt_code,
"Site name"= smpt_name,
"Easting" = easting,
"Northing" = northing)
sites_table <- datatable(table_df, width="100%", height="100%",
caption = htmltools::tags$caption( style = 'caption-side: top;
text-align: center; color:black;
font-size:200% ;','Monitoring sites') )
Then, we search for the closest rainfall monitoring site to each of the water quality sites and put this information in another table.
rainfall_data = data.frame()
for (i in 1:nrow(sites)) {
lat <- sites$lat[i]
lon <- sites$lon[i]
rainfall_st <- read_csv(paste0("https://environment.data.gov.uk/flood-monitoring/id/stations.csv?parameter=rainfall&lat=",lat,"&long=",lon,"&dist=10"))
if (nrow(rainfall_st) > 0) {
rainfall_data <- rbind(rainfall_data, rainfall_st[1,])
}
}
rainfall_data <- as.data.frame(rainfall_data)
# rainfall table
rainfall_table_df <- rainfall_data %>% select(notation, qualifier, easting, northing) %>%
rename("Site code" = notation,
"Site type"= qualifier,
"Easting" = easting,
"Northing" = northing)
rainfall_table <- datatable(rainfall_table_df, width="100%", height="100%",
caption = htmltools::tags$caption( style = 'caption-side: top;
text-align: center; color:black;
font-size:200% ;','Rainfall gauges') )
At this point, we can create a map showing water quality and corresponding rainfall monitoring point.
pal <- colorFactor("Dark2", domain = all_records$smpt_code)
fillColor = ~pal(smpt_code)
# map
m <- leaflet() %>%
fitBounds(min(all_records$lon), min(all_records$lat),
max(all_records$lon), max(all_records$lat),
options = list(padding = c(25,25))) %>%
addTiles() %>%
## add sites
addCircleMarkers(data = sites, ~lon, ~lat,
color = "black",
weight = 0.5,
radius = 6,
layerId = ~smpt_code,
fillOpacity = 0.7,
opacity = 1,
fillColor = fillColor,
popup = ~paste("<strong>Site Id :</strong>", smpt_code, "<br>",
"<strong>Site Name :</strong>", smpt_name, "<br>")
) %>%
## add rainfall sites
addCircleMarkers(data = rainfall_data, ~long, ~lat,
color = "black",
weight = 0.5,
radius = 6,
layerId = ~notation,
fillOpacity = 0.7,
opacity = 1,
fillColor = "blue",
popup = ~paste("<strong>Site Id :</strong>", notation, "<br>",
"<strong>Site type :</strong>", qualifier, "<br>")
)
Then we create an interactive Plotly graph that shows how phosphate changes in time at the two sites we selected.
# time series graph
ts <- ggplot(all_records, aes(sample_datetime, concentration_adj, colour = smpt_name)) +
geom_line() +
geom_point() +
theme_bw() +
labs(title = "Phosphate concentration by site",
colour = "Site name") +
xlab("Sample date") +
ylab("Concentration (mg/l)")
ts1 <- ggplotly(ts)
Finally, we assemble the two tables, the map and the graph and a title in a tagList and then save it to an HTML file.
# assemble final document
doc <- htmltools::tagList(
h1(paste("Phosphate data summary - data extracted on ", now), style="text-align:center;"),
div(m),
div(div(sites_table, style="width:50%;"),
div(rainfall_table, style="width:50%;margin-left:2em;"), style="margin:2em;display:flex;"),
div(ts1)
)
htmltools::save_html(html = doc, file = "Phosphate_report.html")
And here is the final product!
Top comments (0)