a=[-1,1,2;-2,3,1;4,0,9]

for 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，減少了一行的樣本值個數。

>> a=[-1,1,2;-2,3,1;4,0,9]

a =

    -1     1     2

    -2     3     1

     4     0     9

>> cov(a)

ans =

   10.3333   -4.1667   14.0000

   -4.1667    2.3333   -5.5000

   14.0000   -5.5000   19.0000

>> dim1=a(:,1)

dim1 =

    -1

    -2

     4

>> dim2=a(:,2)

dim2 =

     1

     3

     0

>> dim3=a(:,3)

dim3 =

     2

     1

     9

>> dim1-mean(dim1)

ans =

   -1.3333

   -2.3333

    3.6667

>> dim2-mean(dim2)

ans =

   -0.3333

    1.6667

   -1.3333

>> m=dim1-mean(dim1)

m =

   -1.3333

   -2.3333

    3.6667

>> m=dim2-mean(dim2)

m =

   -0.3333

    1.6667

   -1.3333

>> m=dim1-mean(dim1)

m =

   -1.3333

   -2.3333

    3.6667

>> n=dim2-mean(dim2)

n =

   -0.3333

    1.6667

   -1.3333

>> m.*n

ans =

    0.4444

   -3.8889

   -4.8889

>> size(a,1)

ans =

     3

>> sum(m.*n)

ans =

   -8.3333

>> ans/2

ans =

   -4.1667