降维_机器学习

最新推荐文章于 2024-08-24 14:56:37 发布

明朝百晓生

最新推荐文章于 2024-08-24 14:56:37 发布

阅读量473

点赞数

分类专栏：人工智能

4AM_明朝百晓生

本文链接：https://blog.csdn.net/chengxf2/article/details/101207474

版权

人工智能专栏收录该内容

102 篇文章 3 订阅

订阅专栏

前言：

在很多应用中，数据的维度会很高，一方面计算量巨大，另一方面会面临维度灾难。这里结合菜肴推荐系统

，莺尾花数据集，图像压缩介绍PCA, SVD 来简化数据处理

主要知识背景。

线性代数《相似矩阵及二次型》

《向量组的线性相关性》

PCA
SVD
例子

一 PCA(Principal Component Analysis)

主要原理协方差矩阵的特征分解（为社么要用协方差，原理证明概率论中有专门推导过程）

主要流程：

1：输入数据 m*n的矩阵（m 代表样本个数，n代表维度）

2：矩阵转置得到X(n*m)

3: 矩阵每一行减去当前行的均值

4：计算协方差矩阵

5: 计算协方差矩阵对应的特征向量p，以及特征值

6：排序特征值，取其中最大的K个，组成矩阵P(K*n)

计算协方差（为转置矩阵，不同特征值对应的向量之间正交）

计算特征值特征向量

排序特征值

保留最大的N个特征向量

1.2 向量降维：

样本减去平均值向量

左乘投影矩阵，得到降维后的向量

原理推导：

原数据的协方差矩阵为：

其中X是n*m的矩阵，每一行数据代表一个特征，当前维度的平均值。

因为协方差矩阵是转置矩阵：两个定理

1 其不同特征值对应的特征向量必定正交

2 是S的特征方程K重根，矩阵

从而恰巧有k个线性无关的特征向量

设 P为正交矩阵，线性变换 Y = PX

可以看出来，P就是原空间协方差矩阵对应的特征向量，取其中最大的k个，从大到小，必定是正交基

这章原理很简单，只有熟悉线性代数，很容易

问题：

PCA降维时候，有k重根的时候，如何通过施密特正交化，得到标准正交基

二 SVD 分解

优点：

简化数据，去除噪声，提高算法的结果。同时奇异分解最大的优点，使用Eig函数遇到，奇异矩阵（不可逆矩阵）无法进行特征分解。

PCA 是针对协方差矩阵进行特征分解，SVD是直接对A进行特征分解，主要原理跟PCA基本一致，还是用到对称矩阵的对角化原理。模型如下：

基础知识：

1 : 矩阵和它的转置矩阵特征值相同

2 两个矩阵特征值相同

根据矩阵对角化定理五：存在正交矩阵

其中 P,Q 都是正交矩阵，特征值为A的平方，设A 特征值为q

的特征值为

则：

SVD 分解：

def Test():
    
    A = np.mat([[1,2],
                [1,1],
                [0,0]])
    
    s1 = A*A.T
    s2 = A.T*A
    
    a1,b1 = np.linalg.eig(s1)
    a2,b2 = np.linalg.eig(s2)
    
    svdP, svdLambda ,svd = np.linalg.svd(A)
    print("\n A  A^T 特征值 \n ",a1,"\n A^T   A 特征值: \n ",a2, "\n A: 特征值平方:  \n ", np.power(svdLambda,2)) ##求解特征值


输出：
    
 A  A^T 特征值 
         [6.85410197 0.14589803 0.        ] 
 A^T   A 特征值: 
          [0.14589803 6.85410197] 

 A: 特征值平方:  
           [6.85410197 0.14589803]

SVD 降维后，复原数据例子

def LoadExData():
    
    dataList =[[1,1,1,0,0],
              [2,2,2,0,0],
              [1,1,1,0,0],
              [5,5,5,0,0],
              [1,1,0,2,2],
              [0,0,0,3,3],
              [0,0,0,1,1]]
    
    U, sigMa, Q = np.linalg.svd(dataList)
    
    r = 3
    U1 = U[:,:r] ##取对应的列
    Q1 = Q[:r,:]
    
    sigR = np.zeros((r,r)) 
    for i in range(r):
        sigR[i,i]= sigMa[i]
    

    print("\n sigMa: ",np.round(sigR,2))
    A = np.mat(U1)*sigR*np.mat(Q1)
    print("\n A \n",np.round(A,2))

输出：
   sigMa:  [[9.72 0.   0.  ]
 [0.   5.29 0.  ]
 [0.   0.   0.68]]

 A 
 [[ 1.  1.  1. -0. -0.]
 [ 2.  2.  2.  0.  0.]
 [ 1.  1.  1.  0.  0.]
 [ 5.  5.  5. -0. -0.]
 [ 1.  1. -0.  2.  2.]
 [-0.  0. -0.  3.  3.]
 [-0.  0. -0.  1.  1.]]

  可以看出特征值减少后，依然可以复原原来数据

三例子：

3.1： PCA例子

降维之前

降维之后

# -*- coding: utf-8 -*-
"""
Created on Tue Sep 24 16:27:57 2019

@author: chengxf2
"""

import numpy as np
import sys,os
import pandas as pd
from pandas import Series,DataFrame
from numpy import nan as NaN
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler
from sklearn.datasets import load_iris

"""
绘制散点
Args
    dataMat:  数据集
    dataLabel: 标签集
return
    None
"""
def PlotData(dataMat, dataLabel):
    
    m,n = np.shape(dataMat)
    colorList=['r','g','b']
    for i in range(m):
        x1 = dataMat[i,0]
        x2 = dataMat[i,1]
        co = colorList[dataLabel[i]]
        #print("\n x1: ",x1, "\t x2: ",x2, "c ",co)
        
        plt.scatter(x1,x2, c=co)
    plt.show()

"""
数据标准化
Args
    dataList:  数据集
    x-mean/std

return
    stdData: 标准化后的数据
"""
def Standard(dataList):
   
    
    scaler = StandardScaler()
    scaler.fit(dataList)
    stdData = scaler.transform(dataList)
  
    return stdData


"""
加载数据集
Args
  None
return 
    dataList
"""
def loadData():
    iris = load_iris()
    dataList = iris.data
    print("dataList ",dataList)
    dataLabel = iris.target
    return Standard(dataList),dataLabel

"""
PCA降维
Args
    dataList:  数据集
    numFeat: 降维后的维度

return
    RedMat: 降为后的数据集
"""
def Pca(dataList, numFeat):
    

    m,n = np.shape(dataList)
    #print("\n m: ",m ,"n: ",n, "numFeat: ",numFeat)
    mean = np.mean(dataList, axis=0) ##每一列当作一个数据
   # print("\n mean: ",mean)
    dataMat = dataList- mean #remove mean
    #print("\n 减去均值:  ",dataMat)
   # covMat = np.cov(meanRemoved, rowvar=0)
   
    covMat = np.mat(dataMat).T *np.mat(dataMat)/(m-1)
    #covMatAPI = np.cov(dataMat,rowvar=0)  ##0 每一列代表一个观察变量
   
    #print("\n covMat:   ",  covMat)
    #print("\n covMatAPI:  ",covMatAPI)
    
    
    
    eigVals,eigVects = np.linalg.eig(covMat)
    #print("\n eigVals:  \n ",eigVals)  ##矩阵特征值
  #  print("\n eigVals:  \n ",eigVects)  ##矩阵特征向量
    
    eigSort = np.argsort(eigVals)            #从小到大排序，返回对应的Index
    eigRed = eigSort[-1:-(numFeat+1):-1]    #cut off unwanted dimensions
    

    redEigMat = eigVects[:,eigRed]       #reorganize eig vects largest to smallest ##取对应的列
    RedDataMat = dataMat *redEigMat     #transform data into new dimensionsY =PX
    print("dataMat: ",RedDataMat)
    
    return RedDataMat



dataList,dataLabel = loadData()
RedDataMat= Pca(dataList, 3)
PlotData(RedDataMat, dataLabel)

处理缺失数据：

1：使用可用特征的均值来填补缺失值

2：使用特殊值来填补缺失值 -1

3：忽略有缺失值的样本

4：使用相似样本的均值来填补缺失值

5：使用另外的机器学习算法来预测缺失值

3.2 应用例子点菜系统推荐菜肴。

为用户没有评论过的菜肴打分，推荐

# -*- coding: utf-8 -*-
"""
Created on Fri Sep 27 17:01:11 2019

@author: chengxf2
"""
import numpy as np
from numpy import linalg as lin

"""
欧式距离相似度 
Args
    vecA:向量
    vectB： 向量

return
    向量欧式距离倒数
"""
def eculidSim(vecA, vecB):
    
     sim = lin.norm(vecA- vecB)
     
     return 1.0/(1.0+sim)
 

"""
皮尔逊相关系数 相似度
Args
    vecA: 向量A
    vecB: 向量B

return
    协方差，皮尔逊系数
"""
def pearSim(vecA, vecB):
    
    if len(vecA)<3:
        return 1.0
    
    sim = 0.5+0.5*np.corrcoef(vecA, vecB, rowvar =0)[0][1]
    
    return sim


"""
余弦相似度 
Args
    vecA: 列向量 n 行一列
    vecB: 列向量
return
    sim: 相似度
"""
def cosSim(vecA, vecB):
    
   # print("vecA ",vecA)
    m,n = np.shape(vecA)
    
    if n !=1:
        print("error")
        return None
    
    num = vecA.T*vecB
    #print("\n num ",num)
    
    demon = lin.norm(vecA)*lin.norm(vecB)
    
    sim = 0.5+0.5*(num/demon)
    return sim[0,0]
    
"""
找出两列中元素大于0的行，都评分过的行
Args
  dataMat: 矩阵
  colI: I列
  colJ: J列

return
   大于0的列
"""
def GetSubItem(dataMat, colI, colJ):
    
    m,n = np.shape(dataMat)
   # print("m: ",m)
    
    ret =[]
    for i in range(m):
        if dataMat[i, colI]>0 and dataMat[i,colJ]>0:
            ret.append(i)
    return ret
    
"""
给定相似度计算方法条件下，用户对产品评分
Args
   dataMat: 数据矩阵
   user: 用户编号，第几行
   simMeas: 相似度度量
   item: 某一列属性
return
    取相似度平均值
"""
def standEst(dataMat, user, simMeas, item):
   # print("\n ****item****:  \n  ",item)
    n = np.shape(dataMat)[1]
    simTotal = 0.0; ratSimTotal = 0.0
    for j in range(n):
        userRating = dataMat[user,j] ##用户打分
        if userRating == 0: 
            continue
        
        
        overLapK = np.nonzero(np.logical_and(dataMat[:,item].A>0, \
                                      dataMat[:,j].A>0))[0]
        #print("\n j:  ",j)
        overLap = GetSubItem(dataMat, j, item)
       # print("\n j: ",j, "\t  overLap",overLap, "overLapK ",overLapK)
        if len(overLap) == 0: ##没有重叠的数据相似度为0
               similarity = 0
        else: ##计算相似度
            similarity = simMeas(dataMat[overLap,item], dataMat[overLap,j])
        #print 'the %d and %d similarity is: %f' % (item, j, similarity)
        simTotal += similarity
        ratSimTotal += similarity * userRating
    if simTotal == 0: 
        return 0
    else: 
        return ratSimTotal/simTotal
    

"""
推荐系统
Args
   dataMat: 数据集
   user: 第几行
   N: 
   simMeas: 相似度函数
   estMethod: 给定相似度计算方法条件下，用户对产品评分
"""
def recommend(dataMat, user, N=3, simMeas=cosSim, estMethod= standEst):
    #  [0, 0, 0, 0, 4, 0, 0, 1, 0, 4, 0]
    unrateditem = np.nonzero(dataMat[user,:].A==0)[1] ##返回不为0的索引，已经评分过了的菜
    
   
    if 0 == len(unrateditem):
        return "you rated everything"
    
    itemScores =[]
    for uncol in unrateditem:
       # print("\n item ",item)
        estScore = estMethod(dataMat, user, simMeas, uncol)
        itemScores.append((uncol, estScore))
    
    rec = sorted(itemScores, key = lambda j:j[1], reverse=True)[:N] #寻找前N个未评价物品
    
    print("rec: ",rec)
    return rec
        
        
    
    
"""
加载数据集
Args
   m: 数据集个数
   n: 维度
"""    
def loadExData2():
    dataList = [[0, 0, 0, 0, 0, 4, 0, 0, 0, 0, 5],
           [0, 0, 0, 3, 0, 4, 0, 0, 0, 0, 3],
           [0, 0, 0, 0, 4, 0, 0, 1, 0, 4, 0],
           [3, 3, 4, 0, 0, 0, 0, 2, 2, 0, 0],
           [5, 4, 5, 0, 0, 0, 0, 5, 5, 0, 0],
           [0, 0, 0, 0, 5, 0, 1, 0, 0, 5, 0],
           [4, 3, 4, 0, 0, 0, 0, 5, 5, 0, 1],
           [0, 0, 0, 4, 0, 4, 0, 0, 0, 0, 4],
           [0, 0, 0, 2, 0, 2, 5, 0, 0, 1, 2],
           [0, 0, 0, 0, 5, 0, 0, 0, 0, 4, 0],
           [1, 0, 0, 0, 0, 0, 0, 1, 2, 0, 0]]
    
  
    
    return np.mat(dataList)


"""
保留能量的百分百
行列式大小就是特征值乘积
Args
   sigma: 特征值
   percentage 百分比
"""
def cal_eig_k(sigma, percentage):
    sigma2 = sigma ** 2  # 对sigma求平方
    sumsgm2 = sum(sigma2)  # 求所有奇异值sigma的平方和
    sumsgm3 = 0  # sumsgm3是前k个奇异值的平方和
    k = 0
    for i in sigma:
        sumsgm3 += i ** 2
        k += 1
        if sumsgm3 >= sumsgm2 * percentage:
            return k

"""
dataList =[[2,  0,  0,  4,  4],
           [5,  5,  5,  3,  3],
           [2,  4,  2,  1,  1]]
"""

"""
奇异分解评估
Args
   dataMat: 数据集
   user: 对应的输入数据
   simMeas: 相似度度量函数
   item: 比较项目

"""
def svdEst(dataMat, user, simMeas, item):
    m,n = np.shape(dataMat)
   # print("\n m: ",m , "\t n: ",n)
    simTotal = 0
    ratSimTotal = 0
    u,eig,v = lin.svd(dataMat)
    k = cal_eig_k(eig, 0.9)
    #print("\n k:",k)

   # print("\n sigMa ", np.round(sigma,2))
    eigK =np.mat(np.eye(k)*eig[:k])
    
   # print("\n sig4 ",sig4.I)
   # vRed = dataMat.T * u[:,:k] * eigK.I 
    #vRed1 = np.round(dataMat.T * u[:,:k] * eigK.I,4)
    vRed = v[:k,:].T ##A^TA 的特征向量,降维
    
    
   # if False ==bSame:
        #print("vRed1  ",vRed1, " \n vRed ",vRed)
    
  
   # print("\n *******n********")
    #print("\n 计算出来的 \n ",np.round(vRed,2))
    #print("\n SVD 直接分解得到的 \n  ",np.round(vRed1,2))
    
    for j in range(n):
        userRate = dataMat[user, j]  #当前某个属性用户的评分
        if userRate ==0 or j == item:  ##没有评分或者是当前列
            continue

        s1 = simMeas(vRed[item,:].T, vRed[j,:].T)
        #print("\n s1: ",np.round(s1, 2), "sp ",np.shape(vRed[item,:]))
        simTotal += s1
        ratSimTotal += s1*userRate
    if 0== simTotal:
        return 0
    else:
        return ratSimTotal/simTotal
        
    
   
    
    
    

        
def Test(dataMat):
    
   
    
    #overLap =[ 3  ,4 , 6, 10]
    
   # print("\n == standEst 使用原始数据评估 == \n")
    #recommend(dataMat, 1,N=3)#pearSimcosSim
    print("\n ==== Recommd 使用奇异分解后的评估=== \n")
    recommend(dataMat, 2,N=3,estMethod=svdEst)
 
    
   
    
dataMat = loadExData2()
euc =Test(dataMat)

3.3 SVD 图像压缩

# -*- coding: utf-8 -*-
"""
Created on Sun Sep 29 16:30:14 2019

@author: chengxf2
"""

import numpy as np
import matplotlib.image as Img
import matplotlib.pyplot as plt
import matplotlib as mpl
from sklearn.preprocessing import StandardScaler
from PIL import Image


        
"""
根据比例，取k个特征值
Args
    rat: 百分比
    eig: 特征值
"""
def cal_eig_k(rat, eig):
    
    
    #print("\n rat ", rat)
    data = np.reshape(eig,(-1,1))
    
    scaler = StandardScaler()
    scaler.fit(data)
    tran_data = scaler.transform(data)
    
   # print("\n tran_data: \n ",np.shape(tran_data))
    eigNorm =  np.power(tran_data, 2)
    totalEng = np.sum(eigNorm)*rat
    
  
   # print("\n totalEng:  ",totalEng, "\t tran_data: ",tran_data[0:10])
    curEng = 0
    k = 0

    for t in eigNorm:
        k =k+1
        curEng += t
        if curEng> totalEng:
            return k
        
    #print("\n k: ",k)
   
        

"""
奇异分解图片
Args
    n:
    pic: 图片
return
    None
"""
def Image_svd(n, pic):
    
    
    
    a, b, c = np.linalg.svd(pic)
    
    svd = np.zeros((a.shape[0],c.shape[1]))
    for i in range(0, n):
        svd[i, i] = b[i]
    img = np.matmul(a, svd)
    img = np.matmul(img, c)
    img[ img >= 255] = 255
    img[  0 >= img ] = 0
    img = img.astype(np.uint8)
    return img


"""
奇异分解
Args:
    rat: 占比
    r: 对应通道的r,g,b值
    eig: (700,)
"""
def SvdImg(rat, r):
    
    
    U,Eig,V =np.linalg.svd(r)
    m,n = np.shape(U)
    m1,n1 = np.shape(V)
  
    
    #print("\n  sp eig :\t", np.shape(eig))
    k = cal_eig_k(rat, Eig)
    #print("\n**** k:*** \t ",k)
    
    svdEig = np.zeros((k,k))
    
    for i in range(k):
        svdEig[i,i] = Eig[i]
    
    UK = U[:,:k]
    VK = V[:k,:]
    img = np.matmul(UK, svdEig)
    img = np.matmul(img, VK)
    
    img[ img >= 255] = 255
    img[  0 >= img ] = 0
    img = img.astype(np.uint8)
    return img

        
    
    
"""
加载图片,SVD后显示
Args
   path:  图片路径
   img: (980, 700, 3), numpy.ndarray'
   r: [980,700], type: Image.Image
return
    None
"""
def LoadImag():
    path = "D:/0.png"
    
    
    im = Image.open(path)
    r,g,b = im.split()

    
    plt.figure(figsize=(30,30))
    engRatio =[ 0.97, 0.99,0.995]  ##5个元素
    
   
    
    for eng in engRatio:
        
        r_img = SvdImg(eng,r)
        g_img = SvdImg(eng,g)
        b_img = SvdImg(eng,b)
        
        rn = Image.fromarray(r_img)
        gn = Image.fromarray(g_img)
        bn = Image.fromarray(b_img)
        
        #print("\n t: ",type(Image.fromarray(np.uint8(r_img)), "\t tr ",type(r), "\t sp ",np.shape(r_img ),"\t sp2",np.shape(r))
        print("sp ",type(rn))
        im_merge = Image.merge("RGB",[rn,gn,bn])
        Tip  = "能量占比%s:"%str(eng)
        im_merge.show(title=Tip)
        
        
       




load = LoadImag()

不同k值的效果图