In ices-tools-dev/sfdQC: Functions to produce quality control reports and plots for ICES WGSFD VMS data and products
knitr::opts_chunk$set(comment = NA,
warning = FALSE,
message = FALSE,
error = FALSE,
echo = FALSE,
eval = TRUE,
#dev = "my_png",
#fig.ext = "png",
fig.align = 'left',
tab.align = 'left')
#options(knitr.kable.NA = '-')
#options(knitr.table.format = 'latex')
#options(digits = 4)
my_png <- function(file, width, height) {
png(file, width = width, height = height,
res = 800, units = "in", pointsize = 20)
}
library(tidyverse)
library(maps)
library(lubridate)
library(RColorBrewer)
library(viridisLite)
library(sfdQC)
library(vmstools)
grade <- function(x, dx) {
brks <- seq(floor(min(x)), ceiling(max(x)),dx)
ints <- findInterval(x, brks, all.inside = TRUE)
x <- (brks[ints] + brks[ints + 1]) / 2
return(x)
}
#Read in latest submission -->",
ICES_LE <- read.table(params$LE_filename, sep = ',', header = TRUE,
stringsAsFactors = FALSE, na.strings = 'NULL',
colClasses = c('character', 'character', 'integer', 'integer',
'character', 'character', 'character', 'numeric',
'character', 'character',
'numeric', 'numeric', 'numeric'))
ICES_VE <- read.table(params$VE_filename, sep = ',', header = TRUE,
stringsAsFactors = FALSE, na.strings = 'NULL',
colClasses = c('character', 'character', 'integer', 'integer',
'character', 'character', 'character', 'character',
'numeric', 'numeric', 'numeric', 'numeric',
'numeric', 'numeric', 'numeric', 'numeric', 'numeric'))
#Most recent submission -->
country <- ICES_VE$country[1]
ICES_VE <- cbind(ICES_VE, CSquare2LonLat(ICES_VE$c_square, degrees=0.05))
spatBound <- list(xrange = range(ICES_VE$SI_LONG, na.rm=TRUE),
yrange = range(ICES_VE$SI_LATI, na.rm=TRUE))
spatCore <- list(xrange = quantile(ICES_VE$SI_LONG, c(0.025, 0.975)),
yrange = quantile(ICES_VE$SI_LATI, c(0.025, 0.975)))
tempBound <- range(ICES_VE$year, na.rm=TRUE)
#Most recent submission -->
ICES_LE <- cbind(ICES_LE, ICESrectangle2LonLat(ICES_LE$ICES_rectangle, midpoint=TRUE))
spatBoundLog <- list(xrange = range(ICES_LE$SI_LONG, na.rm=TRUE),
yrange = range(ICES_LE$SI_LATI, na.rm=TRUE))
tempBoundLog <- range(ICES_LE$year, na.rm=TRUE)
coord <-
coord_quickmap(xlim = range(spatBound$xrange),
ylim = range(spatBound$yrange))
spatBoundasp <-
coord$aspect(list(x.range = range(spatBound$xrange),
y.range = range(spatBound$yrange)))
VMS
hasvms <- nrow(ICES_VE) > 0
haslogbook <- nrow(ICES_LE) > 0
r if (!hasvms) {"No VMS data available\n\\begin{comment}"}
knitr::opts_chunk$set(eval = hasvms)
Records by year
d <- ICES_VE %>% count(year)
d %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n)) +
theme_icesqc() +
geom_col() +
labs(x = "Year", y = "Number of records")
\newpage
Records by month
d <- ICES_VE %>% count(month, year)
d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n, fill = factor(month))) +
theme_icesqc(legend.position = "right") +
geom_col() +
scale_fill_grey() +
labs(x = "Year", y = "Number of records", fill = "Month")
\newpage
Vessel length categories
d <- ICES_VE %>% count(variable = vessel_length_category, year)
d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n, fill = factor(variable))) +
theme_icesqc(legend.position = "right") +
geom_col() +
scale_fill_grey() +
labs(x = "Year", y = "Number of records", fill = "Length category")
\newpage
Gear code
d <- ICES_VE %>% count(variable = gear_code, year)
d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n, fill = factor(variable))) +
theme_icesqc(legend.position = "right") +
geom_col() +
scale_fill_grey() +
labs(x = "Year", y = "Number of records", fill = "Gear code")
Unique DCF Level 6 codes
ICES_VE %>%
select(variable = LE_MET_level6, year) %>%
distinct() %>%
count(year) %>%
knitr::kable(booktabs = TRUE)
Top 5 DCF Level 6 codes
top5Met6 <- ICES_VE %>% count(LE_MET_level6) %>% slice(1:5) %>% pull(LE_MET_level6)
d <- ICES_VE %>%
filter(LE_MET_level6 %in% top5Met6) %>%
count(variable = LE_MET_level6, year)
d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n, fill = factor(variable))) +
theme_icesqc(legend.position = "right") +
geom_col() +
scale_fill_grey() +
labs(x = "Year", y = "Number of records", fill = "Level6")
\newpage
Average fishing speed
knitr::kable(do.call(rbind, tapply(ICES_VE$avg_fishing_speed, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = avg_fishing_speed)) +
geom_histogram() +
xlab("Average fishing speed") + ylab ("Count") +
facet_wrap( ~ year, ncol = 2) +
theme_icesqc()
\newpage
Average fishing hours
knitr::kable(do.call(rbind, tapply(ICES_VE$fishing_hours, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = fishing_hours)) +
geom_histogram() +
xlab("Average fishing hours") + ylab("Count") +
facet_wrap( ~ year, ncol = 2) +
theme_icesqc()
\newpage
Average length
knitr::kable(do.call(rbind, tapply(ICES_VE$avg_oal, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = avg_oal))+
geom_histogram() +
xlab("Average vessel length") + ylab("Count") +
facet_wrap( ~ year, ncol = 2) +
theme_icesqc()
\newpage
Average kW
knitr::kable(do.call(rbind, tapply(ICES_VE$avg_kw, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = avg_kw))+
geom_histogram() +
xlab("Average kW") + ylab ("Count") +
facet_wrap( ~ year, ncol = 2) +
theme_icesqc()
\newpage
Average kW-hours
knitr::kable(do.call(rbind, tapply(ICES_VE$kw_fishinghours, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = kw_fishinghours))+
geom_histogram() +
xlab("Average kW-hours") + ylab("Count") +
facet_wrap(~ year, ncol = 2) +
theme_icesqc()
\newpage
Landings by gear by year
ps <- gear_splits(ICES_VE$totweight, data = ICES_VE, "kg landed", year_groups = 2, gear_groups = 4, func = sum)
ps$table
for (p in ps$plots) print(p)
\newpage
Mean landing per kW fishing hours
ps <- gear_splits(with(ICES_VE, totweight/kw_fishinghours), data = ICES_VE, "kg/kWh", gear_groups = 4, func = median)
ps$table
for (p in ps$plots) print(p)
\newpage
Value by gear by year
ps <- gear_splits(ICES_VE$totvalue, data = ICES_VE, "EUR landed", gear_groups = 4, func = sum)
ps$table
if(any(!is.na(ICES_VE$totvalue))) for (p in ps$plots) print(p)
\newpage
Median value per KW fishing hours by year
ps <- gear_splits(with(ICES_VE, totvalue/kw_fishinghours), data = ICES_VE, "EUR/kWh", gear_groups = 4, func = median)
ps$table
if(any(!is.na(ICES_VE$totvalue))) for (p in ps$plots) print(p)
\newpage
Average price
ps <- gear_splits(with(ICES_VE, totvalue/totweight), data = ICES_VE, "Mean price (EUR/kg)", gear_groups = 4, func = median)
ps$table
if(any(!is.na(ICES_VE$totvalue))) for (p in ps$plots) print(p)
r if (!hasvms) {"\\end{comment}"}
\newpage
r if (any(is.na(ICES_VE$SI_LONG))) {"\\newpage"}
Invalid Rectangles
if (any(is.na(ICES_VE$SI_LONG))) {
knitr::kable(table(`ICES Rectangle` = ICES_VE$ICES_rectangle[is.na(ICES_VE$SI_LONG)],
Year = ICES_VE$year[is.na(ICES_VE$SI_LONG)]), booktabs = TRUE)
} else {
x <- "There were no invalid Statistical Rectangles reported"
attr(x, "format") <- "markdown"
attr(x, "class") <- "knit_asis"
x
}
\newpage
Spatial extent
Range by year
dat2tab <- as.matrix(aggregate(ICES_VE$SI_LONG,by=list(ICES_VE$year),FUN=range))
dat2tab <- cbind(dat2tab,as.matrix(aggregate(ICES_VE$SI_LATI,by=list(ICES_VE$year),FUN=range)[,-1]))
colnames(dat2tab) <- c("Year", "min lon", "max lon", "min lat", "max lat")
knitr::kable(dat2tab, booktabs = TRUE)
Full area
d <-
ICES_VE %>%
as_tibble() %>%
select(long = SI_LONG, lat = SI_LATI, year) %>%
mutate(long = grade(long, 0.5),
lat = grade(lat, 0.5)) %>%
distinct()
map_data("world",
xlim = spatBound$xrange,
ylim = spatBound$yrange) %>%
ggplot(aes(long, lat)) +
theme_void() +
geom_polygon(aes(group = group), fill = "grey90") +
geom_raster(data = d, fill = "red") +
coord_quickmap(xlim = spatBound$xrange,
ylim = spatBound$yrange) +
geom_path(aes(group = group), colour = "grey", lwd = 0.1) +
facet_wrap(~ year, ncol = 2) +
labs(x = NULL, y = NULL)
Core area
d <-
ICES_VE %>%
as_tibble() %>%
select(long = SI_LONG, lat = SI_LATI, year) %>%
distinct()
map_data("world",
xlim = spatCore$xrange,
ylim = spatCore$yrange) %>%
ggplot(aes(long, lat)) +
theme_void() +
geom_polygon(aes(group = group), fill = "grey90") +
geom_raster(data = d, fill = "red") +
coord_quickmap(xlim = spatCore$xrange,
ylim = spatCore$yrange) +
geom_path(aes(group = group), colour = "grey", lwd = 0.1) +
facet_wrap(~ year, ncol = 2) +
labs(x = NULL, y = NULL)
\newpage
3 most dominant gear
# get 3 main gears
gear_table <- with(ICES_VE, aggregate(fishing_hours, list(gear_code = gear_code), sum, na.rm = TRUE))
gear_table <- gear_table[order(gear_table$x, decreasing = TRUE),]
top3gears <- gear_table$gear_code[1:pmin(3, nrow(gear_table))]
# split by gear
idx <- which(ICES_VE$gear_code %in% top3gears)
coordGrd <- unique(ICES_VE[idx, c("SI_LONG", "SI_LATI", "year", "c_square", "gear_code")])
# create a fortied polygon of csquares : 0.39 secs
polVMS <- make_polVMS(coordGrd)
polVMS$year <- rep(coordGrd$year, each = 5)
polVMS$c_square <- rep(coordGrd$c_square,each=5)
polVMS$gear_code<- rep(coordGrd$gear_code,each=5)
# aggregate fishing days
dat <- with(ICES_VE[idx,],
aggregate(fishing_hours/24,
by = list(year = year, c_square = c_square, gear_code = gear_code),
FUN = sum, na.rm = TRUE))
dat <- dplyr::rename(dat, fishing_days = x)
# legend calculations
steps <- ceiling(max(dat$fishing_days, na.rm = TRUE)/250)
breaks <- unique(c(-1, 0, steps * c(1, 2.5, 5, 10, 25, 50, 100, 250)))
legval <- paste(breaks[-length(breaks)], "<=", breaks[-1L])
legval[1] <- "0"
palette <- c("white", RColorBrewer::brewer.pal(length(breaks)-2, "YlOrRd"))
dat$colgrp <- as.numeric(cut(dat$fishing_days, breaks = breaks))
dat$cols <- palette[dat$colgrp]
# join dat onto polVMS (much faster than left_join or merge)
# make joining keys
a_by <- apply(dat[c("year","c_square","gear_code")], 1, paste, collapse = ".")
b_by <- apply(polVMS[c("year","c_square","gear_code")], 1, paste, collapse = ".")
ib_by <- as.integer(factor(b_by, levels = a_by))
polVMS$cols <- dat$cols[ib_by]
# do plots
for (gear in top3gears) {
p <-
spatialplot(polVMS[polVMS$gear_code == gear,]) +
scale_fill_manual(values = rev(palette), labels = rev(legval)) +
guides(fill = guide_legend(title = paste("Days@Sea", gear))) +
facet_wrap(~ year, ncol = 2) +
theme_icesqc(legend.position = "top")
print(p)
}
top3gears <-
ICES_VE %>%
group_by(gear_code) %>%
summarise(fishing_hours = sum(fishing_hours, na.rm = TRUE)) %>%
arrange(desc(fishing_hours)) %>%
slice(1:3) %>%
pull(gear_code)
d <-
ICES_VE %>%
as_tibble() %>%
filter(gear_code %in% top3gears) %>%
select(year, long = SI_LONG, lat = SI_LATI, gear_code, fishing_hours) %>%
group_by(year, gear_code, long, lat) %>%
summarise(value = sum(fishing_hours, na.rm = TRUE) / 24)
#d %>% group_by(year) %>% summarise(min = min(value), max = max(value))
breaks <- c(0, 1, 2.5, 5, 10, 25, 50, 100, 250)
labels <- cut(0:250, breaks = breaks) %>% unique() %>% as.character()
#x <- cut(d$value, breaks = breaks)
# NEEDS FIXING: The top range is not shown in the legends
for(i in 1:length(top3gears)) {
p <-
map_data("world",
xlim = spatCore$xrange,
ylim = spatCore$yrange) %>%
ggplot(aes(long, lat)) +
theme_void() +
geom_polygon(aes(group = group), fill = "grey90") +
geom_raster(data = d %>% filter(gear_code %in% top3gears[i]),
aes(long, lat, fill = value)) +
scale_fill_gradientn(colours = inferno(256,
# drop the light yellow colour
begin = 0.95,
end = 0),
guide = guide_legend(),
breaks = breaks,
labels = labels) +
coord_quickmap(xlim = spatCore$xrange,
ylim = spatCore$yrange) +
geom_path(aes(group = group), colour = "grey", lwd = 0.1) +
facet_wrap(~ year, ncol = 2) +
labs(x = NULL, y = NULL, fill = paste("Days@Sea", top3gears[i])) +
theme(legend.position = "top")
print(p)
}
\newpage
Spatial distribution of effort
coordGrd <- unique(ICES_VE[,c("SI_LONG","SI_LATI","year","c_square")])
# make fortified DF for csquares
polVMS <- make_polVMS(coordGrd)
polVMS$year <- rep(coordGrd$year, each = 5)
polVMS$c_square <- rep(coordGrd$c_square, each = 5)
# aggregate fishing days
dat <- with(ICES_VE,
aggregate(fishing_hours/24,
by = list(year = year, c_square = c_square),
FUN = sum, na.rm = TRUE))
dat <- dplyr::rename(dat, fishing_days = x)
# legend calculations
steps <- ceiling(max(dat$fishing_days, na.rm = TRUE)/250)
breaks <- unique(c(-1, 0, steps * c(1, 2.5, 5, 10, 25, 50, 100, 250)))
legval <- paste(breaks[-length(breaks)], "<=", breaks[-1L])
legval[1] <- "0"
palette <- c("white", RColorBrewer::brewer.pal(length(breaks)-2, "YlOrRd"))
dat$colgrp <- as.numeric(cut(dat$fishing_days, breaks = breaks))
dat$cols <- palette[dat$colgrp]
# join dat onto polVMS (much faster than left_join or merge)
# make joining keys
a_by <- apply(dat[c("year","c_square")], 1, paste, collapse = ".")
b_by <- apply(polVMS[c("year","c_square")], 1, paste, collapse = ".")
ib_by <- as.integer(factor(b_by, levels = a_by))
polVMS$cols <- dat$cols[ib_by]
# do plots
spatialplot(polVMS) +
scale_fill_manual(values = rev(palette), labels = rev(legval)) +
guides(fill = guide_legend(title = "Days@Sea")) +
facet_wrap(~ year, ncol = 2) +
theme_icesqc(legend.position = "top")
d <-
ICES_VE %>%
as_tibble() %>%
select(year, long = SI_LONG, lat = SI_LATI, fishing_hours) %>%
group_by(year, long, lat) %>%
summarise(value = sum(fishing_hours, na.rm = TRUE) / 24)
#d %>% group_by(year) %>% summarise(min = min(value), max = max(value))
breaks <- c(0, 1, 2.5, 5, 10, 25, 50, 100, 250)
labels <- cut(0:250, breaks = breaks) %>% unique() %>% as.character()
#x <- cut(d$value, breaks = breaks)
# Note: the top range is not shown in the legends
map_data("world",
xlim = spatCore$xrange,
ylim = spatCore$yrange) %>%
ggplot(aes(long, lat)) +
theme_void() +
geom_polygon(aes(group = group), fill = "grey90") +
geom_raster(data = d, aes(long, lat, fill = value)) +
scale_fill_gradientn(colours = inferno(256,
# drop the light yellow colour
begin = 0.95,
end = 0),
guide = guide_legend(),
breaks = breaks,
labels = labels) +
coord_quickmap(xlim = spatCore$xrange,
ylim = spatCore$yrange) +
geom_path(aes(group = group), colour = "grey", lwd = 0.1) +
facet_wrap(~ year, ncol = 2) +
labs(x = NULL, y = NULL, fill = "Days@Sea") +
theme(legend.position = "top")
\newpage
Spatial difference of effort r tempBound[1]:r (tempBound[2]-1) vs r tempBound[2]
base <- with(ICES_VE,
aggregate(fishing_hours,
by = list(c_square = c_square, year = year),
FUN = sum, na.rm = TRUE))
base <- dplyr::rename(base, fishing_hours = x)
# calculate total fishing hours for recent year
recent <- base[base$year == tempBound[2],]
# calculate median of the total fishing hours per square for historical years
base <- with(base[base$year < tempBound[2],],
aggregate(fishing_hours,
by = list(c_square = c_square),
FUN = median, na.rm = TRUE))
base <- dplyr::rename(base, fishing_hours_median = x)
# join
dat2plot <- dplyr::full_join(base,
recent[,c("c_square","fishing_hours")])
# set NAs to zero
dat2plot$fishing_hours_median[is.na(dat2plot$fishing_hours_median)] <- 0
dat2plot$fishing_hours[is.na(dat2plot$fishing_hours)] <- 0
# calculate ratio (with exceptions for zeros)
dat2plot$ratio <- 1/with(dat2plot, pmax(fishing_hours, 1e-9) / pmax(fishing_hours_median, 1e-9))
# add back in lat and long
dat2plot <- cbind(dat2plot,
vmstools::CSquare2LonLat(dat2plot$c_square, degrees = 0.05))
# make 'fortified' data frame
polVMS <- make_polVMS(dat2plot)
polVMS$c_square <- rep(dat2plot$c_square, each = 5)
# set up legend
# This is not what the legend says.... there is no +/-5% break!
breaks <- c(-Inf, 1/2, 1/1.5, 1/1.25, 1/1.05, 1, 1.05, 1.25, 1.5, 2, Inf)
legval <- c("historic >>","historic> +100%","historic> +50%","historic> +25%",
"+/-5%",
"recent> +5%","recent> +25%","recent> +50%","recent> +100%","recent >>")
palette <- RColorBrewer::brewer.pal(length(breaks)-1,"RdYlBu") #colour for fishing intensity
dat2plot$colgrp <- as.numeric(cut(dat2plot$ratio, breaks = breaks))
dat2plot$cols <- palette[dat2plot$colgrp]
# join dat onto polVMS (much faster than left_join or merge)
# make joining keys
a_by <- dat2plot$c_square
b_by <- polVMS$c_square
ib_by <- as.integer(factor(b_by, levels = a_by))
polVMS$cols <- dat2plot$cols[ib_by]
# NOTE dark blue here just means that nothing was observered before, and now there is something...
# - could be 1 hour of fishing
grps <- sort(unique(dat2plot$colgrp))
spatialplot(polVMS) +
guides(fill = guide_legend(title = "Days@Sea")) +
scale_fill_manual(values = rev(palette)[grps], labels = legval[grps]) +
theme_icesqc(legend.position = "right")
\newpage
Spatial difference of effort r tempBound[2]-1 vs r tempBound[2]
base <- with(ICES_VE,
aggregate(fishing_hours,
by = list(c_square = c_square, year = year),
FUN = sum, na.rm = TRUE))
base <- dplyr::rename(base, fishing_hours = x)
# calculate total fishing hours for recent year
recent <- base[base$year == tempBound[2],]
# previous year
base <- base[base$year == tempBound[2]-1,]
base <- dplyr::rename(base, fishing_hours_median = fishing_hours)
# join
dat2plot <- dplyr::full_join(base,
recent[,c("c_square","fishing_hours")])
# set NAs to zero
dat2plot$fishing_hours_median[is.na(dat2plot$fishing_hours_median)] <- 0
dat2plot$fishing_hours[is.na(dat2plot$fishing_hours)] <- 0
# calculate ratio (with exceptions for zeros)
dat2plot$ratio <- 1/with(dat2plot, pmax(fishing_hours, 1e-9) / pmax(fishing_hours_median, 1e-9))
# add back in lat and long
dat2plot <- cbind(dat2plot,
vmstools::CSquare2LonLat(dat2plot$c_square, degrees = 0.05))
# make 'fortified' data frame
polVMS <- make_polVMS(dat2plot)
polVMS$c_square <- rep(dat2plot$c_square, each = 5)
# set up legend
# This is not what the legend says.... there is no +/-5% break!
breaks <- c(-Inf, 1/2, 1/1.5, 1/1.25, 1/1.05, 1, 1.05, 1.25, 1.5, 2, Inf)
legval <- c("historic >>","historic> +100%","historic> +50%","historic> +25%",
"+/-5%",
"recent> +5%","recent> +25%","recent> +50%","recent> +100%","recent >>")
palette <- RColorBrewer::brewer.pal(length(breaks)-1,"RdYlBu") #colour for fishing intensity
dat2plot$colgrp <- as.numeric(cut(dat2plot$ratio, breaks = breaks))
dat2plot$cols <- palette[dat2plot$colgrp]
# join dat onto polVMS (much faster than left_join or merge)
# make joining keys
a_by <- dat2plot$c_square
b_by <- polVMS$c_square
ib_by <- as.integer(factor(b_by, levels = a_by))
polVMS$cols <- dat2plot$cols[ib_by]
# NOTE dark blue here just means that nothing was observered before, and now there is something...
# - could be 1 hour of fishing
grps <- sort(unique(dat2plot$colgrp))
spatialplot(polVMS) +
guides(fill = guide_legend(title = "Days@Sea")) +
scale_fill_manual(values = rev(palette)[grps], labels = legval[grps]) +
theme_icesqc(legend.position = "right")
Logbooks
r if (!haslogbook) {"No Logbook data available\n\\begin{comment}"}
knitr::opts_chunk$set(eval = haslogbook)
Records by year
d <- ICES_LE %>% count(year)
d %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n)) +
theme_icesqc() +
geom_col() +
labs(x = "Year", y = "Number of records")
\newpage
Records by month
d <- ICES_LE %>% count(month, year)
d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n, fill = factor(month))) +
theme_icesqc(legend.position = "right") +
geom_col() +
scale_fill_grey() +
labs(x = "Year", y = "Number of records", fill = "Month")
\newpage
Vessel length categories
d <- ICES_LE %>% count(variable = vessel_length_category, year)
d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n, fill = factor(variable))) +
theme_icesqc(legend.position = "right") +
geom_col() +
scale_fill_grey() +
labs(x = "Year", y = "Number of records", fill = "Length category")
\newpage
Gear code
d <- ICES_LE %>% count(variable = gear_code, year)
d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n, fill = factor(variable))) +
theme_icesqc(legend.position = "right") +
geom_col() +
scale_fill_grey() +
labs(x = "Year", y = "Number of records", fill = "Gear code")
\newpage
Unique DCF Level 6 codes
knitr::kable(plyr::count(unique(ICES_LE[,c("LE_MET_level6","year")]),"year"), booktabs = TRUE)
Top 5 DCF Level 6 codes
top5Met6 <- ICES_LE %>% count(LE_MET_level6) %>% slice(1:5) %>% pull(LE_MET_level6)
d <- ICES_LE %>%
filter(LE_MET_level6 %in% top5Met6) %>%
count(variable = LE_MET_level6, year)
d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE)
d %>%
ggplot(aes(factor(year), n, fill = factor(variable))) +
theme_icesqc(legend.position = "right") +
geom_col() +
scale_fill_grey() +
labs(x = "Year", y = "Number of records", fill = "Level6")
Average fishing days
knitr::kable(do.call(rbind, tapply(ICES_LE$fishing_days, ICES_LE$year, summary)), booktabs = TRUE)
ggplot(ICES_LE, aes(x = fishing_days)) +
geom_histogram() +
xlab("Average fishing days") + ylab("Count") +
facet_wrap( ~ year, ncol = 2) +
theme_icesqc()
\newpage
Average kW-hours
knitr::kable(do.call(rbind, tapply(ICES_LE$kw_fishing_days, ICES_LE$year, summary)), booktabs = TRUE)
ggplot(ICES_LE, aes(x = kw_fishing_days))+
geom_histogram() +
xlab("Average KW-days") + ylab("Count") +
facet_wrap(~ year, ncol = 2) +
theme_icesqc()
\newpage
Landings by gear by year
ps <- gear_splits(ICES_LE$totweight, data = ICES_LE, "kg landed", gear_groups = 4, func = sum)
ps$table
for (p in ps$plots) print(p)
Mean landing per kW fishing days
ps <- gear_splits(with(ICES_LE, totweight/kw_fishing_days), data = ICES_LE, "kg/kWd", gear_groups = 4, func = median)
ps$table
for (p in ps$plots) print(p)
\newpage
Value by gear by year
ps <- gear_splits(ICES_LE$totvalue, data = ICES_LE, "EUR landed", gear_groups = 4, func = sum)
ps$table
if(any(!is.na(ICES_LE$totvalue))) for (p in ps$plots) print(p)
\newpage
Median value per KW fishing days by year
ps <- gear_splits(with(ICES_LE, totvalue/kw_fishing_days), data = ICES_LE, "EUR/kWh", gear_groups = 4, func = median)
ps$table
if(any(!is.na(ICES_LE$totvalue))) for (p in ps$plots) print(p)
\newpage
Average price
ps <- gear_splits(with(ICES_LE, totvalue/totweight), data = ICES_LE, "Mean price (EUR/kg)", gear_groups = 4, func = median)
ps$table
if(any(!is.na(ICES_LE$totvalue))) for (p in ps$plots) print(p)
r if (any(is.na(ICES_LE$SI_LONG))) {"\\newpage"}
Records assigned to invalid Statistical Rectangles
if (any(is.na(ICES_LE$SI_LONG))) {
knitr::kable(table(`ICES Rectangle` = ICES_LE$ICES_rectangle[is.na(ICES_LE$SI_LONG)],
Year = ICES_LE$year[is.na(ICES_LE$SI_LONG)]), booktabs = TRUE)
} else {
x <- "There were no invalid Statistical Rectangles reported"
attr(x, "format") <- "markdown"
attr(x, "class") <- "knit_asis"
x
}
\newpage
Area extent of data submitted by year:
dat2tab <- as.matrix(aggregate(ICES_LE$SI_LONG, by = list(year = ICES_LE$year), FUN=range, na.rm=TRUE))
dat2tab <- cbind(dat2tab,as.matrix(aggregate(ICES_LE$SI_LATI,by=list(ICES_LE$year), FUN=range, na.rm=TRUE)[,-1]))
colnames(dat2tab) <- c("Year","min lon","max lon","min lat","max lat")
knitr::kable(dat2tab, booktabs = TRUE)
Area for which data has been submitted:
d <-
ICES_LE %>%
select(long = SI_LONG, lat = SI_LATI, year) %>%
distinct()
map_data("world",
xlim = spatBoundLog$xrange,
ylim = spatBoundLog$yrange) %>%
ggplot(aes(long, lat)) +
theme_void() +
geom_polygon(aes(group = group), fill = "grey90") +
geom_raster(data = d, fill = "red") +
coord_quickmap(xlim = spatBoundLog$xrange,
ylim = spatBoundLog$yrange) +
geom_path(aes(group = group), colour = "black", lwd = 0.1) +
facet_wrap(~ year, ncol = 2) +
labs(x = NULL, y = NULL)
map_data("world",
xlim = spatCore$xrange,
ylim = spatCore$yrange) %>%
ggplot(aes(long, lat)) +
theme_void() +
geom_polygon(aes(group = group), fill = "grey90") +
geom_raster(data = d, fill = "red") +
coord_quickmap(xlim = spatCore$xrange,
ylim = spatCore$yrange) +
geom_path(aes(group = group), colour = "black", lwd = 0.1) +
facet_wrap(~ year, ncol = 2) +
labs(x = NULL, y = NULL)
\newpage
Frequency of VMS enabled category
knitr::kable(table(ICES_LE$vms_enabled,ICES_LE$vessel_length_category), booktabs = TRUE)
\newpage
Spatial distribution of effort by year:
coordGrd <- unique(ICES_LE[,c("SI_LONG","SI_LATI","year","ICES_rectangle")])
polRect <- make_polVMS(coordGrd, resolution = 1)
polRect$year <- rep(coordGrd$year,each=5)
polRect$ICES_rectangle <- rep(coordGrd$ICES_rectangle, each = 5)
#TIDY ---
dat <- aggregate(ICES_LE$fishing_days,by=list(ICES_LE$year,ICES_LE$ICES_rectangle),FUN=sum,na.rm=T)
steps <- ceiling(max(dat$x,na.rm=T)/250)
cutbreaksval <- unique(c(-1,0,steps*c(1,2.5,5,10,25,50,100,250)))
legval <- outer(ac(cutbreaksval),ac(cutbreaksval),function(x,y){return(paste(x,"<=",y))})
legval <- c("0",diag(legval[-1,-c(1,2)]))
palette <- c("white", brewer.pal(length(cutbreaksval)-2,"YlOrRd"))
cols <- palette[cut(dat$x,breaks=cutbreaksval)]
cols <- cbind(cols,id=1:length(cols),ICES_rectangle=dat$Group.2,year=dat$Group.1)
polRect <- merge(polRect,cols,by=c("ICES_rectangle","year"))
# ----
spatialplot(polRect) +
guides(fill = guide_legend(title = "Days@Sea")) +
scale_fill_manual(values = rev(palette), labels = rev(legval)) +
facet_wrap(~ year, ncol = 2) +
theme_icesqc(legend.position = "top")
d <-
ICES_LE %>%
as_tibble() %>%
select(year, long = SI_LONG, lat = SI_LATI, value = fishing_days) %>%
group_by(year, long, lat) %>%
summarise(value = sum(value, na.rm = TRUE))
steps <- ceiling(max(d$value, na.rm = TRUE) / 250)
breaks <- steps * c(0, 1, 2.5, 5, 10, 25, 50, 100, 250)
labels <- cut(0:(steps * 250), breaks = breaks) %>% unique() %>% as.character()
# Note: the top range is not shown in the legends
map_data("world",
xlim = spatCore$xrange,
ylim = spatCore$yrange) %>%
ggplot(aes(long, lat)) +
theme_void() +
geom_polygon(aes(group = group), fill = "grey90") +
geom_raster(data = d, aes(long, lat, fill = value)) +
scale_fill_gradientn(colours = inferno(256,
# drop the light yellow colour
begin = 0.95,
end = 0),
guide = guide_legend(),
breaks = breaks,
labels = labels) +
coord_quickmap(xlim = spatCore$xrange,
ylim = spatCore$yrange) +
geom_path(aes(group = group), colour = "grey", lwd = 0.1) +
facet_wrap(~ year, ncol = 2) +
labs(x = NULL, y = NULL, fill = "Days@Sea") +
theme(legend.position = "top")
\newpage
Spatial difference of effort r tempBound[1]:r (tempBound[2]-1) vs r tempBound[2]
base <- with(ICES_LE,
aggregate(fishing_days,
by = list(ICES_rectangle = ICES_rectangle, year = year),
FUN = sum, na.rm = TRUE))
base <- dplyr::rename(base, fishing_days = x)
# calculate total fishing hours for recent year
recent <- base[base$year == tempBound[2],]
# calculate median of the total fishing hours per square for historical years
base <- with(base[base$year < tempBound[2],],
aggregate(fishing_days,
by = list(ICES_rectangle = ICES_rectangle),
FUN = median, na.rm = TRUE))
base <- dplyr::rename(base, fishing_days_median = x)
# join
dat2plot <- dplyr::full_join(base,
recent[,c("ICES_rectangle","fishing_days")])
# set NAs to zero
dat2plot$fishing_days_median[is.na(dat2plot$fishing_days_median)] <- 0
dat2plot$fishing_days[is.na(dat2plot$fishing_days)] <- 0
# calculate ratio (with exceptions for zeros)
dat2plot$ratio <- 1/with(dat2plot, pmax(fishing_days, 1e-9) / pmax(fishing_days_median, 1e-9))
# add back in lat and long
dat2plot <- cbind(dat2plot,
vmstools::ICESrectangle2LonLat(dat2plot$ICES_rectangle, midpoint = TRUE))
# make 'fortified' data frame
polRect <- make_polVMS(dat2plot, resolution = 1)
polRect$ICES_rectangle <- rep(dat2plot$ICES_rectangle, each = 5)
## tidy ---
breaks <- rev(c(1e-10,0.5,2/3,0.8,0.952381,1,1.05,1.25,1.5,2,1e10))
legval <- c("historic >>","historic> +100%","historic> +50%","historic> +25%",
"+/-5%",
"recent> +5%","recent> +25%","recent> +50%","recent> +100%","recent >>")
palette <- brewer.pal(length(cutbreaksval)-1,"RdYlBu")
colgrp <- as.numeric(cut(dat2plot$ratio, breaks = breaks))
cols <- cbind(cols = palette[colgrp], ICES_rectangle = dat2plot$ICES_rectangle)
polRect <- merge(polRect, cols, by=c("ICES_rectangle"))
# ---
spatialplot(polRect) +
guides(fill=guide_legend(title="Days@Sea")) +
scale_fill_manual(values = rev(palette), labels = legval) +
theme_icesqc(legend.position = "right")
\newpage
Spatial difference of effort r tempBound[2]-1 vs r tempBound[2]
base <- with(ICES_LE,
aggregate(fishing_days,
by = list(ICES_rectangle = ICES_rectangle, year = year),
FUN = sum, na.rm = TRUE))
base <- dplyr::rename(base, fishing_days = x)
# calculate total fishing hours for recent year
recent <- base[base$year == tempBound[2],]
# previous year
base <- base[base$year == tempBound[2]-1,]
base <- dplyr::rename(base, fishing_days_median = fishing_days)
# join
dat2plot <- dplyr::full_join(base,
recent[,c("ICES_rectangle","fishing_days")])
# set NAs to zero
dat2plot$fishing_days_median[is.na(dat2plot$fishing_days_median)] <- 0
dat2plot$fishing_days[is.na(dat2plot$fishing_days)] <- 0
# calculate ratio (with exceptions for zeros)
dat2plot$ratio <- 1/with(dat2plot, pmax(fishing_days, 1e-9) / pmax(fishing_days_median, 1e-9))
# add back in lat and long
dat2plot <- cbind(dat2plot,
vmstools::ICESrectangle2LonLat(dat2plot$ICES_rectangle, midpoint = TRUE))
# make 'fortified' data frame
polRect <- make_polVMS(dat2plot, resolution = 1)
polRect$ICES_rectangle <- rep(dat2plot$ICES_rectangle, each = 5)
## tidy ---
breaks <- rev(c(1e-10,0.5,2/3,0.8,0.952381,1,1.05,1.25,1.5,2,1e10))
legval <- c("historic >>","historic> +100%","historic> +50%","historic> +25%",
"+/-5%",
"recent> +5%","recent> +25%","recent> +50%","recent> +100%","recent >>")
palette <- brewer.pal(length(cutbreaksval)-1,"RdYlBu")
colgrp <- as.numeric(cut(dat2plot$ratio, breaks = breaks))
cols <- cbind(cols = palette[colgrp], ICES_rectangle = dat2plot$ICES_rectangle)
polRect <- merge(polRect, cols, by=c("ICES_rectangle"))
# ---
spatialplot(polRect) +
guides(fill=guide_legend(title="Days@Sea")) +
scale_fill_manual(values = rev(palette), labels = legval) +
theme_icesqc(legend.position = "right")
\newpage
Relationship fishing days and total weight
ggplot(ICES_LE[ICES_LE$year == tempBoundLog[2],],
aes(x = fishing_days, y = totweight)) +
geom_point() +
facet_wrap(~ gear_code, ncol = 3, scale = "free") +
xlab("Fishing days") + ylab ("Total weight") +
theme_icesqc()
\newpage
VMS vs Logbooks
Comparison of Metier level 6 reporting between logbook and VMS
ledat <- with(ICES_LE, table(LE_MET_level6, year))
vedat <- with(ICES_VE, table(LE_MET_level6, year))
dat2tab <-
rbind(
cbind(as.data.frame.table(ledat), data = "LE (records)"),
cbind(as.data.frame.table(vedat), data = "VE (pings)"))
tab <- with(dat2tab, tapply(Freq, list(LE_MET_level6, data, year), sum))
tab[tab == 0] <- NA
for (i in dimnames(tab)[[3]]) {
x <- tab[,,i]
x <- x[apply(!is.na(x), 1, any),]
cat(knitr::kable(cbind(x, year = i), booktabs = TRUE), sep = "\n")
cat("\n")
}
r if (!haslogbook) {"\\end{comment}"}
ices-tools-dev/sfdQC documentation built on May 30, 2020, 4:57 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
knitr::opts_chunk$set(comment = NA, warning = FALSE, message = FALSE, error = FALSE, echo = FALSE, eval = TRUE, #dev = "my_png", #fig.ext = "png", fig.align = 'left', tab.align = 'left') #options(knitr.kable.NA = '-') #options(knitr.table.format = 'latex') #options(digits = 4)
my_png <- function(file, width, height) { png(file, width = width, height = height, res = 800, units = "in", pointsize = 20) } library(tidyverse) library(maps) library(lubridate) library(RColorBrewer) library(viridisLite) library(sfdQC) library(vmstools) grade <- function(x, dx) { brks <- seq(floor(min(x)), ceiling(max(x)),dx) ints <- findInterval(x, brks, all.inside = TRUE) x <- (brks[ints] + brks[ints + 1]) / 2 return(x) }
#Read in latest submission -->", ICES_LE <- read.table(params$LE_filename, sep = ',', header = TRUE, stringsAsFactors = FALSE, na.strings = 'NULL', colClasses = c('character', 'character', 'integer', 'integer', 'character', 'character', 'character', 'numeric', 'character', 'character', 'numeric', 'numeric', 'numeric')) ICES_VE <- read.table(params$VE_filename, sep = ',', header = TRUE, stringsAsFactors = FALSE, na.strings = 'NULL', colClasses = c('character', 'character', 'integer', 'integer', 'character', 'character', 'character', 'character', 'numeric', 'numeric', 'numeric', 'numeric', 'numeric', 'numeric', 'numeric', 'numeric', 'numeric'))
#Most recent submission --> country <- ICES_VE$country[1] ICES_VE <- cbind(ICES_VE, CSquare2LonLat(ICES_VE$c_square, degrees=0.05)) spatBound <- list(xrange = range(ICES_VE$SI_LONG, na.rm=TRUE), yrange = range(ICES_VE$SI_LATI, na.rm=TRUE)) spatCore <- list(xrange = quantile(ICES_VE$SI_LONG, c(0.025, 0.975)), yrange = quantile(ICES_VE$SI_LATI, c(0.025, 0.975))) tempBound <- range(ICES_VE$year, na.rm=TRUE)
#Most recent submission --> ICES_LE <- cbind(ICES_LE, ICESrectangle2LonLat(ICES_LE$ICES_rectangle, midpoint=TRUE)) spatBoundLog <- list(xrange = range(ICES_LE$SI_LONG, na.rm=TRUE), yrange = range(ICES_LE$SI_LATI, na.rm=TRUE)) tempBoundLog <- range(ICES_LE$year, na.rm=TRUE)
coord <- coord_quickmap(xlim = range(spatBound$xrange), ylim = range(spatBound$yrange)) spatBoundasp <- coord$aspect(list(x.range = range(spatBound$xrange), y.range = range(spatBound$yrange)))
VMS
hasvms <- nrow(ICES_VE) > 0 haslogbook <- nrow(ICES_LE) > 0
r if (!hasvms) {"No VMS data available\n\\begin{comment}"}
knitr::opts_chunk$set(eval = hasvms)
Records by year
d <- ICES_VE %>% count(year) d %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n)) + theme_icesqc() + geom_col() + labs(x = "Year", y = "Number of records")
\newpage
Records by month
d <- ICES_VE %>% count(month, year) d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n, fill = factor(month))) + theme_icesqc(legend.position = "right") + geom_col() + scale_fill_grey() + labs(x = "Year", y = "Number of records", fill = "Month")
\newpage
Vessel length categories
d <- ICES_VE %>% count(variable = vessel_length_category, year) d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n, fill = factor(variable))) + theme_icesqc(legend.position = "right") + geom_col() + scale_fill_grey() + labs(x = "Year", y = "Number of records", fill = "Length category")
\newpage
Gear code
d <- ICES_VE %>% count(variable = gear_code, year) d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n, fill = factor(variable))) + theme_icesqc(legend.position = "right") + geom_col() + scale_fill_grey() + labs(x = "Year", y = "Number of records", fill = "Gear code")
Unique DCF Level 6 codes
ICES_VE %>% select(variable = LE_MET_level6, year) %>% distinct() %>% count(year) %>% knitr::kable(booktabs = TRUE)
Top 5 DCF Level 6 codes
top5Met6 <- ICES_VE %>% count(LE_MET_level6) %>% slice(1:5) %>% pull(LE_MET_level6) d <- ICES_VE %>% filter(LE_MET_level6 %in% top5Met6) %>% count(variable = LE_MET_level6, year) d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n, fill = factor(variable))) + theme_icesqc(legend.position = "right") + geom_col() + scale_fill_grey() + labs(x = "Year", y = "Number of records", fill = "Level6")
\newpage
Average fishing speed
knitr::kable(do.call(rbind, tapply(ICES_VE$avg_fishing_speed, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = avg_fishing_speed)) + geom_histogram() + xlab("Average fishing speed") + ylab ("Count") + facet_wrap( ~ year, ncol = 2) + theme_icesqc()
\newpage
Average fishing hours
knitr::kable(do.call(rbind, tapply(ICES_VE$fishing_hours, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = fishing_hours)) + geom_histogram() + xlab("Average fishing hours") + ylab("Count") + facet_wrap( ~ year, ncol = 2) + theme_icesqc()
\newpage
Average length
knitr::kable(do.call(rbind, tapply(ICES_VE$avg_oal, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = avg_oal))+ geom_histogram() + xlab("Average vessel length") + ylab("Count") + facet_wrap( ~ year, ncol = 2) + theme_icesqc()
\newpage
Average kW
knitr::kable(do.call(rbind, tapply(ICES_VE$avg_kw, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = avg_kw))+ geom_histogram() + xlab("Average kW") + ylab ("Count") + facet_wrap( ~ year, ncol = 2) + theme_icesqc()
\newpage
Average kW-hours
knitr::kable(do.call(rbind, tapply(ICES_VE$kw_fishinghours, ICES_VE$year, summary)), booktabs = TRUE)
ggplot(ICES_VE, aes(x = kw_fishinghours))+ geom_histogram() + xlab("Average kW-hours") + ylab("Count") + facet_wrap(~ year, ncol = 2) + theme_icesqc()
\newpage
Landings by gear by year
ps <- gear_splits(ICES_VE$totweight, data = ICES_VE, "kg landed", year_groups = 2, gear_groups = 4, func = sum) ps$table for (p in ps$plots) print(p)
\newpage
Mean landing per kW fishing hours
ps <- gear_splits(with(ICES_VE, totweight/kw_fishinghours), data = ICES_VE, "kg/kWh", gear_groups = 4, func = median) ps$table for (p in ps$plots) print(p)
\newpage
Value by gear by year
ps <- gear_splits(ICES_VE$totvalue, data = ICES_VE, "EUR landed", gear_groups = 4, func = sum) ps$table if(any(!is.na(ICES_VE$totvalue))) for (p in ps$plots) print(p)
\newpage
Median value per KW fishing hours by year
ps <- gear_splits(with(ICES_VE, totvalue/kw_fishinghours), data = ICES_VE, "EUR/kWh", gear_groups = 4, func = median) ps$table if(any(!is.na(ICES_VE$totvalue))) for (p in ps$plots) print(p)
\newpage
Average price
ps <- gear_splits(with(ICES_VE, totvalue/totweight), data = ICES_VE, "Mean price (EUR/kg)", gear_groups = 4, func = median) ps$table if(any(!is.na(ICES_VE$totvalue))) for (p in ps$plots) print(p)
r if (!hasvms) {"\\end{comment}"}
\newpage
r if (any(is.na(ICES_VE$SI_LONG))) {"\\newpage"}
Invalid Rectangles
if (any(is.na(ICES_VE$SI_LONG))) { knitr::kable(table(`ICES Rectangle` = ICES_VE$ICES_rectangle[is.na(ICES_VE$SI_LONG)], Year = ICES_VE$year[is.na(ICES_VE$SI_LONG)]), booktabs = TRUE) } else { x <- "There were no invalid Statistical Rectangles reported" attr(x, "format") <- "markdown" attr(x, "class") <- "knit_asis" x }
\newpage
Spatial extent
Range by year
dat2tab <- as.matrix(aggregate(ICES_VE$SI_LONG,by=list(ICES_VE$year),FUN=range)) dat2tab <- cbind(dat2tab,as.matrix(aggregate(ICES_VE$SI_LATI,by=list(ICES_VE$year),FUN=range)[,-1])) colnames(dat2tab) <- c("Year", "min lon", "max lon", "min lat", "max lat") knitr::kable(dat2tab, booktabs = TRUE)
Full area
d <- ICES_VE %>% as_tibble() %>% select(long = SI_LONG, lat = SI_LATI, year) %>% mutate(long = grade(long, 0.5), lat = grade(lat, 0.5)) %>% distinct() map_data("world", xlim = spatBound$xrange, ylim = spatBound$yrange) %>% ggplot(aes(long, lat)) + theme_void() + geom_polygon(aes(group = group), fill = "grey90") + geom_raster(data = d, fill = "red") + coord_quickmap(xlim = spatBound$xrange, ylim = spatBound$yrange) + geom_path(aes(group = group), colour = "grey", lwd = 0.1) + facet_wrap(~ year, ncol = 2) + labs(x = NULL, y = NULL)
Core area
d <- ICES_VE %>% as_tibble() %>% select(long = SI_LONG, lat = SI_LATI, year) %>% distinct() map_data("world", xlim = spatCore$xrange, ylim = spatCore$yrange) %>% ggplot(aes(long, lat)) + theme_void() + geom_polygon(aes(group = group), fill = "grey90") + geom_raster(data = d, fill = "red") + coord_quickmap(xlim = spatCore$xrange, ylim = spatCore$yrange) + geom_path(aes(group = group), colour = "grey", lwd = 0.1) + facet_wrap(~ year, ncol = 2) + labs(x = NULL, y = NULL)
\newpage
3 most dominant gear
# get 3 main gears gear_table <- with(ICES_VE, aggregate(fishing_hours, list(gear_code = gear_code), sum, na.rm = TRUE)) gear_table <- gear_table[order(gear_table$x, decreasing = TRUE),] top3gears <- gear_table$gear_code[1:pmin(3, nrow(gear_table))] # split by gear idx <- which(ICES_VE$gear_code %in% top3gears) coordGrd <- unique(ICES_VE[idx, c("SI_LONG", "SI_LATI", "year", "c_square", "gear_code")]) # create a fortied polygon of csquares : 0.39 secs polVMS <- make_polVMS(coordGrd) polVMS$year <- rep(coordGrd$year, each = 5) polVMS$c_square <- rep(coordGrd$c_square,each=5) polVMS$gear_code<- rep(coordGrd$gear_code,each=5) # aggregate fishing days dat <- with(ICES_VE[idx,], aggregate(fishing_hours/24, by = list(year = year, c_square = c_square, gear_code = gear_code), FUN = sum, na.rm = TRUE)) dat <- dplyr::rename(dat, fishing_days = x) # legend calculations steps <- ceiling(max(dat$fishing_days, na.rm = TRUE)/250) breaks <- unique(c(-1, 0, steps * c(1, 2.5, 5, 10, 25, 50, 100, 250))) legval <- paste(breaks[-length(breaks)], "<=", breaks[-1L]) legval[1] <- "0" palette <- c("white", RColorBrewer::brewer.pal(length(breaks)-2, "YlOrRd")) dat$colgrp <- as.numeric(cut(dat$fishing_days, breaks = breaks)) dat$cols <- palette[dat$colgrp] # join dat onto polVMS (much faster than left_join or merge) # make joining keys a_by <- apply(dat[c("year","c_square","gear_code")], 1, paste, collapse = ".") b_by <- apply(polVMS[c("year","c_square","gear_code")], 1, paste, collapse = ".") ib_by <- as.integer(factor(b_by, levels = a_by)) polVMS$cols <- dat$cols[ib_by] # do plots for (gear in top3gears) { p <- spatialplot(polVMS[polVMS$gear_code == gear,]) + scale_fill_manual(values = rev(palette), labels = rev(legval)) + guides(fill = guide_legend(title = paste("Days@Sea", gear))) + facet_wrap(~ year, ncol = 2) + theme_icesqc(legend.position = "top") print(p) }
top3gears <- ICES_VE %>% group_by(gear_code) %>% summarise(fishing_hours = sum(fishing_hours, na.rm = TRUE)) %>% arrange(desc(fishing_hours)) %>% slice(1:3) %>% pull(gear_code) d <- ICES_VE %>% as_tibble() %>% filter(gear_code %in% top3gears) %>% select(year, long = SI_LONG, lat = SI_LATI, gear_code, fishing_hours) %>% group_by(year, gear_code, long, lat) %>% summarise(value = sum(fishing_hours, na.rm = TRUE) / 24) #d %>% group_by(year) %>% summarise(min = min(value), max = max(value)) breaks <- c(0, 1, 2.5, 5, 10, 25, 50, 100, 250) labels <- cut(0:250, breaks = breaks) %>% unique() %>% as.character() #x <- cut(d$value, breaks = breaks) # NEEDS FIXING: The top range is not shown in the legends for(i in 1:length(top3gears)) { p <- map_data("world", xlim = spatCore$xrange, ylim = spatCore$yrange) %>% ggplot(aes(long, lat)) + theme_void() + geom_polygon(aes(group = group), fill = "grey90") + geom_raster(data = d %>% filter(gear_code %in% top3gears[i]), aes(long, lat, fill = value)) + scale_fill_gradientn(colours = inferno(256, # drop the light yellow colour begin = 0.95, end = 0), guide = guide_legend(), breaks = breaks, labels = labels) + coord_quickmap(xlim = spatCore$xrange, ylim = spatCore$yrange) + geom_path(aes(group = group), colour = "grey", lwd = 0.1) + facet_wrap(~ year, ncol = 2) + labs(x = NULL, y = NULL, fill = paste("Days@Sea", top3gears[i])) + theme(legend.position = "top") print(p) }
\newpage
Spatial distribution of effort
coordGrd <- unique(ICES_VE[,c("SI_LONG","SI_LATI","year","c_square")]) # make fortified DF for csquares polVMS <- make_polVMS(coordGrd) polVMS$year <- rep(coordGrd$year, each = 5) polVMS$c_square <- rep(coordGrd$c_square, each = 5) # aggregate fishing days dat <- with(ICES_VE, aggregate(fishing_hours/24, by = list(year = year, c_square = c_square), FUN = sum, na.rm = TRUE)) dat <- dplyr::rename(dat, fishing_days = x) # legend calculations steps <- ceiling(max(dat$fishing_days, na.rm = TRUE)/250) breaks <- unique(c(-1, 0, steps * c(1, 2.5, 5, 10, 25, 50, 100, 250))) legval <- paste(breaks[-length(breaks)], "<=", breaks[-1L]) legval[1] <- "0" palette <- c("white", RColorBrewer::brewer.pal(length(breaks)-2, "YlOrRd")) dat$colgrp <- as.numeric(cut(dat$fishing_days, breaks = breaks)) dat$cols <- palette[dat$colgrp] # join dat onto polVMS (much faster than left_join or merge) # make joining keys a_by <- apply(dat[c("year","c_square")], 1, paste, collapse = ".") b_by <- apply(polVMS[c("year","c_square")], 1, paste, collapse = ".") ib_by <- as.integer(factor(b_by, levels = a_by)) polVMS$cols <- dat$cols[ib_by] # do plots spatialplot(polVMS) + scale_fill_manual(values = rev(palette), labels = rev(legval)) + guides(fill = guide_legend(title = "Days@Sea")) + facet_wrap(~ year, ncol = 2) + theme_icesqc(legend.position = "top")
d <- ICES_VE %>% as_tibble() %>% select(year, long = SI_LONG, lat = SI_LATI, fishing_hours) %>% group_by(year, long, lat) %>% summarise(value = sum(fishing_hours, na.rm = TRUE) / 24) #d %>% group_by(year) %>% summarise(min = min(value), max = max(value)) breaks <- c(0, 1, 2.5, 5, 10, 25, 50, 100, 250) labels <- cut(0:250, breaks = breaks) %>% unique() %>% as.character() #x <- cut(d$value, breaks = breaks) # Note: the top range is not shown in the legends map_data("world", xlim = spatCore$xrange, ylim = spatCore$yrange) %>% ggplot(aes(long, lat)) + theme_void() + geom_polygon(aes(group = group), fill = "grey90") + geom_raster(data = d, aes(long, lat, fill = value)) + scale_fill_gradientn(colours = inferno(256, # drop the light yellow colour begin = 0.95, end = 0), guide = guide_legend(), breaks = breaks, labels = labels) + coord_quickmap(xlim = spatCore$xrange, ylim = spatCore$yrange) + geom_path(aes(group = group), colour = "grey", lwd = 0.1) + facet_wrap(~ year, ncol = 2) + labs(x = NULL, y = NULL, fill = "Days@Sea") + theme(legend.position = "top")
\newpage
Spatial difference of effort r tempBound[1]:r (tempBound[2]-1) vs r tempBound[2]
base <- with(ICES_VE, aggregate(fishing_hours, by = list(c_square = c_square, year = year), FUN = sum, na.rm = TRUE)) base <- dplyr::rename(base, fishing_hours = x) # calculate total fishing hours for recent year recent <- base[base$year == tempBound[2],] # calculate median of the total fishing hours per square for historical years base <- with(base[base$year < tempBound[2],], aggregate(fishing_hours, by = list(c_square = c_square), FUN = median, na.rm = TRUE)) base <- dplyr::rename(base, fishing_hours_median = x) # join dat2plot <- dplyr::full_join(base, recent[,c("c_square","fishing_hours")]) # set NAs to zero dat2plot$fishing_hours_median[is.na(dat2plot$fishing_hours_median)] <- 0 dat2plot$fishing_hours[is.na(dat2plot$fishing_hours)] <- 0 # calculate ratio (with exceptions for zeros) dat2plot$ratio <- 1/with(dat2plot, pmax(fishing_hours, 1e-9) / pmax(fishing_hours_median, 1e-9)) # add back in lat and long dat2plot <- cbind(dat2plot, vmstools::CSquare2LonLat(dat2plot$c_square, degrees = 0.05)) # make 'fortified' data frame polVMS <- make_polVMS(dat2plot) polVMS$c_square <- rep(dat2plot$c_square, each = 5) # set up legend # This is not what the legend says.... there is no +/-5% break! breaks <- c(-Inf, 1/2, 1/1.5, 1/1.25, 1/1.05, 1, 1.05, 1.25, 1.5, 2, Inf) legval <- c("historic >>","historic> +100%","historic> +50%","historic> +25%", "+/-5%", "recent> +5%","recent> +25%","recent> +50%","recent> +100%","recent >>") palette <- RColorBrewer::brewer.pal(length(breaks)-1,"RdYlBu") #colour for fishing intensity dat2plot$colgrp <- as.numeric(cut(dat2plot$ratio, breaks = breaks)) dat2plot$cols <- palette[dat2plot$colgrp] # join dat onto polVMS (much faster than left_join or merge) # make joining keys a_by <- dat2plot$c_square b_by <- polVMS$c_square ib_by <- as.integer(factor(b_by, levels = a_by)) polVMS$cols <- dat2plot$cols[ib_by] # NOTE dark blue here just means that nothing was observered before, and now there is something... # - could be 1 hour of fishing grps <- sort(unique(dat2plot$colgrp)) spatialplot(polVMS) + guides(fill = guide_legend(title = "Days@Sea")) + scale_fill_manual(values = rev(palette)[grps], labels = legval[grps]) + theme_icesqc(legend.position = "right")
\newpage
Spatial difference of effort r tempBound[2]-1 vs r tempBound[2]
base <- with(ICES_VE, aggregate(fishing_hours, by = list(c_square = c_square, year = year), FUN = sum, na.rm = TRUE)) base <- dplyr::rename(base, fishing_hours = x) # calculate total fishing hours for recent year recent <- base[base$year == tempBound[2],] # previous year base <- base[base$year == tempBound[2]-1,] base <- dplyr::rename(base, fishing_hours_median = fishing_hours) # join dat2plot <- dplyr::full_join(base, recent[,c("c_square","fishing_hours")]) # set NAs to zero dat2plot$fishing_hours_median[is.na(dat2plot$fishing_hours_median)] <- 0 dat2plot$fishing_hours[is.na(dat2plot$fishing_hours)] <- 0 # calculate ratio (with exceptions for zeros) dat2plot$ratio <- 1/with(dat2plot, pmax(fishing_hours, 1e-9) / pmax(fishing_hours_median, 1e-9)) # add back in lat and long dat2plot <- cbind(dat2plot, vmstools::CSquare2LonLat(dat2plot$c_square, degrees = 0.05)) # make 'fortified' data frame polVMS <- make_polVMS(dat2plot) polVMS$c_square <- rep(dat2plot$c_square, each = 5) # set up legend # This is not what the legend says.... there is no +/-5% break! breaks <- c(-Inf, 1/2, 1/1.5, 1/1.25, 1/1.05, 1, 1.05, 1.25, 1.5, 2, Inf) legval <- c("historic >>","historic> +100%","historic> +50%","historic> +25%", "+/-5%", "recent> +5%","recent> +25%","recent> +50%","recent> +100%","recent >>") palette <- RColorBrewer::brewer.pal(length(breaks)-1,"RdYlBu") #colour for fishing intensity dat2plot$colgrp <- as.numeric(cut(dat2plot$ratio, breaks = breaks)) dat2plot$cols <- palette[dat2plot$colgrp] # join dat onto polVMS (much faster than left_join or merge) # make joining keys a_by <- dat2plot$c_square b_by <- polVMS$c_square ib_by <- as.integer(factor(b_by, levels = a_by)) polVMS$cols <- dat2plot$cols[ib_by] # NOTE dark blue here just means that nothing was observered before, and now there is something... # - could be 1 hour of fishing grps <- sort(unique(dat2plot$colgrp)) spatialplot(polVMS) + guides(fill = guide_legend(title = "Days@Sea")) + scale_fill_manual(values = rev(palette)[grps], labels = legval[grps]) + theme_icesqc(legend.position = "right")
Logbooks
r if (!haslogbook) {"No Logbook data available\n\\begin{comment}"}
knitr::opts_chunk$set(eval = haslogbook)
Records by year
d <- ICES_LE %>% count(year) d %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n)) + theme_icesqc() + geom_col() + labs(x = "Year", y = "Number of records")
\newpage
Records by month
d <- ICES_LE %>% count(month, year) d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n, fill = factor(month))) + theme_icesqc(legend.position = "right") + geom_col() + scale_fill_grey() + labs(x = "Year", y = "Number of records", fill = "Month")
\newpage
Vessel length categories
d <- ICES_LE %>% count(variable = vessel_length_category, year) d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n, fill = factor(variable))) + theme_icesqc(legend.position = "right") + geom_col() + scale_fill_grey() + labs(x = "Year", y = "Number of records", fill = "Length category")
\newpage
Gear code
d <- ICES_LE %>% count(variable = gear_code, year) d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n, fill = factor(variable))) + theme_icesqc(legend.position = "right") + geom_col() + scale_fill_grey() + labs(x = "Year", y = "Number of records", fill = "Gear code")
\newpage
Unique DCF Level 6 codes
knitr::kable(plyr::count(unique(ICES_LE[,c("LE_MET_level6","year")]),"year"), booktabs = TRUE)
Top 5 DCF Level 6 codes
top5Met6 <- ICES_LE %>% count(LE_MET_level6) %>% slice(1:5) %>% pull(LE_MET_level6) d <- ICES_LE %>% filter(LE_MET_level6 %in% top5Met6) %>% count(variable = LE_MET_level6, year) d %>% spread(year, n) %>% knitr::kable(booktabs = TRUE) d %>% ggplot(aes(factor(year), n, fill = factor(variable))) + theme_icesqc(legend.position = "right") + geom_col() + scale_fill_grey() + labs(x = "Year", y = "Number of records", fill = "Level6")
Average fishing days
knitr::kable(do.call(rbind, tapply(ICES_LE$fishing_days, ICES_LE$year, summary)), booktabs = TRUE)
ggplot(ICES_LE, aes(x = fishing_days)) + geom_histogram() + xlab("Average fishing days") + ylab("Count") + facet_wrap( ~ year, ncol = 2) + theme_icesqc()
\newpage
Average kW-hours
knitr::kable(do.call(rbind, tapply(ICES_LE$kw_fishing_days, ICES_LE$year, summary)), booktabs = TRUE)
ggplot(ICES_LE, aes(x = kw_fishing_days))+ geom_histogram() + xlab("Average KW-days") + ylab("Count") + facet_wrap(~ year, ncol = 2) + theme_icesqc()
\newpage
Landings by gear by year
ps <- gear_splits(ICES_LE$totweight, data = ICES_LE, "kg landed", gear_groups = 4, func = sum) ps$table for (p in ps$plots) print(p)
Mean landing per kW fishing days
ps <- gear_splits(with(ICES_LE, totweight/kw_fishing_days), data = ICES_LE, "kg/kWd", gear_groups = 4, func = median) ps$table for (p in ps$plots) print(p)
\newpage
Value by gear by year
ps <- gear_splits(ICES_LE$totvalue, data = ICES_LE, "EUR landed", gear_groups = 4, func = sum) ps$table if(any(!is.na(ICES_LE$totvalue))) for (p in ps$plots) print(p)
\newpage
Median value per KW fishing days by year
ps <- gear_splits(with(ICES_LE, totvalue/kw_fishing_days), data = ICES_LE, "EUR/kWh", gear_groups = 4, func = median) ps$table if(any(!is.na(ICES_LE$totvalue))) for (p in ps$plots) print(p)
\newpage
Average price
ps <- gear_splits(with(ICES_LE, totvalue/totweight), data = ICES_LE, "Mean price (EUR/kg)", gear_groups = 4, func = median) ps$table if(any(!is.na(ICES_LE$totvalue))) for (p in ps$plots) print(p)
r if (any(is.na(ICES_LE$SI_LONG))) {"\\newpage"}
Records assigned to invalid Statistical Rectangles
if (any(is.na(ICES_LE$SI_LONG))) { knitr::kable(table(`ICES Rectangle` = ICES_LE$ICES_rectangle[is.na(ICES_LE$SI_LONG)], Year = ICES_LE$year[is.na(ICES_LE$SI_LONG)]), booktabs = TRUE) } else { x <- "There were no invalid Statistical Rectangles reported" attr(x, "format") <- "markdown" attr(x, "class") <- "knit_asis" x }
\newpage
Area extent of data submitted by year:
dat2tab <- as.matrix(aggregate(ICES_LE$SI_LONG, by = list(year = ICES_LE$year), FUN=range, na.rm=TRUE)) dat2tab <- cbind(dat2tab,as.matrix(aggregate(ICES_LE$SI_LATI,by=list(ICES_LE$year), FUN=range, na.rm=TRUE)[,-1])) colnames(dat2tab) <- c("Year","min lon","max lon","min lat","max lat") knitr::kable(dat2tab, booktabs = TRUE)
Area for which data has been submitted:
d <- ICES_LE %>% select(long = SI_LONG, lat = SI_LATI, year) %>% distinct() map_data("world", xlim = spatBoundLog$xrange, ylim = spatBoundLog$yrange) %>% ggplot(aes(long, lat)) + theme_void() + geom_polygon(aes(group = group), fill = "grey90") + geom_raster(data = d, fill = "red") + coord_quickmap(xlim = spatBoundLog$xrange, ylim = spatBoundLog$yrange) + geom_path(aes(group = group), colour = "black", lwd = 0.1) + facet_wrap(~ year, ncol = 2) + labs(x = NULL, y = NULL)
map_data("world", xlim = spatCore$xrange, ylim = spatCore$yrange) %>% ggplot(aes(long, lat)) + theme_void() + geom_polygon(aes(group = group), fill = "grey90") + geom_raster(data = d, fill = "red") + coord_quickmap(xlim = spatCore$xrange, ylim = spatCore$yrange) + geom_path(aes(group = group), colour = "black", lwd = 0.1) + facet_wrap(~ year, ncol = 2) + labs(x = NULL, y = NULL)
\newpage
Frequency of VMS enabled category
knitr::kable(table(ICES_LE$vms_enabled,ICES_LE$vessel_length_category), booktabs = TRUE)
\newpage
Spatial distribution of effort by year:
coordGrd <- unique(ICES_LE[,c("SI_LONG","SI_LATI","year","ICES_rectangle")]) polRect <- make_polVMS(coordGrd, resolution = 1) polRect$year <- rep(coordGrd$year,each=5) polRect$ICES_rectangle <- rep(coordGrd$ICES_rectangle, each = 5) #TIDY --- dat <- aggregate(ICES_LE$fishing_days,by=list(ICES_LE$year,ICES_LE$ICES_rectangle),FUN=sum,na.rm=T) steps <- ceiling(max(dat$x,na.rm=T)/250) cutbreaksval <- unique(c(-1,0,steps*c(1,2.5,5,10,25,50,100,250))) legval <- outer(ac(cutbreaksval),ac(cutbreaksval),function(x,y){return(paste(x,"<=",y))}) legval <- c("0",diag(legval[-1,-c(1,2)])) palette <- c("white", brewer.pal(length(cutbreaksval)-2,"YlOrRd")) cols <- palette[cut(dat$x,breaks=cutbreaksval)] cols <- cbind(cols,id=1:length(cols),ICES_rectangle=dat$Group.2,year=dat$Group.1) polRect <- merge(polRect,cols,by=c("ICES_rectangle","year")) # ---- spatialplot(polRect) + guides(fill = guide_legend(title = "Days@Sea")) + scale_fill_manual(values = rev(palette), labels = rev(legval)) + facet_wrap(~ year, ncol = 2) + theme_icesqc(legend.position = "top")
d <- ICES_LE %>% as_tibble() %>% select(year, long = SI_LONG, lat = SI_LATI, value = fishing_days) %>% group_by(year, long, lat) %>% summarise(value = sum(value, na.rm = TRUE)) steps <- ceiling(max(d$value, na.rm = TRUE) / 250) breaks <- steps * c(0, 1, 2.5, 5, 10, 25, 50, 100, 250) labels <- cut(0:(steps * 250), breaks = breaks) %>% unique() %>% as.character() # Note: the top range is not shown in the legends map_data("world", xlim = spatCore$xrange, ylim = spatCore$yrange) %>% ggplot(aes(long, lat)) + theme_void() + geom_polygon(aes(group = group), fill = "grey90") + geom_raster(data = d, aes(long, lat, fill = value)) + scale_fill_gradientn(colours = inferno(256, # drop the light yellow colour begin = 0.95, end = 0), guide = guide_legend(), breaks = breaks, labels = labels) + coord_quickmap(xlim = spatCore$xrange, ylim = spatCore$yrange) + geom_path(aes(group = group), colour = "grey", lwd = 0.1) + facet_wrap(~ year, ncol = 2) + labs(x = NULL, y = NULL, fill = "Days@Sea") + theme(legend.position = "top")
\newpage
Spatial difference of effort r tempBound[1]:r (tempBound[2]-1) vs r tempBound[2]
base <- with(ICES_LE, aggregate(fishing_days, by = list(ICES_rectangle = ICES_rectangle, year = year), FUN = sum, na.rm = TRUE)) base <- dplyr::rename(base, fishing_days = x) # calculate total fishing hours for recent year recent <- base[base$year == tempBound[2],] # calculate median of the total fishing hours per square for historical years base <- with(base[base$year < tempBound[2],], aggregate(fishing_days, by = list(ICES_rectangle = ICES_rectangle), FUN = median, na.rm = TRUE)) base <- dplyr::rename(base, fishing_days_median = x) # join dat2plot <- dplyr::full_join(base, recent[,c("ICES_rectangle","fishing_days")]) # set NAs to zero dat2plot$fishing_days_median[is.na(dat2plot$fishing_days_median)] <- 0 dat2plot$fishing_days[is.na(dat2plot$fishing_days)] <- 0 # calculate ratio (with exceptions for zeros) dat2plot$ratio <- 1/with(dat2plot, pmax(fishing_days, 1e-9) / pmax(fishing_days_median, 1e-9)) # add back in lat and long dat2plot <- cbind(dat2plot, vmstools::ICESrectangle2LonLat(dat2plot$ICES_rectangle, midpoint = TRUE)) # make 'fortified' data frame polRect <- make_polVMS(dat2plot, resolution = 1) polRect$ICES_rectangle <- rep(dat2plot$ICES_rectangle, each = 5) ## tidy --- breaks <- rev(c(1e-10,0.5,2/3,0.8,0.952381,1,1.05,1.25,1.5,2,1e10)) legval <- c("historic >>","historic> +100%","historic> +50%","historic> +25%", "+/-5%", "recent> +5%","recent> +25%","recent> +50%","recent> +100%","recent >>") palette <- brewer.pal(length(cutbreaksval)-1,"RdYlBu") colgrp <- as.numeric(cut(dat2plot$ratio, breaks = breaks)) cols <- cbind(cols = palette[colgrp], ICES_rectangle = dat2plot$ICES_rectangle) polRect <- merge(polRect, cols, by=c("ICES_rectangle")) # --- spatialplot(polRect) + guides(fill=guide_legend(title="Days@Sea")) + scale_fill_manual(values = rev(palette), labels = legval) + theme_icesqc(legend.position = "right")
\newpage
Spatial difference of effort r tempBound[2]-1 vs r tempBound[2]
base <- with(ICES_LE, aggregate(fishing_days, by = list(ICES_rectangle = ICES_rectangle, year = year), FUN = sum, na.rm = TRUE)) base <- dplyr::rename(base, fishing_days = x) # calculate total fishing hours for recent year recent <- base[base$year == tempBound[2],] # previous year base <- base[base$year == tempBound[2]-1,] base <- dplyr::rename(base, fishing_days_median = fishing_days) # join dat2plot <- dplyr::full_join(base, recent[,c("ICES_rectangle","fishing_days")]) # set NAs to zero dat2plot$fishing_days_median[is.na(dat2plot$fishing_days_median)] <- 0 dat2plot$fishing_days[is.na(dat2plot$fishing_days)] <- 0 # calculate ratio (with exceptions for zeros) dat2plot$ratio <- 1/with(dat2plot, pmax(fishing_days, 1e-9) / pmax(fishing_days_median, 1e-9)) # add back in lat and long dat2plot <- cbind(dat2plot, vmstools::ICESrectangle2LonLat(dat2plot$ICES_rectangle, midpoint = TRUE)) # make 'fortified' data frame polRect <- make_polVMS(dat2plot, resolution = 1) polRect$ICES_rectangle <- rep(dat2plot$ICES_rectangle, each = 5) ## tidy --- breaks <- rev(c(1e-10,0.5,2/3,0.8,0.952381,1,1.05,1.25,1.5,2,1e10)) legval <- c("historic >>","historic> +100%","historic> +50%","historic> +25%", "+/-5%", "recent> +5%","recent> +25%","recent> +50%","recent> +100%","recent >>") palette <- brewer.pal(length(cutbreaksval)-1,"RdYlBu") colgrp <- as.numeric(cut(dat2plot$ratio, breaks = breaks)) cols <- cbind(cols = palette[colgrp], ICES_rectangle = dat2plot$ICES_rectangle) polRect <- merge(polRect, cols, by=c("ICES_rectangle")) # --- spatialplot(polRect) + guides(fill=guide_legend(title="Days@Sea")) + scale_fill_manual(values = rev(palette), labels = legval) + theme_icesqc(legend.position = "right")
\newpage
Relationship fishing days and total weight
ggplot(ICES_LE[ICES_LE$year == tempBoundLog[2],], aes(x = fishing_days, y = totweight)) + geom_point() + facet_wrap(~ gear_code, ncol = 3, scale = "free") + xlab("Fishing days") + ylab ("Total weight") + theme_icesqc()
\newpage
VMS vs Logbooks
Comparison of Metier level 6 reporting between logbook and VMS
ledat <- with(ICES_LE, table(LE_MET_level6, year)) vedat <- with(ICES_VE, table(LE_MET_level6, year)) dat2tab <- rbind( cbind(as.data.frame.table(ledat), data = "LE (records)"), cbind(as.data.frame.table(vedat), data = "VE (pings)")) tab <- with(dat2tab, tapply(Freq, list(LE_MET_level6, data, year), sum)) tab[tab == 0] <- NA for (i in dimnames(tab)[[3]]) { x <- tab[,,i] x <- x[apply(!is.na(x), 1, any),] cat(knitr::kable(cbind(x, year = i), booktabs = TRUE), sep = "\n") cat("\n") }
r if (!haslogbook) {"\\end{comment}"}
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.