【故障检测】基于 KPCA 的故障检测(Matlab代码实现)

该博客介绍了如何使用内核PCA(主成分分析)进行故障检测,特别是通过RBF核、多项式核等不同核函数的选择和参数调整来影响决策边界的可视化。文章展示了Matlab代码实现,包括训练、测试数据的处理,以及T2和Q统计指标的计算,用于生成等高线图,以区分正常和错误状态的区域。
摘要由CSDN通过智能技术生成

💥💥💞💞欢迎来到本博客❤️❤️💥💥

🏆博主优势:🌞🌞🌞博客内容尽量做到思维缜密,逻辑清晰,为了方便读者。

⛳️座右铭:行百里者,半于九十。

📋📋📋本文目录如下:🎁🎁🎁

目录

💥1 概述

📚2 运行结果

🎉3 参考文献

🌈4 Matlab代码实现


💥1 概述

数据包含取自模拟流程示例的二维数据集。此数据用于训练和测试内核 PCA 以进行故障检测。训练后,针对输出数据空间中的每个位置计算广泛用于故障检测的T2和Q统计指标,从而生成等高线图。然后将 2% 显著性水平检测限叠加在地图上,作为数据空间的正常(绿色)和错误(洋红色)区域之间的边界。

使用等高线图,可以将各种核类型和参数选择对正常和错误过程状态之间的决策边界的影响可视化。

📚2 运行结果

 

 

 

 

部分代码:

%% Get 2D data
close all; clc; tic;
if nargin == 0
    load dataset.mat p;
    train = p{1}; test = p{2};
    
    % Kernel types and parameters:
    ktype = 'rbf'; kpar = 1;                % RBF kernel
    %ktype = 'rbf'; kpar = 10;               % RBF kernel
    %ktype = 'rbf'; kpar = 0.9;              % RBF kernel
    %ktype = 'rbfpoly'; kpar = [1 1 0.65];   % mixed kernel
    %ktype = 'poly'; kpar = 2;               % polynomial kernel
    %ktype = 'imquad'; kpar = 10;            % inverse multiquadric kernel
    %ktype = 'cauchy'; kpar = 5;             % Cauchy kernel
end

%lax = [-15 15 -15 15];
lax = [-4 10 -3 6];                         % Axes limits
N = length(train); M = length(test);
z0T = train; z1T = test;                    % Training and Test data
[xx,yy] = meshgrid(lax(1):0.05:lax(2),...   % Meshgrid for contours
    lax(3):0.05:lax(4));
z2T = [xx(:) yy(:)]; L = length(z2T);       % Vectorize meshgrid points

K.type = ktype; K.p = kpar;                 % Kernel type and parameters
set(0,'defaultfigurecolor',[1 1 1]);        % Set fig color to w
conf = 0.99;                                % Significance level (*100%)

% Normalize 2D Data
zm = mean(z0T); zs = std(z0T);
z0 = (z0T - zm(ones(N,1),:))./zs(ones(N,1),:);  % Normalize training z
z1 = (z1T - zm(ones(N,1),:))./zs(ones(M,1),:);  % Normalize test z
z2 = (z2T - zm(ones(L,1),:))./zs(ones(L,1),:);  % Normalize surf z

%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  KERNEL PRINCIPAL COMPONENTS ANALYSIS %

[K0c,K0,U0] = kerneltrain(z0,K);            % Populate kernel matrix
K1c = kerneltest(z1,z0,K0,U0,K);            % Project test data to RKHS
K2c = kerneltest(z2,z0,K0,U0,K);            % Project surf data to RKHS

[V,D] = eig(K0c/N);                         % Eigenvalue decomposition
[S,sj] = sort(diag(D),'descend');           % Sort eigenvalues
V = V(:,sj); S = S';                        % Re-arrange eigenvectors
S(S < 1e-7) = [];                           % Remove eigenvalues <= 0
P = V(:,1:length(S))*diag(S.^-0.5);         % Projection matrix

if ~isreal(S)
    disp('Complex eigenvalues detected.');  % Warn about complex eigs
end

%% Perform KPCA Monitoring
CS = cumsum(S)/sum(S)*100;
RP = find(CS >= 99.9,1);                    % Get eigenvalues by %CPV
disp([num2str(RP) ' principal'...
    ' components chosen.']);
t0 = K0c*P(:,1:RP);                         % Kernel subspace (train)
t1 = K1c*P(:,1:RP);                         % Kernel subspace (test)
t2 = K2c*P(:,1:RP);                         % Kernel subspace (surf)
T2 = sum((t0.^2)./S(ones(N,1),1:RP),2);     % T2 statistics (train)
t0n = K0c*P;                                % Full kernel space
Q = abs(sum(t0n.^2,2) - sum(t0.^2,2));      % Q statistics (train)

if strcmp(ktype,'rbf') == 1
    fprintf('\n At infinite fault magnitude:\n');
    U1 = ones(1,N)/N;
    tt = U1*K0*(U0 - eye(N))*P(:,1:RP);
    fprintf('  T2 limit: %.2f\n',...
        sum((tt.^2)./S(1:RP),2));           % Limit of T2 for RBF
    tu = U1*K0*(U0 - eye(N))*P;
    fprintf('  Q limit: %.2f\n\n',...
        abs(sum(tu.^2,2) - sum(tt.^2,2)));  % Limit of Q for RBF
end

T2t = sum((t1.^2)./S(ones(M,1),1:RP),2);    % T2 statistics (test)
T2u = sum((t2.^2)./S(ones(L,1),1:RP),2);    % T2 statistics (surf)
t1n = K1c*P; t2n = K2c*P;
Qt = abs(sum(t1n.^2,2) - sum(t1.^2,2));     % Q statistics (test)
Qu = abs(sum(t2n.^2,2) - sum(t2.^2,2));     % Q statistics (surf)

%% Plot monitoring charts
figure(3); subplot(211);
semilogy(1:N,T2,'b',1:M,T2t,'m','linewidth',1.2);   % T2 monitoring chart
xlabel('Time'); ylabel('T^2'); subplot(212);
semilogy(1:N,Q,'b',1:M,Qt,'m','linewidth',1.2);     % Q monitoring chart
xlabel('Time'); ylabel('Q');

🎉3 参考文献

部分理论来源于网络,如有侵权请联系删除。

[1] K.E.S. Pilario, Y. Cao, and M. Shafiee. Mixed Kernel Canonical Variate Dissimilarity Analysis for Incipient Fault Monitoring in Nonlinear Dynamic Processes. Comput. and Chem. Eng., 123, 143-154. 2019. doi: 10.1016/j.compchemeng.2018.12.027

🌈4 Matlab代码实现

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值