LEM
WORKING PAPER SERIESA Method for Agent-Based Models Validation

Mattia Guerini °Alessio Moneta °°Institute of Economics, Scuola Superiore Sant'Anna, Pisa, Italy

2016/16

April 2016
ISSN(ONLINE) 2284-0400

A Method for Agent-Based Models Validation
Mattia Guerini✯ a and Alessio Moneta❸ aaInstitute of Economics, Scuola Superiore Sant’Anna, Pisa, Italy

April 11, 2016

Abstract
This paper proposes a new method for empirically validate simulation models that generate artificial time series data comparable with real-world data. The approach is basedon comparing structures of vector autoregression models which are estimated from both
artificial and real-world data by means of causal search algorithms. This relatively simple
procedure is able to tackle both the problem of confronting theoretical simulation models
with the data and the problem of comparing different models in terms of their empiricalreliability. Moreover the paper provides an application of the validation procedure to theDosi et al. (2015) macro-model.

Keywords: Models validation; Agent-Based models; Causality; Structural Vector Autoregressions.
JEL classification: C32, C52, E37.
Acknowledgments
We are grateful to the authors of the Dosi et al. (2015) model for providing us the artificial datasets. We
also thank Pietro Battiston, Francesco Lamperti, Matteo Richiardi, Matteo Sostero, Pietro Terna and all theparticipants of the VPDE-BRICK Workshop for useful comments and suggestions. We acknowledge funding
from the European Union Horizon 2020 research and innovation programme under grant agreement No. 649186(ISIGrowth) as well as funding from the Institute for New Economic Thinking under grant agreement INO1500021.
✯ Corresponding Author: m.guerini@sssup.it❸ amoneta@sssup.it

Although demonstrative simulation models are useful,not least at performing “what if ” exercises of exploration of different models,
policy analysis requires validated, descriptive simulation models.
Marks (2013)

1 Introduction
Economics, as any scientific discipline intended to inform policy, has inevitably addressed questions related to identification and measurement of causes and effects. This paper, by identifying
and comparing causal structures, proposes a method which improves the empirical reliabilityof policy-oriented simulation models.
The foundation of the Econometric Society in 1930 paved the way of a rigorous and formal
approach to the analysis of causality, which, as Heckman (2000) points out, constituted the
major contribution of econometrics.1 In the post World War II period causal claims wereintroduced in macroeconomics by means of aggregate, mechanic and dynamic models in which
the ex-ante use of economic theory was pivotal. Under this approach the causal process used tobe partitioned in a deterministic component and a random component. The former was meant
to reflect the causal relations dictated by economic theory. The condition for it to be considered“valid” was to have the random component satisfying the standard Gauss-Markov statisticalproperties. Such a methodology goes under the name of Cowles Commission or Simultaneous
Equations Model (SEM) approach. The most prominent proposers were Haavelmo (1944) andKoopmans (1950).
This approach has been strongly criticized by Lucas (1976) and Sims (1980) on theoretical and methodological grounds respectively: the former insisted on the fact that individuals
endowed with rational expectations would have anticipated the policy interventions supported
by SEMs and their behaviour would have brought to opposite results with respect to the ones
predicted by SEMs; the latter instead stressed the fact that in the Cowles Commission approach the distinction between endogenous and exogenous variable was stated ad-hoc, in orderto ensure system identifiability.
Taking as starting points the Lucas (1976) and Sims (1980) critiques, Kydland and Prescott
1

As Hoover (2004) has shown, however, causal language has not always been explicit in economics and in thesciences in general. In the first half of the twentieth century, under the influence of Karl Pearson, Ernst Machand Bertrand Russell, many research scientists endeavoured to eschew causal concepts in order to privilege talkof functional and statistical dependencies (Illari et al., 2011). Causal talk revived in the second half of the lastcentury (Hoover, 2004). See also Granger (1980)

2

(1982) paved the way toward a new class of models, becoming the founding fathers of thestream of literature that goes under the name of Real Business Cycle (RBC) theory and whichthen evolved in what today is known as the Dynamic Stochastic General Equilibrium (DSGE)
approach. These types of models are nowadays the most widely used to draw and to evaluate
policy claims because they bear the advantage of simultaneously addressing two critical issues
about causal structures. On the one hand, under the acceptance of the rational expectation
hypothesis, the structure modelled by the RBC/DSGE approach remains invariant under policy intervention because it takes into account the forward-looking behaviour of the economic
agents. On the other hand, the theoretical structure has an empirical counterpart in which thedistinction between endogenous and exogenous variables is eschewed. The empirical counterpart is represented by a Structural Vector Autoregressive (SVAR) model.2 But the RBC/DSGE
approach is not exempt from problematic issues: structural stability is grounded on individual
behaviour, but assuming a representative agent which neglects or even denies any form of interaction. Moreover, the identification of the empirical structure in the SVAR model is typicallyachieved by imposing restrictions derived from the theoretical model, which are therefore notsubjected to any severe test. (See Fagiolo and Roventini (2012) for a detailed criticism on
similar issues).
An alternative approach to the problem of representing macroeconomic causal structures, inwhich it is possible to run reliable policy experiments, is to build a class of models which betterreflect the existing economic mechanisms, including the microeconomic interactions. This is theaim of the Agent-Based Model (ABM) approach, also known as the Agent-Based Computational
Economics (ACE) approach, in which the macroeconomic structure is analysed as an emerging
property from the interaction between heterogeneous and bounded rational economic actors.
This modelling strategy has been applied to economic theory only during the last three decades,
but it rapidly gained a significant success and in the very last years has begun to be perceived
as a new valuable paradigm, able to provide a viable alternative to the DSGE framework.3
ABM is a useful and flexible tool for performing rich policy experiments and for evaluating
their implications; among the main advantages of the ABM strategy, stands the possibilityof analysing endogenously generated booms and busts and also to study the reaction of the2See, however, Canova and Sala (2009) and Fukac and Pagan (2006) for cautionary notes about the existenceof the empirical counterpart of a DSGE model.
3
The fast success of ABM might be due both to the huge improvements in computational power and recently,
to their ability of explaining and revealing intriguing aspects of the world financial and economic crisis; DSGE
models instead, has proven to be of few help when facing the crisis. See Trichet (2010) for clarification on thislast point.

3

economy to different stimuli, also applied in periods of distress and not only around a fictitious
locally stable steady state of the economy.
But ABMs pose a serious methodological problem because of their unclear relationshipwith the empirical evidence. This paper aims to address this issue. The difficulties of the ABM
approach, which represent the counterpart of its flexibility, are perceived both in the modeldata confrontation and in the comparison of different models investigating the same piece ofevidence. The quality of ABMs has been up to now evaluated according to the ex-post ability
in reproducing a number of stylized facts even if other validation procedures are available (seeFagiolo et al., 2007). We argue that such an evaluation strategy is not rigorous enough. Indeed
the reproduction, no matter how robust, of a set of statistical properties of the data by a model
is a quite weak form of validation, since, in general, given a set of statistical dependencies
there are possibly many causal structures which may have generated them. Thus models which
incorporate different causal structures, on which diverse and even opposite practical policy
suggestions can be grounded, may well replicate the same empirical facts.4
The present work proposes a procedure to validate a simulation model which proceeds by
first estimating both the causal structure incorporated in the model (using the data artificiallygenerated by the model) and the causal structure underlying the real-world data. Secondly, it
compares the two inferred causal structures. In this manner the proposed procedure offers asolution to both the issue of confronting an ABM to empirical data and the issue of confronting
different simulation models. Indeed causal structures inferred from different simulation data,
generated by different models can be compared in the same way. A good matching between
the causal structure incorporated in the ABM and the causal structure underlying the realworld data provides a more rigorous empirical support to the policy statements drawn from theABM, if compared with the support coming from mere replication of statistical evidence. Other
validation procedures have been recently proposed based on information criteria by Lamperti(2015) and Barde (2015); other researchers such as Grazzini (2012), Lux (2012), Recchioni et al.
(2015), Gilli and Winker (2003) have focused on estimation, or on the analysis of the emergent
properties stemming from ABMs (see Grazzini and Richiardi, 2015); there has also been interestfor parameter space exploration (see Salle and Yildizoglu, 2014). The flourishing of all thesecomplementary approaches devoted to the solution of such interrelated issues can be seen an
indicator of their relevance and a signal of the vitality of the agent-based community.
The paper is organised as follows. Section 2 reviews the different strands of literature on
4

At the root of this underdetermination problem is the fact that statistical relationships are in generalsymmetric, while this is not necessarily the case for causal relationships.

4which our method is built upon; the validation algorithm is presented extensively in Section 3;
Section 4 provides a first application of the method to the “Schumpeter meeting Keynes” model
proposed by Dosi et al. (2015). Section 5 concludes.

2 Background literature
DSGE models are confronted to the data in two ways. The first and traditional approachis through calibration, in which the parameters of the model are chosen from pre-existing
microeconomic studies or in order to replicate the statistical properties of aggregate variables(Kydland and Prescott, 1996). The second approach is through the estimation of a VAR model
built to represent the empirical counterpart of the DSGE model. Having estimated a VAR, one
can identify a SVAR and confront its impulse response functions with the responses to policy
shocks derived from the DSGE model. Alternatively, as proposed by Ireland (2004), one canaugment the DSGE model with the VAR residuals and estimate a hybrid model via maximum
likelihood (for a criticism see Juselius, 2011).
Calibration and replication of statistical properties of data is practised in the ACE community as well. To our knowledge, however, the models by Bianchi et al. (2007) and Bianchi
et al. (2008) are the unique medium-scale agent-based macro-models in which parameters are
estimated ex-ante and calibrated ex-post in order to replicate statistical properties of observed
data.5
Although calibration and replication of statistical properties of data make some step intaking the model to the data, we claim that this is not enough for the reliability of the policy
implications derived from the model, since two models with alternative policy implications may
well be both calibrated in order to replicate certain statistical properties of observed data.
Reliability can be improved only through a validation exercise designed to provide evidence
that the modelled data generating mechanism is an adequate representation of the real-data
generating mechanism.
We are aware that the economic system has the characteristics of a complex system inwhich somehow stable macroscopic aggregate properties emerge from intricate connections at
the microscopic level. But we further believe that representing into a unique model every single
micro-mechanism at work in a complex economy and showing its good match with data at
different levels of aggregation is a very difficult task. A reduction in the complexity of the5

There are, however, several small-scale ACE financial models which are instead calibrated or estimated suchas the ones by Alfarano et al. (2005, 2006, 2007).

5

issue may be necessary, and hence in what follows we will analyse only the relations between
macro variables.6 Our strategy is indeed to focus only on representing causal structures among
aggregate variables of the ABM and test whether they significantly differ from the causal
structures that can be found in the real world from observed aggregate variables, without
further considerations of the micro properties. In other words, we compare a macro-reduced
version of the modelled generated mechanism with a macro-reduced version of the real-data
generating mechanism.7
Our procedure will separately identify the causal structures of the two different data generating processes at their aggregate level, and then will compare the results of the estimations: if thecausal structures are similar, then the model is a good characterization of the causal structureof the real-world data generating process and we will consider it as “valid”. The identification
method is the same for both processes: we will estimate a SVAR model using both observed
and simulated aggregate data. This model, being a model with well known properties, provides
us enough power and flexibility to compare the explanatory performances of ABM with that ofreal-world data. But a crucial feature in the SVAR estimation is the identification procedure,
which we describe in the next subsections.

2.1 SVAR identification: an open issueStarting from a multiple time series dataset composed of K variables collected for T periods
we can denote by Yt = (Y1t , . . . , YKt )′ the values of these variables at a particular time t. A
simple, but useful way of representing the data generating process, is to model the value of eachvariable Ykt as a linear combination of the previous values of all the variables as well as theircontemporaneous values:
Yt = BYt + Γ1 Yt−1 + · · · + Γp Yt−p + εt(1)where the diagonal elements of the matrix B are set equal to zero by definition and where εt
represents a vector of error terms which we will assume to be mutually statistically independent.
Therefore the covariance matrix Σε = E [εt ε′t ] is diagonal. The SVAR form of this model can
6

Cfr. Haldane (2012).
Possible developments of the same method may allow to compare a micro-macro version of the modelledgenerated mechanism with the real-data generating mechanism, using observations at different levels of aggregation.
7

6

also be written as
Γ0 Yt = Γ1 Yt−1 + · · · + Γp Yt−p + εt(2)where Γ0 = I − B. The problem with equations (1) and (2) is that they cannot be directly
estimated without biases, being the contemporaneous variables endogenous. What is typicallydone in the literature is to estimate the reduced form VAR model:−1−1
Yt = Γ−1
0 Γ1 Yt−1 + · · · + Γ0 Γp Yt−p + Γ0 εt(3)= A1 Yt−1 + · · · + Ap Yt−p + ut .
′where ut = Γ−10 εt is a zero-mean white noise process with a covariance matrix Σu = E [ut ut ]that in general is not diagonal.
The problem is that, even if the parameters contained into Ai , for i = 1, . . . , p can be
estimated from equation (3) without incurring in any particular issue, their knowledge is not
sufficient for the recovery of the structural parameter contained in B and in Γi , for i = 1, . . . , pof equation (1), making impossible the inference of any causal and/or policy claim. To do suchclaims we need to recover the matrix Γ0 that contains the contemporaneous causal effects. But
the problem is that any invertible unit-diagonal matrix might be compatible with the coefficients
estimated from the VAR in equation (3).
The problem of finding the appropriate Γ0 (and hence also finding the matrices Γ1 , . . . , Γp )is called the identification problem and it is usually performed by imposing restrictions on the Γ0matrix using a Cholesky factorization of the estimated covariance matrix Σu . But this approachshould only be employed when the recursive ordering implied by the identification scheme isfirmly supported by theoretical consideration. A class of alternative identification procedures
derives from the seminal papers by Bernanke (1986) and Blanchard (1989) and imposes zero
restrictions which are based on economic considerations about contemporaneous interactions.
An alternative identification strategy descends from Shapiro and Watson (1988) and Blanchard
and Quah (1989) by assuming that certain economic shocks (e.g. supply shocks) have long-run
effects on some variables but do not influence in the long-run the level of other variables, whileother shocks (e.g. demand shocks) have only short-run effects on all the variables. Unfortunately
these identification strategies are grounded on some level of theoretical apriorism which does
not completely solves the critique put forward by Sims (1980).
A relatively recent approach for solving the identification issue of a SVAR model in a more
agnostic and data-driven fashion, allowing to avoid as much as possible subjective choices and7theory driven considerations, has been put forward by Swanson and Granger (1997), Bessler
and Lee (2002), Demiralp and Hoover (2003), Moneta (2008) and Moneta et al. (2011) and isbased on graphical causal models (see Pearl, 2000; Spirtes et al., 2000).

2.2 Graphical causal models and SVAR identification
A Causal Graph G is a model that consists of a set V of vertices (nodes) and a set E of edges(links) and might be written concisely as G = hV, Ei. It is aimed at representing and analysing
specific features of the data-generating process underlying the set of observed variables. The
vertices of such a graph correspond to random variables and the edges denote causal relationships among them. In what follows we focus on the simple case of Directed Acyclic Graphs(DAG) in which all the edges are directed and causal loops are not allowed.
The identification procedure based on graphical causal models consists of three steps: (i)estimating the reduced form VAR of equation (3), (ii) applying a search algorithm to theestimated residuals ut to obtain the matrix Γ0 (cfr. equation 2) and (iii) recovering the othermatrices Γi (i = 1, . . . , p) of the SVAR model.
The critical part of the procedure is the second step, in which an algorithm is applied inorder to uncover the causal dependencies among the residuals ut . The literature on causal
search models has developed a plethora of algorithms which differ among each other for theassumptions on which they are based and the computational properties. Assumptions typicallyconcern the form of the causal structure (e.g. cyclic or acyclic), the presence or exclusion oflatent variables (i.e. causally not sufficient or causal sufficient structures), rules of inference(more on that below), and, finally, statistical properties of the residuals (e.g. normality or
linearity) which allow the application of specific tests of conditional independence.
The algorithm presented in Appendix A is the PC algorithm originally developed by Spirtes
et al. (2000). In this algorithm causal loops are not allowed. Indeed it is assumed that thecausal generating mechanism can be modelled by a DAG. In the SVAR framework this amounts
to exclude feedbacks in the contemporaneous causal structure, while feedbacks over time areof course conceivable (e.g. Xt causes Yt+1 which in turn causes Xt+2 ). The PC algorithm also
assumes causal sufficiency, i.e. there is no unmeasured variable which simultaneously affecttwo or more observed variables. Rules of inference are conditions which permit deriving causal
relationships starting from tests of conditional independence. The PC, and similar algorithms

8

of the same class, hinge on two rules of inference (see Spirtes et al., 2000):
Condition 1. (Causal Markov Condition) Any variable in the causal graph G is conditionallyindependent of its graphical nondescendants (i.e. non-effects) — except its graphical parents —given its graphical parents (i.e. direct causes).
Condition 2. (Faithfulness Condition) Let G be a causal graph and P be a probability
distribution associated with the vertices of G. Then every conditional independence relation
true in P is entailed by the Causal Markov Condition applied to G.
The PC algorithm, as many other of the class of constraint-based search algorithms, needs
as input knowledge of the conditional independence relationships among the variables. There
are many possibilities of testing conditional independence which in principle are all compatible
with the PC algorithm. If the probability distribution underlying the data is Gaussian, zeropartial correlation implies conditional independence. Then a typical procedure is to test for
Gaussianity and in case this is not rejected, one can test for zero partial correlations. In manystatistical packages the default option is to test zero partial correlation through the Fisherz-transformation, as proposed by (Spirtes et al., 2000). An alternative option, suited for the
SVAR framework, is to test zero partial correlations among the VAR residuals through a Wald
test that exploits the asymptotic normality of the covariance matrix of the maximum-likelihood
estimated VAR residuals (for details see Moneta, 2008). If Gaussianity is rejected or one is notwilling to make distributional assumptions, one way to proceed is to rely on nonparametric testsof conditional independence, which, however, present the well-known problem of dimensionality(cfr. Chlass and Moneta, 2010).
The PC algorithms follows this scheme:
i. Create a complete graph on the variables (X1 , . . . , Xk );ii. Apply tests for conditional independence in order to prune unnecessary edges;
iii. Apply tests for conditional independence in order to direct remaining edges.
There are other algorithms in the literature that, following a similar scheme, allows for
feedback loops or the possibility of latent variables (e.g. the CCD or FCI algorithm; cfr. Spirtes
et al., 2000).
Usually, the output of the causal search is in general not a unique graph G, but a set ofMarkov equivalent graphs which represent all the possible data generating processes consistent9with the underlying probability P. Hence the information obtained from this approach is
generally not sufficient to provide full identification of the SVAR model requiring again a certain
level of a priori theoretical knowledge. Moreover, if the distribution of the residuals is nonGaussian, it is necessary to apply tests of conditional independence which are different from testsof zero partial correlation. However, Moneta et al. (2013) have shown that if the VAR residuals
are non-Gaussian, one can exploit higher-order statistics of the data and apply Independent
Component Analysis (ICA) (see Comon, 1994; Hyvarinen et al., 2001) in order to fully identify
the SVAR model.

2.3 Independent component analysis and SVAR identification
Let us recall the fact that VAR disturbances ut and structural shocks εt are connected via
ut = Γ−10 εt .

(4)

In this framework the VAR residuals are interpreted as generated by a linear combinationof non-Gaussian and independent structural shocks via the mixing matrix Γ−1
0 . Independent
Component Analysis applied to equation (4) allows the estimation of the mixing matrix Γ−1
0 and
the independent components εt by finding linear combinations of ut whose mutual statisticaldependence is, according to some given measure, minimized. Some points should be noticed: (i)while the assumptions of mutual independence of the structural shocks is usually not necessary
in a SVAR framework (orthogonality is usually sufficient), such an assumption is necessary
to apply ICA; (ii) ICA does not require any specific distribution of the residuals ut but only
requires that they are non-Gaussian (with the possibility of at maximum one Gaussian element);(iii) the ICA-based approach for causal search does not require the faithfulness condition; (iv)in non-Gaussian settings while conditional independence implies zero partial correlation, theconverse does not hold in general.
The application of ICA to the estimated VAR residuals allows identifying the rows of thematrix Γ0 , but not their order, sign and scale (for details see Hyvarinen et al., 2001). In orderto obtain the correct matrix Γ0 , that is the matrix incorporating the contemporaneous causal
structure and such that Γ0 = I − B in equation (1), we further assume that the VAR residuals
can be represented as a Linear Non-Gaussian Acyclic Model (LiNGAM) so that the contemporaneous causal structure can be represented as a DAG. On the basis of this assumption, it
is possible to apply the causal search algorithm presented in Appendix B (VAR-LiNGAM),
10

which draws on the original contributions of Shimizu et al. (2006) and Hyvarinen et al. (2010)(for an application to economics see Moneta et al., 2013). The basic idea by which the VARLiNGAM algorithm solves the order indeterminacy is that if the underlying causal structureis acyclic, there must be only one row-permutation of the ICA-estimated rows of Γ0 such that
all the entries of the main diagonal are different from zero. Hence, the algorithm applies asearch procedure to find such a permutation (for details see Appendix B, step C). The scaling
indeterminacy is solved by normalizing the elements of the ICA-estimated matrix and rightly
row-permuted Γ0 , such that the main diagonal is one (and the main diagonal of B is zero, as
dictated by equation (1).
In the literature there are alternative ICA-based search algorithms, which relax the assumption of acyclicity and causal sufficiency: see for example the algorithms proposed by Lacerda
et al. (2008) and Hoyer et al. (2008), which respectively allow for feedback loops and for latentvariables. However, since this is a first application of this type of framework to validation ofsimulated model, we decided to keep the analysis as simple as possible so that in future works
we might relax assumptions, and understand which are the most critical for validation concerns.

3 The validation method
In this section we describe our validation procedure which is composed of five different steps
as shown in figure (1). In the first step we apply some simple transformations that allow theempirical and the artificial data to be directly comparable; in the second step we analyse theemergent properties of the series produced by the simulated model; in the third step we estimate
the reduced-form VAR model; in the fourth step we identify the structural form of the modelby means of some causal search algorithm; in the last step we compare the two estimated causal
structures according to some distance measure.

3.1 Dataset uniformity
Our method starts by selecting, in the model under validation inquiry, K variables of interest(v1 , . . . , vK ). We then collect a dataset that corresponds to the actual realization of thesevariables in the real world (we call this dataset RW-data) as well as a dataset for the realizationsof M Monte Carlo simulations of the agent-based model (we call this one AB-data). We thus
obtain two preliminary datasets VRW and VAB which might be of different dimensions. In
11

Dataset Uniformity

Analysis of ABM PropertiesVAR Estimation

SVAR Identification

Validation Assessment

Figure 1: The five steps of the validation method.

general we will have
dim(VRW ) = 1 × K × TRW
dim(VAB )= M × K × TAB

meaning that for the real world we observe only 1 realization of the the K variables of interest
for a period of length TRW while for the simulated data (for which we can possibly have an
infinity of observations) we will have M Monte Carlo realizations, of the same K variables, fora period of length TAB ; it often holds true that TAB >> TRW , so that the two datasets are not
perfectly matchable.
The large availability of realizations in the simulated data is in fact an advantage and not
an issue, since this allows a pairwise comparison of each run of the Monte Carlo simulation withthe unique empirical realization. But the presence of different lengths in the time series might
generate issues in two directions: (i) using the whole length of the AB-data time series creates
the risk of capturing the effects present in the transient period, which does not represent thetrue dynamic entailed by the model, but is only due to the choice of the initial conditions, (ii)lag selection might be affected due to unobserved persistence in some of the modelled variables.
Therefore we remove an initial subset of length TAB −TRW from each of the M artificial datasets

12(as shown in figure 2) in order to force each pair of datasets to have the same dimensions:dim(VRW ) = dim(VAB (i)) = 1 × K × TRW

1

for i = 1, . . . , M.

TAB − TRW

TAB

1

TRW

Figure 2: Time window selection. The first periods of the AB-data are canceled with theobjective of homogenising the series.

Moreover the order of magnitude of RW-data and AB-data are typically different; this is not
perceived by the ABM community as an issue, being the concern of a large number of ABMs
the replication of stylized facts (distributions, variations, statistical properties but not levels).
But in our approach this might create comparability issue. We will see that in our application
it is sufficient to take a logarithmic transformation in order to smooth out this scaling issue,
and we speculate that in many applications any monotonic transformation might be applied.

3.2 Analysis of ABM properties
Some considerations about two underlying assumptions are needed. For the model to be agood proxy of the data generating process, we require that it is in a statistical equilibrium
state in which the properties of the analysed series are constant. In particular we require thatthe series, or a transformation of them (e.g. first differences), have distributional properties
that are time-independent; secondly we require that the series are ergodic, meaning that theobserved time series are a random sample of a multivariate stochastic process.
In a context where the series have been generated by a simulation model, which provides
M Monte Carlo realizations, these two assumptions can be tested directly (see Grazzini, 2012).
Indeed if we consider all the M time series realizations of a variable of interest k we will collect ammatrix with dimensions M ×T containing all the generated data Yk,t
, as represented in figure (3).

We call here ensembles the column vectors of such a matrix; therefore each ensemble contains(·)

the M observations Yk,t in which the time dimension is fixed. We instead define samples all themrow vectors of the same matrix, each of which contains the T observations Yk,(·)
in which the13Monte Carlo dimension is fixed. Let’s denote by Ft (Yk ) the empirical cumulative distributionfunction of an ensemble and by Fm (Yk ) the empirical cumulative distribution function of asample. Testing for statistical equilibrium and for ergodicity reduces to test (via KolmogorovSmirnov test) for the following conditions:
Fi (Yk ) = Fj (Yk ), for i, j = 1, . . . , T
Fi (Yk ) = Fj (Yk ), for i = 1, . . . , T

i 6= j(5)

j = 1, . . . , M(6)Therefore we perform two tests as represented in figure (3): we recursively run tests of pairwise
equality of distributions and we present the percentage of non-rejection of such tests. Rejecting
the test would imply that the distribution under investigation are different from each other.
Our two assumptions will result to be supported by the data if we obtain high percentages of
non-rejection.

y1,1 y1,2 . . . y1,T
 y2,1 . . . . . . y2,T

k
Y =
..
..
..
..

.
.
.
.

yM,1 yM,2 . . . yM,T

y1,1 y1,2 . . . y1,T
 y2,1 . . . . . . y2,T

k
Y =
..
..
..
..

.
.
.
.

yM,1 yM,2 . . . yM,T
















Figure 3: The elements of comparison when testing for statistical equilibrium (left) and for
ergodicity (right).

3.3 VAR estimation
Following the Box and Jenkins (1970) methodology, the first task when any time series modelhas to be estimated, is the lag selection: in our case this means choosing the two values pRW
and piAB , for i = 1, . . . , M according to some information criterion like the BIC the HQC or
the AIC. Two cases might anyway emerge from the data:
1. pRW − piAB = 0 which means that our estimations based on the artificial dataset and on
the real-world dataset are perfectly comparable;
2. pRW − piAB 6= 0 which means that one of the two dataset presents at least one effect which14

is not present in the other; we keep this fact into account when computing the similarity
measure.
Once the lag selection has been performed, our procedure estimates the VAR as explicitedin equation (3) via OLS and also in his VECM form, via maximum likelihood estimation (seeLutkepohl, 1993) using the Johansen and Juselius (1990) procedure.

3.4 SVAR identification
In this step we extract the vectors of residuals (u1 , . . . , uK ) from the estimation of the VAR
and analyse their statistical properties and their distributions. We test for normality applying
the Shapiro-Wilk and the Jarque-Bera statistics. Then according to the outcome of the tests,
we select the appropriate causal search algorithm to be adopted for the identification strategy.
Two algorithms, the PC (to be adopted for the Gaussian case) and the VAR-LiNGAM (forthe non-Gaussian case) are presented extensively in the appendices A and B. At the end of theidentification procedure, we have estimated our structural matrices ΓRW
for i = 0, . . . , pRWi, for i = 0, . . . , pAB and for m = 1, . . . , M .
and ΓAB,mi

3.5 Validation assessment
The last step consists in the comparison of the causal effects entailed by the SVAR RW andthe SVAR AB models. This will tell us how many of the real-world estimated causal effects arecaptured also by the agent-based model under validation inquiry.
In order to compare the causal effects we will use the similarity measure Ω, which weconstruct, as already anticipated, starting from the estimates of the SVAR.
RW
AB
the (j, k) element
Let’s denote γi,jk
the (j, k) element of ΓRW
for i = 0, . . . , pRW and γi,jkiof ΓAB,mfor i = 0, . . . , pAB and for m = 1, . . . , M . We define pmax = max {pRW , pAB } andithen we setΓRW = 0
i


ΓAB
=0
i

for pRW < i ≤ pmax

if

pRW < pmax = pAB

for pAB < i ≤ pmaxif

pAB < pmax = pRWThis allows us to penalize the value obtained by the similarity measure for the fact that thecausal effects are completely mismatched after a certain lag. Then we build the indicator15

function:
ωi,jk =


1 if sign(γ RW ) = sign(γ AB )
i,jk
i,jk

(7)


RW
AB
0 if sign(γi,jk) 6= sign(γi,jk)

where i = 0, . . . , pmax is the index for the ith -order matrix containing the causal effects, while j
and k are the row and column indexes of these matrices. The similarity measure is then defined
as
Ω=

P

pmax PK PK
i=1
j=1
k=1 ωi,jk

K 2 pmax.

(8)

Our similarity measure is bounded between [0, 1] allowing us to have an easily interpretableindex that represents the ability of the agent-based model to recover, at the aggregate macrolevel, the same causal relationsips estimated in the RW-data.

4 Application to the “Schumpeter Meeting Keynes” model
The idea of building a simulation macroeconomic laboratory performing policy exercises dates
back to Lucas (1976) and Kydland and Prescott (1982) but one problem of this approach has
always been the external validity. Our methodology can be interpreted as a test for external
validity for any simulation model and it might be of particular interest for researchers engaged
in practical policy matters and policy makers who should take decisions according to models
which must be reliable and valid. We already argued that a policy reliable model is one thatis able not only to replicate a list of stylized facts, but also to represent the real-world causal
structure as much accurately as possible. We want to test here the validity of the agent-basedlaboratory by Dosi et al. (2015), a model that builds upon a series of previous papers (see Dosi
et al., 2013, 2010) and that has attracted great attention in the recent literature.
The Dosi et al. (2015) model aims at investigating the implications of demand and supply
public policies in a model that “bridges Keynesian theories of demand-generation and Schumpeterian theories of technology-fuelled economic growth”. The model by itself is able to reproduce
a long list of stylized facts and in particular is able to reproduce a cross correlation table close
to the one usually computed with the US observed data.

164.1 The dataset
Since the model under validation inquiry is dedicated to the analysis of the real side of an
economic system and to the analysis of fiscal and monetary policies, the K = 6 variables ofmajor interest that we will consider in our system of equation are: aggregate consumption (C),gross private investments (I), unemployment rate (U ), gross domestic product (Y ), current
price index (P ) and effective federal funds rate (R). The RW-data refer to the United States
and are collected from the Federal Reserve Economic Database (FRED); we decided to cover
the period from January 1959 to April 2014 with a quarterly frequency, implying a time series
length T = 222; this is a typical selection when analysing US business cycle. All the variables
are plotted in figure 4.

Real World Time Series

Typical Agent−Based Time Series

C

I

C

I

8
9

24

28
7

6

22
26

I

C

I

C

8

7
20
5

24

6
18
0

50

11

100

150

200

0

50

100

150

Time

Time

U

Y

200

0

50

100

150

200

0

50

100

Time

Time

U

Y

150

200

150

200

150

200

50

40

9

9

28

8

Y

7

U

Y

U

30
26

20
5

7

10

24

0
0

50

100

150

200

0

50

100

150

Time

Time

p

r

200

0

50

100

150

200

0

50

100

Time

Time

p

r

5.5
11

15
5.0

10

r

p

r

p

10
10

4.5

9
4.0

5

5
8

3.5
0

0
0

50

100

Time

150

200

0

50

100

150

200

0

Time

50

100

Time

150

200

0

50

100

Time

Figure 4: Left columns: time series of RW-data. Right columns: time series of a typicalAB-data.

To fulfil the dataset uniformity requirement for the AB-data, we collect the last T time observations, getting rid of possible transients and we consider M = 100 Monte Carlo simulations;17

each of them will be pairwise compared to the unique realization of the RW-data. Finally wetake logs of the C, I, Y, P variables and we uniform U and R by expressing them in percentage
terms.8
We then check whether the assumptions we require for estimation of our AB model are too
stringent or they are supported by the data; we perform the statistical equilibrium and theergodicity tests as described in the section 3.2 and in figure (3); this kind of test is in line alsowith the analysis proposed by Grazzini (2012).
C
I
U
Y
P
R

Equilibrium
0.95380.96340.96080.95320.95600.9609Ergodicity
0.94790.95640.9396
0.95130.90550.9716

Table 1: Percentages of non-rejection of statistical equilibrium and ergodicity.
Table (1) presents the percentage of non-rejection of the Kolmogorov-Smirnov test of eachpairwise comparison, for each stationarised series. The results come from T ×(T2 −1) = 24310 and
T × M = 22100 pairwise comparisons for the statistical equilibrium and for the ergodicity tests
respectively. For all the series we have values higher than 90% and this allows us to conclude
that the assumptions about the model having reached a statistical equilibrium and producing
ergodic series are reasonable.

4.2 Estimation and validation results
The augmented Dickey-Fuller test does not reject the null hypotheses of unit root in all thereal-world time series. For AB-data, the evidence for ubiquity of unit root is weaker since for
I, U and R we can reject at the 5% level the presence of unit root (see table 2). This does notcreate any difficulty to our causal search procedure and it is only a stylized fact not replicatedby the model.
We then estimate the model as a vector error correction model (VECM) with cointegrating
relationships, without taking first difference of any variable, using the Johansen and Juselius(1990) procedure, which is based on a maximum-likelihood estimation with normal errors, but isrobust also to non-Gaussian disturbances. For sake of completeness we also check the robustness
8

Three Monte Carlo simulations, m = {55, 61, 89}, are excluded from the original dataset because in (atleast) one period of the model, investment goes to 0, implying we cannot take logarithm.

18

(a) RW-data

Variable
C
I
U
Y
p
r

ADF p-value for levels
0.99
0.89
0.06
0.99
0.92
0.22

ADF p-value for 1st -differences
0.02
0.01
0.01
0.01
0.19
0.01

Critical level
0.05
0.05
0.05
0.05
0.05
0.05

(b) AB-data

Variable
C
I
U
Y
p
r

ADF p-value for levels
0.43
0.01
0.01
0.19
0.63
0.01

ADF p-value for 1st -differences
0.01
0.01
0.01
0.01
0.01
0.01

Critical level
0.05
0.05
0.05
0.05
0.05
0.05

Table 2: Augmented Dickey-Fuller Test.
of the results by estimating the VAR in level via OLS. We select the number of lags according to
the Bayes-Schwarz Information Criterion (BIC) and the number of cointegrating relationships
following the Johansen procedure. For the real-world dataset the suggestion is that of using 3
lags and 2 cointegrating relationships while in the artificial datasets (typically) we should use
3 lags and 3 cointegrating relationships. This implies that we do not have any comparability
issue for what concerns the number of lags. With respect to the cointegrating relations, it isagain a stylized fact not matched by the model which does not create any estimation issue,
since we are interested in structural form of the model and cointegration is only a reduced form
property.
The empirical distributions of the VAR residuals (uC,t , . . . , uR,t ) are represented in figure (5)both for RW-data and for a typical Monte Carlo realization of the AB-data simulation; moreovertable (3) collects the results of the Shapiro-Wilk and the Jarque-Bera tests for normality; for all
the variables, the residuals from the real-world data and all but one residual (the unemployment)from the artificial data, the tests rejects the null hypothesis of normality.
We conclude that the residuals ut are non-Gaussian and this result leads us toward theidentification of the SVAR via the LiNGAM algorithm.9
After having completed the estimation, we compute the similarity measure as defined in
9

The VAR-LiNGAM algorithm is consistent even if one variable is normally distributed and therefore eventhe unemployment residuals quasi-Gaussianity does not add complications to our identification and structural
estimation procedure.

19

RW−data VECM Residuals Density vs Normal Density
C

Typical AB−data VECM residuals Density vs Normal Density

I

C

15

80

I

12.5
0.9

10.0
7.5

density40

density10

densitydensity

60

5.0

0.6

5
0.3

20

2.5

0

0
−0.03

−0.02

−0.01

0.00

0.01

0.0
−0.10

−0.05

0.00

C

0.05

0.10

0.0

−0.10

−0.05

0.00

I

U

Y

2.0

0.05

0.10

−1

0

C

1

I

U

Y
8

60
0.075

1.5

6

20

0
−0.5

0.0

0.5

0.00

p

0.01

0.02

0.03

4

2

0.000
−0.02 −0.01

U

125

0.050

0.025

0.5

0.0

densitydensitydensitydensity

40
1.0

0
−15

−10

−5

0

Y

U

r

p

5

10

−0.1

0.0

0.1

0.2

Y

r

0.8
0.3

100

densitydensitydensity

50

0.4

0

0.0
−0.03

−0.02

−0.01

0.00

0.01

20

0
−2

p

0

2

4

0.0
−0.04

r

0.2

0.1

10

0.2

25

density

30

0.6
75

−0.02

0.00

p

0.02

−2.5

0.0

2.5

r

Figure 5: Left columns: RW-data VECM residuals distribution (green) and normal distribution (blue). Right columns: typical AB-data VECM residuals distribution (green) andnormal distribution (blue).

equation (8). The results suggest that when we estimate the system using a OLS-VAR strategy
the Schumpeter meeting Keynes model is able to reproduce, on a Monte Carlo average, the78.92% of the causal relations entailed in the real-world dataset (the similarity drops to 64.9%
after accounting only for bootstrapped significant parameters); on the other side, if in the firststep we estimate a VECM by means of maximum likelihood, the similarity measure marks
73.85% (raising to 79.89% when considering only bootstrapped significant parameters). The
results are reported also in table (4), containing not only the means but also standard deviations
across Monte Carlo. Being the dispersion index quite low, we can conclude that neither very
negative nor very positive outliers are present. Therefore a large fraction of simulations entail
the same bulk of causal relations.

20

(a) RW-data
C
I
U
Y
P
R

Shapiro-Wilk Test
0.95
0.96
0.98
0.97
0.88
0.89

Shapiro-Wilk p-value
0.00
0.00
0.02
0.00
0.00
0.00

C
I
U
Y
P
R

Shapiro-Wilk Test
0.99
0.97
0.99
0.98
0.98
0.99

Shapiro-Wilk p-value
0.03
0.00
0.08
0.02
0.00
0.05

Jarque-Bera Test
271.48
54.96
12.16
61.16
1710.58698.46Jarque-Bera p-value
0.00
0.00
0.00
0.00
0.00
0.00

(b) AB-data
Jarque-Bera Test
6.01
29.64
3.23
4.72
10.98
7.30

Jarque-Bera p-value
0.05
0.00
0.20
0.09
0.00
0.03

Table 3: Normality test on the VECM residuals.
Estimation Method
VAR-OLS (all parameters)
VECM-ML (all parameters)
VAR-OLS (significant parameters)
VECM-ML (significant parameters)µ
0.78920.73850.64900.7989σ
0.05170.06280.10300.0689

Table 4: Mean and standard deviation of the similarity measure.

5 Conclusions
In this paper we have presented a new method for validating policy-oriented Agent-Basedmacroeconomic models able to generate artificial time series comparable with the aggregate
time series computed by statistical offices, central banks and institutional organizations. The
approach is based on comparing Structural Vector Autoregressive models which are estimatedfrom both artificial and real-world data by means of causal search algorithms. In the paper wealso have presented a first application of our method to the Dosi et al. (2015) model. We havecalculated that by using the simulated data and according to the proposed similarity measure,the model is able to resemble between 65% and 80% of the causal relations entailed by a SVAR
estimated on real-world data. We guess that this is a positive result for the Schumpeter meeting
Keynes model but in order to reinforce this claim, we should need to compare this results withthose coming from other models. In our opinion, this paper sets a new benchmark on which
members of the agent-based community might build upon. Convinced about the fact that thevalidation issue cannot be settled out in an ultimate manner, other approaches for model validity
can emerge and might bring complementary evidence to ours. Indeed a possible strategy, forresearchers willing to bring their agent-based models to the audience of policymakers, is thatof applying a plurality of methods.
21

References
Alfarano, Simone, Thomas Lux, and Friedrich Wagner (2005). “Estimation of agentbased models: the case of an asymmetric herding model”. In: Computational Economics
26.1.
Alfarano, Simone, Thomas Lux, and Friedrich Wagner (2006). “Estimation of a simple agent-based model of financial markets: An application to australian stock and foreignexchange data”. In: Physica A: Statistical Mechanics and its Applications 370.1, pp. 38–42.
Alfarano, Simone, Thomas Lux, and Friedrich Wagner (2007). “Empirical validationof stochastic models of interacting agents”. In: The European Physical Journal B-Condensed
Matter and Complex Systems 55.2, pp. 183–187.
Barde, Sylvain (2015). “Direct calibration and comparison of agent-based herding models offinancial markets”. In: Studies in Economics, University of Kent.
Bernanke, Ben S. (1986). “Alternative Explanations of the Money-Income Correlation”. In:
NBER Working Paper Series 1842.
Bessler, David A. and Seongpyo Lee (2002). “Money and prices: U.S. Data 1869-1914 (Astudy with directed graphs)”. In: Empirical Economics 27, pp. 427–446.
Bianchi, Carlo, Pasquale Cirillo, Mauro Gallegati, and Pietro Vagliasindi (2007).
“Validating and Calibrating Agent-Based Models: A Case Study”. In: Computational Economics 30, pp. 245–264.
Bianchi, Carlo, Pasquale Cirillo, Mauro Gallegati, and Pietro Vagliasindi (2008).
“Validation in agent-based models: An investigation on the CATS model”. In: Journal of
Economic Behavior & Organization 67, pp. 947–964.
Blanchard, Olivier J. (1989). “A Traditional Interpretation of Macroeconomic Fluctuations”. In: The American Economic Review 79.5, pp. 1146–1164.
Blanchard, Olivier J. and Danny Quah (1989). “The Dynamic Effects of Aggregate Demand and Supply Disturbances”. In: The American Economic Review 79.4, pp. 655–673.
Box, George and Gwilym Jenkins (1970). Time Series Analysis: Forecasting and Control.
Holden-Day.

22

Canova, Fabio and Luca Sala (2009). “Back to Square One: Identification Issues in DSGE
Models”. In: Journal of Monetary Economics 56, pp. 431–449.
Chlass, Nadine and Alessio Moneta (2010). “Can Graphical Causal Inference Be Extended
to Nonlinear Settings? An Assessment of Conditional Independence Tests”. In: EPSA Epistemology and Methodology of Science. Ed. by M. Suarez M. Dorato M. Redei. Springer
Verlag, pp. 63–72.
Comon, Pierre (1994). “Independent component analysis, a new concept?” In: Signal Processing 36.3, pp. 287–314.
Demiralp, Selva and Kevin D. Hoover (2003). “Searching for the Causal Structure of a
Vector Autoregression”. In: Oxford Bullettin of Economics and Statistics 65, pp. 745–767.
Dosi, Giovanni, Giogio Fagiolo, Mauro Napoletano, and Andrea Roventini (2013).
“Income distribution, credit and fiscal policies in an agent-based Keynesian model”. In:
Journal of Economic Dynamics and Control 37.8, pp. 1598–1625.
Dosi, Giovanni, Giorgio Fagiolo, Mauro Napoletano, Andrea Roventini, and Tania
Treibich (2015). “Fiscal and Monetary Policies in Complex Evolving Economies”. In: Journal of Economic Dynamics and Control 52, pp. 166–189.
Dosi, Giovanni, Giorgio Fagiolo, and Andrea Roventini (2010). “Schumpeter meeting
Keynes: A policy-friendly model of endogenous growth and business cycles”. In: Journal of
Economic Dynamics and Control 34.9, pp. 1748–1767.
Fagiolo, Giorgio, Alessio Moneta, and Paul Windrum (2007). “A Critical Guide to
Empirical Validation of Agent-Based Models in Economics and Finance: Methodologies,
Procedures and Open Problems”. In: Computational Economics 30.3, pp. 195–226.
Fagiolo, Giorgio and Andrea Roventini (2012). “Macroeconomic policy in DSGE andagent-based models”. In: Revue de l’OFCE 124, pp. 67–116.
Fukac, Martin and Adrian Pagan (2006). “Issues in Adopting DSGE Models for Use in
the Policy Process”. In: CAMA Working Paper Series 10.
Gilli, Manfred and Peter Winker (2003). “Global optimization heuristic for estimating
agent based models”. In: Computational Statistics & Data Analysis 42.3, pp. 299–312.

23

Granger, Clive W. J. (1980). “Testing for causality: A personal viewpoint”. In: Journal of
Economic Dynamics and Control 2, pp. 329–352.
Grazzini, Jakob (2012). “Analysis of the emergent properties: Stationarity and ergodicity”.
In: Journal of Artificial Societies and Simulation 15.2.
Grazzini, Jakob and Matteo Richiardi (2015). “Estimation of ergodic agent-based models by simulated minimum distance.” In: Journal of Economic Dynamics and Control 51,
pp. 148–165.
Haavelmo, Trygve (1944). “The Probability Approach in Econometrics”. In: Econometrica
12.Supplement.
Haldane, Andrew (2012). The dog and the frisbee. Tech. rep. Speech at the Federal Reserve
Bank of Kansas City 366th economic policy symposium.
Heckman, James (2000). “Causal Parameters and Policy Analysis in Economics: A TwentiethCentury Retrospective”. In: The Quarterly Journal of Economics 115.1, pp. 45–97.
Hoover, Kevin D. (2004). “Lost Causes”. In: Journal of the History of Economic Thought
26.2, pp. 149–164.
Hoyer, Patrik O., Shohei Shimizu, Antti J. Kerminen, and Markus Palviainen(2008). “Estimation of Causal Effects Using Linear Non-Gaussian Causal Models with Hidden Variables”. In: International Journal of Approximate Reasoning 49, pp. 362–278.
Hyvarinen, Aapo, Juha Karhunen, and Erkki Oja (2001). Indepentent Component Analysis. Wiley.
Hyvarinen, Aapo, Kun Zhang, Shohei Shimizu, and Patrick O. Hoyer (2010). “Estimation of a Structural Vector Autoregression Model Using Non-Gaussianity”. In: Journalof Machine Learning Research 11, pp. 1709–1731.
Illari, Phyllis McKay, Federica Russo, and Jon Williamson (2011). “Why look at
Causality in the Sciences? A manifesto”. In: Causality in the Sciences. Oxford University
Press.
Ireland, Peter N (2004). “A method for taking models to the data”. In: Journal of Economic
Dynamics and Control 28.6, pp. 1205–1226.

24

Johansen, Soren and Katerina Juselius (1990). “Maximum Likelihood Estimation and Inference on Cointegration. With Application to the Demand for Money”. In: Oxford Bullettinof Economics and Statistics 52, pp. 169–210.
Juselius, Katarina (2011). “On the role of theory and evidence in macroeconomics”. In: Elgar
Companion to Recent Economic Methodology. Edward Elgar, pp. 404–436.
Koopmans, T. C. (1950). “Statistical Inference in Dynamic Economic Models”. In: Cowles
Commission Monograph. John Wiley & Sons 10.
Kydland, Finn E. and Edward C. Prescott (1982). “Time to Build and Aggregate Fluctuations”. In: Econometrica 50.6, pp. 1345–1370.
Kydland, Finn E and Edward C Prescott (1996). “The computational experiment: an
econometric tool”. In: The Journal of Economic Perspectives, pp. 69–85.
Lacerda, Gustavo, Peter Spirtes, Patrick O. Hoyer, and Joseph Ramsey (2008).
“Discovering Cyclic Causal Models by Independent Components Analysis”. In: Proceedingsof the 24th Conference on Uncertainty in Artificial Intelligence.
Lamperti, Francesco (2015). “An Information Theoretic Criterion for Empirical Validationof Time Series Models”. In: LEM Working Paper.
Lucas, Robert (1976). “Econometric Policy Evaluation: A Critique”. In: Carnegie-Rochester
Conference Series on Public Policy 1.1, pp. 19–46.
Lutkepohl, Helmut (1993). Introduction to Multiple Time Serie Analysis. Springer-Verlag.
Lux, Thomas (2012). “Estimation of an agent-based model of investor sentiment formation in
financial markets”. In: Journal of Economic Dynamics and Control 36.8, pp. 1284–1302.
Marks, Robert (2013).“Validation and Model Selection: three similarity measures compared”.
In: Complexity Economics 2.1.
Moneta, Alessio (2008). “Graphical causal models and VARs: an empirical assessment of thereal business cycles hypothesis”. In: Empirical Economics 35.2, pp. 275–300.
Moneta, Alessio, Nadine Chlass, Doris Entner, and Patrik Hoyer (2011). “Causal
Search in Structural Vector Autoregressive Models”. In: Journal of Machine Learning Research 12, pp. 95–118.

25

Moneta, Alessio, Doris Entner, Patrick Hoyer, and Alex Coad (2013). “Causal Inference by Independent Component Analysis: Theory and Applications”. In: Oxford Bullettinof Economics and Statistics 75.5, pp. 705–730.
Pearl, Judea (2000). Causality: Models, Reasoning and Inference. Cambridge University
Press.
Recchioni, Maria Cristina, Gabriele Tedeschi, and Mauro Gallegati (2015). “A
calibration procedure for analyzing stock price dynamics in an agent-based framework”. In:
Journal of Economic Dynamics and Control 60.1, pp. 1–25.
Salle, Isabelle and Murat Yildizoglu (2014). “Efficient Sampling and Meta-Modelingfor Computational Economic Models”. In: Computational Economics 44.4, pp. 507–536.
Shapiro, Matthew D. and Mark W. Watson (1988). “Sources of Business Cycle Fluctuations”. In: Cowles Foundation Discussion Papers 870.
Shimizu, Shohei, Patrik O. Hoyer, Aapo Hyvarinen, and Antti J. Kerminen (2006).
“A Linear Non-Gaussian Acyclic Model for Causal Discovery”. In: Journal of Machine Learning Research 7, pp. 2003–2030.
Sims, Christopher (1980). “Macroeconomics and Reality”. In: Econometrica 48.1, pp. 1–48.
Spirtes, Peter, Clark Glymur, and Richard Scheines (2000). Causation, Prediction,and Search. MIT Press.
Swanson, Norman R. and Clive W. J. Granger (1997). “Impulse response function basedon a causal approach to residual orthogonalization in vector autoregressions”. In: Journalof the American Statistical Association 92.437, pp. 357–367.
Trichet, Jean Claude (2010). Reflections on the nature of monetary policy non-standard
measures and finance theory. Tech. rep.

26

Appendix A - PC Algorithm
A. Connect everything
Form the complete undirected graph G on the vertex set (u1t , . . . , uKt ) so that each vertex is connected toany other vertex by an undirected edge.
B. Cut some edges
n=0
REPEAT :REPEAT :select an ordered pair of variables uht and uit that are adjacent in G such that the numberof variables adjacent to uht is equal or greater than n + 1. Select a set S of n variables
adjacent to uht suchthat uit ∈/ S. If uht ⊥ uit |S delete edge uht − uit from G.
UNTIL all ordered pairs of adjacent variables uht and uit such that the number of variables
adjacent to uht is equal or greater than n + 1 and all sets S of n variables adjacent to uht
such that uit ∈/ S have been checked to see if uht ⊥ uit |S;n = n + 1;
UNTIL for each ordered pair of adjacent variables uht , uit , the number of adjacent variables to uht is lessthan n + 1.
C. Build colliders
For each triple of vertices uht , uit , ujt such that the pair uht , uit and the pair uit , ujt are each adjacent in G
but the pair uht , ujt is not adjacent in G, orient uht − uit − ujt as uht → uit ← ujt if and only if uit does
not belong to any set of variables S such that uht ⊥ ujt |S.
D. Direct some other edges
REPEAT :if uat → ubt , ubt and uct are adjacent, uat and uct are not adjacent and ubt belongs to every set
S such that uat ⊥ uct |S, then orient ubt − uct as ubt → uct ; if there is a directed path fromuat to ubt and an edge between uat and ubt , then orient uat − ubt as uat → ubt ;UNTIL no more edges can be oriented.

27

Appendix B - VAR-LiNGAM Algorithm
A. Estimate the reduced form VAR model of equation (3) obtaining estimates Âi of the matrices Ai , ∀ i =1, . . . , p. Denote by Û the K × T matrix of the corresponding estimated VAR error terms, that is eachcolumn of Û is ût ≡ (û1t , . . . , ûKt )′ , ∀ t = 1, . . . , T . Check whether the uit (for all rows i) indeed are
non-Gaussian, and proceed only if this is so.
B. Use FastICA or any other suitable ICA algorithm (Hyvarinen et al., 2001) to obtain a decompositionÛ = PE where P is K × K and E is K × T , such that the rows of E are the estimated independent
components of Û. Then validate non-Gaussianity and (at least approximate) statistical independenceof the components before proceeding.
˜ = P−1 . Find Γ̃ , the row-permutated version of Γ̃˜ which minimizes P 1 with respect to
C. Let Γ̃
0
0
0
i |Γ̃
|
0,ii

the permutation. Note that this is a linear matching problem which can be easily solved even for highK (Shimizu et al., 2006).
D. Divide each row of Γ̃0 by its diagonal element, to obtain a matrix Γ̂0 with all ones on the diagonal.
E. Let B̃ = I − Γ̂0 .
F. Find the permutation matrix Z which makes ZB̃ZT as close as possible to lower triangular. Thiscan be formalized as minimizing the sum of squares of the permuted upper-triangular elements, andminimized using a heuristic procedure (Shimizu et al., 2006). Set the upper-triangular elements to zero,and permute back to obtain B̂ which now contains the acyclic contemporaneous structure. (Note that
it is useful to check that ZB̃ZT indeed is close to strictly lower-triangular).
G. B̂ now contains K(K − 1)/2 non-zero elements, some of which may be very small (and statisticallyinsignificant). For improved interpretation and visualization, it may be desired to prune out (set tozero) small elements at this stage, for instance using a bootstrap approach (Shimizu et al., 2006).
H. Finally, calculate estimates of Γi , ∀ i = 1, . . . , p for lagged effects using Γi = (I − B̂)Âi .

28