1. 程式人生 > >高斯分佈中均值,方差,協方差的計算及matlab實現

高斯分佈中均值,方差,協方差的計算及matlab實現

今天看論文的時候又看到了協方差矩陣這個破東西,以前看模式分類的時候就特困擾,沒想到現在還是搞不清楚,索性開始查協方差矩陣的資料,惡補之後決定馬上記錄下來,嘿嘿~本文我將用自認為循序漸進的方式談談協方差矩陣。

統計學的基本概念

學過概率統計的孩子都知道,統計裡最基本的概念就是樣本的均值,方差,或者再加個標準差。首先我們給你一個含有n個樣本的集合,依次給出這些概念的公式描述,這些高中學過數學的孩子都應該知道吧,一帶而過。

協方差矩陣的概念及matlab計算
均值:協方差矩陣的概念及matlab計算


標準差:協方差矩陣的概念及matlab計算


方差:協方差矩陣的概念及matlab計算

很 顯然,均值描述的是樣本集合的中間點,它告訴我們的資訊是很有限的,而標準差給我們描述的則是樣本集合的各個樣本點到均值的距離之平均。以這兩個集合為 例,[0,8,12,20]和[8,9,11,12],兩個集合的均值都是10,但顯然兩個集合差別是很大的,計算兩者的標準差,前者是8.3,後者是 1.8,顯然後者較為集中,故其標準差小一些,標準差描述的就是這種“散佈度”。之所以除以n-1而不是除以n,是因為這樣能使我們以較小的樣本集更好的 逼近總體的標準差,即統計上所謂的“無偏估計”。而方差則僅僅是標準差的平方。

為什麼需要協方差?

上 面幾個統計量看似已經描述的差不多了,但我們應該注意到,標準差和方差一般是用來描述一維資料的,但現實生活我們常常遇到含有多維資料的資料集,最簡單的 大家上學時免不了要統計多個學科的考試成績。面對這樣的資料集,我們當然可以按照每一維獨立的計算其方差,但是通常我們還想了解更多,比如,一個男孩子的 猥瑣程度跟他受女孩子歡迎程度是否存在一些聯絡啊,嘿嘿~協方差就是這樣一種用來度量兩個隨機變數關係的統計量,我們可以仿照方差的定義:

協方差矩陣的概念及matlab計算

來度量各個維度偏離其均值的程度,標準差可以這麼來定義:

協方差矩陣的概念及matlab計算

協方差的結果有什麼意義呢?如果結果為正值,則說明兩者是正相關的(從協方差可以引出“相關係數”的定義),也就是說一個人越猥瑣就越受女孩子歡迎,嘿嘿,那必須的~結果為負值就說明負相關的,越猥瑣女孩子越討厭,可能嗎?如果為0,也是就是統計上說的“相互獨立”。

從協方差的定義上我們也可以看出一些顯而易見的性質,如:

協方差矩陣的概念及matlab計算

協方差多了就是協方差矩陣

上一節提到的猥瑣和受歡迎的問題是典型二維問題,而協方差也只能處理二維問題,那維數多了自然就需要計算

協方差矩陣的概念及matlab計算

多個協方差,比如n維的資料集就需要計算個協方差,那自然而然的我們會想到使用矩陣來組織這些資料。給出協方差矩陣的定義:

協方差矩陣的概念及matlab計算

這個定義還是很容易理解的,我們可以舉一個簡單的三維的例子,假設資料集有三個維度,則協方差矩陣為

協方差矩陣的概念及matlab計算
可見,協方差矩陣是一個對稱的矩陣,而且對角線是各個維度上的方差。

Matlab協方差實戰

上面涉及的內容都比較容易,協方差矩陣似乎也很簡單,但實戰起來就很容易讓人迷茫了。必須要明確一點,協方差矩陣計算的是不同維度之間的協方差,而不是不同樣本之間的。

這個我將結合下面的例子說明,以下的演示將使用Matlab,為了說明計算原理,不直接呼叫Matlab的cov函式(藍色部分為Matlab程式碼)。

首先,隨機產生一個10*3維的整數矩陣作為樣本集,10為樣本的個數,3為樣本的維數。
MySample = fix(rand(10,3)*50)

協方差矩陣的概念及matlab計算

根據公式,計算協方差需要計算均值,那是按行計算均值還是按列呢,我一開始就老是困擾這個問題。前面我們也特別強調了,協方差矩陣是計算不同維度間的協方差,要時刻牢記這一點。樣本矩陣的每行是一個樣本,每列為一個維度,所以我們要按列計算均值。為了描述方便,我們先將三個維度的資料分別賦值:

dim1 = MySample(:,1);
dim2 = MySample(:,2);
dim3 = MySample(:,3);

計算dim1與dim2,dim1與dim3,dim2與dim3的協方差:

sum( (dim1-mean(dim1)) .* (dim2-mean(dim2)) ) / ( size(MySample,1)-1 ) % 得到  74.5333
sum( (dim1-mean(dim1)) .* (dim3-mean(dim3)) ) / ( size(MySample,1)-1 )
 % 得到  -10.0889
sum( (dim2-mean(dim2)) .* (dim3-mean(dim3)) ) / ( size(MySample,1)-1 )
 % 得到  -10***000

搞清楚了這個後面就容易多了,協方差矩陣的對角線就是各個維度上的方差,下面我們依次計算:

std(dim1)^2 % 得到   108.3222
std(dim2)^2
 % 得到   260.6222
std(dim3)^2
 % 得到  94.1778

這樣,我們就得到了計算協方差矩陣所需要的所有資料,呼叫Matlab自帶的cov函式進行驗證:

cov(MySample)

協方差矩陣的概念及matlab計算

把我們計算的資料對號入座,是不是一摸一樣?

Update: 今天突然發現,原來協方差矩陣還可以這樣計算,先讓樣本矩陣中心化,即每一維度減去該維度的均值,使每一維度上的均值為0,然後直接用新的到的樣本矩陣乘 上它的轉置,然後除以(N-1)即可。其實這種方法也是由前面的公式通道而來,只不過理解起來不是很直觀,但在抽象的公式推導時還是很常用的!同樣給出 Matlab程式碼實現:

X = MySample – repmat(mean(MySample),10,1);    % 中心化樣本矩陣,使各維度均值為0
C = (X’*X)./(size(X,1)-1)

總結

理解協方差矩陣的關鍵就在於牢記它計算的是不同維度之間的協方差,而不是不同樣本之間,拿到一個樣本矩陣,我們最先要明確的就是一行是一個樣本還是一個維度,心中明確這個整個計算過程就會順流而下,這麼一來就不會迷茫了~

P.S.寫論文要選Latex,在wordpress裡編輯公式還得用Latex,用Latex還真對得起咱學計算機這張臉~

例2:

=-112-231403for i=1:size(a,2
    
for j=1:size(a,2
        c(i,j)
=sum((a(:,i)-mean(a(:,i))).*(a(:,j)-mean(a(:,j))))/(size(a,1)-1);
    end 
end

c =

   10.3333   -4.1667    3.0000

   -4.1667    2.3333   -1.5000

    3.0000   -1.5000    1.0000 

 c為求得的協方差矩陣,在matlab以矩陣a的每一列為變數,對應的每一行為樣本。這樣在矩陣a中就有3個列變數分別為a(:,1), a(:,2), a(:,3)。

 在協方差矩陣c中,每一個元素c(i,j)為對第i列與第j列的協方差,例如c(1,2) = -4.1667為第一列與第二列的協方差。

 拿c(1,2)的求解過程來說

 c(1,2)=sum((a(:,1)-mean(a(:,1))).*(a(:,2)-mean(a(:,2))))/(size(a,1)-1);

 1. a(:,1)-mean(a(:,1)),第一列的元素減去該列的均值得到

   -1.3333

   -2.3333

    3.6667

2,  a(:,2)-mean(a(:,2)),第二列的元素減去該列的均值得到

   -0.3333

    1.6667

   -1.3333

3, 再將第一步與第二部的結果相乘

   -1.3333        -0.3333           0.4444

   -2.3333  .*     1.6667  =     -3.8889

    3.6667         -1.3333          -4.8889

4, 再將結果求和/size(a,1)-1 得 -4.1667,該值即為c(1,2)的值。

再細看一下是不是與協方差公式:Cov(X,Y) = E{ [ (X-E(X) ] [ (Y-E(Y) ] } 過程基本一致呢,只是在第4步的時候matlab做了稍微的調整,自由度為n-1,減少了一行的樣本值個數。