実データに対する因子分析①
滋賀大学 奥村太一 先生 の上越教育大学での教授資料を学んだ記録です。この場を借りて御礼申し上げます。（先生の滋賀大学へのご転勤にともない、現在、この内容は公開されていません）
library(psych)
library(MASS)
#3D表示
options(rgl.printRglwidget = TRUE)
#WebGLに変換してRStudio, VSCode, ウェブブラウザなど出力
library(rgl)
library(ggplot2)##
## Attaching package: 'ggplot2'## The following objects are masked from 'package:psych':
##
## %+%, alphalibrary(ggExtra)
library("sem")
par(family = "HiraKakuProN-W3") #日本語フォントの指定
library("semPlot")1 データ
data01 <- read.csv("factor.csv",header=T,
row.names=1)
#data01 <- factor
head(data01)| 国語 | 英語 | 社会 | 数学 | 理科 | 技術 |
|---|---|---|---|---|---|
| 4 | 4 | 4 | 4 | 3 | 5 |
| 4 | 4 | 4 | 4 | 5 | 5 |
| 3 | 3 | 4 | 4 | 4 | 4 |
| 3 | 3 | 3 | 3 | 3 | 2 |
| 3 | 3 | 2 | 2 | 3 | 2 |
| 3 | 3 | 2 | 3 | 2 | 2 |
str(data01)## 'data.frame': 200 obs. of 6 variables:
## $ 国語: int 4 4 3 3 3 3 4 3 4 3 ...
## $ 英語: int 4 4 3 3 3 3 3 4 3 3 ...
## $ 社会: int 4 4 4 3 2 2 4 3 4 3 ...
## $ 数学: int 4 4 4 3 2 3 3 2 4 5 ...
## $ 理科: int 3 5 4 3 3 2 3 3 4 3 ...
## $ 技術: int 5 5 4 2 2 2 3 3 4 4 ...cor01 <- cor(data01); cor01## 国語 英語 社会 数学 理科 技術
## 国語 1.0000000 0.4158982 0.4836725 0.2505441 0.2907936 0.2225323
## 英語 0.4158982 1.0000000 0.4980821 0.2649821 0.2946882 0.3051505
## 社会 0.4836725 0.4980821 1.0000000 0.2945884 0.4020571 0.3465945
## 数学 0.2505441 0.2649821 0.2945884 1.0000000 0.5766291 0.5583059
## 理科 0.2907936 0.2946882 0.4020571 0.5766291 1.0000000 0.5253087
## 技術 0.2225323 0.3051505 0.3465945 0.5583059 0.5253087 1.0000000eigen(cor01)## eigen() decomposition
## $values
## [1] 2.9199728 1.1337679 0.5975924 0.5004401 0.4595162 0.3887106
##
## $vectors
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] -0.3590916 -0.4837731 0.6680046 0.37069260 0.1663058 -0.1603623
## [2,] -0.3799172 -0.4216254 -0.6656819 0.36200826 -0.2978789 -0.1223958
## [3,] -0.4219678 -0.3725137 -0.0501626 -0.70076040 0.1594348 0.4051862
## [4,] -0.4192907 0.4426995 0.1286752 0.33581647 -0.2436834 0.6629486
## [5,] -0.4421990 0.3233837 0.1840456 -0.34618027 -0.5034363 -0.5410373
## [6,] -0.4210875 0.3858358 -0.2401897 0.08864763 0.7384859 -0.2508091eigen01 <- eigen(cor01)$values; eigen01## [1] 2.9199728 1.1337679 0.5975924 0.5004401 0.4595162 0.3887106eigen(cor01)$vectors[,1]#第1固有ベクトル## [1] -0.3590916 -0.3799172 -0.4219678 -0.4192907 -0.4421990 -0.4210875cor01 %*% eigen(cor01)$vectors#固有ベクトルが固有値倍に変換される 行列の掛け算## [,1] [,2] [,3] [,4] [,5] [,6]
## 国語 -1.048538 -0.5484864 0.39919446 0.18550943 0.07642021 -0.06233454
## 英語 -1.109348 -0.4780253 -0.39780648 0.18116344 -0.13688017 -0.04757655
## 社会 -1.232135 -0.4223441 -0.02997679 -0.35068859 0.07326286 0.15750017
## 数学 -1.224318 0.5019185 0.07689531 0.16805602 -0.11197649 0.25769515
## 理科 -1.291209 0.3666421 0.10998423 -0.17324248 -0.23133716 -0.21030694
## 技術 -1.229564 0.4374482 -0.14353556 0.04436283 0.33934622 -0.09749215cor01 %*% eigen(cor01)$vectors[,1]#固有ベクトルを変換すると固有値倍になる## [,1]
## 国語 -1.048538
## 英語 -1.109348
## 社会 -1.232135
## 数学 -1.224318
## 理科 -1.291209
## 技術 -1.2295642.9199*eigen(cor01)$vectors#固有ベクトルの固有値倍## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] -1.048512 -1.4125692 1.9505065 1.0823853 0.4855963 -0.4682419
## [2,] -1.109320 -1.2311040 -1.9437247 1.0570279 -0.8697765 -0.3573835
## [3,] -1.232104 -1.0877027 -0.1464698 -2.0461503 0.4655335 1.1831031
## [4,] -1.224287 1.2926383 0.3757186 0.9805505 -0.7115313 1.9357435
## [5,] -1.291177 0.9442481 0.5373946 -1.0108118 -1.4699838 -1.5797748
## [6,] -1.229533 1.1266019 -0.7013300 0.2588422 2.1563048 -0.7323374plot(eigen01, type="b", main="Scree Plot",xlab="Number", ylab="Eigenvalue")
fit00 <- factanal(x=data01, factors=2, rotation="none")#因子分析の初期解
print(fit00, cutoff=0)##
## Call:
## factanal(x = data01, factors = 2, rotation = "none")
##
## Uniquenesses:
## 国語 英語 社会 数学 理科 技術
## 0.609 0.584 0.393 0.372 0.449 0.494
##
## Loadings:
## Factor1 Factor2
## 国語 0.510 0.362
## 英語 0.545 0.344
## 社会 0.656 0.420
## 数学 0.704 -0.363
## 理科 0.715 -0.199
## 技術 0.670 -0.238
##
## Factor1 Factor2
## SS loadings 2.445 0.653
## Proportion Var 0.408 0.109
## Cumulative Var 0.408 0.516
##
## Test of the hypothesis that 2 factors are sufficient.
## The chi square statistic is 2.81 on 4 degrees of freedom.
## The p-value is 0.59print(fit00, cutoff=0.3)##
## Call:
## factanal(x = data01, factors = 2, rotation = "none")
##
## Uniquenesses:
## 国語 英語 社会 数学 理科 技術
## 0.609 0.584 0.393 0.372 0.449 0.494
##
## Loadings:
## Factor1 Factor2
## 国語 0.510 0.362
## 英語 0.545 0.344
## 社会 0.656 0.420
## 数学 0.704 -0.363
## 理科 0.715
## 技術 0.670
##
## Factor1 Factor2
## SS loadings 2.445 0.653
## Proportion Var 0.408 0.109
## Cumulative Var 0.408 0.516
##
## Test of the hypothesis that 2 factors are sufficient.
## The chi square statistic is 2.81 on 4 degrees of freedom.
## The p-value is 0.59psych::alpha(data01)##
## Reliability analysis
## Call: psych::alpha(x = data01)
##
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
## 0.78 0.79 0.78 0.38 3.7 0.024 3 0.5 0.35
##
## 95% confidence boundaries
## lower alpha upper
## Feldt 0.73 0.78 0.83
## Duhachek 0.73 0.78 0.83
##
## Reliability if an item is dropped:
## raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
## 国語 0.77 0.77 0.76 0.41 3.4 0.025 0.015 0.37
## 英語 0.76 0.77 0.75 0.40 3.3 0.027 0.018 0.37
## 社会 0.74 0.75 0.73 0.37 2.9 0.029 0.019 0.30
## 数学 0.74 0.75 0.73 0.38 3.0 0.028 0.010 0.37
## 理科 0.73 0.74 0.73 0.36 2.9 0.030 0.014 0.33
## 技術 0.74 0.75 0.74 0.38 3.0 0.029 0.013 0.35
##
## Item statistics
## n raw.r std.r r.cor r.drop mean sd
## 国語 200 0.66 0.64 0.53 0.45 3 0.83
## 英語 200 0.63 0.67 0.56 0.50 3 0.54
## 社会 200 0.73 0.72 0.65 0.57 3 0.75
## 数学 200 0.71 0.70 0.64 0.55 3 0.77
## 理科 200 0.73 0.74 0.68 0.60 3 0.66
## 技術 200 0.71 0.71 0.64 0.55 3 0.73
##
## Non missing response frequency for each item
## 1 2 3 4 5 miss
## 国語 0.03 0.26 0.44 0.27 0.01 0
## 英語 0.00 0.14 0.72 0.14 0.00 0
## 社会 0.02 0.24 0.50 0.24 0.00 0
## 数学 0.02 0.22 0.54 0.20 0.03 0
## 理科 0.00 0.20 0.60 0.17 0.01 0
## 技術 0.00 0.26 0.51 0.22 0.01 0#biplot(fit00)
#fa.diagram( fit00 )
fit01 <- factanal(x=data01, factors=2, rotation="varimax")#バリマックス回転
print(fit01, cutoff=0)##
## Call:
## factanal(x = data01, factors = 2, rotation = "varimax")
##
## Uniquenesses:
## 国語 英語 社会 数学 理科 技術
## 0.609 0.584 0.393 0.372 0.449 0.494
##
## Loadings:
## Factor1 Factor2
## 国語 0.163 0.604
## 英語 0.202 0.612
## 社会 0.238 0.742
## 数学 0.774 0.169
## 理科 0.678 0.302
## 技術 0.668 0.244
##
## Factor1 Factor2
## SS loadings 1.629 1.469
## Proportion Var 0.272 0.245
## Cumulative Var 0.272 0.516
##
## Test of the hypothesis that 2 factors are sufficient.
## The chi square statistic is 2.81 on 4 degrees of freedom.
## The p-value is 0.59print(fit01, cutoff=0.3)##
## Call:
## factanal(x = data01, factors = 2, rotation = "varimax")
##
## Uniquenesses:
## 国語 英語 社会 数学 理科 技術
## 0.609 0.584 0.393 0.372 0.449 0.494
##
## Loadings:
## Factor1 Factor2
## 国語 0.604
## 英語 0.612
## 社会 0.742
## 数学 0.774
## 理科 0.678 0.302
## 技術 0.668
##
## Factor1 Factor2
## SS loadings 1.629 1.469
## Proportion Var 0.272 0.245
## Cumulative Var 0.272 0.516
##
## Test of the hypothesis that 2 factors are sufficient.
## The chi square statistic is 2.81 on 4 degrees of freedom.
## The p-value is 0.59fit02 <- factanal(x=data01, factors=2, rotation="promax")#プロマックス回転
print(fit02, cutoff=0)##
## Call:
## factanal(x = data01, factors = 2, rotation = "promax")
##
## Uniquenesses:
## 国語 英語 社会 数学 理科 技術
## 0.609 0.584 0.393 0.372 0.449 0.494
##
## Loadings:
## Factor1 Factor2
## 国語 -0.035 0.645
## 英語 0.007 0.641
## 社会 0.001 0.778
## 数学 0.836 -0.081
## 理科 0.675 0.109
## 技術 0.685 0.044
##
## Factor1 Factor2
## SS loadings 1.624 1.452
## Proportion Var 0.271 0.242
## Cumulative Var 0.271 0.513
##
## Factor Correlations:
## Factor1 Factor2
## Factor1 1.000 0.571
## Factor2 0.571 1.000
##
## Test of the hypothesis that 2 factors are sufficient.
## The chi square statistic is 2.81 on 4 degrees of freedom.
## The p-value is 0.59print(fit02, cutoff=0.3)##
## Call:
## factanal(x = data01, factors = 2, rotation = "promax")
##
## Uniquenesses:
## 国語 英語 社会 数学 理科 技術
## 0.609 0.584 0.393 0.372 0.449 0.494
##
## Loadings:
## Factor1 Factor2
## 国語 0.645
## 英語 0.641
## 社会 0.778
## 数学 0.836
## 理科 0.675
## 技術 0.685
##
## Factor1 Factor2
## SS loadings 1.624 1.452
## Proportion Var 0.271 0.242
## Cumulative Var 0.271 0.513
##
## Factor Correlations:
## Factor1 Factor2
## Factor1 1.000 0.571
## Factor2 0.571 1.000
##
## Test of the hypothesis that 2 factors are sufficient.
## The chi square statistic is 2.81 on 4 degrees of freedom.
## The p-value is 0.59par(family = "HiraKakuProN-W3") #日本語フォントの指定
library(semPlot)#パス図を描く
semPaths(fit00,whatLabels="stand",layout="circle",style="lisrel")
semPaths(fit02,whatLabels="stand",layout="circle",style="lisrel")
semPaths(fit01,whatLabels="stand",layout="circle",style="lisrel")
2 パス解析
data02 <- read.csv("regpath.csv",header=T,fileEncoding="CP932")
head(data02)| 生徒 | 学力 | 意欲 | 適応 | 自尊 |
|---|---|---|---|---|
| 1 | 38 | 6 | 14 | 3 |
| 2 | 75 | 16 | 17 | 14 |
| 3 | 44 | 9 | 14 | 11 |
| 4 | 59 | 8 | 13 | 7 |
| 5 | 50 | 13 | 17 | 3 |
| 6 | 75 | 14 | 11 | 16 |
data02 <- data02[,2:5]
#data02
par(family = "HiraKakuProN-W3") #日本語フォントの指定
mcov02 <- cov(data02);mcov02#標準化されていない（因子分析とは異なる）## 学力 意欲 適応 自尊
## 学力 457.00000 51.950000 39.0333333 55.1500000
## 意欲 51.95000 17.026667 6.2933333 6.2950000
## 適応 39.03333 6.293333 14.2266667 0.2033333
## 自尊 55.15000 6.295000 0.2033333 28.6566667cor(data02)## 学力 意欲 適応 自尊
## 学力 1.0000000 0.5889286 0.48408983 0.48191959
## 意欲 0.5889286 1.0000000 0.40435642 0.28498243
## 適応 0.4840898 0.4043564 1.00000000 0.01007034
## 自尊 0.4819196 0.2849824 0.01007034 1.00000000library("sem")
library("semPlot")
fit02 <- factanal(x=data02,factors=1,rotation="none")
#因子分析の初期解
#print(fit02, cutoff=0)
print(fit02, cutoff=0.3)##
## Call:
## factanal(x = data02, factors = 1, rotation = "none")
##
## Uniquenesses:
## 学力 意欲 適応 自尊
## 0.005 0.651 0.765 0.767
##
## Loadings:
## Factor1
## 学力 0.997
## 意欲 0.591
## 適応 0.485
## 自尊 0.482
##
## Factor1
## SS loadings 1.812
## Proportion Var 0.453
##
## Test of the hypothesis that 1 factor is sufficient.
## The chi square statistic is 2.55 on 2 degrees of freedom.
## The p-value is 0.279psych::alpha(data02)## Number of categories should be increased in order to count frequencies.##
## Reliability analysis
## Call: psych::alpha(x = data02)
##
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
## 0.51 0.71 0.71 0.38 2.4 0.062 23 7.2 0.44
##
## 95% confidence boundaries
## lower alpha upper
## Feldt 0.09 0.51 0.76
## Duhachek 0.39 0.51 0.63
##
## Reliability if an item is dropped:
## raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
## 学力 0.45 0.48 0.43 0.23 0.91 0.192 0.0409 0.28
## 意欲 0.41 0.59 0.61 0.33 1.45 0.071 0.0746 0.48
## 適応 0.47 0.71 0.66 0.45 2.47 0.075 0.0238 0.48
## 自尊 0.43 0.74 0.67 0.49 2.91 0.073 0.0086 0.48
##
## Item statistics
## n raw.r std.r r.cor r.drop mean sd
## 学力 25 0.98 0.88 0.86 0.74 56 21.4
## 意欲 25 0.68 0.78 0.67 0.60 12 4.1
## 適応 25 0.55 0.65 0.50 0.45 11 3.8
## 自尊 25 0.58 0.61 0.44 0.44 11 5.4par(family = "HiraKakuProN-W3") #日本語フォントの指定
semPaths(fit02,whatLabels="stand",layout="circle",style="lisrel")