{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "KcHw7S3T6cj_"
   },
   "source": [
    "#Workshop 1:\n",
    "\n",
    "1-1: Attributing emissions under alternative allocation schemes\n",
    "\n",
    "1-2: Avoiding double counting\n",
    "\n",
    "Data: JIOT_2015_390.XLSX"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "2TRhxuwQ6ocB"
   },
   "source": [
    "Preparation for libraries."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "id": "4bJc5Guz6nLI"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import numpy.linalg as la\n",
    "import sys"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "vZmisvEa7Fn9"
   },
   "source": [
    "Data input & labels"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "id": "_eYOQB_i7JlO"
   },
   "outputs": [],
   "source": [
    "IOdata = pd.read_excel('JIOT_2015_390.xlsx','IOdata', header=0,index_col=0)\n",
    "\n",
    "l_col = tuple(IOdata.columns)\n",
    "l_row = tuple(IOdata.index)\n",
    "n_sec = 390 # The number of endogenous sectors\n",
    "\n",
    "l_x = l_col[:n_sec] # List of sectors\n",
    "\n",
    "n_v = l_row.index(\"Gross value added\") # Value added\n",
    "n_x = l_row.index(\"Domestic production\") # Output\n",
    "n_g = l_row.index(\"GHG\") # GHG\n",
    "n_ex = l_col.index(\"Exports\")\n",
    "n_im = l_col.index(\"Imports\")\n",
    "\n",
    "# Convert Pd.DataFrame into NUmpy array\n",
    "Data = IOdata.to_numpy()\n",
    "\n",
    "# Output, ghg, value added\n",
    "x_ = Data[n_x, :n_sec] + sys.float_info.epsilon\n",
    "ghg = Data[ n_g, :n_sec ] \n",
    "v_  = Data[ n_v,:n_sec ] + sys.float_info.epsilon\n",
    "\n",
    "# Final demand\n",
    "n_y = 6 # The number of domestic final demand items, excl. Exports.\n",
    "y_ = Data[ :n_sec, n_sec:n_sec + n_y ]\n",
    "ex_ = Data[ :n_sec, n_ex] # Exports\n",
    "im_ = Data[ :n_sec, n_im] # Imports"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "CkThppCe7uR_"
   },
   "source": [
    "Calculations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "id": "l-Kyx_557zxh"
   },
   "outputs": [],
   "source": [
    "# The input coefficients matrix\n",
    "epsilon = sys.float_info.epsilon\n",
    "A_ = Data[:n_sec, :n_sec] / x_\n",
    "\n",
    "# Leontief inverse\n",
    "L = np.linalg.inv( np.eye(A_.shape[0]) - A_ )\n",
    "\n",
    "# Emission and value added coefficients\n",
    "f_ = ghg / x_  # The emission coefficients\n",
    "pi_ = v_ / x_  # Value added ratios\n",
    "\n",
    "#%% Domestic Leontief Inverse (L_d)\n",
    "# Rate of domestic products in the total supply\n",
    "demand_d = np.sum(Data[:n_sec, :n_sec+n_y], axis=1) # Exports not included in denominator\n",
    "d = 1 + im_/ demand_d\n",
    "\n",
    "# Input coefficients matrix adjusted for imports\n",
    "A_d = np.diag(d) @ A_\n",
    "\n",
    "# Leontief inverse\n",
    "L_d = np.linalg.inv( np.eye(A_.shape[0]) - A_d )\n",
    "\n",
    "# Final demand for domestic products including exports\n",
    "y_d = np.diag(d) @ np.sum(y_,axis=1) + ex_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "AKJVLbPg7_LI"
   },
   "source": [
    "Consistency check"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "id": "VBLe2L_K8Blg"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Inconsistent indices for L_d @ y_d: []\n"
     ]
    }
   ],
   "source": [
    "epsilon = 0.00001\n",
    "indices = np.where(np.abs(x_ - L_d @ y_d) > epsilon)[0]\n",
    "print(\"Inconsistent indices for L_d @ y_d:\", indices)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "nsrqnQ878dP0"
   },
   "source": [
    "#Workshop 1-1 Attributing emissions under alternative allocation schemes\n",
    "\n",
    "Emissions calculations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "id": "n9uxD4wSAzDl"
   },
   "outputs": [],
   "source": [
    "# * Consumption-based emissions \"\"\"\n",
    "fLy_d = f_ @ (L_d @ np.diag(y_d))\n",
    "\n",
    "# * Income based emissions\n",
    "# Ghosh matrix\n",
    "B_ = np.diag(1/x_) @ Data[ :n_sec, :n_sec ]\n",
    "B_sum = B_.sum(axis=1)\n",
    "# Ghosh inverse\n",
    "G_ = np.linalg.inv( np.eye(B_.shape[0]) - B_ )\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "iV0nfFV-BAXX"
   },
   "source": [
    "Consistency check for Ghosh inverse\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "id": "im0Q1GVPBDRf"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Inconsistent indices for v_ @ G_: []\n"
     ]
    }
   ],
   "source": [
    "indices = np.where(np.abs(x_ - v_@G_) > epsilon)[0]\n",
    "print(\"Inconsistent indices for v_ @ G_:\", indices)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "id": "ZlOExy1PBax0"
   },
   "outputs": [],
   "source": [
    "# Income based emissions vGf\n",
    "vGf = f_ @ (G_.T @ np.diag(v_))\n",
    "\n",
    "# * Value-added ratio-based allocation of fLy_d \"\"\"\n",
    "va_based = ( np.diag(pi_) @ L ) @ fLy_d\n",
    "\n",
    "# Outputs \"\"\"\n",
    "Resp = pd.DataFrame( np.column_stack((ghg[:n_sec], fLy_d, vGf, va_based)) ,index = l_x, columns = [\"fx\", \"fly\", \"vGf\", \"vLfLy\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "tS-wL-HSBgm1"
   },
   "source": [
    "Export to an excel file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "id": "7HUnSYPBBkV9"
   },
   "outputs": [],
   "source": [
    "Resp.to_excel('Resp.xlsx', index=l_x, engine='openpyxl')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "Yjnqy20CDX6p"
   },
   "source": [
    "End of WS 1-1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "X-pTBfK0CVEe"
   },
   "source": [
    "# 1-2 Avoiding Double Counting based on the Method of Sebastian Dente\n",
    "Avoiding Double Counting based on the Method of Sebastian Dente"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "id": "DCm0F7-_ChQ4"
   },
   "outputs": [],
   "source": [
    "# Target and Other sectors\n",
    "\n",
    "# Read target sectors from Excel and filter\n",
    "l_t = pd.read_excel('JIOT_2015_390.xlsx', 'Target sectors', index_col=0).query('target == 1').iloc[:, 1].tolist()\n",
    "\n",
    "# List difference to determine other sectors\n",
    "l_o = list(set(l_x) - set(l_t))\n",
    "\n",
    "# Indices for target and other sectors\n",
    "num_t = [l_x.index(i) for i in l_t]\n",
    "num_o = [l_x.index(i) for i in l_o]\n",
    "\n",
    "# Output without double counting for target sectors\n",
    "x_t = x_[num_t]\n",
    "x_t_wdc_d  = la.inv( L_d[np.ix_(num_t,num_t)] ) @ x_t\n",
    "\n",
    "# Output with double counting\n",
    "x_t_dc_d = x_t - x_t_wdc_d\n",
    "\n",
    "# Share of x_t double counted\n",
    "k_d = x_t_dc_d/x_t"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "vSsgTb2lC1qk"
   },
   "source": [
    "Emission calculations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "id": "EjWwjfy1C7Lk"
   },
   "outputs": [],
   "source": [
    "#%% Upstream emissions for target sectors\n",
    "E_t_up1 = f_[num_t] @ np.diag(x_t_wdc_d)\n",
    "E_t_up2 = f_[num_t] @ (L_d[np.ix_(num_t, num_t)] - np.eye(len(l_t)) ) @ np.diag(x_t_wdc_d)\n",
    "E_t_up3 = f_[num_o] @ L_d[np.ix_(num_o, num_t)] @ np.diag(x_t_wdc_d)\n",
    "E_t_up = E_t_up1 + E_t_up2 + E_t_up3\n",
    "\n",
    "\n",
    "#%% Down stream emissions\n",
    "\n",
    "# Amount of emissions related to the production of non-target sectors\n",
    "x_o_up = L_d[np.ix_(num_o,num_t)] @x_t_wdc_d\n",
    "x_o_dn = x_[num_o] - x_o_up\n",
    "E_o_dn = f_[num_o] * x_o_dn\n",
    "\n",
    "# Allocating a share of E_o_dn to target materials\n",
    "tmp = A_d[np.ix_(num_t,num_o)]@ la.inv( np.eye(len(l_o)) - A_d[np.ix_(num_o,num_o)] )\n",
    "E_t_dn = tmp@E_o_dn\n",
    "\n",
    "# Total emissions of tareget materials\n",
    "E_t = E_t_up + E_t_dn\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "v6DPRpnkDH6z"
   },
   "source": [
    "Outputs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "id": "89Mq_uJaDJkH"
   },
   "outputs": [],
   "source": [
    "indexes = [l_x.index(element) for element in l_t]\n",
    "Result = pd.DataFrame( np.column_stack((indexes, k_d,E_t_up1,E_t_up2,E_t_up3,E_t_up,E_t_dn,E_t )) ,index = l_t)\n",
    "Result.columns = [\"sector no.\", \"k_d\",\"E_t_up1\", \"E_t_up2\", \"E_t_up3\", \"E_t_up\", \"E_t_dn\", \"E_t\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "RuNQKxyXDKWH"
   },
   "source": [
    "Export to an excel file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "id": "sLwws_KZDQPQ"
   },
   "outputs": [],
   "source": [
    "Result.to_excel('Result.xlsx', index=l_t, engine='openpyxl')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "jq6IXYK3DTyA"
   },
   "source": [
    "End of WS 1-2"
   ]
  }
 ],
 "metadata": {
  "colab": {
   "provenance": []
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.11.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}