Important Notes

The purpose of this script is to determine species richness by species group, create managed area statistics, generate summary plots, and create reports in pdf form for Coastal Wetlands data.

These scripts were created by J.E. Panzik () for SEACAR. Updated by T.G. Hill ().

All scripts and outputs can be found on the SEACAR GitHub repository:

https://github.com/FloridaSEACAR/SEACAR_Trend_Analyses

This markdown file is designed to be compiled by CW_Plot_Render.R (https://github.com/FloridaSEACAR/SEACAR_Trend_Analyses/blob/main/Coastal_Wetlands/CW_Plot_Render.R).

Coastal Wetlands Species Richness

Apalachicola Bay Aquatic Preserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Apalachicola Bay Aquatic Preserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Mangroves and associates 2 2022 2023 4 1 1 1.0 1.00 0.00 2022 2018
Marsh 10 2014 2023 144 1 5 1.5 2.08 1.21 2022 2018
Marsh succulents 10 2014 2023 56 1 3 3.0 2.20 0.96 2022 2018

Between 2022 and 2023, the median annual number of species for mangroves and associates was 1 based on 4 observations. Between 2014 and 2023, the median annual number of species for marsh was 1.5 based on 144 observations. Between 2014 and 2023, the median annual number of species for marsh succulents was 3 based on 56 observations.

Apalachicola National Estuarine Research Reserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Apalachicola National Estuarine Research Reserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Mangroves and associates 2 2022 2023 4 1 1 1.0 1.00 0.00 2022 2018
Marsh 10 2014 2023 144 1 5 1.5 2.08 1.21 2022 2018
Marsh succulents 10 2014 2023 56 1 3 3.0 2.20 0.96 2022 2018

Between 2022 and 2023, the median annual number of species for mangroves and associates was 1 based on 4 observations. Between 2014 and 2023, the median annual number of species for marsh was 1.5 based on 144 observations. Between 2014 and 2023, the median annual number of species for marsh succulents was 3 based on 56 observations.

Big Bend Seagrasses Aquatic Preserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Big Bend Seagrasses Aquatic Preserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Marsh 2 2013 2016 6 1 3 2 2.17 0.75 2013 2016

Between 2013 and 2016, the median annual number of species for marsh was 2 based on 6 observations.

Cockroach Bay Aquatic Preserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Cockroach Bay Aquatic Preserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Mangroves and associates 1 2015 2015 1 2 2 2 2 - 2015 2015

In the year 2015, 2 species were observed for mangroves and associates based on 1 observation.

Guana River Marsh Aquatic Preserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Guana River Marsh Aquatic Preserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Marsh 12 2012 2023 795 1 3 1 1.0 0.08 2012 2021
Marsh succulents 12 2012 2023 58 1 2 1 1.1 0.31 2012 2021

Between 2012 and 2023, the median annual number of species for marsh was 1 based on 795 observations. Between 2012 and 2023, the median annual number of species for marsh succulents was 1 based on 58 observations.

Guana Tolomato Matanzas National Estuarine Research Reserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Guana Tolomato Matanzas National Estuarine Research Reserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Mangroves and associates 12 2012 2023 387 1 2 1 1.02 0.14 2013 2017
Marsh 12 2012 2023 2568 1 3 1 1.12 0.33 2013 2017
Marsh succulents 12 2012 2023 868 1 4 1 1.38 0.52 2013 2017

Between 2012 and 2023, the median annual number of species for mangroves and associates was 1 based on 387 observations. Between 2012 and 2023, the median annual number of species for marsh was 1 based on 2,568 observations. Between 2012 and 2023, the median annual number of species for marsh succulents was 1 based on 868 observations.

Nature Coast Aquatic Preserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Nature Coast Aquatic Preserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Mangroves and associates 1 2020 2020 8 1 2 1 1.38 0.52 2020 2020
Marsh 1 2020 2020 14 1 3 1 1.29 0.61 2020 2020

In the year 2020, 1 species were observed for mangroves and associates based on 8 observations. In the year 2020, 1 species were observed for marsh based on 14 observations.

Pellicer Creek Aquatic Preserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Pellicer Creek Aquatic Preserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Mangroves and associates 2 2020 2021 2 1 1 1 1.00 0.00 2020 2018
Marsh 12 2012 2023 564 1 3 1 1.41 0.50 2020 2018
Marsh succulents 12 2012 2023 353 1 4 1 1.24 0.45 2020 2018

Between 2020 and 2021, the median annual number of species for mangroves and associates was 1 based on 2 observations. Between 2012 and 2023, the median annual number of species for marsh was 1 based on 564 observations. Between 2012 and 2023, the median annual number of species for marsh succulents was 1 based on 353 observations.

Pine Island Sound Aquatic Preserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Pine Island Sound Aquatic Preserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Marsh succulents 1 2004 2004 17 1 2 1 1.12 0.33 2004 2004

In the year 2004, 1 species were observed for marsh succulents based on 17 observations.

Pinellas County Aquatic Preserve

Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Line graph of annual average coastal wetlands species richness over time for mangroves and associates (triangles), marsh (squares), and marsh succulents (circles). If the time series by species group included more than one year of observations, a line connects data points for visualization.
Coastal Wetlands Species Richness - Pinellas County Aquatic Preserve
SpeciesGroup1 N-Years EarliestYear LatestYear N-Data Min Max Median Mean StDev Year-MinRichness Year-MaxRichness
Mangroves and associates 1 2014 2014 1 3 3 3 3 - 2014 2014
Marsh 1 2014 2014 1 5 5 5 5 - 2014 2014
Marsh succulents 1 2014 2014 1 3 3 3 3 - 2014 2014

In the year 2014, 3 species were observed for mangroves and associates based on 1 observation. In the year 2014, 5 species were observed for marsh based on 1 observation. In the year 2014, 3 species were observed for marsh succulents based on 1 observation.

Libraries and Settings

Loads libraries used in the script. The inclusion of scipen option limits how frequently R defaults to scientific notation. Sets default settings for displaying warning and messages in created document, and sets figure dpi.

library(knitr)
library(data.table)
library(dplyr)
library(lubridate)
library(ggplot2)
library(ggpubr)
library(scales)
library(EnvStats)
library(tidyr)
library(kableExtra)
library(glue)
library(grid)
library(stringr)
options(scipen=999)

File Import

Imports file that is determined in the CW_Plot_Render.R script.

The command fread is used because of its improved speed while handling large data files. Only columns that are used by the script are imported from the file, and are designated in the select input.

The script then gets the name of the parameter as it appears in the data file and units of the parameter.

The latest version of Coastal Wetlands data is available at: https://usf.box.com/s/jpwsi8kram54xt6zyma5wo9mqferddcn

The file being used for the analysis is: All_CW_Parameters-2025-Sep-04.txt

#Import data from coastal wetlands file
data <- fread(file_in, sep="|", header=TRUE, stringsAsFactors=FALSE,
              na.strings=c("NULL"))

file_short <- tail(str_split(file_in, "/")[[1]],1)
cat(paste("The data file used is:", file_short, sep="\n"))

Data Filtering

The processing and filtering that is done to the data is as follows:

  1. Only take data rows that are percent cover
  2. Only keep data rows that are Marsh, Marsh succulents, and Mangroves and associates
  3. Set parameter names to Species Richness
  4. Sets units
  5. Removes rows that contains NA values in ManagedAreaName, GenusName, SpeciesName, Month, Year, SpeciesGroup1, and removes invasive species data
  6. Sets ResultValue to be numeric values and removes rows where percent cover is 0
  7. Removes duplicates (MADup==1)
  8. Combines genus and species names
  9. Corrects some managed area names to match what is being used with other habitats
  10. Creates species richness dataset
    • Grouped based on common ManagedAreaName, ProgramID, ProgramName, ProgramLocationID, and SampleDate
    • SpeciesRichness determined based on the number of unique species (gensp) in each group
  11. Writes to file with “_UsedData” file name to indicate what data was used for species richness.
# Only interested in Percent Cover measurements
data <- data[ParameterName=="Percent Cover", ]
# Only keep data rows that are Marsh, Marsh succulents, and Mangroves and assoc.
keep_spg <- c("Marsh","Marsh succulents","Mangroves and associates")
data <- data[SpeciesGroup1 %in% keep_spg, ]
# Create ParameterName Column
data$ParameterName <- "Species Richness"
parameter <- "Species Richness"

# Sets units for species richness
unit <- "# of species"
data$ParameterUnits <- unit

# Remove rows with missing ManagedAreaName
data <- data[!is.na(data$ManagedAreaName),]
data <- data[data$ManagedAreaName!="NA",]
# Remove rows with missing GenusName
data <- data[!is.na(data$GenusName),]
# Remove rows with missing SpeciesName
data <- data[!is.na(data$SpeciesName),]
# Remove rows with missing Months
data <- data[!is.na(data$Month),]
# Remove rows with missing Years
data <- data[!is.na(data$Year),]
# Set ResultValue to be a number value
data$ResultValue <- as.numeric(data$ResultValue)
# Remove rows where ResultValue is 0
data <- data[data$ResultValue!=0,]
# Remove duplicate rows
data <- data[data$MADup==1,]
# Create variable that combines the genus and species name
data$gensp <- paste(data$GenusName, data$SpeciesName, sep=" ")

# Create Species Richness values for groups of unique combinations of
# ManagedAreaName, ProgramID, ProgramName, ProgramLocationID, and SampleDate.
data <- data %>%
  group_by(ManagedAreaName, ProgramID, ProgramName, ProgramLocationID,
           SampleDate, SpeciesGroup1) %>%
  summarise(ParameterName=parameter,
            Year=unique(Year), Month=unique(Month),
            SpeciesRichness=length(unique(gensp)))

# Adds AreaID for each managed area by combining the MA_All datatable to the
# data based on ManagedAreaName
data <- merge.data.frame(MA_All[,c("AreaID", "ManagedAreaName")],
                         data, by="ManagedAreaName")

# Writes this data that is used by the rest of the script to a text file
fwrite(data, paste0(out_dir,"/CoastalWetlands_", param_file, "_UsedData.txt"),
       sep="|")

# Makes sure SampleDate is being stored as a Date object
data$SampleDate <- as.Date(data$SampleDate)

# Creates a variable with the names of all the managed areas that contain
# species observations
cw_MA_Include <- unique(data$ManagedAreaName[!is.na(data$SpeciesRichness)])

# Puts the managed areas in alphabetical order
cw_MA_Include <- cw_MA_Include[order(cw_MA_Include)]

# Determines the number of managed areas used
n <- length(cw_MA_Include)

Managed Area Statistics

Gets summary statistics for each managed area. Uses piping from dplyr package to feed into subsequent steps. The following steps are performed:

  1. Group data that have the same ManagedAreaName, Year, Month, and SpeciesGroup.
    • Second summary statistics do not use the Month grouping and are only for ManagedAreaName, Year, and SpeciesGroup..
    • Third summary statistics do not use Year grouping and are only for ManagedAreaName, Month, and SpeciesGroup.
    • Fourth summary statistics are only grouped based on ManagedAreaName and SpeciesGroup.
      • Determines the years that the minimum and maximum species richness occurred
  2. For each group, provide the following information: Parameter Name (ParameterName), Number of Entries (N_Data), Lowest Value (Min), Largest Value (Max), Median, Mean, Standard Deviation, and a list of all Programs included in these measurements.
  3. Sort the data in ascending (A to Z and 0 to 9) order based on ManagedAreaName then Year then Month
  4. Write summary stats to a pipe-delimited .txt file in the output directory 
# Create summary statistics for each managed area based on Year and Month
# intervals.
MA_YM_Stats <- data %>%
  group_by(AreaID, ManagedAreaName, Year, Month, SpeciesGroup1) %>%
  summarize(ParameterName=parameter,
            N_Data=length(na.omit(SpeciesRichness)),
            Min=min(SpeciesRichness),
            Max=max(SpeciesRichness),
            Median=median(SpeciesRichness),
            Mean=mean(SpeciesRichness),
            StandardDeviation=sd(SpeciesRichness),
            Programs=paste(sort(unique(ProgramName), decreasing=FALSE),
                           collapse=', '),
            ProgramIDs=paste(sort(unique(ProgramID), decreasing=FALSE),
                             collapse=', '))
# Puts the data in order based on ManagedAreaName, Year, then Month
MA_YM_Stats <- as.data.table(MA_YM_Stats[order(MA_YM_Stats$ManagedAreaName,
                                               MA_YM_Stats$Year,
                                               MA_YM_Stats$Month), ])
# Writes summary statistics to file
fwrite(MA_YM_Stats, paste0(out_dir,"/CoastalWetlands_", param_file,
                           "_MA_MMYY_Stats.txt"), sep="|")
# Removes variable storing data to improve computer memory
rm(MA_YM_Stats)

# Create summary statistics for each managed area based on Year intervals
MA_Y_Stats <- data %>%
  group_by(AreaID, ManagedAreaName, Year, SpeciesGroup1) %>%
  summarize(ParameterName=parameter,
            N_Data=length(na.omit(SpeciesRichness)),
            Min=min(SpeciesRichness),
            Max=max(SpeciesRichness),
            Median=median(SpeciesRichness),
            Mean=mean(SpeciesRichness),
            StandardDeviation=sd(SpeciesRichness),
            Programs=paste(sort(unique(ProgramName), decreasing=FALSE),
                           collapse=', '),
            ProgramIDs=paste(sort(unique(ProgramID), decreasing=FALSE),
                             collapse=', '))
# Puts the data in order based on ManagedAreaName then Year
MA_Y_Stats <- as.data.table(MA_Y_Stats[order(MA_Y_Stats$ManagedAreaName,
                                             MA_Y_Stats$Year), ])
# Writes summary statistics to file
fwrite(MA_Y_Stats, paste0(out_dir,"/CoastalWetlands_", param_file,
                          "_MA_Yr_Stats.txt"), sep="|")

# Create summary statistics for each managed area based on Month intervals.
MA_M_Stats <- data %>%
  group_by(AreaID, ManagedAreaName, Month, SpeciesGroup1) %>%
  summarize(ParameterName=parameter,
            N_Data=length(na.omit(SpeciesRichness)),
            Min=min(SpeciesRichness),
            Max=max(SpeciesRichness),
            Median=median(SpeciesRichness),
            Mean=mean(SpeciesRichness),
            StandardDeviation=sd(SpeciesRichness),
            Programs=paste(sort(unique(ProgramName), decreasing=FALSE),
                           collapse=', '),
            ProgramIDs=paste(sort(unique(ProgramID), decreasing=FALSE),
                             collapse=', '))
# Puts the data in order based on ManagedAreaName then Month
MA_M_Stats <- as.data.table(MA_M_Stats[order(MA_M_Stats$ManagedAreaName,
                                             MA_M_Stats$Month), ])
# Writes summary statistics to file
fwrite(MA_M_Stats, paste0(out_dir,"/CoastalWetlands_", param_file,
                          "_MA_Mo_Stats.txt"), sep="|")
# Removes variable storing data to improve computer memory
rm(MA_M_Stats)

# Create summary overall statistics for each managed area.
MA_Ov_Stats <- data %>%
  group_by(AreaID, ManagedAreaName, SpeciesGroup1) %>%
  summarize(ParameterName=parameter,
            N_Years=length(unique(na.omit(Year))),
            EarliestYear=min(Year),
            LatestYear=max(Year),
            N_Data=length(na.omit(SpeciesRichness)),
            Min=min(SpeciesRichness),
            Max=max(SpeciesRichness),
            Median=median(SpeciesRichness),
            Mean=mean(SpeciesRichness),
            StandardDeviation=sd(SpeciesRichness),
            Programs=paste(sort(unique(ProgramName), decreasing=FALSE),
                           collapse=', '),
            ProgramIDs=paste(sort(unique(ProgramID), decreasing=FALSE),
                             collapse=', '))
# Puts the data in order based on ManagedAreaName
MA_Ov_Stats <- as.data.table(MA_Ov_Stats[order(MA_Ov_Stats$ManagedAreaName), ])
# Creates Year_MinRichness and Year_MaxRichness columns
MA_Ov_Stats$Year_MinRichness <- NA
MA_Ov_Stats$Year_MaxRichness <- NA

# Loops through each ManagedAreaName.
# Determines what year the minimum and maximum species richness occurred
for(m in 1:nrow(MA_Ov_Stats)){
  # Stores ManagedAreaName for this row
  ma <- MA_Ov_Stats$ManagedAreaName[m]
  
  # Skips to next row if there are no data for this combination
  if(MA_Ov_Stats$N_Data[m]==0){
    next
  }
  # Gets subset of data from MA_Y_Stats (yearly summary stats) with this
  # ManagedAreaName
  ds <- MA_Y_Stats[MA_Y_Stats$ManagedAreaName==ma,]
  # Gets the minimum and maximum Mean (yearly averages)
  min <- min(ds$Mean)
  max <- max(ds$Mean)
  #Determines what years those minimum and maximum values occured
  year_min <- ds$Year[ds$Mean==min]
  year_max <- ds$Year[ds$Mean==max]
  # Stores the occurrence years of the minimum and maximum into the overall
  # stats for this row
  MA_Ov_Stats$Year_MinRichness[m] <- year_min
  MA_Ov_Stats$Year_MaxRichness[m] <- year_max
}
# Replaces blank ProgramIDs with NA (missing values)
MA_Ov_Stats$ProgramIDs <- replace(MA_Ov_Stats$ProgramIDs,
                                  MA_Ov_Stats$ProgramIDs=="", NA)
MA_Ov_Stats$Programs <- replace(MA_Ov_Stats$Programs,
                                MA_Ov_Stats$Programs=="", NA)
# Write overall statistics to file
fwrite(MA_Ov_Stats, paste0(out_dir,"/CoastalWetlands_", param_file,
                           "_MA_Overall_Stats.txt"), sep="|")
# Removes entries from the overall statistics that do not have data.
# Based on presence or absence of EarliestYear
MA_Ov_Stats <- MA_Ov_Stats[!is.na(MA_Ov_Stats$EarliestYear), ]

Appendix I: Managed Area Species Richness

The plots shown here are the species richness for each managed area with a yearly average.

  1. Set common plot theme.
  2. Determine the earliest and latest year of the data to create x-axis scale and intervals
  3. Determine the upper and lower limit of the plot for better y-axis labels
  4. Determines what species groups are present and adjusts legend entries
  5. Add the plot line
  6. Set the plot type as a point plot with the size of the points
  7. Create the title, x-axis, y-axis, and color fill labels
  8. Set the y and x limits
  9. Apply common plot theme
  10. Add table with summary statistics below each figure
  • Numerical non-integer values are rounded to 2 decimal places
  • StandardDeviation is renamed StDev for space reasons
  1. Create file name to save figure
  2. Save figure as png file
# Defines standard plot theme: black and white, no major or minor grid lines,
# Arial font. Title is centered, size 12, and blue (hex coded). Subtitle is
# centered, size 10, and blue (hex coded). Legend title is size 10 and the
# legend is left-justified. X-axis title is size 10 and the margins are padded
# at the top and bottom to give more space for angled axis labels. Y-axis title
# is size 10 and margins are padded on the right side to give more space for
# axis labels. Axis labels are size 10 and the x-axis labels are rotated -45
# degrees with a horizontal justification that aligns them with the tick mark
plot_theme <- theme_bw() +
  theme(panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        text=element_text(family="Arial"),
        plot.title=element_text(hjust=0.5, size=12, color="#314963"),
        plot.subtitle=element_text(hjust=0.5, size=10, color="#314963"),
        legend.title=element_text(size=10),
        legend.text = element_text(hjust=0),
        axis.title.x = element_text(size=10, margin = margin(t = 5, r = 0,
                                                             b = 10, l = 0)),
        axis.title.y = element_text(size=10, margin = margin(t = 0, r = 10,
                                                             b = 0, l = 0)),
        axis.text=element_text(size=10),
        axis.text.x=element_text(angle = -45, hjust = 0))

# Color palette for SEACAR
color_palette <- c("#005396", "#0088B1", "#00ADAE", "#65CCB3", "#AEE4C1",
                   "#FDEBA8", "#F8CD6D", "#F5A800", "#F17B00")

# All unique SpeciesGroup1 values get assigned a shape and color
cw_groups <- sort(unique(MA_Y_Stats$SpeciesGroup1), decreasing = T)

group_colors <- color_palette[seq_len(length(cw_groups))]
group_shapes <- c(21,22,24,25)
names(group_colors) <- cw_groups
names(group_shapes) <- cw_groups

# Loop that cycles through each managed area with data
if(n==0){
  # Prints a statement if there are no managed areas with appropriate data
  print("There are no monitoring locations that qualify.")
} else {
  for (i in 1:n) {
    ma_i <- cw_MA_Include[i]
    ma_abrev <- MA_All[ManagedAreaName==ma_i, Abbreviation]
    # Gets data for target managed area
    plot_data <- MA_Y_Stats[MA_Y_Stats$ManagedAreaName==ma_i]
    # Determines most recent year with available data for managed area
    t_max <- max(MA_Ov_Stats$LatestYear[MA_Ov_Stats$ManagedAreaName==
                                          ma_i])
    # Determines earliest recent year with available data for managed area
    t_min <- min(MA_Ov_Stats$EarliestYear[MA_Ov_Stats$ManagedAreaName==
                                            ma_i])
    # Determines how many years of data are present
    t <- t_max-t_min
    
    # Creates break intervals for plots based on number of years of data
    if(t>=30){
      # Set breaks to every 10 years if more than 30 years of data
      brk <- -10
    }else if(t<30 & t>=10){
      # Set breaks to every 5 years if between 30 and 10 years of data
      brk <- -5
    }else if(t<10 & t>=4){
      # Set breaks to every 2 years if between 10 and 4 years of data
      brk <- -2
    }else if(t<4 & t>=2){
      # Set breaks to every year if between 4 and 2 years of data
      brk <- -1
    }else if(t<2){
      # Set breaks to every year if less than 2 years of data
      brk <- -1
      # Sets t_max to be 1 year greater and t_min to be 1 year lower
      # Forces graph to have at least 3 tick marks
      t_max <- t_max+1
      t_min <- t_min-1
    }
    # Determine range of data values for the managed area
    y_range <- max(plot_data$Mean) - min(plot_data$Mean)
    
    # Determines lower bound of y-axis based on data range. Set based on
    # relation of data range to minimum value. Designed to set lower boundary
    # to be 10% of the data range below the minimum value
    y_min <- if(min(plot_data$Mean)-(0.1*y_range)<0){
      # If 10% of the data range below the minimum value is less than 0,
      # set as 0
      y_min <- 0
    } else {
      # Otherwise set minimum bound as 10% data range below minimum value
      y_min <- min(plot_data$Mean)-(0.1*y_range)
    }
    
    # Sets upper bound of y-axis to be 10% of the data range above the
    # maximum value.
    y_max <- max(plot_data$Mean)+(0.1*y_range)
    
    # Determines what combination of groups are present for managed area
    # and subsets color and shape scheme to be used by plots.
    # Used so only group combinations present for managed area appear in
    # the legend.
    group_colors_plot <- group_colors[unique(plot_data$SpeciesGroup1)]
    group_shapes_plot <- group_shapes[unique(plot_data$SpeciesGroup1)]
    
    # Creates plot object using plot_data.
    # Data is plotted as symbols with connected lines.
    p1 <- ggplot(data=plot_data, group=as.factor(SpeciesGroup1)) +
      geom_line(aes(x=Year, y=Mean, color=as.factor(SpeciesGroup1)),
                size=0.75, alpha=1) +
      geom_point(aes(x=Year, y=Mean, fill=as.factor(SpeciesGroup1),
                     shape=as.factor(SpeciesGroup1)), size=2,
                 color="#333333", alpha=1) +
      labs(title="Coastal Wetlands Species Richness",
           subtitle=ma_i,
           x="Year", y="Richness (# of species)",
           fill="Species group", color="Species group",
           shape="Species group") +
      scale_x_continuous(limits=c(t_min-0.25, t_max+0.25),
                         breaks=seq(t_max, t_min, brk)) +
      scale_y_continuous(limits=c(y_min, y_max),
                         breaks=pretty_breaks(n=5)) +
      scale_fill_manual(values=group_colors_plot) +
      scale_color_manual(values=group_colors_plot) +
      scale_shape_manual(values=group_shapes_plot) +
      plot_theme
    # Sets file name of plot created
    outname <- paste0("CoastalWetlands_", param_file, "_", ma_abrev, ".png")
    # Saves plot as a png image
    png(paste0(out_dir, "/Figures/", outname),
        width = 8,
        height = 4,
        units = "in",
        res = 200)
    print(p1)
    dev.off()
  }
}