By

文章目录
  1. 主成分和因子分析
    1. 主成分分析
      1. 1、主成分个数
      2. 2、提取主成分
      3. 3、主成分旋转
      4. 4、获取主成分得分
    2. 探索性因子分析
      1. 1、公因子数
      2. 2、提取公共因子
      3. 3、因子旋转
      4. 4、因子得分

主成分和因子分析

主成分分析

1、主成分个数

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
> ##主成分和因子分析
> 
> #主成分分析
> #PCA的目标是用一组较少的不相关变量代替大量相关变量，同时尽可能保留初始变量的信
> #息，这些推导所得的变量称为主成分，它们是观测变量的线性组合。如第一成分
> # PC1 = a1X1 + a2X2 + ... + akXk
> #它是k个观测变量的加权组合，对初始变量集的方差解释性最大。第二主成分也是初始变量
> #的线性组合，对方差的解释性排第二，同时与第一主成分正交（不相关）。
> #后面每一个主成分都最大化它对方差的解释程度，同时与之前所有的主成分都正交。
> #你可以选取与变量数相同的主成分，但从实用的角度来看，我们都希望能用较少的主成分来近似全变量集。
> library(psych)
> fa.parallel(USJudgeRatings[,-1], fa="pc", n.iter=100, show.legend = FALSE, main="Scree plot with parallel analysis")
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
Parallel analysis suggests that the number of factors =  NA  and the number of components =  1 
> #碎石图（直线与x符号）、特征值大于1准则（水平线）和100次模拟的平行分析（虚线）都表明保留一个主成分即可

2、提取主成分

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
> #提取主成分principal
> #factors设定主成分数（默认为1）；
> #rotate指定旋转的方法[默认最大方差旋转（varimax)]。
> #scores设定是否需要计算主成分得分（默认不需要）。
> pc <- principal(USJudgeRatings[,-1], nfactors = 1)
> pc
Principal Components Analysis
Call: principal(r = USJudgeRatings[, -1], nfactors = 1)
Standardized loadings (pattern matrix) based upon correlation matrix
      PC1   h2     u2 com
INTG 0.92 0.84 0.1565   1
DMNR 0.91 0.83 0.1663   1
DILG 0.97 0.94 0.0613   1
CFMG 0.96 0.93 0.0720   1
DECI 0.96 0.92 0.0763   1
PREP 0.98 0.97 0.0299   1
FAMI 0.98 0.95 0.0469   1
ORAL 1.00 0.99 0.0091   1
WRIT 0.99 0.98 0.0196   1
PHYS 0.89 0.80 0.2013   1
RTEN 0.99 0.97 0.0275   1
                 PC1
SS loadings    10.13
Proportion Var  0.92
Mean item complexity =  1
Test of the hypothesis that 1 component is sufficient.
The root mean square of the residuals (RMSR) is  0.04 
 with the empirical chi square  6.2  with prob <  1 
Fit based upon off diagonal values = 1
> #PC1栏包含了成分载荷,指观测变量与主成分的相关系数,第一主成分（PC1）与每个变量都高度相关，也就是说，它是一个可用来进行一般性评价
> #的维度
> #h2栏指成分公因子方差——主成分对每个变量的方差解释度。 u2栏指成分唯一性——方差无法被主成分解释的比例（1-h2）。
> #例如，体能（PHYS） 80%的方差都可用第一主成分来解释，20%不能。相比而言， PHYS是用第一主成分表示性最差的变量。
> 
> #SS loadings行包含了与主成分相关联的特征值，指的是与特定主成分相关联的标准化后的
> #方差值（本例中，第一主成分的值为10）。最后， Proportion Var行表示的是每个主成分对整
> #个数据集的解释程度。此处可以看到，第一主成分解释了11个变量92%的方差
> 
> #多个主成分例子
> library(psych)
> fa.parallel(Harman23.cor$cov, n.obs = 302, fa="pc", n.iter=100, 
+             show.legend = FALSE, main="Scree plot with parallel analysis")
Parallel analysis suggests that the number of factors =  NA  and the number of components =  2

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
> pc <- principal(Harman23.cor$cov, nfactors = 2, rotate="none")
> pc
Principal Components Analysis
Call: principal(r = Harman23.cor$cov, nfactors = 2, rotate = "none")
Standardized loadings (pattern matrix) based upon correlation matrix
                PC1   PC2   h2    u2 com
height         0.86 -0.37 0.88 0.123 1.4
arm.span       0.84 -0.44 0.90 0.097 1.5
forearm        0.81 -0.46 0.87 0.128 1.6
lower.leg      0.84 -0.40 0.86 0.139 1.4
weight         0.76  0.52 0.85 0.150 1.8
bitro.diameter 0.67  0.53 0.74 0.261 1.9
chest.girth    0.62  0.58 0.72 0.283 2.0
chest.width    0.67  0.42 0.62 0.375 1.7
                       PC1  PC2
SS loadings           4.67 1.77
Proportion Var        0.58 0.22
Cumulative Var        0.58 0.81
Proportion Explained  0.73 0.27
Cumulative Proportion 0.73 1.00
Mean item complexity =  1.7
Test of the hypothesis that 2 components are sufficient.
The root mean square of the residuals (RMSR) is  0.05 
Fit based upon off diagonal values = 0.99
> #第一主成分与每个身体测量指标都正相关，看起来似乎是
> #一个一般性的衡量因子；第二主成分与前四个变量（height、arm.span、forearm和lower.leg）
> #负相关，与后四个变量（weight、 bitro.diameter、 chest.girth和chest.width）正相关，
> #因此它看起来似乎是一个长度—容量因子。但理念上的东西都不容易构建，当提取了多个成分时，
> #对它们进行旋转可使结果更具解释性，接下来我们便讨论该问题。
>

3、主成分旋转

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
> #旋转
> #旋转是一系列将成分载荷阵变得更容易解释的数学方法，它们尽可能地对成分去噪。旋转方
> #法有两种：使选择的成分保持不相关（正交旋转），和让它们变得相关（斜交旋转）
> #最流行的正交旋转是方差极大旋转
> pc <- principal(Harman23.cor$cov, nfactors = 2, rotate="varimax")
> pc
Principal Components Analysis
Call: principal(r = Harman23.cor$cov, nfactors = 2, rotate = "varimax")
Standardized loadings (pattern matrix) based upon correlation matrix
                RC1  RC2   h2    u2 com
height         0.90 0.25 0.88 0.123 1.2
arm.span       0.93 0.19 0.90 0.097 1.1
forearm        0.92 0.16 0.87 0.128 1.1
lower.leg      0.90 0.22 0.86 0.139 1.1
weight         0.26 0.88 0.85 0.150 1.2
bitro.diameter 0.19 0.84 0.74 0.261 1.1
chest.girth    0.11 0.84 0.72 0.283 1.0
chest.width    0.26 0.75 0.62 0.375 1.2
                       RC1  RC2
SS loadings           3.52 2.92
Proportion Var        0.44 0.37
Cumulative Var        0.44 0.81
Proportion Explained  0.55 0.45
Cumulative Proportion 0.55 1.00
Mean item complexity =  1.1
Test of the hypothesis that 2 components are sufficient.
The root mean square of the residuals (RMSR) is  0.05 
Fit based upon off diagonal values = 0.99
> #列的名字都从PC变成了RC，以表示成分被旋转。观察RC1栏的载荷，你可以发现第一主成分
> #主要由前四个变量来解释（长度变量）。 RC2栏的载荷表示第二主成分主要由变量5到变量8来解
> #释（容量变量）
> #两个主成分旋转后的累积方差解释性没有变化（81%），变的只是各个主成分对
> #方差的解释度（成分1从58%变为44%，成分2从22%变为37%）。各成分的方差解释度趋同，准确
> #来说，此时应该称它们为成分而不是主成分（因为单个主成分方差最大化性质没有保留）。

4、获取主成分得分

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
> library(psych)
> pc <- principal(USJudgeRatings[,-1], nfactors = 1, scores = TRUE)
> head(pc$scores)
                 PC1
AARONSON,L.H.  -0.19
ALEXANDER,J.M.  0.75
ARMENTANO,A.J.  0.07
BERDON,R.I.     1.14
BRACKEN,J.J.   -2.16
BURNS,E.B.      0.77
> 
> #当主成分分析基于相关系数矩阵时，原始数据便不可用了，也不可能获取每个观测的主成分
> #得分，但是你可以得到用来计算主成分得分的系数
> 
> rc <- principal(Harman23.cor$cov, nfactors = 2, rotate="varimax")
> round(unclass(rc$weight), 2)
                 RC1   RC2
height          0.28 -0.05
arm.span        0.30 -0.08
forearm         0.30 -0.09
lower.leg       0.28 -0.06
weight         -0.06  0.33
bitro.diameter -0.08  0.32
chest.girth    -0.10  0.34
chest.width    -0.04  0.27
> 
> #则各记录的得分PC1 = 0.28*height +...-0.04 * chest.width
> #得分PC2 = -0.05*height +...+0.27 * chest.width

探索性因子分析

1、公因子数

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
> ##探索性因子分析
> #EFA的目标是通过发掘隐藏在数据下的一组较少的、更为基本的无法观测的变量，来解释一
> #组可观测变量的相关性。这些虚拟的、无法观测的变量称作因子。每个因子被认为可解释多个
> #观测变量间共有的方差，因此准确来说，它们应该称作公共因子。
> #112个人参与了六个测验，包括非语言的普通智力测验（general）、画图测验（picture）、积木图案测验（blocks）、迷津测验（maze）、
> #阅读测验（reading）和词汇测验（vocab）。我们如何用一组较少的、潜在的心理学因素来解
> #释参与者的测验得分呢？
> options(digits=2)
> covariances <- ability.cov$cov
> correlations <- cov2cor(covariances)
> 
> fa.parallel(correlations, n.obs = 112, fa="both", n.iter=100, main="Scree plots with parallel analysis")
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
The estimated weights for the factor scores are probably incorrect.  Try a different factor extraction method.
Parallel analysis suggests that the number of factors =  2  and the number of components =  1 
There were 18 warnings (use warnings() to see them)
> #图中同时展示了PCA和EFA的结果。PCA结果建议提取一个或者两个成分， EFA建议提取两个因子
> #对于EFA， Kaiser-Harris准则的特征值数大于0，而不是1。
>

2、提取公共因子

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
> #提取公共因子
> # X1 = a1F1 + a2F2 + ... +anFn
> #fa(r, nfactors=, n.objs=, rotate=, scores=, fm=)
> #r是相关系数矩阵或者原始数据矩阵；
> # nfactors设定提取的因子数（默认为1）；
> # n.obs是观测数（输入相关系数矩阵时需要填写）；
> # rotate设定旋转的方法（默认互变异数最小法）；
> # scores设定是否计算因子得分（默认不计算）；
> # fm设定因子化方法（默认极小残差法）
> 
> #提取公共因子的方法很多，包括最大似然法（ml）、主轴迭代法（pa）、加权
> #最小二乘法（wls）、广义加权最小二乘法（gls）和最小残差法（minres）
> #统计学家青睐使用最大似然法，因为它有良好的统计性质。不过有时候最大似然法不会收敛，此时使用主轴迭代法
> #效果会很好。
> 
> fa <- fa(correlations, nfactors = 2, rotate = "none", fm="pa")
> fa
Factor Analysis using method =  pa
Call: fa(r = correlations, nfactors = 2, rotate = "none", fm = "pa")
Standardized loadings (pattern matrix) based upon correlation matrix
         PA1   PA2   h2    u2 com
general 0.75  0.07 0.57 0.432 1.0
picture 0.52  0.32 0.38 0.623 1.7
blocks  0.75  0.52 0.83 0.166 1.8
maze    0.39  0.22 0.20 0.798 1.6
reading 0.81 -0.51 0.91 0.089 1.7
vocab   0.73 -0.39 0.69 0.313 1.5
                       PA1  PA2
SS loadings           2.75 0.83
Proportion Var        0.46 0.14
Cumulative Var        0.46 0.60
Proportion Explained  0.77 0.23
Cumulative Proportion 0.77 1.00
Mean item complexity =  1.5
Test of the hypothesis that 2 factors are sufficient.
The degrees of freedom for the null model are  15  and the objective function was  2.5
The degrees of freedom for the model are 4  and the objective function was  0.07 
The root mean square of the residuals (RMSR) is  0.03 
The df corrected root mean square of the residuals is  0.06 
Fit based upon off diagonal values = 0.99
Measures of factor score adequacy             
                                                   PA1  PA2
Correlation of (regression) scores with factors   0.96 0.92
Multiple R square of scores with factors          0.93 0.84
Minimum correlation of possible factor scores     0.86 0.68
> 
> #两个因子解释了六个心理学测验60%的方差。不过因子载荷阵的意义并不太好解
> #释，此时使用因子旋转将有助于因子的解释。

3、因子旋转

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
> #正交旋转
> fa.varimax <- fa(correlations, nfactors=2, rotate="varimax", fm="pa")
> fa.varimax
Factor Analysis using method =  pa
Call: fa(r = correlations, nfactors = 2, rotate = "varimax", fm = "pa")
Standardized loadings (pattern matrix) based upon correlation matrix
         PA1  PA2   h2    u2 com
general 0.49 0.57 0.57 0.432 2.0
picture 0.16 0.59 0.38 0.623 1.1
blocks  0.18 0.89 0.83 0.166 1.1
maze    0.13 0.43 0.20 0.798 1.2
reading 0.93 0.20 0.91 0.089 1.1
vocab   0.80 0.23 0.69 0.313 1.2
                       PA1  PA2
SS loadings           1.83 1.75
Proportion Var        0.30 0.29
Cumulative Var        0.30 0.60
Proportion Explained  0.51 0.49
Cumulative Proportion 0.51 1.00
Mean item complexity =  1.3
Test of the hypothesis that 2 factors are sufficient.
The degrees of freedom for the null model are  15  and the objective function was  2.5
The degrees of freedom for the model are 4  and the objective function was  0.07 
The root mean square of the residuals (RMSR) is  0.03 
The df corrected root mean square of the residuals is  0.06 
Fit based upon off diagonal values = 0.99
Measures of factor score adequacy             
                                                   PA1  PA2
Correlation of (regression) scores with factors   0.96 0.92
Multiple R square of scores with factors          0.91 0.85
Minimum correlation of possible factor scores     0.82 0.71
> 
> #阅读和词汇在第一因子上载荷较大，画图、积木图案和迷宫
> #在第二因子上载荷较大，非语言的普通智力测量在两个因子上载荷较为平均，这表明存在一个语
> #言智力因子和一个非语言智力因子。
> 
> #使用正交旋转将人为地强制两个因子不相关。如果想允许两个因子相关该怎么办呢？此时可
> #以使用斜交转轴法，比如promax
> 
> library(GPArotation)
> fa.promax <- fa(correlations, nfactors=2, rotate="promax", fm="pa")
Warning message:
In fac(r = r, nfactors = nfactors, n.obs = n.obs, rotate = rotate,  :
   A loading greater than abs(1) was detected.  Examine the loadings carefully.
> #因子斜交的结果图形
> factor.plot(fa.promax, labels=rownames(fa.promax$loadings))


1
2
3
> #若使simple = TRUE，那么将仅显示每个因子下最大的载荷，以及因子
> #间的相关系数。这类图形在有多个因子时十分实用。
> fa.diagram(fa.promax, simple = FALSE)

4、因子得分

1
2
3
4
5
6
7
8
9
10
11
12
13
14
> #因子得分
> fa.promax <- fa(correlations, nfactors=2, rotate="promax", fm="pa", score=TRUE)
Warning message:
In fac(r = r, nfactors = nfactors, n.obs = n.obs, rotate = rotate,  :
   A loading greater than abs(1) was detected.  Examine the loadings carefully.
> #获得用来计算因子得分的权重
> fa.promax$weights
          PA1   PA2
general 0.078 0.211
picture 0.020 0.090
blocks  0.037 0.702
maze    0.027 0.035
reading 0.743 0.030
vocab   0.177 0.036
文章目录
  1. 主成分和因子分析
    1. 主成分分析
      1. 1、主成分个数
      2. 2、提取主成分
      3. 3、主成分旋转
      4. 4、获取主成分得分
    2. 探索性因子分析
      1. 1、公因子数
      2. 2、提取公共因子
      3. 3、因子旋转
      4. 4、因子得分