GitLab at IIASA

capri_sankey.R 20.8 KiB
Newer Older
Xinxin Yang's avatar
Xinxin Yang committed
# Script: filter different markets( just for cake oil market and sugar, cereals...)

#' function reads balance detailed, split it into "biofuels" and "non-biofuels"
#'
#' @param balance_detailed A data frame.
#' @param p_biofuels boolean.
#' @return links
#' @export
#'
prelinks <- function(balance_detailed, p_biofuels = TRUE){
  # balance_detailed <- balance_detailed_Scenario

  # add variable key, itemName and color from products
  balance_detailed <- left_join(balance_detailed, products[,c(2,3,9)], by = c(".i4" = "key"))

  # rename new variable
  colnames(balance_detailed)[13] <- "Commodities"
  colnames(balance_detailed)[1] <- "key"

  balance_detailed$Commodities <- sub("\\[.*?\\]", "", balance_detailed$Commodities)


  library(magrittr)
  # prod_list_fullname <- balance_detailed %>% select(Commodities) %>% extract2(1)

  # remove character begin with [ and end with ]
  # prod_list_fullname <- gsub("\\[.*?\\]", "", prod_list_fullname)
  # p_biofuels = T

  if (p_biofuels){
    selected_list = c("Rape seed",
                      "Soya seed",
                      "Sunflower seed",
                      "Rape seed oil",
                      "Soya oil",
Xinxin Yang's avatar
Xinxin Yang committed
                      "Sunflower seed oil",
Xinxin Yang's avatar
Xinxin Yang committed
                      "Soya cake",
                      "Sunflowe seed cake",
                      "Rape seed cake",
                      "Destilled dried grains from bio-ethanol processing",
                      "Pulses",
                      "Bio ethanol")
  }else{
    selected_list = c("Wheat","Grain maize","Barley", "Other cereals","Poultry meat",
                      "Sugar", "Fish and other acquatic products", "Rice milled")
  }

  #
  balance_detailed_ <- balance_detailed %>% as.data.frame()
  balance_detailed_selected <-  filter(balance_detailed_, balance_detailed_[,13] %in% selected_list)

  as_tibble(balance_detailed_selected)

  #rename "Destilled dried grains from bio-ethanol processing" -> "DDGS"
  # balance_detailed_selected <- balance_detailed_selected %>%
  #   mutate(Commodities = ifelse(Commodities!="Destilled dried grains from bio-ethanol processing",Commodities,"DDGS"))


  # filter
  balance_detailed_all <- balance_detailed_selected %>%
    select(-.i1)  %>% # remove region
    select(-.i5)  %>% # remove year: 2030
    select(-interv_ch)  %>%
    group_by(Commodities) %>%
    summarise_at(vars(supply:exports),sum,na.rm=TRUE)

  library(magrittr)
  prod_list_fullname <- balance_detailed_selected %>% select(Commodities) %>% extract2(1)

  #print(prod_list_fullname)
  vars_to_run <- prod_list_fullname


  left_side <- c("supply","imports")

  links_list <- list()
  #nodes_list <- list()

  vars_to_run <- vars_to_run[!vars_to_run %in%
                               c("Soya oil",
                                 "Soya cake",
                                 "Rape seed oil",
                                 "Rape seed cake",
Xinxin Yang's avatar
Xinxin Yang committed
                                 "Sunflower seed oil",
Xinxin Yang's avatar
Xinxin Yang committed
                                 "Sunflowe seed cake",
                                 "Bio ethanol",
                                 "Destilled dried grains from bio-ethanol processing"
                               )]
  # "Rape seed"      "Sunflower seed" "Soya seed"      "Pulses"
  for (i in vars_to_run){

    data_df <- reshape2::melt(balance_detailed_all)

    data_df_1 <- data_df %>% filter(Commodities == i)
    data_df_1_l <- data_df_1  %>% filter(variable %in% left_side )
    data_df_1_r <- data_df_1  %>% filter(!variable %in% left_side)
    data_df_1_l_agg <- data_df_1_l %>% group_by(Commodities) %>% summarise(value=sum(value))
    data_df_1_l_agg$variable <- paste0(i)


    first <- merge(data_df_1_l,data_df_1_l_agg,by="Commodities",all=TRUE)

    first <- first[,c("variable.x","variable.y","value.x")]
    colnames(first)[3] <- "value"


    second <- merge(data_df_1_l_agg,data_df_1_r,by="Commodities",all=TRUE)
    second <- second[,c("variable.x","variable.y","value.y")]
    colnames(second)[3] <- "value"

    if (p_biofuels){
      if (i != "Pulses"){
        if (i == "Soya seed") {produced  = c("Soya oil","Soya cake")}
        else if (i == "Rape seed") {produced  = c("Rape seed oil","Rape seed cake") }
Xinxin Yang's avatar
Xinxin Yang committed
        else if (i == "Sunflower seed") {produced  = c("Sunflower seed oil","Sunflowe seed cake") }
Xinxin Yang's avatar
Xinxin Yang committed
        else {print(i);break}

        data_df_proc_help <- balance_detailed_all %>% filter(Commodities %in% produced)
        data_df_proc_help <- reshape2::melt(data_df_proc_help)
        data_df_proc_help <- data_df_proc_help %>% filter(variable== "supply" )
        data_df_proc_help$variable <- paste(data_df_proc_help$Commodities)

        data_df_proc_help$Commodities <- i
        colnames(data_df_proc_help) <- c("variable.x","variable.y","value")

        second <- rbind(second,data_df_proc_help)
        second <- second %>% filter(variable.y!="processing")

        data_df_proc_help <- balance_detailed_all %>% filter(Commodities %in% produced)
        data_df_proc_help <- reshape2::melt(data_df_proc_help)
        data_df_proc_help <- data_df_proc_help %>% filter(variable!="supply")
        colnames(data_df_proc_help) <- c("variable.x","variable.y","value")
        second <- rbind(second,data_df_proc_help)
        second_imports <- second %>% filter(variable.y=="imports")
        second <- second %>% filter(!variable.y=="imports")
        colnames(second_imports) <- c("variable.y","variable.x","value")
        second <- rbind(second,second_imports)
      }

    }

    combined <- rbind(first,second)
    links_list[[i]] <- combined


    links <- rbindlist(links_list)
    links <- unique(links)

  }

  if(p_biofuels){

    # links[grepl("Bio",variable.x) | grepl("Bio",variable.y)]
    setnames(links,old=c("variable.x","variable.y"),new=c("source","target"))
    # links[grepl("Bio",source) | grepl("Bio",target)]

    i <- c("Bio ethanol","Destilled dried grains from bio-ethanol processing")

    data_df_proc_help <- balance_detailed_all %>% filter(Commodities %in% i)
    data_df_proc_help <- reshape2::melt(data_df_proc_help)
    data_df_proc_help <- data_df_proc_help %>% filter(variable=="supply")
    data_df_proc_help$variable <- paste(data_df_proc_help$Commodities)

    data_df_proc_help$Commodities <- "biofuels"
    colnames(data_df_proc_help) <- c("variable.x","variable.y","value")

    data_df <- reshape2::melt(balance_detailed_all)

    data_df_1 <- data_df[data_df$Commodities %in% i,]
    data_df_1_l <- data_df_1[data_df_1$variable %in% left_side,]
    data_df_1_l$variable <- as.character(data_df_1_l$variable)
    data_df_1_l <- data_df_1_l %>% filter(variable!="supply")
    data_df_1_r <- data_df_1[!data_df_1$variable %in% left_side,]
    data_df_1_r$variable <- as.character(data_df_1_r$variable)
    data_df_1_l_agg <- data_df_1_l%>%group_by(Commodities) %>% summarise(value=sum(value))
    data_df_1_l_agg$variable <- paste0(i)
    first <- merge(data_df_1_l,data_df_1_l_agg,by="Commodities",all=TRUE)
    first <- first[,c("variable.x","variable.y","value.x")]
    colnames(first)[3] <- "value"
    second <- merge(data_df_1_l_agg,data_df_1_r,by="Commodities",all=TRUE)
    second <- second[,c("variable.x","variable.y","value.y")]
    colnames(second)[3] <- "value"
    second <- rbind(second,data_df_proc_help)


    second <- second %>% filter(!variable.y %in% c("processing","biofuels"))


    #nodes <- rbind(nodes,data.frame("name"=i))
    combined <- rbind(first,second)
    setnames(combined,old=c("variable.x","variable.y"),new=c("source","target"))

    links <- rbind(links,combined)

  }
  else{
    setnames(links,old=c("variable.x","variable.y"),new=c("source","target"))

  }

  # remove 0 value and <=900
Xinxin Yang's avatar
Xinxin Yang committed
  links[, value := ifelse(value==0,NA,value)]
Xinxin Yang's avatar
Xinxin Yang committed
  # links[,value:=ifelse(value<=50,NA,value)]
  links<-na.omit(links)

  return(links)
}

#' function calulates sum of links
#'
#' @param balance_detailed_links links.
#' @return sum
#'
linksum <- function(balance_detailed_links){

  # the total value of different products from source
  linkcopysum<- aggregate(balance_detailed_links$value,
                          by=list(Category=balance_detailed_links$source),
                          FUN=sum)

  # add a column: the total value of different prodcuts from target
  linkcopysum <- merge(linkcopysum,
                       aggregate(balance_detailed_links$value,
                                 by=list(Category=balance_detailed_links$target),
                                 FUN=sum),
                       by="Category")

  return(linkcopysum)

}

#'  function Replaces specific characters within strings
#' @param old old column names.
#' @param new new column names.
#' @param df replace data frame
#' @return df
#'
#'
#'
replacestr <- function(old, new, df){
  df$label <- gsub(old, new, df$label)
  return(df)
}

#'  function get links and nodes for drawing a sankey diagram
#' @param baseline Baseline balance market.
#' @param scenario Scenario balance market.
#' @param p_baseline boolean.
#' @param fixedNodePosition boolean.
Xinxin Yang's avatar
Xinxin Yang committed
#'
#' @return links and nodes
#' @export
#'
#'
links_nodes <- function(baseline, scenario, p_baseline, fixedNodePosition = TRUE){
Xinxin Yang's avatar
Xinxin Yang committed

  linksum_baseline = linksum(baseline)
  linksum_scenario = linksum(scenario)

  # percentage diferences between scenario and baseline
  linksum_scenario$percent_source = (linksum_scenario$x.x*100/linksum_baseline$x.x)-100
  linksum_scenario$percent_target = (linksum_scenario$x.y*100/linksum_baseline$x.y)-100

  #links = baseline
  if (p_baseline){
    links = baseline
  }else{links= scenario}


  nodes <- data.frame("name" = c(unique(c(as.character(links$source),links$target)) ))

  # product: names -> number
  links$source <- sapply(links$source,function(x)which(nodes$name %in% x))
  links$source  <- links$source -1
  links$target <- sapply(links$target,function(x)which(nodes$name %in% x))
  links$target <- links$target-1

  if (p_baseline){
    links$group <- paste0(gsub(" ", "_", baseline$source),"->", gsub(" ", "_", baseline$target))
  }else{
    links$group <- paste0(gsub(" ", "_", scenario$source),"->", gsub(" ", "_", scenario$target))
  }

  source_baseline <- baseline %>% as.data.frame() %>%
    filter(source %in% c("imports", "supply")) %>%
    group_by(source) %>%
    summarise(sum(value))%>%
    rename(name=source, sum = `sum(value)`)

  source_scenario <- scenario %>% as.data.frame() %>%
    filter(source %in% c("imports", "supply")) %>%
    group_by(source) %>%
    summarise(sum(value))%>%
    rename(name=source, sum = `sum(value)`)

  # calculate percent
  source_scenario$percent<- (source_scenario$sum/source_baseline$sum)*100-100

  target_baseline <- baseline %>% as.data.frame() %>%
    filter(!(target %in% c("imports", "supply"))) %>%
    group_by(target) %>%
    summarise(sum(value)) %>%
    rename(name=target, sum = `sum(value)`)


  target_scenario <- scenario %>% as.data.frame() %>%
    filter(!(target %in% c("imports", "supply"))) %>%
    group_by(target) %>%
    summarise(sum(value)) %>%
    rename(name=target, sum = `sum(value)`)

  # calculate percent
  target_scenario$percent<- (target_scenario$sum/target_baseline$sum)*100-100

  all_baseline <- rbind(source_baseline, target_baseline)
  all_scenario <- rbind(source_scenario, target_scenario)

  # if p_basline != TRUE -> add percent for nodes
  if(!p_baseline){
    # nodes$label <- ifelse(nodes$percent_target=="",
    #                       as.character(nodes$label),
    #                       paste0(nodes$label, " (", nodes$percent_target,"%", ")"))

    all_scenario$label <- paste0(all_scenario$name, ": ",
                                 format(round((all_scenario$sum)/1000,1), decimal.mark = ","),
                                 "K (",round(all_scenario$percent,1), "%)" )
    all <- all_scenario
    #nodes$levels <- 1:3
  }else{
    all_baseline$label <- paste0(all_baseline$name, ": ",
                                 format(round((all_baseline$sum)/1000,1), decimal.mark = "," ), "K" )
    all <- all_baseline }


  nodes_ <- left_join(nodes, all)

  rm(nodes)
  nodes <- nodes_

  nodes$rank=c(1:nrow(nodes))

  nodes = merge(nodes, linksum_scenario[,c("Category","percent_target")],
                by.x="name",
                by.y="Category",
                all.x=TRUE)

  nodes = nodes[order(nodes$rank,decreasing = FALSE),]

  nodes$percent_target = round(nodes$percent_target,0)


  nodes$percent_target <- ifelse(nodes$name=="biofuels","",nodes$percent_target)
  nodes[is.na(nodes)]<-""

  selectedcolor <- bind_rows(products %>% select(c(itemName,color)),
                             dim5s %>% select(c(itemName,color)))

  selectedcolor$itemName <- sub("\\[.*?\\]", "", selectedcolor$itemName)
  nodes <-  left_join(nodes, selectedcolor, by = c("name" = "itemName"))


  #print(nodes)
  nodes <- nodes %>% select(c(rank,name,label,color))


  # add links color

  links_temp <- links
  links_temp$color <- ""
  colorlink <- c("Rape", "Sunflowe", "Soya", "Pulses",
                 "Wheat","Grain","Barley",
                 "Other","Poultry",
                 "Sugar", "Fish",
                 "Rice")


  # get colors from nodes

  color_nodes <- c("Wheat|Grain|Barley|Other|Poultry|Sugar|Fish|Rice|Rape seed$|Sunflowe seed|Soya seed$|Pulses")

  getcolor <-  nodes[grep(color_nodes, nodes$name), ] %>% select(c(name, color))


  library(colorspace)
  for (colorcol in colorlink){
    temp <- getcolor %>%
      filter(grepl(colorcol, getcolor$name))
Xinxin Yang's avatar
Xinxin Yang committed
    # message(colorcol)
    # print(links_temp[grepl(colorcol, group), ])
    # print(links_temp[grep(colorcol, group), ])
Xinxin Yang's avatar
Xinxin Yang committed
    links_temp[grep(colorcol, group), ]$color <- paste0(colorspace::lighten(temp$color, 0.05),"80")

  }
Xinxin Yang's avatar
Xinxin Yang committed
  # print(links_temp)
Xinxin Yang's avatar
Xinxin Yang committed


  links <- links_temp %>%
    mutate(color = ifelse(color=="","rgba(123,123,123,0.1)",color))

  # Replace specific characters within strings
  old=c("supply",
        "human_cons",
        "processing",
        "biofuels",
        "feed",
        "imports",
        "exports",
        "Destilled dried grains from bio-ethanol Processing")
  new=c("Production",
        "Human Consumption",
        "Processing",
        "Biofuels",
        "Feed",
        "Imports",
        "Exports",
        "DDG")
  for (i in 1:length(old)){
    nodes <- replacestr(old[i],new[i], nodes)
  }

  # add nodes position

  if (fixedNodePosition) {nodes_lv <- nodes_position(nodes);nodes <- nodes %>% left_join(nodes_lv)}


Xinxin Yang's avatar
Xinxin Yang committed
  return(list(links,nodes))
}
#' sets the position of nodes
#' @param nodes nodes for drawing a sankey diagram.
#' @return nodes with fixed position
#'
nodes_position <- function(nodes){
  # nodes$level_x = ""
  # nodes$level_y = ""

  # if( "Pulses" %in% nodes$name){
  #   # for oil cake
  #   nodes <- nodes %>%
  #     mutate(level_x = ifelse(name %in% c("supply", "imports"), 0.01,
  #                             ifelse(name %in% c("human_cons", "feed", "exports", "processing"),0.8,
  #                                    ifelse(name %in% c("Rape seed", "Soya seed", "Sunflower seed","Pulses"), 0.2,
  #                                           ifelse(name %in% c("Rape seed oil", "Soya oil", "Sunflower seed oil"),0.35,
  #                                                  ifelse(name %in% c("Rape seed cake", "Soya cake", "Sunflowe seed cake", "biofuels"), 0.5, 0.65))) ))) %>%
  #     mutate(level_y = ifelse(name %in% c("supply", "imports"), c(0.01,0.8),
  #                             ifelse(name %in% c("human_cons","Bio ethanol", "feed", "exports", "processing"),c(0.01,0.3,0.5,0.7,0.9),
  #                                    ifelse(name %in% c("Rape seed", "Soya seed","Sunflower seed","Pulses"), c(0.1,0.3,0.6,0.9),
  #                                           ifelse(name %in% c("Rape seed oil", "Soya oil", "Sunflower seed oil"),c(0.1, 0.4, 0.7),
  #                                                  ifelse(name %in% c("Rape seed cake", "Soya cake", "Sunflowe seed cake", "biofuels"), c(0.1, 0.3, 0.5, 0.7), 0.5))) )))
  #
  #
  # }else{
  #   nodes <- nodes %>%
  #     # for rest
  #     mutate(level_x = ifelse(name %in% c("supply", "imports"),0.01, ifelse(name %in% c("human_cons","biofuels", "feed", "exports", "processing"),0.8,0.4)))%>%
  #     mutate(level_y = ifelse(name %in% c("supply", "imports"),c(0.1,0.8), ifelse(name %in% c("human_cons","biofuels", "feed", "exports", "processing"),c(0.01,0.3,0.5,0.7,0.9),c(0.01,0.1,0.2,0.3,0.5,0.6,0.7,0.8))))
  #
  # }
Xinxin Yang's avatar
Xinxin Yang committed

  if( "Pulses" %in% nodes$name){

    nodes_level <- read.table(text="name	level_x	level_y
                              'supply' 0.01 0.1
                              'imports' 0.01 0.5
                              'Rape seed'	0.2	0.01
Xinxin Yang's avatar
Xinxin Yang committed
                              'Rape seed oil'	0.35	0.2
                              'Sunflowe seed cake'	0.5	0.5
                              'Soya seed'	0.2	0.3
Xinxin Yang's avatar
Xinxin Yang committed
                              'Soya oil'	0.35	0.3
                              'Bio ethanol'	0.65	0.14
                              'biofuels'	0.5	0.7
                              'feed'	0.8	0.3
                              'processing'	0.8	0.8
                              'Rape seed cake'	0.5	0.1
Xinxin Yang's avatar
Xinxin Yang committed
                              'Sunflower seed'	0.2	0.45
                              'Sunflower seed oil'	0.35	0.55
                              'Soya cake'	0.5	0.3
Xinxin Yang's avatar
Xinxin Yang committed
                              'Pulses'	0.2	0.6
                              'Destilled dried grains from bio-ethanol processing'	0.65	0.6
                              'human_cons'	0.8	0.01
                              'exports'	0.8	0.5", header=TRUE)}
  else{nodes_level <- read.table(text = "name	level_x	level_y
                                 'supply' 0.01 0.2
                                 'imports' 0.01 0.5
                                 'biofuels'	0.8	0.7
                                 'feed'	0.8	0.3
                                 'processing'	0.8	0.8
                                 'human_cons'	0.8	0.01
                                 'exports'	0.8	0.5
                                 'Grain maize' 0.4 0.01
                                 'Barley' 0.4 0.15
Xinxin Yang's avatar
Xinxin Yang committed
                                 'Other cereals' 0.4 0.25
                                 'Wheat' 0.4 0.35
                                 'Poultry meat' 0.4 0.5
                                 'Rice milled' 0.4 0.6
                                 'Sugar' 0.4 0.7
                                 'Fish and other acquatic products' 0.4 0.8"
                            , header=TRUE)}




  return(nodes_level)
Xinxin Yang's avatar
Xinxin Yang committed
}

#' function draws sankey diagram and saves html pages for sankey
#'
#'
#' @param data a data frame object has two lists, which contains the links between the nodes and the nodes, the nodes has node id and properties of the nodes.links should have include the Source and Target for each link. An optional Value variable can be included to specify how close the nodes are to one another.
#'   If no ID is specified then the nodes must be in the same order as the Source variable column in the Links data frame. Currently only grouping variable is allowed.
Xinxin Yang's avatar
Xinxin Yang committed
#' @param p_baseline boolean.
#' @param png boolean. if TRUE, the sankey chart will be saved.
#'
#' @return sankey diagram.
#' @export
Xinxin Yang's avatar
Xinxin Yang committed
plot_sankey <- function(data, p_baseline, png, outdata.dir){

  links <- data[[1]]
  nodes <- data[[2]]
Xinxin Yang's avatar
Xinxin Yang committed
  N_ethanol <- nodes %>%
                  mutate_all(list(~ . %in% "Bio ethanol")) %>%
                  as.matrix %>%
                  which(., arr.ind = TRUE) %>%
                  nrow()

  N_cereals <- nodes %>%
    mutate_all(list(~ . %in% "Other cereals")) %>%
    as.matrix %>%
    which(., arr.ind = TRUE) %>%
    nrow()

  if (N_ethanol==2) NAME <- "comparison_cake_oil"
  else if(N_cereals == 2 ) NAME <- "comparison_cereals"
  else if(p_baseline){if("Bio ethanol" %in% nodes$name){ NAME <- "baseline"; print("Biofuels")}else{ NAME <- "baseline_Rest"}}
  else{ if("Bio ethanol" %in% nodes$name){ NAME <- "scenario"; print("Biofuels")} else{ NAME <- "scenario_Rest"}}
Xinxin Yang's avatar
Xinxin Yang committed

  fig <- plot_ly(
Xinxin Yang's avatar
Xinxin Yang committed
    width = 1000,
Xinxin Yang's avatar
Xinxin Yang committed
    type = "sankey",
Xinxin Yang's avatar
Xinxin Yang committed
    # domain = list(
    #   x =  c(0,1),
    #   y =  c(0,1)
    # ),
Xinxin Yang's avatar
Xinxin Yang committed
    orientation = "h",
    arrangement = 'snap',
    # arrangement = 'perpendicular',
Xinxin Yang's avatar
Xinxin Yang committed
    valueformat = ".0f",
    valuesuffix = " tonnes",

    node = list(
      label = nodes$label,
      color = nodes$color,
      x = nodes$level_x,
      y= nodes$level_y,
      pad = 20,
      thickness = 15,
      line = list(
        color = "blau",
        width = 1
      )
    ),

    link = list(
      source = links$source,
      target = links$target,
      value =  links$value,
      color = links$color
      #label =  links$group
    )
  )  %>% layout(
    autosize = F,
    margin = list(l=10, r=10, b=50, t=50, pad=40),
    title = "Market Balance in Tonnens",
    font = list(size = 13, family ="Arial"),
    xaxis = list(showgrid = F, zeroline = F),
    yaxis = list(showgrid = F, zeroline = F)
  )

  #
  library(htmlwidgets)
  print(paste0("Export sankey diagram: ",paste0(outdata.dir,"/",NAME,"_fig.html")))
  saveWidget(fig, paste0(outdata.dir,"/",NAME,"_fig.html"), selfcontained = F, libdir = "lib")
  library(webshot)
  #you convert it as png
  if(png){
Xinxin Yang's avatar
Xinxin Yang committed
  webshot(paste0(outdata.dir,"/",NAME,"_fig.html"),paste0(outdata.dir,"/",NAME,"_fig.png"),
Xinxin Yang's avatar
Xinxin Yang committed
  vwidth = 1000, vheight = 600)}
Xinxin Yang's avatar
Xinxin Yang committed

  return(fig)