1. 程式人生 > 其它 >數模-灰色關聯分析

數模-灰色關聯分析




水質問題的資料概覽

程式碼

code.m

clear;clc
load data_water_quality.mat
% 不會匯入資料的同學可以看看第二講topsis模型,我們也可以自己在工作區新建變數,並把Excel的資料貼上過來
% 注意Matlab的當前資料夾一定要切換到有資料檔案的這個資料夾內

%%  判斷是否需要正向化
[n,m] = size(X);
disp(['共有' num2str(n) '個評價物件, ' num2str(m) '個評價指標']) 
Judge = input(['這' num2str(m) '個指標是否需要經過正向化處理,需要請輸入1 ,不需要輸入0:  ']);   %1

if Judge == 1
    Position = input('請輸入需要正向化處理的指標所在的列,例如第2、3、6三列需要處理,那麼你需要輸入[2,3,6]: '); %[2,3,4]
    disp('請輸入需要處理的這些列的指  標型別(1:極小型, 2:中間型, 3:區間型) ')
    Type = input('例如:第2列是極小型,第3列是區間型,第6列是中間型,就輸入[1,3,2]:  '); %[2,1,3]
    % 注意,Position和Type是兩個同維度的行向量
    for i = 1 : size(Position,2)  %這裡需要對這些列分別處理,因此我們需要知道一共要處理的次數,即迴圈的次數
        X(:,Position(i)) = Positivization(X(:,Position(i)),Type(i),Position(i));
    % Positivization是我們自己定義的函式,其作用是進行正向化,其一共接收三個引數
    % 第一個引數是要正向化處理的那一列向量 X(:,Position(i))   回顧上一講的知識,X(:,n)表示取第n列的全部元素
    % 第二個引數是對應的這一列的指標型別(1:極小型, 2:中間型, 3:區間型)
    % 第三個引數是告訴函式我們正在處理的是原始矩陣中的哪一列
    % 該函式有一個返回值,它返回正向化之後的指標,我們可以將其直接賦值給我們原始要處理的那一列向量
    end
    disp('正向化後的矩陣 X =  ')
    disp(X)
end

%% 對正向化後的矩陣進行預處理
Mean = mean(X);  % 求出每一列的均值以供後續的資料預處理
Z = X ./ repmat(Mean,size(X,1),1);  
disp('預處理後的矩陣為:'); disp(Z)

%% 構造母序列和子序列
Y = max(Z,[],2);  % 母序列為虛擬的,用每一行的最大值構成的列向量表示母序列
X = Z; % 子序列就是預處理後的資料矩陣

%% 計算得分
absX0_Xi = abs(X - repmat(Y,1,size(X,2)))  % 計算|X0-Xi|矩陣
a = min(min(absX0_Xi))    % 計算兩級最小差a
b = max(max(absX0_Xi))  % 計算兩級最大差b
rho = 0.5; % 分辨係數取0.5
gamma = (a+rho*b) ./ (absX0_Xi  + rho*b)  % 計運算元序列中各個指標與母序列的關聯絡數
weight = mean(gamma) / sum(mean(gamma));  % 利用子序列中各個指標的灰色關聯度計算權重
score = sum(X .* repmat(weight,size(X,1),1),2);   % 未歸一化的得分
stand_S = score / sum(score);   % 歸一化後的得分
[sorted_S,index] = sort(stand_S ,'descend') % 進行排序

Inter2Max.m

function [posit_x] = Inter2Max(x,a,b)
    r_x = size(x,1);  % row of x 
    M = max([a-min(x),max(x)-b]);
    posit_x = zeros(r_x,1);   %zeros函式用法: zeros(3)  zeros(3,1)  ones(3)
    % 初始化posit_x全為0  初始化的目的是節省處理時間
    for i = 1: r_x
        if x(i) < a
           posit_x(i) = 1-(a-x(i))/M;
        elseif x(i) > b
           posit_x(i) = 1-(x(i)-b)/M;
        else
           posit_x(i) = 1;
        end
    end
end

Mid2Max.m

function [posit_x] = Mid2Max(x,best)
    M = max(abs(x-best));
    posit_x = 1 - abs(x-best) / M;
end

Min2Max.m

function [posit_x] = Min2Max(x)
    posit_x = max(x) - x;
     %posit_x = 1 ./ x;    %如果x全部都大於0,也可以這樣正向化
end

Positivization.m

% function [輸出變數] = 函式名稱(輸入變數)  
% 函式的中間部分都是函式體
% 函式的最後要用end結尾
% 輸出變數和輸入變數可以有多個,用逗號隔開
% function [a,b,c]=test(d,e,f)
%     a=d+e;
%     b=e+f;
%     c=f+d;
% end
% 自定義的函式要單獨放在一個m檔案中,不可以直接放在主函式裡面(和其他大多數語言不同)

function [posit_x] = Positivization(x,type,i)
% 輸入變數有三個:
% x:需要正向化處理的指標對應的原始列向量
% type: 指標的型別(1:極小型, 2:中間型, 3:區間型)
% i: 正在處理的是原始矩陣中的哪一列
% 輸出變數posit_x表示:正向化後的列向量
    if type == 1  %極小型
        disp(['第' num2str(i) '列是極小型,正在正向化'] )
        posit_x = Min2Max(x);  %呼叫Min2Max函式來正向化
        disp(['第' num2str(i) '列極小型正向化處理完成'] )
        disp('~~~~~~~~~~~~~~~~~~~~分界線~~~~~~~~~~~~~~~~~~~~')
    elseif type == 2  %中間型
        disp(['第' num2str(i) '列是中間型'] )
        best = input('請輸入最佳的那一個值: ');
        posit_x = Mid2Max(x,best);
        disp(['第' num2str(i) '列中間型正向化處理完成'] )
        disp('~~~~~~~~~~~~~~~~~~~~分界線~~~~~~~~~~~~~~~~~~~~')
    elseif type == 3  %區間型
        disp(['第' num2str(i) '列是區間型'] )
        a = input('請輸入區間的下界: ');
        b = input('請輸入區間的上界: '); 
        posit_x = Inter2Max(x,a,b);
        disp(['第' num2str(i) '列區間型正向化處理完成'] )
        disp('~~~~~~~~~~~~~~~~~~~~分界線~~~~~~~~~~~~~~~~~~~~')
    else
        disp('沒有這種型別的指標,請檢查Type向量中是否有除了1、2、3之外的其他值')
    end
end

讀入的資料

第一列越大越好,第二列中間型,第三列極小型,第四列區間型

code執行

最終結果

如果只是需要係數的話sorted_S就是需要的係數,後面的index是對應的河流