library(haven)
library(psych)
library(car)
library(apaTables)
library(tidyverse)Faktory související s ochotou vyhledat
odbornou psychologickou pomoc
Balíčky
V rámci analýzy budou použity následující balíčky: haven (Wickham et al., 2022) pro import dat, psych (2023) pro výpočet odhadů reliability a pro korelační analýzu, tidyverse (Wickham et al., 2019) pro transformaci dat a tvorbu grafů, apaTables (Stanley, 2021) pro tvorbu regresních tabulek ve formátu APA a balíčky car (Fox & Weisberg, 2019) a broom (Robinson et al., 2023) pro diagnostiku regresních modelů.
Úvod
V následující práci se budeme faktory, které souvisejí s ochotou vyhledat odbornou psychologickou pomoc. Konkrétně se zaměříme na pohlaví, věk a vyjadřování emocí (emoční expresivitu). Předpokládáme přitom následující:
Muži mají menší ochotu vyhledat odbornou psychologickou pomoc než ženy (H1), protože mají sklon méně emoce vyjadřovat navenek a vyhledání této pomoci mohou častěji chápat jako známku selhání (Pawar & Varma, 2023).
Mladší lidé mají větší ochotu vyhledat odbornou psychologickou pomoc než starší lidé (H2), a to kvůli menší míře stigmatizace duševních poruch u mladší generace (Chen et al., 2020).
Emoční otevřenost kladně souvisí s ochotou vyhledat odbornou psychologickou pomoc (H3), protože řešení duševních potíží vyžaduje jistou otevřenost, takže lidé, kteří neradi projevují své emoci, se budou zdráhat služby psychologa využít (Komiya et al., 2000).
Emoční otevřenost vysvětluje dodatečný rozptyl v ochotě vyhledat odbornou psychologickou pomoc
i po kontrole efektu pohlaví a věku (H4).
Představení dat a plán analýzy
Sběr dat
Data byla získána v rámci předmětu Statistická analýza dat v psychologii. Sběr dat probíhal online metodou příležitostného výběru během listopadu 2022. Respondenti byli oslovování na sociální síti Facebook s prosbou o vyplnění dotazníku zaměřeného na ochotu vyhledat odbornou psychologickou pomoc a faktory, které s ní mohou souviset.
Pokud se respondenti rozhodli kliknout na odkaz pro vyplnění dotazníku, byli informováni o tom, že získaná data budou použita pro seminární práci a mohou být dále sdílena pro edukační i vědecké účely. Dále byli respondenti ujištěni o anonymitě výzkumu a požádáni o udělení informovaného souhlasu s účastní na výzkumu a se sběrem a použitím dat.
Respondenti, kteří informovaný souhlas udělili, pak mohli pokračovat ve vyplňování samotného dotazníku. První sekce dotazníku tvořily otázky týkající se základních demografických údajů (pohlaví, věk, zda mají partnera/partnerku, zda již mají vysokoškolský titul, zda v současnosti studují na vysoké škole a jak hodnotí svou finanční situaci); druhou sekci dotazníku tvořily položky Škály postojů k vyhledávání odborné psychologické pomoci; a třetí a poslední sekci tvořily položky i Škály emoční otevřenosti
Nástroje měření
Škála postojů k vyhledávání odborné psychologické pomoci (Attitudes Towards Seeking Professional Psychological Help – ATSPPH; Picco et al. (2016)) měří to, jak respondenti vnímají odbornou psychologickou pomoc (jako je psychologické poradenství či psychoterapie), tj. zda ji považují za přínosnou a užitečnou a byli by ochotni ji využít, nebo nikoli. Obsahuje 10 položek (např. “Kdybych prožíval/a vážnou emoční krizi, jsem si jistý/á, že by mi psychoterapie pomohla”), které jsou hodnoceny odpovědní škála od 1 = “Nesouhlasím” do 4 = “Souhlasím”.
Škála emoční otevřenosti (Emotional Expressivity Scale – EES; Kring et al. (1994)) měří sklon vyjadřovat emoce navenek, nebo naopak “držet v sobě”. Obsahuje celkem 17 položek (např. “Své emoce projevuji před ostatními lidmi.”), které jsou hodnoceny na škále od 1 = “Nikdy pravda” po 7 = “Vždy pravda”.
Import a základní přehled dat
Data si nejprve data importujeme do R prostřednictvím balíčku haven, neboť jsou ve formátu, který používá SPSS (.sav). Importovaný dataset si nazveme df.
# Import dat
df <- haven::read_sav("data/help_seeking.sav")
glimpse(df)
#> Rows: 172
#> Columns: 33
#> $ gender <dbl+lbl> 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 1, 1, 1, 1, 1, 1…
#> $ partner <dbl+lbl> 1, 1, 1, 2, 2, 2, 1, 1, 1, 1, 1, 1, 2, 2, 1, 2…
#> $ age <dbl> 22, 21, 20, 21, 21, 20, 21, 23, 21, 23, 22, 22, 21…
#> $ college_title <dbl+lbl> 1, 1, 1, 1, 1, 1, 1, 2, 1, 1, 1, 1, 1, 1, 1, 1…
#> $ college_student <dbl+lbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2…
#> $ finantial_situation <dbl+lbl> 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 3, 2, 2, 2, 3, 2…
#> $ atspph01 <dbl+lbl> 4, 2, 3, 3, 1, 2, 3, 2, 2, 1, 2, 2, 3, 2, 3, 3…
#> $ atspph02 <dbl+lbl> 1, 3, 1, 1, 1, 1, 1, 1, 1, 2, 1, 1, 1, 1, 1, 1…
#> $ atspph03 <dbl+lbl> 4, 2, 3, 3, 2, 3, 3, 4, 3, 2, 3, 3, 4, 3, 3, 4…
#> $ atspph04 <dbl+lbl> 3, 4, 3, 2, 4, 3, 2, 4, 3, 4, 2, 4, 4, 2, 3, 3…
#> $ atspph05 <dbl+lbl> 4, 2, 3, 3, 3, 3, 3, 4, 4, 3, 3, 4, 4, 3, 4, 4…
#> $ atspph06 <dbl+lbl> 4, 2, 3, 2, 3, 4, 4, 4, 3, 3, 1, 3, 3, 3, 4, 2…
#> $ atspph07 <dbl+lbl> 4, 3, 3, 3, 3, 3, 3, 4, 3, 3, 3, 3, 3, 4, 3, 3…
#> $ atspph08 <dbl+lbl> 1, 3, 2, 2, 3, 2, 1, 3, 2, 2, 4, 1, 1, 2, 1, 2…
#> $ atspph09 <dbl+lbl> 1, 2, 2, 1, 2, 2, 1, 1, 1, 3, 2, 1, 1, 3, 1, 1…
#> $ atspph10 <dbl+lbl> 1, 2, 2, 2, 2, 2, 1, 3, 2, 2, 2, 2, 2, 2, 1, 1…
#> $ ees01 <dbl+lbl> 2, 2, 5, 3, 2, 4, 2, 3, 2, 3, 5, 3, 2, 3, 2, 1…
#> $ ees02 <dbl+lbl> 2, 2, 5, 3, 5, 3, 4, 2, 3, 3, 4, 2, 3, 3, 3, 2…
#> $ ees03 <dbl+lbl> 5, 3, 3, 5, 5, 5, 4, 4, 3, 4, 2, 5, 5, 4, 6, 5…
#> $ ees04 <dbl+lbl> 5, 2, 2, 4, 2, 3, 4, 5, 4, 4, 2, 4, 4, 3, 4, 5…
#> $ ees05 <dbl+lbl> 2, 5, 5, 3, 5, 3, 4, 3, 3, 3, 5, 3, 3, 4, 2, 1…
#> $ ees06 <dbl+lbl> 3, 4, 5, 3, 4, 4, 2, 2, 3, 2, 4, 2, 3, 3, 3, 2…
#> $ ees07 <dbl+lbl> 5, 3, 2, 4, 2, 2, 3, 5, 4, 3, 2, 4, 5, 4, 5, 3…
#> $ ees08 <dbl+lbl> 1, 3, 4, 1, 1, 3, 2, 2, 2, 1, 3, 2, 2, 1, 1, 1…
#> $ ees09 <dbl+lbl> 4, 5, 6, 3, 6, 5, 3, 2, 3, 3, 4, 3, 3, 3, 3, 2…
#> $ ees10 <dbl+lbl> 4, 2, 2, 3, 2, 5, 3, 5, 4, 5, 2, 2, 4, 5, 5, 2…
#> $ ees11 <dbl+lbl> 2, 3, 5, 3, 5, 2, 3, 2, 3, 3, 4, 2, 2, 2, 2, 1…
#> $ ees12 <dbl+lbl> 1, 5, 3, 3, 2, 2, 2, 2, 4, 2, 2, 2, 2, 2, 3, 1…
#> $ ees13 <dbl+lbl> 4, 2, 1, 3, 1, 1, 2, 2, 3, 2, 1, 3, 4, 3, 4, 6…
#> $ ees14 <dbl+lbl> 4, 4, 5, 3, 5, 2, 3, 4, 4, 5, 4, 3, 3, 2, 3, 4…
#> $ ees15 <dbl+lbl> 6, 2, 2, 4, 4, 4, 3, 4, 4, 6, 3, 2, 5, 3, 5, 6…
#> $ ees16 <dbl+lbl> 3, 4, 5, 3, 4, 3, 3, 2, 3, 5, 4, 2, 2, 3, 2, 2…
#> $ ees17 <dbl+lbl> 2, 5, 5, 3, 5, 3, 3, 3, 3, 3, 4, 3, 2, 2, 3, 1…Jak můžeme vidět data obsahují celkem 172 řádků a 33 sloupců, výzkumný soubor tudíž tvořilo N = 172 osob. Podíváme se nyní na to, zda se v datech vyskytují nějaké chybějící hodnoty.
df %>%
is.na() %>%
sum()
#> [1] 0Jak můžeme vidět, v datech se nenacházejí žádné chybějící hodnoty. Tento výsledek je téměř jistě zapříčiněn tím, že data byla sbírána skrze Google Forms a ty zaznamenávají pouze finalizované dotazníky (když někdo přeruší vyplňování dotazníku a neklikne na tlačítko “odeslat”, odpověď se nezaznamená).
Popis výzkumného souboru
Aby se nám s kategorickými proměnnými lépe pracovalo, převedeme si je na faktory.
df <- df %>%
mutate(
across(c(gender, partner, starts_with("college"),
finantial_situation),
as_factor)
)Abychom omezili opakování kódu, vytvoříme si funkci freq() pro výpočet absolutních a relativních četností, funkci ci_mean() pro výpočet statistik týkajících se průměrů a funkci quartiles()pro výpočet kvartilů.
# Funkce pro výpočet četností
freq <- function(...) {
df %>%
count(...) %>%
mutate(prc = n/sum(n)*100)
}
# Funkce pro výpočet průměrů a dalších statistik pro numerické proměnné
ci_mean <- function(x, level = .95, na.rm = TRUE) {
if (na.rm == TRUE) {
x <- x[!is.na(x)]
}
n <- length(x)
df <- n - 1
s <- sd(x)
m <- mean(x)
se <- s/sqrt(n)
alpha <- 1 - .95
t <- qt(1 - alpha/2, df = df)
me <- t*se
ci_lower <- m - me
ci_upper <- m + me
tibble(m = m,
ci_lower = ci_lower,
ci_upper = ci_upper,
sd = s,
se = se,
n = n,
skew = psych::skew(x),
kurt = psych::kurtosi(x))
}
# Funkce pro výpočet kvartilů
quartiles <- function(x,
na.rm = TRUE) {
tibble(
stat = c("Min", "q25", "Mdn", "q75", "Max"),
value = quantile(x,
probs = c(0, .25, .5, .75, 1),
na.rm = na.rm)
) %>%
pivot_wider(names_from = stat, values_from = value)
}Jak můžeme vidět níže, ženy tvořilo téměř 80 % respondentů a pouze tři osoby se označily za nebinární. Co se týče partnerských vztahů, asi 62 % osob uvedlo, že má v současné době partnera
freq(gender)
#> # A tibble: 3 × 3
#> gender n prc
#> <fct> <int> <dbl>
#> 1 Žena 134 77.9
#> 2 Muž 35 20.3
#> 3 Jiné 3 1.74
freq(partner)
#> # A tibble: 2 × 3
#> partner n prc
#> <fct> <int> <dbl>
#> 1 Ne 66 38.4
#> 2 Ano 106 61.6Věk respondentů se pohyboval od 16 do 78 let (M = 30,6, SD = 13,3 let), ale rozdělení věku bylo silně pravostranně zešikmené, takže přibližně dvě třetiny respondentů bylo ve věku 16–30 let.
df %>%
summarise(ci_mean(age))
#> # A tibble: 1 × 8
#> m ci_lower ci_upper sd se n skew kurt
#> <dbl> <dbl> <dbl> <dbl> <dbl> <int> <dbl> <dbl>
#> 1 30.6 28.6 32.6 13.3 1.01 172 1.13 0.330
df %>%
summarise(quartiles(age))
#> # A tibble: 1 × 5
#> Min q25 Mdn q75 Max
#> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 16 21 24 40 78
df %>%
ggplot(aes(age, y = after_stat(density)*2)) +
geom_histogram(binwidth = 2, center = 20,
fill = "grey", color = "black") +
geom_density(linewidth = 1, color = "red") +
labs(x = "Věk v letech", y = "Relativní četnost") +
scale_y_continuous(labels = scales::percent,
breaks = seq(0, 1, by = 0.05)) +
scale_x_continuous(breaks = seq(20, 80, by = 10))
mean(df$age <= 30)
#> [1] 0.6453488
Pokud jde o vzdělání, jak můžeme vidět níže, pouze asi 31 % osob uvedlo, že má vysokoškolský titul a přibližně 37 % uvedlo, že v současnosti studuje vysokou školu.
# Vysokoškolský titul
freq(college_title)
#> # A tibble: 2 × 3
#> college_title n prc
#> <fct> <int> <dbl>
#> 1 Ne 118 68.6
#> 2 Ano 54 31.4
# Studium VŠ v současnosti
freq(college_student)
#> # A tibble: 2 × 3
#> college_student n prc
#> <fct> <int> <dbl>
#> 1 Ne 109 63.4
#> 2 Ano 63 36.6Konečně pokud jde o finanční situaci (viz tabulka a graf níže), většina respondentů (69 %) ji hodnotila jako přijatelnou.
freq(finantial_situation)
#> # A tibble: 3 × 3
#> finantial_situation n prc
#> <fct> <int> <dbl>
#> 1 Špatná 14 8.14
#> 2 Přijatelná 118 68.6
#> 3 Dobrá 40 23.3
df %>%
count(finantial_situation) %>%
mutate(p = n/sum(n)) %>%
ggplot(aes(finantial_situation, p)) +
geom_col() +
scale_y_continuous(labels = scales::percent) +
labs(x = "Finanční situace", y = "Relativní četnost")
Plán analýzy
Nejprve provedeme základní položkovou analýzu ATSPPH a EES a odhadneme reliabilitu (vnitřní konzistenci) obou škál pomocí Cronbachova koeficientu alfa a McDonaldova koeficientu omega. Poté vypočteme celkové skóry obou škál jako průměr položek. Rozdíly mezi muži a ženami v ATSPPH (H1) ověříme pomocí Welchova t-testu pro nezávislé výběru. Bivariační vztah mezi věkem a ATSPPH ověříme pomocí Spearmanova korelačního koeficientu a bivariační vztah mezi EES a ATSPPH pomocí Pearsonova korelačního koeficientu. Nakonec pro ověření hypotézy H4 odhadneme lineárně regresní model s ATSPPH jako závislou proměnnou. Prediktory do něj zadáme ve dvou krocích: v prvních kroku do modelu jako prediktory pohlaví a věk, ve druhém kroku EES.
Výsledky
Položková analýza
Nejdříve si připravíme dva textové vektory s názvy reverzních položek.
rev_help <- str_c("atspph",
str_pad(c(2, 4, 8, 9, 10),
pad = 0,
width = 2))
rev_emotion <- str_c("ees",
str_pad(c(1, 2, 5, 6, 8, 9, 11, 12, 14, 16, 17),
pad = 0,
width = 2))Poté skóry reverzních položek obrátíme.
df <- df %>%
mutate(
across(all_of(rev_help), ~5 - .x),
across(all_of(rev_emotion), ~7 - .x)
)Po rekódování se můžeme podívat na reliabilitu obou škál.
df %>%
select(starts_with("atspph")) %>%
psych::alpha()
#>
#> Reliability analysis
#> Call: psych::alpha(x = .)
#>
#> raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
#> 0.82 0.82 0.84 0.32 4.7 0.02 3 0.54 0.31
#>
#> 95% confidence boundaries
#> lower alpha upper
#> Feldt 0.78 0.82 0.86
#> Duhachek 0.78 0.82 0.86
#>
#> Reliability if an item is dropped:
#> raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
#> atspph01 0.80 0.80 0.81 0.31 4.1 0.022 0.023 0.31
#> atspph02 0.80 0.81 0.82 0.31 4.1 0.022 0.025 0.35
#> atspph03 0.79 0.79 0.80 0.30 3.8 0.023 0.019 0.30
#> atspph04 0.84 0.84 0.84 0.37 5.3 0.018 0.012 0.39
#> atspph05 0.79 0.79 0.80 0.29 3.8 0.024 0.018 0.30
#> atspph06 0.81 0.81 0.81 0.31 4.1 0.022 0.024 0.35
#> atspph07 0.81 0.81 0.82 0.33 4.4 0.021 0.025 0.39
#> atspph08 0.80 0.80 0.82 0.31 4.0 0.022 0.024 0.31
#> atspph09 0.79 0.80 0.81 0.30 3.9 0.023 0.026 0.30
#> atspph10 0.82 0.82 0.83 0.33 4.5 0.021 0.025 0.39
#>
#> Item statistics
#> n raw.r std.r r.cor r.drop mean sd
#> atspph01 172 0.65 0.65 0.61 0.54 2.5 0.94
#> atspph02 172 0.64 0.64 0.58 0.53 3.5 0.84
#> atspph03 172 0.73 0.73 0.71 0.64 3.1 0.84
#> atspph04 172 0.32 0.32 0.20 0.16 2.0 0.90
#> atspph05 172 0.76 0.76 0.75 0.68 3.3 0.87
#> atspph06 172 0.65 0.64 0.59 0.52 3.0 0.96
#> atspph07 172 0.56 0.57 0.50 0.45 3.3 0.78
#> atspph08 172 0.67 0.66 0.61 0.56 2.8 0.93
#> atspph09 172 0.72 0.72 0.68 0.62 3.3 0.91
#> atspph10 172 0.51 0.52 0.44 0.39 3.4 0.77
#>
#> Non missing response frequency for each item
#> 1 2 3 4 miss
#> atspph01 0.15 0.35 0.33 0.17 0
#> atspph02 0.04 0.10 0.19 0.66 0
#> atspph03 0.06 0.14 0.47 0.33 0
#> atspph04 0.36 0.36 0.22 0.06 0
#> atspph05 0.06 0.08 0.37 0.49 0
#> atspph06 0.08 0.23 0.32 0.37 0
#> atspph07 0.02 0.13 0.39 0.45 0
#> atspph08 0.10 0.25 0.40 0.24 0
#> atspph09 0.06 0.12 0.30 0.52 0
#> atspph10 0.03 0.09 0.33 0.56 0
df %>%
select(starts_with("ees")) %>%
psych::alpha()
#>
#> Reliability analysis
#> Call: psych::alpha(x = .)
#>
#> raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
#> 0.92 0.92 0.93 0.4 11 0.0088 3.6 0.89 0.4
#>
#> 95% confidence boundaries
#> lower alpha upper
#> Feldt 0.9 0.92 0.94
#> Duhachek 0.9 0.92 0.94
#>
#> Reliability if an item is dropped:
#> raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
#> ees01 0.92 0.92 0.93 0.40 11 0.0092 0.019 0.40
#> ees02 0.91 0.91 0.93 0.40 11 0.0093 0.018 0.40
#> ees03 0.92 0.92 0.93 0.42 12 0.0086 0.016 0.42
#> ees04 0.91 0.91 0.93 0.40 11 0.0093 0.019 0.40
#> ees05 0.91 0.91 0.92 0.38 10 0.0100 0.016 0.39
#> ees06 0.91 0.92 0.93 0.40 11 0.0093 0.019 0.40
#> ees07 0.91 0.91 0.93 0.39 10 0.0098 0.018 0.38
#> ees08 0.92 0.92 0.93 0.40 11 0.0092 0.020 0.40
#> ees09 0.91 0.91 0.93 0.39 10 0.0097 0.018 0.39
#> ees10 0.92 0.92 0.93 0.41 11 0.0090 0.018 0.41
#> ees11 0.91 0.91 0.93 0.39 10 0.0098 0.018 0.38
#> ees12 0.92 0.92 0.93 0.41 11 0.0090 0.018 0.41
#> ees13 0.92 0.92 0.93 0.41 11 0.0088 0.018 0.41
#> ees14 0.91 0.91 0.93 0.39 10 0.0097 0.018 0.39
#> ees15 0.91 0.91 0.93 0.39 10 0.0096 0.019 0.39
#> ees16 0.92 0.92 0.93 0.41 11 0.0090 0.018 0.41
#> ees17 0.91 0.91 0.92 0.39 10 0.0100 0.017 0.38
#>
#> Item statistics
#> n raw.r std.r r.cor r.drop mean sd
#> ees01 172 0.61 0.61 0.57 0.55 3.6 1.3
#> ees02 172 0.65 0.65 0.62 0.59 3.7 1.4
#> ees03 172 0.43 0.44 0.39 0.36 3.9 1.4
#> ees04 172 0.63 0.64 0.61 0.58 3.0 1.3
#> ees05 172 0.84 0.83 0.84 0.80 3.4 1.4
#> ees06 172 0.62 0.63 0.60 0.57 3.7 1.2
#> ees07 172 0.79 0.79 0.78 0.76 3.4 1.3
#> ees08 172 0.62 0.64 0.61 0.57 5.1 1.0
#> ees09 172 0.74 0.73 0.72 0.69 3.2 1.5
#> ees10 172 0.57 0.57 0.53 0.50 3.2 1.4
#> ees11 172 0.77 0.77 0.76 0.73 3.8 1.4
#> ees12 172 0.52 0.53 0.50 0.46 4.5 1.2
#> ees13 172 0.56 0.54 0.51 0.48 2.8 1.6
#> ees14 172 0.75 0.75 0.74 0.71 3.6 1.4
#> ees15 172 0.74 0.74 0.73 0.69 3.7 1.4
#> ees16 172 0.55 0.55 0.52 0.48 3.5 1.3
#> ees17 172 0.82 0.81 0.81 0.78 3.2 1.5
#>
#> Non missing response frequency for each item
#> 1 2 3 4 5 6 miss
#> ees01 0.05 0.17 0.24 0.27 0.20 0.06 0
#> ees02 0.06 0.16 0.20 0.24 0.25 0.08 0
#> ees03 0.03 0.13 0.26 0.19 0.26 0.13 0
#> ees04 0.12 0.26 0.28 0.19 0.13 0.03 0
#> ees05 0.08 0.25 0.20 0.21 0.20 0.06 0
#> ees06 0.03 0.20 0.16 0.34 0.23 0.04 0
#> ees07 0.04 0.26 0.26 0.23 0.16 0.05 0
#> ees08 0.00 0.01 0.08 0.22 0.23 0.47 0
#> ees09 0.13 0.25 0.18 0.22 0.16 0.06 0
#> ees10 0.06 0.32 0.24 0.17 0.13 0.08 0
#> ees11 0.06 0.17 0.17 0.26 0.22 0.12 0
#> ees12 0.01 0.05 0.15 0.21 0.36 0.22 0
#> ees13 0.28 0.25 0.16 0.11 0.13 0.07 0
#> ees14 0.07 0.19 0.20 0.28 0.20 0.07 0
#> ees15 0.03 0.18 0.26 0.25 0.16 0.13 0
#> ees16 0.08 0.17 0.20 0.28 0.22 0.05 0
#> ees17 0.12 0.27 0.18 0.19 0.17 0.08 0Jak můžeme vidět, reliabilita škály ATSPPH byla solidní, Cronbachovo \(\alpha\) = 0,82 (95% CI [0,78; 0,86]). Všechny položky kladně korelovaly s celkovým skórem škály, ačkoli položka 4 (“Obdivuji lidi, kteří jsou ochotni zvládnout své problémy a starosti, aniž by vyhledali odbornou pomoc.”) poměrně slabě (r = 0,20).
Reliabilita škály EES byla také dobrá, Cronbachovo \(\alpha\) = 0,92 (95% CI [0,90; 0,94]) a všechny položky korelovaly s celkovým skórem škály kladně.
Výpočet a přehled celkových skórů
Jak bylo řečeno v plánu analýzy, celkové skóry obou použitých škál vypočteme jako průměr položek.
df <- df %>%
mutate(help = pick(starts_with("atspph")) %>%
rowMeans(),
emotion = pick(starts_with("ees")) %>%
rowMeans())Když je máme vypočteny, můžeme se podívat na jejich deskriptivní statistiky a rozdělení.
df %>%
pivot_longer(c(help, emotion), names_to = "variable") %>%
group_by(variable) %>%
summarise(ci_mean(value))
#> # A tibble: 2 × 9
#> variable m ci_lower ci_upper sd se n skew kurt
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <int> <dbl> <dbl>
#> 1 emotion 3.61 3.47 3.74 0.894 0.0681 172 -0.00940 -0.605
#> 2 help 3.00 2.92 3.09 0.545 0.0415 172 -0.823 0.644
df %>%
pivot_longer(c(help, emotion), names_to = "variable") %>%
group_by(variable) %>%
summarise(quartiles(value))
#> # A tibble: 2 × 6
#> variable Min q25 Mdn q75 Max
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 emotion 1.24 2.94 3.65 4.29 5.65
#> 2 help 1.1 2.8 3.1 3.32 3.9
df %>%
ggplot(aes(emotion, y = after_stat(density))) +
geom_histogram(binwidth = 0.2, center = 4) +
geom_density(color = "red", linewidth = 1) +
labs(x = "Emoční otevřenost")
df %>%
ggplot(aes(help, y = after_stat(density))) +
geom_histogram(binwidth = 0.1, center = 2) +
geom_density(color = "red", linewidth = 1) +
labs(x = "Ochota vyhledat odbnornou psychologickou pomoc")
Jak můžeme vidět výše, skóry emoční otevřenosti (EES) byly přibližně normálně rozděleny, zatímco skóry ochoty vyhledat odbnornou psychologickou pomoc byly mírně levostranně zešikmeny, takže se jejich vrchol jejich rozdělení nacházel nad pomyslným středem odpovědní škály (2,5), a to někde okolo hodnoty 3,0.
Testy hypotéz
Hypotéza 1
Nejprve se zaměříme na ochotu můžu a žen k vyhledání odborné psychologické pomoci. Vytvoříme si proměnnou gender_binary, která bude zahrnovat pouze muže a ženy, protože výzkumný soubor zahrnuje pouze tři nebinární osoby a naše hypotéza se týká pouze mužů a žen.
df <- df %>%
mutate(gender_binary = case_when(
gender == "Žena" ~ 1L,
gender == "Muž" ~ 2L,
TRUE ~ NA_integer_) %>%
factor(labels = c("Žena", "Muž")))Poté se podíváme, jak vypadá rozdělení naší závislé proměnné v obou skupinách.
df %>%
drop_na(gender_binary) %>%
ggplot(aes(help, y = after_stat(density))) +
geom_histogram(binwidth = 0.2) +
geom_density(color = "red", linewidth = 1) +
facet_wrap(~gender_binary)
df %>%
drop_na(gender_binary) %>%
ggplot(aes(sample = help)) +
geom_qq_line(color = "red", linewidth = 1) +
geom_qq() +
facet_wrap(~gender_binary)
Jak můžeme vidět, v obou skupinách je rozdělení ATSPPH podobně negativně zešikmené.
Nyní provedeme samotný nezávislý t-test a vypočteme si deskriptivní statistiky.
df %>%
drop_na(gender_binary) %>%
summarise(ci_mean(help),
.by = gender_binary)
#> # A tibble: 2 × 9
#> gender_binary m ci_lower ci_upper sd se n skew kurt
#> <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <int> <dbl> <dbl>
#> 1 Žena 3.08 2.99 3.16 0.506 0.0437 134 -0.857 1.15
#> 2 Muž 2.71 2.50 2.91 0.609 0.103 35 -0.443 -0.765
t.test(help ~ gender_binary, data = df)
#>
#> Welch Two Sample t-test
#>
#> data: help by gender_binary
#> t = 3.333, df = 47.002, p-value = 0.001681
#> alternative hypothesis: true difference in means between group Žena and group Muž is not equal to 0
#> 95 percent confidence interval:
#> 0.1477254 0.5975624
#> sample estimates:
#> mean in group Žena mean in group Muž
#> 3.078358 2.705714
effectsize::cohens_d(help ~ gender_binary, data = df,
pooled_sd = FALSE)
#> Cohen's d | 95% CI
#> ------------------------
#> 0.67 | [0.25, 1.08]
#>
#> - Estimated using un-pooled SD.Jak můžeme vidět výše, ženy (M = 3,08, SD = 2,99) vykazovaly v průměru vyšší skór ATSPPH než muži
(M = 3,08, SD = 2,50) a tento rozdíl byl statisticky významný, t(47) = 3,33, p = 0,002. Rozdíl mezi průměry činil \(\Delta\)M = 0,37 (95% CI [0,14; 0,59]) a Cohenovo d = 0.67 (95% CI [0,25; 1,08]). Výsledky tedy podporují hypotézu, že ženy mají větší ochotu vyhledat odbornou psychologickou pomoc než muži.
Hypotéza 2
Dále se podíváme na vztah mezi věkem a ATSPPH. Nejprve jej prozkoumáme pomocí scatterplotu.
df %>%
ggplot(aes(age, help)) +
geom_point(position = position_jitter(seed = 123),
alpha = .5, size = 2) +
geom_smooth(fill = "red", color = "red", alpha = 0.1) +
geom_smooth(fill = "blue", color = "blue", alpha = 0.1,
method = "lm", se = F)
Jak můžeme vidět, výsledky v tomto případě neodpovídají naší hypotéze, že ochota vyhledat odbnornou psychologickou pomoc je vyšší u mladších lidí. Zdá se naopak, že u starších lidí (nad 25 let) se skór ATSPPH příliš nemění zatímco u osob do 25 let skór ATSPPH s věkem roste.
df %>%
select(age, help) %>%
psych::cor.ci(plot = FALSE) %>%
print(digits = 3)
#> Call:corCi(x = x, keys = keys, n.iter = n.iter, p = p, overlap = overlap,
#> poly = poly, method = method, plot = plot, minlength = minlength,
#> n = n)
#>
#> Coefficients and bootstrapped confidence intervals
#> age help
#> age 1.00
#> help 0.03 1.00
#>
#> scale correlations and bootstrapped confidence intervals
#> lower.emp lower.norm estimate upper.norm upper.emp p
#> age-help -0.139 -0.129 0.028 0.16 0.138 0.827Nicméně Pearsonův korelační koeficient by nasvědčoval pouze slabému a nesignifikantnímu vztahu, r(170) = 0,03 (95% CI [–0,15; 0,18]), p = 0,887. Výsledky tedy nepodporují hypotézu, že starší lidé jsou méně ochotni vyhledat odbornou psychologickou pomoc než lidé mladší.
Hypotéza 3
Vztah mezi EES a ATSPPH opět prozkoumáme nejdříve pomocí scatterplotu.
df %>%
ggplot(aes(emotion, help)) +
geom_point(position = position_jitter(seed = 123),
alpha = .5, size = 2) +
geom_smooth(fill = "red", color = "red", alpha = 0.1) +
geom_smooth(fill = "blue", color = "blue", alpha = 0.1,
method = "lm", se = F)
Jak můžeme vidět, bivariační vztah mezi těmito proměnnými byl kladný a lineární. Pearsonův korelační koeficient by jej tedy měl dobře popsat.
df %>%
select(help, emotion) %>%
psych::cor.ci(plot = FALSE) %>%
print(digits = 3)
#> Call:corCi(x = x, keys = keys, n.iter = n.iter, p = p, overlap = overlap,
#> poly = poly, method = method, plot = plot, minlength = minlength,
#> n = n)
#>
#> Coefficients and bootstrapped confidence intervals
#> help emotn
#> help 1.00
#> emotion 0.33 1.00
#>
#> scale correlations and bootstrapped confidence intervals
#> lower.emp lower.norm estimate upper.norm upper.emp p
#> help-emotn 0.187 0.188 0.326 0.453 0.443 0Korelace mezi EES a ATSPPH byla středně silná a signifikantní, r(170) = 0,33 (95% CI [0,20; 0,44]), p < 0,001. Výsledky tedy podporují hypotézu, že lidé s otevřenějším projevem emocí jsou ochotnější vyhledat odbornou psychologickou pomoc.
Hypotéza 4
K ověření H4 nejdříve odhadneme model s pohlavím (ženy jsou referenční kategorií) a věkem jako prediktory ATSPPH.
fit1 <- lm(help ~ gender + age, data = df)
broom::glance(fit1)
#> # A tibble: 1 × 12
#> r.squared adj.r.squared sigma statistic p.value df logLik AIC BIC
#> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 0.0777 0.0613 0.528 4.72 0.00345 3 -132. 274. 290.
#> # ℹ 3 more variables: deviance <dbl>, df.residual <int>, nobs <int>
broom::tidy(fit1, conf.int = TRUE) %>%
mutate(across(where(is.numeric), ~round(.x, digits = 3)))
#> # A tibble: 4 × 7
#> term estimate std.error statistic p.value conf.low conf.high
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 (Intercept) 3.06 0.106 28.9 0 2.85 3.27
#> 2 genderMuž -0.371 0.1 -3.7 0 -0.57 -0.173
#> 3 genderJiné 0.094 0.31 0.304 0.761 -0.517 0.705
#> 4 age 0.001 0.003 0.195 0.845 -0.005 0.007Jak můžeme vidět výše, tento model byl signifikantní, F(3; 168) = 4,72, p = 0,003 a vysvětloval přibližně 8 % rozptylu. Jediným významným prediktorem zde bylo pohlaví: mužům model predikoval o 0,37 (95% CI [0,17; 0,57]) nižší hodnotu ATSPPH než ženám.
V druhém kroku jsme jako prediktor přidali skóry EES.
fit2 <- lm(help ~ gender + age + emotion, data = df)
broom::glance(fit2)
#> # A tibble: 1 × 12
#> r.squared adj.r.squared sigma statistic p.value df logLik AIC BIC
#> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 0.161 0.140 0.505 7.99 0.00000647 4 -124. 260. 279.
#> # ℹ 3 more variables: deviance <dbl>, df.residual <int>, nobs <int>
broom::tidy(fit2, conf.int = TRUE) %>%
mutate(across(where(is.numeric), ~round(.x, digits = 3)))
#> # A tibble: 5 × 7
#> term estimate std.error statistic p.value conf.low conf.high
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 (Intercept) 2.44 0.184 13.2 0 2.07 2.80
#> 2 genderMuž -0.319 0.097 -3.29 0.001 -0.51 -0.128
#> 3 genderJiné -0.007 0.297 -0.023 0.981 -0.594 0.58
#> 4 age 0 0.003 -0.102 0.919 -0.006 0.005
#> 5 emotion 0.178 0.044 4.06 0 0.092 0.265
anova(fit1, fit2)
#> Analysis of Variance Table
#>
#> Model 1: help ~ gender + age
#> Model 2: help ~ gender + age + emotion
#> Res.Df RSS Df Sum of Sq F Pr(>F)
#> 1 168 46.802
#> 2 167 42.598 1 4.204 16.481 7.537e-05 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Jak můžeme vidět výše, po přidání EES jako prediktoru se podíl vysvětleného rozptylu zvýšil přibližně na 16 % a tento přírůstek byl statisticky významný, F(1, 167) = 16,48, p < 0,001. V souladu s očekáváním emoční otevřenost kladně souvisela s ochotou vyhledat odbornou psychologickou pomoc. Zvýšení EES o jeden bod bylo spojeno s průměrným přírůstkem ATSPPH přibližně o 0,18 bodu (95% CI [0,09; 0,26]).
Tabulka níže pro lepší přehled zobrazuje b-koeficienty a semiparciální korelace koeficientů obou modelů, ale také podíly vysvětleného rozptylu a přírůstek vysvětleného rozptylu společně s 95% intervaly spolehlivosti.
apaTables::apa.reg.table(fit1, fit2)
#>
#>
#> Regression results using help as the criterion
#>
#>
#> Predictor b b_95%_CI sr2 sr2_95%_CI Fit
#> (Intercept) 3.06** [2.85, 3.27]
#> genderMuž -0.37** [-0.57, -0.17] .08 [-.00, .15]
#> genderJiné 0.09 [-0.52, 0.71] .00 [-.01, .01]
#> age 0.00 [-0.01, 0.01] .00 [-.00, .00]
#> R2 = .078**
#> 95% CI[.01,.15]
#>
#> (Intercept) 2.44** [2.07, 2.80]
#> genderMuž -0.32** [-0.51, -0.13] .05 [-.01, .12]
#> genderJiné -0.01 [-0.59, 0.58] .00 [-.00, .00]
#> age -0.00 [-0.01, 0.01] .00 [-.00, .00]
#> emotion 0.18** [0.09, 0.26] .08 [.01, .16]
#> R2 = .161**
#> 95% CI[.06,.24]
#>
#> Difference
#>
#>
#>
#>
#>
#>
#>
#>
#>
#>
#>
#>
#> Delta R2 = .083**
#> 95% CI[.01, .16]
#>
#>
#> Note. A significant b-weight indicates the semi-partial correlation is also significant.
#> b represents unstandardized regression weights.
#> sr2 represents the semi-partial correlation squared.
#> Square brackets are used to enclose the lower and upper limits of a confidence interval.
#> * indicates p < .05. ** indicates p < .01.
#> K ověření předpokladu linearity můžeme použít tzv. ceres grafy. Jak můžeme vidět, vztah mezi EES a ATSPPH byl lineární, zatímco linearitu vztahu mezi věkem a ATSPPH mírně narušovalo to, čeho už jsme si všimli dříve: u starších lidí (nad 25 let) se skór ATSPPH příliš neměnil v závislosti na věku, zatímco u osob do 25 let skór ATSPPH s věkem mírně vzrůstal.
car::ceresPlots(fit2, layout = c(2, 1))
S multikolinearitou náš model problém neměl, soudě podle hodnoty Variance Inflation Factor, které se odchylují od 1 jen zanedbatelně, jak můžeme vidět níže.
car::vif(fit2)
#> GVIF Df GVIF^(1/(2*Df))
#> gender 1.039952 2 1.009842
#> age 1.018452 1 1.009184
#> emotion 1.032649 1 1.016193Nezávislost mezi sousedními rezidui jsme ověřili pomocí Durbin-Watsonova testu. Jak můžeme vidět níže, korelace mezi sousedními rezidui byla mizivá.
car::durbinWatsonTest(fit2)
#> lag Autocorrelation D-W Statistic p-value
#> 1 -0.08299901 2.145973 0.354
#> Alternative hypothesis: rho != 0Co se týče normality reziduí, jak ukazuje histogram a Q-Q graf níže, jejich rozdělení bylo jen mírně negativně zešikmené.
fit2 %>%
resid() %>%
scale() %>%
hist(main = NULL)
car::qqPlot(fit2)
#> [1] 111 155
Graf reziduí a predikovaných hodnot (viz níže) pak nenasvědčoval tomu, že rozptyl reziduí je konstantní a nezávisí na predikované hodnotě.
plot(fit2, 1)
Jak ukazuje graf níže, některé případy vykazovaly vysoké reziduum a některé vysokou hodnotu leverage, ale žádné nevykazovaly kombinaci vysokého rezidua a vysoké hodnoty leverage. Případy napravo (s čísly 49, 50 a 26) jsou jedinými zástupci své kategorie (tři osoby s nebinárním pohlavím), a proto mají zákonitě vysokou hodnotu leverage.
df_augment <- broom::augment(fit2) %>%
mutate(id = row_number())
df_augment %>%
ggplot(aes(.hat, .std.resid, size = .cooksd)) +
geom_text(aes(label = id))
I kdybychom vyřadili extrémní a vlivné případy (26, 49, 50 a 155), regresní analýza by vedla ke stejným výsledkům, jak můžeme vidět níže.
df_augment %>%
filter(.hat < 0.3, abs(.std.resid) < 3) %>%
lm(help ~ gender + age + emotion, data = .) %>%
summary()
#>
#> Call:
#> lm(formula = help ~ gender + age + emotion, data = .)
#>
#> Residuals:
#> Min 1Q Median 3Q Max
#> -1.32810 -0.24443 0.03459 0.36339 1.06286
#>
#> Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) 2.505161 0.178617 14.025 < 2e-16 ***
#> genderMuž -0.337308 0.093610 -3.603 0.000416 ***
#> age -0.001137 0.002840 -0.400 0.689329
#> emotion 0.170223 0.042437 4.011 9.16e-05 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> Residual standard error: 0.4872 on 164 degrees of freedom
#> Multiple R-squared: 0.1702, Adjusted R-squared: 0.155
#> F-statistic: 11.21 on 3 and 164 DF, p-value: 9.891e-07Závěry
Výsledky podpořily naše hypotézy, že ženy vykazují vyšší ochotu vyhledat odbornou psychologickou pomoc než muži a že emoční otevřenost kladně souvisí s ochotou vyhledat tuto pomoc, a to i po kontrole efektu věku a pohlaví. Hypotéza, že starší lidé mají nižší ochotu vyhledat odbornou psychologickou pomoc, ovšem podpořena nebyla. Ve skutečnosti jsme pozorovali určitý trend, kdy ochota vyhledat odbornou psychologickou pomoc mírně rostla s věkem u osob přibližně do 25 roku věku, ale u osob nad 25 let byl vztah mezi věkem a ochotou vyhledat psychologickou pomoc zanedbatelný.
Co se týče statistických limitů, při ověřování rozdílu mezi muži a ženami pomocí nezávislého t-testu bylo rozdělení závislé proměnné u mužů a žen mírně negativně zešikmené. Ale vzhledem k tomu, že tvar rozdělení byl u osobu skupin podobný, měly být být výsledky t-testu robustní.
Ačkoli rozdělení reziduí z lineárně regresního modelu bylo mírně zešikmené a v datech se vyskytlo několik extrémních či vlivných případů, výsledky rovněž považujeme za robustní, protože (a) odchylky od normality byly poměrně malé; (b) vzorek zahrnoval více než 100 osob a nároky na normalitu s velikostí vzorku klesají; (c) i po vyřazení extrémních a vlivných případů by regresní analýza vedla k totožným závěrům.
Session info
Pro reprodukovatelnost uvádíme tzv. session info, obsahující informace o R, operačním systému a všech načtených balíčcích.
sessionInfo()
#> R version 4.2.1 (2022-06-23 ucrt)
#> Platform: x86_64-w64-mingw32/x64 (64-bit)
#> Running under: Windows 10 x64 (build 19045)
#>
#> Matrix products: default
#>
#> locale:
#> [1] LC_COLLATE=Czech_Czechia.utf8 LC_CTYPE=Czech_Czechia.utf8
#> [3] LC_MONETARY=Czech_Czechia.utf8 LC_NUMERIC=C
#> [5] LC_TIME=Czech_Czechia.utf8
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] lubridate_1.9.3 forcats_1.0.0 stringr_1.5.1 dplyr_1.1.4
#> [5] purrr_1.0.2 readr_2.1.4 tidyr_1.3.0 tibble_3.2.1
#> [9] ggplot2_3.4.4 tidyverse_2.0.0 apaTables_2.0.8 car_3.1-2
#> [13] carData_3.0-5 psych_2.3.9 haven_2.5.1
#>
#> loaded via a namespace (and not attached):
#> [1] mvtnorm_1.1-3 lattice_0.20-45 zoo_1.8-10 digest_0.6.29
#> [5] utf8_1.2.2 R6_2.5.1 backports_1.4.1 coda_0.19-4
#> [9] evaluate_0.16 pillar_1.9.0 rlang_1.1.1 multcomp_1.4-20
#> [13] rstudioapi_0.14 Matrix_1.6-4 effectsize_0.8.6 rmarkdown_2.23
#> [17] splines_4.2.1 labeling_0.4.2 htmlwidgets_1.6.2 munsell_0.5.0
#> [21] broom_1.0.5 compiler_4.2.1 xfun_0.39 pkgconfig_2.0.3
#> [25] parameters_0.21.1 mnormt_2.1.0 mgcv_1.8-42 htmltools_0.5.6
#> [29] insight_0.19.6 tidyselect_1.2.0 codetools_0.2-18 fansi_1.0.3
#> [33] tzdb_0.3.0 withr_2.5.0 MASS_7.3-60 grid_4.2.1
#> [37] xtable_1.8-4 nlme_3.1-157 jsonlite_1.8.7 gtable_0.3.1
#> [41] lifecycle_1.0.3 magrittr_2.0.3 bayestestR_0.13.1 scales_1.2.1
#> [45] datawizard_0.9.0 estimability_1.4.1 cli_3.6.1 stringi_1.7.8
#> [49] farver_2.1.1 ellipsis_0.3.2 generics_0.1.3 vctrs_0.6.5
#> [53] sandwich_3.0-2 TH.data_1.1-1 tools_4.2.1 glue_1.6.2
#> [57] hms_1.1.2 emmeans_1.8.9 survival_3.3-1 abind_1.4-5
#> [61] parallel_4.2.1 fastmap_1.1.0 yaml_2.3.5 timechange_0.2.0
#> [65] colorspace_2.0-3 MBESS_4.9.2 knitr_1.43Literatura
Fox, J., & Weisberg, S. (2019). An R Companion to Applied Regression. https://socialsciences.mcmaster.ca/jfox/Books/Companion/
Chen, P., Liu, X. J., Wang, X. Q., Yang, B. X., Ruan, J., & Liu, Z. (2020). Attitude toward seeking professional psychological help among community-dwelling population in China. Frontiers in Psychiatry, 11, 417. https://doi.org/10.3389/fpsyt.2020.00417
Komiya, N., Good, G. E., & Sherrod, N. B. (2000). Emotional openness as a predictor of college students’ attitudes toward seeking psychological help. Journal of Counseling Psychology, 47(1), 138–143. https://doi.org/10.1037/0022-0167.47.1.138
Kring, A. M., Smith, D. A., & Neale, J. M. (1994). Individual differences in dispositional expressiveness: Development and validation of the Emotional Expressivity Scale. Journal of Personality and Social Psychology, 66(5), 934–949. https://doi.org/10.1037/0022-3514.66.5.934
Pawar, N., & Varma, D. P. (2023). Gender-role identification, emotional expressivity and attitude towards seeking professional psychological help in middle adulthood. International Journal of Indian Psychology, 11(3). https://doi.org/10.25215/1103.430
Picco, L., Abdin, E., Chong, S. A., Pang, S., Shafie, S., Chua, B. Y., Vaingankar, J. A., Ong, L. P., Tay, J., & Subramaniam, M. (2016). Attitudes toward seeking professional psychological help: Factor structure and socio-demographic predictors. Frontiers in Psychology, 7, 547. https://doi.org/10.3389/fpsyg.2016.00547
Robinson, D., Hayes, A., & Couch, S. (2023). broom: Convert Statistical Objects into Tidy Tibbles. https://CRAN.R-project.org/package=broom
Stanley, D. (2021). apaTables: Create American Psychological Association (APA) Style Tables. https://CRAN.R-project.org/package=apaTables
Wickham, H., Averick, M., Bryan, J., Chang, W., McGowan, L. D., François, R., Grolemund, G., Hayes, A., Henry, L., Hester, J., Kuhn, M., Pedersen, T. L., Miller, E., Bache, S. M., Müller, K., Ooms, J., Robinson, D., Seidel, D. P., Spinu, V., … Yutani, H. (2019). Welcome to the tidyverse. 4, 1686. https://doi.org/10.21105/joss.01686
Wickham, H., Miller, E., & Smith, D. (2022). haven: Import and Export ’SPSS’, ’Stata’ and ’SAS’ Files. https://CRAN.R-project.org/package=haven
William Revelle. (2023). psych: Procedures for Psychological, Psychometric, and Personality Research. https://CRAN.R-project.org/package=psych