PCA故障診斷步驟

1、基于PCA算法的故障診斷步驟離線PCA監(jiān)測模型的計算步驟：（1） 選擇監(jiān)控變量，收集正常工況下的各變量的樣本，記為訓練樣本數(shù)據(jù)X_train和檢驗數(shù)據(jù)X_test；X_train =x11 x12x1mx21 x22x2mxn1 xn2xnmX_train為nm矩陣，即n個樣本，m個觀測變量（即以列向量來看的話，為一個觀測變量各個采樣點的值）對樣本數(shù)據(jù)X_train和檢驗數(shù)據(jù)X_test進行中心化和標準化處理得到和Xstest；中心化處理：按列對X_train減去觀測變量的均值xj=xj-xj 觀測變量某一采樣點的值減去這一觀測變量所有采樣點的平均值xj=1ni=1nxij 求取一列（即某一觀

2、測變量）的平均值標準化處理：對X_train除以觀測變量的標準差（按列（觀測變量）進行）xs j=xjsjsj=i=1n（xij-xj）2n-1 標準差求出標準化矩陣的協(xié)方差矩陣；的協(xié)方差矩陣對為：=1n-1XstainT*Xstrain（2） 對進行特征分解,求得特征值1,2,m（12m）及其對應的特征向量p1,p2,pm（負荷向量）；（3） 確定主元個數(shù), 確定了主元個數(shù)k,就得到了k個特征值12k,及其對應的特征向量p1,p2,pk；A:累計貢獻率法：前k個主元的累積方差貢獻率為： i=1kii=1mi當前k個主元的累積方差貢獻率達到85%，則主元個數(shù)取k值B:交叉檢驗估計法：將采集到的

3、數(shù)據(jù)分成k個部分，1部分數(shù)據(jù)用來建立主元模型,剩下的k-1部分用來作為檢驗數(shù)據(jù)去檢驗所建的模型。如此,建立若干個不同主元個數(shù)的模型,并測試所建立的模型,從中選取一個通過檢驗后誤差最小的模型的主元個數(shù)作為系統(tǒng)主元個數(shù)。（4） 建立PCA主元模型，并進行交叉驗證以確定誤差最小按照ti=XstrainPi,求出第i個主元，并依據(jù)Xstrain=t1p1T+t2p2T+tkpkT =Xstrainp1p1T+Xstrainp2p2T+XstrainpkpkT求出其主元模型用Xstest帶入得到另一主元模型Xstest，依據(jù)E=X-X,求出模型誤差，確定模型誤差最小的那個模型即為主元模型。（5） 計算T

4、2統(tǒng)計量控制限和SPE統(tǒng)計量控制限；對于樣本個數(shù)為n,主元個數(shù)為k的過程變量X_train, T2統(tǒng)計量服從自由度為k和n一k的F分布,則置信度為的T2統(tǒng)計量控制上限為:T2UCL=k(n-1)n-kF(k,n-k)或T2UCL=k(n-1)(n+1)n(n-k)F(k,n-k) 檢驗水平為的SPE統(tǒng)計量控制上限為：QUCL=11+h0C221+2h0(h0-1)121h01=i=k+1mi, 2=i=k+1mi2,3=i=k+1mi3 , h0=1-223322, C是與（1-）分位點對應的標準差在線過程監(jiān)測與故障診斷步驟：（1） 采集第i時刻的在線實時數(shù)據(jù)Xi（Xi為1m矩陣），并進行中心

5、化和標準化處理得到Xis；（2） 按照t(i)(1k)=Xis(1m)pk(mk),求出Xis的得分向量，依據(jù)Xis=tipkT(km),求出PCA模型估計量Xis,這里pk=p1k,p2k,pkk；（3） 計算Xis的T2統(tǒng)計量和SPE統(tǒng)計量，并畫出T2統(tǒng)計量和SPE統(tǒng)計量的控制圖；T(i)2=t(i)k-1t(i)T=Xispkk-1pkTXisTQ(i)=Xis(I-pkpkT)XisT（4） 將上述計算結(jié)果與T2統(tǒng)計量控制限和SPE統(tǒng)計量控制限比較,以檢測過程運行有無異常，當有異常狀態(tài)發(fā)生時,繪制貢獻圖，找出與故障相關的系統(tǒng)變量：1) 檢查每個觀測值x的標準化得分tii, 并確定造成失

6、控狀態(tài)的r(ra)個得分: tiiT2UCLa；2) 計算每個變量xi相對于失控得分ti的貢獻率是:conti,j=tiipi,jxj3) 當conti,j是負時,設它為零；4) 計算第j個過程變量xj的總貢獻率：CONTj=i=1r(conti,j)5) 把所有m個過程變量xj的CONTj畫在一個曲線圖上。PCA_TE仿真程序：%TE過程的傳統(tǒng)主元分析在Matlab中的仿真程序%建立模型：%載入模型數(shù)據(jù)，以故障11為例Xtrain = pkx101; Xtrain = double(Xtrain);%載入測試數(shù)據(jù)Xtest = pkx102;Xtest = double(Xtest);%標準

7、化處理：X_mean = mean(Xtrain); %按列求Xtrain平均值 X_std = std(Xtrain); %求標準差 X_row,X_col = size(Xtrain); %求Xtrain行、列數(shù) % for i = 1:X_col%Xtrain(:,i)=(Xtrain(:,i) - X_mean(i)./X_std(i);%Xtest(:,i) = (Xtest(:,i) - X_mean(i)./X_std(i);% end Xtrain=(Xtrain-repmat(X_mean,X_row,1)./repmat(X_std,X_row,1);%求協(xié)方差矩陣sigm

8、aXtrain = cov(Xtrain);%對協(xié)方差矩陣進行特征分解，lamda為特征值構(gòu)成的對角陣，T的列為單位特征向量，且與lamda中的特征值一一對應：T,lamda = eig(sigmaXtrain); % disp(特征根（由小到大）);% disp(lamda);% disp(特征向量：);% disp(T); %取對角元素(結(jié)果為一列向量)，即lamda值，并上下反轉(zhuǎn)使其從大到小排列，主元個數(shù)初值為1，若累計貢獻率小于90%則增加主元個數(shù)D = flipud(diag(lamda); num_pc = 1; while sum(D(1:num_pc)/sum(D) T2UCL/num_pc r = cat(2,r,i); endend%2.計算每個變量相對于上述失控得分的貢獻cont = zeros(length(r),20);for i = length(r) for j = 1:20 cont(i,j) = abs(S(i)/D(i)*P(j,i)*Xtest(400,j); endend%3.計算每個變量的總貢獻CONTJ = zeros(20,1);for j = 1:20 CONTJ(j) = sum(cont(:,j);end%4.計算每個變量對Q的貢獻e = Xtest(400,:)*(I - P*P);contq =