机器学习实战教程(四):从特征分解到协方差矩阵:详细剖析和实现PCA算法

news2024/7/4 5:40:50

1. 协方差

概念

方差和标准差的原理和实例演示,请参考

方差

方差(Variance)是度量一组数据的分散程度。方差是各个样本与样本均值的差的平方和的均值:
在这里插入图片描述

标准差

标准差是数值分散的测量。
标准差的符号是 σ (希腊语字母 西格马,英语 sigma)
公式很简单:方差的平方根。
在这里插入图片描述

协方差

通俗理解

可以通俗的理解为:两个变量在变化过程中是同方向变化?还是反方向变化?同向或反向程度如何?
你变大,同时我也变大,说明两个变量是同向变化的,这时协方差就是正的。
你变大,同时我变小,说明两个变量是反向变化的,这时协方差就是负的。
从数值来看,协方差的数值越大,两个变量同向程度也就越大。反之亦然。
通俗易懂的理解看知乎文章 或者 gitlab转载

协方差矩阵

协方差(Covariance)在概率论和统计学中用于衡量两个变量的总体误差。而方差是协方差的一种特殊情况,即当两个变量是相同的情况。 这个解释摘自维基百科,看起来很是抽象,不好理解。其实简单来讲,协方差就是衡量两个变量相关性的变量。当协方差为正时,两个变量呈正相关关系(同增同减);当协方差为负时,两个变量呈负相关关系(一增一减)。而协方差矩阵,只是将所有变量的协方差关系用矩阵的形式表现出来而已。通过矩阵这一工具,可以更方便地进行数学运算。
数学定义
回想概率统计里面关于方差的数学定义:
在这里插入图片描述
协方差的数学定义异曲同工:
在这里插入图片描述
这里的 x和y表示两个变量空间。用机器学习的话讲,就是样本有 x和 y两种特征,
而 X 就是包含所有样本的 x特征的集合,
Y就是包含所有样本的 y特征的集合。
用一个例子来解释会更加形象。
在这里插入图片描述
用一个矩阵表示为:
在这里插入图片描述
现在,我们用两个变量空间X ,Y 来表示这两个特征:
在这里插入图片描述
由于协方差反应的是两个变量之间的相关性,因此,协方差矩阵表示的是所有变量之间两两相关的关系,具体来讲,一个包含两个特征的矩阵,其协方差矩阵应该有 2*2 大小:
在这里插入图片描述
接下来,就来逐一计算 Cov(Z)的值。 首先,我们需要先计算出 X,Y 两个特征空间的平均值:
AVG(X)=3.25,AVG(Y)=3 , 。 然后,根据协方差的数学定义,计算协方差矩阵的每个元素:
在这里插入图片描述
所以协方差矩阵
在这里插入图片描述
好了,虽然这只是一个二维特征的例子,但我们已经可以从中总结出协方差矩阵
的「计算套路」:
在这里插入图片描述
python协方差原理

# 协方差主要是多个特征
pa=np.array([
     [1,2] ,
     [3,6] ,
     [4,2] ,
     [5,2] 
   ])
'''
 x应该是个二维矩阵表示
 [
     [1] ,
     [3] ,
     [4] ,
     [5] 
   ]
'''
x=np.array([pa[:,0]]).reshape((4,1))
'''
 y应该是个二维矩阵表示
 [
     [2] ,
     [6] ,
     [2] ,
     [2] 
   ]
'''
y=np.array([pa[:,1]]).reshape((4,1))
print("分别获取X和Y:",x,y)
x_mean=np.mean(x)
y_mean=np.mean(y)
print("x和y特征的均值",x_mean,y_mean)
'''
这里只是求第一个特征x(第一列)和第二个特征(第二列)的方差cov(x,y),
,实际还有cov(x,x),cov(y,x),cov(y,y)
x-x_mean转置T就变成了
 [
     [1-xmean,3-xmean,4-xmean,5-xmean] 
 ]
y-ymean是
 [
     [2-ymean] ,
     [6-ymean] ,
     [2-ymean] ,
     [2-ymean] 
 ]
(x-x_mean).T.dot(y-y_mean)就变成了矩阵乘法了

(1-xmean)*(2-ymean)+(3-xmean)*(6-ymean)+(4-xmean)*(2-ymean)+(5-xmean)*(2-ymean)
然后除以n-1就是协方差了cov(x,y)
'''
print("cov(x,y)=",np.sum((x-x_mean).T.dot(y-y_mean))/(len(pa)-1))

#数学表示横表示列,竖表示行,默认行横表示特征
'''
[[1 3 4 5]
 [2 6 2 2]]
'''
print(pa.T)
print("conv",np.cov(pa.T))
# 使用rowvar=False,表示列是变量是特征
print("conv",np.cov(pa,rowvar=False))

输出结果可查看:github

2. PCA

概念

主成分分析(Principal Component Analysis):

  • 一个非监督的机器学习算法
  • 主要用于数据的降维
  • 通过降维,可以发现更便于人类理解的特征
  • 其他应用:可视化,降噪

假设存在一根直线,将所有的点都映射在该条指直线上,这样的话点的整体分布和原来的点的分布就没有很大的差异(点和点的距离比映射到x轴或者映射到y轴都要大,区分度就更加明显),与此同时所有的点都在一个轴上(理解成一个维度),虽然这个轴是斜着的。用这种方式将二维降到了一维度

公式推导

那么如何找到这个让样本间距最大的轴?

如何定义样本间间距? 使用方差(Variance)
在这里插入图片描述
方差越大代表样本之间越稀疏,方差越小代表样本之间越紧密。
移动坐标轴,使得样本在每一个维度均值都为0:

创建一个3x+4线性附近20个随机样本,样例和结果:https://github.com/lzeqian/machinelearn/blob/master/sklean_pca/demean.ipynb

import numpy as np;
import matplotlib.pyplot as plot
from commons.common import setXY
#生成一个 3x+4附近的点
np.random.seed(100);
# 获取randc个随机点
randc=20
x1=np.random.rand(randc);
x2=x1*3+4+np.random.rand(randc);
plot.plot(x1,x2,"o");

在这里插入图片描述
转换为矩阵表示

#使用矩阵数组表示,x1,x2  x1和x2是两个特征
'''
[   [x1,x2]
    [1,2],
    [3,4]
]
X=x1.reshape(-1,1); 等价于x1.reshape(len(x1),1)
[1,2]转换为
x1=[
  [1],
  [2]
]
x2=[
  [4],
  [5]
]
x1.hstack(x2)
[
  [1,4],
  [2,5]
]
'''
X=np.hstack((x1.reshape(len(x1),1),x2.reshape(len(x2),1)))
print(X)
[[0.54340494 6.06191901]
 [0.27836939 5.77513797]
 [0.42451759 6.09120215]
 [0.84477613 6.87044035]
 [0.00471886 4.18956702]
 [0.12156912 4.73753941]
 [0.67074908 6.01793576]
 [0.82585276 6.72998462]
 [0.13670659 5.20578228]
 [0.57509333 5.74053496]
 [0.89132195 7.27280924]
 [0.20920212 5.23141091]
 [0.18532822 4.66113234]
 [0.10837689 4.70707412]
 [0.21969749 4.69556853]
 [0.97862378 7.82628292]
 [0.81168315 7.4159703 ]
 [0.17194101 4.57576503]
 [0.81622475 7.33922019]
 [0.27407375 5.39912274]]

demean,均值归0处理

np.set_printoptions(suppress=True) #不使用科学计数法
setXY()
def demean(X):
    return X-np.mean(X,axis=0) #取对应列的均值
#均值归零的算法是x1-xmean,x2-x2mean
X_demean=demean(X)
plot.plot(X_demean[:,0],X_demean[:,1],"o");

在这里插入图片描述
demean之后的方差最大其实就是求映射后每个点到(0,0)的距离最大再求和,假设降维后轴的方向是w=(w1, w2)
在这里插入图片描述
,Xi是映射前的向量,Xi(project)是映射后的向量,这里注意w向量是单位向量 |w|=1
在这里插入图片描述
以上是推导公式
目标函数是:
在这里插入图片描述
通过公司可知
在这里插入图片描述

注意i是样本索引,下标1,2是特征数,x1是特征1,x2是特征2 xj是特征j,每个特征xj对应一个维度wj

目标函数即:
在这里插入图片描述
对目标函数求梯度:
在这里插入图片描述
转化为:
在这里插入图片描述
由于最终转换的结果是一个1行m列的矩阵,而我们想要得到一个n行1列的矩阵,所以还要进行一次转置

在这里插入图片描述

编程实现(第一主成分)

产生一个 3x+4附近的点

import numpy as np;
import matplotlib.pyplot as plot
from commons.common import setXY
#生成一个 3x+4附近的点
#np.random.seed(100);
# 获取randc个随机点
randc=100
# 0-1的数*多少倍,注意太小的样本点,abs(f(w=w, X=X) - f(w=last_w, X=X)) 的值增量过小可能导致循环次数后,还没有到<epsilon导致拟合不准确
# 如果是1,建议循环次数加大100000
# 如果是100 可以设置为100
blow=100
x1=np.random.rand(randc)*blow;
x2=x1*3+4+np.random.rand(randc)*blow;
plot.plot(x1,x2,"o");
X=np.hstack((x1.reshape(len(x1),1),x2.reshape(len(x2),1)))
np.set_printoptions(suppress=True) #不使用科学计数法

定义目标函数

# 这是目标函数 np.sum((X*W)**2)/M
# 注意目标函数是传入X是已知的数据样本,w是个2个特征向量 ,f(w1,w2)是个三位空间

def f(X,w):
    return np.sum((X.dot(w))**2)/len(X)

获取w在各个特征的导数

# 获取各个维度的导数    
def df_w(X,w):
    return X.T.dot(X.dot(w))*2/len(X)

也可以用这个通用的方法求导数

'''
通用计算某个点的斜率的方法
为了验证我们的这个是正确的,使用这个df_debug这个函数,
和线性下降法一样,使两个点之间连成的直线不断的靠近应得的直线,
使其斜率相当,注意的是,这里的epsilon取值比较小,是因为在PCA的梯度上升法中,
w是一个方向向量,其模为1,所以w的每一个维度其实都很小,那么为了适应,相应的epsilon也要小一些
'''
def df_debug(w,X,epsilon=0.0001):
      res = np.empty(len(w))
      for i in range(len(w)):
          w_1 = w.copy()
          w_1[i] += epsilon
          w_2 = w.copy()
          w_2[i] -= epsilon
          res[i] = (f(w=w_1,X=X) - f(w=w_2,X=X)) / (2*epsilon)

将w向量转换为单位向量

# 将任意向量转换为单位向量 np.linalg.norm(w)是 x**2+x1**2开根号
# (3,4)/5=(3/5,4/5)就是单位向量,模是1
def direction(w):
    return w / np.linalg.norm(w)

梯度上升测试

wapp=np.array([])
def gradient_ascent(df, X, initial_w, eta, n_iters=1e4, epsilon=1e-8):
    w = direction(initial_w)
    i_iter = 1
    global wapp
    wapp=np.append(wapp,w).reshape((1,len(w)))
    while i_iter < n_iters:
        gradient = df(w=w, X=X)
        last_w = w
        # gradient是对每个维度求偏导得到的列表,如果偏导数为负则w的这个维度加上一个负值,降维后的方差趋于变大
        # 如果偏导数为正,则w的这个维度加上一个正值,降维后的方差趋于变大,因此w加上导数值,降维后的方差趋于变大
        # 在eta合适的情况下,随着循环进行,导数值逐渐趋近0,eta是常数,降维后的方差的变化量会越来越小
        w = w + eta * gradient
        w = direction(w)  # 注意1,每次求一个单位向量
        wapp=np.vstack((wapp,np.array([w])))
        # abs求绝对值
        if (abs(f(w=w, X=X) - f(w=last_w, X=X)) < epsilon):
            print("精度:",abs(f(w=w, X=X) - f(w=last_w, X=X)),epsilon)
            print("梯度",gradient)
            
            break
        i_iter += 1
    return w


initial_w = np.random.random(X.shape[1])  # 注意2:不能用0向量开始
eta = 0.0001
# print(gradient_ascent(df_debug, X=X_demean, initial_w=initial_w, eta=eta))
w=(gradient_ascent(df_w, X=X_demean, initial_w=initial_w, eta=eta, n_iters=100))
setXY()
plot.plot(X_demean[:,0],X_demean[:,1],"o");
print(w)
# 单位向量乘同一个,方向是相同的
plot.plot([0,w[0]*(blow)],[0,w[1]*blow])
#plot.plot(X_demean[:,0],w[1]/w[0]*X_demean[:,0])
#%%
fig = plot.figure()
#创建梯度上升的过程
ax = plot.axes(projection='3d')
wappy=[f(w=w_t, X=X) for w_t in wapp]
ax.plot3D(wapp[:,0],wapp[:,1],wappy, 'red')
print("w1值",wapp[:,0])
print("w2值",wapp[:,1])
print("方差:",wappy)

ax.set_title('3D line plot')
plot.show()

拟合的直线
在这里插入图片描述
梯度上升过程
在这里插入图片描述

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/196264.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

【redis6】第十五章(应用问题解决)

缓存穿透 问题描述 key对应的数据在数据源并不存在&#xff0c;每次针对此key的请求从缓存获取不到&#xff0c;请求都会压到数据源&#xff0c;从而可能压垮数据源。比如用一个不存在的用户id获取用户信息&#xff0c;不论缓存还是数据库都没有&#xff0c;若黑客利用此漏洞…

Yolov5环境部署步骤记录

目录1.Anaconda1.1 虚拟环境1.1.1 虚拟环境手动配置Pytorch库2.Pycharm社区版2.1 Yolov5源码下载2.2 Pycharm设置3. Yolov53.1 安装所需的插件3.2 运行detect.py1.Anaconda 安装&#xff0c;Anaconda3-5.3.1-Windows-X86_64.exe&#xff0c;装好之后&#xff1b; 配置环境变量…

Code:美团代码托管平台的演进与实践

美团代码托管平台经过长期的打磨&#xff0c;完成了分布式架构的改造落地&#xff0c;托管数以万计的仓库&#xff0c;日均Git相关请求达到千万级别。本文主要介绍了美团代码托管平台在迭代演进过程中面临的挑战及解决思路&#xff0c;希望对大家有所帮助或启发。 1. 引言 2. …

六: 数 组(eclipse的使用)

目录复习隐藏eclipse中的工程3.2 一维数组的使用&#xff1a;数组元素的引用p103.2 一维数组的使用&#xff1a;数组元素的默认初始化值p123.2 一维数组的使用p213.2 一维数组的使用p223.3 多维数组的使用p243.3 多维数组的使用p25引用类型的变量&#xff0c;保存的要么是地址值…

国际标准下载的几个网站 请点赞收藏

名称以ITU开头的国际标准下载网站名称以ITU开头的国际标准&#xff0c;是国际电信联盟组织制定的国际标准。国际电信联盟组织简称国际电联。它的官网是http://www.itu.int。 通过在这个网站下输入标准的部门名称就可以搜索下载电信标准。已ITU Q.763为例打开官网输入标准名称点…

C++11 类的新功能

作者&#xff1a;小萌新 专栏&#xff1a;C进阶 作者简介&#xff1a;大二学生 希望能和大家一起进步&#xff01; 本篇博客简介&#xff1a;介绍C11类的新功能和一些关键字 类的新功能类的新功能默认成员函数类成员变量的初始化C11新关键字defaultdeletefinaloverride类的新功…

Java poi之Excel文本图片内容提取

目录结构前言文档准备引入Maven依赖代码块提取结果验证excel03.xls 提取结果excel07.xlsx 提取结果前言 应公司需求&#xff0c;需实现以下功能 Excel文本内容的替换&#xff1b;Excel文本内容的提取&#xff1b;Excel中图片的提取存放 此文章将使用Apache POI实现Excel文件…

我问 ChatGPT:怎样成为优秀的架构师?看它怎么回答的……

要成为一名优秀的架构师,需要以下几个方面的努力: 1. 系统的学习计算机科学和工程相关的知识,如计算机网络,数据结构,算法,操作系统等。 2. 实践和经验积累。参与许多实际的项目,不断积累经验,提高解决问题的能力。 3. 持续学习和追求卓越。保持对新技术和趋势的敏锐观…

Docker容器 01

前言 1.1 从环境配置说起 环境配置是软件开发的一大难题。开发、测试及运维人员需要相同的代码运行环境&#xff0c;如此一来就需要多次搭建环境&#xff0c;想想就觉得麻烦&#xff0c;实际上&#xff0c;在不了解docker等容器技术以前&#xff0c;还真就是这么干的&#xff…

IDEA 中动态web 工程的操作

目录a)IDEA 中如何创建动态web 工程1、创建一个新工程exer&#xff1a;2、在exer下创建module&#xff1a;test3、动态web工程创建成功 如下图b)Web 工程的目录介绍c)如何给动态 web 工程添加额外jar 包1 添加lib目录2 将jar包复制到lib目录中3 将jar包添加到工程4 可以打开项目…

西湖论剑2022部分misc

文章目录签到题喵take_the_zip_easymp3机你太美签到题喵 把文件尾的16进制复制出来&#xff0c;再转换字符串 私信后台即可获得flag take_the_zip_easy 明文攻击 echo -n “dasflow.pcapng” > 1.txt time ./bkcrack -C zipeasy.zip -c dasflow.zip -p 1.txt -o 30 -x 0…

六类网线为啥那么受欢迎,网络工程师必知

目前&#xff0c;国内千兆网络已大规模普及&#xff0c;从前的“销冠”百兆超五类网线已经渐渐淡出了人们的视野&#xff0c;已然被千兆的六类网线取代成为现代布线入门级主力军。万兆超六类网线也同时是城市智能化5G、万物互联时代入门首选&#xff0c;各自顺应时代成为不同领…

APP在Google Play上架被拒的原因

即便了解了Google Play商店的相关政策和应用指南&#xff0c;我们也不能避免应用在上架时或者是应用在更新时被拒的情况发生。那今天我们就展开讲讲Google Play商店被拒的原因及解决方案。 出现不当言论或内容&#xff08;比如&#xff0c;色情内容&#xff0c;带有种族歧视和…

基于Springboot搭建java项目(三十五)—— Ngnix配置的使用

Ngnix配置的使用 一、Nginx配置文件(nginx.conf) 1、配置文件的层级 配置文件目前分为三大部分&#xff0c;全局块、event块和http块&#xff0c;下面是具体的结构 2、配置文件概览 # 全局快 ---------------------------------------------------------------------------…

Bean实例化的基本流程

Bean实例化的基本流程 Bean实例化的基本流程-BeanDefinition Spring容器在进行初始化时&#xff0c;会将xml配置的的信息封装成一个BeanDefinition对象&#xff0c;所有的BeanDefinition存储到一个名为beanDefinitionMap的Map集合中去&#xff0c;Spring框架在对该Map进行遍历…

[NeurIPS 2018] Hyperbolic neural networks

ContentsIntroductionThe Geometry of the Poincar BallHyperbolic space: the Poincar ballGyrovector spaces (陀螺矢量空间)Mbius additionMbius scalar multiplicationDistanceHyperbolic trigonometryConnecting Gyrovector spaces and Riemannian geometry of the Poinca…

如何用提取网页内容的工具快速提取网站内容

随着社会的不断的进步&#xff0c;我们已经进入一个效率时代&#xff0c;相信每个人在互联网上下载或者复制粘贴过内容。特别是整理行业的数据&#xff0c;以及收集资料。今天小编就教大家如何用提取网页内容的工具快速提取到你想要的信息&#xff0c;只需要点几下鼠标就能提取…

GitHub2022年十大热门编程语言榜单

全球知名代码托管平台 GitHub发布的2022年GitHub Octoverse年度报告公布了全球最流行的十大编程语言&#xff0c;其中JavaScript蝉联第一&#xff0c;Python位列次席。 编程是技术革新的核心&#xff0c;对于所有的编程开发人员来说&#xff0c;对世界范围内编程语言发展和趋势…

磨金石教育摄影技能干货分享|人物系列摄影作品欣赏

人间烟火气&#xff0c;最能抚人心。生活中一些平平静静的瞬间&#xff0c;聊天、走路、欢笑&#xff0c;构成了人生当中闪闪的光。今天我们来欣赏一组充满烟火气的人物摄影。没有刻意的姿势&#xff0c;没有华丽的造景&#xff0c;有的就是真实与自然。《放学路上》小时候最欢…

linux 中的压缩和解压操作

1、压缩/解压操作 在开发中&#xff0c;很多时候会遇到某些文件要进行压缩的操作&#xff0c;比如文件较大不方便传输的时候&#xff0c;可能会考虑对文件进行压缩&#xff0c;以减少文件传输的时间。 比如在网络中传输文件的时候&#xff0c;就会考虑先将文件进行压缩&#xf…