www.gusucode.com > mbcmodels 工具箱 matlab 源码程序 > mbcmodels/@xregtwostage/pevinit.m
function [TS,OK] =pevinit(TS,Xgc,Yrf,W0,univariate,isCoded); %PEVINIT initializes a two-stage model ready for evaluating pev % Copyright 2000-2015 The MathWorks, Inc. and Ford Global Technologies, Inc. if nargin<5 univariate=0; end if nargin<6 | ~isCoded [Xgc,Yrf,W0,OkGroups,lossdf]= checkdata(TS,Xgc,Yrf,W0); else lossdf=0; end if length(TS.covmodel)==0 TS= covinit(TS,Xgc,Yrf); end nf= size(Yrf,2); nc= length(TS.covmodel); % degrees of freedom % nf*ns - p - numDispersionParams % df= numel(Yrf)-numParams(TS)- nc - length(covmodel(TS.Local)) - lossdf; % local dispersion parameters df= numel(Yrf)-numParams(TS)- lossdf; % local dispersion parameters OK= canMLE(TS) && df>0; if OK if univariate % make covariance diagonal sig= unstruct(TS.covmodel); sig= diag(diag(sig)); TS.covmodel= covupdate(TS.covmodel,sig); lam= RidgeMatrix(TS); else lam= []; end % scale data [Ts2,Xs,Ys,W0s,S]= mle_scale(TS,Xgc,Yrf,W0); Wci= choltinv(Ts2.covmodel,W0s); % make sure we are doing qr on a full matrix if nnz(lam) Xa= [Wci*Xs;lam]; else Xa= Wci*Xs; end [Xs,S]= xregprecond(Xa); R= qr(Xs,0); % check rank of R rd= abs(diag(R)); tol= length(Xa)*eps*max(rd); if tol==0 tol= eps; end OK= ~(size(Xa,2)>size(R,1) | any(rd<tol)); if OK if issparse(R) ri= S/R; else n= size(R,2); ri= S/(triu(R(1:n,:))); end % store results TS= var(TS,ri,[],df); end end