[問題] 照著論文做狀態分析但就是無法得出一樣

[dennisshyu12 (路人)]

想請教版上的大神，是否是小弟程式有誤導至結果不同，或者是公式的部分打錯了...
請大神幫幫忙

[狀況說明]
四種狀態數據分數階導數後，做可拓的物質元素分析 ，皆照著步驟與公式做，
但結果就是與論文不同，達不到百分之百，
論文名:Fractional-Order Chaotic Self-Synchronization Based Tracking Faults Diagnosis
of  Ball Bearing System
以下是我的程式，附上步驟與說明

%-------------------------------------------------------------------------------
% % 步驟
% training
% -------------------------------
% 1.從測試集中，隨機取一段距離(一千筆數據)，將4種狀態的數據代到分數階混沌方程
，並把結果作歐式(norm)存到一個矩陣，重複作直到迴圈結束(測試也一樣，不過只要代
入一種數據即可)。
%
% 2.尋找各種狀態得最大最小值，儲存為經典域，並各往外擴長10%
% testing
% ----------------------------
% 3.將每筆測試數據，計算出與可拓域與經典域之距離
%
% 4.計算每筆測試數據與4種狀態得關聯度，並進行正規化
%
% 5.比較結果，最大值即判斷為該種狀態。

%4種DATA合併
S(:,1) = X097_DE_time(1:240000); %normal
S(:,2) = X122_DE_time(1:240000); %ball
S(:,3) = X109_DE_time(1:240000); %inner
S(:,4) = X135_DE_time(1:240000); %outer

%參數固定
a = 5;
b = -10;
c = -3.8;
alpha = 0.3;

% % trainning
%
-----------------------------------------------------------------------------------------------------------------
R=randi([48000,134000],1,50); % 隨機幾段 %幾次
for s=1:4
    if s==1  % s 是狀態
        for i = 1:50
            for T = R(i) : (R(i)+1000)
                %以下是 分數階混沌方程 公式
                E(1,s) = S(T,s);
                E(2,s) = S(T+1,s);
                E(3,s) = S(T+2,s);
                g1= [ a , -E(3,s) , 0 ; E(3,s) , b , 0 ; 1/3*E(2,s) , 0 , c ];
                h1= ( gamma(2) / gamma( 2 - alpha )) ;
                h2= [ E(1,s)^(1- alpha ) ; E(2,s)^(1-alpha) ;
E(3,s)^(1-alpha) ];
                e1(:,i)=g1*h2.*h1;
            end
        O1(1,i) = norm(e1(1,:));
        O1(2,i) = norm(e1(2,:));
        end
    end
    if s==2
        for i = 1:50
            for T = R(i) : (R(i)+1000)
                E(1,s) = S(T,s);
                E(2,s) = S(T+1,s);
                E(3,s) = S(T+2,s);
                g1= [ a , -E(3,s) , 0 ; E(3,s) , b , 0 ; 1/3*E(2,s) , 0 , c ];
                h1= ( gamma(2) / gamma( 2 - alpha )) ;
                h2= [ E(1,s)^(1- alpha ) ; E(2,s)^(1-alpha) ;
E(3,s)^(1-alpha) ];
                e2(:,i)=g1*h2.*h1;
            end
        O2(1,i) = norm(e2(1,:));
        O2(2,i) = norm(e2(2,:));
        end

    end
    if s==3
        for i = 1:50
            for T = R(i) : (R(i)+1000)
                E(1,s) = S(T,s);
                E(2,s) = S(T+1,s);
                E(3,s) = S(T+2,s);
                g1= [ a , -E(3,s) , 0 ; E(3,s) , b , 0 ; 1/3*E(2,s) , 0 , c ];
                h1= ( gamma(2) / gamma( 2 - alpha )) ;
                h2= [ E(1,s)^(1- alpha ) ; E(2,s)^(1-alpha) ;
E(3,s)^(1-alpha) ];
                e3(:,i)=g1*h2.*h1;
            end
        O3(1,i) = norm(e3(1,:));
        O3(2,i) = norm(e3(2,:));
        end
    end
    if s==4
        for i = 1:50
            for T = R(i) : (R(i)+1000)
                E(1,s) = S(T,s);
                E(2,s) = S(T+1,s);
                E(3,s) = S(T+2,s);
                g1= [ a , -E(3,s) , 0 ; E(3,s) , b , 0 ; 1/3*E(2,s) , 0 , c ];
                h1= ( gamma(2) / gamma( 2 - alpha )) ;
                h2= [ E(1,s)^(1- alpha ) ; E(2,s)^(1-alpha) ;
E(3,s)^(1-alpha) ];
                e4(:,i)=g1*h2.*h1;
            end
        O4(1,i) = norm(e4(1,:));
        O4(2,i) = norm(e4(2,:));
        end
    end
end

% % 找各種狀態大小 存為經典域 並擴大10%當做可拓
%
-----------------------------------------------------------------------------------------------
% % 經典
 N1(1,:) = [min(O1(1,:)), max(O1(1,:))];
 N1(2,:) = [min(O1(2,:))  max(O1(2,:))];
 N2(1,:) = [min(O2(1,:)), max(O2(1,:))];
 N2(2,:) = [min(O2(2,:))  max(O2(2,:))];
 N3(1,:) = [min(O3(1,:)), max(O3(1,:))];
 N3(2,:) = [min(O3(2,:))  max(O3(2,:))];
 N4(1,:) = [min(O4(1,:)), max(O4(1,:))];
 N4(2,:) = [min(O4(2,:))  max(O4(2,:))];
 % % 可拓
 EX1(:,1) = N1(:,1).*0.9;
 EX1(:,2) = N1(:,2).*1.1;

 EX2(:,1) = N2(:,1).*0.9;
 EX2(:,2) = N2(:,2).*1.1;

 EX3(:,1) = N3(:,1).*0.9;
 EX3(:,2) = N3(:,2).*1.1;

 EX4(:,1) = N4(:,1).*0.9;
 EX4(:,2) = N4(:,2).*1.1;

%
----------------------------------------------------------------------------------
%每筆測試數據
m=10;
R=randi([48000,134000],1, m);
s=3%第幾種狀態
if s==3
    for i = 1: m
            for T = R(i) : (R(i)+1000)
                E(1,s) = S(T,s);
                E(2,s) = S(T+1,s);
                E(3,s) = S(T+2,s);
                g1= [ a , -E(3,s) , 0 ; E(3,s) , b , 0 ; 1/3*E(2,s) , 0 , c ];
                h1= ( gamma(2) / gamma( 2 - alpha )) ;
                h2= [ E(1,s)^(1- alpha ) ; E(2,s)^(1-alpha) ;
E(3,s)^(1-alpha) ];
                e1(:,i)=g1*h2.*h1;
            end
        O5(1,i) = norm(e1(1,:));
        O5(2,i) = norm(e1(2,:));
        end
end

% % --------------------------------------------------
%計算出待測物(O5)與可拓經典域的距離
for i = 1 : m
    %經典
    D1(1,i) = abs( O5(1,i)-(N1(1,1)+N1(1,2))/2 ) - ((N1(1,2)-N1(1,1))/2) ;
    D1(2,i) = abs( O5(2,i)-(N1(2,1)+N1(2,2))/2 ) - ((N1(2,2)-N1(2,1))/2) ;

    D2(1,i) = abs( O5(1,i)-(N2(1,1)+N2(1,2))/2 ) - ((N2(1,2)-N2(1,1))/2) ;
    D2(2,i) = abs( O5(2,i)-(N2(2,1)+N2(2,2))/2 ) - ((N2(2,2)-N2(2,1))/2) ;

    D3(1,i) = abs( O5(1,i)-(N3(1,1)+N3(1,2))/2 ) - ((N3(1,2)-N3(1,1))/2) ;
    D3(2,i) = abs( O5(2,i)-(N3(2,1)+N3(2,2))/2 ) - ((N3(2,2)-N3(2,1))/2) ;

    D4(1,i) = abs( O5(1,i)-(N4(1,1)+N4(1,2))/2 ) - ((N4(1,2)-N4(1,1))/2) ;
    D4(2,i) = abs( O5(2,i)-(N4(2,1)+N4(2,2))/2 ) - ((N4(2,2)-N4(2,1))/2) ;

    % 可拓
    E1(1,i) = abs( O5(1,i) - (EX1(1,1)+EX1(1,2))/2 ) -
((EX1(1,2)-EX1(1,1))/2) ;
    E1(2,i) = abs( O5(2,i) - (EX1(2,1)+EX1(2,2))/2 ) -
((EX1(2,2)-EX1(2,1))/2) ;

    E2(1,i) = abs( O5(1,i) - (EX2(1,1)+EX2(1,2))/2 ) -
((EX2(1,2)-EX2(1,1))/2) ;
    E2(2,i) = abs( O5(2,i) - (EX2(2,1)+EX2(2,2))/2 ) -
((EX2(2,2)-EX2(2,1))/2) ;

    E3(1,i) = abs( O5(1,i) - (EX3(1,1)+EX3(1,2))/2 ) -
((EX3(1,2)-EX3(1,1))/2) ;
    E3(2,i) = abs( O5(2,i) - (EX3(2,1)+EX3(2,2))/2 ) -
((EX3(2,2)-EX3(2,1))/2) ;

    E4(1,i) = abs( O5(1,i) - (EX4(1,1)+EX4(1,2))/2 ) -
((EX4(1,2)-EX4(1,1))/2) ;
    E4(2,i) = abs( O5(2,i) - (EX4(2,1)+EX4(2,2))/2 ) -
((EX4(2,2)-EX4(2,1))/2) ;

end

% %
%
---------------------------------------------------------------------------------------
%關聯函式 判斷待測目標(O5)與四種狀態的關係
% 狀態一
for i = 1 : m
    %狀態1 e1
    if  O5(1,i) >=N1(1,1)  && O5(1,i) <=N1(1,2)
        K1(1,i)  = -D1(1,i) / abs( N1(1,1)-N1(1,2) );
    else
        K1(1,i) = D1(1,i) / ( E1(1,i) - D1(1,i));
    end
    %狀態1 e1
    if  O5(2,i) >=N1(2,1)  && O5(2,i) <=N1(2,2)
        K1(2,i)  = -D1(2,i) / abs(N1(2,1)-N1(2,2));
    else
        K1(2,i) = D1(2,i) / ( E1(2,i) - D1(2,i));
    end
    Kk1(1,i) = K1(1,i)*0.5+ K1(2,i)*0.5;
end
Kk1;
% 狀態二
for i = 1 : m
    %狀態2 e1
    if  O5(1,i) >=N2(1,1)  && O5(1,i) <=N2(1,2)
        K2(1,i)  = -D2(1,i) / abs( N2(1,1)-N2(1,2) );
    else
        K2(1,i) = D2(1,i) / ( E2(1,i) - D2(1,i));
    end
    %狀態2 e1
    if  O5(2,i) >=N2(2,1)  && O5(2,i) <=N2(2,2)
        K2(2,i)  = -D2(2,i) / abs(N2(2,1)-N2(2,2));
    else
        K2(2,i) = D2(2,i) / ( E2(2,i) - D2(2,i));
    end
    Kk2(1,i) = K2(1,i)*0.5+ K2(2,i)*0.5;
end
Kk2;
% 狀態三
for i = 1 : m
    %狀態3 e1
    if  O5(1,i) >=N3(1,1)  && O5(1,i) <=N3(1,2)
        K3(1,i)  = -D3(1,i) / abs( N3(1,1)-N3(1,2) );
    else
        K3(1,i) = D3(1,i) / ( E3(1,i) - D3(1,i));
    end
    %狀態3 e2
    if  O5(2,i) >=N3(2,1)  && O5(2,i) <=N3(2,2)
        K3(2,i)  = -D3(2,i) / abs(N3(2,1)-N3(2,2));
    else
        K3(2,i) = D3(2,i) / ( E3(2,i) - D3(2,i));
    end
    Kk3(1,i) = K3(1,i)*0.5+ K3(2,i)*0.5;
end
Kk3;
for i = 1 : m
    %狀態4 e1
    if  O5(1,i) >=N4(1,1)  && O5(1,i) <=N4(1,2)
        K4(1,i)  = -D4(1,i) / abs( N4(1,1)-N4(1,2) );
    else
        K4(1,i) = D4(1,i) / ( E4(1,i) - D4(1,i));
    end
    %狀態4 e2
    if  O5(2,i) >=N4(2,1)  && O5(2,i) <=N4(2,2)
        K4(2,i)  = -D4(2,i) / abs(N4(2,1)-N4(2,2));
    else
        K4(2,i) = D4(2,i) / ( E4(2,i) - D4(2,i));
    end
    Kk4(1,i) = K4(1,i)*0.5+ K4(2,i)*0.5;
end
Kk4;
F=[Kk1;Kk2;Kk3;Kk4]';

%正規化

for  j= 1:m
    for i=1:4

Q(j,i)=(2*F(j,i)-max(F(j,1:4))-min(F(j,1:4)))/(max(F(j,1:4))-min(F(j,1:4)));
    end
end
Q

--
※ 發信站: 批踢踢實業坊(ptt.cc), 來自: 111.250.212.32 (臺灣)
※ 文章網址: https://www.ptt.cc/bbs/MATLAB/M.1613924874.A.2E1.html
推 book5601: 對於這種請別人debug 的文章，好像都乏人問津。勸你各位
→ book5601: 做論文還是要多主動思考把大問題分割成數個小問題，確認
→ book5601: 每個環節都在自己的掌握之後再連結各個小問題，問題應該
→ book5601: 就解決了。但不能排除有許多條件是論文沒有提及的，最快
→ book5601: 且最不出錯的方法是問本人。
噓 jatj: ......
→ guteres: 推一樓，我想作者都有留連絡方式，我們以前都是直接問
推 el3qu4: 偷偷跟你說 論文都會暗藏小地方不可能全部都給你看....問
→ el3qu4: 原作者可能好一點