Валидация целесообразности использования более продвинутых моделей описания процентного риска на развивающихся рынках

29.11.2022

PLBPS0000040

CAPITAL PARK SA 2014 F/R 06/17 B

13.06.2017

PLCPPRK00052

ECHO INVESTMENT SA 2014 F/R 06/16 A

26.06.2016

PLECHPS00142

ECHO INVESTMENT SA 2014 F/R 07/16 B

02.07.2016

PLECHPS00167

GETIN NOBLE BANK SA 2012 F/R 09/19 PP2-II

26.09.2019

PLGETBK00111

GHELAMCO INVEST SP 2014 F/R 01/18 PA

29.01.2018

PLGHLMC00099

KREDYT INKASO SA 2012 F/R 09/16 U03

05.09.2016

PLKRINK00121

AOW FAKTORING SA 2014 F/R 05/17 E

28.05.2017

PLAOWFK00050

CASUS FINANSE SA 2013 F/R 12/16 C

28.12.2016

PLCSSFN00026

FAST FINANCE SA 2013 F/R 11/16 S-G

15.11.2016

PLFSTFC00061

GETIN NOBLE BANK SA 2013 F/R 04/20 PP2-IX

26.04.2020

PLGETBK00178

GETIN NOBLE BANK SA 2013 F/R 10/20 PP3-V

07.10.2020

PLGETBK00210

KRUK SA 2012 CMX.CPN. 12/16 O2

06.12.2016

PLKRK0000200

KRUK SA 2013 F/R 08/18 R3

13.08.2018

PLKRK0000275

PCC AUTOCHEM SP 2014 6.8% 05/02/17 S-A

05.02.2017

PLPCCTH00011

BANCO PENTA SA 2013 3.6% 25/03/18 INDXLK.

25.03.2018

CL0001976266

2C PARTNERS SA 2013 10% 31/03/16 D

31.03.2016

PL2CPRT00048

PODKARPACKI BANK 2013 F/R 05/23 PBS0523

28.05.2023

PLPBS0000032

ROBYG SA 2013 F/R 06/16 I

23.06.2016

PLROBYG00123

ROBYG SA 2014 F/R 02/18 K2

26.02.2018

PLROBYG00164

OBRGS IT ZHILSHCHNY 2012 11% 07/02/17 02

07.02.2017

RU000A0JS3Q8

ABS INVESTMENT SA 2014 8 1/2% 24/07/17 S-A

24.07.2017

PLABSIN00020

GETIN NOBLE BANK SA 2014 F/R 02/21 PP3-IX

15.02.2021

PLGETBK00251

PPHU MEXPOL SP 2014 F/R 03/16 A

31.03.2016

PLMXPOL00019

WIERZYCIEL SA 2014 F/R 05/16 A

13.05.2016

PLWRZCL00048

SPETSSTROYFINANS 2010 8 1/2% 30/05/16 01

30.05.2016

RU000A0JQWE7

BANK ZACHODNI WBK 2013 F/R 12/16 A

19.12.2016

PLBZ00000150

FAST FINANCE SA 2013 F/R 09/16 F

27.09.2016

PLFSTFC00053

GETIN NOBLE BANK SA 2013 F/R 11/20 PP3-VI

09.11.2020

PLGETBK00228

KUJAWSKI BK.SPLDY.W 2012 F/R 05/18 KBS0518

30.05.2018

PLKBSAK00014

KRUK SA 2013 F/R 05/17 P4

20.05.2017

PLKRK0000242

M.W. TRADE SA 2013 F/R 04/16 G2013

26.04.2016

PLMWTRD00112

UNIDEVELOPMENT SA 2014 F/R 03/17 A

07.03.2017

PLUNDVL00014

VTB LEASING FINANS 2009 G/R 30/11/16 07

30.11.2016

RU000A0JQMB4

GLAVNAYA DOROGA 2010 G/R 30/10/28 03

30.10.2028

RU000A0JR4U9

BANCO DE CHILE 2011 3% 01/05/17 INDXLK.

01.05.2017

CL0001790253

BANCO VE POR MAS 2010 F/R 08/20 10

06.08.2020

MX0QBC070001

AOW FAKTORING SA 2013 F/R 04/16 C

26.04.2016

PLAOWFK00035

MBANK SA 2013 F/R 12/23 BREO201

20.12.2023

PLBRE0005177

FERRATUM CAP.PLD. 2014 F/R 05/17 B2

23.05.2017

PLFRRTM00067

GETIN NOBLE BANK SA 2013 F/R 08/20 PP3-III

28.08.2020

PLGETBK00202

GETIN NOBLE BANK SA 2014 F/R 03/21 PP3-X

12.03.2021

PLGETBK00269

GETIN NOBLE BANK SA 2014 F/R 04/21 PP3-XI

07.04.2021

PLGETBK00277

GETIN NOBLE BANK SA 2014 F/R 06/17 GNB1400

30.06.2017

PLGETBK00285

GHELAMCO INVEST SP 2014 6.27% 09/07/18 I

09.07.2018

PLGHLMC00123

GHELAMCO INVEST SP 2014 F/R 11/17 PD

27.11.2017

PLGHLMC00164

KRUK SA 2013 F/R 12/17 T1

06.12.2017

PLKRK0000317

MCI MANAGEMENT SA 2013 F/R 12/16 H3

19.12.2016

PLMCIMG00160

M.W. TRADE SA 2014 F/R 04/17 A2014

17.04.2017

PLMWTRD00138

BBI DEVELOPMENT SA 2014 F/R 02/17 BBI0217

07.02.2017

PLNFI1200158

GETIN NOBLE BANK SA 2011 F/R 09/17 D

20.09.2017

PLNOBLE00066

POLSKI BK.SPLDY. 2012 F/R 09/27 BSC0927

14.09.2027

PLPBSC000032

POLNORD SA 2014 F/R 02/17 NS1

11.02.2017

PLPOLND00134

WARSZAWSKI BANK 2011 F/R 05/21 WBS0521

27.05.2021

PLWAWBS00014

WLODARZEWSKA SA 2013 F/R 05/16 C

15.05.2016

PLWDRZW00043

VTB LEASING FINANS 2010 G/R 01/08/17 08

01.08.2017

RU000A0JQZ67

VTB LEASING FINANS 2010 G/R 02/08/17 09

02.08.2017

RU000A0JQZ75

BANK URALSIB PAO 2012 G/R 25/08/17 05

25.08.2017

RU000A0JSW84

BANCO SECURITY SA 2009 3% 01/07/19 INDXLK.

01.07.2019

CL0001516567

CORPBANCA 2010 3% 01/07/16 INDXLK.

01.07.2016

CL0001745026

BANCO DE CHILE 2011 3.4% 01/09/21 INDXLK.

01.09.2021

CL0001841130

CORPBANCA 2010 3% 01/07/17 INDXLK.

01.07.2017

CL0001847715

DIAKHITEL KOZPONT 2013 5 1/2% 22/12/16

22.12.2016

HU0000352448

FHB FOLDHITEL 2011 8 1/2% 24/03/17

24.03.2017

HU0000652508

FHB FOLDHITEL 2012 8 3/4% 14/12/17 FJ17NF0

14.12.2017

HU0000652730

BANCO MERCANTIL 2012 F/R 05/22 12

27.05.2022

MX0QBA070060

BBVA BANCO CONT. 2010 7.1875% 13/08/20 4 A

13.08.2020

PEP11600M145

AOW FAKTORING SA 2014 F/R 03/17 D

28.03.2017

PLAOWFK00043

BEST SA 2012 F/R 05/16 G

28.05.2016

PLBEST000101

BEST SA 2014 F/R 04/18 K1

30.04.2018

PLBEST000127

BANK PLSK.SPLDZ. 2010 7.27% 12/07/20 BPS0720

12.07.2020

PLBPS0000024

BANK PLSK.SPLDZ. 2011 F/R 07/18 BPS0718

15.07.2018

PLBPS0000032

BANK SPLDY.W BLJ. 2010 F/R 04/25 BRB0425

09.04.2025

PLBSBR000023

BK.SPLDY.W OOW.MZJ. 2011 F/R 02/21 BOM0221

11.02.2021

PLBSPOM00010

BK.SPLDY.W SKIERWH. 2011 F/R 04/18 BSS0418

01.04.2018

PLBSSKR00011

COPERNICUS SECS. 2014 F/R 04/16 S-A

10.04.2016

PLCPRSC00020

CZERWONA TOREBKA SA 2013 F/R 04/16 S-A

15.04.2016

PLCRWTR00048

ECHO INVESTMENT SA 2014 F/R 03/18 C

04.03.2018

PLECHPS00175

ED INVEST SA 2014 F/R 11/16 A

16.11.2016

PLEDINV00048

EUROPEJSKI FUNDUSZ 2013 F/R 12/16 E

20.12.2016

PLERPFM00030

GETIN NOBLE BANK SA 2012 F/R 05/18 PP-V

23.05.2018

PLGETBK00020

GETIN NOBLE BANK SA 2012 F/R 10/19 PP2-III

19.10.2019

PLGETBK00129

GETIN NOBLE BANK SA 2013 F/R 02/20 PP2-VII

28.02.2020

PLGETBK00152

GETIN NOBLE BANK SA 2013 F/R 03/20 PP2-VII

30.03.2020

PLGETBK00160

GETIN NOBLE BANK SA 2013 4.77% 10/07/20 PP3-II

10.07.2020

PLGETBK00194

GHELAMCO INVEST SP 2014 F/R 07/18 J

04.07.2018

PLGHLMC00131

GHELAMCO INVEST SP 2014 6.27% 09/07/18 K

09.07.2018

PLGHLMC00149

KRAKOWSKI BK.SPLDY. 2010 F/R 10/20 A

29.10.2020

PLKRBSP00017

KRUK SA 2013 F/R 06/18 R2

25.06.2018

PLKRK0000267

GETIN NOBLE BANK SA 2011 F/R 10/17 E

17.10.2017

PLNOBLE00074

GETIN NOBLE BANK SA 2012 F/R 02/18 PP-I PP

23.02.2018

PLNOBLE00082

GETIN NOBLE BANK SA 2012 F/R 03/18 PP-III

23.03.2018

PLNOBLE00090

OT LOGISTICS SA 2014 F/R 02/17 B

17.02.2017

PLODRTS00058

PODKARPACKI BANK 2010 F/R 07/20 PBS0720

02.07.2020

PLPBS0000016

PODKARPACKI BANK 2011 F/R 10/21 PBS1021

06.10.2021

PLPBS0000024

POLSKI BK.SPLDY. 2011 F/R 06/21 BSC0621

22.06.2021

PLPBSC000024

PRAGMA INWESTYCJE 2013 F/R 10/16 S-D

09.10.2016

PLPRINW00042

RANK PROGRESS SA 2013 F/R 06/16D

14.06.2016

PLRNKPR00048

ROBYG SA 2014 F/R 02/18 K1

26.02.2018

PLROBYG00156

RAIFFEISEN BANK SA 2013 5 1/2% 20/07/16

20.07.2016

RORFZBDBC028

VTB LEASING FINANS 2009 G/R 07/06/16 3

07.06.2016

RU000A0JPZ84

SEVERO ZAPADNAYA 2011 G/R 26/09/31 03

26.09.2031

RU000A0JRU20

ALFA BANK AO 2011 G/R 02/12/16 02

02.12.2016

RU000A0JRZB4

AKB SVYAZ-BANK PAO 2011 G/R 09/12/21 03

09.12.2021

RU000A0JS1F5

IAV 1 ZAO 2011 G/R 10/08/44 A

10.08.2044

RU000A0JS215

ROSNANO OAO 2012 8.6%(F/R) 04/19 04

10.04.2019

RU000A0JS6S7

ROSNANO OAO 2012 8.6%(F/R) 04/19 05

10.04.2019

RU000A0JS6T5

DOMASHNIYE DEN GI 2012 19%(G/R) 04/17 01

27.04.2017

RU000A0JS7J4

IPOAG.TKB 1 ZAO 2013 8% 14/10/45 A

14.10.2045

RU000A0JUDD1



ПРИЛОЖЕНИЕ 2. КОД VBA ДЛЯ ОБРАБОТКИ МНОЖЕСТВЕННЫХ ФАЙЛОВ С ВРЕМЕННЫМИ РЯДАМИ ДНЕВНЫХ КОТИРОВОК С САЙТА CBONDS.RU

Option Base 1

Sub Filling()

Dim wb As Workbook

Dim myPath As String

Dim myFile As String

Dim myExtension As String

Dim sh1 As Worksheet

Dim sh2 As Worksheet

Dim rng As Range

'Optimize Macro Speed

Application.ScreenUpdating = False

Application.EnableEvents = False

Application.Calculation = xlCalculationManual

Set sh1 = ActiveWorkbook.Sheets("Summary")

Set sh2 = ActiveWorkbook.Sheets("Time_series")

sh1.Activate

anc1 = sh1.Rows(1).Find("ISIN").Address

a = sh1.Range(anc1).Address

Do While Not IsEmpty(sh1.Range(a).Value)

a = Range(a).Offset(1, 0).Address

Loop

iLastRow = Range(a).Row() - 1

colISIN = Range(anc1).Column()

a = sh2.Range("A1").Address

Do While Not IsEmpty(sh2.Range(a).Value)

a = Range(a).Offset(1, 0).Address

Loop

rows_final = Range(a).Row - 2

a = sh2.Range("A1").Address

Do While Not IsEmpty(sh2.Range(a).Value)

a = Range(a).Offset(0, 1).Address

Loop

cols_final = Range(a).Column - 1

For k = 2 To iLastRow

fname = sh1.Cells(k, colISIN).Value

myExtension = ".xls"

myPath = "B:\Bonds\Series\"

myFile = Dir(myPath & fname & myExtension)

Set wb = Workbooks.Open(Filename:=myPath & myFile)

anc_temp = wb.Worksheets(1).Rows(1).Find("Дата торгов").Address

a = wb.Worksheets(1).Range(anc_temp).Address

Do While Not IsEmpty(wb.Worksheets(1).Range(a).Value)

a = Range(a).Offset(1, 0).Address

Loop

rownum_temp = Range(a).Row() - 2

rng1 = Range(anc_temp).Offset(1, 0).Address

rng2 = Range(anc_temp).Offset(rownum_temp, 0).Address

For Each c In wb.Worksheets(1).Range(rng1, rng2).Cells

c.Value = CDate(c.Value)

Next c

ReDim mas_temp(rownum_temp, 2)

For j = 1 To 2

For i = 1 To rownum_temp

mas_temp(i, j) = wb.Worksheets(1).Range(anc_temp).Offset(i, j - 1).Value

Next i

Next j

wb.Close SaveChanges:=False

col_temp = sh2.Rows(1).Find(fname).Column

For i = 1 To rows_final

For j = 1 To rownum_temp

If sh2.Cells(i + 1, 1).Value = mas_temp(j, 1) Then sh2.Cells(i + 1, col_temp).Value = mas_temp(j, 2)

Next j

Next i

a = 0

Next k

sh2.Activate

Set rng = sh2.Range(Cells(1, 1), Cells(rows_final, cols_final))

For i = rows_final To 1 Step (-1)

If WorksheetFunction.CountA(rng.Rows(i)) = 1 Then rng.Rows(i).Delete

Next i

MsgBox "Task Complete!"

Application.EnableEvents = True

Application.Calculation = xlCalculationAutomatic

Application.ScreenUpdating = True

End Sub



ПРИЛОЖЕНИЕ 3. ОБЩИЙ КОД ПРОЦЕДУРЫ MATLAB

%% Loading data into script

clear;clc;

tic

gen=readtable('maturity.xlsx');

gen.Maturity = datetime(gen.Maturity,'ConvertFrom','excel');

date=readtable('dt.xlsx');

date.Date=datetime(date.Date,'ConvertFrom','excel');

date=date.Date;

data=importdata('data.xlsx');

ytm=data.data;

cl=readtable('cl.xlsx');

cl.min_date=datetime(cl.min_date,'ConvertFrom','excel');

cl.max_date=datetime(cl.max_date,'ConvertFrom','excel');

%load('cir_par.mat','cir_par');

%load('vas_par.mat','vas_par');

%load('hw_par.mat','hw_par');

%load('bekk_par.mat','bekk_par');

%% General data transformation

[n, m]=size(ytm);

j=sum(data.data(end,:)>100);

par_value=zeros(j,1);

par_value_2=zeros(m-j,1);

k=0; a=0;

for i=1:m

if data.data(end,i)>100

k=k+1;

par_value(k)=i;

else a=a+1;

par_value_2(a)=i;

end

end

for j = 1:m

for i=1:n

if not(isnan(data.data(i,j)))

if data.data(i,j)>100

ytm(i,j)=-log(100/data.data(i,j))/...

years(gen.Maturity(j,1)-date(i));

else ytm(i,j)=-log(data.data(i,j)/100)/...

years(gen.Maturity(j,1)-date(i));

end

end

end

end

alpha=99;

%% GARCH-BEKK and GARCH-RM estimation

cluster=nan(m,2);

for i=1:m

for j=1:length(unique(cl.cluster));

if gen.Maturity(i)>=cl.min_date(j) && gen.Maturity(i)<=...

cl.max_date(j)

cluster(i,1)=i; % Order of ISIN

cluster(i,2)=cl.cluster(j); % Cluster of ISIN

continue

end

end

end

cl=unique(cl.cluster);

logr=nan(n-1,m);

for j=1:m

for i=1:n-1

if not(isnan(data.data(i+1,j))) && not(isnan(data.data(i,j)))

logr(i,j)=log(data.data(i+1,j)/data.data(i,j));

end

end

end

% additional check on length of data

cl_vol=zeros(length(cl),3);

for u=1:length(cl)

cl_vol(u,1)=cl(u);

cl_vol(u,2)=sum(cluster(:,2)==cl(u));

end

cl_vol(:,3)=cl_vol(:,2).^2;

rm_var=nan(sum(cl_vol(:,3)),250);

rm_2_var=nan(sum(cl_vol(:,3)),250);

bekk_var=nan(sum(cl_vol(:,3)),250);

cl_dim=zeros(length(cl),1);

for i=1:length(cl);

cl_dim(i)=cl_vol(i,2)*(cl_vol(i,2)+1)/2+2;

end

bekk_mu=zeros(250,length(cl));

backdates=zeros(260,m);

%{

bekk_par=zeros(sum(cl_dim),1);

%}

for t=1:250

h=waitbar(0,horzcat('GARCH estimation process for period t = ',...

num2str(t),' out of 250'));

for u=1:length(cl)

waitbar(u/length(cl));

cl_temp = sum(cluster(:,2)==cl(u));

cl_temp=nan(cl_temp,1);

a=0;

for i=1:m

if cluster(i,2)==cl(u)

a=a+1;

cl_temp(a)=cluster(i,1);

end

end

gens=logr(:,cl_temp);

gens2=[gens,(2:n)'];

gens=gens2;

gens=gens(~any(isnan(gens),2),:);

if t==1

for a=1:length(cl_temp)

backdates(:,cl_temp(a))=gens(end-259:end,end);

end

end

gens=gens(1:end-259+t,1:end-1);

bekk_mu_temp=mean(mean(gens));

gens=gens-bekk_mu_temp;

a=cumsum(cl_dim(1:u));

%{

if t==1

if u==1

bekk_par(1:a(end))= bekk(gens,[],1,0,1);

else

bekk_par(1+a(end-1):a(end))= bekk(gens,[],1,0,1);

end

end

%}

if u==1

[~,bekk_var_temp]=bekk_simulate(1,length(cl_temp),...

bekk_par(1:a(end)),1,0,1,'Scalar');

else

[~,bekk_var_temp]=bekk_simulate(1,length(cl_temp),...

bekk_par(1+a(end-1):a(end)),1,0,1,'Scalar');

end

if u==1

bekk_var(1:cl_vol(1,3),t)=reshape(bekk_var_temp,...

cl_vol(u,3),1);

else bekk_var(1+sum(cl_vol(1:u-1,3)):sum(cl_vol(1:u,3)),t)=...

reshape(bekk_var_temp,cl_vol(u,3),1);

end

rm_var_temp=riskmetrics(gens,0.94);

rm_var_temp=rm_var_temp(:,:,end);

rm_2_var_temp=riskmetrics2006(gens);

rm_2_var_temp=rm_2_var_temp(:,:,end);

bekk_mu(t,u)=bekk_mu_temp;

if u==1

rm_var(1:cl_vol(1,3),t)=reshape(rm_var_temp,cl_vol(u,3),1);

else rm_var(1+sum(cl_vol(1:u-1,3)):sum(cl_vol(1:u,3)),t)=...

reshape(rm_var_temp,cl_vol(u,3),1);

end

if u==1

rm_2_var(1:cl_vol(1,3),t)=reshape(rm_2_var_temp,cl_vol(u,3),1);

else rm_2_var(1+sum(cl_vol(1:u-1,3)):sum(cl_vol(1:u,3)),t)=...

reshape(rm_2_var_temp,cl_vol(u,3),1);

end

end

close(h);

end

lags=[1,5,10,20];

%------------------------------------------------------------------------%

% Portfolio creation

portfolio=nan(m,1);

for u=1:length(cl)

cl_temp = sum(cluster(:,2)==cl(u));

cl_temp=nan(cl_temp,1);

a=0;

for i=1:m

if cluster(i,2)==cl(u)

a=a+1;

cl_temp(a)=cluster(i,1);

end

end

portfolio(cl_temp)=1;

end

%-----------------------------------------------------------------------%

% GARCH VaR calculation

bekk_std=zeros(250,m);

rm_std=zeros(250,m);

rm_2_std=zeros(250,m);

garch_mu=zeros(250,m);

for t=1:250

for u=1:length(cl)

cl_temp = sum(cluster(:,2)==cl(u));

port_temp=nan(cl_temp,1);

cl_temp=nan(cl_temp,1);

a=0;

for i=1:m

if cluster(i,2)==cl(u)

a=a+1;

cl_temp(a)=cluster(i,1);

port_temp(a)=portfolio(i,1);

end

end

garch_mu(t,cl_temp)=bekk_mu(t,u);

if u==1

Ht=reshape(bekk_var(1:cl_vol(1,3),t),cl_vol(1,2),cl_vol(1,2));

else Ht=reshape(bekk_var(1+sum(cl_vol(1:u-1,3)):...

sum(cl_vol(1:u,3)),t), cl_vol(u,2),cl_vol(u,2));

end

bekk_std(t,cl_temp)=(port_temp'*Ht*port_temp)^0.5;

if u==1

Ht=reshape(rm_var(1:cl_vol(1,3),t),cl_vol(1,2),cl_vol(1,2));

else Ht=reshape(rm_var(1+sum(cl_vol(1:u-1,3)):...

sum(cl_vol(1:u,3)),t), cl_vol(u,2),cl_vol(u,2));

end

rm_std(t,cl_temp)=(port_temp'*Ht*port_temp)^0.5;

if u==1

Ht=reshape(rm_2_var(1:cl_vol(1,3),t),cl_vol(1,2),cl_vol(1,2));

else Ht=reshape(rm_2_var(1+sum(cl_vol(1:u-1,3)):...

sum(cl_vol(1:u,3)),t), cl_vol(u,2),cl_vol(u,2));

end

rm_2_std(t,cl_temp)=(port_temp'*Ht*port_temp)^0.5;

end

end

garch_back=nan(size(backdates));

for i=1:length(backdates)

for j=1:size(backdates,2)

if not(isnan(backdates(i,j))) && backdates(i,j)~=0

garch_back(i,j)=data.data(backdates(i,j),j);

end

end

end

bekk_risk=zeros(250,length(lags));

rm_risk=zeros(250,length(lags));

rm_2_risk=zeros(250,length(lags));

for i=1:length(lags)

bekk_risk_rt=garch_mu+norminv(alpha/100,0,1).*(lags(i)^0.5).*bekk_std;

rm_risk_rt=garch_mu+norminv(alpha/100,0,1).*(lags(i)^0.5).*rm_std;

rm_2_risk_rt=garch_mu+norminv(alpha/100,0,1).*(lags(i)^0.5).*rm_2_std;

for t=1:250

rm_risk(t,i)=sum(garch_back(t,:).*exp(rm_risk_rt(t,:)));

bekk_risk(t,i)=sum(garch_back(t,:).*exp(bekk_risk_rt(t,:)));

rm_2_risk(t,i)=sum(garch_back(t,:).*exp(rm_2_risk_rt(t,:)));

end

end

%-----------------------------------------------------------------------%

%% Backtesting portfolio creation

back_matrix=garch_back(2:end,:);

y=sum(back_matrix,2);

back_val=zeros(250,length(lags)*m);

for i=1:250

for j=1:length(lags)

if lags(j)<=10

[~,b]=min(y(i:i+lags(j)-1));

back_val(i,1+m*(j-1):m*j)=back_matrix(i+b-1,:);

else

for k=1:m

a=data.data(:,k);

a_2=a(~any(isnan(a),2));

a=a_2(end-249-lags(j):end);

back_val(i,k+m*(j-1))=min(a(i:i+lags(j)-1));

end

end

end

end

%% Nelson-Sigel

y_temp=zeros(m,2);

ns_risk=zeros(250,length(lags));

ns_port=zeros(250,length(lags)*m);

h=waitbar(0,horzcat('Nelson-Siegel estimation process'));

for t=1:250

waitbar(t/250);

ns_vol=zeros(m,1);

if t==1

ns_beta=zeros(n-260+t,3);

ns_tau=zeros(n-260+t,1);

for i=1:m

y=ytm(1:end-260+t,i);

ns_vol(i)=std(y(~any(isnan(y),2)));

end

for i=1:n-260+t

y_temp(:,1)=ytm(i,:)';

y_temp(:,2)=(1:m)';

y=y_temp(~any(isnan(y_temp),2),:);

x=years(gen.Maturity(y(:,2))-date(i));

y=y(:,1);

[ns_par] = nelsonfit(x,y);

ns_beta(i,:)=ns_par.beta(:)';

ns_tau(i)=ns_par.tau;

end

ns_beta_f=ns_beta(~any(abs(ns_beta)>1,2),:);

ns_tau_f=mean(ns_tau(~any(ns_tau>9.9,2),:));

a=ns_betas(ns_beta_f,lags);

end

x=years(gen.Maturity(:)-date(end-260+t));

for i=1:length(lags)

ns_ytm=nelsonfun(x,a(i,:),ns_tau_f)+ns_vol*(lags(i)^0.5)*...

norminv(alpha/100);

for j=1:length(par_value)

ns_ytm(par_value(j))=100*exp(-x(par_value(j))...

*ns_ytm(par_value(j)));

end

for j=1:length(par_value_2)

ns_ytm(par_value_2(j))=100/exp(-x(par_value_2(j))...

*ns_ytm(par_value_2(j)));

end

ns_risk(t,i)=sum(ns_ytm);

ns_port(t,1+m*(i-1):m*i)=ns_ytm;

end

end

close(h);

%% Vasicek, CIR and Hull-White models

occ=isnan(ytm);

occ=sum(occ)';

occ=[occ,(1:m)'];

occ2=occ(occ(:,1)<=130,:);

ytm_2=[ytm,(1:n)'];

matur=gen.Maturity(occ2(:,2))';

ytm_2=ytm_2(:,[occ2(:,2);end]);

y=ytm_2(~any(isnan(ytm_2),2),:);

periods=date(y(:,end));

tau=zeros(length(periods),length(matur));

for i=1:length(periods)

for j=1:length(matur)

tau(i,j)=years(matur(j)-periods(i));

end

end

y=y(:,1:end-1);

dt=years(diff(periods));

dt=[1.5/365;dt];

r0=zeros(250,m);

for i=1:250

for j=1:m

r0(i,j)=ytm(backdates(i,j),j);

end

end

cir_ytm=zeros(250,m);

cir_risk=zeros(250,length(lags));

vas_ytm=zeros(250,m);

vas_risk=zeros(250,length(lags));

hw_ytm=zeros(250,m);

hw_risk=zeros(250,length(lags));

cir_port=zeros(250,length(lags)*m);

vas_port=zeros(250,length(lags)*m);

hw_port=zeros(250,length(lags)*m);

h=waitbar(0,horzcat('Affine models estimation process'));

%{

cir_par=zeros(4,1);

vas_par=zeros(4,1);

hw_par=zeros(2,1);

%}

for t=1:250

waitbar(t/250);

y_est=y(1:end-260+t,:);

tau_est=tau(1:end-260+t,:);

dt_est=dt(1:end-260+t,:);

r0_est=r0(t,:);

x=years(gen.Maturity(:)-date(end-260+t));

%{

if t==1

cir_par=cir_estim(y_est,tau_est,dt_est);

vas_par=vas_estim(y_est,tau_est,dt_est);

[hw_par,b]=hw_estim(y_est,tau_est,dt_est);

end

%}

for i=1:length(lags)

cir_ytm(t,:)=cirpath(lags(i),cir_par,r0_est,alpha)';

vas_ytm(t,:)=vaspath(lags(i),vas_par,r0_est,alpha)';

hw_ytm(t,:)=hwpath(lags(i),hw_par,b,r0_est,alpha)';

for j=1:length(par_value)

cir_ytm(t,par_value(j))=100/exp(-x(par_value(j))...

*cir_ytm(t,par_value(j)));

hw_ytm(t,par_value(j))=100/exp(-x(par_value(j))...

*hw_ytm(t,par_value(j)));

vas_ytm(t,par_value(j))=100/exp(-x(par_value(j))...

*vas_ytm(t,par_value(j)));

end

for j=1:length(par_value_2)

vas_ytm(t,par_value_2(j))=exp(-x(par_value_2(j))...

*vas_ytm(t,par_value_2(j)))*100;

cir_ytm(t,par_value_2(j))=exp(-x(par_value_2(j))...

*cir_ytm(t,par_value_2(j)))*100;

hw_ytm(t,par_value_2(j))=exp(-x(par_value_2(j))...

*hw_ytm(t,par_value_2(j)))*100;

end

cir_risk(t,i)=sum(cir_ytm(t,:));

cir_port(t,1+m*(i-1):m*i)=cir_ytm(t,:);

vas_risk(t,i)=sum(vas_ytm(t,:));

vas_port(t,1+m*(i-1):m*i)=vas_ytm(t,:);

hw_risk(t,i)=sum(hw_ytm(t,:));

hw_port(t,1+m*(i-1):m*i)=hw_ytm(t,:);

end

end

close(h)

%% Backtesting

back_val_garch=zeros(250,length(lags));

for i=1:length(lags)

back_val_garch(:,i)=sum(back_val(:,1+(i-1)*m:i*m),2);

end

%save('cir_par.mat','cir_par');

%save('vas_par.mat','vas_par');

%save('hw_par.mat','hw_par');

%save('bekk_par.mat','bekk_par');

bekk_res=zeros(9*length(lags),3);

rm_res=zeros(9*length(lags),3);

rm_2_res=zeros(9*length(lags),3);

for i=1:length(lags)

bekk_res(1+9*(i-1):9*i,:)=VaRLR(back_val_garch(:,i),...

bekk_risk(:,i), alpha);

rm_res(1+9*(i-1):9*i,:)=VaRLR(back_val_garch(:,i),...

rm_risk(:,i), alpha);

rm_2_res(1+9*(i-1):9*i,:)=VaRLR(back_val_garch(:,i),...

rm_2_risk(:,i), alpha);

end

for i=1:length(lags)

ns_res=VaRLR(back_val(:,1+(i-1)*m:i*m), ns_port(:,1+(i-1)*m:i*m),...

alpha, portfolio);

vas_res=VaRLR(back_val(:,1+(i-1)*m:i*m), vas_port(:,1+(i-1)*m:i*m),...

alpha, portfolio);

cir_res=VaRLR(back_val(:,1+(i-1)*m:i*m), cir_port(:,1+(i-1)*m:i*m),...

alpha, portfolio);

hw_res=VaRLR(back_val(:,1+(i-1)*m:i*m), hw_port(:,1+(i-1)*m:i*m),...

alpha, portfolio);

end

time=toc;



ПРИЛОЖЕНИЕ 4. КОДЫ ФУНКЦИЙ MATLAB

function Ht = riskmetrics(data,lambda,backCast)

switch nargin

case 2

backCast = [];

case 3

% nothing

otherwise

error('2 or 3 inputs required.')

end

if ndims(data)==2

[T,k] = size(data);

temp = zeros(k,k,T);

for t=1:T

temp(:,:,t) = data(t,:)'*data(t,:);

end

data = temp;

end

if lambda<=0 || lambda>=1

error('LAMBDA must be between 0 and 1.')

end

if isempty(backCast)

endPoint = max(min(floor(log(.01)/log(lambda)),T),k);

weights = (1-lambda).*lambda.^(0:endPoint-1);

weights = weights/sum(weights);

backCast = zeros(k);

for i=1:endPoint

backCast = backCast + weights(i)*data(:,:,i);

end

end

backCast = (backCast+backCast)/2;

if min(eig(backCast))<0

error('BACKCAST must be positive semidefinite if provided.')

end

Ht = zeros(k,k,T);

Ht(:,:,1) = backCast;

for i=2:T

Ht(:,:,i) = (1-lambda)*data(:,:,i-1) + lambda * Ht(:,:,i-1);

End

function [Ht,weights] = riskmetrics2006(data,tau0,tau1,kmax,rho)

switch nargin

case 1

tau0=[];

tau1=[];

kmax=[];

rho =[];

case 2

tau1=[];

kmax=[];

rho =[];

case 3

kmax=[];

rho =[];

case 4

rho =[];

case 5

% Nothing

otherwise

error('1 to 5 inputs required.')

end

if isempty(tau0)

tau0 = 1560;

end

if isempty(tau1)

tau1 = 4;

end

if isempty(kmax)

kmax = 14;

end

if isempty(rho)

rho = sqrt(2);

end

if tau1*rho^(kmax-1)>tau0

error('The inputs must satisfy: TAU1*RHO^(KMAX-1)>TAU0')

end

if tau1<0

error('TAU1 must be positive')

end

if tau0<0

error('TAU0 must be positive')

end

if kmax<1 || floor(kmax)~=kmax

error('KMAX must be an integer (weakly) larger than 1.')

end

if rho<0

error('RHO must be positive')

end

if ndims(data)==2

[T,K] = size(data);

temp = zeros(K,K,T);

for t=1:T

temp(:,:,t) = data(t,:)'*data(t,:);

end

data = temp;

else

K = size(data,1);

T = size(data,3);

end

tauks = tau1*rho.^((1:kmax)-1);

w= 1-log(tauks)/log(tau0);

w=w/sum(w);

Ht = zeros(K,K,T);

Httilde = zeros(K,K,T);

for k=1:kmax

tauk = tauks(k);

mu = exp(-1/tauk);

% back casting

endPoint = max(min(floor(log(.01)/log(mu)),T),k);

weights = (1-mu).*mu.^(0:endPoint-1);

weights = weights/sum(weights);

backCast = zeros(K);

for i=1:endPoint

backCast = backCast + weights(i)*data(:,:,i);

end

Httilde(:,:,1) = backCast;

for t=2:T

Httilde(:,:,t) = mu*Httilde(:,:,t-1) + (1-mu)*data(:,:,t-1);

end

Ht = Ht + w(k) * Httilde;

end

if nargout>1

weights = zeros(T,T);

for k=1:kmax

tauk = tauks(k);

mu = exp(-1/tauk);

weightMatrix = [mu.^(0:T-1)' zeros(T,T-1)];

weightMatrix(T,T) = (1-mu);

for j=1:(T-2);

weightMatrix(T,T-j) = mu * weightMatrix(T,T-j+1);

end

for t=2:T-1

weightMatrix(t,2:t) = weightMatrix(T,T+(-t+2:0));

end

weights = weights + w(k) * weightMatrix;

end

end

function[par] = nelsonfit(x,y)

par.tau = fminbnd(@(tau) nelsonsse(tau),0,10);

par.beta = lsbetas(par.tau);

function[f] = nelsonsse(tau)

[b,f] = lsbetas(tau); %#ok

end

function[b,varargout] = lsbetas(tau)

i = x(:)/tau;

j = 1-exp(-i);

n = length(x);

z = [ones(n,1) j./i (j./i)+j-1];

b = (z'*z)\(z'*y(:));

e = y(:) - z*b;

varargout(1) = {e'*e};

end

end

function[y] = nelsonfun(x,beta,tau)

i = x(:)/tau;

j = 1-exp(-i);

y = beta(1) + beta(2)*j./i + beta(3)*((j./i)+j-1);

end

function [est_beta_t]=ns_betas(betas,lags)

est_beta_t=zeros(length(lags),3);

for j=1:3

for i=1:length(lags)

y=betas(1:end-lags(i),j);

x=ones(length(betas)-lags(i),2);

x(:,2)=betas(1+lags(i):end,j);

[b]=regress(y,x);

est_beta_t(i,j)=b(1)+b(2)*betas(end,1);

end

end

end

function [para]=vas_estim(Y, tau, dt)

% S(t+1) = mu + F S(t) + noise(Q)

% Y(t) = A + H S(t) + noise(R)

[nrow, ncol] = size(Y);

para0=[0.06,0.05,0.025,-0.2,0.1*rand(1,ncol).*ones(1,ncol)];

lb=[0.0001,0.0001,0.0001,-1, 0.00001*ones(1,ncol)];

ub=ones(1,length(para0));

% Minimizer

options = optimset('fmincon');

options.Display = 'iter';

options.Diagnostics = 'on';

options.Algorithm = 'interior-point';

options.MaxFunEvals=19000;

options.StepTolerance=1e-4;

options.FunctionTolerance=1e-4;

[x]=fmincon(@LLoneVASICEK,para0,[],[],[],[],lb,ub,[],options,Y,...

tau, dt, nrow, ncol);

para=x(1:4);

end

function sumll=LLoneVASICEK(para,Y, tau, dt, nrow, ncol)

% initialize the parameters for VASICEK model

theta=para(1); kappa=para(2); sigma=para(3); lambda=para(4);

sigmai=para(5:end);

R=eye(ncol);

for i=1:ncol

R(i,i)=sigmai(i)^2;

end

%dt=1/12;

initx=theta;

initV=sigma^2/(2*kappa);

% System Matrices Initialization

C=nan(nrow,1);

F=nan(nrow,1);

for i=1:nrow

C(i)=theta*(1-exp(-kappa*dt(i)));

F(i) = exp(-kappa*dt(i));

end

A = zeros(1, ncol);

H = A;

% Create A and B

for j = 1:ncol

AffineG=kappa^2*(theta-(sigma*lambda)/kappa)-sigma^2/2; % eqn a.6

AffineB=1/kappa*(1-exp(-kappa*tau(j))); % eqn a.5

AffineA=AffineG*(AffineB-tau(j))/kappa^2-(sigma^2*AffineB^2)/...

(4*kappa); % eqn a.4

A(j)=-AffineA/tau(j); % eqn b.2

H(j)=AffineB/tau(j); % eqn b.2

end

%% Kalman Filter

% Starting values

AdjS=initx;

VarS=initV;

LL=zeros(nrow,1); % log-likelihood vector initialization

for i=1:nrow

PredS=C(i)+F(i)*AdjS; % eqn c.10

Q=theta*sigma*sigma*(1-exp(-kappa*dt(i)))^2/(2*kappa)+sigma*...

sigma/kappa*(exp(-kappa*dt(i))-exp(-2*kappa*dt(i)))*AdjS; % eqn b.4

VarS=F(i)*VarS*F(i)'+Q; % eqn c.11

% Step 2

PredY=A+H*PredS; % eqn c.4

VarY=H'*VarS*H+R; % eqn c.5

% Step 3

PredError=Y(i,:)-PredY; % eqn c.6

KalmanGain=VarS*H*inv(VarY); % eqn c.8

AdjS=PredS+KalmanGain*PredError'; % eqn c.7

VarS=VarS*(1-KalmanGain*H'); % eqn c.9

% Step 5 Construct the likelihood function

DetY=det(VarY);

LL(i)=-(ncol/2)*log(2*pi)-0.5*PredError*...

inv(VarY)*PredError' -0.5*log(DetY);

end

sumll = -sum(LL);

end

function [rt]=vaspath(t,para,r0,alpha)

dt2 = 1/365;

theta=para(1); kappa=para(2); sigma=para(3); lambda=para(4);

% Short Rate Dynamics

rt=nan(length(r0),1);

for k=1:length(r0)

r = [r0(k); nan*t];

trial=zeros(10000,1);

for tr=1:10000

for i=1:t

r(i+1)=r(i)*exp(-(kappa+lambda*sigma)*dt2)+theta*(1-...

exp(-kappa*dt2))+sigma*((1-exp(-2*kappa*dt2))/...

2*kappa)^0.5*randn(1);

end

trial(tr)=r(end);

end

rt(k)=prctile(trial,alpha);

end

end

function [para] = cir_estim(Y,tau,dt)

% S(t+1) = mu + F S(t) + noise(Q)

% Y(t) = A + H S(t) + noise(R)

[nrow, ncol] = size(Y);

para0 = [0.05, 0.1, 0.1, 0.01*rand(1,ncol)];

options = optimset('fmincon');

options.Display = 'iter';

options.Diagnostics = 'on';

options.Algorithm = 'interior-point';

options.MaxFunEvals=19000;

options.StepTolerance=1e-4;

options.FunctionTolerance=1e-4;

[x] = fmincon(@loglik, para0,[],[],[],[],[0.0001,0.0001,0.0001,...

0.00001*ones(1,ncol)],ones(1,length(para0)),[],options,Y, tau, dt, nrow, ncol);

para = x(1:4);

end

function sumll = loglik(para, Y, tau, dt, nrow, ncol)

% initialize the parameter for CIR model

theta = para(1); kappa = para(2); sigma = para(3);

% sigmai = 0.005*ones(1,ncol);%volatility of measurement error

sigmai = para(4:end);

R = eye(ncol);

for i = 1:ncol

R(i,i) = sigmai(i)^2;

end

%dt = 1/365; %monthly Y

initx = theta;

initV = sigma^2*theta/(2*kappa);

% parameter setting for transition equation

mu=nan(1,nrow);

F=nan(1, nrow);

for i=1:nrow

mu(i) = theta*(1-exp(-kappa*dt(i)));

F(i) = exp(-kappa*dt(i));

end

% parameter setting for measurement equation

A = zeros(1, ncol);

H = A;

for i = 1:ncol

AffineGamma = sqrt(kappa^2+2*sigma^2);

AffineBeta = 2*(exp(AffineGamma*tau(1,i))-1)/((AffineGamma+kappa)*...

(exp(AffineGamma*tau(1,i))-1)+2*AffineGamma);

AffineAlpha = 2*kappa*theta/(sigma^2)*log(2*AffineGamma*...

exp((AffineGamma+kappa)*tau(1,i)/2)/((AffineGamma+kappa)*...

(exp(AffineGamma*tau(1,i))-1)+2*AffineGamma));

A(i) = -AffineAlpha/tau(1,i);

H(i) = AffineBeta/tau(1,i);

end

%now recursive steps

AdjS = initx;

VarS = initV;

ll = zeros(nrow,1); %log-likelihood

for i = 1:nrow

PredS = mu(i)+F(i)*AdjS; %predict values for S and Y

Q = theta*sigma*sigma*(1-exp(-kappa*dt(i)))^2/(2*kappa)+...

sigma*sigma/kappa*(exp(-kappa*dt(i))-exp(-2*kappa*dt(i)))*AdjS;

VarS = F*VarS*F'+Q;

PredY = A+H*PredS;

PredError = Y(i,:)-PredY;

VarY = H'*VarS*H+R;

DetY=det(VarY);

InvVarY = inv(VarY);

%updating

KalmanGain = VarS*H*InvVarY;

AdjS = PredS+KalmanGain*PredError';

VarS = VarS*(1-KalmanGain*H');

ll(i) = -(ncol/2)*log(2*pi)-0.5*PredError*InvVarY*PredError'-0.5*log(DetY);

end

sumll = -sum(ll);

end

function [rt] = cirpath(t,para,r0,alpha)

a=para(1); b=para(2); s=para(3);

t=(1:t)/365;

dt2 = diff(t(:));

n = length(t);

v = s^2;

d = 4*a*b/v;

e = exp(-a*dt2);

c = v.*(1-e)/(4*a);

rt=nan(length(r0),1);

for k=1:length(r0)

r = [r0(k); nan*dt2];

trial=zeros(10000,1);

for tr=1:10000

for i = 1:(n-1)

l = r(i)*e(i)/c(i);

r(i+1) = c(i)*ncx2rnd(d,l);

end

trial(tr)=r(end);

end

rt(k)=prctile(trial,alpha);

end

end

function [para, b]=hw_estim(Y, tau, dt)

% S(t+1) = mu + F S(t) + noise(Q)

% Y(t) = A + H S(t) + noise(R)

[nrow, ncol] = size(Y);

para0=[0.06,0.05,0.025,0.1*rand(1,ncol).*ones(1,ncol)];

lb=[0.0001,0.0001,0.0001, 0.00001*ones(1,ncol)];

ub=ones(1,length(para0));

options = optimset('fmincon');

options.Display = 'iter';

options.Diagnostics = 'on';

options.Algorithm = 'interior-point';

options.MaxFunEvals=19000;

options.StepTolerance=1e-4;

options.FunctionTolerance=1e-4;

% Minimizer

[x]=fmincon(@LLHW,para0,[],[],[],[],lb,ub,[],options,Y, tau, dt, nrow, ncol);

para=x(2:3);

[b]=ThetaHW(Y, para, dt);

end

function [b]=ThetaHW(Y, para, dt)

kappa=para(1); sigma=para(2);

y_theta=zeros(length(Y)-1,size(Y,2));

for i=2:length(Y)

for j=1:size(Y,2)

y_theta(i-1,j)=(Y(i,j)-Y(i-1,j).*exp(-kappa*dt(i))-sigma*((1-exp(-2*...

kappa*dt(i)))/(2*kappa))^0.5*randn())/(1-exp(-kappa*dt(i)));

end

end

b=mean(y_theta,2);

end

function sumll=LLHW(para,Y, tau, dt, nrow, ncol)

% initialize the parameters for VASICEK model

theta=para(1); kappa=para(2); sigma=para(3);

sigmai=para(4:end);

R=eye(ncol);

for i=1:ncol

R(i,i)=sigmai(i)^2;

end

%dt=1/12;

initx=theta;

initV=sigma^2/(2*kappa);

% System Matrices Initialization

C=nan(nrow,1);

F=nan(nrow,1);

for i=1:nrow

C(i)=theta*(1-exp(-kappa*dt(i)));

F(i) = exp(-kappa*dt(i));

end

A = zeros(1, ncol);

H = A;

% Create A and B

for j = 1:ncol

AffineG=kappa^2*theta-sigma^2/2; % eqn a.6

AffineB=1/kappa*(1-exp(-kappa*tau(j))); % eqn a.5

AffineA=AffineG*(AffineB-tau(j))/kappa^2-(sigma^2*AffineB^2)/...

(4*kappa); % eqn a.4

A(j)=-AffineA/tau(j); % eqn b.2

H(j)=AffineB/tau(j); % eqn b.2

end

AdjS=initx;

VarS=initV;

LL=zeros(nrow,1); % log-likelihood vector initialization

for i=1:nrow

PredS=C(i)+F(i)*AdjS; % eqn c.10

Q=theta*sigma*sigma*(1-exp(-kappa*dt(i)))^2/(2*kappa)+sigma*sigma/...

kappa*(exp(-kappa*dt(i))-exp(-2*kappa*dt(i)))*AdjS; % eqn b.4

VarS=F(i)*VarS*F(i)'+Q; % eqn c.11

% Step 2

PredY=A+H*PredS; % eqn c.4

VarY=H'*VarS*H+R; % eqn c.5

% Step 3

PredError=Y(i,:)-PredY; % eqn c.6

KalmanGain=VarS*H*inv(VarY); % eqn c.8

AdjS=PredS+KalmanGain*PredError'; % eqn c.7

VarS=VarS*(1-KalmanGain*H'); % eqn c.9

% Step 5 Construct the likelihood function

DetY=det(VarY);

LL(i)=-(ncol/2)*log(2*pi)-0.5*PredError*...

inv(VarY)*PredError'-0.5*log(DetY);

end

sumll = -sum(LL);

end

function [rt]=hwpath(t,para,b,r0,alpha)

dt2 = 1/365;

kappa=para(1); sigma=para(2);

theta=zeros(t,1);

for i=1:t

theta(i)=hw_thetas(b,i);

end

rt=nan(length(r0),1);

for k=1:length(r0)

r = [r0(k); nan(t,1)];

trial=zeros(10000,1);

for tr=1:10000

for i=1:t

r(i+1)=r(i)*exp(-kappa*dt2)+theta(i)*(1-exp(-kappa*dt2))+sigma*((1-exp(-2*kappa*dt2))/2*kappa)^0.5*randn(1);

end

trial(tr)=r(end);

end

rt(k)=prctile(trial,alpha);

end

end

function [est_beta_t]=hw_thetas(betas,lags)

est_beta_t=zeros(length(lags),1);

for i=1:length(lags)

y=betas(1:end-lags(i));

x=ones(length(betas)-lags(i),2);

x(:,2)=betas(1+lags(i):end);

[b]=regress(y,x);

est_beta_t(i)=b(1)+b(2)*betas(end,1);

end

end



ПРИЛОЖЕНИЕ 5. КОД ФУНКЦИИ ОБРАТНОГО ТЕСТИРОВАНИЯ MATLAB

function results = VaRLR(Back_data, VaR_data, alpha, position)

%{

-----------------------------------------------------------------------

USAGE:

results = VaRLR(fdata, VaR, alpha, position)

results: ( 9 x 3 ) matrix

INPUTS:

fdata: ( t x m ) vector of the out-of-sample data,

VaR: ( t x m ) vector of VaR estimates

alpha: a%

position: ( t x 1 ) volume of short/long position for each instrument

-----------------------------------------------------------------------

OUTPUTS:

results: PF: Percentage of Failures

TUFF: Time Until First Failure

LRTUFF: Likelihood Ratio of Time Until First Failure

LRUC: Likelihood Ratio Unconditional Coverage

LRIND: Likelihood Ratio Independence Coverage

LRCC: Likelihood Ratio Conditional Coverage

Basel: Basel II Accord

LR_MIX: Interval Forecast Test

LRDW: Haas Test

L: Lopez sum

%}

%% General actions

alpha=1-alpha/100;

if nargin==3

position=1;

end

Back_data=Back_data*position;

VaR_data=VaR_data*position;

indicator=Back_data<VaR_data;

T = size(indicator,1);

n1 = sum(indicator); % Number of Violations

p = n1/T;

%% Percentage of Failures

PF = p;

% Kupiec (1995) POF Likelihood Ratio Test

LRPOF =-2*log(((1-alpha)./(1-p)).^(T-n1)*(alpha/p).^n1);

%% Kupiec (1995) Time Until First Failure

tt = find(indicator(:)== 1,1); % Find the First Failure

if isempty(tt)

LRTUFF=0;

else

LRTUFF = -2*log((alpha*(1-alpha)^(tt-1)))+ 2*log((1/tt)*(1-(1/tt))^(tt-1));

end

clear tt

%% Christoffersen (2003) Tests

% Independence Coverage

n00=0;

n01=0;

n10=0;

n11=0;

for i = 1:T-1

if indicator(i) == 0 && indicator(i+1) == 0 % 0 followed by 0

n00 = n00 + 1;

elseif indicator(i) == 0 && indicator(i+1) == 1 % 0 followed by 1

n01 = n01 + 1;

elseif indicator(i) == 1 && indicator(i+1) == 0 % 1 followed by 0

n10 = n10+1;

elseif indicator(i) == 1 && indicator(i+1) == 1 % 1 followed by 1

n11 = n11+1;

end % if

end

p01 = n01/(n00+n01);

p11 = n11/(n10+n11);

p2 = (n01+n11)/(n00+n01+n10+n11);

if n11 == 0

LRIND = ((1-p01).^n00).*(p01.^n01);

else

LRIND = -2*log((((1-p2).^(n00+n10)).*(p2.^(n01+n11)))/...

(((1-p01).^n00)*(p01.^n01).*((1-p11).^n10)*(p11.^n11)));

end

% Conditional Coverage

LRCC = LRPOF + LRIND;

%% BASEL II Accord

limits = binocdf(1:50,T,0.01)';

green = length(limits(limits<0.90));

yellow = green + length(limits(limits>0.90 & limits<0.9999));

if n1 >= yellow

BASEL = -1; % Red

elseif n1 <= yellow && n1 > green

BASEL = 0; % Yellow

else

BASEL = 1; % Green

end

%% Mixed Kupiec test

tt=sum(indicator);

if tt==0

d=T;

else d=zeros(tt,1); % durations of non-violations, min=1;

end

a=0; b=1;

for i=1:T

if indicator(i)==1

a=a+1;

d(a)=b;

b=1;

else b=b+1;

end

end

LRIND_MIX=zeros(tt,1);

for i=1:tt

LRIND_MIX(i)=-2*log((alpha*(1-alpha)^(d(i)-1)))+ 2*log((1/d(i))*...

(1-1/d(i))^(d(i)-1));

end

LR_MIX=sum(LRIND_MIX);

clear i a b LRIND_MIX;

%% Haas Test

paraur = [1, 1];

if d==T

LL_UR=0;

LL_R=0;

else

[~, LL_UR] = fmincon(@HaasLL_UR, paraur,[],[],[],[],[],[],[],[],d);

LL_UR=-LL_UR;

parar=1;

[~, LL_R] = fmincon(@HaasLL_R, parar,[],[],[],[],[],[],[],[],d);

LL_R=-LL_R;

end

LRDW=-2*(LL_UR-LL_R);

%% Lopez Test

L=1+sum((VaR_data(indicator)-Back_data(indicator)).^2)/tt;

%% Adjust the results for the non-violation models

if n1 == 0 && p ==0

PF = 0;

LRTUFF = 0;

LRPOF = 0;

LRIND = 0;

LRCC = 0;

BASEL = 1;

LR_MIX = 0;

end

% Save results

results = [PF, 0, 0;...

LRTUFF, chi2inv(1-alpha,1), 1-chi2cdf(LRTUFF,1);...

LRPOF, chi2inv(1-alpha,1), 1-chi2cdf(LRPOF,1);...

LRIND, chi2inv(1-alpha,1), 1-chi2cdf(LRIND,1);...

LRCC, chi2inv(1-alpha,2), 1-chi2cdf(LRCC,2);...

BASEL, 0, 0;...

LR_MIX, chi2inv(1-alpha,tt+1), 1-chi2cdf(LR_MIX,tt+1);...

LRDW, chi2inv(1-alpha,1), 1-chi2cdf(LRDW,1);...

L, 0, 0];

end

%% Discrete Weibull estimation

function LL_UR=HaasLL_UR(d, paraur)

a=paraur(1); b=paraur(2);

tt=length(d);

LL=length(d);

for i=1:tt

LL(i)=log(exp(-a^b*(d(i)-1)^b)-exp(-a^b*d(i)^b));

end

LL_UR=-sum(LL);

end

function LL_R=HaasLL_R(d, parar)

a=parar(1);

tt=length(d);

LL=length(d);

for i=1:tt

LL(i)=log(exp(-a*(d(i)-1))-exp(-a*d(i)));

end

LL_R=-sum(LL);

end

Размещено на Allbest.ru