R 中的資料探勘演算法/降維/奇異值分解

在本章中，我們將探討奇異值分解（SVD），這是一種矩陣分解方法，利用線性代數知識來進行這樣的分解。

奇異值分解

在資料探勘中，此演算法可用於透過顯示重要維度的數量來更好地理解資料庫，並透過減少資料探勘過程中使用的屬性數量來簡化資料庫。這種減少消除了從線性代數角度來看線性相關的非必要資料。例如，假設一個數據庫包含一個欄位，該欄位儲存多個樣本的水溫，另一個欄位儲存其狀態（固體、液體或氣體）。很容易看出，第二個欄位依賴於第一個欄位，因此，SVD 可以輕鬆地向我們表明它對於分析並不重要。

主成分分析 (PCA) 是 SVD 的一個特例。

演算法

SVD 是對具有$n$ 行和 $d$ 列的實數或複數矩陣 $X$ 的因式分解。

$X=UAV^{\top }$

其中 $U$ 是一個維度為 $n\times n$ 的矩陣，$V$ 是另一個維度為 $d\times d$ 的矩陣，$A$ 是一個維度為 $n\times d$ 的矩陣，與 $X$ 的維度相同。

此外

$U^{\top }U=I_{n}$ 和 $V^{\top }V=I_{d}$

其中 $I_{n}$ 和 $I_{d}$ 是單位矩陣，大小分別為 $n$ 和 $d$.

矩陣 $U$ 的列是矩陣 $X$ 的 **左奇異向量**，矩陣 $V$ （或 $V^{\top }$ 的行）的列是 **右奇異向量**。

矩陣 $A$ 是對角矩陣，其對角線上的值是矩陣 $X$ 的 **奇異值**。矩陣 $A$ 中每行上的奇異值永遠不會小於其下方行的值。所有奇異值都大於 $0$.

計算 SVD 就是找到 $XX^{\top }$ 和 $X^{\top }X$ 的特徵值和特徵向量。 $X^{\top }X$ 的特徵向量是 $V$ 的列，而 $XX^{\top }$ 的特徵向量是 $U$ 的列。矩陣 $X$ 的奇異值（在矩陣 $A$ 的對角線上）是 $XX^{\top }$ 和 $X^{\top }X$ 共有的正特徵值的平方根。

如果 $XX^{\top }$ 和 $X^{\top }X$ 具有相同數量的特徵值，則 $X$ 是一個方陣；否則，特徵值較少的矩陣的特徵值是特徵值較多的矩陣的特徵值。因此，$X$ 的奇異值是矩陣的特徵值，位於 $XX^{\top }$ 和 $X^{\top }X$ 中，特徵值較少的矩陣中。

矩陣的奇異值數量等於該矩陣的秩，即矩陣中線性無關的列或行的數量。秩不大於 $\min(n,d)$，因為這是矩陣對角線上元素的數量。奇異值是矩陣 $A$ 對角線上的元素。正奇異值的數量等於矩陣的秩。

因此，演算法如下

1. 計算 $XX^{\top }$ 通常。
2. 計算 $XX^{\top }$ 的特徵值和特徵向量通常。
3. 計算 $X^{\top }X$ 。
4. 計算 $X^{\top }X$ 的特徵值和特徵向量。
5. 計算 $XX^{\top }$ 和 $X^{\top }X$ 的公共正特徵值的平方根。
6. 最後，將計算出的值分配給 $U$，$V$ 和 $A$ 。

示例

X = ${\begin{bmatrix}1&1&0&2\\1&1&2&0\\2&0&1&1\end{bmatrix}}$

XX^T = ${\begin{bmatrix}6&2&4\\2&6&4\\4&4&6\end{bmatrix}}$

XX^T 的特徵值為

12.744563, 4.000000 和 1.255437

XX^T 的特徵向量為

1. 0.5417743, 0.5417743, 0.6426206
2. 0.7071068, -0.7071068, -2.144010 x 10^-16
3. 0.4544013, 0.4544013, -0.7661846

X^TX = ${\begin{bmatrix}6&2&4&4\\2&2&2&2\\4&2&5&1\\4&2&1&5\end{bmatrix}}$

X^TX 的特徵值為

12.744563, 4.000000, 1.255437 和 (5.940557 x 10^-18)

X^TX 的特徵向量為

1. 0.6635353, 0.3035190, 0.4835272, 0.4835272
2. 0.0000000, 0.5181041 x 10^-16, -0.7071068, 0.7071068
3. -0.5565257, 0.8110957, 0.1272850, 0.1272850
4. 0.5000000, 0.5000000, -0.5000000, -0.5000000

X 的奇異值為 X^TX 和 XX^T 的共同特徵值的平方根

${\sqrt {12.744563}}=3.569953$

${\sqrt {4.000000}}=2.000000$

${\sqrt {1.255437}}=1.120463$

因此

A = ${\begin{bmatrix}3.569953&0&0&0\\0&2.000000&0&0\\0&0&1.120463&0\end{bmatrix}}$

最後，X 被分解為

X = UAV^T = ${\begin{bmatrix}0.5417743&0.7071068&0.4544013\\0.5417743&-0.7071068&0.4544013\\0.6426206&-2.144010\cdot 10^{-16}&-0.7661846\end{bmatrix}}{\begin{bmatrix}3.569953&0&0&0\\0&2.000000&0&0\\0&0&1.120463&0\end{bmatrix}}{\begin{bmatrix}0.6635353&0.000000&-0.5565257&0.5\\0.3035190&5.181041\cdot 10^{-16}&0.8110957&0.5\\0.4835272&-0.7071068&0.1272850&-0.5\\0.4835272&0.7071068&0.1272850&-0.5\end{bmatrix}}$

實現

R 內建了一個名為'svd()'的函式用於計算 SVD。

預設情況下，它接收 R 的 原生矩陣 作為引數，並返回一個包含 U、A 和 V 的 資料框。

引數

如果不需要計算所有奇異值和向量，可以告訴 svd() 需要多少個元素。

這可以透過將這些值分別分配給 **nu** 和 **nv** 來實現，它們分別代表所需的左奇異向量和右奇異向量數量。

例如，為了只計算一半的向量，可以執行以下操作

svd(X, nu = min(nrow(X), ncol(X)) / 2, nv = min(nrow(X), ncol(X))/2)

返回的物件

s = svd(X)

考慮到

X = UAV^T

返回的物件是一個包含三個欄位的資料結構

s$d 是包含 X 的奇異值的向量，它來自矩陣 A 的對角線。

s$u 是一個矩陣，它的列包含 X 的左奇異向量。它的行數與 X 的行數相同，它的列數與傳遞給引數 nu 的數字相同。注意，如果 nu 為 0，則不會建立此矩陣。

s$v 是一個矩陣，它的列包含 X 的右奇異向量。它的行數與 X 的列數相同，它的列數與傳遞給引數 nv 的數字相同。注意，如果 nv 為 0，則不會建立此矩陣。

例子

在 R 終端或 R 程式中執行以下命令序列，並檢視結果

示例 1

dat <- seq(1,240,2)

X <- matrix(dat,ncol=12)

s <- svd(X)

A <- diag(s$d)

s$u %*% A %*% t(s$v) #  X = U A V'

示例 2

dat <- seq(1,240,2)

X <- matrix(dat,ncol=12)

s <- svd(X, nu = nrow(X), nv = ncol(X))

A <- diag(s$d)

A <- cbind(A, 0)  # Add two columns with zero, in order to A have the same dimensions of X.

A <- cbind(A, 0)

s$u %*% A %*% t(s$v) #  X = U A V'

案例研究

場景

為了更好地視覺化 SVD 的結果，在這個案例研究中，我們將使用一個表示影像的矩陣，這樣任何對矩陣的更改都可以輕鬆觀察到。

輸入資料

要使用 R 處理影像，應安裝 rimage 包

> install.packages("rimage")

然後將影像載入到 R 中，將其轉換為灰度影像。這樣，我們將得到一個二進位制矩陣，其中 1 表示白色，0 表示黑色。

library(rimage)
tux <- read.jpeg('tux.jpeg')
tux <- imagematrix(tux,type='grey')

我們可以看到此匯入的結果

plot(tux)

為了幫助我們進行降維，讓我們建立一個小的輔助函式，該函式將接收我們的 tux 和我們想要的維度數，並返回我們的新的 tux。

reduce <- function(A,dim) {
    #Calculates the SVD
    sing <- svd(A)

    #Approximate each result of SVD with the given dimension
    u<-as.matrix(sing$u[, 1:dim])
    v<-as.matrix(sing$v[, 1:dim])
    d<-as.matrix(diag(sing$d)[1:dim, 1:dim])

    #Create the new approximated matrix
    return(imagematrix(u%*%d%*%t(v),type='grey'))
}

執行和輸出

現在我們有了矩陣，讓我們看看它有多少奇異值

tux_d <- svd(tux)
length(tux_d$d)
[1] 335

因此，這個矩陣有 335 個奇異值。讓我們首先嚐試將其減少到只有一個奇異值

plot(reduce(tux,1))

正如我們所看到的，這種近似值從我們的矩陣中刪除了許多有用的資訊。如果它是資料庫，我們肯定已經丟失了重要的資料。

讓我們嘗試使用 10% (35) 個奇異值

plot(reduce(tux,35))

分析

僅使用 10% 的真實資料，我們就能創建出真實資料的非常好的近似值。此外，透過這種方法，我們只需使用最重要的奇異值就可以去除噪聲和線性相關元素。這在資料探勘中非常有用，因為很難確定資料庫是否“乾淨”，如果不是，哪些元素對我們的分析有用或無用。

參考文獻

• 資料探勘演算法基礎 - Mohammed J. Zaki，Wagner Meira Jr. - 第 7.4 章
• 奇異值分解
• http://www.ats.ucla.edu/stat/r/code/svd_demos.htm
• http://stat.ethz.ch/R-manual/R-patched/library/base/html/svd.html