stepan
/
dynare
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge remote-tracking branch 'Dynare/master' into enterprise

time-shift
Stéphane Adjemian (Odysseus) 2021-01-09 17:56:38 +01:00 committed by Sébastien Villemot
parent e066556d7a c3e7ef42af
commit 0e1c9312e2
No known key found for this signature in database
GPG Key ID: 2CECE9350ECEBE4A
199 changed files with 8055 additions and 3293 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

5
.gitignore vendored
View File

@ -70,12 +70,9 @@ doc/internals/*.html
doc/internals/ltxpng
# MATLAB dir
/matlab/preprocessor*
/matlab/preprocessor64/
/matlab/dynare_version.m
# JULIA dir
/julia/preprocessor*
# DLL rules
*.mex
*.oct

2
.gitlab-ci.yml
View File

@ -33,7 +33,7 @@ build_binaries:
- make -j $(nproc) LN_S="cp -p"
artifacts:
paths:
- matlab/preprocessor*/*
- preprocessor/dynare-preprocessor
- mex/octave/
- mex/matlab/
- dynare++/parser/cc/*_tab.cc

43
Makefile.am
View File

@ -27,25 +27,29 @@ EXTRA_DIST = \
scripts \
.dir-locals.el
all-local: preprocessor/src/dynare_m$(EXEEXT)
if file preprocessor/src/dynare_m$(EXEEXT) | grep -q x86.64; then \
ARCH="64"; \
else \
ARCH="32"; \
fi; \
mkdir -p $(abs_srcdir)/matlab/preprocessor$$ARCH && \
$(LN_S) -f $(abs_builddir)/preprocessor/src/dynare_m$(EXEEXT) $(abs_srcdir)/matlab/preprocessor$$ARCH
all-local:
# Create top-level preprocessor symlink needed by matlab/dynare.m (when Dynare is run
# from the locally-built copy)
$(LN_S) -f $(abs_builddir)/preprocessor/src/dynare-preprocessor$(EXEEXT) $(abs_builddir)/preprocessor/dynare-preprocessor$(EXEEXT)
# Create backward-compatibility symlink for old location of preprocessor
$(MKDIR_P) matlab/preprocessor64/
$(LN_S) -f $(abs_builddir)/preprocessor/src/dynare-preprocessor$(EXEEXT) $(abs_srcdir)/matlab/preprocessor64/dynare_m$(EXEEXT)
clean-local:
rm -rf $(abs_srcdir)/matlab/preprocessor32 $(abs_srcdir)/matlab/preprocessor64
rm -f preprocessor/dynare-preprocessor$(EXEEXT)
rm -rf matlab/preprocessor64/
dist-hook:
rm -rf `find $(distdir)/matlab $(distdir)/examples -name *~`
rm -rf $(distdir)/matlab/preprocessor* $(distdir)/matlab/dynare_version.m
rm -f $(distdir)/matlab/dynare_version.m
rm -rf $(distdir)/matlab/preprocessor64/
$(MKDIR_P) $(distdir)/mex/matlab $(distdir)/mex/octave
rm -rf `find $(distdir) -name '.git*'`
install-exec-local:
$(MKDIR_P) $(DESTDIR)$(pkglibdir)/preprocessor
# The following will break under Windows, but we don’t use the install rule there
$(LN_S) -f $(bindir)/dynare-preprocessor$(EXEEXT) $(DESTDIR)$(pkglibdir)/preprocessor/dynare-preprocessor$(EXEEXT)
$(MKDIR_P) $(DESTDIR)$(pkglibdir)/contrib/ms-sbvar/TZcode
cp -r contrib/ms-sbvar/TZcode/MatlabFiles $(DESTDIR)$(pkglibdir)/contrib/ms-sbvar/TZcode
$(MKDIR_P) $(DESTDIR)$(pkglibdir)/contrib/jsonlab
@ -53,14 +57,17 @@ install-exec-local:
cp -r examples $(DESTDIR)$(pkglibdir)
cp -r matlab $(DESTDIR)$(pkglibdir)
find $(DESTDIR)$(pkglibdir) -name LICENSE.md -delete
rm -rf $(DESTDIR)$(pkglibdir)/matlab/preprocessor*
if file preprocessor/src/dynare_m | grep -q x86.64; then \
ARCH="64"; \
else \
ARCH="32"; \
fi; \
mkdir -p $(DESTDIR)$(pkglibdir)/matlab/preprocessor$$ARCH && \
cp preprocessor/src/dynare_m $(DESTDIR)$(pkglibdir)/matlab/preprocessor$$ARCH
# Recreate backward-compatibility symlink
rm -f $(DESTDIR)$(pkglibdir)/matlab/preprocessor64/dynare_m$(EXEEXT)
$(LN_S) -f $(bindir)/dynare-preprocessor$(EXEEXT) $(DESTDIR)$(pkglibdir)/matlab/preprocessor64/dynare_m$(EXEEXT)
# We don’t fail over doc install rules, since the user may deliberately not have compiled them
# (e.g. with the “nodoc” option under Debian)
$(MKDIR_P) $(docdir)
-cp doc/*.pdf doc/gsa/gsa.pdf doc/parallel/parallel.pdf doc/dseries-and-reporting/dseriesReporting.pdf preprocessor/doc/preprocessor/preprocessor.pdf preprocessor/doc/macroprocessor/macroprocessor.pdf doc/manual/build/latex/dynare-manual.pdf $(docdir)
-cp -r doc/manual/build/html $(docdir)/dynare-manual.html
$(MKDIR_P) $(docdir)/dynare++
-cp dynare++/doc/*.pdf $(docdir)/dynare++
uninstall-local:
rm -rf $(DESTDIR)$(pkglibdir)
rm -rf $(DESTDIR)$(docdir)

41
NEWS.md
View File

@ -1,3 +1,44 @@
Announcement for Dynare 4.6.3 (on 2020-11-23)
=============================================
We are pleased to announce the release of Dynare 4.6.3.
This maintenance release fixes various bugs.
The Windows, macOS and source packages are already available for download at
[the Dynare website](https://www.dynare.org/download/).
All users are strongly encouraged to upgrade.
This release is compatible with MATLAB versions ranging from 7.9 (R2009b) to
9.9 (R2020b), and with GNU Octave versions 5.2.0 (under Windows) and 4.4.1
(under macOS).
Here is a list of the problems identified in version 4.6.2 and that have been
fixed in version 4.6.3:
* Using an unknown symbol in `irf_shocks` option of `stoch_simul` would lead to
a crash of the preprocessor
* `identification` would crash for purely forward-looking models
* The `endogenous_prior` option did not properly handle missing observations
* The auxiliary particle filter with pruning and resampling would crash
* Initialization of the state variance for particle filters was buggy
* An `@#else` clause after an `@#ifndef` was not correctly interpreted
* An `@#elseif` clause after an `@#ifdef` or an `@#ifndef` was not correctly
interpreted
* Perfect foresight simulations of models with a single equation would crash
when using either the `lmmcp` option or the `linear_approximation`
* Inequality constraints on endogenous variables (when using the `lmmcp`
option) were not enforced on purely backward or purely forward models
* Perfect foresight simulations with `bytecode` and `block` options could crash
if there was a purely forward variable whose value in all periods could be
evaluated backward (typically a process of the form `y=a*y(+1)+e`)
* `extended_path` was broken with `bytecode`
* Under Windows, with Octave, the k-order perturbation and MS-SBVAR MEX files
could not be loaded
* On Fedora (and possibly other GNU/Linux distributions), compilation from
source would fail against Octave 5
Announcement for Dynare 4.6.2 (on 2020-09-07)
=============================================

296
README.md
View File

@ -56,6 +56,7 @@ a 32-bit Octave.
1. [**General Instructions**](#general-instructions)
1. [**Debian or Ubuntu**](#debian-or-ubuntu)
1. [**Fedora, CentOS or RHEL**](#fedora-centos-or-rhel)
1. [**Windows**](#windows)
1. [**macOS**](#macos)
@ -68,9 +69,12 @@ A number of tools and libraries are needed in order to recompile everything. You
- A POSIX compliant shell and an implementation of Make (mandatory)
- The [GNU Compiler Collection](http://gcc.gnu.org/), version 8 or later, with
gcc, g++ and gfortran (mandatory)
- MATLAB (if you want to compile the MEX for MATLAB)
- [GNU Octave](http://www.octave.org), with the development headers (if you
want to compile the MEX for Octave)
- [MATLAB](https://mathworks.com) (if you want to compile the MEX for MATLAB)
- [GNU Octave](http://www.octave.org) with
- the development headers (if you want to compile the MEX for Octave)
- the development libraries corresponding to the [UMFPACK](https://people.engr.tamu.edu/davis/suitesparse.html) packaged with Octave
- Optionally, the [Control](https://wiki.octave.org/Control_package), [IO](https://wiki.octave.org/IO_package), [Optimization](https://wiki.octave.org/Optimization_package) and [Statistics](https://wiki.octave.org/Statistics_package) package either installed via your package manager or through [Octave Forge](https://wiki.octave.org/Category:Octave_Forge).
- [Boost libraries](http://www.boost.org), version 1.36 or later
- [Bison](http://www.gnu.org/software/bison/), version 3.2 or later (only if you get the source through Git)
- [Flex](http://flex.sourceforge.net/), version 2.5.4 or later (only if you get the source through Git)
@ -86,48 +90,35 @@ A number of tools and libraries are needed in order to recompile everything. You
- [Sphinx](http://www.sphinx-doc.org/)
- [MathJax](https://www.mathjax.org/)
- [Doxygen](http://www.stack.nl/%7Edimitri/doxygen/) (if you want to build Dynare preprocessor source documentation)
- For Octave, the development libraries corresponding to the UMFPACK packaged with Octave
- [X-13ARIMA-SEATS Seasonal Adjustment Program](https://www.census.gov/srd/www/x13as/)
### Preparing the sources
If you have downloaded the sources from an official source archive or the source snapshot, just unpack it.
If you want to use Git, do the following from a terminal:
git clone --recurse-submodules https://git.dynare.org/Dynare/dynare.git
cd dynare
autoreconf -si
The last line runs Autoconf and Automake in order to prepare the build environment (this is not necessary if you got the sources from an official source archive or the source snapshot).
```sh
git clone --recurse-submodules https://git.dynare.org/Dynare/dynare.git
cd dynare
autoreconf -si
```
The last line runs Autoconf and Automake in order to prepare the build environment (this is not necessary if you got the sources from an official source archive or the source snapshot). If you want a certain version (e.g. 4.6) , then add `--single-branch --branch 4.6` to the git clone command.
### Configuring the build tree
Simply launch the configure script from a terminal:
```sh
./configure --with-matlab=<…>
```
./configure
```
If you have MATLAB, you need to indicate both the MATLAB location. For example, on GNU/Linux:
```
./configure --with-matlab=/usr/local/MATLAB/R2019a
```
where the path to MATLAB is specified.
Alternatively, you can disable the compilation of MEX files for MATLAB with the `--disable-matlab` flag, and MEX files for Octave with `--disable-octave`.
Some important options:
You may need to specify additional options to the configure script, see the
output of the `--help` option, and also the platform specific instructions
below.
- `--disable-matlab`: skip the compilation of MEX files for MATLAB
- `--disable-octave`: skip the compilation of MEX files for Octave
- `--disable-doc`: skip the compilation of the documentation (PDF and HTML)
Note that if you don't want to compile the C/C++ programs with debugging information, you can specify the `CFLAGS` and `CXXFLAGS` variables to the configure script, such as:
```
./configure CFLAGS="-O3" CXXFLAGS="-O3"
```
To remove debugging information for MATLAB MEX functions, the analogous call would be:
```
./configure MATLAB_MEX_CFLAGS="-O3" MATLAB_MEX_CXXFLAGS="-O3"
```
If the configuration goes well, the script will tell you which components are
correctly configured and will be built.
You may need to specify additional options to the configure script, see the output of the `--help` option, and also the platform specific instructions below. If the configuration goes well, the script will tell you which components are correctly configured and will be built.
Note that it is possible that some MEX files cannot be compiled, due to missing
build dependencies. If you find no way of installing the missing dependencies,
@ -143,70 +134,62 @@ addpath <DYNARE_ROOT>/matlab/missing/mex/gensylv
### Building
Binaries and Info documentation are built with:
```
Binaries are built with:
```sh
make
```
PDF and HTML documentation are respectively built with:
```
```sh
make pdf
make html
```
The testsuites can be run with:
```
make check
```
Note that running the testsuite with Octave requires the additional packages
`pstoedit`, `epstool`, `xfig`, and `gnuplot`.
### Check
The Git source comes with unit tests (in the MATLAB functions) and integration tests (under the `tests` subfolder). All the tests can be run with:
```
```sh
make check
```
In the `tests` subfolder. If Dynare has been compiled against MATLAB and Octave, the tests will be run with MATLAB and Octave. Depending on
the performance of your machine, this can take several hours. It is possible to run the tests only with MATLAB:
```
in the `tests` subfolder. If Dynare has been compiled against MATLAB and Octave, the tests will be run with both MATLAB and Octave. Depending on the performance of your machine, this can take several hours. It is possible to run the tests only with MATLAB:
```sh
make check-matlab
```
or only with Octave:
```
```sh
make check-octave
```
A summary of the results is available in `tests/run_test_matlab_output.txt` or `tests/run_test_octave_output.txt`. Often, it does not make sense
to run the complete testsuite. For instance, if you modify codes only related to the perfect foresight model solver, you can decide to run only a
subset of the integration tests, with:
```
Note that running the testsuite with Octave requires the additional packages `pstoedit`, `epstool`, `xfig`, and `gnuplot`.
A summary of the results is available in `tests/run_test_matlab_output.txt` or `tests/run_test_octave_output.txt`. Often, it does not make sense to run the complete testsuite. For instance, if you modify codes only related to the perfect foresight model solver, you can decide to run only a subset of the integration tests, with:
```sh
make deterministic_simulations
```
This will run all the integration tests in `tests/deterministic_simulations` with MATLAB and Octave. Again, it is possible to do this only with MATLAB:
```
```sh
make m/deterministic_simulations
```
or with Octave:
```
```sh
make o/deterministic_simulations
```
Finally if you want to run a single integration test, e.g. `deterministic_simulations/lbj/rbc.mod` with MATLAB:
```
```sh
make deterministic_simulations/lbj/rbc.m.trs
```
or with Octave:
```
```sh
make deterministic_simulations/lbj/rbc.o.trs
```
The result of the test (`PASSED` or `FAILED`) will be printed in the terminal, the produced log can be displayed with:
```
```sh
make deterministic_simulations/lbj/rbc.m.drs
```
or
```
```sh
make deterministic_simulations/lbj/rbc.o.drs
```
Note that only tests will be executed where the `m.trs/o.trs` does not yet exist. You can run
```
```sh
make clean
```
in the `tests` folder to delete files that were created by the run of the testsuite. You can also manually delete the desired `m.trs/o.trs` file(s).
@ -239,11 +222,120 @@ All the prerequisites are packaged:
- `latexmk`
- `libjs-mathjax`
- `doxygen`
- `x13as`
You can install them all at once with:
```sh
apt install build-essential gfortran liboctave-dev libboost-graph-dev libgsl-dev libmatio-dev libslicot-dev libslicot-pic libsuitesparse-dev flex bison autoconf automake texlive texlive-publishers texlive-latex-extra texlive-fonts-extra texlive-latex-recommended texlive-science texlive-plain-generic lmodern python3-sphinx latexmk libjs-mathjax doxygen x13as
```
apt install build-essential gfortran liboctave-dev libboost-graph-dev libgsl-dev libmatio-dev libslicot-dev libslicot-pic libsuitesparse-dev flex bison autoconf automake texlive texlive-publishers texlive-latex-extra texlive-fonts-extra texlive-latex-recommended texlive-science texlive-plain-generic lmodern python3-sphinx latexmk libjs-mathjax doxygen
If you use MATLAB, we strongly advise to also `apt install matlab-support` and confirm to rename the GCC libraries shipped with MATLAB to avoid possible conflicts with GCC libraries shipped by your distribution.
Tested on
- Debian 10
- Ubuntu 20.04
- Ubuntu 20.10
## Fedora, CentOS or RHEL
Almost all prerequisites are packaged:
- `gcc`, `gcc-c++`, `make`
- `gcc-gfortran`
- `lapack` and `lapack-devel`
- `openblas` and `openblas-devel`
- `boost-devel`
- `gsl-devel`
- `matio-devel`
- `suitesparse-devel`
- `flex`
- `bison`
- `autoconf`
- `automake`
- `redhat-rpm-config`
- `octave`, `octave-devel`, `octave-statistics`, `octave-io`, `octave-optim`, `octave-control`
- `texlive-scheme-minimal`, `texlive-collection-publishers`, `texlive-collection-latexextra`, `texlive-collection-fontsextra`, `texlive-collection-latexrecommended`, `texlive-collection-science`, `texlive-collection-plaingeneric`, `texlive-lm`
- `python3-sphinx`
- `latexmk`
- `mathjax`
- `doxygen`
You can install them all at once with:
```sh
# Minimal packages (use --disable-doc and --disable-octave flags)
dnf install -y gcc gcc-c++ make gcc-gfortran lapack lapack-devel openblas openblas-devel boost-devel gsl-devel matio-devel suitesparse-devel flex bison autoconf automake redhat-rpm-config
# Octave packages (use --disable-doc flag)
dnf install octave octave-devel octave-statistics octave-io octave-optim octave-control
# Documentation packages
dnf install texlive-scheme-minimal texlive-collection-publishers texlive-collection-latexextra texlive-collection-fontsextra texlive-collection-latexrecommended texlive-collection-science texlive-collection-plaingeneric texlive-lm python3-sphinx latexmk mathjax doxygen
```
In Fedora these are available from the default repositories; whereas for CentOS and RHEL you need to enable the [Extra Packages for Enterprise Linux (EPEL)](https://fedoraproject.org/wiki/EPEL) repository and either the PowerTools repository for CentOS or the CodeReady Linux Builder repository for RHEL:
```sh
yum install https://dl.fedoraproject.org/pub/epel/epel-release-latest-8.noarch.rpm
# CentOS 8
dnf config-manager --set-enabled PowerTools
# RHEL 8
ARCH=$( /bin/arch )
subscription-manager repos --enable "codeready-builder-for-rhel-8-${ARCH}-rpms"
```
The documentation packages have slightly different names in CentOS and RHEL, you can also choose to pass the `--disable-doc` flag to your configure script to skip these dependencies.
`Slicot` and `x13as` need to be compiled from source:
```sh
# compile slicot from source and put it into /home/$USER/dynare/slicot/lib/
mkdir -p /home/$USER/dynare/slicot
cd /home/$USER/dynare/slicot
wget https://deb.debian.org/debian/pool/main/s/slicot/slicot_5.0+20101122.orig.tar.gz
tar xf slicot_5.0+20101122.orig.tar.gz
cd slicot-5.0+20101122
make FORTRAN=gfortran OPTS="-O2 -fPIC -fdefault-integer-8" LOADER=gfortran lib
mkdir -p /home/$USER/dynare/slicot/lib
cp slicot.a /home/$USER/dynare/slicot/lib/libslicot64_pic.a #for matlab
cp slicot.a /home/$USER/dynare/slicot/lib/libslicot_pic.a #for octave
# compile x13as from source and put it into /usr/bin/
mkdir -p /home/$USER/dynare/x13as
cd /home/$USER/dynare/x13as
wget https://www.census.gov/ts/x13as/unix/x13assrc_V1.1_B39.tar.gz
tar xf x13assrc_V1.1_B39.tar.gz
sed -i "s|-static| |" makefile.gf # this removes '-static' in the makefile.gf
make -f makefile.gf FFLAGS="-O2 -std=legacy" PROGRAM=x13as
sudo cp x13as /usr/bin/
```
If you use MATLAB, we strongly advise to also rename or exclude the GCC libraries shipped with MATLAB to avoid possible conflicts with GCC libraries shipped by Fedora, see e.g. [Matlab on Fedora 33](https://mutschler.eu/linux/install-guides/fedora-post-install/#matlab) or [MATLAB-ArchWiki](https://wiki.archlinux.org/index.php/MATLAB) for instructions.
Keep in mind to use the `--with-slicot` option to the configure command, e.g.:
```sh
cd /home/$USER/dynare
git clone --recurse-submodules https://git.dynare.org/dynare/dynare.git unstable
cd unstable
autoreconf -si
./configure --with-slicot=/home/$USER/dynare/slicot --with-matlab=/usr/local/MATLAB/R2020b
make -j$(($(nproc)+1)) #rule of thumb: one more than CPUs as shown by e.g. lscpu
```
If your distribution ships an older version of `bison`, compile it from source and append it *temporarily* to your path before calling the configure script:
```sh
bison --version # bison (GNU Bison) 3.0.4
mkdir -p /home/$USER/dynare/bison
cd /home/$USER/dynare/bison
wget http://ftp.gnu.org/gnu/bison/bison-3.6.4.tar.gz #change the version number accordingly
tar xf bison-3.6.4.tar.gz
cd bison-3.6.4
./configure --prefix=/home/$USER/dynare/bison
make
make install
export PATH=/home/$USER/dynare/bison/bin:$PATH
bison --version # bison (GNU Bison) 3.6.4
```
Now configure dynare as above.
Tested on
- CentOS 8
- Fedora Workstation 32
- Fedora Workstation 33
- Red Hat Enterprise Linux 8
## Windows
@ -251,18 +343,18 @@ apt install build-essential gfortran liboctave-dev libboost-graph-dev libgsl-dev
have a 32-bit Windows, in which case see below)
- Run a MSYS MinGW 64-bit shell
- Update the system:
```
```sh
pacman -Syu
```
You may be asked to close the window at the end of the
first upgrade batch, in which case you should rerun the upgrade in a new
window to complete the upgrade.
- Install all needed dependencies:
```
```sh
pacman -S git autoconf automake-wrapper bison flex make tar texinfo mingw-w64-x86_64-gcc mingw-w64-x86_64-gcc-fortran mingw-w64-x86_64-boost mingw-w64-x86_64-gsl mingw-w64-x86_64-matio mingw-w64-x86_64-openblas
```
- Compile and install SLICOT, needed for the `kalman_steady_state` MEX file
```
```sh
wget https://deb.debian.org/debian/pool/main/s/slicot/slicot_5.0+20101122.orig.tar.gz
tar xf slicot_5.0+20101122.orig.tar.gz
cd slicot-5.0+20101122
@ -272,13 +364,13 @@ cp slicot.a /usr/local/lib/libslicot64_pic.a
cd ..
```
- Prepare the Dynare sources, either by unpacking the source tarball, or with:
```
```sh
git clone --recurse-submodules https://git.dynare.org/Dynare/dynare.git
cd dynare
autoreconf -si
```
- Configure Dynare from the source directory:
```
```sh
./configure --with-slicot=/usr/local --with-matlab=<…> --disable-octave --disable-doc
```
where the path of MATLAB is specified. Note that you should use
@ -288,11 +380,11 @@ does not have short filenames (8dot3), then you can run `mkdir -p
/usr/local/MATLAB && mount c:/Program\ Files/MATLAB /usr/local/MATLAB`, and
then pass `/usr/local/MATLAB/…` as MATLAB path to the configure script.
- Compile:
```
```sh
make
```
- Run the testsuite:
```
```sh
make -C tests check-matlab
```
@ -309,43 +401,53 @@ currently not supported.
## macOS
To simply use a snapshot of Dynare, you have two choices. On MATLAB, you can
use the [snapshot build](https://www.dynare.org/snapshot/macos/) provided by
Dynare. On Octave, you can simply install [Homebrew](https://brew.sh/) and run
```brew install dynare --HEAD``` (See the Install Dynare (unstable) section of
[this webpage](https://archives.dynare.org/DynareWiki/InstallOnMacOSX) for more
details).
If you do not wish to use the snapshots provided by Dynare or Homebrew, follow
the directions below to build Dynare on your local machine.
Preparatory work:
- Install the Xcode Command Line Tools:
- Open Terminal.app and type `xcode-select --install`
- Install the Xcode Command Line Tools. Open Terminal.app and type:
```sh
xcode-select --install
```
- Install [Homebrew](https://brew.sh/) by following the instructions on their website
The following commands will install the programs that Dynare needs to
compile. They should be entered at the command prompt in Terminal.app.
- `brew install automake bison flex boost fftw gcc gsl hdf5 libmatio metis veclibfort`
- **(Optional)** To compile Dynare mex files for use on Octave:
- `brew install octave`
- **(Optional)** To compile Dynare documentation
- Install the latest version of [MacTeX](http://www.tug.org/mactex/), deselecting the option to install Ghostscript
- `brew install doxygen latex2html`
- Install [MacTeX](http://www.tug.org/mactex/index.html). Alternatively, if you
don’t want to install MacTeX, you should pass the `--disable-doc` flag to the
`configure` command below.
- Install required Homebrew packages. Open Terminal.app and type:
```sh
brew install automake bison flex boost gcc gsl libmatio veclibfort octave sphinx-doc wget
brew link --force sphinx-doc
```
- Compile and install SLICOT, needed for the `kalman_steady_state` MEX file.
Still from Terminal.app:
```sh
wget https://deb.debian.org/debian/pool/main/s/slicot/slicot_5.0+20101122.orig.tar.gz
tar xf slicot_5.0+20101122.orig.tar.gz
cd slicot-5.0+20101122
make -j$(nproc) FORTRAN=gfortran OPTS="-O2" LOADER=gfortran lib
cp slicot.a /usr/local/lib/libslicot_pic.a
make clean
make -j$(nproc) FORTRAN=gfortran OPTS="-O2 -fdefault-integer-8" LOADER=gfortran lib
cp slicot.a /usr/local/lib/libslicot64_pic.a
cd ..
```
The following commands will download the Dynare source code and compile
it. They should be entered at the command prompt in Terminal.app from the
folder where you want Dynare installed.
- `git clone --recurse-submodules https://git.dynare.org/Dynare/dynare.git`
- `cd dynare`
- `PATH="/usr/local/opt/bison/bin:/usr/local/opt/flex/bin:$PATH"`
- `autoreconf -si`
- `CC=gcc-9 CXX=g++-9 ./configure --disable-octave --with-matlab=/Applications/MATLAB_R2019a.app --with-matio=/usr/local --with-gsl=/usr/local --with-slicot=/usr/local`, adjusting the MATLAB path to accord with your local installation. If you don't have MATLAB, simply remove `--with-matlab=/Applications/MATLAB_R2019a.app` from the above command
- `make -j`
- **(Optional)** To then build mex files for Octave, run
- `cd mex/build/octave`
- `CC=gcc-9 CXX=g++-9 ./configure --with-matio=/usr/local --with-gsl=/usr/local --with-slicot=/usr/local LDFLAGS=-L/usr/local/lib`
- `make -j`
- Prepare the Dynare sources:
```sh
git clone --recurse-submodules https://git.dynare.org/Dynare/dynare.git`
cd dynare
autoreconf -si
```
- Configure Dynare from the source directory:
```sh
./configure --with-matlab=<…> CC=gcc-10 CXX=g++-10 CPPFLAGS=-I/usr/local/include LDFLAGS=-L/usr/local/lib LEX=/usr/local/opt/flex/bin/flex YACC=/usr/local/opt/bison/bin/bison
```
where the path to MATLAB is specified, typically of the form
`/Applications/MATLAB_R2020b.app`. If you don’t have MATLAB, simply replace
`--with-matlab=<…>` in the above command by `--disable-matlab`.
- Compile:
```sh
make -j$(nproc)
```

1
doc/manual/source/bibliography.rst
View File

@ -71,4 +71,5 @@ Bibliography
* Smets, Frank and Rafael Wouters (2003): “An Estimated Dynamic Stochastic General Equilibrium Model of the Euro Area,” *Journal of the European Economic Association*, 1(5), 1123–1175.
* Stock, James H. and Mark W. Watson (1999). “Forecasting Inflation,”, *Journal of Monetary Economics*, 44(2), 293–335.
* Uhlig, Harald (2001): “A Toolkit for Analysing Nonlinear Dynamic Stochastic Models Easily,” in *Computational Methods for the Study of Dynamic Economies*, Eds. Ramon Marimon and Andrew Scott, Oxford University Press, 30–61.
* U.S. Census Bureau (2017): “X-13 ARIMA-SEATSReference Manual”.
* Villemot, Sébastien (2011): “Solving rational expectations models at first order: what Dynare does,” *Dynare Working Papers*, 2, CEPREMAP.

8
doc/manual/source/examples.rst
View File

@ -61,4 +61,10 @@ description, please refer to the comments inside the files themselves.
File demonstrating how to conduct optimal policy experiments in a
simple New Keynesian model either under commitment (Ramsey) or using
optimal simple rules (OSR)
optimal simple rules (OSR)
``Ramsey_steady_file.mod``
File demonstrating how to conduct optimal policy experiments in a
simple New Keynesian model under commitment (Ramsey) with a user-defined
conditional steady state file

10
doc/manual/source/index.rst
View File

@ -3,12 +3,11 @@ The Dynare Reference Manual, version |version|
Currently the development team of Dynare is composed of:
* Stéphane Adjemian (Université du Maine, Gains)
* Stéphane Adjemian (Le Mans Université, Gains)
* Houtan Bastani
* Michel Juillard (Banque de France)
* Sumudu Kankanamge (Toulouse School of Economics)
* Frédéric Karamé (Université du Maine, Gains and CEPREMAP)
* Dóra Kocsis (CEPREMAP)
* Sumudu Kankanamge (Toulouse School of Economics and CEPREMAP)
* Frédéric Karamé (Le Mans Université, Gains and CEPREMAP)
* Junior Maih (Norges Bank)
* Ferhat Mihoubi (Université Paris-Est Créteil, Érudite and CEPREMAP)
* Willi Mutschler (University of Münster)
@ -22,10 +21,11 @@ The following people used to be members of the team:
* Alejandro Buesa
* Fabrice Collard
* Assia Ezzeroug
* Dóra Kocsis
* Stéphane Lhuissier
* George Perendia
Copyright © 1996-2020, Dynare Team.
Copyright © 1996-2021, Dynare Team.
Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.3 or any later version published by the Free Software Foundation; with no Invariant Sections, no Front-Cover Texts, and no Back-Cover Texts.

8
doc/manual/source/introduction.rst
View File

@ -95,8 +95,8 @@ You should cite Dynare if you use it in your research. The
recommended way todo this is to cite the present manual, as:
Stéphane Adjemian, Houtan Bastani, Michel Juillard, Frédéric
Karamé, Junior Maih, Ferhat Mihoubi, George Perendia, Johannes Pfeifer, Marco
Ratto and Sébastien Villemot (2011), “Dynare: Reference Manual,
Karamé, Junior Maih, Ferhat Mihoubi, Willi Mutschler, George Perendia, Johannes Pfeifer,
Marco Ratto and Sébastien Villemot (2011), “Dynare: Reference Manual,
Version 4,” *Dynare Working Papers*, 1, CEPREMAP
For convenience, you can copy and paste the following into your BibTeX file:
@ -106,8 +106,8 @@ For convenience, you can copy and paste the following into your BibTeX file:
@TechReport{Adjemianetal2011,
author = {Adjemian, St\'ephane and Bastani, Houtan and
Juillard, Michel and Karam\'e, Fr\'ederic and
Maih, Junior and Mihoubi, Ferhat and
Perendia, George and Pfeifer, Johannes and
Maih, Junior and Mihoubi, Ferhat and Mutschler, Willi
and Perendia, George and Pfeifer, Johannes and
Ratto, Marco and Villemot, S\'ebastien},
title = {Dynare: Reference Manual Version 4},
year = {2011},

5
doc/manual/source/running-dynare.rst
View File

@ -534,7 +534,10 @@ by the ``dynare`` command.
``M_``,``oo_``, and ``options_`` structures are saved in a file
called ``FILENAME_results.mat``. If they exist, ``estim_params_``,
``bayestopt_``, ``dataset_``, ``oo_recursive_`` and
``estimation_info`` are saved in the same file.
``estimation_info`` are saved in the same file. Note that Matlab
by default only allows ``.mat``-files up to 2GB. You can lift this
restriction by enabling the ``save -v7.3``-option in
``Preferences -> General -> MAT-Files``.
.. matvar:: M_

54
doc/manual/source/the-model-file.rst
View File

@ -1137,9 +1137,8 @@ equations using the ``write_latex_static_model`` command.
* Endogenous variables with leads or lags greater or equal
than two will have been removed, replaced by new auxiliary
variables and equations,
* For a stochastic model, exogenous variables with leads or
lags will also have been replaced by new auxiliary variables
and equations.
* Exogenous variables with leads or lags will also have been replaced
by new auxiliary variables and equations.
For the required LaTeX packages, see
:comm:`write_latex_original_model`.
@ -1221,9 +1220,8 @@ auxiliary endogenous variables—along with corresponding auxiliary
equations—which will appear in the final output.
The main transformation concerns leads and lags. Dynare will perform a
transformation of the model so that there is only one lead and one lag
on endogenous variables and, in the case of a stochastic model, no
leads/lags on exogenous variables.
transformation of the model so that there is only one lead and one lag on
endogenous variables and no leads/lags on exogenous variables.
This transformation is achieved by the creation of auxiliary variables
and corresponding equations. For example, if ``x(+2)`` exists in the
@ -3842,7 +3840,7 @@ Computing the stochastic solution
:mvar:`oo_.conditional_variance_decomposition_ME`). The
variance decomposition is only conducted, if theoretical
moments are requested, *i.e.* using the ``periods=0``-option.
Only available at ``order<3``. In case of ``order=2``,
Only available at ``order<3`` and without ``pruning''. In case of ``order=2``,
Dynare provides a second-order accurate
approximation to the true second moments based on the linear
terms of the second-order solution (see *Kim, Kim,
@ -3859,7 +3857,11 @@ Computing the stochastic solution
algorithm of *Kim, Kim, Schaumburg and Sims (2008)*, while at
third order its generalization by *Andreasen,
Fernández-Villaverde and Rubio-Ramírez (2018)* is used.
Not available above third order.
Not available above third order. When specified, theoretical moments
are based on the pruned state space, i.e. the computation of second moments
uses all terms as in *Andreasen, Fernández-Villaverde and Rubio-Ramírez (2018), page 10*
as opposed to simply providing a second-order accurate result based on the
linear solution as in *Kim, Kim, Schaumburg and Sims (2008)*.
.. option:: partial_information
@ -4261,7 +4263,9 @@ which is described below.
the endogenous variables are generated by assuming that the
agents believe that there will no more shocks after period
:math:`t+S`. This is an experimental feature and can be quite
slow. Default: ``0``.
slow. A non-zero value is not compatible with either the
``bytecode`` or the ``block`` option of the ``model`` block.
Default: ``0``.
.. option:: hybrid
@ -5114,9 +5118,9 @@ block decomposition of the model (see :opt:`block`).
.. option:: mh_replic = INTEGER
Number of replications for Metropolis-Hastings algorithm. For
the time being, ``mh_replic`` should be larger
than 1200. Default: ``20000``.
Number of replications for each chain of the Metropolis-Hastings algorithm.
The number of draws should be sufficient to achieve convergence of the MCMC and
to meaningfully compute posterior objects. Default: ``20000``.
.. option:: sub_draws = INTEGER
@ -5357,6 +5361,14 @@ block decomposition of the model (see :opt:`block`).
(available under MATLAB if the Global Optimization
Toolbox is installed; not available under Octave).
``13``
Uses the ``lsqnonlin`` non-linear least squares optimization
routine (available under MATLAB if the Optimization Toolbox is
installed; available under Octave if the `optim
<https://octave.sourceforge.io/optim/>`__ package from
Octave-Forge is installed). ``101``
``101``
Uses the SolveOpt algorithm for local nonlinear
@ -8950,11 +8962,6 @@ Optimal policy under commitment (Ramsey)
optimal policy for the first time and committing not to
re-optimize in the future.
Because it entails computing at least a second order approximation, the
computation of the planner objective value is skipped with a message when
the model is too large (more than 180 state variables, including lagged
Lagrange multipliers).
.. command:: ramsey_policy [VARIABLE_NAME...];
ramsey_policy (OPTIONS...) [VARIABLE_NAME...];
@ -9039,8 +9046,9 @@ Optimal policy under discretion
under discretion. The algorithm implemented is essentially an LQ
solver, and is described by *Dennis (2007)*.
You should ensure that your model is linear and your objective is
quadratic. Also, you should set the ``linear`` option of the
You must ensure that your objective is quadratic. Regarding the model, it must
either be linear or solved at first order with an analytical steady state provided.
In the first case, you should set the ``linear`` option of the
``model`` block.
It is possible to use the :comm:`estimation` command after the
@ -9842,6 +9850,12 @@ Performing identification analysis
Tolerance level for selecting non-zero singular values.
Default: ``1.e-3``.
.. option:: schur_vec_tol = DOUBLE
Tolerance level used to find nonstationary variables in Schur decomposition
of the transition matrix.
Default: ``1.e-11``.
*Identification Strength Options*
.. option:: no_identification_strength
@ -12212,7 +12226,7 @@ Misc commands
and ``morefloats``.
.. _Dynare wiki: https://archives.dynare.org/DynareWiki/EquationsTags
.. _Dynare wiki: https://archives.dynare.org/DynareWiki/
.. _io: https://octave.sourceforge.io/io/
.. _AIM website: https://www.federalreserve.gov/econres/ama-index.htm

467
doc/manual/source/time-series.rst
View File

@ -10,10 +10,12 @@ Time Series
Dynare provides a MATLAB/Octave class for handling time series data,
which is based on a class for handling dates. Dynare also provides a
new type for dates, so that the basic user does not have to worry
about class and methods for dates. Below, you will first find the
class and methods used for creating and dealing with dates and then
the class used for using time series.
new type for dates, so that the user does not have to worry about
class and methods for dates. Below, you will first find the class and
methods used for creating and dealing with dates and then the class
used for using time series. Dynare also provides an interface to the
X-13 ARIMA-SEATS seasonal adjustment program produced, distributed, and
maintained by the US Census Bureau.
Dates
@ -85,44 +87,42 @@ below. Basic operations can be performed on dates:
Tests if two ``dates`` objects are equal. ``+1950Q1==1950Q2``
returns ``true``, ``1950Q1==1950Q2`` returns ``false``. If the compared
objects have both ``n>1`` elements, the ``eq`` operator returns a
column vector, ``n`` by ``1``, of zeros and ones.
column vector, ``n`` by ``1``, of logicals.
**ne operator (not equal, ~=)**
Tests if two ``dates`` objects are not equal. ``+1950Q1~=``
returns ``false`` while ``1950Q1~=1950Q2`` returns ``true``. If the
compared objects both have ``n>1`` elements, the ``ne`` operator
returns an ``n`` by ``1`` column vector of zeros and ones.
returns an ``n`` by ``1`` column vector of logicals.
**lt operator (less than, <)**
Tests if a ``dates`` object preceeds another ``dates`` object. For
instance, ``1950Q1<1950Q3`` returns ``true``. If the compared objects
have both ``n>1`` elements, the ``lt`` operator returns a column
vector, ``n`` by ``1``, of zeros and ones.
vector, ``n`` by ``1``, of logicals.
**gt operator (greater than, >)**
Tests if a ``dates`` object follows another ``dates`` object. For
instance, ``1950Q1>1950Q3`` returns ``false``. If the compared objects
have both ``n>1`` elements, the ``gt`` operator returns a column
vector, ``n`` by ``1``, of zeros and ones.
vector, ``n`` by ``1``, of logicals.
**le operator (less or equal, <=)**
Tests if a ``dates`` object preceeds another ``dates`` object or
is equal to this object. For instance, ``1950Q1<=1950Q3`` returns
``true``. If the compared objects have both ``n>1`` elements, the
``le`` operator returns a column vector, ``n`` by ``1``, of zeros
and ones.
``le`` operator returns a column vector, ``n`` by ``1``, of logicals.
**ge operator (greater or equal, >=)**
Tests if a ``dates`` object follows another ``dates`` object or is
equal to this object. For instance, ``1950Q1>=1950Q3`` returns
``false``. If the compared objects have both ``n>1`` elements, the
``ge`` operator returns a column vector, ``n`` by ``1``, of zeros
and ones.
``ge`` operator returns a column vector, ``n`` by ``1``, of logicals.
One can select an element, or some elements, in a ``dates`` object as
he would extract some elements from a vector in MATLAB/Octave. Let ``a
@ -164,13 +164,15 @@ The dates class
.. class:: dates
:arg int freq: equal to 1, 4, or 12 (resp. for annual,
quarterly, or monthly dates).
:arg int ndat: the number of declared dates in the object.
:arg int time: a ``ndat*2`` array, the years are stored in the
first column, the subperiods (1 for annual dates,
1-4 for quarterly dates, and 1-12 for monthly
dates) are stored in the second column.
:arg freq: equal to 1, 4, 12 or 365 (resp. for annual, quarterly,
monthly, or daily dates).
:arg time: a ``n*2`` array of integers. If `freq` is equal to 1,
4, or 12, the years are stored in the first column, the
subperiods (1 for annual dates, 1-4 for quarterly
dates, and 1-12 for monthly dates) are stored in the
second column. If `freq` is equal to 365, the first
column stores the number of days since the first day of
year 0, the second column is not used.
Each member is private, one can display the content of a member
but cannot change its value directly. Note that it is not possible
@ -185,14 +187,14 @@ The dates class
|br| Returns an empty ``dates`` object with a given frequency
(if the constructor is called with one input
argument). ``FREQ`` is a character equal to ’Y’ or ’A’ for
annual dates, ’Q’ for quarterly dates, or ’M’ for monthly
dates. Note that ``FREQ`` is not case sensitive, so that, for
instance, ’q’ is also allowed for quarterly dates. The
frequency can also be set with an integer scalar equal to 1
(annual), 4 (quarterly), or 12 (monthly). The instantiation of
empty objects can be used to rename the ``dates`` class. For
instance, if one only works with quarterly dates, object
``qq`` can be created as::
annual dates, ’Q’ for quarterly dates, ’M’ for monthly dates,
or ’D’ for daily dates. Note that ``FREQ`` is not case
sensitive, so that, for instance, ’q’ is also allowed for
quarterly dates. The frequency can also be set with an integer
scalar equal to 1 (annual), 4 (quarterly), 12
(monthly), or 365 (daily). The instantiation of empty objects can be used to
rename the ``dates`` class. For instance, if one only works
with quarterly dates, object ``qq`` can be created as::
qq = dates('Q')
@ -201,7 +203,14 @@ The dates class
d0 = qq(2009,2);
which is much simpler if ``dates`` objects have to be defined
programmatically.
programmatically. For daily dates, we would instantiate an
empty daily dates object as::
dd = dates('D')
and a ``dates`` object holding the date ``2020-12-31``::
d1 = dd(2020,12,31);
.. construct:: dates(STRING)
@ -211,11 +220,11 @@ The dates class
given by the string ``STRING``. This string has to be
interpretable as a date (only strings of the following forms
are admitted: ``'1990Y'``, ``'1990A'``, ``'1990Q1'``,
``'1990M2'``), the routine ``isdate`` can be used to test if a
string is interpretable as a date. If more than one argument
is provided, they should all be dates represented as strings,
the resulting ``dates`` object contains as many elements as
arguments to the constructor.
``'1990M2'``, or ``'2020-12-31'``), the routine ``isdate`` can
be used to test if a string is interpretable as a date. If
more than one argument is provided, they should all be dates
represented as strings, the resulting ``dates`` object
contains as many elements as arguments to the constructor.
.. construct:: dates(DATES)
@ -229,15 +238,16 @@ The dates class
constructor.
.. construct:: dates (FREQ, YEAR, SUBPERIOD)
.. construct:: dates (FREQ, YEAR, SUBPERIOD[, S])
|br| where ``FREQ`` is a single character (’Y’, ’A’, ’Q’, ’M’)
or integer (1, 4, or 12) specifying the frequency, ``YEAR``
and ``SUBPERIOD`` are ``n*1`` vectors of integers. Returns a
``dates`` object with ``n`` elements. If ``FREQ`` is equal to
``'Y'``, ``'A'`` or ``1``, the third argument is not needed
(because ``SUBPERIOD`` is necessarily a vector of ones in this
case).
|br| where ``FREQ`` is a single character (’Y’, ’A’, ’Q’, ’M’,
’D’) or integer (1, 4, 12, or 365) specifying the frequency,
``YEAR`` and ``SUBPERIOD`` and ``S`` are ``n*1`` vectors of
integers. Returns a ``dates`` object with ``n`` elements. The
last argument, ``S``, is only to be used for daily
frequency. If ``FREQ`` is equal to ``'Y'``, ``'A'`` or ``1``,
the third argument is not needed (because ``SUBPERIOD`` is
necessarily a vector of ones in this case).
*Example*
@ -248,6 +258,7 @@ The dates class
do2 = dates('1950Q2','1950Q3');
do3 = dates(do1,do2);
do4 = dates('Q',1950, 1);
do5 = dates('D',1973, 1, 25);
A list of the available methods, by alphabetical order, is given
@ -338,11 +349,11 @@ The dates class
'1950Q1'
.. datesmethod:: C = colon (A, B)
.. datesmethod:: C = colon (A, B)
C = colon (A, i, B)
|br| Overloads the MATLAB/Octave colon (``:``) operator. A and B
are ``dates`` objects. The optional increment ``i`` is a
|br| Overloads the MATLAB/Octave colon (``:``) operator. A and
B are ``dates`` objects. The optional increment ``i`` is a
scalar integer (default value is ``i=1``). This method returns
a ``dates`` object and can be used to create ranges of dates.
@ -366,12 +377,12 @@ The dates class
|br| Returns a copy of a ``dates`` object.
.. datesmethod:: disp (A)
.. datesmethod:: disp (A)
|br| Overloads the MATLAB/Octave disp function for ``dates`` object.
.. datesmethod:: display (A)
.. datesmethod:: display (A)
|br| Overloads the MATLAB/Octave display function for ``dates`` object.
@ -555,7 +566,7 @@ The dates class
1
.. datesmethod:: C = isequal (A, B)
.. datesmethod:: C = isequal (A, B)
|br| Overloads the MATLAB/Octave ``isequal`` function.
@ -616,7 +627,7 @@ The dates class
203
.. datesmethod:: C = lt (A, B)
.. datesmethod:: C = lt (A, B)
|br| Overloads the MATLAB/Octave ``lt`` (less than,
``<``) operator. ``dates`` objects ``A`` and ``B`` must have
@ -861,7 +872,7 @@ The dates class
>> A = dates('1950Q1');
>> A = A:A+1;
>> strings(A)
>> A.strings()
ans =
@ -873,7 +884,8 @@ The dates class
.. datesmethod:: B = subperiod (A)
|br| Returns the subperiod of a date (an integer scalar
between 1 and ``A.freq``).
between 1 and ``A.freq``). This method is not implemented for
daily dates.
*Example*
@ -987,16 +999,14 @@ The dseries class
data. As any MATLAB/Octave statements, this class can be used in a
Dynare’s mod file. A ``dseries`` object has six members:
:arg name: A ``nobs*1`` cell of strings or a ``nobs*p`` character
array, the names of the variables.
:arg tex: A ``nobs*1`` cell of strings or a ``nobs*p`` character
array, the tex names of the variables.
:arg name: A ``vobs*1`` cell of strings or a ``vobs*p`` character array, the names of the variables.
:arg tex: A ``vobs*1`` cell of strings or a ``vobs*p`` character array, the tex names of the variables.
:arg dates dates: An object with ``nobs`` elements, the dates of the sample.
:arg double data: A ``nobs`` by ``vobs`` array, the data.
:arg ops: The history of operations on the variables.
:arg tags: The user-defined tags on the variables.
``data``, ``name``, ``tex`` are private members. The following
``data``, ``name``, ``tex``, and ``ops`` are private members. The following
constructors are available:
.. construct:: dseries ()
@ -1021,8 +1031,6 @@ The dseries class
INIT__ = '1994Q3';
NAMES__ = {'azert';'yuiop'};
TEX__ = {'azert';'yuiop'};
TAGS__ = struct()
DATA__ = {}
azert = randn(100,1);
yuiop = randn(100,1);
@ -1367,7 +1375,7 @@ The dseries class
5Y | 0.66997
.. dseriesmethod:: B = cumprod (A[, d[, v]])
.. dseriesmethod:: B = cumprod (A[, d[, v]])
cumprod_ (A[, d[, v]])
|br| Overloads the MATLAB/Octave ``cumprod`` function for
@ -1510,6 +1518,12 @@ The dseries class
detrended with a different polynomial.
.. dseriesmethod:: B = dgrowth (A)
dgrowth_ (A)
|br| Computes daily growth rates.
.. dseriesmethod:: B = diff (A)
diff_ (A)
@ -1659,6 +1673,13 @@ The dseries class
|br| Returns the first period where all the variables in ``dseries`` object ``A`` are observed (non NaN).
.. dseriesmethod:: B = flip (A)
flip_ (A)
|br| Flips the rows in the data member (without changing the
periods order).
.. dseriesmethod:: f = frequency (B)
|br| Returns the frequency of the variables in ``dseries`` object ``B``.
@ -2028,8 +2049,10 @@ The dseries class
.. dseriesmethod:: B = mdiff (A)
mdiff_ (A)
B = mgrowth (A)
mgrowth_ (A)
|br| Computes monthly growth rates of variables in
|br| Computes monthly differences or growth rates of variables in
``dseries`` object ``A``.
@ -2299,6 +2322,16 @@ The dseries class
otherwise (default) the arithmetic mean is reported.
.. dseriesmethod:: B = nanstd (A[, geometric])
|br| Overloads the MATLAB/Octave ``nanstd`` function for
``dseries`` objects. Returns the standard deviation of each
variable in ``dseries`` object ``A`` ignoring the NaN
values. If the second argument is ``true`` the geometric std
is computed, default value of the second argument is
``false``.
.. dseriesmethod:: C = ne (A, B)
|br| Overloads the MATLAB/Octave ``ne`` (not equal, ``~=``)
@ -2478,6 +2511,39 @@ The dseries class
3Y | 1 | 1
.. dseriesmethod:: A = projection (A, info, periods)
|br| Projects variables in dseries object ``A``. ``info`` is
is a :math:`n \times 3` cell array. Each row provides
informations necessary to project a variable. The first column
contains the name of variable (row char array). the second
column contains the name of the method used to project the
associated variable (row char array), possible values are
``'Trend'``, ``'Constant'``, and ``'AR'``. Last column
provides quantitative information about the projection. If the
second column value is ``'Trend'``, the third column value is
the growth factor of the (exponential) trend. If the second
column value is ``'Constant'``, the third column value is the
level of the variable. If the second column value is ``'AR'``,
the third column value is the autoregressive parameter. The
variables can be projected with an AR(p) model, if the third
column contains a 1×p vector of doubles. The stationarity of
the AR(p) model is not tested. The case of the constant
projection, using the last value of the variable, is covered
with 'Trend' and a growth factor equal to 1, or 'AR' with an
autoregressive parameter equal to one (random walk). This
projection routine only deals with exponential trends.
*Example*
::
>> data = ones(10,4);
>> ts = dseries(data, '1990Q1', {'A1', 'A2', 'A3', 'A4'});
>> info = {'A1', 'Trend', 1.2; 'A2', 'Constant', 0.0; 'A3', 'AR', .5; 'A4', 'AR', [.4, -.2]};
>> ts.projection(info, 10);
.. dseriesmethod:: B = qdiff (A)
B = qgrowth (A)
qdiff_ (A)
@ -2592,6 +2658,44 @@ The dseries class
2Y | 1 | 1 | 1
.. dseriesmethod:: A = resetops (A, ops)
|br| Redefine ``ops`` member.
.. dseriesmethod:: A = resetags (A, ops)
|br| Redefine ``tags`` member.
.. dseriesmethod:: B = round (A[, n])
round_ (A[, n])
|br| Rounds to the nearest decimal or integer. ``n`` is the
precision parameter (number of decimals), default value is 0
meaning that that by default the method rounds to the nearest
integer.
*Example*
::
>> ts = dseries(pi)
ts is a dseries object:
| Variable_1
1Y | 3.1416
>> ts.round_();
>> ts
ts is a dseries object:
| Variable_1
1Y | 3
.. dseriesmethod:: save (A, basename[, format])
|br| Overloads the MATLAB/Octave ``save`` function and saves
@ -2697,6 +2801,28 @@ The dseries class
(default value of the second argument is ``false``).
.. dseriesmethod:: B = subsample (A, d1, d2)
|br| Returns a subsample, for periods between ``dates`` ``d1``
and ``d2``. The same can be achieved by indexing a
``dseries`` object with a ``dates`` object, but the
``subsample`` method is easier to use programmatically.
*Example*
::
>> o = dseries(transpose(1:5));
>> o.subsample(dates('2y'),dates('4y'))
ans is a dseries object:
| Variable_1
2Y | 2
3Y | 3
4Y | 4
.. dseriesmethod:: A = tag (A, a[, b, c])
|br| Add a tag to a variable in ``dseries`` object ``A``.
@ -2795,9 +2921,224 @@ The dseries class
2
.. dseriesmethod:: B = ydiff (A)
B = ygrowth (A)
ydiff_ (A)
ygrowth_ (A)
.. dseriesmethod:: B = ydiff (A)
B = ygrowth (A)
ydiff_ (A)
ygrowth_ (A)
|br| Computes yearly differences or growth rates.
.. _x13-members:
X-13 ARIMA-SEATS interface
==========================
.. class:: x13
|br| The x13 class provides a method for each X-13 command as
documented in the X-13 ARIMA-SEATS reference manual (`x11`,
`automdl`, `estimate`, ...), options can then be passed by
key/value pairs. The ``x13`` class has 22 members:
:arg y: ``dseries`` object with a single variable.
:arg x: ``dseries`` object with an arbitrary number of variables (to be used in the REGRESSION block).
:arg arima: structure containing the options of the ARIMA model command.
:arg automdl: structure containing the options of the ARIMA model selection command.
:arg regression: structure containing the options of the Regression command.
:arg estimate: structure containing the options of the estimation command.
:arg transform: structure containing the options of the transform command.
:arg outlier: structure containing the options of the outlier command.
:arg forecast: structure containing the options of the forecast command.
:arg check: structure containing the options of the check command.
:arg x11: structure containing the options of the X11 command.
:arg force: structure containing the options of the force command.
:arg history: structure containing the options of the history command.
:arg metadata: structure containing the options of the metadata command.
:arg identify: structure containing the options of the identify command.
:arg pickmdl: structure containing the options of the pickmdl command.
:arg seats: structure containing the options of the seats command.
:arg slidingspans: structure containing the options of the slidingspans command.
:arg spectrum: structure containing the options of the spectrum command.
:arg x11regression: structure containing the options of the x11Regression command.
:arg results: structure containing the results returned by x13.
:arg commands: cell array containing the list of commands.
All these members are private. The following constructors are available:
.. construct:: x13 (y)
|br| Instantiates an ``x13`` object with `dseries` object
``y``. The ``dseries`` object passed as an argument must
contain only one variable, the one we need to pass to X-13.
.. construct:: x13 (y, x)
|br| Instantiates an ``x13`` object with `dseries` objects
``y`` and ``x``. The first ``dseries`` object passed as an
argument must contain only one variable, the second
``dseries`` object contains the exogenous variables used by
some of the X-13 commands. Both objects must be defined on the
same time span.
The Following methods allow to set sequence of X-13 commands, write an `.spc` file and run the X-13 binary:
.. x13method:: A = arima (A, key, value[, key, value[, [...]]])
Interface to the ``arima`` command, see the X-13 ARIMA-SEATS
reference manual. All the options must be passed by key/value
pairs.
.. x13method:: A = automdl (A, key, value[, key, value[, [...]]])
Interface to the ``automdl`` command, see the X-13 ARIMA-SEATS
reference manual. All the options must be passed by key/value
pairs.
.. x13method:: A = regression (A, key, value[, key, value[, [...]]])
Interface to the ``regression`` command, see the X-13
ARIMA-SEATS reference manual. All the options must be passed
by key/value pairs.
.. x13method:: A = estimate (A, key, value[, key, value[, [...]]])
Interface to the ``estimate`` command, see the X-13
ARIMA-SEATS reference manual. All the options must be passed
by key/value pairs.
.. x13method:: A = transform (A, key, value[, key, value[, [...]]])
Interface to the ``transform`` command, see the X-13
ARIMA-SEATS reference manual. All the options must be passed
by key/value pairs.
.. x13method:: A = outlier (A, key, value[, key, value[, [...]]])
Interface to the ``outlier`` command, see the X-13 ARIMA-SEATS
reference manual. All the options must be passed by key/value
pairs.
.. x13method:: A = forecast (A, key, value[, key, value[, [...]]])
Interface to the ``forecast`` command, see the X-13
ARIMA-SEATS reference manual. All the options must be passed
by key/value pairs.
.. x13method:: A = check (A, key, value[, key, value[, [...]]])
Interface to the ``check`` command, see the X-13 ARIMA-SEATS
reference manual. All the options must be passed by key/value
pairs.
.. x13method:: A = x11 (A, key, value[, key, value[, [...]]])
Interface to the ``x11`` command, see the X-13 ARIMA-SEATS
reference manual. All the options must be passed by key/value
pairs.
.. x13method:: A = force (A, key, value[, key, value[, [...]]])
Interface to the ``force`` command, see the X-13 ARIMA-SEATS
reference manual. All the options must be passed by key/value
pairs.
.. x13method:: A = history (A, key, value[, key, value[, [...]]])
Interface to the ``history`` command, see the X-13 ARIMA-SEATS
reference manual. All the options must be passed by key/value
pairs.
.. x13method:: A = metadata (A, key, value[, key, value[, [...]]])
Interface to the ``metadata`` command, see the X-13
ARIMA-SEATS reference manual. All the options must be passed
by key/value pairs.
.. x13method:: A = identify (A, key, value[, key, value[, [...]]])
Interface to the ``identify`` command, see the X-13
ARIMA-SEATS reference manual. All the options must be passed
by key/value pairs.
.. x13method:: A = pickmdl (A, key, value[, key, value[, [...]]])
Interface to the ``pickmdl`` command, see the X-13 ARIMA-SEATS
reference manual. All the options must be passed by key/value
pairs.
.. x13method:: A = seats (A, key, value[, key, value[, [...]]])
Interface to the ``seats`` command, see the X-13 ARIMA-SEATS
reference manual. All the options must be passed by key/value
pairs.
.. x13method:: A = slidingspans (A, key, value[, key, value[, [...]]])
Interface to the ``slidingspans`` command, see the X-13
ARIMA-SEATS reference manual. All the options must be passed
by key/value pairs.
.. x13method:: A = spectrum (A, key, value[, key, value[, [...]]])
Interface to the ``spectrum`` command, see the X-13
ARIMA-SEATS reference manual. All the options must be passed
by key/value pairs.
.. x13method:: A = x11regression (A, key, value[, key, value[, [...]]])
Interface to the ``x11regression`` command, see the X-13
ARIMA-SEATS reference manual. All the options must be passed
by key/value pairs.
.. x13method:: print (A[, basefilename])
Prints an ``.spc`` file with all the X-13 commands. The
optional second argument is a row char array specifying the
name (without extension) of the file.
.. x13method:: run (A)
Calls the X-13 binary and run the previously defined
commands. All the results are stored in the structure
``A.results``. When it makes sense these results are saved in
``dseries`` objects (*e.g.* for forecasts or filtered variables).
*Example*
::
>> ts = dseries(rand(100,1),'1999M1');
>> o = x13(ts);
>> o.x11('save','(d11)');
>> o.automdl('savelog','amd','mixed','no');
>> o.outlier('types','all','save','(fts)');
>> o.check('maxlag',24,'save','(acf pcf)');
>> o.estimate('save','(mdl est)');
>> o.forecast('maxlead',18,'probability',0.95,'save','(fct fvr)');
>> o.run();

11
doc/manual/utils/dynare_dom.py
View File

@ -121,6 +121,8 @@ class DynObject(ObjectDescription):
return _('%s (dates method)') % name
elif self.objtype == 'dseriesmethod':
return _('%s (dseries method)') % name
elif self.objtype == 'x13method':
return _('%s (x13 method)') % name
elif self.objtype == 'reportingmethod':
return _('%s (reporting method)') % name
elif self.objtype == 'matcomm':
@ -176,6 +178,10 @@ class DseriesMethod(DynCallable):
display_prefix = 'Method: '
allow_nesting = True
class X13Method(DynCallable):
display_prefix = 'Method: '
allow_nesting = True
class ReportingMethod(DynCallable):
display_prefix = 'Method: '
allow_nesting = True
@ -294,7 +300,8 @@ class DynareDomain(Domain):
'function': ObjType(_('function'), 'func'),
'datesmethod': ObjType(_('method'), 'datmeth'),
'dseriesmethod': ObjType(_('method'), 'dsermeth'),
'reportingmethod': ObjType(_('method'), 'repmeth'),
'x13method': ObjType(_('method'), 'x13meth'),
'reportingmethod': ObjType(_('method'), 'repmeth'),
'matcomm': ObjType(_('matlab command'), 'mcomm'),
'command': ObjType(_('command'), 'comm'),
'class': ObjType(_('class'), 'class'),
@ -312,6 +319,7 @@ class DynareDomain(Domain):
'function': DynFunction,
'datesmethod': DatesMethod,
'dseriesmethod': DseriesMethod,
'x13method': X13Method,
'reportingmethod': ReportingMethod,
'matcomm': MatComm,
'command': DynComm,
@ -330,6 +338,7 @@ class DynareDomain(Domain):
'func': DynareXRefRole(),
'datmeth': DynareXRefRole(),
'dsermeth': DynareXRefRole(),
'x13meth': DynareXRefRole(),
'repmeth': DynareXRefRole(),
'mcomm': DynareXRefRole(),
'comm': DynareXRefRole(),

2
doc/parallel/parallel.tex
View File

@ -248,7 +248,7 @@ The configuration file is designed as follows:
The list of slave options includes:
\begin{description}
\item[Name]: name of the node;
\item[CPUnbr]: this is the number of CPU's to be used on that computer; if \verb"CPUnbr" is a vector of integers, the syntax is \verb"[s:d]", with \verb"d>=s" (\verb"d, s" are integer); the first core has number 1 so that, on a quad-core, use \verb"4" to use all cores, but use \verb [3:4] to specify just the last two cores (this is particularly relevant for Windows where it is possible to assign jobs to specific processors);
\item[CPUnbr]: this is the number of CPU's to be used on that computer; if \verb"CPUnbr" is a vector of integers, the syntax is \verb"[s:d]", with \verb"d>=s" (\verb"d, s" are integer); the first core has number 1 so that, on a quad-core, use \verb"4" to use all cores, but use \verb"[3:4]" to specify just the last two cores (this is particularly relevant for Windows where it is possible to assign jobs to specific processors);
\item[ComputerName]: Computer name on the network or IP address; use the NETBIOS name under Windows\footnote{In Windows XP it is possible find this name in 'My Computer' $->$ mouse right click $->$ 'Property' $->$ 'Computer Name'.}, or the DNS name under Unix.;
\item[UserName]: required for remote login; in order to assure proper communications between the master and the slave threads, it must be the same user name actually logged on the `master' machine. On a Windows network, this is in the form \verb"DOMAIN\username", like \verb"DEPT\JohnSmith", i.e. user JohnSmith in windows group DEPT;
\item[Password]: required for remote login (only under Windows): it is the user password on \verb"DOMAIN" and \verb"ComputerName";

16
dynare++/kord/journal.cc
View File

@ -1,6 +1,6 @@
/*
* Copyright © 2004-2011 Ondra Kamenik
* Copyright © 2019 Dynare Team
* Copyright © 2019-2020 Dynare Team
*
* This file is part of Dynare.
*
@ -203,8 +203,8 @@ Journal::printHeader()
*this << "Dynare++ v. " << VERSION << '\n'
<< '\n'
<< u8"Copyright © 2004-2011 Ondra Kamenik\n"
<< u8"Copyright © 2019 Dynare Team\n"
<< "Dynare++ comes with ABSOLUTELY NO WARRANTY and is distributed under the GNU GPL,"
<< u8"Copyright © 2019-2020 Dynare Team\n"
<< "Dynare++ comes with ABSOLUTELY NO WARRANTY and is distributed under the GNU GPL,\n"
<< "version 3 or later (see https://www.gnu.org/licenses/gpl.html)\n"
<< "\n\n"
<< "System info: ";
@ -219,7 +219,15 @@ Journal::printHeader()
*this << ", processors online: " << std::thread::hardware_concurrency()
<< "\n\nStart time: ";
std::time_t t = std::time(nullptr);
*this << std::put_time(std::localtime(&t), "%c %Z")
// NB: in the date/time string, we avoid using locale-specific strings (#1751)
*this << std::put_time(std::localtime(&t),
#ifndef _WIN32
"%F %T %z"
#else
// Specifiers introduced in C++11 don’t work under Windows…
"%Y-%m-%d %H:%M:%S"
#endif
)
<< "\n\n"
<< u8" ┌────╼ elapsed time (seconds) \n"
<< u8" │ ┌────╼ record unique identifier \n"

12
dynare++/kord/journal.hh
View File

@ -1,6 +1,6 @@
/*
* Copyright © 2004 Ondra Kamenik
* Copyright © 2019 Dynare Team
* Copyright © 2019-2020 Dynare Team
*
* This file is part of Dynare.
*
@ -57,14 +57,18 @@ struct SystemResources
class Journal : public std::ofstream
{
int ord;
int depth;
int ord{0};
int depth{0};
public:
explicit Journal(const std::string &fname)
: std::ofstream(fname), ord(0), depth(0)
: std::ofstream(fname)
{
printHeader();
}
/* Constructor that does not initialize the std::ofstream. To be used when an
on-disk journal is not wanted. */
Journal() = default;
Journal &operator=(Journal &&) = default;
~Journal() override
{
flush();

4
dynare++/src/main.cc
View File

@ -1,6 +1,6 @@
/*
* Copyright © 2004-2011 Ondra Kamenik
* Copyright © 2019 Dynare Team
* Copyright © 2019-2020 Dynare Team
*
* This file is part of Dynare.
*
@ -47,7 +47,7 @@ main(int argc, char **argv)
std::cout << "Dynare++ v. " << VERSION << '\n'
<< '\n'
<< u8"Copyright © 2004-2011 Ondra Kamenik\n"
<< u8"Copyright © 2019 Dynare Team\n"
<< u8"Copyright © 2019-2020 Dynare Team\n"
<< "Dynare++ comes with ABSOLUTELY NO WARRANTY and is distributed under the GNU GPL,\n"
<< "version 3 or later (see https://www.gnu.org/licenses/gpl.html)\n";
return EXIT_SUCCESS;

119
examples/Ramsey_steady_file.mod Normal file
View File

@ -0,0 +1,119 @@
/*
* This file replicates the model studied in:
* Lawrence J. Christiano, Roberto Motto and Massimo Rostagno (2007):
* "Notes on Ramsey-Optimal Monetary Policy", Section 2
* The paper is available at http://faculty.wcas.northwestern.edu/~lchrist/d16/d1606/ramsey.pdf
*
* Notes:
* - This mod-files allows to simulate a simple New Keynesian Model with Rotemberg price
* adjustment costs under fully optimal monetary under commitment (Ramsey)
*
* - This files shows how to use a userd-defined conditional steady state file in the Ramsey case. It takes
* the value of the defined instrument R as given and then computes the rest of the steady
* state, including the steady state inflation rate, based on this value. The initial value
* of the instrument for steady state search must then be defined in an initval-block.
*
* This implementation was written by Johannes Pfeifer.
*
* If you spot mistakes, email me at jpfeifer@gmx.de
*
* Please note that the following copyright notice only applies to this Dynare
* implementation of the model.
*/
/*
* Copyright (C) 2019 Dynare Team
*
* This is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* It is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* For a copy of the GNU General Public License,
* see <http://www.gnu.org/licenses/>.
*/
var C $C$ (long_name='Consumption')
pi $\pi$ (long_name='Gross inflation')
h $h$ (long_name='hours worked')
Z $Z$ (long_name='TFP')
R $R$ (long_name='Net nominal interest rate')
log_C ${\ln C}$ (long_name='Log Consumption')
log_h ${\ln h}$ (long_name='Log hours worked')
pi_ann ${\pi^{ann}}$ (long_name='Annualized net inflation')
R_ann ${R^{ann}}$ (long_name='Annualized net nominal interest rate')
r_real ${r^{ann,real}}$ (long_name='Annualized net real interest rate')
y_nat ${y^{nat}}$ (long_name='Natural (flex price) output')
y_gap ${r^{gap}}$ (long_name='Output gap')
;
varexo epsilon ${\varepsilon}$ (long_name='TFP shock')
;
parameters beta ${\beta}$ (long_name='discount factor')
theta ${\theta}$ (long_name='substitution elasticity')
tau ${\tau}$ (long_name='labor subsidy')
chi ${\chi}$ (long_name='labor disutility')
phi ${\phi}$ (long_name='price adjustment costs')
rho ${\rho}$ (long_name='TFP autocorrelation')
;
beta=0.99;
theta=5;
phi=100;
rho=0.9;
tau=0;
chi=1;
model;
[name='Euler equation']
1/(1+R)=beta*C/(C(+1)*pi(+1));
[name='Firm FOC']
(tau-1/(theta-1))*(1-theta)+theta*(chi*h*C/(exp(Z))-1)=phi*(pi-1)*pi-beta*phi*(pi(+1)-1)*pi(+1);
[name='Resource constraint']
C*(1+phi/2*(pi-1)^2)=exp(Z)*h;
[name='TFP process']
Z=rho*Z(-1)+epsilon;
[name='Definition log consumption']
log_C=log(C);
[name='Definition log hours worked']
log_h=log(h);
[name='Definition annualized inflation rate']
pi_ann=4*log(pi);
[name='Definition annualized nominal interest rate']
R_ann=4*R;
[name='Definition annualized real interest rate']
r_real=4*log((1+R)/pi(+1));
[name='Definition natural output']
y_nat=exp(Z)*sqrt((theta-1)/theta*(1+tau)/chi);
[name='output gap']
y_gap=log_C-log(y_nat);
end;
initval;
R=1/beta-1;
end;
shocks;
var epsilon = 0.01^2;
end;
//use Ramsey optimal policy
//define planner objective, which corresponds to utility function of agents
planner_objective log(C)-chi/2*h^2;
//set up Ramsey optimal policy problem with interest rate R as the instrument,...
// defining the discount factor in the planner objective to be the one of private agents
ramsey_model(instruments=(R),planner_discount=beta,planner_discount_latex_name=$\beta$);
//conduct stochastic simulations of the Ramsey problem
stoch_simul(order=1,irf=20,periods=500) pi_ann log_h R_ann log_C Z r_real;
evaluate_planner_objective;

79
examples/Ramsey_steady_file_steadystate.m Normal file
View File

@ -0,0 +1,79 @@
function [ys,params,check] = Ramsey_steady_file_steadystate(ys,exo,M_,options_)
% function [ys,params,check] = Ramsey_steady_file_steadystate(ys,exo,M_,options_)
% computes the steady state for the Ramsey_steady_file.mod, conditional on
% the instrument value provided
%
% Inputs:
% - ys [vector] vector of initial values for the steady state of
% the endogenous variables
% - exo [vector] vector of values for the exogenous variables
% - M_ [structure] Dynare model structure
% - options [structure] Dynare options structure
%
% Output:
% - ys [vector] vector of steady state values for the the endogenous variables
% - params [vector] vector of parameter values
% - check [scalar] set to 0 if steady state computation worked and to
% 1 of not (allows to impose restrictions on parameters)
% Copyright (C) 2020 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
% read out parameters to access them with their name
beta=NaN; %make parameter known to Matlab function, prevents crashes due to Matlab function with same name;
%will be overwritten next
NumberOfParameters = M_.param_nbr;
for ii = 1:NumberOfParameters
paramname = M_.param_names{ii};
eval([ paramname ' = M_.params(' int2str(ii) ');']);
end
% read in instrument values
for ii = 1:size(options_.instruments,1)
eval([options_.instruments{ii} ' = ys(strmatch(options_.instruments{ii},M_.endo_names,''exact'')) ;']);
end
% initialize indicator
check = 0;
%% Enter model equations here
Z=0;
pi=(R+1)*beta;
C=sqrt((1+1/theta*((1-beta)*(pi-1)*pi-(tau-1/(theta-1))*(1-theta)))/(chi*(1+phi/2*(pi-1)^2)));
h=C*(1+phi/2*(pi-1)^2);
log_C=log(C);
log_h=log(h);
pi_ann=4*log(pi);
R_ann=4*R;
r_real=4*log((1+R)/pi);
y_nat=sqrt((theta-1)/theta*(1+tau)/chi);
y_gap=log_C-log(y_nat);
%% end own model equations
params=NaN(NumberOfParameters,1);
for iter = 1:length(M_.params) %update parameters set in the file
eval([ 'params(' num2str(iter) ') = ' M_.param_names{iter} ';' ])
end
NumberOfEndogenousVariables = M_.orig_endo_nbr; %auxiliary variables are set automatically
for ii = 1:NumberOfEndogenousVariables
varname = M_.endo_names{ii};
eval(['ys(' int2str(ii) ') = ' varname ';']);
end

6
license.txt
View File

@ -1,11 +1,10 @@
Format: https://www.debian.org/doc/packaging-manuals/copyright-format/1.0/
Upstream-Name: Dynare
Upstream-Contact: Dynare Team, whose members in 2020 are:
Upstream-Contact: Dynare Team, whose members in 2021 are:
- Stéphane Adjemian <stephane.adjemian@univ-lemans.fr>
- Houtan Bastani <houtan@dynare.org>
- Michel Juillard <michel.juillard@mjui.fr>
- Frédéric Karamé <frederic.karame@univ-lemans.fr>
- Dóra Kocsis <dora@dynare.org>
- Junior Maih <junior.maih@gmail.com>
- Ferhat Mihoubi <fmihoubi@univ-evry.fr>
- Willi Mutschler <willi@mutschler.eu>
@ -17,12 +16,13 @@ Upstream-Contact: Dynare Team, whose members in 2020 are:
- Alejandro Buesa
- Fabrice Collard
- Assia Ezzeroug
- Dóra Kocsis
- Stéphane Lhuissier
- George Perendia
Source: https://www.dynare.org
Files: *
Copyright: 1996-2020 Dynare Team
Copyright: 1996-2021 Dynare Team
License: GPL-3+
Files: matlab/AIM/SP*

7
m4/ax_mexopts.m4
View File

@ -64,7 +64,7 @@ case ${MATLAB_ARCH} in
OPENMP_LDFLAGS="-Wl,-Bstatic,--whole-archive -lgomp -Wl,-Bdynamic,--no-whole-archive"
ax_mexopts_ok="yes"
;;
maci | maci64)
maci64)
MACOSX_DEPLOYMENT_TARGET='10.9'
MATLAB_DEFS="$MATLAB_DEFS -DNDEBUG"
MATLAB_CFLAGS="-fno-common -arch x86_64 -mmacosx-version-min=$MACOSX_DEPLOYMENT_TARGET -fexceptions"
@ -72,7 +72,10 @@ case ${MATLAB_ARCH} in
MATLAB_FCFLAGS="-g -O2 -fexceptions -fbackslash -arch x86_64"
MATLAB_LDFLAGS_NOMAP="-Wl,-twolevel_namespace -undefined error -arch x86_64 -mmacosx-version-min=$MACOSX_DEPLOYMENT_TARGET -bundle"
MATLAB_LDFLAGS="$MATLAB_LDFLAGS_NOMAP -Wl,-exported_symbols_list,\$(abs_top_srcdir)/mexFunction-MacOSX.map"
MATLAB_LIBS="-L$MATLAB/bin/maci64 -lmx -lmex -lmat -lmwlapack -lmwblas"
# This -L flag is put here, hence later on the linker command line, so as
# to avoid linking against the HDF5 shipped by MATLAB (which would
# otherwise override the HDF5 from Homebrew)
MATLAB_LIBS="-L$MATLAB/bin/maci64 -lmx -lmex -lmat -lmwlapack -lmwblas -lstdc++"
ax_mexopts_ok="yes"
;;
*)

2
m4/ax_slicot.m4
View File

@ -43,7 +43,7 @@ AC_DEFUN([AX_SLICOT],
AC_FC_FUNC(sb02od)
if test "$1" = matlab; then
LDFLAGS="$MATLAB_LDFLAGS_NOMAP $LDFLAGS_SLICOT"
LDFLAGS="$LDFLAGS $MATLAB_LDFLAGS_NOMAP $LDFLAGS_SLICOT"
case ${MATLAB_ARCH} in
glnxa64 | win64 | maci64)

29
macOS/build.sh
View File

@ -1,6 +1,6 @@
#!/usr/bin/env bash
# Copyright © 2019-2020 Dynare Team
# Copyright © 2019-2021 Dynare Team
#
# This file is part of Dynare.
#
@ -64,11 +64,15 @@ fi
##
cd "$ROOTDIR"
[[ -f configure ]] || autoreconf -si
CC=$CC CXX=$CXX ./configure \
./configure \
PACKAGE_VERSION="$VERSION" \
PACKAGE_STRING="dynare $VERSION" \
CXXFLAGS=-I/usr/local/include \
CC=$CC \
CXX=$CXX \
CPPFLAGS=-I/usr/local/include \
LDFLAGS=-static-libgcc \
LEX=/usr/local/opt/flex/bin/flex \
YACC=/usr/local/opt/bison/bin/bison \
--with-gsl="$LIB64"/gsl \
--with-matio="$LIB64"/matio \
--with-slicot="$LIB64"/Slicot/with-underscore \
@ -88,6 +92,7 @@ make -j"$NTHREADS"
NAME=dynare-"$VERSION"
PKGFILES="$ROOTDIR"/macOS/pkg/"$NAME"
mkdir -p \
"$PKGFILES"/preprocessor \
"$PKGFILES"/mex/matlab/maci64-8.3-9.3 \
"$PKGFILES"/mex/matlab/maci64-9.4-9.9 \
"$PKGFILES"/mex/octave \
@ -107,6 +112,12 @@ cp -p "$ROOTDIR"/license.txt "$PKGFILES"
cp -pr "$ROOTDIR"/matlab "$PKGFILES"
cp -pr "$ROOTDIR"/examples "$PKGFILES"
cp -p "$ROOTDIR"/preprocessor/src/dynare-preprocessor "$PKGFILES"/preprocessor
# Recreate backward-compatibility symlink
rm -f "$ROOTDIR"/matlab/preprocessor64/dynare_m
ln -sf ../../preprocessor/dynare-preprocessor "$PKGFILES"/matlab/preprocessor64/dynare_m
cp -L "$ROOTDIR"/mex/matlab/* "$PKGFILES"/mex/matlab/maci64-8.3-9.3
cp -p "$ROOTDIR"/scripts/dynare.el "$PKGFILES"/scripts
@ -136,10 +147,12 @@ cp -p "$ROOTDIR"/macOS/deps/lib64/x13as/x13as "$PKGFILES"
##
cd "$ROOTDIR"/mex/build/matlab
make clean
CC=$CC CXX=$CXX ./configure \
./configure \
PACKAGE_VERSION="$VERSION" \
PACKAGE_STRING="dynare $VERSION" \
CXXFLAGS=-I/usr/local/include \
CC=$CC \
CXX=$CXX \
CPPFLAGS=-I/usr/local/include \
LDFLAGS=-static-libgcc \
--with-gsl="$LIB64"/gsl \
--with-matio="$LIB64"/matio \
@ -156,10 +169,12 @@ cd "$ROOTDIR"/mex/build/octave
OCTAVE_VERSION=$(grep OCTAVE_VERSION "$ROOTDIR"/macOS/deps/versions.mk | cut -d'=' -f2 | sed -e 's/^[[:space:]]*//' -e 's/[[:space:]]*$//')
OCTAVE_USR_DIR="/Applications/Octave-$OCTAVE_VERSION.app/Contents/Resources/usr"
OCTAVE_BIN_DIR="$OCTAVE_USR_DIR/Cellar/octave-octave-app@$OCTAVE_VERSION/$OCTAVE_VERSION/bin"
PATH="$OCTAVE_BIN_DIR:$PATH" CC=$CC CXX=$CXX ./configure \
PATH="$OCTAVE_BIN_DIR:$PATH" ./configure \
PACKAGE_VERSION="$VERSION" \
PACKAGE_STRING="dynare $VERSION" \
CXXFLAGS=-I/usr/local/include \
CC=$CC \
CXX=$CXX \
CPPFLAGS=-I/usr/local/include \
LDFLAGS="-static-libgcc -L$OCTAVE_USR_DIR/lib " \
--with-gsl="$LIB64"/gsl \
--with-matio="$LIB64"/matio \

12
matlab/UnivariateSpectralDensity.m
View File

@ -98,7 +98,7 @@ end
[vx, u] = lyapunov_symm(A,B*M_.Sigma_e*B',options_.lyapunov_fixed_point_tol,options_.qz_criterium,options_.lyapunov_complex_threshold,[],options_.debug);
iky = iv(ivar);
if ~isempty(u)
iky = iky(find(any(abs(ghx(iky,:)*u) < options_.Schur_vec_tol,2)));
iky = iky(find(any(abs(ghx(iky,:)*u) < options_.schur_vec_tol,2)));
ivar = oo_.dr.order_var(iky);
end
@ -151,11 +151,11 @@ else
end
if ~options_.nograph
if ~exist(M_.fname, 'dir')
mkdir('.',M_.fname);
if ~exist(M_.dname, 'dir')
mkdir('.',M_.dname);
end
if ~exist([M_.fname '/graphs'],'dir')
mkdir(M_.fname,'graphs');
if ~exist([M_.dname '/graphs'],'dir')
mkdir(M_.dname,'graphs');
end
for i= 1:nvar
@ -165,6 +165,6 @@ if ~options_.nograph
ylabel('f(\omega)')
box on
axis tight
dyn_saveas(hh,[M_.fname ,filesep,'graphs', filesep, 'SpectralDensity_' M_.endo_names{ivar(i)}],options_.nodisplay,options_.graph_format)
dyn_saveas(hh,[M_.dname ,filesep,'graphs', filesep, 'SpectralDensity_' M_.endo_names{ivar(i)}],options_.nodisplay,options_.graph_format)
end
end

74
matlab/accessors/get_posterior_covariance.m Normal file
View File

@ -0,0 +1,74 @@
function r = get_posterior_covariance(endo1, endo2)
% Copyright (C) 2020 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
global M_ options_
id1 = find(strcmp(endo1, M_.endo_names));
if nargin < 2
id2 = id1;
else
id2 = find(strcmp(endo2, M_.endo_names));
end
ivar = get_variables_list(options_,M_);
if ~ismember(id1, ivar)
error('Posterior moments have not been computed for %s', endo1)
end
if ~ismember(id2, ivar)
error('Posterior moments have not been computed for %s', endo2)
end
r = struct('draw', [], 'moment', []);
% Fetch parameter draws
NumberDrawsFiles = length(dir([M_.dname '/metropolis/' M_.fname '_posterior_draws*' ]));
if NumberDrawsFiles == 0
error('Can''t find posterior draws file(s)')
end
idx = 1;
for file = 1:NumberDrawsFiles
load([M_.dname '/metropolis/' M_.fname '_posterior_draws' int2str(file) ],'pdraws');
for i = 1:size(pdraws, 1)
r(idx).draw = pdraws{i, 1};
idx = idx+1;
end
end
% Fetch moments
symidx = symmetric_matrix_index(find(ivar==id1), find(ivar==id2), length(ivar));
NumberMomentsFiles = length(dir([M_.dname '/metropolis/' M_.fname '_Posterior2ndOrderMoments*']));
if NumberMomentsFiles == 0
error('Can''t find posterior 2nd order moments file(s)')
end
idx = 1;
for file = 1:NumberMomentsFiles
load([M_.dname '/metropolis/' M_.fname '_Posterior2ndOrderMoments' int2str(file) ],'Covariance_matrix');
for i = 1:size(Covariance_matrix, 1)
r(idx).moment = Covariance_matrix(i, symidx);
idx = idx+1;
end
end

74
matlab/accessors/get_posterior_irf.m Normal file
View File

@ -0,0 +1,74 @@
function r = get_posterior_irf(endo, exo)
% Copyright (C) 2020 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
global M_ options_
endo_id = find(strcmp(endo, M_.endo_names));
exo_id = find(strcmp(exo, M_.exo_names));
if isempty(exo_id)
error('Unknown exogenous: %s', exo)
end
ivar = get_variables_list(options_,M_);
if ~ismember(endo_id, ivar)
error('Posterior IRFs have not been computed for %s', endo)
end
r = struct('draw', [], 'irf', []);
% Fetch parameter draws
NumberDrawsFiles = length(dir([M_.dname '/metropolis/' M_.fname '_param_irf*' ]));
if NumberDrawsFiles == 0
error('Can''t find posterior draws file(s)')
end
idx = 1;
for file = 1:NumberDrawsFiles
load([M_.dname '/metropolis/' M_.fname '_param_irf' int2str(file) ],'stock');
for i = 1:size(stock, 1)
r(idx).draw = stock(i, :);
idx = idx+1;
end
end
% Fetch IRFs
% We use the raw files, i.e. those that are *not* the output of
% ReshapeMatFiles.m. Their filename contains 'irf_dsge' in *lowercase*.
filesUpperAndLower = dir([M_.dname '/metropolis/' M_.fname '_irf_dsge*']);
filesLower = cellfun(@any, regexp({filesUpperAndLower.name}, [M_.fname '_irf_dsge\d\.mat']));
NumberIRFsFiles = sum(filesLower);
if NumberIRFsFiles == 0
error('Can''t find posterior IRFs file(s)')
end
endo_id_varlist = find(ivar == endo_id);
idx = 1;
for file = 1:NumberIRFsFiles
load([M_.dname '/metropolis/' M_.fname '_irf_dsge' int2str(file) ],'stock_irf_dsge');
for i = 1:size(stock_irf_dsge, 4)
r(idx).irf = stock_irf_dsge(:, endo_id_varlist, exo_id, i);
idx = idx+1;
end
end

110
matlab/accessors/get_posterior_timeseries.m Normal file
View File

@ -0,0 +1,110 @@
function r = get_posterior_timeseries(type, endo)
% Returns the posterior distribution of either
% smoothed/updated/filtered/forecast timeseries for a given endogenous variable.
% "type" must be either 'smoothed', 'updated', 'filtered' or 'forecast'.
%
% For filtered variables, returns a matrix with step-ahead in lines and periods
% in columns.
% For forecasts, returns the "point" forecast (i.e. with uncertainty about both
% parameters and shocks).
% Copyright (C) 2020 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
global M_
endo_id = find(strcmp(endo, M_.endo_names));
if isempty(endo_id)
error('Unknown endogenous variable %s', endo)
end
r = struct('draw', [], type, []);
% Fetch parameter draws
DrawsFiles = dir([M_.dname '/metropolis/' M_.fname '_param*' ]);
NumberDrawsFiles = 0;
for i=1:length(DrawsFiles)
% We need to filter out the _param_irf* files
if ~isempty(regexp(DrawsFiles(i).name, '.*_param[0-9]+\.mat'))
NumberDrawsFiles = NumberDrawsFiles + 1;
end
end
if NumberDrawsFiles == 0
error('Can''t find posterior draws file(s)')
end
idx = 1;
for file = 1:NumberDrawsFiles
load([M_.dname '/metropolis/' M_.fname '_param' int2str(file) ],'stock');
for i = 1:size(stock, 1)
r(idx).draw = stock(i, :);
idx = idx+1;
end
end
% Fetch timeseries
switch type
case 'smoothed'
basename = 'smooth';
case 'updated'
basename = 'update';
case 'filtered'
basename = 'filter_step_ahead';
case 'forecast'
basename = 'forc_point';
otherwise
error('Unknown type requested. Should be one of: smoothed, updated, filtered, forecast')
end
if strcmp(type, 'smoothed')
SmoothedFiles = dir([M_.dname '/metropolis/' M_.fname '_smooth*']);
NumberTimeseriesFiles = 0;
for i=1:length(SmoothedFiles)
%% We need to filter out the _smoothed_{constant,trend}* files
if ~isempty(regexp(SmoothedFiles(i).name, '.*_smooth[0-9]+\.mat'))
NumberTimeseriesFiles = NumberTimeseriesFiles + 1;
end
end
else
NumberTimeseriesFiles = length(dir([M_.dname '/metropolis/' M_.fname '_' ...
basename '*']));
end
if NumberTimeseriesFiles == 0
error('Can''t find file(s) with posterior timeseries of requested type')
end
idx = 1;
for file = 1:NumberTimeseriesFiles
load([M_.dname '/metropolis/' M_.fname '_' basename int2str(file) ],'stock');
if strcmp(type, 'filtered')
for i = 1:size(stock, 4)
r(idx).(type) = reshape(stock(:, endo_id, :, i), size(stock,1), size(stock,3));
idx = idx+1;
end
else
for i = 1:size(stock, 3)
r(idx).(type) = stock(endo_id, :, i)';
idx = idx+1;
end
end
end

8
matlab/bvar_forecast.m
View File

@ -114,7 +114,7 @@ sims_with_shocks_sort = sort(sims_with_shocks, 3);
sims_with_shocks_down_conf = sims_with_shocks_sort(:, :, sort_idx(1));
sims_with_shocks_up_conf = sims_with_shocks_sort(:, :, sort_idx(2));
OutputDirectoryName = CheckPath('graphs',M_.fname);
OutputDirectoryName = CheckPath('graphs',M_.dname);
dyn_graph=dynare_graph_init(sprintf('BVAR forecasts (nlags = %d)', nlags), ny, {'b-' 'g-' 'g-' 'r-' 'r-'});
@ -153,10 +153,10 @@ end
% Store results
DirectoryName = [ M_.fname '/bvar_forecast' ];
DirectoryName = [ M_.dname '/bvar_forecast' ];
if ~isdir(DirectoryName)
if ~isdir(M_.fname)
mkdir(M_.fname);
if ~isdir(M_.dname)
mkdir(M_.dname);
end
mkdir(DirectoryName);
end

2
matlab/bvar_irf.m
View File

@ -126,7 +126,7 @@ for shock=1:ny
end
% Save intermediate results
DirectoryName = [ M_.fname '/bvar_irf' ];
DirectoryName = [ M_.dname '/bvar_irf' ];
if ~isdir(DirectoryName)
mkdir('.',DirectoryName);
end

4
matlab/check_posterior_sampler_options.m
View File

@ -378,11 +378,13 @@ end
% here are all samplers requiring a proposal distribution
if ~strcmp(posterior_sampler_options.posterior_sampling_method,'slice')
if ~options_.cova_compute && ~(options_.load_mh_file && posterior_sampler_options.use_mh_covariance_matrix)
if ~options_.cova_compute && ~(options_.load_mh_file && posterior_sampler_options.use_mh_covariance_matrix)
if strcmp('hessian',options_.MCMC_jumping_covariance)
skipline()
disp('check_posterior_sampler_options:: I cannot start the MCMC because the Hessian of the posterior kernel at the mode was not computed')
disp('check_posterior_sampler_options:: or there is no previous MCMC to load ')
error('check_posterior_sampler_options:: MCMC cannot start')
end
end
end

2
matlab/compute_moments_varendo.m
View File

@ -151,7 +151,7 @@ if M_.exo_nbr > 1
end
skipline();
if ~all(diag(M_.H)==0)
if isoctave
if isoctave && octave_ver_less_than('6')
[observable_name_requested_vars, varlist_pos] = intersect_stable(var_list_, options_.varobs);
else
[observable_name_requested_vars, varlist_pos] = intersect(var_list_, options_.varobs, 'stable');

2
matlab/conditional_variance_decomposition.m
View File

@ -89,7 +89,7 @@ end
% Measurement error
if ~all(diag(StateSpaceModel.measurement_error)==0)
if isoctave
if isoctave && octave_ver_less_than('6')
[observable_pos,index_subset,index_observables]=intersect_stable(SubsetOfVariables,StateSpaceModel.observable_pos);
else
[observable_pos,index_subset,index_observables]=intersect(SubsetOfVariables,StateSpaceModel.observable_pos,'stable');

2
matlab/conditional_variance_decomposition_ME_mc_analysis.m
View File

@ -63,7 +63,7 @@ if isempty(exogenous_variable_index)
end
end
if isoctave
if isoctave && octave_ver_less_than('6')
[observable_pos_requested_vars,index_subset,index_observables]=intersect_stable(var_list,options_.varobs);
else
[observable_pos_requested_vars,index_subset,index_observables]=intersect(var_list,options_.varobs,'stable');

2
matlab/default_option_values.m
View File

@ -46,7 +46,7 @@ options_.gstep(2) = 1.0;
options_.scalv = 1;
options_.debug = false;
options_.initval_file = false;
options_.Schur_vec_tol = 1e-11; % used to find nonstationary variables in Schur decomposition of the
options_.schur_vec_tol = 1e-11; % used to find nonstationary variables in Schur decomposition of the
% transition matrix
options_.qz_criterium = [];
options_.qz_zero_threshold = 1e-6;

5
matlab/discretionary_policy/discretionary_policy.m
View File

@ -36,15 +36,16 @@ options_.discretionary_policy = 1;
options_.order = 1;
[info, oo_, options_, M_] = stoch_simul(M_, options_, oo_, var_list);
oo_.steady_state = oo_.dr.ys;
if ~options_.noprint
disp_steady_state(M_,oo_)
disp_steady_state(M_,oo_,options_)
for i=M_.orig_endo_nbr:M_.endo_nbr
if strmatch('mult_', M_.endo_names{i})
fprintf('%s \t\t %g\n', M_.endo_names{i}, oo_.dr.ys(i));
end
end
end
oo_.planner_objective_value = evaluate_planner_objective(M_,options_,oo_);
options_.order = origorder;

26
matlab/discretionary_policy/discretionary_policy_1.m
View File

@ -41,11 +41,13 @@ beta = get_optimal_policy_discount_factor(M_.params, M_.param_names);
%call steady_state_file if present to update parameters
if options_.steadystate_flag
% explicit steady state file
[~,M_.params,info] = evaluate_steady_state_file(oo_.steady_state,[oo_.exo_steady_state; oo_.exo_det_steady_state],M_, ...
[ys,M_.params,info] = evaluate_steady_state_file(oo_.steady_state,[oo_.exo_steady_state; oo_.exo_det_steady_state],M_, ...
options_,false);
if info(1)
return;
end
else
ys=zeros(M_.endo_nbr,1);
end
[U,Uy,W] = feval([M_.fname,'.objective.static'],zeros(M_.endo_nbr,1),[], M_.params);
if any(any(isnan(Uy)))
@ -73,8 +75,10 @@ W=reshape(W,M_.endo_nbr,M_.endo_nbr);
klen = M_.maximum_lag + M_.maximum_lead + 1;
iyv=M_.lead_lag_incidence';
% Find the jacobian
z = repmat(zeros(M_.endo_nbr,1),1,klen);
z = z(nonzeros(iyv)) ;
z = repmat(ys,1,klen);
iyr0 = find(iyv(:)) ;
z = z(iyr0);
it_ = M_.maximum_lag + 1 ;
if M_.exo_nbr == 0
@ -82,16 +86,16 @@ if M_.exo_nbr == 0
end
[junk,jacobia_] = feval([M_.fname '.dynamic'],z, [zeros(size(oo_.exo_simul)) ...
oo_.exo_det_simul], M_.params, zeros(M_.endo_nbr,1), it_);
if any(junk~=0)
info = 65; %the model must be written in deviation form and not have constant terms
return;
oo_.exo_det_simul], M_.params, ys, it_);
if max(abs(junk))>options_.solve_tolf
info = 65; %the model must be written in deviation form and not have constant terms or have a steady state provided
return;
end
Indices={'lag','contemp','lead'};
iter=1;
for j=1:numel(Indices)
A.(Indices{j})=zeros(M_.orig_eq_nbr,M_.endo_nbr);
A.(Indices{j})=zeros(M_.eq_nbr,M_.endo_nbr);
if strcmp(Indices{j},'contemp')||(strcmp(Indices{j},'lag') && M_.maximum_lag)||(strcmp(Indices{j},'lead') && M_.maximum_lead)
[~,row,col]=find(M_.lead_lag_incidence(iter,:));
A.(Indices{j})(:,row)=jacobia_(:,col);
@ -116,10 +120,12 @@ else
end
%write back solution to dr
dr.ys =zeros(M_.endo_nbr,1);
dr.ys =ys;
dr=set_state_space(dr,M_,options_);
T=H(dr.order_var,dr.order_var);
dr.ghu=G(dr.order_var,:);
Selection=M_.lead_lag_incidence(1,dr.order_var)>0;%select state variables
if M_.maximum_endo_lag
Selection=M_.lead_lag_incidence(1,dr.order_var)>0;%select state variables
end
dr.ghx=T(:,Selection);
oo_.dr = dr;

22
matlab/disp_dr.m
View File

@ -67,6 +67,9 @@ header_label_length=16; %default
for ii=1:length(ivar)
header_label_length = max(header_label_length,length(M_.endo_names{k1(ivar(ii))})+2);
end
if options_.loglinear
header_label_length=header_label_length+5;
end
header_label_format = sprintf('%%%ds',header_label_length);
value_format_float = sprintf('%%%d.6f',header_label_length);
value_format_zero = sprintf('%%%dd',header_label_length);
@ -83,8 +86,11 @@ else
aux_var_additional_characters=0;
end
var_name_width = max([cellofchararraymaxlength(M_.endo_names(k1(ivar))), cellofchararraymaxlength(M_.exo_names)]);
if options_.loglinear
var_name_width = max([cellofchararraymaxlength(M_.endo_names(k1(ivar)))+5, cellofchararraymaxlength(M_.exo_names)]);
else
var_name_width = max([cellofchararraymaxlength(M_.endo_names(k1(ivar))), cellofchararraymaxlength(M_.exo_names)]);
end
%deal with covariances
if order > 1
var_name_width=max(2*(var_name_width+aux_var_additional_characters)+2,20); %account for covariances, separated by comma
@ -99,7 +105,11 @@ disp('POLICY AND TRANSITION FUNCTIONS')
% variable names
str = char(32*ones(1,var_name_width));
for i=1:nvar
str = [str sprintf(header_label_format, M_.endo_names{k1(ivar(i))})];
if options_.loglinear
str = [str sprintf(header_label_format, ['log(',M_.endo_names{k1(ivar(i))},')'])];
else
str = [str sprintf(header_label_format, M_.endo_names{k1(ivar(i))})];
end
end
disp(str);
%
@ -138,7 +148,11 @@ for k=1:nx
else
str1 = subst_auxvar(k1(klag(k,1)),klag(k,2)-M_.maximum_lag-2);
end
str = sprintf(label_format,str1);
if options_.loglinear
str = sprintf(label_format,['log(',str1,')']);
else
str = sprintf(label_format,str1);
end
for i=1:nvar
x = dr.ghx(ivar(i),k);
[str,flag]=get_print_string(str,x,value_format_zero,value_format_float,flag,options_);

15
matlab/disp_moments.m
View File

@ -50,7 +50,7 @@ y = y(ivar,options_.drop+1:end)';
ME_present=0;
if ~all(M_.H==0)
if isoctave
if isoctave && octave_ver_less_than('6')
[observable_pos_requested_vars, index_subset, index_observables] = intersect_stable(ivar, options_.varobs_id);
else
[observable_pos_requested_vars, index_subset, index_observables] = intersect(ivar, options_.varobs_id, 'stable');
@ -94,8 +94,8 @@ oo_.var(zero_variance_var_index,:)=0;
oo_.var(:,zero_variance_var_index)=0;
labels = M_.endo_names(ivar);
labels_TeX = M_.endo_names_tex(ivar);
labels=get_labels_transformed_vars(M_.endo_names,ivar,options_,false);
labels_TeX=get_labels_transformed_vars(M_.endo_names_tex,ivar,options_,true);
if ~options_.nomoments
z = [ m' s' s2' oo_.skewness oo_.kurtosis ];
@ -204,10 +204,11 @@ if ~options_.nodecomposition
headers = vertcat(' ', headers);
lh = cellofchararraymaxlength(M_.endo_names(ivar))+2;
dyntable(options_, title, vertcat(headers, 'Tot. lin. contr.'), ...
M_.endo_names(ivar), [oo_.variance_decomposition sum(oo_.variance_decomposition,2)], lh, 8, 2);
labels, [oo_.variance_decomposition sum(oo_.variance_decomposition,2)], lh, 8, 2);
if ME_present
headers_ME = vertcat(headers, 'ME');
dyntable(options_, [title,' WITH MEASUREMENT ERROR'], vertcat(headers_ME, 'Tot. lin. contr.'), M_.endo_names(ivar(index_subset)), ...
labels_obs=get_labels_transformed_vars(M_.endo_names,ivar(index_subset),options_,false);
dyntable(options_, [title,' WITH MEASUREMENT ERROR'], vertcat(headers_ME, 'Tot. lin. contr.'), labels_obs, ...
[oo_.variance_decomposition_ME sum(oo_.variance_decomposition_ME, 2)], lh, 8, 2);
end
if options_.TeX
@ -219,9 +220,10 @@ if ~options_.nodecomposition
labels_TeX, [oo_.variance_decomposition sum(oo_.variance_decomposition, 2)], lh, 8, 2);
if ME_present
headers_ME = vertcat(headers, 'ME');
labels_obs_TeX=get_labels_transformed_vars(M_.endo_names_tex,ivar(index_subset),options_,true);
dyn_latex_table(M_, options_, [title, ' WITH MEASUREMENT ERROR'], 'sim_var_decomp_ME', ...
vertcat(headers_ME, 'Tot. lin. contr.'), ...
labels_TeX(ivar(index_subset)), ...
labels_obs_TeX, ...
[oo_.variance_decomposition_ME sum(oo_.variance_decomposition_ME, 2)], lh, 8, 2);
end
end
@ -256,3 +258,4 @@ else
end
end

13
matlab/disp_steady_state.m
View File

@ -1,10 +1,11 @@
function disp_steady_state(M,oo)
% function disp_steady_state(M,oo)
function disp_steady_state(M,oo,options)
% function disp_steady_state(M,oo,options)
% computes and prints the steady state calculations
%
% INPUTS
% M structure of parameters
% oo structure of results
% options structure of options
%
% OUTPUTS
% none
@ -12,7 +13,7 @@ function disp_steady_state(M,oo)
% SPECIAL REQUIREMENTS
% none
% Copyright (C) 2001-2018 Dynare Team
% Copyright (C) 2001-2020 Dynare Team
%
% This file is part of Dynare.
%
@ -30,7 +31,11 @@ function disp_steady_state(M,oo)
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
skipline()
disp('STEADY-STATE RESULTS:')
if options.loglinear
disp('STEADY-STATE RESULTS FOR THE UNLOGGED VARIABLES:')
else
disp('STEADY-STATE RESULTS:')
end
skipline()
endo_names = char(M.endo_names);
steady_state = oo.steady_state;

29
matlab/disp_th_moments.m
View File

@ -55,7 +55,7 @@ oo_.var = oo_.gamma_y{1};
ME_present=0;
if ~all(diag(M_.H)==0)
if isoctave
if isoctave && octave_ver_less_than('6')
[observable_pos_requested_vars,index_subset,index_observables]=intersect_stable(ivar,options_.varobs_id);
else
[observable_pos_requested_vars,index_subset,index_observables]=intersect(ivar,options_.varobs_id,'stable');
@ -82,11 +82,11 @@ if size(stationary_vars, 1) > 0
end
title = add_filter_subtitle(title, options_);
headers = {'VARIABLE';'MEAN';'STD. DEV.';'VARIANCE'};
labels = M_.endo_names(ivar);
labels=get_labels_transformed_vars(M_.endo_names,ivar,options_,false);
lh = cellofchararraymaxlength(labels)+2;
dyntable(options_, title, headers, labels, z, lh, 11, 4);
if options_.TeX
labels = M_.endo_names_tex(ivar);
labels=get_labels_transformed_vars(M_.endo_names_tex,ivar,options_,true);
lh = cellofchararraymaxlength(labels)+2;
dyn_latex_table(M_, options_, title, 'th_moments', headers, labels, z, lh, 11, 4);
end
@ -102,26 +102,29 @@ if size(stationary_vars, 1) > 0
headers = M_.exo_names;
headers(M_.exo_names_orig_ord) = headers;
headers = vertcat(' ', headers);
lh = cellofchararraymaxlength(M_.endo_names(ivar(stationary_vars)))+2;
dyntable(options_, title, headers, M_.endo_names(ivar(stationary_vars)), 100*oo_.gamma_y{options_.ar+2}(stationary_vars,:), lh, 8, 2);
labels=get_labels_transformed_vars(M_.endo_names,ivar(stationary_vars),options_,false);
lh = cellofchararraymaxlength(labels)+2;
dyntable(options_, title, headers, labels, 100*oo_.gamma_y{options_.ar+2}(stationary_vars,:), lh, 8, 2);
if ME_present
if isoctave
if isoctave && octave_ver_less_than('6')
[stationary_observables, pos_index_subset] = intersect_stable(index_subset, stationary_vars);
else
[stationary_observables, pos_index_subset] = intersect(index_subset, stationary_vars, 'stable');
end
headers_ME = vertcat(headers, 'ME');
dyntable(options_, [title,' WITH MEASUREMENT ERROR'], headers_ME, M_.endo_names(ivar(stationary_observables)), ...
labels=get_labels_transformed_vars(M_.endo_names,ivar(stationary_observables),options_,false);
dyntable(options_, [title,' WITH MEASUREMENT ERROR'], headers_ME, labels, ...
oo_.variance_decomposition_ME(pos_index_subset,:), lh, 8, 2);
end
if options_.TeX
headers = M_.exo_names_tex;
headers = vertcat(' ', headers);
labels = M_.endo_names_tex(ivar(stationary_vars));
labels=get_labels_transformed_vars(M_.endo_names_tex,ivar(stationary_vars),options_,true);
lh = cellofchararraymaxlength(labels)+2;
dyn_latex_table(M_, options_, title, 'th_var_decomp_uncond', headers, labels, 100*oo_.gamma_y{options_.ar+2}(stationary_vars,:), lh, 8, 2);
if ME_present
headers_ME = vertcat(headers, 'ME');
labels=get_labels_transformed_vars(M_.endo_names_tex,ivar(stationary_observables),options_,true);
dyn_latex_table(M_, options_, [title,' WITH MEASUREMENT ERROR'], ...
'th_var_decomp_uncond_ME', headers_ME, labels, oo_.variance_decomposition_ME(pos_index_subset,:), lh, 8, 2);
end
@ -172,12 +175,12 @@ if ~options_.nocorr && size(stationary_vars, 1)>0
title = 'MATRIX OF CORRELATIONS';
end
title = add_filter_subtitle(title, options_);
labels = M_.endo_names(ivar(i1));
labels=get_labels_transformed_vars(M_.endo_names,ivar(i1),options_,false);
headers = vertcat('Variables', labels);
lh = cellofchararraymaxlength(labels)+2;
dyntable(options_, title, headers, labels, corr(i1,i1), lh, 8, 4);
if options_.TeX
labels = M_.endo_names_tex(ivar(i1));
labels=get_labels_transformed_vars(M_.endo_names_tex,ivar(i1),options_,true);
headers = vertcat('Variables', labels);
lh = cellofchararraymaxlength(labels)+2;
dyn_latex_table(M_, options_, title, 'th_corr_matrix', headers, labels, corr(i1,i1), lh, 8, 4);
@ -199,15 +202,15 @@ if options_.ar > 0 && size(stationary_vars, 1) > 0
title = 'COEFFICIENTS OF AUTOCORRELATION';
end
title = add_filter_subtitle(title, options_);
labels = M_.endo_names(ivar(i1));
labels=get_labels_transformed_vars(M_.endo_names,ivar(i1),options_,false);
headers = vertcat('Order ', cellstr(int2str([1:options_.ar]')));
lh = cellofchararraymaxlength(labels)+2;
dyntable(options_, title, headers, labels, z, lh, 8, 4);
if options_.TeX
labels = M_.endo_names_tex(ivar(i1));
labels=get_labels_transformed_vars(M_.endo_names_tex,ivar(i1),options_,true);
headers = vertcat('Order ', cellstr(int2str([1:options_.ar]')));
lh = cellofchararraymaxlength(labels)+2;
dyn_latex_table(M_, options_, title, 'th_autocorr_matrix', headers, labels, z, lh, 8, 4);
end
end
end
end

0
matlab/disp_th_moments_order3.m → matlab/disp_th_moments_pruned_state_space.m
View File

9
matlab/display_conditional_variance_decomposition.m
View File

@ -63,9 +63,12 @@ else
end
headers = vertcat(' ', headers);
lh = cellofchararraymaxlength(M_.endo_names(SubsetOfVariables))+2;
labels=get_labels_transformed_vars(M_.endo_names,SubsetOfVariables,options_,false);
lh = cellofchararraymaxlength(labels)+2;
if options_.TeX
labels_TeX = M_.endo_names_tex(SubsetOfVariables);
labels_TeX=get_labels_transformed_vars(M_.endo_names,SubsetOfVariables,options_,true);
lh = cellofchararraymaxlength(labels_TeX)+2;
end
@ -76,7 +79,7 @@ for i=1:length(Steps)
for j=1:shock_number
vardec_i(:,j) = 100*conditional_decomposition_array(:,i,j);
end
dyntable(options_, '', headers, M_.endo_names(SubsetOfVariables), vardec_i, lh, 8, 2);
dyntable(options_, '', headers, labels, vardec_i, lh, 8, 2);
if options_.TeX
dyn_latex_table(M_, options_, [title, '; Period ' int2str(Steps(i))], ['th_var_decomp_cond_h', int2str(Steps(i))], headers_TeX, labels_TeX, vardec_i, lh, 8, 2);
end

2
matlab/dsge_likelihood.m
View File

@ -791,7 +791,7 @@ if DynareOptions.endogenous_prior==1
if DynareOptions.lik_init==2 || DynareOptions.lik_init==3
error('Endogenous prior not supported with non-stationary models')
else
[lnpriormom] = endogenous_prior(Y,Pstar,BayesInfo,H);
[lnpriormom] = endogenous_prior(Y,DatasetInfo,Pstar,BayesInfo,H);
fval = (likelihood-lnprior-lnpriormom);
end
else

2
matlab/dsge_simulated_theoretical_conditional_variance_decomposition.m
View File

@ -83,7 +83,7 @@ MaXNumberOfConditionalDecompLines = ceil(options_.MaxNumberOfBytes/NumberOfSaved
ME_present=0;
if ~all(diag(M_.H)==0)
if isoctave
if isoctave && octave_ver_less_than('6')
[observable_pos_requested_vars,index_subset,index_observables]=intersect_stable(ivar,options_.varobs_id);
else
[observable_pos_requested_vars,index_subset,index_observables]=intersect(ivar,options_.varobs_id,'stable');

2
matlab/dsge_simulated_theoretical_variance_decomposition.m
View File

@ -86,7 +86,7 @@ MaXNumberOfDecompLines = ceil(options_.MaxNumberOfBytes/NumberOfSavedElementsPer
ME_present=0;
if ~all(diag(M_.H)==0)
if isoctave
if isoctave && octave_ver_less_than('6')
[observable_pos_requested_vars,index_subset,index_observables]=intersect_stable(ivar,options_.varobs_id);
else
[observable_pos_requested_vars,index_subset,index_observables]=intersect(ivar,options_.varobs_id,'stable');

2
matlab/dyn_second_order_solver.m
View File

@ -62,7 +62,7 @@ k1 = nonzeros(M.lead_lag_incidence(:,dr.order_var)');
kk1 = [k1; length(k1)+(1:M.exo_nbr+M.exo_det_nbr)'];
nk = size(kk1,1);
kk2 = reshape(1:nk^2,nk,nk);
ic = [ M.nstatic+(1:M.nspred) M.endo_nbr+(1:size(dr.ghx,2)-M.nspred) ]';
ic = [ M.nstatic+(1:M.nspred) ]';
klag = M.lead_lag_incidence(1,dr.order_var); %columns are in DR order
kcurr = M.lead_lag_incidence(2,dr.order_var); %columns are in DR order

40
matlab/dynare.m
View File

@ -16,7 +16,7 @@ function dynare(fname, varargin)
% SPECIAL REQUIREMENTS
% none
% Copyright (C) 2001-2020 Dynare Team
% Copyright (C) 2001-2021 Dynare Team
%
% This file is part of Dynare.
%
@ -184,16 +184,14 @@ nolog = ismember('nolog', varargin) || ismember('nolog', file_opts);
onlymacro = ismember('onlymacro', varargin) || ismember('onlymacro', file_opts);
onlyjson = ismember('onlyjson', varargin) || ismember('onlyjson', file_opts);
if ispc
arch = getenv('PROCESSOR_ARCHITECTURE');
else
[~, arch] = system('uname -m');
end
if isempty(strfind(arch, '64'))
arch_ext = '32';
else
arch_ext = '64';
% Start journal
diary off
if ~nolog
logfile = [ fname(1:end-4) '.log' ];
if exist(logfile, 'file')
delete(logfile)
end
diary(logfile)
end
if preprocessoroutput
@ -206,7 +204,7 @@ if preprocessoroutput
end
end
command = ['"' dynareroot 'preprocessor' arch_ext filesep 'dynare_m" ' fname] ;
command = ['"' dynareroot '..' filesep 'preprocessor' filesep 'dynare-preprocessor" ' fname] ;
command = [ command ' mexext=' mexext ' "matlabroot=' matlabroot '"'];
% Properly quote arguments before passing them to the shell
if ~isempty(varargin)
@ -265,15 +263,8 @@ if exist(fname(1:end-4),'dir') && exist([fname(1:end-4) filesep 'hooks'],'dir')
run([fname(1:end-4) filesep 'hooks/postprocessing'])
end
% Save preprocessor result in logfile (if `no_log' option not present)
if ~nolog
logname = [fname(1:end-4) '.log'];
fid = fopen(logname, 'w');
fprintf(fid, '%s', result);
fclose(fid);
end
if status
diary off
% Should not use "error(result)" since message will be truncated if too long
error('Dynare: preprocessing failed')
end
@ -286,7 +277,14 @@ end
% within the driver will clean the rest)
clear(['+' fname '/driver'])
evalin('base',[fname '.driver']) ;
try
evalin('base',[fname '.driver']) ;
catch ME
diary off
rethrow(ME)
end
diary off
end

4
matlab/dynare_config.m
View File

@ -90,8 +90,8 @@ if isoctave && octave_ver_less_than('5')
p{end+1} = '/missing/ordeig';
end
%% intersect(…, 'stable') doesn't exist in Octave
if isoctave
%% intersect(…, 'stable') doesn't exist in Octave < 6
if isoctave && octave_ver_less_than('6')
p{end+1} = '/missing/intersect_stable';
end

4
matlab/dynare_estimation.m
View File

@ -147,7 +147,7 @@ if nnobs > 1 && horizon > 0
[~,nr,nc,~,~,nstar] = pltorg(nvar);
m = 1;
plot_index=0;
OutputDirectoryName = CheckPath('graphs',M_.fname);
OutputDirectoryName = CheckPath('graphs',M_.dname);
for i = 1:length(var_list)
if mod(i,nstar) == 1
plot_index=plot_index+1;
@ -208,7 +208,7 @@ if nnobs > 1 && horizon > 0
xlim([nobs(1)-offsetx nobs(end)+horizon])
m = m + 1;
if mod(i+1,nstar) == 1 || i==length(var_list)
dyn_saveas(hfig,[M_.fname,filesep,'graphs',filesep M_.fname '_RecursiveForecasts_' int2str(plot_index)],options_.nodisplay,options_.graph_format);
dyn_saveas(hfig,[M_.dname,filesep,'graphs',filesep M_.fname '_RecursiveForecasts_' int2str(plot_index)],options_.nodisplay,options_.graph_format);
end
end
end

22
matlab/dynare_estimation_1.m
View File

@ -557,11 +557,11 @@ if (~((any(bayestopt_.pshape > 0) && options_.mh_replic) || (any(bayestopt_.psha
if ~options_.nograph
[nbplt,nr,nc,lr,lc,nstar] = pltorg(M_.exo_nbr);
if ~exist([M_.fname '/graphs'],'dir')
mkdir(M_.fname,'graphs');
if ~exist([M_.dname '/graphs'],'dir')
mkdir(M_.dname,'graphs');
end
if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
fidTeX = fopen([M_.fname, '/graphs/' M_.fname '_SmoothedShocks.tex'],'w');
fidTeX = fopen([M_.dname, '/graphs/' M_.fname '_SmoothedShocks.tex'],'w');
fprintf(fidTeX,'%% TeX eps-loader file generated by dynare_estimation_1.m (Dynare).\n');
fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
fprintf(fidTeX,' \n');
@ -599,11 +599,11 @@ if (~((any(bayestopt_.pshape > 0) && options_.mh_replic) || (any(bayestopt_.psha
title(name,'Interpreter','none')
end
end
dyn_saveas(fh,[M_.fname, '/graphs/' M_.fname '_SmoothedShocks' int2str(plt)],options_.nodisplay,options_.graph_format);
dyn_saveas(fh,[M_.dname, '/graphs/' M_.fname '_SmoothedShocks' int2str(plt)],options_.nodisplay,options_.graph_format);
if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
fprintf(fidTeX,'\\begin{figure}[H]\n');
fprintf(fidTeX,'\\centering \n');
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_SmoothedShocks%s}\n',options_.figures.textwidth*min(i/nc,1),[M_.fname, '/graphs/' M_.fname],int2str(plt));
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_SmoothedShocks%s}\n',options_.figures.textwidth*min(i/nc,1),[M_.dname, '/graphs/' M_.fname],int2str(plt));
fprintf(fidTeX,'\\caption{Smoothed shocks.}');
fprintf(fidTeX,'\\label{Fig:SmoothedShocks:%s}\n',int2str(plt));
fprintf(fidTeX,'\\end{figure}\n');
@ -628,7 +628,7 @@ if (~((any(bayestopt_.pshape > 0) && options_.mh_replic) || (any(bayestopt_.psha
if ~options_.nograph
[nbplt,nr,nc,lr,lc,nstar] = pltorg(number_of_plots_to_draw);
if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
fidTeX = fopen([M_.fname, '/graphs/' M_.fname '_SmoothedObservationErrors.tex'],'w');
fidTeX = fopen([M_.dname, '/graphs/' M_.fname '_SmoothedObservationErrors.tex'],'w');
fprintf(fidTeX,'%% TeX eps-loader file generated by dynare_estimation_1.m (Dynare).\n');
fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
fprintf(fidTeX,' \n');
@ -665,11 +665,11 @@ if (~((any(bayestopt_.pshape > 0) && options_.mh_replic) || (any(bayestopt_.psha
title(name,'Interpreter','none')
end
end
dyn_saveas(fh,[M_.fname, '/graphs/' M_.fname '_SmoothedObservationErrors' int2str(plt)],options_.nodisplay,options_.graph_format);
dyn_saveas(fh,[M_.dname, '/graphs/' M_.fname '_SmoothedObservationErrors' int2str(plt)],options_.nodisplay,options_.graph_format);
if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
fprintf(fidTeX,'\\begin{figure}[H]\n');
fprintf(fidTeX,'\\centering \n');
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_SmoothedObservationErrors%s}\n',options_.figures.textwidth*min(i/nc,1),[M_.fname, '/graphs/' M_.fname],int2str(plt));
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_SmoothedObservationErrors%s}\n',options_.figures.textwidth*min(i/nc,1),[M_.dname, '/graphs/' M_.fname],int2str(plt));
fprintf(fidTeX,'\\caption{Smoothed observation errors.}');
fprintf(fidTeX,'\\label{Fig:SmoothedObservationErrors:%s}\n',int2str(plt));
fprintf(fidTeX,'\\end{figure}\n');
@ -689,7 +689,7 @@ if (~((any(bayestopt_.pshape > 0) && options_.mh_replic) || (any(bayestopt_.psha
if ~options_.nograph
[nbplt,nr,nc,lr,lc,nstar] = pltorg(n_varobs);
if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
fidTeX = fopen([M_.fname, '/graphs/' M_.fname '_HistoricalAndSmoothedVariables.tex'],'w');
fidTeX = fopen([M_.dname, '/graphs/' M_.fname '_HistoricalAndSmoothedVariables.tex'],'w');
fprintf(fidTeX,'%% TeX eps-loader file generated by dynare_estimation_1.m (Dynare).\n');
fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
fprintf(fidTeX,' \n');
@ -727,11 +727,11 @@ if (~((any(bayestopt_.pshape > 0) && options_.mh_replic) || (any(bayestopt_.psha
title(name,'Interpreter','none')
end
end
dyn_saveas(fh,[M_.fname, '/graphs/' M_.fname '_HistoricalAndSmoothedVariables' int2str(plt)],options_.nodisplay,options_.graph_format);
dyn_saveas(fh,[M_.dname, '/graphs/' M_.fname '_HistoricalAndSmoothedVariables' int2str(plt)],options_.nodisplay,options_.graph_format);
if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
fprintf(fidTeX,'\\begin{figure}[H]\n');
fprintf(fidTeX,'\\centering \n');
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_HistoricalAndSmoothedVariables%s}\n',options_.figures.textwidth*min(i/nc,1),[M_.fname, '/graphs/' M_.fname],int2str(plt));
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_HistoricalAndSmoothedVariables%s}\n',options_.figures.textwidth*min(i/nc,1),[M_.dname, '/graphs/' M_.fname],int2str(plt));
fprintf(fidTeX,'\\caption{Historical and smoothed variables.}');
fprintf(fidTeX,'\\label{Fig:HistoricalAndSmoothedVariables:%s}\n',int2str(plt));
fprintf(fidTeX,'\\end{figure}\n');

12
matlab/dynare_estimation_init.m
View File

@ -555,6 +555,18 @@ end
%set options for old interface from the ones for new interface
if ~isempty(dataset_)
options_.nobs = dataset_.nobs;
if options_.endogenous_prior
if dataset_info.missing.no_more_missing_observations<dataset_.nobs-10
fprintf('\ndynare_estimation_init: There are missing observations in the data.\n')
fprintf('dynare_estimation_init: I am computing the moments for the endogenous prior only\n')
fprintf('dynare_estimation_init: on the observations after the last missing one, i.e. %u.\n',dataset_info.missing.no_more_missing_observations)
else
fprintf('\ndynare_estimation_init: There are too many missing observations in the data.\n')
fprintf('dynare_estimation_init: The endogenous_prior-option needs a consistent sample of \n')
fprintf('dynare_estimation_init: at least 10 full observations at the end.\n')
error('The endogenous_prior-option does not support your missing data.')
end
end
end
% setting steadystate_check_flag option

4
matlab/dynare_identification.m
View File

@ -139,6 +139,8 @@ options_ident = set_default_option(options_ident,'tol_deriv',1.e-8);
% tolerance level for selecting columns of non-zero derivatives
options_ident = set_default_option(options_ident,'tol_sv',1.e-3);
% tolerance level for selecting non-zero singular values in identification_checks.m
options_ident = set_default_option(options_ident,'schur_vec_tol',1e-11);
% tolerance level used to find nonstationary variables in Schur decomposition of the transition matrix.
%check whether to compute identification strength based on information matrix
if ~isfield(options_ident,'no_identification_strength')
@ -289,7 +291,7 @@ if options_ident.order == 3
end
options_.ar = options_ident.ar;
options_.prior_mc = options_ident.prior_mc;
options_.Schur_vec_tol = 1.e-8;
options_.schur_vec_tol = options_ident.schur_vec_tol;
options_.nomoments = 0;
options_.analytic_derivation=options_ident.analytic_derivation;
% 1: analytic derivation of gradient and hessian of likelihood in dsge_likelihood.m, only works for stationary models, i.e. kalman_algo<3

10
matlab/endogenous_prior.m
View File

@ -1,8 +1,9 @@
function [lnpriormom] = endogenous_prior(data,Pstar,BayesInfo,H)
function [lnpriormom] = endogenous_prior(data,dataset_info, Pstar,BayesInfo,H)
% Computes the endogenous log prior addition to the initial prior
%
% INPUTS
% data [double] n*T vector of data observations
% dataset_info [structure] various information about the dataset
% Pstar [double] k*k matrix of
% BayesInfo [structure]
%
@ -11,7 +12,7 @@ function [lnpriormom] = endogenous_prior(data,Pstar,BayesInfo,H)
% Code to implement notes on endogenous priors by Lawrence Christiano,
% specified in the appendix of:
% ’Introducing Financial Frictions and Unemployment into a Small Open Economy Model’
% Introducing Financial Frictions and Unemployment into a Small Open Economy Model
% by Lawrence J. Christiano, Mathias Trabandt and Karl Walentin (2011), Journal of Economic Dynamics and Control
% this is the 'mother' of the priors on the model parameters.
% the priors include a metric across some choosen moments of the (supposedly
@ -41,7 +42,8 @@ function [lnpriormom] = endogenous_prior(data,Pstar,BayesInfo,H)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
Y=data';
Y=data(:,dataset_info.missing.no_more_missing_observations:end)';
[Tsamp,n]=size(Y); % sample length and number of matched moments (here set equal to nr of observables)
hmat=zeros(n,Tsamp);
@ -85,4 +87,4 @@ Z=II(mf,:);
Ftheta=diag(Z*Pstar(:,mf)+H);
% below commented out line is for Del Negro Schorfheide style priors:
% lnpriormom=-.5*n*TT*log(2*pi)-.5*TT*log(det(sigma))-.5*TT*trace(inv(sigma)*(gamyy-2*phi'*gamxy+phi'*gamxx*phi));
lnpriormom=.5*n*log(Tsamp/(2*pi))-.5*log(det(Shat))-.5*Tsamp*(Fhat-Ftheta)'/Shat*(Fhat-Ftheta);
lnpriormom=.5*n*log(Tsamp/(2*pi))-.5*log(det(Shat))-.5*Tsamp*(Fhat-Ftheta)'/Shat*(Fhat-Ftheta);

2
matlab/ep/ep_notes.org
View File

@ -37,7 +37,6 @@ ep:
ortpol: 'hermite'
nodes: 5
pruned: [1x1 struct]
use_bytecode: 0
verbosity: 0
initial_period: NaN
lmmcp:
@ -59,7 +58,6 @@ ut: (unscented free parameter)
pfm.stochastic_order = DynareOptions.ep.stochastic.order;
pfm.periods = DynareOptions.ep.periods;
pfm.verbose = DynareOptions.ep.verbosity;
pfm.use_bytecode = DynareOptions.ep.use_bytecode ;
* in extended_path_core, one passes options.ep and individual options

8
matlab/ep/euler_equation_error.m
View File

@ -1,6 +1,6 @@
function e = euler_equation_error(y0,x,innovations,M,options,oo,pfm,nodes,weights)
% Copyright (C) 2016-2017 Dynare Team
% Copyright (C) 2016-2020 Dynare Team
%
% This file is part of Dynare.
%
@ -23,7 +23,7 @@ ep = options.ep;
M.endo_nbr, M.exo_nbr, ...
innovations.positive_var_indx, ...
x, ep.init, y0, oo.steady_state, ...
0, ep.use_bytecode, ...
0, ...
ep.stochastic.order, M, ...
pfm, ep.stochastic.algo, ...
ep.solve_algo, ...
@ -39,7 +39,7 @@ for i=1:length(nodes)
[y2, info_convergence, endogenousvariablespaths] = ...
extended_path_core(ep.periods, M.endo_nbr, M.exo_nbr, ...
innovations.positive_var_indx, x2, ep.init, ...
y1, oo.steady_state, 0, ep.use_bytecode, ...
y1, oo.steady_state, 0, ...
ep.stochastic.order, M, pfm, ep.stochastic.algo, ...
ep.solve_algo, ep.stack_solve_algo, options.lmmcp, ...
options, oo, []);
@ -47,4 +47,4 @@ for i=1:length(nodes)
z = [y0(i_pred); y1; y2(i_fwrd)];
res(:,i) = dynamic_model(z,x,M.params,oo.steady_state,2);
end
e = res*weights;
e = res*weights;

6
matlab/ep/extended_path.m
View File

@ -19,7 +19,7 @@ function [ts, DynareResults] = extended_path(initialconditions, samplesize, exog
%
% SPECIAL REQUIREMENTS
% Copyright (C) 2009-2018 Dynare Team
% Copyright (C) 2009-2020 Dynare Team
%
% This file is part of Dynare.
%
@ -70,7 +70,7 @@ while (t <= samplesize)
[endogenous_variables_paths(:,t), info_convergence, endogenousvariablespaths] = extended_path_core(ep.periods, DynareModel.endo_nbr, DynareModel.exo_nbr, innovations.positive_var_indx, ...
spfm_exo_simul, ep.init, endogenous_variables_paths(:,t-1), ...
DynareResults.steady_state, ...
verbosity, ep.use_bytecode, ep.stochastic.order, ...
verbosity, ep.stochastic.order, ...
DynareModel, pfm, ep.stochastic.algo, ep.solve_algo, ep.stack_solve_algo, ...
DynareOptions.lmmcp, ...
DynareOptions, ...
@ -105,4 +105,4 @@ assignin('base', 'Simulated_time_series', ts);
if ~nargout || nargout<2
assignin('base', 'oo_', DynareResults);
end
end

66
matlab/ep/extended_path_core.m
View File

@ -1,7 +1,7 @@
function [y, info_convergence, endogenousvariablespaths] = extended_path_core(periods,endo_nbr,exo_nbr,positive_var_indx, ...
exo_simul,init,initial_conditions,...
steady_state, ...
debug,bytecode_flag,order,M,pfm,algo,solve_algo,stack_solve_algo,...
debug,order,M,pfm,algo,solve_algo,stack_solve_algo,...
olmmcp,options,oo,initialguess)
% Copyright (C) 2016-2020 Dynare Team
@ -40,46 +40,38 @@ if debug
save ep_test_1.mat endo_simul exo_simul
end
if bytecode_flag && ~ep.stochastic.order
try
tmp = bytecode('dynamic', endo_simul, exo_simul, M.params, endo_simul, periods);
flag = false;
catch ME
disp(ME.message);
flag = true;
end
else
flag = true;
if options.bytecode && order > 0
error('Option order > 0 of extended_path command is not compatible with bytecode option.')
end
if options.block && order > 0
error('Option order > 0 of extended_path command is not compatible with block option.')
end
if flag
if order == 0
options.periods = periods;
options.block = pfm.block;
oo.endo_simul = endo_simul;
oo.exo_simul = exo_simul;
oo.steady_state = steady_state;
options.bytecode = bytecode_flag;
options.lmmcp = olmmcp;
options.solve_algo = solve_algo;
options.stack_solve_algo = stack_solve_algo;
tmp = perfect_foresight_solver_core(M, options, oo);
if ~tmp.deterministic_simulation.status
info_convergence = false;
else
info_convergence = true;
end
if order == 0
options.periods = periods;
options.block = pfm.block;
oo.endo_simul = endo_simul;
oo.exo_simul = exo_simul;
oo.steady_state = steady_state;
options.lmmcp = olmmcp;
options.solve_algo = solve_algo;
options.stack_solve_algo = stack_solve_algo;
tmp = perfect_foresight_solver_core(M, options, oo);
if ~tmp.deterministic_simulation.status
info_convergence = false;
else
switch(algo)
case 0
[flag, tmp.endo_simul] = ...
solve_stochastic_perfect_foresight_model(endo_simul, exo_simul, pfm, ep.stochastic.quadrature.nodes, ep.stochastic.order);
case 1
[flag, tmp.endo_simul] = ...
solve_stochastic_perfect_foresight_model_1(endo_simul, exo_simul, options, pfm, ep.stochastic.order);
end
info_convergence = ~flag;
info_convergence = true;
end
else
switch(algo)
case 0
[flag, tmp.endo_simul] = ...
solve_stochastic_perfect_foresight_model(endo_simul, exo_simul, pfm, ep.stochastic.quadrature.nodes, ep.stochastic.order);
case 1
[flag, tmp.endo_simul] = ...
solve_stochastic_perfect_foresight_model_1(endo_simul, exo_simul, options, pfm, ep.stochastic.order);
end
info_convergence = ~flag;
end
if ~info_convergence && ~options.no_homotopy

3
matlab/ep/extended_path_initialization.m
View File

@ -35,9 +35,6 @@ function [initial_conditions, innovations, pfm, ep, verbosity, DynareOptions, Dy
ep = DynareOptions.ep;
% Set bytecode flag (see #1717)
ep.use_bytecode = DynareOptions.bytecode;
% Set verbosity levels.
DynareOptions.verbosity = ep.verbosity;
verbosity = ep.verbosity+ep.debug;

6
matlab/ep/extended_path_mc.m
View File

@ -19,7 +19,7 @@ function Simulations = extended_path_mc(initialconditions, samplesize, replic, e
%
% SPECIAL REQUIREMENTS
% Copyright (C) 2016-2017 Dynare Team
% Copyright (C) 2016-2020 Dynare Team
%
% This file is part of Dynare.
%
@ -83,7 +83,7 @@ if ep.parallel
[endogenous_variables_paths(:,t), info_convergence] = extended_path_core(ep.periods, DynareModel.endo_nbr, DynareModel.exo_nbr, innovations_.positive_var_indx, ...
spfm_exo_simul, ep.init, endogenous_variables_paths(:,t-1), ...
DynareResults_.steady_state, ...
ep.verbosity, ep.use_bytecode, ep.stochastic.order, ...
ep.verbosity, ep.stochastic.order, ...
DynareModel, pfm,ep.stochastic.algo, ep.solve_algo, ep.stack_solve_algo, ...
DynareOptions.lmmcp, DynareOptions, DynareResults_);
if ~info_convergence
@ -112,7 +112,7 @@ else
[endogenous_variables_paths(:,t), info_convergence] = extended_path_core(ep.periods, DynareModel.endo_nbr, DynareModel.exo_nbr, innovations.positive_var_indx, ...
spfm_exo_simul, ep.init, endogenous_variables_paths(:,t-1), ...
DynareResults.steady_state, ...
ep.verbosity, ep.use_bytecode, ep.stochastic.order, ...
ep.verbosity, ep.stochastic.order, ...
DynareModel, pfm,ep.stochastic.algo, ep.solve_algo, ep.stack_solve_algo, ...
DynareOptions.lmmcp, DynareOptions, DynareResults);
if ~info_convergence

7
matlab/ep/setup_stochastic_perfect_foresight_model_solver.m
View File

@ -1,6 +1,6 @@
function pfm = setup_stochastic_perfect_foresight_model_solver(DynareModel,DynareOptions,DynareOutput)
% Copyright (C) 2013-2017 Dynare Team
% Copyright (C) 2013-2020 Dynare Team
%
% This file is part of Dynare.
%
@ -66,8 +66,9 @@ else
end
pfm.i_cols_j = 1:pfm.nd;
pfm.i_upd = pfm.ny+(1:pfm.periods*pfm.ny);
pfm.dynamic_model = str2func([DynareModel.fname,'.dynamic']);
if ~DynareOptions.bytecode
pfm.dynamic_model = str2func([DynareModel.fname,'.dynamic']);
end
pfm.verbose = DynareOptions.ep.verbosity;
pfm.maxit_ = DynareOptions.simul.maxit;
pfm.tolerance = DynareOptions.dynatol.f;
pfm.use_bytecode = DynareOptions.ep.use_bytecode ;

13
matlab/evaluate_planner_objective.m
View File

@ -37,11 +37,6 @@ dr = oo.dr;
exo_nbr = M.exo_nbr;
nstatic = M.nstatic;
nspred = M.nspred;
if nspred > 180
fprintf('\nevaluate_planner_objective: model too large, can''t evaluate planner objective\n')
planner_objective_value = NaN;
return
end
beta = get_optimal_policy_discount_factor(M.params, M.param_names);
Gy = dr.ghx(nstatic+(1:nspred),:);
@ -64,7 +59,13 @@ Uyygygu = A_times_B_kronecker_C(Uyy,gy,gu);
Wbar =U/(1-beta); %steady state welfare
Wy = Uy*gy/(eye(nspred)-beta*Gy);
Wu = Uy*gu+beta*Wy*Gu;
Wyy = Uyygygy/(eye(nspred*nspred)-beta*kron(Gy,Gy));
% Wyy = Uyygygy/(eye(nspred*nspred)-beta*kron(Gy,Gy)); %solve Wyy=Uyy*kron(gy,gy)+beta*Wyy*kron(Gy,Gy)
if isempty(options.qz_criterium)
options.qz_criterium = 1+1e-6;
end
%solve Lyapunuv equation Wyy=gy'*Uyy*gy+beta*Gy'Wyy*Gy
Wyy = reshape(lyapunov_symm(sqrt(beta)*Gy',reshape(Uyygygy,nspred,nspred),options.lyapunov_fixed_point_tol,options.qz_criterium,options.lyapunov_complex_threshold, 3, options.debug),1,nspred*nspred);
Wyygugu = A_times_B_kronecker_C(Wyy,Gu,Gu);
Wyygygu = A_times_B_kronecker_C(Wyy,Gy,Gu);
Wuu = Uyygugu+beta*Wyygugu;

20
matlab/evaluate_steady_state.m
View File

@ -309,12 +309,22 @@ if M.static_and_dynamic_models_differ
if options.bytecode
[r, ~]= bytecode('dynamic','evaluate', z, zx, params, ys, 1);
elseif options.block
[r, oo.dr] = feval([M.fname '.dynamic'], z', zx, params, ys, M.maximum_lag+1, oo.dr);
T=NaN(M.block_structure.dyn_tmp_nbr, 1);
for i = 1:length(M.block_structure.block)
[rr, yy, T, g] = feval([M.fname '.dynamic'], i, ...
dynvars_from_endo_simul(z, M.maximum_lag+1, M), ...
zx, params, ys, T, M.maximum_lag+1, false);
if M.block_structure.block(i).Simulation_Type == 1 || ... % evaluateForward
M.block_structure.block(i).Simulation_Type == 2 % evaluateBackward
vidx = M.block_structure.block(i).variable;
rr = yy(M.lead_lag_incidence(M.maximum_endo_lag+1, vidx)) - oo.steady_state(vidx);
end
idx = M.block_structure.block(i).equation;
r(idx) = rr;
end
else
iyv = M.lead_lag_incidence';
iyr0 = find(iyv(:));
xys = z(iyr0);
r = feval([M.fname '.dynamic'], z(iyr0), zx, params, ys, M.maximum_lag + 1);
r = feval([M.fname '.dynamic'], dynvars_from_endo_simul(z, M.maximum_lag+1, M), ...
zx, params, ys, M.maximum_lag + 1);
end
% Fail if residual greater than tolerance
if max(abs(r)) > options.solve_tolf

2
matlab/get_error_message.m
View File

@ -122,7 +122,7 @@ switch info(1)
case 64
message = 'discretionary_policy: the derivatives of the objective function contain NaN.';
case 65
message = 'discretionary_policy: the model must be written in deviation form and not have constant terms.';
message = 'discretionary_policy: the model must be written in deviation form and not have constant terms or an analytical steady state meeds to be provided.';
case 66
message = 'discretionary_policy: the objective function must have zero first order derivatives.';
case 71

2
matlab/get_variance_of_endogenous_variables.m
View File

@ -48,7 +48,7 @@ n = length(i_var);
[vx,u] = lyapunov_symm(A,B*Sigma_e*B',options_.lyapunov_fixed_point_tol,options_.qz_criterium,options_.lyapunov_complex_threshold, [], options_.debug);
if size(u,2) > 0
i_stat = find(any(abs(ghx*u) < options_.Schur_vec_tol,2)); %only set those variances of objective function for which variance is finite
i_stat = find(any(abs(ghx*u) < options_.schur_vec_tol,2)); %only set those variances of objective function for which variance is finite
ghx = ghx(i_stat,:);
ghu = ghu(i_stat,:);
else

5
matlab/global_initialization.m
View File

@ -144,11 +144,6 @@ dyn_first_order_solver();
% Set dynare random generator and seed.
set_dynare_seed('default');
% Create directories
[~,~]=mkdir(M_.fname);
[~,~]=mkdir([M_.fname filesep 'Output']);
% Load user configuration file.
if isfield(options_, 'global_init_file')
if isfile(options_.global_init_file)

2
matlab/histvalf.m
View File

@ -31,8 +31,6 @@ function [endo_histval, exo_histval, exo_det_histval] = histvalf(M, options)
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
series = histvalf_initvalf('HISTVAL', M, options);
% capture the difference between stochastic and
% perfect foresight setup
k = M.orig_maximum_lag - M.maximum_lag + 1;
endo_histval = series{M.endo_names{:}}.data(k:end, :)';

27
matlab/identification_analysis.m
View File

@ -147,18 +147,23 @@ if info(1) == 0 %no errors in solution
if init
%check stationarity
if ~no_identification_moments
ind_dMOMENTS = (find(max(abs(dMOMENTS'),[],1) > tol_deriv)); %index for non-zero rows
if isempty(ind_dMOMENTS) && any(any(isnan(dMOMENTS)))
error('There are NaN in the dMOMENTS matrix. Please check whether your model has units roots and you forgot to set diffuse_filter=1.' )
end
if any(any(isnan(MOMENTS)))
error('There are NaN''s in the theoretical moments: make sure that for non-stationary models stationary transformations of non-stationary observables are used for checking identification. [TIP: use first differences].')
if options_.diffuse_filter == 1 % use options_ as it inherits diffuse_filter from options_ident if set by user
error('There are NaN''s in the theoretical moments. Make sure that for non-stationary models stationary transformations of non-stationary observables are used for checking identification. [TIP: use first differences].')
else
error('There are NaN''s in the theoretical moments. Please check whether your model has units roots, and you forgot to set diffuse_filter=1.' )
end
end
ind_dMOMENTS = (find(max(abs(dMOMENTS'),[],1) > tol_deriv)); %index for non-zero rows
if isempty(ind_dMOMENTS) && any(any(isnan(dMOMENTS)))
error('There are NaN in the dMOMENTS matrix.' )
end
end
if ~no_identification_spectrum
ind_dSPECTRUM = (find(max(abs(dSPECTRUM'),[],1) > tol_deriv)); %index for non-zero rows
if isempty(ind_dSPECTRUM) && any(any(isnan(dSPECTRUM)))
warning_SPECTRUM = 'WARNING: There are NaN in the dSPECTRUM matrix. Please check whether your model has units roots and your forgot to set diffuse_filter=1.\n';
if any(any(isnan(dSPECTRUM)))
warning_SPECTRUM = 'WARNING: There are NaN in the dSPECTRUM matrix. Note that identification based on spectrum does not support non-stationary models (yet).\n';
warning_SPECTRUM = [warning_SPECTRUM ' Skip identification analysis based on spectrum.\n'];
fprintf(warning_SPECTRUM);
%reset options to neither display nor plot dSPECTRUM anymore
@ -168,8 +173,8 @@ if info(1) == 0 %no errors in solution
end
if ~no_identification_minimal
ind_dMINIMAL = (find(max(abs(dMINIMAL'),[],1) > tol_deriv)); %index for non-zero rows
if isempty(ind_dMINIMAL) && any(any(isnan(dMINIMAL)))
warning_MINIMAL = 'WARNING: There are NaN in the dMINIMAL matrix. Please check whether your model has units roots and you forgot to set diffuse_filter=1.\n';
if any(any(isnan(dMINIMAL)))
warning_MINIMAL = 'WARNING: There are NaN in the dMINIMAL matrix. Note that identification based on minimal system does not support non-stationary models (yet).\n';
warning_MINIMAL = [warning_MINIMAL ' Skip identification analysis based on minimal system.\n'];
fprintf(warning_MINIMAL);
%reset options to neither display nor plot dMINIMAL anymore
@ -178,8 +183,8 @@ if info(1) == 0 %no errors in solution
end
end
if no_identification_moments && no_identification_minimal && no_identification_spectrum
%display error if all three criteria fail
error('identification_analyis: Stationarity condition(s) failed and/or diffuse_filter option missing');
%display error if all three criteria fail
error(sprintf('identification_analyis: Stationarity condition(s) failed and/or diffuse_filter option missing.\nMake sure that for non-stationary models stationary transformations of non-stationary observables are used for checking identification.\n[TIP: use first differences].'));
end
% Check order conditions

4
matlab/lmmcp/catstruct.m
View File

@ -48,7 +48,7 @@ function A = catstruct(varargin)
% Copyright (C) 2005 Jos van der Geest <jos@jasen.nl>
% Copyright (C) 2013 Christophe Gouel
% Copyright (C) 2016-2017 Dynare Team
% Copyright (C) 2016-2020 Dynare Team
%
% Redistribution and use in source and binary forms, with or without
% modification, are permitted provided that the following conditions are
@ -150,7 +150,7 @@ else
FN = squeeze(FN) ;
VAL = squeeze(VAL) ;
MatlabVersion = version;
if isoctave || str2double(MatlabVersion(end-5:end-2))<2013 % Equivalent to, but faster than if verLessThan('matlab','8.1')
if isoctave && octave_ver_less_than('6')
[UFN,ind] = unique(FN) ;
else
[UFN,ind] = unique(FN,'legacy') ;

20
matlab/lmmcp/get_complementarity_conditions.m
View File

@ -1,5 +1,5 @@
function [lb,ub,eq_index] = get_complementarity_conditions(M,ramsey_policy)
% [lb,ub,eq_index] = get_complementarity_conditions(M,ramsey_policy)
function [lb,ub,eq_index] = get_complementarity_conditions(M_,ramsey_policy)
% [lb,ub,eq_index] = get_complementarity_conditions(M_,ramsey_policy)
% INPUTS
% - M [struct] contains a description of the model.
% - ramsey_policy [boolean] indicator whether a Ramsey problem is considered
@ -29,12 +29,12 @@ function [lb,ub,eq_index] = get_complementarity_conditions(M,ramsey_policy)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
ub = inf(M.endo_nbr,1);
ub = inf(M_.endo_nbr,1);
lb = -ub;
eq_index = (1:M.endo_nbr)';
eq_index = (1:M_.endo_nbr)';
if ramsey_policy
if isfield(M,'ramsey_model_constraints')
rc = M.ramsey_model_constraints;
if isfield(M_,'ramsey_model_constraints')
rc = M_.ramsey_model_constraints;
for i = 1:length(rc)
switch rc{i}{2}
case {'>','>='}
@ -48,17 +48,17 @@ if ramsey_policy
end
end
etags = M.equations_tags;
etags = M_.equations_tags;
for i=1:size(etags,1)
if strcmp(etags{i,2},'mcp')
eq_nbr = etags{i,1};
if ramsey_policy
eq_nbr = eq_nbr + M.ramsey_eq_nbr;
eq_nbr = eq_nbr + M_.ramsey_eq_nbr;
end
str = etags{i,3};
kop = strfind(etags{i,3},'<');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)), M.endo_names)); %get variable index with restriction
k = find(strcmp(strtrim(str(1:kop-1)), M_.endo_names)); %get variable index with restriction
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized'],etags{i,3},strtrim(str(1:kop-1))))
@ -69,7 +69,7 @@ for i=1:size(etags,1)
else
kop = strfind(etags{i,3},'>');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)), M.endo_names)); %get variable index with restriction
k = find(strcmp(strtrim(str(1:kop-1)), M_.endo_names)); %get variable index with restriction
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized'],etags{i,3},strtrim(str(1:kop-1))))

206
matlab/method_of_moments/method_of_moments.m
View File

@ -49,8 +49,8 @@ function [oo_, options_mom_, M_] = method_of_moments(bayestopt_, options_, oo_,
% o get_all_parameters.m
% o get_matrix_entries_for_psd_check.m
% o makedataset.m
% o method_of_moments_check_plot.m
% o method_of_moments_data_moments.m
% o method_of_moments_mode_check.m
% o method_of_moments_objective_function.m
% o method_of_moments_optimal_weighting_matrix
% o method_of_moments_standard_errors
@ -63,7 +63,7 @@ function [oo_, options_mom_, M_] = method_of_moments(bayestopt_, options_, oo_,
% o set_all_parameters.m
% o test_for_deep_parameters_calibration.m
% =========================================================================
% Copyright (C) 2020 Dynare Team
% Copyright (C) 2020-2021 Dynare Team
%
% This file is part of Dynare.
%
@ -86,21 +86,13 @@ function [oo_, options_mom_, M_] = method_of_moments(bayestopt_, options_, oo_,
% =========================================================================
%% TO DO LIST
% - [ ] why does lsqnonlin take less time in Andreasen toolbox?
% - [ ] test user-specified weightning matrix
% - [ ] which qz_criterium value?
% - [ ] document that in method_of_moments_data_moments.m NaN are replaced by mean of moment
% - [ ] add IRF matching
% - [ ] test estimated_params_bounds block
% - [ ] test what happens if all parameters will be estimated but some/all are not calibrated
% - [ ] speed up lyapunov equation by using doubling with old initial values
% - [ ] check smm at order > 3 without pruning
% - [ ] provide option to use analytical derivatives to compute std errors (similar to what we already do in identification)
% - [ ] add Bayesian GMM/SMM estimation
% - [ ] useautocorr
% - [ ] do we need dirname?
% - [ ] decide on default weighting matrix scheme, I would propose 2 stage with Diagonal of optimal matrix
% - [ ] check smm with product moments greater than 2
% - [ ] speed up pruned_state_space_system (by using doubling with old initial values, hardcoding zeros, other "tricks" used in e.g. nlma)
% - [ ] add option to use autocorrelations (we have useautocorr in identification toolbox already)
% - [ ] SMM with extended path
% - [ ] deal with measurement errors (once @wmutschl has implemented this in identification toolbox)
% - [ ] improve check for duplicate moments by using the cellfun and unique functions
% - [ ] dirname option to save output to different directory not yet implemented
% -------------------------------------------------------------------------
% Step 0: Check if required structures and options exist
% -------------------------------------------------------------------------
@ -125,7 +117,7 @@ else
options_mom_.loglinear = false;
end
fprintf('\n==== Method of Moments (%s) Estimation ====\n\n',options_mom_.mom.mom_method)
fprintf('\n==== Method of Moments Estimation (%s) ====\n\n',options_mom_.mom.mom_method)
% -------------------------------------------------------------------------
% Step 1a: Prepare options_mom_ structure
@ -142,15 +134,15 @@ if strcmp(options_mom_.mom.mom_method,'GMM') || strcmp(options_mom_.mom.mom_meth
options_mom_.mom = set_default_option(options_mom_.mom,'bartlett_kernel_lag',20); % bandwith in optimal weighting matrix
options_mom_.mom = set_default_option(options_mom_.mom,'penalized_estimator',false); % include deviation from prior mean as additional moment restriction and use prior precision as weight
options_mom_.mom = set_default_option(options_mom_.mom,'verbose',false); % display and store intermediate estimation results
options_mom_.mom = set_default_option(options_mom_.mom,'weighting_matrix',{'DIAGONAL'; 'DIAGONAL'}); % weighting matrix in moments distance objective function at each iteration of estimation; cell of strings with
options_mom_.mom = set_default_option(options_mom_.mom,'weighting_matrix',{'DIAGONAL'; 'DIAGONAL'}); % weighting matrix in moments distance objective function at each iteration of estimation;
% possible values are 'OPTIMAL', 'IDENTITY_MATRIX' ,'DIAGONAL' or a filename. Size of cell determines stages in iterated estimation.
options_mom_.mom = set_default_option(options_mom_.mom,'weighting_matrix_scaling_factor',1); % scaling of weighting matrix
options_mom_.mom = set_default_option(options_mom_.mom,'weighting_matrix_scaling_factor',1); % scaling of weighting matrix in objective function
options_mom_.mom = set_default_option(options_mom_.mom,'se_tolx',1e-5); % step size for numerical computation of standard errors
options_mom_ = set_default_option(options_mom_,'order',1); % order of Taylor approximation in perturbation
options_mom_ = set_default_option(options_mom_,'pruning',true); % use pruned state space system at higher-order
options_mom_ = set_default_option(options_mom_,'pruning',false); % use pruned state space system at higher-order
% Checks for perturbation order
if options_mom_.order < 1
error('method_of_moments:: The order of the Taylor approximation cannot be 0!')
error('method_of_moments: The order of the Taylor approximation cannot be 0!')
end
end
if strcmp(options_mom_.mom.mom_method,'SMM')
@ -169,11 +161,16 @@ if strcmp(options_mom_.mom.mom_method,'GMM')
fprintf('GMM at higher order only works with pruning, so we set pruning option to 1.\n');
options_mom_.pruning = true;
end
if options_mom_.order > 3
error('method_of_moments: perturbation orders higher than 3 are not implemented for GMM estimation, try using SMM.\n');
end
options_mom_.mom = set_default_option(options_mom_.mom,'analytic_standard_errors',false); % compute standard errors numerically (0) or analytically (1). Analytical derivatives are only available for GMM.
end
options_mom_.mom.compute_derivs = false;% flag to compute derivs in objective function (needed for analytic standard errors with GMM)
% General options that can be set by the user in the mod file, otherwise default values are provided
options_mom_ = set_default_option(options_mom_,'dirname',M_.fname); % directory in which to store estimation output
options_mom_ = set_default_option(options_mom_,'dirname',M_.dname); % specify directory in which to store estimation output [not yet working]
options_mom_ = set_default_option(options_mom_,'graph_format','eps'); % specify the file format(s) for graphs saved to disk
options_mom_ = set_default_option(options_mom_,'nodisplay',false); % do not display the graphs, but still save them to disk
options_mom_ = set_default_option(options_mom_,'nograph',false); % do not create graphs (which implies that they are not saved to the disk nor displayed)
@ -203,12 +200,12 @@ options_mom_ = set_default_option(options_mom_,'mode_compute',13);
options_mom_ = set_default_option(options_mom_,'additional_optimizer_steps',[]); % vector of additional mode-finders run after mode_compute
options_mom_ = set_default_option(options_mom_,'optim_opt',[]); % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute
options_mom_ = set_default_option(options_mom_,'silent_optimizer',false); % run minimization of moments distance silently without displaying results or saving files in between
% Mode_check plot options that can be set by the user in the mod file, otherwise default values are provided
% Check plot options that can be set by the user in the mod file, otherwise default values are provided
options_mom_.mode_check.nolik = false; % we don't do likelihood (also this initializes mode_check substructure)
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'status',false); % plot the target function for values around the computed mode for each estimated parameter in turn. This is helpful to diagnose problems with the optimizer.
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'neighbourhood_size',.5); % width of the window around the mode to be displayed on the diagnostic plots. This width is expressed in percentage deviation. The Inf value is allowed, and will trigger a plot over the entire domain
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'symmetric_plots',true); % ensure that the check plots are symmetric around the mode. A value of 0 allows to have asymmetric plots, which can be useful if the posterior mode is close to a domain boundary, or in conjunction with mode_check_neighbourhood_size = Inf when the domain is not the entire real line
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'number_of_points',20); % number of points around the mode where the target function is evaluated (for each parameter)
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'status',false); % plot the target function for values around the computed minimum for each estimated parameter in turn. This is helpful to diagnose problems with the optimizer.
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'neighbourhood_size',.5); % width of the window around the computed minimum to be displayed on the diagnostic plots. This width is expressed in percentage deviation. The Inf value is allowed, and will trigger a plot over the entire domain
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'symmetric_plots',true); % ensure that the check plots are symmetric around the minimum. A value of 0 allows to have asymmetric plots, which can be useful if the minimum is close to a domain boundary, or in conjunction with neighbourhood_size = Inf when the domain is not the entire real line
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'number_of_points',20); % number of points around the minimum where the target function is evaluated (for each parameter)
% Numerical algorithms options that can be set by the user in the mod file, otherwise default values are provided
options_mom_ = set_default_option(options_mom_,'aim_solver',false); % use AIM algorithm to compute perturbation approximation instead of mjdgges
@ -226,8 +223,11 @@ options_mom_ = set_default_option(options_mom_,'lyapunov_fixed_point_tol',1e-10)
options_mom_ = set_default_option(options_mom_,'lyapunov_doubling_tol',1e-16); % convergence criterion used in the doubling algorithm
options_mom_ = set_default_option(options_mom_,'sylvester_fp',false); % determines whether to use fixed point algorihtm to solve Sylvester equation (gensylv_fp), faster for large scale models
options_mom_ = set_default_option(options_mom_,'sylvester_fixed_point_tol',1e-12); % convergence criterion used in the fixed point Sylvester solver
options_mom_ = set_default_option(options_mom_,'qz_criterium',1-1e-6); % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems [IS THIS CORRET @wmutschl]
options_mom_ = set_default_option(options_mom_,'qz_criterium',1-1e-6); % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems
% if there are no unit roots one can use 1.0 (or slightly below) which we set as default; if they are possible, you may have have multiple unit roots and the accuracy decreases when computing the eigenvalues in lyapunov_symm
% Note that unit roots are only possible at first-order, at higher order we set it to 1 in pruned_state_space_system and focus only on stationary observables.
options_mom_ = set_default_option(options_mom_,'qz_zero_threshold',1e-6); % value used to test if a generalized eigenvalue is 0/0 in the generalized Schur decomposition
options_mom_ = set_default_option(options_mom_,'schur_vec_tol',1e-11); % tolerance level used to find nonstationary variables in Schur decomposition of the transition matrix.
if options_mom_.order > 2
fprintf('Dynare will use ''k_order_solver'' as the order>2\n');
options_mom_.k_order_solver = true;
@ -326,6 +326,7 @@ options_mom_.solveopt = options_.solveopt;
options_mom_.gradient_method = options_.gradient_method;
options_mom_.gradient_epsilon = options_.gradient_epsilon;
options_mom_.analytic_derivation = 0;
options_mom_.analytic_derivation_mode = 0; % needed by get_perturbation_params_derivs.m, ie use efficient sylvester equation method to compute analytical derivatives as in Ratto & Iskrev (2012)
options_mom_.vector_output= false; % specifies whether the objective function returns a vector
@ -362,8 +363,6 @@ end
% -------------------------------------------------------------------------
% Step 2: Checks and transformations for matched moments structure (preliminary)
% -------------------------------------------------------------------------
% Note that we do not have a preprocessor interface yet for this, so this
% will need much improvement later on. @wmutschl
% Initialize indices
options_mom_.mom.index.E_y = false(options_mom_.obs_nbr,1); %unconditional first order product moments
@ -423,8 +422,6 @@ for jm=1:size(M_.matched_moments,1)
end
end
% @wmutschl: add check for duplicate moments by using the cellfun and unique functions
%Remove duplicate elements
UniqueMomIdx = [nonzeros(options_mom_.mom.index.E_y_pos); nonzeros(tril(options_mom_.mom.index.E_yy_pos)); nonzeros(options_mom_.mom.index.E_yyt_pos)];
DuplicateMoms = setdiff(1:size(M_.matched_moments,1),UniqueMomIdx);
@ -614,16 +611,14 @@ if ~isempty(dataset_)
options_mom_.nobs = dataset_.nobs;
end
% provide info on missing observations
if any(any(isnan(dataset_.data)))
fprintf('missing observations will be replaced by the sample mean of the corresponding moment')
end
% Check length of data for estimation of second moments
if options_mom_.ar > options_mom_.nobs+1
error('method_of_moments: Data set is too short to compute second moments');
end
% Provide info on data moments handling
fprintf('Computing data moments. Note that NaN values in the moments (due to leads and lags or missing data) are replaced by the mean of the corresponding moment\n');
% Get data moments for the method of moments
[oo_.mom.data_moments, oo_.mom.m_data] = method_of_moments_data_moments(dataset_.data, oo_, M_.matched_moments, options_mom_);
@ -632,9 +627,17 @@ if strcmp(options_mom_.mom.mom_method,'SMM')
options_mom_.mom.long = round(options_mom_.mom.simulation_multiple*options_mom_.nobs);
options_mom_.mom.variance_correction_factor = (1+1/options_mom_.mom.simulation_multiple);
% draw shocks for SMM
smmstream = RandStream('mt19937ar','Seed',options_mom_.mom.seed);
temp_shocks = randn(smmstream,options_mom_.mom.long+options_mom_.mom.burnin,M_.exo_nbr);
temp_shocks_ME = randn(smmstream,options_mom_.mom.long,length(M_.H));
if ~isoctave
smmstream = RandStream('mt19937ar','Seed',options_mom_.mom.seed);
temp_shocks = randn(smmstream,options_mom_.mom.long+options_mom_.mom.burnin,M_.exo_nbr);
temp_shocks_ME = randn(smmstream,options_mom_.mom.long,length(M_.H));
else
[state_u,state_n] = get_dynare_random_generator_state; %get state for later resetting
set_dynare_random_generator_state(options_mom_.mom.seed,options_mom_.mom.seed);
temp_shocks = randn(options_mom_.mom.long+options_mom_.mom.burnin,M_.exo_nbr);
temp_shocks_ME = randn(options_mom_.mom.long,length(M_.H));
set_dynare_random_generator_state(state_u,state_n); %reset state for later resetting
end
if options_mom_.mom.bounded_shock_support == 1
temp_shocks(temp_shocks>2) = 2;
temp_shocks(temp_shocks<-2) = -2;
@ -643,6 +646,9 @@ if strcmp(options_mom_.mom.mom_method,'SMM')
end
options_mom_.mom.shock_series = temp_shocks;
options_mom_.mom.ME_shock_series = temp_shocks_ME;
if options_mom_.k_order_solver && ~options_mom_.pruning % dynare++ routines will be called in simult_.m, store some additional stuff
options_mom_.DynareRandomStreams.seed = options_mom_.mom.seed;
end
end
% -------------------------------------------------------------------------
@ -680,6 +686,11 @@ if isfield(estim_params_,'param_vals') && ~isempty(estim_params_.param_vals)
fprintf('This will override parameter values and may lead to wrong results.\n')
fprintf('Check whether this is really intended.\n')
warning('The steady state file internally changes the values of the estimated parameters.')
if strcmp(options_mom_.mom.mom_method,'GMM') && options_mom_.mom.analytic_standard_errors
fprintf('For analytical standard errors, the parameter-Jacobians of the dynamic model and of the steady-state will be computed numerically\n'),
fprintf('(re-set options_mom_.analytic_derivation_mode= -2)'),
options_mom_.analytic_derivation_mode= -2;
end
end
end
@ -729,11 +740,6 @@ catch last_error% if check fails, provide info on using calibration if present
rethrow(last_error);
end
if options_mom_.mode_compute == 0 %We only report value of moments distance at initial value of the parameters
fprintf('No minimization of moments distance due to ''mode_compute=0''\n')
return
end
% -------------------------------------------------------------------------
% Step 7a: Method of moments estimation: print some info
% -------------------------------------------------------------------------
@ -751,34 +757,69 @@ end
if options_mom_.mom.penalized_estimator
fprintf('\n - penalized estimation using deviation from prior mean and weighted with prior precision');
end
if options_mom_.mode_compute == 1; fprintf('\n - optimizer (mode_compute=1): fmincon');
elseif options_mom_.mode_compute == 2; fprintf('\n - optimizer (mode_compute=2): continuous simulated annealing');
elseif options_mom_.mode_compute == 3; fprintf('\n - optimizer (mode_compute=3): fminunc');
elseif options_mom_.mode_compute == 4; fprintf('\n - optimizer (mode_compute=4): csminwel');
elseif options_mom_.mode_compute == 5; fprintf('\n - optimizer (mode_compute=5): newrat');
elseif options_mom_.mode_compute == 6; fprintf('\n - optimizer (mode_compute=6): gmhmaxlik');
elseif options_mom_.mode_compute == 7; fprintf('\n - optimizer (mode_compute=7): fminsearch');
elseif options_mom_.mode_compute == 8; fprintf('\n - optimizer (mode_compute=8): Dynare Nelder-Mead simplex');
elseif options_mom_.mode_compute == 9; fprintf('\n - optimizer (mode_compute=9): CMA-ES');
elseif options_mom_.mode_compute == 10; fprintf('\n - optimizer (mode_compute=10): simpsa');
elseif options_mom_.mode_compute == 11; fprintf('\n - optimizer (mode_compute=11): online_auxiliary_filter');
elseif options_mom_.mode_compute == 12; fprintf('\n - optimizer (mode_compute=12): particleswarm');
elseif options_mom_.mode_compute == 101; fprintf('\n - optimizer (mode_compute=101): SolveOpt');
elseif options_mom_.mode_compute == 102; fprintf('\n - optimizer (mode_compute=102): simulannealbnd');
elseif options_mom_.mode_compute == 13; fprintf('\n - optimizer (mode_compute=13): lsqnonlin');
elseif ischar(minimizer_algorithm); fprintf(['\n - user-defined optimizer: ' minimizer_algorithm]);
else
error('method_of_moments: Unknown optimizer, please contact the developers ')
end
if options_mom_.silent_optimizer
fprintf(' (silent)');
optimizer_vec=[options_mom_.mode_compute;num2cell(options_mom_.additional_optimizer_steps)]; % at each stage one can possibly use different optimizers sequentially
for i = 1:length(optimizer_vec)
if i == 1
str = '- optimizer (mode_compute';
else
str = ' (additional_optimizer_steps';
end
switch optimizer_vec{i}
case 0
fprintf('\n %s=0): no minimization',str);
case 1
fprintf('\n %s=1): fmincon',str);
case 2
fprintf('\n %s=2): continuous simulated annealing',str);
case 3
fprintf('\n %s=3): fminunc',str);
case 4
fprintf('\n %s=4): csminwel',str);
case 5
fprintf('\n %s=5): newrat',str);
case 6
fprintf('\n %s=6): gmhmaxlik',str);
case 7
fprintf('\n %s=7): fminsearch',str);
case 8
fprintf('\n %s=8): Dynare Nelder-Mead simplex',str);
case 9
fprintf('\n %s=9): CMA-ES',str);
case 10
fprintf('\n %s=10): simpsa',str);
case 11
fprintf('\n %s=11): online_auxiliary_filter',str);
case 12
fprintf('\n %s=12): particleswarm',str);
case 101
fprintf('\n %s=101): SolveOpt',str);
case 102
fprintf('\n %s=102): simulannealbnd',str);
case 13
fprintf('\n %s=13): lsqnonlin',str);
otherwise
if ischar(optimizer_vec{i})
fprintf('\n %s=%s): user-defined',str,optimizer_vec{i});
else
error('method_of_moments: Unknown optimizer, please contact the developers ')
end
end
if options_mom_.silent_optimizer
fprintf(' (silent)');
end
end
fprintf('\n - perturbation order: %d', options_mom_.order)
if options_mom_.order > 1 && options_mom_.pruning
fprintf(' (with pruning)')
end
if strcmp(options_mom_.mom.mom_method,'GMM') && options_mom_.mom.analytic_standard_errors
fprintf('\n - standard errors: analytic derivatives');
else
fprintf('\n - standard errors: numerical derivatives');
end
fprintf('\n - number of matched moments: %d', options_mom_.mom.mom_nbr);
fprintf('\n - number of parameters: %d\n\n', length(xparam0));
fprintf('\n - number of parameters: %d', length(xparam0));
fprintf('\n\n');
% -------------------------------------------------------------------------
% Step 7b: Iterated method of moments estimation
@ -787,8 +828,6 @@ if size(options_mom_.mom.weighting_matrix,1)>1 && ~(any(strcmpi('diagonal',optio
fprintf('\nYou did not specify the use of an optimal or diagonal weighting matrix. There is no point in running an iterated method of moments.\n')
end
optimizer_vec=[options_mom_.mode_compute,options_mom_.additional_optimizer_steps]; % at each stage one can possibly use different optimizers sequentially
for stage_iter=1:size(options_mom_.mom.weighting_matrix,1)
fprintf('Estimation stage %u\n',stage_iter);
Woptflag = false;
@ -834,15 +873,20 @@ for stage_iter=1:size(options_mom_.mom.weighting_matrix,1)
end
for optim_iter= 1:length(optimizer_vec)
if optimizer_vec(optim_iter)==13
options_mom_.vector_output = true;
if optimizer_vec{optim_iter}==0
xparam1=xparam0; %no minimization, evaluate objective at current values
fval = feval(objective_function, xparam1, Bounds, oo_, estim_params_, M_, options_mom_);
else
options_mom_.vector_output = false;
end
[xparam1, fval, exitflag] = dynare_minimize_objective(objective_function, xparam0, optimizer_vec(optim_iter), options_mom_, [Bounds.lb Bounds.ub], bayestopt_laplace.name, bayestopt_laplace, [],...
Bounds, oo_, estim_params_, M_, options_mom_);
if options_mom_.vector_output
fval = fval'*fval;
if optimizer_vec{optim_iter}==13
options_mom_.vector_output = true;
else
options_mom_.vector_output = false;
end
[xparam1, fval, exitflag] = dynare_minimize_objective(objective_function, xparam0, optimizer_vec{optim_iter}, options_mom_, [Bounds.lb Bounds.ub], bayestopt_laplace.name, bayestopt_laplace, [],...
Bounds, oo_, estim_params_, M_, options_mom_);
if options_mom_.vector_output
fval = fval'*fval;
end
end
fprintf('\nStage %d Iteration %d: value of minimized moment distance objective function: %12.10f.\n',stage_iter,optim_iter,fval)
if options_mom_.mom.verbose
@ -853,8 +897,12 @@ for stage_iter=1:size(options_mom_.mom.weighting_matrix,1)
options_mom_.vector_output = false;
% Update M_ and DynareResults (in particular to get oo_.mom.model_moments)
M_ = set_all_parameters(xparam1,estim_params_,M_);
if strcmp(options_mom_.mom.mom_method,'GMM') && options_mom_.mom.analytic_standard_errors
options_mom_.mom.compute_derivs = true; % for GMM we compute derivatives analytically in the objective function with this flag
end
[fval, ~, ~,~,~, oo_] = feval(objective_function, xparam1, Bounds, oo_, estim_params_, M_, options_mom_);
% Compute Standard errors
options_mom_.mom.compute_derivs = false; % reset to not compute derivatives in objective function during optimization
SE = method_of_moments_standard_errors(xparam1, objective_function, Bounds, oo_, estim_params_, M_, options_mom_, Woptflag);
% Store results in output structure
@ -901,7 +949,7 @@ if options_mom_.TeX
end
if options_mom_.mode_check.status
method_of_moments_mode_check(objective_function,xparam1,SE,options_mom_,M_,estim_params_,Bounds,bayestopt_laplace,...
method_of_moments_check_plot(objective_function,xparam1,SE,options_mom_,M_,estim_params_,Bounds,bayestopt_laplace,...
Bounds, oo_, estim_params_, M_, options_mom_)
end

35
matlab/method_of_moments/method_of_moments_mode_check.m → matlab/method_of_moments/method_of_moments_check_plot.m
View File

@ -1,8 +1,8 @@
function method_of_moments_mode_check(fun,xparam,SE_vec,options_,M_,estim_params_,Bounds,bayestopt_,varargin)
% Checks the estimated ML mode or Posterior mode.
function method_of_moments_check_plot(fun,xparam,SE_vec,options_,M_,estim_params_,Bounds,bayestopt_,varargin)
% Checks the estimated local minimum of the moment's distance objective
% Copyright (C) 2020 Dynare Team
% Copyright (C) 2020-2021 Dynare Team
%
% This file is part of Dynare.
%
@ -28,7 +28,7 @@ fval = feval(fun,xparam,varargin{:});
if ~isempty(SE_vec)
skipline()
disp('MODE CHECK')
disp('LOCAL MINIMUM CHECK')
skipline()
fprintf('Fval obtained by the minimization routine: %f', fval);
skipline()
@ -39,12 +39,12 @@ end
[nbplt,nr,nc,lr,lc,nstar] = pltorg(length(xparam));
if ~exist([M_.fname filesep 'graphs'],'dir')
mkdir(M_.fname,'graphs');
if ~exist([M_.dname filesep 'graphs'],'dir')
mkdir(M_.dname,'graphs');
end
if TeX && any(strcmp('eps',cellstr(options_.graph_format)))
fidTeX = fopen([M_.fname, '/graphs/', M_.fname '_MoMCheckPlots.tex'],'w');
fprintf(fidTeX,'%% TeX eps-loader file generated by method_of_moments_mode_check.m (Dynare).\n');
fidTeX = fopen([M_.dname, '/graphs/', M_.fname '_MoMCheckPlots.tex'],'w');
fprintf(fidTeX,'%% TeX eps-loader file generated by method_of_moments_check_plot.m (Dynare).\n');
fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
fprintf(fidTeX,' \n');
end
@ -54,14 +54,14 @@ if isinf(ll)
options_.mode_check.symmetric_plots = false;
end
mcheck = struct('cross',struct(),'emode',struct());
mcheck = struct('cross',struct(),'emin',struct());
for plt = 1:nbplt
if TeX
NAMES = [];
TeXNAMES = [];
end
hh = dyn_figure(options_.nodisplay,'Name','Mode check plots');
hh = dyn_figure(options_.nodisplay,'Name','Minimum check plots');
for k=1:min(nstar,length(xparam)-(plt-1)*nstar)
subplot(nr,nc,k)
kk = (plt-1)*nstar+k;
@ -119,7 +119,7 @@ for plt = 1:nbplt
else
y(i,1) = NaN;
if options_.debug
fprintf('mode_check:: could not solve model for parameter %s at value %4.3f, error code: %u\n',name,z(i),info(1))
fprintf('method_of_moments_check_plot:: could not solve model for parameter %s at value %4.3f, error code: %u\n',name,z(i),info(1))
end
end
if options_.mom.penalized_estimator
@ -127,9 +127,9 @@ for plt = 1:nbplt
y(i,2) = (y(i,1)+prior-dy);
end
end
mcheck.cross = setfield(mcheck.cross, name, [transpose(z), -y]);
mcheck.emode = setfield(mcheck.emode, name, xparam(kk));
fighandle=plot(z,-y);
mcheck.cross = setfield(mcheck.cross, name, [transpose(z), y]);
mcheck.emin = setfield(mcheck.emin, name, xparam(kk));
fighandle=plot(z,y);
hold on
yl=get(gca,'ylim');
plot( [xparam(kk) xparam(kk)], yl, 'c', 'LineWidth', 1)
@ -165,12 +165,12 @@ for plt = 1:nbplt
text(0.25,0.5,'log-post')
text(0.69,0.5,'log-lik kernel')
end
dyn_saveas(hh,[M_.fname, '/graphs/', M_.fname '_MoMCheckPlots' int2str(plt) ],options_.nodisplay,options_.graph_format);
dyn_saveas(hh,[M_.dname, '/graphs/', M_.fname '_MoMCheckPlots' int2str(plt) ],options_.nodisplay,options_.graph_format);
if TeX && any(strcmp('eps',cellstr(options_.graph_format)))
% TeX eps loader file
fprintf(fidTeX,'\\begin{figure}[H]\n');
fprintf(fidTeX,'\\centering \n');
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%_MoMCheckPlots%s}\n',options_.figures.textwidth*min(k/nc,1),[M_.fname, '/graphs/',M_.fname],int2str(plt));
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_MoMCheckPlots%s}\n',options_.figures.textwidth*min(k/nc,1),[M_.dname, '/graphs/',M_.fname],int2str(plt));
fprintf(fidTeX,'\\caption{Method of Moments check plots.}');
fprintf(fidTeX,'\\label{Fig:MoMCheckPlots:%s}\n',int2str(plt));
fprintf(fidTeX,'\\end{figure}\n');
@ -181,5 +181,4 @@ if TeX && any(strcmp('eps',cellstr(options_.graph_format)))
fclose(fidTeX);
end
OutputDirectoryName = CheckPath('modecheck',M_.dname);
save([OutputDirectoryName '/MoM_check_plot_data.mat'],'mcheck');
save([M_.dname filesep 'graphs' filesep M_.fname '_MoMCheckPlots_data.mat'],'mcheck');

8
matlab/method_of_moments/method_of_moments_data_moments.m
View File

@ -16,7 +16,7 @@ function [dataMoments, m_data] = method_of_moments_data_moments(data, oo_, match
% o method_of_moments.m
% o method_of_moments_objective_function.m
% =========================================================================
% Copyright (C) 2020 Dynare Team
% Copyright (C) 2020-2021 Dynare Team
%
% This file is part of Dynare.
%
@ -59,7 +59,11 @@ for jm = 1:options_mom_.mom.mom_nbr
end
end
% We replace NaN (due to leads and lags and missing values) with the corresponding mean
dataMoments(jm,1) = mean(m_data_tmp,'omitnan');
if isoctave
dataMoments(jm,1) = nanmean(m_data_tmp);
else
dataMoments(jm,1) = mean(m_data_tmp,'omitnan');
end
m_data_tmp(isnan(m_data_tmp)) = dataMoments(jm,1);
m_data(:,jm) = m_data_tmp;
end

59
matlab/method_of_moments/method_of_moments_objective_function.m
View File

@ -31,7 +31,7 @@ function [fval, info, exit_flag, junk1, junk2, oo_, M_, options_mom_] = method_o
% o resol
% o set_all_parameters
% =========================================================================
% Copyright (C) 2020 Dynare Team
% Copyright (C) 2020-2021 Dynare Team
%
% This file is part of Dynare.
%
@ -109,7 +109,32 @@ if strcmp(options_mom_.mom.mom_method,'GMM')
%--------------------------------------------------------------------------
% 3. Set up pruned state-space system and compute model moments
%--------------------------------------------------------------------------
pruned_state_space = pruned_state_space_system(M_, options_mom_, dr, oo_.dr.obs_var, options_mom_.ar, 0, 0);
if options_mom_.mom.compute_derivs && options_mom_.mom.analytic_standard_errors
indpmodel = []; %initialize index for model parameters
if ~isempty(estim_params_.param_vals)
indpmodel = estim_params_.param_vals(:,1); %values correspond to parameters declaration order, row number corresponds to order in estimated_params
end
indpstderr=[]; %initialize index for stderr parameters
if ~isempty(estim_params_.var_exo)
indpstderr = estim_params_.var_exo(:,1); %values correspond to varexo declaration order, row number corresponds to order in estimated_params
end
indpcorr=[]; %initialize matrix for corr paramters
if ~isempty(estim_params_.corrx)
indpcorr = estim_params_.corrx(:,1:2); %values correspond to varexo declaration order, row number corresponds to order in estimated_params
end
if estim_params_.nvn || estim_params_.ncn %nvn is number of stderr parameters and ncn is number of corr parameters of measurement innovations as declared in estimated_params
error('Analytic computation of standard errrors does not (yet) support measurement errors.\nInstead, define them explicitly as varexo and provide measurement equations in the model definition.\nAlternatively, use numerical standard errors.')
end
modparam_nbr = estim_params_.np; % number of model parameters as declared in estimated_params
stderrparam_nbr = estim_params_.nvx; % number of stderr parameters
corrparam_nbr = estim_params_.ncx; % number of corr parameters
totparam_nbr = stderrparam_nbr+corrparam_nbr+modparam_nbr;
dr.derivs = get_perturbation_params_derivs(M_, options_mom_, estim_params_, oo_, indpmodel, indpstderr, indpcorr, 0); %analytic derivatives of perturbation matrices
oo_.mom.model_moments_params_derivs = NaN(options_mom_.mom.mom_nbr,totparam_nbr);
pruned_state_space = pruned_state_space_system(M_, options_mom_, dr, oo_.dr.obs_var, options_mom_.ar, 0, 1);
else
pruned_state_space = pruned_state_space_system(M_, options_mom_, dr, oo_.dr.obs_var, options_mom_.ar, 0, 0);
end
oo_.mom.model_moments = NaN(options_mom_.mom.mom_nbr,1);
offset = 0;
@ -118,6 +143,9 @@ if strcmp(options_mom_.mom.mom_method,'GMM')
E_y = pruned_state_space.E_y;
E_y_nbr = nnz(options_mom_.mom.index.E_y);
oo_.mom.model_moments(offset+1:E_y_nbr,1) = E_y(options_mom_.mom.index.E_y);
if options_mom_.mom.compute_derivs && options_mom_.mom.analytic_standard_errors
oo_.mom.model_moments_params_derivs(offset+1:E_y_nbr,:) = pruned_state_space.dE_y(options_mom_.mom.index.E_y,:);
end
offset = offset + E_y_nbr;
end
% Second moments
@ -125,22 +153,47 @@ if strcmp(options_mom_.mom.mom_method,'GMM')
if isfield(options_mom_.mom.index,'E_yy') && nnz(options_mom_.mom.index.E_yy) > 0
if options_mom_.prefilter
E_yy = pruned_state_space.Var_y;
if options_mom_.mom.compute_derivs && options_mom_.mom.analytic_standard_errors
dE_yy = pruned_state_space.dVar_y;
end
else
E_yy = pruned_state_space.Var_y + pruned_state_space.E_y*pruned_state_space.E_y';
if options_mom_.mom.compute_derivs && options_mom_.mom.analytic_standard_errors
dE_yy = pruned_state_space.dVar_y;
for jp=1:totparam_nbr
dE_yy(:,:,jp) = dE_yy(:,:,jp) + pruned_state_space.dE_y(:,jp)*pruned_state_space.E_y' + pruned_state_space.E_y*pruned_state_space.dE_y(:,jp)';
end
end
end
E_yy_nbr = nnz(tril(options_mom_.mom.index.E_yy));
oo_.mom.model_moments(offset+(1:E_yy_nbr),1) = E_yy(tril(options_mom_.mom.index.E_yy));
if options_mom_.mom.compute_derivs && options_mom_.mom.analytic_standard_errors
oo_.mom.model_moments_params_derivs(offset+(1:E_yy_nbr),:) = reshape(dE_yy(repmat(tril(options_mom_.mom.index.E_yy),[1 1 totparam_nbr])),E_yy_nbr,totparam_nbr);
end
offset = offset + E_yy_nbr;
end
% Lead/lags covariance
if isfield(options_mom_.mom.index,'E_yyt') && nnz(options_mom_.mom.index.E_yyt) > 0
if options_mom_.prefilter
E_yyt = pruned_state_space.Var_yi;
if options_mom_.mom.compute_derivs && options_mom_.mom.analytic_standard_errors
dE_yyt = pruned_state_space.dVar_yi;
end
else
E_yyt = pruned_state_space.Var_yi + repmat(pruned_state_space.E_y*pruned_state_space.E_y',[1 1 size(pruned_state_space.Var_yi,3)]);
if options_mom_.mom.compute_derivs && options_mom_.mom.analytic_standard_errors
dE_yyt = pruned_state_space.dVar_yi;
for jp=1:totparam_nbr
dE_yyt(:,:,:,jp) = dE_yyt(:,:,:,jp) + repmat(pruned_state_space.dE_y(:,jp)*pruned_state_space.E_y',[1 1 size(pruned_state_space.Var_yi,3)])...
+ repmat(pruned_state_space.E_y*pruned_state_space.dE_y(:,jp)',[1 1 size(pruned_state_space.Var_yi,3)]);
end
end
end
E_yyt_nbr = nnz(options_mom_.mom.index.E_yyt);
oo_.mom.model_moments(offset+(1:E_yyt_nbr),1) = E_yyt(options_mom_.mom.index.E_yyt);
if options_mom_.mom.compute_derivs && options_mom_.mom.analytic_standard_errors
oo_.mom.model_moments_params_derivs(offset+(1:E_yyt_nbr),:) = reshape(dE_yyt(repmat(options_mom_.mom.index.E_yyt,[1 1 1 totparam_nbr])),E_yyt_nbr,totparam_nbr);
end
end
elseif strcmp(options_mom_.mom.mom_method,'SMM')
@ -179,7 +232,7 @@ elseif strcmp(options_mom_.mom.mom_method,'SMM')
i_ME = setdiff([1:size(M_.H,1)],find(diag(M_.H) == 0)); % find ME with 0 variance
chol_S = chol(M_.H(i_ME,i_ME)); %decompose rest
shock_mat=zeros(size(options_mom_.mom.ME_shock_series)); %initialize
shock_mat(:,i_ME)=options_mom_.mom.ME_shock_series(:,i_exo_var)*chol_S;
shock_mat(:,i_ME)=options_mom_.mom.ME_shock_series(:,i_ME)*chol_S;
y_sim = y_sim+shock_mat;
end

2
matlab/method_of_moments/method_of_moments_optimal_weighting_matrix.m
View File

@ -19,7 +19,7 @@ function W_opt = method_of_moments_optimal_weighting_matrix(m_data, moments, q_l
% This function calls:
% o CorrMatrix (embedded)
% =========================================================================
% Copyright (C) 2020 Dynare Team
% Copyright (C) 2020-2021 Dynare Team
%
% This file is part of Dynare.
%

58
matlab/method_of_moments/method_of_moments_standard_errors.m
View File

@ -29,7 +29,7 @@ function [SE_values, Asympt_Var] = method_of_moments_standard_errors(xparam, obj
% o SMM_objective_function.m
% o method_of_moments_optimal_weighting_matrix
% =========================================================================
% Copyright (C) 2020 Dynare Team
% Copyright (C) 2020-2021 Dynare Team
%
% This file is part of Dynare.
%
@ -57,30 +57,46 @@ dim_params = size(xparam,1);
D = zeros(num_mom,dim_params);
eps_value = options_mom_.mom.se_tolx;
for i=1:dim_params
%Positive step
xparam_eps_p = xparam;
xparam_eps_p(i,1) = xparam_eps_p(i) + eps_value;
[~, info_p, ~, ~,~, oo__p] = feval(objective_function, xparam_eps_p, Bounds, oo_, estim_params_, M_, options_mom_);
% Negative step
xparam_eps_m = xparam;
xparam_eps_m(i,1) = xparam_eps_m(i) - eps_value;
[~, info_m, ~, ~,~, oo__m] = feval(objective_function, xparam_eps_m, Bounds, oo_, estim_params_, M_, options_mom_);
% The Jacobian:
if nnz(info_p)==0 && nnz(info_m)==0
D(:,i) = (oo__p.mom.model_moments - oo__m.mom.model_moments)/(2*eps_value);
else
problpar = get_the_name(i,options_mom_.TeX, M_, estim_params_, options_mom_);
message_p = get_error_message(info_p, options_mom_);
message_m = get_error_message(info_m, options_mom_);
warning('method_of_moments:info','Cannot compute the Jacobian for parameter %s - no standard errors available\n %s %s\nCheck your bounds and/or priors, or use a different optimizer.\n',problpar, message_p, message_m)
if strcmp(options_mom_.mom.mom_method,'GMM') && options_mom_.mom.analytic_standard_errors
fprintf('\nComputing standard errors using analytical derivatives of moments\n');
D = oo_.mom.model_moments_params_derivs; %already computed in objective function via get_perturbation_params.m
idx_nan = find(any(isnan(D)));
if any(idx_nan)
for i = idx_nan
fprintf('No standard errors available for parameter %s\n',get_the_name(i,options_mom_.TeX, M_, estim_params_, options_mom_))
end
warning('There are NaN in the analytical Jacobian of Moments. Check your bounds and/or priors, or use a different optimizer.')
Asympt_Var = NaN(length(xparam),length(xparam));
SE_values = NaN(length(xparam),1);
return
end
else
fprintf('\nComputing standard errors using numerical derivatives of moments\n');
for i=1:dim_params
%Positive step
xparam_eps_p = xparam;
xparam_eps_p(i,1) = xparam_eps_p(i) + eps_value;
[~, info_p, ~, ~,~, oo__p] = feval(objective_function, xparam_eps_p, Bounds, oo_, estim_params_, M_, options_mom_);
% Negative step
xparam_eps_m = xparam;
xparam_eps_m(i,1) = xparam_eps_m(i) - eps_value;
[~, info_m, ~, ~,~, oo__m] = feval(objective_function, xparam_eps_m, Bounds, oo_, estim_params_, M_, options_mom_);
% The Jacobian:
if nnz(info_p)==0 && nnz(info_m)==0
D(:,i) = (oo__p.mom.model_moments - oo__m.mom.model_moments)/(2*eps_value);
else
problpar = get_the_name(i,options_mom_.TeX, M_, estim_params_, options_mom_);
message_p = get_error_message(info_p, options_mom_);
message_m = get_error_message(info_m, options_mom_);
warning('method_of_moments:info','Cannot compute the Jacobian for parameter %s - no standard errors available\n %s %s\nCheck your bounds and/or priors, or use a different optimizer.\n',problpar, message_p, message_m)
Asympt_Var = NaN(length(xparam),length(xparam));
SE_values = NaN(length(xparam),1);
return
end
end
end
T = options_mom_.nobs; %Number of observations

12
matlab/mh_autocorrelation_function.m
View File

@ -95,13 +95,13 @@ hh=dyn_figure(options_.nodisplay,'Name',FigureName);
bar(0:options_.mh_autocorrelation_function_size,autocor,'k');
axis tight
% create subdirectory <fname>/graphs if it doesn't exist
if ~exist(M_.fname, 'dir')
mkdir('.',M_.fname);
% create subdirectory <dname>/graphs if it doesn't exist
if ~exist(M_.dname, 'dir')
mkdir('.',M_.dname);
end
if ~exist([M_.fname filesep 'graphs'])
mkdir(M_.fname,'graphs');
if ~exist([M_.dname filesep 'graphs'])
mkdir(M_.dname,'graphs');
end
plot_name=get_the_name(column,0,M_,estim_params_,options_);
dyn_saveas(hh,[M_.fname, filesep, 'graphs', filesep, 'MH_Autocorrelation_' plot_name],options_.nodisplay,options_.graph_format)
dyn_saveas(hh,[M_.dname, filesep, 'graphs', filesep, 'MH_Autocorrelation_' plot_name],options_.nodisplay,options_.graph_format)

18
matlab/missing/mex/mjdgges/mjdgges.m
View File

@ -1,4 +1,4 @@
function [ss, tt, zz, sdim, eigval, info] = mjdgges(e, d, qz_criterium, zhreshold)
function [ss, tt, zz, sdim, eigval, info] = mjdgges(e, d, qz_criterium, zhreshold) % --*-- Unitary tests --*--
%
% INPUTS
% e [double] real square (n*n) matrix.
@ -65,3 +65,19 @@ try
catch
info = 1; % Not as precise as lapack's info!
end
%@test:1
%$ try
%$ E =[0,0,0,0,0,0,0,-1,0,0;0,0,0,0,0,0,0,0,-1,0;0,0,0,0,0,0,0,0,0,-0.990099009900990;0,0,0,0,0,0,0,0,0,0.0990099009900990;0,0,0,-1.01010101010101,0,0.0427672955974843,0,0,0,0;0,0,0,0,0,0.128301886792453,-1,0,0,0;0.800000000000000,0,0,0,0,0,0,0,0,0;0,1,0,0,0,0,1,0,0,0;0,0,0.900000000000000,0,0,0,0,0,0,0;0,0,0,-1.01010101010101,-1,0,2,0,-1,0];
%$ D=[0,0,0,0,-1,0,0,-0.792000000000000,0,0;0,0,0,0,0,0,0,0,-0.990000000000000,0;0,0,0,-0.000493818030899887,0,0,0,0,0,-0.882178217821782;0,0,0,-1.00493818030900,0,0,0,0,0,0.0882178217821782;0,0,0,-1,0.128301886792453,0,0,0,0,0;-1,0,0,0,0,0,-0.990000000000000,0,0,0;1,0,0,0,0,0,0,0,0,0;0,1,0,0,0,0,0,0,0,0;0,0,1,0,0,0,0,0,0,0;0,0,0,0,-1,0,0,0,0,0];
%$ [ss, tt, w, sdim, dr.eigval, info1]=mjdgges(E, D, 1.000001, 1e-06);
%$ if sdim==5
%$ t(1) = 1;
%$ else
%$ t(1) = 0;
%$ end
%$ catch
%$ t(1) = 0;
%$ end
%$ T = all(t);
%@eof:1

29
matlab/missing/mex/num_procs/num_procs.m Normal file
View File

@ -0,0 +1,29 @@
function num_logical_cores=num_procs
% function num_logical_cores=num_procs
% Dummy function providing fake number of cores as m-file do not use
% this functionality
%
% OUTPUTS
% num_logical_cores [double] number of logical cores
%
% SPECIAL REQUIREMENTS
% none
% Copyright (C) 2020 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
num_logical_cores=1;

10
matlab/mode_check.m
View File

@ -77,11 +77,11 @@ end
[nbplt,nr,nc,lr,lc,nstar] = pltorg(length(x));
if ~exist([Model.fname filesep 'graphs'],'dir')
mkdir(Model.fname,'graphs');
if ~exist([Model.dname filesep 'graphs'],'dir')
mkdir(Model.dname,'graphs');
end
if TeX && any(strcmp('eps',cellstr(DynareOptions.graph_format)))
fidTeX = fopen([Model.fname, '/graphs/', Model.fname '_CheckPlots.tex'],'w');
fidTeX = fopen([Model.dname, '/graphs/', Model.fname '_CheckPlots.tex'],'w');
fprintf(fidTeX,'%% TeX eps-loader file generated by mode_check.m (Dynare).\n');
fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
fprintf(fidTeX,' \n');
@ -201,12 +201,12 @@ for plt = 1:nbplt
text(0.25,0.5,'log-post')
text(0.69,0.5,'log-lik kernel')
end
dyn_saveas(hh,[Model.fname, '/graphs/', Model.fname '_CheckPlots' int2str(plt) ],DynareOptions.nodisplay,DynareOptions.graph_format);
dyn_saveas(hh,[Model.dname, '/graphs/', Model.fname '_CheckPlots' int2str(plt) ],DynareOptions.nodisplay,DynareOptions.graph_format);
if TeX && any(strcmp('eps',cellstr(DynareOptions.graph_format)))
% TeX eps loader file
fprintf(fidTeX,'\\begin{figure}[H]\n');
fprintf(fidTeX,'\\centering \n');
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_CheckPlots%s}\n',DynareOptions.figures.textwidth*min(k/nc,1),[Model.fname, '/graphs/',Model.fname],int2str(plt));
fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_CheckPlots%s}\n',DynareOptions.figures.textwidth*min(k/nc,1),[Model.dname, '/graphs/',Model.fname],int2str(plt));
fprintf(fidTeX,'\\caption{Check plots.}');
fprintf(fidTeX,'\\label{Fig:CheckPlots:%s}\n',int2str(plt));
fprintf(fidTeX,'\\end{figure}\n');

2
matlab/modules/reporting

@ -1 +1 @@
Subproject commit 13c98b256d8d42f62692aa9d48d379e897f2ca5e
Subproject commit 57a4d6679636e0caa4f98e5a82fb4e33781643b4

4
matlab/ms-sbvar/ms_write_markov_file.m
View File

@ -52,6 +52,10 @@ for i_chain = 1:n_chains
%//=====================================================//
fprintf(fh,'//== Number of states for state_variable[%d] ==//\n', ...
i_chain);
if ~isfield(options.ms.ms_chain(i_chain),'regime') ...
|| isempty(options.ms.ms_chain(i_chain).regime)
error('Chain %u does not feature any regime states. Did you forget to specify them?',i_chain)
end
n_states = length(options.ms.ms_chain(i_chain).regime);
fprintf(fh,'%d\n\n',n_states);

6
matlab/non_linear_dsge_likelihood.m
View File

@ -135,14 +135,16 @@ end
switch DynareOptions.particle.initialization
case 1% Initial state vector covariance is the ergodic variance associated to the first order Taylor-approximation of the model.
StateVectorMean = ReducedForm.constant(mf0);
StateVectorVariance = lyapunov_symm(dr.ghx(mf0,:), dr.ghu(mf0,:)*Q*dr.ghu(mf0,:)', DynareOptions.lyapunov_fixed_point_tol, ...
[A,B] = kalman_transition_matrix(dr,dr.restrict_var_list,dr.restrict_columns,Model.exo_nbr);
StateVectorVariance = lyapunov_symm(A, B*Q*B', DynareOptions.lyapunov_fixed_point_tol, ...
DynareOptions.qz_criterium, DynareOptions.lyapunov_complex_threshold, [], DynareOptions.debug);
StateVectorVariance = StateVectorVariance(mf0,mf0);
case 2% Initial state vector covariance is a monte-carlo based estimate of the ergodic variance (consistent with a k-order Taylor-approximation of the model).
StateVectorMean = ReducedForm.constant(mf0);
old_DynareOptionsperiods = DynareOptions.periods;
DynareOptions.periods = 5000;
y_ = simult(DynareResults.steady_state, dr,Model,DynareOptions,DynareResults);
y_ = y_(state_variables_idx,2001:5000);
y_ = y_(dr.order_var(state_variables_idx),2001:5000); %state_variables_idx is in dr-order while simult_ is in declaration order
StateVectorVariance = cov(y_');
DynareOptions.periods = old_DynareOptionsperiods;
clear('old_DynareOptionsperiods','y_');

21
matlab/occbin/call_solve_one_constraint.m
View File

@ -1,21 +0,0 @@
% Solve model, generate model IRFs
[zdatalinear, zdatapiecewise, zdatass, oobase_, Mbase_ ] = ...
solve_one_constraint(modnam,modnamstar,...
constraint, constraint_relax,...
shockssequence,irfshock,nperiods,maxiter);
% unpack the IRFs
for i=1:Mbase_.endo_nbr
eval([Mbase_.endo_names{i}, '_uncdifference=zdatalinear(:,i);']);
eval([Mbase_.endo_names{i}, '_difference=zdatapiecewise(:,i);']);
eval([Mbase_.endo_names{i}, '_ss=zdatass(i);']);
end
nparams = length(Mbase_.param_names);
for i = 1:nparams
eval([Mbase_.param_names{i}, '= Mbase_.params(i);']);
end

21
matlab/occbin/call_solve_two_constraints.m
View File

@ -1,21 +0,0 @@
[zdatalinear, zdatapiecewise, zdatass, oobase_, Mbase_] = solve_two_constraints(...
modnam_00,modnam_10,modnam_01,modnam_11,...
constraint1, constraint2,...
constraint_relax1, constraint_relax2,...
scalefactormod,irfshock,nperiods,curb_retrench,maxiter);
for i=1:Mbase_.endo_nbr
eval([Mbase_.endo_names{i}, '_uncdifference=zdatalinear(:,i);']);
eval([Mbase_.endo_names{i}, '_difference=zdatapiecewise(:,i);']);
eval([Mbase_.endo_names{i}, '_ss=zdatass(i);']);
end
constraint1_difference = process_constraint(constraint1, '_difference', Mbase_.endo_names,0);
constraint2_difference = process_constraint(constraint2, '_difference', Mbase_.endo_names,0);
nparams = length(Mbase_.param_names);
for i = 1:nparams
eval([Mbase_.param_names{i}, '= Mbase_.params(i);']);
end

19
matlab/occbin/evaluate_model.m
View File

@ -1,19 +0,0 @@
function [r,g1,g2,g3] = evaluate_model(z,x,M,ss)
ll = M.lead_lag_incidence';
y = z(find(ll(:)));
switch nargout
case 1
r = feval([M.fname '.dynamic'],y,x, ...
M.params, ss, 1);
case 2
[r,g1] = feval([M.fname '.dynamic'],y,x, ...
M.params, ss, 1);
case 3
[r,g1,g2] = feval([M.fname '.dynamic'],y,x, ...
M.params, ss, 1);
case 4
[r,g1,g2,g3] = feval([M.fname '.dynamic'],y,x, ...
M.params, ss, 1);
end

26
matlab/occbin/get_coef.m
View File

@ -1,26 +0,0 @@
function [coef_y,coef_u] = get_coef(jacobian,M)
ll = M.lead_lag_incidence;
endo_nbr = M.endo_nbr;
coef_y = zeros(endo_nbr,3*endo_nbr);
coef_u = zeros(endo_nbr,M.exo_nbr);
if M.maximum_lag > 0
[~,c1,c2] = find(ll(1,:));
coef_y(:,c1) = jacobian(:,c2);
[~,c1,c2] = find(ll(2,:));
coef_y(:,c1+endo_nbr) = jacobian(:,c2);
if M.maximum_lead > 0
[~,c1,c2] = find(ll(3,:));
coef_y(:,c1+2*endo_nbr) = jacobian(:,c2);
end
else
[~,c1,c2] = find(ll(1,:));
coef_y(:,c1+endo_nbr) = jacobian(:,c2);
if M.maximum_lead > 0
[~,c1,c2] = find(ll(2,:));
coef_y(:,c1+2*endo_nbr) = jacobian(:,c2);
end
end
coef_u = jacobian(:,max(ll(end,:))+1:end);

71
matlab/occbin/get_complementarity_conditions.m
View File

@ -1,71 +0,0 @@
function [lb,ub,eq_index] = get_complementarity_conditions(M,ramsey_policy)
% Copyright (C) 2014-2018 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
ub = inf(M.endo_nbr,1);
lb = -ub;
eq_index = (1:M.endo_nbr)';
if ramsey_policy
if isfield(M,'ramsey_model_constraints')
rc = M.ramsey_model_constraints;
for i = 1:length(rc)
switch rc{i}{2}
case {'>','>='}
lb(rc{i}{1}) = eval(rc{i}{3});
case {'<','<='}
ub(rc{i}{1}) = eval(rc{i}{3});
otherwise
error('Wrong operator in get_complementarity_conditions')
end
eq_index(i) = 1;
end
end
end
etags = M.equations_tags;
for i=1:size(etags,1)
if strcmp(etags{i,2},'mcp')
str = etags{i,3};
kop = strfind(etags{i,3},'<');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)), M.endo_names));
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized',etags{i,3},b{1}]))
end
ub(k) = str2num(str(kop+1:end));
eq_index(etags{i,1}) = k;
eq_index(k) = etags{i,1};
else
kop = strfind(etags{i,3},'>');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)), M.endo_names));
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized',etags{i},b{1}]))
end
lb(k) = str2num(str(kop+1:end));
eq_index(etags{i,1}) = k;
eq_index(k) = etags{i,1};
else
error(sprintf(['Complementarity condition %s can''t be ' ...
'parsed'],etags{i,3}))
end
end
end
end

82
matlab/occbin/get_deriv.m
View File

@ -1,82 +0,0 @@
function [hm1,h,hl1,j,resid] = get_deriv(M_,ys_)
iy_ = M_.lead_lag_incidence;
it_ = 1;
x = zeros(1,M_.exo_nbr);
% For most models, there are leads, lags and current values of variables
if size(iy_,1)==3
% find non-zero columns of hm1
lag_cols = find(iy_(1,:)~=0);
% find non-zero columns of h
con_cols = find(iy_(2,:));
% find non-zero columns of hl1
lea_cols = find(iy_(3,:));
% If models either lacks leads or lags, iy_ will have two rows
% In this case, we guess that the row with more nonzeros is the row with current variables
elseif size(iy_,1)==2
% if first row has more nonzero entries than the second, assume model lacks lagged variables
if length(find(iy_(1,:)))>length(find(iy_(2,:)))
warning('Model does not have lagged endogenous variables')
con_cols = find(iy_(1,:));
lea_cols = find(iy_(2,:));
lag_cols = [];
else
warning('Model does not have expected future endogenous variables')
lag_cols = find(iy_(1,:));
con_cols = find(iy_(2,:));
lea_cols = [];
end
end
% find number of entries for y vector
ny = length(find(iy_~=0));
% build steady state y
y = ys_(lag_cols);
y = [y;ys_(con_cols)];
y = [y;ys_(lea_cols)];
if ismac
eval(['[resid,g1]=',M_.fname,'.dynamic(y,x, M_.params, ys_, it_);']);
% Older versions of DYNARE for Mac did not include ys_ in the call structure
%eval(['[resid,g1]=',M_.fname,'.dynamic(y,x, M_.params, it_);']);
else
eval(['[resid,g1]=',M_.fname,'.dynamic(y,x, M_.params, ys_, it_);']);
end
hm1=zeros(M_.endo_nbr);
h = hm1;
hl1 = hm1;
j = zeros(M_.endo_nbr,M_.exo_nbr);
% build hm1
nlag_cols = length(lag_cols);
for i=1:nlag_cols
hm1(:,lag_cols(i)) = g1(:,i);
end
% build h
ncon_cols = length(con_cols);
for i=1:ncon_cols
h(:,con_cols(i)) = g1(:,i+nlag_cols);
end
% build hl1
nlea_cols = length(lea_cols);
for i=1:nlea_cols
hl1(:,lea_cols(i)) = g1(:,i+nlag_cols+ncon_cols);
end
for i = 1:M_.exo_nbr
j(:,i) =g1(:,i+ny);
end

74
matlab/occbin/get_occbin_complementarity_conditions.m
View File

@ -1,74 +0,0 @@
function [ivar,ieq,lb,ub] = get_occbin_complementarity_conditions(M,ramsey_policy)
% Copyright (C) 2015-2018 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
nrow = 1;
if ramsey_policy
if isfield(M,'ramsey_model_constraints')
rc = M.ramsey_model_constraints;
for i = 1:length(rc)
switch rc{i}{2}
case {'>','>='}
ivar(nrow) = rc{i}{1};
ieq(nrow) = rc{i}{1};
lb(nrow) = eval(rc{i}{3});
case {'<','<='}
ivar(nrow) = rc{i}{1};
ieq(nrow) = rc{i}{1};
ub(nrow) = eval(rc{i}{3});
otherwise
error('Wrong operator in get_complementarity_conditions')
end
nrow = nrow + 1;
end
end
end
etags = M.equations_tags;
for i=1:size(etags,1)
if strcmp(etags{i,2},'mcp')
str = etags{i,3};
kop = strfind(etags{i,3},'<');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)), M.endo_names));
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized',etags{i,3},b{1}]))
end
ivar(nrow) = k;
ieq(nrow) = etags{i,1};
ub(nrow) = eval(str(kop+1:end));
else
kop = strfind(etags{i,3},'>');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)), M.endo_names));
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized',etags{i},b{1}]))
end
ivar(nrow) = k;
ieq(nrow) = etags{i,1};
lb(k) = eval(str(kop+1:end));
else
error(sprintf(['Complementarity condition %s can''t be ' ...
'parsed'],etags{i,3}))
end
end
nrow = nrow + 1;
end
end

88
matlab/occbin/get_occbin_constraints.m
View File

@ -1,88 +0,0 @@
function [i_base,i_alt,c_base,c_alt] = get_occbin_constraints(M,steady_state,ramsey_policy)
% Copyright (C) 2015 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
nrow = 1;
if ramsey_policy
if isfield(M,'ramsey_model_constraints')
rc = M.ramsey_model_constraints;
for i = 1:length(rc)
switch rc{i}{2}
case {'>','>='}
ivar(nrow) = rc{i}{1};
ieq(nrow) = rc{i}{1};
lb(nrow) = eval(rc{i}{3});
case {'<','<='}
ivar(nrow) = rc{i}{1};
ieq(nrow) = rc{i}{1};
ub(nrow) = eval(rc{i}{3});
otherwise
error('Wrong operator in get_complementarity_conditions')
end
nrow = nrow + 1;
end
end
end
i_base = {};
i_alt = {};
etags = M.equations_tags;
m = 1;
base = true;
for i=1:size(etags,1)
[iv,boundary,operator] = parse_constraint(etags{i,3},M.endo_names,M.params,M.param_names);
if strcmp(etags{i,2},'OCCBIN')
if base
i_alt{m} = 1:M.eq_nbr;
i_alt{m}(etags{i,1}) = [];
c_base{m,1} = etags{i,1};
c_base{m,2} = iv;
c_base{m,3} = boundary - steady_state(iv);
c_base{m,4} = operator;
base = false;
else
i_base{m} = 1:M.eq_nbr;
i_base{m}(etags{i,1}) = [];
c_alt{m,1} = etags{i,1};
c_alt{m,2} = iv;
c_alt{m,3} = boundary - steady_state(iv);
c_alt{m,4} = operator;
base = true;
m = m + 1;
end
end
end
if ~base
error('OCCBIN: constraints must come by pair')
end
function [iv,boundary,operator] = parse_constraint(str,endo_names,params,param_names)
delim = {'<=','>=','<','>'};
[c,operator] = strsplit(str,delim);
operator = operator{1};
iv = strmatch(strtrim(c{1}),endo_names);
% try for a number
boundary = str2num(strtrim(c{2}));
% if not a number try for a parameter name
if isempty(boundary)
k = strmatch(strtrim(c{2}),param_names);
if isempty(k)
error(['OCCBIN: illegal constraint ' str]);
end
boundary = params(k);
end

28
matlab/occbin/get_pq.m
View File

@ -1,28 +0,0 @@
function [p,q]=get_pq(dr_,nstatic,nfwrd)
nvars = size(dr_.ghx,1);
nshocks = size(dr_.ghu,2);
statevar_pos = (nstatic +1):(nvars-nfwrd);
p = zeros(nvars);
% interlace matrix
nnotzero = length(statevar_pos);
for i=1:nnotzero
p(:,statevar_pos(i)) = dr_.ghx(:,i);
end
% reorder p matrix according to order in lgy_
inverse_order = zeros(nvars,1);
for i=1:nvars
inverse_order(i) = find(i==dr_.order_var);
end
p_reordered = zeros(nvars);
q = zeros(nvars,nshocks);
for i=1:nvars
for j=1:nvars
p_reordered(i,j)=p(inverse_order(i),inverse_order(j));
end
q(i,:)=dr_.ghu(inverse_order(i),:);
end
p=p_reordered;

10
matlab/occbin/get_residuals.m
View File

@ -1,10 +0,0 @@
function r = get_residuals(ivar,lb,ub,M,oo)
ss = oo.steady_state;
for i = 1:length(ivar)
% only one is different from zero
ss(ivar(i)) = lb(i) + ub(i);
end
oo.steady_state = ss;
r = evaluate_model(M,oo);

81
matlab/occbin/makechart.m
View File

@ -1,81 +0,0 @@
function makechart(titlelist,legendlist,figlabel,ylabels,zdata1,zdata2,zdata3)
figure
titlelist = char(strrep(cellstr(titlelist),'_','.'));
ndsets=3; % default, changed below as applicable
if nargin==5
zdata2=nan*zdata1;
zdata3=nan*zdata1;
ndsets =1;
elseif nargin == 6
zdata3 =nan*zdata1;
ndsets=2;
elseif ((nargin>8) | (nargin <=4))
error ('makechart takes 5 to 6 arguments')
end
nobs = size(zdata1,1);
xvalues = (1:nobs)';
nvars = size(titlelist,1);
if nvars==1
nrows=1;
ncols = 1;
elseif nvars==2
nrows =2;
ncols = 1;
elseif (nvars == 3 | nvars ==4)
nrows = 2;
ncols =2;
elseif (nvars==5 |nvars ==6)
nrows = 3;
ncols = 2;
elseif (nvars==7 | nvars==8)
nrows = 4;
ncols = 2;
elseif (nvars==9 | nvars==10)
nrows = 5;
ncols = 2;
else
error('too many variables (makechart)')
end
for i = 1:nvars
subplot(nrows,ncols,i)
h1=plot(xvalues,zdata1(:,i),'b-','linewidth',2); hold on
h1=plot(xvalues,zdata2(:,i),'r--','linewidth',2); hold on
h2=plot(xvalues,zdata3(:,i),'b-','LineWidth',3);
[x0, x1, y10, y11] = pickaxes(xvalues,zdata1(:,i));
[x0, x1, y20, y21] = pickaxes(xvalues,zdata2(:,i));
[x0, x1, y30, y31] = pickaxes(xvalues,zdata3(:,i));
y0 = min([y10,y20,y30]);
y1 = max([y11,y21,y31]);
if y0==y1
y1=y0+1;
end
axis([x0 x1 y0 y1])
set(h1);
if i==1 && isempty(legendlist)==0
legend(legendlist)
text('String',figlabel,'Units','normalized','Position',[1.2 1.24],...
'FontSize',14,'FontWeight','bold','HorizontalAlignment','center');
end
if i==nvars | i==nvars-1
xlabel('Time');
end
% set(gca,'XTick',xtick)
% set(gca,'XTickLabel',xticklabel)
title([num2str(i),'. ',titlelist(i,:)]);
ylabel(ylabels(i,:))
end
% sets printing preferences
%printpref

136
matlab/occbin/makechart9.m
View File

@ -1,136 +0,0 @@
function makechart9(titlelist,legendlist,figlabel,yearshock,ylabels,...
zdata1,zdata2,zdata3,zdata4,zdata5,zdata6,zdata7)
figure
titlelist = char(strrep(cellstr(titlelist),'_','.'));
ndsets=7; % default, changed below as applicable
if nargin==6
zdata2=nan*zdata1;
zdata3=nan*zdata1;
zdata4=nan*zdata1;
zdata5=nan*zdata1;
zdata6=nan*zdata1;
zdata7=nan*zdata1;
ndsets =1;
elseif nargin==7
zdata3=nan*zdata1;
zdata4=nan*zdata1;
zdata5=nan*zdata1;
zdata6=nan*zdata1;
zdata7=nan*zdata1;
ndsets =2;
elseif nargin == 8
zdata4 =nan*zdata1;
zdata5 =nan*zdata1;
zdata6=nan*zdata1;
zdata7=nan*zdata1;
ndsets=3;
elseif nargin == 9
zdata5 =nan*zdata1;
zdata6=nan*zdata1;
zdata7=nan*zdata1;
ndsets=4;
elseif nargin == 10
zdata6 =nan*zdata1;
zdata7=nan*zdata1;
ndsets=5;
elseif nargin == 11
zdata7=nan*zdata1;
ndsets=6;
elseif ((nargin>=13) | (nargin <=3))
error ('makechart takes 4 to 10 arguments')
end
nobs = size(zdata1,1);
if yearshock>-100
xvalues = yearshock+(0:nobs-1)'/4; % Matteo plot year on x axis
else
xvalues = (1:nobs)'; % Matteo plot year on x axis
end
nvars = size(titlelist,1);
if nvars==1
nrows=1;
ncols = 1;
elseif nvars==2
nrows =2;
ncols = 1;
elseif nvars == 3
nrows = 3;
ncols = 1;
elseif nvars==4
nrows = 2;
ncols = 2;
elseif (nvars==5 | nvars ==6)
nrows = 3;
ncols = 2;
elseif (nvars==7 | nvars==8)
nrows = 4;
ncols = 2;
elseif nvars>8 & nvars<=12
nrows = 3;
ncols = 4;
elseif nvars>12 & nvars<=15
nrows = 5;
ncols = 3;
else
error('too many variables (makechart)')
end
for i = 1:nvars
subplot(nrows,ncols,i)
h1=plot(xvalues,zdata1(:,i),'k',...
xvalues,zdata2(:,i),'r',...
xvalues,zdata3(:,i),'b',...
xvalues,zdata4(:,i),'g',...
xvalues,zdata5(:,i),'g',...
xvalues,zdata6(:,i),'c',...
xvalues,zdata7(:,i),'y');
[x0, x1, y10, y11] = pickaxes(xvalues,zdata1(:,i));
[x0, x1, y20, y21] = pickaxes(xvalues,zdata2(:,i));
[x0, x1, y30, y31] = pickaxes(xvalues,zdata3(:,i));
[x0, x1, y40, y41] = pickaxes(xvalues,zdata4(:,i));
[x0, x1, y50, y51] = pickaxes(xvalues,zdata5(:,i));
[x0, x1, y60, y61] = pickaxes(xvalues,zdata6(:,i));
[x0, x1, y70, y71] = pickaxes(xvalues,zdata7(:,i));
grid on
y0 = min([y10,y20,y30,y40,y50,y60,y70]);
y1 = max([y11,y21,y31,y41,y51,y61,y71]);
if y0==y1
y1=y0+1;
end
axis([x0 x1 y0 y1])
set(h1,'linewidth',2);
if i==1
if numel(strvcat(legendlist(1,:)))
h=legend(legendlist,'Location','Northwest');
set(h,'Fontsize',8)
end
end
if i==1
if nvars>3
text('String',figlabel,'Units','normalized','Position',[1.2 1.21],...
'FontSize',13,'FontWeight','bold','HorizontalAlignment','center');
else
text('String',figlabel,'Units','normalized','Position',[0.4 1.24],...
'FontSize',13,'FontWeight','bold','HorizontalAlignment','center');
end
end
%set(gca,'XTick',xtick)
%set(gca,'XTickLabel',xticklabel)
title(titlelist(i,:),'Fontsize',11);
ylabel(ylabels(i,:))
end
% sets printing preferences
%printpref

24
matlab/occbin/map_regime.m
View File

@ -1,24 +0,0 @@
function [regime, regimestart]=map_regime(violvecbool)
nperiods = length(violvecbool)-1;
% analyse violvec and isolate contiguous periods in the other regime.
regime(1) = violvecbool(1);
regimeindx = 1;
regimestart(1) = 1;
for i=2:nperiods
if violvecbool(i)~=regime(regimeindx)
regimeindx=regimeindx+1;
regime(regimeindx) = violvecbool(i);
regimestart(regimeindx)=i;
end
end
if (regime(1) == 1 & length(regimestart)==1)
warning('Increase nperiods');
end
if (regime(end)==1)
warning('Increase nperiods');
end

63
matlab/occbin/mkdata.m
View File

@ -1,63 +0,0 @@
function [zdata]=mkdata(nperiods,decrulea,decruleb,endog_,exog_,wishlist,irfshock,scalefactormod,init)
%[nsim, ksim, ysim, isim, csim] = mkdata(nperiods,cofb,endog_)
% given decision rule
neqs = length(endog_);
if nargin<9
init = zeros(neqs,1);
end
if nargin<8
scalefactormod=1;
end
if nargin<7
error('Not enough inputs')
end
history = zeros(neqs,nperiods+1);
nshocks = size(irfshock,1);
for i = 1:nshocks
shockpos = strmatch(irfshock(i,:), exog_,'exact');
if ~isempty(shockpos)
irfshockpos(i) = shockpos;
else
error(['Shock ',irfshock(i,:),' is not in the model']);
end
end
% generate data
% history will contain data, the state vector at each period in time will
% be stored columnwise.
history = zeros(neqs,nperiods);
history(:,1)= init;
lengthshock = size(scalefactormod,1);
errvec = zeros(length(exog_), 1);
for i = 2:nperiods+1
if i<=(lengthshock+1)
for j = 1:nshocks
errvec(irfshockpos(j)) = scalefactormod(i-1,j);
end
history(:,i) = decrulea * history(:,i-1)+decruleb*errvec;
else
% update endogenous variables
history(:,i) = decrulea * history(:,i-1);
end
end
% extract desired variables
nwish=size(wishlist,1);
wishpos = zeros(nwish,1);
history=history';
for i=1:nwish
wishpos(i) = strmatch(wishlist(i,:), endog_, 'exact');
end
zdata = history(2:end,wishpos);

124
matlab/occbin/mkdatap_anticipated.m
View File

@ -1,124 +0,0 @@
function [zdata]=mkdatap_anticipated(nperiods,decrulea,decruleb,...
cof,Jbarmat,cofstar,Jstarbarmat,Dstarbarmat,...
regime,regimestart,violvecbool,...
endog_,exog_,irfshock,scalefactormod,init)
nvars = length(endog_);
if nargin<16
init=zeros(nvars,1);
end
if nargin<15
scalefactormod=1;
end
nshocks = size(irfshock,1);
for i = 1:nshocks
shockpos = strmatch(irfshock(i,:), exog_, 'exact');
if ~isempty(shockpos)
irfshockpos(i) = shockpos;
else
error(['Shock ', irfshock(i,:), ' is not in the model']);
end
end
nregimes = length(regime);
Cbarmat = cof(:,1:nvars);
Bbarmat = cof(:,nvars+1:2*nvars);
Abarmat = cof(:,2*nvars+1:3*nvars);
% cofstar contains the system for the model when the constraint binds
Cstarbarmat = cofstar(:,1:nvars);
Bstarbarmat = cofstar(:,nvars+1:2*nvars);
Astarbarmat = cofstar(:,2*nvars+1:3*nvars);
% get the time-dependent decision rules
Tmax = regimestart(nregimes)-1; % Tmax is the position of the last period
% when the constraint binds
if Tmax > 0
P = zeros(nvars,nvars,Tmax);
D = zeros(nvars,Tmax);
invmat = inv((Astarbarmat*decrulea+Bstarbarmat));
P(:,:,Tmax) = -invmat*Cstarbarmat;
D(:,Tmax) = -invmat*Dstarbarmat;
% equivalent to pre-multiplying by the inverse above if the target
% matrix is invertible. Otherwise it yields the minimum state solution
%P(:,:,Tmax) = -(Astarbarmat*decrulea+Bstarbarmat)\Cstarbarmat;
%D(:,Tmax) = -(Astarbarmat*decrulea+Bstarbarmat)\Dstarbarmat;
for i = Tmax-1:-1:1
if violvecbool(i)
invmat = inv(Bstarbarmat+Astarbarmat*P(:,:,i+1));
P(:,:,i)=-invmat*Cstarbarmat;
D(:,i) = -invmat*(Astarbarmat*D(:,i+1)+Dstarbarmat);
else
invmat = inv(Bbarmat+Abarmat*P(:,:,i+1));
P(:,:,i)=-invmat*Cbarmat;
D(:,i) = -invmat*(Abarmat*D(:,i+1));
end
end
if Tmax > 1
if violvecbool(1)
E = -invmat*Jstarbarmat;
else
E = -invmat*Jbarmat;
end
else
invmat = inv(Astarbarmat*decrulea+Bstarbarmat);
E = -invmat*Jstarbarmat;
end
end
% generate data
% history will contain data, the state vector at each period in time will
% be stored columnwise.
history = zeros(nvars,nperiods+1);
history(:,1) = init;
errvec = zeros(length(exog_), 1);
% deal with predetermined conditions
for i = 1:nshocks
errvec(irfshockpos(i)) = scalefactormod(i);
end
% deal with shocks
irfpos =1;
if irfpos <=Tmax
history(:,irfpos+1) = P(:,:,irfpos)* history(:,irfpos)+...
D(:,irfpos) + E*errvec;
else
history(:,irfpos+1) = decrulea*history(:,irfpos)+decruleb*errvec;
end
% all other periods
for irfpos=2:nperiods+1
if irfpos <=Tmax
history(:,irfpos+1) = P(:,:,irfpos)* history(:,irfpos)+...
D(:,irfpos);
else
history(:,irfpos+1) = decrulea*history(:,irfpos);
end
end
history=history';
zdata = history(2:end,:);

180
matlab/occbin/mkdatap_anticipated_2constraints.m
View File

@ -1,180 +0,0 @@
function [zdata]=mkdatap_anticipated_2constraints(nperiods,decrulea,decruleb,...
cof,Jbarmat,...
cof10,Jbarmat10,Dbarmat10,...
cof01,Jbarmat01,Dbarmat01,...
cof11,Jbarmat11,Dbarmat11,...
regime1,regimestart1,...
regime2,regimestart2,...
violvecbool,endog_,exog_,...
irfshock,scalefactormod,init)
nvars = size(endog_,1);
if nargin<16
init=zeros(nvars,1);
end
if nargin<15
scalefactormod=1;
end
nshocks = size(irfshock,1);
for i = 1:nshocks
shockpos = strmatch(irfshock(i,:),exog_,'exact');
if ~isempty(shockpos)
irfshockpos(i) = shockpos;
else
error(['Shock ',irfshock(i,:),' is not in the model']);
end
end
Cbarmat = cof(:,1:nvars);
Bbarmat = cof(:,nvars+1:2*nvars);
Abarmat = cof(:,2*nvars+1:3*nvars);
% cofstar contains the system for the model when the constraint binds
Cbarmat10 = cof10(:,1:nvars);
Bbarmat10 = cof10(:,nvars+1:2*nvars);
Abarmat10 = cof10(:,2*nvars+1:3*nvars);
Cbarmat01 = cof01(:,1:nvars);
Bbarmat01 = cof01(:,nvars+1:2*nvars);
Abarmat01 = cof01(:,2*nvars+1:3*nvars);
Cbarmat11 = cof11(:,1:nvars);
Bbarmat11 = cof11(:,nvars+1:2*nvars);
Abarmat11 = cof11(:,2*nvars+1:3*nvars);
% get the time-dependent decision rules
nregimes1 = length(regime1);
nregimes2 = length(regime2);
Tmax = max([regimestart1(nregimes1) regimestart2(nregimes2)])-1; % Tmax is the position of the last period
% when the constraint binds
if Tmax > 0
P = zeros(nvars,nvars,Tmax);
D = zeros(nvars,Tmax);
% invmat = inv((Astarbarmat*decrulea+Bstarbarmat));
% P(:,:,Tmax) = -invmat*Cstarbarmat;
% D(:,Tmax) = -invmat*Dstarbarmat;
if (violvecbool(Tmax,1) & ~violvecbool(Tmax,2))
%XXX fix next three lines
invmat = inv((Abarmat10*decrulea+Bbarmat10));
P(:,:,Tmax) = -invmat*Cbarmat10;
D(:,Tmax) = -invmat*Dbarmat10;
elseif (violvecbool(Tmax,1) & violvecbool(Tmax,2))
invmat = inv((Abarmat11*decrulea+Bbarmat11));
P(:,:,Tmax) = -invmat*Cbarmat11;
D(:,Tmax) = -invmat*Dbarmat11;
else
invmat = inv((Abarmat01*decrulea+Bbarmat01));
P(:,:,Tmax) = -invmat*Cbarmat01;
D(:,Tmax) = -invmat*Dbarmat01;
end
for i = Tmax-1:-1:1
if (violvecbool(i,1) & ~violvecbool(i,2))
invmat = inv(Bbarmat10+Abarmat10*P(:,:,i+1));
P(:,:,i)=-invmat*Cbarmat10;
D(:,i) = -invmat*(Abarmat10*D(:,i+1)+Dbarmat10);
elseif (~violvecbool(i,1) & violvecbool(i,2))
invmat = inv(Bbarmat01+Abarmat01*P(:,:,i+1));
P(:,:,i)=-invmat*Cbarmat01;
D(:,i) = -invmat*(Abarmat01*D(:,i+1)+Dbarmat01);
elseif (violvecbool(i,1) & violvecbool(i,2))
invmat = inv(Bbarmat11+Abarmat11*P(:,:,i+1));
P(:,:,i)=-invmat*Cbarmat11;
D(:,i) = -invmat*(Abarmat11*D(:,i+1)+Dbarmat11);
else
invmat = inv(Bbarmat+Abarmat*P(:,:,i+1));
P(:,:,i)=-invmat*Cbarmat;
D(:,i) = -invmat*(Abarmat*D(:,i+1));
end
end
% Double check the appropriate invmat in each case
% right now -- inherited from previous loop
if Tmax > 1
if ( ~violvecbool(1,1) & violvecbool(1,2) )
E = -invmat*Jbarmat01;
elseif ( violvecbool(1,1) & ~violvecbool(1,2) )
E = -invmat*Jbarmat10;
elseif ( violvecbool(1,1) & violvecbool(1,2) )
E = -invmat*Jbarmat11;
else
E = -invmat*Jbarmat;
end
else % Tmax is equal to 1
% invmat = inv((Astarbarmat*decrulea+Bstarbarmat));
% E = -invmat*Jstarbarmat;
if ( ~violvecbool(1,1) & violvecbool(1,2) )
invmat = inv((Abarmat01*decrulea+Bbarmat01));
E = -invmat*Jbarmat01;
elseif ( violvecbool(1,1) & violvecbool(1,2) )
invmat = inv((Abarmat11*decrulea+Bbarmat11));
E = -invmat*Jbarmat11;
else
invmat = inv((Abarmat10*decrulea+Bbarmat10));
E = -invmat*Jbarmat10;
end
end
end
% generate data
% history will contain data, the state vector at each period in time will
% be stored columnwise.
history = zeros(nvars,nperiods+1);
history(:,1) = init;
errvec = zeros(size(exog_,1),1);
for i = 1:nshocks
errvec(irfshockpos(i)) = scalefactormod(i);
end
% deal with shocks
irfpos =1;
if irfpos <=Tmax
history(:,irfpos+1) = P(:,:,irfpos)* history(:,irfpos)+...
D(:,irfpos) + E*errvec;
else
history(:,irfpos+1) = decrulea*history(:,irfpos)+decruleb*errvec;
end
% all other periods
for irfpos=2:nperiods+1
if irfpos <=Tmax
history(:,irfpos+1) = P(:,:,irfpos)* history(:,irfpos)+...
D(:,irfpos);
else
history(:,irfpos+1) = decrulea*history(:,irfpos);
end
end
history=history';
zdata = history(2:end,:);

19
matlab/occbin/occbin.m
View File

@ -1,19 +0,0 @@
function [endo_simul,endo_simul_no_constraint,status] = occbin(M,oo,options)
% function oo=occbin(M,oo,options) solves linear models with occasionally
% binding constraints using OCCBIN by L. Guerrieri
status = 1;
constraint_nbr = sum(strcmp(upper(M.equations_tags(:,2)),'OCCBIN'))/2;
switch(constraint_nbr)
case 1
[zdatalinear_ zdatapiecewise_ zdatass_ oobase_ ] = ...
solve_one_constraint(M,oo,options);
case 2
[zdatalinear_ zdatapiecewise_ zdatass_ oobase_ ] = ...
solve_two_constraints(M,oo,options);
otherwise
error('OCCBIN can only handle two constraints in a model')
end
endo_simul = zdatapiecewise_';
endo_simul_no_constraint = zdatalinear_';

16
matlab/occbin/pickaxes.m
View File

@ -1,16 +0,0 @@
function [x0,x1,y0,y1] = pickaxes(xvalues,yvalues)
x0=xvalues(1);
nobs = length(xvalues);
x1=xvalues(nobs);
maxy = max(yvalues);
miny = min(yvalues);
y0 = miny - .05*abs(miny);
if (miny>0 & y0<0)
y0 = 0;
end
y1 = maxy + .05*abs(maxy);

Some files were not shown because too many files have changed in this diff Show More