Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Updates to Ghosh implementation #146

Open
wants to merge 6 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
23 changes: 14 additions & 9 deletions doc/source/math.rst
Original file line number Diff line number Diff line change
Expand Up @@ -181,13 +181,15 @@ Upstream and Downstream scope 3
------------

In the context of impact analyses the factors of production are often categorized into scope 1, 2 and 3, with scope 3

sub-divided into upstream and downstream.
For a MRIO the scope 1 is the direct impact of the industries. The factors of production scope 1 associated
with some product or service in sector 'i' of monetary value 'r' is given by :math:`S e_i r',
where :math:`e_i' is the 'i^{th}' unit vector.
Scope 2 is the indirect impact through directly consumed energy (mostly electricity). The precise defintion of scope 2 in an MRIO depends
on the list of MRIO sectors that are classified as scope 2 energy suppliers. Scope 2 is therefore included in the
upstream scope 3, which we refer to as upstream indirect impacts. The upstream multipliers are
Scope 2 is the indirect impact through directly consumed energy (mostly electricity). The precise definition of scope 2
in an MRIO depends on the list of MRIO sectors that are classified as scope 2 energy suppliers. Scope 2 is therefore
included in the upstream scope 3, which we refer to as upstream indirect impacts. The upstream multipliers are

.. math::

\begin{equation}
Expand All @@ -201,25 +203,28 @@ The downstream attribution according to Ghosh is done by the input share matrix
.. math::

\begin{equation}
A^{*} = Z^{T} \hat{x}^{-1}
B = \hat{x}^{-1} Z
\end{equation}

Note that we have defined this matrix in analogy with :math:`A`, meaning that the factors of production coefficient
are applied from the right-hand side. The full downstream multiplier (including scope 1) is given
by :math:`S G` where
are given by :math:`S G^{T}` where

.. math::

\begin{equation}
G = (\mathrm{I}- Z^{T}\hat{x}^{-1})^{-1}
G = (\mathrm{I}- B)^{-1} = \hat{x} L \hat{x}^{-1}
\end{equation}
is the transpose of the Ghosh inverse matrix.

is the Ghosh inverse matrix.
The pure downstream multiplier (excluding scope 1) is given by


.. math::

\begin{equation}
M_{down} = S((\mathrm{I}- Z^{T}\hat{x}^{-1})^{-1} -I) = S(G-I) = S(\hat{x}^{-1} L^{T} \hat{x} - I)

M_{down} = S(G^{T}-I)

\end{equation}

The sector's total impact multiplier is simply the sum of :math:`M_{up}`, :math:`S` and :math:`M_{down}`.
Expand Down
12 changes: 6 additions & 6 deletions pymrio/core/mriosystem.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,7 @@
from pymrio.tools.iomath import (
calc_A,
calc_accounts,
calc_As,
calc_B,
calc_F,
calc_F_Y,
calc_G,
Expand Down Expand Up @@ -1720,7 +1720,7 @@ def __init__(
Z=None,
Y=None,
A=None,
As=None,
B=None,
x=None,
L=None,
G=None,
Expand All @@ -1740,7 +1740,7 @@ def __init__(
self.Y = Y
self.x = x
self.A = A
self.As = As
self.B = B
self.L = L
self.G = G
self.unit = unit
Expand Down Expand Up @@ -1867,16 +1867,16 @@ def calc_system(self):
self.A = calc_A(self.Z, self.x)
self.meta._add_modify("Coefficient matrix A calculated")

if self.As is None:
self.As = calc_As(self.Z, self.x)
if self.B is None:
self.B = calc_B(self.Z, self.x)
self.meta._add_modify("Coefficient matrix As calculated")

if self.L is None:
self.L = calc_L(self.A)
self.meta._add_modify("Leontief matrix L calculated")

if self.G is None:
self.G = calc_G(self.As)
self.G = calc_G(self.B)
self.meta._add_modify("Ghosh matrix G calculated")

return self
Expand Down
54 changes: 26 additions & 28 deletions pymrio/tools/iomath.py
Original file line number Diff line number Diff line change
Expand Up @@ -152,10 +152,10 @@ def calc_A(Z, x):
return Z * recix


def calc_As(Z, x):
"""Calculate the As matrix (coefficients) from Z and x
def calc_B(Z, x):
"""Calculate the B matrix (coefficients) from Z and x

As is a normalized version of the industrial flows of the transpose of Z, which quantifies the input to downstream
B is a normalized version of the industrial flows of the transpose of Z, which quantifies the input to downstream
sectors.

Parameters
Expand Down Expand Up @@ -187,13 +187,15 @@ def calc_As(Z, x):
recix = recix.reshape((1, -1))
# use numpy broadcasting - factor ten faster
# Mathematical form - slow
# return Z.dot(np.diagflat(recix))
# return np.diagflat(recix).dot(Z)
if type(Z) is pd.DataFrame:
return pd.DataFrame(
np.transpose(Z.values) * recix, index=Z.index, columns=Z.columns
np.transpose(np.transpose(Z.values) * recix),
index=Z.index,
columns=Z.columns,
)
else:
return np.transpose(Z) * recix
return np.transpose(np.transpose(Z) * recix)


def calc_L(A):
Expand Down Expand Up @@ -232,30 +234,26 @@ def calc_L(A):
return np.linalg.inv(I - A)


def calc_G(As, L=None, x=None):
"""Calculate the Ghosh inverse matrix G either from As (high computation effor) or from Leontief matrix L and x
def calc_G(B, L=None, x=None):
"""Calculate the Ghosh inverse matrix G either from B (high computation effort) or from Leontief matrix L and x
(low computation effort)

G = inverse matrix of (I - As) = hat(x)^{-1} * L^T * hat(x)
G = inverse matrix of (I - B) = hat(x) * L * hat(x)^{-1}

where I is an identity matrix of same shape as As, and hat(x) is the diagonal matrix with values of x on the
where I is an identity matrix of same shape as B, and hat(x) is the diagonal matrix with values of x on the
diagonal.

Note that we define G as the transpose of the Ghosh inverse matrix, so that we can apply the factors of
production intensities from the left-hand-side for both Leontief and Ghosh attribution. In this way the
multipliers have the same (vector) dimensions and can be added.

Parameters
----------
As : pandas.DataFrame or numpy.array
B : pandas.DataFrame or numpy.array
Symmetric input output table (coefficients)

Returns
-------
pandas.DataFrame or numpy.array
Ghosh input output table G
The type is determined by the type of As.
If DataFrame index/columns as As
The type is determined by the type of B.
If DataFrame index/columns as B

"""
# if L has already been calculated, then G can be derived from it with low computational cost.
Expand All @@ -275,20 +273,20 @@ def calc_G(As, L=None, x=None):

if type(L) is pd.DataFrame:
return pd.DataFrame(
recix * np.transpose(L.values * x), index=Z.index, columns=Z.columns
np.transpose(recix * np.transpose(L.values * x)),
index=Z.index,
columns=Z.columns,
)
else:
# G = hat(x)^{-1} * L^T * hat(x) in mathematical form np.linalg.inv(hatx).dot(L.transpose()).dot(hatx).
# G = hat(x) * L * hat(x)^{-1} in mathematical form hatx.dot(L.transpose()).dot(np.linalg.inv(hatx)).
# it is computationally much faster to multiply element-wise because hatx is a diagonal matrix.
return recix * np.transpose(L * x)
return np.transpose(recix * np.transpose(L * x))
else: # calculation of the inverse of I-As has a high computational cost.
I = np.eye(As.shape[0]) # noqa
if type(As) is pd.DataFrame:
return pd.DataFrame(
np.linalg.inv(I - As), index=As.index, columns=As.columns
)
I = np.eye(B.shape[0]) # noqa
if type(B) is pd.DataFrame:
return pd.DataFrame(np.linalg.inv(I - B), index=B.index, columns=B.columns)
else:
return np.linalg.inv(I - As) # G = inverse matrix of (I - As)
return np.linalg.inv(I - B) # G = inverse matrix of (I - B)


def calc_S(F, x):
Expand Down Expand Up @@ -421,7 +419,7 @@ def calc_M(S, L):
def calc_M_down(S, G):
"""Calculate downstream multipliers of the extensions

M_down = S * ( G - I )
M_down = S * ( G^T - I )

Where I is an identity matrix of same shape as G

Expand All @@ -440,7 +438,7 @@ def calc_M_down(S, G):
If DataFrame index/columns as S

"""
return S.dot(G - np.eye(G.shape[0]))
return S.dot(np.transpose(G) - np.eye(G.shape[0]))


def calc_e(M, Y):
Expand Down
2 changes: 1 addition & 1 deletion pymrio/tools/ioparser.py
Original file line number Diff line number Diff line change
Expand Up @@ -941,7 +941,7 @@ def parse_wiod(path, year=None, names=("isic", "c_codes"), popvector=None):
)

# remove meta data, empty rows, total column
wiot_data.iloc[0 : wiot_meta["end_row"], wiot_meta["col"]] = np.NaN
wiot_data.iloc[0 : wiot_meta["end_row"], wiot_meta["col"]] = np.nan
wiot_data.drop(wiot_empty_top_rows, axis=0, inplace=True)
wiot_data.drop(wiot_data.columns[wiot_marks["total_column"]], axis=1, inplace=True)
# at this stage row and column header should have the same size but
Expand Down
103 changes: 55 additions & 48 deletions tests/test_math.py
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@
# the function which should be tested here
from pymrio.tools.iomath import calc_A # noqa
from pymrio.tools.iomath import calc_accounts # noqa
from pymrio.tools.iomath import calc_As # noqa
from pymrio.tools.iomath import calc_B # noqa
from pymrio.tools.iomath import calc_e # noqa
from pymrio.tools.iomath import calc_F # noqa
from pymrio.tools.iomath import calc_F_Y # noqa
Expand Down Expand Up @@ -235,21 +235,28 @@ class IO_Data:
columns=_Z_multiindex,
)

As = pd.DataFrame(
B = pd.DataFrame(
data=[
[0.19607843, 0.0, 0.17241379, 0.1, 0.0, 0.12820513], # noqa
[0.09803922, 0.13333333, 0.05172414, 0.0, 0.19230769, 0.0], # noqa
[0.01960784, 0.0, 0.34482759, 0.0, 0.01923077, 0.0], # noqa
[0.11764706, 0.0, 0.06896552, 0.02, 0.0, 0.02564103], # noqa
[
0.19607843,
0.09803922,
0.33333333,
0.01960784,
0.11764706,
0.09803922,
0.1372549,
], # noqa
[0.0, 0.13333333, 0.0, 0.0, 0.33333333, 0.2], # noqa
[
0.17241379,
0.05172414,
0.34482759,
0.06896552,
0.03448276,
0.2,
0.38461538,
0.25641026,
0.0,
], # noqa
[0.1372549, 0.2, 0.0, 0.18, 0.01923077, 0.25641026], # noqa
[0.1, 0.0, 0.0, 0.02, 0.2, 0.18], # noqa
[0.0, 0.19230769, 0.01923077, 0.0, 0.38461538, 0.01923077], # noqa
[0.12820513, 0.0, 0.0, 0.02564103, 0.25641026, 0.25641026],
],
index=_Z_multiindex,
columns=_Z_multiindex,
Expand Down Expand Up @@ -314,52 +321,52 @@ class IO_Data:
data=[
[
1.33871463,
0.07482651,
0.38065717,
0.19175489,
0.04316348,
0.25230906,
0.28499301,
0.05727924,
0.1747153,
0.58648041,
0.38121958,
], # noqa
[
0.28499301,
0.07482651,
1.37164729,
0.22311379,
0.15732254,
0.44208094,
0.20700334,
0.02969614,
0.02185805,
0.93542302,
0.41222074,
], # noqa
[
0.05727924,
0.02969614,
0.38065717,
0.22311379,
1.54929428,
0.02349465,
0.05866155,
0.03091401,
0.15941084,
0.39473223,
0.17907022,
], # noqa
[
0.1747153,
0.02185805,
0.15941084,
0.19175489,
0.15732254,
0.02349465,
1.05420074,
0.01404088,
0.07131676,
0.60492361,
0.34854312,
], # noqa
[
0.58648041,
0.93542302,
0.39473223,
0.60492361,
0.04316348,
0.44208094,
0.05866155,
0.01404088,
1.95452858,
0.79595212,
0.18081884,
], # noqa
[
0.38121958,
0.41222074,
0.17907022,
0.34854312,
0.18081884,
0.25230906,
0.20700334,
0.03091401,
0.07131676,
0.79595212,
1.48492515,
], # noqa
],
],
index=_Z_multiindex,
columns=_Z_multiindex,
Expand Down Expand Up @@ -637,15 +644,15 @@ def test_calc_A_MRIO(td_small_MRIO):
pdt.assert_frame_equal(td_small_MRIO.A, calc_A(td_small_MRIO.Z, x_Tvalues))


def test_calc_As_MRIO(td_small_MRIO):
pdt.assert_frame_equal(td_small_MRIO.As, calc_As(td_small_MRIO.Z, td_small_MRIO.x))
def test_calc_B_MRIO(td_small_MRIO):
pdt.assert_frame_equal(td_small_MRIO.B, calc_B(td_small_MRIO.Z, td_small_MRIO.x))
# we also test the different methods to provide x:
x_values = td_small_MRIO.x.values
x_Tvalues = td_small_MRIO.x.T.values
x_series = pd.Series(td_small_MRIO.x.iloc[:, 0])
pdt.assert_frame_equal(td_small_MRIO.As, calc_As(td_small_MRIO.Z, x_series))
pdt.assert_frame_equal(td_small_MRIO.As, calc_As(td_small_MRIO.Z, x_values))
pdt.assert_frame_equal(td_small_MRIO.As, calc_As(td_small_MRIO.Z, x_Tvalues))
pdt.assert_frame_equal(td_small_MRIO.B, calc_B(td_small_MRIO.Z, x_series))
pdt.assert_frame_equal(td_small_MRIO.B, calc_B(td_small_MRIO.Z, x_values))
pdt.assert_frame_equal(td_small_MRIO.B, calc_B(td_small_MRIO.Z, x_Tvalues))


def test_calc_Z_MRIO(td_small_MRIO):
Expand All @@ -664,7 +671,7 @@ def test_calc_L_MRIO(td_small_MRIO):


def test_calc_G_MRIO(td_small_MRIO):
pdt.assert_frame_equal(td_small_MRIO.G, calc_G(td_small_MRIO.As))
pdt.assert_frame_equal(td_small_MRIO.G, calc_G(td_small_MRIO.B))


def test_calc_S_MRIO(td_small_MRIO):
Expand Down
Loading
Loading