Scholars wave 1 - partial data

Initial Analysis

1 Define functions

## Define function that recodes to numeric, but watches out to coercion to not introduce NAs
colstonumeric <- function(df){
  tryCatch({
    df_num <- as.data.frame(
      lapply(df,
             function(x) { as.numeric(as.character(x))})) 
  },warning = function(stop_on_warning) {
    message("Stoped the execution of numeric conversion: ", conditionMessage(stop_on_warning))
  }) 
}
##
## Define function that reverse codes items
ReverseCode <- function(df, tonumeric = FALSE, min = NULL, max = NULL) {
  if(tonumeric) df <- colstonumeric(df)
  df <- (max + min) - df
}
##
## Define function that scores only rows with less than 10% NAs (returns NA if all or above threshold percentage of rows are NA); can reverse code if vector of column indexes and min, max are provided.
ScoreLikert <- function(df, napercent = .1, tonumeric = FALSE, reversecols = NULL, min = NULL, max = NULL, engine = "sum") {
  reverse_list <- list(reversecols = reversecols, min = min, max = max)
  reverse_check <- !sapply(reverse_list, is.null)
  
  # Recode to numeric, but watch out to coercion to not introduce NAs
  colstonumeric <- function(df){
    tryCatch({
      df_num <- as.data.frame(
        lapply(df,
               function(x) { as.numeric(as.character(x))})) 
    },warning = function(stop_on_warning) {
      message("Stoped the execution of numeric conversion: ", conditionMessage(stop_on_warning))
    }) 
  }
  
  if(tonumeric) df <- colstonumeric(df)
  
  if(all(reverse_check)){
    df[ ,reversecols] <- (max + min) - df[ ,reversecols]
  }else if(any(reverse_check)){
    stop("Insuficient info for reversing. Please provide: ", paste(names(reverse_list)[!reverse_check], collapse = ", "))
  }
  
  if(engine == "sum") {
    return(
      ifelse(rowSums(is.na(df)) > ncol(df) * napercent,
             NA,
             rowSums(df, na.rm = TRUE) * NA ^ (rowSums(!is.na(df)) == 0)
      )
    )  
  }
  
  if(engine == "mean") {
    return(
      ifelse(rowMeans(is.na(df)) > ncol(df) * napercent,
             NA,
             rowMeans(df, na.rm = TRUE) * NA ^ (rowSums(!is.na(df)) == 0)
      )       
    )
  }
  
    if(engine == "mean_na") {
      df[is.na(df)] <- 0
      rowMeans(df)
    }
}

quick_score_hist <- function(x, xlab = "X", na.rm = TRUE, cut_into_bins = TRUE) {
  if(na.rm) {x = x[!is.na(x)]}
  brx <- pretty(range(x, na.rm = TRUE), 
    n = nclass.Sturges(x), min.n = 1)
  df <- data.frame(x = x)
  if(cut_into_bins) {df <- data.frame(x = cut(x, brx, ordered_result = TRUE))}
  df %>% 
  # mutate(x = as.factor(as.character(x))) %>%
  count(x) %>%
  mutate(pct = prop.table(n),
         Percent = paste0(round(pct * 100, 2), " %")) %>% 
  ggplot(aes(x = x, y = pct, label = scales::percent(pct))) + 
     geom_bar(stat = "identity") +
     geom_text(aes(label = Percent), vjust = -0.25) +
     scale_y_continuous(labels = scales::percent) +
     ylab("Percentage %") +  xlab(xlab)    
}

2 Read survey structure

folder <- here::here("Rsyntax&data_650")
data_name <- "survey_686732_R_data_file.csv"
script_name <- "survey_686732_R_syntax_file.R"
  
# Check most recent .csv file
last_csv_file <- 
  dir(folder, pattern = ".*csv", full.names = TRUE) %>% 
  file.info() %>%
  dplyr::arrange(dplyr::desc(ctime)) %>%
  dplyr::slice(1) %>%
  row.names()
if(identical(last_csv_file, file.path(folder, data_name))) {
  cat("Most recent .csv is used.")
} else {
  cat("NOT using the most recent .csv!")
}

# -------------------------------------------------------------------------
# Read data
library(limonaid)
library(sticky)  # need this for sticky labels

df <- limonaid::ls_import_data(
  datafile = file.path(folder, data_name),
  scriptfile = file.path(folder, script_name),
  massConvertToNumeric = FALSE
)

df_compl <-
  df %>%
  filter(lastpage == 17)

# -------------------------------------------------------------------------
# Labels to factor levels levels ("label" = question text; "labels" = response options text)
# library(labelled)
# library(sjlabelled)
# sjlabelled::get_labels(df$G01Q59_SQ008, attr.only = TRUE, values = "as.prefix")
# sjlabelled::get_values(df$G01Q59_SQ008)
# sjlabelled::as_label(df$G01Q59_SQ008, prefix = TRUE, keep.labels = TRUE) 
# sjlabelled::as_character(df$G01Q59_SQ008, prefix = TRUE, keep.labels = TRUE)
# labelled::var_label(df$G01Q59_SQ008)
# labelled::to_factor(df$G01Q59_SQ008, levels = "values")

lime_label_recode <- function (x, prefix = FALSE) {
  labels <- attr(x, "labels", exact = TRUE)
  if (is.null(labels)) {
    x
  } else {
    labels <- unname(labels)
    values <- names(attr(x, "labels", exact = TRUE))
    if (prefix) {
      labels <- sprintf("[%s] %s", values, labels)
    }
    # No recoding solution preserve attributes, even with sticky
      x_rec <- c(labels, x)[match(x, c(values, x))]
    attributes(x_rec) <- attributes(x)  # reattach attributes
    x_rec
  }
}
# test_df <- cbind(df$G02Q02_SQ021, lime_label_recode(df$G02Q02_SQ021))
# lime_label_recode(df$G01Q59_SQ008)
# lime_label_recode(df$G04Q05_SQ001)

# -------------------------------------------------------------------------
# Recode using labels
# cols_to_recode <- lapply(df, function(x) {!is.null(attr(x, "labels", exact = TRUE))})
# cols_to_recode <- which(unlist(cols_to_recode))

# df_recoded <- df
# list_recoded <- lapply(df_recoded[, cols_to_recode], lime_label_recode)
# df_recoded[, cols_to_recode] <- as.data.frame(do.call(cbind, list_recoded))

# df_recoded <-
#   df %>%
#   mutate(across(all_of(cols_to_recode), lime_label_recode)) 

df_recoded <-
  df %>%
  mutate(across(everything(), lime_label_recode)) %>%   # some values have same labels: df$G01Q60_SQ006
  mutate(across(where(is.character), function(col) iconv(col, to="UTF-8")))  # encoding: df_recoded$G01Q56

3 Score 3 Questionnaires

# ------------------------------------------------------------------------------
# Define 3 scales
# ------------------------------------------------------------------------------
# ATSPPH - 10 items (likert 0-3) total sum
atspph_idx <- 184:193  # grep("G06Q13", names(df));  df[, grep("G06Q13", names(df), value = TRUE)]
atspph_labs <- unique(lapply(df[, atspph_idx], attr, "labels"))
atspph_rev <- c(2, 4, 8, 9, 10)

atspph_recode <- function(df, rev) {
  df %>%
    mutate(
      across(everything(),
        ~ case_when(
          . == "AO02" ~ 0,
          . == "AO03" ~ 1,
          . == "AO04" ~ 2,
          . == "AO05" ~ 3
        )
      )       
    ) %>%
    mutate(   # here reverse code
      across(rev,
      ~ 3 - .x 
      )
    )
}  # atspph_recode(df_compl[, atspph_idx], atspph_rev)

# FSozU - 6 items (likert 1-5) total mean
fsozu_idx <- 222:227 # grep("G12Q45", names(df)); df[, grep("G12Q45", names(df), value = TRUE)]  
fsozu_labs <- unique(lapply(df[, fsozu_idx], attr, "labels"))  

fsozu_recode <- function(df) {
  df %>%
    mutate(
      across(everything(),
        ~ case_when(
          . == "AO01" ~ 1,
          . == "AO02" ~ 2,
          . == "AO03" ~ 3,
          . == "AO04" ~ 4,
          . == "AO05" ~ 5
        )
      )       
    ) 
}  # fsozu_recode(df_compl[, fsozu_idx])

# PMHSS - 24 items (likert 1-5) subscale sum
pmhss_idx <- 228:251   # grep("G13Q46", names(df)); df[, grep("G13Q46", names(df), value = TRUE)]
pmhss_labs <- unique(lapply(df[, pmhss_idx], attr, "labels"))

pmhss_aware <- c(2, 4, 5, 6, 8, 10, 11, 12)  
pmhss_agree <- c(14, 16, 17, 18, 20, 22, 23, 24) 
pmhss_posit <- c(1, 3, 7, 9, 13, 15, 19, 21)

pmhss_recode <- function(df) {
  df %>%
    mutate(
      across(everything(),
        ~ case_when(
          . == "AO01" ~ 1,
          . == "AO02" ~ 2,
          . == "AO03" ~ 3,
          . == "AO04" ~ 4,
          . == "AO05" ~ 5
        )
      )       
    )
}  # pmhss_recode(df_compl[, pmhss_idx])

# ------------------------------------------------------------------------------
# Recode & Score
df_compl[, atspph_idx] <- atspph_recode(df_compl[, atspph_idx], atspph_rev)
df_compl[, fsozu_idx] <- fsozu_recode(df_compl[, fsozu_idx])
df_compl[, pmhss_idx] <- pmhss_recode(df_compl[, pmhss_idx])

df_compl$help_seek <- ScoreLikert(df_compl[, atspph_idx], napercent = .5, engine = "sum")
df_compl$soc_supp <- ScoreLikert(df_compl[, fsozu_idx], napercent = .5, engine = "mean")

df_compl$aware <- ScoreLikert(df_compl[, pmhss_idx][pmhss_aware], napercent = .5, engine = "sum") 
df_compl$agree <- ScoreLikert(df_compl[, pmhss_idx][pmhss_agree], napercent = .5, engine = "sum")
df_compl$posit <- ScoreLikert(df_compl[, pmhss_idx][pmhss_posit], napercent = .5, engine = "sum")

4 Some analyses on 3 Questionnaires

4.0.1 Just checks

# find_mod(df_scales)
# moderation_model_list #1,2,3,6,7,10,11,12

mod_synth <-
  moderation_model_list %>%
  purrr::pluck("Syntax") %>%
  stringr::str_match("# Regressions\\\n(.*?)\\\n\\\n#") %>%   # string between "# Regressions\n" and "\n\n#"
  as.data.frame() %>%
  dplyr::pull(2) %>% 
  stringr::str_remove_all(fixed("b0*1 + ")) 

mod_tabl <- 
  moderation_model_list %>%
  purrr::pluck("Model")

for(i in seq_len(length(mod_tabl))) {print(mod_synth[i]); print(mod_tabl[[i]])}

4.0.2 Med - just check

# find_med(df_scales)
# mediation_model_list

for(i in seq_len(length(mediation_model_list$MedEs))) {print(mediation_model_list$MedEs[i]); print(mediation_model_list$PathEs[[i]])}

–>

4.0.3 Odd stigma patterns

ggplot(df_scales, aes(aware, agree, color = posit)) +
  geom_smooth(method = "loess", formula = y ~ x, se = TRUE, alpha = 0.1, color = "red", fill = "red") +
  geom_point() +
  scale_colour_distiller(palette = "Blues", direction = 1)

df_scales %>%
  mutate(posit_cat = cut(posit,
    breaks = c(5, 10, 20, 30, 40))
  ) %>%
  ggplot(aes(aware, agree, color = posit_cat)) +
  geom_point() -> plot_stigma1
plotly::ggplotly(plot_stigma1)

ggplot(df_scales, aes(posit, agree, color = aware)) +
  geom_smooth(method = "loess", formula = y ~ x, se = TRUE, alpha = 0.1, color = "red", fill = "red") +
  geom_point() +
  scale_colour_distiller(palette = "Blues", direction = 1)

df_scales %>%
  mutate(aware_cat = cut(posit,
    breaks = c(5, 10, 20, 30, 40))
  ) %>%
  ggplot(aes(posit, agree, color = aware_cat)) +
  geom_point() -> plot_stigma2
plotly::ggplotly(plot_stigma2)

coplot(agree ~ posit | aware, overlap = 0, data = df_scales,
  panel = function(x, y, ...) {
          points(x, y, ...)
          abline(lm(y ~ x), col = "red")}
)

 Missing rows: 23, 26, 28, 29, 30, 42, 76, 80, 134, 136, 137, 138, 139, 142, 182, 185, 188, 224, 230, 237, 271, 284, 285, 333, 378, 403, 405, 408, 417, 419, 429, 431, 432, 467, 476, 579, 580, 587, 612, 613, 615, 632, 640, 653 

4.0.4 Partial correlations stigma (partial everything from everything)

psych::lowerMat(psych::partial.r(df_scales[, c("agree", "aware", "posit")]))

      agree aware posit
agree 1.00             
aware 0.57  1.00       
posit 0.03  0.39  1.00 

4.0.5 Interaction stigma

df_agree <- na.omit(df_scales[, c("agree", "sex", "age", "resid", "aware", "soc_supp", "posit")])

mod_stigma_interac <- lm(agree ~ aware * posit, data = df_agree)
interactions::interact_plot(mod_stigma_interac, pred = posit, modx = aware)

# interactions::sim_slopes(mod_stigma_interac, pred = posit, modx = aware)

4.0.6 Gender diff stigma

ggstatsplot::ggbetweenstats(df_agree, x = sex, y = agree)

ggstatsplot::ggbetweenstats(df_agree, x = sex, y = aware)

ggstatsplot::ggbetweenstats(df_agree, x = sex, y = posit)

4.0.7 Silly model that works smh (0 m, 1 fem)

df_scales %>% 
  mutate(sex = as.numeric(as.factor(sex)) - 1) %>% 
  psych::mediate(posit ~ sex + aware:agree + (aware), data = .)

Mediation/Moderation Analysis 
Call: psych::mediate(y = posit ~ sex + aware:agree + (aware), data = .)

The DV (Y) was  posit . The IV (X) was  sex agree aware*agree . The mediating variable(s) =  aware .

Total effect(c) of  sex  on  posit  =  -2.35   S.E. =  0.5  t  =  -4.71  df=  653   with p =  0.0000031
Direct effect (c') of  sex  on  posit  removing  aware  =  -1.66   S.E. =  0.49  t  =  -3.42  df=  652   with p =  0.00067
Indirect effect (ab) of  sex  on  posit  through  aware   =  -0.69 
Mean bootstrapped indirect effect =  -0.68  with standard error =  0.18  Lower CI =  -1.06    Upper CI =  -0.35

Total effect(c) of  agree  on  posit  =  0.33   S.E. =  0.04  t  =  8.85  df=  653   with p =  0.0000000000000000082
Direct effect (c') of  agree  on  posit  removing  aware  =  0.09   S.E. =  0.05  t  =  1.91  df=  652   with p =  0.056
Indirect effect (ab) of  agree  on  posit  through  aware   =  0.24 
Mean bootstrapped indirect effect =  0.24  with standard error =  0.04  Lower CI =  0.16    Upper CI =  0.32

Total effect(c) of  aware*agree  on  posit  =  -0.03   S.E. =  0  t  =  -6.89  df=  653   with p =  0.000000000013
Direct effect (c') of  aware*agree  on  posit  removing  aware  =  -0.02   S.E. =  0  t  =  -3.6  df=  652   with p =  0.00034
Indirect effect (ab) of  aware*agree  on  posit  through  aware   =  -0.01 
Mean bootstrapped indirect effect =  -0.01  with standard error =  0  Lower CI =  -0.02    Upper CI =  -0.01
R = 0.52 R2 = 0.27   F = 59.98 on 4 and 652 DF   p-value:  0.00000000000000000000000000000000000000000000000000213 

 To see the longer output, specify short = FALSE in the print statement or ask for the summary

# psych::mediate(agree ~ posit + (aware), data = df_scales) # silly but works
# psych::mediate(agree ~ posit * aware, data = df_scales)

4.0.8 Reg step - everything about stigma is wacky

mod_agree <- lm(agree ~ sex + age + resid + aware + soc_supp, data = df_agree)
best_mod_agree <- step(mod_agree, scope = help_seek ~ .^2, direction = "both", data = mod_agree$model, trace = 0) # BIC with k = log(nrow(mod_agree$model))
summary(best_mod_agree)

Call:
lm(formula = agree ~ sex + aware + soc_supp + sex:aware, data = df_agree)

Residuals:
     Min       1Q   Median       3Q      Max 
-22.8090  -2.9789   0.4495   3.1011  15.5731 

Coefficients:
                  Estimate Std. Error t value             Pr(>|t|)    
(Intercept)        6.81884    1.14539   5.953        0.00000000446 ***
sexMasculin       -3.60276    1.32908  -2.711               0.0069 ** 
aware              0.48422    0.03922  12.346 < 0.0000000000000002 ***
soc_supp           0.34104    0.18238   1.870               0.0620 .  
sexMasculin:aware  0.23489    0.05816   4.039        0.00006065122 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.129 on 604 degrees of freedom
Multiple R-squared:  0.4256,    Adjusted R-squared:  0.4218 
F-statistic: 111.9 on 4 and 604 DF,  p-value: < 0.00000000000000022

df_helpseek <- na.omit(df_scales[, c("help_seek", "sex", "age", "agree", "aware", "soc_supp", "agree", "posit")])
mod_helpseek <- lm(help_seek ~ sex + age + agree  + aware + soc_supp, data = df_helpseek)
best_mod_helpseek <- step(mod_helpseek, scope = help_seek ~ .^2, direction = "both", data = mod_helpseek$model, trace = 0) # BIC with k = log(nrow(mod_helpseek$model))
summary(best_mod_helpseek)


summary(lm(help_seek ~ sex + age * aware, data = df_scales)) 

4.1 Demographics

demo_sex <- 
  df_scales %>% 
  mutate(sex = as.factor(as.character(sex))) %>%
  group_by(sex) %>%
  dplyr::summarise(counts = n()) %>%
  mutate(prop = round(counts*100/sum(counts), 1),
         lab.ypos = cumsum(prop) - .5*prop,
         Percent = paste0(prop, " %")) 
demo_sex

demo_sex %>%
  ggpubr::ggpie(x = "prop", label = "Percent",
                fill = "sex", color = "white", 
                lab.pos = "in", lab.font = list(color = "white"),
                palette = "grey")

demo_resid <- 
  df_scales %>% 
  mutate(sex = as.factor(as.character(resid))) %>%
  group_by(resid) %>%
  dplyr::summarise(counts = n()) %>%
  mutate(prop = round(counts*100/sum(counts), 1),
         lab.ypos = cumsum(prop) - .5*prop,
         Percent = paste0(prop, " %")) 
demo_resid

demo_resid %>%
  ggpubr::ggpie(x = "prop", label = "Percent",
                fill = "resid", color = "white", 
                lab.pos = "in", lab.font = list(color = "white"),
                palette = "grey")

demo_residsex <- 
  df_scales %>% 
  mutate(sex = as.factor(as.character(sex))) %>%
  mutate(resid = as.factor(as.character(resid))) %>%
  dplyr::count(resid, sex, .drop = FALSE) %>%         # Group by, then count number in each group (dont drop 0 counts)
  mutate(pct = prop.table(n),                            # Calculate percent within each var
         Percent = paste0(round(pct * 100, 2), " %"))
demo_residsex

df_scales %>% 
  mutate(sex = as.factor(as.character(sex))) %>%
  mutate(resid = as.factor(as.character(resid))) %>%
  ggstatsplot::grouped_ggpiestats(
    x = sex,
    grouping.var = resid,
    package = "RColorBrewer",
    palette = "Greys",
    bf.message = FALSE,
    ggplot.component = list(scale_fill_grey())
  )

demo_age <- 
  df_scales %>% 
  mutate(age = as.factor(as.character(age))) %>%
  count(age) %>%
  mutate(pct = prop.table(n),
         Percent = paste0(round(pct * 100, 2), " %"))
demo_age

demo_age %>%
  ggplot(aes(x = age, y = pct, label = scales::percent(pct))) + 
     geom_bar(stat = "identity") +
     geom_text(aes(label = Percent), vjust = -0.25) +
     scale_y_continuous(labels = scales::percent) +
     ylab("Percentage %") +  xlab("")  

NA

demo_agesex <-
  df_scales %>% 
  mutate(sex = as.factor(as.character(sex))) %>%
  mutate(age = as.factor(as.character(age))) %>%
  dplyr::count(age, sex, .drop = FALSE) %>%         # Group by, then count number in each group (dont drop 0 counts)
  mutate(pct = prop.table(n),                            # Calculate percent within each var
         Percent = paste0(round(pct * 100, 2), " %"))
demo_agesex

demo_agesex %>%     
  ggplot(aes(x = age, y = pct, fill = sex, label = scales::percent(pct))) + 
      geom_col(position = position_dodge(preserve = "single"), stat = "identity",) +    # Don't drop zero count
      geom_text(position = position_dodge(width = .9),      # move to center of bars
                vjust = -0.5,                               # nudge above top of bar
                size = 3) + 
      scale_y_continuous(labels = scales::percent) +
      ggtitle("") +
      xlab("Varsta") + ylab("Percentage %") + 
      guides(fill = guide_legend(title = "Gen", ncol = 1)) + 
      scale_fill_grey(start = 0.8, end = 0.2, na.value = "red", aesthetics = "fill") +
      theme(legend.position = "right", legend.direction = "vertical", 
            legend.justification = c(0, 1), panel.border = element_rect(fill = NA, colour = "black"))

4.2 Mental Health scales scoring

# ------------------------------------------------------------------------------
# Mental Health scales
# ------------------------------------------------------------------------------

# Screening  DSM-5-TR 11-17 (22 items likert 0-4; 3 items 1/0; 1 item 1/0)
# https://www.psychiatry.org/getmedia/9352851c-d69f-411a-8933-3212e8c29063/APA-DSM5TR-Level1MeasureChildAge11To17.pdf
screen_1_idx <- 59:80 # grep("G02Q02", names(df)); df[, grep("G02Q02", names(df), value = TRUE)]  
screen_1_labs <- unique(lapply(df[, screen_1_idx], attr, "labels"))  
screen_2_idx <- 81:84 # grep("G02Q47|G02Q48", names(df)); df[, grep("G02Q47|G02Q48", names(df), value = TRUE)];   
screen_2_labs <- unique(lapply(df[, screen_2_idx], attr, "labels"))   # 84 is suicide item

# df[, screen_1_idx] %>% map(~ attr(.x, "label"))  # get item text

screen_1_recode <- function(df) {
  df %>%
    mutate(
      across(everything(),
             ~ case_when(
               . == "AO01" ~ 0,
               . == "AO02" ~ 1,
               . == "AO03" ~ 2,
               . == "AO04" ~ 3,
               . == "AO05" ~ 4
             )
      )       
    ) 
}  # screen_1_recode(df_compl[, screen_1_idx])

screen_2_recode <- function(df) {
  df %>%
    mutate(
      across(everything(),
             ~ case_when(
               . == "AO01" ~ 1,
               . == "AO02" ~ 0
             )
      )       
    ) 
}  # screen_2_recode(df_compl[, screen_2_idx])



# ------------------------------------------------------------------------------
# ADHD - 8 items (likert 5-1) total sum;
adhd_idx <- 114:120 # grep("G03Q03", names(df)); df[, grep("G03Q03", names(df), value = TRUE)]  
adhd_labs <- unique(lapply(df[, adhd_idx], attr, "labels"))  

adhd_recode <- function(df) {
  df %>%
    mutate(
      across(everything(),
             ~ case_when(
               . == "AO01" ~ 5,
               . == "AO02" ~ 4,
               . == "AO03" ~ 3,
               . == "AO04" ~ 2,
               . == "AO05" ~ 1
             )
      )       
    ) 
}  # adhd_recode(df_compl[, adhd_idx])

# ------------------------------------------------------------------------------
# Depresie (PHQ-9) – v. versiunea modificata pentru adolescenți (likert 0-3) total sum; cutoff 11
# https://www.childrenshospital.org/sites/default/files/2022-03/PHQ%20Form.pdf
phq_idx <- 85:94 # grep("G09Q39", names(df)); df[, grep("G09Q39", names(df), value = TRUE)]  
phq_labs <- unique(lapply(df[, phq_idx], attr, "labels"))  

phq_recode <- function(df) {
  df %>%
    mutate(
      across(everything(),
             ~ case_when(
               . == "AO02" ~ 0,
               . == "AO03" ~ 1,
               . == "AO04" ~ 2,
               . == "AO05" ~ 3
             )
      )       
    ) 
}  # phq_recode(df_compl[, phq_idx])

# ------------------------------------------------------------------------------
# Anxietate (GAD-7) (likert 0-3) total sum; cutoff 10
gad_idx <- 96:102 # grep("G10Q41", names(df)); df[, grep("G10Q41", names(df), value = TRUE)]
gad_labs <- unique(lapply(df[, gad_idx], attr, "labels"))  

gad_recode <- function(df) {
  df %>%
    mutate(
      across(everything(),
             ~ case_when(
               . == "AO02" ~ 0,
               . == "AO03" ~ 1,
               . == "AO04" ~ 2,
               . == "AO05" ~ 3
             )
      )       
    ) 
}  # gad_recode(df_compl[, gad_idx])

# ------------------------------------------------------------------------------
# Anxietate socială SMSAD 11-17 (likert 0-4) total sum; cutoff 20
# https://www.psychiatry.org/File%20Library/Psychiatrists/Practice/DSM/APA_DSM5_Severity-Measure-For-Social-Anxiety-Disorder-Child-Age-11-to-17.pdf
smsad_idx <- 104:113 # grep("G11Q42", names(df)); df[, grep("G11Q42", names(df), value = TRUE)]; 
smsad_labs <- unique(lapply(df[, smsad_idx], attr, "labels"))  

smsad_recode <- function(df) {
  df %>%
    mutate(
      across(everything(),
             ~ case_when(
               . == "AO02" ~ 0,
               . == "AO03" ~ 1,
               . == "AO04" ~ 2,
               . == "AO05" ~ 3,
               . == "AO06" ~ 4
             )
      )       
    ) 
}  # smsad_recode(df_compl[, smsad_idx])

# ------------------------------------------------------------------------------
# Tulburări de alimentație (NEDA) (slider 0-100) total sum - also has likert items but didnt consider them here 
neda_part1_idx <- 201:205 # grep("G08Q36", names(df)); df[, grep("G08Q36", names(df), value = TRUE)]; 
neda_part1_labs <- unique(lapply(df[, neda_part1_idx], attr, "labels"))  
# --- complicated scoring

# ------------------------------------------------------------------------------
# ACE (likert ) part 1: 10 items, part 2: 9 items; total sum
ace_part1_idx <- 163:172 # grep("G04Q05", names(df)); df[, grep("G04Q05", names(df), value = TRUE)]; 
ace_part2_idx <- 174:182 # grep("G04Q07", names(df)); df[, grep("G04Q07", names(df), value = TRUE)];
ace_idx <- c(ace_part1_idx, ace_part2_idx)
ace_labs <- unique(lapply(df[, ace_idx], attr, "labels")) 




# ------------------------------------------------------------------------------
# Recode & Score
df_compl[, screen_1_idx] <- screen_1_recode(df_compl[, screen_1_idx])
df_compl[, screen_2_idx] <- screen_2_recode(df_compl[, screen_2_idx])

df_compl[, adhd_idx] <- adhd_recode(df_compl[, adhd_idx])
df_compl$adhd <- ScoreLikert(df_compl[, adhd_idx], napercent = .5, engine = "sum")

df_compl[, phq_idx] <- phq_recode(df_compl[, phq_idx])
df_compl$phq <- ScoreLikert(df_compl[, phq_idx], napercent = .5, engine = "sum")

df_compl[, gad_idx] <- gad_recode(df_compl[, gad_idx])
df_compl$gad <- ScoreLikert(df_compl[, gad_idx], napercent = .5, engine = "sum")

df_compl[, smsad_idx] <- smsad_recode(df_compl[, smsad_idx])
df_compl$smsad <- ScoreLikert(df_compl[, smsad_idx], napercent = .5, engine = "sum")

df_compl$ace_part1 <- ScoreLikert(df_compl[, ace_part1_idx], napercent = .5, engine = "sum")
df_compl$ace_part2 <- ScoreLikert(df_compl[, ace_part2_idx], napercent = .5, engine = "sum")
df_compl$ace <- ScoreLikert(df_compl[, ace_idx], napercent = .5, engine = "sum")

# df_compl$neda_part1 <- ScoreLikert(df_compl[, neda_part1_idx], napercent = .5, engine = "sum")

4.3 Mental Health descriptives

4.3.1 Screening

df[, screen_1_idx] %>%     # map(~ attr(.x, "label"))  # get item text
  setNames(
    purrr::map(., ~ attr(.x, "label")) %>% 
      str_replace("În ultimele DOUĂ.*$", "") %>%  
      str_wrap(60)
  ) %>% 
  screen_1_recode() %>% 
  mutate_all(as.factor) %>% 
  likert() %>% 
  plot() +
  theme(# text = element_text(size = 20),
        axis.text.y = element_text(size = 10),
        legend.position = "none") 

df[, screen_2_idx] %>%     # map(~ attr(.x, "label"))  # get item text
  setNames(
    purrr::map(., ~ attr(.x, "label")) %>% 
      str_replace("În ultimele DOUĂ.*$", "") %>%  
      str_wrap(60)
  ) %>% 
  screen_2_recode() %>% 
  mutate_all(as.factor) %>% 
  likert() %>% 
  plot() +
  theme(# text = element_text(size = 20),
        axis.text.y = element_text(size = 10),
        legend.position = "none") 

4.3.2 Same plot without Suicide

df[, screen_2_idx] %>%     # map(~ attr(.x, "label"))  # get item text
  select(-G02Q48) %>%  # remove the suicide item
  setNames(
    purrr::map(., ~ attr(.x, "label")) %>% 
      str_replace("În ultimele DOUĂ.*$", "") %>%  
      str_wrap(60)
  ) %>% 
  screen_2_recode() %>% 
  mutate_all(as.factor) %>% 
  likert() %>% 
  plot() +
  theme(# text = element_text(size = 20),
        axis.text.y = element_text(size = 10),
        legend.position = "none") 

4.3.3 ADHD

df_compl$adhd %>% 
  quick_score_hist("ADHD")

4.3.4 PHQ-9

df_compl$phq %>% 
  quick_score_hist(xlab = "Depresie", na.rm = TRUE) + 
  geom_vline(xintercept = "(10,12]", color = "red")  # cutoff 11

4.3.5 GAD-7

df_compl$gad %>% 
  quick_score_hist(xlab = "Depresie", na.rm = TRUE) + 
  geom_vline(xintercept = "(8,10]", color = "red")  # cutoff 10

4.3.6 SMSAD

df_compl$smsad %>% 
  quick_score_hist(xlab = "Depresie", na.rm = TRUE) + 
  geom_vline(xintercept = "(20,25]", color = "red")  # cutoff 20

5 Session Info

R version 4.3.2 (2023-10-31 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 11 x64 (build 22621)

Matrix products: default


locale:
[1] LC_COLLATE=English_United Kingdom.utf8  LC_CTYPE=English_United Kingdom.utf8    LC_MONETARY=English_United Kingdom.utf8
[4] LC_NUMERIC=C                            LC_TIME=English_United Kingdom.utf8    

time zone: Europe/Bucharest
tzcode source: internal

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] likert_1.3.5       xtable_1.8-4       sticky_0.5.6.1     limonaid_0.1.5     gtsummary_1.7.2    report_0.5.7       scales_1.3.0      
 [8] ggstatsplot_0.12.5 psych_2.3.9        rio_1.0.1          conflicted_1.2.0   fs_1.6.3           papaja_0.1.2       tinylabels_0.2.4  
[15] lubridate_1.9.3    forcats_1.0.0      stringr_1.5.1      dplyr_1.1.4        purrr_1.0.2        readr_2.1.4        tidyr_1.3.0       
[22] tibble_3.2.1       ggplot2_3.5.1      tidyverse_2.0.0    pacman_0.5.1      

loaded via a namespace (and not attached):
  [1] RColorBrewer_1.1-3         rstudioapi_0.15.0          jsonlite_1.8.7             datawizard_0.13.0          correlation_0.8.6         
  [6] magrittr_2.0.3             TH.data_1.1-2              estimability_1.4.1         rmarkdown_2.25             farver_2.1.1              
 [11] vctrs_0.6.4                memoise_2.0.1              paletteer_1.6.0            base64enc_0.1-3            effectsize_0.8.9          
 [16] rstatix_0.7.2              htmltools_0.5.7            broom_1.0.5                Formula_1.2-5              sass_0.4.7                
 [21] parallelly_1.37.1          bslib_0.5.1                htmlwidgets_1.6.2          interactions_1.2.0         plyr_1.8.9                
 [26] sandwich_3.1-0             emmeans_1.9.0              plotly_4.10.4              zoo_1.8-12                 cachem_1.0.8              
 [31] gt_0.10.0                  lifecycle_1.0.4            pkgconfig_2.0.3            Matrix_1.6-5               R6_2.5.1                  
 [36] fastmap_1.1.1              future_1.33.1              BayesFactor_0.9.12-4.7     digest_0.6.33              colorspace_2.1-0          
 [41] GGally_2.2.1               furrr_0.3.1                rematch2_2.1.2             patchwork_1.3.0            rprojroot_2.0.4           
 [46] prismatic_1.1.1            Hmisc_5.1-1                crosstalk_1.2.0            ggpubr_0.6.0               labeling_0.4.3            
 [51] fansi_1.0.5                timechange_0.2.0           abind_1.4-5                httr_1.4.7                 mgcv_1.9-0                
 [56] compiler_4.3.2             here_1.0.1                 withr_3.0.2                pander_0.6.5               htmlTable_2.4.2           
 [61] backports_1.4.1            carData_3.0-5              ggstats_0.6.0              R.utils_2.12.3             ggsignif_0.6.4            
 [66] broom.mixed_0.2.9.4        MASS_7.3-60                PerformanceAnalytics_2.0.4 tools_4.3.2                foreign_0.8-85            
 [71] statsExpressions_1.6.1     nnet_7.3-19                R.oo_1.25.0                glue_1.8.0                 quadprog_1.5-8            
 [76] nlme_3.1-163               grid_4.3.2                 checkmate_2.3.1            cluster_2.1.4              reshape2_1.4.4            
 [81] generics_0.1.3             gtable_0.3.4               tzdb_0.4.0                 R.methodsS3_1.8.2          data.table_1.14.8         
 [86] hms_1.1.3                  car_3.1-2                  xml2_1.3.5                 utf8_1.2.4                 ggrepel_0.9.6             
 [91] pillar_1.9.0               splines_4.3.2              lattice_0.21-9             survival_3.5-7             tidyselect_1.2.0          
 [96] pbapply_1.7-2              knitr_1.45                 gridExtra_2.3              xfun_0.41                  stringi_1.8.1             
[101] lazyeval_0.2.2             yaml_2.3.7                 evaluate_0.23              codetools_0.2-19           cli_3.6.1                 
[106] rpart_4.1.21               parameters_0.23.0          jquerylib_0.1.4            munsell_0.5.0              Rcpp_1.0.11               
[111] globals_0.16.2             zeallot_0.1.0              coda_0.19-4                parallel_4.3.2             MatrixModels_0.5-3        
[116] ellipsis_0.3.2             bayestestR_0.15.0          listenv_0.9.1              broom.helpers_1.14.0       viridisLite_0.4.2         
[121] mvtnorm_1.2-3              xts_0.13.1                 insight_0.20.5             crayon_1.5.2               rlang_1.1.4               
[126] multcomp_1.4-25            mnormt_2.1.1               jtools_2.3.0              

 

---

A work by Claudiu Papasteri

claudiu.papasteri@gmail.com