OMP求解稀疏表示，matlab和opencv2

其中OMP算法的步骤如下

本文采用OMP算法来求解稀疏系数。首先随机生成字典数据和待测试数据

字典数据：

　　这是一个5*10的矩阵，行数代表维度，列数代表样本数。列数在字典中也叫字典原子，此处有10个原子，原子数大于维数，符合过完备要求。

信号数据：

　　为了简便，只模拟了一个信号数据，是一个5*1的矩阵，如果有多个数据，则应该是5*n的矩阵。求解的时候，可用循环求解。

一、在matlab中实现稀疏表示，求解稀疏系数

clc;close all;clear all;
dic =[
6,
7,
9,
9,
7,
0,
6,
3,
6,
9;
1,
8,
7,
8,
5,
3,
8,
1,
7,
3;
3,
3,
5,
4,
8,
2,
6,
1,
2,
2;
6,
1,
0,
7,
3,
5,
0,
6,
3,
3;
7,
5,
0,
5,
3,
0,
2,
7,
1,
7];
%字典
signal=[
9;
8;
8;
3;
9];
%原始信号
dic=dic*diag(1./sqrt(sum(dic.^2)));
%字典原子单位化，即每列的norm为1
signal=signal/norm(signal);
%信号单位化
[A,res]=OMP(dic,signal,6);
%稀疏度设定为6，即非零元素最多为6个
A
%输出系数
res
%输出残差
epsilon=norm(signal-dic*A)
%验证残差
||Y-Dx||2

其中OMP算法：

　　结果：

A =
0.1450
0.9391
0
0
0.4210
0.1049
0
0
-0.5503
0
res =
3.1402e-16
epsilon =
3.1402e-16

从系数中可以看出，从10个原子共选出了5个原子进行表示，最后的残差非常小，说明稀疏表示的结果和原数据非常接近。

二、在opencv2中实现稀疏表示 

void getData(Mat &data, Mat &signal);
int main(int argc, char* argv[])
{
Mat dic, signal;
getData(dic, signal);
//获取模拟数据
Mat temp(1, dic.cols, CV_32F);
//用一个矩阵保存每个原子的模长
for (int i = 0; i<dic.cols; i++)
{
temp.col(i) = norm(dic.col(i));
//每个原子的模长
}
divide(dic, repeat(temp, dic.rows, 1), dic); //字典原子单位化
signal = signal / norm(signal);
//信号单位化
Mat A=src.OMP(dic, signal, 8);
//调用OPM求解
float res =(float)norm(signal - dic*A); //计算残差
cout << "系数：" <<endl<< A << endl;
cout<<endl<<"残差:"<< endl<<res << endl; //输出残差
waitKey(0);
return 0;
}
void getData(Mat &dic, Mat &signal)
{
dic = (Mat_<float>(5, 10) <<
6, 7, 9, 9, 7, 0, 6, 3, 6, 9,
1, 8, 7, 8, 5, 3, 8, 1, 7, 3,
3, 3, 5, 4, 8, 2, 6, 1, 2, 2,
6, 1, 0, 7, 3, 5, 0, 6, 3, 3,
7, 5, 0, 5, 3, 0, 2, 7, 1, 7);
signal = (Mat_<float>(5, 1) << 9, 8, 8, 3, 9);
}

其中，OMP函数为：

Mat SRC::OMP(Mat& dic, Mat& signal,int sparsity)
{
if (signal.cols>1)
{
cout << "wrong signal" << endl;
exit(-1);
}
vector<int> selectedAtomOrder;
//保存所有选出的字典原子序号
Mat coef(dic.cols, 1, CV_32F, Scalar::all(0)); //需要返回的系数
Mat residual = signal.clone();
//初始化残差
Mat indx(0, 1, CV_32F);//初始化临时系数
Mat phi;
//保存已选出的原子向量
float max_coefficient;
unsigned int atomOrder;
//每次所选择的原子的序号
for (;;)
{
max_coefficient = 0;
//取出内积最大列
for (int i = 0; i <dic.cols; i++)
{
float coefficient = (float)dic.col(i).dot(residual);
if (abs(coefficient) > abs(max_coefficient))
{
max_coefficient = coefficient;
atomOrder = i;
}
}
selectedAtomOrder.push_back(atomOrder); //添加选出的原子序号
Mat& temp_atom = dic.col(atomOrder); //取出该原子
if (phi.cols == 0)
phi = temp_atom;
else
hconcat(phi, temp_atom, phi); //将新原子合并到原子集合中（都是列向量）
indx.push_back(0.0f);
//对系数矩阵新加一项
solve(phi, signal, indx, DECOMP_SVD);
//求解最小二乘问题
residual = signal - phi*indx;
//更新残差
float res_norm = (float)norm(residual);
if (indx.rows >= sparsity || res_norm <= 1e-6) //如果残差小于阈值或达到要求的稀疏度，就返回
{
for (int k = 0; k < selectedAtomOrder.size(); k++)
{
coef.row(selectedAtomOrder[k]).setTo(indx.row(k));
//得到最终的系数
}
return coef;
}
}
}

最终输出结果为：

系数：
[0.14503297;
0.9391216;
0;
0;
0.42096639;
0.1048916;
0;
0;
-0.55029994;
0]
残差:
1.70999e-007

看以看出，opencv得到的系数和matlab得到的系数基本是一样，只是小数点后保留的位数区别。因为小数位数不相同，所以最后残差有点不同，但不影响最终结果，我们只需要系数相同即可。

最后

以上就是沉默芝麻最近收集整理的关于OMP求解稀疏表示，matlab和opencv2的全部内容，更多相关OMP求解稀疏表示内容请搜索靠谱客的其他文章。