clcclear
y=[15.19	15.91	16.88	18	16.85	17.12];?%數據1行6列的矩陣
n=length（y）;?%?計算出數據y長度
disp（‘級比檢驗‘）
lambda=y（1:end-1）./y（2:end）;
range=minmax（lambda）
range1=[exp（-2/（n+1））exp（2/（n+1））]
yy=ones（n1）;??%?n行1列的1矩陣
for?i=2:n????%for循環?累加數據
yy（i）=yy（i-1）+y（i）;
end???%原始數據累加
B=ones（n-12）;??%?n-1行2列的1矩陣
for?i=1:（n-1）??%for循環?計算B矩陣?詳見灰色預測模型ppt
B（i1）=-（yy（i）+yy（i+1））/2;
B（i2）=1;
end
BT=B‘;??%?‘表示轉置?行列互換
for?j=1:（n-1）???%相當于取y矩陣的后面九個數，其實可以寫的更簡單，YN=y（:2:10）;
YN（j）=y（j+1）;
end
YN=YN‘;?%轉置
A=inv（BT*B）*BT*YN;?%ppt中的公式
a=A（1）;??%?上面計算出來的a是個矩陣，取A矩陣的第一個數
u=A（2）;%?取第二個
t=u/a;??
t_test=input（‘輸入預測個數‘）;%預測后續數據
i=1:t_test+n;??%?相當于i是1到（t_test+n）的數列，間隔為1
yys（i+1）=（y（1）-t）.*exp（-a.*i）+t;?%數學公式詳見ppt??yys為預測數據
yys（1）=y（1）;?%將y矩陣的第一個數字賦給yys矩陣的第一個位置
for?j=n+t_test:-1:2??%for循環?將預測數據yys還原（預測出來的是累加數據）
ys（j）=yys（j）-yys（j-1）;
end
x=1:n;?%x為1到n的數列
xs=2:n+t_test;?%?xs為2到n+t_text的數列
yn=ys（2:n+t_test）;?%將ys矩陣第2到n+t_text的數賦值給yn
xs=[1xs];