Vignette
The EEAaq package allows users to retrieve air
quality data for multiple geographical zones, pollutants, and time
periods in a single request. Queries are submitted as lists, which
enables flexibility in specifying combinations of parameters.
zone_name:
Accepts only NUTS_ID and LATN_NAME for Nomenclature of Territorial Units
for Statistics (NUTS) and LAU_ID and LAU_NAME for Local Administrative
Units (LAU).
LAU_ISO:
For LAU nomenclature, you must specify the ISO code of the country
(e.g., IT for Italy) as LAU_ID values are not unique across
Europe.
pollutant:
You can specify the pollutant using its name (e.g., “NO2”, “PM10”) or
the corresponding Code from the pollutant dataframe (available via
EEAaq_get_dataframe()). Using the Code is recommended to
ensure accurate measurement units, especially for pollutants with
multiple unit-specific codes.
EEAaq_get_data
Below we demonstrate the use of query by using LAU_ID in one and
NUTS_ID in the other:
- `LAU_ID of LAU:
#Query NO2 data for a specific LAU zone
data_lau<- EEAaq::EEAaq_get_data(
zone_name = "15146", # LAU zone code
NUTS_level = "LAU", # NUTS level
LAU_ISO = "IT", # Country code for Italy
pollutants = "PM10", # Pollutant
from = "2022-01-01", # Start date
to = "2023-12-31", # End date
verbose = FALSE # Print detailed progress
)
# Preview the first few rows of the dataset
head(data_lau)
## # A tibble: 6 × 6
## AirQualityStationEoICode AirQualityStationName PM10 AveragingTime
## <chr> <chr> <dbl> <chr>
## 1 IT0477A MILANO - V.LE MARCHE 135. day
## 2 IT0477A MILANO - V.LE MARCHE 77.6 day
## 3 IT0477A MILANO - V.LE MARCHE 44.9 day
## 4 IT0477A MILANO - V.LE MARCHE 51.5 day
## 5 IT0477A MILANO - V.LE MARCHE 23.6 day
## 6 IT0477A MILANO - V.LE MARCHE 21.0 day
## # ℹ 2 more variables: DatetimeBegin <dttm>, DatetimeEnd <dttm>
#Identify the names of the areas from which to download the data
zones <- c("Région de Bruxelles-Capitale/Brussels Hoofdstedelijk Gewest","Vlaams Gewest","West-Nederland","Zuid-Nederland")
data <- EEAaq::EEAaq_get_data(
zone_name = zones, # LAU zone code
NUTS_level = "NUTS1", # NUTS level
pollutants = c("NO2", "PM10"), # Pollutant
from = "2023-01-01", # Start date
to = "2024-08-29", # End date
verbose = FALSE # Print detailed progress
)
base::unique(data$AirQualityStationEoICode)
## [1] "BELAL01" "BELAT83" "BELHB23" "BETB001" "BETB004" "BETB006" "BETB008"
## [8] "BETB011" "BETBUL1" "BETCHA1" "BETE013" "BETE714" "BETE716" "BETM802"
## [15] "BETMEU1" "BETN043" "BETR001" "BETR002" "BETR012" "BETR701" "BETR702"
## [22] "BETR721" "BETR740" "BETR801" "BETR802" "BETR803" "BETR804" "BETR805"
## [29] "BETR806" "BETR817" "BETR818" "BETR822" "BETR831" "BETR842" "BETR891"
## [36] "BETR897" "BETREG1" "BETVBX1" "BETVBX2" "BETVBX3" "NL00136" "NL00138"
## [43] "NL00236" "NL00237" "NL00240" "NL00241" "NL00247" "NL00546" "NL00551"
## [50] "NL00553" "NL00556" "NL00570" "NL00572" "NL00573" "NL00701" "NL00704"
If the query’s zone_name parameter corresponds to a valid
CITY_NAME (i.e., not NULL in the dataset), the function will
return the corresponding data. If no valid CITY_NAME is associated with
the zone_name, the function attempts to retrieve all available data for
the entire country and subsequently filter for the specified
zone_name.
Note: For very small towns or certain countries such
as Turkey or Albania, data may not currently be available in the
dataset.This limitation reflects the data unavailability at the EEA
Air Quality Viewer.
If the parameters used in the query include polygon or
quadrant, the function outputs an EEAaq_df_sfc object.
Otherwise, it returns an EEAaq_df object, which is a tibble
dataframe.
EEAaq map stations
EEAaq map stations and accepts as input either an object of the EEAaq
df class (output of the EEAaq_get_data function), or all
other parameters specifying the area and the pollutants.
#Map the stations using as EEAaq_df object, the output of the previous function
EEAaq::EEAaq_map_stations(data = data, bounds_level = "NUTS3", color = FALSE, dynamic = FALSE)
## Simple feature collection with 4 features and 8 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 2.546088 ymin: 50.688 xmax: 6.225231 ymax: 53.18511
## Geodetic CRS: WGS 84
## NUTS_ID LEVL_CODE CNTR_CODE
## 1 BE1 1 BE
## 2 BE2 1 BE
## 3 NL3 1 NL
## 4 NL4 1 NL
## NAME_LATN
## 1 Région de Bruxelles-Capitale/Brussels Hoofdstedelijk Gewest
## 2 Vlaams Gewest
## 3 West-Nederland
## 4 Zuid-Nederland
## NUTS_NAME MOUNT_TYPE
## 1 Région de Bruxelles-Capitale/Brussels Hoofdstedelijk Gewest NA
## 2 Vlaams Gewest NA
## 3 West-Nederland NA
## 4 Zuid-Nederland NA
## URBN_TYPE COAST_TYPE geometry
## 1 NA NA MULTIPOLYGON (((4.415738 50...
## 2 NA NA MULTIPOLYGON (((5.776583 50...
## 3 NA NA MULTIPOLYGON (((5.171192 52...
## 4 NA NA MULTIPOLYGON (((5.518671 51...
## points Country ISO AirQualityStationEoICode AirQualityStationNatCode AirQualityStationName Altitude NUTS1 NUTS1_ID NUTS2 NUTS2_ID NUTS3 NUTS3_ID LAU_NAME LAU_ID AirPollutant geometry n
Using the parameter bounds_level = “NUTS3”, the map is generated with
internal boundaries corresponding to the NUTS3 level.
The same output could be obtained specifying explicitly the zone
information.
EEAaq::EEAaq_map_stations(zone_name = zones, NUTS_level = "NUTS1",
pollutant = c("NO2", "PM10"), bounds_level = "NUTS3", color = FALSE, dynamic = FALSE)
## Simple feature collection with 4 features and 8 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 2.546088 ymin: 50.688 xmax: 6.225231 ymax: 53.18511
## Geodetic CRS: WGS 84
## NUTS_ID LEVL_CODE CNTR_CODE
## 1 BE1 1 BE
## 2 BE2 1 BE
## 3 NL3 1 NL
## 4 NL4 1 NL
## NAME_LATN
## 1 Région de Bruxelles-Capitale/Brussels Hoofdstedelijk Gewest
## 2 Vlaams Gewest
## 3 West-Nederland
## 4 Zuid-Nederland
## NUTS_NAME MOUNT_TYPE
## 1 Région de Bruxelles-Capitale/Brussels Hoofdstedelijk Gewest NA
## 2 Vlaams Gewest NA
## 3 West-Nederland NA
## 4 Zuid-Nederland NA
## URBN_TYPE COAST_TYPE geometry
## 1 NA NA MULTIPOLYGON (((4.415738 50...
## 2 NA NA MULTIPOLYGON (((5.776583 50...
## 3 NA NA MULTIPOLYGON (((5.171192 52...
## 4 NA NA MULTIPOLYGON (((5.518671 51...
## points Country ISO AirQualityStationEoICode AirQualityStationNatCode AirQualityStationName Altitude NUTS1 NUTS1_ID NUTS2 NUTS2_ID NUTS3 NUTS3_ID LAU_NAME LAU_ID AirPollutant geometry n
while LAU represents the finest level of granularity compared to
NUTS3 which shows only the boundaries of the corresponding
municipality.
EEAaq::EEAaq_map_stations(data =data_lau, color =
TRUE, dynamic=TRUE, ID= TRUE )
## Simple feature collection with 1 feature and 9 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 9.042173 ymin: 45.38829 xmax: 9.277104 ymax: 45.53535
## Geodetic CRS: WGS 84
## GISCO_ID ISO LAU_ID LAU_NAME POP_2021 POP_DENS_2 AREA_KM2 YEAR NUTS3_ID
## 1 IT_015146 IT 15146 Milano 1397715 NA 181.8155 2021 ITC4C
## geometry
## 1 MULTIPOLYGON (((9.120405 45...
## points Country ISO AirQualityStationEoICode AirQualityStationNatCode AirQualityStationName Altitude NUTS1 NUTS1_ID NUTS2 NUTS2_ID NUTS3 NUTS3_ID LAU_NAME LAU_ID AirPollutant geometry n
EEAaq summary
This function aims to describe the dataset that has been previously
imported, both at a global level, which means considering it as a whole,
and at the station level, where summary statistics and information are
grouped by monitoring station. In addition to these two levels, the
function also provides information about the gap length and the
correlation between pollutants if at least two pollutants are considered
simultaneously. The EEAaq summary function receives as
input an EEAaq df object, i.e. the output of the EEAaq get data
function.
summ <- EEAaq::EEAaq_summary(data = data)
## The dataset contains:
## ** 1188623 total observations
## ** 56 stations
## ** 14520 time stamps: from 2023-01-01 01:00:00 to 2024-08-29 01:00:00
## # A tibble: 2 × 8
## Pollutant NA_count NA_perc negative_count min mean max sd
## <chr> <int> <dbl> <int> <dbl> <dbl> <dbl> <dbl>
## 1 PM10 684273 57.6 5877 0 17.6 1565. 12.0
## 2 NO2 504350 42.4 76 0 19.6 299 14.1
summ$Summary_byStat$Mean_byStat
## # A tibble: 56 × 4
## AirQualityStationEoICode AirQualityStationName PM10 NO2
## <chr> <chr> <dbl> <dbl>
## 1 BELAL01 40AL01 - ANTWERPEN 17.6 19.6
## 2 BELAT83 40AT83 - BERENDRECHT 17.6 19.6
## 3 BELHB23 40HB23 - HOBOKEN 17.6 19.6
## 4 BETB001 41B001 - BRUSSEL (Kunst-Wet) 17.6 19.6
## 5 BETB004 41B004 - STE.CATHERI 17.6 19.6
## 6 BETB006 41B006 - PARL.EUROPE 17.6 19.6
## 7 BETB008 41B008 - Brussel (Beliardstraat) 17.6 19.6
## 8 BETB011 41B011 - BERCHEM S.A 17.6 19.6
## 9 BETBUL1 41BUL1 - BRUXELLES 17.6 19.6
## 10 BETCHA1 41CHA1 - GANSHOREN 17.6 19.6
## # ℹ 46 more rows
## # A tibble: 56 × 4
## AirQualityStationEoICode AirQualityStationName PM10_NO2 NO2_PM10
## <chr> <chr> <dbl> <dbl>
## 1 BELAL01 40AL01 - ANTWERPEN NA NA
## 2 BELAT83 40AT83 - BERENDRECHT NA NA
## 3 BELHB23 40HB23 - HOBOKEN NA NA
## 4 BETB001 41B001 - BRUSSEL (Kunst-Wet) NA NA
## 5 BETB004 41B004 - STE.CATHERI NA NA
## 6 BETB006 41B006 - PARL.EUROPE NA NA
## 7 BETB008 41B008 - Brussel (Beliardstraat) NA NA
## 8 BETB011 41B011 - BERCHEM S.A NA NA
## 9 BETBUL1 41BUL1 - BRUXELLES NA NA
## 10 BETCHA1 41CHA1 - GANSHOREN NA NA
## # ℹ 46 more rows
EEAaq time aggregate
Most pollutants are monitored hourly or daily, posing challenges for
interpretation and representation. The EEAaq_time_aggregate
function simplifies this by aggregating data into annual, monthly,
weekly, daily, or hourly intervals, generating summary statistics for
each station in an EEAaq_taggr_df object.
#Get the monthly minimum, maximum, mean and median values of the pollutant concentrations
t_aggr <- EEAaq::EEAaq_time_aggregate(data = data, frequency =
"monthly", aggr_fun = c("min", "max", "mean", "median" ))
## # A tibble: 1,093 × 11
## AirQualityStationEoICode AirQualityStationName Date PM10_min PM10_max
## <chr> <chr> <date> <dbl> <dbl>
## 1 BELAL01 40AL01 - ANTWERPEN 2023-01-01 3.9 77.4
## 2 BELAL01 40AL01 - ANTWERPEN 2023-02-01 5.4 76.4
## 3 BELAL01 40AL01 - ANTWERPEN 2023-03-01 4.4 87.4
## 4 BELAL01 40AL01 - ANTWERPEN 2023-04-01 3.9 82.9
## 5 BELAL01 40AL01 - ANTWERPEN 2023-05-01 8.4 165.
## 6 BELAL01 40AL01 - ANTWERPEN 2023-06-01 6.9 65.9
## 7 BELAL01 40AL01 - ANTWERPEN 2023-07-01 4.9 49.9
## 8 BELAL01 40AL01 - ANTWERPEN 2023-08-01 5.4 54.4
## 9 BELAL01 40AL01 - ANTWERPEN 2023-09-01 5.4 115.
## 10 BELAL01 40AL01 - ANTWERPEN 2023-10-01 4.9 51.4
## # ℹ 1,083 more rows
## # ℹ 6 more variables: PM10_mean <dbl>, PM10_median <dbl>, NO2_min <dbl>,
## # NO2_max <dbl>, NO2_mean <dbl>, NO2_median <dbl>
t_aggr$TimeAggr_byPollutant$PM10
## # A tibble: 1,093 × 7
## AirQualityStationEoICode AirQualityStationName Date min max mean
## <chr> <chr> <date> <dbl> <dbl> <dbl>
## 1 BELAL01 40AL01 - ANTWERPEN 2023-01-01 3.9 77.4 19.1
## 2 BELAL01 40AL01 - ANTWERPEN 2023-02-01 5.4 76.4 27.3
## 3 BELAL01 40AL01 - ANTWERPEN 2023-03-01 4.4 87.4 16.1
## 4 BELAL01 40AL01 - ANTWERPEN 2023-04-01 3.9 82.9 20.1
## 5 BELAL01 40AL01 - ANTWERPEN 2023-05-01 8.4 165. 24.5
## 6 BELAL01 40AL01 - ANTWERPEN 2023-06-01 6.9 65.9 24.9
## 7 BELAL01 40AL01 - ANTWERPEN 2023-07-01 4.9 49.9 14.8
## 8 BELAL01 40AL01 - ANTWERPEN 2023-08-01 5.4 54.4 15.3
## 9 BELAL01 40AL01 - ANTWERPEN 2023-09-01 5.4 115. 20.6
## 10 BELAL01 40AL01 - ANTWERPEN 2023-10-01 4.9 51.4 17.3
## # ℹ 1,083 more rows
## # ℹ 1 more variable: median <dbl>
EEAaq_idw_map
To enable quick and intuitive visual analysis, a function has been
developed to create spatial interpolation maps using the Inverse
Distance Weighting (IDW) method, introduced by D. Shepard in 1968. This
technique estimates the value of a variable at unknown locations by
calculating a weighted average of known values, with weights inversely
proportional to the distance from known points. Closer points contribute
more heavily to the estimate, making it a practical approach for
interpolating geolocated air quality data.
EEAaq::EEAaq_idw_map(data = t_aggr, pollutant = "NO2", aggr_fun =
"mean", distinct = TRUE, gradient = TRUE, idp = 2)
## Map initialization started at 2024-12-23 14:18:32.401687
## Map initialization ended at 2024-12-23 14:18:36.064331
## Computing IDW interpolation started at 2024-12-23 14:18:36.064457
## Computing IDW interpolation for: 2023-01-01, 1 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-02-01, 2 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-03-01, 3 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-04-01, 4 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-05-01, 5 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-06-01, 6 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-07-01, 7 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-08-01, 8 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-09-01, 9 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-10-01, 10 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-11-01, 11 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-12-01, 12 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2024-01-01, 13 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2024-02-01, 14 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2024-03-01, 15 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2024-04-01, 16 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2024-05-01, 17 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2024-06-01, 18 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2024-07-01, 19 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2024-08-01, 20 of 20
## [inverse distance weighted interpolation]
## Computing IDW interpolation ended at 2024-12-23 14:21:03.585806
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
##
## [[5]]
##
## [[6]]
##
## [[7]]
##
## [[8]]
##
## [[9]]
##
## [[10]]
##
## [[11]]
##
## [[12]]
##
## [[13]]
##
## [[14]]
##
## [[15]]
##
## [[16]]
##
## [[17]]
##
## [[18]]
##
## [[19]]
##
## [[20]]
t_aggr_1 <- EEAaq::EEAaq_time_aggregate(data = data_lau, frequency =
"monthly", aggr_fun = c("min", "max", "mean", "median" ))
EEAaq::EEAaq_idw_map(data = t_aggr_1, pollutant = "NO2", aggr_fun = "mean",
distinct = TRUE, gradient = FALSE, dynamic = TRUE, fill_NUTS_level = "LAU")
## Map initialization started at 2024-12-23 14:21:12.38424
## Map initialization ended at 2024-12-23 14:21:18.306087
## Computing IDW interpolation started at 2024-12-23 14:21:18.306213
## Computing IDW interpolation for: 2022-01-01, 1 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-02-01, 2 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-03-01, 3 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-04-01, 4 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-05-01, 5 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-06-01, 6 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-07-01, 7 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-08-01, 8 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-09-01, 9 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-10-01, 10 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-11-01, 11 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2022-12-01, 12 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-01-01, 13 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-02-01, 14 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-03-01, 15 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-04-01, 16 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-05-01, 17 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-06-01, 18 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-07-01, 19 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-08-01, 20 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-09-01, 21 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-10-01, 22 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-11-01, 23 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation for: 2023-12-01, 24 out of 24
## [inverse distance weighted interpolation]
## Computing IDW interpolation ended at 2024-12-23 14:23:41.563519