time-shift
Johannes Pfeifer 2013-06-20 18:45:23 +02:00
commit 7a7ef657fa
35 changed files with 670 additions and 214 deletions

View File

@ -23,6 +23,7 @@ EXTRA_DIST = \
contrib \
NEWS \
license.txt \
README.md \
windows/dynare.nsi \
windows/mexopts-win32.bat \
windows/mexopts-win64.bat \

277
README.md Normal file
View File

@ -0,0 +1,277 @@
# Dynare
Described on the homepage: <http://www.dynare.org/>
Most users should use the precompiled package available for your OS, also
available via the Dynare homepage: <http://www.dynare.org/download/dynare-stable>.
# License
Most of the source files are covered by the GNU General Public Licence version
3 or later (there are some exceptions to this, see [license.txt](license.txt) in
Dynare distribution for specifics).
# Building Dynare From Source
Here, we explain how to build from source:
- Dynare, including preprocessor and MEX files for MATLAB and Octave
- Dynare++
- all the associated documentation (PDF and HTML)
This source can be retrieved in three forms:
- via git, at <https://github.com/DynareTeam/dynare.git>
- using the stable source archive of the latest Dynare version (currently 4.3) from <http://www.dynare.org/download/dynare-4.3/source>
- using a source snapshot of the unstable version, from <http://www.dynare.org/download/dynare-unstable/source-snapshot>
Note that if you obtain the source code via git, you will need to install more tools (see below).
The first section of this page gives general instructions, which apply to all platforms. Then some specific platforms are discussed.
**NB**: Here, when we refer to 32-bit or 64-bit, we refer to the type of MATLAB installation, not the type of Windows installation. It is perfectly possible to run a 32-bit MATLAB on a 64-bit Windows: in that case, instructions for Windows 32-bit should be followed. To determine the type of your MATLAB installation, type:
```matlab
>> computer
```
at the MATLAB prompt: if it returns `PCWIN`, then you have a 32-bit MATLAB; if it returns `PCWIN64`, then you have a 64-bit MATLAB.
**Contents**
1. [**General Instructions**](#general-instructions)
1. [**Debian or Ubuntu**](#debian-or-ubuntu)
1. [**Fedora**](#fedora)
1. [**Windows**](#windows)
1. [**Mac OS X**](#mac-os-x)
## General Instructions
### Prerequisites
A number of tools and libraries are needed in order to recompile everything. You don't necessarily need to install everything, depending on what you want to compile.
- A POSIX compliant shell and an implementation of Make (mandatory)
- The [GNU Compiler Collection](http://gcc.gnu.org/), with gcc, g++ and gfortran (mandatory)
- [MATLAB](http://www.dynare.org/DynareWiki/BuildingDynareFromSource?action=AttachFile&do=view&target=dynare-mingw64-libs.zip) (if you want to compile MEX for MATLAB)
- [GNU Octave](http://www.octave.org), with the development headers (if you want to compile MEX for Octave)
- [Boost libraries](http://www.boost.org), version 1.36 or later
- [Bison](http://www.gnu.org/software/bison/), version 2.3 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)
- [Autoconf](http://www.gnu.org/software/autoconf/), version 2.62 or later (only if you get the source through Git) (see [Installing an updated version of Autoconf in your own directory, in GNU/Linux](http://www.dynare.org/DynareWiki/AutoMake))
- [Automake](http://www.gnu.org/software/automake/), version 1.11.2 or later (only if you get the source through Git) (see [Installing an updated version of AutoMake in your own directory, in GNU/Linux](http://www.dynare.org/DynareWiki/AutoMake))
- [CWEB](http://www-cs-faculty.stanford.edu/%7Eknuth/cweb.html), with its tools `ctangle` and `cweave` (only if you want to build Dynare++ and get the source through Git)
- An implementation of BLAS and LAPACK: either [ATLAS](http://math-atlas.sourceforge.net/), [OpenBLAS](http://xianyi.github.com/OpenBLAS/), Netlib ([BLAS](http://www.netlib.org/blas/), [LAPACK](http://www.netlib.org/lapack/)) or [MKL](http://software.intel.com/en-us/intel-mkl/) (only if you want to build Dynare++)
- An implementation of [POSIX Threads](http://en.wikipedia.org/wiki/POSIX_Threads) (optional, for taking advantage of multi-core)
- [MAT File I/O library](http://sourceforge.net/projects/matio/) (if you want to compile Markov-Switching code, the estimation DLL, k-order DLL and Dynare++ in unstable)
- [SLICOT](http://www.slicot.org) (if you want to compile the Kalman steady state DLL)
- [GSL library](http://www.gnu.org/software/gsl/) (if you want to compile Markov-Switching code)
- A decent LaTeX distribution (if you want to compile PDF documentation). The following extra components may be needed:
- The Econometrica bibliography style: you need [harvard](http://www.ctan.org/tex-archive/macros/latex/contrib/harvard/) and [economic](http://www.ctan.org/tex-archive/biblio/bibtex/contrib/economic/) packages from CTAN (only if you want to build Dynare user guide, no more needed with Dynare unstable)
- [Eplain](http://www.tug.org/eplain/) TeX macros (only if you want to build Dynare++ source documentation)
- [Beamer](http://latex-beamer.sourceforge.net/) (for some PDF presentations)
- For building the reference manual:
- [GNU Texinfo](http://www.gnu.org/software/texinfo/)
- [Texi2HTML](http://www.nongnu.org/texi2html) and [Latex2HTML](http://www.latex2html.org), if you want nice mathematical formulas in HTML output
- [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 (only in unstable)
### 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 http://github.com/DynareTeam/dynare.git
cd dynare
git submodule update --init
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).
### Configuring the build tree
Simply launch the configure script from a terminal:
```
./configure
```
If you have MATLAB, you need to indicate both the MATLAB location and version. For example, on GNU/Linux:
```
./configure --with-matlab=/usr/local/MATLAB/R2013a MATLAB_VERSION=8.1
```
Note that the MATLAB version can also be specified via the MATLAB family product release (R2009a, R2008b, ...).
**NB**: For MATLAB versions strictly older than 7.1, you need to explicitly give the MEX extension, via `MEXEXT` variable of the configure script (for example, `MEXEXT=dll` for Windows with MATLAB \< 7.1).
Alternatively, you can disable the compilation of MEX files for MATLAB with the `--disable-matlab` flag, and MEX files for Octave with `--disable-octave`.
You may need to specify additional options to the configure script, see the platform specific instructions below.
Note that if you don't want to compile with debugging information, you can specify the `CFLAGS` and `CXXFLAGS` variables to configure, such as:
```
./configure CFLAGS="-O3" CXXFLAGS="-O3"
```
If you want to give a try to the parallelized versions of some mex files (`A_times_B_kronecker_C` and `sparse_hessian_times_B_kronecker_C` used to get the reduced form of the second order approximation of the model) you can add the `--enable-openmp` flag, for instance:
```
./configure --with-matlab=/usr/local/matlab78 MATLAB_VERSION=7.8 --enable-openmp
```
If the configuration goes well, the script will tell you which components are correctly configured and will be built.
### Bulding
Binaries and Info documentation are built with:
```
make
```
PDF and HTML documentation are respectively built with:
```
make pdf
make html
```
The testsuites can be run with:
```
make check
```
## Debian or Ubuntu
All the prerequisites are packaged.
The easiest way to install the pre-requisites in Debian is to use Debian's dynare package and do:
```
apt-get build-dep dynare
```
Alternatively, if you want to build everything, manually install the following packages:
- `build-essential` (for gcc, g++ and make)
- `octave3.2-headers` or `liboctave-dev` (will install ATLAS)
- `libboost-graph-dev`
- `libgsl0-dev`
- `libmatio-dev`
- `libslicot-dev` and `libslicot-pic`
- `libsuitesparse-dev` (only for Unstable)
- `flex`
- `bison`
- `autoconf`
- `automake`
- `texlive`
- `texlive-publishers` (for Econometrica bibliographic style)
- `texlive-extra-utils` (for CWEB)
- `texlive-formats-extra` (for Eplain)
- `texlive-latex-extra` (for fullpage.sty)
- `latex-beamer`
- `texinfo`
- `texi2html`, `latex2html`
- `doxygen`
## Fedora
**NB**: Documentation still in progress…
- `octave-devel`
- `boost-devel`
- `gsl-devel`
- `matio-devel`
- `flex`
- `bison`
- `autoconf`
- `automake`
- `texlive`
- `texinfo`
- `texi2html`, `latex2html`
- `doxygen`
## Windows
The following instructions are compatible with MATLAB or with Octave/MinGW (as downloadable [here](http://www.dynare.org/download/octave)).
### Setting up the Compilation Environment
- First, you need to setup a Cygwin environment, following the instructions at <http://www.cygwin.com>. You need the following packages:
- `make`
- `bison`
- `flex`
- `autoconf` and `autoconf2.5`
- `automake` and `automake1.11`
- `texlive`, `texlive-collection-latexextra`, `texlive-collection-formatsextra`, `texlive-collection-publishers`
- `texinfo`
- `doxygen`
- `mingw64-i686-gcc`, `mingw64-i686-gcc-g++`, `mingw64-i686-gcc-fortran` (if you have Octave/MinGW or if you have MATLAB 32-bit)
- `mingw64-x86_64-gcc`, `mingw64-x86_64-gcc-g++`, `mingw64-x86_64-gcc-fortran` (if you have MATLAB 64-bit)
- Second, install precompiled librairies for BLAS, LAPACK, Boost and GSL:
- If you have Octave or MATLAB 32-bit, download [dynare-mingw32-libs.zip](http://www.dynare.org/DynareWiki/BuildingDynareFromSource?action=AttachFile&do=view&target=dynare-mingw32-libs.zip), and uncompress it in `c:\cygwin\usr\local\lib\mingw32`
- If you have MATLAB 64-bit, download [dynare-mingw64-libs.zip](http://www.dynare.org/DynareWiki/BuildingDynareFromSource?action=AttachFile&do=view&target=dynare-mingw64-libs.zip), and uncompress it in `c:\cygwin\usr\local\lib\mingw64`
### Compiling the preprocessor, Dynare++, the MEX for MATLAB and the documentation
Download and uncompress the Dynare source tree, lets say in `c:\cygwin\home\user\dynare`.
Launch a Cygwin shell, and enter the Dynare source tree:
```
cd dynare
```
If you retrieved the source from Git, don't forget to do:
```
autoreconf -i -s
```
Then, configure the package.
- If your MATLAB is 32-bit, let's say version R2008b installed in `c:\Program Files\MATLAB\R2008b`
```
./configure --host=i686-w64-mingw32 --with-boost=/usr/local/lib/mingw32/boost --with-blas=/usr/local/lib/mingw32/blas/libopenblas.a --with-lapack=/usr/local/lib/mingw32/lapack/liblapack.a --with-gsl=/usr/local/lib/mingw32/gsl --with-matio=/usr/local/lib/mingw32/matio --with-slicot=/usr/local/lib/mingw32/slicot --with-matlab=/cygdrive/c/Progra~1/MATLAB/R2008b MATLAB_VERSION=R2008b --disable-octave
```
- If your MATLAB is 64-bit:
```
./configure --host=x86_64-w64-mingw32 --with-boost=/usr/local/lib/mingw64/boost --with-blas=/usr/local/lib/mingw64/blas/libopenblas.a --with-lapack=/usr/local/lib/mingw64/lapack/liblapack.a --with-gsl=/usr/local/lib/mingw64/gsl --with-matio=/usr/local/lib/mingw64/matio --with-slicot=/usr/local/lib/mingw64/slicot --with-matlab=/cygdrive/c/Progra~1/MATLAB/R2008b MATLAB_VERSION=R2008b --disable-octave
```
A few remarks:
- Note that here we use `Progra~1` (the 8.3 filename) instead of `Program Files`. This is because spaces in filenames confuse the configuration scripts.
- If you dont have MATLAB, then drop the `--with-matlab` and `MATLAB_VERSION` options
- If your MATLAB is 32-bit and your Windows is 64-bit, you need to explicitly give the MEX extension, with `MEXEXT=mexw32`
Then compile everything with:
```
make all pdf html
```
This should build:
- Dynare preprocessor
- Dynare MEX files for MATLAB (provided you gave the MATLAB path to configure)
- Dynare++
- Part of the documentation
### Compiling the MEX for Octave (MinGW package)
Launch a Cygwin shell, and enter the Dynare source tree for Octave MEX:
cd dynare/mex/build/octave
Configure and make:
./configure MKOCTFILE=/usr/local/lib/mingw32/mkoctfile-win --with-boost=/usr/local/lib/mingw32/boost --with-gsl=/usr/local/lib/mingw32/gsl --with-matio=/usr/local/lib/mingw32/matio --with-slicot=/usr/local/lib/mingw32/slicot-underscore
make
## Mac OS X
- Install the Xcode Common Tools:
- Install [Xcode](http://developer.apple.com/xcode/) from the App Store
- Open Xcode
- Go to `Xcode->Preferences...`
- In the window that opens, click on the `Downloads` tab
- In the tab that appears, click on the `Components` button
- Next to `Command Line Tools`, click on `Install`
- Download [MacOSX10.6.sdk.zip](http://www.jamesgeorge.org/uploads/MacOSX10.6.sdk.zip) and unzip it in `/Applications/Xcode.app/Contents/Developer/Platforms/MacOSX.platform/Developer/SDKs`. Change the owner to be `root` and the group to be `wheel`
- Install [Homebrew](http://mxcl.github.io/homebrew/)
- Install the following brews:
```
brew install automake
brew install gsl
brew install boost
brew install gfortran
brew install matlab2tikz --HEAD
brew install libmatio --with-hdf5
brew install slicot --with-default-integer-8
```
- **(Optional)** To compile Dynare mex files for use on Octave, first install Octave following the [Simple Installation Instructions](http://wiki.octave.org/Octave_for_MacOS_X#Simple_Installation_Instructions_3). Then, you will probably also want to install graphicsmagick via Homebrew with `brew install graphicsmagick`.
- **(Optional)** To compile Dynare's documentation, first install the latest version of [MacTeX](http://www.tug.org/mactex/). Then install `doxygen` and `latex2html` via Homebrew with the following commands:
- **(On OS X 10.7 Only)** Copy [FlexLexer.h](http://www.dynare.org/DynareWiki/BuildingDynareFromSource?action=AttachFile&do=view&target=FlexLexer.h) into the `preprocessor` directory (there was an error in the `FlexLexer.h` file distributed with 10.7)
- Finally, switch to the root dynare directory. Ensure your path contains `/usr/bin:/bin:/usr/sbin:/sbin:/usr/local/bin:/opt/X11/bin:/usr/texbin:/usr/local/sbin`. Run:
- `autoconf -si`
- `./configure --with-matlab=/Applications/MATLAB_R2013a.app MATLAB_VERSION=8.1` for builds with Matlab or `./configure` for builds just using Octave
- `make`

View File

@ -212,7 +212,7 @@ if test "x$enable_octave" = "xyes"; then
AC_CHECK_PROG([OCTAVE], [octave], [octave])
fi
AM_CONDITIONAL([ENABLE_OCTAVE], [test "x$enable_octave" = "xyes"])
AM_CONDITIONAL([HAVE_OCTAVE], [test "x$enable_octave" = "xyes" -a test "x$OCTAVE" != "x"])
AM_CONDITIONAL([HAVE_OCTAVE], [test "x$enable_octave" = "xyes" -a "x$OCTAVE" != "x"])
# Enable exporting of Org files
# The clean way would be to test for Emacs, Org-mode, latex, dvipng...

@ -1 +1 @@
Subproject commit 458d9a3f23927247a7b4c3967c4dceca108c40a8
Subproject commit ed0942d17510fcc561e45d16c0e7a59f5646e597

View File

@ -160,10 +160,12 @@ Configuration
Running Dynare
* Dynare invocation::
* Dynare hooks::
* Understanding Preprocessor Error Messages::
Dynare invocation
* Dynare hooks
* Understanding Preprocessor Error Messages::
The Model file
@ -519,8 +521,7 @@ You need to download Dynare source code from the
Then you will need to recompile the pre-processor and the dynamic
loadable libraries. Please refer to
@uref{http://www.dynare.org/DynareWiki/BuildingDynareFromSource,Dynare
Wiki}.
@uref{https://github.com/DynareTeam/dynare/blob/master/README.md,README.md}.
@node Configuration
@section Configuration
@ -638,6 +639,7 @@ required by the user. Its contents is described in @ref{The Model file}.
@menu
* Dynare invocation::
* Dynare hooks::
* Understanding Preprocessor Error Messages::
@end menu
@ -814,6 +816,20 @@ during the computation.
Structure containing the various results of the computations.
@end defvr
@node Dynare hooks
@section Dynare hooks
It is possible to call pre and post dynare preprocessor hooks written as matlab scripts.
The script @file{@var{FILENAME}_pre_dynare_preprocessor_hook.m} is executed before the
call to Dynare's preprocessor, and can be used to programatically transform the mod file
that will be read by the preprocessor. The script @file{@var{FILENAME}_post_dynare_preprocessor_hook.m}
is executed just after the call to Dynare's preprocessor, and can be used to programatically
transform the files generated by Dynare's preprocessor before actual computations start. The
pre and/or post dynare preprocessor hooks are executed if and only if the aforementioned scripts
are detected in the same folder as the the model file, @file{@var{FILENAME}.mod}.
@node Understanding Preprocessor Error Messages
@section Understanding Preprocessor Error Messages
@ -4192,7 +4208,7 @@ the convergence diagnostics are based on comparing pooled and within
MCMC moments (Dynare displays the second and third order moments, and
the length of the Highest Probability Density interval covering 80% of
the posterior distribution). Due to computational reasons, the
multivariate convergence diagnostic does not follow of @cite{Brooks and
multivariate convergence diagnostic does not follow @cite{Brooks and
Gelman (1998)} strictly, but rather applies their idea for univariate
convergence diagnostics to the range of the posterior likelihood
function instead of the individual parameters. The posterior kernel is
@ -4446,6 +4462,9 @@ which can be useful if the posterior mode is close to a domain
boundary, or in conjunction with @code{mode_check_neighbourhood_size =
Inf} when the domain in not the entire real line. Default: @code{1}.
@item mode_check_number_of_points = @var{INTEGER}
Number of points around the posterior mode where the posterior kernel is evaluated (for each parameter). Default is @code{20}
@item prior_trunc = @var{DOUBLE}
@anchor{prior_trunc} Probability of extreme values of the prior
density that is ignored when computing bounds for the

View File

@ -67,10 +67,10 @@ switch a.freq
c.time(1) = a.time(1) + b - 1;
case {4,12,52}
c = a;
n1 = b;
n2 = floor(n1/a.freq);
n3 = mod(n1,a.freq);
c.time(2) = c.time(2)+n3-1;
n1 = b + a.time(2);
n2 = floor((n1 - 1)/a.freq);
n3 = mod(n1 - 1,a.freq) + 1;
c.time(2) = n3;
c.time(1) = c.time(1)+n2;
otherwise
error('dynDate::plus: Unknown frequency!')
@ -94,9 +94,9 @@ end
%$
%$ % Expected results.
%$ e1 = dynDate(1952);
%$ e2 = dynDate('1951Q4');
%$ e3 = dynDate('2001M5');
%$ e4 = dynDate('2000M12');
%$ e2 = dynDate('1952Q1');
%$ e3 = dynDate('2001M6');
%$ e4 = dynDate('2001M1');
%$
%$ % Check the results.
%$ t(1) = dyn_assert(e1.time,d1.time);

View File

@ -56,7 +56,6 @@ time_range_of_b = b.init:b.init+b.nobs;
last_a = time_range_of_a(a.nobs);
last_b = time_range_of_b(b.nobs);
common_time_range = intersect(time_range_of_a,time_range_of_b);
if isempty(common_time_range)
@ -85,9 +84,10 @@ if last_a>last_b
end
if last_a<last_b
n = last_a-last_b;
a.data = [a.data; NaN(n,b.vobs)];
n = last_b-last_a;
a.data = [a.data; NaN(n,a.vobs)];
a.nobs = a.nobs+n;
return
end
%@test:1
@ -120,4 +120,68 @@ end
%$ t(5) = dyn_assert(ts2.data,[B; NaN(4,2)], 1e-15);
%$ end
%$ T = all(t);
%@eof:1
%@eof:1
%@test:2
%$ % Define a datasets.
%$ A = rand(8,3); B = rand(7,2);
%$
%$ % Define names
%$ A_name = {'A1';'A2';'A3'};
%$ B_name = {'B1';'B2'};
%$
%$ % Define initial dates
%$ A_init = '1990Q1';
%$ B_init = '1990Q1';
%$
%$ % Instantiate two dynSeries objects
%$ ts1 = dynSeries(A,A_init,A_name);
%$ ts2 = dynSeries(B,B_init,B_name);
%$
%$ try
%$ [ts1, ts2] = align(ts1, ts2);
%$ t(1) = 1;
%$ catch
%$ t(1) = 0;
%$ end
%$
%$ if t(1)
%$ t(2) = dyn_assert(ts1.nobs,ts2.nobs);
%$ t(3) = dyn_assert(ts1.init==ts2.init,1);
%$ t(4) = dyn_assert(ts1.data,A, 1e-15);
%$ t(5) = dyn_assert(ts2.data,[B; NaN(1,2)], 1e-15);
%$ end
%$ T = all(t);
%@eof:2
%@test:3
%$ % Define a datasets.
%$ A = rand(8,3); B = rand(7,2);
%$
%$ % Define names
%$ A_name = {'A1';'A2';'A3'};
%$ B_name = {'B1';'B2'};
%$
%$ % Define initial dates
%$ A_init = '1990Q1';
%$ B_init = '1990Q1';
%$
%$ % Instantiate two dynSeries objects
%$ ts1 = dynSeries(A,A_init,A_name);
%$ ts2 = dynSeries(B,B_init,B_name);
%$
%$ try
%$ [ts2, ts1] = align(ts2, ts1);
%$ t(1) = 1;
%$ catch
%$ t(1) = 0;
%$ end
%$
%$ if t(1)
%$ t(2) = dyn_assert(ts1.nobs,ts2.nobs);
%$ t(3) = dyn_assert(ts1.init==ts2.init,1);
%$ t(4) = dyn_assert(ts1.data,A, 1e-15);
%$ t(5) = dyn_assert(ts2.data,[B; NaN(1,2)], 1e-15);
%$ end
%$ T = all(t);
%@eof:3

View File

@ -106,6 +106,8 @@ for i=1:nargin-1
error('dynSeries::extract: Cannot handle more than two regular expressions!')
end
end
% Remove trailing white spaces if any
VariableName_ = strtrim(VariableName_);
end
% Get indices of the selected variables

View File

@ -40,7 +40,7 @@ function A = minus(B,C)
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
if ~isequal(B.vobs,C.vobs) && ~(isequal(B.vobs,1) || isequal(C.vobs,1))
error(['dynSeries::plus: Cannot substract ' inputname(1) ' and ' inputname(2) ' (wrong number of variables)!'])
error(['dynSeries::minus: Cannot substract ' inputname(1) ' and ' inputname(2) ' (wrong number of variables)!'])
else
if B.vobs>C.vobs
idB = 1:B.vobs;

View File

@ -39,8 +39,6 @@ function A = plus(B,C)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
% AUTHOR(S) stephane DOT adjemian AT univ DASH lemans DOT fr
if ~isequal(B.vobs,C.vobs) && ~(isequal(B.vobs,1) || isequal(C.vobs,1))
error(['dynSeries::plus: Cannot add ' inputname(1) ' and ' inputname(2) ' (wrong number of variables)!'])
else
@ -173,3 +171,51 @@ A.data = bsxfun(@plus,B.data,C.data);
%$ end
%$ T = all(t);
%@eof:3
%@test:4
%$ t = zeros(7,1);
%$
%$ try
%$ ts = dynSeries(transpose(1:5),'1950q1',{'Output'}, {'Y_t'});
%$ us = dynSeries(transpose(1:5),'1949q4',{'Consumption'}, {'C_t'});
%$ vs = ts+us;
%$ t(1) = 1;
%$ catch
%$ t = 0;
%$ end
%$
%$ if length(t)>1
%$ t(2) = dyn_assert(ts.freq,4);
%$ t(3) = dyn_assert(us.freq,4);
%$ t(4) = dyn_assert(ts.init.time,[1950, 1]);
%$ t(5) = dyn_assert(us.init.time,[1949, 4]);
%$ t(6) = dyn_assert(vs.init.time,[1949, 4]);
%$ t(7) = dyn_assert(vs.nobs,6);
%$ end
%$
%$ T = all(t);
%@eof:4
%@test:5
%$ t = zeros(7,1);
%$
%$ try
%$ ts = dynSeries(transpose(1:5),'1950q1',{'Output'}, {'Y_t'});
%$ us = dynSeries(transpose(1:7),'1950q1',{'Consumption'}, {'C_t'});
%$ vs = ts+us;
%$ t(1) = 1;
%$ catch
%$ t = 0;
%$ end
%$
%$ if length(t)>1
%$ t(2) = dyn_assert(ts.freq,4);
%$ t(3) = dyn_assert(us.freq,4);
%$ t(4) = dyn_assert(ts.init.time,[1950, 1]);
%$ t(5) = dyn_assert(us.init.time,[1950, 1]);
%$ t(6) = dyn_assert(vs.init.time,[1950, 1]);
%$ t(7) = dyn_assert(vs.nobs,7);
%$ end
%$
%$ T = all(t);
%@eof:5

View File

@ -62,9 +62,9 @@ elseif TotalNumberOfMhFiles == 2 && FirstMhFile > 1
record.KeepedDraws.Distribution = [MAX_nruns-FirstLine+1 ; record.MhDraws(end,3)];
end
save([DirectoryName '/' M_.fname '_mh_history.mat'],'record');
fprintf('MH: Total number of Mh draws: %d.\n',TotalNumberOfMhDraws);
fprintf('MH: Total number of generated Mh files: %d.\n',TotalNumberOfMhFiles);
fprintf('MH: Total number of MH draws: %d.\n',TotalNumberOfMhDraws);
fprintf('MH: Total number of generated MH files: %d.\n',TotalNumberOfMhFiles);
fprintf('MH: I''ll use mh-files %d to %d.\n',FirstMhFile,TotalNumberOfMhFiles);
fprintf('MH: In mh-file number %d i''ll start at line %d.\n',FirstMhFile,FirstLine);
fprintf('MH: In MH-file number %d I''ll start at line %d.\n',FirstMhFile,FirstLine);
fprintf('MH: Finally I keep %d draws.\n',TotalNumberOfMhDraws-FirstDraw);
disp(' ');

View File

@ -224,10 +224,9 @@ if ~isinf(dsge_prior_weight)% Evaluation of the likelihood of the dsge-var model
tmp1 = dsge_prior_weight*DynareDataset.info.ntobs*GYX + mYX;
tmp2 = inv(dsge_prior_weight*DynareDataset.info.ntobs*GXX+mXX);
SIGMAu = tmp0 - tmp1*tmp2*tmp1'; clear('tmp0');
if ~ispd(SIGMAu)
v = diag(SIGMAu);
k = find(v<0);
fval = objective_function_penalty_base + sum(v(k).^2);
[SIGMAu_is_positive_definite, penalty] = ispd(SIGMAu)
if ~SIGMAu_is_positive_definite
fval = objective_function_penalty_base + penalty;
info = 52;
exit_flag = 0;
return;

View File

@ -1,4 +1,4 @@
function datatomfile (s,var_list)
function datatomfile (s,var_list, names)
% function datatomfile (s,var_list)
% This optional command saves the simulation results in a text file. The name of each
% variable preceeds the corresponding results. This command must follow SIMUL.
@ -6,6 +6,7 @@ function datatomfile (s,var_list)
% INPUTS
% s: data file name
% var_list: vector of selected endogenous variables
% names: vector of strings (alternative names for the endogenous variables in the data file)
%
% OUTPUTS
% none
@ -32,14 +33,26 @@ function datatomfile (s,var_list)
global M_ oo_
% Open the data file.
sm=[s,'.m'];
fid=fopen(sm,'w') ;
if nargin < 2 || size(var_list,1) == 0
var_list = M_.endo_names(1:M_.orig_endo_nbr,:);
end
n = size(var_list,1);
ivar=zeros(n,1);
if nargin==3
if ~isequal(size(var_list,1),n)
error('datatomfile:: Second and third arguments must have the same number of rows (variables)!')
end
else
names = var_list;
end
% Get indices for the endogenous variables.
for i=1:n
i_tmp = strmatch(var_list(i,:),M_.endo_names,'exact');
if isempty(i_tmp)
@ -49,15 +62,15 @@ for i=1:n
end
end
% Save the selected data.
for i = 1:n
fprintf(fid,[M_.endo_names(ivar(i),:), '=['],'\n') ;
fprintf(fid,[strtrim(names(i,:)), ' = ['],'\n') ;
fprintf(fid,'\n') ;
fprintf(fid,'%15.8g\n',oo_.endo_simul(ivar(i),:)') ;
fprintf(fid,'\n') ;
fprintf(fid,'];\n') ;
fprintf(fid,'\n') ;
end
% Close the data file.
fclose(fid) ;

View File

@ -206,36 +206,24 @@ Q = Model.Sigma_e;
H = Model.H;
% Test if Q is positive definite.
if EstimatedParameters.ncx
% Try to compute the cholesky decomposition of Q (possible iff Q is positive definite)
[CholQ,testQ] = chol(Q);
if testQ
% The variance-covariance matrix of the structural innovations is not definite positive. We have to compute the eigenvalues of this matrix in order to build the endogenous penalty.
a = diag(eig(Q));
k = find(a < 0);
if k > 0
fval = objective_function_penalty_base+sum(-a(k));
exit_flag = 0;
info = 43;
return
end
if ~issquare(Q) && EstimatedParameters.ncx
[Q_is_positive_definite, penalty] = ispd(Q);
if ~Q_is_positive_definite
fval = objective_function_penalty_base+penalty;
exit_flag = 0;
info = 43;
return
end
end
% Test if H is positive definite.
if EstimatedParameters.ncn
% Try to compute the cholesky decomposition of H (possible iff H is positive definite)
[CholH,testH] = chol(H);
if testH
% The variance-covariance matrix of the measurement errors is not definite positive. We have to compute the eigenvalues of this matrix in order to build the endogenous penalty.
a = diag(eig(H));
k = find(a < 0);
if k > 0
fval = objective_function_penalty_base+sum(-a(k));
exit_flag = 0;
info = 44;
return
end
if ~issquare(H) && EstimatedParameters.ncn
[H_is_positive_definite, penalty] = ispd(H);
if ~H_is_positive_definite
fval = objective_function_penalty_base+penalty;
exit_flag = 0;
info = 44;
return
end
end

View File

@ -96,15 +96,24 @@ if length(d) == 0
error(['DYNARE: can''t open ' fname])
end
% pre-dynare-preprocessor-hook
if exist([fname(1:end-4) '_pre_dynare_preprocessor_hook.m'],'file')
eval([fname(1:end-4) '_pre_dynare_preprocessor_hook'])
end
command = ['"' dynareroot 'dynare_m" ' fname] ;
for i=2:nargin
command = [command ' ' varargin{i-1}];
end
[status, result] = system(command);
disp(result)
% post-dynare-prerocessor-hook
if exist([fname(1:end-4) '_post_dynare_preprocessor_hook.m'],'file')
eval([fname(1:end-4) '_post_dynare_preprocessor_hook'])
end
% Save preprocessor result in logfile (if `no_log' option not present)
no_log = 0;
for i=2:nargin

View File

@ -399,6 +399,9 @@ if ~isequal(options_.mode_compute,0) && ~options_.mh_posterior_mode_estimation
[x, fval, COUNTEVAL, STOPFLAG, OUT, BESTEVER] = cmaes(func2str(objective_function),xparam1,H0,options_.cmaes,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
xparam1=BESTEVER.x;
disp(sprintf('\n Objective function at mode: %f',fval))
case 10
options_.cova_compute = 0 ;
[xparam1,stdh,lb_95,ub_95,med_param] = online_auxiliary_filter(xparam1,dataset_,options_,M_,estim_params_,bayestopt_,oo_) ;
case 101
myoptions=soptions;
myoptions(2)=1e-6; %accuracy of argument
@ -492,7 +495,7 @@ if ~options_.mh_posterior_mode_estimation && options_.cova_compute
end
end
if options_.mode_check == 1 && ~options_.mh_posterior_mode_estimation
if options_.mode_check.status == 1 && ~options_.mh_posterior_mode_estimation
ana_deriv = options_.analytic_derivation;
options_.analytic_derivation = 0;
mode_check(objective_function,xparam1,hh,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
@ -572,7 +575,7 @@ if (any(bayestopt_.pshape >0 ) && options_.mh_replic) || ...
if options_.cova_compute
oo_.MC_record=feval(options_.posterior_sampling_method,objective_function,options_.proposal_distribution,xparam1,invhess,bounds,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
else
error('I Cannot start the MCMC because the hessian of the posterior kernel at the mode was not computed.')
error('I Cannot start the MCMC because the Hessian of the posterior kernel at the mode was not computed.')
end
options_.analytic_derivation = ana_deriv;
end
@ -580,7 +583,7 @@ if (any(bayestopt_.pshape >0 ) && options_.mh_replic) || ...
CutSample(M_, options_, estim_params_);
return
else
if ~options_.nodiagnostic && options_.mh_replic > 1000 && options_.mh_nblck > 1
if ~options_.nodiagnostic && options_.mh_replic > 2000 && options_.mh_nblck > 1
McMCDiagnostics(options_, estim_params_, M_);
end
%% Here i discard first half of the draws:

View File

@ -57,8 +57,11 @@ options_.solve_tolf = eps^(1/3);
options_.solve_tolx = eps^(2/3);
options_.solve_maxit = 500;
options_.mode_check_neighbourhood_size = 0.5;
options_.mode_check_symmetric_plots = 1;
options_.mode_check.status = 0;
options_.mode_check.neighbourhood_size = .5;
options_.mode_check.symmetric_plots = 1;
options_.mode_check.number_of_points = 20;
options_.mode_check.nolik = 0;
% Default number of threads for parallelized mex files.
options_.threads.kronecker.A_times_B_kronecker_C = 1;
@ -237,6 +240,9 @@ particle.resampling.number_of_partitions = 200;
particle.mixture_state_variables = 5 ;
particle.mixture_structural_shocks = 1 ;
particle.mixture_measurement_shocks = 1 ;
% Online approach
particle.liu_west_delta = 0.99 ;
particle.liu_west_chol_sigma_bar = .01 ;
% Copy ep structure in options_ global structure
options_.particle = particle;
@ -368,8 +374,7 @@ options_.mh_nblck = 2;
options_.mh_recover = 0;
options_.mh_replic = 20000;
options_.recursive_estimation_restart = 0;
options_.mode_check = 0;
options_.mode_check_nolik = 0;
options_.mode_compute = 4;
options_.mode_file = '';
options_.moments_varendo = 0;

View File

@ -1,18 +1,30 @@
function test = ispd(A)
function [test, penalty] = ispd(A)
% function test = ispd(A)
% Tests if a square matrix is positive definite.
%
% INPUTS
% o A [double] a square matrix.
%
% OUTPUTS
% o test [integer] is equal to one if A is pd, 0 otherwise.
%
% SPECIAL REQUIREMENTS
% None.
%@info:
%! @deftypefn {Function File} {[@var{test}, @var{penalty} =} ispd (@var{A})
%! @anchor{ispd}
%! @sp 1
%! Tests if the square matrix @var{A} is positive definite.
%! @sp 2
%! @strong{Inputs}
%! @sp 1
%! @table @ @var
%! @item A
%! A square matrix.
%! @end table
%! @sp 2
%! @strong{Outputs}
%! @sp 1
%! @table @ @var
%! @item test
%! Integer scalar equal to 1 if @var{A} is a positive definite sqquare matrix, 0 otherwise.
%! @item penalty
%! Absolute value of the uum of the negative eigenvalues of A. This output argument is optional.
%! @end table
%! @end deftypefn
%@eod:
% Copyright (C) 2007-2009 Dynare Team
% Copyright (C) 2007-2009, 2013 Dynare Team
%
% This file is part of Dynare.
%
@ -29,12 +41,20 @@ function test = ispd(A)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
m = length(A);% I do not test for a square matrix...
test = 1;
if ~isquare(A)
error(['ispd:: Input argument ' inputname(1) ' has to be a square matrix!'])
end
[cholA, info] = chol(A);
test = ~info;
for i=1:m
if ( det( A(1:i, 1:i) ) < 2.0*eps )
test = 0;
break
if nargout>1
penalty = 0;
if info
a = diag(eig(A));
k = find(a<0);
if k > 0
penalty = sum(-a(k));
end
end
end

45
matlab/isquare.m Normal file
View File

@ -0,0 +1,45 @@
function info = isquare(A)
%@info:
%! @deftypefn {Function File} {@var{info} =} isquare (@var{A})
%! @anchor{isquare}
%! @sp 1
%! Tests if @var{A} is square matrix.
%! @sp 2
%! @strong{Inputs}
%! @sp 1
%! @table @ @var
%! @item A
%! Matlab's array.
%! @end table
%! @sp 2
%! @strong{Outputs}
%! @sp 1
%! @table @ @var
%! @item info
%! Integer scalar equal to 1 if @var{A} is a square matrix, 0 otherwise.
%! @end table
%! @end deftypefn
%@eod:
% Copyright (C) 2013 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/>.
info = 0;
if isequal(ndims(A),2) && isequal(size(A,1),size(A,2))
info = 1;
end

View File

@ -119,6 +119,7 @@ if ~exist('OCTAVE_VERSION')
error('load_csv_file_data:: Shouldn''t arrive here');
end
freq = init.freq;
varlist = transpose(varlist);
return
end

View File

@ -62,7 +62,7 @@ fprintf(' Done!\n');
save([M_.fname '_mean.mat'],'xparam1','hh','SIGMA');
disp(' ');
disp('MH: I''m computing the posterior log marginale density (modified harmonic mean)... ');
disp('MH: I''m computing the posterior log marginal density (modified harmonic mean)... ');
detSIGMA = det(SIGMA);
invSIGMA = inv(SIGMA);
marginal = zeros(9,2);

View File

@ -84,9 +84,9 @@ if TeX
fprintf(fidTeX,' \n');
end
ll = DynareOptions.mode_check_neighbourhood_size;
ll = DynareOptions.mode_check.neighbourhood_size;
if isinf(ll),
DynareOptions.mode_check_symmetric_plots = 0;
DynareOptions.mode_check.symmetric_plots = 0;
end
for plt = 1:nbplt,
@ -111,7 +111,7 @@ for plt = 1:nbplt,
xx = x;
l1 = max(BayesInfo.lb(kk),(1-sign(x(kk))*ll)*x(kk)); m1 = 0;
l2 = min(BayesInfo.ub(kk),(1+sign(x(kk))*ll)*x(kk));
if DynareOptions.mode_check_symmetric_plots,
if DynareOptions.mode_check.symmetric_plots,
if l2<(1+ll)*x(kk)
l1 = x(kk) - (l2-x(kk));
m1 = 1;
@ -120,10 +120,10 @@ for plt = 1:nbplt,
l2 = x(kk) + (x(kk)-l1);
end
end
z1 = l1:((x(kk)-l1)/10):x(kk);
z2 = x(kk):((l2-x(kk))/10):l2;
z1 = l1:((x(kk)-l1)/(DynareOptions.mode_check.number_of_points/2)):x(kk);
z2 = x(kk):((l2-x(kk))/(DynareOptions.mode_check.number_of_points/2)):l2;
z = union(z1,z2);
if DynareOptions.mode_check_nolik==0,
if DynareOptions.mode_check.nolik==0,
y = zeros(length(z),2);
dy = priordens(xx,BayesInfo.pshape,BayesInfo.p6,BayesInfo.p7,BayesInfo.p3,BayesInfo.p4);
end
@ -135,7 +135,7 @@ for plt = 1:nbplt,
else
y(i,1) = NaN;
end
if DynareOptions.mode_check_nolik==0
if DynareOptions.mode_check.nolik==0
lnprior = priordens(xx,BayesInfo.pshape,BayesInfo.p6,BayesInfo.p7,BayesInfo.p3,BayesInfo.p4);
y(i,2) = (y(i,1)+lnprior-dy);
end
@ -153,7 +153,7 @@ for plt = 1:nbplt,
axis tight
drawnow
end
if DynareOptions.mode_check_nolik==0,
if DynareOptions.mode_check.nolik==0,
if exist('OCTAVE_VERSION'),
axes('outerposition',[0.3 0.93 0.42 0.07],'box','on'),
else

View File

@ -158,81 +158,31 @@ if (DynareOptions.mode_compute~=1) && any(xparam1>BayesInfo.ub)
return
end
% Get the diagonal elements of the covariance matrices for the structural innovations (Q) and the measurement error (H).
Model = set_all_parameters(xparam1,EstimatedParameters,Model);
Q = Model.Sigma_e;
H = Model.H;
for i=1:EstimatedParameters.nvx
k =EstimatedParameters.var_exo(i,1);
Q(k,k) = xparam1(i)*xparam1(i);
end
offset = EstimatedParameters.nvx;
if EstimatedParameters.nvn
for i=1:EstimatedParameters.nvn
k = EstimatedParameters.var_endo(i,1);
H(k,k) = xparam1(i+offset)*xparam1(i+offset);
if ~isscalar(Q) && EstimatedParameters.ncx
[Q_is_positive_definite, penalty] = ispd(Q);
if ~Q_is_positive_definite
fval = objective_function_penalty_base+penalty;
exit_flag = 0;
info = 43;
return
end
offset = offset+EstimatedParameters.nvn;
else
H = zeros(nvobs);
end
% Get the off-diagonal elements of the covariance matrix for the structural innovations. Test if Q is positive definite.
if EstimatedParameters.ncx
for i=1:EstimatedParameters.ncx
k1 =EstimatedParameters.corrx(i,1);
k2 =EstimatedParameters.corrx(i,2);
Q(k1,k2) = xparam1(i+offset)*sqrt(Q(k1,k1)*Q(k2,k2));
Q(k2,k1) = Q(k1,k2);
if ~isscalar(H) && EstimatedParameters.ncn
[H_is_positive_definite, penalty] = ispd(H);
if ~H_is_positive_definite
fval = objective_function_penalty_base+penalty;
exit_flag = 0;
info = 44;
return
end
% Try to compute the cholesky decomposition of Q (possible iff Q is positive definite)
[CholQ,testQ] = chol(Q);
if testQ
% The variance-covariance matrix of the structural innovations is not definite positive. We have to compute the eigenvalues of this matrix in order to build the endogenous penalty.
a = diag(eig(Q));
k = find(a < 0);
if k > 0
fval = objective_function_penalty_base+sum(-a(k));
exit_flag = 0;
info = 43;
return
end
end
offset = offset+EstimatedParameters.ncx;
end
% Get the off-diagonal elements of the covariance matrix for the measurement errors. Test if H is positive definite.
if EstimatedParameters.ncn
for i=1:EstimatedParameters.ncn
k1 = DynareOptions.lgyidx2varobs(EstimatedParameters.corrn(i,1));
k2 = DynareOptions.lgyidx2varobs(EstimatedParameters.corrn(i,2));
H(k1,k2) = xparam1(i+offset)*sqrt(H(k1,k1)*H(k2,k2));
H(k2,k1) = H(k1,k2);
end
% Try to compute the cholesky decomposition of H (possible iff H is positive definite)
[CholH,testH] = chol(H);
if testH
% The variance-covariance matrix of the measurement errors is not definite positive. We have to compute the eigenvalues of this matrix in order to build the endogenous penalty.
a = diag(eig(H));
k = find(a < 0);
if k > 0
fval = objective_function_penalty_base+sum(-a(k));
exit_flag = 0;
info = 44;
return
end
end
offset = offset+EstimatedParameters.ncn;
end
% Update estimated structural parameters in Mode.params.
if EstimatedParameters.np > 0
Model.params(EstimatedParameters.param_vals(:,1)) = xparam1(offset+1:end);
end
% Update Model.Sigma_e and Model.H.
Model.Sigma_e = Q;
Model.H = H;
%------------------------------------------------------------------------------
% 2. call model setup & reduction program
%------------------------------------------------------------------------------
@ -345,11 +295,9 @@ ReducedForm.StateVectorVariance = StateVectorVariance;
% 4. Likelihood evaluation
%------------------------------------------------------------------------------
DynareOptions.warning_for_steadystate = 0;
seed_norm_old = randn('state') ;
seed_uniform_old = rand('state') ;
[s1,s2] = get_dynare_random_generator_state();
LIK = feval(DynareOptions.particle.algorithm,ReducedForm,Y,[],DynareOptions);
randn('state',seed_norm_old);
rand('state',seed_uniform_old);
set_dynare_random_generator_state(s1,s2);
if imag(LIK)
likelihood = objective_function_penalty_base;
exit_flag = 0;

View File

@ -4,22 +4,22 @@ function [fOutVar,nBlockPerCPU, totCPU] = masterParallel(Parallel,fBlock,nBlock,
% computing.
% It is the top-level function called on the master computer when parallelizing a task.
% This function have two main computational startegy for manage the matlab worker (slave process).
% This function has two main computational strategies for managing the matlab worker (slave process).
% 0 Simple Close/Open Stategy:
% In this case the new matlab istances (slave process) are open when
% In this case the new Matlab instances (slave process) are open when
% necessary and then closed. This can happen many times during the
% simulation of a model.
% 1 Always Open Strategy:
% In this case we have a more sophisticated management of slave processes,
% which are no longer closed at the end of each job. The slave processes
% waits for a new job (if exist). If a slave do not receives a new job after a
% wait for a new job (if it exists). If a slave does not receive a new job after a
% fixed time it is destroyed. This solution removes the computational
% time necessary to Open/Close new matlab istances.
% time necessary to Open/Close new Matlab instances.
% The first (point 0) is the default Strategy
% i.e.(Parallel_info.leaveSlaveOpen=0). This value can be changed by the
% user in xxx.mod file or it is changed by the programmer if it necessary to
% user in xxx.mod file or it is changed by the programmer if it is necessary to
% reduce the overall computational time. See for example the
% prior_posterior_statistics.m.
@ -30,7 +30,7 @@ function [fOutVar,nBlockPerCPU, totCPU] = masterParallel(Parallel,fBlock,nBlock,
% o fBlock [int] index number of the first thread
% (between 1 and nBlock)
% o nBlock [int] index number of the last thread
% o NamFileInput [cell array] containins the list of input files to be
% o NamFileInput [cell array] contains the list of input files to be
% copied in the working directory of remote slaves
% 2 columns, as many lines as there are files
% - first column contains directory paths
@ -49,8 +49,8 @@ function [fOutVar,nBlockPerCPU, totCPU] = masterParallel(Parallel,fBlock,nBlock,
% struct per thread
% o nBlockPerCPU [int vector] for each CPU used, indicates the number of
% threads run on that CPU
% o totCPU [int] total number of CPU used (can be lower than
% the number of CPU declared in "Parallel", if
% o totCPU [int] total number of CPUs used (can be lower than
% the number of CPUs declared in "Parallel", if
% the number of required threads is lower)
% Copyright (C) 2009-2013 Dynare Team
@ -106,8 +106,8 @@ if isfield(Parallel_info,'local_files')
end
end
% Deactivate some 'Parallel/Warning' message in Octave!
% Comment the line 'warning('off');' in order to view the warning message
% Deactivate some 'Parallel/Warning' messages in Octave!
% Comment the line 'warning('off');' in order to view the warning messages
% in Octave!
if exist('OCTAVE_VERSION'),

View File

@ -37,7 +37,7 @@ function [LIK,lik] = conditional_particle_filter(ReducedForm,Y,start,DynareOptio
%
% NOTES
% The vector "lik" is used to evaluate the jacobian of the likelihood.
% Copyright (C) 2009-2013 Dynare Team
% Copyright (C) 2009-2010 Dynare Team
%
% This file is part of Dynare.
%
@ -116,7 +116,7 @@ for t=1:sample_size
if (strcmp(DynareOptions.particle.resampling.status,'generic') && neff(SampleWeights)<DynareOptions.particle.resampling.neff_threshold*sample_size ) || ...
strcmp(DynareOptions.particle.resampling.status,'systematic')
ks = ks + 1 ;
StateParticles = resample(StateParticles',SampleWeights,DynareOptions)';
StateParticles = resample(StateParticles',SampleWeights',DynareOptions)';
SampleWeights = ones(1,number_of_particles)/number_of_particles ;
end
end

View File

@ -19,7 +19,7 @@ function IncrementalWeights = gaussian_densities(obs,mut_t,sqr_Pss_t_t,st_t_1,sq
%
% NOTES
% The vector "lik" is used to evaluate the jacobian of the likelihood.
% Copyright (C) 2009-2012 Dynare Team
% Copyright (C) 2009-2010 Dynare Team
%
% This file is part of Dynare.
%

View File

@ -59,7 +59,6 @@ function new_particles = multivariate_smooth_resampling(particles,weights)
% stephane DOT adjemian AT univ DASH lemans DOT fr
number_of_particles = length(weights);
weights = weights' ;
number_of_states = size(particles,2);
[P,D] = eig(particles'*(bsxfun(@times,1/number_of_particles,particles))) ;
D = diag(D) ;

View File

@ -41,7 +41,6 @@ persistent Y init_flag mf0 mf1 bounds number_of_particles number_of_parameters l
persistent start_param sample_size number_of_observed_variables number_of_structural_innovations
% Set seed for randn().
%set_dynare_seed('mt19937ar',1234) ;
set_dynare_seed('default') ;
pruning = DynareOptions.particle.pruning;
second_resample = 1 ;
@ -62,11 +61,10 @@ if isempty(init_flag)
number_of_observed_variables = length(mf1);
number_of_structural_innovations = length(ReducedForm.Q);
liu_west_delta = DynareOptions.particle.liu_west_delta ;
liu_west_chol_sigma_bar = DynareOptions.particle.liu_west_chol_sigma_bar*eye(number_of_parameters) ;
%start_param = xparam1 ;
% Conditions initiales
%liu_west_chol_sigma_bar = bsxfun(@times,eye(number_of_parameters),BayesInfo.p2) ;
start_param = BayesInfo.p1 ;
%liu_west_chol_sigma_bar = DynareOptions.particle.liu_west_chol_sigma_bar*eye(number_of_parameters) ;
start_param = xparam1 ;
%liu_west_chol_sigma_bar = sqrt(bsxfun(@times,eye(number_of_parameters),BayesInfo.p2)) ;
%start_param = BayesInfo.p1 ;
bounds = [BayesInfo.lb BayesInfo.ub] ;
init_flag = 1;
end
@ -87,11 +85,12 @@ small_a = sqrt(1-h_square) ;
% Initialization of parameter particles
xparam = zeros(number_of_parameters,number_of_particles) ;
stderr = sqrt(bsxfun(@power,bounds(:,2)+bounds(:,1),2)/12)/1000 ;
stderr = sqrt(bsxfun(@power,bounds(:,2)+bounds(:,1),2)/12)/100 ;
stderr = sqrt(bsxfun(@power,bounds(:,2)+bounds(:,1),2)/12)/50 ;
i = 1 ;
while i<=number_of_particles
candidate = start_param + 10*liu_west_chol_sigma_bar*randn(number_of_parameters,1) ;
%candidate = start_param + bsxfun(@times,stderr,randn(number_of_parameters,1)) ;
%candidate = start_param + .001*liu_west_chol_sigma_bar*randn(number_of_parameters,1) ;
candidate = start_param + bsxfun(@times,stderr,randn(number_of_parameters,1)) ;
if all(candidate(:) >= bounds(:,1)) && all(candidate(:) <= bounds(:,2))
xparam(:,i) = candidate(:) ;
i = i+1 ;
@ -113,7 +112,7 @@ ub95_xparam = zeros(number_of_parameters,sample_size) ;
%% The Online filter
for t=1:sample_size
disp(t)
disp(t)
% Moments of parameters particles distribution
m_bar = xparam*(weights') ;
temp = bsxfun(@minus,xparam,m_bar) ;
@ -210,7 +209,7 @@ for t=1:sample_size
if (variance_update==1) && (neff(weights)<DynareOptions.particle.resampling.neff_threshold*sample_size)
variance_update = 0 ;
end
% final resampling (advised)
% final resampling (advised)
if second_resample==1
indx = index_resample(0,weights,DynareOptions);
StateVectors = StateVectors(:,indx) ;
@ -259,7 +258,8 @@ for t=1:sample_size
end
disp([lb95_xparam(:,t) mean_xparam(:,t) ub95_xparam(:,t)])
end
xparam = mean_xparam(:,sample_size) ;
distrib_param = xparam ;
xparam = mean_xparam(:,sample_size) ;
std_param = std_xparam(:,sample_size) ;
lb_95 = lb95_xparam(:,sample_size) ;
ub_95 = ub95_xparam(:,sample_size) ;
@ -345,8 +345,8 @@ for plt = 1:nbplt,
TeXNAMES = char(TeXNAMES,texname);
end
end
optimal_bandwidth = mh_optimal_bandwidth(xparam(kk,:)',number_of_particles,bandwidth,kernel_function);
[density(:,1),density(:,2)] = kernel_density_estimate(xparam(kk,:)',number_of_grid_points,...
optimal_bandwidth = mh_optimal_bandwidth(distrib_param(kk,:)',number_of_particles,bandwidth,kernel_function);
[density(:,1),density(:,2)] = kernel_density_estimate(distrib_param(kk,:)',number_of_grid_points,...
number_of_particles,optimal_bandwidth,kernel_function);
plot(density(:,1),density(:,2));
hold on
@ -370,3 +370,4 @@ for plt = 1:nbplt,
fprintf(fidTeX,' \n');
end
end

View File

@ -54,6 +54,7 @@ function resampled_particles = resample(particles,weights,DynareOptions)
% AUTHOR(S) frederic DOT karame AT univ DASH evry DOT fr
% stephane DOT adjemian AT univ DASH lemans DOT fr
switch DynareOptions.particle.resampling.method1
case 'residual'
if strcmpi(DynareOptions.particle.resampling.method2,'kitagawa')

View File

@ -101,7 +101,7 @@ ghxu = ReducedForm.ghxu;
% Get covariance matrices.
Q = ReducedForm.Q; % Covariance matrix of the structural innovations.
H = ReducedForm.H; % COvariance matrix of the measurement errors.
H = ReducedForm.H; % Covariance matrix of the measurement errors.
if isempty(H)
H = 0;
end
@ -123,7 +123,6 @@ state_variance_rank = size(StateVectorVarianceSquareRoot,2);
% Factorize the covariance matrix of the structural innovations
Q_lower_triangular_cholesky = chol(Q)';
StateVectors = bsxfun(@plus,StateVectorVarianceSquareRoot*randn(state_variance_rank,number_of_particles),StateVectorMean) ;
% Set seed for randn().
set_dynare_seed('default');
@ -161,7 +160,7 @@ for t=1:sample_size
StateVectors = temp(:,1:number_of_state_variables)' ;
StateVectors_ = temp(:,number_of_state_variables+1:2*number_of_state_variables)';
else
StateVectors = resample(tmp(mf0,:)',weights,DynareOptions)';
StateVectors = resample(tmp(mf0,:)',weights',DynareOptions)';
end
weights = ones(1,number_of_particles)/number_of_particles;
elseif strcmp(DynareOptions.particle.resampling.status,'none')

View File

@ -168,8 +168,7 @@ end
offset = EstimatedParameters.nvx;
if EstimatedParameters.nvn
for i=1:EstimatedParameters.nvn
k = EstimatedParameters.var_endo(i,1);
H(k,k) = xparam1(i+offset)*xparam1(i+offset);
H(i,i) = xparam1(i+offset)*xparam1(i+offset);
end
offset = offset+EstimatedParameters.nvn;
else

View File

@ -57,7 +57,7 @@ if dataset_.info.nvobs-size(rawdata,2)
end
if size(rawdata,1)~=dataset_.info.ntobs
fprintf('Restricting the sample to observations %d to %d. Using in total %d observations. \n',first,dataset_.info.ntobs,dataset_.info.ntobs-first+1)
fprintf('Restricting the sample to observations %d to %d. Using in total %d observations. \n',first,first+dataset_.info.ntobs-1,dataset_.info.ntobs)
end
rawdata = rawdata(first:(first+dataset_.info.ntobs-1),:);

View File

@ -112,7 +112,7 @@ class ParsingDriver;
%token LABELS LAPLACE LIK_ALGO LIK_INIT LINEAR LOAD_IDENT_FILES LOAD_MH_FILE LOAD_PARAMS_AND_STEADY_STATE LOGLINEAR LYAPUNOV
%token LYAPUNOV_FIXED_POINT_TOL LYAPUNOV_DOUBLING_TOL LYAPUNOV_SQUARE_ROOT_SOLVER_TOL LOG_DEFLATOR LOG_TREND_VAR LOG_GROWTH_FACTOR MARKOWITZ MARGINAL_DENSITY MAX MAXIT
%token MFS MH_DROP MH_INIT_SCALE MH_JSCALE MH_MODE MH_NBLOCKS MH_REPLIC MH_RECOVER MIN MINIMAL_SOLVING_PERIODS SOLVE_MAXIT
%token MODE_CHECK MODE_CHECK_NEIGHBOURHOOD_SIZE MODE_CHECK_SYMMETRIC_PLOTS MODE_COMPUTE MODE_FILE MODEL MODEL_COMPARISON MODEL_INFO MSHOCKS ABS SIGN
%token MODE_CHECK MODE_CHECK_NEIGHBOURHOOD_SIZE MODE_CHECK_SYMMETRIC_PLOTS MODE_CHECK_NUMBER_OF_POINTS MODE_COMPUTE MODE_FILE MODEL MODEL_COMPARISON MODEL_INFO MSHOCKS ABS SIGN
%token MODEL_DIAGNOSTICS MODIFIEDHARMONICMEAN MOMENTS_VARENDO DIFFUSE_FILTER SUB_DRAWS
%token <string_val> NAME
%token NAN_CONSTANT NO_STATIC NOBS NOCONSTANT NODISPLAY NOCORR NODIAGNOSTIC NOFUNCTIONS
@ -1519,6 +1519,7 @@ estimation_options : o_datafile
| o_mode_check
| o_mode_check_neighbourhood_size
| o_mode_check_symmetric_plots
| o_mode_check_number_of_points
| o_prior_trunc
| o_mh_mode
| o_mh_nblocks
@ -2354,9 +2355,10 @@ o_mh_init_scale : MH_INIT_SCALE EQUAL non_negative_number { driver.option_num("m
o_mode_file : MODE_FILE EQUAL filename { driver.option_str("mode_file", $3); };
o_mode_compute : MODE_COMPUTE EQUAL INT_NUMBER { driver.option_num("mode_compute", $3); };
| MODE_COMPUTE EQUAL symbol { driver.option_str("mode_compute", $3); };
o_mode_check : MODE_CHECK { driver.option_num("mode_check", "1"); };
o_mode_check_neighbourhood_size : MODE_CHECK_NEIGHBOURHOOD_SIZE EQUAL signed_number_w_inf { driver.option_num("mode_check_neighbourhood_size", $3); };
o_mode_check_symmetric_plots : MODE_CHECK_SYMMETRIC_PLOTS EQUAL INT_NUMBER { driver.option_num("mode_check_symmetric_plots", $3); };
o_mode_check : MODE_CHECK { driver.option_num("mode_check.status", "1"); };
o_mode_check_neighbourhood_size : MODE_CHECK_NEIGHBOURHOOD_SIZE EQUAL signed_number_w_inf { driver.option_num("mode_check.neighbourhood_size", $3); };
o_mode_check_number_of_points : MODE_CHECK_NUMBER_OF_POINTS EQUAL INT_NUMBER { driver.option_num("mode_check.number_of_points", $3); };
o_mode_check_symmetric_plots : MODE_CHECK_SYMMETRIC_PLOTS EQUAL INT_NUMBER { driver.option_num("mode_check.symmetric_plots", $3); };
o_prior_trunc : PRIOR_TRUNC EQUAL non_negative_number { driver.option_num("prior_trunc", $3); };
o_mh_mode : MH_MODE EQUAL INT_NUMBER { driver.option_num("mh_mode", $3); };
o_mh_nblocks : MH_NBLOCKS EQUAL INT_NUMBER { driver.option_num("mh_nblck", $3); };

View File

@ -241,6 +241,7 @@ string eofbuff;
<DYNARE_STATEMENT>mode_check {return token::MODE_CHECK;}
<DYNARE_STATEMENT>mode_check_neighbourhood_size {return token::MODE_CHECK_NEIGHBOURHOOD_SIZE;}
<DYNARE_STATEMENT>mode_check_symmetric_plots {return token::MODE_CHECK_SYMMETRIC_PLOTS;}
<DYNARE_STATEMENT>mode_check_number_of_points {return token::MODE_CHECK_NUMBER_OF_POINTS;}
<DYNARE_STATEMENT>prior_trunc {return token::PRIOR_TRUNC;}
<DYNARE_STATEMENT>mh_mode {return token::MH_MODE;}
<DYNARE_STATEMENT>mh_nblocks {return token::MH_NBLOCKS;}

View File

@ -35,12 +35,12 @@ steady;
//disp(oo_.mean) ;
estimated_params;
alp, uniform_pdf,,, 0.0001, 1;
bet, uniform_pdf,,, 0.75, 0.999;
alp, uniform_pdf,,, 0.0001, 0.5;
bet, uniform_pdf,,, 0.0001, 0.99;
tet, uniform_pdf,,, 0.0001, 1;
tau, uniform_pdf,,, 0.0001, 100;
delt, uniform_pdf,,, 0.0001, 0.05;
rho, uniform_pdf,,, 0.0001, 0.999;
rho, uniform_pdf,,, 0.8, 0.99;
stderr e_a, uniform_pdf,,, 0.00001, 0.1;
stderr y, uniform_pdf,,, 0.00001, 0.1;
stderr l, uniform_pdf,,, 0.00001, 0.1;
@ -49,11 +49,11 @@ end;
estimated_params_init;
alp, 0.4;
bet, 0.99;
bet, 0.97;
tet, 0.357 ;
tau, 50;
delt, 0.02;
rho, 0.95;
rho, 0.9 ;
stderr e_a, .035;
stderr y, .0175;//.00158;
stderr l, .00312;//.0011;
@ -66,28 +66,42 @@ varobs y l i ;
//options_.gstep(2) = .1;
options_.particle.status = 1;
options_.particle.algorithm = 'sequential_importance_particle_filter';
options_.particle.initialization = 1;
options_.particle.pruning = 1;
options_.particle.number_of_particles = 2000;
options_.particle.pruning = 0;
options_.particle.number_of_particles = 1000 ;
options_.particle.resampling.status = 'systematic';
options_.particle.resampling.neff_threshold = .1;
//options_.particle.resampling.method1 = 'traditional' ;
//options_.particle.resampling.method1 = 'residual' ;
options_.particle.resampling.method1 = 'smooth' ;
options_.particle.reampling.method2 = 'kitagawa' ;
//options_.particle.resampling.method2 = 'stratified' ;
options_.particle.resampling.number_of_partitions = 1;
options_.particle.resampling.neff_threshold = .5;
set_dynare_threads('local_state_space_iteration_2',3);
options_.particle.algorithm = 'sequential_importance_particle_filter';
//options_.particle.algorithm = 'auxiliary_particle_filter';
//options_.particle.algorithm = 'gaussian_mixture_filter';
//options_.particle.algorithm = 'each_gaussian_filter';
//options_.particle.algorithm = 'conditional_particle_filter';
//options_.particle.algorithm = 'gaussian_filter';
//options_.particle.IS_approximation_method = 'quadrature' ;
//options_.particle.IS_approximation_method = 'cubature' ;
options_.particle.IS_approximation_method = 'cubature' ;
//options_.particle.IS_approximation_method = 'unscented' ;
//options_.particle.approximation_method = 'quadrature' ;
//options_.particle.approximation_method = 'cubature' ;
options_.particle.approximation_method = 'cubature' ;
//options_.particle.approximation_method = 'unscented' ;
//options_.particle.approximation_method = 'MonteCarlo' ;
//options_.mh_posterior_mode_estimation=1 ;
// online
options_.particle.liu_west_delta = 0.9 ;
options_.mode_check_node_number = 250 ;
estimation(datafile=data_risky_perturb2,nograph,order=2,nobs=100,mh_replic=0,mode_compute=7,mode_check);