K-SVD字典学习算法

1.提出问题:什么是稀疏表示
假设我们用一个MN的矩阵表示数据集Y,每一行代表一个样本,每一列代表样本的一个属性,一般而言,该矩阵是稠密的,即大多数元素不为0。
稀疏表示的含义是,寻找一个系数矩阵X(K
N)以及一个字典矩阵D(MK),使得DX尽可能的还原Y,且X尽可能的稀疏。X便是Y的稀疏表示。

算法思想

算法求解思路为交替迭代的进行稀疏编码和字典更新两个步骤. K-SVD在构建字典步骤中,K-SVD不仅仅将原子依次更新,对于原子对应的稀疏矩阵中行向量也依次进行了修正. 不像MOP,K-SVD不需要对矩阵求逆,而是利用SVD数学分析方法得到了一个新的原子和修正的系数向量.

固定系数矩阵X和字典矩阵D,字典的第k个原子为dk,同时dk对应的稀固定系数矩阵X和字典矩阵D,字典的第k个原子为d_k,同时d_k对应的稀XDkdkdk
疏矩阵为X中的第k个行向量xTk.假设当前更新进行到原子dk,样本矩阵和疏矩阵为X中的第k个行向量x^k_T. 假设当前更新进行到原子d_k,样本矩阵和XkxTk.dk字典逼近的误差为:
在这里插入图片描述
在得到当前误差矩阵Ek后,需要调整dk和XTk,使其乘积与Ek的误差尽可能的小.在得到当前误差矩阵E_k后,需要调整d_k和X^k_T,使其乘积与E_k的误差尽可能的小.EkdkXTk使Ek.
如果直接对dk和XTk进行更新,可能导致xTk不稀疏.所以可以先把原有向量xTk中零如果直接对d_k和X^k_T进行更新,可能导致x^k_T不稀疏. 所以可以先把原有向量x^k_T中零dkXTkxTk.xTk
元素去除,保留非零项,构成向量xRk,然后从误差矩阵Ek中取出相应的列向量,元素去除,保留非零项,构成向量x^k_R,然后从误差矩阵E_k中取出相应的列向量,xRkEk
构成矩阵EkR.对EkR进行SVD(SingularValueDecomposition)分解,有EkR=构成矩阵E^R_k. 对E^R_k进行SVD(Singular Value Decomposition)分解,有E^R_k=EkR.EkRSVDSingularValueDecompositionEkR=
UΔVT,由U的第一列更新dk,由V的第一列乘以Δ(1,1)所得结果更新xRk.UΔV^T,由U的第一列更新d_k,由V的第一列乘以Δ(1,1)所得结果更新x^k_R.UΔVTUdkVΔ(1,1)xRk.

代码

import numpy as np
import matplotlib.pyplot as plt
from numpy import *
import scipy.io as sio
import random
from sklearn import linear_model
import scipy.misc
from PIL import Image
def esErrDic(data,recons):
    m,n=data.shape
    esErr=0
    for i in range(m):
        for j in range(n):
            esErr+=(data[i][j]-recons[i][j])**2
    return esErr/(m*n) 
class KSVD(object):
    def __init__(self, n_components, max_iter=100, tol=1e-6,n_nonzero_coefs=None):
        """
        稀疏模型Y = DX,Y为样本矩阵,使用KSVD动态更新字典矩阵D和稀疏矩阵X
        :param n_components: 字典所含原子个数(字典的列数)
        :param max_iter: 最大迭代次数
        :param tol: 稀疏表示结果的容差
        :param n_nonzero_coefs: 稀疏度
        """
        self.dictionary = None
        self.sparsecode = None 
        self.max_iter = max_iter 
        self.tol = tol 
        self.n_components = n_components 
        self.n_nonzero_coefs = n_nonzero_coefs
    def _initialize(self, y):
        """
        初始化字典矩阵
        """
        u, s, v = np.linalg.svd(y)
        self.dictionary = u[:, :self.n_components]
    def _update_dict(self, y, d, x):
        """
        使用KSVD更新字典的过程
        """
        for i in range(self.n_components):
            index = np.nonzero(x[i, :])[0]#选出Xk中非零的元素下标
            if len(index) == 0:
                continue
            d[:, i] = 0
            r = (y - np.dot(d, x))[:, index] 
            u, s, v = np.linalg.svd(r, full_matrices=False) 
            d[:, i] = u[:, 0].T 
            x[i, index] = s[0] * v[0, :]
        return d, x
    def fit(self, y):
        """
        KSVD迭代过程
        """
        self._initialize(y)
        for i in range(self.max_iter):
            x = linear_model.orthogonal_mp(self.dictionary, y, n_nonzero_coefs=self.n_nonzero_coefs)
            e = np.linalg.norm(y - np.dot(self.dictionary, x)) 
            if e < self.tol:
                break
            self._update_dict(y, self.dictionary, x)
        self.sparsecode = linear_model.orthogonal_mp(self.dictionary, y, n_nonzero_coefs=self.n_nonzero_coefs) 
        return self.dictionary, self.sparsecode
if __name__ == '__main__':
    file='G:/lecture of grade one/pattern recognition/trial_two/train_data2_807802844.mat'
    Img=loadData(file)
    #字典学习部分代码,其中KSVD的参数就是原子个数,可修改
    ksvd = KSVD(100)
    dictionary, sparsecode = ksvd.fit(Img['Data'])
    recons=dictionary.dot(sparsecode)
    err=esErrDic(Img['Data'],recons)
	'''
    #测试KSVD代码的调试代码,选择一张图片,用自己编写的KSVD对其进行字典学习
    image =Image.open('/home/swh/Downloads/scene.jpeg')
    image = np.array(image)
    image=image[:,:,0]
    im_ascent = image.astype('float32')
    #im_ascent = scipy.misc.ascent().astype(np.float)
    ksvd = KSVD(30) 
    dictionary, sparsecode = ksvd.fit(im_ascent) 
    plt.figure() 
    plt.subplot(1, 2, 1) 
    plt.imshow(im_ascent) 
    plt.subplot(1, 2, 2) 
    recon=dictionary.dot(sparsecode)
    plt.imshow(recon) 
    plt.show()
    '''
        
评论 3
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值