function [p,sumar]=arlagselect(x,P,XX,CC,F)
%
% Luca Benati
% Bank of England
% Monetary Assessment and Strategy Division
% February 2003
%
% function p=arlagselect(x,P,XX,CC,F)
% This program selects the optimal lag order for an AR(p) process based on AIC or SIC.
%                             Input of the program is:
% x     = a Tx1 vector
% P     = the maximum order of the autoregression the algorithm has to consider
% XX    = must be equal to 'SA', if the series is seasonally adjusted (so no seasonal dummies), or equal to 'NA', if the series is
%         seasonally unadjusted (so we need seasonal dummies, and the program takes care automatically of generating the dummies).
% CC    = the criterion for selecting the lag order: must be set equal to 'AIC' or 'SIC'
% F     = the frequency, 'M' or 'Q'
%
%                             Output of the program is:
% p     = the optimal lag order
% sumar = the sum of the ar coefficients
%
if nargin<5
    F='Q';
end
%
N=length(x);
AIC=[];
SIC=[];
for j=1:P
    Y=x(j+1:N);
    NY=length(Y);
    X(:,1)=ones(NY,1);
    for h=1:j
	    X(:,h+1)=x(j+1-h:N-h);
    end
    if XX=='SA'
        X=X;
    else
        X=seasadd(X,F);
    end
    bhat=(X'*X)\(X'*Y);
    uhat=Y-X*bhat;
	s2=(uhat'*uhat)/(NY-length(bhat));
    S=s2*ones(NY,1);
	L=-0.5*log(S)-0.5*((uhat.^2)./S);
	aic=NY*log(s2)+2*j;
	%aic=(2/NY)*(length(bhat)-sum(L));
%	sic=(1/NY)*(length(bhat)*log(NY)-2*sum(L));
	sic=log(s2)+log(NY-P)/(NY-P);
    AIC=[AIC' aic]';
    SIC=[SIC' sic]';
    clear Y
    clear X
    clear L
end
%
if CC=='AIC'
    p=(1:1:P)*(AIC==min(AIC));
else
    p=(1:1:P)*(SIC==min(SIC));
end

% now re-estimate the model with optimal p and get the sum of coefficients
    j=p;
    Y=x(j+1:N);
    NY=length(Y);
    X(:,1)=ones(NY,1);
    for h=1:j
	    X(:,h+1)=x(j+1-h:N-h);
    end
    if XX=='SA'
        X=X;
        bhat=(X'*X)\(X'*Y);
        sumar = sum(bhat(2:end));
    else
        X=seasadd(X,F);
        bhat=(X'*X)\(X'*Y);
        sumar = sum(bhat(2:end-11));      
    end
