File: xydata.m

package info (click to toggle)
dynare 4.6.3-4
  • links: PTS, VCS
  • area: main
  • in suites: bullseye
  • size: 74,896 kB
  • sloc: cpp: 98,057; ansic: 28,929; pascal: 13,844; sh: 5,947; objc: 4,236; yacc: 4,215; makefile: 2,583; lex: 1,534; fortran: 877; python: 647; ruby: 291; lisp: 152; xml: 22
file content (188 lines) | stat: -rw-r--r-- 8,393 bytes parent folder | download | duplicates (8)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
function [phi,y,fss,xtx,xty,yty,Bols] = xydata(xdgel,lags,mu)
% [phi,y,xtx,xty,yty,Bols,fss] = xydata(xdgel,lags,mu)
%  Construct matrices X, Y, X'X, etc. so that y = X*B + U, where y is T-by-1, X T-by-ncoef,
%     and B is ncoef-by-1.
%
% xdgel: the general matrix of the original data (no manipulation involved)
%                     with sample size including lags
% lags:  the lag length in the AR(p) process
% mu: 2-by-1 vector of hyperparameters for dummy observations (FRA's parameter settings)
%       mu(1): weight on nvar sums of coeffs dummy observations (unit roots);  (5)
%       mu(2): weight on single dummy initial observation including constant
%               (cointegration, unit roots, and stationarity);  (5)
%---------------
% phi:  X: T-by-ncoef where ncoef=nvar*lags + constant and T=fss
% y:  Y: T-by-nvar where T=fss
% fss: effective sample size (including dummies if mu is provided) exclusing all lags
% xtx:  X'X (Data)
% xty:  X'Y (Data)
% yty:  Y'Y (Data)
% Bols:  ncoef-by-1: OLS estimate of B
%
%function [Gb,Sbd,Bh,SpH,fss,ndobs,phi,y,nvar,ncoef,xxhpc,a0indx,na0p,...
%             idmat0,idmatpp] = szasbvar(idfile,q_m,lags,xdgel,mu)
%  Now version: [Gb,Sbd,Bh,SpH,fss,ndobs,phi,y,nvar,ncoef,xxhpc,a0indx,na0p,...
%             idmat0,idmatpp] = szasbvar(idfile,q_m,lags,xdgel,mu)
%
%    Estimating the Bayesian VAR of Sims and Zha with asymmetric prior (as)
%
% idfile:  Identification filename with rows corresponding to equations.
%              But all the output will be transposed to columns
%              corresponding to equations.
% q_m:  quarter or month
% lags: the maximum length of lag
% nhp:  number of haperparameters
% xdgel: the general matrix of the original data (no manipulation involved)
%                    with sample size including lags
% mu: 6-by-1 vector of hyperparameters (the following numbers for Atlanta Fed's forecast)
%       mu(1): overall tightness and also for A0;  (0.57)
%       mu(2): relative tightness for A+;  (0.13)
%       mu(3): relative tightness for the constant term;  (0.1)
%       mu(4): tightness on lag decay;  (1)
%       mu(5): weight on nvar sums of coeffs dummy observations (unit roots);  (5)
%       mu(6): weight on single dummy initial observation including constant
%               (cointegration, unit roots, and stationarity);  (5)
% --------------------
% --------------------
% Gb: cell(nvar,1). Each cell, postmultiplied by A0(:,i), is used to compute A+(:,i)
%               where A+ is k-by-m, where k--ncnoef, m--nvar, if asymmetric prior
% Sbd: cell(nvar,1). Sbd=diag(S(1), ..., S(m)).  Already divided by 'fss' for the
%        posterior of a0 or A0(:) in Waggoner and Zha when one has asymmetric prior.
%        Note,"bd" stands for block diagonal.
% Bh: reduced form B, k(mp+1)-by-m, column--equation. Bh=NaN if asymmetric prior.
%     In the form of y=X*B+U where B=[B1|B2| ...|Bp|C]
%       where Bi: m-by-m (i=1,..,p--lag length), C: 1-by-ml, constant.
% SpH:  divided by T, the final S in the Waggoner-Zha exponential part of p(A0|y)
%               SpH=NaN if asymmetric prior
% fss:  in-sample-size (for forecasting). Including dummy observations
% ndobs: number of dummy observations, ndobs=nvar+1
% phi:  X in the form of y = X*B+U. Row: nSmaple-lags+ndobs. Column: ncoef
% y:    y in the form y = X*B+U. Including the dummy observations too,
%         T-lags+ndobs-by-nvar.
% nvar: number of variables
% ncoef: number of coefficients in *each* equation. RHS coefficients only, nvar*lags+1
% xxhpc: chol(X'X+inv(H_p_tilde)): upper triangular but its transpose
%                                      is lower triagular Choleski
% a0indx: the index number for free parameters in A0.  Column meaning equation
% na0p:   number of free parameters in A0
% idmat0:  identification matrix for A0 with asymmetric prior;  column -- equation.
% idmatpp: asymmetric prior variance for A+ without constant; column -- equation.
%
% Revisions by CAS 9/27/96:  If we transmit idmat, rather than a0indx, we can scale
%   the elements of idmat0 or idmatp so that it carries information about relative
%   prior variances.
% Revsions by TZ, 10/13/96:  Efficiency improvement by streamlining the previous
%   code according to the general setup in Sims and Zha's IER and according
%   to Zha's notes Forecast (0) (p.3) and Forecast (1) (p.9).
%
%
%  NOTE: "nSample" is something I deleted as an input recently, so it may not
%    be compatible with old programs, 10/15/98 by TAZ.
%
%  NOTE2: I put "mu" as an input arg and delete "nhp" as an input arg, so it may
%    not be compatible with old programs, 03/06/99 by TAZ
%
%
% Copyright (C) 1997-2012 Lutz Kilian and Tao Zha
%
% This 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.
%
% If you did not received a copy of the GNU General Public License
% with this software, see <http://www.gnu.org/licenses/>.
%


[nSample,nvar]=size(xdgel);   % sample size including lags and number of variables
ncoef=nvar*lags+1;  % number of coefficients in *each* equation, RHS coefficients only,
                    % with only constant
ndobs=nvar+1;    % number of dummy observations
sb = lags+1;   % the beginning of original sample without dummies

if nargin==3     % activate the option of dummy observations
   fss = nSample+ndobs-lags;
                    % Effective sample size including dummies but exclusing all lags

   %
   % **** nvar prior dummy observations with the sum of coefficients
   % ** construct X for Y = X*B + U where phi = X: (T-lags)*k, Y: (T-lags)*nvar
   % **    columns: k = # of [nvar for 1st lag, ..., nvar for last lag, const]
   % **    Now, T=T+ndobs -- added with "ndobs" dummy observations
   %
   phi = zeros(fss,ncoef);
   const = ones(fss,1);
   const(1:nvar) = zeros(nvar,1);
            % the first nvar periods: no or zero constant (designed for dummies)!
   phi(:,ncoef) = const;
   %
   xdgelint = mean(xdgel(1:lags,:),1); % mean of the first lags initial conditions
   %**** Dummies for phi
   for k=1:nvar
      for m=1:lags
         phi(ndobs,nvar*(m-1)+k) = xdgelint(k);
         phi(k,nvar*(m-1)+k) = xdgelint(k);
         % <<>> multiply hyperparameter later
      end
   end
   %* True data for phi
   for k=1:lags
      phi(ndobs+1:fss,nvar*(k-1)+1:nvar*k) = xdgel(sb-k:nSample-k,:);
      % row: T-lags; column: [nvar for 1st lag, ..., nvar for last lag, const]
      %                     Thus, # of columns is nvar*lags+1 = ncoef.
   end
   %
   % ** Y with "ndobs" dummies added
   y = zeros(fss,nvar);
   %* Dummies for y
   for k=1:nvar
      y(ndobs,k) = xdgelint(k);
      y(k,k) = xdgelint(k);
      % multiply hyperparameter later
   end
   %* True data for y
   y(ndobs+1:fss,:) = xdgel(sb:nSample,:);
   %
   %*** Dummies once again (finfal version)
   phi(1:nvar,:) = 1*mu(1)*phi(1:nvar,:);    % standard Sims and Zha prior
   y(1:nvar,:) = mu(1)*y(1:nvar,:);      % standard Sims and Zha prior
   %
   phi(nvar+1,:) = mu(2)*phi(nvar+1,:);
   y(nvar+1,:) = mu(2)*y(nvar+1,:);
else      % dummy observations
   fss = nSample-lags;
            % Effective sample size with no dummies but exclusing all lags

   %
   % **** nvar prior dummy observations with the sum of coefficients
   % ** construct X for Y = X*B + U where phi = X: (T-lags)*k, Y: (T-lags)*nvar
   % **    columns: k = # of [nvar for 1st lag, ..., nvar for last lag, const]
   % **    Now, T=T+ndobs -- added with "ndobs" dummy observations
   %
   phi = zeros(fss,ncoef);
   const = ones(fss,1);
   phi(:,ncoef) = const;
   %* True data for phi
   for k=1:lags
      phi(1:fss,nvar*(k-1)+1:nvar*k) = xdgel(sb-k:nSample-k,:);
      % row: T-lags; column: [nvar for 1st lag, ..., nvar for last lag, const]
      %                     Thus, # of columns is nvar*lags+1 = ncoef.
   end
   %
   %* True data for y
   y(1:fss,:) = xdgel(sb:nSample,:);
end
%*** data input for the posterior density
[xq,xr]=qr(phi,0);
xtx=xr'*xr;   % X'X
xty=phi'*y;  % X'Y
yty = y'*y;   % Y'Y

Bols = xr\(xr'\xty);  % inv(X'X)*X'Y