Responses on a 7-pt scale (1=strongly disagree; 7=strongly agree) are recorded for above variables for 30 respondents. Conduct Factor Analysis/PCA on the data. Analyze and interpret the factors.
Exploratory data analysis:
dim(toothpaste)
## [1] 30 7
summary(toothpaste)
## RESPONDENT.NUMBER V1 V2 V3
## Min. : 1.00 Min. :1.000 Min. :2.0 Min. :1.0
## 1st Qu.: 8.25 1st Qu.:2.000 1st Qu.:3.0 1st Qu.:2.0
## Median :15.50 Median :4.000 Median :4.0 Median :4.0
## Mean :15.50 Mean :3.933 Mean :3.9 Mean :4.1
## 3rd Qu.:22.75 3rd Qu.:6.000 3rd Qu.:5.0 3rd Qu.:6.0
## Max. :30.00 Max. :7.000 Max. :7.0 Max. :7.0
## V4 V5 V6
## Min. :2.0 Min. :1.0 Min. :2.000
## 1st Qu.:3.0 1st Qu.:2.0 1st Qu.:3.000
## Median :4.0 Median :3.5 Median :4.000
## Mean :4.1 Mean :3.5 Mean :4.167
## 3rd Qu.:5.0 3rd Qu.:5.0 3rd Qu.:4.750
## Max. :7.0 Max. :7.0 Max. :7.000
nm_fa <-toothpaste[,]
PreventCav = nm_fa$V1
ShinyTeeth = nm_fa$V2
StrengthGum = nm_fa$V3
Fresh = nm_fa$V4
Decay = 8-nm_fa$V5
Attractive = nm_fa$V6
mydata<- data.frame(PreventCav,ShinyTeeth,StrengthGum,Fresh,Decay,Attractive)
mydata
## PreventCav ShinyTeeth StrengthGum Fresh Decay Attractive
## 1 7 3 6 4 6 4
## 2 1 3 2 4 3 4
## 3 6 2 7 4 7 3
## 4 4 5 4 6 6 5
## 5 1 2 2 3 2 2
## 6 6 3 6 4 6 4
## 7 5 3 6 3 4 3
## 8 6 4 7 4 7 4
## 9 3 4 2 3 2 3
## 10 2 6 2 6 1 6
## 11 6 4 7 3 6 3
## 12 2 3 1 4 3 4
## 13 7 2 6 4 7 3
## 14 4 6 4 5 5 6
## 15 1 3 2 2 2 4
## 16 6 4 6 3 5 4
## 17 5 3 6 3 5 4
## 18 7 3 7 4 7 4
## 19 2 4 3 3 2 3
## 20 3 5 3 6 4 6
## 21 1 3 2 3 3 3
## 22 5 4 5 4 6 4
## 23 2 2 1 5 4 4
## 24 4 6 4 6 4 7
## 25 6 5 4 2 7 4
## 26 3 5 4 6 4 7
## 27 4 4 7 2 6 5
## 28 3 7 2 6 4 3
## 29 4 6 3 7 6 7
## 30 2 3 2 4 1 2
attach(mydata)
## The following objects are masked _by_ .GlobalEnv:
##
## Attractive, Decay, Fresh, PreventCav, ShinyTeeth, StrengthGum
Load libraries required for analysis
#install.packages("nFactors")
library(nFactors)
## Warning: package 'nFactors' was built under R version 3.5.1
## Loading required package: MASS
## Loading required package: psych
## Warning: package 'psych' was built under R version 3.5.1
## Loading required package: boot
##
## Attaching package: 'boot'
## The following object is masked from 'package:psych':
##
## logit
## Loading required package: lattice
##
## Attaching package: 'lattice'
## The following object is masked from 'package:boot':
##
## melanoma
##
## Attaching package: 'nFactors'
## The following object is masked from 'package:lattice':
##
## parallel
Determine correlation matrix, Eigen values an Eigen vectors
cor(mydata)
## PreventCav ShinyTeeth StrengthGum Fresh Decay
## PreventCav 1.000000000 -0.05321785 0.87309020 -0.086162233 0.857636627
## ShinyTeeth -0.053217850 1.00000000 -0.15502002 0.572212066 -0.019745647
## StrengthGum 0.873090198 -0.15502002 1.00000000 -0.247787899 0.777848036
## Fresh -0.086162233 0.57221207 -0.24778790 1.000000000 0.006581882
## Decay 0.857636627 -0.01974565 0.77784804 0.006581882 1.000000000
## Attractive 0.004168129 0.64046495 -0.01806881 0.640464946 0.136402944
## Attractive
## PreventCav 0.004168129
## ShinyTeeth 0.640464946
## StrengthGum -0.018068814
## Fresh 0.640464946
## Decay 0.136402944
## Attractive 1.000000000
ev = eigen(cor(mydata))
ev
## eigen() decomposition
## $values
## [1] 2.73118833 2.21811927 0.44159791 0.34125765 0.18262823 0.08520861
##
## $vectors
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 0.5617362 -0.17002785 -0.01160665 0.24440063 0.17120132 0.7525734
## [2,] -0.1818493 -0.53396282 0.69061335 0.43660211 -0.06215148 -0.1018996
## [3,] 0.5664794 -0.08788446 0.15766817 -0.15898911 0.58984309 -0.5228057
## [4,] -0.2066896 -0.52974393 -0.68162756 0.34707290 0.26361452 -0.1486015
## [5,] 0.5256802 -0.23553722 -0.17876526 0.03986449 -0.74090048 -0.2927504
## [6,] -0.1068835 -0.58493089 0.03854977 -0.77609317 -0.02205218 0.2052773
EigenValue = ev$values
EigenValue
## [1] 2.73118833 2.21811927 0.44159791 0.34125765 0.18262823 0.08520861
Decide the number of factors using Scree plot
Factor = c(1,2,3,4,5,6)
Scree = data.frame(Factor, EigenValue)
plot(Scree, main="Scree Plot", col="Blue")
lines(Scree, col="Red")
Determine the Principal components without rotation
library(psych)
Unrotate = principal(mydata, nfactors = 3, rotate = "none")
print(Unrotate, digits =4)
## Principal Components Analysis
## Call: principal(r = mydata, nfactors = 3, rotate = "none")
## Standardized loadings (pattern matrix) based upon correlation matrix
## PC1 PC2 PC3 h2 u2 com
## PreventCav 0.9283 0.2532 -0.0077 0.9260 0.07400 1.148
## ShinyTeeth -0.3005 0.7952 0.4589 0.9334 0.06664 1.925
## StrengthGum 0.9362 0.1309 0.1048 0.9045 0.09545 1.065
## Fresh -0.3416 0.7890 -0.4530 0.9443 0.05568 2.012
## Decay 0.8688 0.3508 -0.1188 0.8919 0.10810 1.360
## Attractive -0.1766 0.8712 0.0256 0.7908 0.20923 1.084
##
## PC1 PC2 PC3
## SS loadings 2.7312 2.2181 0.4416
## Proportion Var 0.4552 0.3697 0.0736
## Cumulative Var 0.4552 0.8249 0.8985
## Proportion Explained 0.5066 0.4115 0.0819
## Cumulative Proportion 0.5066 0.9181 1.0000
##
## Mean item complexity = 1.4
## Test of the hypothesis that 3 components are sufficient.
##
## The root mean square of the residuals (RMSR) is 0.0515
## with the empirical chi square 2.3845 with prob < NA
##
## Fit based upon off diagonal values = 0.9882
UnrotatedProfile = plot(Unrotate,row.names(Unrotate$loadings))
Determine the Principal components with rotation
library(psych)
Rotate = principal(mydata, nfactors = 3, rotate = "varimax")
print(Rotate, digits =4)
## Principal Components Analysis
## Call: principal(r = mydata, nfactors = 3, rotate = "varimax")
## Standardized loadings (pattern matrix) based upon correlation matrix
## RC1 RC2 RC3 h2 u2 com
## PreventCav 0.9617 -0.0127 -0.0304 0.9260 0.07400 1.002
## ShinyTeeth -0.0648 0.9369 0.2268 0.9334 0.06664 1.127
## StrengthGum 0.9308 -0.0267 -0.1933 0.9045 0.09545 1.088
## Fresh -0.0852 0.3313 0.9096 0.9443 0.05568 1.280
## Decay 0.9359 -0.0050 0.1267 0.8919 0.10810 1.037
## Attractive 0.0856 0.6764 0.5709 0.7908 0.20923 1.981
##
## RC1 RC2 RC3
## SS loadings 2.6860 1.4459 1.2590
## Proportion Var 0.4477 0.2410 0.2098
## Cumulative Var 0.4477 0.6887 0.8985
## Proportion Explained 0.4983 0.2682 0.2335
## Cumulative Proportion 0.4983 0.7665 1.0000
##
## Mean item complexity = 1.3
## Test of the hypothesis that 3 components are sufficient.
##
## The root mean square of the residuals (RMSR) is 0.0515
## with the empirical chi square 2.3845 with prob < NA
##
## Fit based upon off diagonal values = 0.9882
RotatedProfile = plot(Rotate,nfactors = 3, row.names(Rotate$loadings), cex=1.0)
## Warning in plot.window(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy, type, ...): "nfactors" is not a graphical parameter
## Warning in title(...): "nfactors" is not a graphical parameter
## Warning in plot.window(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy, type, ...): "nfactors" is not a graphical parameter
## Warning in title(...): "nfactors" is not a graphical parameter
## Warning in axis(side = side, at = at, labels = labels, ...): "nfactors" is
## not a graphical parameter
## Warning in plot.xy(xy.coords(x, y), type = type, ...): "nfactors" is not a
## graphical parameter
## Warning in text.default(x, y, vnames, ...): "nfactors" is not a graphical
## parameter
## Warning in plot.window(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy, type, ...): "nfactors" is not a graphical parameter
## Warning in title(...): "nfactors" is not a graphical parameter
## Warning in axis(side = side, at = at, labels = labels, ...): "nfactors" is
## not a graphical parameter
## Warning in plot.xy(xy.coords(x, y), type = type, ...): "nfactors" is not a
## graphical parameter
## Warning in text.default(x, y, vnames, ...): "nfactors" is not a graphical
## parameter
## Warning in plot.window(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy, type, ...): "nfactors" is not a graphical parameter
## Warning in title(...): "nfactors" is not a graphical parameter
## Warning in axis(side = side, at = at, labels = labels, ...): "nfactors" is
## not a graphical parameter
## Warning in plot.xy(xy.coords(x, y), type = type, ...): "nfactors" is not a
## graphical parameter
## Warning in text.default(x, y, vnames, ...): "nfactors" is not a graphical
## parameter
## Warning in plot.window(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy, type, ...): "nfactors" is not a graphical parameter
## Warning in title(...): "nfactors" is not a graphical parameter
## Warning in plot.window(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy, type, ...): "nfactors" is not a graphical parameter
## Warning in title(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy.coords(x, y), type = type, ...): "nfactors" is not a
## graphical parameter
## Warning in text.default(x, y, vnames, ...): "nfactors" is not a graphical
## parameter
## Warning in plot.window(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy, type, ...): "nfactors" is not a graphical parameter
## Warning in title(...): "nfactors" is not a graphical parameter
## Warning in axis(side = side, at = at, labels = labels, ...): "nfactors" is
## not a graphical parameter
## Warning in plot.xy(xy.coords(x, y), type = type, ...): "nfactors" is not a
## graphical parameter
## Warning in text.default(x, y, vnames, ...): "nfactors" is not a graphical
## parameter
## Warning in plot.window(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy, type, ...): "nfactors" is not a graphical parameter
## Warning in title(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy.coords(x, y), type = type, ...): "nfactors" is not a
## graphical parameter
## Warning in text.default(x, y, vnames, ...): "nfactors" is not a graphical
## parameter
## Warning in plot.window(...): "nfactors" is not a graphical parameter
## Warning in plot.xy(xy, type, ...): "nfactors" is not a graphical parameter
## Warning in title(...): "nfactors" is not a graphical parameter
## Warning in axis(side = side, at = at, labels = labels, ...): "nfactors" is
## not a graphical parameter
## Warning in axis(side = side, at = at, labels = labels, ...): "nfactors" is
## not a graphical parameter
Rotate$scores
## RC1 RC2 RC3
## [1,] 1.15800103 -0.62899119 0.25025824
## [2,] -1.14924111 -0.68020893 0.25503587
## [3,] 1.30900482 -1.58258582 0.48502410
## [4,] 0.29880448 0.30842958 1.31118408
## [5,] -1.39907016 -1.50338963 -0.54150990
## [6,] 0.97695875 -0.62945273 0.24351854
## [7,] 0.37911545 -0.43622841 -0.90416167
## [8,] 1.31159396 0.07424164 -0.11737071
## [9,] -1.03439646 0.15436092 -1.13042533
## [10,] -1.30798968 1.47222011 0.67063385
## [11,] 1.08169029 0.19909900 -1.13619761
## [12,] -1.13327293 -0.74348087 0.36466473
## [13,] 1.32497300 -1.64585776 0.59465296
## [14,] 0.12608783 1.62057368 0.22790864
## [15,] -1.34762781 -0.02082536 -1.34483916
## [16,] 0.76698047 0.44344269 -0.99758948
## [17,] 0.60327444 -0.26580825 -0.62421548
## [18,] 1.50997248 -0.63408620 0.26948268
## [19,] -1.05036464 0.21763287 -1.24005418
## [20,] -0.38384498 0.62114023 1.32093890
## [21,] -1.19224741 -0.62418008 -0.64167741
## [22,] 0.61350612 0.01514161 -0.04044569
## [23,] -0.92329465 -1.81637630 1.60577275
## [24,] -0.01780322 1.63337334 1.00250832
## [25,] 0.78754607 1.18238560 -1.66657835
## [26,] -0.18150925 0.92412240 1.37588232
## [27,] 0.78838433 0.97195078 -1.57289256
## [28,] -0.69537648 1.25723949 0.19010294
## [29,] 0.19666207 1.13670530 2.08850541
## [30,] -1.41651682 -1.02058771 -0.29811680
factor.scores(mydata,f=Rotate$loadings,method = "Harman")
## $scores
## RC1 RC2 RC3
## [1,] 1.171514272 -0.59420838 0.20862734
## [2,] -1.168175555 -0.68600418 0.11829400
## [3,] 1.248929558 -1.40032983 0.42532777
## [4,] 0.262618106 0.46376080 1.26834763
## [5,] -1.364876363 -1.37893438 -0.52933987
## [6,] 0.957763245 -0.59957888 0.19347464
## [7,] 0.407087612 -0.46163357 -0.75294626
## [8,] 1.251796391 0.11864798 -0.05623194
## [9,] -0.934507478 0.15078442 -0.98059419
## [10,] -1.211544257 1.32815219 0.78258098
## [11,] 1.066009654 0.23449467 -0.91865281
## [12,] -1.105367042 -0.71612873 0.18522790
## [13,] 1.311738071 -1.43045438 0.49226167
## [14,] 0.106706752 1.48696033 0.24312546
## [15,] -1.350688534 -0.28462801 -1.43372901
## [16,] 0.782890581 0.32504470 -0.88373323
## [17,] 0.577186794 -0.39052672 -0.63211322
## [18,] 1.473594656 -0.58618245 0.22569346
## [19,] -0.997315990 0.18090898 -1.04752809
## [20,] -0.385389865 0.57640945 1.21926718
## [21,] -1.202824420 -0.59762661 -0.67527954
## [22,] 0.585022466 0.06975651 -0.03666955
## [23,] -0.930494377 -1.64075230 1.26243577
## [24,] -0.009782253 1.42605188 0.96785166
## [25,] 0.759546501 1.09627082 -1.65083668
## [26,] -0.215486040 0.71048046 1.21947167
## [27,] 0.660742670 0.60785913 -1.58657463
## [28,] -0.609310791 1.66558839 0.58154588
## [29,] 0.157238529 1.14866507 1.89891537
## [30,] -1.294622893 -0.82284737 -0.10821938
##
## $weights
## RC1 RC2 RC3
## PreventCav 0.42355276 0.01064176 0.03002543
## ShinyTeeth -0.01105202 0.97538461 -0.36638362
## StrengthGum 0.31048088 0.07301146 -0.10651125
## Fresh 0.02154517 -0.25676057 1.01849115
## Decay 0.28827072 -0.05239284 0.13131501
## Attractive 0.02638813 0.13718648 0.07234762
##
## $r.scores
## RC1 RC2 RC3
## RC1 1.000000000 0.001997905 -0.004965744
## RC2 0.001997905 1.000000000 0.109062480
## RC3 -0.004965744 0.109062480 1.000000000
##
## $missing
## [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
##
## $R2
## [1] 0.9672721 0.9197127 0.9208932