%	portfolio_tangency.m
%   Produces Figure 5.8

bmu = [.08;.03;.05] ;
bOmega = [ .3 .02 .01 ;
	  .02 .15 .03 ;
	  .01 .03 .18 ] ;

muf = .02 ;

bone = ones(length(bmu),1) ;


muP = linspace(min(bmu),max(bmu),50) ;
sigmaP = zeros(1,50) ;

ibOmega = inv(bOmega) ;
A = bone'*ibOmega*bmu ;
B = bmu'*ibOmega*bmu ;
C = bone'*ibOmega*bone ;
D = B*C - A^2 ;
bg = (B*ibOmega*bone - A*ibOmega*bmu)/D ;
bh = (C*ibOmega*bmu - A*ibOmega*bone)/D ;

for i=1:50 ;
	omegaP = bg + muP(i)*bh ;
	sigmaP(i) = sqrt(omegaP'*bOmega*omegaP) ;
end ;


%	Compute minimum expected return and minimum return SD
gg = bg'*bOmega*bg ;
hh = bh'*bOmega*bh ;
gh = bg'*bOmega*bh ;
mumin = - gh/hh ;
sdmin = sqrt( gg * ( 1 - gh^2/(gg*hh)) ) ;


bomegabar = ibOmega*(bmu - muf*bone) ;
bomegaT = bomegabar/(bone'*bomegabar) ;
sigmaT = sqrt(bomegaT'*bOmega*bomegaT) ;
muT = bmu'*bomegaT ;

fsize = 16 ;
fsize2 = 35 ;
bomegaP2 = [0;.3;.7] ;
sigmaP2 = sqrt(bomegaP2'*bOmega*bomegaP2) ;
muP2 = bmu'*bomegaP2 ;
ind = (muP > mumin) ;	%	Indicates efficient horizon
ind2 = (muP < mumin) ;	%	Indicates locus below efficient horizon

figure(1)
clf ;
p1 = plot(sigmaP(ind),muP(ind),'-',sigmaP(ind2),muP(ind2), ... 
	'--' ,sdmin,mumin,'.') ;
l1 = line([0,sigmaT],[muf,muT]) ;
t1= text(sigmaP2,muP2,'* P','fontsize',fsize2) ;
t2= text(sigmaT,muT,'* T','fontsize',fsize2) ;
t3=text(0,muf+.006,'F','fontsize',fsize2) ;
t3B= text(0,muf,'*','fontsize',fsize2) ;

set(p1(1:2),'linewidth',4) ;
set(p1(3),'markersize',40) ;
set(p1(3),'color','red') ;
set(l1,'linewidth',4) ;
set(l1,'linestyle','--') ;


xlabel('standard deviation of return','fontsize',fsize) ;
ylabel('expected return','fontsize',fsize) ;

orient portrait ;
print portfolio04.ps -deps ;
