function?varargout?=?arburgwithcriterion（?xcriterion?）
%ARBURGWITHCRITERION???AR?parameter?estimation?via?Burg?method.
%???A?=?ARBURGWITHCRITERION（X）?returns?the?polynomial?A?corresponding?to?the?AR
%???parametric?signal?model?estimate?of?vector?X?using?Burg‘s?method.
%???模型階次選擇默認(rèn)為最小預(yù)測(cè)定誤差階準(zhǔn)則（Final?Prediction?Error?CriterionFPE）
%???可選用的模型定階準(zhǔn)則有FPE?AIC?MDL?CAT定階準(zhǔn)則
%
%???[AE]?=?ARBURGWITHCRITERION（...）?returns?the?final?prediction?error?E?（the?variance
%???estimate?of?the?white?noise?input?to?the?AR?model）.
%
%???[AEK]?=?ARBURGWITHCRITERION（...）?returns?the?vector?K?of?reflection?
%???coefficients?（parcor?coefficients）.
%
%???[AEKP]?=?ARBURGWITHCRITERION（...）returns?the?order?P?of?the?AR
%???model.

%???Ref:?祈才君，數(shù)字信號(hào)處理技術(shù)的算法分析與應(yīng)用
%??????????????機(jī)械工業(yè)出版社?2005?Chapter?8

%???最初版本：V1.0???????時(shí)間：2009.03.31

error（nargchk（12nargin））

if?nargin?==?1
????criterion?=?‘FPE‘;
end

[mxnx]?=?size（x）;
if?isempty（x）?||?min（mxnx）?>?1
???error（‘X?must?be?a?vector?with?length?greater?than?twice?the?model?order.‘）;
elseif?~（strcmp（criterion‘FPE‘）?||?strcmp（criterion‘AIC‘）?||?strcmp（criterion‘MDL‘）?||?strcmp（criterion‘CAT‘））
???error（‘不正確的定階準(zhǔn)則。‘）
end
if?issparse（x）
???error（‘Input?signal?cannot?be?sparse.‘）
end

x??=?x（:）;
N??=?length（x）;

%?Initialization
orderpridict?=?1;
ef?=?x;
eb?=?x;
a?=?1;

%?Initial?error
E?=?x‘*x./N;
if?strcmp（criterion‘FPE‘）
????objectfun?=?（N+1）/（N-1）*E;
elseif?strcmp（criterion‘AIC‘）
????objectfun?=?N*log（E）+2*1;
elseif?strcmp（criterion‘MDL‘）
????objectfun?=?N*log（E）+1*log（N）;
elseif?strcmp（criterion‘CAT‘）
????temp?=?（N-1）/（N*E）;
????objectfun?=?1/N*（N-1）/（N*E）-（N-1）/（N*E）;
end

%?Preallocate?‘k‘?for?speed.
%?一般AR模型階次Pk?=?zeros（1?N）;

for?m?=?1:N-1
???%?Calculate?the?next?order?reflection?（parcor）?coefficient
???efp?=?ef（2:end）;
???ebp?=?eb（1:end-1）;
???num?=?-2.*ebp‘*efp;
???den?=?efp‘*efp+ebp‘*ebp;
???
???k（m）?=?num?./?den;
???
???%?Update?the?forward?and?backward?prediction?errors
???ef?=?efp?+?k（m）*ebp;
???eb?=?ebp?+?k（m）‘*efp;
???
???%?Update?the?AR?coeff.
???a=[a;0]?+?k（m）*[0;conj（flipud（a））];
???
???%?Update?the?prediction?error
???E（m+1）?=?（1?-?k（m）‘*k（m））*E（m）;
???
???%?判斷是否達(dá)到所選定階準(zhǔn)則的要求
???if?strcmp（criterion‘FPE‘）
???????objectfun（m+1）?=?（N+（m+1））/（N-（m+1））*E（m+1）;
???elseif?strcmp（criterion‘AIC‘）
???????objectfun（m+1）?=?N*log（E（m+1））+2*（m+1）;
???elseif?strcmp（criterion‘MDL‘）
???????objectfun（m+1）?=?N*log（E（m+1））+（m+1）*log（N）;
???elseif?strcmp（criterion‘CAT‘）
???????for?index?=?1:m+1
???????????temp?=?temp+（N-index）/（N*E（index））;
???????end
???????objectfun（m+1）?=?1/N*temp-（N-（m+1））/（N*E（m+1））;
???end
???
???if?objectfun（m+1）?>=?objectfun（m）?%?該階預(yù)測(cè)誤差大于上一階預(yù)測(cè)誤差，則階次定為上一階
???????orderpredict?=?m;
???????break;
???end
end

%?實(shí)際k矩陣
k?=?k（1:orderpredict）;

a?=?a（:）.‘;?%?By?convention?all?polynomials?are?row?vectors
varargout{1}?=?a（1:orderpredict）;
if?nargout?>=?2
????varargout{2}?=?E（ord

?屬性????????????大小?????日期????時(shí)間???名稱
-----------?---------??----------?-----??----

?????文件???????3290??2009-08-28?21:23??arburgwithcriterion.m

?????文件???????5152??2009-08-28?21:23??arspectrawithcriterion.m

?????文件???????1652??2009-08-28?21:23??criterion.m

?????文件???????3243??2009-08-28?21:23??pburgwithcriterion.m

-----------?---------??----------?-----??----

????????????????13337????????????????????4