function PlotMargLikeModifiedHarmonic(PostSample,LogPost,NumBurnin,CurrChain,DSGEModel,CurrINI)
% PlotMargLikeModifiedHarmonic: Calculates sequentially the modified harmonic mean estimator
%                               of the log marginal likelihood.
%
% USAGE:
%
%       PlotMargLikeModifiedHarmonic(PostSample,NumBurnin,CurrChain,DSGEModel,CurrINI)
%
% REQUIRED INPUT: PostSample (matrix) of dimension NumIterxNumParam with the posterior sample draws
%                            of the phi parameters.
%
%                 LogPost (vector) of dimension (NumIterx1) with all the log posterior values for the
%                            non-discarded posterior draws. These posterior values are calculated
%                            using the "logPosteriorPhiDSGE" function.
%
%                 NumBurnin (integer), with the number of burn-in draws.
%
%                 CurrChain (integer) with the current MCMC chain number.
%
%                 DSGEModel (structure) whose fields contain all required information
%                            about the model, such as input files, variable and sample
%                            data. Actual observations are not stored in this file.
%
%                 CurrINI (structure) with initialization information.
%
%
%                       Written by: Anders Warne
%                                   New Area Wide Model Project
%                                   DG-R/EMO
%                                   European Central Bank (ECB)
%                                   Email: anders.warne@ecb.europa.eu
%                                   Copyright © 2006-2008 European Central Bank.
%
%                       First version: December 12, 2006.
%                        This version: December 19, 2008.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% LICENSE INFORMATION:
%
%      YADA is 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.
%
%      This program 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.
%
%      You should have received a copy of the GNU General Public License
%      along with this program.  If not, see <http://www.gnu.org/licenses/>.
%
%      YADA is released under the GNU General Public License, Version 3,
%      29 June 2007 <http://www.gnu.org/licenses/>. The current release of
%      the program was last modified by the ECB on the "This version" date
%      above.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% CHANGELOG:
%
% * 13-12-2006: Finished work on the function.
%
% * 14-12-2006: Made sure that the wait dialog displays a cancel button.
%
% * 27-2-2007:  Added the call "SetAxesFonts(CurrINI,legHandle)" after the legend call.
%
% * 10-4-2007:  Took the name change RWMH to RWM into account.
%
% * 18-4-2007:  Added the average line. Also, the line widths have been made smaller and the
%               markers have been removed.
%
% * 9-11-2007:  Updated the documentation.
%
% * 10-1-2008:  Changed the last drawnow call to pause(0.02).
%
% * 10-3-2008:  Added the "CurrChain" input variable.
%
% * 11-3-2008:  Added code for writing output to file.
%
% * 23-5-2008:  Updated the documentation.
%
% * 28-5-2008:  Changed name text on figure window.
%
% * 19-12-2008: Updated the documentation.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%

%
% locate the controls structure
%
maingui = findobj('Type','figure','Tag','YADA');
controls = get(maingui,'UserData');
%
% Get the base values number of draws, etc.
%
% 1. Number of posterior draws
%
DrawsStr = get(controls.posterior.draws,'String');
NumPostDraws = str2num(StringTrim(strrep(DrawsStr(DSGEModel.PosteriorDrawsValue,:),',','')));
%
% Coverage probabilities for the Modified Harmonic Mean estimator of
% the marginal likelihood
%
CovStartStr = get(controls.dsge.coveragestart,'String');
DSGEModel.CovStartValue = str2num(StringTrim(CovStartStr(DSGEModel.CoverageStart,:)));
CovIncStr = get(controls.dsge.coverageincrement,'String');
DSGEModel.CovIncValue = str2num(StringTrim(CovIncStr(DSGEModel.CoverageIncrement,:)));
CovEndStr = get(controls.dsge.coverageend,'String');
DSGEModel.CovEndValue = str2num(StringTrim(CovEndStr(DSGEModel.CoverageEnd,:)));
%
% sequential marginal likelihood start value and step length
%
SeqStartStr = get(controls.dsge.sequentialstart,'String');
DSGEModel.SequentialStartIterationValue = str2num(strrep(StringTrim(SeqStartStr(DSGEModel.SequentialStartIteration,:)),',',''));
SeqStepStr = get(controls.dsge.sequentialstep,'String');
DSGEModel.SequentialStepLengthValue = str2num(strrep(StringTrim(SeqStepStr(DSGEModel.SequentialStepLength,:)),',',''));
%
% Open up a wait dialog.
%
txt = ['Please wait while YADA computes Geweke''s modified harmonic mean estimator of the log marginal likelihood. Estimation started at: ' StringTrim(datestr(now,14)) '.'];
WaitHandle = WaitDLG(txt,'information','Sequential Marginal Likelihood',500,150,CurrINI,1);
WaitControls = get(WaitHandle,'UserData');
%
% estimate the marginal likelihood
%
[LogMargsMH,CovProb] = MargLikeModifiedHarmonic(PostSample(NumBurnin+1:NumPostDraws,:),LogPost(NumBurnin+1:NumPostDraws,:),1,[],DSGEModel,WaitControls.cancel);
%
% shut down the wait dialog
%
set(WaitControls.text,'UserData','done');
delete(WaitHandle);
drawnow;
pause(0.02);
%
if isempty(LogMargsMH)==0;
   %
   % ask if we write to file
   %
   txt = 'Would you like to save the sequential estimation data to file?';
   answer = Query(txt,'question',200,'Save - Modified Harmonic Mean',500,'yes',CurrINI);
   if strcmp(lower(answer),'yes')==1;
      %
      Xvalues = LogMargsMH(:,1)+NumBurnin;
      LogMargLikeMH = LogMargsMH(:,2:size(LogMargsMH,2));
      %
      file = [DSGEModel.OutputDirectory '\SeqModifiedHarmonic-' DSGEModel.NameOfModel '-' int2str(NumPostDraws) '-' int2str(CurrChain) '.mat'];
      save(file,'CovProb','Xvalues','LogMargLikeMH');
      %
      mfile = [DSGEModel.OutputDirectory '\SeqMH' strrep(DSGEModel.NameOfModel,'-','') int2str(NumPostDraws) int2str(CurrChain) '.m'];
      fid = fopen(mfile,'wt');
      fprintf(fid,['%%\n%% load the data in ' GetFilename(file) '\n%%\n']);
      fprintf(fid,'DataStructure = load(''%s'');\n',GetFilename(file));
      fprintf(fid,'%%\n%% Number of posterior draws used for sequential estimation (N x 1) vector\n%%\n');
      fprintf(fid,'Xvalues = DataStructure.Xvalues;\n');
      fprintf(fid,'%%\n%% coverage probabilities used (1 x p vector)\n%%\n');
      fprintf(fid,'CovProb = DataStructure.CovProb;\n');
      fprintf(fid,'%%\n%% N x p matrix with log marginal likelihood values using\n%% the modified harmonic mean estimator\n%%\n');
      fprintf(fid,'LogMargLikeMH = DataStructure.LogMargLikeMH;\n');
      fprintf(fid,'%%\n%% Add your own commands below\n%%\n\n\n');
      %
      % end the script file
      %
      fprintf(fid,'%%\n%% Created by YADA on %s\n%%\n',datestr(now,0));
      fclose(fid);
      %
      txt = ['The data on sequential estimation of the log marginal likelihood has been saved to the file "' GetFilename(file) '" in the directory "' GetPath(file) '". The file contains 3 entries: Xvalues (Nx1 vector with the number of posteror draws used for the sequential estimation), CovProb (1xp vector with coverage probabilities), and LogMargLikeMH (Nxp matrix with estimates of the log marginal likelihood with the modified harmonic mean). To access this data you may run the matlab script file "' GetFilename(mfile) '".'];
      About(txt,'information','Modified Harmonic Mean',200,500,CurrINI);
   end;
   %
   % prepare the graph
   %
   GraphDirectory = [DSGEModel.OutputDirectory '\graphs'];
   %
   FigHandle =figure('Visible','off', ...
          'Position',[(CurrINI.scrsz(3)-700)/2 (CurrINI.scrsz(4)-500)/2 700 500], ...
          'Tag','MargLikeMHM', ...
          'Name','Sequential Marginal Likelihood - Modified Harmonic Mean');
   %
   %
   SetFigureProperties(FigHandle,CurrINI);
   %
   Xvalues = LogMargsMH(:,1)+NumBurnin;
   Yvalues = LogMargsMH(:,2:size(LogMargsMH,2));
   CurrMap = colormap;
   NumColor = 1;
   NumStyle = 1;
   NumMarker = 1;
   ColorOrder = [(1:8:57) (2:8:58) (3:8:59) (4:8:60) (5:8:61) (6:8:62) (7:8:63) (8:8:64)];
   LegendStr = '';
   for i=1:size(Yvalues,2);
      if NumColor<=length(CurrMap);
         CurrColor = CurrMap(ColorOrder(i),:);
      else;
         CurrColor = CurrMap(ColorOrder(1),:);
         NumColor = 1;
      end;
      if NumStyle==1;
         CurrStyle = '-';
      elseif NumStyle==2;
         CurrStyle = '--';
      elseif NumStyle==3;
         CurrStyle = ':';
      elseif NumStyle==4;
         NumStyle = 0;
         CurrStyle = '-.';
      end;
      plot(Xvalues,Yvalues(:,i),'LineStyle',CurrStyle,'Color',CurrColor,'LineWidth',1);
      hold('on');
      NumColor = NumColor+1;
      NumStyle = NumStyle+1;
      LegendStr = strvcat(LegendStr,[num2str(CovProb(i)*100,'%0.0f') ' percent coverage prob']);
   end;
   %
   % plot the average
   %
   plot(Xvalues,(1/size(Yvalues,2))*sum(Yvalues,2),'LineStyle','-','Color',[0 0 0],'LineWidth',1);
   LegendStr = strvcat(LegendStr,'Overall mean');
   hold('off');
   SetAxesFonts(CurrINI,'Title','Modified Harmonic Mean Estimates','XLabel',['RWM draw number (with ' num2str(NumBurnin,'%0.0f') ' burn-in draws)'],'YLabel','Log marginal likelihood');
   axis('tight');
   %
   % check the YLim
   %
   ylim = get(gca,'YLim');
   add_on = 0.1*(ylim(2)-ylim(1));
   ylim = [ylim(1)-add_on ylim(2)+add_on];
   set(gca,'YLim',ylim);
   %
   legHandle = legend(LegendStr,0);
   SetAxesFonts(CurrINI,legHandle);
   %
   % save the graph
   %
   if DirExist(GraphDirectory);
      file = [GraphDirectory '\LogMarginalLikelihood-ModifiedHarmonic-' num2str(NumPostDraws,'%0.0f') '-' int2str(CurrChain) '.fig'];
      saveas(FigHandle,file);
   end;
   %
   % display the graph
   %
   SetToolbarImages(FigHandle,CurrINI.images);
   set(FigHandle,'Visible','on','CloseRequestFcn','delete(gcf); drawnow; pause(0.02);');
   drawnow;
   pause(0.02);
end;

%
% end of PlotMargLikeModifiedHarmonic.m
%
