Fish Biomass Report

---

This code creates fish biomass report card visualizations from MERMAID fishbelt data. It is designed to be adapted for any MERMAID project by changing the configuration variables below. These visualizations are based on ideas that were co-developed with WCS Solomon Islands.

---

Configuration

Set your project tag, biomass thresholds, and optional site grouping here. These are the only variables you need to change to adapt this report to your project.

#### PROJECT CONFIGURATION ####
# Change this tag to match your MERMAID project(s)
project_name <- "SERF2_AlanFriedlander_Niue"

#### BIOMASS THRESHOLDS ####
# Stoplight thresholds for fish biomass (kg/ha)
# Yellow threshold: biomass below this is considered "poor" (red)
# Green threshold: biomass above this is considered "good" (green)
# Biomass between the two thresholds is considered "moderate" (yellow)
yellow_threshold <- 500
green_threshold <- 1000

#### SITE GROUPING (OPTIONAL) ####
# Set to TRUE to group sites by a common prefix (e.g., community name)
# Sites are assumed to be named as "[GroupName] [Number]" (e.g., "Haleta 1")
use_site_grouping <- TRUE

# If you want to manually set the order of groups (e.g., west to east),
# provide a character vector below. Set to NULL to order by mean biomass.
# Example: custom_group_order <- c("GroupA", "GroupB", "GroupC")
custom_group_order <- NULL

# If certain group names should be merged, provide a named list below.
# Each element maps original names to the merged name.
# Example: group_merges <- list("MergedName" = c("OriginalA", "OriginalB"))
# Set to NULL if no merging is needed.
group_merges <- NULL

Getting fishbelt data from MERMAID

Loads the necessary R packages and gets fishbelt sample event data from MERMAID for projects matching the configured tag.

Note — the data will only export if you are a user in a project with the specified tag. Replace the tag in the Configuration section above to access other projects for which you are listed as a user.

options(scipen = 999)

#### Load packages ####
# install.packages("remotes")
# remotes::install_github("data-mermaid/mermaidr")

library(mermaidr)
library(tidyverse)
library(plotly)
library(htmlwidgets)
library(DT)
library(ggplot2)
library(ggpubr)

#### Get data from MERMAID ####

target_project <- mermaid_search_projects(name = project_name)

mermaidr::mermaid_auth()

fish_tbl <- mermaid_get_project_data(project = target_project$id,
                                     method = "fishbelt",
                                     data = "sampleevents")

# Create a unique label for each sample event. If a site was surveyed
# more than once, append the date so each event gets its own bar.
# If there are still duplicates (same site + date but different depths),
# also append the depth.
fish_tbl <- fish_tbl %>%
  add_count(site, name = "n_per_site") %>%
  mutate(site_label = ifelse(
    test = n_per_site > 1,
    yes = paste0(site, " (", sample_date, ")"),
    no = as.character(site))) %>%
  add_count(site_label, name = "n_per_label") %>%
  mutate(site_label = ifelse(
    test = n_per_label > 1,
    yes = paste0(site, " (", sample_date, ") [", depth_avg, "m]"),
    no = site_label)) %>%
  select(-n_per_site, -n_per_label)

---

View the data

An interactive table showing the site, management rules, and average fish biomass for each sample event.

datatable(fish_tbl %>%
            select(site_label, management_rules, biomass_kgha_avg))

---

Fish biomass histogram with stoplight thresholds

An interactive histogram of site-level mean fish biomass, color-coded by the stoplight thresholds set in the Configuration section. Dotted vertical lines mark the threshold boundaries.

This is based on code used for the MERMAID Explore visualizations.

biomass <- fish_tbl$biomass_kgha_avg

# Truncate high biomass values so the histogram x-axis stays readable.
# Any survey above the cutoff is placed into a single "cutoff+" bin.
biomass_cutoff <- 5000
biomass_truncated <- ifelse(biomass > biomass_cutoff, biomass_cutoff + 1, biomass)

# Build x-axis tick labels: 0, 1000, 2000, ... up to "cutoff+"
tick_vals <- seq(0, biomass_cutoff, by = 1000)
tick_labels <- as.character(tick_vals)
tick_labels[length(tick_labels)] <- paste0(biomass_cutoff, "+")

# Split biomass into the three stoplight bins
biomass_red <- biomass_truncated[biomass_truncated < yellow_threshold]
biomass_yellow <- biomass_truncated[biomass_truncated >= yellow_threshold &
                                    biomass_truncated <= green_threshold]
biomass_green <- biomass_truncated[biomass_truncated > green_threshold]

# Build the histogram, only adding traces for non-empty bins
fish_histogram <- plot_ly()

if (length(biomass_red) > 0) {
  fish_histogram <- fish_histogram %>%
    add_histogram(
      x = biomass_red,
      xbins = list(start = 0, size = 100),
      marker = list(color = "#d13823"),
      name = paste0("< ", yellow_threshold),
      showlegend = FALSE
    )
}

if (length(biomass_yellow) > 0) {
  fish_histogram <- fish_histogram %>%
    add_histogram(
      x = biomass_yellow,
      xbins = list(start = 0, size = 100),
      marker = list(color = "#f3a224"),
      name = paste0(yellow_threshold, "–", green_threshold),
      showlegend = FALSE
    )
}

if (length(biomass_green) > 0) {
  fish_histogram <- fish_histogram %>%
    add_histogram(
      x = biomass_green,
      xbins = list(start = 0, size = 100),
      marker = list(color = "#277d1d"),
      name = paste0("> ", green_threshold),
      showlegend = FALSE
    )
}

fish_histogram <- fish_histogram %>%
  config(
    displayModeBar = TRUE,
    modeBarButtonsToRemove = c(
      "zoom", "pan", "select", "zoomIn", "zoomOut",
      "autoScale", "resetScale", "lasso2d",
      "hoverClosestCartesian", "hoverCompareCartesian"
    ),
    displaylogo = FALSE
  ) %>%
  layout(
    barmode = "overlay",
    bargap = 0.1,
    shapes = list(
      list(
        type = "line",
        x0 = yellow_threshold, x1 = yellow_threshold,
        y0 = 0, y1 = 1, yref = "paper",
        line = list(color = "black", dash = "dot")
      ),
      list(
        type = "line",
        x0 = green_threshold, x1 = green_threshold,
        y0 = 0, y1 = 1, yref = "paper",
        line = list(color = "black", dash = "dot")
      )
    ),
    xaxis = list(
      title = "Fish biomass (kg/ha)",
      linecolor = "black",
      linewidth = 2,
      range = c(0, biomass_cutoff + 200),
      tickmode = "array",
      tickvals = tick_vals,
      ticktext = tick_labels
    ),
    yaxis = list(
      title = "Number of surveys",
      linecolor = "black",
      linewidth = 2
    ),
    annotations = list(
      list(
        x = 0, y = 1.15,
        #This next line is the title and can be changed if desired
        text = paste0(project_name, " Fish Biomass"), 
        showarrow = FALSE,
        xref = "paper", yref = "paper",
        xanchor = "left", yanchor = "top",
        font = list(size = 20)
      ),
      list(
        x = 0, y = 1.08,
        text = paste0(length(biomass), " Surveys"),
        showarrow = FALSE,
        xref = "paper", yref = "paper",
        xanchor = "left", yanchor = "top",
        font = list(size = 12)
      )
    ),
    margin = list(t = 80, b = 75)
  ) %>%
  style(hovertemplate = "Bin: %{x}<br>%{y} surveys<extra></extra>")

fish_histogram

---

Fish biomass by trophic group

Stacked horizontal bar plot

A horizontal stacked bar plot showing fish biomass at each site, broken down by trophic group. Colors match the MERMAID dashboard. Horizontal lines mark the stoplight thresholds.

This is based on code originally created for a Madagascar report card.

#### Apply site grouping ####
if (use_site_grouping) {
  fish_tbl <- fish_tbl %>%
    mutate(community = gsub("\\s*\\d+$", "", site))
  # Note: community is derived from original site name, not site_label

  # Apply group merges if specified
  if (!is.null(group_merges)) {
    for (merged_name in names(group_merges)) {
      fish_tbl <- fish_tbl %>%
        mutate(community = if_else(
          community %in% group_merges[[merged_name]],
          merged_name,
          community
        ))
    }
  }

  # Determine community ordering
  if (!is.null(custom_group_order)) {
    community_levels <- custom_group_order
  } else {
    # Order by mean biomass (descending)
    comm_order <- fish_tbl %>%
      group_by(community) %>%
      summarise(comm_biomass = mean(biomass_kgha_avg, na.rm = TRUE),
                .groups = "drop") %>%
      arrange(desc(comm_biomass))
    community_levels <- comm_order$community
  }

  # Check if grouping is meaningful: at least one community must contain

  # more than one site. If every site is its own community, faceting just
  # adds clutter with truncated strip labels.
  sites_per_community <- fish_tbl %>%
    group_by(community) %>%
    summarise(n = n(), .groups = "drop")

  grouping_useful <- any(sites_per_community$n > 1)

  if (!grouping_useful) {
    # Fall back to no grouping
    fish_tbl <- fish_tbl %>%
      mutate(community = "All sites")
    community_levels <- "All sites"
    site_levels <- fish_tbl %>%
      distinct(site_label, .keep_all = TRUE) %>%
      arrange(desc(biomass_kgha_avg)) %>%
      pull(site_label)
  } else {
    # Order sites within communities by biomass
    site_order <- fish_tbl %>%
      distinct(site_label, .keep_all = TRUE) %>%
      mutate(community = factor(community, levels = community_levels)) %>%
      arrange(community, desc(biomass_kgha_avg))

    site_levels <- site_order$site_label
  }

} else {
  # No grouping — order all sites by biomass
  fish_tbl <- fish_tbl %>%
    mutate(community = "All sites")

  community_levels <- "All sites"
  site_levels <- fish_tbl %>%
    distinct(site_label, .keep_all = TRUE) %>%
    arrange(desc(biomass_kgha_avg)) %>%
    pull(site_label)
}

#### Reshape data to long format by trophic group ####
fish_trophic_tbl <- fish_tbl %>%
  select(community, site_label, latitude, longitude,
         management_rules, biomass_kgha_avg,
         starts_with("biomass_kgha_trophic_group_avg_")) %>%
  pivot_longer(
    starts_with("biomass_kgha_trophic_group_avg_"),
    names_prefix = "biomass_kgha_trophic_group_avg_",
    names_to = "trophic_group",
    values_to = "biomass"
  ) %>%
  mutate(
    biomass = ifelse(is.na(biomass), 0, biomass),
    community = factor(community, levels = community_levels),
    site_label = factor(site_label, levels = rev(site_levels)),
    trophic_group = case_when(
      trophic_group == "planktivore" ~ "Planktivore",
      trophic_group == "herbivore_macroalgae" ~ "Herbivore (macroalgae)",
      trophic_group == "herbivore_detritivore" ~ "Herbivore (detritivore)",
      trophic_group == "invertivore_sessile" ~ "Invertivore (sessile)",
      trophic_group == "invertivore_mobile" ~ "Invertivore (mobile)",
      trophic_group == "omnivore" ~ "Omnivore",
      trophic_group == "piscivore" ~ "Piscivore",
      trophic_group == "other" ~ "Other",
      .default = trophic_group
    )
  ) %>%
  mutate(trophic_group = factor(trophic_group, levels = c(
    "Planktivore",
    "Herbivore (macroalgae)",
    "Herbivore (detritivore)",
    "Invertivore (sessile)",
    "Invertivore (mobile)",
    "Omnivore",
    "Piscivore",
    "Other"
  )))

#### Trophic group colors (matching MERMAID dashboard) ####
tg_colors <- c(
  "Planktivore" = "#bebad8",
  "Herbivore (macroalgae)" = "#659034",
  "Herbivore (detritivore)" = "#ddee96",
  "Invertivore (sessile)" = "#f1da54",
  "Invertivore (mobile)" = "#e7b16d",
  "Omnivore" = "#9ccbc1",
  "Piscivore" = "#597db4",
  "Other" = "grey"
)

max_biomass <- max(fish_trophic_tbl$biomass_kgha_avg, na.rm = TRUE)

# Dynamic axis: choose a break interval that avoids overlapping labels
axis_break <- if (max_biomass > 10000) {
  2000
} else if (max_biomass > 5000) {
  1000
} else {
  500
}

# Add enough room beyond the longest bar for the biomass text label
label_padding <- max_biomass * 0.08
axis_max <- ceiling((max_biomass + label_padding) / axis_break) * axis_break

#### Create the plot ####
trophic_biomass_plot <- ggplot(
  data = fish_trophic_tbl,
  aes(x = site_label, y = biomass, fill = trophic_group)
) +
  geom_bar(position = "stack", stat = "identity", alpha = 0.9,
           color = "black", linewidth = 0.25) +
  geom_hline(yintercept = yellow_threshold, colour = "#f3a224",
             linewidth = 1) +
  geom_hline(yintercept = green_threshold, colour = "#277d1d",
             linewidth = 1) +
  labs(x = "", y = "Fish Biomass (kg/ha)", fill = "Trophic group",
       title = paste0(project_name, " — Fish Biomass by Trophic Group")) +
  coord_flip() +
  scale_fill_manual(values = tg_colors) +
  geom_text(data = fish_trophic_tbl,
            aes(x = site_label, y = biomass_kgha_avg + (max_biomass * 0.01),
                label = paste0(round(biomass_kgha_avg, 0))),
            size = 3, hjust = 0, vjust = 0.5) +
  scale_y_continuous(
    expand = c(0, 0),
    limits = c(0, axis_max),
    breaks = seq(0, axis_max, by = axis_break),
    labels = seq(0, axis_max, by = axis_break)
  ) +
  theme_classic() +
  theme(
    axis.text = element_text(size = 11, colour = "black"),
    axis.line = element_line(colour = "black"),
    axis.ticks = element_line(colour = "black"),
    axis.title = element_text(size = 12, colour = "black"),
    plot.subtitle = element_text(colour = "black", size = 11, hjust = 0.5),
    legend.background = element_rect(fill = "white", color = NA),
    legend.position = "right",
    plot.title = element_text(colour = "black", size = 14,
                              hjust = 0.5, face = "bold"),
    legend.box.background = element_blank(),
    legend.key = element_rect(color = "black", linewidth = 0.25),
    legend.title = element_text(colour = "black", face = "bold"),
    plot.margin = unit(c(0.2, 1, 0, 0.2), "cm"),
    axis.ticks.y = element_blank(),
    axis.line.x = element_line(color = "black"),
    axis.line.y = element_blank(),
    axis.ticks.x = element_line(color = "black"),
    axis.text.y = element_text(hjust = 0.5, size = 8),
    panel.border = element_blank(),
    strip.background = element_rect(fill = "grey90", colour = NA),
    strip.text.y = element_text(face = "bold")
  )

# Only add faceting if there are multiple meaningful communities
if (length(community_levels) > 1 &&
    !identical(community_levels, "All sites")) {
  trophic_biomass_plot <- trophic_biomass_plot +
    facet_grid(
      community ~ .,
      scales = "free_y",
      space = "free_y"
    )
}

# Save the plot
ggsave("../plots/fishTrophicBiomassPlot.svg",
       plot = trophic_biomass_plot,
       width = 9, height = 10)

trophic_biomass_plot

---

Fish biomass by target fish family (optional)

This section creates a stacked bar plot by fish family instead of trophic group. Because which families are of interest will depend on the question or region, this section requires customization. Define your target fish families in the target_families vector below and provide pretty labels in family_labels.

If you do not need this plot, you can remove this section entirely.

library(colorspace)

#### CUSTOMIZE: Define your target fish families ####
# List the family names (lowercase) that are of particular interest in your region.
# For example, food fish families.
# Any families not listed here will be grouped into "Other families".
target_families <- c(
  "acanthuridae",
  "balistidae",
  "caesionidae",
  "carangidae",
  "ephippidae",
  "haemulidae",
  "holocentridae",
  "kyphosidae",
  "labridae",
  "lethrinidae",
  "lutjanidae",
  "mullidae",
  "nemipteridae",
  "pomacanthidae",
  "scaridae",
  "scombridae",
  "siganidae",
  "sphyraenidae"
)

#### CUSTOMIZE: Pretty labels for each family ####
family_labels <- c(
  acanthuridae  = "Acanthuridae (surgeonfish)",
  balistidae    = "Balistidae (triggerfish)",
  caesionidae   = "Caesionidae (fusiliers)",
  carangidae    = "Carangidae (jacks)",
  ephippidae    = "Ephippidae (spade/batfish)",
  haemulidae    = "Haemulidae (sweetlips)",
  holocentridae = "Holocentridae (soldier/squirrelfish)",
  kyphosidae    = "Kyphosidae (sea chubs)",
  labridae      = "Labridae (wrasse)",
  lethrinidae   = "Lethrinidae (emperors)",
  lutjanidae    = "Lutjanidae (snappers)",
  mullidae      = "Mullidae (goatfish)",
  nemipteridae  = "Nemipteridae (breams)",
  pomacanthidae = "Pomacanthidae (angelfish)",
  scaridae      = "Scaridae (parrotfish)",
  scombridae    = "Scombridae (mackerels/tunas)",
  siganidae     = "Siganidae (rabbitfish)",
  sphyraenidae  = "Sphyraenidae (barracudas)",
  other         = "Other families"
)

#### Reshape data by fish family ####
fish_family_tbl <- fish_tbl %>%
  select(
    community, site_label, latitude, longitude,
    biomass_kgha_avg,
    starts_with("biomass_kgha_fish_family_avg_")
  ) %>%
  pivot_longer(
    starts_with("biomass_kgha_fish_family_avg_"),
    names_to = "family",
    values_to = "biomass",
    names_prefix = "biomass_kgha_fish_family_avg_"
  ) %>%
  mutate(
    biomass = ifelse(is.na(biomass), 0, biomass),
    family = if_else(family %in% target_families, family, "other")
  ) %>%
  group_by(community, site_label, latitude, longitude,
           biomass_kgha_avg, family) %>%
  summarise(biomass = sum(biomass), .groups = "drop") %>%
  mutate(
    family = recode(family, !!!family_labels),
    community = factor(community, levels = community_levels),
    site_label = factor(site_label, levels = rev(site_levels))
  ) %>%
  mutate(
    family = factor(family),
    family = factor(
      family,
      levels = c(
        setdiff(levels(family), "Other families"),
        "Other families"
      )
    )
  )

max_biomass_fam <- max(fish_family_tbl$biomass_kgha_avg, na.rm = TRUE)

# Dynamic axis (same logic as trophic group plot)
axis_break_fam <- if (max_biomass_fam > 10000) {
  2000
} else if (max_biomass_fam > 5000) {
  1000
} else {
  500
}
label_padding_fam <- max_biomass_fam * 0.08
axis_max_fam <- ceiling((max_biomass_fam + label_padding_fam) / axis_break_fam) * axis_break_fam

#### Colors ####
family_levels_vec <- levels(fish_family_tbl$family)
n_fam <- length(family_levels_vec)
n_main <- n_fam - 1

main_cols <- qualitative_hcl(n_main, palette = "Dark 3")
family_colors <- setNames(c(main_cols, "grey40"), family_levels_vec)

#### Create the plot ####
family_biomass_plot <- ggplot(
  data = fish_family_tbl,
  aes(x = site_label, y = biomass, fill = family)
) +
  geom_bar(position = "stack", stat = "identity", alpha = 0.9,
           color = "black", linewidth = 0.25) +
  geom_hline(yintercept = yellow_threshold, colour = "#f3a224",
             linewidth = 1) +
  geom_hline(yintercept = green_threshold, colour = "#277d1d",
             linewidth = 1) +
  labs(x = "", y = "Fish Biomass (kg/ha)", fill = "Fish family",
       title = paste0(project_name, " — Fish Biomass by Fish Family")) +
  coord_flip() +
  scale_fill_manual(values = family_colors) +
  geom_text(
    aes(x = site_label, y = biomass_kgha_avg + (max_biomass_fam * 0.01),
        label = round(biomass_kgha_avg, 0)),
    size = 3, hjust = 0, vjust = 0.5
  ) +
  scale_y_continuous(
    expand = c(0, 0),
    limits = c(0, axis_max_fam),
    breaks = seq(0, axis_max_fam, by = axis_break_fam),
    labels = seq(0, axis_max_fam, by = axis_break_fam)
  ) +
  theme_classic() +
  theme(
    axis.text = element_text(size = 11, colour = "black"),
    axis.line = element_line(colour = "black"),
    axis.ticks = element_line(colour = "black"),
    axis.title = element_text(size = 12, colour = "black"),
    plot.subtitle = element_text(colour = "black", size = 11, hjust = 0.5),
    legend.background = element_rect(fill = "white", color = NA),
    legend.position = "right",
    plot.title = element_text(colour = "black", size = 14,
                              hjust = 0.5, face = "bold"),
    legend.box.background = element_blank(),
    legend.key = element_rect(color = "black", linewidth = 0.25),
    legend.title = element_text(colour = "black", face = "bold"),
    plot.margin = unit(c(0.2, 1, 0, 0.2), "cm"),
    axis.ticks.y = element_blank(),
    axis.line.x = element_line(color = "black"),
    axis.line.y = element_blank(),
    axis.ticks.x = element_line(color = "black"),
    axis.text.y = element_text(hjust = 0.5, size = 8),
    panel.border = element_blank(),
    strip.background = element_rect(fill = "grey90", colour = NA),
    strip.text.y = element_text(face = "bold")
  )

# Only add faceting if there are multiple meaningful communities
if (length(community_levels) > 1 &&
    !identical(community_levels, "All sites")) {
  family_biomass_plot <- family_biomass_plot +
    facet_grid(
      community ~ .,
      scales = "free_y",
      space = "free_y"
    )
}

# Save the plot
ggsave("../plots/fishFamilyBiomassPlot.svg",
       plot = family_biomass_plot,
       width = 9, height = 10)

family_biomass_plot