function output = var(y,nlag,const);

% PURPOSE: performs vector autogressive estimation
%--------------------------------------------------------
% USAGE:  result = var(y,nlag,const) 
% where:    y    = an (nobs x neqs) matrix of y-vectors
%           nlag = the lag length
%           const = binary: 1 = include constant 0 = no constant
%--------------------------------------------------------
% RETURNS a structure
% output.meth = 'vector autoregression'
% output.nobs = nobs, # of observations
% output.neqs = neqs, # of equations
% output.nlags = nlag, # of lags
% output.beta  = bhat 
% output.tstat = t-statistics 
% output.stder = standard errors
% output.resid = residuals 
% output.yhat  = predicted values 
% output.y     = actual values 
% output.omega  = covariance matrix
% output.rsqr  = r-squared
% output.rbar  = r-squared adjusted
%-------------------------------------------------------
% written by: Gregory Givens
%             3 March 2003
%-------------------------------------------------------

if (nargin ~= 3);
    error('Wrong # of arguments to var');
    else
 [nobs neqs] = size(y);
end;

output.meth = 'vector autoregression';

% adjust nobs to feed the lags
nobse = nobs - nlag;
output.nobs = nobse; output.nlags = nlag; output.neqs = neqs;

%create X matrix
for j = 1:nlag;
    X0(:,(j-1)*neqs+1:j*neqs) = lag(y,j);
end;
if const == 0;
    X = X0;
end;
if const == 1;
    X = [ones(nobs,1) X0];
end;
X = X(nlag+1:nobs,:);

% equation by equation least squares
for j=1:neqs;
    yvec = y(nlag+1:nobs,j);
    res = ols(yvec,X);
    beta(j,:)  = res.beta';      % bhats
    tstat(j,:) = res.tstat';     % t-stats
    stder(j,:) = res.stder';     % std errors
    resid(:,j) = res.resid;      % resids 
    yhat(:,j) = res.yhat;        % yhats
    yy(:,j) = yvec;              % actual y
    rsqr(j,1) = res.rsqr;        % r-squared
    rbar(j,1) = res.rbar;        % r-adjusted
end;
omega = (1/(nobse-size(X,2)))*(resid'*resid);
output.y = y; output.beta = beta; output.tstat = tstat;
output.stder = stder; output.resid = resid; output.yhat = yhat;
output.rsqr = rsqr; output.rbar = rbar; output.omega = omega;