library(ltm)
library(psych)
library(MPsychoR)
library(dplyr)
library(tidyr)
library(knitr)
#Bu kod baya işlevsel oldu
setwd(dirname(rstudioapi::getActiveDocumentContext()$path))

1 Giriş

Bu günlükte:

  • Tek boyutluluk
  • Model-veri uyumu
  • Yen Q3 istatistiği
  • Model seçimi ve karşılaştırma
  • RSM
  • PCM
  • GPCM
  • GRM kestirimleri yapılmıştır.

Veri seti LISS 2024 veri setinde yer alan PANAS ölçeğine ait Positive Affect boyutudur.

data <- haven::read_sav("data.sav")
data <- data %>%
  select(cp24p146:cp24p165) %>% na.omit()
data <- data[,c(1,3,5,8,10,12,14,16,17,19)]
#cp24p020 - cp24p069: BIG-V (IPIP, Goldberg)
#For more information: IPIP Home (ori.org)
kable(psych::describe(data, fast = T),format = "html",digits = 2)
vars n mean sd median min max range skew kurtosis se
cp24p146 1 5421 5.17 1.28 5 1 7 6 -0.88 0.94 0.02
cp24p148 2 5421 3.07 1.53 3 1 7 6 0.26 -0.78 0.02
cp24p150 3 5421 4.66 1.31 5 1 7 6 -0.53 0.22 0.02
cp24p153 4 5421 1.76 1.19 1 1 7 6 1.87 3.36 0.02
cp24p155 5 5421 4.61 1.51 5 1 7 6 -0.54 -0.15 0.02
cp24p157 6 5421 4.52 1.49 5 1 7 6 -0.52 -0.14 0.02
cp24p159 7 5421 3.97 1.54 4 1 7 6 -0.26 -0.52 0.02
cp24p161 8 5421 4.41 1.52 5 1 7 6 -0.51 -0.20 0.02
cp24p162 9 5421 4.74 1.43 5 1 7 6 -0.65 0.19 0.02
cp24p164 10 5421 4.42 1.52 5 1 7 6 -0.39 -0.41 0.02

2 Tek Boyutluluk

library(Gifi)
## Warning: package 'Gifi' was built under R version 4.4.3
## 
## Attaching package: 'Gifi'
## The following object is masked from 'package:MASS':
## 
##     mammals
pc <- Gifi::princals(as.data.frame(data))
plot(pc)

148 ve 153 kodlu maddeler tek boyutluluğu bozuyor gibi duruyor.

library(EGAnet)
## Warning: package 'EGAnet' was built under R version 4.4.3
## 
## EGAnet (version 2.3.0) 
## 
## For help getting started, see <https://r-ega.net> 
## 
## For bugs and errors, submit an issue to <https://github.com/hfgolino/EGAnet/issues>
## 
## Attaching package: 'EGAnet'
## The following object is masked from 'package:ltm':
## 
##     information
UVA(data)
## Variable pairs with wTO > 0.30 (large-to-very large redundancy)
## 
##    node_i   node_j   wto
##  cp24p157 cp24p162 0.370
##  cp24p148 cp24p153 0.337
## 
## ----
## 
## Variable pairs with wTO > 0.25 (moderate-to-large redundancy)
## 
## ----
## 
## Variable pairs with wTO > 0.20 (small-to-moderate redundancy)
## 
##    node_i   node_j   wto
##  cp24p150 cp24p155 0.230
##  cp24p161 cp24p162 0.207

148-153 ve 157-162 arasında oldukça yüksek reduncancy kestirildiği için bu maddeleri çiftlerinden 148 ve 153’ü çıkaralım.

data <- data[,-c(2,4)]
pc2 <- princals(as.data.frame(data))
plot(pc2)

psych::fa.parallel(data)

## Parallel analysis suggests that the number of factors =  4  and the number of components =  1

PA 1 bileşen önerdi.

EGA(data)

## Model: GLASSO (EBIC with gamma = 0.5)
## Correlations: auto
## Lambda: 0.0699399504202639 (n = 100, ratio = 0.1)
## 
## Number of nodes: 8
## Number of edges: 25
## Edge density: 0.893
## 
## Non-zero edge weights: 
##      M    SD   Min   Max
##  0.134 0.104 0.015 0.458
## 
## ----
## 
## Algorithm:  Louvain
## 
## Number of communities:  1
## 
## cp24p146 cp24p150 cp24p155 cp24p157 cp24p159 cp24p161 cp24p162 cp24p164 
##        1        1        1        1        1        1        1        1 
## 
## ----
## 
## Unidimensional Method: Louvain
## Unidimensional: Yes
## 
## ----
## 
## TEFI: 0

EGA da 1 boyut önerdiği için kestirimlere devam edebiliriz.

fa1 <- psych::fa(data,cor = "poly")
kable(fa1$Vaccounted,digits = 2)
MR1
SS loadings 4.05
Proportion Var 0.51 
kable(fa1$loadings,format = "html",digits = 2)
MR1
cp24p146 0.57
cp24p150 0.71
cp24p155 0.69
cp24p157 0.68
cp24p159 0.71
cp24p161 0.78
cp24p162 0.76
cp24p164 0.76

3 Yerel Bağımsızlık

library(mirt)
## Warning: package 'mirt' was built under R version 4.4.2
## Loading required package: stats4
## Loading required package: lattice
## 
## Attaching package: 'mirt'
## The following object is masked from 'package:ltm':
## 
##     Science
est1 <- mirt(data, 1, verbose = FALSE)
Q3 <- residuals(est1, type = "Q3")
## Q3 summary statistics:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  -0.346  -0.175  -0.128  -0.118  -0.050   0.251 
## 
##          cp24p146 cp24p150 cp24p155 cp24p157 cp24p159 cp24p161 cp24p162
## cp24p146    1.000   -0.005   -0.042   -0.118   -0.036   -0.167   -0.137
## cp24p150   -0.005    1.000    0.144   -0.237   -0.119   -0.171   -0.331
## cp24p155   -0.042    0.144    1.000   -0.184   -0.028   -0.116   -0.346
## cp24p157   -0.118   -0.237   -0.184    1.000   -0.149   -0.225    0.251
## cp24p159   -0.036   -0.119   -0.028   -0.149    1.000   -0.082   -0.240
## cp24p161   -0.167   -0.171   -0.116   -0.225   -0.082    1.000   -0.063
## cp24p162   -0.137   -0.331   -0.346    0.251   -0.240   -0.063    1.000
## cp24p164   -0.051   -0.094   -0.158   -0.151   -0.047   -0.240   -0.149
##          cp24p164
## cp24p146   -0.051
## cp24p150   -0.094
## cp24p155   -0.158
## cp24p157   -0.151
## cp24p159   -0.047
## cp24p161   -0.240
## cp24p162   -0.149
## cp24p164    1.000
abs(Q3) >= 0.2
##          cp24p146 cp24p150 cp24p155 cp24p157 cp24p159 cp24p161 cp24p162
## cp24p146     TRUE    FALSE    FALSE    FALSE    FALSE    FALSE    FALSE
## cp24p150    FALSE     TRUE    FALSE     TRUE    FALSE    FALSE     TRUE
## cp24p155    FALSE    FALSE     TRUE    FALSE    FALSE    FALSE     TRUE
## cp24p157    FALSE     TRUE    FALSE     TRUE    FALSE     TRUE     TRUE
## cp24p159    FALSE    FALSE    FALSE    FALSE     TRUE    FALSE     TRUE
## cp24p161    FALSE    FALSE    FALSE     TRUE    FALSE     TRUE    FALSE
## cp24p162    FALSE     TRUE     TRUE     TRUE     TRUE    FALSE     TRUE
## cp24p164    FALSE    FALSE    FALSE    FALSE    FALSE     TRUE    FALSE
##          cp24p164
## cp24p146    FALSE
## cp24p150    FALSE
## cp24p155    FALSE
## cp24p157    FALSE
## cp24p159    FALSE
## cp24p161     TRUE
## cp24p162    FALSE
## cp24p164     TRUE

Yerel bağımsızlık varsayımının ihlal edildiği açıkça görülmektedir. (\(Q\leq 0.20\)) Testlet modelin nasıl dahil edileceğini bilemedim o yüzden bu hali ile kestirimlere devam edeceğim.

4 Çok Kategorili

eRm paketi için veti seti 0’dan başlamalıydı.

library(eRm)
data <- data - 1
fitRSM <- RSM(data)
fitRSM
## 
## Results of RSM estimation: 
## 
## Call:  RSM(X = data) 
## 
## Conditional log-likelihood: -47621.98 
## Number of iterations: 45 
## Number of parameters: 12 
## 
## Item (Category) Difficulty Parameters (eta):
##             cp24p150    cp24p155   cp24p157   cp24p159   cp24p161    cp24p162
## Estimate -0.07499758 -0.02853705 0.04783006 0.49860646 0.14246751 -0.14563089
## Std.Err   0.01210244  0.01200758 0.01187818 0.01188461 0.01176669  0.01226858
##          cp24p164      Cat 2      Cat 3     Cat 4     Cat 5     Cat 6
## Estimate 0.136547 0.04367397 0.04956308 1.1872902 3.1716794 6.7035835
## Std.Err  0.011772 0.04936818 0.07934171 0.1127944 0.1480665 0.1862036
ppar <- person.parameter(fitRSM)
eRm::itemfit(ppar)
## 
## Itemfit Statistics: 
##             Chisq   df p-value Outfit MSQ Infit MSQ Outfit t Infit t Discrim
## cp24p146 6334.600 5365   0.000      1.181     1.146    8.348   6.641   0.512
## cp24p150 4296.121 5365   1.000      0.801     0.776  -10.723 -12.156   0.659
## cp24p155 5387.215 5365   0.413      1.004     0.981    0.205  -0.974   0.656
## cp24p157 5424.623 5365   0.281      1.011     0.999    0.556  -0.041   0.632
## cp24p159 4764.963 5365   1.000      0.888     0.845   -6.020  -8.806   0.669
## cp24p161 4160.161 5365   1.000      0.775     0.770  -12.390 -12.858   0.737
## cp24p162 4219.714 5365   1.000      0.786     0.807  -11.505 -10.262   0.725
## cp24p164 4533.414 5365   1.000      0.845     0.806   -8.319 -10.694   0.717

Infit ve Outfit 0.5 - 1.5 arasında değişir (De ayala, 2005). 0.7 - 1.3 arasında olabileceğine de değinmiştik. Bu biçimiyle Outfit ve Infıt değerleri kabul edilebilir aralıktadır.

5 GPCM Kestirimleri

gpcmfit <- gpcm(data)
gpcmfit
## 
## Call:
## gpcm(data = data)
## 
## Coefficients:
##           Catgr.1  Catgr.2  Catgr.3  Catgr.4  Catgr.5  Catgr.6  Dscrmn
## cp24p146   -1.950   -2.037   -2.683   -1.274   -0.073    2.176   0.617
## cp24p150   -1.987   -1.815   -1.671   -0.311    0.727    2.413   0.968
## cp24p155   -1.396   -1.610   -1.527   -0.274    0.546    2.034   0.822
## cp24p157   -1.406   -1.496   -1.564   -0.232    0.755    2.242   0.807
## cp24p159   -1.098   -1.046   -1.097    0.463    1.344    2.676   0.895
## cp24p161   -1.305   -1.212   -1.251   -0.102    0.775    2.043   1.277
## cp24p162   -1.608   -1.519   -1.485   -0.451    0.504    1.843   1.209
## cp24p164   -1.680   -1.291   -1.063   -0.110    0.722    1.995   1.114
## 
## Log.Lik: -66635.88
plot(gpcmfit)

mirt::plot(gpcmfit, type = "IIC")

mirt::plot(gpcmfit, type = "OCCu")

5.1 Graded Response - Dereceli Tepki

gradedfit <- mirt(data, 1, itemtype = "graded",verbose = F)
plot(gradedfit, type = "trace")

plot(gradedfit, type = "info")

# Model Karşılaştırması

rasch <- mirt(data, 1, itemtype = "Rasch", verbose = F, SE = T)
rsm <- mirt(data, 1, itemtype = "rsm", verbose = F, SE = T)
grsm <- mirt(data, 1, itemtype = "grsm", verbose = F, SE = T)
graded <- mirt(data, 1, itemtype = "graded", verbose = F, SE = T)
gpcm <- mirt(data, 1, itemtype = "gpcm", verbose = F, SE = T)
anova(rasch,rsm, grsm, graded, gpcm)
##             AIC    SABIC       HQ      BIC    logLik        X2  df   p
## rasch  133799.6 133967.2 133912.4 134122.9 -66850.78                  
## rsm    134033.2 134081.1 134065.5 134125.6 -67002.61  -303.646 -35 NaN
## grsm   131910.6 131982.5 131959.0 132049.2 -65934.32  2136.585   7   0
## graded 131460.7 131652.2 131589.7 131830.2 -65674.35   519.938  35   0
## gpcm   133384.7 133576.2 133513.7 133754.2 -66636.34 -1923.995   0 NaN

AIC, saBIC ve BIC değerleri ile Loglikelihood değerleri karşılaştırıldığında en iyi uyum “Graded” modelden elde edilmiştir.

irtparsgrm <- coef(graded,IRTpars = T,simplify = T)
kable(irtparsgrm$items, digits = 2)
a b1 b2 b3 b4 b5 b6
cp24p146 1.31 -3.66 -2.89 -2.20 -1.06 0.20 1.92
cp24p150 1.88 -2.72 -2.02 -1.33 -0.28 0.79 2.18
cp24p155 1.80 -2.35 -1.76 -1.14 -0.23 0.67 1.90
cp24p157 1.80 -2.31 -1.71 -1.12 -0.19 0.80 2.06
cp24p159 1.91 -1.82 -1.21 -0.63 0.35 1.30 2.50
cp24p161 2.49 -1.87 -1.38 -0.91 -0.08 0.79 1.90
cp24p162 2.33 -2.23 -1.71 -1.20 -0.37 0.56 1.73
cp24p164 2.21 -2.12 -1.45 -0.87 -0.07 0.79 1.89

GRM modelinden elde edilen madde ayırt edicilikleri (a) incelendiğinde en yüksek ayırt edicilik madde 161 (a = 2.4884744), en düşük ayırt edicilik ise madde 146’den (a = 1.3074061) elde edilmiştir. Madde güçlüğü için eşik değerleri incelendiğinde ise 7 kategori için 6 eşik değer elde edilmiştir ve b1’den b6’ya doğru bir artış görülmektedir.

plot(graded,type = "trace")

Madde eşik değerlerine ilişkin bilgiler görseldeki gibidir.

graded <- ltm::grm(data)
plot(graded, type = "IIC")

Her bir eşik için madde bilgi fonksiyonları görseldeki gibidir.

plot(gradedfit, type = "info")

Testin en yüksek bilgi verdiği nokta yaklaşık -2 \(\theta\) düzeyindedir. Uç yetenek düzeylerinde bilgi düzeyinin arttığı görülmektedir.

plot(gradedfit, type = "SE")

Yüksek yetenek düzeylerinde teste ilişkin standart hatanın ise arttığı belirlenmiştir.

6 Theta’ların Dağılımı

theta <- fscores(gradedfit,method = "EAP",full.scores.SE = T)
hist(theta[,1])

6.1 LR testi ile ölçme değişmezliği yap

n <- nrow(data)
group_codes <- rep(1:2, length.out = n)
data$GROUP <- as.factor(group_codes)
group <- as.factor(data$GROUP) 
resp <- data[,1:8]
conf <- mirt::multipleGroup(resp,model = 1,itemtype = "graded", group = group, verbose = F)
metric <- mirt::multipleGroup(resp,model = 1,itemtype = "graded", group = group, invariance = c("slopes"), verbose = F)
scalar <- mirt::multipleGroup(resp,model = 1,itemtype = "graded", group = group, invariance = c("slopes","intercepts"), verbose = F)
anova(conf,metric,scalar)
##             AIC    SABIC       HQ      BIC    logLik      X2  df   p
## conf   131502.5 131885.6 131760.5 132241.5 -65639.27                
## metric 131491.0 131846.7 131730.5 132177.2 -65641.51  -4.474  -8 NaN
## scalar 131460.7 131652.2 131589.7 131830.2 -65674.35 -65.671 -48 NaN

Ölçme değişmezliği sağlanmadı olarak algıladım. (LRtest için medyan vb. ayırma vardı ancak mirt’teki bu ÇGDFA’ya daha çok benzediği için kolay geldi buradan yaptım.)

DIF(conf, which.par = c("d1","d2","d3","d4","d5","d6"),p.adjust = "BH",verbose = F)
## NOTE: No hyper-parameters were estimated in the DIF model. 
##       For effective DIF testing, freeing the focal group hyper-parameters is recommended.
##          groups converged    AIC  SABIC     HQ    BIC     X2 df     p adj_p
## cp24p146    1,2      TRUE -7.818 12.704  6.000 31.770 19.818  6 0.003 0.024
## cp24p150    1,2      TRUE  9.268 29.790 23.086 48.856  2.732  6 0.842 0.842
## cp24p155    1,2      TRUE  2.764 23.286 16.582 42.352  9.236  6 0.161 0.643
## cp24p157    1,2      TRUE  8.862 29.385 22.681 48.451  3.138  6 0.791 0.842
## cp24p159    1,2      TRUE  6.409 26.931 20.227 45.997  5.591  6 0.471 0.842
## cp24p161    1,2      TRUE  5.461 25.984 19.280 45.050  6.539  6 0.366 0.842
## cp24p162    1,2      TRUE  7.103 27.625 20.921 46.691  4.897  6 0.557 0.842
## cp24p164    1,2      TRUE  9.220 29.742 23.038 48.808   2.78  6 0.836 0.842
plot(conf, which.item = 1,type = "trace")

Madde 1’in 1. grup için yanlı olduğu belirlenmiştir (p<0.05)