matlab实现主成分分析方法图像压缩和传输重建

原创风一样的航哥航哥小站 2024年11月12日 15:23 江苏

为了研究图像的渐进式传输技术，前文提到过小波变换，但是发现小波变换非常适合传输缩略图，实现渐进式传输每次传输的数据量不一样，这是因为每次变换之后低频成分大约是上一次的1/4，这样导致重建最小的图之后，继续重建上一层，传输参数会增加为4倍，在带宽有限的情况下，采用分段传输就会让用户等待的时间变得越来越长，这个某些应用场景是不太能接受的，结合压缩算法可以缓解一部分这个问题，但是这不算完美的解决方案。于是我继续研究，发现可以采用主成分分析来完成我想要的功能。

主成分分析（Principal Component Analysis, PCA）是一种常用的数据降维技术，主要用于减少数据集的维度，同时尽可能保留数据的主要特征和结构。PCA通过将原始数据转换到一个新的坐标系中，使得新坐标系的第一个轴（第一主成分）最大化数据的方差，第二个轴（第二主成分）在与第一个轴正交的条件下最大化剩余数据的方差，以此类推。

什么是PCA，对于科普类的知识可以直接请教AI，快速了解大概。

PCA的基本原理：

数据标准化：PCA通常要求数据是标准化的，即每个特征的均值为0，方差为1。这是因为PCA对不同尺度的特征敏感，标准化可以避免某些特征因尺度较大而占据主导地位。
协方差矩阵：计算原始数据的协方差矩阵。协方差矩阵描述了各个特征之间的线性关系和数据的分布情况。
特征值和特征向量：求解协方差矩阵的特征值和特征向量。特征值表示对应特征向量方向上的方差大小，特征向量表示新的坐标轴方向。
选择主成分：根据特征值的大小选择前k个最大的特征值对应的特征向量，这些特征向量就是新的坐标轴，称为主成分。通常选择的主成分数量k使得累计方差贡献率达到一定的阈值（如95%）。
数据投影：将原始数据投影到新的坐标系中，得到降维后的数据。

PCA的应用：

数据可视化：通过将高维数据降到二维或三维，便于在图表中展示和分析。
特征提取：提取数据的主要特征，去除噪声和冗余信息。
数据压缩：减少数据存储和传输的成本。
机器学习：作为预处理步骤，提高模型的训练效率和预测性能。

PCA的优缺点：

优点：

降维：有效减少数据维度，简化模型复杂度。
去噪：去除数据中的噪声，提高数据质量。
可视化：将高维数据降到低维，便于可视化和理解。

缺点：

信息损失：降维过程中会丢失部分信息，特别是在选择较少主成分时。
线性假设：PCA假设数据之间的关系是线性的，对于非线性关系的数据效果不佳。
解释性：降维后的主成分通常难以直接解释其物理意义。

有了上述的基本了解，再查查别的资料，重要的不是手写PCA（当然手写也不是很困难），而是用PCA，先用matlab仿真试一下，编写pcasample函数能够实现基于样本（变量）的主成分分析：

本文处理的原图：

function [coeff,score,rate]=pcasample(X,p)% X:样本矩阵% p：提取前p个主成分% coeff：特征向量矩阵（系数矩阵）% score：得分向量% reta：贡献率% % 将样本归一化% X=zscore(X); %这里好像不用归一化% 计算样本方差的特征向量[V,D]=eig(X'*X);% 将特征向量中的最大值置为正数for i=1:size(V,2)[~,idx]=max(abs(V(:,i))); V(:,i)=V(:,i)*sign(V(idx,i));end% 将特征根按照从大到小的顺序排列[lambda,locs]=sort(diag(D),'descend');V=V(:,locs);% 只提取前p个主成分coeff=V(:,1:p);% 计算得分矩阵score=X*V(:,1:p);% 计算贡献率rate=sum(lambda(1:p))/sum(lambda);end

找个图片运行一下，然后发现主成分个数每次都是1？？？说好的有多个主成分，按照重要性排列呢？怎么只能有一个，那不就变成纯纯的压缩了嘛？没办法，继续研究，在《计算机视觉和深度学习实战》这本书里面，终于明白了怎么回事。

PCA主要是降维，默认输入一张图片就是一整块矩阵，最多算是2维数据。所以要处理图像，生成多个主成分，需要预处理一下。

在一般情况下，数字图像矩阵可以被视为二维数组，为了将图像数组转换为样本矩阵，需要首先对图像进行子块划分，然后将每个子块都拉伸成一维的，最后将所有子块都组合成一个样本矩阵。其中，MATLAB自带的im2col函数可以实现二维数组的分块及向量整合。

继续阅读：“主成分分析（PCA）计算协方差矩阵的特征值和特征向量，并选择少数几个主分量代表多变量的方差（即协方差）结构，是一种有效的特征提取方法。数字图像是二维矩阵，对其通过PCA处理来提取特征，可以在一定比例上保留原始图像的特征信息，并且能够大大减少计算量。因此，PCA图像压缩处理属于一种降维方法，它通过对高维图像块向量空间进行降维处理，将多变量的图像块数据表进行最佳综合、简化，导出少数几个主分量，进而实现在一定比例上保留原始图像信息，又能保持图像块之间的不相关性，进而保证图像压缩的有效性。”

参考例程：

先写一个例程进行PCA分析和重建。

clc;clear;close all;% 读取图像image = imread('pic1.png');%%k=1;for p=1:5:20[Ipca,ratio,contribution]=pcaimage(image,p,[24 24]); subplot(2,2,k);imshow(Ipca)title(['主成分个数=',num2str(p),'压缩比=',num2str(ratio),'贡献率=',num2str(contribution)]);k=k+1;endfunction [Ipca,ratio,contribution] =pcaimage(I,pset,block)%pcaimage 使用主成分实现图像的压缩% 此处提供详细说明% I：进行压缩的图像% pset：主成分个数% Ipca：主成分分析重构图像% ratio：压缩比% contribution：贡献率if nargin<1disp('argument is too few.')endif nargin<2pset=3;endif nargin<3block=[16 16];end% 将彩色图像转换为灰度图if ndims(I)==3I=rgb2gray(I);end% 将图像数组转换为样本矩阵X=im2col(double(I),block,'distinct')';% 样本和变量个数[n,p]=size(X);% 主成分个数不能超过变量个数m=min(pset,p);% 提取前p个主成分，在压缩之后只需要保存coeff和score[coeff,score,contribution]=pcasample(X,m);% 根据系数矩阵重建X=score*coeff';% 将样本矩阵转换为图像数组Ipca=cast(col2im(X',block,size(I),'distinct'),class(I));% 计算压缩比ratio=n*p/(n*m+p*m);endfunction [coeff,score,rate]=pcasample(X,p)% X:样本矩阵% p：提取前p个主成分% coeff：特征向量矩阵（系数矩阵）% score：得分向量% reta：贡献率% % 将样本归一化% X=zscore(X);% 计算样本方差的特征向量[V,D]=eig(X'*X);% 将特征向量中的最大值置为正数for i=1:size(V,2)[~,idx]=max(abs(V(:,i)));V(:,i)=V(:,i)*sign(V(idx,i));end% 将特征根按照从大到小的顺序排列[lambda,locs]=sort(diag(D),'descend');V=V(:,locs);% 只提取前p个主成分coeff=V(:,1:p);% 计算得分矩阵score=X*V(:,1:p);% 计算贡献率rate=sum(lambda(1:p))/sum(lambda);end

得到如下的结果：

看结果，嗯，恢复得还行。

回到渐进式传输的实现，主成分是可以生成多个的，那么可以按照第1主成分、第2主成分、第3主成分……直到传输到够清晰为止，看实际应用，我觉得20已经差不多了。然后写代码来模拟这个过程：

clc;clear;close all;row=4;column=5;% 读取图像image = imread('pic1.png');I=rgb2gray(image);pset=row*column;block=[30 30]; %需要选择合适的参数，让系数最少，初步测试30比较少% 将图像数组转换为样本矩阵X=im2col(double(I),block,'distinct')';% 样本和变量个数[n,p]=size(X);% 主成分个数不能超过变量个数m=min(pset,p);% 提取前p个主成分，在压缩之后只需要保存coeff和score[coeff,score,contribution]=pcasample(X,m);% 根据系数矩阵重建X=score*coeff';% 将样本矩阵转换为图像数组Ipca=cast(col2im(X',block,size(I),'distinct'),class(I));% 计算压缩比ratio=n*p/(n*m+p*m);figuresubplot(1,2,1);imshow(I)title('原图')subplot(1,2,2);imshow(Ipca);title(['主成分个数=',num2str(m),'压缩比=',num2str(ratio),'贡献率=',num2str(contribution)]);figure% 模拟传输过程，每次传输一点，重建一点gScore=[];gcoeff=[];for k=1:psetgScore=[gScore score(:,k)];gcoeff=[gcoeff;coeff(:,k)']; %根据系数矩阵重建X=gScore*gcoeff; %将样本矩阵转换为图像数组Ipca=cast(col2im(X',block,size(I),'distinct'),class(I));subplot(row,column,k);imshow(Ipca)end