reconciliation_USGS_BACI.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### **Mapping Copper Flows And Environmental Impacts For The Dutch Energy Transition**\n",
    "\n",
    "This notebook has been created to reconcile production data for copper provided the U.S. Geological Survey as well as global trade data and link it to impact factors derived from Life Cycle Assessments. This was done as part of a study project in the course \"Sustainability Challenge\" of the Master of Science \"Industrial Ecology\" at the Leiden University (LU) and Delft University of Technology (TUD). The project was commissioned by PBL Netherlands Environmental Assessment Agency (PBL). PBL, LU and TUD are not responsible for the content of this notebook. \n",
    "\n",
    "**Disclaimer**: \n",
    "- The code is not cleaned and optimised. There are a lot of small issues that are unresolved. Please consider it work in progress. Currently, we do not have time to pursue it further, so we would be happy about anyone continuing our work.  \n",
    "- In principal, the script can be adapted to other metals and countries with minimal changes. However, currently the code is not self-explanatory! If you are interested in  using it, I (Aaron) would be happy to answer any questions that might help with adapting the script. \n",
    "- The limitations of the approach are described in more detail in the report that can also be found in the github repository. \n",
    "- The Sankey currently show the impact \"flows\" for each stage separately. In the end, we decided to adapt them for the report so that the aggregate impacts along the supply chain are shown. The code for this change is not part of this notebook, but can be made available. \n",
    "\n",
    "**Created by**: Aaron Paris, Christina Drotenko, Alessia Linares Capurro, Kevin Heideman, Adrien Perello-y-bestard\n",
    "\n",
    "**Contact**: a.paris@umail.leidenuniv.nl"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "### **1. Import of modules**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import plotly.graph_objects as go"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "### **2. User input / decisions**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Inputs for the Material Flow Analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Shares of NL\n",
    "NL_share = 0.00808 # share of global refined copper going into NL ultimately (also in form of intermediates etc.) \n",
    "NL_energy_share = 0.2506 # share of copper in NL used for the energy sector\n",
    "\n",
    "# Year\n",
    "year = 2019"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Visualisation choices"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Threshold for display (does not reduce the total flows shown, but the aggregation of countries in \"Other\")\n",
    "relevance_threshold = 0.04 # threshold in %, this is the actual threshold for the USGS values after the impacts\n",
    "\n",
    "# Scaling\n",
    "scale_to_NL = True # if True than all values will be scaled for NL and the energy sector, if false the global values are shown\n",
    "\n",
    "# Choice of impact category to display\n",
    "# display = 'copper'\n",
    "display = 'ghg'\n",
    "# display = 'land_use'\n",
    "# display = 'water'\n",
    "\n",
    "# Sankey style\n",
    "highlight_LMIC = True # True means that the sankey will have colours grouped by high income and LMIC, false that colours are unique for each country. \n",
    "sankey_labels_absolute = True # if True, absolute numbers are shown in the sankey, otherwise percent"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "### **3. Import and preparation of trade data**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Trade data import (Source: http://www.cepii.fr/CEPII/en/bdd_modele/bdd_modele_item.asp?id=37)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# trade data\n",
    "trade_data = pd.read_csv(\"data/BACI_HS92_V202301/BACI_HS92_Y\"+str(year)+\"_V202301.csv\")\n",
    "country_codes = pd.read_csv(\"data/BACI_HS92_V202301/country_codes_V202301.csv\")\n",
    "product_codes = pd.read_csv(\"data/BACI_HS92_V202301/product_codes_HS92_V202301.csv\")\n",
    "\n",
    "# renaming columns\n",
    "trade_data.rename(columns={\n",
    "                    't': 'year', \n",
    "                    'i': 'exporter_code', \n",
    "                    'j': 'importer_code', \n",
    "                    'k': 'hscode', \n",
    "                    'v': 'kUSD_current', \n",
    "                    'q': 'tonnes'},\n",
    "                    inplace = True)\n",
    "\n",
    "# Cleaning product dataframe for merging\n",
    "product_codes = product_codes.drop(product_codes[product_codes['code'] == '9999AA'].index)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Defining product codes to include and filtering the trade data accordingly"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Manual input of stages, hscodes to filter for and copper content of the products\n",
    "copper_stages = pd.DataFrame({\n",
    "                            'hscode': ['2603', '2620', '7401', '7402'],\n",
    "                            'stage': ['mined', 'mined', 'smelted', 'smelted'],\n",
    "                            'content': [0.3, 0.1, 0.6, 0.95]\n",
    "                            })\n",
    "\n",
    "# Changing the hscodes for merging later\n",
    "trade_data['hscode'] = trade_data['hscode'].astype(str)\n",
    "trade_data['hscode_short'] = trade_data['hscode'].str[:4]\n",
    "copper_stages['hscode'] = copper_stages['hscode'].astype(str)\n",
    "\n",
    "#%% Filtering for hscodes\n",
    "copper_trade = trade_data[trade_data['hscode_short'].isin(copper_stages['hscode'])]\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Cleaning the dataframe, adding country and product names"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\2495299667.py:3: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  copper_trade.drop('year', axis = 1, inplace = True)\n"
     ]
    }
   ],
   "source": [
    "#%% Cleaning the resulting dataframe\n",
    "copper_trade.reset_index(drop=True, inplace = True)\n",
    "copper_trade.drop('year', axis = 1, inplace = True)\n",
    "copper_trade = copper_trade[copper_trade['tonnes'] != '           NA']\n",
    "\n",
    "#%% Adding country names, ISO3 codes and product names\n",
    "copper_trade = pd.merge(copper_trade, country_codes[['country_code', 'country_name_full', 'iso_3digit_alpha']], left_on='exporter_code', right_on='country_code', how='left')\n",
    "copper_trade = pd.merge(copper_trade, country_codes[['country_code', 'country_name_full', 'iso_3digit_alpha']], left_on='importer_code', right_on='country_code', how='left')\n",
    "copper_trade = pd.merge(copper_trade, product_codes, left_on='hscode', right_on='code', how='left')\n",
    "\n",
    "# Cleaning and renaming\n",
    "copper_trade.drop(['exporter_code', 'importer_code', 'country_code_x', 'country_code_y', 'code'], axis = 1, inplace = True)\n",
    "\n",
    "copper_trade.rename(columns={                \n",
    "                    'iso_3digit_alpha_x': 'exporter_ISO3', \n",
    "                     'country_name_full_x': 'exporter',\n",
    "                     'iso_3digit_alpha_y': 'importer_ISO3', \n",
    "                     'country_name_full_y': 'importer',\n",
    "                     'description': 'product'},\n",
    "                     inplace = True)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculating the copper content of the trade flows"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "#%% Adding stages and copper content\n",
    "copper_trade = pd.merge(copper_trade, copper_stages, left_on='hscode_short', right_on='hscode', how='left')\n",
    "\n",
    "# Calculating the copper content of each flow\n",
    "copper_trade['tonnes_copper'] = copper_trade['tonnes'].astype(float) * copper_trade['content']\n",
    "\n",
    "# Cleaning and renaming\n",
    "copper_trade.drop(['hscode_y', 'content'], axis = 1, inplace = True)\n",
    "copper_trade.rename(columns={'hscode_x': 'hscode'}, inplace = True)\n",
    "\n",
    "# Changing order of columns\n",
    "copper_trade = copper_trade[['exporter_ISO3', 'exporter', 'importer_ISO3', 'importer', 'stage', 'hscode', 'hscode_short', 'product', 'kUSD_current', 'tonnes', 'tonnes_copper']]\n",
    "\n",
    "# Changing names for countries in the copper trade data\n",
    "copper_trade = copper_trade.replace({'Democratic Republic of the Congo': 'DRC', \n",
    "                                        'Russian Federation': 'Russia', \n",
    "                                        'USA, Puerto Rico and US Virgin Islands': 'USA'})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "### **4. Import and preparation of USGS data**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Import of USGS production data (Source: https://www.usgs.gov/centers/national-minerals-information-center/copper-statistics-and-information). Only primary production is included. For those values for which there is no differentiation in the USGS data, it is assumed to be primary. The USGS data was cleaned by hand in the Excel so that it is machine readable. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "usgs_mining = pd.read_excel(\"data/usgs_data.xlsx\", sheet_name = 'mining')[['country',year]]\n",
    "usgs_mining.rename(columns={year:'mining_tonnes'}, inplace=True)\n",
    "\n",
    "usgs_smelting = pd.read_excel(\"data/usgs_data.xlsx\", sheet_name = 'smelting')[['country',year]]\n",
    "usgs_smelting.rename(columns={year:'smelting_tonnes'}, inplace=True)\n",
    "\n",
    "usgs_refining = pd.read_excel(\"data/usgs_data.xlsx\", sheet_name = 'refining')[['country',year]]\n",
    "usgs_refining.rename(columns={year:'refining_tonnes'}, inplace=True)\n",
    "\n",
    "# Merging all stages\n",
    "usgs_data = pd.merge(usgs_mining, usgs_smelting, on='country', how='outer').merge(usgs_refining, on='country', how='outer')\n",
    "usgs_data = usgs_data.sort_values(by='country', ascending= True)\n",
    "usgs_data = usgs_data.fillna(0)\n",
    "\n",
    "# Calculating shares of countries in each stage\n",
    "usgs_data['mining_share'] = usgs_data['mining_tonnes'] / usgs_data['mining_tonnes'].sum()\n",
    "usgs_data['smelting_share'] = usgs_data['smelting_tonnes'] / usgs_data['smelting_tonnes'].sum()\n",
    "usgs_data['refining_share'] = usgs_data['refining_tonnes'] / usgs_data['refining_tonnes'].sum()\n",
    "\n",
    "# Changing names for countries in the USGS data\n",
    "usgs_data = usgs_data.replace({'Democratic Republic of the Congo': 'DRC', \n",
    "                                     'Korea, Republic of': 'Republic of Korea', \n",
    "                                     'United States': 'USA'})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculating the domestic flows from USGS data. The underlying assumption is that countries process as much of the former stage as possible. If 700,000 tonnes are mined, for example, and 400,000 tonnes smelted, the flow from mining to smelting will be 400,000 tonnes and 300,000 tonnes of mined metal need to be exported. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Applying threshold from above to filter the USGS data for relevance\n",
    "usgs_size = usgs_data[usgs_data[['mining_share', 'smelting_share', 'refining_share']].max(axis=1) >= 0.001]\n",
    "\n",
    "# Aggregating countries below the threshold to \"Other\"\n",
    "usgs_size_below = usgs_data[usgs_data[['mining_share', 'smelting_share', 'refining_share']].max(axis=1) < 0.001]\n",
    "usgs_size_below = usgs_size_below.sum(axis=0)  # axis=0 sums along columns, axis=1 would sum along rows\n",
    "usgs_size_below = usgs_size_below.to_frame().T\n",
    "usgs_size_below.iloc[0,0] = 'a_Other' # named with a_ so that it is sorted after all other values but above the black boxes\n",
    "\n",
    "# Adding the \"Other\" category to usgs_size\n",
    "usgs_size = pd.concat([usgs_size,usgs_size_below])\n",
    "usgs_size.reset_index(drop=True, inplace=True)\n",
    "\n",
    "# building dataframe in sankey format for  domestic flows\n",
    "flows_domestic = usgs_size.copy()\n",
    "flows_domestic['mining_to_smelting_domestic'] = flows_domestic[['mining_tonnes', 'smelting_tonnes']].min(axis=1)\n",
    "flows_domestic['smelting_to_refining_domestic'] = flows_domestic[['smelting_tonnes', 'refining_tonnes']].min(axis=1)\n",
    "flows_domestic = flows_domestic[['country', 'mining_to_smelting_domestic', 'smelting_to_refining_domestic']]\n",
    "\n",
    "# Mined copper\n",
    "flows_domestic_mined = flows_domestic[['country', 'mining_to_smelting_domestic']]\n",
    "flows_domestic_mined['target_country'] = flows_domestic_mined['country']\n",
    "flows_domestic_mined['source_stage'] = 'mining'\n",
    "flows_domestic_mined['target_stage'] = 'smelting'\n",
    "flows_domestic_mined.rename(columns={'country': 'source_country', 'mining_to_smelting_domestic': 'tonnes_copper'}, inplace=True)\n",
    "flows_domestic_mined = flows_domestic_mined[['source_stage', 'source_country', 'target_stage', 'target_country', 'tonnes_copper']]\n",
    "\n",
    "# Smelted copper\n",
    "flows_domestic_smelted = flows_domestic[['country', 'smelting_to_refining_domestic']]\n",
    "flows_domestic_smelted['target_country'] = flows_domestic_smelted['country']\n",
    "flows_domestic_smelted['source_stage'] = 'smelting'\n",
    "flows_domestic_smelted['target_stage'] = 'refining'\n",
    "flows_domestic_smelted.rename(columns={'country': 'source_country', 'smelting_to_refining_domestic': 'tonnes_copper'}, inplace=True)\n",
    "flows_domestic_smelted = flows_domestic_smelted[['source_stage', 'source_country', 'target_stage', 'target_country', 'tonnes_copper']]\n",
    "\n",
    "# Mined and smelted copper\n",
    "flows_domestic = pd.concat([flows_domestic_mined,flows_domestic_smelted])\n",
    "flows_domestic = flows_domestic[(flows_domestic['tonnes_copper'] != 0)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculating the resulting foreign flows (imports/exports)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\1749641001.py:9: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  flows_foreign_mined.rename(columns={'mining_to_smelting_balance':'balance'}, inplace=True)\n",
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\1749641001.py:15: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  flows_foreign_smelted.rename(columns={'smelting_to_refining_balance':'balance'}, inplace=True)\n"
     ]
    }
   ],
   "source": [
    "#%% Calculating the necessary foreign flows\n",
    "flows_foreign = usgs_size.copy()\n",
    "flows_foreign['mining_to_smelting_balance'] = flows_foreign['smelting_tonnes'] - flows_foreign['mining_tonnes']\n",
    "flows_foreign['smelting_to_refining_balance'] = flows_foreign['refining_tonnes'] - flows_foreign['smelting_tonnes']\n",
    "flows_foreign = flows_foreign[['country', 'mining_to_smelting_balance', 'smelting_to_refining_balance']]\n",
    "\n",
    "# Mined copper\n",
    "flows_foreign_mined = flows_foreign[['country', 'mining_to_smelting_balance']]\n",
    "flows_foreign_mined.rename(columns={'mining_to_smelting_balance':'balance'}, inplace=True)\n",
    "flows_foreign_mined['source_stage'] = 'mining'\n",
    "flows_foreign_mined['target_stage'] = 'smelting'\n",
    "\n",
    "# Smelted copper\n",
    "flows_foreign_smelted = flows_foreign[['country', 'smelting_to_refining_balance']]\n",
    "flows_foreign_smelted.rename(columns={'smelting_to_refining_balance':'balance'}, inplace=True)\n",
    "flows_foreign_smelted['source_stage'] = 'smelting'\n",
    "flows_foreign_smelted['target_stage'] = 'refining'\n",
    "\n",
    "# Mined and smelted copper\n",
    "flows_foreign = pd.concat([flows_foreign_mined,flows_foreign_smelted])\n",
    "flows_foreign.reset_index(drop=True, inplace=True)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "### **5. Reconciliation of USGS and trade data**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculating necessary import / export balancing flows.  \n",
    "\n",
    "4 cases are possible: \n",
    "1. import of mining necessary because not enough in smelting - positive value\n",
    "2. export of mining necessary because too much in smelting - negative value\n",
    "3. import of smelted necessary because not enough in refining - positive value\n",
    "4. export of smelted necessary because too much in refining - negative value"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\1150928826.py:5: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  copper_trade_cleaned.loc[copper_trade_cleaned[\"stage\"] == \"mined\", \"exporter_stage\"] = 'mining'\n",
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\1150928826.py:8: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  copper_trade_cleaned.loc[copper_trade_cleaned[\"stage\"] == \"mined\", \"importer_stage\"] = 'smelting'\n",
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\1150928826.py:12: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  copper_trade_cleaned.rename(columns={'exporter': 'source_country',\n"
     ]
    }
   ],
   "source": [
    "# Initiating dataframe for Sankey with relevant columns from the trade data\n",
    "copper_trade_cleaned = copper_trade[['stage', 'exporter', 'importer', 'tonnes_copper']]\n",
    "\n",
    "# Adding stages\n",
    "copper_trade_cleaned.loc[copper_trade_cleaned[\"stage\"] == \"mined\", \"exporter_stage\"] = 'mining'\n",
    "copper_trade_cleaned.loc[copper_trade_cleaned[\"stage\"] == \"smelted\", \"exporter_stage\"] = 'smelting'\n",
    "\n",
    "copper_trade_cleaned.loc[copper_trade_cleaned[\"stage\"] == \"mined\", \"importer_stage\"] = 'smelting'\n",
    "copper_trade_cleaned.loc[copper_trade_cleaned[\"stage\"] == \"smelted\", \"importer_stage\"] = 'refining'\n",
    "\n",
    "# Renaming columns\n",
    "copper_trade_cleaned.rename(columns={'exporter': 'source_country', \n",
    "                                     'importer': 'target_country', \n",
    "                                     'exporter_stage': 'source_stage',\n",
    "                                     'importer_stage': 'target_stage'}, inplace=True)\n",
    "copper_trade_cleaned = copper_trade_cleaned[['source_stage', 'source_country', 'target_stage', 'target_country', 'tonnes_copper']]\n",
    "\n",
    "# Getting the share countries have of imports/exports of other countries\n",
    "copper_trade_cleaned['share_exports_by_stage'] = copper_trade_cleaned['tonnes_copper'] / copper_trade_cleaned.groupby(['source_stage', 'source_country'])['tonnes_copper'].transform('sum')\n",
    "copper_trade_cleaned['share_imports_by_stage'] = copper_trade_cleaned['tonnes_copper'] / copper_trade_cleaned.groupby(['target_stage', 'target_country'])['tonnes_copper'].transform('sum')\n",
    "\n",
    "# Removing 0 values (nothing to allocate) and differentiating positive (more import necessary) and negative (more export necessary) flows to allocated\n",
    "flows_foreign_pos = flows_foreign[flows_foreign['balance'] > 0]\n",
    "flows_foreign_pos.reset_index(drop=True, inplace=True)\n",
    "\n",
    "flows_foreign_neg = flows_foreign[flows_foreign['balance'] < 0]\n",
    "flows_foreign_neg.reset_index(drop=True, inplace=True)\n",
    "\n",
    "# Merging dataframe to have the necessary allocation amounts in a column\n",
    "copper_trade_allocated = pd.merge(copper_trade_cleaned, flows_foreign_pos,\n",
    "                                  left_on=['target_stage', 'target_country'],\n",
    "                                  right_on=['target_stage', 'country'], \n",
    "                                  how='left', suffixes=('', '_y'))\n",
    "\n",
    "copper_trade_allocated = pd.merge(copper_trade_allocated, flows_foreign_neg,\n",
    "                                  left_on=['source_stage', 'source_country'],\n",
    "                                  right_on=['source_stage', 'country'], \n",
    "                                  how='left', suffixes=('', '_y'))\n",
    "\n",
    "\n",
    "copper_trade_allocated = copper_trade_allocated[['source_stage', 'source_country', 'target_stage', 'target_country', \n",
    "                                                 'tonnes_copper', 'share_exports_by_stage', 'share_imports_by_stage', \n",
    "                                                 'balance', 'balance_y']]\n",
    "copper_trade_allocated.rename(columns={'balance':'imports_to_allocate', 'balance_y':'exports_to_allocate'}, inplace=True)\n",
    "\n",
    "# Filling NaNs (nothing to allocate there) with 0\n",
    "copper_trade_allocated[['imports_to_allocate', 'exports_to_allocate']] = copper_trade_allocated[['imports_to_allocate', 'exports_to_allocate']].fillna(value=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Allocation of the values by multiplying with import and export shares and calculation of flows that need to go to / come from black boxes. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\1729633393.py:19: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  copper_trade_no_black_boxes.rename(columns={'tonnes_copper_allocated':'tonnes_copper'}, inplace=True)\n"
     ]
    }
   ],
   "source": [
    "# Allocating\n",
    "copper_trade_allocated['tonnes_copper_allocated_by_imports'] = copper_trade_allocated['imports_to_allocate'] * copper_trade_allocated['share_imports_by_stage']\n",
    "copper_trade_allocated['tonnes_copper_allocated_by_exports'] = copper_trade_allocated['exports_to_allocate'] * -1 * copper_trade_allocated['share_exports_by_stage']\n",
    "copper_trade_allocated[['tonnes_copper_allocated_by_imports', 'tonnes_copper_allocated_by_exports']] = copper_trade_allocated[['tonnes_copper_allocated_by_imports', 'tonnes_copper_allocated_by_exports']].fillna(value=0)\n",
    "\n",
    "# Calculating what volume of flows that need to go to black box\n",
    "copper_trade_allocated['to_black_box_source'] = copper_trade_allocated['tonnes_copper_allocated_by_exports'] - copper_trade_allocated['tonnes_copper_allocated_by_imports']\n",
    "copper_trade_allocated['to_black_box_source'] = copper_trade_allocated['to_black_box_source'].apply(lambda x: max(0, x))\n",
    "\n",
    "copper_trade_allocated['from_black_box_target'] = copper_trade_allocated['tonnes_copper_allocated_by_imports'] - copper_trade_allocated['tonnes_copper_allocated_by_exports']\n",
    "copper_trade_allocated['from_black_box_target'] = copper_trade_allocated['from_black_box_target'].apply(lambda x: max(0, x))\n",
    "\n",
    "# Calculating the volume of flows that has been successfully allocated with the trade data\n",
    "copper_trade_allocated['tonnes_copper_allocated'] = copper_trade_allocated[['tonnes_copper_allocated_by_exports', 'tonnes_copper_allocated_by_imports']].min(axis=1)\n",
    "\n",
    "# Setting up the dataframes for flows that do not need black boxes\n",
    "copper_trade_no_black_boxes = copper_trade_allocated[['source_stage', 'source_country', 'target_stage', 'target_country', \n",
    "                                                 'tonnes_copper_allocated']]\n",
    "copper_trade_no_black_boxes.rename(columns={'tonnes_copper_allocated':'tonnes_copper'}, inplace=True)\n",
    "\n",
    "copper_trade_no_black_boxes = copper_trade_no_black_boxes[copper_trade_no_black_boxes['tonnes_copper'] != 0]\n",
    "copper_trade_no_black_boxes.reset_index(inplace=True, drop=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Setting up the dataframes for flows that need black boxes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\2272345297.py:13: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  copper_trade_black_boxes_from.rename(columns={'from_black_box_target': 'tonnes_copper'}, inplace=True)\n",
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\2272345297.py:15: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  copper_trade_black_boxes_to.rename(columns={'to_black_box_source':'tonnes_copper'}, inplace=True)\n"
     ]
    }
   ],
   "source": [
    "# Setting up the dataframes for flows that need black boxes\n",
    "copper_trade_black_boxes = copper_trade_allocated[['source_stage', 'source_country', 'target_stage', 'target_country', 'from_black_box_target', 'to_black_box_source']]\n",
    "copper_trade_black_boxes = copper_trade_black_boxes[(copper_trade_black_boxes['from_black_box_target'] != 0) | (copper_trade_black_boxes['to_black_box_source'] != 0)]\n",
    "\n",
    "# Renaming source and target columns to black boxes (4 cases)\n",
    "copper_trade_black_boxes.loc[(copper_trade_black_boxes['to_black_box_source'] > 0) & (copper_trade_black_boxes['source_stage'] == 'mining'), 'target_country'] = 'balance_smelting'\n",
    "copper_trade_black_boxes.loc[(copper_trade_black_boxes['from_black_box_target'] > 0) & (copper_trade_black_boxes['target_stage'] == 'smelting'), 'source_country'] = 'balance_mining'\n",
    "copper_trade_black_boxes.loc[(copper_trade_black_boxes['to_black_box_source'] > 0) & (copper_trade_black_boxes['source_stage'] == 'smelting'), 'target_country'] = 'balance_refining'\n",
    "copper_trade_black_boxes.loc[(copper_trade_black_boxes['from_black_box_target'] > 0) & (copper_trade_black_boxes['target_stage'] == 'refining'), 'source_country'] = 'balance_smelting'\n",
    "\n",
    "# Transforming to dataframe in sankey format\n",
    "copper_trade_black_boxes_from = copper_trade_black_boxes[['source_stage', 'source_country', 'target_stage', 'target_country', 'from_black_box_target']]\n",
    "copper_trade_black_boxes_from.rename(columns={'from_black_box_target': 'tonnes_copper'}, inplace=True)\n",
    "copper_trade_black_boxes_to = copper_trade_black_boxes[['source_stage', 'source_country', 'target_stage', 'target_country', 'to_black_box_source']]\n",
    "copper_trade_black_boxes_to.rename(columns={'to_black_box_source':'tonnes_copper'}, inplace=True)\n",
    "copper_trade_black_boxes = pd.concat([copper_trade_black_boxes_from, copper_trade_black_boxes_to])\n",
    "copper_trade_black_boxes = copper_trade_black_boxes[(copper_trade_black_boxes['tonnes_copper'] != 0)]\n",
    "copper_trade_black_boxes = copper_trade_black_boxes.sort_values(by=['source_stage', 'source_country', 'target_country']).reset_index(drop=True)\n",
    "\n",
    "# Summing up values that are double because of renaming to black boxes ()\n",
    "copper_trade_black_boxes = copper_trade_black_boxes.groupby(['source_stage', 'source_country', 'target_stage', 'target_country'], as_index=False)['tonnes_copper'].sum()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Allocating the remainder of the flows that have not been allocated so far. These are flows for which there are differences between the USGS production stages, but no flows are in the trade data. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Getting all country/stage combinations in the countries considered after filtering the USGS data\n",
    "countries_stages = usgs_size.copy()\n",
    "countries_stages['mining_tonnes'] = 'mining'\n",
    "countries_stages['smelting_tonnes'] = 'smelting'\n",
    "countries_stages['refining_tonnes'] = 'refining'\n",
    "countries_stages = pd.melt(countries_stages, id_vars='country')\n",
    "countries_stages.rename(columns={'value':'stage'}, inplace = True)\n",
    "countries_stages.drop(['variable'], axis=1, inplace=True)\n",
    "countries_stages.reset_index(drop=True, inplace=True)\n",
    "\n",
    "# Getting all country/stage combinations that are in the trade data for the allocation\n",
    "countries_stages_allocated1 = copper_trade_allocated[['source_stage', 'source_country']].drop_duplicates()\n",
    "countries_stages_allocated1.rename(columns={'source_stage':'stage', 'source_country':'country'}, inplace = True)\n",
    "countries_stages_allocated2 = copper_trade_allocated[['target_stage', 'target_country']].drop_duplicates()\n",
    "countries_stages_allocated2.rename(columns={'target_stage':'stage', 'target_country':'country'}, inplace = True)\n",
    "countries_stages_allocated = pd.concat([countries_stages_allocated1, countries_stages_allocated2], ignore_index=True).drop_duplicates()\n",
    "\n",
    "# Getting the difference to know what country / stage combinations were not available for allocation\n",
    "df1_common = pd.merge(countries_stages, countries_stages_allocated, how='inner')\n",
    "countries_stages_not_allocated = pd.concat([countries_stages, df1_common, df1_common]).drop_duplicates(keep=False)\n",
    "countries_stages_not_allocated.reset_index(inplace=True)\n",
    "\n",
    "# Creating dummy combinations for those country/stage combinations missing (the stages must be created precisely for what is missing, otherwise double values are generated)\n",
    "new_rows = []\n",
    "\n",
    "for index in countries_stages_not_allocated.index:\n",
    "    \n",
    "    # Mining stage missing --> flow from mining stage of missing country to smelting black box must be added\n",
    "    if countries_stages_not_allocated.loc[index, 'stage'] == 'mining':\n",
    "        new_rows.extend([['mining', countries_stages_not_allocated.loc[index, 'country'], 'smelting', 'balance_smelting']])\n",
    "    \n",
    "    # Smelting stage missing --> flow from mining stage black box to smelting stage (missing country) and to refining stage (black box) must be added\n",
    "    elif countries_stages_not_allocated.loc[index, 'stage'] == 'smelting':\n",
    "        new_rows.extend([\n",
    "            ['mining', 'balance_mining', 'smelting', countries_stages_not_allocated.loc[index, 'country']],\n",
    "            ['smelting', countries_stages_not_allocated.loc[index, 'country'], 'refining', 'balance_refining'],\n",
    "        ])\n",
    "    \n",
    "    # Refining stage missing --> flow from smelting (black box) to refining stage of missing country must be added\n",
    "    elif countries_stages_not_allocated.loc[index, 'stage'] == 'refining':\n",
    "        new_rows.extend([['smelting', 'balance_smelting', 'refining', countries_stages_not_allocated.loc[index, 'country']]])\n",
    "\n",
    "USGS_values_not_allocated = pd.DataFrame(new_rows, columns=['source_stage', 'source_country', 'target_stage', 'target_country'])\n",
    "\n",
    "# Allocating the remainder of flows calculated from USGS\n",
    "USGS_values_not_allocated = pd.merge(USGS_values_not_allocated, flows_foreign_pos,\n",
    "                                  left_on=['target_stage', 'target_country'],\n",
    "                                  right_on=['target_stage', 'country'], \n",
    "                                  how='left', suffixes=('', '_y'))\n",
    "\n",
    "USGS_values_not_allocated = pd.merge(USGS_values_not_allocated, flows_foreign_neg,\n",
    "                                  left_on=['source_stage', 'source_country'],\n",
    "                                  right_on=['source_stage', 'country'], \n",
    "                                  how='left', suffixes=('', '_y'))\n",
    "\n",
    "# Cleaning the dataframe from those rows without any allocated copper flows\n",
    "USGS_values_not_allocated.drop(['country', 'source_stage_y', 'country_y', 'target_stage_y'], axis=1, inplace=True)\n",
    "USGS_values_not_allocated.rename(columns={'balance':'imports_to_allocate', 'balance_y':'exports_to_allocate'}, inplace=True)\n",
    "USGS_values_not_allocated = USGS_values_not_allocated.dropna(subset=['imports_to_allocate', 'exports_to_allocate'], how='all')\n",
    "USGS_values_not_allocated = USGS_values_not_allocated.drop_duplicates()\n",
    "USGS_values_not_allocated = USGS_values_not_allocated.fillna(0)\n",
    "USGS_values_not_allocated.reset_index(drop=True, inplace=True)\n",
    "\n",
    "# Merging the two columns describing exports to allocated and imports to allocate as no further allocation is required\n",
    "# Due to structure of dummy dataframe these already represent flows\n",
    "USGS_values_not_allocated['exports_to_allocate'] = USGS_values_not_allocated['exports_to_allocate'] * -1\n",
    "USGS_values_not_allocated['tonnes_copper'] = USGS_values_not_allocated[['exports_to_allocate', 'imports_to_allocate']].max(axis=1)\n",
    "USGS_values_not_allocated.drop(['exports_to_allocate', 'imports_to_allocate'], axis=1, inplace=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Building the final dataframe with all flows for the Sankey"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [],
   "source": [
    "copper_flows_sankey = pd.concat([flows_domestic, copper_trade_no_black_boxes, copper_trade_black_boxes, USGS_values_not_allocated])\n",
    "copper_flows_sankey = copper_flows_sankey.sort_values(by=['source_stage', 'source_country', 'target_country']).reset_index(drop=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Scaling the flows to NL"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [],
   "source": [
    "if scale_to_NL == True: \n",
    "    copper_flows_sankey['tonnes_copper'] = copper_flows_sankey['tonnes_copper'] * NL_share\n",
    "\n",
    "    refining_df = copper_flows_sankey[copper_flows_sankey['target_stage'] == 'refining'].drop_duplicates()\n",
    "    refining_df = refining_df.groupby(['target_stage', 'target_country'])['tonnes_copper'].sum().reset_index()\n",
    "    \n",
    "    use_rows = []\n",
    "\n",
    "    for index in refining_df.index:\n",
    "    \n",
    "        use_rows.extend([\n",
    "            ['refining', refining_df.loc[index, 'target_country'], 'use', 'Energy', refining_df.loc[index, 'tonnes_copper'] * NL_energy_share],\n",
    "            ['refining', refining_df.loc[index, 'target_country'], 'use', 'Other', refining_df.loc[index, 'tonnes_copper'] * (1 - NL_energy_share)],\n",
    "        ])\n",
    "    \n",
    "    use_rows = pd.DataFrame(use_rows, columns=['source_stage', 'source_country', 'target_stage', 'target_country', 'tonnes_copper'])\n",
    "\n",
    "    copper_flows_sankey = pd.concat([copper_flows_sankey, use_rows], ignore_index=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "### **6. Importing LCA data and calculation of impact flows**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Data import"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [],
   "source": [
    "lca_data = pd.read_excel(\"data/compilation_ecoinvent.xlsx\", sheet_name= \"Impacts\")\n",
    "lca_data.drop([\"Region\", \"climate change, specific for country\", \"water depletion, specific for country\", \"natural land transformation, specific for country\"], axis=1, inplace=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Mering with the MFA and calculating absolute impacts of the flows (scaling LCA impacts intensities with the MFA flow volumes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Merging with copper flows\n",
    "impact_flows = pd.merge(copper_flows_sankey, lca_data, left_on=['source_stage', 'source_country'], right_on=['Stage', 'Geography'], how='left').fillna(\"RoW\")\n",
    "\n",
    "# Transforming units to impact/tonne not per kg as in the LCA data\n",
    "impact_flows['GHG_emissions'] = impact_flows['GHG_emissions'] * 1000\n",
    "impact_flows['water_depletion'] = impact_flows['water_depletion'] * 1000\n",
    "impact_flows['natural_land_transformation'] = impact_flows['natural_land_transformation'] * 1000\n",
    "\n",
    "# Calculating absolute impacts per flow\n",
    "impact_flows['GHG (tCO₂-eq.)'] = impact_flows['tonnes_copper'] * impact_flows['GHG_emissions'] / 1000 # conversion from kg\n",
    "impact_flows['Water depletion (m3)'] = impact_flows['tonnes_copper'] * impact_flows['water_depletion']\n",
    "impact_flows['Natural land transformation (m2)'] = impact_flows['tonnes_copper'] * impact_flows['natural_land_transformation'] \n",
    "\n",
    "impact_flows = impact_flows.drop(['Stage', 'Geography', 'GHG_emissions', 'water_depletion', 'natural_land_transformation'], axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Grouping countries that are under the relevance threshold as \"Other\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Filtering the dataframe by countries that have a share of at least what is defined under relevance_threshold\n",
    "filter_DF = usgs_data[usgs_data[['mining_share', 'smelting_share', 'refining_share']].max(axis=1) >= relevance_threshold]\n",
    "filter_countries = filter_DF['country'].to_list()\n",
    "filter_countries.append('balance_mining')\n",
    "filter_countries.append('balance_smelting')\n",
    "filter_countries.append('balance_refining')\n",
    "filter_countries.append('Energy')\n",
    "filter_countries.append('Other')\n",
    "\n",
    "# Grouping other countries as \"Other\"\n",
    "impact_flows_final = impact_flows.copy()\n",
    "impact_flows_final['source_country'] = impact_flows_final['source_country'].apply(lambda x: x if x in filter_countries else 'a_Other')\n",
    "impact_flows_final['target_country'] = impact_flows_final['target_country'].apply(lambda x: x if x in filter_countries else 'a_Other')\n",
    "\n",
    "impact_flows_final = impact_flows_final.groupby(['source_stage', 'source_country', 'target_stage', 'target_country']).agg({\n",
    "    'tonnes_copper': 'sum',\n",
    "    'GHG (tCO₂-eq.)': 'sum',\n",
    "    'Natural land transformation (m2)': 'sum',\n",
    "    'Water depletion (m3)': 'sum'\n",
    "}).reset_index()\n",
    "\n",
    "# Sorting by stage and country\n",
    "sort1 = impact_flows_final[(impact_flows_final['source_stage'] == 'mining') | (impact_flows_final['source_stage'] == 'smelting')]\n",
    "sort1 = sort1.sort_values(by=['source_stage', 'source_country'])\n",
    "\n",
    "sort2 = impact_flows_final[(impact_flows_final['source_stage'] == 'refining')]\n",
    "sort2 = sort2.sort_values(by=['source_country'])\n",
    "\n",
    "impact_flows_final = pd.concat([sort1, sort2])\n",
    "impact_flows_final.reset_index(inplace=True, drop=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Setting impact to show in the Sankey (based on user decision under 2.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\aaron\\AppData\\Local\\Temp\\ipykernel_24192\\2754204386.py:8: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  sankey_flows.rename(columns={'GHG (tCO₂-eq.)':'value'}, inplace=True)\n"
     ]
    }
   ],
   "source": [
    "#%% Deciding which impact to show and changing the title\n",
    "if display == 'copper': \n",
    "    sankey_flows = impact_flows_final[['source_stage', 'source_country', 'target_stage', 'target_country', 'tonnes_copper']]\n",
    "    sankey_flows.rename(columns={'tonnes_copper':'value'}, inplace=True)\n",
    "    \n",
    "elif display == 'ghg': \n",
    "    sankey_flows = impact_flows_final[['source_stage', 'source_country', 'target_stage', 'target_country', 'GHG (tCO₂-eq.)']]\n",
    "    sankey_flows.rename(columns={'GHG (tCO₂-eq.)':'value'}, inplace=True)\n",
    "    \n",
    "elif display == 'land_use':\n",
    "    sankey_flows = impact_flows_final[['source_stage', 'source_country', 'target_stage', 'target_country', 'Natural land transformation (m2)']]\n",
    "    sankey_flows.rename(columns={'Natural land transformation (m2)':'value'}, inplace=True)\n",
    "    \n",
    "elif display == 'water':    \n",
    "    sankey_flows = impact_flows_final[['source_stage', 'source_country', 'target_stage', 'target_country', 'Water depletion (m3)']]  \t \n",
    "    sankey_flows.rename(columns={'Water depletion (m3)':'value'}, inplace=True)\n",
    "    \n",
    "else:\n",
    "    print('The choice of what to display in the Sankey is invalid')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "----\n",
    "### **7. Visualisation**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "#### 7.1 **Preparing dataframe for the sankey**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### 7.1.1 Node numbers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Creating node numbers that are unique for the combination of country and stage (e.g. smelting & Australia)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Mining stage\n",
    "mining_codes = sankey_flows[['source_stage', 'source_country']].drop_duplicates()\n",
    "mining_codes = mining_codes[mining_codes['source_stage'] == 'mining']\n",
    "mining_codes.rename(columns={'source_stage': 'stage', 'source_country': 'country'}, inplace = True)\n",
    "mining_codes.reset_index(drop=True, inplace=True)\n",
    "\n",
    "# Smelting stage\n",
    "smelting_codes_source = sankey_flows[['source_stage', 'source_country']].drop_duplicates()\n",
    "smelting_codes_source = smelting_codes_source[smelting_codes_source['source_stage'] == 'smelting']\n",
    "smelting_codes_source.rename(columns={'source_stage': 'stage', 'source_country': 'country'}, inplace = True)\n",
    "\n",
    "smelting_codes_target = sankey_flows[['target_stage', 'target_country']].drop_duplicates()\n",
    "smelting_codes_target = smelting_codes_target[smelting_codes_target['target_stage'] == 'smelting']\n",
    "smelting_codes_target.rename(columns={'target_stage': 'stage', 'target_country': 'country'}, inplace = True)\n",
    "\n",
    "smelting_codes = pd.concat([smelting_codes_source, smelting_codes_target]).drop_duplicates()\n",
    "smelting_codes = smelting_codes.sort_values(by='country', ascending= True)\n",
    "smelting_codes.reset_index(drop=True, inplace=True)\n",
    "\n",
    "# Refining stage\n",
    "refining_codes_source = sankey_flows[[\"source_stage\", \"source_country\"]].drop_duplicates()\n",
    "refining_codes_source = refining_codes_source[refining_codes_source[\"source_stage\"] == 'refining']\n",
    "refining_codes_source.rename(\n",
    "    columns={\"source_stage\": \"stage\", \"source_country\": \"country\"}, inplace=True\n",
    ")\n",
    "\n",
    "refining_codes_target = sankey_flows[[\"target_stage\", \"target_country\"]].drop_duplicates()\n",
    "refining_codes_target = refining_codes_target[refining_codes_target[\"target_stage\"] == 'refining']\n",
    "refining_codes_target.rename(\n",
    "    columns={\"target_stage\": \"stage\", \"target_country\": \"country\"}, inplace=True\n",
    ")\n",
    "\n",
    "refining_codes = pd.concat([refining_codes_source, refining_codes_target]).drop_duplicates()\n",
    "refining_codes = refining_codes.sort_values(by=\"country\", ascending=True)\n",
    "refining_codes.reset_index(drop=True, inplace=True)\n",
    "\n",
    "if scale_to_NL == True:\n",
    "    # Use stage\n",
    "    use_codes = sankey_flows[['target_stage', 'target_country']].drop_duplicates()\n",
    "    use_codes = use_codes[use_codes['target_stage'] == 'use']\n",
    "    use_codes.rename(columns={'target_stage': 'stage', 'target_country': 'country'}, inplace = True)\n",
    "    use_codes = use_codes.sort_values(by='country', ascending= True)\n",
    "    use_codes.reset_index(drop=True, inplace=True)\n",
    "\n",
    "# Building complete dataframe and resetting index to get a column with a unique number per node, sorted by stage and then country\n",
    "if scale_to_NL == True:\n",
    "    sankey_nodes = pd.concat([mining_codes, smelting_codes, refining_codes, use_codes])\n",
    "else:\n",
    "    sankey_nodes = pd.concat([mining_codes, smelting_codes, refining_codes])\n",
    "sankey_nodes = sankey_nodes.reset_index(drop=True)\n",
    "sankey_nodes = sankey_nodes.reset_index(drop=False)\n",
    "sankey_nodes.rename(columns={'index': 'node_number'}, inplace = True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Adding the node numbers to the flows"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Adding exporter nodes\n",
    "sankey_flows_final = pd.merge(sankey_flows, sankey_nodes, \n",
    "                                     left_on=['source_stage', 'source_country'], \n",
    "                                     right_on=['stage', 'country'], \n",
    "                                     how='left')\n",
    "\n",
    "sankey_flows_final.drop(['stage', 'country'], axis = 1, inplace = True)\n",
    "\n",
    "# Adding importer nodes\n",
    "sankey_flows_final = pd.merge(sankey_flows_final, sankey_nodes, \n",
    "                                     left_on=['target_stage', 'target_country'], \n",
    "                                     right_on=['stage', 'country'], \n",
    "                                     how='left')\n",
    "\n",
    "sankey_flows_final.drop(['stage', 'country'], axis = 1, inplace = True)\n",
    "\n",
    "# Renaming columns\n",
    "sankey_flows_final.rename(columns={'node_number_x': 'source', \n",
    "                                'node_number_y': 'target'}, \n",
    "                                 inplace = True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### 7.2.2 Node positions\n",
    "\n",
    "This step is necessary because plotly does not allow to sort the nodes, e.g. by countries. Instead, the position of each node has to be defined using x and y values in the figure space. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculating the total size of each stage as a reference for the diagram (max stage size <-> 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {},
   "outputs": [],
   "source": [
    "export_sums = sankey_flows_final.groupby(['source_stage', 'source_country'])['value'].sum().reset_index()\n",
    "export_sums.rename(columns={'source_stage':'stage','source_country':'country'}, inplace=True)\n",
    "import_sums = sankey_flows_final.groupby(['target_stage', 'target_country'])['value'].sum().reset_index()\n",
    "import_sums.rename(columns={'target_stage':'stage','target_country':'country'}, inplace=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculatingn node sizes from export/import values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {},
   "outputs": [],
   "source": [
    "# mining\n",
    "sizes_mining = pd.merge(sankey_nodes[sankey_nodes['stage'] == 'mining'], export_sums[export_sums['stage'] == 'mining'], \n",
    "                     on='country', how='outer', suffixes=('', '_y'))\n",
    "sizes_mining = sizes_mining.drop(['node_number', 'stage_y'], axis = 1,)\n",
    "\n",
    "# smelting\n",
    "sizes_smelting_exports = pd.merge(sankey_nodes[sankey_nodes['stage'] == 'smelting'], export_sums[export_sums['stage'] == 'smelting'], \n",
    "                     on='country', how='outer', suffixes=('', '_y'))\n",
    "sizes_smelting_exports = sizes_smelting_exports.drop(['node_number', 'stage_y'], axis = 1,)\n",
    "\n",
    "sizes_smelting_imports = pd.merge(sankey_nodes[sankey_nodes['stage'] == 'smelting'], import_sums[import_sums['stage'] == 'smelting'], \n",
    "                     on='country', how='outer', suffixes=('', '_y'))\n",
    "sizes_smelting_imports = sizes_smelting_imports.drop(['node_number', 'stage_y'], axis = 1,)\n",
    "\n",
    "# when calculating the impacts, the inflows and outflows are not always the same --> need to take the maximum of both to calculate the node size for the positions\n",
    "sizes_smelting = pd.merge(sizes_smelting_exports, sizes_smelting_imports, on=['stage', 'country']).assign(value=lambda x: x[['value_x', 'value_y']].max(axis=1))\n",
    "sizes_smelting = sizes_smelting[['stage', 'country', 'value']]\n",
    "\n",
    "if scale_to_NL == True:\n",
    "\n",
    "    sizes_refining_exports = pd.merge(sankey_nodes[sankey_nodes['stage'] == 'refining'], export_sums[export_sums['stage'] == 'refining'], \n",
    "                          on='country', how='outer', suffixes=('', '_y'))\n",
    "    sizes_refining_exports = sizes_refining_exports.drop(['node_number', 'stage_y'], axis = 1,)\n",
    "\n",
    "    sizes_refining_imports = pd.merge(sankey_nodes[sankey_nodes['stage'] == 'refining'], import_sums[import_sums['stage'] == 'refining'], \n",
    "                          on='country', how='outer', suffixes=('', '_y'))\n",
    "    sizes_refining_imports = sizes_refining_imports.drop(['node_number', 'stage_y'], axis = 1,)\n",
    "\n",
    "    # when calculating the impacts, the inflows and outflows are not always the same --> need to take the maximum of both to calculate the node size for the positions\n",
    "    sizes_refining = pd.merge(sizes_refining_exports, sizes_refining_imports, on=['stage', 'country']).assign(value=lambda x: x[['value_x', 'value_y']].max(axis=1))\n",
    "    sizes_refining = sizes_refining[['stage', 'country', 'value']]\n",
    "    \n",
    "    sizes_use = pd.merge(sankey_nodes[sankey_nodes['stage'] == 'use'], import_sums[import_sums['stage'] == 'use'], \n",
    "                          on='country', how='outer', suffixes=('', '_y'))\n",
    "    sizes_use = sizes_use.drop(['node_number', 'stage_y'], axis = 1,)\n",
    "else:    \n",
    "    # refining only has exports if allocating to NL/use\n",
    "    sizes_refining = pd.merge(sankey_nodes[sankey_nodes['stage'] == 'refining'], import_sums[import_sums['stage'] == 'refining'], \n",
    "                          on='country', how='outer', suffixes=('', '_y'))\n",
    "    sizes_refining = sizes_refining.drop(['node_number', 'stage_y'], axis = 1,)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculating the node positions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [],
   "source": [
    "# calculating scale factor to normalise node size to scale 0 to 1\n",
    "if scale_to_NL == True:\n",
    "    total_scale = 1 / max(sizes_mining['value'].sum(), sizes_smelting['value'].sum(), sizes_refining['value'].sum(), sizes_use['value'].sum())\n",
    "else: \n",
    "    total_scale = 1 / max(sizes_mining['value'].sum(), sizes_smelting['value'].sum(), sizes_refining['value'].sum())\n",
    "\n",
    "# Mining nodes\n",
    "sizes_mining['x_pos'] = 0.005\n",
    "sizes_mining['y_size'] = sizes_mining['value'] * total_scale\n",
    "sizes_mining.loc[0, 'y_pos'] = sizes_mining.loc[0, 'y_size']/2\n",
    "for i in range(1, len(sizes_mining)):\n",
    "    sizes_mining.loc[i, 'y_pos'] = sizes_mining.loc[i, 'y_size']/2 + sizes_mining.loc[i-1, 'y_size']/2 + sizes_mining.loc[i-1, 'y_pos']\n",
    "\n",
    "# Smelting nodes\n",
    "if scale_to_NL == True:\n",
    "    sizes_smelting['x_pos'] = 0.335\n",
    "else:\n",
    "    sizes_smelting['x_pos'] = 0.5\n",
    "sizes_smelting['y_size'] = sizes_smelting['value'] * total_scale\n",
    "sizes_smelting.loc[0, 'y_pos'] = sizes_smelting.loc[0, 'y_size']/2\n",
    "\n",
    "for i in range(1, len(sizes_smelting)):\n",
    "    sizes_smelting.loc[i, 'y_pos'] = sizes_smelting.loc[i, 'y_size']/2 + sizes_smelting.loc[i-1, 'y_size']/2 + sizes_smelting.loc[i-1, 'y_pos']\n",
    "\n",
    "# Refining nodes\n",
    "if scale_to_NL == True:\n",
    "    sizes_refining['x_pos'] = 0.665\n",
    "else:\n",
    "    sizes_refining['x_pos'] = 0.995\n",
    "sizes_refining['y_size'] = sizes_refining['value'] * total_scale\n",
    "sizes_refining.loc[0, 'y_pos'] = sizes_refining.loc[0, 'y_size']/2\n",
    "\n",
    "for i in range(1, len(sizes_refining)):\n",
    "    sizes_refining.loc[i, 'y_pos'] = sizes_refining.loc[i, 'y_size']/2 + sizes_refining.loc[i-1, 'y_size']/2 + sizes_refining.loc[i-1, 'y_pos']\n",
    "\n",
    "# sizes use\n",
    "if scale_to_NL == True:\n",
    "    sizes_use['x_pos'] = 0.995\n",
    "    sizes_use['y_size'] = sizes_use['value'] * total_scale\n",
    "    sizes_use.loc[0, 'y_pos'] = sizes_use.loc[0, 'y_size']/2\n",
    "\n",
    "    for i in range(1, len(sizes_use)):\n",
    "        sizes_use.loc[i, 'y_pos'] = sizes_use.loc[i, 'y_size']/2 + sizes_use.loc[i-1, 'y_size']/2 + sizes_use.loc[i-1, 'y_pos']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Assembling all sizes and positions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {},
   "outputs": [],
   "source": [
    "# assembling all sizes and positions\n",
    "if scale_to_NL == True:\n",
    "    sizes = pd.concat([sizes_mining, sizes_smelting, sizes_refining, sizes_use]).reset_index(drop = True)\n",
    "else: \n",
    "    sizes = pd.concat([sizes_mining, sizes_smelting, sizes_refining]).reset_index(drop = True)\n",
    "sizes.reset_index(drop =True, inplace=True)\n",
    "\n",
    "# Getting the share of the country for the specific stage\n",
    "sizes['share_by_stage'] = sizes['value'] / sizes.groupby('stage')['value'].transform('sum')\n",
    "\n",
    "# extracing positional values for sankey\n",
    "y_position = sizes['y_pos'].tolist()\n",
    "x_position = sizes['x_pos'].tolist()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### 7.2.3 Node and flow colours"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Node colours (imported from external excel)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Importing node colours\n",
    "if highlight_LMIC == True: \n",
    "    country_colours = pd.read_excel(\"data/colours.xlsx\", sheet_name='grouped') # grouped by industrial, LMIC and black boxes\n",
    "else:\n",
    "    country_colours = pd.read_excel(\"data/colours.xlsx\", sheet_name='unique') # unique for each country\n",
    "\n",
    "sizes = pd.merge(left=sizes, right=country_colours, on='country', how='left')\n",
    "\n",
    "# Extracting the list of node colours for the Sankey\n",
    "colours_nodes = sizes['colour'].tolist()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Flow colours (imported from external excel)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Assigning flow colours so that the colour is the same as for the exporter node\n",
    "flow_colours = pd.merge(sankey_flows_final[['source_country']], country_colours, \n",
    "                        left_on='source_country', right_on='country', how='left')\n",
    "\n",
    "# Extracting the list of flow colours for the Sankey\n",
    "colours_links = flow_colours['colour'].tolist() \n",
    "\n",
    "# Changing opacity of flow colours\n",
    "opacity_links = '0.5)'\n",
    "for i in range(len(colours_links)):\n",
    "    colours_links[i] = colours_links[i].replace(\"1)\", opacity_links)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### 7.2.4 Assembling label for the Sankey and converting for format for plotly"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Country names\n",
    "names = sankey_nodes['country'].tolist()\n",
    "\n",
    "# Changing the names for the black boxes and \"Other\" category for the visualization\n",
    "for i in range(len(names)):\n",
    "    if names[i] == 'balance_mining':\n",
    "        names[i] = 'Black Box'\n",
    "    elif names[i] == 'balance_smelting':\n",
    "        names[i] = 'Black Box'\n",
    "    elif names[i] == 'balance_refining':\n",
    "        names[i] = 'Black Box'\n",
    "    elif names[i] == 'a_Other':\n",
    "        names[i] = 'Other'\n",
    "    elif names[i] == 'Energy':\n",
    "        names[i] = 'Renewables' \n",
    "        \n",
    "if sankey_labels_absolute == True: \n",
    "    numbers = sizes['value'].tolist()\n",
    "    formatted_numbers = [\"{:,.0f}\".format(num) for num in numbers]\n",
    "    labels = [str(x) + ' ' + str(y) for x, y in zip(names, formatted_numbers)]\n",
    "else: \n",
    "    numbers = sizes['share_by_stage'].tolist()\n",
    "    formatted_percentages = [\"{:.0%}\".format(num) for num in numbers]\n",
    "    labels = [str(x) + ' ' + str(y) for x, y in zip(names, formatted_percentages)]    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### 7.2.5 Assembling label for the Sankey and converting for format for plotly"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Extract source, target, and value columns\n",
    "source = sankey_flows_final['source'].tolist()\n",
    "target = sankey_flows_final['target'].tolist()\n",
    "value = sankey_flows_final['value'].tolist()\n",
    "\n",
    "# Define link and node dictionaries for the Sankey diagram\n",
    "link = dict(\n",
    "    source=source, \n",
    "    target=target, \n",
    "    value=value, \n",
    "    color=colours_links\n",
    "    )\n",
    "\n",
    "node = {\"label\": labels, \n",
    "        'pad': 0, \n",
    "        'thickness': 30,\n",
    "        \"x\": x_position,  \n",
    "        \"y\": y_position, \n",
    "        \"color\": colours_nodes, \n",
    "        }"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### **7.2 Plotting**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Figure size\n",
    "width = (24.62 / 2.54) * 96\n",
    "height = (13.95 / 2.54) * 96\n",
    "\n",
    "# Creating the figure\n",
    "fig = go.Figure(go.Sankey(\n",
    "    arrangement='perpendicular',\n",
    "    link=link, \n",
    "    node=node\n",
    "    )) \n",
    "\n",
    "# Updating the layout\n",
    "fig.update_layout(\n",
    "    hovermode='x',\n",
    "    font_family=\"Georgia\",\n",
    "    font=dict(size=12), \n",
    "    margin=dict(l=15, r=15, b=5, t=20),\n",
    "    height = height,    \n",
    "    width = width\n",
    ")\n",
    "\n",
    "# Adding headers for mining, smelting and refining (and use)\n",
    "if scale_to_NL == True: \n",
    "    annotations = [\"Mining\", \"Smelting\", \"Refining\", \"Use\"]\n",
    "    annotation_positions = [-0.01, 0.305, 0.6925, 1.005]\n",
    "else: \n",
    "    annotations = [\"Mining\", \"Smelting\", \"Refining\"]\n",
    "    annotation_positions = [-0.01, 0.5, 1.01]\n",
    "\n",
    "# Adding annotations\n",
    "for i in range(len(annotations)): \n",
    "    fig.add_annotation(x=annotation_positions[i], \n",
    "                   y=1.04, \n",
    "                   text=annotations[i],\n",
    "                   showarrow=False,\n",
    "                   font=dict(size=12, family='Georgia'),\n",
    "                   )\n",
    "\n",
    "# Showing and saving\n",
    "fig.show(renderer=\"browser\")\n",
    "fig.write_image(\"sankey.svg\", engine=\"kaleido\")\n",
    "fig.write_html(\"sankey.html\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "thesis",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}