#' Making an overview barplot for total landings per scenario for a selected set of stocks
#'
#' This function processes further the lst_loglike objects
#'
#' @param fname First name
#' @param lname Last name
#' @export
#' @examples
#'
#' \dontrun{
#' general <- setGeneralOverallVariable (pathToRawInputs =file.path("C:", "Users", "fbas",
#' "Documents", "GitHub", paste0("DISPLACE_input_gis_",
#' "DanishFleet")),
#' pathToDisplaceInputs = file.path("C:", "Users", "fbas",
#' "Documents", "GitHub", paste0("DISPLACE_input_", "DanishFleet")),
#' pathToOutputs =file.path("C:","DISPLACE_outputs"),
#' caseStudy="DanishFleet",
#' iGraph=41,
#' iYear="2015",
#' iCountry="DEN",
#' nbPops=39,
#' nbSzgroup=14,
#' theScenarios= c("svana_baseline",
#' "svana_sub1mx20",
#' "svana_sub4mx20",
#' "svana_sub4mx5ns20bt",
#' "svana_sub4mx20ns5bt",
#' "svana_sub4mx5ns5bt" ),
#' nbSimus=20,
#' useSQLite=FALSE
#' )
#'
#'
#'
#'
#'
#' loadLoglikeFiles(general=general, use_port_info=FALSE)
#'
#'
#'
#' barplotTotLandingsPerSce (general=general,
#' selected="_selected_set1_",
#' type_of_column="pop",
#' selected_pops=c(0:38),
#' group1= c(0, 1, 2, 3, 11, 23, 24, 26, 30, 31, 32),
#' selected_scenarios=general$namefolderoutput,
#' scenarios_names=general$namefolderoutput,
#' nby=5, firsty="2015", lasty="2019",
#' a_width=3500, a_height=2000, black_and_white = FALSE, ylims = c(0,305))
#'
#' Note that if black_and_white at true then the barplot is split into 2 groups
#' (first group is group1 stocks, second is all other ones)
#'
#'
#' barplotTotLandingsPerSce (general=general,
#' type_of_column="disc",
#' selected="_selected_set1_",
#' selected_pops=c(0:38),
#' group1= c(0, 1, 2, 3, 11, 23, 24, 26, 30, 31, 32),
#' selected_scenarios=general$namefolderoutput,
#' scenarios_names=general$namefolderoutput,
#' nby=5, firsty="2015", lasty="2019",
#' a_width=3500, a_height=2000, black_and_white = FALSE, ylims = c(0,305))
#'
#'
#' }
barplotTotLandingsPerSce <- function(general=general,
type_of_column="pop", # or "disc"
selected="_selected_set1_",
selected_pops=c(0, 1, 2, 3, 11, 23, 24, 26, 30, 31, 32),
group1= c(0, 1, 2, 3, 11, 23, 24, 26, 30, 31, 32),
selected_scenarios=general$namefolderoutput,
scenarios_names=general$namefolderoutput,
nby=5, firsty="2015", lasty="2019",
a_width=3500, a_height=2000, black_and_white = TRUE, ylims = c(0,10)){
c.listquote <- function (...)
{
args <- as.list(match.call()[-1])
lstquote <- list(as.symbol("list"))
for (i in args) {
if (class(i) == "name" || (class(i) == "call" && i[[1]] !=
"list")) {
i <- eval(substitute(i), sys.frame(sys.parent()))
}
if (class(i) == "call" && i[[1]] == "list") {
lstquote <- c(lstquote, as.list(i)[-1])
}
else if (class(i) == "character") {
for (chr in i) {
lstquote <- c(lstquote, list(parse(text = chr)[[1]]))
}
}
else stop(paste("[", deparse(substitute(i)), "] Unknown class [",
class(i), "] or is not a list()", sep = ""))
}
return(as.call(lstquote))
}
#------------------------
reshape_per_pop_for_one_sce <- function(lst_loglike1=lst_loglike,
namesimu=paste("simu", 1:5, sep=''),
selected_pops=explicit_pops,
type_of_column="pop",
nby=5){
# filter lst_loglike to trash away the failed (i.e. non-complete) simus:
# detection according to the number of rows...
dd <- table(unlist(lapply(lst_loglike1, nrow)))
namesimu <- names(unlist(lapply(lst_loglike1, function(x) nrow(x)>=nby)))
print(namesimu)
# subset
refsimu <- namesimu[1]
lst_loglike <- lst_loglike1[ namesimu ]
# explicit pops??
pops <- selected_pops
idx_col <- grep(type_of_column, colnames(lst_loglike[[refsimu]]))
# check
#print("check if same numbers of row across the simus")
#print( lapply(lst_loglike1[ namesimu ], nrow) )
# first, reshape for pop
#browser()
for (i in 1:length(lst_loglike)){
lst_loglike[[ i ]] <- lst_loglike[[ i ]][, colnames(lst_loglike[[ i ]]) %in% c("year", paste0(type_of_column,".",pops)) ] [1:nby, (1:(length(selected_pops)+1))]
lst_loglike[[ i ]] <- cbind( reshape(lst_loglike[[ i ]], idvar="year", varying=list(2:ncol(lst_loglike[[ i ]])),
v.names=type_of_column, direction="long"), a_pop = rep(pops, each=nrow(lst_loglike[[ i ]])) )
}
# reshape...
res <- NULL
for(pop in pops){ # for each (explicit) pop
mat.sim1 <- matrix(unlist(lapply(lst_loglike[ namesimu ], function(x){
res <- try(x[x$a_pop==pop, type_of_column], silent=TRUE); if(class(res)=="try-error") res <- rep(NA, ncol(lst_loglike[[refsimu]])); res
})), nrow=nrow(lst_loglike[[refsimu]][lst_loglike[[refsimu]]$a_pop==pop,]) , byrow=FALSE)
colnames(mat.sim1) <- c(paste(type_of_column,"_", namesimu , sep=''))
mat.sim1 <- replace(mat.sim1, is.na(mat.sim1), 0)
mat.sim1 <- cbind.data.frame(lst_loglike[[refsimu]] [ lst_loglike[[1]]$a_pop==pop, c("year","a_pop") ], mat.sim1)
if(!is.null(res)) res <- rbind.data.frame(res, mat.sim1) else res <- mat.sim1
}
res <- res[!is.na(res$year),]
# get the median....
res <- cbind (res,
median=apply(res[,-c(1,2)], 1, quantile, probs=0.5),
mean= apply(res[,-c(1,2)], 1, mean)
)
return(res)
}
### 1- calls for reshaping the lst_popdyn1----------------
count <- 0
for (sce in selected_scenarios){
count <- count+1
lst_loglike <- get(paste("lst_loglike_agg_weight", selected, sce, sep=''), env=.GlobalEnv)
print(sce)
# aggregate landings weight PER YEAR
for (i in 1:length(lst_loglike)){
library(data.table)
loglike <- lst_loglike[[i]] # input
if(selected=="_met_") loglike <- loglike[loglike$metier==met,]
loglike$year <- unlist(lapply(strsplit(as.character(loglike$year.month), split="\\."), function(x) x[1]))
nm <- colnames(loglike)
idx.col <- grep(paste0(type_of_column,'.'), nm) # 40: DEBUG
DT <- data.table(loglike)
eq1 <- c.listquote( paste ("sum(",nm[idx.col],",na.rm=TRUE)",sep="") )
loglike.agg <- DT[,eval(eq1),by=list(year)]
loglike.agg <- data.frame( loglike.agg)
colnames(loglike.agg) <- c("year", paste(paste0(type_of_column,'.'), 0:(general$nbpops-1), sep=''))
loglike.agg <- loglike.agg[order(loglike.agg$year),] # order
lst_loglike[[i]] <- loglike.agg # output
}
# keep complete simu i.e. 5 years
lst_loglike <- lst_loglike [ names(lst_loglike) [lapply(lst_loglike, nrow)>=nby] ]
loglike_reshaped <- reshape_per_pop_for_one_sce(
lst_loglike1=lst_loglike,
namesimu=general$namesimu[[sce]],
selected_pops=selected_pops,
type_of_column=type_of_column,
nby=nby)
assign(paste("loglike_", sce, sep=''), loglike_reshaped)
}
### the code for the plot
all_sces_first_y <- data.frame(NULL)
all_sces_last_y <- data.frame(NULL)
count <- 0
for (sce in selected_scenarios){
count <- count+1
loglike <- get(paste("loglike_", sce, sep=''))
loglike_last_tstep <- loglike[loglike$year %in% lasty ,]
loglike_first_tstep <- loglike[loglike$year %in% firsty ,]
pop_names <- read.table(file.path(general$main.path.ibm, paste("pop_names_", general$namefolderinput, ".txt", sep='')), header=TRUE)
if(count==1) group1 <- as.character(pop_names[ pop_names[,1] %in% group1 ,2])
loglike_last_tstep$a_pop <- factor(loglike_last_tstep$a_pop)
levels(loglike_last_tstep$a_pop) <- pop_names[ pop_names[,1] %in% levels(loglike_last_tstep$a_pop) ,2]
loglike_first_tstep$a_pop <- factor(loglike_first_tstep$a_pop)
levels(loglike_first_tstep$a_pop) <- pop_names[ pop_names[,1] %in% levels(loglike_first_tstep$a_pop) ,2]
if(count==1) all_sces_first_y <- cbind.data.frame(a_pop=as.character(loglike_first_tstep[,"a_pop"]))
all_sces_first_y <- cbind.data.frame(all_sces_first_y, loglike_first_tstep[,"mean"])
colnames(all_sces_first_y)[ncol(all_sces_first_y)] <- scenarios_names[count]
if(count==1) all_sces_last_y <- cbind.data.frame(a_pop=as.character(loglike_last_tstep[,"a_pop"]))
all_sces_last_y <- cbind.data.frame(all_sces_last_y, loglike_last_tstep[,"mean"])
colnames(all_sces_last_y)[ncol(all_sces_last_y)] <- scenarios_names[count]
} # end sce
if(black_and_white){
# split into two groups...
levels(all_sces_first_y$a_pop) [levels(all_sces_first_y$a_pop) %in% group1] <- "GROUP1"
levels(all_sces_first_y$a_pop) [!levels(all_sces_first_y$a_pop) %in% c("GROUP1",group1)] <- "GROUP2"
levels(all_sces_last_y$a_pop) [levels(all_sces_last_y$a_pop) %in% group1] <- "GROUP1"
levels(all_sces_last_y$a_pop) [!levels(all_sces_last_y$a_pop) %in% c("GROUP1",group1)] <- "GROUP2"
all_sces_first_y <- aggregate(all_sces_first_y[,-1], list(all_sces_first_y[,1]), sum)
colnames(all_sces_first_y)[1] <- "a_pop"
all_sces_last_y <- aggregate(all_sces_last_y[,-1], list(all_sces_last_y[,1]), sum)
colnames(all_sces_last_y)[1] <- "a_pop"
# black and white version
some_colors <- c(grey(0.2),grey(0.5),grey(1))
the_density <- rep(500, length(some_colors))
} else{
# reorder stocks
all_sces_first_y$a_pop <- factor(all_sces_first_y$a_pop)
#species_order <- c(
# grep("nsea", all_sces_first_y$a_pop),
# grep("kask", all_sces_first_y$a_pop ),
# grep("kat", all_sces_first_y$a_pop),
# grep("3a22", all_sces_first_y$a_pop ),
# grep("3a2223", all_sces_first_y$a_pop ),
# grep("2224", all_sces_first_y$a_pop),
# grep("2532", all_sces_first_y$a_pop ),
# grep("2232", all_sces_first_y$a_pop )
# )
#species_order <- unique(species_order)
#all_sces_first_y <- all_sces_first_y[rev(species_order),]
#all_sces_last_y <- all_sces_last_y[rev(species_order),]
}
some_colors <- c("#a6cee3","#1f78b4","red","#b2df8a","green","#33a02c", "#fb9a99",grey(0.5),"#fdbf6f","#ff7f00","#cab2d6","#6a3d9a","#ffff99","#b15928", "black")
the_density <- rep(500, length(some_colors))
##-------A PLOT-----------
graphics.off()
if(type_of_column=="pop") namefileplot <- paste("landings_per_stock_per_scenario", selected, sep='')
if(type_of_column=="disc") namefileplot <- paste("discards_per_stock_per_scenario", selected, sep='')
tiff(file=file.path(general$main.path, general$namefolderinput, paste(namefileplot, ".tiff", sep="")),
width = a_width, height = a_height, compression="lzw",
units = "px", pointsize = 12, res=300)
par(mfrow=c(1,2))
par(mar=c(7,3,2,2))
par(oma=c(3,2,1,1))
rownames(all_sces_first_y) <- all_sces_first_y[,1]
rownames(all_sces_last_y) <- all_sces_last_y[,1]
mp_first_y <- barplot(as.matrix(all_sces_first_y[,-1])/1e6, las=2, ylim=ylims, xlab="", ylab="",
col =some_colors , density=the_density, legend=TRUE, axes = FALSE,axisnames = FALSE,
args.legend = list(x = "topright", bty = "n", ncol=2, cex=0.9))
text(mp_first_y, par("usr")[3], labels = colnames(all_sces_first_y[-1]), srt = 45, adj = 1, xpd = TRUE, cex = 1)
axis(2, las=2)
title ("First year", adj=0)
mtext("Scenario", 1, line=0, cex=1.5, outer=TRUE)
if(type_of_column=="pop") mtext(side=2, "Annual landings [*000 tons]", line=0, cex=1., outer=TRUE)
if(type_of_column=="disc") mtext(side=2, "Annual discards [*000 tons]", line=0, cex=1., outer=TRUE)
mp_last_y <- barplot(as.matrix(all_sces_last_y[,-1])/1e6, las=2, , ylim=ylims, xlab="", ylab="",
col =some_colors , density=the_density, legend=FALSE, axes = FALSE,axisnames = FALSE,
args.legend = list(x = "topright", bty = "n", ncol=2))
text(mp_last_y, par("usr")[3], labels = colnames(all_sces_last_y[-1]), srt = 45, adj = 1, xpd = TRUE, cex = 1)
axis(2, las=2)
title ("Last year", adj=0)
mtext("Scenario", 1, line=0, cex=1.5, outer=TRUE)
if(type_of_column=="pop") mtext(side=2,"Annual landings [*000 tons]",line=0, cex=1., outer=TRUE)
if(type_of_column=="disc") mtext(side=2,"Annual discards [*000 tons]",line=0, cex=1., outer=TRUE)
dev.off()
return()
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.