Bu hafta online ders sonrası HSB veri setiyle ders notlarındaki anlatıma paralel olacak şekilde farklı değişkenlerle uygulama yaptım.
Öncelikle boş model kurarak okullar arası ve öğrenci düzeyi varyansları inceledim, ardından rastgele katsayı modeli ile cinsiyetin matematik başarısı üzerindeki etkisinin okullara göre değişip değişmediğini test ettim. Sonrasında kesişim ve eğim katsayılarının Düzey-2 değişkenlerle açıklandığı modeli kurarak okulun ortalama sosyoekonomik düzeyinin hem başarı üzerindeki etkisini hem de cinsiyet farkını açıklayıp açıklamadığını incelemeye çalıştım.
Geleneksel istatistiksel yöntemler veriyi ya birey düzeyinde ya da grup düzeyinde toplulaştırır. Bu durum; bağımsızlık varsayımı ihlaline yol açar. Çünkü aynı gruptaki bireyler (örneğin aynı okuldan öğrenciler) birbirine benzer. Ancak klasik regresyon bunları bağımsız kabul eder . Buna bağlı olarak gözlemlerin bağımsızlığı varsayımını ihlal ederek hatalı standart hatalara ve bilgi kaybına neden olur.
Hiyerarşik doğrusal modellerin (HLM) temel mantığı, iç içe geçmiş (nested) veri yapılarını uygun şekilde analiz etmektir. Bu yapıya göre, bireyler (Düzey 1) daha büyük birimlerin (Düzey 2), bu birimler de bazen daha da büyük organizasyonların (Düzey 3) içinde yer alır.
HLM, her bir grup için modele benzersiz bir rassal etki (random effect) dahil ederek bu bağımlılığı istatistiksel olarak hesaba katar ve böylece hatalı (genellikle olması gerekenden küçük çıkan) standart hataların önüne geçer.
Varyans Ayrıştırma: Varyansın ne kadarının birey düzeyinden (Düzey-1), ne kadarının grup düzeyinden (Düzey-2) kaynaklandığını açıklar.
Eğimlerin Değişkenliği: Gruplar arasındaki ilişkilerin (eğimlerin) farklılaşmasına izin vererek daha gerçekçi tahminler sunar.
Rastgele Etkiler ANOVA Modeli (Random Effects ANOVA Model), en basit HLM modelidir. Sadece-Kesişim Modeli (Intercept-only model) veya Tamamen Koşulsuz Model (Fully Unconditional Model) olarak da adlandırılır. Her iki düzeyde de yordayıcı değişken bulunmaz.
library(tidyverse)## Warning: package 'ggplot2' was built under R version 4.5.3## Warning: package 'dplyr' was built under R version 4.5.3## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.0 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.2
## ✔ ggplot2 4.0.2 ✔ tibble 3.3.0
## ✔ lubridate 1.9.4 ✔ tidyr 1.3.1
## ✔ purrr 1.1.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(sjmisc) # frekan tabloları## Warning: package 'sjmisc' was built under R version 4.5.3##
## Attaching package: 'sjmisc'
##
## The following object is masked from 'package:purrr':
##
## is_empty
##
## The following object is masked from 'package:tidyr':
##
## replace_na
##
## The following object is masked from 'package:tibble':
##
## add_caselibrary(car) #özet istatitik## Warning: package 'car' was built under R version 4.5.3## Zorunlu paket yükleniyor: carData
##
## Attaching package: 'car'
##
## The following object is masked from 'package:dplyr':
##
## recode
##
## The following object is masked from 'package:purrr':
##
## somelibrary(lme4) # lmer() functions## Warning: package 'lme4' was built under R version 4.5.2## Zorunlu paket yükleniyor: Matrix
##
## Attaching package: 'Matrix'
##
## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpack
##
## Registered S3 method overwritten by 'lme4':
## method from
## na.action.merMod car# library(nlme) #
library(interactions) # grafikler## Warning: package 'interactions' was built under R version 4.5.2load("HSB.RData")
colnames(HSB)## [1] "school" "minority" "female" "ses" "mathach" "size"
## [7] "sector" "pracad" "disclim" "himinty" "meanses"model_null <- lmer(mathach ~ 1 + (1 | school), data = HSB)
summary(model_null)## Linear mixed model fit by REML ['lmerMod']
## Formula: mathach ~ 1 + (1 | school)
## Data: HSB
##
## REML criterion at convergence: 47116.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.0631 -0.7539 0.0267 0.7606 2.7426
##
## Random effects:
## Groups Name Variance Std.Dev.
## school (Intercept) 8.614 2.935
## Residual 39.148 6.257
## Number of obs: 7185, groups: school, 160
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 12.6370 0.2444 51.71modelin matematiksel hali
mathachij =β0+uj+eij
Modelin anlamı: Sadece ortalama var (bağımsız değişken yok) Okulların ortalamaları birbirinden farklı olabilir Bu farklar rastgele etki (uj) ile modellenir
1.Sabit etki (Genel ortalama) :öğrencilerin genel ortalama matematik puanı yaklaşık 12.64
okullar arası varyans 8.614(Okulların ortalamaları birbirinden ne kadar farklı?) öğrenci düzeyi varyans 39.148 (Aynı okul içindeki öğrenciler ne kadar farklı?)
Genel ortalama için %95 güven aralığı
95%CI(γ00)=γ00±1.96(SE) 95%CI(γ00)=12.637 ±1.96(0.2444)=(12.158,13.116)
Yorum: Genel ortalamanın gerçek değeri %95 olasılıkla 12.158 ile 13.116 arasında olacaktır.
Olası Değerler Aralığı (Plausible Values Range): Sınıf ortalamalarının evrende nasıl dağıldığını görmek için grup-arası varyans (τ00=0.87 ) kullanılır:
12.637 ± 1.96 √8,64=(6.884, 18.390)
Yorum: Bu sonuç okulların ortalama matematik başarılarının evrende yaklaşık olarak 6.88 ile 18.39 arasında dağılabileceğini gösterir.
Rastgele etkiler ANOVA modeli bize şu kritik varyans bileşenlerini sağlar:
σ2 : Grup-içi değişkenlik (Düzey 1 hata varyansı)
τ00: Gruplar-arası değişkenlik (Düzey 2 hata varyansı)
Bu değerler kullanılarak varyansın yüzde kaçının Düzey 2’den kaynaklandığını gösteren ICC (ρ) hesaplanır:
8.614 / (8.614 + 39.148)## [1] 0.1803526Matematik başarısındaki toplam varyansın yaklaşık %18’i okullar arası farklardan, %82’si öğrenciler arası farklardan kaynaklanmaktadır.
Yani aynı okul içindeki öğrenciler birbirine tamamen bağımsız değil, Okul etkisi var.
Denklemler
Level-1 (öğrenci düzeyi) mathachij =β0+eij
mathachij → i. öğrencinin j. okuldaki matematik puanı
β0j → j. okulun ortalaması
eij → öğrenciye ait hata (bireysel sapma)Level-2 (okul düzeyi) β0j=γ00+uj
γ00 → tüm okulların genel ortalaması (12.637) uj→ j. okulun genel ortalamadan sapması (rastgele etki)
*Birleştirilmiş model mathachij = γ00+uj+eij mathachij = 12,637+uj+eij
sadece ortalama değil, bir değişkenin etkisinin (eğim) okullara göre değişip değişmediğine bakıyoruz.
Soru: Cinsiyet farkı her okulda aynı mı?
model_random <- lmer(mathach ~ female + (1 + female | school), data = HSB)
summary(model_random)## Linear mixed model fit by REML ['lmerMod']
## Formula: mathach ~ female + (1 + female | school)
## Data: HSB
##
## REML criterion at convergence: 47053.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2087 -0.7512 0.0364 0.7710 2.6480
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## school (Intercept) 8.6604 2.9429
## femaleFemale 0.6842 0.8272 -0.30
## Residual 38.7331 6.2236
## Number of obs: 7185, groups: school, 160
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.3454 0.2619 50.950
## femaleFemale -1.3817 0.1867 -7.402
##
## Correlation of Fixed Effects:
## (Intr)
## femaleFemal -0.422VarCorr(model_random)## Groups Name Std.Dev. Corr
## school (Intercept) 2.94286
## femaleFemale 0.82715 -0.304
## Residual 6.22359MODELİN DENKLEMİ
Yij=γ00+γ10Xij+u0j+u1jXij+rij Xij=female γ00: Erkek öğrenciler için genel ortalama başarı (13,3454) γ10: Kız ve erkek öğrenciler arasındaki ortalama fark (-1,3817) u0j: Okulların ortalamalarının farklı olması u1j: female etkisinin okullara göre değişmesi Intercept (13.345): Erkek öğrenciler için ortalama matematik puanı
(female=0 , referans grup erkek)
Kız öğrencilerin puanı erkeklere göre 1.38 puan daha düşük
kesişim varyansı τ00=8.6604 eğim varyansı τ11=0.6842
boş model artık varyansı σ2=39.148 rastgele katsayılı model artık varyansı σ2=38.7331
D1 kesişimi
β0j için %95 olası değerler aralığı 13.3454 ±1.96 √8,6604= (7.58, 19.11)
Yorum: Sınıflardaki/okullardaki erkek öğrencilerin ortalama matematik başarılarının yaklaşık %95’i 7.58 ile 19.11 aralığında yer almaktadır.
β1j için %95 olası değerler aralığı -1.3817 ±1.96 √0,6842= (−3.00, 0.24)
Yorum: Okullardaki cinsiyet etkilerinin yaklaşık %95’i -3.00 ile 0.24 aralığında yer almaktadır. Bu sonuç, çoğu okulda kız öğrencilerin matematik puanlarının erkek öğrencilere göre daha düşük olma eğiliminde olduğunu, ancak bazı okullarda bu farkın oldukça küçüldüğünü, hatta sıfıra yaklaşabildiğini göstermektedir. Aralığın 0’ı içermesi, cinsiyet etkisinin her okulda aynı yön ve büyüklükte olmadığını düşündürür.
# Boş modelin öğrenci düzeyi varyansı
var_null <- 39.148
# Rastgele katsayılı modelin öğrenci düzeyi varyansı
var_random <- 38.7331
# Açıklanan varyans oranı
R1_2 <- (var_null - var_random) / var_null
R1_2## [1] 0.01059824Yorum: Cinsiyet değişkeni, matematik başarısındaki öğrenci düzeyi varyansın yaklaşık %1.1’ini açıklamaktadır.
# Boş model (ML)
model_null_ml <- lmer(mathach ~ 1 + (1 | school), data = HSB, REML = FALSE)
# Rastgele katsayı modeli (ML)
model_random_ml <- lmer(mathach ~ female + (1 + female | school), data = HSB, REML = FALSE)## boundary (singular) fit: see help('isSingular')deviance(model_null_ml)## [1] 47115.81deviance(model_random_ml)## [1] 47052.94anova(model_null_ml, model_random_ml) ## Data: HSB
## Models:
## model_null_ml: mathach ~ 1 + (1 | school)
## model_random_ml: mathach ~ female + (1 + female | school)
## npar AIC BIC logLik -2*log(L) Chisq Df Pr(>Chisq)
## model_null_ml 3 47122 47142 -23558 47116
## model_random_ml 6 47065 47106 -23527 47053 62.868 3 1.433e-13 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Null model: 47116 Random slope model: 47053
47053 < 47116 Rastgele katsayılı model, boş modele göre veriye anlamlı derecede daha iyi uyum sağlamaktadır
Bu modelde;
öğrencinin cinsiyetinin matematik başarısı üzerindeki etkisinin, okulun ortalama sosyoekonomik düzeyine göre değişip değişmediği incelenmektedir.
Öncesinde, Okulun ortalama sosyoekonomik düzeyi, modelde katsayıların anlamlı ve yorumlanabilir hale gelmesi, özellikle kesişim teriminin gerçekçi bir referans noktası (ortalama SES düzeyi) üzerinden yorumlanabilmesi için merkezlenmiştir.
HSB$meanses_c <- scale(HSB$meanses, center = TRUE, scale = FALSE) #merkezlemelibrary(lmerTest)## Warning: package 'lmerTest' was built under R version 4.5.3##
## Attaching package: 'lmerTest'## The following object is masked from 'package:lme4':
##
## lmer## The following object is masked from 'package:stats':
##
## stepmodel_cross <- lmer(
mathach ~ female * meanses_c + (1 + female | school),
data = HSB
)## Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model failed to converge with max|grad| = 0.00201792 (tol = 0.002, component 1)
## See ?lme4::convergence and ?lme4::troubleshooting.summary(model_cross)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: mathach ~ female * meanses_c + (1 + female | school)
## Data: HSB
##
## REML criterion at convergence: 46893
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2609 -0.7466 0.0348 0.7688 2.6899
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## school (Intercept) 3.1668 1.7796
## femaleFemale 0.8102 0.9001 -0.61
## Residual 38.7274 6.2231
## Number of obs: 7185, groups: school, 160
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.38304 0.18446 139.44750 72.554 < 2e-16 ***
## femaleFemale -1.36244 0.18341 124.62059 -7.428 1.52e-11 ***
## meanses_c 5.68247 0.45550 142.53472 12.475 < 2e-16 ***
## femaleFemale:meanses_c 0.07055 0.45942 135.87948 0.154 0.878
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) fmlFml mnss_c
## femaleFemal -0.635
## meanses_c -0.019 0.018
## fmlFml:mns_ 0.020 0.018 -0.653
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00201792 (tol = 0.002, component 1)
## See ?lme4::convergence and ?lme4::troubleshooting.VarCorr(model_cross)## Groups Name Std.Dev. Corr
## school (Intercept) 1.77956
## femaleFemale 0.90013 -0.615
## Residual 6.22314Modelde yakınsama (convergence) sorunu gözlendi.
Sabit Etkiler
Intercept = 13.383 Ortalama sosyoekonomik düzeye sahip bir okulda, erkek öğrencilerin beklenen matematik başarı puanı 13.383’tür.
female = -1.362 Ortalama sosyoekonomik düzeyde, kız öğrencilerin matematik başarı puanları erkek öğrencilere göre ortalama 1.362 puan daha düşüktür.
meanses_c = 5.682 Okulun ortalama sosyoekonomik düzeyi bir birim arttığında, erkek öğrencilerin matematik başarı puanları ortalama 5.682 puan artmaktadır.
interaction = 0.0706 Okulun sosyoekonomik düzeyindeki artışın, kız ve erkek öğrenciler arasındaki başarı farkı üzerindeki etkisi çok zayıf olup anlamlı değildir.
Rastgele Etkiler
(Intercept) varyansı = 3.1668
Okulların ortalama matematik başarıları birbirinden farklıdır ve bu farklılıkların varyansı 3.1668’dir.
female (eğim) varyansı = 0.8102
Cinsiyetin matematik başarısı üzerindeki etkisi okullar arasında değişmektedir ve bu değişimin varyansı 0.8102’dir.
Kesişim (Intercept) Önceki model (random slope): 8.6604 Bu model: 3.1668
var_int_1 <- 8.6604
var_int_2 <- 3.1668
R2_intercept <- (var_int_1 - var_int_2) / var_int_1
R2_intercept## [1] 0.6343356Yorum: Okul SES değişkeni, okullar arası ortalama farkların yaklaşık %63’ünü açıklamaktadır.
var_slope_1 <- 0.6842
var_slope_2 <- 0.8102
R2_slope <- (var_slope_1 - var_slope_2) / var_slope_1
R2_slope## [1] -0.1841567Yorum: Cinsiyet etkisindeki değişkenlik, okul SES ile açıklanamamıştır.
deviance(model_null)## Warning in deviance.merMod(model_null): deviance() is deprecated for REML fits;
## use REMLcrit for the REML criterion or deviance(.,REML=FALSE) for deviance
## calculated at the REML fit## [1] 47116.79deviance(model_random)## Warning in deviance.merMod(model_random): deviance() is deprecated for REML
## fits; use REMLcrit for the REML criterion or deviance(.,REML=FALSE) for
## deviance calculated at the REML fit## [1] 47053.86deviance(model_cross)## Warning in deviance.merMod(model_cross): deviance() is deprecated for REML
## fits; use REMLcrit for the REML criterion or deviance(.,REML=FALSE) for
## deviance calculated at the REML fit## [1] 46893.03Boş model: 47116.79 Rastgele katsayı modeli: 47053.86 Kesişim ve eğim katsayılarının bağımlı değişken olduğu model: 46893.03
Yorum: Kesişim ve eğim katsayılarının Düzey-2 değişkenleriyle modellendiği bu son model, sapma değerinin en düşük olması nedeniyle veriye diğer iki modelden çok daha iyi uyum sağlamaktadır.