% Summary Calculations 
% May 17, 2015
%


clear all;
small = 1.0e-10;
big = 1.0e+6;                  
rng(63761);
 
  % -- File Directories  
  outdir = '/Users/mwatson/Dropbox/p_comp/ddisk/matlab/out/';
  figdir = '/Users/mwatson/Dropbox/p_comp/ddisk/matlab/fig/';
  matdir = '/Users/mwatson/Dropbox/p_comp/ddisk/matlab/mat/';
  
  outfile_name = [outdir 'Table_2.out'];
  fileID = fopen(outfile_name,'w');
  
  % -- Read in Data --- 
  load_data = 1;  % 1 if reloading data from Excel, etc 
  pcomp_data_calendar_mtoq;
  
  % Dates
  first_date = [1960 1];
  last_date = [2015 1];
  ismpl_60_end = smpl(calvec,first_date,last_date,4);
  cal_60_end = calvec(ismpl_60_end==1);
  dnobs_60_end = size(cal_60_end,1);
  
  dp_agg_60_end = dp_agg(ismpl_60_end==1,:);
  dp_agg_xe_60_end = dp_agg_xe(ismpl_60_end==1,:);
  dp_agg_xfe_60_end = dp_agg_xfe(ismpl_60_end==1,:);
  dp_disagg_60_end = dp_disagg(ismpl_60_end==1,:);
  share_avg_60_end = share_avg(ismpl_60_end==1,:);
  
  tmp = dp_agg;
  tmp1 = NaN*zeros(dnobs,1);
  tmp1(4:end) = (tmp(4:end)+tmp(3:end-1)+tmp(2:end-2)+tmp(1:end-3))/4;
  dp4_agg_60_end = tmp1(ismpl_60_end==1);
  
  tmp = dp_agg_xe;
  tmp1 = NaN*zeros(dnobs,1);
  tmp1(4:end) = (tmp(4:end)+tmp(3:end-1)+tmp(2:end-2)+tmp(1:end-3))/4;
  dp4_agg_xe_60_end = tmp1(ismpl_60_end==1);
  
  tmp = dp_agg_xfe;
  tmp1 = NaN*zeros(dnobs,1);
  tmp1(4:end) = (tmp(4:end)+tmp(3:end-1)+tmp(2:end-2)+tmp(1:end-3))/4;
  dp4_agg_xfe_60_end = tmp1(ismpl_60_end==1);
  
  tmp = dp_disagg;
  tmp1 = NaN*zeros(dnobs,17);
  tmp1(4:end,:) = (tmp(4:end,:)+tmp(3:end-1,:)+tmp(2:end-2,:)+tmp(1:end-3,:))/4;
  dp4_disagg_60_end = tmp1(ismpl_60_end==1,:);
  
  % -- Some results for the 17-component model
  mlabel = '_mucsv_tvp_17c';
  str_tmp = [matdir 'calvec_ismpl' mlabel]; load(str_tmp);
  % check dates
  tmp = calvec_ismpl-cal_60_end;
  if sum(tmp) ~= 0;
    error('Dates not synchronized');
  end;  
  % Components from 17
  s = share_avg(ismpl_60_end==1,:);
  ifood = 5;
  ienergy = [7 10];
  icore = [1:4 6 8:9 11:17];
  s_energy = s(:,ienergy);
  s_energy = s_energy./repmat(sum(s_energy,2),1,2);
  s_core = s(:,icore);
  s_core = s_core./repmat(sum(s_core,2),1,14);
  % -- Mean Values of Parameters
  str_tmp = [matdir 'alpha_eps_mean' mlabel]; load(str_tmp);
  str_tmp = [matdir 'alpha_tau_mean' mlabel]; load(str_tmp);
  str_tmp = [matdir 'var_eps_common_mean' mlabel]; load(str_tmp);
  str_tmp = [matdir 'var_eps_common_total_mean' mlabel]; load(str_tmp);
  str_tmp = [matdir 'var_eps_unique_mean' mlabel]; load(str_tmp);
  str_tmp = [matdir 'var_eps_unique_total_mean' mlabel]; load(str_tmp); 
  str_tmp = [matdir 'var_dtau_common_mean' mlabel]; load(str_tmp); 
  str_tmp = [matdir 'var_dtau_unique_mean' mlabel]; load(str_tmp);

  alpha_eps_17c = alpha_eps_mean;
  alpha_tau_17c = alpha_tau_mean;
  
  var_eps_common_17c = var_eps_common_mean;
  var_eps_common_total_17c = var_eps_common_total_mean;
  var_eps_unique_17c = var_eps_unique_mean;
  var_eps_unique_total_17c = var_eps_unique_total_mean;
  var_dtau_common_17c = var_dtau_common_mean;
  var_dtau_unique_17c = var_dtau_unique_mean;
 
  
% -- Carry out filtering and so forth to determine mean and SD from various
% tau(t/t) estimators

% Exercise 1: Use system parameters from 17-component model
alpha_eps = alpha_eps_17c;
alpha_tau = alpha_tau_17c;

sigma_eps_common = sqrt(var_eps_common_total_17c);
sigma_dtau_common = sqrt(var_dtau_common_17c);
sigma_eps_unique = sqrt(var_eps_unique_total_17c);
sigma_dtau_unique = sqrt(var_dtau_unique_17c);

n_y = 17;
share_h = share_avg_60_end;
share_xe = share_h;
share_xe(:,[ienergy]) = 0;
share_xfe = share_h;
share_xfe(:,[ifood ienergy]) = 0;
share_xe = share_xe./repmat(sum(share_xe,2),1,n_y);
share_xfe = share_xfe./repmat(sum(share_xfe,2),1,n_y);


% Set up State Vector
  % --- State Vector
  %     (1) eps(t)
  %     (2) tau(t)
  %     (3) tau_u(t)
  %     (4) eps_u(t)
  ns = 2+2*n_y;      % size of state
  Q = zeros(ns,ns);
  F = zeros(ns,ns);
  F(2,2) = 1;
  F(3:2+n_y,3:2+n_y) = eye(n_y);
  H = zeros(ns,n_y);
  H(3:2+n_y,1:n_y) = eye(n_y);
  H(2+n_y+1:2+2*n_y,1:n_y) = eye(n_y);
  
  % Set up KF to run
  % Initial conditions
  X1_init = zeros(ns,1);
  P1_init = zeros(ns,ns);
  var_dtau_unique = big*ones(n_y,1);
  P1_init(3:2+n_y,3:2+n_y) = diag(var_dtau_unique);
  
  tau_17c_mean = zeros(dnobs_60_end,3);
  tau_17c_var = zeros(dnobs_60_end,3);
  X1 = X1_init;
  P1 = P1_init;
  for t = 1:dnobs_60_end;
    H_t = H;
    H_t(1,:) = alpha_eps(t,:);
    H_t(2,:) = alpha_tau(t,:);
    Q_t = diag([(sigma_eps_common(t,1).^2) (sigma_dtau_common(t,1).^2) (sigma_dtau_unique(t,:).^2) (sigma_eps_unique(t,:).^2)]');
    S_Y_DP = eye(n_y);
    yt = S_Y_DP*dp_disagg_60_end(t,:)';
    Rt = zeros(size(yt,1),size(yt,1));
    Ht = H_t*S_Y_DP';
    yt = dp_disagg_60_end(t,:)';
    Rt = zeros(size(yt,1),size(yt,1));
    Ht = H_t;
    [X1,P1,X2,P2,e,ht,llf] = kfilt_llf(yt,X1,P1,Ht,F,Rt,Q_t);
    % Parameters mapping state into tau-headline 
    S = zeros(ns,1);
    S(2) = share_h(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_h(t,:)';
    tauh = S'*X1;
    varh = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xe(t,:)';
    tauxe = S'*X1;
    varxe = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xfe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xfe(t,:)';
    tauxfe = S'*X1;
    varxfe = S'*P1*S;
    
    tau_17c_mean(t,:) = [tauh tauxe tauxfe];
    tau_17c_var(t,:) = [varh varxe varxfe];
  
  end;
  
  tau_3c_mean = zeros(dnobs_60_end,3);
  tau_3c_var = zeros(dnobs_60_end,3);
  X1 = X1_init;
  P1 = P1_init;
  for t = 1:dnobs_60_end;
    H_t = H;
    H_t(1,:) = alpha_eps(t,:);
    H_t(2,:) = alpha_tau(t,:);
    Q_t = diag([(sigma_eps_common(t,1).^2) (sigma_dtau_common(t,1).^2) (sigma_dtau_unique(t,:).^2) (sigma_eps_unique(t,:).^2)]');
    % Mapping from dp_disagg to observable
    S_Y_DP = zeros(3,n_y);
    S_Y_DP(1,icore) = share_h(t,icore);
    S_Y_DP(2,ifood) = share_h(t,ifood);
    S_Y_DP(3,ienergy) = share_h(t,ienergy);
    yt = S_Y_DP*dp_disagg_60_end(t,:)';
    Rt = zeros(size(yt,1),size(yt,1));
    Ht = H_t*S_Y_DP';
    [X1,P1,X2,P2,e,ht,llf] = kfilt_llf(yt,X1,P1,Ht,F,Rt,Q_t);
    % Parameters mapping state into tau-headline 
    S = zeros(ns,1);
    S(2) = share_h(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_h(t,:)';
    tauh = S'*X1;
    varh = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xe(t,:)';
    tauxe = S'*X1;
    varxe = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xfe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xfe(t,:)';
    tauxfe = S'*X1;
    varxfe = S'*P1*S;
    
    tau_3c_mean(t,:) = [tauh tauxe tauxfe];
    tau_3c_var(t,:) = [varh varxe varxfe];
  
  end;
  
  tau_h_mean = zeros(dnobs_60_end,3);
  tau_h_var = zeros(dnobs_60_end,3);
  X1 = X1_init;
  P1 = P1_init;
  for t = 1:dnobs_60_end;
    H_t = H;
    H_t(1,:) = alpha_eps(t,:);
    H_t(2,:) = alpha_tau(t,:);
    Q_t = diag([(sigma_eps_common(t,1).^2) (sigma_dtau_common(t,1).^2) (sigma_dtau_unique(t,:).^2) (sigma_eps_unique(t,:).^2)]');
    % Mapping from dp_disagg to observable
    S_Y_DP = zeros(1,n_y);
    S_Y_DP(1,:) = share_h(t,:);
    yt = S_Y_DP*dp_disagg_60_end(t,:)';
    Rt = zeros(size(yt,1),size(yt,1));
    Ht = H_t*S_Y_DP';
    [X1,P1,X2,P2,e,ht,llf] = kfilt_llf(yt,X1,P1,Ht,F,Rt,Q_t);
    % Parameters mapping state into tau-headline 
    S = zeros(ns,1);
    S(2) = share_h(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_h(t,:)';
    tauh = S'*X1;
    varh = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xe(t,:)';
    tauxe = S'*X1;
    varxe = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xfe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xfe(t,:)';
    tauxfe = S'*X1;
    varxfe = S'*P1*S;
    
    tau_h_mean(t,:) = [tauh NaN NaN];
    tau_h_var(t,:) = [varh NaN NaN];
  
  end;
  
  tau_xe_mean = zeros(dnobs_60_end,3);
  tau_xe_var = zeros(dnobs_60_end,3);
  X1 = X1_init;
  P1 = P1_init;
  for t = 1:dnobs_60_end;
    H_t = H;
    H_t(1,:) = alpha_eps(t,:);
    H_t(2,:) = alpha_tau(t,:);
    Q_t = diag([(sigma_eps_common(t,1).^2) (sigma_dtau_common(t,1).^2) (sigma_dtau_unique(t,:).^2) (sigma_eps_unique(t,:).^2)]');
    % Mapping from dp_disagg to observable
    S_Y_DP = zeros(1,n_y);
    S_Y_DP(1,:) = share_xe(t,:);
    yt = S_Y_DP*dp_disagg_60_end(t,:)';
    Rt = zeros(size(yt,1),size(yt,1));
    Ht = H_t*S_Y_DP';
    [X1,P1,X2,P2,e,ht,llf] = kfilt_llf(yt,X1,P1,Ht,F,Rt,Q_t);
    % Parameters mapping state into tau-headline 
    S = zeros(ns,1);
    S(2) = share_h(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_h(t,:)';
    tauh = S'*X1;
    varh = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xe(t,:)';
    tauxe = S'*X1;
    varxe = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xfe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xfe(t,:)';
    tauxfe = S'*X1;
    varxfe = S'*P1*S;
    
    tau_xe_mean(t,:) = [NaN tauxe NaN];
    tau_xe_var(t,:) = [NaN varxe NaN];
  
  end;
  
  tau_xfe_mean = zeros(dnobs_60_end,3);
  tau_xfe_var = zeros(dnobs_60_end,3);
  X1 = X1_init;
  P1 = P1_init;
  for t = 1:dnobs_60_end;
    H_t = H;
    H_t(1,:) = alpha_eps(t,:);
    H_t(2,:) = alpha_tau(t,:);
    Q_t = diag([(sigma_eps_common(t,1).^2) (sigma_dtau_common(t,1).^2) (sigma_dtau_unique(t,:).^2) (sigma_eps_unique(t,:).^2)]');
    % Mapping from dp_disagg to observable
    S_Y_DP = zeros(1,n_y);
    S_Y_DP(1,:) = share_xfe(t,:);
    yt = S_Y_DP*dp_disagg_60_end(t,:)';
    Rt = zeros(size(yt,1),size(yt,1));
    Ht = H_t*S_Y_DP';
    [X1,P1,X2,P2,e,ht,llf] = kfilt_llf(yt,X1,P1,Ht,F,Rt,Q_t);
    % Parameters mapping state into tau-headline 
    S = zeros(ns,1);
    S(2) = share_h(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_h(t,:)';
    tauh = S'*X1;
    varh = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xe(t,:)';
    tauxe = S'*X1;
    varxe = S'*P1*S;
    
    S = zeros(ns,1);
    S(2) = share_xfe(t,:)*alpha_tau(t,:)';
    S(3:2+n_y) = share_xfe(t,:)';
    tauxfe = S'*X1;
    varxfe = S'*P1*S;
    
    tau_xfe_mean(t,:) = [NaN NaN tauxfe];
    tau_xfe_var(t,:) = [NaN NaN varxfe];
  
  end;
  
  h_mean = NaN*zeros(dnobs_60_end,3);
  h_var = NaN*zeros(dnobs_60_end,3);
  xe_mean = NaN*zeros(dnobs_60_end,3);
  xe_var = NaN*zeros(dnobs_60_end,3);
  xfe_mean = NaN*zeros(dnobs_60_end,3);
  xfe_var = NaN*zeros(dnobs_60_end,3);
  for t = 1:dnobs_60_end;
    S_Y_DP = zeros(1,n_y);
    S_Y_DP = share_h(t,:);
    h_mean(t,1) = S_Y_DP*dp_disagg_60_end(t,:)';
    tmp1 = ((S_Y_DP*alpha_eps(t,:)')*sigma_eps_common(t,1))^2;
    tmp2 = sum((S_Y_DP.*sigma_eps_unique(t,:)).^2,2);
    h_var(t,1) = tmp1 + tmp2;
    
    S_Y_DP = zeros(1,n_y);
    S_Y_DP = share_xe(t,:);
    xe_mean(t,2) = S_Y_DP*dp_disagg_60_end(t,:)';
    tmp1 = ((S_Y_DP*alpha_eps(t,:)')*sigma_eps_common(t,1))^2;
    tmp2 = sum((S_Y_DP.*sigma_eps_unique(t,:)).^2,2);
    xe_var(t,2) = tmp1 + tmp2;
    
    S_Y_DP = zeros(1,n_y);
    S_Y_DP = share_xfe(t,:);
    xfe_mean(t,3) = S_Y_DP*dp_disagg_60_end(t,:)';
    tmp1 = ((S_Y_DP*alpha_eps(t,:)')*sigma_eps_common(t,1))^2;
    tmp2 = sum((S_Y_DP.*sigma_eps_unique(t,:)).^2,2);
    xfe_var(t,3) = tmp1 + tmp2;
    
  end;
  
  
  ao_h_mean = NaN*zeros(dnobs_60_end,3);
  ao_h_var = NaN*zeros(dnobs_60_end,3);
  ao_xe_mean = NaN*zeros(dnobs_60_end,3);
  ao_xe_var = NaN*zeros(dnobs_60_end,3);
  ao_xfe_mean = NaN*zeros(dnobs_60_end,3);
  ao_xfe_var = NaN*zeros(dnobs_60_end,3);
  for t = 1:dnobs_60_end;
    % Use approximation, locally constant parameter values 
    S_Y_DP = zeros(1,n_y);
    S_Y_DP = share_h(t,:);
    ao_h_mean(t,1) = S_Y_DP*dp4_disagg_60_end(t,:)';
    tmp1 = ((S_Y_DP*alpha_eps(t,:)')*sigma_eps_common(t,1))^2;
    tmp2 = sum((S_Y_DP.*sigma_eps_unique(t,:)).^2,2);
    var_eps_comp = (tmp1 + tmp2)/4;
    tmp1 = ((S_Y_DP*alpha_tau(t,:)')*sigma_dtau_common(t,1))^2;
    tmp2 = sum((S_Y_DP.*sigma_dtau_unique(t,:)).^2,2);
    var_tau_comp = (tmp1+tmp2)*14/16;
    ao_h_var(t,1) = var_eps_comp+var_tau_comp; 
    
    S_Y_DP = zeros(1,n_y);
    S_Y_DP = share_xe(t,:);
    ao_xe_mean(t,2) = S_Y_DP*dp4_disagg_60_end(t,:)';
    tmp1 = ((S_Y_DP*alpha_eps(t,:)')*sigma_eps_common(t,1))^2;
    tmp2 = sum((S_Y_DP.*sigma_eps_unique(t,:)).^2,2);
    var_eps_comp = (tmp1 + tmp2)/4;
    tmp1 = ((S_Y_DP*alpha_tau(t,:)')*sigma_dtau_common(t,1))^2;
    tmp2 = sum((S_Y_DP.*sigma_dtau_unique(t,:)).^2,2);
    var_tau_comp = (tmp1+tmp2)*14/16;
    ao_xe_var(t,2) = var_eps_comp+var_tau_comp; 
    
    S_Y_DP = zeros(1,n_y);
    S_Y_DP = share_xfe(t,:);
    ao_xfe_mean(t,3) = S_Y_DP*dp4_disagg_60_end(t,:)';
    tmp1 = ((S_Y_DP*alpha_eps(t,:)')*sigma_eps_common(t,1))^2;
    tmp2 = sum((S_Y_DP.*sigma_eps_unique(t,:)).^2,2);
    var_eps_comp = (tmp1 + tmp2)/4;
    tmp1 = ((S_Y_DP*alpha_tau(t,:)')*sigma_dtau_common(t,1))^2;
    tmp2 = sum((S_Y_DP.*sigma_dtau_unique(t,:)).^2,2);
    var_tau_comp = (tmp1+tmp2)*14/16;
    ao_xfe_var(t,3) = var_eps_comp+var_tau_comp; 
    
  end;
  
  
  i_fcst = smpl(cal_60_end,[1965 1],[2015 1],4);
  tmp = cal_60_end(i_fcst==1);
  fprintf(fileID,'Results for Sample Period: %6.2f - %6.2f \n',[tmp(1) tmp(end)]);
  fprintf(fileID,'MUCSV_17c,');
  tmp = sqrt(mean(tau_17c_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'MUCSV_3c,');
  tmp = sqrt(mean(tau_3c_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'UCSV_All,');
  tmp = sqrt(mean(tau_h_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'UCSV_xE,');
  tmp = sqrt(mean(tau_xe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'UCSV_xFE,');
  tmp = sqrt(mean(tau_xfe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'PCE_All,');
  tmp = sqrt(mean(h_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'PCE_xE,');
  tmp = sqrt(mean(xe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'PCE_xFE,');
  tmp = sqrt(mean(xfe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'4QMA_PCE_All,');
  tmp = sqrt(mean(ao_h_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'4QMA_PCE_xE,');
  tmp = sqrt(mean(ao_xe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'4QMA_PCE_xFE,');
  tmp = sqrt(mean(ao_xfe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');  
  
 
  i_fcst = smpl(cal_60_end,[1965 1],[1989 4 1],4);
  tmp = cal_60_end(i_fcst==1);
  fprintf(fileID,'\n \n Results for Sample Period: %6.2f - %6.2f \n',[tmp(1) tmp(end)]);
  fprintf(fileID,'MUCSV_17c,');
  tmp = sqrt(mean(tau_17c_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'MUCSV_3c,');
  tmp = sqrt(mean(tau_3c_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'UCSV_All,');
  tmp = sqrt(mean(tau_h_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'UCSV_xE,');
  tmp = sqrt(mean(tau_xe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'UCSV_xFE,');
  tmp = sqrt(mean(tau_xfe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'PCE_All,');
  tmp = sqrt(mean(h_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'PCE_xE,');
  tmp = sqrt(mean(xe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'PCE_xFE,');
  tmp = sqrt(mean(xfe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'4QMA_PCE_All,');
  tmp = sqrt(mean(ao_h_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'4QMA_PCE_xE,');
  tmp = sqrt(mean(ao_xe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'4QMA_PCE_xFE,');
  tmp = sqrt(mean(ao_xfe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');  
  
  
  i_fcst = smpl(cal_60_end,[1990 1],[2015 1 1],4);
  tmp = cal_60_end(i_fcst==1);
  fprintf(fileID,'\n \n Results for Sample Period: %6.2f - %6.2f \n',[tmp(1) tmp(end)]);
  fprintf(fileID,'MUCSV_17c,');
  tmp = sqrt(mean(tau_17c_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'MUCSV_3c,');
  tmp = sqrt(mean(tau_3c_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'UCSV_All,');
  tmp = sqrt(mean(tau_h_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'UCSV_xE,');
  tmp = sqrt(mean(tau_xe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'UCSV_xFE,');
  tmp = sqrt(mean(tau_xfe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'PCE_All,');
  tmp = sqrt(mean(h_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'PCE_xE,');
  tmp = sqrt(mean(xe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'PCE_xFE,');
  tmp = sqrt(mean(xfe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'4QMA_PCE_All,');
  tmp = sqrt(mean(ao_h_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'4QMA_PCE_xE,');
  tmp = sqrt(mean(ao_xe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');
  fprintf(fileID,'4QMA_PCE_xFE,');
  tmp = sqrt(mean(ao_xfe_var(i_fcst==1,:),1));
  prtmat_comma(tmp,fileID,'%5.2f ','\n');  
  