computeParametricVaR(returns,confidence_level,plot_flag)

function VaR = computeParametricVaR(returns,confidence_level,plot_flag)
% Inputs:
% returns               Vector of daily returns
% confidence_level      Confidence level (default 0.95)
% plot_flag             if true, visualize result (default is true)
%
% Outputs:
% VaR                   Value at Risk

% handle inputs
if nargin < 3
    plot_flag = true;
end
if nargin < 2
    confidence_level = 0.95;
end

% Calculate mean and standard deviation of returns
[mu,sigma] = normfit(returns);

% Calculate VaR Estimate using Normal Distribution Fit
VaR = sigma*norminv(1-confidence_level) + mu;

% Plot results if requested
if plot_flag
    % Construct domain of probabilities to graph distribution.
    % 100 equally spaced points between min and max return
    x_min = min(returns);
    x_max = max(returns);
    x_full = linspace(x_min,x_max,100);
    x_partial = x_full(x_full < VaR);
    y_full = normpdf(x_full,mu,sigma);
    y_partial = normpdf(x_partial,mu,sigma);
    area(x_full,y_full,'FaceColor','b');
    hold on;
    area(x_partial,y_partial,'FaceColor','r');
    grid on;
    % Histogram data
    [count,bins] = hist(returns,30);
    % Scale bins
    num_returns = numel(returns);
    scale = (bins(2)-bins(1))*num_returns;
    % Plot full data set
    a = bar(bins,count/scale,'w');
    set(get(a,'Children'),'FaceAlpha',0.8)
    hold off;
    title(['Best Normal Fit of Returns. Red Indicates Returns Below VaR: ',num2str(VaR)],'FontWeight','bold');
end