Citizen Science vs Professional Data comparisons

Citizen Science

Task Load CS dataset and clean data

install.packages(c("tidyverse", "skimr", "janitor", "corrplot"))

trying URL 'https://cran.rstudio.com/bin/macosx/big-sur-x86_64/contrib/4.5/tidyverse_2.0.0.tgz'
trying URL 'https://cran.rstudio.com/bin/macosx/big-sur-x86_64/contrib/4.5/skimr_2.2.2.tgz'
trying URL 'https://cran.rstudio.com/bin/macosx/big-sur-x86_64/contrib/4.5/janitor_2.2.1.tgz'
trying URL 'https://cran.rstudio.com/bin/macosx/big-sur-x86_64/contrib/4.5/corrplot_0.95.tgz'

The downloaded binary packages are in
    /var/folders/5_/qmm6kpd90mggm8b4nhx7s09c0000gn/T//RtmpR5woME/downloaded_packages

library(tidyverse)
library(skimr)
library(janitor)
library(corrplot)
library(tidyr)

#load CS dataset from CSV
cs_raw <- read.csv("BA_records.csv",
                       header = TRUE, stringsAsFactors = FALSE)

#omit "incorrect" verification substatus

clean_data <- subset(cs_raw, verification_substatus_text !="incorrect")

#clean latitude values to only include entries in ireland

clean_latdata <- subset(clean_data, lat >= 51 & lat <= 55.5)

cs_data <- subset(clean_latdata, lon >= -10.25 & lon <= -5.5)

#remove "notes" column (irrelevant data)
cs_data$notes <- NULL

#remove "location" column (irrelevant data)
cs_data$location <- NULL

##Marine Data

Task Load professional datasets, merge marine datasets, and clean data

#load marine datasets

#1997 to 2016 data
marine1_data <- read.csv("HABs_raw_CoCliME.csv",
                       header = TRUE, stringsAsFactors = FALSE)


#2011 to 2019 data
marine2_data <- read.csv("HABs_Biotoxin_weekly_max_2011_2019.csv",
                         header = TRUE, stringsAsFactors = FALSE)

#clean marine2 dataset to match columns

marine2_data$species <- NULL
marine2_data$HY_OA <- NULL
marine2_data$HY_DTX1 <- NULL
marine2_data$HY_DTX2 <- NULL
marine2_data$PTX1 <- NULL
marine2_data$PTX2 <- NULL
marine2_data$D_spp <- NULL

marine2_rdata <- rename(marine2_data, "sample_date" = "End_of_Week")

#merge marine datasets
marine12_data <- rbind(marine1_data, marine2_rdata)

#remove duplicate points
marine_cdata <- marine12_data[!duplicated(marine_data),]

#remove rows that sum to zeros (no species recorded)
marine_data <- marine_cdata[rowSums(marine_cdata[, c(5, 6, 7, 8, 9, 10, 11, 12)] > 0) != 0, ]

##Freshwater Data

Task Load professional dataset and clean data

#load freshwater dataset

fw_raw <- read.csv("PhytoplanktonSurveyData_2007_2025.csv",)

#remove extra column
fw_raw$Local.Authority <- NULL
fw_raw$WFD.Lake.Waterbody.Code <- NULL

#change "Not found" coordinate entries to NA
fw_raw$Easting[fw_raw$Easting == "Not found"] <- NA
fw_raw$Northing[fw_raw$Northing == "Not found"] <- NA

#drop any rows with NA
fw_data <- drop_na(fw_raw)

#reformat dates
fw_data$Sample.date <- as.Date(fw_data$Sample.date, format = "%m/%d/%y")

#reformat northings and eastings to Irish Grid coordinates

install.packages("sf")

trying URL 'https://cran.rstudio.com/bin/macosx/big-sur-x86_64/contrib/4.5/sf_1.1-2.tgz'
Content type 'application/x-gzip' length 31655373 bytes (30.2 MB)
==================================================
downloaded 30.2 MB

The downloaded binary packages are in
    /var/folders/5_/qmm6kpd90mggm8b4nhx7s09c0000gn/T//RtmpR5woME/downloaded_packages

library(sf)

# Convert Easting/Northing to an sf object
sites_sf <- st_as_sf(fw_data,
                     coords = c("Easting", "Northing"),
                     crs = 29903)

# Transform to WGS84 (latitude/longitude)
sites_ll <- st_transform(sites_sf, 4326)

# Extract coordinates
coords <- st_coordinates(sites_ll)

fw_data$Longitude <- coords[, 1]
fw_data$Latitude <- coords[, 2]

#remove "Easting" and "Northing" columns
fw_data$Easting <- NULL
fw_data$Northing <- NULL

#reorder columns
fw_data <- fw_data[, c(3, 2, 1, 14:15, 4:13:ncol(fw_data))]

#remove repeat columns
fw_data$Longitude.1 <- NULL
fw_data$Latitude.1 <- NULL

#remove rows that sum to zeros (no species recorded)
fw_data <- fw_data[rowSums(fw_data[, c(7, 8, 9, 10, 11, 12, 13, 14, 15)] > 0) != 0, ]

#Time Series Trends

#yearly counts

#load packages
library(dplyr)
library(lubridate)

#CS yearly
cs_yearly <- cs_data %>%
  mutate(Year = year(sample_date)) %>%
  group_by(Year) %>%
  summarise(Observations = n()) %>%
  mutate(Dataset = "Blooming Algae CS Data")

#marine yearly
marine_yearly <- marine_data %>%
  mutate(Year = year(sample_date)) %>%
  group_by(Year) %>%
  summarise(Observations = n()) %>%
  mutate(Dataset = "Marine Institute Coastal Data")

#FW yearly
fw_yearly <- fw_data %>%
  mutate(Year = year(Sample.date)) %>%
  group_by(Year) %>%
  summarise(Observations = n()) %>%
  mutate(Dataset = "EPA Freshwater Data")

#combine yearly counts
annual_counts <- bind_rows(
  cs_yearly,
  marine_yearly,
  fw_yearly
)

#remove any observation counts that do not have a year
annual_counts <- drop_na(annual_counts)

#plot
library(ggplot2)

ggplot(annual_counts,
       aes(x = Year,
           y = Observations,
           colour = Dataset)) +

  geom_line(linewidth = 1.2) +
  geom_point(size = 2) +

  geom_smooth(
    method = "lm",
    se = FALSE,
    linetype = "dashed"
  ) +

  labs(
    title = "Annual Harmful Algal Bloom Observations",
    x = "Year",
    y = "Number of observations",
    colour = "Dataset"
  ) +

  theme_classic(base_size = 14)

#standardize datasets
annual_standardized <- annual_counts %>%
  group_by(Dataset) %>%
  mutate(
    Std_Obs = (Observations - min(Observations)) /
              (max(Observations) - min(Observations))
  ) %>%
  ungroup()

library(ggplot2)

#plot
ggplot(
  annual_standardized,
  aes(
    x = Year,
    y = Std_Obs,
    colour = Dataset
  )
) +

  geom_line(linewidth = 1.2) +

  geom_point(size = 2) +

  geom_smooth(
    method = "lm",
    se = FALSE,
    linetype = "dashed"
  ) +

  labs(
    title = "Standardized Annual HAB Observation Trends",
    subtitle = "Comparison of temporal patterns between monitoring programmes",
    x = "Year",
    y = "Standardized Observation Frequency",
    colour = "Dataset"
  ) +

  theme_classic(base_size = 14)

NA
NA
NA

#percentage of maximum
annual_standardized2 <- annual_counts %>%
  group_by(Dataset) %>%
  mutate(
    Percent_Max = Observations /
                  max(Observations) * 100
  )

library(ggplot2)

#plot
ggplot(
  annual_standardized2,
  aes(
    x = Year,
    y = Percent_Max,
    colour = Dataset
  )
) +

  geom_line(linewidth = 1.2) +

  geom_point(size = 0) +

  geom_smooth(
    method = "lm",
    se = FALSE,
    linetype = "dashed"
  ) +

  labs(
    title = "Standardized Annual HAB Observation Trends",
    subtitle = "Comparison of temporal patterns between monitoring programmes",
    x = "Year",
    y = "Standardized Observation Frequency",
    colour = "Dataset"
  ) +

  theme_classic(base_size = 14)

NA
NA
NA

#Z-score standardization
annual_standardized <- annual_counts %>%
  group_by(Dataset) %>%
  mutate(
    Zscore = scale(Observations)[,1]
  )

library(ggplot2)

#plot
ggplot(
  annual_standardized,
  aes(
    x = Year,
    y = Zscore,
    colour = Dataset
  )
) +

  geom_line(linewidth = 1.2) +

  geom_point(size = 2) +

  geom_smooth(
    method = "lm",
    se = FALSE,
    linetype = "dashed"
  ) +

  labs(
    title = "Standardized Annual HAB Observation Trends",
    subtitle = "Comparison of temporal patterns between monitoring programmes",
    x = "Year",
    y = "Standardized Observation Frequency",
    colour = "Dataset"
  ) +

  theme_classic(base_size = 14)

NA
NA
NA

#Monthly Observations

cs_monthly <- cs_data %>%
  mutate(Month = month(sample_date,
                       label = TRUE,
                       abbr = TRUE)) %>%
  group_by(Month) %>%
  summarise(Observations = n(), .groups = "drop") %>%
  mutate(Dataset = "Blooming Algae CS Data")

marine_monthly <- marine_data %>%
  mutate(Month = month(sample_date,
                       label = TRUE,
                       abbr = TRUE)) %>%
  group_by(Month) %>%
  summarise(Observations = n(), .groups = "drop") %>%
  mutate(Dataset = "Marine Institute Coastal Data")

fw_monthly <- fw_data %>%
  mutate(Month = month(Sample.date,
                       label = TRUE,
                       abbr = TRUE)) %>%
  group_by(Month) %>%
  summarise(Observations = n(), .groups = "drop") %>%
  mutate(Dataset = "EPA Freshwater Data")

#combine rows
monthly_counts <- bind_rows(
  cs_monthly,
  marine_monthly,
  fw_monthly
)

#remove any observation counts that do not have a month
monthly_counts <- drop_na(monthly_counts)

#plot
ggplot(monthly_counts,
       aes(x = Month,
           y = Observations,
           colour = Dataset,
           group = Dataset)) +

  geom_line(linewidth = 1.2) +
  geom_point(size = 3) +

  labs(
    title = "Seasonal Distribution of HAB Observations",
    x = "Month",
    y = "Number of observations"
  ) +

  theme_classic(base_size = 14)

#Compare Datasets

#calculate unique locations for each dataset

#marine
length(unique(marine_data[["location"]]))

[1] 184

#fresh water
length(unique(fw_data[["Site.Code"]]))

[1] 310

#cs data
n_distinct(cs_data$lat, cs_data$lon)

[1] 627

##EXPORT DATASETS

#cs data
write.csv(cs_data, file = "cs_data.csv")

#marine data
write.csv(marine_data, file = "marine_data.csv")

#fw data
write.csv(fw_data, file = "fw_data.csv")