協方差矩陣(協方差矩陣的意義)

這篇文章的目標是從數據預處理的基礎知識到深度學習中使用的現代技術。我的觀點是我們可以使用代碼（例如Python / NumPy）來更好地理解抽象的數學概念。

我們將從數據科學和機器學習/深度學習中的基本但非常有用的概念開始，如方差和協方差矩陣。我們將進一步介紹一些用于將圖像輸入神經網絡的預處理技術。我們將嘗試使用Python代碼獲得更具體的見解，以實際了解每個等式的作用。

預處理是指在將原始數據饋送到機器學習或深度學習算法之前對原始數據進行的所有轉換。例如，在原始圖像上訓練卷積神經網絡可能會導致糟糕的分類性能。預處理對于加速訓練也很重要。

背景

A.方差和協方差

方差是在概率論和統計方差衡量隨機變量或一組數據時離散程度的度量。協方差表示的是兩個變量的總體的誤差，這與只表示一個變量誤差的方差不同。 如果兩個變量的變化趨勢一致，也就是說如果其中一個大于自身的期望值，另外一個也大于自身的期望值，那么兩個變量之間的協方差就是正值。 如果兩個變量的變化趨勢相反，即其中一個大于自身的期望值，另外一個卻小于自身的期望值，那么兩個變量之間的協方差就是負值。

協方差值取決于變量的大小，因此很難對其進行分析。可以使用更容易解釋的相關系數。相關系數就是歸一化協方差。

上圖中：正協方差意味著一個變量的大值與另一個變量的大值相關聯（左）。負協方差意味著一個變量的大值與另一個變量的小值相關聯（右）。

協方差矩陣是一個矩陣，它總結了一組向量的方差和協方差，它可以告訴你關于變量的很多事情。對角線對應于每個向量的方差：

矩陣 A及其協方差矩陣。對角線對應于每個列向量的方差

我們來看看方差的公式：

例如，A的第一列向量的方差是：

這是我們協方差矩陣的第一個單元格。對角線上的第二個元素對應于第二列向量與A的方差，依此類推。

注:從矩陣A中提取的向量對應于矩陣A的列。

其他單元格對應的是來自a的兩個列向量之間的協方差，例如，第一列和第三列之間的協方差位于協方差矩陣中，即第1列和第3行(或第3列和第1行)。

讓我們檢查A的第一和第三列向量之間的協方差是否等于-2.67。兩個變量X和Y之間的協方差公式為：

變量X和Y是上面例子中的第一和第三列向量。讓我們分開這個公式，以確保它非常清晰：

1、求和符號(Σ)意味著我們將遍歷元素的向量。我們從第一個元素(i=1)開始計算X的第一個元素減去向量X的均值。

2.將結果乘以Y的第一個元素減去向量Y的平均值。

3.重復向量的每個元素的過程并計算所有結果的總和。

4.除以向量中的元素數。

例1

讓我們從矩陣A開始：

我們將計算第一個和第三個列向量之間的協方差：

和

x? = 3，?= 4，n = 3，所以我們有：

就是協方差矩陣的值。

簡單的方法。使用NumPy，可以使np.cov用函數計算協方差矩陣。

值得注意的是，如果您希望NumPy將列用作向量，則rowvar=Fal，必須使用該參數。此外，bias=True 除以n 而不是n-1。

讓我們先用Python創建數組：

# First things first: let's import some librariesimport numpy as npimport matplotlib.pyplot as pltimport aborn as sns A = np.array([[1, 3, 5], [5, 4, 1], [3, 8, 6]])

現在我們將使用NumPy函數計算協方差：

np.cov(A, rowvar=Fal, bias=True)

array([[ 2.66666667, 0.66666667, -2.66666667],

[ 0.66666667, 4.66666667, 2.33333333],

[-2.66666667, 2.33333333, 4.66666667]])

使用點積查找協方差矩陣

還有另一種計算A的協方差矩陣的方法。你可以將A居中于0左右。從向量的每個元素中減去向量的平均值，得到一個平均值等于0的向量。它乘以自己的轉置，并除以觀察的數量。 讓我們從一個Python實現開始：

def calculateCovariance(X): meanX = np.mean(X, axis = 0) lenX = X.shape[0] X = X - meanX covariance = X.T.dot(X)/lenXreturn covariance

讓我們在矩陣A上測試它：

calculateCovariance(A)

array([[ 2.66666667, 0.66666667, -2.66666667],

[ 0.66666667, 4.66666667, 2.33333333],

[-2.66666667, 2.33333333, 4.66666667]])

我們最終得到與以前相同的結果。

解釋很簡單。兩個向量之間的點積可以表示為：

它是向量的每個元素的乘積之和：

點積對應于矢量的每個元素的乘積之和

如果n是我們向量中的元素數，并且我們除以n：

你可以注意到，這與我們之前看到的協方差公式相差不遠：

唯一的區別是，在協方差公式中，我們從每個元素中減去向量的均值。這就是為什么我們需要在做點積之前對數據進行center處理。

現在，如果我們有一個矩陣A，A和它的轉置之間的點積將給你一個新的矩陣：

如果從zero-centered矩陣開始，此矩陣與其轉置之間的點積將給出每個向量的方差和它們之間的協方差，即協方差矩陣。

B.可視化數據和協方差矩陣

為了獲得關于協方差矩陣以及它如何有用的更多見解，我們將創建一個函數來將其與2D數據一起可視化。

該函數將如上所述計算協方差矩陣。它將創建兩個子圖 - 一個用于協方差矩陣，另一個用于數據。Seaborn的heatmap()函數用于創建顏色漸變 - 小值將以淺綠色著色，大值以深藍色著色。我們選擇了一種調色板顏色，但您可能更喜歡其他顏色。數據表示為散點圖。

def plotDataAndCov(data): ACov = np.cov(data, rowvar=Fal, bias=True) print 'Covariance matrix: ', ACov fig, ax = plt.subplots(nrows=1, ncols=2) fig.t_size_inches(10, 10) ax0 = plt.subplot(2, 2, 1) # Choosing the colors cmap = sns.color_palette("GnBu", 10) sns.heatmap(ACov, cmap=cmap, vmin=0) ax1 = plt.subplot(2, 2, 2) # data can include the colors if data.shape[1]==3: c=data[:,2] el: c="#0A98BE" ax1.scatter(data[:,0], data[:,1], c=c, s=40) # Remove the top and right axes from the data plot ax1.spines['right'].t_visible(Fal)ax1.spines['top'].t_visible(Fal)

c.模擬數據

不相關的數據

現在我們有了繪圖函數，我們將生成一些隨機數據來可視化協方差矩陣能告訴我們什么。我們將從使用NumPy函數np.random.normal()從正態分布中提取的一些數據開始。

使用NumPy從正態分布中繪制樣本

該函數需要平均值，標準偏差和分布的觀察數作為輸入。我們將創建300個觀測值的兩個隨機變量，標準差為1.第一個平均值為1，第二個平均值為2.如果我們從正態分布中隨機抽取兩組300個觀測值，則兩個矢量將是不相關的。

np.random.ed(1234)a1 = np.random.normal(2, 1, 300)a2 = np.random.normal(1, 1, 300)A = np.array([a1, a2]).TA.shape

(300,2)

注意1:我們用.T轉置數據，因為原始的shape是(2,300)，我們想要觀察到的數據的行數(shape(300,2))。

注意2:我們使用np.random.ed函數進行再現性。下次運行單元格時將使用相同的隨機數。

讓我們看看數據是怎樣的:

A[:10,:]

array([[ 2.47143516, 1.52704645],

[ 0.80902431, 1.7111124 ],

[ 3.43270697, 0.78245452],

[ 1.6873481 , 3.63779121],

[ 1.27941127, -0.74213763],

[ 2.88716294, 0.90556519],

[ 2.85958841, 2.43118375],

[ 1.3634765 , 1.59275845],

[ 2.01569637, 1.1702969 ],

[-0.24268495, -0.75170595]])

現在，我們可以檢查分布是否正常：

sns.distplot(A[:,0], color="#53BB04")sns.distplot(A[:,1], color="#0A98BE")plt.show()plt.clo()

我們可以看到分布具有相同的標準偏差但是具有不同的均值（1和2）。這正是我們所要求的。

現在我們可以用我們的函數繪制數據集及其協方差矩陣：

plotDataAndCov(A)plt.show()plt.clo()

Covariance matrix:

[[ 0.95171641 -0.0447816 ]

[-0.0447816 0.87959853]]

我們可以在散點圖上看到兩個維度是不相關的。注意，我們有一個平均值為1（y軸），另一個平均值為2（x軸）。

協方差矩陣表明每個變量的方差非常大(約1)，列1和列2的協方差非常小(約0)，因為我們保證了兩個向量是獨立的，所以這是一致的。反之則未必正確:協方差為0并不保證獨立性

相關數據

現在，讓我們通過從另一列中指定一列來構造依賴數據。Python代碼如下：

np.random.ed(1234)b1 = np.random.normal(3, 1, 300)b2 = b1 + np.random.normal(7, 1, 300)/2.B = np.array([b1, b2]).TplotDataAndCov(B)plt.show()plt.clo()

Covariance matrix:

[[ 0.95171641 0.92932561]

[ 0.92932561 1.12683445]]

在散點圖上可以看到兩個維度之間的相關性。我們可以看到可以繪制一條線并用于從x預測y，反之亦然。協方差矩陣不是對角線的（對角線外有非零單元）。這意味著維度之間的協方差不為零。

2.預處理

A.平均歸一化

平均歸一化只是去除每次觀察的平均值。

其中X'是歸一化數據集，X是原始數據集，x?是X的平均值。

均值歸一化具有將數據居中于0的效果。我們將創建執行此操作的函數center()：

def center(X): newX = X - np.mean(X, axis = 0) return newX

讓我們試試我們之前創建的矩陣B：

BCentered = center(B) print 'Before: ' plotDataAndCov(B)plt.show()plt.clo() print 'After: ' plotDataAndCov(BCentered)plt.show()plt.clo()

Before:

Covariance matrix:

[[ 0.95171641 0.92932561]

[ 0.92932561 1.12683445]]

After:

Covariance matrix:

[[ 0.95171641 0.92932561]

[ 0.92932561 1.12683445]]

第一個圖再次顯示原始數據B，第二個圖顯示中心數據(look at the scale)。

B.標準化or歸一化

標準化用于將所有特征放在相同的比例上。每個零中心尺寸除以其標準偏差。

其中X'是標準化的數據集，X是原始數據集，x? 是的平均值X， 并且σ是標準偏差X。

def standardize(X): newX = center(X)/np.std(X, axis = 0) return newX

讓我們創建另一個具有不同scale 的數據集來檢查它是否工作。

np.random.ed(1234)c1 = np.random.normal(3, 1, 300)c2 = c1 + np.random.normal(7, 5, 300)/2.C = np.array([c1, c2]).T plotDataAndCov(C)plt.xlim(0, 15)plt.ylim(0, 15)plt.show()plt.clo()

Covariance matrix:

[[ 0.95171641 0.83976242]

[ 0.83976242 6.22529922]]

我們可以看到x和y的尺度不同。另請注意，由于scale 差異，相關性似乎較小。現在讓我們標準化它：

CStandardized = standardize(C) plotDataAndCov(CStandardized)plt.show()plt.clo()

Covariance matrix:

[[ 1. 0.34500274]

[ 0.34500274 1. ]]

看起來不錯。您可以看到scales 相同，并且數據集根據兩個軸以零為中心。

現在，看一下協方差矩陣。您可以看到每個坐標的方差 - 左上角單元格和右下角單元格 - 等于1。

這個新的協方差矩陣實際上是相關矩陣。兩個變量（c1和c2）之間的Pearson相關系數是0.54220151。

C.白化

白化，或球化，數據意味著我們要轉換它，使它有一個協方差矩陣它是單位矩陣對角線上是1，其他單元格是0。它被稱為白噪聲。

白化比其他預處理要復雜一些，但我們現在擁有了所需的所有工具。它涉及以下步驟：

Zero-center the data Decorrelate the data Rescale the data

讓我們再次采取C并嘗試執行這些步驟。

1.Zero-centering

這指的是均值歸一化（2. A）。請查看有關該center()函數的詳細信息。

CCentered = center(C) plotDataAndCov(CCentered)plt.show()plt.clo()

Covariance matrix:

[[ 0.95171641 0.83976242]

[ 0.83976242 6.22529922]]

2. Decorrelate

此時，我們需要對數據進行去相關。直覺上，這意味著我們想要旋轉數據直到不再存在相關性。看下面的圖片，

左邊的圖顯示了相關的數據。例如，如果取一個數據點的x值很大，那么相關聯的y也會很大。

現在取所有數據點，做一個旋轉(大概是逆時針方向45度左右)。右邊繪制的新數據不再相關。可以看到，大小y值與相同類型的x值相關。

問題是:我們如何找到正確的旋轉來得到不相關的數據?

實際上，它正是協方差矩陣的特征向量所做的。它們表示數據擴散最大的方向:

協方差矩陣的特征向量為您提供最大化方差的方向。綠線的方向是方差最大的地方。只要看看這條線上預測的最小和最大點 - 差距很大。將其與橙色線上的projection相比- spread 非常小。

因此，我們可以通過使用特征向量projection 來解相關數據。這將具有應用所需旋轉并移除尺寸之間的相關性的效果。以下是步驟：

計算協方差矩陣計算協方差矩陣的特征向量將特征向量矩陣應用于數據 - 這將應用旋轉

讓我們將其打包到函數中：

def decorrelate(X): newX = center(X) cov = X.T.dot(X)/float(X.shape[0]) # Calculate the eigenvalues and eigenvectors of the covariance matrix eigVals, eigVecs = np.linalg.eig(cov) # Apply the eigenvectors to X decorrelated = X.dot(eigVecs) return decorrelated

讓我們試著去掉我們的zero-centered的矩陣C來看它的實際效果：

plotDataAndCov(C)plt.show()plt.clo() CDecorrelated = decorrelate(CCentered)plotDataAndCov(CDecorrelated)plt.xlim(-5,5)plt.ylim(-5,5)plt.show()plt.clo()

Covariance matrix:

[[ 0.95171641 0.83976242]

[ 0.83976242 6.22529922]]

Covariance matrix:

[[ 5.96126981e-01 -1.48029737e-16]

[ -1.48029737e-16 3.15205774e+00]]

這很有效。

我們可以看到相關性不再存在了。協方差矩陣，現在是對角矩陣，證實了兩個維度之間的協方差等于0。

3.Rescale the data

下一步是縮放不相關矩陣以獲得對應于單位矩陣的協方差矩陣。為此，我們通過將每個維度除以其對應特征值的平方根來縮放我們的去相關數據。

def whiten(X): newX = center(X) cov = X.T.dot(X)/float(X.shape[0]) # Calculate the eigenvalues and eigenvectors of the covariance matrix eigVals, eigVecs = np.linalg.eig(cov) # Apply the eigenvectors to X decorrelated = X.dot(eigVecs) # Rescale the decorrelated data whitened = decorrelated / np.sqrt(eigVals + 1e-5)return whitened

注意：我們添加一個小值（此處為10 ^ -5）以避免除以0。

CWhitened = whiten(CCentered) plotDataAndCov(CWhitened)plt.xlim(-5,5)plt.ylim(-5,5)plt.show()plt.clo()

Covariance matrix:

[[ 9.99983225e-01 -1.06581410e-16]

[ -1.06581410e-16 9.99996827e-01]]

我們可以看到，使用協方差矩陣，這一切都很好。我們的東西看起來像一個單位矩陣 - 對角線上1，其他地方0。

3.圖像白化

我們將看到如何應用白化來預處理圖像數據集。為此，我們將使用Pal＆Sudeep（2016）的論文(https://ieeexplore.ieee.org/document/7808140/)，其中提供了有關該過程的一些細節。這種預處理技術稱為零分量分析（ZCA）。

查看論文，但這是他們得到的結果。顯示原始圖像（左）和ZCA（右）后的圖像。

從CIFAR10數據集中白化圖像。Pal＆Sudeep（2016）的論文結果 。

首先要做的事情。我們將從CIFAR數據集加載圖像。此數據集可從Keras獲得，您也可以在此處下載(https://www.cs.toronto.edu/~kriz/cifar.html)。

from keras.datats import cifar10 (X_train, y_train), (X_test, y_test) = cifar10.load_data() X_train.shape

(50000, 32, 32, 3)

CIFAR10數據集的訓練集包含50000個圖像。X_train是(50000, 32, 32, 3)。每幅圖像為32像素×32像素，每個像素包含3個維度（R，G，B）。每個值是0到255之間相應顏色的亮度。

我們將從僅選擇圖像的一個子集開始，假設為1000：

X = X_train[:1000]print X.shape

(1000,32,32,3)

現在我們將對陣列進行reshape ，以獲得每行一個圖像的平面圖像數據。每個圖像都是(1, 3072)因為32 x 32 x 3 = 3072.因此，包含所有圖像的數組將是(1000, 3072)：

X = X.reshape(X.shape[0], X.shape[1]*X.shape[2]*X.shape[3])print X.shape

(1000,3072)

下一步是能夠看到圖像。Matplotlib中的函數imshow()可用于顯示圖像。它需要具有shape的圖像(M x N x 3)，所以讓我們創建一個函數來reshape圖像并能夠從shape(1,3072)中可視化它們。

def plotImage(X): plt.figure(figsize=(1.5, 1.5)) plt.imshow(X.reshape(32,32,3)) plt.show() plt.clo()

例如，讓我們繪制我們加載的圖像之一：

plotImage(X[12, :])

我們現在可以實現圖像的白化。Pal＆Sudeep（2016）描述了這個過程：

1.第一步是通過除以255（像素的最大值）來重新縮放圖像以獲得范圍[0,1]。

回想一下，獲得范圍[0,1]的公式是：

但是，這里，最小值為0，因此這導致：

X_norm = X / 255.print 'X.min()', X_norm.min()print 'X.max()', X_norm.max()

X.min() 0.0

X.max() 1.0

2.從所有圖像中減去均值。

一種方法是拍攝每張圖像并從每個像素中刪除該圖像的平均值（Jarrett等，2009）。這個過程背后的直覺是它將每個圖像的像素集中在0左右。

另一種方法是為每個圖像拍攝我們擁有的3072個像素（RGB，32乘32像素），并在所有圖像中減去該像素的平均值。根據所有圖像，每個像素將以0為中心。當您使用圖像為網絡提供信息時，每個像素都被視為不同的特征。通過每像素平均減法，我們將每個特征（像素）集中在0左右。這種技術是常用的。

我們現在將從1000張圖像中進行每像素平均減法。我們的數據按這些維度進行組織(images, pixels)。這是(1000, 3072)因為有1000張圖像，32 x 32 x 3 = 3072像素。因此，平均每像素可以從第一軸得到:

X_norm.mean（axis = 0）.shape

(3072,)

這就得到3072個值，也就是均值的數量每像素1。

X_norm.mean(axis=0)

array([ 0.5234 , 0.54323137, 0.5274 , …, 0.50369804,

0.50011765, 0.45227451])

這接近0.5，因為我們已經歸一化到范圍[0,1]。但是，我們仍然需要從每個像素中刪除均值：

X_norm = X_norm - X_norm.mean(axis=0)

為了說服自己它有效，我們將計算第一個像素的平均值。我們希望它是0。

X_norm.mean(axis=0)

array([ -5.30575583e-16, -5.98021632e-16, -4.23439062e-16, …,

-1.81965554e-16, -2.49800181e-16, 3.98570066e-17])

這不完全是0，但它足夠小，我們可以認為它有效！

現在我們想要計算零中心數據的協方差矩陣。就像我們在上面看到的那樣，我們可以使用NumPy中的np.cov()函數來計算它。

請注意，我們的變量是不同的圖像。這意味著變量是矩陣x的行。為了明確一點，我們將用參數rowvar=TRUE將這個信息告訴NumPy，即使它默認為TRUE：

cov = np.cov(X_norm, rowvar=True)

現在是神奇的部分 - 我們將計算協方差矩陣的奇異值和向量，并使用它們來旋轉我們的數據集。

注意：大量圖像可能需要一些時間，這就是我們僅使用1000的原因。在論文中，他們使用了10000張圖像。您可以根據您使用的圖像數量來比較結果：

U,S,V = np.linalg.svd(cov)

在論文中，他們使用以下等式：

其中U為左奇異向量，S為圖像初始歸一化數據集協方差的奇異值，X為歸一化數據集。?hyper-parameter叫做白化系數。diag(a)對應于一個矩陣，向量a是對角線，其他單元格中為0。

我們將嘗試實現這個方程。讓我們從檢查SVD的維數開始:

print U.shape, S.shape

(1000, 1000) (1000,)

S是一個包含1000個元素(奇異值)的向量。因此，diag(S)的shape為(1000,1000)，以S為對角線:

print np.diag(S)print ' shape:', np.diag(S).shape

[[ 8.15846654e+00 0.00000000e+00 0.00000000e+00 …, 0.00000000e+00

0.00000000e+00 0.00000000e+00]

[ 0.00000000e+00 4.68234845e+00 0.00000000e+00 …, 0.00000000e+00

0.00000000e+00 0.00000000e+00]

[ 0.00000000e+00 0.00000000e+00 2.41075267e+00 …, 0.00000000e+00

0.00000000e+00 0.00000000e+00]

…,

[ 0.00000000e+00 0.00000000e+00 0.00000000e+00 …, 3.92727365e-05

0.00000000e+00 0.00000000e+00]

[ 0.00000000e+00 0.00000000e+00 0.00000000e+00 …, 0.00000000e+00

3.52614473e-05 0.00000000e+00]

[ 0.00000000e+00 0.00000000e+00 0.00000000e+00 …, 0.00000000e+00

0.00000000e+00 1.35907202e-15]]

shape: (1000, 1000)

檢查這部分：

shape(1000、1000) 。我們也看到X的shape是(1000 3072)X_ZCA的形狀為:

它對應于初始數據集的shape。

epsilon = 0.1X_ZCA = U.dot(np.diag(1.0/np.sqrt(S + epsilon))).dot(U.T).dot(X_norm) plotImage(X[12, :])plotImage(X_ZCA[12, :])

令人失望！如果你看一下這篇論文，這不是他們所展示的那種結果。實際上，這是因為我們沒有重新調整像素，并且存在負值。為此，我們可以使用與上述相同的技術將其放回[0,1]范圍內：

X_ZCA_rescaled = (X_ZCA - X_ZCA.min()) / (X_ZCA.max() - X_ZCA.min())print 'min:', X_ZCA_rescaled.min()print 'max:', X_ZCA_rescaled.max()

min: 0.0

max: 1.0

plotImage(X[12, :])plotImage(X_ZCA_rescaled[12, :])

很棒！它看起來像論文上的圖像。如前所述，他們使用了10000張圖像而不是像我們這樣的1000張。

要根據您使用的圖像數量和超參數ε的效果來查看結果的差異，以下是不同值的結果：

根據我們使用的圖像的數量和超參數ε的值，白化的結果是不同的。左側的圖像是原始圖像。在論文中，Pal＆Sudeep（2016）使用了10000個圖像，epsilon = 0.1。這對應于左下圖像。