R语言与主成分分析
學習筆記;
學習書目:《統計學:從數據到結論》–吳喜之;
具體原理部分請移步:主成分分析
主成分分析
主成分分析是利用降維的思想,在損失很少信息的前提下,把多個指標轉化為幾個綜合指標的多元統計方法。實際上,主成分分析可以說是因子分析的一個特例。這兩種方法都是尋找眾多相關變量的少數代表,這些代表變量、又稱為成分或因子,都是原先變量的線性組合。由于代表變量的數目顯著地小于原先變量的數目,數據的維數也就降低了。
主成分分析數學較簡單,發展也較早,因子分析需要的數學假定較多,理論稍微有些復雜但結果可能會比主成分分析更理想。
為了直觀地描述主成分分析降維的過程,我們先假定原先數據只是兩個變量的觀測值,即二維數據。如果這兩個變量分別由橫軸和縱軸所代表,每個觀測值都有相應于這兩個坐標軸的兩個坐標值,也就是這個二維坐標系中的一個點。如果這些數據點形成一個有橢圓形輪廓的點陣,那么這個橢圓有一個長軸和一個短軸,稱為主軸。主軸之間是互相垂直的,在短軸方向上數據變化較小,在長軸方向上數據變化較大。如果兩個坐標軸和橢圓的長短軸平行,那么代表長軸的變量就描述了數據的主要變化,而代表短軸的變量就描述了數據的次要變化。
但是,坐標軸通常并不和橢圓的長短軸平行。因此并進行變換,使得新變量和橢圓的長短軸平行。如果長軸變量代表了數據包含的大部分信息,就用該變量代替原先的兩個變量,并舍去次要的短軸變量,降維就完成了。在極端的情況,短軸如果退化成一點,那只用長軸變量就能夠完全解釋這些點的變化。這樣,由二維到一維的降維就自然完成了。橢圓的長短軸相差得越大,降維也越有道理。
多維變量的情況和二維類似。也有高維的橢球,只不過無法直觀地看見罷了。首先,把高維橢球的各個主軸找出來,再用代表大多數數據信息的最長的幾個軸作為新變量。這樣,主成分分析就基本完成了。注意,和二維的情況類似,高維橢球的主軸也是互相垂直的,這些互相正交的新變量是原先變量的線性組合,叫做主成分。
那么我們該如何選取多少個主成分呢?
我們用主抽長度比例,來說明這件事。我們選取的主成分所代表的主軸的長度之和應該占主軸長度總和的大部分。那么這個大部分到底是多大一部分呢?有的文獻建議,所選的主軸總長度占所有主軸長度之和大約80%(也有的建議85%)即可。但我們具體要選取幾個,還是要看實際情祝而定。如果所有涉及的變量都不那么相關,就很難降維,不相關的變量就只有自己能代表自己。
R語言實現
現在我有全國31個城市的廢水排放總量(FSPF)、一般工業固體廢物產生量(GYGTFW)和廢氣排放量(FQPF)數據,我想對其進行主成分分析。
首先,我們計算數據集的相關矩陣的特征值和特征向量:
> cor_city <- cor(city_data) > cor_cityFSPF GYGTFW FQPF FSPF 1.0000000 0.1233397 0.4804970 GYGTFW 0.1233397 1.0000000 0.8076458 FQPF 0.4804970 0.8076458 1.0000000計算相關矩陣的特征值以及特征向量:
> (b <- eigen(cor_city)) eigen() decomposition $values [1] 1.9971866 0.8920420 0.1107714$vectors[,1] [,2] [,3] [1,] -0.4050807 0.86138173 0.3064818 [2,] -0.6051700 -0.50388584 0.6163346 [3,] -0.6853312 -0.06419163 -0.7253968列出特征值、貢獻率以及累計貢獻率:
> data.frame(Eigenvalue = b$values, Contribution = b$values/sum(b$values),CC = cumsum(b$values)/sum(b$values))Eigenvalue Contribution CC 1 1.9971866 0.66572887 0.6657289 2 0.8920420 0.29734734 0.9630762 3 0.1107714 0.03692379 1.0000000碎石圖和累計貢獻率圖: par(mfrow = c(1, 2)) plot(b$values, type = 'o', ylab = 'Eigenvalue') plot(cumsum(b$values)/sum(b$values), type = 'o', ylab = 'CC')
因為我們的變量較少,所以圖形呈現出來的效果非常不咋地,但我還是貼一下:
因子載荷矩陣:
> (loadings = sweep(b$vectors, 2, sqrt(b$values), '*'))[,1] [,2] [,3] [1,] -0.5724676 0.81355762 0.1020043 [2,] -0.8552375 -0.47590998 0.2051305 [3,] -0.9685228 -0.06062769 -0.2414290各個省市在第一和第二主成分上的主成分得分圖(原始數據需要標準化):
> city_data_scale <-as.matrix(scale(city_data)) > plot(city_data_scale %*% b$vectors[,1:2], type='n', xlab = '第一主成分', ylab = '第二主成分') > text(city_data_scale %*% b$vectors[,1:2], row.names(city_data))圖像:
總結
- 上一篇: 腾达 AC18 无线路由器打印机服务使用
- 下一篇: 小白的算法初识课堂(part1)--二分