This page describes how to use Python libraries for macroeconomic modelling, allowing to extend Dynare functionalities.

Python ressources

• Formal calculus

  • sympy python library to perform formal calculus.

  • sage (includes sympy) python wrapper around many opensource libraries which converts all of them into a single webbased interface. Is also able to delegate calculus to Maple or Mathematica.

• Numerical computing

  • numpy : contains an efficient implementation of matrix processing whose performance is close to matlab

  • scipy : set of scientific routines including optimization routines

  • pyrex and psyco (included in scipy) can speedup computations. The first allows to write C extensions in a python-like syntax, while the second performs just-in-time compilation of python functions.

• Other

  • epydoc generates online documentation from inline documentation

  • mlabwrap uses the Matlab Engine to interact with Matlab. It permits to send and retrieve data, as well as to call Matlab function from Python.

  • rpy : low-level access to statistical program R permitting to call every R functions as python functions.

  • py2exe packs a python program with all its libraries in a single windows executable which includes the interpreter.

  • pymaclab another attempt to use python estimate DSGE models. It uses its own Modfile format which is then parsed by Python.

  • qz.py : a wrapper around LAPACK's qz routine.

  • eclipse and pydev provide a good development environment for Python.

  • IPython is a python console which may prove useful when learning new commands. In particular command? returns incode documentation for command command.

Writing/extending modfiles in Python

There is a small prototype in python, called DareDare] which can be called from the command line for quick processing or conversion (not functional) or can be imported as python library. Its initial aim was to extend Dynare to deal with the formal representation of objects.

Here is the code of ramsey.py, which translates ramst.mod into python code. It is run by the standard Python interpreter. Line-by-line explication follows.

   1 ####################
   2 # Model Definition #
   3 ####################
   4 
   5 from daredare import *
   6 
   7 set_variables("k c")
   8 set_exovariables("z")
   9 set_shocks("x")
  10 set_parameters("alpha beta gamma delta aa")
  11 
  12 equations = [
  13     Equation(c + k , (1 - delta)*k(-1) + k(-1)**alpha*(1 + z)*aa),
  14     Equation(0 , -1/(1 + beta)*c(1)**(-gamma)*(1 - delta + k**(-1 + alpha)*(1 + z(1))*aa*alpha) + c**(-gamma)),
  15     Equation(z(1),x(1))
  16 ]
  17 
  18 parameters_values = {
  19     alpha : 0.5,
  20     gamma : 1*alpha,
  21     delta : 0.02,
  22     beta : 0.05,
  23     aa : 0.5
  24 }
  25 
  26 covariances = matrix([[0.01]])
  27 
  28 init_values = {
  29     k:(1/aa/alpha*(beta + delta))**(1/(-1 + alpha)),
  30     c:k**alpha*aa - delta*k
  31 }
  32 
  33 ####################
  34 # Model processing #
  35 ####################
  36 
  37 model = Model("ramsey",lookup=True)
  38 model.check_all('uhlig')
  39 
  40 solver = Solver(model,lookup=True)
  41 
  42 ss = solver.find_steady_state()
  43 print(ss)
  44 
  45 res = solver.solve_with_uhligs_toolkit()
  46 
  47 
  48 print('PP',res.PP)
  49 print('QQ',res.QQ)

Language extension (for model definition)

On the very first line of this code from daredare import * imports the library. The second block

   1 set_variables("k c")
   2 set_exovariables("z")
   3 set_shocks("x")
   4 set_parameters("alpha beta gamma delta aa")

creates formal objects "k,c,z" which all derive of the same class Variable. This class is a subclass of the standard Sympy class Symbol, which only adds the possibility to take lags or leads with a function call : for instance k(1) returns another Variable whose formal name is  k_{+1} .  k.lag  returns the lag so that  k(-k.lag)  or equivalently k.p() will always return a variable with lag 0.

Parameters are of the class Parameter which derives also from Symbol but cannot be called.

It should be noted that, a priori, the formal object, and the Python reference name which is used by the compiler are to different things. For instance,

   1 somestrangename = Variable('x')
   2 anotherstrangename = somestrangename
   3 y = Variable('y')
   4 eq = somestrangename * anotherstrangename / y
   5 print(eq)

will print : x**2/y. In this respect the functions like set_variables are notable because instead of returning the created objects, they also add references to the interpreter context, with the same names as the variables. They also add the python names variables,shocks,exovariables,parameters. The next block

   1 equations = [
   2     Equation(c + k , (1 - delta)*k(-1) + k(-1)**alpha*(1 + z)*aa),
   3     Equation(0 , -1/(1 + beta)*c(1)**(-gamma)*(1 - delta + k**(-1 + alpha)*(1 + z(1))*aa*alpha) + c**(-gamma)),
   4     Equation(z(1),x)
   5 ]

defines equations. The object Equation is a child of Sympy's Equality object. It contains the expression as well as some attributes as (optionally) equation name, expectational nature, and other information that could be computed later. Three remarks :

• unfortunately we cannot write x+y == z because in python == must always return true or false (see Sympy mailing list).

• to convert from/to matlab, you only have to convert * from/to ^. For instance, if eq is a python formal expression, str(eq).replace('==','=').replace('**','^') will give the valid Dynare expression.

• every sympy operations on equalities are available

  • printing.latex(eq) (or printing.mathml(eq) will return the latex (or mathml) representation of the equality

  • preview(eq) will launch an viewer to preview the latexified output

  • eq.rhs.diff(c) will return the right side of the equation differentiated with respect to c

  • eq.rhs.subs(c(1),c) will replace c(1) by c in the right side

To define initial values, we use dictionaries whose whose keys are the symbolic and whose values are either numeric or symbolic. Dictionaries values are accesses using brackets.

   1 init_values[k]

will return the symbolic expression k:(1/aa/alpha*(beta + delta))**(1/(-1 + alpha))

Processing model(s)

In the next lines we create the model and check its consistency :

   1 model = Model("ramsey",lookup=True)
   2 model.check_all('uhlig')

Note the option lookup=True. Normally, in the scope of a function, you can only access local variables and arguments. model function uses python's introspection routines to lookup in the stack for the lists variables,exovariables,shocks,parameters and equations and to attach them to the object model. We have already used this magic in the functions set_variables. Then, equations can be accessed with model.equations.

The second line checks for model consistency. So far the formal objects we have created are totally agnostic to modelling conventions which vary across software and users. The instance function check_all makes some basic checks to see if the model is consistent with one formulation. As for now, it accepts two arguments 'uhlig' and 'dynare'. Note that exovariables as it is defined in Uhlig's conventions don't exist in Dynare. shocks denote independent random variables and correspond to the varexo statement in Dynare.

   1 solver = Solver(model,lookup=True)
   2 res = solver.solve_with_uhligs_toolkit()

A solver object depends on a model and maps a set of parameters to a decision rule. Here, we have only one value for this set, so we can get them from the initial declaration using lookup=True. The second statement is useless when Dynare is the destination. It returns a dictionary containing the steady state. Many steps are involved :

• Initial values and initial parameters are computed from the formal declarations. This involve, solving a triangular system, which is easy.
• Formal equations are converted into a numeric python function which is added to interpreting context at runtime.
• This function is optimized to return the steady state

   1 res = solver.solve_with_uhligs_toolkit()
   2 print('PP',res.PP)
   3 print('QQ',res.QQ)

The last lines compute the decision rule. In the case of Uhlig, it returns two matrix representing the decision rules.

If the model were written following Dynare convention, we would have called :

   1 res = send_to_octave(solver,interactive=False,append="steady;\ncheck;\nstoch_simul(nograph);")
   2 M = res['M_']
   3 oo = res['oo_']
   4 options = res['options_']

to get the results from octave or

   1 res = send_to_matlab(solver,interactive=False,append="steady;\ncheck;\nstoch_simul(nograph);")
   2 print(res)
   3 M = res['M_']
   4 oo = res['oo_']
   5 options = res['options_']

to get the results from Matlab. For Octave, we write a modfile and send it to a new session while we use the Matlab engine to send it to Matlab with the added benefit that we can do many calculations without clearing the Matlab workspace each time.

Import and export

Modfile import

Instead of rewriting a model that is already in a modfile, we can import it to Python with the command import_modfile and process the model as seen previously. The command adds all variables of the model block that are founds in the modfile to python interpreter. It doesn't even try to understand other elements of modfile syntax such as macrocommands and additional processing commands.

Given that a ramst.mod file exists, the following code does just that and replaces all occurrences of c by log(c/css) where css is a parameter

   1 import_modfile('ramst.mod')
   2 
   3 # following two lines could be collapsed in a new function add_parameter('css') which would do :
   4 css = Parameter('css') 
   5 parameters.append( css )
   6 
   7 for i in range(equations):
   8     equations[i] = equations[i].replace(c,log(c/css))
   9 
  10 ###########################
  11 # End of model definition #
  12 ###########################

Modfile export

A solver object can be exported to a modfile and written to a file using :

   1 modfile = solver.export_to_modfile(append="stoch_simul(nograph);")
   2 f = file('output')
   3 f.write(modfile)
   4 f.close()

XML serialization

A solver object can be serialized to XML using

   1 xmltext = solver.toxml()

The resulting string is a quick prototype where stricly needed informations are stored (equations are stored unparsed in Dynare format) it could be adapted easily to communicate with Dynare or another program using a simple standardized file.