beta_estimation.m

Implements Maximum likelihood estimation of beta and other parameters for model of stock portfolio v

%Implements Maximum likelihood estimation of beta and other parameters
%for model of stock portfolio vs. index as described in page 5 of paper
%Estimating Value at Risk with the Kalman Filter
%http://www.icer.it/workshop/Berardi_Corradin_Sommacampagna.pdf

clear all;
global Rmt; %observable returns of market index
global betait; %unobservable beta process
global Rit; %observable stock returns process
global xt_last; %variable to store final beta value
Rmt=ones(1000,1);
function [betait,Rit]=SimulateData
	global Rmt;
	n=size(Rmt,1);
	betait=zeros(n,1);
	Rit=zeros(n,1);
	
	disp('\nActual parameters:');
	%AR coefficient of beta process
	T=.7 
	
	%initial starting value of beta
	beta0=.5 
	
	%covariance of unobservable process
	Q=0.1 
	
	%constant of Rit process
	alphai=0.2 
	
	%covariance of observable process
	H=0.1 
	
	shock_beta=normal_rnd(0,Q,n,1);
	shock_Rit=normal_rnd(0,H,n,1);
	
	betait(1)=beta0;
	Rit(1)=alphai+betait(1)*Rmt(1)+shock_Rit(1);
	for i=2:n,
		betait(i)=T*betait(i-1)+shock_beta(i);
		Rit(i)=alphai+betait(i)*Rmt(i)+shock_Rit(i);
	end
endfunction

[betait,Rit]=SimulateData;

%%%%% Evaluate process parameters using Maximum Likelihood Estimation

%Kalman filter equations are described in 
function ret = loglikelihoodfn(p)
	T=p(1);
	beta0=p(2);
	Q=p(3)^2;
	alpha=p(4);
	H=p(5)^2;

	global Rit;
	global Rmt;
	global xt_last;
	yt=Rit;
	nsamples=size(Rit,1);
	n=2;
	c=0;
	tmpsum=-0.5*nsamples*log(2*pi);
	xt=beta0;
	Pt=0.05;
	d=alpha;
	
	for i=2:nsamples,
		Zt=Rmt(i);
		vt=yt(i)-Zt*xt-d; %eqn (11)
		Ft=Zt*Pt*Zt'+H; %eqn(12)
		xt=T*xt+T*Pt*Zt'*vt/Ft;  %eqn 9 and 13                  
		Pt=T*Pt*T'-T*Pt*Zt'*Zt*Pt*T/Ft+Q; %eq 10 & 14
		tmpsum=tmpsum-0.5*log(Ft)-0.5*vt'*vt/Ft; %based on eqn(8)
	end
	xt_last=xt;
	ret=-tmpsum;
endfunction

%Supply initial guess values
T=.1;
beta0=0.2;
Q=0.01;
alphai=0.02;
H=0.01;
	
p=zeros(5,1);
p(1)=T;
p(2)=beta0;
p(3)=Q^0.5;
p(4)=alphai;
p(5)=H^0.5;
pnew=fmins('loglikelihoodfn',p);


disp('\nCalibrated parameters:');
T=pnew(1)
beta0=pnew(2)
Q=pnew(3)^2
alphai=pnew(4)
H=pnew(5)^2
disp('\n');
beta_projected=xt_last
n=size(Rit,1);
beta_actual=betait(n)
beta_ols=ols(Rit,Rmt)