diff --git a/indicator_processing/fishery_dependent/.Rhistory b/indicator_processing/fishery_dependent/.Rhistory deleted file mode 100644 index 5f540ee..0000000 --- a/indicator_processing/fishery_dependent/.Rhistory +++ /dev/null @@ -1,112 +0,0 @@ -rm(list=ls()) -if (!"maps" %in% installed.packages()) install.packages("maps", repos='http://cran.us.r-project.org') -library(maps) -source("C:/Users/mandy.karnauskas/Desktop/RS_FATEproject/MASTER_codes/findMinDepth.R") -load("C:/Users/mandy.karnauskas/Desktop/RS_FATEproject/MASTER_codes/GOMreleaseForScaling.RData") -rm(list=ls()) -if (!"maps" %in% installed.packages()) install.packages("maps", repos='http://cran.us.r-project.org') -library(maps) -source("C:/Users/mandy.karnauskas/Desktop/completed_manuscripts/RS_FATEproject/MASTER_codes/findMinDepth.R") -load("C:/Users/mandy.karnauskas/Desktop/completed_manuscripts/RS_FATEproject/MASTER_codes/GOMreleaseForScaling.RData") -GOM <- matfinGOM -load("C:/Users/mandy.karnauskas/Desktop/completed_manuscripts/RS_FATEproject/MASTER_codes/ATLreleaseForScaling.RData") -ATL <- matS -sum(GOM$V5)/10^5 -sum(ATL$V5)/10^5 -areaGOM <- length(unique(paste(GOM$V2, GOM$V3))) -areaATL <- length(unique(paste(ATL$V2, ATL$V3))) -areaGOM -areaATL -arearat <- areaGOM/areaATL -arearat -# calculate constants for scaling ------------------------------------- -const_scaler <- sum(GOM$V5)/sum(ATL$V5) -const_scaler -unscaled <- rbind(GOM, ATL) -dim(unscaled) -GOMsc <- GOM -ATLsc <- ATL -GOMsc$V5 <- round(GOMsc$V5 * 10) -ATLsc$V5 <- round(ATLsc$V5 * 10 * const_scaler / arearat) -scaled <- rbind(GOMsc, ATLsc) -dim(scaled) -which(scaled$V5==0) -par(mfrow=c(1,2), mex=0.8) -barplot(c(sum(GOM$V5), sum(ATL$V5)), names.arg=c("GoM", "Atl"), main="unscaled", ylab=""); abline(0,0) -barplot(c(sum(GOMsc$V5), sum(ATLsc$V5)), names.arg=c("GoM", "Atl"), main="scaled", ylab=""); abline(0,0) -# check ratios ----------------------------------------- -(sum(GOMsc$V5)) / (sum(ATLsc$V5)) -arearat -(sum(GOM$V5)) / (sum(ATL$V5)) -dev.off() -par(mar = c(5, 7, 1, 1)) -sc2 <- scaled[which(scaled$V6==2013 & scaled$V7==6 & scaled$V8==27),] -x <- sc2$V5 / 130 -pos <- c(0.005, 0.05, 0.1, 0.2, 0.5, 1, 2, 4, 5, 10, 20, 50, 100, 1000) -a <- floor(min(x)) -b <- max((x)-a)*1.03 -pind <- round((x-a)/b*100+1); print(min(pind)); print(max(pind)) -cols <- c(rainbow(30, start=0.82, end=0.99), rainbow(70, start=0.01, end=0.17))[100:1] -map('state', fill = 1, interior=F, col = gray(0.85), ylim=c(22.5, 35), xlim=c(-88,-76)) -#mtext(side = 1, line = 2.5, "longitude") -#mtext(side = 2, line = 2.5, "latitude") -points(sc2$V2, sc2$V3, col=cols[pind], pch=15, cex=0.5) -#box(); axis(1); axis(2, las = 2) -xloc <- seq(-86, -78, length.out=100) -for (j in 1:100) { polygon(c(xloc[j], xloc[j+1],xloc[j+1], xloc[j]), c(23.0,23.0,23.4,23.4), col=cols[j], border=NA) } -w <- which.min(abs(((max(x)-min(x))/6) - pos)) -if(-pos[w]min(x)] } else { xx <- c(seq(-pos[w], min(x), -pos[w]), seq(0, max(x), pos[w])) } -text(xloc[round((xx-a)/b*100+1)], y=22.75, xx, pos=2) -text(-82, 24.25, "relative index of spawning output", pos = 1) -degs = seq(88, 76, -2) -a = sapply(degs, function(x) bquote(.(x)*degree ~ W)) -rm(list = ls()) -library(rfishbase) -setwd("C:/Users/mandy.karnauskas/Desktop/CONFIDENTIAL/CaribbeanData/Jun2022/") -rm(list = ls()) -setwd("C:/Users/mandy.karnauskas/Desktop/Caribbean-ESR/indicator_processing/fishery_dependent") -sc <- read.csv("C:/Users/mandy.karnauskas/Desktop/CONFIDENTIAL/CaribbeanData/Jun2022/shellcatch_pr_data_req_02152023_C.csv") # original shellcatch data -sc$POUNDS_LANDED2 <- as.numeric(as.vector(sc$POUNDS_LANDED)) -sc$ADJUSTED_POUNDS2 <- as.numeric(as.vector(sc$ADJUSTED_POUNDS)) -sc$CORRECTION_FACTOR <- as.numeric(as.vector(sc$CORRECTION_FACTOR)) -sc$PRICE <- as.numeric(as.vector(sc$PRICE)) -summary(sc$POUNDS_LANDED == sc$POUNDS_LANDED2) -summary(sc$ADJUSTED_POUNDS == sc$ADJUSTED_POUNDS2) -apply(sc[1:4], 2, table, useNA = "always") -table(sc$LANDING_AREA_COUNTY_OR_MUNICIPALITY) -table(sc$COAST) -table(sc$GEAR_NAME) -hist(sc$POUNDS_LANDED2) -hist(sc$ADJUSTED_POUNDS2) -table(sc$CORRECTION_FACTOR) -co <- sc$POUNDS_LANDED2 / sc$CORRECTION_FACTOR -table(round(co- sc$ADJUSTED_POUNDS2)) -hist(sc$PRICE) -table(sc$AREA_CD1) -sc$SPECIES_NM <- as.character(sc$SPECIES_NM) -# replace duplicate names ----------------------- -sc$SPECIES_NM[which(sc$SPECIES_NM == "CROAKERS")] <- "CROAKER" -sc$SPECIES_NM[which(sc$SPECIES_NM == "HERRING")] <- "THREAD HERRING" -sc$SPECIES_NM[which(sc$SPECIES_NM == "SMOOTHTAIL SPINY LOBSTER")] <- "SPINY LOBSTER" -sc$SPECIES_NM[which(sc$SPECIES_NM == "KING MACKAREL, KINGFISH")] <- "KINGFISH MACKEREL" -sc$SPECIES_NM[which(sc$SPECIES_NM == "WENCHMAN")] <- "CARDINAL" -sc$SPECIES_NM[which(sc$SPECIES_NM == "SPOTTED TRUNKFISH")] <- "TRUNKFISH" -sc$SPECIES_NM[which(sc$SPECIES_NM == "SMOOTH TRUNKFISH")] <- "TRUNKFISH" -sc$SPECIES_NM[which(sc$SPECIES_NM == "SEA BASSES")] <- "GROUPERS" -nam <- read.csv("name_matches_edited.csv") -nam -head(nam) -dim(nam) -unique(nam$shellcatch) -length(unique(nam$shellcatch)) -length(unique(nam$logbook)) -length(unique(nam$logbook)) -unique(nam$second_name) -tail(nam) -unique(nam$third_name) -length(unique(nam$shellcatch)) -unique(nam$shellcatch) -unique(nam$logbook) -135/300 -dat -nam diff --git a/indicator_processing/fishery_dependent/INDICATOR_total_landings.R b/indicator_processing/fishery_dependent/INDICATOR_total_landings.R index 33bcd36..8e44911 100644 --- a/indicator_processing/fishery_dependent/INDICATOR_total_landings.R +++ b/indicator_processing/fishery_dependent/INDICATOR_total_landings.R @@ -5,17 +5,13 @@ rm(list = ls()) # input data for Puerto Rico --------------------------- setwd("C:/Users/mandy.karnauskas/Desktop/Caribbean-ESR/indicator_processing/fishery_dependent/") -dat <- read.csv("C:/Users/mandy.karnauskas/Desktop/CONFIDENTIAL/CaribbeanData/Jun2022/PR_landings_83_20_wSC.csv") +dat <- read.csv("C:/Users/mandy.karnauskas/Desktop/CONFIDENTIAL/CaribbeanData/Jun2022/PR_landings_83_20_wSC_2005cor.csv") table(dat$YEAR_LANDED) # check multiplier adjustments ------------------------- -dat$xADJ <- dat$POUNDS_LANDED * 1/dat$CORRECTION_FACTOR -hist(dat$ADJUSTED_POUNDS - dat$xADJ) -max(abs(dat$ADJUSTED_POUNDS - dat$xADJ), na.rm = T) -table(round(dat$ADJUSTED_POUNDS) - round(dat$xADJ)) -#dat$ADJUSTED_POUNDS <- dat$xADJ +summary(dat$xADJ == dat$POUNDS_LANDED * 1/dat$CORRECTION_FACTOR) # define start and end years --------------------------- styear <- 1990 @@ -35,20 +31,25 @@ sort(table(d$ITIS_COMMON_NAME[grep("LOBSTER", d$ITIS_COMMON_NAME)])) sort(table(d$ITIS_COMMON_NAME[grep("CONCH", d$ITIS_COMMON_NAME)])) d$sppgrp <- "other" -#d$sppgrp[grep("LOBSTER", d$ITIS_COMMON_NAME)] <- "lobster" -d$sppgrp[grep("LOBSTERS,SPINY", d$ITIS_COMMON_NAME)] <- "lobster" +d$sppgrp[grep("LOBSTER", d$ITIS_COMMON_NAME)] <- "lobster" +#d$sppgrp[grep("LOBSTERS,SPINY", d$ITIS_COMMON_NAME)] <- "lobster" # to compare directly to only spiny lobster d$sppgrp[grep("CONCH", d$ITIS_COMMON_NAME)] <- "conch" table(d$ITIS_COMMON_NAME, d$sppgrp) +# check for rule of 3 ------------------------ + +table(d$sppgrp, d$YEAR_LANDED) +table(table(d$sppgrp, d$YEAR_LANDED) < 3) + # sum landings by year ---------------------- totland_pr <- tapply(d$ADJUSTED_POUNDS, list(d$YEAR_LANDED, d$sppgrp), sum, na.rm = T) / 10^3 dim(totland_pr) totland_pr -matplot(totland_pr, type = "l", lty = 1, lwd = 2) +matplot(1990:2020, totland_pr, type = "l", lty = 1, lwd = 2) -ls() +rm(list = ls()[-match(c("totland_pr", "styear", "enyear"), ls())]) # calculate for STT -------------------------------------- @@ -76,28 +77,81 @@ d$sppgrp[grep("CONCH", d$SPECIES_NM)] <- "conch" table(d$SPECIES_NM, d$sppgrp) +# check for rule of 3 ------------------------ + +table(d$TRIP_YEAR, d$sppgrp) +table(table(d$TRIP_YEAR, d$sppgrp) <= 3) + # sum landings by year ---------------------- totland_st <- tapply(d$POUNDS_LANDED, list(d$TRIP_YEAR, d$sppgrp), sum, na.rm = T) / 10^3 +totland_st[3, 1] <- NA # fix confidentiality issue dim(totland_st) totland_st totland_st[is.na(totland_st)] <- 0 matplot(totland_st, type = "l") +# calculate for STX -------------------------------------- + +rm(list = ls()[-match(c("totland_st", "totland_pr", "styear", "enyear"), ls())]) + +dat <- read.csv("C:/Users/mandy.karnauskas/Desktop/CONFIDENTIAL/CaribbeanData/STX_072011_present_LANDINGS_trip_2021-03-11.csv") +head(dat) + +# take a look at data fields ---------------------------- + +table(dat$LANDING_AREA_NAME) +table(dat$TRIP_YEAR) +d <- dat[which(dat$TRIP_YEAR >= styear & dat$TRIP_YEAR <= enyear), ] +d$TRIP_YEAR <- factor(d$TRIP_YEAR, levels = c(styear: enyear)) +table(d$TRIP_YEAR) + +# check field codes ----------------- + +table(d$TRIP_YEAR, useNA = "always") + +# look at top landings ----------------------- + +sort(tapply(d$POUNDS_LANDED, d$SPECIES_NM, sum, na.rm = T)) +sort(table(d$SPECIES_NM[grep("LOBSTER", d$SPECIES_NM)])) +sort(table(d$SPECIES_NM[grep("CONCH", d$SPECIES_NM)])) + +d$sppgrp <- "other" +d$sppgrp[grep("LOBSTER", d$SPECIES_NM)] <- "lobster" +d$sppgrp[grep("CONCH", d$SPECIES_NM)] <- "conch" + +table(d$SPECIES_NM, d$sppgrp) + +# check for rule of 3 ------------------------ + +table(d$TRIP_YEAR, d$sppgrp) +table(table(d$TRIP_YEAR, d$sppgrp) <= 3) + +# sum landings by year ---------------------- + +totland_sx <- tapply(d$POUNDS_LANDED, list(d$TRIP_YEAR, d$sppgrp), sum, na.rm = T) / 10^3 +dim(totland_sx) +totland_sx +#totland_sx[is.na(totland_sx)] <- 0 +matplot(totland_sx, type = "l") + # summarize and plot --------------------- ls()[grep("totland", ls())] datdata <- styear:enyear -inddata <- data.frame(cbind(totland_pr[, 2], totland_st[, 2], - totland_pr[, 1], totland_st[, 1], - totland_pr[, 3], totland_st[, 3])) +inddata <- data.frame(cbind(totland_pr[, 2], totland_st[, 2], totland_sx[, 2], + totland_pr[, 1], totland_st[, 1], totland_sx[, 1], + totland_pr[, 3], totland_st[, 3], totland_sx[, 3])) labs <- c("Lobster landings", "thousands of pounds", "Puerto Rico", "Lobster landings", "thousands of pounds", "St. Thomas and St. John", + "Lobster landings", "thousands of pounds", "St. Croix", "Conch landings", "thousands of pounds", "Puerto Rico", "Conch landings", "thousands of pounds", "St. Thomas and St. John", + "Conch landings", "thousands of pounds", "St. Croix", "Landings of all other species", "thousands of pounds", "Puerto Rico", - "Landings of all other species", "thousands of pounds", "St. Thomas and St. John") + "Landings of all other species", "thousands of pounds", "St. Thomas and St. John", + "Landings of all other species", "thousands of pounds", "St. Croix") indnames <- data.frame(matrix(labs, nrow = 3, byrow = F)) s <- list(labels = indnames, indicators = inddata, datelist = datdata) #, ulim = ulidata, llim = llidata) @@ -105,7 +159,7 @@ class(s) <- "indicatordata" setwd("C:/Users/mandy.karnauskas/Desktop/Caribbean-ESR/indicator_plots/") -plotIndicatorTimeSeries(s, coltoplot = 1:6, plotrownum = 2, sublabel = T, sameYscale = F, +plotIndicatorTimeSeries(s, coltoplot = 1:9, plotrownum = 3, plotcolnum = 3, sublabel = T, sameYscale = F, widadj = 0.8, hgtadj = 0.7, trendAnalysis = F) # outtype = "png") #inddata <- s diff --git a/indicator_processing/fishery_dependent/PREPROCESS_build_reference_file.R b/indicator_processing/fishery_dependent/PREPROCESS_build_reference_file.R index 5e0e65e..052df17 100644 --- a/indicator_processing/fishery_dependent/PREPROCESS_build_reference_file.R +++ b/indicator_processing/fishery_dependent/PREPROCESS_build_reference_file.R @@ -166,7 +166,76 @@ dim(dat2) #write.table(dat2, file = "C:/Users/mandy.karnauskas/Desktop/CONFIDENTIAL/CaribbeanData/Jun2022/PR_landings_83_20_wSC.csv", # sep = ",", col.names = T, row.names = F) -d <- dat2 +# fix error with 2005 correction factor ---------------- + +rm(list = ls()) +dat <- read.csv("C:/Users/mandy.karnauskas/Desktop/CONFIDENTIAL/CaribbeanData/Jun2022/PR_landings_83_20_wSC.csv") + +table(dat$YEAR_LANDED) + +# check multiplier adjustments ------------------------- + +dat$xADJ <- dat$POUNDS_LANDED * 1/dat$CORRECTION_FACTOR +hist(dat$ADJUSTED_POUNDS - dat$xADJ) +max(abs(dat$ADJUSTED_POUNDS - dat$xADJ), na.rm = T) +table(round(dat$ADJUSTED_POUNDS) - round(dat$xADJ)) +#dat$ADJUSTED_POUNDS <- dat$xADJ + +# fix 2005 correction factor --------------------------- + +d5 <- dat[which(dat$YEAR_LANDED == 2005), ] +table(d5$CORRECTION_FACTOR, d5$MUNICIPALITY) +plot(d5$CORRECTION_FACTOR ~ d5$MUNICIPALITY) +tapply(d5$CORRECTION_FACTOR, d5$MUNICIPALITY, sd) + +d4 <- dat[which(dat$YEAR_LANDED == 2004), ] +d6 <- dat[which(dat$YEAR_LANDED == 2006), ] + +tapply(d4$CORRECTION_FACTOR, d4$MUNICIPALITY, sd, na.rm = T) +tapply(d6$CORRECTION_FACTOR, d6$MUNICIPALITY, sd, na.rm = T) + +c4 <- tapply(d4$CORRECTION_FACTOR, d4$MUNICIPALITY, mean, na.rm = T) +c6 <- tapply(d6$CORRECTION_FACTOR, d6$MUNICIPALITY, mean, na.rm = T) +summary(names(c4) == names(c6)) + +c5 <- (c4 + c6) / 2 +c5[10] +c6[10] +table(c4) +table(c6) +table(c5) + +table(d5$CORRECTION_FACTOR) +lis265 <- names(which(c5 == 0.265)) +lis465 <- names(which(c5 == 0.465)) +lis475 <- names(which(c5 == 0.475)) +lis835 <- names(which(c5 == 0.835)) + +dat$CORRECTION_FACTOR[which(dat$YEAR_LANDED == 2005 & dat$MUNICIPALITY %in% lis265)] <- 0.265 +dat$CORRECTION_FACTOR[which(dat$YEAR_LANDED == 2005 & dat$MUNICIPALITY %in% lis465)] <- 0.465 +dat$CORRECTION_FACTOR[which(dat$YEAR_LANDED == 2005 & dat$MUNICIPALITY %in% lis475)] <- 0.475 +dat$CORRECTION_FACTOR[which(dat$YEAR_LANDED == 2005 & dat$MUNICIPALITY %in% lis835)] <- 0.835 + +table(dat$CORRECTION_FACTOR[which(dat$YEAR_LANDED == 2005)]) +dat[which(dat$CORRECTION_FACTOR > 0.9 & dat$YEAR_LANDED == 2005), ] +c6[10] +dat$CORRECTION_FACTOR[which(dat$CORRECTION_FACTOR > 0.9 & dat$YEAR_LANDED == 2005)] <- c6[10] +table(dat$CORRECTION_FACTOR[which(dat$YEAR_LANDED == 2005)]) + + +# check +dat$xADJ <- dat$POUNDS_LANDED * 1/dat$CORRECTION_FACTOR +hist(dat$ADJUSTED_POUNDS - dat$xADJ) +max(abs(dat$ADJUSTED_POUNDS - dat$xADJ), na.rm = T) +table(round(dat$ADJUSTED_POUNDS) - round(dat$xADJ)) +dat$ADJUSTED_POUNDS <- dat$xADJ + +#write.table(dat, file = "C:/Users/mandy.karnauskas/Desktop/CONFIDENTIAL/CaribbeanData/Jun2022/PR_landings_83_20_wSC_2005cor.csv", +# sep = ",", col.names = T, row.names = F) + +################################################################# + +d <- dat # take a look at data fields ---------------------------- table(d$VESSEL, useNA = "always")