Öncelikle, dersteki kodları tekrar çalıştırdım bu kısımda yapılanları tekrar etmiş oldum.
# 1. aşama madde parametreleri
set.seed(37)
madde <- 8
# Sırasıyla a,b,c parametreleri
maddepar <- cbind(
rnorm( madde, mean = 0,sd =0.15),
rnorm( madde, mean = 0.40, sd =0.20),
rnorm( madde, mean = 0.10, sd =0.05)
)
maddepar ## [,1] [,2] [,3]
## [1,] 0.01871310 0.57119088 -0.03524248
## [2,] 0.05731119 0.44319910 0.01315379
## [3,] 0.08688642 0.32445958 0.12015376
## [4,] -0.04406222 0.40773871 0.06339636
## [5,] -0.12425237 0.68496301 0.08608195
## [6,] -0.04990704 0.59646198 0.16421153
## [7,] -0.02882393 0.46209289 0.14120147
## [8,] 0.20444741 0.06649433 -0.01922178# yetenek paremetreleri
set.seed(37)
birey <- 1000
yetenek <- rnorm( birey,0,1)
head(yetenek)## [1] 0.1247540 0.3820746 0.5792428 -0.2937481 -0.8283492 -0.3327136# irtoys sim fonksiyonu her bir bireyin her bir maddeyi doğru cevaplama olasılığını seçilen modele göre hesaplar.
# 1000 kişiye ait simülatif veri
veri <- irtoys::sim(ip=maddepar,x = yetenek)
head(veri)## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,] 0 0 1 1 0 1 1 0
## [2,] 1 0 0 1 1 1 1 1
## [3,] 0 0 1 1 1 1 1 1
## [4,] 1 1 0 1 0 1 0 0
## [5,] 1 0 1 0 0 0 1 1
## [6,] 0 1 0 0 1 0 1 0colnames(veri) <- paste0("madde",1:madde)
head(veri)## madde1 madde2 madde3 madde4 madde5 madde6 madde7 madde8
## [1,] 0 0 1 1 0 1 1 0
## [2,] 1 0 0 1 1 1 1 1
## [3,] 0 0 1 1 1 1 1 1
## [4,] 1 1 0 1 0 1 0 0
## [5,] 1 0 1 0 0 0 1 1
## [6,] 0 1 0 0 1 0 1 0library(irtoys)
ip=maddepar
x = yetenek
i = irf(ip = ip, x = x)
head(i$x) ## [1] 0.1247540 0.3820746 0.5792428 -0.2937481 -0.8283492 -0.3327136head(i$f) ## [,1] [,2] [,3] [,4] [,5] [,6] [,7]
## [1,] 0.4802166 0.5020744 0.5562603 0.5346178 0.5589384 0.5870245 0.5726883
## [2,] 0.4814628 0.5057126 0.5611780 0.5319630 0.5516387 0.5843414 0.5710959
## [3,] 0.4824178 0.5085005 0.5649460 0.5299287 0.5460422 0.5822853 0.5698758
## [4,] 0.4781898 0.4961585 0.5482647 0.5389350 0.5707915 0.5913873 0.5752781
## [5,] 0.4756010 0.4886060 0.5380632 0.5444480 0.5858764 0.5969573 0.5785857
## [6,] 0.4780011 0.4956078 0.5475206 0.5393369 0.5718935 0.5917934 0.5755192
## [,8]
## [1,] 0.4934241
## [2,] 0.5068233
## [3,] 0.5170760
## [4,] 0.4716310
## [5,] 0.4439024
## [6,] 0.4696042# f: her bir bireyin bir maddeyi cevaplama olasılığı d = dim(i$f)
u = runif(d[1] * d[2])
dim(u) = d
head(ifelse(i$f > u, 1, 0))## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,] 0 1 0 1 1 1 1 0
## [2,] 0 0 0 1 0 0 0 1
## [3,] 0 0 0 0 1 1 1 1
## [4,] 0 1 0 1 1 1 0 0
## [5,] 0 0 0 1 1 1 1 1
## [6,] 1 1 1 1 1 1 0 1veri_uretimi <- function(maddesay, bireysay, seed){
# seed fonksiyonu
set.seed(seed)
# madde parametresi üretimi
maddepar <- cbind(
rnorm( madde, mean = 1.13 , sd =0.25)*1.702,
rnorm( madde, mean = 0.30 , sd =0.50)*1.702,
rnorm( madde, mean = 0.16, sd =0.05))
# yetenek paremetresi
yetenek <- rnorm(bireysay,0,1)
# 3pl modele göre veri üretimi
cevaplar <- irtoys::sim(maddepar, yetenek)
colnames(cevaplar) <- paste0("madde",1:madde)
veri <- list(maddepar = maddepar,
yetenek = yetenek,
seed=seed,
cevaplar=cevaplar)
return(veri)
}veri_1 <- veri_uretimi(maddesay=8,bireysay = 1000, seed=37)## kestirilen par isimli fonksiyon, veri_uretimi fonksiyonunun cevaplar bileşeni girdi olarak alır.
kestirilen_par <- function(veri, par=3){
if(par==3){
model <- mirt::mirt(veri,
1,
itemtype = "3PL",
verbose=FALSE)
}else if(par ==2){
model <- mirt::mirt(veri,
1,
itemtype = "2PL",
verbose=FALSE)
}else{
model <- mirt::mirt(veri,
"1PL",
itemtype = "Rasch",
verbose=FALSE)
}
kestirim <- mirt::coef(model,simplify=TRUE, IRTpars=TRUE)$item[,1:3]
kestirim
}kestirilen_par(veri_1$cevaplar)## a b g
## madde1 1.854323 1.1405108 2.419344e-03
## madde2 1.649962 0.6696361 1.820476e-02
## madde3 2.958612 0.2888204 2.293011e-01
## madde4 1.685040 0.3195506 8.491345e-03
## madde5 0.643828 1.9113242 3.182818e-05
## madde6 2.131575 1.3631337 2.610244e-01
## madde7 2.658378 0.9415949 2.743031e-01
## madde8 3.014089 -0.3691584 3.834051e-01veri_1$maddepar## [,1] [,2] [,3]
## [1,] 1.976343 1.2390172 0.02475752
## [2,] 2.085833 0.6944122 0.07315379
## [3,] 2.169728 0.1891755 0.18015376
## [4,] 1.798270 0.5435282 0.12339636
## [5,] 1.570797 1.7231176 0.14608195
## [6,] 1.781690 1.3465457 0.22421153
## [7,] 1.841496 0.7748053 0.20120147
## [8,] 2.503209 -0.9084666 0.04077822# kestirdiğimiz parametreler, istediğimiz parametreler arasında çok fark var.# Hata hesaplama
# Tek tekrar yeterli değil, hata yüksek oluyor. Tekrar sayısı arttıkça SE ve RMSE değerleri düşer ve tahminler kararlı hale gelir.
# RMSE^2 = BIAS^2 + SE^2
hata <- function(kestirilen, gercek){
result <- data.frame(parametreler = c("a","b","c"),
rmse=sapply(1:3, function(i) sqrt(mean((kestirilen[,i] - gercek[,i])^2) )),
bias=sapply(1:3, function(i) mean(kestirilen[,i]-gercek[,i])),
se= sapply(1:3, function(i) sd(kestirilen[,i])))
result
}
hata(kestirilen_par(veri_1$cevaplar),veri_1$maddepar )# TEKRAR VERİ ÜRETİMİ
veri_uretimi <- function(maddesay, bireysay, seed = NULL){
if(!is.null(seed)){
set.seed(seed)
}else{
seed <-sample.int(10000,1)
set.seed(seed)
}
# seed fonksiyonu
cat("atanan seed", seed, "\n")
# madde parametresi üretimi
print("madde parametreleri üretildi")
maddepar <- cbind(
rlnorm( maddesay, meanlog = 0 , sdlog =0.5)*1.702,
rnorm( maddesay, mean = 0.30 , sd =0.50)*1.702,
rnorm( maddesay, mean = 0.16, sd =0.05)
)
# yetenek parametresi üretimi
yetenek <- rnorm(bireysay,0,1)
print("yetenek parametreleri üretildi")
# 3pl modele göre veri üretimi
cevaplar <- irtoys::sim(maddepar, yetenek)
colnames(cevaplar) <- paste0("madde",1:maddesay)
print("bileşenler birleştiriliyor")
veri <- list(maddepar = maddepar,
yetenek = yetenek,
seed=seed,
cevaplar=cevaplar)
return(veri)
}veri2 <- veri_uretimi(10,1000)## atanan seed 21
## [1] "madde parametreleri üretildi"
## [1] "yetenek parametreleri üretildi"
## [1] "bileşenler birleştiriliyor"veriler <- replicate(5L,veri_uretimi(10,1000)) ## atanan seed 1098
## [1] "madde parametreleri üretildi"
## [1] "yetenek parametreleri üretildi"
## [1] "bileşenler birleştiriliyor"
## atanan seed 2149
## [1] "madde parametreleri üretildi"
## [1] "yetenek parametreleri üretildi"
## [1] "bileşenler birleştiriliyor"
## atanan seed 1606
## [1] "madde parametreleri üretildi"
## [1] "yetenek parametreleri üretildi"
## [1] "bileşenler birleştiriliyor"
## atanan seed 4814
## [1] "madde parametreleri üretildi"
## [1] "yetenek parametreleri üretildi"
## [1] "bileşenler birleştiriliyor"
## atanan seed 1493
## [1] "madde parametreleri üretildi"
## [1] "yetenek parametreleri üretildi"
## [1] "bileşenler birleştiriliyor"# tekrar sayısı 5
lapply(1:5, function(i) veriler[,i][1])## [[1]]
## [[1]]$maddepar
## [,1] [,2] [,3]
## [1,] 6.7490143 0.7312388 0.16687171
## [2,] 1.9471021 0.9890589 0.15315606
## [3,] 1.3498383 -1.0891373 0.07978747
## [4,] 1.0360599 0.3043721 0.18651368
## [5,] 1.7002150 0.7604314 0.08946298
## [6,] 0.9238229 2.1551889 0.10173746
## [7,] 1.0383398 -0.1562188 0.12774338
## [8,] 1.4036868 0.1050944 0.10998328
## [9,] 2.8123108 2.1047701 0.26053763
## [10,] 1.3944989 0.5365599 0.19725134
##
##
## [[2]]
## [[2]]$maddepar
## [,1] [,2] [,3]
## [1,] 1.0732286 -0.1133397 0.1581915
## [2,] 1.7810067 1.1078248 0.2082098
## [3,] 1.4059870 -0.5673548 0.1387106
## [4,] 3.0259620 0.1979059 0.1551955
## [5,] 0.9806083 -0.7460035 0.1438454
## [6,] 0.8312758 0.6716454 0.0972568
## [7,] 4.0087481 -0.9634126 0.1074729
## [8,] 2.2957808 0.1739008 0.2053455
## [9,] 0.7202705 -1.1924169 0.1621419
## [10,] 1.4411788 0.3788448 0.1811188
##
##
## [[3]]
## [[3]]$maddepar
## [,1] [,2] [,3]
## [1,] 1.725826 1.27133802 0.14675258
## [2,] 2.259488 -0.13124641 0.13906841
## [3,] 1.025015 1.47948210 0.18939314
## [4,] 2.010506 -0.52046096 0.22421886
## [5,] 1.063474 0.26097513 0.20923462
## [6,] 4.416694 1.40542189 0.15876380
## [7,] 0.946570 0.77246813 0.07682463
## [8,] 1.145450 0.46902732 0.20864395
## [9,] 4.701380 1.97989268 0.11663079
## [10,] 2.983169 0.09222526 0.23412805
##
##
## [[4]]
## [[4]]$maddepar
## [,1] [,2] [,3]
## [1,] 4.7981807 1.0320508 0.09801986
## [2,] 2.5896632 1.4692000 0.23225315
## [3,] 0.8596521 0.5792187 0.15377904
## [4,] 1.6188236 0.8638407 0.08980358
## [5,] 1.3292055 1.9717870 0.14933491
## [6,] 1.4198357 1.7814724 0.14058574
## [7,] 1.5736118 0.5559212 0.16916567
## [8,] 1.6863692 -0.3992936 0.18468751
## [9,] 1.2843819 -0.2638391 0.18542958
## [10,] 0.8954598 1.1705315 0.18677579
##
##
## [[5]]
## [[5]]$maddepar
## [,1] [,2] [,3]
## [1,] 2.6196720 0.02649782 0.13641306
## [2,] 2.4017650 0.18059299 0.14735100
## [3,] 0.9307662 -0.78936881 0.13818842
## [4,] 1.1846370 -0.38372662 0.13047796
## [5,] 1.8614056 0.72589475 0.11345486
## [6,] 1.4014350 0.23534865 0.08837373
## [7,] 1.6727941 0.42374940 0.19319209
## [8,] 1.4292553 0.92290279 0.13807274
## [9,] 3.4048163 0.60729769 0.28504693
## [10,] 0.9748528 0.19234520 0.10530914sapply(1:5, function(i) veriler[,i][1])## $maddepar
## [,1] [,2] [,3]
## [1,] 6.7490143 0.7312388 0.16687171
## [2,] 1.9471021 0.9890589 0.15315606
## [3,] 1.3498383 -1.0891373 0.07978747
## [4,] 1.0360599 0.3043721 0.18651368
## [5,] 1.7002150 0.7604314 0.08946298
## [6,] 0.9238229 2.1551889 0.10173746
## [7,] 1.0383398 -0.1562188 0.12774338
## [8,] 1.4036868 0.1050944 0.10998328
## [9,] 2.8123108 2.1047701 0.26053763
## [10,] 1.3944989 0.5365599 0.19725134
##
## $maddepar
## [,1] [,2] [,3]
## [1,] 1.0732286 -0.1133397 0.1581915
## [2,] 1.7810067 1.1078248 0.2082098
## [3,] 1.4059870 -0.5673548 0.1387106
## [4,] 3.0259620 0.1979059 0.1551955
## [5,] 0.9806083 -0.7460035 0.1438454
## [6,] 0.8312758 0.6716454 0.0972568
## [7,] 4.0087481 -0.9634126 0.1074729
## [8,] 2.2957808 0.1739008 0.2053455
## [9,] 0.7202705 -1.1924169 0.1621419
## [10,] 1.4411788 0.3788448 0.1811188
##
## $maddepar
## [,1] [,2] [,3]
## [1,] 1.725826 1.27133802 0.14675258
## [2,] 2.259488 -0.13124641 0.13906841
## [3,] 1.025015 1.47948210 0.18939314
## [4,] 2.010506 -0.52046096 0.22421886
## [5,] 1.063474 0.26097513 0.20923462
## [6,] 4.416694 1.40542189 0.15876380
## [7,] 0.946570 0.77246813 0.07682463
## [8,] 1.145450 0.46902732 0.20864395
## [9,] 4.701380 1.97989268 0.11663079
## [10,] 2.983169 0.09222526 0.23412805
##
## $maddepar
## [,1] [,2] [,3]
## [1,] 4.7981807 1.0320508 0.09801986
## [2,] 2.5896632 1.4692000 0.23225315
## [3,] 0.8596521 0.5792187 0.15377904
## [4,] 1.6188236 0.8638407 0.08980358
## [5,] 1.3292055 1.9717870 0.14933491
## [6,] 1.4198357 1.7814724 0.14058574
## [7,] 1.5736118 0.5559212 0.16916567
## [8,] 1.6863692 -0.3992936 0.18468751
## [9,] 1.2843819 -0.2638391 0.18542958
## [10,] 0.8954598 1.1705315 0.18677579
##
## $maddepar
## [,1] [,2] [,3]
## [1,] 2.6196720 0.02649782 0.13641306
## [2,] 2.4017650 0.18059299 0.14735100
## [3,] 0.9307662 -0.78936881 0.13818842
## [4,] 1.1846370 -0.38372662 0.13047796
## [5,] 1.8614056 0.72589475 0.11345486
## [6,] 1.4014350 0.23534865 0.08837373
## [7,] 1.6727941 0.42374940 0.19319209
## [8,] 1.4292553 0.92290279 0.13807274
## [9,] 3.4048163 0.60729769 0.28504693
## [10,] 0.9748528 0.19234520 0.10530914library(doParallel)
detectCores()## [1] 4cl <- makeCluster(2)
registerDoParallel(cl)
foreach(i=1:4) %dopar% sqrt(i)## [[1]]
## [1] 1
##
## [[2]]
## [1] 1.414214
##
## [[3]]
## [1] 1.732051
##
## [[4]]
## [1] 2foreach(i=1:4,j=1:4) %dopar% { sqrt(i * j) }## [[1]]
## [1] 1
##
## [[2]]
## [1] 2
##
## [[3]]
## [1] 3
##
## [[4]]
## [1] 4stopCluster(cl)cl <- makeCluster(2)
registerDoParallel(cl)
foreach(i=1:5) %dopar% {sum(rnorm(1e8))}## [[1]]
## [1] 6965.615
##
## [[2]]
## [1] -7357.249
##
## [[3]]
## [1] 9084.347
##
## [[4]]
## [1] 2789.959
##
## [[5]]
## [1] 12209.85stopCluster(cl)cl <- makeCluster(2)
registerDoParallel(cl)
tekrar=100L
maddesay=10
bireysay=1000
simulasyon <- foreach(i=1:tekrar,
.packages = c("mirt","irtoys","doParallel"),
.combine = rbind) %dopar% {
adim1 <-veri_uretimi(maddesay = maddesay,
bireysay = bireysay)
adim2 <- kestirilen_par(adim1$cevaplar)
hata(adim2, adim1$maddepar)
}stopCluster(cl)library(tidyverse)
simulasyon_v1 <- simulasyon %>%
group_by(parametreler) %>%
summarise(rmse = round(mean(rmse),2),
bias = round(mean(bias),2),
se = round(mean(se),2)) %>%
mutate(maddesay= maddesay,
bireysay= bireysay) %>%
as.data.frame()
simulasyon_v1Gerçek Veriden Parametre Kestirip Yeni Veri Üreten Fonksiyon
veri_uret<- function(veri, bireysay, seed = NULL){
if(!is.null(seed)){
set.seed(seed)
} else {
seed <- sample.int(10000,1)
set.seed(seed)
}
cat("atanan seed:", seed, "\n")
maddesay <- ncol(veri)
model <- mirt::mirt(veri, 1, itemtype="3PL", verbose=FALSE)
maddepar <- mirt::coef(model, IRTpars=TRUE, simplify=TRUE)$item[,1:3]
yetenek <- rnorm(bireysay, 0, 1)
cevaplar <- irtoys::sim(maddepar, yetenek)
colnames(cevaplar) <- paste0("madde", 1:maddesay)
return(list(
orijinal_veri = veri,
maddepar = maddepar,
yetenek = yetenek,
yeni_veri = cevaplar
))
}DENEME
library(mirt)
data <- expand.table(LSAT7)
head(data)sonuc <- veri_uret(data, bireysay = 1000)## atanan seed: 9994head(sonuc$yeni_veri)## madde1 madde2 madde3 madde4 madde5
## [1,] 1 1 1 0 1
## [2,] 1 1 1 0 0
## [3,] 0 1 1 0 1
## [4,] 1 1 1 1 1
## [5,] 1 1 1 1 1
## [6,] 0 0 0 1 0IRF fonksiyonu yerine yeni fonksiyon yazma:
Madde Tepki Kuramı’nda IRF (Item Response Function), bir kişinin belirli bir yetenek düzeyinde (θ) bir maddeyi doğru cevaplama olasılığını gösteren fonksiyondur.
1PL modelinde (Rasch modeli) IRF:
\[ P(\theta)=\frac{1}{1+e^{-(\theta-b)}} \]
Burada: θ: yetenek düzeyi b: madde güçlüğü e: üssel dönüşüm
# 1PL IRF fonksiyonu
irf_1pl <- function(theta, b) {
pay <- 1
payda <- 1 + exp(-(theta - b))
p <- pay / payda
return(p)
}2PL için formül ve fonksiyon \[ P(\theta)=\frac{1}{1+e^{-a(\theta-b)}} \]
irf_2pl <- function(theta, a, b) {
pay <- 1
payda <- 1 + exp(-a * (theta - b))
p <- pay / payda
return(p)
}3PL için formül ve fonksiyon
\[ P(\theta)=c+\frac{1-c}{1+e^{-a(\theta-b)}} \] a:ayırt edicilik b:güçlük c:tahmin/şans parametresi
irf_3pl <- function(theta, a, b, c) {
pay <- 1 - c
payda <- 1 + exp(-a * (theta - b))
p <- c + (pay / payda)
return(p)
}Ders videosunu tekrar en baştan izledim, bu şekilde simülasyon mantığını biraz daha iyi anladım. sim fonksiyonu ile maddelerin ayırt edicilik, güçlük, şans parametrelerini ve bireylerin yetenek düzeylerini bir araya getirerek 0-1 matrisi ürettik. Daha sonra üretilen bu yapay verilerden mirt paketi aracılığıyla parametre kestirimi yaptık. Gerçek değerler ile kestirilen değerler arasındaki farka ilişkin RMSE, BIAS ve SE değerleri hesaplayan bir hata fonksiyonu yazdık. Aslında bu kısımda, tek bir simülasyon tekrarının yüksek hatalar verebileceğini, tahminlerin kararlı hale gelmesi için tekrar sayısının artırılması gerektiğini somut bir şekilde görmüş olduk.doParallel ile işlem süresini ciddi oranda kısaltabileceğimizi öğrendik. Simülasyon işi ilginç ve eğlenceli geldi, özellikle gerçek hayatta veri toplamanın güçlüğünü düşününce oturduğun yerden veri üretmek güzel işmiş dedim kendi kendime:) Tabi bu simülasyon verisi ile yapılan çalışmalar ne kadar kabul görüyor literatürde orası aklımda bir soru işareti yarattı. Sanırım bu son öğrenme günlüğümüz, dokuz ayda ne kadar yol katettiğimi düşününce gerçekten kendimle gurur duyuyorum, emekleriniz için çok teşekkürler hocam…