PortVaRmc(nsim)

function var = PortVaRmc(nsim)
% portfolio VaR using Monte Carlo
% 
% nosim is the number of simulations of the possible future price
% of the 9 assets in the Equity.xls file.
%
% Select the Equity.xls file when prompted
%
% Example:
%       var = PortVaRmc(1000)
%       var =
%           7.9407e+005
%
%
% Note:  Comment out lines 51 through 58 to suppress the plotting
% of the output of the nosim x 9 lines, which shows the output of each
% simulation

[fname, pname] = uigetfile('*.xls','Select a file with equity data');
if fname==0
	return
end

fname= [pname, fname];
range='a1:j789';

% get equity data
edat = excelget(fname, 1, range);

eprice = cell2mat(edat(2:end, 2:end)); % price data
eret = log(eprice(2:end,:)./eprice(1:end-1,:));

% positions in number of shares
positions = [905   569   632   234 ...
    549   932   335   656   392];

% positions in monetary value
posvalues = positions.*eprice(end,:);

% confidence level
alpha = norminv(0.95,0,1);

time = 0:20;

% Call Monte Carlo Function
aVal = multimc(eprice, time, nsim);
%aVal time-by-sim-by-assets

% Plot output of the simulations over 9 assets
figure;
for k = 1:9
    subplot(3,3,k)
    plot(aVal(:,:,k));
    set(gca,'xlim',[0 time(end)+1])
    grid on
    title(['Asset ',num2str(k),': ',sprintf('%6.2f',mean(mean(aVal(:,:,k))))])
end

% calculate returns at the end of the MC simulations
mcRets = log(aVal(end,:,:)./aVal(1,:,:));

row = sum(~isnan(mcRets(1,:,1)));
nbins = ceil(sqrt(row));

figure;
[n,xbin]=hist(mcRets(1,:,1),nbins);
hh = bar(xbin,n,1); % Plots the histogram.  No gap between bars.
%hist(mcRets(1,:,1));
title('Distribution of Returns');
hold on

xd = get(hh,'Xdata'); % Gets the x-data of the bins.

rangex = max(xd(:)) - min(xd(:)); % Finds the range of this data.
binwidth = rangex/nbins;    % Finds the width of each bin.

mn = mean(mcRets(1,:,1));
st = std(mcRets(1,:,1));
xdata = linspace(min(xbin),max(xbin),100);% Evenly spaced samples of the expected data range.
cdf = normcdf(xdata,mn,st);
cdf = cdf*max(n);
plot(xdata,cdf)
set(gca,'xlim',[min(xbin),max(xbin)])

% standard deviation of assets returns
stdRets = std(squeeze(mcRets));
% value at risk as an absolute monetary loss
% var = conf_level * position * std_returns

var = alpha * posvalues * stdRets';
var = sum(var);