概述
其中OMP算法的步骤如下
本文采用OMP算法来求解稀疏系数。首先随机生成字典数据和待测试数据
字典数据:
1 2 3 4 5 6 |
|
这是一个5*10的矩阵,行数代表维度,列数代表样本数。列数在字典中也叫字典原子,此处有10个原子,原子数大于维数,符合过完备要求。
信号数据:
1 |
|
为了简便,只模拟了一个信号数据,是一个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算法:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 |
|
结果:
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求解稀疏表示,matlab和opencv2所遇到的程序开发问题。
如果觉得靠谱客网站的内容还不错,欢迎将靠谱客网站推荐给程序员好友。
发表评论 取消回复