Wczytanie danych z ankiet
# KROK 2: Wczytanie danych ------------------------------------------------
plik <- "data/wyniki_ankiety.xlsx"
# Sprawdzenie nazw arkuszy
excel_sheets(plik)
## [1] "Liczba odpowiedzi 1"
# Wczytanie pierwszego arkusza
dane_raw <- read_excel(plik, sheet = 1)
# Podstawowa kontrola
dim(dane_raw)
## [1] 123 49
names(dane_raw)
## [1] "Sygnatura czasowa"
## [2] "Wyrażam świadomą zgodę na udział w badaniu naukowym."
## [3] "Płeć"
## [4] "Wiek w latach"
## [5] "Wykształcenie"
## [6] "Miejsce zamieszkania"
## [7] "Stan cywilny"
## [8] "1. Dobrze się czuję z ludźmi wokół mnie."
## [9] "2. Brakuje mi towarzystwa."
## [10] "3. Nie mam nikogo, do kogo mógłbym/mogłabym się zwrócić."
## [11] "4. Nie czuję się samotny."
## [12] "5. Czuję się częścią grupy przyjaciół."
## [13] "6. Mam wiele wspólnego z ludźmi wokół mnie."
## [14] "7. Już nie jestem dla nikogo bliską osobą."
## [15] "8. Moje zainteresowania nie są podzielane przez ludzi wokół mnie."
## [16] "9. Jestem osobą towarzyską."
## [17] "10. Są ludzie, którzy są mi bliscy."
## [18] "11. Czuję się opuszczony/a."
## [19] "12. Moje relacje społeczne są powierzchowne."
## [20] "13. Nikt tak naprawdę nie zna mnie dobrze."
## [21] "14. Czuję się odizolowany/a od innych."
## [22] "15. Potrafię znaleźć towarzystwo, kiedy tylko chcę."
## [23] "16. Są ludzie, którzy naprawdę mnie rozumieją."
## [24] "17. Jestem nieszczęśliwy/a z powodu bycia zamkniętym w sobie."
## [25] "18. Ludzie są wokół mnie, ale nie ze mną."
## [26] "19. Są ludzie, z którymi mogę porozmawiać."
## [27] "20. Są ludzie, do których mogę się zwrócić."
## [28] "1. Uważam, że jestem osobą wartościową przynajmniej w takim samym stopniu, co inni."
## [29] "2. Uważam, że posiadam wiele pozytywnych cech."
## [30] "3. Ogólnie biorąc jestem skłonny(a) sądzić, że nie wiedzie mi się."
## [31] "4. Potrafię robić różne rzeczy tak dobrze, jak większość innych ludzi."
## [32] "5. Uważam, że nie mam wielu powodów, aby być z siebie dumn(ą)ym."
## [33] "6. Lubię siebie."
## [34] "7. Ogólnie rzecz biorąc, jestem z siebie zadowolon(a)y."
## [35] "8. Chciał(a)bym mieć więcej szacunku dla samego siebie."
## [36] "9. Czasami czuję się bezużyteczn(a)y."
## [37] "10. Niekiedy uważam, że jestem do niczego."
## [38] "1. Niewielkie zainteresowanie lub odczuwanie przyjemności\nz wykonywania czynności."
## [39] "2. Uczucie smutku, przygnębienia lub beznadziejności."
## [40] "3. Kłopoty z zaśnięciem lub przerywany sen, albo zbyt długi sen."
## [41] "4. Uczucie zmęczenia lub brak energii."
## [42] "5. Brak apetytu lub przejadanie się."
## [43] "6. Poczucie niezadowolenia z siebie - lub uczucie, że jest\nsię do niczego, albo że zawiódł/zawiodła Pan/Pani siebie lub rodzinę."
## [44] "7. Problemy ze skupieniem się na przykład przy czytaniu\ngazety lub oglądaniu telewizji."
## [45] "8. Poruszanie się lub mówienie tak wolno, że inni mogliby to zauważyć? Albo wręcz przeciwnie - niemożność\nusiedzenia w miejscu lub podenerwowanie powodujące\nruchliwość znacznie większą niż zwykle"
## [46] "9. Myśli, że lepiej byłoby umrzeć, albo chęć zrobienia sobie\njakiejś krzywdy."
## [47] "10. Jeżeli zaznaczył/-a Pan/Pani którekolwiek z problemów, jak bardzo utrudniły one Panu/Pani wykonywanie pracy, zajmowanie się domem lub relacje z innymi ludźmi?"
## [48] "Komentarze"
## [49] "W tym miejscu można wpisać uwagi/sugestie"
# KROK 3: Uporządkowanie nazw kolumn -------------------------------------
dane <- dane_raw
# Pierwsze 7 kolumn to metryczka
names(dane)[1:7] <- c(
"timestamp",
"zgoda",
"plec",
"wiek",
"wyksztalcenie",
"miejsce_zamieszkania",
"stan_cywilny"
)
# Kolejne kolumny to pozycje skal
names(dane)[8:27] <- paste0("UCLA", 1:20)
names(dane)[28:37] <- paste0("SES", 1:10)
names(dane)[38:46] <- paste0("PHQ", 1:9)
# Kolumna 47 to dodatkowe pytanie PHQ o utrudnienie funkcjonowania
names(dane)[47] <- "PHQ_utrudnienie"
# Kolumny 48-49 to komentarze
names(dane)[48:49] <- c("komentarze", "uwagi")
# Kontrola nazw
names(dane)
## [1] "timestamp" "zgoda" "plec"
## [4] "wiek" "wyksztalcenie" "miejsce_zamieszkania"
## [7] "stan_cywilny" "UCLA1" "UCLA2"
## [10] "UCLA3" "UCLA4" "UCLA5"
## [13] "UCLA6" "UCLA7" "UCLA8"
## [16] "UCLA9" "UCLA10" "UCLA11"
## [19] "UCLA12" "UCLA13" "UCLA14"
## [22] "UCLA15" "UCLA16" "UCLA17"
## [25] "UCLA18" "UCLA19" "UCLA20"
## [28] "SES1" "SES2" "SES3"
## [31] "SES4" "SES5" "SES6"
## [34] "SES7" "SES8" "SES9"
## [37] "SES10" "PHQ1" "PHQ2"
## [40] "PHQ3" "PHQ4" "PHQ5"
## [43] "PHQ6" "PHQ7" "PHQ8"
## [46] "PHQ9" "PHQ_utrudnienie" "komentarze"
## [49] "uwagi"
# KROK 4: Kontrola danych -------------------------------------------------
# Liczba odpowiedzi
nrow(dane)
## [1] 123
# Zgoda na udział
table(dane$zgoda, useNA = "ifany")
##
## Tak
## 123
# Wiek
summary(dane$wiek)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 18.00 24.00 27.00 26.82 29.50 35.00
# Zakres wieku
range(dane$wiek, na.rm = TRUE)
## [1] 18 35
# Osoby spoza zakresu 18-35 lat
dane %>%
filter(wiek < 18 | wiek > 35)
## # A tibble: 0 × 49
## # ℹ 49 variables: timestamp <dttm>, zgoda <chr>, plec <chr>, wiek <dbl>,
## # wyksztalcenie <chr>, miejsce_zamieszkania <chr>, stan_cywilny <chr>,
## # UCLA1 <chr>, UCLA2 <chr>, UCLA3 <chr>, UCLA4 <chr>, UCLA5 <chr>,
## # UCLA6 <chr>, UCLA7 <chr>, UCLA8 <chr>, UCLA9 <chr>, UCLA10 <chr>,
## # UCLA11 <chr>, UCLA12 <chr>, UCLA13 <chr>, UCLA14 <chr>, UCLA15 <chr>,
## # UCLA16 <chr>, UCLA17 <chr>, UCLA18 <chr>, UCLA19 <chr>, UCLA20 <chr>,
## # SES1 <chr>, SES2 <chr>, SES3 <chr>, SES4 <chr>, SES5 <chr>, SES6 <chr>, …
# Braki danych w głównych kolumnach
colSums(is.na(dane[, 1:47]))
## timestamp zgoda plec
## 0 0 0
## wiek wyksztalcenie miejsce_zamieszkania
## 0 0 0
## stan_cywilny UCLA1 UCLA2
## 0 0 0
## UCLA3 UCLA4 UCLA5
## 0 0 0
## UCLA6 UCLA7 UCLA8
## 0 0 0
## UCLA9 UCLA10 UCLA11
## 0 0 0
## UCLA12 UCLA13 UCLA14
## 0 0 0
## UCLA15 UCLA16 UCLA17
## 0 0 0
## UCLA18 UCLA19 UCLA20
## 0 0 0
## SES1 SES2 SES3
## 0 0 0
## SES4 SES5 SES6
## 0 0 0
## SES7 SES8 SES9
## 0 0 0
## SES10 PHQ1 PHQ2
## 0 0 0
## PHQ3 PHQ4 PHQ5
## 0 0 0
## PHQ6 PHQ7 PHQ8
## 0 0 0
## PHQ9 PHQ_utrudnienie
## 0 0
# KROK 5: Duplikaty -------------------------------------------------------
# Tworzymy klucz duplikatu bez kolumny timestamp
dane_z_kluczem <- dane %>%
mutate(
klucz_duplikatu = apply(
select(., -timestamp),
1,
paste,
collapse = " | "
)
)
# Sprawdzenie liczby duplikatów
sum(duplicated(dane_z_kluczem$klucz_duplikatu))
## [1] 1
# Podgląd duplikatów
dane_z_kluczem %>%
filter(duplicated(klucz_duplikatu) | duplicated(klucz_duplikatu, fromLast = TRUE)) %>%
select(timestamp, plec, wiek, wyksztalcenie, miejsce_zamieszkania, stan_cywilny)
## # A tibble: 2 × 6
## timestamp plec wiek wyksztalcenie miejsce_zamieszkania
## <dttm> <chr> <dbl> <chr> <chr>
## 1 2026-04-25 12:29:19 Kobieta 35 Podstawowe Wieś/obszar wiejski
## 2 2026-04-25 12:29:22 Kobieta 35 Podstawowe Wieś/obszar wiejski
## # ℹ 1 more variable: stan_cywilny <chr>
# Usunięcie duplikatu - zostawiamy pierwszą odpowiedź
dane_clean <- dane_z_kluczem %>%
filter(!duplicated(klucz_duplikatu)) %>%
select(-klucz_duplikatu)
# Kontrola liczby osób po usunięciu duplikatu
nrow(dane_clean)
## [1] 122
# KROK 6: Charakterystyka próby ------------------------------------------
# Płeć
table(dane_clean$plec)
##
## Inna Kobieta Mężczyzna
## 2 69 51
round(prop.table(table(dane_clean$plec)) * 100, 2)
##
## Inna Kobieta Mężczyzna
## 1.64 56.56 41.80
# Wykształcenie
table(dane_clean$wyksztalcenie)
##
## Podstawowe Średnie Wyższe Zawodowe
## 18 50 30 24
round(prop.table(table(dane_clean$wyksztalcenie)) * 100, 2)
##
## Podstawowe Średnie Wyższe Zawodowe
## 14.75 40.98 24.59 19.67
# Miejsce zamieszkania
table(dane_clean$miejsce_zamieszkania)
##
## Miasto od 10 do 100 tys. mieszkańców Miasto od 100 do 500 tys. mieszkańców
## 37 31
## Miasto powyżej 500 tys. mieszkańców Wieś/obszar wiejski
## 27 27
round(prop.table(table(dane_clean$miejsce_zamieszkania)) * 100, 2)
##
## Miasto od 10 do 100 tys. mieszkańców Miasto od 100 do 500 tys. mieszkańców
## 30.33 25.41
## Miasto powyżej 500 tys. mieszkańców Wieś/obszar wiejski
## 22.13 22.13
# Stan cywilny
table(dane_clean$stan_cywilny)
##
## W związku małżeńskim W związku nieformalnym Wolny (singiel/singielka)
## 33 47 42
round(prop.table(table(dane_clean$stan_cywilny)) * 100, 2)
##
## W związku małżeńskim W związku nieformalnym Wolny (singiel/singielka)
## 27.05 38.52 34.43
# Wiek
summary(dane_clean$wiek)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 18.00 24.00 27.00 26.75 29.00 35.00
mean(dane_clean$wiek)
## [1] 26.7541
sd(dane_clean$wiek)
## [1] 4.388775
# KROK 7: Funkcja kodowania odpowiedzi -----------------------------------
wyciagnij_liczbe <- function(x) {
as.integer(str_extract(as.character(x), "^\\d+"))
}
# KROK 8: Kodowanie UCLA --------------------------------------------------
ucla_items <- paste0("UCLA", 1:20)
ucla_reverse <- c(
"UCLA1", "UCLA4", "UCLA5", "UCLA6", "UCLA9",
"UCLA10", "UCLA15", "UCLA16", "UCLA19", "UCLA20"
)
# Zamiana odpowiedzi tekstowych na liczby
dane_scores <- dane_clean %>%
mutate(
across(all_of(ucla_items), wyciagnij_liczbe)
)
# Odwracanie wybranych pozycji UCLA
dane_scores <- dane_scores %>%
mutate(
across(all_of(ucla_reverse), ~ 5 - .x)
)
# Wynik ogólny UCLA
dane_scores <- dane_scores %>%
mutate(
UCLA_wynik = rowSums(select(., all_of(ucla_items)), na.rm = FALSE)
)
# Kontrola zakresu
summary(dane_scores$UCLA_wynik)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 25.00 46.00 50.00 48.83 53.75 72.00
range(dane_scores$UCLA_wynik, na.rm = TRUE)
## [1] 25 72
# KROK 9: Kodowanie SES ---------------------------------------------------
ses_items <- paste0("SES", 1:10)
ses_reverse <- c(
"SES1", "SES2", "SES4", "SES6", "SES7"
)
# Zamiana odpowiedzi tekstowych na liczby
dane_scores <- dane_scores %>%
mutate(
across(all_of(ses_items), wyciagnij_liczbe)
)
# Odwracanie pozycji pozytywnych SES
dane_scores <- dane_scores %>%
mutate(
across(all_of(ses_reverse), ~ 5 - .x)
)
# Wynik ogólny SES
dane_scores <- dane_scores %>%
mutate(
SES_wynik = rowSums(select(., all_of(ses_items)), na.rm = FALSE)
)
# Kontrola zakresu
summary(dane_scores$SES_wynik)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 12.00 24.00 26.00 26.02 28.00 40.00
range(dane_scores$SES_wynik, na.rm = TRUE)
## [1] 12 40
# KROK 10: Kodowanie PHQ-9 ------------------------------------------------
phq_items <- paste0("PHQ", 1:9)
# Zamiana odpowiedzi tekstowych na liczby
dane_scores <- dane_scores %>%
mutate(
across(all_of(phq_items), wyciagnij_liczbe)
)
# Wynik ogólny PHQ-9
dane_scores <- dane_scores %>%
mutate(
PHQ9_wynik = rowSums(select(., all_of(phq_items)), na.rm = FALSE)
)
# Kontrola zakresu
summary(dane_scores$PHQ9_wynik)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 11.00 13.00 12.82 15.00 25.00
range(dane_scores$PHQ9_wynik, na.rm = TRUE)
## [1] 0 25
# KROK 11: Statystyki opisowe wyników skal -------------------------------
statystyki_skal <- dane_scores %>%
summarise(
N = n(),
UCLA_min = min(UCLA_wynik, na.rm = TRUE),
UCLA_max = max(UCLA_wynik, na.rm = TRUE),
UCLA_M = mean(UCLA_wynik, na.rm = TRUE),
UCLA_SD = sd(UCLA_wynik, na.rm = TRUE),
SES_min = min(SES_wynik, na.rm = TRUE),
SES_max = max(SES_wynik, na.rm = TRUE),
SES_M = mean(SES_wynik, na.rm = TRUE),
SES_SD = sd(SES_wynik, na.rm = TRUE),
PHQ9_min = min(PHQ9_wynik, na.rm = TRUE),
PHQ9_max = max(PHQ9_wynik, na.rm = TRUE),
PHQ9_M = mean(PHQ9_wynik, na.rm = TRUE),
PHQ9_SD = sd(PHQ9_wynik, na.rm = TRUE)
)
statystyki_skal
## # A tibble: 1 × 13
## N UCLA_min UCLA_max UCLA_M UCLA_SD SES_min SES_max SES_M SES_SD PHQ9_min
## <int> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 122 25 72 48.8 8.62 12 40 26.0 4.68 0
## # ℹ 3 more variables: PHQ9_max <dbl>, PHQ9_M <dbl>, PHQ9_SD <dbl>
# KROK 12: Zapis danych ---------------------------------------------------
# Utworzenie folderu outputs, jeśli nie istnieje
if (!dir.exists("outputs")) {
dir.create("outputs")
}
# Zapis pełnych danych po kodowaniu
write.csv(
dane_scores,
file = "outputs/dane_po_kodowaniu.csv",
row.names = FALSE
)
# Zapis jako plik RDS - wygodny do dalszej pracy w R
saveRDS(
dane_scores,
file = "outputs/dane_po_kodowaniu.rds"
)
# Zapis podstawowych statystyk
write.csv(
statystyki_skal,
file = "outputs/statystyki_skal_podstawowe.csv",
row.names = FALSE
)
# KROK 13: Tabele pomocnicze ---------------------------------------------
tabela_plec <- dane_scores %>%
count(plec) %>%
mutate(
procent = round(n / sum(n) * 100, 2)
)
tabela_wyksztalcenie <- dane_scores %>%
count(wyksztalcenie) %>%
mutate(
procent = round(n / sum(n) * 100, 2)
)
tabela_miejsce <- dane_scores %>%
count(miejsce_zamieszkania) %>%
mutate(
procent = round(n / sum(n) * 100, 2)
)
tabela_stan_cywilny <- dane_scores %>%
count(stan_cywilny) %>%
mutate(
procent = round(n / sum(n) * 100, 2)
)
# Zapis do jednego pliku Excel z kilkoma arkuszami
write_xlsx(
list(
"plec" = tabela_plec,
"wyksztalcenie" = tabela_wyksztalcenie,
"miejsce_zamieszkania" = tabela_miejsce,
"stan_cywilny" = tabela_stan_cywilny,
"statystyki_skal" = statystyki_skal
),
path = "outputs/tabele_krok_1_2.xlsx"
)
# KROK 14: Wczytanie danych po kodowaniu ---------------------------------
dane_scores <- readRDS("outputs/dane_po_kodowaniu.rds")
# Kontrola
nrow(dane_scores)
## [1] 122
names(dane_scores)
## [1] "timestamp" "zgoda" "plec"
## [4] "wiek" "wyksztalcenie" "miejsce_zamieszkania"
## [7] "stan_cywilny" "UCLA1" "UCLA2"
## [10] "UCLA3" "UCLA4" "UCLA5"
## [13] "UCLA6" "UCLA7" "UCLA8"
## [16] "UCLA9" "UCLA10" "UCLA11"
## [19] "UCLA12" "UCLA13" "UCLA14"
## [22] "UCLA15" "UCLA16" "UCLA17"
## [25] "UCLA18" "UCLA19" "UCLA20"
## [28] "SES1" "SES2" "SES3"
## [31] "SES4" "SES5" "SES6"
## [34] "SES7" "SES8" "SES9"
## [37] "SES10" "PHQ1" "PHQ2"
## [40] "PHQ3" "PHQ4" "PHQ5"
## [43] "PHQ6" "PHQ7" "PHQ8"
## [46] "PHQ9" "PHQ_utrudnienie" "komentarze"
## [49] "uwagi" "UCLA_wynik" "SES_wynik"
## [52] "PHQ9_wynik"
# KROK 15: Rzetelność skal ------------------------------------------------
# Listy pozycji
ucla_items <- paste0("UCLA", 1:20)
ses_items <- paste0("SES", 1:10)
phq_items <- paste0("PHQ", 1:9)
# Alfa Cronbacha dla UCLA
alpha_ucla <- psych::alpha(
dane_scores[, ucla_items],
check.keys = FALSE
)
# Alfa Cronbacha dla SES
alpha_ses <- psych::alpha(
dane_scores[, ses_items],
check.keys = FALSE
)
# Alfa Cronbacha dla PHQ-9
alpha_phq <- psych::alpha(
dane_scores[, phq_items],
check.keys = FALSE
)
# Podgląd wyników
alpha_ucla$total
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd
## 0.8246629 0.8246351 0.8552314 0.190362 4.702396 0.02283698 2.441393 0.431157
## median_r
## 0.1979458
alpha_ses$total
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd
## 0.7397391 0.7403599 0.7495592 0.2218798 2.851485 0.03510122 2.602459 0.468418
## median_r
## 0.2055471
alpha_phq$total
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd
## 0.7315132 0.7317413 0.7429231 0.232589 2.727745 0.03639686 1.424408 0.5314353
## median_r
## 0.234392
# KROK 16: Tabela rzetelności --------------------------------------------
tabela_rzetelnosc <- data.frame(
Narzedzie = c(
"Skala Poczucia Osamotnienia UCLA",
"Skala Samooceny Rosenberga SES",
"Kwestionariusz PHQ-9"
),
Liczba_pozycji = c(20, 10, 9),
Alfa_Cronbacha = c(
alpha_ucla$total$raw_alpha,
alpha_ses$total$raw_alpha,
alpha_phq$total$raw_alpha
)
)
tabela_rzetelnosc <- tabela_rzetelnosc %>%
mutate(
Alfa_Cronbacha = round(Alfa_Cronbacha, 3)
)
tabela_rzetelnosc
## Narzedzie Liczba_pozycji Alfa_Cronbacha
## 1 Skala Poczucia Osamotnienia UCLA 20 0.825
## 2 Skala Samooceny Rosenberga SES 10 0.740
## 3 Kwestionariusz PHQ-9 9 0.732
# KROK 17: Statystyki opisowe zmiennych ilościowych -----------------------
zmienne_glowne <- dane_scores %>%
select(
wiek,
UCLA_wynik,
SES_wynik,
PHQ9_wynik
)
statystyki_opisowe <- zmienne_glowne %>%
psych::describe() %>%
as.data.frame()
statystyki_opisowe <- statystyki_opisowe %>%
select(
n,
mean,
sd,
median,
min,
max,
skew,
kurtosis
)
statystyki_opisowe <- round(statystyki_opisowe, 2)
statystyki_opisowe
## n mean sd median min max skew kurtosis
## wiek 122 26.75 4.39 27 18 35 0.07 -0.71
## UCLA_wynik 122 48.83 8.62 50 25 72 -0.57 0.80
## SES_wynik 122 26.02 4.68 26 12 40 0.04 1.34
## PHQ9_wynik 122 12.82 4.78 13 0 25 -0.28 0.90
# KROK 18: Tabela statystyk opisowych do pracy ----------------------------
tabela_statystyki_opisowe <- statystyki_opisowe
tabela_statystyki_opisowe$Zmienna <- rownames(tabela_statystyki_opisowe)
tabela_statystyki_opisowe <- tabela_statystyki_opisowe %>%
select(
Zmienna,
n,
mean,
sd,
median,
min,
max,
skew,
kurtosis
) %>%
mutate(
Zmienna = recode(
Zmienna,
"wiek" = "Wiek",
"UCLA_wynik" = "Poczucie osamotnienia",
"SES_wynik" = "Samoocena",
"PHQ9_wynik" = "Symptomy depresyjne"
)
)
tabela_statystyki_opisowe
## Zmienna n mean sd median min max skew kurtosis
## wiek Wiek 122 26.75 4.39 27 18 35 0.07 -0.71
## UCLA_wynik Poczucie osamotnienia 122 48.83 8.62 50 25 72 -0.57 0.80
## SES_wynik Samoocena 122 26.02 4.68 26 12 40 0.04 1.34
## PHQ9_wynik Symptomy depresyjne 122 12.82 4.78 13 0 25 -0.28 0.90
# KROK 19: Charakterystyka płci ------------------------------------------
tabela_plec <- dane_scores %>%
count(plec) %>%
mutate(
procent = round(n / sum(n) * 100, 2)
)
tabela_plec
## # A tibble: 3 × 3
## plec n procent
## <chr> <int> <dbl>
## 1 Inna 2 1.64
## 2 Kobieta 69 56.6
## 3 Mężczyzna 51 41.8
# KROK 20: Charakterystyka wieku -----------------------------------------
wiek_opis <- dane_scores %>%
summarise(
N = n(),
M = round(mean(wiek, na.rm = TRUE), 2),
SD = round(sd(wiek, na.rm = TRUE), 2),
Me = round(median(wiek, na.rm = TRUE), 2),
Min = min(wiek, na.rm = TRUE),
Max = max(wiek, na.rm = TRUE)
)
wiek_opis
## # A tibble: 1 × 6
## N M SD Me Min Max
## <int> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 122 26.8 4.39 27 18 35
# KROK 21: Zapis tabel z rzetelnością i statystykami ----------------------
write_xlsx(
list(
"rzetelnosc_skal" = tabela_rzetelnosc,
"statystyki_opisowe" = tabela_statystyki_opisowe,
"plec" = tabela_plec,
"wiek" = wiek_opis
),
path = "outputs/tabele_rzetelnosc_i_opisowe.xlsx"
)
# KROK 22: Kontrola końcowa ----------------------------------------------
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
tabela_rzetelnosc
## Narzedzie Liczba_pozycji Alfa_Cronbacha
## 1 Skala Poczucia Osamotnienia UCLA 20 0.825
## 2 Skala Samooceny Rosenberga SES 10 0.740
## 3 Kwestionariusz PHQ-9 9 0.732
tabela_statystyki_opisowe
## Zmienna n mean sd median min max skew kurtosis
## wiek Wiek 122 26.75 4.39 27 18 35 0.07 -0.71
## UCLA_wynik Poczucie osamotnienia 122 48.83 8.62 50 25 72 -0.57 0.80
## SES_wynik Samoocena 122 26.02 4.68 26 12 40 0.04 1.34
## PHQ9_wynik Symptomy depresyjne 122 12.82 4.78 13 0 25 -0.28 0.90
tabela_plec
## # A tibble: 3 × 3
## plec n procent
## <chr> <int> <dbl>
## 1 Inna 2 1.64
## 2 Kobieta 69 56.6
## 3 Mężczyzna 51 41.8
wiek_opis
## # A tibble: 1 × 6
## N M SD Me Min Max
## <int> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 122 26.8 4.39 27 18 35
# KROK 23: Przygotowanie zmiennych do analizy rozkładów -------------------
library(dplyr)
library(tidyr)
library(ggplot2)
library(writexl)
# Jeżeli dane_scores nie jest aktywne w środowisku, wczytujemy je ponownie
dane_scores <- readRDS("outputs/dane_po_kodowaniu.rds")
# Wybór głównych zmiennych ilościowych
zmienne_glowne <- dane_scores %>%
select(
wiek,
UCLA_wynik,
SES_wynik,
PHQ9_wynik
)
# Nazwy do tabel i wykresów
nazwy_zmiennych <- c(
"wiek" = "Wiek",
"UCLA_wynik" = "Poczucie osamotnienia",
"SES_wynik" = "Samoocena",
"PHQ9_wynik" = "Symptomy depresyjne"
)
# Kontrola
head(zmienne_glowne)
## # A tibble: 6 × 4
## wiek UCLA_wynik SES_wynik PHQ9_wynik
## <dbl> <dbl> <dbl> <dbl>
## 1 25 43 29 13
## 2 22 62 15 17
## 3 18 37 28 14
## 4 20 59 12 25
## 5 29 45 26 15
## 6 23 62 29 19
# KROK 24: Test normalności rozkładów -------------------------------------
tabela_normalnosc <- zmienne_glowne %>%
pivot_longer(
cols = everything(),
names_to = "Zmienna",
values_to = "Wartosc"
) %>%
group_by(Zmienna) %>%
summarise(
Shapiro_W = round(shapiro.test(Wartosc)$statistic, 3),
Shapiro_p = round(shapiro.test(Wartosc)$p.value, 4),
.groups = "drop"
) %>%
mutate(
Zmienna = nazwy_zmiennych[Zmienna]
)
tabela_normalnosc
## # A tibble: 4 × 3
## Zmienna Shapiro_W Shapiro_p
## <chr> <dbl> <dbl>
## 1 Symptomy depresyjne 0.955 0.0005
## 2 Samoocena 0.96 0.001
## 3 Poczucie osamotnienia 0.95 0.0002
## 4 Wiek 0.975 0.023
# KROK 25: Utworzenie folderu na wykresy ---------------------------------
dir.create("outputs/wykresy", showWarnings = FALSE)
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
# KROK 26: Dane do wykresów -----------------------------------------------
wykres_dane <- zmienne_glowne %>%
pivot_longer(
cols = everything(),
names_to = "Zmienna",
values_to = "Wartosc"
) %>%
mutate(
Zmienna = nazwy_zmiennych[Zmienna]
)
head(wykres_dane)
## # A tibble: 6 × 2
## Zmienna Wartosc
## <chr> <dbl>
## 1 Wiek 25
## 2 Poczucie osamotnienia 43
## 3 Samoocena 29
## 4 Symptomy depresyjne 13
## 5 Wiek 22
## 6 Poczucie osamotnienia 62
# KROK 27: Zbiorczy histogram rozkładów -----------------------------------
p_hist_all <- ggplot(wykres_dane, aes(x = Wartosc)) +
geom_histogram(
aes(y = after_stat(density)),
bins = 12,
fill = "grey85",
color = "grey30"
) +
geom_density(linewidth = 1) +
facet_wrap(~ Zmienna, scales = "free", ncol = 2) +
labs(
title = "Rozkłady wyników analizowanych zmiennych",
x = "Wynik",
y = "Gęstość"
) +
theme_minimal(base_size = 12)
p_hist_all
ggsave(
filename = "outputs/wykresy/histogramy_zmiennych_glownych.png",
plot = p_hist_all,
width = 9,
height = 7,
dpi = 300
)
# KROK 28: Wykresy pudełkowe ---------------------------------------------
p_box_all <- ggplot(wykres_dane, aes(x = "", y = Wartosc)) +
geom_boxplot(
fill = "grey85",
color = "grey30"
) +
facet_wrap(~ Zmienna, scales = "free_y", ncol = 2) +
labs(
title = "Wykresy pudełkowe analizowanych zmiennych",
x = "",
y = "Wynik"
) +
theme_minimal(base_size = 12)
p_box_all
ggsave(
filename = "outputs/wykresy/boxploty_zmiennych_glownych.png",
plot = p_box_all,
width = 9,
height = 7,
dpi = 300
)
# KROK 29: Wykresy Q-Q ----------------------------------------------------
p_qq_all <- ggplot(wykres_dane, aes(sample = Wartosc)) +
stat_qq(color = "grey30") +
stat_qq_line() +
facet_wrap(~ Zmienna, scales = "free", ncol = 2) +
labs(
title = "Wykresy Q-Q analizowanych zmiennych",
x = "Kwantyle teoretyczne",
y = "Kwantyle empiryczne"
) +
theme_minimal(base_size = 12)
p_qq_all
ggsave(
filename = "outputs/wykresy/qqploty_zmiennych_glownych.png",
plot = p_qq_all,
width = 9,
height = 7,
dpi = 300
)
# KROK 30: Zapis tabel do pliku Excel -------------------------------------
write_xlsx(
list(
"rzetelnosc_skal" = tabela_rzetelnosc,
"statystyki_opisowe" = tabela_statystyki_opisowe,
"normalnosc" = tabela_normalnosc,
"plec" = tabela_plec,
"wiek" = wiek_opis
),
path = "outputs/tabele_rzetelnosc_opisowe_normalnosc.xlsx"
)
# KROK 31: Kontrola końcowa po wykresach ----------------------------------
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
list.files("outputs/wykresy")
## [1] "boxploty_zmiennych_glownych.png"
## [2] "H1_osamotnienie_samoocena.png"
## [3] "H2_osamotnienie_depresja.png"
## [4] "H3_samoocena_depresja.png"
## [5] "H5_poczucie_osamotnienia_plec.png"
## [6] "H5_roznice_plec_boxplot.png"
## [7] "H5_samoocena_plec.png"
## [8] "H5_symptomy_depresyjne_plec.png"
## [9] "H6_model_mediacji.png"
## [10] "H7_moderacja_samoocena.png"
## [11] "histogramy_zmiennych_glownych.png"
## [12] "qqploty_zmiennych_glownych.png"
## [13] "regresja_obserwowane_przewidywane_PHQ9.png"
## [14] "regresja_qq_reszty.png"
## [15] "regresja_reszty_wartosci_przewidywane.png"
tabela_rzetelnosc
## Narzedzie Liczba_pozycji Alfa_Cronbacha
## 1 Skala Poczucia Osamotnienia UCLA 20 0.825
## 2 Skala Samooceny Rosenberga SES 10 0.740
## 3 Kwestionariusz PHQ-9 9 0.732
tabela_statystyki_opisowe
## Zmienna n mean sd median min max skew kurtosis
## wiek Wiek 122 26.75 4.39 27 18 35 0.07 -0.71
## UCLA_wynik Poczucie osamotnienia 122 48.83 8.62 50 25 72 -0.57 0.80
## SES_wynik Samoocena 122 26.02 4.68 26 12 40 0.04 1.34
## PHQ9_wynik Symptomy depresyjne 122 12.82 4.78 13 0 25 -0.28 0.90
tabela_normalnosc
## # A tibble: 4 × 3
## Zmienna Shapiro_W Shapiro_p
## <chr> <dbl> <dbl>
## 1 Symptomy depresyjne 0.955 0.0005
## 2 Samoocena 0.96 0.001
## 3 Poczucie osamotnienia 0.95 0.0002
## 4 Wiek 0.975 0.023
tabela_plec
## # A tibble: 3 × 3
## plec n procent
## <chr> <int> <dbl>
## 1 Inna 2 1.64
## 2 Kobieta 69 56.6
## 3 Mężczyzna 51 41.8
wiek_opis
## # A tibble: 1 × 6
## N M SD Me Min Max
## <int> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 122 26.8 4.39 27 18 35
# KROK 32: Przygotowanie danych do analizy korelacji ----------------------
library(dplyr)
library(writexl)
# Wczytanie danych po kodowaniu
dane_scores <- readRDS("outputs/dane_po_kodowaniu.rds")
# Wybór zmiennych do korelacji
dane_korelacje <- dane_scores %>%
select(
wiek,
UCLA_wynik,
SES_wynik,
PHQ9_wynik
)
# Kontrola danych
head(dane_korelacje)
## # A tibble: 6 × 4
## wiek UCLA_wynik SES_wynik PHQ9_wynik
## <dbl> <dbl> <dbl> <dbl>
## 1 25 43 29 13
## 2 22 62 15 17
## 3 18 37 28 14
## 4 20 59 12 25
## 5 29 45 26 15
## 6 23 62 29 19
summary(dane_korelacje)
## wiek UCLA_wynik SES_wynik PHQ9_wynik
## Min. :18.00 Min. :25.00 Min. :12.00 Min. : 0.00
## 1st Qu.:24.00 1st Qu.:46.00 1st Qu.:24.00 1st Qu.:11.00
## Median :27.00 Median :50.00 Median :26.00 Median :13.00
## Mean :26.75 Mean :48.83 Mean :26.02 Mean :12.82
## 3rd Qu.:29.00 3rd Qu.:53.75 3rd Qu.:28.00 3rd Qu.:15.00
## Max. :35.00 Max. :72.00 Max. :40.00 Max. :25.00
# KROK 33: Funkcja do korelacji Spearmana ---------------------------------
korelacja_spearman <- function(data, x, y, opis) {
test <- cor.test(
data[[x]],
data[[y]],
method = "spearman",
exact = FALSE,
use = "complete.obs"
)
data.frame(
Hipoteza = opis,
Zmienna_1 = x,
Zmienna_2 = y,
rho_Spearmana = round(as.numeric(test$estimate), 3),
p = round(test$p.value, 4)
)
}
# KROK 34: Korelacje Spearmana dla hipotez H1-H4 --------------------------
tabela_korelacje <- bind_rows(
korelacja_spearman(
dane_korelacje,
"UCLA_wynik",
"SES_wynik",
"H1: poczucie osamotnienia a samoocena"
),
korelacja_spearman(
dane_korelacje,
"UCLA_wynik",
"PHQ9_wynik",
"H2: poczucie osamotnienia a symptomy depresyjne"
),
korelacja_spearman(
dane_korelacje,
"SES_wynik",
"PHQ9_wynik",
"H3: samoocena a symptomy depresyjne"
),
korelacja_spearman(
dane_korelacje,
"wiek",
"UCLA_wynik",
"H4: wiek a poczucie osamotnienia"
),
korelacja_spearman(
dane_korelacje,
"wiek",
"SES_wynik",
"H4: wiek a samoocena"
),
korelacja_spearman(
dane_korelacje,
"wiek",
"PHQ9_wynik",
"H4: wiek a symptomy depresyjne"
)
)
tabela_korelacje
## Hipoteza Zmienna_1 Zmienna_2
## 1 H1: poczucie osamotnienia a samoocena UCLA_wynik SES_wynik
## 2 H2: poczucie osamotnienia a symptomy depresyjne UCLA_wynik PHQ9_wynik
## 3 H3: samoocena a symptomy depresyjne SES_wynik PHQ9_wynik
## 4 H4: wiek a poczucie osamotnienia wiek UCLA_wynik
## 5 H4: wiek a samoocena wiek SES_wynik
## 6 H4: wiek a symptomy depresyjne wiek PHQ9_wynik
## rho_Spearmana p
## 1 -0.574 0.0000
## 2 0.476 0.0000
## 3 -0.414 0.0000
## 4 -0.283 0.0016
## 5 0.078 0.3951
## 6 -0.176 0.0521
# KROK 35: Interpretacja istotności korelacji -----------------------------
tabela_korelacje <- tabela_korelacje %>%
mutate(
Istotnosc = case_when(
p < 0.001 ~ "p < 0,001",
p < 0.01 ~ "p < 0,01",
p < 0.05 ~ "p < 0,05",
TRUE ~ "nieistotne statystycznie"
),
Kierunek = case_when(
rho_Spearmana > 0 ~ "dodatni",
rho_Spearmana < 0 ~ "ujemny",
TRUE ~ "brak"
),
Sila_zwiazku = case_when(
abs(rho_Spearmana) < 0.10 ~ "bardzo słaby",
abs(rho_Spearmana) < 0.30 ~ "słaby",
abs(rho_Spearmana) < 0.50 ~ "umiarkowany",
abs(rho_Spearmana) < 0.70 ~ "silny",
TRUE ~ "bardzo silny"
)
)
tabela_korelacje
## Hipoteza Zmienna_1 Zmienna_2
## 1 H1: poczucie osamotnienia a samoocena UCLA_wynik SES_wynik
## 2 H2: poczucie osamotnienia a symptomy depresyjne UCLA_wynik PHQ9_wynik
## 3 H3: samoocena a symptomy depresyjne SES_wynik PHQ9_wynik
## 4 H4: wiek a poczucie osamotnienia wiek UCLA_wynik
## 5 H4: wiek a samoocena wiek SES_wynik
## 6 H4: wiek a symptomy depresyjne wiek PHQ9_wynik
## rho_Spearmana p Istotnosc Kierunek Sila_zwiazku
## 1 -0.574 0.0000 p < 0,001 ujemny silny
## 2 0.476 0.0000 p < 0,001 dodatni umiarkowany
## 3 -0.414 0.0000 p < 0,001 ujemny umiarkowany
## 4 -0.283 0.0016 p < 0,01 ujemny słaby
## 5 0.078 0.3951 nieistotne statystycznie dodatni bardzo słaby
## 6 -0.176 0.0521 nieistotne statystycznie ujemny słaby
# KROK 36: Pełna macierz korelacji Spearmana ------------------------------
macierz_korelacji <- cor(
dane_korelacje,
method = "spearman",
use = "complete.obs"
)
macierz_korelacji <- round(macierz_korelacji, 3)
macierz_korelacji
## wiek UCLA_wynik SES_wynik PHQ9_wynik
## wiek 1.000 -0.283 0.078 -0.176
## UCLA_wynik -0.283 1.000 -0.574 0.476
## SES_wynik 0.078 -0.574 1.000 -0.414
## PHQ9_wynik -0.176 0.476 -0.414 1.000
# KROK 37: Zapis wyników korelacji ----------------------------------------
write_xlsx(
list(
"korelacje_H1_H4" = tabela_korelacje,
"macierz_korelacji" = as.data.frame(macierz_korelacji)
),
path = "outputs/tabele_korelacje_H1_H4.xlsx"
)
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
# KROK 38: Wykresy zależności dla głównych hipotez ------------------------
library(ggplot2)
library(tidyr)
dir.create("outputs/wykresy", showWarnings = FALSE)
# H1: poczucie osamotnienia a samoocena
p_H1 <- ggplot(dane_korelacje, aes(x = UCLA_wynik, y = SES_wynik)) +
geom_point(alpha = 0.7) +
geom_smooth(method = "loess", se = TRUE) +
labs(
title = "Związek między poczuciem osamotnienia a samooceną",
x = "Poczucie osamotnienia",
y = "Samoocena"
) +
theme_minimal(base_size = 12)
p_H1
ggsave(
filename = "outputs/wykresy/H1_osamotnienie_samoocena.png",
plot = p_H1,
width = 7,
height = 5,
dpi = 300
)
# H2: poczucie osamotnienia a symptomy depresyjne
p_H2 <- ggplot(dane_korelacje, aes(x = UCLA_wynik, y = PHQ9_wynik)) +
geom_point(alpha = 0.7) +
geom_smooth(method = "loess", se = TRUE) +
labs(
title = "Związek między poczuciem osamotnienia a symptomami depresyjnymi",
x = "Poczucie osamotnienia",
y = "Symptomy depresyjne"
) +
theme_minimal(base_size = 12)
p_H2
ggsave(
filename = "outputs/wykresy/H2_osamotnienie_depresja.png",
plot = p_H2,
width = 7,
height = 5,
dpi = 300
)
# H3: samoocena a symptomy depresyjne
p_H3 <- ggplot(dane_korelacje, aes(x = SES_wynik, y = PHQ9_wynik)) +
geom_point(alpha = 0.7) +
geom_smooth(method = "loess", se = TRUE) +
labs(
title = "Związek między samooceną a symptomami depresyjnymi",
x = "Samoocena",
y = "Symptomy depresyjne"
) +
theme_minimal(base_size = 12)
p_H3
ggsave(
filename = "outputs/wykresy/H3_samoocena_depresja.png",
plot = p_H3,
width = 7,
height = 5,
dpi = 300
)
# KROK 39: Kontrola końcowa po korelacjach --------------------------------
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
list.files("outputs/wykresy")
## [1] "boxploty_zmiennych_glownych.png"
## [2] "H1_osamotnienie_samoocena.png"
## [3] "H2_osamotnienie_depresja.png"
## [4] "H3_samoocena_depresja.png"
## [5] "H5_poczucie_osamotnienia_plec.png"
## [6] "H5_roznice_plec_boxplot.png"
## [7] "H5_samoocena_plec.png"
## [8] "H5_symptomy_depresyjne_plec.png"
## [9] "H6_model_mediacji.png"
## [10] "H7_moderacja_samoocena.png"
## [11] "histogramy_zmiennych_glownych.png"
## [12] "qqploty_zmiennych_glownych.png"
## [13] "regresja_obserwowane_przewidywane_PHQ9.png"
## [14] "regresja_qq_reszty.png"
## [15] "regresja_reszty_wartosci_przewidywane.png"
tabela_korelacje
## Hipoteza Zmienna_1 Zmienna_2
## 1 H1: poczucie osamotnienia a samoocena UCLA_wynik SES_wynik
## 2 H2: poczucie osamotnienia a symptomy depresyjne UCLA_wynik PHQ9_wynik
## 3 H3: samoocena a symptomy depresyjne SES_wynik PHQ9_wynik
## 4 H4: wiek a poczucie osamotnienia wiek UCLA_wynik
## 5 H4: wiek a samoocena wiek SES_wynik
## 6 H4: wiek a symptomy depresyjne wiek PHQ9_wynik
## rho_Spearmana p Istotnosc Kierunek Sila_zwiazku
## 1 -0.574 0.0000 p < 0,001 ujemny silny
## 2 0.476 0.0000 p < 0,001 dodatni umiarkowany
## 3 -0.414 0.0000 p < 0,001 ujemny umiarkowany
## 4 -0.283 0.0016 p < 0,01 ujemny słaby
## 5 0.078 0.3951 nieistotne statystycznie dodatni bardzo słaby
## 6 -0.176 0.0521 nieistotne statystycznie ujemny słaby
macierz_korelacji
## wiek UCLA_wynik SES_wynik PHQ9_wynik
## wiek 1.000 -0.283 0.078 -0.176
## UCLA_wynik -0.283 1.000 -0.574 0.476
## SES_wynik 0.078 -0.574 1.000 -0.414
## PHQ9_wynik -0.176 0.476 -0.414 1.000
# KROK 40: Przygotowanie danych do analizy H5 -----------------------------
library(dplyr)
library(tidyr)
library(ggplot2)
library(writexl)
# Wczytanie danych po kodowaniu
dane_scores <- readRDS("outputs/dane_po_kodowaniu.rds")
# Do analizy H5 zostawiamy tylko kobiety i mężczyzn
dane_H5 <- dane_scores %>%
filter(plec %in% c("Kobieta", "Mężczyzna")) %>%
mutate(
plec = factor(plec, levels = c("Kobieta", "Mężczyzna"))
)
# Kontrola liczebności
table(dane_H5$plec)
##
## Kobieta Mężczyzna
## 69 51
# Kontrola liczby osób
nrow(dane_H5)
## [1] 120
# KROK 41: Statystyki opisowe według płci ---------------------------------
tabela_H5_opis <- dane_H5 %>%
group_by(plec) %>%
summarise(
n = n(),
UCLA_M = round(mean(UCLA_wynik, na.rm = TRUE), 2),
UCLA_SD = round(sd(UCLA_wynik, na.rm = TRUE), 2),
UCLA_Me = median(UCLA_wynik, na.rm = TRUE),
UCLA_Min = min(UCLA_wynik, na.rm = TRUE),
UCLA_Max = max(UCLA_wynik, na.rm = TRUE),
SES_M = round(mean(SES_wynik, na.rm = TRUE), 2),
SES_SD = round(sd(SES_wynik, na.rm = TRUE), 2),
SES_Me = median(SES_wynik, na.rm = TRUE),
SES_Min = min(SES_wynik, na.rm = TRUE),
SES_Max = max(SES_wynik, na.rm = TRUE),
PHQ9_M = round(mean(PHQ9_wynik, na.rm = TRUE), 2),
PHQ9_SD = round(sd(PHQ9_wynik, na.rm = TRUE), 2),
PHQ9_Me = median(PHQ9_wynik, na.rm = TRUE),
PHQ9_Min = min(PHQ9_wynik, na.rm = TRUE),
PHQ9_Max = max(PHQ9_wynik, na.rm = TRUE),
.groups = "drop"
)
tabela_H5_opis
## # A tibble: 2 × 17
## plec n UCLA_M UCLA_SD UCLA_Me UCLA_Min UCLA_Max SES_M SES_SD SES_Me
## <fct> <int> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Kobieta 69 47.5 9.09 49 25 68 26.2 4.68 26
## 2 Mężczyzna 51 50.9 7.04 51 28 72 25.6 4.55 26
## # ℹ 7 more variables: SES_Min <dbl>, SES_Max <dbl>, PHQ9_M <dbl>,
## # PHQ9_SD <dbl>, PHQ9_Me <dbl>, PHQ9_Min <dbl>, PHQ9_Max <dbl>
# KROK 42: Tabela opisowa H5 w układzie do pracy --------------------------
tabela_H5_opis_long <- dane_H5 %>%
select(
plec,
UCLA_wynik,
SES_wynik,
PHQ9_wynik
) %>%
pivot_longer(
cols = c(UCLA_wynik, SES_wynik, PHQ9_wynik),
names_to = "Zmienna",
values_to = "Wynik"
) %>%
mutate(
Zmienna = recode(
Zmienna,
"UCLA_wynik" = "Poczucie osamotnienia",
"SES_wynik" = "Samoocena",
"PHQ9_wynik" = "Symptomy depresyjne"
)
) %>%
group_by(Zmienna, plec) %>%
summarise(
n = n(),
M = round(mean(Wynik, na.rm = TRUE), 2),
SD = round(sd(Wynik, na.rm = TRUE), 2),
Me = median(Wynik, na.rm = TRUE),
Min = min(Wynik, na.rm = TRUE),
Max = max(Wynik, na.rm = TRUE),
.groups = "drop"
)
tabela_H5_opis_long
## # A tibble: 6 × 8
## Zmienna plec n M SD Me Min Max
## <chr> <fct> <int> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Poczucie osamotnienia Kobieta 69 47.5 9.09 49 25 68
## 2 Poczucie osamotnienia Mężczyzna 51 50.9 7.04 51 28 72
## 3 Samoocena Kobieta 69 26.2 4.68 26 14 40
## 4 Samoocena Mężczyzna 51 25.6 4.55 26 12 39
## 5 Symptomy depresyjne Kobieta 69 12.4 5.3 13 0 25
## 6 Symptomy depresyjne Mężczyzna 51 13.6 3.78 14 0 25
# KROK 43: Test U Manna-Whitneya dla H5 -----------------------------------
test_mann_whitney <- function(data, zmienna, opis_zmiennej) {
test <- wilcox.test(
as.formula(paste(zmienna, "~ plec")),
data = data,
exact = FALSE
)
# Liczebności grup
n_k <- sum(data$plec == "Kobieta")
n_m <- sum(data$plec == "Mężczyzna")
# Statystyka W z R; przy dwóch grupach odpowiada statystyce rangowej
W <- as.numeric(test$statistic)
# Efekt rangowo-dwuseryjny.
# Wartość dodatnia oznacza wyższe wyniki w grupie kobiet,
# wartość ujemna oznacza wyższe wyniki w grupie mężczyzn.
r_rb <- (2 * W) / (n_k * n_m) - 1
data.frame(
Zmienna = opis_zmiennej,
W = round(W, 2),
p = round(test$p.value, 4),
r_rangowo_dwuseryjne = round(r_rb, 3)
)
}
tabela_H5_testy <- bind_rows(
test_mann_whitney(
dane_H5,
"UCLA_wynik",
"Poczucie osamotnienia"
),
test_mann_whitney(
dane_H5,
"SES_wynik",
"Samoocena"
),
test_mann_whitney(
dane_H5,
"PHQ9_wynik",
"Symptomy depresyjne"
)
)
tabela_H5_testy
## Zmienna W p r_rangowo_dwuseryjne
## 1 Poczucie osamotnienia 1470.0 0.1243 -0.165
## 2 Samoocena 1783.5 0.9003 0.014
## 3 Symptomy depresyjne 1429.5 0.0788 -0.188
# KROK 44: Interpretacja wyników H5 ---------------------------------------
tabela_H5_testy <- tabela_H5_testy %>%
mutate(
Istotnosc = case_when(
p < 0.001 ~ "p < 0,001",
p < 0.01 ~ "p < 0,01",
p < 0.05 ~ "p < 0,05",
TRUE ~ "nieistotne statystycznie"
),
Sila_efektu = case_when(
abs(r_rangowo_dwuseryjne) < 0.10 ~ "bardzo słaba",
abs(r_rangowo_dwuseryjne) < 0.30 ~ "słaba",
abs(r_rangowo_dwuseryjne) < 0.50 ~ "umiarkowana",
TRUE ~ "silna"
),
Kierunek = case_when(
r_rangowo_dwuseryjne > 0 ~ "wyższe wyniki u kobiet",
r_rangowo_dwuseryjne < 0 ~ "wyższe wyniki u mężczyzn",
TRUE ~ "brak różnicy kierunkowej"
)
)
tabela_H5_testy
## Zmienna W p r_rangowo_dwuseryjne
## 1 Poczucie osamotnienia 1470.0 0.1243 -0.165
## 2 Samoocena 1783.5 0.9003 0.014
## 3 Symptomy depresyjne 1429.5 0.0788 -0.188
## Istotnosc Sila_efektu Kierunek
## 1 nieistotne statystycznie słaba wyższe wyniki u mężczyzn
## 2 nieistotne statystycznie bardzo słaba wyższe wyniki u kobiet
## 3 nieistotne statystycznie słaba wyższe wyniki u mężczyzn
# KROK 45: Wykresy różnic między kobietami i mężczyznami ------------------
dir.create("outputs/wykresy", showWarnings = FALSE)
dane_H5_wykres <- dane_H5 %>%
select(
plec,
UCLA_wynik,
SES_wynik,
PHQ9_wynik
) %>%
pivot_longer(
cols = c(UCLA_wynik, SES_wynik, PHQ9_wynik),
names_to = "Zmienna",
values_to = "Wynik"
) %>%
mutate(
Zmienna = recode(
Zmienna,
"UCLA_wynik" = "Poczucie osamotnienia",
"SES_wynik" = "Samoocena",
"PHQ9_wynik" = "Symptomy depresyjne"
)
)
p_H5_box <- ggplot(dane_H5_wykres, aes(x = plec, y = Wynik)) +
geom_boxplot(fill = "grey85", color = "grey30", outlier.shape = NA) +
geom_jitter(width = 0.15, alpha = 0.5) +
facet_wrap(~ Zmienna, scales = "free_y", ncol = 3) +
labs(
title = "Różnice między kobietami i mężczyznami w zakresie badanych zmiennych",
x = "Płeć",
y = "Wynik"
) +
theme_minimal(base_size = 12)
p_H5_box
ggsave(
filename = "outputs/wykresy/H5_roznice_plec_boxplot.png",
plot = p_H5_box,
width = 11,
height = 5,
dpi = 300
)
# KROK 46: Osobne wykresy dla H5 ------------------------------------------
zmienne_H5 <- data.frame(
zmienna = c("UCLA_wynik", "SES_wynik", "PHQ9_wynik"),
etykieta = c(
"Poczucie osamotnienia",
"Samoocena",
"Symptomy depresyjne"
),
plik = c(
"H5_poczucie_osamotnienia_plec.png",
"H5_samoocena_plec.png",
"H5_symptomy_depresyjne_plec.png"
)
)
for (i in 1:nrow(zmienne_H5)) {
var <- zmienne_H5$zmienna[i]
label <- zmienne_H5$etykieta[i]
file <- zmienne_H5$plik[i]
p <- ggplot(dane_H5, aes(x = plec, y = .data[[var]])) +
geom_boxplot(fill = "grey85", color = "grey30", outlier.shape = NA) +
geom_jitter(width = 0.15, alpha = 0.5) +
labs(
title = paste("Różnice ze względu na płeć:", label),
x = "Płeć",
y = "Wynik"
) +
theme_minimal(base_size = 12)
ggsave(
filename = file.path("outputs/wykresy", file),
plot = p,
width = 6,
height = 4,
dpi = 300
)
}
# KROK 46: Osobne wykresy dla H5 ------------------------------------------
zmienne_H5 <- data.frame(
zmienna = c("UCLA_wynik", "SES_wynik", "PHQ9_wynik"),
etykieta = c(
"Poczucie osamotnienia",
"Samoocena",
"Symptomy depresyjne"
),
plik = c(
"H5_poczucie_osamotnienia_plec.png",
"H5_samoocena_plec.png",
"H5_symptomy_depresyjne_plec.png"
)
)
for (i in 1:nrow(zmienne_H5)) {
var <- zmienne_H5$zmienna[i]
label <- zmienne_H5$etykieta[i]
file <- zmienne_H5$plik[i]
p <- ggplot(dane_H5, aes(x = plec, y = .data[[var]])) +
geom_boxplot(fill = "grey85", color = "grey30", outlier.shape = NA) +
geom_jitter(width = 0.15, alpha = 0.5) +
labs(
title = paste("Różnice ze względu na płeć:", label),
x = "Płeć",
y = "Wynik"
) +
theme_minimal(base_size = 12)
ggsave(
filename = file.path("outputs/wykresy", file),
plot = p,
width = 6,
height = 4,
dpi = 300
)
}
# KROK 48: Kontrola końcowa po analizie H5 --------------------------------
tabela_H5_opis
## # A tibble: 2 × 17
## plec n UCLA_M UCLA_SD UCLA_Me UCLA_Min UCLA_Max SES_M SES_SD SES_Me
## <fct> <int> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Kobieta 69 47.5 9.09 49 25 68 26.2 4.68 26
## 2 Mężczyzna 51 50.9 7.04 51 28 72 25.6 4.55 26
## # ℹ 7 more variables: SES_Min <dbl>, SES_Max <dbl>, PHQ9_M <dbl>,
## # PHQ9_SD <dbl>, PHQ9_Me <dbl>, PHQ9_Min <dbl>, PHQ9_Max <dbl>
tabela_H5_opis_long
## # A tibble: 6 × 8
## Zmienna plec n M SD Me Min Max
## <chr> <fct> <int> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Poczucie osamotnienia Kobieta 69 47.5 9.09 49 25 68
## 2 Poczucie osamotnienia Mężczyzna 51 50.9 7.04 51 28 72
## 3 Samoocena Kobieta 69 26.2 4.68 26 14 40
## 4 Samoocena Mężczyzna 51 25.6 4.55 26 12 39
## 5 Symptomy depresyjne Kobieta 69 12.4 5.3 13 0 25
## 6 Symptomy depresyjne Mężczyzna 51 13.6 3.78 14 0 25
tabela_H5_testy
## Zmienna W p r_rangowo_dwuseryjne
## 1 Poczucie osamotnienia 1470.0 0.1243 -0.165
## 2 Samoocena 1783.5 0.9003 0.014
## 3 Symptomy depresyjne 1429.5 0.0788 -0.188
## Istotnosc Sila_efektu Kierunek
## 1 nieistotne statystycznie słaba wyższe wyniki u mężczyzn
## 2 nieistotne statystycznie bardzo słaba wyższe wyniki u kobiet
## 3 nieistotne statystycznie słaba wyższe wyniki u mężczyzn
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
list.files("outputs/wykresy")
## [1] "boxploty_zmiennych_glownych.png"
## [2] "H1_osamotnienie_samoocena.png"
## [3] "H2_osamotnienie_depresja.png"
## [4] "H3_samoocena_depresja.png"
## [5] "H5_poczucie_osamotnienia_plec.png"
## [6] "H5_roznice_plec_boxplot.png"
## [7] "H5_samoocena_plec.png"
## [8] "H5_symptomy_depresyjne_plec.png"
## [9] "H6_model_mediacji.png"
## [10] "H7_moderacja_samoocena.png"
## [11] "histogramy_zmiennych_glownych.png"
## [12] "qqploty_zmiennych_glownych.png"
## [13] "regresja_obserwowane_przewidywane_PHQ9.png"
## [14] "regresja_qq_reszty.png"
## [15] "regresja_reszty_wartosci_przewidywane.png"
# KROK 49: Pakiety do analizy regresji ------------------------------------
library(dplyr)
library(ggplot2)
library(writexl)
# Pakiety pomocnicze do regresji
# install.packages(c("broom", "car", "lmtest", "sandwich"))
library(broom)
library(car)
library(lmtest)
library(sandwich)
# KROK 50: Przygotowanie danych do regresji -------------------------------
dane_scores <- readRDS("outputs/dane_po_kodowaniu.rds")
dane_reg <- dane_scores %>%
select(
wiek,
plec,
UCLA_wynik,
SES_wynik,
PHQ9_wynik
) %>%
mutate(
UCLA_wynik = as.numeric(UCLA_wynik),
SES_wynik = as.numeric(SES_wynik),
PHQ9_wynik = as.numeric(PHQ9_wynik),
wiek = as.numeric(wiek)
)
# Kontrola danych
summary(dane_reg)
## wiek plec UCLA_wynik SES_wynik
## Min. :18.00 Length:122 Min. :25.00 Min. :12.00
## 1st Qu.:24.00 Class :character 1st Qu.:46.00 1st Qu.:24.00
## Median :27.00 Mode :character Median :50.00 Median :26.00
## Mean :26.75 Mean :48.83 Mean :26.02
## 3rd Qu.:29.00 3rd Qu.:53.75 3rd Qu.:28.00
## Max. :35.00 Max. :72.00 Max. :40.00
## PHQ9_wynik
## Min. : 0.00
## 1st Qu.:11.00
## Median :13.00
## Mean :12.82
## 3rd Qu.:15.00
## Max. :25.00
# Braki danych
colSums(is.na(dane_reg))
## wiek plec UCLA_wynik SES_wynik PHQ9_wynik
## 0 0 0 0 0
# KROK 51: Model regresji liniowej ----------------------------------------
model_reg <- lm(
PHQ9_wynik ~ UCLA_wynik + SES_wynik,
data = dane_reg
)
summary(model_reg)
##
## Call:
## lm(formula = PHQ9_wynik ~ UCLA_wynik + SES_wynik, data = dane_reg)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.0383 -2.1490 0.2275 2.2744 11.3786
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 12.09659 4.81654 2.511 0.013364 *
## UCLA_wynik 0.20241 0.05387 3.757 0.000267 ***
## SES_wynik -0.35199 0.09917 -3.549 0.000554 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.647 on 119 degrees of freedom
## Multiple R-squared: 0.4282, Adjusted R-squared: 0.4186
## F-statistic: 44.57 on 2 and 119 DF, p-value: 3.582e-15
# KROK 52: Standaryzowane współczynniki beta ------------------------------
dane_reg_z <- dane_reg %>%
mutate(
UCLA_z = as.numeric(scale(UCLA_wynik)),
SES_z = as.numeric(scale(SES_wynik)),
PHQ9_z = as.numeric(scale(PHQ9_wynik))
)
model_reg_std <- lm(
PHQ9_z ~ UCLA_z + SES_z,
data = dane_reg_z
)
summary(model_reg_std)
##
## Call:
## lm(formula = PHQ9_z ~ UCLA_z + SES_z, data = dane_reg_z)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.09879 -0.44931 0.04757 0.47552 2.37900
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.550e-17 6.903e-02 0.000 1.000000
## UCLA_z 3.649e-01 9.712e-02 3.757 0.000267 ***
## SES_z -3.447e-01 9.712e-02 -3.549 0.000554 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7625 on 119 degrees of freedom
## Multiple R-squared: 0.4282, Adjusted R-squared: 0.4186
## F-statistic: 44.57 on 2 and 119 DF, p-value: 3.582e-15
# KROK 53: Tabela współczynników regresji ---------------------------------
tabela_regresja <- broom::tidy(
model_reg,
conf.int = TRUE
) %>%
mutate(
beta_standaryzowany = case_when(
term == "UCLA_wynik" ~ coef(model_reg_std)["UCLA_z"],
term == "SES_wynik" ~ coef(model_reg_std)["SES_z"],
TRUE ~ NA_real_
),
term = dplyr::recode(
term,
"(Intercept)" = "Stała",
"UCLA_wynik" = "Poczucie osamotnienia",
"SES_wynik" = "Samoocena"
),
across(
c(estimate, std.error, statistic, p.value, conf.low, conf.high, beta_standaryzowany),
~ round(.x, 3)
)
) %>%
rename(
Predyktor = term,
B = estimate,
SE = std.error,
t = statistic,
p = p.value,
CI_95_dolny = conf.low,
CI_95_gorny = conf.high,
Beta = beta_standaryzowany
)
tabela_regresja
## # A tibble: 3 × 8
## Predyktor B SE t p CI_95_dolny CI_95_gorny Beta
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Stała 12.1 4.82 2.51 0.013 2.56 21.6 NA
## 2 Poczucie osamotnienia 0.202 0.054 3.76 0 0.096 0.309 0.365
## 3 Samoocena -0.352 0.099 -3.55 0.001 -0.548 -0.156 -0.345
# KROK 54: Dopasowanie modelu regresji ------------------------------------
tabela_dopasowanie_modelu <- broom::glance(model_reg) %>%
select(
r.squared,
adj.r.squared,
statistic,
p.value,
df,
df.residual,
AIC,
BIC
) %>%
mutate(
across(where(is.numeric), ~ round(.x, 3))
) %>%
rename(
R2 = r.squared,
R2_skorygowane = adj.r.squared,
F = statistic,
p_modelu = p.value,
df_modelu = df,
df_reszt = df.residual
)
tabela_dopasowanie_modelu
## # A tibble: 1 × 8
## R2 R2_skorygowane F p_modelu df_modelu df_reszt AIC BIC
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0.428 0.419 44.6 0 2 119 667. 678.
# KROK 55: Współliniowość predyktorów -------------------------------------
vif_model <- car::vif(model_reg)
tabela_vif <- data.frame(
Predyktor = names(vif_model),
VIF = round(as.numeric(vif_model), 3)
) %>%
mutate(
Predyktor = dplyr::recode(
Predyktor,
"UCLA_wynik" = "Poczucie osamotnienia",
"SES_wynik" = "Samoocena"
)
)
tabela_vif
## Predyktor VIF
## 1 Poczucie osamotnienia 1.963
## 2 Samoocena 1.963
# KROK 56: Sprawdzenie założeń regresji -----------------------------------
# Reszty modelu
reszty_modelu <- residuals(model_reg)
# Test Shapiro-Wilka dla reszt
shapiro_reszty <- shapiro.test(reszty_modelu)
# Test Breuscha-Pagana dla homoskedastyczności
bp_test <- lmtest::bptest(model_reg)
# Obserwacje wpływowe — Cook's distance
cook <- cooks.distance(model_reg)
prog_cook <- 4 / nrow(dane_reg)
liczba_wplywowych <- sum(cook > prog_cook)
tabela_zalozenia_regresji <- data.frame(
Element = c(
"Normalność reszt",
"Homoskedastyczność",
"Obserwacje wpływowe"
),
Test = c(
"Shapiro-Wilk",
"Breusch-Pagan",
"Cook's distance > 4/n"
),
Statystyka = c(
round(as.numeric(shapiro_reszty$statistic), 3),
round(as.numeric(bp_test$statistic), 3),
liczba_wplywowych
),
p = c(
round(shapiro_reszty$p.value, 4),
round(bp_test$p.value, 4),
NA
)
)
tabela_zalozenia_regresji
## Element Test Statystyka p
## 1 Normalność reszt Shapiro-Wilk 0.988 0.3767
## 2 Homoskedastyczność Breusch-Pagan 0.465 0.7924
## 3 Obserwacje wpływowe Cook's distance > 4/n 13.000 NA
# KROK 57: Regresja z odpornymi błędami standardowymi HC3 -----------------
regresja_HC3 <- lmtest::coeftest(
model_reg,
vcov. = sandwich::vcovHC(model_reg, type = "HC3")
)
tabela_regresja_HC3 <- data.frame(
Predyktor = rownames(regresja_HC3),
B = regresja_HC3[, "Estimate"],
SE_HC3 = regresja_HC3[, "Std. Error"],
t = regresja_HC3[, "t value"],
p = regresja_HC3[, "Pr(>|t|)"],
row.names = NULL
) %>%
mutate(
Predyktor = dplyr::recode(
Predyktor,
"(Intercept)" = "Stała",
"UCLA_wynik" = "Poczucie osamotnienia",
"SES_wynik" = "Samoocena"
),
across(
c(B, SE_HC3, t, p),
~ round(.x, 3)
)
)
tabela_regresja_HC3
## Predyktor B SE_HC3 t p
## 1 Stała 12.097 6.029 2.006 0.047
## 2 Poczucie osamotnienia 0.202 0.070 2.892 0.005
## 3 Samoocena -0.352 0.117 -3.004 0.003
# KROK 58: Wykres reszt względem wartości przewidywanych ------------------
dir.create("outputs/wykresy", showWarnings = FALSE)
dane_reg_diag <- broom::augment(model_reg, data = dane_reg)
p_reszty <- ggplot(dane_reg_diag, aes(x = .fitted, y = .resid)) +
geom_point(alpha = 0.7) +
geom_hline(yintercept = 0, linetype = "dashed") +
geom_smooth(method = "loess", se = TRUE) +
labs(
title = "Wykres reszt względem wartości przewidywanych",
x = "Wartości przewidywane",
y = "Reszty"
) +
theme_minimal(base_size = 12)
p_reszty
ggsave(
filename = "outputs/wykresy/regresja_reszty_wartosci_przewidywane.png",
plot = p_reszty,
width = 7,
height = 5,
dpi = 300
)
# KROK 59: Wykres Q-Q reszt modelu ----------------------------------------
p_qq_reszty <- ggplot(dane_reg_diag, aes(sample = .std.resid)) +
stat_qq(color = "grey30") +
stat_qq_line() +
labs(
title = "Wykres Q-Q reszt modelu regresji",
x = "Kwantyle teoretyczne",
y = "Standaryzowane reszty"
) +
theme_minimal(base_size = 12)
p_qq_reszty
ggsave(
filename = "outputs/wykresy/regresja_qq_reszty.png",
plot = p_qq_reszty,
width = 7,
height = 5,
dpi = 300
)
# KROK 60: Wykres wartości obserwowanych i przewidywanych -----------------
p_pred <- ggplot(dane_reg_diag, aes(x = .fitted, y = PHQ9_wynik)) +
geom_point(alpha = 0.7) +
geom_smooth(method = "lm", se = TRUE) +
labs(
title = "Wartości obserwowane i przewidywane symptomów depresyjnych",
x = "Wartości przewidywane przez model",
y = "Wartości obserwowane PHQ-9"
) +
theme_minimal(base_size = 12)
p_pred
ggsave(
filename = "outputs/wykresy/regresja_obserwowane_przewidywane_PHQ9.png",
plot = p_pred,
width = 7,
height = 5,
dpi = 300
)
# KROK 61: Zapis wyników regresji -----------------------------------------
write_xlsx(
list(
"regresja_wspolczynniki" = tabela_regresja,
"regresja_dopasowanie" = tabela_dopasowanie_modelu,
"regresja_VIF" = tabela_vif,
"regresja_zalozenia" = tabela_zalozenia_regresji,
"regresja_HC3" = tabela_regresja_HC3
),
path = "outputs/tabele_regresja_liniowa.xlsx"
)
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
list.files("outputs/wykresy")
## [1] "boxploty_zmiennych_glownych.png"
## [2] "H1_osamotnienie_samoocena.png"
## [3] "H2_osamotnienie_depresja.png"
## [4] "H3_samoocena_depresja.png"
## [5] "H5_poczucie_osamotnienia_plec.png"
## [6] "H5_roznice_plec_boxplot.png"
## [7] "H5_samoocena_plec.png"
## [8] "H5_symptomy_depresyjne_plec.png"
## [9] "H6_model_mediacji.png"
## [10] "H7_moderacja_samoocena.png"
## [11] "histogramy_zmiennych_glownych.png"
## [12] "qqploty_zmiennych_glownych.png"
## [13] "regresja_obserwowane_przewidywane_PHQ9.png"
## [14] "regresja_qq_reszty.png"
## [15] "regresja_reszty_wartosci_przewidywane.png"
# KROK 62: Kontrola końcowa po regresji -----------------------------------
tabela_regresja
## # A tibble: 3 × 8
## Predyktor B SE t p CI_95_dolny CI_95_gorny Beta
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Stała 12.1 4.82 2.51 0.013 2.56 21.6 NA
## 2 Poczucie osamotnienia 0.202 0.054 3.76 0 0.096 0.309 0.365
## 3 Samoocena -0.352 0.099 -3.55 0.001 -0.548 -0.156 -0.345
tabela_dopasowanie_modelu
## # A tibble: 1 × 8
## R2 R2_skorygowane F p_modelu df_modelu df_reszt AIC BIC
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0.428 0.419 44.6 0 2 119 667. 678.
tabela_vif
## Predyktor VIF
## 1 Poczucie osamotnienia 1.963
## 2 Samoocena 1.963
tabela_zalozenia_regresji
## Element Test Statystyka p
## 1 Normalność reszt Shapiro-Wilk 0.988 0.3767
## 2 Homoskedastyczność Breusch-Pagan 0.465 0.7924
## 3 Obserwacje wpływowe Cook's distance > 4/n 13.000 NA
tabela_regresja_HC3
## Predyktor B SE_HC3 t p
## 1 Stała 12.097 6.029 2.006 0.047
## 2 Poczucie osamotnienia 0.202 0.070 2.892 0.005
## 3 Samoocena -0.352 0.117 -3.004 0.003
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
list.files("outputs/wykresy")
## [1] "boxploty_zmiennych_glownych.png"
## [2] "H1_osamotnienie_samoocena.png"
## [3] "H2_osamotnienie_depresja.png"
## [4] "H3_samoocena_depresja.png"
## [5] "H5_poczucie_osamotnienia_plec.png"
## [6] "H5_roznice_plec_boxplot.png"
## [7] "H5_samoocena_plec.png"
## [8] "H5_symptomy_depresyjne_plec.png"
## [9] "H6_model_mediacji.png"
## [10] "H7_moderacja_samoocena.png"
## [11] "histogramy_zmiennych_glownych.png"
## [12] "qqploty_zmiennych_glownych.png"
## [13] "regresja_obserwowane_przewidywane_PHQ9.png"
## [14] "regresja_qq_reszty.png"
## [15] "regresja_reszty_wartosci_przewidywane.png"
# KROK 63: Przygotowanie danych do analizy mediacji -----------------------
library(dplyr)
library(tidyr)
library(ggplot2)
library(writexl)
library(broom)
library(lmtest)
library(grid)
# Wczytanie danych po kodowaniu
dane_scores <- readRDS("outputs/dane_po_kodowaniu.rds")
# Dane do mediacji
dane_med <- dane_scores %>%
select(
UCLA_wynik,
SES_wynik,
PHQ9_wynik
) %>%
drop_na()
# Kontrola liczebności i zakresów
nrow(dane_med)
## [1] 122
summary(dane_med)
## UCLA_wynik SES_wynik PHQ9_wynik
## Min. :25.00 Min. :12.00 Min. : 0.00
## 1st Qu.:46.00 1st Qu.:24.00 1st Qu.:11.00
## Median :50.00 Median :26.00 Median :13.00
## Mean :48.83 Mean :26.02 Mean :12.82
## 3rd Qu.:53.75 3rd Qu.:28.00 3rd Qu.:15.00
## Max. :72.00 Max. :40.00 Max. :25.00
# KROK 64: Modele regresji dla mediacji -----------------------------------
# Ścieżka a: poczucie osamotnienia -> samoocena
model_a <- lm(
SES_wynik ~ UCLA_wynik,
data = dane_med
)
# Ścieżki b i c': poczucie osamotnienia + samoocena -> symptomy depresyjne
model_b_cprime <- lm(
PHQ9_wynik ~ UCLA_wynik + SES_wynik,
data = dane_med
)
# Ścieżka c: efekt całkowity poczucia osamotnienia -> symptomy depresyjne
model_c_total <- lm(
PHQ9_wynik ~ UCLA_wynik,
data = dane_med
)
# Podgląd wyników
summary(model_a)
##
## Call:
## lm(formula = SES_wynik ~ UCLA_wynik, data = dane_med)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.1541 -2.0322 0.2531 1.7995 10.3679
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 44.60348 1.75457 25.42 <2e-16 ***
## UCLA_wynik -0.38050 0.03539 -10.75 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.357 on 120 degrees of freedom
## Multiple R-squared: 0.4906, Adjusted R-squared: 0.4864
## F-statistic: 115.6 on 1 and 120 DF, p-value: < 2.2e-16
summary(model_b_cprime)
##
## Call:
## lm(formula = PHQ9_wynik ~ UCLA_wynik + SES_wynik, data = dane_med)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.0383 -2.1490 0.2275 2.2744 11.3786
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 12.09659 4.81654 2.511 0.013364 *
## UCLA_wynik 0.20241 0.05387 3.757 0.000267 ***
## SES_wynik -0.35199 0.09917 -3.549 0.000554 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.647 on 119 degrees of freedom
## Multiple R-squared: 0.4282, Adjusted R-squared: 0.4186
## F-statistic: 44.57 on 2 and 119 DF, p-value: 3.582e-15
summary(model_c_total)
##
## Call:
## lm(formula = PHQ9_wynik ~ UCLA_wynik, data = dane_med)
##
## Residuals:
## Min 1Q Median 3Q Max
## -12.5773 -2.2162 0.2816 1.8222 11.8312
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.60327 1.99607 -1.805 0.0736 .
## UCLA_wynik 0.33634 0.04026 8.354 1.34e-13 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.819 on 120 degrees of freedom
## Multiple R-squared: 0.3677, Adjusted R-squared: 0.3624
## F-statistic: 69.79 on 1 and 120 DF, p-value: 1.341e-13
# KROK 65: Obliczenie efektu pośredniego, bezpośredniego i całkowitego ----
a <- coef(model_a)["UCLA_wynik"]
b <- coef(model_b_cprime)["SES_wynik"]
c_prime <- coef(model_b_cprime)["UCLA_wynik"]
c_total <- coef(model_c_total)["UCLA_wynik"]
efekt_posredni <- a * b
efekt_bezposredni <- c_prime
efekt_calkowity <- c_total
proporcja_mediacji <- efekt_posredni / efekt_calkowity
efekty_mediacji <- data.frame(
Efekt = c(
"a: poczucie osamotnienia -> samoocena",
"b: samoocena -> symptomy depresyjne",
"c: efekt całkowity",
"c': efekt bezpośredni",
"a*b: efekt pośredni",
"proporcja mediacji"
),
Wartosc = c(
a,
b,
efekt_calkowity,
efekt_bezposredni,
efekt_posredni,
proporcja_mediacji
)
) %>%
mutate(
Wartosc = round(Wartosc, 3)
)
efekty_mediacji
## Efekt Wartosc
## 1 a: poczucie osamotnienia -> samoocena -0.380
## 2 b: samoocena -> symptomy depresyjne -0.352
## 3 c: efekt całkowity 0.336
## 4 c': efekt bezpośredni 0.202
## 5 a*b: efekt pośredni 0.134
## 6 proporcja mediacji 0.398
# KROK 66: Tabela ścieżek mediacyjnych ------------------------------------
sciezka_a <- broom::tidy(model_a, conf.int = TRUE) %>%
filter(term == "UCLA_wynik") %>%
mutate(Sciezka = "a: poczucie osamotnienia -> samoocena")
sciezka_b <- broom::tidy(model_b_cprime, conf.int = TRUE) %>%
filter(term == "SES_wynik") %>%
mutate(Sciezka = "b: samoocena -> symptomy depresyjne")
sciezka_cprime <- broom::tidy(model_b_cprime, conf.int = TRUE) %>%
filter(term == "UCLA_wynik") %>%
mutate(Sciezka = "c': poczucie osamotnienia -> symptomy depresyjne")
sciezka_c <- broom::tidy(model_c_total, conf.int = TRUE) %>%
filter(term == "UCLA_wynik") %>%
mutate(Sciezka = "c: efekt całkowity poczucia osamotnienia")
tabela_mediacja_sciezki <- bind_rows(
sciezka_a,
sciezka_b,
sciezka_c,
sciezka_cprime
) %>%
select(
Sciezka,
estimate,
std.error,
statistic,
p.value,
conf.low,
conf.high
) %>%
mutate(
across(
c(estimate, std.error, statistic, p.value, conf.low, conf.high),
~ round(.x, 4)
)
) %>%
rename(
B = estimate,
SE = std.error,
t = statistic,
p = p.value,
CI_95_dolny = conf.low,
CI_95_gorny = conf.high
)
tabela_mediacja_sciezki
## # A tibble: 4 × 7
## Sciezka B SE t p CI_95_dolny CI_95_gorny
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 a: poczucie osamotnienia -… -0.380 0.0354 -10.8 0 -0.451 -0.310
## 2 b: samoocena -> symptomy d… -0.352 0.0992 -3.55 6e-4 -0.548 -0.156
## 3 c: efekt całkowity poczuci… 0.336 0.0403 8.35 0 0.257 0.416
## 4 c': poczucie osamotnienia … 0.202 0.0539 3.76 3e-4 0.0957 0.309
# KROK 67: Bootstrap efektu pośredniego -----------------------------------
set.seed(123)
n_boot <- 5000
n <- nrow(dane_med)
boot_efekty <- data.frame(
efekt_posredni = rep(NA, n_boot),
efekt_bezposredni = rep(NA, n_boot),
efekt_calkowity = rep(NA, n_boot),
proporcja_mediacji = rep(NA, n_boot)
)
for (i in 1:n_boot) {
indeksy <- sample(1:n, size = n, replace = TRUE)
dane_b <- dane_med[indeksy, ]
model_a_b <- lm(
SES_wynik ~ UCLA_wynik,
data = dane_b
)
model_b_cprime_b <- lm(
PHQ9_wynik ~ UCLA_wynik + SES_wynik,
data = dane_b
)
model_c_total_b <- lm(
PHQ9_wynik ~ UCLA_wynik,
data = dane_b
)
a_b <- coef(model_a_b)["UCLA_wynik"]
b_b <- coef(model_b_cprime_b)["SES_wynik"]
c_prime_b <- coef(model_b_cprime_b)["UCLA_wynik"]
c_total_b <- coef(model_c_total_b)["UCLA_wynik"]
boot_efekty$efekt_posredni[i] <- a_b * b_b
boot_efekty$efekt_bezposredni[i] <- c_prime_b
boot_efekty$efekt_calkowity[i] <- c_total_b
boot_efekty$proporcja_mediacji[i] <- (a_b * b_b) / c_total_b
}
# Przedziały ufności bootstrap
CI_posredni <- quantile(
boot_efekty$efekt_posredni,
probs = c(0.025, 0.975),
na.rm = TRUE
)
CI_bezposredni <- quantile(
boot_efekty$efekt_bezposredni,
probs = c(0.025, 0.975),
na.rm = TRUE
)
CI_calkowity <- quantile(
boot_efekty$efekt_calkowity,
probs = c(0.025, 0.975),
na.rm = TRUE
)
CI_proporcja <- quantile(
boot_efekty$proporcja_mediacji,
probs = c(0.025, 0.975),
na.rm = TRUE
)
# KROK 68: Tabela efektów mediacyjnych ------------------------------------
tabela_mediacja_efekty <- data.frame(
Efekt = c(
"Efekt pośredni a*b",
"Efekt bezpośredni c'",
"Efekt całkowity c",
"Proporcja mediacji"
),
Oszacowanie = c(
efekt_posredni,
efekt_bezposredni,
efekt_calkowity,
proporcja_mediacji
),
CI_95_dolny = c(
CI_posredni[1],
CI_bezposredni[1],
CI_calkowity[1],
CI_proporcja[1]
),
CI_95_gorny = c(
CI_posredni[2],
CI_bezposredni[2],
CI_calkowity[2],
CI_proporcja[2]
)
) %>%
mutate(
across(
c(Oszacowanie, CI_95_dolny, CI_95_gorny),
~ round(.x, 4)
),
Istotnosc = ifelse(
CI_95_dolny > 0 | CI_95_gorny < 0,
"istotny efekt",
"nieistotny efekt"
)
)
tabela_mediacja_efekty
## Efekt Oszacowanie CI_95_dolny CI_95_gorny Istotnosc
## 1 Efekt pośredni a*b 0.1339 0.0475 0.2382 istotny efekt
## 2 Efekt bezpośredni c' 0.2024 0.0645 0.3267 istotny efekt
## 3 Efekt całkowity c 0.3363 0.2347 0.4292 istotny efekt
## 4 Proporcja mediacji 0.3982 0.1462 0.7678 istotny efekt
# KROK 69: Standaryzowane współczynniki mediacji --------------------------
dane_med_z <- dane_med %>%
mutate(
UCLA_z = as.numeric(scale(UCLA_wynik)),
SES_z = as.numeric(scale(SES_wynik)),
PHQ9_z = as.numeric(scale(PHQ9_wynik))
)
model_a_z <- lm(
SES_z ~ UCLA_z,
data = dane_med_z
)
model_b_cprime_z <- lm(
PHQ9_z ~ UCLA_z + SES_z,
data = dane_med_z
)
model_c_total_z <- lm(
PHQ9_z ~ UCLA_z,
data = dane_med_z
)
a_z <- coef(model_a_z)["UCLA_z"]
b_z <- coef(model_b_cprime_z)["SES_z"]
c_prime_z <- coef(model_b_cprime_z)["UCLA_z"]
c_total_z <- coef(model_c_total_z)["UCLA_z"]
tabela_mediacja_standaryzowana <- data.frame(
Sciezka = c(
"a: poczucie osamotnienia -> samoocena",
"b: samoocena -> symptomy depresyjne",
"c: efekt całkowity",
"c': efekt bezpośredni",
"a*b: efekt pośredni standaryzowany"
),
Beta = c(
a_z,
b_z,
c_total_z,
c_prime_z,
a_z * b_z
)
) %>%
mutate(
Beta = round(Beta, 3)
)
tabela_mediacja_standaryzowana
## Sciezka Beta
## 1 a: poczucie osamotnienia -> samoocena -0.700
## 2 b: samoocena -> symptomy depresyjne -0.345
## 3 c: efekt całkowity 0.606
## 4 c': efekt bezpośredni 0.365
## 5 a*b: efekt pośredni standaryzowany 0.241
# KROK 70: Diagnostyka modeli mediacyjnych --------------------------------
# Reszty modelu a
reszty_a <- residuals(model_a)
# Reszty modelu b/c'
reszty_b_cprime <- residuals(model_b_cprime)
# Testy normalności reszt
shapiro_a <- shapiro.test(reszty_a)
shapiro_b_cprime <- shapiro.test(reszty_b_cprime)
# Testy homoskedastyczności
bp_a <- lmtest::bptest(model_a)
bp_b_cprime <- lmtest::bptest(model_b_cprime)
tabela_mediacja_zalozenia <- data.frame(
Model = c(
"Model a: samoocena ~ poczucie osamotnienia",
"Model b/c': symptomy depresyjne ~ poczucie osamotnienia + samoocena",
"Model a: samoocena ~ poczucie osamotnienia",
"Model b/c': symptomy depresyjne ~ poczucie osamotnienia + samoocena"
),
Zalozenie = c(
"normalność reszt",
"normalność reszt",
"homoskedastyczność",
"homoskedastyczność"
),
Test = c(
"Shapiro-Wilk",
"Shapiro-Wilk",
"Breusch-Pagan",
"Breusch-Pagan"
),
Statystyka = c(
round(as.numeric(shapiro_a$statistic), 3),
round(as.numeric(shapiro_b_cprime$statistic), 3),
round(as.numeric(bp_a$statistic), 3),
round(as.numeric(bp_b_cprime$statistic), 3)
),
p = c(
round(shapiro_a$p.value, 4),
round(shapiro_b_cprime$p.value, 4),
round(bp_a$p.value, 4),
round(bp_b_cprime$p.value, 4)
)
)
tabela_mediacja_zalozenia
## Model
## 1 Model a: samoocena ~ poczucie osamotnienia
## 2 Model b/c': symptomy depresyjne ~ poczucie osamotnienia + samoocena
## 3 Model a: samoocena ~ poczucie osamotnienia
## 4 Model b/c': symptomy depresyjne ~ poczucie osamotnienia + samoocena
## Zalozenie Test Statystyka p
## 1 normalność reszt Shapiro-Wilk 0.984 0.1741
## 2 normalność reszt Shapiro-Wilk 0.988 0.3767
## 3 homoskedastyczność Breusch-Pagan 5.307 0.0212
## 4 homoskedastyczność Breusch-Pagan 0.465 0.7924
# KROK 71: Wykres modelu mediacji -----------------------------------------
dir.create("outputs/wykresy", showWarnings = FALSE)
label_a <- paste0("a = ", round(a, 3))
label_b <- paste0("b = ", round(b, 3))
label_cprime <- paste0("c' = ", round(c_prime, 3))
label_indirect <- paste0("a*b = ", round(efekt_posredni, 3))
p_mediacja <- ggplot() +
xlim(0, 10) +
ylim(0, 6) +
annotate(
"rect",
xmin = 0.7, xmax = 3.2,
ymin = 2.4, ymax = 3.6,
fill = "grey85",
color = "grey30"
) +
annotate(
"text",
x = 1.95,
y = 3,
label = "Poczucie\nosamotnienia\n(X)",
size = 4
) +
annotate(
"rect",
xmin = 3.9, xmax = 6.1,
ymin = 4.4, ymax = 5.4,
fill = "grey85",
color = "grey30"
) +
annotate(
"text",
x = 5,
y = 4.9,
label = "Samoocena\n(M)",
size = 4
) +
annotate(
"rect",
xmin = 6.8, xmax = 9.4,
ymin = 2.4, ymax = 3.6,
fill = "grey85",
color = "grey30"
) +
annotate(
"text",
x = 8.1,
y = 3,
label = "Symptomy\ndepresyjne\n(Y)",
size = 4
) +
annotate(
"segment",
x = 3.2, xend = 4.0,
y = 3.6, yend = 4.5,
arrow = arrow(length = unit(0.2, "cm")),
color = "grey30"
) +
annotate(
"text",
x = 3.45,
y = 4.25,
label = label_a,
size = 4
) +
annotate(
"segment",
x = 6.1, xend = 6.8,
y = 4.5, yend = 3.6,
arrow = arrow(length = unit(0.2, "cm")),
color = "grey30"
) +
annotate(
"text",
x = 6.55,
y = 4.25,
label = label_b,
size = 4
) +
annotate(
"segment",
x = 3.2, xend = 6.8,
y = 3.0, yend = 3.0,
arrow = arrow(length = unit(0.2, "cm")),
color = "grey30"
) +
annotate(
"text",
x = 5,
y = 2.65,
label = label_cprime,
size = 4
) +
annotate(
"text",
x = 5,
y = 1.4,
label = label_indirect,
size = 4
) +
labs(
title = "Model mediacji: samoocena jako mediator relacji między osamotnieniem a symptomami depresyjnymi"
) +
theme_void() +
theme(
plot.title = element_text(hjust = 0.5, size = 12)
)
p_mediacja
ggsave(
filename = "outputs/wykresy/H6_model_mediacji.png",
plot = p_mediacja,
width = 9,
height = 5,
dpi = 300
)
# KROK 72: Zapis wyników mediacji -----------------------------------------
write_xlsx(
list(
"mediacja_sciezki" = tabela_mediacja_sciezki,
"mediacja_efekty_bootstrap" = tabela_mediacja_efekty,
"mediacja_standaryzowana" = tabela_mediacja_standaryzowana,
"mediacja_zalozenia" = tabela_mediacja_zalozenia,
"bootstrap_surowe" = boot_efekty
),
path = "outputs/tabele_H6_mediacja.xlsx"
)
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
list.files("outputs/wykresy")
## [1] "boxploty_zmiennych_glownych.png"
## [2] "H1_osamotnienie_samoocena.png"
## [3] "H2_osamotnienie_depresja.png"
## [4] "H3_samoocena_depresja.png"
## [5] "H5_poczucie_osamotnienia_plec.png"
## [6] "H5_roznice_plec_boxplot.png"
## [7] "H5_samoocena_plec.png"
## [8] "H5_symptomy_depresyjne_plec.png"
## [9] "H6_model_mediacji.png"
## [10] "H7_moderacja_samoocena.png"
## [11] "histogramy_zmiennych_glownych.png"
## [12] "qqploty_zmiennych_glownych.png"
## [13] "regresja_obserwowane_przewidywane_PHQ9.png"
## [14] "regresja_qq_reszty.png"
## [15] "regresja_reszty_wartosci_przewidywane.png"
# KROK 73: Kontrola końcowa po mediacji -----------------------------------
efekty_mediacji
## Efekt Wartosc
## 1 a: poczucie osamotnienia -> samoocena -0.380
## 2 b: samoocena -> symptomy depresyjne -0.352
## 3 c: efekt całkowity 0.336
## 4 c': efekt bezpośredni 0.202
## 5 a*b: efekt pośredni 0.134
## 6 proporcja mediacji 0.398
tabela_mediacja_sciezki
## # A tibble: 4 × 7
## Sciezka B SE t p CI_95_dolny CI_95_gorny
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 a: poczucie osamotnienia -… -0.380 0.0354 -10.8 0 -0.451 -0.310
## 2 b: samoocena -> symptomy d… -0.352 0.0992 -3.55 6e-4 -0.548 -0.156
## 3 c: efekt całkowity poczuci… 0.336 0.0403 8.35 0 0.257 0.416
## 4 c': poczucie osamotnienia … 0.202 0.0539 3.76 3e-4 0.0957 0.309
tabela_mediacja_efekty
## Efekt Oszacowanie CI_95_dolny CI_95_gorny Istotnosc
## 1 Efekt pośredni a*b 0.1339 0.0475 0.2382 istotny efekt
## 2 Efekt bezpośredni c' 0.2024 0.0645 0.3267 istotny efekt
## 3 Efekt całkowity c 0.3363 0.2347 0.4292 istotny efekt
## 4 Proporcja mediacji 0.3982 0.1462 0.7678 istotny efekt
tabela_mediacja_standaryzowana
## Sciezka Beta
## 1 a: poczucie osamotnienia -> samoocena -0.700
## 2 b: samoocena -> symptomy depresyjne -0.345
## 3 c: efekt całkowity 0.606
## 4 c': efekt bezpośredni 0.365
## 5 a*b: efekt pośredni standaryzowany 0.241
tabela_mediacja_zalozenia
## Model
## 1 Model a: samoocena ~ poczucie osamotnienia
## 2 Model b/c': symptomy depresyjne ~ poczucie osamotnienia + samoocena
## 3 Model a: samoocena ~ poczucie osamotnienia
## 4 Model b/c': symptomy depresyjne ~ poczucie osamotnienia + samoocena
## Zalozenie Test Statystyka p
## 1 normalność reszt Shapiro-Wilk 0.984 0.1741
## 2 normalność reszt Shapiro-Wilk 0.988 0.3767
## 3 homoskedastyczność Breusch-Pagan 5.307 0.0212
## 4 homoskedastyczność Breusch-Pagan 0.465 0.7924
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
list.files("outputs/wykresy")
## [1] "boxploty_zmiennych_glownych.png"
## [2] "H1_osamotnienie_samoocena.png"
## [3] "H2_osamotnienie_depresja.png"
## [4] "H3_samoocena_depresja.png"
## [5] "H5_poczucie_osamotnienia_plec.png"
## [6] "H5_roznice_plec_boxplot.png"
## [7] "H5_samoocena_plec.png"
## [8] "H5_symptomy_depresyjne_plec.png"
## [9] "H6_model_mediacji.png"
## [10] "H7_moderacja_samoocena.png"
## [11] "histogramy_zmiennych_glownych.png"
## [12] "qqploty_zmiennych_glownych.png"
## [13] "regresja_obserwowane_przewidywane_PHQ9.png"
## [14] "regresja_qq_reszty.png"
## [15] "regresja_reszty_wartosci_przewidywane.png"
# KROK 74: Przygotowanie danych do analizy moderacji ----------------------
library(dplyr)
library(tidyr)
library(ggplot2)
library(writexl)
library(broom)
library(lmtest)
library(sandwich)
# Wczytanie danych po kodowaniu
dane_scores <- readRDS("outputs/dane_po_kodowaniu.rds")
# Dane do moderacji
dane_mod <- dane_scores %>%
select(
UCLA_wynik,
SES_wynik,
PHQ9_wynik
) %>%
drop_na()
# Kontrola
nrow(dane_mod)
## [1] 122
summary(dane_mod)
## UCLA_wynik SES_wynik PHQ9_wynik
## Min. :25.00 Min. :12.00 Min. : 0.00
## 1st Qu.:46.00 1st Qu.:24.00 1st Qu.:11.00
## Median :50.00 Median :26.00 Median :13.00
## Mean :48.83 Mean :26.02 Mean :12.82
## 3rd Qu.:53.75 3rd Qu.:28.00 3rd Qu.:15.00
## Max. :72.00 Max. :40.00 Max. :25.00
# KROK 75: Centrowanie zmiennych i interakcja -----------------------------
srednia_UCLA <- mean(dane_mod$UCLA_wynik, na.rm = TRUE)
srednia_SES <- mean(dane_mod$SES_wynik, na.rm = TRUE)
sd_UCLA <- sd(dane_mod$UCLA_wynik, na.rm = TRUE)
sd_SES <- sd(dane_mod$SES_wynik, na.rm = TRUE)
dane_mod <- dane_mod %>%
mutate(
UCLA_c = UCLA_wynik - srednia_UCLA,
SES_c = SES_wynik - srednia_SES,
interakcja_UCLA_SES = UCLA_c * SES_c
)
# Kontrola
summary(dane_mod[, c("UCLA_c", "SES_c", "interakcja_UCLA_SES")])
## UCLA_c SES_c interakcja_UCLA_SES
## Min. :-23.828 Min. :-14.02459 Min. :-270.25013
## 1st Qu.: -2.828 1st Qu.: -2.02459 1st Qu.: -24.37923
## Median : 1.172 Median : -0.02459 Median : -4.19686
## Mean : 0.000 Mean : 0.00000 Mean : -28.06134
## 3rd Qu.: 4.922 3rd Qu.: 1.97541 3rd Qu.: 0.08183
## Max. : 23.172 Max. : 13.97541 Max. : 70.51216
# KROK 76: Modele moderacji -----------------------------------------------
# Model bez interakcji
model_mod_glowne <- lm(
PHQ9_wynik ~ UCLA_c + SES_c,
data = dane_mod
)
# Model z interakcją
model_mod <- lm(
PHQ9_wynik ~ UCLA_c * SES_c,
data = dane_mod
)
summary(model_mod_glowne)
##
## Call:
## lm(formula = PHQ9_wynik ~ UCLA_c + SES_c, data = dane_mod)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.0383 -2.1490 0.2275 2.2744 11.3786
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 12.81967 0.33017 38.827 < 2e-16 ***
## UCLA_c 0.20241 0.05387 3.757 0.000267 ***
## SES_c -0.35199 0.09917 -3.549 0.000554 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.647 on 119 degrees of freedom
## Multiple R-squared: 0.4282, Adjusted R-squared: 0.4186
## F-statistic: 44.57 on 2 and 119 DF, p-value: 3.582e-15
summary(model_mod)
##
## Call:
## lm(formula = PHQ9_wynik ~ UCLA_c * SES_c, data = dane_mod)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.2060 -2.1457 0.0847 2.2093 11.2995
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 13.023949 0.363196 35.859 < 2e-16 ***
## UCLA_c 0.180117 0.056254 3.202 0.001755 **
## SES_c -0.360364 0.099052 -3.638 0.000409 ***
## UCLA_c:SES_c 0.007280 0.005474 1.330 0.186154
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.635 on 118 degrees of freedom
## Multiple R-squared: 0.4367, Adjusted R-squared: 0.4224
## F-statistic: 30.49 on 3 and 118 DF, p-value: 1.147e-14
# Porównanie modeli
porownanie_modeli_mod <- anova(model_mod_glowne, model_mod)
porownanie_modeli_mod
## Analysis of Variance Table
##
## Model 1: PHQ9_wynik ~ UCLA_c + SES_c
## Model 2: PHQ9_wynik ~ UCLA_c * SES_c
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 119 1582.6
## 2 118 1559.3 1 23.367 1.7683 0.1862
# KROK 77: Standaryzowany model moderacji ---------------------------------
dane_mod_z <- dane_mod %>%
mutate(
UCLA_z = as.numeric(scale(UCLA_wynik)),
SES_z = as.numeric(scale(SES_wynik)),
PHQ9_z = as.numeric(scale(PHQ9_wynik)),
interakcja_z = UCLA_z * SES_z
)
model_mod_std <- lm(
PHQ9_z ~ UCLA_z + SES_z + interakcja_z,
data = dane_mod_z
)
summary(model_mod_std)
##
## Call:
## lm(formula = PHQ9_z ~ UCLA_z + SES_z + interakcja_z, data = dane_mod_z)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.13384 -0.44861 0.01771 0.46191 2.36246
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.04271 0.07594 0.562 0.574880
## UCLA_z 0.32473 0.10142 3.202 0.001755 **
## SES_z -0.35292 0.09701 -3.638 0.000409 ***
## interakcja_z 0.06148 0.04623 1.330 0.186154
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.76 on 118 degrees of freedom
## Multiple R-squared: 0.4367, Adjusted R-squared: 0.4224
## F-statistic: 30.49 on 3 and 118 DF, p-value: 1.147e-14
# KROK 78: Tabela współczynników moderacji --------------------------------
tabela_moderacja <- broom::tidy(
model_mod,
conf.int = TRUE
) %>%
mutate(
Beta = case_when(
term == "UCLA_c" ~ coef(model_mod_std)["UCLA_z"],
term == "SES_c" ~ coef(model_mod_std)["SES_z"],
term == "UCLA_c:SES_c" ~ coef(model_mod_std)["interakcja_z"],
TRUE ~ NA_real_
),
term = dplyr::recode(
term,
"(Intercept)" = "Stała",
"UCLA_c" = "Poczucie osamotnienia",
"SES_c" = "Samoocena",
"UCLA_c:SES_c" = "Poczucie osamotnienia × samoocena"
),
across(
c(estimate, std.error, statistic, p.value, conf.low, conf.high, Beta),
~ round(.x, 4)
)
) %>%
rename(
Predyktor = term,
B = estimate,
SE = std.error,
t = statistic,
p = p.value,
CI_95_dolny = conf.low,
CI_95_gorny = conf.high
)
tabela_moderacja
## # A tibble: 4 × 8
## Predyktor B SE t p CI_95_dolny CI_95_gorny Beta
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Stała 13.0 0.363 35.9 0 12.3 13.7 NA
## 2 Poczucie osamotni… 0.180 0.0563 3.20 0.0018 0.0687 0.292 0.325
## 3 Samoocena -0.360 0.0991 -3.64 0.0004 -0.556 -0.164 -0.353
## 4 Poczucie osamotni… 0.0073 0.0055 1.33 0.186 -0.0036 0.0181 0.0615
# KROK 79: Dopasowanie modelu moderacji -----------------------------------
dopasowanie_mod <- broom::glance(model_mod)
R2_glowne <- summary(model_mod_glowne)$r.squared
R2_interakcja <- summary(model_mod)$r.squared
delta_R2 <- R2_interakcja - R2_glowne
tabela_moderacja_dopasowanie <- data.frame(
Model = "Model z interakcją",
R2 = R2_interakcja,
R2_skorygowane = summary(model_mod)$adj.r.squared,
F_modelu = dopasowanie_mod$statistic,
df_modelu = dopasowanie_mod$df,
df_reszt = dopasowanie_mod$df.residual,
p_modelu = dopasowanie_mod$p.value,
Delta_R2_po_dodaniu_interakcji = delta_R2,
F_zmiany = porownanie_modeli_mod$F[2],
df_zmiany = porownanie_modeli_mod$Df[2],
p_zmiany = porownanie_modeli_mod$`Pr(>F)`[2]
) %>%
mutate(
across(where(is.numeric), ~ round(.x, 4))
)
tabela_moderacja_dopasowanie
## Model R2 R2_skorygowane F_modelu df_modelu df_reszt
## value Model z interakcją 0.4367 0.4224 30.4914 3 118
## p_modelu Delta_R2_po_dodaniu_interakcji F_zmiany df_zmiany p_zmiany
## value 0 0.0084 1.7683 1 0.1862
# KROK 80: Współliniowość w modelu moderacji ------------------------------
predyktory_vif <- c("UCLA_c", "SES_c", "interakcja_UCLA_SES")
oblicz_vif <- function(data, zmienna, predyktory) {
pozostale <- setdiff(predyktory, zmienna)
formula_vif <- as.formula(
paste(zmienna, "~", paste(pozostale, collapse = " + "))
)
model_vif <- lm(formula_vif, data = data)
r2 <- summary(model_vif)$r.squared
vif <- 1 / (1 - r2)
return(vif)
}
vif_moderacja <- sapply(
predyktory_vif,
function(x) oblicz_vif(dane_mod, x, predyktory_vif)
)
tabela_vif_moderacja <- data.frame(
Predyktor = names(vif_moderacja),
VIF = round(as.numeric(vif_moderacja), 3)
) %>%
mutate(
Predyktor = dplyr::recode(
Predyktor,
"UCLA_c" = "Poczucie osamotnienia",
"SES_c" = "Samoocena",
"interakcja_UCLA_SES" = "Poczucie osamotnienia × samoocena"
)
)
tabela_vif_moderacja
## Predyktor VIF
## 1 Poczucie osamotnienia 2.155
## 2 Samoocena 1.971
## 3 Poczucie osamotnienia × samoocena 1.145
# KROK 81: Diagnostyka modelu moderacji -----------------------------------
reszty_moderacja <- residuals(model_mod)
# Normalność reszt
shapiro_moderacja <- shapiro.test(reszty_moderacja)
# Homoskedastyczność
bp_moderacja <- lmtest::bptest(model_mod)
# Obserwacje wpływowe
cook_moderacja <- cooks.distance(model_mod)
prog_cook_mod <- 4 / nrow(dane_mod)
liczba_wplywowych_mod <- sum(cook_moderacja > prog_cook_mod)
tabela_moderacja_zalozenia <- data.frame(
Element = c(
"Normalność reszt",
"Homoskedastyczność",
"Obserwacje wpływowe"
),
Test = c(
"Shapiro-Wilk",
"Breusch-Pagan",
"Cook's distance > 4/n"
),
Statystyka = c(
round(as.numeric(shapiro_moderacja$statistic), 4),
round(as.numeric(bp_moderacja$statistic), 4),
liczba_wplywowych_mod
),
p = c(
round(shapiro_moderacja$p.value, 4),
round(bp_moderacja$p.value, 4),
NA
)
)
tabela_moderacja_zalozenia
## Element Test Statystyka p
## 1 Normalność reszt Shapiro-Wilk 0.9867 0.2770
## 2 Homoskedastyczność Breusch-Pagan 0.9327 0.8175
## 3 Obserwacje wpływowe Cook's distance > 4/n 7.0000 NA
# KROK 82: Moderacja z odpornymi błędami standardowymi HC3 ----------------
moderacja_HC3 <- lmtest::coeftest(
model_mod,
vcov. = sandwich::vcovHC(model_mod, type = "HC3")
)
tabela_moderacja_HC3 <- data.frame(
Predyktor = rownames(moderacja_HC3),
B = moderacja_HC3[, "Estimate"],
SE_HC3 = moderacja_HC3[, "Std. Error"],
t = moderacja_HC3[, "t value"],
p = moderacja_HC3[, "Pr(>|t|)"],
row.names = NULL
) %>%
mutate(
Predyktor = dplyr::recode(
Predyktor,
"(Intercept)" = "Stała",
"UCLA_c" = "Poczucie osamotnienia",
"SES_c" = "Samoocena",
"UCLA_c:SES_c" = "Poczucie osamotnienia × samoocena"
),
across(
c(B, SE_HC3, t, p),
~ round(.x, 4)
)
)
tabela_moderacja_HC3
## Predyktor B SE_HC3 t p
## 1 Stała 13.0239 0.3610 36.0776 0.0000
## 2 Poczucie osamotnienia 0.1801 0.0757 2.3807 0.0189
## 3 Samoocena -0.3604 0.1175 -3.0666 0.0027
## 4 Poczucie osamotnienia × samoocena 0.0073 0.0067 1.0838 0.2806
# KROK 83: Proste nachylenia ----------------------------------------------
vcov_HC3_mod <- sandwich::vcovHC(model_mod, type = "HC3")
oblicz_simple_slope <- function(model, vcov_mat, moderator_value, label, ses_raw) {
coefs <- coef(model)
b1 <- coefs["UCLA_c"]
b3 <- coefs["UCLA_c:SES_c"]
slope <- b1 + b3 * moderator_value
var_slope <- vcov_mat["UCLA_c", "UCLA_c"] +
moderator_value^2 * vcov_mat["UCLA_c:SES_c", "UCLA_c:SES_c"] +
2 * moderator_value * vcov_mat["UCLA_c", "UCLA_c:SES_c"]
se <- sqrt(var_slope)
t_value <- slope / se
df <- df.residual(model)
p_value <- 2 * pt(abs(t_value), df = df, lower.tail = FALSE)
ci_low <- slope - qt(0.975, df) * se
ci_high <- slope + qt(0.975, df) * se
data.frame(
Poziom_samooceny = label,
Wartosc_SES = round(ses_raw, 2),
B = round(slope, 4),
SE_HC3 = round(se, 4),
t = round(t_value, 4),
p = round(p_value, 4),
CI_95_dolny = round(ci_low, 4),
CI_95_gorny = round(ci_high, 4)
)
}
ses_low <- srednia_SES - sd_SES
ses_mean <- srednia_SES
ses_high <- srednia_SES + sd_SES
tabela_simple_slopes <- bind_rows(
oblicz_simple_slope(
model_mod,
vcov_HC3_mod,
moderator_value = -sd_SES,
label = "Niska samoocena (M - 1 SD)",
ses_raw = ses_low
),
oblicz_simple_slope(
model_mod,
vcov_HC3_mod,
moderator_value = 0,
label = "Średnia samoocena (M)",
ses_raw = ses_mean
),
oblicz_simple_slope(
model_mod,
vcov_HC3_mod,
moderator_value = sd_SES,
label = "Wysoka samoocena (M + 1 SD)",
ses_raw = ses_high
)
)
tabela_simple_slopes
## Poziom_samooceny Wartosc_SES B SE_HC3 t p
## UCLA_c...1 Niska samoocena (M - 1 SD) 21.34 0.1460 0.0934 1.5633 0.1207
## UCLA_c...2 Średnia samoocena (M) 26.02 0.1801 0.0757 2.3807 0.0189
## UCLA_c...3 Wysoka samoocena (M + 1 SD) 30.71 0.2142 0.0686 3.1235 0.0022
## CI_95_dolny CI_95_gorny
## UCLA_c...1 -0.0389 0.3310
## UCLA_c...2 0.0303 0.3299
## UCLA_c...3 0.0784 0.3500
# KROK 84: Wykres moderacji -----------------------------------------------
dir.create("outputs/wykresy", showWarnings = FALSE)
ucla_seq <- seq(
from = min(dane_mod$UCLA_wynik, na.rm = TRUE),
to = max(dane_mod$UCLA_wynik, na.rm = TRUE),
length.out = 100
)
poziomy_ses <- data.frame(
Poziom_samooceny = c(
"Niska samoocena (M - 1 SD)",
"Średnia samoocena (M)",
"Wysoka samoocena (M + 1 SD)"
),
SES_wynik = c(
ses_low,
ses_mean,
ses_high
)
)
nowe_dane_mod <- tidyr::crossing(
UCLA_wynik = ucla_seq,
poziomy_ses
) %>%
mutate(
UCLA_c = UCLA_wynik - srednia_UCLA,
SES_c = SES_wynik - srednia_SES
)
predykcje_mod <- predict(
model_mod,
newdata = nowe_dane_mod,
interval = "confidence"
)
nowe_dane_mod <- cbind(
nowe_dane_mod,
as.data.frame(predykcje_mod)
)
p_moderacja <- ggplot() +
geom_point(
data = dane_mod,
aes(x = UCLA_wynik, y = PHQ9_wynik),
alpha = 0.35
) +
geom_line(
data = nowe_dane_mod,
aes(
x = UCLA_wynik,
y = fit,
linetype = Poziom_samooceny
),
linewidth = 1
) +
geom_ribbon(
data = nowe_dane_mod,
aes(
x = UCLA_wynik,
ymin = lwr,
ymax = upr,
group = Poziom_samooceny
),
alpha = 0.12
) +
labs(
title = "Moderacyjna rola samooceny w relacji między osamotnieniem a symptomami depresyjnymi",
x = "Poczucie osamotnienia",
y = "Symptomy depresyjne",
linetype = "Poziom samooceny"
) +
theme_minimal(base_size = 12)
p_moderacja
ggsave(
filename = "outputs/wykresy/H7_moderacja_samoocena.png",
plot = p_moderacja,
width = 9,
height = 6,
dpi = 300
)
# KROK 85: Automatyczny wniosek dla H7 ------------------------------------
wiersz_interakcja <- tabela_moderacja_HC3 %>%
filter(Predyktor == "Poczucie osamotnienia × samoocena")
wniosek_H7 <- ifelse(
wiersz_interakcja$p < 0.05 & wiersz_interakcja$B < 0,
"H7 potwierdzona: samoocena istotnie osłabia związek między poczuciem osamotnienia a symptomami depresyjnymi.",
ifelse(
wiersz_interakcja$p < 0.05 & wiersz_interakcja$B > 0,
"Interakcja istotna, ale kierunek przeciwny do H7.",
"H7 niepotwierdzona: interakcja nie jest istotna statystycznie."
)
)
wniosek_H7
## [1] "H7 niepotwierdzona: interakcja nie jest istotna statystycznie."
# KROK 86: Zapis wyników moderacji ----------------------------------------
write_xlsx(
list(
"moderacja_wspolczynniki" = tabela_moderacja,
"moderacja_HC3" = tabela_moderacja_HC3,
"moderacja_dopasowanie" = tabela_moderacja_dopasowanie,
"moderacja_VIF" = tabela_vif_moderacja,
"moderacja_zalozenia" = tabela_moderacja_zalozenia,
"simple_slopes" = tabela_simple_slopes,
"wniosek_H7" = data.frame(Wniosek = wniosek_H7)
),
path = "outputs/tabele_H7_moderacja.xlsx"
)
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
list.files("outputs/wykresy")
## [1] "boxploty_zmiennych_glownych.png"
## [2] "H1_osamotnienie_samoocena.png"
## [3] "H2_osamotnienie_depresja.png"
## [4] "H3_samoocena_depresja.png"
## [5] "H5_poczucie_osamotnienia_plec.png"
## [6] "H5_roznice_plec_boxplot.png"
## [7] "H5_samoocena_plec.png"
## [8] "H5_symptomy_depresyjne_plec.png"
## [9] "H6_model_mediacji.png"
## [10] "H7_moderacja_samoocena.png"
## [11] "histogramy_zmiennych_glownych.png"
## [12] "qqploty_zmiennych_glownych.png"
## [13] "regresja_obserwowane_przewidywane_PHQ9.png"
## [14] "regresja_qq_reszty.png"
## [15] "regresja_reszty_wartosci_przewidywane.png"
# KROK 87: Kontrola końcowa po moderacji ----------------------------------
tabela_moderacja
## # A tibble: 4 × 8
## Predyktor B SE t p CI_95_dolny CI_95_gorny Beta
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Stała 13.0 0.363 35.9 0 12.3 13.7 NA
## 2 Poczucie osamotni… 0.180 0.0563 3.20 0.0018 0.0687 0.292 0.325
## 3 Samoocena -0.360 0.0991 -3.64 0.0004 -0.556 -0.164 -0.353
## 4 Poczucie osamotni… 0.0073 0.0055 1.33 0.186 -0.0036 0.0181 0.0615
tabela_moderacja_HC3
## Predyktor B SE_HC3 t p
## 1 Stała 13.0239 0.3610 36.0776 0.0000
## 2 Poczucie osamotnienia 0.1801 0.0757 2.3807 0.0189
## 3 Samoocena -0.3604 0.1175 -3.0666 0.0027
## 4 Poczucie osamotnienia × samoocena 0.0073 0.0067 1.0838 0.2806
tabela_moderacja_dopasowanie
## Model R2 R2_skorygowane F_modelu df_modelu df_reszt
## value Model z interakcją 0.4367 0.4224 30.4914 3 118
## p_modelu Delta_R2_po_dodaniu_interakcji F_zmiany df_zmiany p_zmiany
## value 0 0.0084 1.7683 1 0.1862
tabela_vif_moderacja
## Predyktor VIF
## 1 Poczucie osamotnienia 2.155
## 2 Samoocena 1.971
## 3 Poczucie osamotnienia × samoocena 1.145
tabela_moderacja_zalozenia
## Element Test Statystyka p
## 1 Normalność reszt Shapiro-Wilk 0.9867 0.2770
## 2 Homoskedastyczność Breusch-Pagan 0.9327 0.8175
## 3 Obserwacje wpływowe Cook's distance > 4/n 7.0000 NA
tabela_simple_slopes
## Poziom_samooceny Wartosc_SES B SE_HC3 t p
## UCLA_c...1 Niska samoocena (M - 1 SD) 21.34 0.1460 0.0934 1.5633 0.1207
## UCLA_c...2 Średnia samoocena (M) 26.02 0.1801 0.0757 2.3807 0.0189
## UCLA_c...3 Wysoka samoocena (M + 1 SD) 30.71 0.2142 0.0686 3.1235 0.0022
## CI_95_dolny CI_95_gorny
## UCLA_c...1 -0.0389 0.3310
## UCLA_c...2 0.0303 0.3299
## UCLA_c...3 0.0784 0.3500
wniosek_H7
## [1] "H7 niepotwierdzona: interakcja nie jest istotna statystycznie."
list.files("outputs")
## [1] "dane_po_kodowaniu.csv"
## [2] "dane_po_kodowaniu.rds"
## [3] "statystyki_skal_podstawowe.csv"
## [4] "tabele_H6_mediacja.xlsx"
## [5] "tabele_H7_moderacja.xlsx"
## [6] "tabele_korelacje_H1_H4.xlsx"
## [7] "tabele_krok_1_2.xlsx"
## [8] "tabele_regresja_liniowa.xlsx"
## [9] "tabele_rzetelnosc_i_opisowe.xlsx"
## [10] "tabele_rzetelnosc_opisowe_normalnosc.xlsx"
## [11] "wykresy"
list.files("outputs/wykresy")
## [1] "boxploty_zmiennych_glownych.png"
## [2] "H1_osamotnienie_samoocena.png"
## [3] "H2_osamotnienie_depresja.png"
## [4] "H3_samoocena_depresja.png"
## [5] "H5_poczucie_osamotnienia_plec.png"
## [6] "H5_roznice_plec_boxplot.png"
## [7] "H5_samoocena_plec.png"
## [8] "H5_symptomy_depresyjne_plec.png"
## [9] "H6_model_mediacji.png"
## [10] "H7_moderacja_samoocena.png"
## [11] "histogramy_zmiennych_glownych.png"
## [12] "qqploty_zmiennych_glownych.png"
## [13] "regresja_obserwowane_przewidywane_PHQ9.png"
## [14] "regresja_qq_reszty.png"
## [15] "regresja_reszty_wartosci_przewidywane.png"