clear;?

P?=?sin（0.1:0.1:9.6）;
F?=?sin（0.1:0.1:9）;
?
%----------------------由于時間序列有不平穩趨勢，進行兩次差分運算，消除趨勢性----------------------%?
for?i=2:96?
????Yt（i）=P（i）-P（i-1）;?
end?
for?i=3:96?
????L（i）=Yt（i）-Yt（i-1）;?
end?
L=L（3:96）;?
Y=L（1:88）;?

%畫圖
figure;?
plot（P）;?
title（‘原數據序列圖‘）;?
hold?on;?
pause??
plot（Y‘r‘）;?
title（‘兩次差分后的序列圖和原數對比圖‘）;?
pause???

%--------------------------------------對數據標準化處理----------------------------------------------%?
%處理的算法?:?（data?-?期望）/方差
Ux=sum（Y）/88???????????????????????????%?求序列均值?
yt=Y-Ux;?
b=0;?
for?i=1:88?
???b=yt（i）^2/88+b;?
end?
v=sqrt（b）??????????????????????????????%?求序列方差?
Y=yt/v;????????????????????????????%?標準化處理公式?
f=F（1:88）;?
t=1:88;?

%畫圖
figure;?
plot（tftY‘r‘）?
title（‘原始數據和標準化處理后對比圖‘）;?
xlabel（‘時間t‘）ylabel（‘油價y‘）;?
legend（‘原始數據?F?‘‘標準化后數據Y?‘）;?
pause???
%--------------------------------------對數據標準化處理----------------------------------------------%?
?
?
%------------------------檢驗預處理后的數據是否符合AR建模要求，計算自相關和偏相關系數---------------%?

%---------------------------------------計算自相關系數-----------------------------------%?
R0=0;
for?i=1:88??
?????R0=Y（i）^2/88+R0;???%標準化處理后的數據的方差
end?

for?k=1:20?
????
????%R??協方差???
????R（k）=0;?
????for?i=k+1:88
????????R（k）=Y（i）*Y（i-k）/88+R（k）;???
????end?
end?
x=R/R0??????????????????????%自相關系數x?=?協方差/方差

%畫圖
figure;?
plot（x）?
title（‘自相關系數分析圖‘）;?
pause???
%-----------------------------------計算自相關系數-------------------------------------%?
?
%-----------------------解Y-W方程，其系數矩陣是Toeplitz矩陣（多普里茲矩陣）。求得偏相關函數X-------------------

X1=x（1）;?
X11=x（1）;?
B=[x（1）?x（2）]‘;?
x2=[1?x（1）];?
A=toeplitz（x2）;???????????????????????
X2=A\B??????????????????????????%x=a\b是方程a*x?=b的解
X22=X2（2）?
?
B=[x（1）?x（2）?x（3）]‘;?
x3=[1?x（1）?x（2）];?
A=toeplitz（x3）;???????????????????????
X3=A\B?
X33=X3（3）?
?
B=[x（1）?x（2）?x（3）?x（4）]‘;?
x4=[1?x（1）?x（2）?x（3）];?
A=toeplitz（x4）;???????????????????????
X4=A\B?
X44=X4（4）?
?
B=[x（1）?x（2）?x（3）?x（4）?x（5）]‘;?
x5=[1?x（1）?x（2）?x（3）?x（4）];?
A=toeplitz（x5）;???????????????????????
X5=A\B?
X55=X5（5）?
?
B=[x（1）?x（2）?x（3）?x（4）?x（5）?x（6）]‘;?
x6=[1?x（1）?x（2）?x（3）?x（4）?x（5）];?
A=toeplitz（x6）;???????????????????????
X6=A\B?
X66=X6（6）?
?
B=[x（1）?x（2）?x（3）?x（4）?x（5）?x（6）?x（7）]‘;?
x7=[1?x（1）?x（2）?x（3）?x（4）?x（5）?x（6）];?
A=toeplitz（x7）;???????????????????????
X7=A\B?
X77=X7（7）?
?
B=[x（1）?x（2）?x（3）?x（4）?x（5）?x（6）?x（7）?x（8）]‘;?
x8=[1?x（1）?x（2）?x（3）?x（4）?x（5）?x（6）?x（7）];?
A=toeplitz（x8）;???????????????????????
X8=A\B?
X88=X8（8）?
?
B=[x（1）?x（2）?x（3）?x（4）?x（5）?x（6）?x（7）?x（8）?x（9）]‘;?
x9=[1?x（1）?x（2）?x（3）?x（4）?x（5）?x（6）?x（7）?x（8）];?
A=toeplitz（x9）;???????????????????????
X9=A\B?
X99=X9（9）?
?
B=[x（1）?x（2）?x（3）?x（4）?x（5）?x（6）?x（7）?x（8）?x（9）?x（10）]‘;?
x10=[1?x（1）?x（2）?x（3）?x（4）?x（5）?x（6）?x（7）?x（8）?x（9）];?
A=toeplitz（x10）;???????????????????????
X10=A\B????
X1010=X10（10）?
??????
B=[x（1）?x（2）?x（3）?x（4）?x（5）?x（6）?x（7）?x（8）?x（9）?x（10）?x（11）]‘;?
x11=[1?x（1）?x（2）?x（3）?x（4）?x（5）?