本次的操作只要的基于《Exploratory.Multivariate.Analysis.by.Example.Using.R.2nd.Edition.2017》。本书的第一张是使用PCA对多维度的数据进行数据分析。基本的PCA的原理。网上自己搜一下就行。这里主要还是联系书中提到的例子 ### import necessary data 首先需要加载需要的包和数据
if (!require("FactoMineR")){ install.packages("FactoMineR")}
## Loading required package: FactoMineR
library(FactoMineR)
该数据集包括了雅典奥运会上十项全能比赛的结果。 对于这两项比赛,每个运动员都可获得以下信息:10场比赛的每场比赛的表现,总分数(每场比赛,运动员根据表现获得积分;此处得分的总和),以及最终比赛排名(见表1.13)。事件按以下顺序进行:100米,跳远,铅球,跳高,400米(第一天)和110米栏,铁饼,撑杆跳,标枪,1500米(第二天) 在此表中,行对应于运动员在比赛时的表现情况,每列对应于描述运动员在会议期间表现的变量。有12个定量变量(10个事件的结果,运动员的排名,以及获得的总分数)和一个分类变量(运动员参与的比赛)。
data("decathlon")
dim(decathlon)
## [1] 41 13
head(decathlon)
## 100m Long.jump Shot.put High.jump 400m 110m.hurdle Discus
## SEBRLE 11.04 7.58 14.83 2.07 49.81 14.69 43.75
## CLAY 10.76 7.40 14.26 1.86 49.37 14.05 50.72
## KARPOV 11.02 7.30 14.77 2.04 48.37 14.09 48.95
## BERNARD 11.02 7.23 14.25 1.92 48.93 14.99 40.87
## YURKOV 11.34 7.09 15.19 2.10 50.42 15.31 46.26
## WARNERS 11.11 7.60 14.31 1.98 48.68 14.23 41.10
## Pole.vault Javeline 1500m Rank Points Competition
## SEBRLE 5.02 63.19 291.7 1 8217 Decastar
## CLAY 4.92 60.15 301.5 2 8122 Decastar
## KARPOV 4.92 50.31 300.2 3 8099 Decastar
## BERNARD 5.32 62.77 280.1 4 8067 Decastar
## YURKOV 4.72 63.44 276.4 5 8036 Decastar
## WARNERS 4.92 51.77 278.1 6 8030 Decastar
由于我们要研究不同变量之运动员之间的相似性,因此我们需要纳入的变量都是连续性的变量。至于其他的比如等级变量,分类变量等都会被当作补充材料。
对于不同的标量之间是否应该进行标准化处理。答案是肯定的。因为如果不进行标准化的处理。不同的变量之间的由于单位的不同,其标准误也会产生很大的区别。
这里我们使用FactoMineR包的PCA函数进行分析。 这个函数主要的参数就是数据集,其他一些包括数据集是否需要标准化scale.unit(默认为T);数据集中是否作为补充的的列quanti.sup以及数据集中是分类变量的列quali.sup graph可以制定在分析的时候是否出图。 - 我们对数据进行PCA分析
res.pca <- PCA(decathlon, quanti.sup = 11:12, quali.sup = 13)
PCA函数运行后,我们还可以得到两个图,一个是各个个体的分布图,一个是各个变量的分布图。
summary.PCA(res.pca)
##
## Call:
## PCA(X = decathlon, quanti.sup = 11:12, quali.sup = 13)
##
##
## Eigenvalues
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5 Dim.6
## Variance 3.272 1.737 1.405 1.057 0.685 0.599
## % of var. 32.719 17.371 14.049 10.569 6.848 5.993
## Cumulative % of var. 32.719 50.090 64.140 74.708 81.556 87.548
## Dim.7 Dim.8 Dim.9 Dim.10
## Variance 0.451 0.397 0.215 0.182
## % of var. 4.512 3.969 2.148 1.822
## Cumulative % of var. 92.061 96.030 98.178 100.000
##
## Individuals (the 10 first)
## Dist Dim.1 ctr cos2 Dim.2 ctr cos2
## SEBRLE | 2.369 | 0.792 0.467 0.112 | 0.772 0.836 0.106 |
## CLAY | 3.507 | 1.235 1.137 0.124 | 0.575 0.464 0.027 |
## KARPOV | 3.396 | 1.358 1.375 0.160 | 0.484 0.329 0.020 |
## BERNARD | 2.763 | -0.610 0.277 0.049 | -0.875 1.074 0.100 |
## YURKOV | 3.018 | -0.586 0.256 0.038 | 2.131 6.376 0.499 |
## WARNERS | 2.428 | 0.357 0.095 0.022 | -1.685 3.986 0.482 |
## ZSIVOCZKY | 2.563 | 0.272 0.055 0.011 | -1.094 1.680 0.182 |
## McMULLEN | 2.561 | 0.588 0.257 0.053 | 0.231 0.075 0.008 |
## MARTINEAU | 3.742 | -1.995 2.968 0.284 | 0.561 0.442 0.022 |
## HERNU | 2.794 | -1.546 1.782 0.306 | 0.488 0.335 0.031 |
## Dim.3 ctr cos2
## SEBRLE 0.827 1.187 0.122 |
## CLAY 2.141 7.960 0.373 |
## KARPOV 1.956 6.644 0.332 |
## BERNARD 0.890 1.375 0.104 |
## YURKOV -1.225 2.606 0.165 |
## WARNERS 0.767 1.020 0.100 |
## ZSIVOCZKY -1.283 2.857 0.250 |
## McMULLEN -0.418 0.303 0.027 |
## MARTINEAU -0.730 0.925 0.038 |
## HERNU 0.841 1.227 0.091 |
##
## Variables
## Dim.1 ctr cos2 Dim.2 ctr cos2 Dim.3 ctr
## 100m | -0.775 18.344 0.600 | 0.187 2.016 0.035 | -0.184 2.420
## Long.jump | 0.742 16.822 0.550 | -0.345 6.869 0.119 | 0.182 2.363
## Shot.put | 0.623 11.844 0.388 | 0.598 20.607 0.358 | -0.023 0.039
## High.jump | 0.572 9.998 0.327 | 0.350 7.064 0.123 | -0.260 4.794
## 400m | -0.680 14.116 0.462 | 0.569 18.666 0.324 | 0.131 1.230
## 110m.hurdle | -0.746 17.020 0.557 | 0.229 3.013 0.052 | -0.093 0.611
## Discus | 0.552 9.328 0.305 | 0.606 21.162 0.368 | 0.043 0.131
## Pole.vault | 0.050 0.077 0.003 | -0.180 1.873 0.033 | 0.692 34.061
## Javeline | 0.277 2.347 0.077 | 0.317 5.784 0.100 | -0.390 10.807
## 1500m | -0.058 0.103 0.003 | 0.474 12.946 0.225 | 0.782 43.543
## cos2
## 100m 0.034 |
## Long.jump 0.033 |
## Shot.put 0.001 |
## High.jump 0.067 |
## 400m 0.017 |
## 110m.hurdle 0.009 |
## Discus 0.002 |
## Pole.vault 0.479 |
## Javeline 0.152 |
## 1500m 0.612 |
##
## Supplementary continuous variables
## Dim.1 cos2 Dim.2 cos2 Dim.3 cos2
## Rank | -0.671 0.450 | 0.051 0.003 | -0.058 0.003 |
## Points | 0.956 0.914 | -0.017 0.000 | -0.066 0.004 |
##
## Supplementary categories
## Dist Dim.1 cos2 v.test Dim.2 cos2 v.test
## Decastar | 0.946 | -0.600 0.403 -1.430 | -0.038 0.002 -0.123 |
## OlympicG | 0.439 | 0.279 0.403 1.430 | 0.017 0.002 0.123 |
## Dim.3 cos2 v.test
## Decastar 0.289 0.093 1.050 |
## OlympicG -0.134 0.093 -1.050 |
PCA分析的一个方面是我们可以看到变量之间是否是结构化的(变量之间是否存在关系)。同时还可以了解每个组分的解释程度。 res.pca$eig对应每个组分的相关联的特征值。处理单纯的观察。我们可以通过图形来观察。
head(res.pca$eig)
## eigenvalue percentage of variance cumulative percentage of variance
## comp 1 3.2719055 32.719055 32.71906
## comp 2 1.7371310 17.371310 50.09037
## comp 3 1.4049167 14.049167 64.13953
## comp 4 1.0568504 10.568504 74.70804
## comp 5 0.6847735 6.847735 81.55577
## comp 6 0.5992687 5.992687 87.54846
barplot(res.pca$eig[,1], names.arg=paste("dim",1:nrow(res.pca$eig)))
个体变异的图,在PCA函数时候的时候会输出出来。在图中我们可以看到作为补充的列也会画上去。如果要查看详细结果的可以可以在 res.pca$ind里面看到分析的结果。可以在res.pca$quali.sup里面看到作为补充的列的结果。 如果要绘制个体在其他成分之间的图的话。可以使用 plot或者plog.PCA函数。在函数当中,choix="ind"需要制定。axes = 3:4代表想要绘制哪两个成分
plot(res.pca, choix = "ind", axes = c(1,3))
我们可以根据分类变量来绘制PCA个体图上增加不同的颜色。这个使用 habillage参数制定即可。
plot(res.pca,choix="ind",habillage=13,cex=0.7)
我们也可以绘制PCA个体图上的分类变量的椭圆。
plotellipses(res.pca,cex=0.8)
变量之间的相关性。在PCA函数当中也是默认出图的。其中活动变量是实线表示,补充变量是用虚线表示 图的横纵坐标代表和成分的相关性系数 整个图分为四个象限,其中都在同一个象限的代表变量之间存在相关性。在相对象限的是负相关。 这些结果同样的可以在 summary.PCA中看到。如果要做图的话。也是用plot函数,使用的是 choix="var"参数以及使用axes制定绘制的成分
plot(res.pca,choix="var")
我们可以使用dimdesc函数来绘制不同变量和各个成分之间的相关性。这个函数存在两个参数。axes制定返回的主成分的结果。proba具有差异的P值。默认的是返回1:3成分中P小于0.05的结果。 如果将 proba制定为0.2的话。则返回分类变量的相关性。
dimdesc(res.pca)
## $Dim.1
## $Dim.1$quanti
## correlation p.value
## Points 0.9561543 2.099191e-22
## Long.jump 0.7418997 2.849886e-08
## Shot.put 0.6225026 1.388321e-05
## High.jump 0.5719453 9.362285e-05
## Discus 0.5524665 1.802220e-04
## Rank -0.6705104 1.616348e-06
## 400m -0.6796099 1.028175e-06
## 110m.hurdle -0.7462453 2.136962e-08
## 100m -0.7747198 2.778467e-09
##
##
## $Dim.2
## $Dim.2$quanti
## correlation p.value
## Discus 0.6063134 2.650745e-05
## Shot.put 0.5983033 3.603567e-05
## 400m 0.5694378 1.020941e-04
## 1500m 0.4742238 1.734405e-03
## High.jump 0.3502936 2.475025e-02
## Javeline 0.3169891 4.344974e-02
## Long.jump -0.3454213 2.696969e-02
##
##
## $Dim.3
## $Dim.3$quanti
## correlation p.value
## 1500m 0.7821428 1.554450e-09
## Pole.vault 0.6917567 5.480172e-07
## Javeline -0.3896554 1.179331e-02
dimdesc(res.pca, proba = 0.2)
## $Dim.1
## $Dim.1$quanti
## correlation p.value
## Points 0.9561543 2.099191e-22
## Long.jump 0.7418997 2.849886e-08
## Shot.put 0.6225026 1.388321e-05
## High.jump 0.5719453 9.362285e-05
## Discus 0.5524665 1.802220e-04
## Javeline 0.2771108 7.942460e-02
## Rank -0.6705104 1.616348e-06
## 400m -0.6796099 1.028175e-06
## 110m.hurdle -0.7462453 2.136962e-08
## 100m -0.7747198 2.778467e-09
##
## $Dim.1$quali
## R2 p.value
## Competition 0.05110487 0.1552515
##
## $Dim.1$category
## Estimate p.value
## OlympicG 0.4393744 0.1552515
## Decastar -0.4393744 0.1552515
##
##
## $Dim.2
## $Dim.2$quanti
## correlation p.value
## Discus 0.6063134 2.650745e-05
## Shot.put 0.5983033 3.603567e-05
## 400m 0.5694378 1.020941e-04
## 1500m 0.4742238 1.734405e-03
## High.jump 0.3502936 2.475025e-02
## Javeline 0.3169891 4.344974e-02
## 110m.hurdle 0.2287933 1.501925e-01
## Long.jump -0.3454213 2.696969e-02
##
##
## $Dim.3
## $Dim.3$quanti
## correlation p.value
## 1500m 0.7821428 1.554450e-09
## Pole.vault 0.6917567 5.480172e-07
## High.jump -0.2595119 1.013160e-01
## Javeline -0.3896554 1.179331e-02
我们通过观察分类变量的相关性。可以说明分类变量之间在某一成分是否是正相关或者负相关的
我们可以比较不同成分之间在个体和变量的不同,来进行汇总分析。 需要注意的是其实个体和变量之间的关系。一方面只是数据的相似性,另外一方面是变量之间的相关性分析。因此有必要通过回顾性的数据来得到支持。 标准化的数据以及相关矩阵的均值和标准差可以通过下面的代码实现。
res.pca$call$centre ###均值
## [1] 10.998049 7.260000 14.477073 1.976829 49.616341 14.605854
## [7] 44.325610 4.762439 58.316585 279.024878
res.pca$call$ecart.type ###标准差
## [1] 0.25979560 0.31251927 0.81431175 0.08785906 1.13929751
## [6] 0.46599998 3.33639725 0.27458865 4.76759315 11.53001177
round(scale(decathlon[,1:12]),2)[1:5,1:5] ###标准化后的数据
## 100m Long.jump Shot.put High.jump 400m
## SEBRLE 0.16 1.01 0.43 1.05 0.17
## CLAY -0.91 0.44 -0.26 -1.31 -0.21
## KARPOV 0.08 0.13 0.36 0.71 -1.08
## BERNARD 0.08 -0.09 -0.28 -0.64 -0.60
## YURKOV 1.30 -0.54 0.86 1.38 0.70
round(cor(decathlon[,1:12]),2)[1:5,1:5] ###数据的相关分析
## 100m Long.jump Shot.put High.jump 400m
## 100m 1.00 -0.60 -0.36 -0.25 0.52
## Long.jump -0.60 1.00 0.18 0.29 -0.60
## Shot.put -0.36 0.18 1.00 0.49 -0.14
## High.jump -0.25 0.29 0.49 1.00 -0.19
## 400m 0.52 -0.60 -0.14 -0.19 1.00
###观察不同变量之间的相关性
pairs(decathlon[,c(1,2,6,10)])
在这个例子当中我们分析了欧洲主要国家的温度变化。数据中包括每个月主要城市的温度以及包括两个分类变量:经度和纬度。以及一个作为区域的变量。
城市的选择:我们在这个数据当中想要分析温度在各个国家当中的变化。每个国家使用他们的首都当作分析变量。对于不是首都的城市,就可以补充变量进行分析了。通过这样的分析,我可以知道 哪些国家之间的差距最大
变量的选择: 我们选择每个月温度的变量作为研究变量。
数据的标准化只有在数据是通过不同的标准检测的才需要进行处理。但是相同测量单位的数据并不代表不需要标准化。一月份的1摄氏度与七月份相同吗?未能使数据标准化意味着为每个变量赋予与其方差成比例的权重。必须注意的是,从一个月到另一个月,标准偏差很小(最多它们加倍)。因此,认为标准化对分析结果没有实际影响是合理的。从另一个角度来看,在计算城镇之间的距离时,未能使数据标准化意味着将相同的影响归因于1度的差异,无论月份如何。当标准化时,这种差异被放大,因为在一个城镇之间的温度变化很小的月份,这种差异变得更加明显。对于这个例子,决定标准化数据,从而将相同的权重归因于每个月。
基本的分析代码如下
temperature <- read.table("http://factominer.free.fr/bookV2/temperature.csv",
header=TRUE,sep=";",dec=".",row.names=1) ###下载数据
res<-PCA(temperature,ind.sup=24:35,quanti.sup=13:16,quali.sup=17) ###PCA分析
plot.PCA(res,choix="ind",habillage=17) ###绘制不同区域之间PCA的图
summary(res,nb.dec=2)
##
## Call:
## PCA(X = temperature, ind.sup = 24:35, quanti.sup = 13:16, quali.sup = 17)
##
##
## Eigenvalues
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5 Dim.6 Dim.7
## Variance 9.95 1.85 0.13 0.04 0.02 0.01 0.01
## % of var. 82.90 15.40 1.05 0.32 0.14 0.11 0.05
## Cumulative % of var. 82.90 98.29 99.35 99.67 99.81 99.91 99.96
## Dim.8 Dim.9 Dim.10 Dim.11 Dim.12
## Variance 0.00 0.00 0.00 0.00 0.00
## % of var. 0.02 0.01 0.01 0.00 0.00
## Cumulative % of var. 99.98 99.99 99.99 100.00 100.00
##
## Individuals (the 10 first)
## Dist Dim.1 ctr cos2 Dim.2 ctr cos2 Dim.3 ctr
## Amsterdam | 1.44 | 0.23 0.02 0.02 | -1.37 4.43 0.90 | -0.10 0.38
## Athens | 7.68 | 7.60 25.25 0.98 | 0.93 2.04 0.01 | 0.56 10.85
## Berlin | 0.50 | -0.29 0.04 0.33 | 0.02 0.00 0.00 | -0.29 2.91
## Brussels | 1.36 | 0.63 0.17 0.22 | -1.18 3.26 0.75 | -0.15 0.80
## Budapest | 2.45 | 1.67 1.22 0.46 | 1.71 6.90 0.49 | -0.50 8.57
## Copenhagen | 1.63 | -1.46 0.93 0.81 | -0.49 0.57 0.09 | 0.44 6.68
## Dublin | 2.74 | -0.51 0.11 0.03 | -2.67 16.82 0.96 | -0.18 1.10
## Elsinki | 4.13 | -4.04 7.12 0.96 | 0.46 0.50 0.01 | 0.59 12.12
## Kiev | 2.65 | -1.71 1.28 0.42 | 2.01 9.48 0.57 | -0.17 1.00
## Krakow | 1.57 | -1.26 0.69 0.65 | 0.87 1.80 0.31 | -0.27 2.58
## cos2
## Amsterdam 0.01 |
## Athens 0.01 |
## Berlin 0.33 |
## Brussels 0.01 |
## Budapest 0.04 |
## Copenhagen 0.07 |
## Dublin 0.00 |
## Elsinki 0.02 |
## Kiev 0.00 |
## Krakow 0.03 |
##
## Supplementary individuals (the 10 first)
## Dist Dim.1 cos2 Dim.2 cos2 Dim.3 cos2
## Antwerp | 1.33 | 0.59 0.20 | -1.16 0.76 | -0.13 0.01 |
## Barcelona | 6.01 | 5.99 0.99 | -0.38 0.00 | 0.25 0.00 |
## Bordeaux | 2.97 | 2.85 0.92 | -0.73 0.06 | -0.27 0.01 |
## Edinburgh | 2.69 | -1.29 0.23 | -2.34 0.76 | -0.18 0.00 |
## Frankfurt | 0.68 | 0.39 0.33 | 0.15 0.05 | -0.44 0.43 |
## Geneva | 0.54 | 0.32 0.36 | 0.17 0.10 | -0.26 0.22 |
## Genoa | 5.87 | 5.85 0.99 | -0.13 0.00 | 0.34 0.00 |
## Milan | 3.55 | 3.16 0.79 | 1.56 0.19 | -0.25 0.00 |
## Palermo | 7.41 | 7.29 0.97 | -0.21 0.00 | -0.59 0.01 |
## Seville | 7.87 | 7.86 1.00 | 0.26 0.00 | 0.22 0.00 |
##
## Variables (the 10 first)
## Dim.1 ctr cos2 Dim.2 ctr cos2 Dim.3 ctr cos2
## January | 0.84 7.13 0.71 | -0.53 15.28 0.28 | 0.07 3.64 0.00 |
## February | 0.88 7.86 0.78 | -0.46 11.25 0.21 | 0.00 0.01 0.00 |
## March | 0.95 8.98 0.89 | -0.29 4.47 0.08 | -0.12 11.56 0.01 |
## April | 0.97 9.53 0.95 | 0.10 0.54 0.01 | -0.20 31.13 0.04 |
## May | 0.87 7.61 0.76 | 0.46 11.34 0.21 | -0.16 19.30 0.02 |
## June | 0.83 6.98 0.69 | 0.55 16.09 0.30 | 0.05 1.94 0.00 |
## July | 0.84 7.16 0.71 | 0.51 14.00 0.26 | 0.15 18.71 0.02 |
## August | 0.91 8.31 0.83 | 0.40 8.74 0.16 | 0.09 6.06 0.01 |
## September | 0.99 9.77 0.97 | 0.15 1.26 0.02 | 0.02 0.41 0.00 |
## October | 0.99 9.88 0.98 | -0.08 0.39 0.01 | 0.00 0.00 0.00 |
##
## Supplementary continuous variables
## Dim.1 cos2 Dim.2 cos2 Dim.3 cos2
## Annual | 1.00 1.00 | -0.07 0.00 | 0.00 0.00 |
## Amplitude | -0.31 0.10 | 0.94 0.89 | 0.04 0.00 |
## Latitude | -0.91 0.83 | -0.22 0.05 | 0.18 0.03 |
## Longitude | -0.36 0.13 | 0.64 0.42 | -0.04 0.00 |
##
## Supplementary categories
## Dist Dim.1 cos2 v.test Dim.2 cos2 v.test Dim.3
## East | 1.78 | -1.10 0.38 -1.08 | 1.38 0.60 3.15 | -0.26
## North | 2.65 | -2.44 0.85 -2.40 | -0.99 0.14 -2.26 | 0.27
## South | 4.57 | 4.56 1.00 3.58 | 0.14 0.00 0.25 | 0.12
## West | 1.03 | 0.50 0.23 0.34 | -0.86 0.69 -1.36 | -0.16
## cos2 v.test
## East 0.02 -2.25 |
## North 0.01 2.35 |
## South 0.00 0.80 |
## West 0.03 -1.00 |
通过总结我们发现,第一组分可以解释82.9%的变异。第二成分可以解释15.4%的变异。通过两者成分我们就看可以解释变量之间的差异了。
dimdesc(res)
## $Dim.1
## $Dim.1$quanti
## correlation p.value
## Annual 0.9975483 9.575809e-26
## October 0.9916246 3.734359e-20
## September 0.9856254 1.056964e-17
## April 0.9738876 5.295583e-15
## November 0.9523567 2.659701e-12
## March 0.9450521 1.151422e-11
## August 0.9092443 1.904695e-09
## February 0.8842848 2.180620e-08
## December 0.8731191 5.450923e-08
## May 0.8698517 7.013036e-08
## July 0.8441626 4.129012e-07
## January 0.8424506 4.594275e-07
## June 0.8333141 7.957991e-07
## Latitude -0.9099106 1.768178e-09
##
## $Dim.1$quali
## R2 p.value
## Area 0.6787608 6.282082e-05
##
## $Dim.1$category
## Estimate p.value
## South 4.182744 2.359137e-05
## North -2.822709 1.241529e-02
##
##
## $Dim.2
## $Dim.2$quanti
## correlation p.value
## Amplitude 0.9444140 1.296159e-11
## Longitude 0.6449726 8.912273e-04
## June 0.5453220 7.119942e-03
## July 0.5086619 1.319151e-02
## May 0.4578116 2.804268e-02
## February -0.4558325 2.881655e-02
## December -0.4728656 2.268397e-02
## January -0.5313576 9.076631e-03
##
## $Dim.2$quali
## R2 p.value
## Area 0.5461533 0.001525752
##
## $Dim.2$category
## Estimate p.value
## East 1.4622991 0.0004492948
## North -0.9063514 0.0201613275
##
##
## $Dim.3
## $Dim.3$quali
## R2 p.value
## Area 0.3936374 0.02091359
##
## $Dim.3$category
## Estimate p.value
## North 0.2779913 0.01510325
## East -0.2478747 0.02080116