The judge cracked with the back of an axe the shinbone on an antelope and the hot marrow dripped smoking on the stones. They watched him. The subject was war.

The good book says that he that lives by the sword shall perish by the sword, said the black.

The judge smiled, his face shining with grease. What right man would have it any other way? he said.

The good book does indeed count war an evil, said Irving. Yet there's many a bloody tale of war inside it.

It makes no difference what men think of war, said the judge. War endures. As well ask men what they think of stone. War was always here. Before man was, war waited for him. The ultimate trade awaiting its ultimate practitioner. That is the way it was and will be. That way and not some other way.

addpath(genpath('F:/RussiaCEE/tools/matlab/v3dot0'))

clear all; clear session; close all; clc
warning off all

% Load data
[xlsdata, xlstext] = xlsread("databq.xlsx",'Sheet1');

X = xlsdata;
dates = xlstext(2:end,1);
vnames_long = xlstext(1,2:end);
vnames = xlstext(1,2:end);
nvar = length(vnames);
data   = Num2NaN(xlsdata);
% Store variables in the structure DATA
for ii=1:length(vnames)
    DATA.(vnames{ii}) = data(:,ii);
end
% Convert the first date to numeric
year = str2double(xlstext{3,1}(1:4));
quarter = str2double(xlstext{3,1}(6));
% Observations
nobs = size(data,1);

% Define number of variables and of observations
[nobs, nvar] = size(X);
% Set deterministics for the VAR
det = 1;
% Set number of nlags
nlags = 12;
% Estimate VAR
[VAR, VARopt] = VARmodel(X,nlags,det);

VARopt.nsteps = 40;
VARopt.ident = 'long';
VARopt.vnames = vnames_long;
VARopt.FigSize = [26,12];

% Compute IR
[IR, VAR] = VARir(VAR,VARopt);

[HD, VAR] = VARhd(VAR,VARopt)

tsio = HD.init(:,1) + HD.const(:,1) + HD.shock(:,1,1);

dfoi = TimeSeriesImport["databq.xlsx"]

val = Values[dfoi];
data = Transpose[#@"Values" & /@ val];
dates = First[val]["Dates"];

varlag = 12
varm = VARModelFit[data, varlag];

svar = VARImpulseResponse[varm, Length@data - varm["nlags"],   "BlanchardQuah", 0];
irf = svar["IRF"];
irfs = Chop@irf[[ ;; ]];
structuralResiduals = 
  varm["FitResiduals"] . Inverse[svar["invA"]]\[Transpose];

empirf = irfs[[All, 1, All]];
histDecomp = ((Plus @@ (Reverse[empirf[[;; #]]]*
         structuralResiduals[[;; #]]))) & /@ 
   Range[(Length@data - varm["nlags"])];
resid = Total /@ histDecomp;


bchart = (Join[
    histDecomp, {yd[[-(Length@data - varm["nlags"]) ;;]] - 
       resid}\[Transpose], 2]);

trendbq = 
 TimeSeries[{dates[[varlag +1 ;;]], 
    Total /@ (bchart[[All, {1, 5}]])}\[Transpose]]
tsoi = TimeSeriesIntegrate[trendbq]