过程监控中会用到很多中方法，如主成分分析（PCA）、慢特征分析（SFA）、概率MVA方法或独立成分分析（ICA）等为主流算法。
其中PCA主要多用于降维及特征提取，且只对正太分布（高斯分布）数据样本有效；SFA被用来学习过程监控的时间相关表示，SFA不仅可以通过监测稳态分布来检测与运行条件的偏差，还可以根据时间分布来识别过程的动态异常，多用于分类分析；概率MVA方法，多以解决动力学、时变、非线性等问题。
今天要介绍的是独立成分分析（ICA），由浅入深，细细道来。
此外文末还附有ICA可实现的代码哟~不要错过
独立成分分析（Independent Component Analysis，ICA）
基本原理
在信号处理中，独立成分分析（ICA）是一种用于将多元信号分离为加性子分量的计算方法。这是通过假设子分量是非高斯信号，并且在统计上彼此独立来完成的。ICA是盲源分离的特例。一个常见的示例应用程序是在嘈杂的房间中聆听一个人的语音的“ 鸡尾酒会问题 ”。
首先，引入一下经典的鸡尾酒宴会问题（Cocktail Party Problem）。
假设我们令一个未知的混合系数矩阵（mixing coefficient matrix）为A ，用来组合叠加信号S ，
ICA的不确定性（ICA ambiguities）
由于h 和s 都不确定，那么在没有先验知识的情况下，无法同时确定这两个相关参数。
ICA 算法
下面直接上ICA算法。
独立成分分析 ICA(Independent Component Correlation Algorithm)是一种函数，X为n维观测信号矢量，S为独立的m（m<=n)维未知源信号矢量，矩阵A被称为混合矩阵。ICA的目的就是寻找解混矩阵W（A的逆矩阵），然后对X进行线性变换，得到输出向量U。
此公式代表一个假设前提：每个人发出的声音信号各自独立。
求导可得，
举个Paper的栗子
下面为我们观测到的信号：
然后，再通过ICA还原后的信号为：
MATLAB代码实现
MATLAB代码：
Fast ICA
% Input：X 行变量维数，列采样个数；需要对原始矩阵转置 % Output：Sources重构的原信号, Q白化矩阵, P白化信号解混矩阵 function [Sources, Q, P] = FastICA(X, P) % 白化处理 [dim, numSample] = size(X); Xcov = cov(X'); [U, lambda] = eig(Xcov); Q = lambda^(-1/2)*U'; Z = Q*X; % FastICA maxiteration = 10000; %最大迭代次数 error = 1e-5; % 收敛误差 % P = randn(dim,dim); % 随机初始化P，并按照列更新 for k = 1:dim Pk = P(:,k); Pk = Pk./norm(Pk); % 向量归一化 lastPk = zeros(dim,1); % 0不需要再归一化 count = 0; while abs(Pk - lastPk)&abs(Pk + lastPk) > error count = count + 1; lastPk = Pk; g = tanh(lastPk'*Z); % g(y)函数 dg = 1 - g.^2; % g(y)的一阶导函数 %-------------------------------核心公式------------------------------------ Pk = mean(Z.*repmat(g,dim,1), 2) - repmat(mean(dg),dim,1).*lastPk; Pk = Pk - sum(repmat(Pk'*P(:,1:k-1),dim,1).*P(:,1:k-1),2); Pk = Pk./norm(Pk); %-------------------------------------------------------------------------- if count == maxiteration fprintf('第%d个分量在%d次迭代内不收敛！\n',k,maxiteration); break; end end P(:,k) = Pk; end Sources = P'*Z; % end
此外还有基于故障诊断的ICA算法代码实现->在这里
下面给出部分代码：
% 基于ICA的故障诊断 clear;clc;close all; load('MPD2000.mat'); Xnormal = MPD0'; % 数据归一化 [dim, numSample] = size(Xnormal); XnormalMean = mean(Xnormal, 2); XnormalStd = std(Xnormal, 0, 2); XnormalNorm = normalization(Xnormal, XnormalMean, XnormalStd); % 正常数据计算 解混矩阵W % P = rand(dim,dim)*100; load('P.mat'); [S, Q, P] = FastICA(XnormalNorm, P); W = P'*Q; % ------------------------利用2范数大小对W的行重新排列--------------------- Wnorm = zeros(dim,1); for k = 1:dim Wnorm(k) = norm(W(k,:)); end [Wnorm, indices] = sort(Wnorm, 'descend'); % -------------------------确定主导成分Sd与参与成分Se---------------------- threshold = 0.80; percentage = cumsum(Wnorm)./sum(Wnorm); for k = 1:dim if percentage(k) > threshold break; end end
效果如下：
