In adhil0/alohaal:

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

library(alohaal) #TODO: Change
library(tidyverse)

According to the EPA, lake ice duration can be an indicator of climate change. This is because lake ice is dependent on several environmental factors, so changes in these factors will influence the formation of ice on top of lakes. As a result, the study and analysis of lake ice formation can inform scientists about how quickly the climate is changing, and are critical to minimizing disruptions to lake ecosystems. We can examine the ice duration of Lake Mendota, Lake Monona, and Lake Wingra, three lakes in the Madison, WI area.

ntl_ice

lake_ice <-
  ggplot(data = ntl_ice %>% filter(!is.na(ice_duration)), aes(x = lakeid, y = ice_duration)) +
  geom_boxplot(aes(color = lakeid, shape = lakeid),
               alpha = 0.8,
               width = 0.5) +
  theme_minimal() +
  labs(title = "Ice Duration of Lakes in the Madison Area",
       y = "Ice Duration (Days)", x = "Lake")
lake_ice

Time Series Plot:

# Add title, labels, update legend, add trendlines
ice_duration <-
  ggplot(ntl_ice %>% filter(!is.na(ice_duration)), aes(x = year)) +
  geom_line(aes(y = ice_duration, color = lakeid), alpha = 0.7) +
  theme_minimal() +
  labs(x = "Year", y = "Ice Duration (Days)", title = "Ice Duration of Lakes in the Madison Area")
ice_duration

Since the box plots and time series plots are similar, we can average the ice duration across all three lakes:

# Average Ice duration across all three lakes
avg_ice <- ntl_ice %>%
  group_by(year) %>%
  summarise(ice_duration = mean(ice_duration)) %>%
  rename(avg_ice_duration = ice_duration)
avg_ice

Time Series Plot:

avg_ice_duration_plot <-
  ggplot(avg_ice %>% filter(!is.na(avg_ice_duration)), aes(x = year, y = avg_ice_duration)) +
  geom_line(alpha = 0.7) +
  theme_minimal() +
  geom_smooth(
    method = lm,
    se = FALSE,
    color = "blue",
    linetype = "dashed"
  ) + 
  labs(y = "Ice Duration (Days)", x = "Year", title = "Average Ice Duration")
avg_ice_duration_plot

The average ice duration in Madison is trending downwards. What environmental factors are influencing this change?

Research suggests that mean annual temperature is one of the primary factors that alter lake ice formation. We can look at the temperature data of Madison, WI found in ntl_temps to see if there is a corresponding change in climate that may have influenced the change in ice duration.

(Note that according to the original metadata: "Daily temperature data prior to 1884 were estimated from 3 times per day sampling and biases are expected and should not be comparable with data after that time.")

#ntl_temps

Compute the mean annual temperature:

#ntl_annual_avg <- ntl_temps %>% 
#  group_by(year) %>% 
#  summarise(ave_air_temp_adjusted = mean(ave_air_temp_adjusted))
#ntl_annual_avg

Time Series Plot:

#temp_plot <-
#  ggplot(ntl_annual_avg %>% filter(!is.na(ave_air_temp_adjusted)),
#         aes(x = year, y = ave_air_temp_adjusted)) +
#  geom_line(alpha = 0.7) +
#  theme_minimal() +
#  geom_smooth(
#    method = lm,
#    se = FALSE,
#    color = "blue",
#    linetype = "dashed"
#  ) +
#  labs(y = "Temperature", x = "Year", title = "Mean Annual Temperature in Madison, WI")
#temp_plot

There seems to be a general upward trend in the mean annual temperature. As a result, there may be a relationship between temperature and ice duration.

To compare the ice cover and temperature data directly, the ntl_temps and ntl_ice datasets can be joined together. A common variable needs to be identified to do this. In this case, the datasets can be joined on the year of observation. Since there is one temperature for each year in ntl_annual_avg, each temperature will be mapped to the corresponding year in the joined table. For example, a ave_air_temp_adjusted value of 6.416 degrees Celsius will be added to each row that has a year value of 1869.

#ntl_full_join <- ntl_annual_avg %>% 
#  full_join(ntl_ice)
#ntl_full_join

Plot the temperature against the ice duration of each of the three lakes:

#scatter <-
# ggplot(data = ntl_full_join %>%
#           filter(!is.na(ave_air_temp_adjusted) &
#                    !is.na(ice_duration)),
#         aes(y = ice_duration, x = ave_air_temp_adjusted)) + geom_point(alpha = 0.8) +
#  theme_minimal() +
#  labs(title = "Air Temperature and Ice Duration of Lakes in Madison, WI",
#       y = "Ice Duration (Days)", x = "Temperature (Celsius)") + geom_smooth(
#         method = "lm",
#         color = "black",
#         se = FALSE,
#         size = 0.3
#       )
#scatter

It appears that there is a negative correlation between the mean annual temperature and the amount of time that each lake was frozen. This means that as temperature increases, the ice duration decreases.

To compare the temperature with the ice duration on a time-series plot, calculate the percent change in ice duration and temperature from the earliest data point for both variables. The earliest year of observation for ice duration was 1851, and the earliest year of observation for temperature was 1869. These values can be determined from the avg_ice and ntl_annual_avg dataframes that were calculated earlier. First, join the full dataset with the average ice duration data that we calculated earlier, and then add two columns with the percent changes.

#pct_change <-
 # ntl_full_join %>% full_join(avg_ice) %>% mutate(
 #   pct_change_temp = (ave_air_temp_adjusted - ntl_annual_avg$ave_air_temp_adjusted[1]) *
 #     100 / ntl_annual_avg$ave_air_temp_adjusted[1],
 #   pct_change_ice = (avg_ice_duration - avg_ice$avg_ice_duration[1]) * 100 / avg_ice$avg_ice_duration[1]
#  ) %>% select(year, pct_change_temp, pct_change_ice)

#pct_change <- unique(arrange(pct_change, year))

Then plot a time series to compare the magnitude in changes in temperature and ice duration. Use the pivot longer function so that both lines can easily be plotted on the same axes.

#comparison_df <-
#  pct_change  %>% pivot_longer(cols = c(pct_change_ice, pct_change_temp))

Time Series Plot:

#comparison_plot <-
#  ggplot(comparison_df %>% filter((name == "pct_change_ice" &
#                                     !is.na(value)) |
#                                    name == "pct_change_temp"),
#         aes(x = year)) +
#  geom_line(aes(y = value, color = name)) +
#  theme_minimal() +
#  labs(x = "Year", y = "Percent Change", title = "Percent Change in Ice Duration and Temperature")
#comparison_plot

It appears that the change in ice somewhat mirrors the change in temperature each year, as several spikes in temperature have a corresponding dip in ice duration. This may indicate a relationship between the two variables.

Citation

Anderson, L. and D. Robertson. 2020. Madison Wisconsin Daily Meteorological Data 1869 - current ver 32. Environmental Data Initiative. https://doi.org/10.6073/pasta/e3ff85971d817e9898bb8a83fb4c3a8b (Accessed 2021-03-08).

Magnuson, J.J., S.R. Carpenter, and E.H. Stanley. 2021. North Temperate Lakes LTER: Ice Duration - Madison Lakes Area 1853 - current ver 35. Environmental Data Initiative. https://doi.org/10.6073/pasta/ab31f2489ee436beb73fc8f1d0213d97 (Accessed 2021-03-08).