Note: the estimation process can be time consuming depending on the computing power. You can same some time by reducing the length of the chains.
library(LAWBL)
dat <- sim18cfa0$dat
J <- ncol(dat) # no. of items
K <- 3 # no. of factors
qlam <- sim18cfa0$qlam
qlam
R> [,1] [,2] [,3]
R> [1,] 0.7 0.0 0.0
R> [2,] 0.7 0.0 0.0
R> [3,] 0.7 0.0 0.0
R> [4,] 0.7 0.0 0.0
R> [5,] 0.7 0.4 0.0
R> [6,] 0.7 0.4 0.0
R> [7,] 0.0 0.7 0.0
R> [8,] 0.0 0.7 0.0
R> [9,] 0.0 0.7 0.0
R> [10,] 0.0 0.7 0.0
R> [11,] 0.0 0.7 0.4
R> [12,] 0.0 0.7 0.4
R> [13,] 0.0 0.0 0.7
R> [14,] 0.0 0.0 0.7
R> [15,] 0.0 0.0 0.7
R> [16,] 0.0 0.0 0.7
R> [17,] 0.4 0.0 0.7
R> [18,] 0.4 0.0 0.7
Q<-matrix(-1,J,K); # -1 for unspecified items
Q[1:2,1]<-Q[7:8,2]<-Q[13:14,3]<-1 # 1 for specified items
Q
R> [,1] [,2] [,3]
R> [1,] 1 -1 -1
R> [2,] 1 -1 -1
R> [3,] -1 -1 -1
R> [4,] -1 -1 -1
R> [5,] -1 -1 -1
R> [6,] -1 -1 -1
R> [7,] -1 1 -1
R> [8,] -1 1 -1
R> [9,] -1 -1 -1
R> [10,] -1 -1 -1
R> [11,] -1 -1 -1
R> [12,] -1 -1 -1
R> [13,] -1 -1 1
R> [14,] -1 -1 1
R> [15,] -1 -1 -1
R> [16,] -1 -1 -1
R> [17,] -1 -1 -1
R> [18,] -1 -1 -1
m0 <- pcfa(dat = dat, Q = Q,LD = FALSE, burn = 2000, iter = 2000,verbose = TRUE)
R>
R> Tot. Iter = 1000
R> [,1] [,2] [,3]
R> Feigen 4.744 4.269 2.380
R> NLA_lg3 8.000 8.000 6.000
R> Shrink 3.041 3.041 3.041
R>
R> Tot. Iter = 2000
R> [,1] [,2] [,3]
R> Feigen 4.079 4.030 3.224
R> NLA_lg3 8.000 8.000 6.000
R> Shrink 2.734 2.734 2.734
R>
R> Tot. Iter = 3000
R> [,1] [,2] [,3]
R> Feigen 3.258 3.460 3.418
R> NLA_lg3 8.000 8.000 8.000
R> Shrink 3.809 3.809 3.809
R> Adj PSR 2.313 1.231 1.517
R>
R> Tot. Iter = 4000
R> [,1] [,2] [,3]
R> Feigen 3.196 3.430 2.891
R> NLA_lg3 8.000 8.000 8.000
R> Shrink 3.527 3.527 3.527
R> Adj PSR 1.21 1.34 1.121
R>
R> #Sign change: 0 0 0 0 0 0
R> user system elapsed
R> 15.35 0.02 15.35
# summarize basic information
summary(m0)
R> $N
R> [1] 1000
R>
R> $J
R> [1] 18
R>
R> $K
R> [1] 3
R>
R> $`Miss%`
R> [1] 0
R>
R> $`LD enabled`
R> [1] FALSE
R>
R> $`Burn in`
R> [1] 2000
R>
R> $Iteration
R> [1] 2000
R>
R> $`No. of sig lambda`
R> [1] 24
R>
R> $`Adj. PSR`
R> Point est. Upper C.I.
R> F1 1.210460 1.752360
R> F2 1.340018 2.156112
R> F3 1.120779 1.193162
#summarize significant loadings in pattern/Q-matrix format
summary(m0, what = 'qlambda')
R> [,1] [,2] [,3]
R> I1 0.6943921 0.0000000 0.0000000
R> I2 0.7018397 0.0000000 0.0000000
R> I3 0.6940404 0.0000000 0.0000000
R> I4 0.6821522 0.0000000 0.0000000
R> I5 0.6970800 0.4219813 0.0000000
R> I6 0.6922918 0.4145146 0.0000000
R> I7 0.0000000 0.7633059 0.0000000
R> I8 0.0000000 0.6836233 0.0000000
R> I9 0.0000000 0.7454941 0.0000000
R> I10 0.0000000 0.6788261 0.0000000
R> I11 0.0000000 0.7309471 0.3569716
R> I12 0.0000000 0.7170643 0.3551296
R> I13 0.0000000 0.0000000 0.6787302
R> I14 0.0000000 0.0000000 0.6770629
R> I15 0.0000000 0.0000000 0.7500404
R> I16 0.0000000 0.0000000 0.6942952
R> I17 0.4276174 0.0000000 0.7058218
R> I18 0.4191974 0.0000000 0.7049021
#factorial eigenvalue
summary(m0,what='eigen')
R> est sd lower upper sig
R> F1 3.308617 0.4444696 2.363488 4.291951 1
R> F2 3.530416 0.4479498 2.519580 4.262045 1
R> F3 3.310640 0.4478500 2.599306 4.240286 1
#plotting factorial eigenvalue
plot_eigen(m0) # trace
Q<-matrix(-1,J,K);
tmp<-summary(m0, what="qlambda")
cind<-apply(tmp,1,which.max)
Q[cbind(c(1:J),cind)]<-1
#alternatively
#Q[1:6,1]<-Q[7:12,2]<-Q[13:18,3]<-1 # 1 for specified items
m1 <- pcfa(dat = dat, Q = Q, burn = 2000, iter = 2000,verbose = TRUE)
summary(m1)
summary(m1, what = 'qlambda')
summary(m1, what = 'offpsx') #summarize significant LD terms
summary(m1,what='eigen')
#plotting factorial eigenvalue
# par(mar = rep(2, 4))
plot_eigen(m1) # trace
plot_eigen(m1, what='density') #density
plot_eigen(m1, what='APSR') #adj, PSRF
dat <- sim18cfa1$dat
J <- ncol(dat) # no. of items
K <- 3 # no. of factors
sim18cfa1$qlam
sim18cfa1$LD # effect size = .3
Q<-matrix(-1,J,K); # -1 for unspecified items
Q[1:2,1]<-Q[7:8,2]<-Q[13:14,3]<-1 # 1 for specified items
Q
m0 <- pcfa(dat = dat, Q = Q,LD = FALSE, burn = 4000, iter = 4000,verbose = TRUE)
summary(m0)
summary(m0, what = 'qlambda')
summary(m0,what='eigen')
plot_eigen(m0) # trace
plot_eigen(m0, what='APSR')
Q<-matrix(-1,J,K);
tmp<-summary(m0, what="qlambda")
cind<-apply(tmp,1,which.max)
Q[cbind(c(1:J),cind)]<-1
Q
m1 <- pcfa(dat = dat, Q = Q,burn = 4000, iter = 4000,verbose = TRUE)
summary(m1)
summary(m1, what = 'qlambda')
summary(m1,what='eigen')
summary(m1, what = 'offpsx')