【故障诊断】用于故障检测的内核 PCA 等值线图（Matlab实现）

荧光Matlab

于 2024-10-10 17:12:44 发布

阅读量241

点赞数 2

文章标签： matlab 人工智能算法

本文链接：https://blog.csdn.net/Matlab88888/article/details/142826960

版权

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

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

⛳️座右铭：行百里者，半于九十。

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

💥1 概述

用于故障检测的内核 PCA（Kernel Principal Component Analysis，核主成分分析）等值线图具有重要的应用价值。内核 PCA 是一种强大的数据分析技术，它能够有效地处理非线性数据。在故障检测中，内核 PCA 可以将高维的故障数据映射到低维空间，从而揭示数据中的潜在结构和模式。等值线图是内核 PCA 在故障检测中的一种可视化工具。通过绘制等值线图，可以直观地展示数据在低维空间中的分布情况。在故障检测中，等值线图可以帮助工程师快速识别出异常数据点，从而及时发现故障。具体来说，用于故障检测的内核 PCA 等值线图具有以下几个特点： 1. 非线性处理能力：内核 PCA 可以处理非线性数据，能够更好地适应复杂的故障检测场景。 2. 降维可视化：等值线图将高维数据映射到低维空间，使得数据的分布情况更加直观，便于工程师进行故障分析。 3. 异常检测：通过观察等值线图中的异常数据点，可以快速发现故障，提高故障检测的效率。 4. 多变量分析：内核 PCA 可以同时处理多个变量，等值线图可以展示多个变量之间的关系，为故障检测提供更全面的信息。总之，用于故障检测的内核 PCA 等值线图是一种有效的故障检测工具，它能够帮助工程师快速发现故障，提高故障检测的准确性和效率。

📚2 运行结果

部分函数代码：


%% Get 2D data
close all; clc;
if nargin == 0
    load dataset.mat p;
    train = p{1}; test = p{2};
end
N = length(train); M = length(test);
z0T = train; z1T = test;                    % Training and Test data

%lax = [-10 15 -10 10];
lax = [-4 10 -3 6];                         % Axes limits
[xx,yy] = meshgrid(lax(1):0.02:lax(2),...   % Meshgrid for contours
    lax(3):0.02:lax(4));
z2T = [xx(:) yy(:)]; L = length(z2T);       % Vectorize meshgrid points
set(0,'defaultfigurecolor',[1 1 1]);        % Set fig color to w
conf = 0.99;                                % Significance level (*100%)

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

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

%% PCA Transformation
[V,D,~] = eig(z0'*z0/(N-1));                % Eigenvalue decomposition
[S,sj] = sort(diag(D),'descend');           % Sort eigenvalues
V = V(:,sj); S = S';                        % Re-arrange eigenvectors

P = V*diag(S.^-0.5);                        % Projection matrix
p0 = z0*P; p1 = z1*P; p2 = z2*P;            % Principal vectors
P = eye(size(P)) - P*P';                    % Residual projection matrix
r0 = z0*P; r1 = z1*P; r2 = z2*P;            % Residual vectors

%% Computation of T2 and Q UCLs using KDE
T2 = sum(p0.^2,2);                          % Compute T2 statistics
Q = sum(r0.^2,2);                           % Compute Q statistics
[Ap1.f,Ap1.x] = ksdensity(T2);              % Estimate PDF of T2
Af1 = cumsum([0;diff(Ap1.x(:))].*Ap1.f(:)); % Numerical integration on PDF
[Ap2.f,Ap2.x] = ksdensity(Q);               % Estimate PDF of Q
Af2 = cumsum([0;diff(Ap2.x(:))].*Ap2.f(:)); % Numerical integration on PDF

T2ucl = max(Ap1.x(Af1<=conf));              % Get T2 upper control limit
Qucl = max(Ap2.x(Af2<=conf));               % Get Q upper control limit

T2t = sum(p1.^2,2); Qt = sum(r1.^2,2);      % Statistics of Test Data
T2u = sum(p2.^2,2); Qu = sum(r2.^2,2);      % Statistics of Surf Data

%% Plotting of PCA Monitoring Charts
figure(3); subplot(211);
semilogy(1:N,T2,'b',1:M,T2t,'m');           % T2 monitoring chart
xlabel('Time'); ylabel('T^2');
figure(3); subplot(212);
semilogy(1:N,Q,'b',1:M,Qt,'m');             % Q monitoring chart
xlabel('Time'); ylabel('Q');

% Plot UCLs on the monitoring charts
figure(3); subplot(211); hold on;
semilogy([1 max(N,M)],T2ucl*[1 1],'r--');
legend({'Training data','Test data','T^2 UCL'},...
    'Location','southwest');
hold off; title('PCA Monitoring Charts');
subplot(212); hold on;
semilogy([1 max(N,M)],Qucl*[1 1],'r-.'); 
legend({'Training data','Test data','Q UCL'},...
    'Location','southwest'); hold off;

% Detection Time
fprintf('Detection Time for T2: Sample %d\n',detecTime(T2ucl,T2t));
fprintf('Detection Time for Q:  Sample %d\n',detecTime(Qucl,Qt));

%% Plot T2 contour map
T2u = log10(T2u); T2ucl = log10(T2ucl);     % Put in log scale
clin = [linspace(min(T2u),max(T2u),10),...
    T2ucl+(-0.01:0.002:0.01)];
clin = sort(clin);

figure(1); ZZ = reshape(T2u,size(xx));
contourf(xx,yy,ZZ,clin); hold on; colormap(greenmag(T2ucl)); 
c = colorbar; title('T^2 Contour Map'); c.Label.String = 'log(T^2)';