数值线性代数:奇异值分解SVD

本文记录计算矩阵奇异值分解SVD的原理与流程。

注1:限于研究水平,分析难免不当,欢迎批评指正。

零、预修

0.1 矩阵的奇异值

A\in \mathbb{R}^{m\times n}列满秩矩阵,若\boldsymbol{A}^{T}\boldsymbol{A}的特征值为\lambda_{1}\geq \lambda_{2}\geq \cdots \lambda_{r}>0,则称\sigma_{i}=\sqrt{\lambda_{i}}为矩阵\boldsymbol{A}的奇异值。

0.2 SVD(分解)定理

\boldsymbol{A}\in \mathbb{R}^{m\times n},则存在正交矩阵\boldsymbol{U}\in \mathbb{R}^{m\times m}\boldsymbol{V}\in \mathbb{R}^{n\times n},使得

\boldsymbol{U}^{T}\boldsymbol{A}\boldsymbol{V}=\begin{pmatrix} \sum_{r} & 0\\ 0& 0 \end{pmatrix}

其中,\sum_{r}=diag\left ( \sigma_{1},\sigma_{2},\cdots ,\sigma_{r} \right )\sigma_{1}\geq \sigma_{2}\geq \cdots \sigma_{r}>0\sigma_{i}即为矩阵\boldsymbol{A}的奇异值。

考虑下述两种情形:

  • 情形1:\tilde{\boldsymbol{A}}=\boldsymbol{A}^{T}\boldsymbol{A}

\left (\boldsymbol{U}^{T}\boldsymbol{A}\boldsymbol{V} \right )^{T}\left ( \boldsymbol{U}^{T}\boldsymbol{A}\boldsymbol{V} \right )=\boldsymbol{V}^{T}\left (\boldsymbol{A}^{T}\boldsymbol{A} \right )\boldsymbol{V}=\begin{pmatrix} \sum_{r}^{2} &\boldsymbol{0} \\ \boldsymbol{0}& \boldsymbol{0} \end{pmatrix}

其中,\sum_{r}^{2}=diag\left ( \sigma_{1}^{2},\sigma_{2}^{2},\cdots ,\sigma_{r}^{2} \right )=diag\left ( \lambda_{1},\lambda_{2},\cdots ,\lambda_{r} \right )

由此可以看出,\tilde{\boldsymbol{A}}=\boldsymbol{A}^{T}\boldsymbol{A}\in \mathbb{R}^{n\times n},通过计算矩阵\boldsymbol{A}的奇异值\sigma_{i},便可矩阵\tilde{\boldsymbol{A}}的特征值\lambda_{i}=\sigma_{i}^{2},而矩阵\boldsymbol{V}即为矩阵\tilde{\boldsymbol{A}}的特征向量

  • 情形2:\boldsymbol{A}^{T}=\boldsymbol{A}, m=n

\boldsymbol{A}\boldsymbol{x}=\sigma \mathbf{x},则\boldsymbol{A}^{T}\boldsymbol{A}\boldsymbol{x}=\sigma \boldsymbol{A}^{T} \boldsymbol{x}=\sigma \boldsymbol{A}\boldsymbol{x}=\sigma ^{2}\boldsymbol{x},也就是说,\sigma^{2}\boldsymbol{A}^{T}\boldsymbol{A}的特征值,\boldsymbol{x}也是\boldsymbol{A}^{T}\boldsymbol{A}的特征向量。同时考虑到实对称矩阵\boldsymbol{A}的秩为n,所以\boldsymbol{A}^{T}\boldsymbol{A}的特征值/特征向量也是\boldsymbol{A}的特征值/特征向量。

0.3 Householder变换

\boldsymbol{w}\in \mathbb{R}^{n},且\boldsymbol{w}^{T}\boldsymbol{w}=1,定义\boldsymbol{H}=I-2\boldsymbol{w}\boldsymbol{w}^{T}为Householder变换。

对于非零向量\boldsymbol{x}\in \mathbb{R}^{n},可构造\boldsymbol{H}=I-2\boldsymbol{w}\boldsymbol{w}^{T},使得

Hx=\alpha \boldsymbol{e}_{1}

其中,\alpha =\pm \left \| x \right \|_{2}\boldsymbol{w}=\frac{x-\alpha e_{1}}{\left \| x-\alpha e_{1} \right \|_{2}}\boldsymbol{e}_{1}^{T}=\begin{pmatrix} 1 &0 &\cdots &0 \end{pmatrix}\in \mathbb{R}^{1\times n}

\boldsymbol{x}=\begin{pmatrix} \boldsymbol{x}^{\left ( 1 \right )}\\ \boldsymbol{x}^{\left ( 2 \right )} \end{pmatrix}\boldsymbol{x}^{\left ( 1 \right )}=\boldsymbol{x}\left ( 1:k-1 \right )\boldsymbol{x}^{\left ( 2 \right )}=\boldsymbol{x}\left ( k:n \right ),对于\boldsymbol{H}_{k}=\begin{pmatrix} \boldsymbol{I}_{k-1} & \boldsymbol{0}\\ \boldsymbol{0} & \mathbf{H}^{\left (k \right )} \end{pmatrix}

\boldsymbol{H}_{k}\boldsymbol{x}=\begin{pmatrix} \boldsymbol{I}_{k-1} & \boldsymbol{0}\\ \boldsymbol{0} & \mathbf{H}^{\left (k \right )} \end{pmatrix}\begin{pmatrix} x^{\left ( 1 \right )}\\ x^{\left ( 2 \right )} \end{pmatrix}=\begin{pmatrix} x^{\left ( 1 \right )}\\ \boldsymbol{H}^{\left ( k \right )}x^{\left ( 2 \right )}\end{pmatrix}

根据上述结论可知,可以构造\boldsymbol{H}^{\left (k \right )},使得\boldsymbol{H}^{\left ( k \right )}x^{\left ( 2 \right )}=\begin{pmatrix} -M\sigma _{k} & 0& 0& 0 \end{pmatrix}^{T}\in \mathbb{R}^{n-k+1}

具体来说,可按照下述流程进行操作:

M=max\left ( \left | x_{k} \right |,\left | x_{k+1} \right |,\cdots ,\left | x_{n} \right | \right )\\ x_{i}\leftarrow \frac{x_{i}}{M}\ \\ \sigma_{k}=sign\left ( x_{k} \right )\sum_{i=k}^{i=n}x_{i}^{2} \\ \boldsymbol{U}_{k}=\left ( \sigma_{k}+x_{k},x_{k+1},\cdots ,x_{n} \right )\in \mathbb{R}^{n-k+1} \\ \rho_{k}=\frac{1}{2}\boldsymbol{U}_{k}^{T}\boldsymbol{U}=\sigma_{k}\left ( \sigma_{k}+x_{k} \right )\\ \boldsymbol{H}^{\left ( k \right )}=\boldsymbol{I}-\frac{1}{\rho_{k}}\boldsymbol{U}_{k}\boldsymbol{U}_{k}^{T}=\boldsymbol{I}-\frac{\sigma_{k}+x_{k}}{\sigma_{k}}\frac{\boldsymbol{U}_{k}}{\sigma_{k}+x_{k}}\frac{\boldsymbol{U}_{k}^{T}}{\sigma_{k}+x_{k}}\in \mathbb{R}^{\left (n-k+1 \right )\times \left ( n-k+1 \right )}

由此,通过Householder变换,可以将某一列向量的部分连续元素约化为0。

0.4 Givens变换

\boldsymbol{e}_{1},\boldsymbol{e}_{2},\cdots ,\boldsymbol{e}_{n}是n维Euclid空间\mathbb{R}^{n}中的一组标准正交基,\forall x, y \in \mathbb{R}^{n},则在平面\left ( \boldsymbol{e}_{i}, \boldsymbol{e}_{j}\right )中存在旋转变换矩阵\boldsymbol{G}\left ( i,j,\theta \right ),满足

\boldsymbol{y}=\mathbf{G}\mathbf{x}

其中,

x=\begin{pmatrix} x\left ( 1:i-1 \right )\\ x_{i}\\ x\left ( i+1:j-1 \right )\\ x_{j}\\ x\left ( j+1:n \right ) \end{pmatrix}, y=\begin{pmatrix} x\left ( 1:i-1 \right )\\ 0\\ x\left ( i+1:j-1 \right )\\ \sqrt{x_{i}^{2}+x_{j}^{2}}\\ x\left ( j+1:n \right ) \end{pmatrix}

\boldsymbol{G}=\begin{pmatrix} \boldsymbol{I}_{i-1,i-1}& & & & \\ & cos\theta & & -sin\theta & \\ & & \boldsymbol{I}_{j-i-1,j-i-1}& & \\ & sin\theta & & cos\theta & \\ & & & & \boldsymbol{I}_{n-j,n-j} \end{pmatrix}

 cos\theta =\frac{x_{i}}{\sqrt{x_{i}^{2}+x_{j}^{2}}},sin\theta =\frac{-x_{j}}{\sqrt{x_{i}^{2}+x_{j}^{2}}},r=\sqrt{x_{i}^{2}+x_{j}^{2}}

由此可以看出,Givens变换可以将向量的某个元素约化为0。

一、大型矩阵特征值/特征向量的求解思路

大型矩阵\boldsymbol{A}特征值/特征向量求解一般按照以下流程进行

  1. 将矩阵\boldsymbol{A}约化为特征值/特征向量容易求解的矩阵\boldsymbol{H}
  2. 求解矩阵\boldsymbol{H}的特征值/特征向量;
  3. 将矩阵\boldsymbol{H}的特征值/特征向量转化成矩阵\boldsymbol{A}约化为特征值/特征向量;

二、隐式QR计算矩阵奇异值分解

参考书籍

Golub G H , Loan C F V .Matrix Computations.Johns Hopkins University Press,1996.

Ford W .Numerical Linear Algebra with Applications using MATLAB. 2014.

徐树方. 数值线性代数(第二版).  北京大学出版社, 2010.

参考文献

Golub G. and Kahan W.. Calculating the Singular Values and Pseudo-Inverse of a Matrix. Journal of the Society for Industrial and Applied Mathematics: Series B, Numerical Analysis, 1965, 2(2) : 205-224.

Demmel J., Kahan W..Accurate Singular Values of Bidiagonal Matrices. SIAM Journal on Scientific and StatisticalComputing, 1990, 11(5):873-912.

P. A. Businger,G. H. Golub. Singular value decomposition of a complex matrix. communications of the acm, 1969.

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

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

相关文章

神经网络简单介绍

人工神经网络(artififial neural network) 简称神经网络,它是一种模仿生物神经网络结构和功能的非线性数学模型。 神经网络通过输入层接受原始特征信息,再通过隐藏层进行特征信息的加工和提取,最后通过输出层输出结果。 根据需要神经网络可以…

RTPSv2.2(中文版)

实时发布订阅协议 (RTPS) DDS互操作性 有线协议规范 V2.2 (2014-09-01正式发布) https://www.omg.org/spec/DDSI-RTPS/2.2/PDF 目 录 1 范围Scope 9 2 一致性Conformance 9 3 参考文献References 9 4 术语和定义Terms a…

centos(linux)手动配置ip地址

查看系统 查看ip 进入网卡配置的目录 [root178-119-30-16 ~]# cd /etc/sysconfig/network-scripts/ [root178-119-30-16 network-scripts]# ls ifcfg-ens192 ifdown ifdown-ippp ifdown-post ifdown-sit ifdown-tunnel ifup-bnep ifup-ipv6 ifup-plus…

【Mysql】万字长文带你快速掌握数据库基础概念及SQL基本操作

文章目录 前言一、数据库相关概念1. 什么是数据库2. 数据库的种类3. Mysql 简介4. SQL简介5. 数据库中常见的数据类型 二、SQL基础1. SQL通用语法2. SQL的主要分类3. DDL(数据库,表,索引,视图)4. DML(数据的…

14:00面试,14:06就出来了,问的问题有点变态。。。

从小厂出来,没想到在另一家公司又寄了。 到这家公司开始上班,加班是每天必不可少的,看在钱给的比较多的份上,就不太计较了。没想到5月一纸通知,所有人不准加班,加班费不仅没有了,薪资还要降40%,…

【每日一题】—— A - 1-2-4 Test (AtCoder Beginner Contest 270)

🌏博客主页:PH_modest的博客主页 🚩当前专栏:每日一题 💌其他专栏: 🔴 每日反刍 🟡 C跬步积累 🟢 C语言跬步积累 🌈座右铭:广积粮,缓称…

HiveSQL SparkSQL中常用知识点记录

目录 0. 相关文章链接 1. hive中多表full join主键重复问题 2. Hive中选出最新一个分区中新增和变化的数据 3. Hive中使用sort_array函数解决collet_list列表排序混乱问题 4. SQL中对小数位数很多的数值转换成文本的时候不使用科学计数法 5. HiveSQL & SparkSQL中炸裂…

小程序动态隐藏分享按钮

// 禁用分享 wx.hideShareMenu({menus: [shareAppMessage, shareTimeline] })// 显示分享 wx.showShareMenu({withShareTicket: true,menus: [shareAppMessage, shareTimeline] })//私密消息 wx.updateShareMenu({isPrivateMessage: true, })

list与erase()

运行代码: //list与erase() #include"std_lib_facilities.h" //声明Item类 struct Item {string name;int iid;double value;Item():name(" "),iid(0),value(0.0){}Item(string ss,int ii,double vv):name(ss),iid(ii),value(vv){}friend istr…

JMM是如何保证原子性呢?

3.3.1加锁 锁是一种通用技术,比如Java提供的Synchronized关键字就是锁的一种实现,Synchronized是排他锁/独占锁,就是有你没我的意思,只要其他线程到来访问,发现锁还未释放,就要在外面等待,因为S…

Redis—相关背景

Redis—相关背景 🔎Redis—特性In-memory data structures—在内存中存储数据Programmability—可编程性Extensibility—可扩展性Persistence—持久化Clustering—集群High availability—高可用 🔎Redis 为什么快🔎Redis 的使用场景Real-tim…

【TiDB理论知识06】PD架构与作用

目录 一 PD的架构与功能 PD架构 PD作用 名词解释 路由功能 二 TSO的分配 概念 分配过程 性能问题 高可用问题 三 PD的调度原理 总流程 1 信息收集 2 生成调度 3 执行调度 四 Label的作用 Label的配置 给TiKV打标签 PD配置 一 PD的架构与功能 PD架构 PD集群…

学习Maven Web 应用

Maven Web 应用 本章节我们将学习如何使用版本控制系统 Maven 来管理一个基于 web 的项目,如何创建、构建、部署已经运行一个 web 应用。 创建 Web 应用 我们可以使用 maven-archetype-webapp 插件来创建一个简单的 Java web 应用。 打开命令控制台,…

数据库应用:Redis安装部署

目录 一、理论 1.缓存 2.关系型数据库与非关系型数据库 3.Redis 4.Redis安装部署 5.Redis命令工具 6.Redis数据库常用命令 7.Redis多数据库操作 二、实验 1.Redis安装部署 2.Redis命令工具 3.Redis数据库命令 4.Redis多数据库操作 三、问题 1.RESP连接CentOS 7 R…

【单片机】温控系统参数辨识及单片机PID控制

温控系统参数辨识及单片机PID控制 1. 温控系统组成2. matlab辨识系统参数2.1 采集阶跃响应信号导入matlab系统辨识模块 PID控制 1. 温控系统组成 半导体制冷片正向通电制冷,反向通电制热。系统采用半导体制冷片(帕尔贴)作为执行单元&#xf…

微信小程序 样式和全局配置

WXSS wxss 把屏幕分为750个物理像素,大屏大,小屏小,随着设备不一致自动适配 推荐使用iPhone6作为标准,1个rpx 0.5个px,把px乘以2就是rpx的参数 import 导入外部样式表 import /common/common.wxss 样式 权重一…

Vue2 第四节 计算属性,监视属性

1.计算属性 2.监视属性 3.计算属性与监视属性之间的关系 一.计算属性 定义:要用的属性不存在,要通过已有属性计算得来原理:底层借助了Object.defineproperty方法提供的getter和setterget函数什么时候会执行:初次读取的时候会执…

day47-SSM分页

SSM分页(增删改查登录注册) applicationContext.xml中加入mybatis-config.xml路径 mybatis-config.xml Mapper接口 Service接口及其实现类 Mapper.xml page.jsp personDetail.jsp addPerson.jsp updatePerson.jsp login.jsp regist…

搭建自己的Git服务器

环境 服务端:Ubuntu 22.04 客户端:Win11_x64 前提条件:需要确保在Windows机器上能够ping通Ubuntu服务器, 并且服务端与客户端均已安装了Git软件 服务端上的配置操作 以Ubuntu服务器作为Git服务端的运行环境,并方便后期免密推…

如何做好IT类的技术面试

目录 一、IT行业的招聘渠道 二、如何做好技术面试官 三、谈谈IT行业如何做好招聘工作 四、面试IT公司的小技巧 五、面试有哪些常见的问题 六、关于面试的一些建议 面试可能是我们每个人都必须会遇到的事情,而技术面试更具有专业性,以下会从几个方面…