求助，关于GARCH-BEEK模型的MATLAB编程

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

谁有相关程序哈,小弟编不出来.希望那位大狭能帮助下,油箱 sugang880@163.com,万分感谢。

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享1 收藏3 回帖

关键词：MATLAB编程 MATLAB GARCH matla atlab MATLAB 模型

function [parameters, loglikelihood, Ht, likelihoods, stdresid, stderrors, A, B, scores]  = full_bekk_mvgarch(data,p,q, BEKKoptions);
% PURPOSE:
%      To Estimate a full BEKK multivariate GARCH model.  
%
% USAGE:
%      [parameters, loglikelihood, Ht, likelihoods, stdresid, stderrors, A, B, scores]  = full_bekk_mvgarch(data,p,q,options);
%
%
% INPUTS:
%      data          - A t by k matrix of zero mean residuals
%      p             - The lag length of the innovation process
%      q             - The lag length of the AR process
%      options       - (optional) Options for the optimization(fminunc)
%
% OUTPUTS:
%      parameters    - A (k*(k+1))/2+p*k^2+q*k^2 vector of estimated parameteters. F
%                         or any k^2 set of Innovation or AR parameters X,
%                         reshape(X,k,k) will give the correct matrix
%                         To recover C, use ivech(parmaeters(1:(k*(k+1))/2)
%      loglikelihood - The loglikelihood of the function at the optimum
%      Ht            - A k x k x t 3 dimension matrix of conditional covariances
%      likelihoods   - A t by 1 vector of individual likelihoods
%      stdresid      - A t by k matrix of multivariate standardized residuals
%      stderrors     - A numParams^2 square matrix of robust Standad Errors(A^(-1)*B*A^(-1)*t^(-1))
%      A             - The estimated inverse of the non-robust Standard errors
%      B             - The estimated covariance of teh scores
%      scores        - A t by numParams matrix of individual scores
%
%
% COMMENTS:
%    You should multiply the data by a constant so that the min std(data) is at least 10.  This will help estimation
%
%      ***************************************************************************************
%      *  THIS FUNCTION INVOLVES ESTIMATING MANY PARAMETERS.  THE EXACT NUMBER OF PARAMETERS
%      *  NEEDING TO BE ESTIMATED IS (k*(k+1))/2+pk^2+qk^2.  FOR A 5 VARIATE (1,1) MODEL THIS
%      *  65 PARAMETERS.  ESTIMATION CAN TAKE A VERY LONG TIME.  A 10 ASSET MODEL TOOK 12
%      *  HOURS ON A PIII-700.
%      ***************************************************************************************
%
% 

 

% need to try and get some smart startgin values
if size(data,2) > size(data,1)
    data=data';
end

[t k]=size(data);
k2=k*(k+1)/2;

scalaropt=optimset('fminunc');
scalaropt=optimset(scalaropt,'TolFun',1e-1,'Display','iter','Diagnostics','on','DiffMaxChange',1e-2);
startingparameters=scalar_bekk_mvgarch(data,p,q,scalaropt);
CChol=startingparameters(1:(k*(k+1))/2);
C=ivech(startingparameters(1:(k*(k+1))/2))*ivech(startingparameters(1:(k*(k+1))/2))';
newA=[];
newB=[];
for i=1:p
    newA=[newA diag(ones(k,1))*startingparameters(((k*(k+1))/2)+i)];
end
for i=1:q
    newB=[newB diag(ones(k,1))*startingparameters(((k*(k+1))/2)+i+p)];
end
newA=reshape(newA,k*k*p,1);
newB=reshape(newB,k*k*q,1);
startingparameters=[CChol;newA;newB];

if nargin<=3 | isempty(BEKKoptions)
    options=optimset('fminunc');
    options.Display='iter';
    options.Diagnostics='on';
    options.TolX=1e-4;
    options.TolFun=1e-4;
    options.MaxFunEvals=5000*length(startingparameters);
    options.MaxIter=5000*length(startingparameters);  
else
    options=BEKKoptions;
end
parameters=fminunc('full_bekk_mvgarch_likelihood',startingparameters,options,data,p,q,k,k2,t);

[loglikelihood,likelihoods,Ht]=full_bekk_mvgarch_likelihood(parameters,data,p,q,k,k2,t);
loglikelihood=-loglikelihood;
likelihoods=-likelihoods;

% Standardized residuals
stdresid=zeros(size(data));
for i=1:t
    stdresid(i,:)=data(i,:)*Ht(:,:,i)^(-0.5);
end

%Std Errors
if nargout>=6
    A=hessian_2sided('full_bekk_mvgarch_likelihood',parameters,data,p,q,k,k2,t);
   
    h=max(abs(parameters/2),1e-2)*eps^(1/3);
    hplus=parameters+h;
    hminus=parameters-h;
    likelihoodsplus=zeros(t,length(parameters));
    likelihoodsminus=zeros(t,length(parameters));
    for i=1:length(parameters)
        hparameters=parameters;
        hparameters(i)=hplus(i);
        [HOLDER, indivlike] = full_bekk_mvgarch_likelihood(hparameters,data,p,q,k,k2,t);
        likelihoodsplus(:,i)=indivlike;
    end
    for i=1:length(parameters)
        hparameters=parameters;
        hparameters(i)=hminus(i);
        [HOLDER, indivlike] = full_bekk_mvgarch_likelihood(hparameters,data,p,q,k,k2,t);
        likelihoodsminus(:,i)=indivlike;
    end
    scores=(likelihoodsplus-likelihoodsminus)./(2*repmat(h',t,1));
    B=cov(scores);
    A=A/t;
    stderrors=A^(-1)*B*A^(-1)*t^(-1);
end

[此贴子已经被作者于2009-4-20 14:38:28编辑过]

专业啊

谢谢达人哈

请问这程序怎么分上下两部分，哪部分才是啊？有人试过吗？我试过，怎么不好用啊？

请问scalar_bekk_mvgarch是怎么定义的？

看不懂啊……难道是直接就可以套用？

顶，先收藏了。

对6楼的问答: 下面这个函数是scalar_bekk_mvgarch
function [parameters, loglikelihood, Ht, likelihoods, stdresid, stderrors, A, B, scores] = scalar_bekk_mvgarch(data,p,q,BEKKoptions);
% PURPOSE:
%      To Estimate a scalar BEKK multivariate GARCH model.  ****SEE WARNING AT END OF HELP FILE****
%
%
% USAGE:
%      [parameters, loglikelihood, Ht, likelihoods, stdresid, stderrors, A, B, scores]  = scalar_bekk_mvgarch(data,p,q,options);
%
%
% INPUTS:
%      data          - A t by k matrix of zero mean residuals
%      p             - The lag length of the innovation process
%      q             - The lag length of the AR process
%      options       - (optional) Options for the optimization(fminunc)
%
%
% OUTPUTS:
%      parameters    - A (k*(k+1))/2+p+q vector of estimated parameteters.
%                         To recover C, use ivech(parmaeters(1:(k*(k+1))/2)
%      loglikelihood - The loglikelihood of the function at the optimum
%      Ht            - A k x k x t 3 dimension matrix of conditional covariances
%      likelihoods   - A t by 1 vector of individual likelihoods
%      stdresid      - A t by k matrix of multivariate standardized residuals
%      stderrors     - A numParams^2 square matrix of robust Standad Errors(A^(-1)*B*A^(-1)*t^(-1))
%      A             - The estimated inverse of the non-robust Standard errors
%      B             - The estimated covariance of the scores
%      scores        - A t by numParams matrix of individual scores
%
%
% COMMENTS:
%      ***************************************************************************************
%      *  THIS FUNCTION INVOLVES ESTIMATING QUITE A FEW PARAMETERS WHEN K IS LARGE.  THE EXACT
%      *  NUMBER OF  PARAMETERS NEEDING TO BE ESTIMATED IS (k*(k+1))/2+p+q.  FOR A 10 VARIATE
%      *  (1,1) MODEL THIS INVLOVES 56 PARAMETERS.  
%      ***************************************************************************************
%
%
% Author: Kevin Sheppard
% kksheppard@ucsd.edu
% Revision: 2    Date: 12/31/2001


% need to try and get some smart startgin values
if size(data,2) > size(data,1)
    data=data';
end

[t k]=size(data);
garchmat=zeros(k,1+p+q);
options=optimset('fmincon');
options=optimset(options,'TolCon',1e-3,'Display','off','Diagnostics','off','LargeScale','off','LevenbergMarquardt','on');
for i=1:k
    temparam=fattailed_garch(data(:,i),p,q,'NORMAL',[],options);
    garchmat(i,:)=temparam';
end

A=mean(garchmat(:,2:p+1));
B=mean(garchmat(:,p+2:p+q+1));


C=cov(data);
alpha0=sqrt(A);
beta0=sqrt(B);

StartC=C*(1-sum(alpha0.^2)-sum(beta0.^2));
CChol=chol(StartC)';
warning off
startingparameters=[vech(CChol);alpha0;beta0];

k2=k*(k+1)/2;
if nargin<=3 | isempty(BEKKoptions)
    options=optimset('fminunc');
    options.Display='iter';
    options.Diagnostics='on';
else
    options=BEKKoptions;
end

parameters=fminunc('scalar_bekk_mvgarch_likelihood',startingparameters,options,data,p,q,k,k2,t);
[loglikelihood,likelihoods,Ht]=scalar_bekk_mvgarch_likelihood(parameters,data,p,q,k,k2,t);
loglikelihood=-loglikelihood;
likelihoods=-likelihoods;

stdresid=zeros(size(data));
for i=1:t
    stdresid(i,:)=data(i,:)*Ht(:,:,i)^(-0.5);
end


%Std Errors
if nargout>=6
    %Std Errors
    A=hessian('scalar_bekk_mvgarch_likelihood',parameters,data,p,q,k,k2,t);
   
    h=max(abs(parameters/2),1e-2)*eps^(1/3);
    hplus=parameters+h;
    hminus=parameters-h;
    likelihoodsplus=zeros(t,length(parameters));
    likelihoodsminus=zeros(t,length(parameters));
    for i=1:length(parameters)
        hparameters=parameters;
        hparameters(i)=hplus(i);
        [HOLDER, indivlike] = scalar_bekk_mvgarch_likelihood(hparameters,data,p,q,k,k2,t);
        likelihoodsplus(:,i)=indivlike;
    end
    for i=1:length(parameters)
        hparameters=parameters;
        hparameters(i)=hminus(i);
        [HOLDER, indivlike] = scalar_bekk_mvgarch_likelihood(hparameters,data,p,q,k,k2,t);
        likelihoodsminus(:,i)=indivlike;
    end
    scores=(likelihoodsplus-likelihoodsminus)./(2*repmat(h',t,1));
    B=cov(scores);
    A=A/t;
    stderrors=A^(-1)*B*A^(-1)*t^(-1);
end

其中还要用到scalar_bekk_mvgarch_likelihood计算似然函数:
function [LLF,likelihoods,Ht]=scalar_bekk_likelihood(parameters,errors,p,q,k,k2,t);
% PURPOSE:
%      To Estimate a scalar BEKK multivariate GARCH likelihood.
%
%
% USAGE:
%      [LLF,likelihoods,Ht]=scalar_bekk_mvgarch_likelihood(parameters,errors,p,q,k,k2,t);
%
%
% INPUTS:
%      parameters - a k*(k+1)/2 + p +q vector of model parameters of the form
%                   [ivech(C);(A(1));...;(A(p));(B(1)); ...(B(q))]
%      errors     - A zeromean t by k martix of residuals
%      p          - The lag length of the innovation process
%      q          - The lag length of the AR process
%      k          - The number of data series
%      k2         - k*k
%      t          - the length of the data series
%
%
% OUTPUTS:
%      LLF           - The loglikelihood of the function at the optimum
%      Ht            - A k x k x t 3 dimension matrix of conditional covariances
%      likelihoods   - A t by 1 vector of individual likelihoods
%
%
% COMMENTS:
%
%
% Author: Kevin Sheppard
% kksheppard@ucsd.edu
% Revision: 2    Date: 12/31/2001

%The first k(k+1)/2 parameters are C, the next p are A, and the next q are B
C=parameters(1:(k2));
A=parameters(k2+1:k2+p);
B=parameters(k2+p+1:k2+p+q);

C=ivech(C);
C=tril(C);
const=C*C';

uncond=cov(errors);
% for starting up, both ee' and H have expectation uncond.  We cna leverage thsi to help the loops.
m=max(p,q);
eeprime=zeros(k,k,t+m);
Ht=zeros(k,k,t+m);
for i=1:m
    eeprime(:,:,i)=uncond;
    Ht(:,:,i)=uncond;
end
LLF=0;
errors=[repmat(diag(uncond)',m,1);errors];
likelihoods=zeros(t+m,1);
for i=m+1:t+m;
    Ht(:,:,i)=const;
    for j=1:p
         Ht(:,:,i)=Ht(:,:,i)+A(j)*(errors(i-j,:))'*(errors(i-j,:))*A(j);
    end
    for j=1:q
         Ht(:,:,i)=Ht(:,:,i)+B(j)*Ht(:,:,i-j)*B(j);
    end
    likelihoods(i)=k*log(2*pi)+(log(det(Ht(:,:,i)))+errors(i,:)*Ht(:,:,i)^(-1)*errors(i,:)');
    LLF=LLF+likelihoods(i);
end
LLF=0.5*(LLF);
likelihoods=0.5*likelihoods(m+1:t+m);
Ht=Ht(:,:,m+1:t+m);
if isnan(LLF)
    LLF=1e6;
end