Economic resilience from input-output susceptibility improves predictions of economic growth and recovery
Peter Klimek Sebastian Poledna Stefan Thurner







^1Section for Science of Complex Systems, CeMSIIS, Medical University of Vienna, Spitalgasse 23, A-1090, Austria

^2Complexity Science Hub Vienna, Josefstädter Strasse 39, A-1080 Vienna, Austria

^3IIASA, Schlossplatz 1, A-2361 Laxenburg, Austria

^4Santa Fe Institute, 1399 Hyde Park Road, Santa Fe, NM 85701, USA

^*To whom correspondence should be addressed. E-mail: stefan.thurner@muv.ac.at




Modern macroeconomic theories were unable to foresee the last Great Recession and could neither predict its prolonged duration nor the recovery rate. They are based on supply-demand equilibria that do not exist during recessionary shocks. Here we focus on resilience as a nonequilibrium property of networked production systems and develop a linear response theory for input-output economics. By calibrating the framework to data from 56 industrial sectors in 43 countries between 2000 and 2014, we find that the susceptibility of individual industrial sectors to economic shocks varies greatly across countries, sectors, and time. We show that susceptibility-based predictions that take sector- and country-specific recovery into account, outperform—by far—standard econometric growth-models. Our results are analytically rigorous, empirically testable, and flexible enough to address policy-relevant scenarios. We illustrate the latter by estimating the impact of recently imposed tariffs on US imports (steel and aluminum) on specific sectors across European countries.


 
 
Economic resilience from input–output susceptibility improves predictions of economic growth and recovery
    Peter Klimek^1,2, Sebastian Poledna^2,3, Stefan Thurner^1,2,3,4,*
    December 30, 2023
=========================================================================================================



In 2008 several advanced economies were hit by the largest recessionary shock in history <cit.>.
This Great Recession was followed by an even more remarkable event, namely a puzzlingly slow rate of economic recovery <cit.>.
Economists not only got the likelihood of a crisis of this severity wrong, as Paul Krugman famously noted, but also how fast we would recover from it <cit.>.
Efforts to understand the origin of this blind spot in economic theory and its failure to predict systemic events have fueled the interest in how economic systems absorb shocks and how they recover <cit.>.

Our lack of understanding of economic resilience <cit.> has been explained by a fundamental mismatch between macroeconomic theories and the reality of how markets work, especially in the presence of extreme events <cit.>.
General equilibrium theory holds that economic growth is characterized by a balance of demand and supply which results in prices that signal an overall equilibrium <cit.>.
However, the crisis was  a “story of contagion, interdependence, interaction, networks, and trust” <cit.> that led these equilibrium assumptions ad absurdum.
The inappropriate use of equilibrium concepts in economics 
in the context of extreme events was pointed out some time ago <cit.>.
So far, the only way to address and study economic non-equilibrium 
are highly stylized statistical models of money exchange <cit.> and computer simulations <cit.>.

In physics, non-equilibrium systems are equally hard to understand and control. 
Aside from some seminal contributions <cit.>, a unified framework for out-of-equilibrium phenomena has yet to be found <cit.>.
However, to understand how systems in equilibrium behave in response to shocks has been successfully addressed within the framework of linear response theory (LRT). 

According to LRT, an external force, X(t), acting on a system induces a proportional flux, J(t)=ρ X(t).
The proportionality is given by transport coefficients or susceptibilities, ρ, 
that are formally related to the decay of the system's equilibrium autocorrelation functions—the so-called Green–Kubo relations <cit.>.
LRT provides the theoretical basis for many linear phenomenological laws that constitute the core of each high school physics curriculum, such as Ohm's law, Newtonian viscosity, or magnetic and electric susceptibilities, see Supplementary Material (SM) Note S1 and Table S1 <ref>.

In the following we give an intuitive account of LRT.
Imagine someone hands you a serving tray with an elaborate house of cards on it.
Which card will fall first and cause the collapse of the house?
You try to answer this question by doing minimal damage: you slightly nudge the tray and observe how the cards respond.
If a tiny nudge moves certain cards, you might conclude that those are the first to fall if the tray was pushed harder.
The first cards that you observe to move you call the most “susceptible” to the shock (nudge).
If you observe no movements of cards whatsoever, you might be tempted to apply a stronger kick; the cards could be glued together.
A similar way of reasoning underlies the theory of linear response.
In the language of statistical mechanics, the tiny nudge that you initially apply plays the role of equilibrium fluctuations. 
These fluctuations may or may not move certain cards as a response—they induce a flux.
One then assumes that this response is proportional to the magnitude of the nudge, the proportionality being described by transport coefficient.
We will show that the same rationale can be used to study the response of Leontief IO economies <cit.> to large shocks.
The resulting framework is applied to IO data from 56 industrial sectors in 43 countries between 2000 and 2014 <cit.>; see SM Note S2.
We identify economic sectors with high susceptibility to production shocks in other parts of the economy. 
In the picture above, highly susceptible sectors correspond to cards that move first.
We show that the lack of recovery after the Great Recession can be related to the susceptibility of individual sectors.

How Leontief economies respond to shocks has been studied in a number of works briefly summarized as follows <cit.>.
Considering a (variant of a) Leontief IO economy, a shock is specified using varying degrees of external assumptions.
It is then studied how the economy relaxes to the old or new equilibrium configuration, unless yet another shock is assumed.
Our approach follows a completely different strategy.
We consider shocks that drive the economy away from its equilibrium into a non-equilibrium stationary state.
This stationary state is different from the original equilibrium state and the equilibrium state implied by the perturbed productivity or technology.
Instead, the new state is characterized by the system trying to achieve a balance between two opposing forces, namely (i) a relaxation to the (unaltered) equilibrium state and (ii) the direct and indirect influences of the external shock that actively drives the system away from equilibrium.
We call an economy in such a state a driven economy.
Market participants (sectors) in a driven IO economy incorporate the external shock in their production functions without altering their demand or required input from other sectors.
The perturbed output of these sectors then propagates along the IO network to other sectors, thereby driving the economy into a new non-equilibrium stationary state.

In this work we develop an analytic and empirically testable framework for the non-equilibrium response and recovery of severely disrupted economies. 
For the first time we formulate a theory of linear response for input–output (IO) economics <cit.>.
As in statistical physics, in the economic context LRT serves as a firm analytic link between the microscopic equilibrium fluctuations of a system and its macroscopic out-of-equilibrium response to large shocks.



The underlying ideas of LRT have been exploited in other contexts, such as the theory of linear time-invariant systems, 
with applications in signal processing and control theory <cit.>.
In econometrics, impulse response functions describe how external shocks drive macroeconomic variables such as output, consumption, or employment in vector autoregressive (VAR) models <cit.>.
Instead of studying the relaxation dynamics of macroeconomic variables within highly stylized VAR models,  
we focus on structural characteristics of dynamical IO matrices that capture the interactions between economic sectors.

 




§ RESULTS

   
   

 §.§ Obtaining economic susceptibilities from input–output data
  


Our formalism provides a quantitative and data-driven method to benchmark individual countries and 
production sectors in terms of their economic susceptibilities to shocks; see Methods.
To illustrate the method, we measure country and sector level economic susceptibilities respectively, 
by using the world input–output database (WIOD) <cit.>.
We consider data for 56 sectors in 43 countries between 2000 and 2014.
For each country, c, and year, t, we extract demands D_i(t), technical coefficients A_ij(t), and outputs Y_i(t), where subscripts refer to sectors.
Our aim is to compute the economic susceptibility matrix for a country and year, ρ^c_ij(t).
Here t denotes the year the data was taken to compute ρ^c_ij(t).
Based on data from t, we model output changes forward in time on a scale denoted by t'>t.
We numerically integrate the stochastic differential equation for a Leontief economy, Eq. (<ref>), in the absence of an external shock (X(t')=0 for all t' ≥ t). Now the time-lagged equilibrium correlation functions between two sectors in Eq. (<ref>) can be computed.
The entries in ρ^c_ij(t) correspond to the area under the curve of these correlation functions when plotted as a function of the time lag in Eq. (<ref>).
Susceptibilities of individual sectors ρ^c_i(t) are the column-wise sums of matrix ρ^c_ij(t).

Response curves of individual sectors are obtained by integrating the correlation functions under specific shocks.
In Fig. <ref> we assume an impulse demand shock of unit size applied at time t'=t, X_i(t') = δ(t'-t) in each sector i (A), 
leading to different response curves for each sector in the USA in 2014 (B).
The shock is applied at t=2014 and results in the same large decrease of output in all sectors immediately after t.
For t'>t, there appear substantial differences between sectors.
For some sectors the shock is amplified, such as for public administration, real estate activities, health, or wholesale trade.
Other sectors immediately start to rebound from the shock, for instance, the various manufacturing sectors.
The fastest rebound is observed for the construction sector, where production even exceeds the equilibrium level (0) for an extended period of time.
Overall, it can take up to six to ten years for each sector to return to its equilibrium state (sectoral recovery time).
Whether a shock is amplified or suppressed in a sector depends on the structure of the susceptibility matrix ρ^c_ij(t), see Fig. <ref>(C).
There we show the backbone of ρ^c_ij(t) (obtained after applying the disparity filter with p=0.05 <cit.>) as a directed weighted network.
Blue (red) links show positive (negative) susceptibilities. Node colors indicate groups of similar sectors, thickness of the node border is proportional to the sum of the weights of all incoming links, see SM Table S2; 
node sizes are inversely proportional to the values of the response functions in Fig. <ref>(B) at a particular point in time.
We show three snapshots of this network at the time when the initial shock is applied (t'=t), and one (t'=t+1), and six (t'=t+6) years afterwards. 
Figure <ref>(C) shows that some but not all of the sectors with a particularly strong shock amplification tend to be among those with a large 
number of incoming links (and weigths thereof), compare for instance the administration (large shock, many incoming links with strong weights) and construction (almost negative shock amplification, small number of incoming links) sectors.  

We apply the  above procedure for every year t (where the shock is applied), every country c, and every sector i, 
to compute a susceptibility value, ρ^c_i(t) .
The average country susceptibilities, ρ^c=⟨ρ^c_i(t) ⟩_i, are obtained by averaging ρ^c_i(t) over all sectors i and years t, see SM Fig. S1. 
The higher the values of ρ^c, the higher is the chance that any sector i in c will be impacted by a shock in any other sector j.

We find similar levels of susceptibility in a large number countries across Europe, North America, and China.
Substantially smaller susceptibilities are found for Croatia, Greece, Malta, and Luxembourg.
For those countries, our findings suggest a higher production concentration in a smaller number of sectors and consequently a smaller exposure to cascading impacts between different sectors (within the country).
At the other end of the spectrum, it is striking to see that four out of the five BRICS countries appear as the most susceptible countries, namely Russia, China, India, and Brazil; data for South Africa is not included in the WIOD due to the lack of available data with sufficient quality <cit.>.
This suggests that the sustained above-average growth of these countries in the last ten to twenty years
did not go along with the formation of resilient economic production structures.

In SM Fig. S2 we show the output-weighted average sector susceptibility, ρ_i=⟨ρ^c_i(t) ⟩_c, see also Table S2. 
Sectors with the highest susceptibilities include wholesale trade, administrative services, electricity, and financial service activities.
This means that if a country experiences an economic shock, those sectors are most likely to be affected by shocks in other sectors.
In contrast, we find that sectors like scientific research, activities of extraterritorial organizations, 
manufacture of transport equipment, or air and water transport are relatively immune to cascading events.




 §.§ Empirical validation of linear response theory for input–output economics
  

We now show to what extent the economic susceptibilbricity matrix ρ_ij is predictive of the size and direction of sectors' future output-changes.
First, it can be shown that the average size of sectoral output-changes can be predicted (out-of-sample) by means of sector-size-dependent random fluctuations.
To evaluate the linear response relation ⟨Δ Y_k ⟩_X = ρ_ki X_i it is necessary to specify the shock, X_i.
A particularly simple assumption is that X_i is itself 
noise with a magnitude proportional to the output of sector i, X_i = η_i Y_i, where η_i has the same expectation value in each sector, 
⟨η_i ⟩_i = η.
The hypothesis is that if ρ indeed captures structural characteristics of economies that relate to their recovery from shocks, 
one should be able to extract how violently Y_k fluctuates in the future, based on its current susceptibility.
To test this, for every sector k in every country c we consider its 
annual absolute output change, Y_k^c(t+1)-Y_k^c(t),  time-averaged over the range t=2000,…,2013, 
⟨Δ Y_k^c  ⟩_t = (1/13)∑_t=2000^2013( Y_k^c(t+1)-Y_k^c(t) ). 
According to the above hypothesis, ⟨Δ Y_k^c  ⟩_t should be a function of susceptibility, ρ^c_ki(t_0), and output, Y_i(t_0), in the year t_0=2000.   
We therefore test the quality of the out-of-sample prediction given by 

    ⟨Δ Y_k^c ⟩_t = η∑_i ρ^c_ki(t_0) Y_i(t_0)  .

Figure <ref> shows that this relation indeed holds (Pearson's correlation coefficient of r=0.83).
This correlation is substantially stronger than the correlation between output change ⟨Δ Y_k^c ⟩_t and output size Y_k(t_0) alone (r=0.56).
Performing a linear regression of ⟨Δ Y_k^c ⟩_t on ∑_i ρ^c_ki(t_0) Y_i(t_0) and Y_k(t_0) indeed yields a similar correlation as Eq. (<ref>) alone
 (giving r=0.83, with a regression coefficient of -0.000(2) for Y_k(t_0)).
Therefore, Eq. (<ref>) adequately captures output fluctuations that go beyond trivial sector size effects.
This confirms that the notion of economic susceptibility—the matrix ρ^c_ki(t_0)—coincides with (and is actually predictive of) the 
intuitive understanding that sectors with high susceptibility are those that are “more easily moved” by external events than low-susceptibility sectors.


  

We now show how the framework can be used to boost the quality of predictions of econometric timeseries models by extracting “implied shocks” from economic data.
Finally, we illustrate potential applications of our results by discussing estimates for economic impacts of recent tariffs imposed on US–EU trades in steel and aluminium. 



  §.§.§ Output predictions based on implied shocks


The linear response formalism requires the specification of a demand shock in one or several sectors.
Such shocks, however, can rarely be observed directly in the data.
If a step demand shock occurs at the beginning of year t, the data from t will not only contain the shock itself, 
but also of how the shock was “digested” by the economy during the year. 
As we have seen, recovery typically takes several years, see Fig. <ref>.

However, one can compute implied shocks from the data as follows (for clarity we omit the country index c from now on).
Consider the “truncated” susceptibility matrix ρ_ik(t,T), given by the area under the curve of the response function of i to a shock in k, evaluated until T years after the shock was applied, see Methods.
Assume that changes in output between year t and t+1 are due to a step demand shock X̃_i(t') = θ(t'-t)X̃_i, with θ the Heaviside step function, see Methods.
The size of this shock as implied by the output data from years t and t+1 can be estimated by using Eq. (<ref>),

    X̃_i = (ρ(t,T=1))_ik^-1( Y_k(t+1) -Y_k(t) )  .

We refer to X̃_i(t) as the implied shock at year t.
Positive (negative) output changes typically coincide with implied shocks that are of even larger (smaller) value, though some sectors defy these general trends, see SM Fig. S3. 

To test the validity of predictions of the linear response formalism, one can now take the implied shock from year t and estimate the 
output in year t+2 using Eq. (<ref>). Note that, by construction, the output in year t+1 is identical in the model and the data.
This yields a LRT timeseries model for individual countries with a driven economy, 

    ⟨ Y_k^LRT(t+2) ⟩_X =  Y_k(t) + ∫_0^2 (σ^-1)_ij⟨ Y_k(t+τ) Y_j(t) ⟩_0 X̃(t) dτ  .
 
The predictions of the LRT timeseries model are compared with expectations from econometric timeseries forecasting methods, 
in particular to results from autoregressive integrated moving average (ARIMA) models <cit.>; see SM Note S3 for a brief introduction.

The respective performance of the ARIMA and LRT model is evaluated by Pearson's correlation coefficient 
between the actual (empirically observed) and predicted output changes. 
For each year t and country c, we compute the correlation coefficient r^LRT(c,t), 
between the empirical output, Y_k(t), and the predictions from the LRT model (Y_k^LRT(t)) in Eq. (<ref>).

Similarly, we compute the the correlation coefficients r^ARIMA(c,t) for predictions from the ARIMA model 
(correlation of Y_k^ARIMA(t) with Y_k(t)). 
Values for r^LRT(c,t) and r^ARIMA(c,t) for are shown in SM Figs. S4 and S5.

The differences between the correlation coefficients of two different models for the same country and year, is referred to as the “predictability gain”,  PG(c,t) = r^LRT(c,t)-r^ARIMA(c,t), see Fig. <ref>(A).
Red (blue) values indicate that for the given country and year the LRT model performs better (poorer) than the ARIMA model.
For every year, we perform a t-test to reject the null hypothesis that the true mean of PG(c,t), taken over all countries, is zero (p<0.05).
The right panel in  Fig. <ref>(A) shows the PG(c,t) averages over all countries taken at each year with a 95% confidence interval (significant values are shown in black, non-significant in grey).
The bottom panels show the results for every country (significant vales are highlighted in black), and the histogram of PG(c,t) taken over all years and countries.
The LRT model performs significantly better than the ARIMA model in almost each year and country.
We find predictability gains of up to 100% and a p-value of p<10^-46 
to reject the null hypothesis that the true mean of the distribution of PG(c,t) is zero in this timespan.
Most intriguingly, for predictions from 2009 to 2010 (two years after the crisis occurred) the LRT model shows by far the largest predictability gains.
This result suggests that the LRT formalism works particularly well to describe the slow economic recovery during the Great Recession.

We design a further test, where it becomes harder for the LRT model to outcompete the ARIMA model, by comparing the 
out-of-sample predictions of the LRT model with the in-sample predictions of the ARIMA model. 
For this, we estimate the parameters of the sectoral ARIMA models over the entire timespan, from 2000 to 2014.
This should clearly stack the deck against the LRT model, as the ARIMA model is now calibrated using full timeseries information, 
in particular on the speed of economic recovery after the crisis.
Results are shown in Fig. <ref>(B). 
Overall, the LRT model again performs significantly better than the ARIMA model (p<10^-12).
The only exception is the prediction for 2009 where there is not clear which model is superior.
In this case the ARIMA model had the chance to “learn” the speed of the autoregression directly after the crisis.
In the following year, however, the LRT model shows the largest predictability gains, 
which again confirms that the LRT formalism is particularly useful to understand economic recovery.
Given that the ARIMA model has access to the full information of the timeseries, 
whereas the predictions of the LRT model are always taken entirely out-of-sample, 
these results once more confirm the superiority of the LRT formalism in describing the response of economies to large recessionary shocks.

Here we showed results for the ARIMA(1,1,1) model. However, qualitatively the same results 
are obtained (in many cases with even stronger relative performance of the LRT model) for other types of model. 
In particular, in SM we show results for the predictability gains PG(c,t) for an ARIMA(1,1,0) model (differenced first order autoregressive model), an ARIMA(0,1,1) model (exponential smoothing), and an ARIMA(1,0,1) model (first order autoregressive moving average model) 
in Figs. S6, S7, and S8.
We also confirmed that the LRT model performs vastly superior to a sectoral VAR model, see SI Note S4 and Fig. S9.


 




  §.§.§ Indirect effects of the US–EU trade war


Finally we show how the LRT model can be used to estimate the economic impact of instances such as the currently escalating trade war 
between the EU and US <cit.>.
Starting from June 1, 2018, the US imposes a 25% tariff on steel and a 10% tariff on aluminum imports from member countries of the EU.
These tariffs are expected to lead to direct negative effects on EU steel and aluminum producers, 
which could be further amplified by other countries that redirect their exports from the US to the EU.
The indirect effects of these tariffs, however, are not so clear.
Increased supply of steel and aluminum in the EU might lead to a decrease in price with positive effects 
on industries that require those metals as inputs.
In the LRT model, the US tariffs impose a negative export demand shock on the manufacturing sector of basic metals (ISIC Code C24) on EU countries and a positive demand shock on the US.
We assume that US demand in this category will reduce by 100% for European countries (and US domestic final consumption will increase accordingly) and estimate the resulting changes to the sectoral outputs using the linear relationship in Eq. (<ref>).
Note that the impacts of shocks with an arbitrary size of x% of current export demand can simply be estimated by multiplying these results by x/100.
Results for ⟨Δ Y_k ⟩_X obtained from Eq. (<ref>) using the most recent data available in WIOD (t=2014) are shown in Fig. (<ref>)(A) for the 25 largest sectors.
In general, output changes fall in the range between ±0.5%.
In European countries, positive effects are particularly strong in the manufacturing sectors (motor vehicles, computers, electronics, machinery, or electrical equipment), whereas there are consistently negative indirect effects for the energy sector.
These findings are consistent with an expectation of positive effects further down the supply chain of steel and aluminum (due to price decreases).
Decreases in the output of steel and aluminium production on the other hand coincide with a decrease in energy consumption.
It is also apparent that the indirect effects are distributed unequally across countries.
Manufacturing activities in Germany, Greece, or Ireland show consistently increased levels of output.
Indirect effects in the US often show opposite signs compared to the impact on European countries.
We find that negative indirect effects prevail for fabricated metal products and motor vehicles while the electricity sector, land transport, and wholesale trade experience positive effects. 
By summing the expected output-changes (in USD terms) over each sector in a country, we obtain the aggregated indirect effects; Fig. (<ref>)(B).
Overall, almost all countries experience positive indirect effects with output increases of up to several billion dollars; the exceptions being Spain, Finland, Italy, and Romania.
Our framework suggests that these countries might either (i) depend to a higher extent on sectors that provide input to the manufacture of basic metals (such as electricity), (ii) lack sectors that can profit from an increased supply of basic metals, or (iii) both of the former might be the case.  
Also, note that for European countries with positive aggregated indirect effects, these effects are typically outweighed by negative direct effects from the tariffs. 
Figure (<ref>)(C) shows the temporal impact (response curves) for Germany, for the a step demand shock for aluminum and steel. 








§ DISCUSSION


We developed the theory of linear response for IO economies to quantify the resilience of national economies to production shocks.
We established an analytic link between stationary output fluctuations and their out-of-equilibrium behavior.
In particular, we derived the Green–Kubo relations for Leontief IO models, in full analogy to a wide range of physical phenomena, 
ranging from electrical and magnetic susceptibilities to shear viscosity and electrical resistance.
Our framework can be applied to other types of IO model, as long as they are (i) linear and (ii) permit a stationary solution.
This includes IO models that use a higher geographic resolution (i.e. regional IO models), 
but also several of their generalizations, such as environmentally extended IO models <cit.>, 
or commodity-by-industry IO models <cit.>.

The central result of our work is a linear “phenomenological law” between demand shock, X_i, and the induced output change Δ Y_k, namely that 
⟨Δ Y_k ⟩_X = ρ_ki X_i, with ρ being a sector-by-sector matrix of economic susceptibilities.
The output change ⟨Δ Y_k ⟩_X  characterizes a driven economy in a nonequilibrium stationary state.
The original equilibrium state, (𝕀 -𝐀)^-1𝐃, is recovered for X_i(t)=0 for all i and t.
The LRT solution ⟨Δ Y_k ⟩_X is also fundamentally different from the “perturbed equilibrium state” implied by a step demand shock of the form 𝐃_𝐏 = 𝐃+𝐗 with 𝐗(t) =𝐗θ(t), 
namely the perturbed equilibrium state  (𝕀 -𝐀)^-1𝐃_𝐏.
To see the difference, note that in LRT the expectation values are taken over the probability density function of the stationary solution of Ẏ = (𝐀-𝕀)𝐘 + 𝐃 + 𝐅(t) (from which the non-equilibrium expectation value ⟨Δ Y_k ⟩_X is estimated), 
whereas the perturbed equilibrium state would be given by expectation values using the stationary solution of Ẏ = (𝐀-𝕀)𝐘 + 𝐃_𝐏 + 𝐅(t) as probability measure; see also SM Note S6 and Fig. S10.

We demonstrated that our measures for economic susceptibility that can be derived from data, 
are indeed predictive of future output fluctuations, even when no knowledge of future shocks is available.
This finding corroborates that sectors with high susceptibility are indeed those that tend to be 
“more easily movable” by external shocks than low-susceptibility sectors.
We showed that out-of-sample predictions from the LRT model consistently outperform standard econometric forecasting methods, 
such as different types of ARIMA model.
Predictions of the LRT model work particularly well in the years that followed the recent financial crisis.
This suggests that the LRT formalism allows us to get an analytic and quantitative understanding 
of the slow economic recovery of certain countries in the wake of the Great Recession.
Because of the versatility and conceptual simplicity of input–output models, our framework can lead to more accurate quantitative estimates for the impact of disruptive events in various applications and scales, ranging from global recessions to regional, critical infrastructure systems.
We illustrate the practical usefulness of our approach in providing concrete estimates 
for the indirect effects of the currently escalating US–EU trade war.
In particular, we considered a negative export demand shock on the manufacturing sector of basic metals on EU countries and a corresponding positive demand shock on the US.
We find that in European countries there is a trend toward positive indirect effects for manufacturing sectors further down the supply chain from basic metals, whereas electricity outputs show negative indirect effects.
In the US we find similar results with reversed signs; positive (negative) effects moving further up (down) along the supply chain.

A limitation of the Leontief IO model that extends to our work is that prices play no role in the model.
Firms in real economies can respond to shocks by adjusting produced quantities as well as prices.
It therefore remains to be seen how prices can be incorporated in the LRT framework, i.e., within a linear time-invariant formulation of the underlying microscopic dynamics.
Besides linearity and time-invariance, our approach also assumes an external shock that may depend only on time and for which we only consider first order correction with respect to the unperturbed state.

In summary, in this work we extended current mainstream economic theories to out-of-equilibrium situations 
in a way that is analytically rigorous, empirically testable, and flexible enough to immediately address a wide range of scenarios with a direct political relevance, such as identifying those parts of a country's economy that are particularly vulnerable in a trade war.      






§ METHODS




 §.§ Linear response theory of input-output economics
  

Consider an economy with N sectors, each sector producing Y_i units of a single homogeneous good.
Assume that sector j requires A_ij units from sector i as input to produce one unit itself, giving the so-called technical coefficients A_ij.
Each sector sells some of its output to consumers, the demand D_i.
The open Leontief IO model, the standard model in economics to depict and analyze inter-sectoral relationships, 
assumes linear production functions given by 𝐘 = 𝐀𝐘+𝐃 (matrix notation).  
The stationary (equilibrium) state of this economy is given by 𝐘^0=(𝕀 -𝐀)^-1𝐃 (𝕀 being the N-dimensional identity matrix).
For the time evolution of an economy in its stationary state, 
assuming that differences in dynamic demand 𝐀𝐘+𝐃 and dynamic production 𝐘 are compensated by 
production changes, this model gives the differential equation <cit.>

    Ẏ = (𝐀-𝕀)𝐘 + 𝐃  .

We assume (i) that each sector i experiences a time-dependent demand shock, X_i(t), and (ii) the presence of multivariate white noise, i.e. a stochastic force, F_i(t).
In the picture of the example of the house of cards given above, the noise F_i(t) represents the “tiny nudge” that we apply to the serving tray to understand if the house would survive a much larger shock, X_i(t).
More formally, the nudge consists of noise with mean value ⟨ F_i(t) ⟩_0=0, 
and covariance ⟨ F_i(t) F_j(s) ⟩_0=ν_ijδ(t-s).
Here, δ(x) denotes the Dirac-delta function, ν is a matrix of constants, and ⟨𝐱(t) ⟩_0 = ∫ d^N𝐘𝐱(𝐘) f_0(𝐘)
 is the equilibrium expectation value of the function 𝐱(t), evaluated in the absence of an external force (𝐗(t)=0), with f_0(𝐘) being the probability distribution to find a given value of 𝐘 under noise F_i(t).
This leads to the stochastic differential equation, 

    Ẏ = (𝐀-𝕀)𝐘 + 𝐃 + 𝐗(t) + 𝐅(t)   .

From the central limit theorem it follows immediately that the stationary or equilibrium solution f_0(𝐘) in the absence of external shocks (𝐗(t)=0) of Eq. (<ref>) is given by a multivariate normal distribution with covariance σ_ij = lim_t→∞⟨ Y_i(t) Y_j(t) ⟩_0.
In the presence of external shocks, i.e., for 𝐗(t) ≠ 0, a solution of Eq. (<ref>) with first order corrections from the shock can be obtained using LRT, see SM Note S5.
We denote the expectation value for the output change of sector k with nonzero shock 𝐗(t) by ⟨Δ Y_k(t) ⟩_X ≡⟨ Y_k(t) ⟩_X-Y_k^0.
That is, averages with a subscript 0 refer to values taken at equilibrium, 
whereas averages with a subscript X refer to out-of-equilibrium properties.
Following LRT <cit.>, we get the general solution for the time evolution of the output changes, ⟨Δ Y_k(t) ⟩_X,

    ⟨Δ Y_k(t) ⟩_X = ∫_-∞^t (σ^-1)_ij⟨ Y_k(τ) Y_j(0) ⟩_0 X_i(τ) dτ  .

Remarkably, we have related the out-of-equilibrium response of the sectoral outputs, 
⟨Δ Y_k(t) ⟩_X, to their correlation functions taken at equilibrium.
Equation <ref> characterizes the state of a driven economy.

For certain types of demand shock, the resulting output change takes a particularly simple form.
For an impulse demand shock, X_i(t) = δ(t) X_i, we get 

    ⟨Δ Y_k^pulse(t) ⟩_X = (σ^-1)_ij⟨ Y_k(t) Y_j(0) ⟩_0 X_i   .

For a step demand shock, X_i(t) = θ(t) X_i with the Heaviside step function θ(t ≥ 0)=1 and θ(t<0)=0, we get 

    ⟨Δ Y_k^step(t) ⟩_X = ∫_0^t (σ^-1)_ij⟨ Y_k(τ) Y_j(0) ⟩_0 X_i dτ  .

For t≫0 we obtain the linear relation

    ⟨Δ Y_k  ⟩_X = ρ_ki X_i,  with ρ_ki = ∫_0^∞ (σ^-1)_ij⟨ Y_k(τ) Y_j(0) ⟩_0 dτ  ,

where we introduced the economic susceptibility ρ,  
in full analogy to the derivation of electric or magnetic susceptibilities in statistical mechanics, see SM Table 1.
The economic susceptibility ρ_ki has the precise meaning of output-change in sector k, given that a step demand shock of unit size occurs in sector i.  

In this paper we encounter different types of susceptibility, depending on how averages are taken. 
In particular we will use the following definitions: 

The N × N susceptibility matrix of a country c at year t is defined by 

    ρ_ij^c(t)=  ∫_0^∞ (σ^-1)_ij⟨ Y_k^c(t+τ) Y_j^c(t) ⟩_0 dτ  ,
 
where the output Y_k^c(t), technical coefficients A_ij^c(t), and the demand D_i^c(t) of a particular country c, 
are read off the WIOD <cit.>.
A truncated version ρ_ij^c(t,T) of this susceptibility matrix is obtained by taking t+T, T>0, as the upper boundary of the integration range in Eq. (<ref>).
The susceptibility of a sector i in country c at year t, ρ_i^c(t), is defined as the corresponding column sum of the susceptibility matrix, ρ_i^c(t) = ∑_j ρ_ij^c(t).

We define the averaged country susceptibility as the average of the sector susceptibility taken over all N sectors and N_t=15 years, 
ρ^c = (N_tN)^-1∑_i,tρ_i^c(t). 
The output-weighted average sector susceptibility, ρ_i is defined as, 
ρ_i = (N_t N_c ∑_t,c Y^c_i(t))^-1∑_t,c Y^c_i(t) ρ_i^c(t), where N_c=43 is the number of countries in the data. 




 §.§ Acknowledgments
  
We thank M. Miess and A. Pichler for helpful discussions, Johannes Sorger for help with the visualizations, and acknowledge support from the European Commission, H2020 SmartResilience No. 700621, FFG Project 857136, and OeNB  Jubiläumsfond project 17795.



 §.§ Author contributions

PK and ST designed research, PK performed research and analyzed data, SP contributed timeseries models, PK and ST wrote the paper.



 §.§ Conflict of interest

The authors declare no conflict of interest.



 §.§ Data availability

The study is based on the 2016 release of the World Input–Output Tables <cit.> (see also http://www.wiod.orghttp://www.wiod.org, accessed 15 January 2019) .



 §.§ Code availability

Code is available upon request directly from the authors.


99

IMF09 International Monetary Fund, World Economic Outlook 2009, Crisis and Recovery (International Monetary Fund, Washington, DC, USA, 2009).

Barnett14 Barnett A, Batten S, Chiu A, Franklin J, Sebastiá-Barriel, The UK productivity puzzle. Bank of England Quarterly Bulletin 54(2), 114–28, 2014.

Bryson14 Bryson A, Forth J,  Askenazy P, Productivity Puzzles in Europe: a Comparison of the UK, France, Germany and Spain. NIESR Discussion Paper 448, 2014.

Schweitzer09 Schweitzer F, Fagiolo G, Sornette D, Vega-Redondo F, Vespignani A, White D, Economic networks: the new challenges. Science 325, 422–25, 2009.

Garas10 Garas A, Argyrakis P, Rozenblat C, Tomassini M, Havlin S, Worldwide spreading of economic crisis. New Journal of Physics 12, 113043, 2010.

Haldane11 Haldane A, May R, Systemic risk in banking ecosystems. Nature 469, 351–5, 2011.

Farmer12 Farmer JD, Gallegati M, Hommes C, Kirman A, Ormerod P, Cincotti S, Sanchez A, Helbing D, A complex systems approach to constructing better models for managing financial markets and the economy. Eur. Phys. J. Special Topics 214, 295–324 (2012).

Martin11 Martin R, Regional economic resilience, hysteresis and recessionary shocks. Journal of economic geography 12(1), 1–32, 2011.

Kirman10 Kirman A, The economic crisis is a crisis for economic theory. CESifo Economic Studies 56(4), 498–35, 2010.

Arrow54 Arrow KJ, Debreu G, Existence of an equilibrium for a competitive economy. Econometrica 22(3), 265–90, 1954.

Smets07 Smets F, Wouters R, Shocks and frictions in US business cycles: A Bayesian DSGE approach. The American Economic Review 97(3), 586–606, 2007.

Mitra08 Mitra-Kahn BH, Debunking the myths of computable general equilibrium models. SCEPA Working Paper 2008-1, 2008.

Farmer05 Farmer JD, Smith DE, Shubik M, Is economics the next physical science? Physics Today 58(9), 37–42, 2005.

Onsager31 Onsager L, Reciprocal relations in irreversible processes. I. Physical Review 37(4), 405, 1931.

Kubo66 Kubo R, The fluctuation-dissipation theorem. Reports on progress in physics 29(1), 255, 1966.

Nicolis77 Nicolis G, Prigogine I, Self-organization in nonequilibrium systems (Wiley, New York, USA, 1977).

Hanel14 Hanel R, Thurner S, Gell-Mann M, How multiplicity determines entropy and the derivation of the maximum entropy principle for complex systems. Proceedings of the National Academy of Sciences USA 111(19), 6905–10, 2014.

Jarzynski15 Jarzynski C, Diverse phenomena, common themes. Nature Physics 11(2), 105–7, 2015.

Gallegatti06 Gallegatti M, Keen S, Lux T, Ormerod P, Worrying trends in econophysics. Physica A 370(1), 1–6, 2006.

Farmer09 Farmer JD, Foley D, The economy needs agent-based modelling. Nature 460, 685–6, 2009.

Klimek14 Klimek P, Poledna S, Farmer JD, Thurner S, To bail-out or to bail-in? Answers from an agent-based model. Journal of Economic Dynamics and Control 50, 144–54, 2014.

Evans90 Evans DJ, Morriss GP, Statistical mechanics of nonequilibrium liquids (Theoretical Chemistry Monograph Series, Academic Press, London, UK, 1990).

Leontief86 Leontief W, Input-output economic (Oxford University Press, New York, USA, 1986).

Green54 Green MS, Markoff random processes and the statistical mechanics of time-dependent phenomena. II. Irreversible processes in fluids.  J. Chem. Phys 22, 398–413, 1954.

Kubo57 Kubo R, Statistical-mechanical theory of irreversible processes. I.
General theory and simple applications to magnetic and conduction problems. J. Phys. Soc. Jpn. 12, 570–586, 1957. 


Timmer15 Timmer MP, Dietzenbacher E, Los B, Stehrer R, de Vries GJ, An illustrated user guide to the world input-output database: the case of global automotive production. Review of International Economics 23, 575–605, 2015.

Long83 Long JB, Plosser CI, Real Business Cycles. Journal of Political Economy 91(1), 39–69, 1983.

Haimes01 Haimes YY, Jiang P, Leontief-based model of risk in complex interconnected infrastructures. Journal of Infrastructure Systems 7(1), 1–12, 2001.

Contreras14 Contreras MGA, Fagiolo G, Propagation of economic shocks in input–output networks: A cross-country analysis. Physical Review E 90(6), 062812, 2014.

Acemoglu16 Acemoglue D, Akcigit U, Kerr W, Networks and the macroeconomy: An empirical exploration. NBER Macroeconomics Annual 30(1), 273–335, 2016.

He17 He P, Ng TS, Su B, Energy-economic recovery resilience with input–output linear programming models. Energy Economics 68, 177–191, 2017.

Phillips13 Phillips CL, Parr JM, Riskin EA, Signals, systems, and transforms (Prentice Hall, 2013). 

Pesaran98 Pesaran HH, Shin Y, Generalized impulse response analysis in linear multivariate models. Economic letters 58(1), 17–29, 1998.

Lutkepohl08 Lütkepohl M, Impulse response function. In: Palgrave Macmillan (eds) The New Palgrave Dictionary of Economics (Palgrave Macmillan, London, 2008).

Serrano09 Angeles Serrano M, Boguna M, Vespignani A, Extracting the multiscale backbone of complex weighted networks. Proceedings of the National Academy of Sciences USA 106(16), 6483-8, 2009.

Timmer16 Timmer M, Los B, Stehrer R, de Vries G, An anatomy of the global trade slowdown based on the WIOD 2016 release (No. GD-162) (Groningen Growth and Development Centre, University of Groningen, 2016). 

Box70 Box GE, Jenkins GM, Reinsel GC, Ljung GM, Time series analysis: forecasting and control (John Wiley & Sons, New Jersey, USA, 1970).

Chu18 Chu B, How economically damaging will Trump's steel tariffs be? The Independent, 31 May 2018. (https://www.independent.co.uk/news/business/analysis-and-features/trump-us-steel-tariff-eu-damaging-explained-trade-war-uk-jobs-a8377871.htmlhttps://www.independent.co.uk, accessed 18 July 2018).

Wiedmann13 Wiedmann TO, Schandl H, Lenzen M, Moran D, Suh S, West J, Kanemoto K, The material footprint of nations. Proceedings of the National Academy of Sciences USA 112(20), 6271-6, 2013.

Miller09 Miller RE, Blair PD, Input-output analysis: foundations and extensions (Cambridge University Press, Cambridge, UK, 2009).
















§ SUPPLEMENTARY MATERIAL





 §.§ Note S1:  General information on linear response theory


The basic idea of LRT is to compute the explicit time evolution of an external shock to a system in equilibrium. 
Originally, LRT has been formulated with thermal and mechanical transport processes in mind.
The idea is that a thermodynamic system (e.g. a canister of gas, a block of metal, a piece of lumber, ...) can be prevented from relaxing to an equilibrium configuration by an external field.
This field acts either directly on the system's microscopic constituents (think of an electric or magnetic field influencing the electrons in a block of metal) or that acts on the system boundaries (e.g. by generating a temperature gradient or by deforming the container of a system).
In both situations the external field drives the system from equilibrium by means of induced transport processes.
An electric field induces a current (e.g. transport of electrons) while a temperature gradient induces a heat flow (transport of particles with high kinetic energy).
The central result of LRT is that a specific current (or flux) can be associated with such an external field.
Furthermore, the magnitude of this flux is proportional to the magnitude of the external field.

Melville S. Green <cit.> and Ryogbo Kubo <cit.> independently found the exact mathematical expression that relates the field to its induced flux, the so-called transport coefficients ρ.
These expressions are now called Green–Kubo relations. 
Note that ρ can be a scalar, a matrix, or a higher dimensional tensor, see Table S1.
Their point of departure was that thermodynamic systems undergo a certain amount of fluctuations even at equilibrium.
Given a perturbation in one of the microscopic degrees of freedom of the system (e.g. a fluctuation in a particle's velocity along a certain direction), the microscopic equations of motion allow to compute how strong and at which time other degrees of freedom (other particles) will “feel” this fluctuation.
LRT assumes that other, larger perturbations will be felt in a similar way, i.e. with a similar delay and attenuation.
Note that in many cases it is not even remotely possible to get a firm analytical understanding of the thermodynamics of a system once it is removed from equilibrium, given that it makes sense to talk about thermodynamics in such regimes at all.
LRT provides one of the few analytic links that allows us to extrapolate out-of-equilibrium properties of a system from its behavior in equilibrium and therefore one of the paramount achievements of 20th century statistical mechanics.




 §.§ Note S2: Description of the World Input–Output Database


In this work we use the World Input–Output Database (WIOD) release from November 2016 (http://www.wiod.org/release16).
This release consists of several data tables that cover 28 EU countries and 15 other major countries (Australia, Brazil, Canada, China, India, Indonesia, Japan, Mexico, Norway, Russia, South Korea, Switzerland, Taiwan, Turkey, and the United States) from 2000 to 2014.
The 15 non-European countries were chosen in order to cover all major parts of the world economy while ensuring that the data is available at a sufficient quality <cit.>.
The resulting 43 countries account for more than 85% of world GDP <cit.>.
For each year and country, WIOD contains data on 56 sectors according to the 2-digit ISIC revision 4 level.

The basic construction steps of the WIOD are as follows <cit.>.
First, timeseries of national supply and use tables are constructed for each country based solely on national account statistics.
With the use of international trade databases, these timeseries are disaggregated into imports by country of origin and use category.
Finally, the resulting data items are integrated into a single input–output table.
The construction of WIOD uses solely national account statistics as input for the construction of national supply and use tables.
This ensures that (i) WIOD is fully consistent with national accounts and (ii) the use of a consistent methodology in categorizing products and services across time.
This comes at the cost of a smaller coverage of countries and years with respect to other comparable databases that do not ensure such consistencies, see also <cit.> and <cit.>.





 §.§ Note S3: A brief introduction to ARIMA models


ARIMA models are a quite general class of models for stationary timeseries, that is for timeseries that are characterized by statistical properties that do not change over time <cit.> (we only consider nonseasonable models here).
Stationary timeseries show variations around their means that are of constant amplitude and that look similar over time.
That is, given a snippet from a timeseries it is impossible to deduce when exactly the snippet was observed, e.g. after a couple of time steps or after observing the process for a million years.
ARIMA models serve the purpose of forecasting timeseries based on observations of previous values that were collected within a specific time interval.
It turns out that stationary timeseries have a quite limited set of properties that fully characterize them.  
These properties are encoded in the acronym ARIMA, which stands for auto-regressive integrated moving average.
Auto-regression (AR) means that the next value of a timeseries can be predicted as a multiple of its prior values (maybe plus a constant term).
Integration (I) means the timeseries should be forecasted by considering (also) differences between prior values, rather than (only) the values themselves.
Moving average (MA) finally indicates that the timeseries is best predicted by considering not the last observed value (as in random walks without memory), but rather the average over several previously observed values.
Each of these three properties, AR, I, and MA, can be present up to a specific order.
For AR this order is the number p of auto-regressive terms, for I the number d of differences that are relevant, and for MA it is the number q of past observations that are included in the average.
An ARIMA model is specified by the choice of these three numbers as an ARIMA(p,d,q) model. 
If the number of observations available to calibrate the ARIMA model is low (as it is the case for our work), it is often not meaningful to consider orders of parameters higher than one or two, as the corresponding coefficients for the higher order correction terms can not reliably be estimated.




 §.§ Note S4: Comparison of LRT and VAR models

Vector autoregression (VAR) models can be used to describe stochastic processes driven by linear dependencies on multiple other variables.
In our case we can consider the sectoral outputs Y_k^c(t) of country c as a stochastic process that may depend on the output of all other N sectors in the same country at the previous timestep, the outputs Y_j^c(t-1) with j∈{1,…,N}.
Note that ARIMA models assume that output changes in each sector can be predicted based on past values of the same sector only.
By benchmarking the LRT model against a VAR model we can go beyond this limitation and test the LRT framework against a regression model that captures the structure of inter-industry dependencies in a more comprehensively way.
The drawbacks of VAR models are that they do not scale well with system size.
For a first order VAR model with N=56 sectors one needs to estimate no less than N^2=3,136 different parameters.
It is therefore completely hopeless to specify a sectoral VAR model from data given that we have only 15 observations for each sector.
However, we can use a strategy to calibrate a sectoral VAR model similar to how we measured the response functions in the LRT framework.
By assuming that the economy can be represented by the stochastic differential equation,

    Ẏ = (A-𝕀)Y + D  + F(t)  ,

we generated 10,000 synthetic observations for the output of each sector.
Using supercomputing resources (the Vienna Scientific Cluster 3, one of the hundred fastest computers worldwide (http://vsc.ac.at/systems/vsc-3/, accessed Sep 14, 2018.) we were then able to estimate the parameters of a first order VAR model using data from one year (t=2000) for each country.
By construction, the model has no linear trend.
It remains to estimate the entries of the autoregression matrix AR^c for each country c, given by Y_k^c(t+1) = AR^c Y_k^c(t) + e_k^c (e_k^c being the intercepts).
We then compare the forecasts of this VAR models with predictions of the LRT model, similar to the comparison with ARIMA models.
The results are shown in Fig. S7.
For almost all years and countries, the LRT model performs vastly superior than the VAR model (p<10^-117).
The advantage of the LRT model is least pronounced in the years 2009 and 2010 (where LRT still performs significantly better), whereas other years show predictability gains in the range between one and two.
Note that these results are based on a VAR model calibrated in 2000 only due to the costly requirement for supercomputing resources.
For a single country (Germany), we evaluated the VAR model also for several years afterwards and in particular for the crisis year 2008.
Still, we found qualitatively the same result as shown in Fig. S7 with the majority of predictability gains lying in the range between one and two.




 §.§ Note S5: Derivation of Green–Kubo relations for input–output economics
 

Let y ≡{ y_i(t) }, i=1,…,N be the stochastic variables describing the time evolution of a dynamical system.
In the following we assume the system to be linear and time invariant.
Linearity means that the time evolution of the system is governed by a linear operator.
Time invariance means that the response of the system does not depend on the time at which we apply the shock.
We describe the state of the system by its probability density function, f(y,t).
The time evolution of the system is described by a linear time operator, L_0(y), for instance a Fokker-Planck operator.
In general, L_0(y) can be any linear operator with a stationary solution for f(y,t) for the LRT framework to be applicable.
Further, the system is perturbed by a time-depending external field X(t), i.e. we allow the external field to vary arbitrarily over time but not as a function of the stochastic variables y..
Under perturbation, the time evolution of f(y,t) is then described by a perturbed operator Λ_X(y,t), as

    Λ_X(y,t)     =     L_0(y) + L_X (y,t)   , 
    
    L_X(y,t)     =     L_X(y) X(t)   ,

such that 

    ∂ f/∂ t = L_0(y,t) f(y,t)   .

Let us denote the stationary solution for L_0 by f_0(y), given by L_0(y) f_0(y) = 0.




  §.§.§ Response Function


Consider a small perturbation of the stationary solution of the form f(y,t) = f_0(y) + Δ f(y,t).
For the time evolution we have then

    ḟ(y,t)     =    Δḟ(y,t) = Λ_X(y,t) f(y,t) 
        =    ( L_0(y) +L_X(y,t) ) ( f_0(y) + Δ f(y,t) )  
        =     L_0(y) Δ f(y,t) + L_X(y,t) f_0(y)   .

The above equation can be solved using well-known properties of the Laplace transformation, given by f̃ (s) = ∫_0^t dt e^-stf(t).
For the l.h.s. of the equation  ḟ(y,t) = L_0(y) Δ f(y,t) + L_X(y,t) f_0(y) we get,

    ḟ(y,t) = Δḟ(y,t) = s Δf̃(s) - Δ f(0)   ,

where we can use that Δ f(0) =0.
The r.h.s. simply transforms to L_0(y) Δ f(y,t) + L_X(y,t) f_0(y) = L_0(y) Δf̃(s) + f_0(y) L̃_X.
We obtain, 

    Δf̃ (s) = L̃_̃X̃ f_0(y)/s-L_0(y)  .

 Now we apply an inverse Laplace transformation to the l.h.s. and r.h.s. of Eq. <ref>.
 Therefore, we make use of the Laplace transform properties that (i) for a product of functions we have f̃(s) g̃(s) =∫ dτ f(τ) g(t-τ) and (ii) that the exponential function, e^at transforms to 1/(s-a).
 Using these properties, Eq. <ref> becomes

    Δ f(y,t) = ∫_-∞^t dτ e^L_0(y) · (t-τ) L_X(τ) f_0(y)   .


We are now interested in the expectation value of any dynamical variable B(t) under the considered perturbations.
Formally, we have

    ⟨ B(t) ⟩    =    ∫ d^N y B(y) f(y,t) 
        =    ⟨ B ⟩_0 + ∫ d^N y B(y) Δ f(y,t)   .

By plugging in the solution for Δ f(y,t) obtained in Eq. <ref> we get

    ⟨ B(t) ⟩ =  ⟨ B ⟩_0 + ∫ _-∞^t R_B,X (t-τ) X(τ) d τ  ,

where R_B,X(t) is the response function

    R_B,X(t)     =    ∫ d^N y  B(y) e^L_0(y) t L_X(y) f_0(y)     t≥ 0, 
    
    R_B,X(t)     =     0   for  t< 0.

R_B,X(t) describes the response of the variable or observable B(t) to the external force X(t).



  §.§.§ Correlation Functions


The response function can be expressed in terms of equilibrium correlation functions c_A,B(τ) between an observable A(y(t)) and B(y(t)) with lag τ.
First, observe that

    c_A,B(τ)     =    ⟨ A(y(τ)) B(y(0)) ⟩
        =    ∫ d^Ny  d^N y_0  A(y) B(y_0) P(y,τ; q_0,0)   ,

where P(y,τ; y_0,0) is the joint probability that the system is in configuration y_0 at time t=0 and in configuration y at t.
This can be written as 

    P(y,τ; y_0,0) = P(y,τ| y_0,0) P(y_0,0)   .

Here, P(y,τ| y_0,0) is the conditional probability that the system evolved from state y_0 at t=0 to state y at t.
Formally, this conditional probability is equal to the propagator e^L_0(y)·τδ(y-y_0).
Note that P(y_0,0) is the stationary solution f_0(y_0).
For the equilibrium correlation function this means that

    c_A,B(τ) = ∫ d^Ny  A(y) e^L_0(y)·τ B(y) f_0(y)   .

Define the function A(y) as A(y) ≡ f_0^-1(y) L_X(y) f_0(y) and introduce the generalized potential ϕ(y) through the relation

    f_0(y) ≡𝒩e^- ϕ(y)  ,

with 𝒩 as normalization constant.
With these definitions we re-write the correlation function as

    c_B,A(τ) =  ∫ d^N y  B(y) e^L_0(y) t L_X(y) f_0(y)   ,

and obtain the result that

    c_B,A(τ) =  R_B,X (τ)   .




  §.§.§ Application to stochastic processes


Consider a stochastic dynamical system with an external force X_i(t) acting on component i and the stochastic force F_R,i(t) with the properties 
⟨ F_R,i(t) ⟩ = 0 and ⟨ F_R,i(t)· F_R,j(t') ⟩ = ϵ_ijδ(t-t').
In particular, we consider linear systems of the type,

    ẏ_i + ∑_j=1^N γ_ijy_j = F_R,i(t)+X_i(t)   .

Clearly, the stochastic IO dynamics that we consider in the main text is of the above type.
To this end, note that one can always use a variable transformation of the form Y → y = Y - (𝕀 - A)^-1D to rewrite the IO model into a homogeneous differential equation.
From now on, we will use the summation convention.
The covariance in the stationary regime (in the absence of external shocks X(t)) between y_i(t) and y_j(t) is denoted by σ_ij(t)= ⟨ ( y_i(t) -⟨ y_i ⟩ ) ( y_j(t) - ⟨ y_j ⟩ ) ⟩, and σ_ij(t →∞) ≡σ_ij.
By the central limit theorem, the stationary solution for this process is immediately given by

    f_0(y) = 1/√((2π)^N |σ|)exp(-1/2 (σ^-1)_ij y_i y_j  )   .

The perturbation is of the form L_X(y,t) = L_X(y_i) X_i(t).
We can identify the time evolution operator, L_X(y_i), and the generalized potential, ϕ, as

    L_X(y_i) = ∂/∂ y_i ;   ∂ϕ/∂ y_i = (σ^-1)_ij y_j   .

The response of component y_k to an external force acting on component i is then

    R_y_k,X_i(t) = (σ^-1)_ij⟨ y_k(t) y_j(0) ⟩· X_i(t)   .


We derived these results making the following three assumptions, namely that our system is (i) linear and (ii) time invariant and that (iii) the external field only depends on time.
In addition, we only considered correction terms of first order in Eq. <ref> and therefore neglect potential non-linear higher order effects.
The applicability of the central limit theorem to the stationary solution is another consequence of these assumptions.
Note that the fact that many economic time series follow fat-tailed or power law distributions is not at variance with the property that the stationay solution is given by a multivariate normal distribution.
After all, there is no reason to assume that the stationary solution f_0(y) can at any point in time be actually observed in the data due to the ceaseless impacts of direct and indirect shocks; f_0(y) is a computational crutch.




 §.§ Note S6: On the difference between the susceptibility matrix and the Leontief inverse matrix


In this note we clarify the difference between the susceptibilty matrix ρ and the Leontief inverse matrix.
In particular, for a step demand shock in LRT we find a non-equilibrium stationary state given by the linear relation ⟨Δ𝐘⟩_X = ρ 𝐗,
whereas in the “standard” Leontief IO economy we would expect a perturbed equilibrium state given by Δ𝐘 = (𝕀-𝐀)^-1𝐗.
Wherein lies the fundamental difference between these two expressions?

To obtain the LRT solution for the output change, we describe the Leontief IO model in terms of a stochastic differential equation, see Eq. (5) in the main text,

    Y = (𝐀-𝕀)𝐘 + 𝐃 + 𝐗(t) + 𝐅(t)   .

The formal and stationary solution of this equation in the absence of a shock (𝐗(t)=0!) is given in the Supporting Note S5, Eq.(32); it is a probability density function (where we changed variables to have a homogeneous equation, 𝐘→𝐲 = 𝐘 - (𝕀 - 𝐀)^-1𝐃),

    𝐟_0(𝐲) = 1/√((2π)^N |σ|)exp(-1/2 (𝐲σ^-1) 𝐲)   .

The LRT prediction is an expectation value computed using this probability density, which is a multivariate normal distribution centred on the solution of the unperturbed Leontief IO model, (𝕀 - 𝐀)^-1𝐃, when expressed in terms of the variable 𝐘.
For a step demand shock we get Eq. (9) as stationary solution (t →∞),

    ⟨Δ Y_k  ⟩_X = ρ_ki X_i,  with ρ_ki = ∫_0^∞ (σ^-1)_ij⟨ Y_k(τ) Y_j(0) ⟩_0 dτ  ,

where ⟨·⟩_0 means that the expectation value is computed assuming 𝐟_0 as underlying probability density function.
We emphasize again that Eq. (<ref>) shows how a non-equilibrium expectation value (⟨Δ Y_k  ⟩_X for 𝐗≠0) can be defined in terms of a known equilibrium expectation value, namely ρ_ki in Eq. (<ref>).

The computation of output changes in a perturbed Leontief IO model would proceed along a different route.
Let us again consider a step demand shock 𝐗(t) =𝐗θ(t).
We are now interested in the equilibrium state of an economy under such a demand shock.
Therefore, we introduce the perturbed demand, 𝐃_𝐏 as 𝐃_𝐏 = 𝐃+𝐗 and ask for the stationary solution of the stochastic differential equation,
 
    Y = (𝐀-𝕀)𝐘 + 𝐃_𝐏 + 𝐅(t)   .

Formally, Eqs. (<ref>) and  (<ref>) are identical.
However, the perturbed equilibrium perspective assumes a different stationary solution, namely

    𝐟_𝐏(𝐳) = 1/√((2π)^N |σ_𝐏|)exp(-1/2  (𝐳σ_𝐏^-1) 𝐳)   ,

where  𝐳 = 𝐘 - (𝕀 - 𝐀)^-1𝐃_𝐏≠𝐘 - (𝕀 - 𝐀)^-1𝐃 and (σ_P)_ij(t)= ⟨ ( z_i(t) -⟨ z_i ⟩ ) ( z_j(t) - ⟨ z_j ⟩ ) ⟩ with σ_𝐏(t →∞) ≡σ_𝐏.
As could be expected, the stationary state of Eq.(<ref>), i.e. the distribution of values of Δ𝐘 for t →∞, is now a multivariate normal distribution around the perturbed equilibrium state (𝕀-𝐀)^-1𝐃_𝐏.
It is clear that this solution coincides with the LRT solution for 𝐃_𝐏=𝐃, i.e. in the absence of a shock, 𝐗=0.
Also, only for 𝐗=0 would the correct expectation value to compute the susceptibility matrix and output change in LRT be given by the distribution 𝐟_𝐏(𝐳).

Within this note we have now encountered three different types of expectation value, namely (i) the equilibrium expectation value using measure 𝐟_0,  ⟨·⟩_0, (ii) the perturbed equilibrium expectation value given by 𝐟_𝐏,  call it ⟨·⟩_P, and (iii) the nonequilibrium expectation value ⟨·⟩_X.
These three expectation values coincide only in the absence of shocks, 𝐗=0.
From a physical point of view, they describe three different types of system, namely a system relaxing to the state (i) 𝐟_0, (ii) 𝐟_𝐏, or (iii) 𝐟_0 while responding to an external driving force, 𝐗(t) (a physicist would say that the external force 𝐗 “does work on the system”).
In physics, the latter class of systems are closely related to “dissipative structures”, i.e. systems in a steady non-equilibrium state driven by an exchange of energy and/or matter with the environment.

In brief, while both the LRT and perturbed IO approach start from the same stochastic differential equation, Eqs. (<ref>) and  (<ref>), they fundamentally differ in their definitions of expectation values.
In constrast to expected values in the Leontief IO model in the perturbed equilibrium approach, the LRT approach assumes that the stationary solution of the system does  not change after the application of an external shock.
It might also be instructive to compare predictions of the LRT model with those from the perturbed Leontief IO model, i.e. the prediction ⟨Δ Y(t+Δ t) ⟩_P = (𝕀-𝐀)^-1𝐗̃, where 𝐗̃ is the implied shock from Eq. (<ref>).
We have to emphasize that such a comparison is problematic, as the use of the implied shock is only properly defined within the LRT framework, whereas the forecast using the perturbed IO model assumes a different type of dynamical system, as described above.
In this sense, these two models are incompatible and should not be treated on equal footing as we do it here.
Nevertheless, using the same evaluation strategy as we did for comparisons with time series models, we find that the LRT predictions significantly outperform predictions from the perturbed IO model (p<10^-90); see also the Supplementary Figure 10.

 

 


 


 

 


 

 

 


 

 









