当前位置：文档视界 › 清华大学_李津老师_数值分析第二次实验作业

清华大学_李津老师_数值分析第二次实验作业

第二次计算实验 SVD 及其应用

一、方法

从已有函数库或软件包(比如Matlab 库函数及其工具箱，Lapack ，Linpack ，Propack 等)选择至少两个用不同方法计算SVD 的程序，解释算法的思想和主要步骤。

这里，我们主要分析Matlab 中的svd()函数和Propack 中的lansvd()函数。(lansvd()函数程序见附件lansvd.m ，svd()由于是Matlab 中封装好函数无法提取，所以这里不在附附件)。

通过阅读相关文献及程序，我们发现这两种算法都是通过将原矩阵上双对角化，再对此上双对角阵进行奇异值分解，他们主要的区别在于如何将原矩阵上双对角化。

我们令带处理矩阵设为A ，则奇异值分解可描述为下式

000T

A V U ∑??=????

1.Matlab 中的svd()算法描述：

(1)先将m n A ?(m>n)上双对角化，这里使用Householder 方法，其具体过程如下所示将m n A ? 分块为

()11n m A v A -???=??

其中1v 为m 阶列向量，先计算m 阶Householder 变换使得

111111(,)

Pv e R e R δδ=∈∈

并且形成：

()()11

11111

T n m n u A P A A m δ--?-??==??

-???

再计算n －1阶Householder 变换1H 使得

1112121(,)n n H u e R e R γγ--=∈∈

并形成：

()()1

21111,01000T m n A A H δγ-?-??

???=?????????

然后对()()11m n A -?-再循环以上运算，共需循环n-1次进行后，即可得到上双对角阵B ，以及U,V

1223

0000

000

0000

n n B δγδ

γγδ?????

???=?

???????

121...n U P P P =，1211...n V H H H -= 1

i i i I

P P -??=???? ，1

00i i T i I V H -??

=????

(2)对B 进行奇异值分解，我们令

T B B =

对T 进行带Wilkinson 位移的对称QR 迭代即可得到B 的奇异值分解，进而得到A 的奇异值分解。

2.Propack 中lansvd()算法

(1)使用lanczos 方法将m n A ?上双对角化 1. 初始化向量0m

p ∈，10p β=，101u p β=，00v =

2. for i=1,2,…

111

,,T i i i i i i i i i

i i i i

i i i i i r A u v r v r p Av u p u p βαααββ-+++=-===-==

end

这样迭代k 步以后，A 有如下分解形式：

1n n n AV U B +=

其中k V 和1k U +为算法B.1中的Lanczos 向量组合而成，每一列都正交，且已归一化。而

n B 为如下形式的二对角矩阵：

122

31n n n B αβαβαβ+??

? ? ?= ? ? ??

再对n B 进行奇异值分解，其方法与svd()中对上双对角阵的奇异值分解方法相同，这里

就不再赘述了。

二、实际计算

生成十个不同的(最好属不同类型或有不同性质的)的m n ? 矩阵，这里,100m n > ，用你选择的算法对其做SVD,比较不同方法的效果(比如计算小气一直和对应左右奇异向量的误差，效率等)，计算时间和所需存储量等，根据结果提出对算法的认识。

1.误差

在实验中，我们取m=200，n=100，利用orth()函数生成了正交矩阵U 、V ，再生成了不同奇异值分布的奇异值矩阵Σ，再通过A =UΣV ，计算出不同的待分解矩阵。

i is

u u -

i is

v v -

通过对比上面两个误差表，我们发现两种方法在对奇异向量矩阵的估计误差相近，并没有明显的优异性，在估计奇异值上svd 方法略强于lansvd 方法。

他们的运行时间以及使用内存情况如下所示， svd 方法

lansvd 方法

由以上两数据对比可知，svd 方法效率明显优于lansvd 方法，而二者的内存使用情况相近。

综上，我们不难得到如下结论，两种方法在内存使用以及估计精度上没有明显的差别，而在效率上svd 方法远远优于lansvd 方法。

三、图像压缩

与传统图像压缩的概念相同，SVD 方法也是舍弃部分图像数据从而达到压缩的效果，这里我们舍弃的数据是奇异值较小奇异向量所包含的数据，具体方法如下所示:

1.将图像A 做奇异值分解，即

000A U V ∑??=????

2.将∑中较小的奇异值舍去，并且抛弃其所对应的奇异向量舍弃，并将剩余的数据

存储，形成压缩的图片信息

我们选择的5幅数字图像如下所示： (1)黑白网格图

(2)横向渐变灰度图

(3)斜向渐变灰度图

(4)Lena灰度图

(5)白噪声灰度图

我们选择抛弃小于最大奇异值的1%所对应的的奇异值与奇异向量如下所示：(1)黑白网格图

(2)横向渐变灰度图

(4)Lena灰度图

通过对比以上五张数字图像，直观上我们发现，形状越复杂的图像压缩后的效果越差，我们利用峰值信噪比来评价图片压缩后的质量

通过上表中的数据对比，我们就可以发现越复杂的图像越难压缩，且压缩后的质量很差。

比较SVD方法与传统的JPEG方法可知，再存储简单图像时SVD方法可能会获得更高的压缩率，因为他不需要想JPEG记录图像上每块儿的信息，而是直接对整体图像进行变换存储，而对于较复杂的图像时，由于SVD方法需要记录很多的奇异向量所以他的压缩率可能会低于有固定正交向量矩阵的JPEG方法。

附：

generator.m 为第二部分中生成计算矩阵的程序

my_svd.m 为利用不同房法求解矩阵奇异值分解的程序fig_process. 为处理图像的程序

lansvd.m 为Propack中的奇异值分解程序