OHI Science | Citation policy
The commercial fishing pressure layers are created from catch data provided by the Sea Around Us Project and net primary production data from the Vertically Generalized Production Model (VGPM) as described in Behrenfeld and Falkowski (1997)
Two layers are created in this analysis, commercial fishing pressure from high bycatch gear and low bycatch gear. Quantifying the amount of catch caught by these different categories of fishing gear is done by combining two different pieces of information.
The amount of fish catch caught per gear type was calculated for the first global Cumulative Human Impact analysis (Halpern 2008).
Total fish catch (tons) per half degree cell globally
The amount of catch caught per gear type has not been recalculated since 2008. The total proportion of catch caught per gear type is calculated from this original data, and then applied to the more recent catch data to get an updated estimate of catch per gear type. This assumes that the proportion of total catch caught by high and low bycatch has remained consistent since 2008.
The fisheries catch data has been updated by the Sea Around Us Project and is provided at a spatial resolution of 0.5 degree. Previously, only aggregate catch at the country level was used.
Reference: Pauly D. and Zeller D. (Editors), 2015. Sea Around Us Concepts, Design and Data (seaaroundus.org)
Downloaded: June 27, 2016 (sent by email from Ar’ash Tavakolie)
Description: Catch per half degree cell (tons)
Native data resolution: 0.5 degree
Time range: 1950 - 2010
Format: Tabular
Load all relevant libraries including parallel processing packages.
#set options for all chunks in code
knitr::opts_chunk$set(warning=FALSE, message=FALSE,fig.width = 6, fig.height = 4, fig.path = 'figs/')
library(parallel)
library(foreach)
library(doParallel)
library(raster)
library(rasterVis)
library(RColorBrewer)
options(scipen=999)
cols = rev(colorRampPalette(brewer.pal(11, 'Spectral'))(255)) # rainbow color scheme
mytheme=rasterTheme(region=cols)
source("~/github/ohiprep/src/R/common.R")
#ocean raster at 1km
ocean = raster(file.path(dir_M, 'model/GL-NCEAS-Halpern2008/tmp/ocean.tif'))
#set mollweide projection CRS
mollCRS=crs('+proj=moll +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs')
#paths
dir_git <- file.path(dir_M,'git-annex/globalprep/prs_fish/v2016')These two rasters are at a resolution of 1km2 and were created in this dataprep.R script
gear_hb = raster(file.path(dir_M,'git-annex/globalprep/prs_fish/v2015/gear_prop_hb_gcs.tif'))
gear_lb = raster(file.path(dir_M,'git-annex/globalprep/prs_fish/v2015/gear_prop_lb_gcs.tif'))
stack(gear_hb,gear_lb)%>%
`names<-`(c("High Bycatch Gear","Low Bycatch Gear"))%>%
levelplot(par.settings = mytheme, main = "Proportion of total catch caught by gear type")
The Net Primary Production data was prepared in npp.Rmd.
npp <- list.files(file.path(dir_M,'git-annex/globalprep/prs_fish/v2016/VGPM_primary_productivity/int/annual_npp'),pattern = 'npp_2',full.names=T)
plot(raster(npp[13]),col=cols,axes=F,main = 'Net Primary Production (mg C/m2/day) \n 2015')
We only use catch data from 2003 - 2010 to calculate these layers since those are the time periods where we also have primary productivity data. And we don’t need to go further back than 2003 since all OHI assessments will use this time period.
catch_files <- list.files('int/annual_industrial_rasters',full.names = T)[54:61]
plot(raster(catch_files[8]),col=cols,axes=F,box=F,main = "Catch (tons) in 2010")
fish_fun is a function that takes each annual catch raster and splits it into two rasters according to the proportion of catch in each cell caught by high and low bycatch gear.
# function that takes each annual catch raster, reprojects and resamples, splits it according to high and low bycatch and then standardizes by npp for that year.
fish_fun <- function(file){
yr <- substr(file,31,34)
r <- raster(file)
hb <- overlay(r,gear_hb,fun=function(x,y){x*y},filename = paste0(file.path(dir_git,'int/annual_hb_catch/annual_hb_catch_'),yr,'.tif',sep=''),overwrite=T)
lb <- overlay(r,gear_lb,fun=function(x,y){x*y},filename = paste0(file.path(dir_git,'int/annual_lb_catch/annual_lb_catch_'),yr,'.tif',sep=''),overwrite=T)
}
mclapply(catch_files,fish_fun, mc.cores = 6)Total catch per cell is standardized by the NPP values. This is done because the same fishing pressure can have different impacts depending on the productivity in the region.
Before standardizing, the NPP data is aggregated to the same spatial resolution as the catch data, 0.5 degree cells, to accurately calculate catch in tons/km2 standardized by mg C/m2/day.
hb_files <- list.files(file.path(dir_git,'int/annual_hb_catch'),full.names=T)
lb_files <- list.files(file.path(dir_git,'int/annual_lb_catch'),full.names=T)
registerDoParallel(4)
foreach (yr = c(2003:2010)) %dopar%{
#1. get net primary production for that year
n <- npp[substr(npp,111,114)==yr]%>%
raster()%>%
projectRaster(gear_hb)%>%
resample(.,gear_hb)
#2. Get bycatch rasters for the yr
cat_hb <- hb_files[substr(hb_files,90,93) == yr]%>%
raster()
cat_lb <- lb_files[substr(lb_files,90,93) == yr]%>%
raster()
#3. Divide catch by npp and save to git-annex
overlay(cat_hb,n,fun=function(x,y){x/y},filename = paste0(dir_git,'/int/annual_catch_npp/hb_catch_npp_',yr,'.tif'), overwrite=T)
overlay(cat_lb,n,fun=function(x,y){x/y},filename = paste0(dir_git,'/int/annual_catch_npp/lb_catch_npp_',yr,'.tif'), overwrite=T)
}Mean catch per cell is calculated over a rolling window of 5 years to account for interannual variability. The data is then log transformed.
lb_npp <- list.files(file.path(dir_git,'/int/annual_catch_npp'),full.names=T, pattern = 'lb')
hb_npp <- list.files(file.path(dir_git,'/int/annual_catch_npp'),full.names=T, pattern = 'hb')
registerDoParallel(4)
foreach (i = 2003:2006) %dopar%{
yrs <- c(i:(i+4))
out_lb <- lb_npp[which(substr(lb_npp,89,92) %in% yrs)]%>%
stack()%>%
calc(fun=function(x){mean(x,na.rm=T)}, filename = paste0(dir_git,'/int/mean_catch_5yr/mean_lb_catch_',yrs[1],'_',yrs[5],'.tif'),overwrite=T)%>%
calc(fun=function(x){log(x+1)})
writeRaster(out_lb,filename = paste0(dir_git,'/int/mean_catch_5yr/mean_lb_catch_log',yrs[1],'_',yrs[5],'.tif'),overwrite=T)
out_hb <- hb_npp[which(substr(hb_npp,89,92) %in% yrs)]%>%
stack()%>%
calc(fun=function(x){mean(x,na.rm=T)}, filename = paste0(dir_git,'/int/mean_catch_5yr/mean_hb_catch_',yrs[1],'_',yrs[5],'.tif'),overwrite=T)%>%
calc(fun=function(x){log(x+1)})
writeRaster(out_hb,filename = paste0(dir_git,'/int/mean_catch_5yr/mean_hb_catch_log',yrs[1],'_',yrs[5],'.tif'),overwrite=T)
}Look at all catch data standardized by NPP from 2003 - 2010 and use the 99.99th quantile as the reference point.
lb_npp <- list.files(file.path(dir_git,'/int/annual_catch_npp'),full.names=T, pattern = 'lb')
hb_npp <- list.files(file.path(dir_git,'/int/annual_catch_npp'),full.names=T, pattern = 'hb')
# HIGH BYCATCH
#get data across all years
vals <- c()
for(i in 2003:2010){
#print(i)
m <- hb_npp[which(substr(hb_npp,89,92) == i)]%>%
raster()%>%
getValues()
n <- m[!is.na(m)]
vals <- c(vals,n)
}
min_v_hb <- min(vals,na.rm=T) # 0
max_v_hb <- max(vals,na.rm=T) # 11.28059
ref_hb <- log(quantile(vals,prob = 0.9999,na.rm=T)+1) #10.29217
#LOW BYCATCH
#get data across all years
lb_vals <- c()
for(i in 2003:2010){
m <- lb_npp[which(substr(lb_npp,89,92) == i)]%>%
raster()%>%
getValues()
n <- m[!is.na(m)]
lb_vals <- c(lb_vals,n)
}
min_v_lb <- min(lb_vals,na.rm=T) # 0
max_v_lb <- max(lb_vals,na.rm=T) # 12.25
ref_lb <- log(quantile(lb_vals,prob = 0.9999,na.rm=T)+1) #11.13The reference point is 10.2921663 (0 - 79266.8067432) for high bycatch and 11.1288256 (0 - 208704.3484551) for low bycatch.
mean_hb <- list.files(file.path(dir_git,'int/mean_catch_5yr'),full.names=T, pattern = 'hb_catch_log')
mean_lb <- list.files(file.path(dir_git,'int/mean_catch_5yr'),full.names=T, pattern = 'lb_catch_log')
registerDoParallel(4)
foreach (i = 2003:2006) %dopar%{
yrs <- c(i:(i+4))
out_hb <- mean_hb[which(substr(mean_hb,90,93) == i)]%>%
raster()%>%
calc(fun=function(x){ifelse(x>ref_hb,1,x/ref_hb)})%>%
projectRaster(crs = mollCRS,over=T)%>%
resample(ocean,method = 'ngb',filename = paste0(file.path(dir_git),'/out/high_bycatch/hb_fish_pressure_',yrs[1],'-',yrs[5],'.tif'),overwrite=T)
out_lb <- mean_lb[which(substr(mean_lb,90,93) == i)]%>%
raster()%>%
calc(fun=function(x){ifelse(x>ref_lb,1,x/ref_lb)})%>%
projectRaster(crs = mollCRS, over=T)%>%
resample(ocean,method = 'ngb',filename = paste0(file.path(dir_git),'/out/low_bycatch/lb_fish_pressure_',yrs[1],'-',yrs[5],'.tif'),overwrite=T)
}hb <- raster(file.path(dir_git,'out/high_bycatch/hb_fish_pressure_2006-2010.tif'))
lb <- raster(file.path(dir_git,'out/low_bycatch/lb_fish_pressure_2006-2010.tif'))
plot(hb,col=cols,main = 'High Bycatch Fishing Pressure Layer \n OHI 2016', axes = F)
plot(lb, col=cols, main = 'Low Bycatch Fishing Pressure Layer \n OHI 2016', axes = F)
setwd('globalprep/prs_fish/v2016')
source("~/github/ohiprep//src/R/common.R")
#paths
dir_git <- file.path(dir_M,'git-annex/globalprep/prs_fish/v2016')
library(raster)
library(rgdal)
library(dplyr)
# raster/zonal data
zones <- raster(file.path(dir_M, "git-annex/globalprep/spatial/d2014/data/rgn_mol_raster_1km/sp_mol_raster_1km.tif")) # raster data
rgn_data <- read.csv(file.path(dir_M, "git-annex/globalprep/spatial/d2014/data/rgn_mol_raster_1km/regionData.csv")) # data for sp_id's used in raster
### Low bycatch data first
# get raster data:
rasts <- list.files(file.path(dir_M,'git-annex/globalprep/prs_fish/v2016/out/low_bycatch'))
pressure_stack <- stack()
for(raster in rasts){ # raster = "lb_fish_pressure_2003-2007.tif"
tmp <- raster(file.path(dir_M,'git-annex/globalprep/prs_fish/v2016/out/low_bycatch', raster))
pressure_stack <- stack(pressure_stack, tmp)
}
# extract data for each region:
regions_stats <- zonal(pressure_stack, zones, fun="mean", na.rm=TRUE, progress="text")
regions_stats2 <- data.frame(regions_stats)
setdiff(regions_stats2$zone, rgn_data$ant_id) # antarctica regions are in there, makes sense....no land
setdiff(rgn_data$ant_id, regions_stats2$zone) # 213 is in there, that makes sense (Antarctica)
data <- merge(rgn_data, regions_stats, all.y=TRUE, by.x="rgn_id", by.y="zone") %>%
gather("year", "pressure_score", starts_with("lb"))
lb_data <- data %>%
mutate(year = substring(year, 23, 26)) %>%
mutate(year = as.numeric(year)) %>%
mutate(year = year + 6) %>%
filter(ant_typ == "eez") %>%
dplyr::select(rgn_id, rgn_nam, year, pressure_score)
write.csv(lb_data, "int/lb.csv", row.names=FALSE)
## save toolbox data for different years/regions
# function to extract data more easily
saveData <- function(newYear){ # newYear= 2012
criteria_year <- ~year == newYear
tmp <- lb_data %>%
filter_(criteria_year) %>%
dplyr::select(rgn_id, pressure_score) %>%
arrange(rgn_id)
write.csv(tmp, sprintf('output/comm_fish_lb_%s.csv', newYear), row.names=FALSE)
}
### extract data
for(newYear in (max(lb_data$year) - 3):(max(lb_data$year))){
saveData(newYear)
}
### try visualizing the data using googleVis plot
library(googleVis)
plotData <- lb_data %>%
dplyr::select(rgn_nam, year, pressure_score)
Motion=gvisMotionChart(plotData,
idvar="rgn_nam",
timevar="year")
plot(Motion)
print(Motion, file='lb.html')
### Get the high bycatch data
rasts <- list.files(file.path(dir_M,'git-annex/globalprep/prs_fish/v2016/out/high_bycatch'))
pressure_stack <- stack()
for(raster in rasts){ # raster = "lb_fish_pressure_rescaled_2003-2007.tif"
tmp <- raster(file.path(dir_M,'git-annex/globalprep/prs_fish/v2016/out/high_bycatch', raster))
pressure_stack <- stack(pressure_stack, tmp)
}
# extract data for each region:
regions_stats <- zonal(pressure_stack, zones, fun="mean", na.rm=TRUE, progress="text")
regions_stats2 <- data.frame(regions_stats)
setdiff(regions_stats2$zone, rgn_data$ant_id) # antarctica regions are in there, makes sense....no land
setdiff(rgn_data$ant_id, regions_stats2$zone) # 213 is in there, that makes sense (Antarctica)
data <- merge(rgn_data, regions_stats, all.y=TRUE, by.x="rgn_id", by.y="zone") %>%
gather("year", "pressure_score", starts_with("hb"))
hb_data <- data %>%
mutate(year = substring(year, 23, 26)) %>%
mutate(year = as.numeric(year)) %>%
mutate(year = year + 6) %>%
filter(ant_typ == "eez") %>%
dplyr::select(rgn_id, rgn_nam, year, pressure_score)
write.csv(hb_data, "int/hb.csv", row.names=FALSE)
## save toolbox data for different years/regions
# function to extract data more easily
saveData <- function(newYear){ # newYear= 2012
criteria_year <- ~year == newYear
tmp <- lb_data %>%
filter_(criteria_year) %>%
dplyr::select(rgn_id, pressure_score) %>%
arrange(rgn_id)
write.csv(tmp, sprintf('output/comm_fish_hb_%s.csv', newYear), row.names=FALSE)
}
### extract data
for(newYear in (max(hb_data$year) - 3):(max(hb_data$year))){
saveData(newYear)
}
### try visualizing the data using googleVis plot
library(googleVis)
plotData <- hb_data %>%
dplyr::select(rgn_nam, year, pressure_score)
Motion=gvisMotionChart(plotData,
idvar="rgn_nam",
timevar="year")
plot(Motion)
print(Motion, file='hb.html')Pauly D. and Zeller D. (Editors), 2015. Sea Around Us Concepts, Design and Data (seaaroundus.org)