If you have any questions or queries, please reach me out at luisfca@puc-rio.br
last updated: 02 April, 2022
pacman::p_load(tidyverse, janitor, arsenal, DT, DataExplorer,summarytools, psych, lavaan, mirt)load("C:/Users/luisf/Dropbox/Puc-Rio/Projeto - ASQ 4 2021/Base - ASQ 4.RData")Data name: data_final_merged
ds_final_merged %>%
mutate(race = case_when(
race == 2 ~ "White",
race == 5 ~ "Black / African America",
TRUE ~ "Other races"
),
gender = if_else(gender == 3, NA_character_,gender)
) %>%
tableby(quest~gender + race + age_9,
data = .,
test=FALSE) %>% summary() ds_final_merged %>%
filter(quest == 4) %>%
select(com_a3_1, com_a3_2, com_a3_3, com_a3_4, com_a3_5, com_a3_6) %>%
alpha()ds_final_merged %>%
filter(quest == 4) %>%
select(com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6) %>%
alpha()cronbach_asq_3 <- ds_final_merged %>%
select(quest, com_a3_1, com_a3_2, com_a3_3, com_a3_4, com_a3_5, com_a3_6,
gm_a3_1, gm_a3_2, gm_a3_3, gm_a3_4, gm_a3_5, gm_a3_6,
fm_a3_1, fm_a3_2, fm_a3_3, fm_a3_4, fm_a3_5, fm_a3_6,
cg_a3_1, cg_a3_2, cg_a3_3, cg_a3_4, cg_a3_5, cg_a3_6,
ps_a3_1, ps_a3_2, ps_a3_3, ps_a3_4, ps_a3_5, ps_a3_6) %>%
pivot_longer(cols = -quest,
names_to = c("scale", ".value"),
names_pattern = "(\\w+_\\w+_)(.)") %>%
nest_by(quest, scale) %>%
mutate(alpha(data)$total)cronbach_asq_4 <- ds_final_merged %>%
select(quest, com_a4_1, com_a4_2,com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2,gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2,fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6) %>%
pivot_longer(cols = -quest,
names_to = c("scale", ".value"),
names_pattern = "(\\w+_\\w+_)(.)") %>%
nest_by(quest, scale) %>%
mutate(alpha(data)$total)ds_final_merged %>%
filter(quest == 2) %>%
select(ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6) %>%
alpha(.)ds_final_merged %>%
filter(quest == 8) %>%
select(ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6) %>%
alpha(.)ds_final_merged %>%
filter(quest == 12) %>%
select(com_a3_1, com_a3_2, com_a3_3, com_a3_4, com_a3_5, com_a3_6) %>%
alpha(.)cronbach_asq_4 %>%
group_by(quest, scale) %>%
summarise(alpha=mean(raw_alpha), cor=average_r)
left_join( #just to add sample size
bind_cols(
cronbach_asq_3 %>%
select(quest, scale, raw_alpha, average_r), #get results of asq3
cronbach_asq_4 %>%
select(quest, scale, raw_alpha, average_r) #get results of asq4
) %>%
as.data.frame() %>% #transform into dataframe
janitor::clean_names() %>%
rename(asq3 = raw_alpha_3,
asq4 = raw_alpha_7,
cor3 = average_r_4,
cor4 = average_r_8) %>% #rename to make easier
select(-scale_6, -quest_5) %>%
pivot_longer(cols = -c(quest_1, scale_2)) %>% #transpose
mutate(name = case_when(
name == "asq3" ~ "alpha_3",
name == "cor3" ~ "cor_3",
name == "asq4" ~ "alpha_4",
name == "cor4" ~ "cor_4")) %>% #rename to make easier
mutate(scale_2 = case_when(
scale_2 == "com_a3_" ~ "Communication",
scale_2 == "gm_a3_" ~ "Gross Motor",
scale_2 == "fm_a3_" ~ "Fine Motor",
scale_2 == "cg_a3_" ~ "Problem Solving",
scale_2 == "ps_a3_" ~ "Personal-Social")) %>%
pivot_wider(names_from = name, values_from = value)#inverse tranpose
#pivot_wider(name, scale_2, values_fn = mean)
,
ds_final_merged %>% count(quest) %>% rename(quest_1=quest) #rename to make easier
) %>%
select(quest_1, scale_2, alpha_3, alpha_4, cor_3, cor_4) %>% #order
mutate(delta_percent = (alpha_4-alpha_3)/alpha_3*100) %>%
mutate_if(is.numeric, round, 2) %>% #round
arrange(desc(delta_percent))bind_cols(
cronbach_asq_3 %>%
select(quest, scale, raw_alpha),
cronbach_asq_4 %>%
select(quest, scale, raw_alpha)) %>%
as.data.frame() %>%
janitor::clean_names() %>%
rename(asq3 = raw_alpha_3,
asq4 = raw_alpha_6) %>%
select(-scale_5, -quest_4) %>%
pivot_longer(cols = -c(quest_1, scale_2))%>%
mutate(scale_2 = case_when(
scale_2 == "com_a3_" ~ "Communication",
scale_2 == "gm_a3_" ~ "Gross Motor",
scale_2 == "fm_a3_" ~ "Fine Motor",
scale_2 == "cg_a3_" ~ "Problem Solving",
scale_2 == "ps_a3_" ~ "Personal-Social")) %>%
ggplot(., aes(x=factor(scale_2, level = c('Communication', 'Gross Motor', 'Fine Motor',
"Problem Solving", "Personal-Social")), y = value, fill = name)) +
geom_col(stat = "summary", position="dodge") +
facet_wrap(~quest_1) +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "Domain", y = "Questionnaire")
ASQ-4
ds_final_merged %>%
filter(quest==9) %>%
select(ps_a4_1,ps_a4_2,ps_a4_3,ps_a4_4,ps_a4_5,ps_a4_6) %>%
alpha() ds_final_merged %>%
filter(quest == 9) %>%
{cor(.$ps_a3_4, .$ps_a4_5)}ASQ-3
ds_final_merged %>%
filter(quest==9) %>%
select(ps_a3_1,ps_a3_2,ps_a3_3,ps_a3_4,ps_a3_5,ps_a3_6) %>%
alpha()ds_final_merged %>%
select(quest, c_sum_a3, gm_sum_a3, fm_sum_a3, cg_sum_a3, ps_sum_a3,
c_sum, gm_sum, fm_sum, cg_sum, ps_sum) %>%
group_by(quest) %>%
summarise(across(contains("sum"),
~mean(.)),
sample = n()) %>%
ungroup() %>%
mutate_if(is.numeric,round,2) %>%
select(quest, sample, order(colnames(.)))ds_final_merged %>%
select(quest, c_sum_a3, gm_sum_a3, fm_sum_a3, cg_sum_a3, ps_sum_a3,
c_sum, gm_sum, fm_sum, cg_sum, ps_sum) %>%
tableone::CreateTableOne(vars = names(.), strata = c("quest"), data = .)
#describeBy(., group = "quest")ds_final_merged %>%
select(quest, contains("sum")) %>%
nest_by(quest) %>%
summarise(categ = combn(names(data), 2, paste, collapse="-"),
pval = combn(data, 2, function(x)
t.test(x[[1]], x[[2]], paired = TRUE)$p.value), .groups = 'drop') %>%
as.data.frame() %>%
separate(categ, into = c("a3","a4"), sep = "-") %>%
mutate(pval = round(pval, 2)) %>%
mutate_at(vars(a3, a4), ~str_remove_all(.,"_a3")) %>%
filter(a3 == a4) %>%
filter(pval < 0.05) %>%
arrange(quest)Looping the ds to perform all pairwise comparisons
ds_final_merged %>%
select(quest, contains("sum")) %>%
pivot_longer(cols = -quest) %>%
group_by(quest) %>%
summarise(pout = list(broom::tidy(pairwise.t.test(value, name,
p.adjust.method = "none", paired = TRUE)))) %>%
unnest(pout) %>%
as.data.frame() %>%
mutate_at(vars(group1, group2), ~str_remove_all(.,"_a3")) %>%
filter(group1 == group2) %>%
arrange(quest) %>%
mutate(p.value = round(p.value,2))ds_final_merged %>%
filter(quest==16) %>%
{t.test(.$gm_sum, .$gm_sum_a3, paired=T, data =.)}ds_final_merged %>%
filter(quest==12) %>%
{t.test(.$ps_sum, .$ps_sum_a3, paired=T, data =.)}library(correlation)
ds_final_merged %>%
group_by(quest) %>%
select(quest, c_sum_a3, gm_sum_a3, fm_sum_a3, cg_sum_a3, ps_sum_a3,
c_sum, gm_sum, fm_sum, cg_sum, ps_sum) %>%
correlation() %>%
as.data.frame() %>%
mutate_at(vars(Parameter1, Parameter2), ~str_remove_all(.,"_a3")) %>%
filter(Parameter1 == Parameter2) ds_final_merged %>%
select(quest, c_sum_a3, gm_sum_a3, fm_sum_a3, cg_sum_a3, ps_sum_a3,
c_sum, gm_sum, fm_sum, cg_sum, ps_sum) %>%
group_by(quest) %>%
nest() %>%
mutate(
correlations = map(data, corrr::correlate)
) %>%
unnest(correlations) %>%
select(-c(contains("_a3"))) %>%
filter(str_detect(term, 'a3')) %>%
select(-data) %>%
pivot_longer(-c(quest,term)) %>%
filter(term == paste0(name,"_a3")) ds %>%
#filter(base == "sup") %>%
select(quest, base,
c_sum, gm_sum, fm_sum, cg_sum, ps_sum) %>%
tableby(interaction(quest, base) ~ ., data = .) %>%
summary()
#tableone::CreateTableOne(vars = names(.[-1]), strata = c("quest","base"), data = .) %>%
#print()ds %>%
select(quest, base, c_sum:ps_sum) %>%
pivot_longer(-c(quest,base)) %>%
mutate(name = case_when(
name == "c_sum" ~ "Communication",
name == "gm_sum" ~ "Gross Motor",
name == "fm_sum" ~ "Fine Motor",
name == "cg_sum" ~ "Problem Solving",
name == "ps_sum" ~ "Personal-Social")) %>%
ggplot(., aes(x = name, y = value, fill = base)) +
geom_bar(stat = "summary", fun = "mean", position = "dodge") +
facet_wrap(~quest) +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1, size=11,face="bold"))ds %>%
filter(quest !=9 & quest != 72) %>%
select(quest, base, c_sum:ps_sum) %>%
pivot_longer(cols = -c(quest, base)) %>%
group_by(quest, name) %>%
nest() %>%
mutate(mean1 = map_dbl(data, ~mean(.x$value[.x$base == "base1"]))) %>%
mutate(mean2 = map_dbl(data, ~mean(.x$value[.x$base == "sup"]))) %>%
mutate(pval = map(data, ~t.test(.x$value ~ .x$base)$p.value)) %>%
unnest(pval) %>%
mutate(pval = round(pval, 2)) %>%
filter(pval > 0.05) %>%
arrange(quest) %>%
select(-data)asq_3_table <- read_csv("C:/Users/luisf/Downloads/asq_3_table.csv")
asq_3_table <- clean_names(asq_3_table)
asq_3_table <- remove_empty(asq_3_table)
asq_3_table <- asq_3_table %>% mutate(quest = as.factor(quest))
asq_3_tableleft_join(
ds %>% #get ds
select(quest, ends_with("sum")) %>% #select the focus variables
pivot_longer(-quest) %>% #transpose to long formata
nest_by(quest, name) %>% #group
mutate(mean = list(map_dbl(data, ~mean(.))),
sd = list(map_dbl(data, ~sd(.))),
m_1sd = mean-sd,
m_1_half_sd = mean-1.5*sd,
m_2sd = mean-2*sd) %>%
unnest(-data) %>% #unnest
pivot_wider(id_cols = quest, names_from = name, values_from = mean:m_2sd) %>%
mutate_if(is.numeric, round, 2) %>%
ungroup()
,
asq_3_table,
by = "quest") %>%
pivot_longer(-quest) %>%
arrange(quest,name) %>%
pivot_wider(id_cols = quest, names_from = name, values_from = value) %>%
rename_all(.,~stringr::str_replace_all(., 'x', 'asq4')) %>%
rename_all(.,~stringr::str_replace_all(., 'y', 'asq3')) %>%
select(quest,contains("m_2sd"))ds %>%
tabyl(quest,website) %>%
adorn_totals(where = c("row","col")) %>%
adorn_percentages(denominator = "col") %>%
adorn_pct_formatting(digits = 0) %>%
adorn_ns(position = "front")
ds %>%
tableby(gender~quest,
data = .,
test=FALSE) %>% summary()ds %>% count(gender) %>% adorn_totals()ds %>%
group_by(quest) %>%
summarise(pval = chisq.test(table(gender))$p.value) %>%
arrange(pval)ds %>% tableby(base~momage_numeric,data = .) %>% summary() ds %>%
tableby(~momed,data = ., test=FALSE) %>% summary()
#%>% xlsx::write.xlsx(., file = "raw_results.xlsx", sheetName="momed", append=TRUE)On Dec 28, 2021
ds %>%
mutate(income = ifelse(as.integer(income)<7, income,NA)) %>% #I changed the order of income in ds sup to compute the risk factor, but I did not changed this variable by itself
tableby(~as.factor(income),data = ., test=FALSE) %>%
summary() On Dec 28, 2021
ds %>%
tableby(~summative_risk,data = ., test=FALSE) %>% summary() ds %>%
select(quest,summative_risk) %>%
mutate(summative_risk_sup = as_factor(summative_risk),
summative_risk_sup = fct_inseq(summative_risk)) %>%
tableone::CreateTableOne(vars = "summative_risk", strata = "quest", data = .) %>%
print(.) %>% t(.) %>% as.data.frame() %>% rownames_to_column("quest")ds %>%
filter(base=="base1") %>%
{gmodels::CrossTable(.$atrisk, .$summative_risk, chisq = T)} ds %>%
tableby(~race,data = ., test=FALSE) %>% summary() ds %>%
filter(!is.na(quest)) %>%
select(quest, ends_with("sum")) %>%
tableby(quest ~ ., control = tableby.control(numeric.stats=c("mean", "sd")), data = .) %>% summary(. , digits = 2)
#tableone::CreateTableOne(vars = names(.), strata = c("quest"), data = .) %>% transpose()
#table1::table1(~ .| quest,
# transpose = TRUE,
# data = .) ds %>%
filter(!is.na(quest)) %>%
select(quest, ends_with("sum")) %>%
pivot_longer(-quest) %>%
mutate_at(vars(name), ~case_when(
. == "ps_sum" ~ "Personal & Social",
. == "gm_sum" ~ "Gross motor",
. == "fm_sum" ~ "Fine motor",
. == "c_sum" ~ "Communication",
. == "cg_sum" ~ "Problem solving",
)) %>%
ggplot(., aes(x = value, y = name, fill = name)) +
ggridges::geom_density_ridges(rel_min_height = 0.01) +
facet_wrap(~quest) +
ggridges::theme_ridges(grid = FALSE, center_axis_labels = TRUE) +
theme(legend.position = "hide") + labs(y="")
reg_com <- "^com_a4_.*"
reg_fm <- "^fm_a4_.*"
reg_gm <- "^gm_a4_.*"
reg_cg <- "^cg_a4_.*"
reg_ps <- "^ps_a4_.*"
regs <- c(reg_fm, reg_com, reg_gm, reg_cg, reg_ps) %>%
set_names(c("fm_a4_", "com_a4_", "gm_a4_", "cg_a4_",
"ps_a4_"))
cronbachs_alpha <-
map_df(regs, ~
ds %>%
select(dplyr::matches(.x)) %>%
psych::alpha(check.keys = TRUE) %>% .$total %>%
tibble::rownames_to_column()
,.id = "scale"
)apply_alpha <- function(data, nest_contains) {
x<-data %>%
select(quest, contains(nest_contains)) %>%
group_by(quest) %>%
do(alpha(.[-1])$total) #compute alpha
y<-data %>%
count(quest) #get n
z <- left_join(x,y, by = "quest")
z <- z %>% mutate_if(is.numeric, round,2)
return(z)
}apply_alpha(ds, 'com_a4_')left_join(
ds %>%
#select items
select(quest, com_a4_1:com_a4_6, gm_a4_1:gm_a4_6, fm_a4_1:fm_a4_6, cg_a4_1:cg_a4_6, ps_a4_1:ps_a4_6) %>%
pivot_longer(cols = -quest,
names_to = c("scale", ".value"),
names_pattern = "(\\w+_\\w+_)(.)") %>%
nest_by(quest, scale) %>%
#Compute cronbach´s alpha for all questionnaires and domains
mutate(alpha(data)$total) %>%
select(quest, scale, std.alpha, average_r)
,
ds %>% count(quest)
) %>%
mutate_if(is.numeric, round,2) %>%
pivot_wider(quest, names_from = scale, values_from = std.alpha:n)#https://stackoverflow.com/questions/69302457/using-dplyr-to-nest-or-group-two-variables-then-perform-the-cronbachs-alpha-fu/69303641#69303641
ds %>%
select(quest, com_a4_1:com_a4_6, gm_a4_1:gm_a4_6, fm_a4_1:fm_a4_6, cg_a4_1:cg_a4_6, ps_a4_1:ps_a4_6) %>%
mutate(id = 1:n()) %>%
pivot_longer(cols = c(-id, -quest)) %>%
separate(col = name,
into = c("scale", "item"),
sep = "_",
extra = "merge") %>%
pivot_wider(names_from = item) %>%
select(-id) %>%
group_by(quest, scale) %>%
nest() %>%
mutate(alpha_results = map(data, ~alpha(.)$total)) %>%
unnest_wider(alpha_results) %>% #get alpha results
select(quest, scale, std.alpha, average_r) %>% #what I want to get
arrange(quest, scale) %>%
pivot_wider(names_from = scale, values_from = std.alpha:average_r) %>%
mutate_if(is.numeric, round, 2) ds %>%
select(quest, com_a4_1:com_a4_6, gm_a4_1:gm_a4_6, fm_a4_1:fm_a4_6, cg_a4_1:cg_a4_6, ps_a4_1:ps_a4_6) %>%
group_by(quest) %>%
do(alpha(.[-1])$total) %>%
select(quest, std.alpha) %>%
mutate_if(is.numeric, round, 2)ds_1 %>%
filter(quest == 16) %>%
select( com_a4_1:com_a4_6, gm_a4_1:gm_a4_6, fm_a4_1:fm_a4_6, cg_a4_1:cg_a4_6, ps_a4_1:ps_a4_6) %>%
psych::alpha(.)#make data
ds_1 %>%
select(quest, com_a4_1:com_a4_6, gm_a4_1:gm_a4_6, fm_a4_1:fm_a4_6, cg_a4_1:cg_a4_6, ps_a4_1:ps_a4_6) %>%
pivot_longer(cols = -quest,
names_to = c("scale", ".value"),
names_pattern = "(\\w+_\\w+_)(.)") %>%
nest_by(quest, scale)%>%
mutate(alpha=alpha(data)$total$raw_alpha) ->x
x %>%
mutate_at(vars(scale), ~case_when(
. == "cg_a4_" ~ "Problem Solving",
. == "com_a4_" ~ "Communication",
. == "fm_a4_" ~ "Fine Motor",
. == "gm_a4_" ~ "Gross Motor",
. == "ps_a4_" ~ "Personal-Social",
)) %>%
ggplot(., aes(x=quest, y = alpha, color = scale)) +
stat_summary(geom = "line", size = 1) +
stat_summary(geom = "point") +
labs(y="Cronbach's Alpha", x="Age-interval") +
theme_bw() +
theme(legend.position = "bottom")ds %>%
select(quest, com_a4_1:com_a4_6, gm_a4_1:gm_a4_6, fm_a4_1:fm_a4_6, cg_a4_1:cg_a4_6, ps_a4_1:ps_a4_6) %>%
pivot_longer(cols = -quest,
names_to = c("scale", ".value"),
names_pattern = "(\\w+_\\w+_)(.)") %>%
nest_by(quest, scale) %>%
mutate(cor=mean(alpha(data)$item.stats$r.cor)) %>%
mutate_if(is.numeric, round,2) %>%
pivot_wider(quest, names_from = scale, values_from = cor)ds_1 %>%
select(quest, fm_a4_1:fm_a4_6) %>%
filter(quest == 4) %>%
{mean(alpha(.)$item.stats$r.cor)}ds_1 %>%
select(quest, com_a4_1:com_a4_6, gm_a4_1:gm_a4_6, fm_a4_1:fm_a4_6, cg_a4_1:cg_a4_6, ps_a4_1:ps_a4_6) %>%
pivot_longer(cols = -quest,
names_to = c("scale", ".value"),
names_pattern = "(\\w+_\\w+_)(.)") %>%
nest_by(quest, scale) %>%
mutate(cor=mean(alpha(data)$item.stats$r.cor))->x
x %>%
mutate_at(vars(scale), ~case_when(
. == "cg_a4_" ~ "Problem Solving",
. == "com_a4_" ~ "Communication",
. == "fm_a4_" ~ "Fine Motor",
. == "gm_a4_" ~ "Gross Motor",
. == "ps_a4_" ~ "Personal-Social",
)) %>%
ggplot(., aes(x=quest, y = cor, color = scale)) +
stat_summary(geom = "line", size = 1) +
stat_summary(geom = "point") +
labs(y="Correlation coefficient", x="Age-interval") +
scale_x_continuous(breaks = seq(2, 72, by = 2)) +
theme_bw() +
theme(legend.position = "bottom")Asked by Jane and Kimberly on March 3, 2022
compare_communication <- function(quest, items) {
quest <- enquo(quest)
items <- enquo(items)
ds %>% #get data
filter(quest == !!quest) %>% #select which questionnaire will be used
mutate(com_exp = rowSums(select(., !!items), na.rm=T)) %>% #create a summative score for expressive items
mutate(com_rec = c_sum-com_exp) %>%
select(com_exp, com_rec, c_sum) %>%
summarise(n=n(),
mean(com_exp),
mean(com_rec),
mean(c_sum),
p = t.test(com_exp, com_rec, alternative = "two.sided", paired = T)$p.value) %>%
t() #compare scores
}
list(
compare_communication(quest = 2, items = c(com_a4_1,com_a4_2,com_a4_5)),
compare_communication(quest = 4 , items = c( com_a4_1,com_a4_4,com_a4_6 )),
compare_communication(quest = 6 , items = c( com_a4_1,com_a4_2,com_a4_5 )),
compare_communication(quest = 8 , items = c( com_a4_3,com_a4_4,com_a4_6 )),
compare_communication(quest = 10 , items = c( com_a4_1,com_a4_3,com_a4_6 )),
compare_communication(quest = 12 , items = c( com_a4_1,com_a4_4,com_a4_6 )),
compare_communication(quest = 14 , items = c( com_a4_1,com_a4_2,com_a4_5 )),
compare_communication(quest = 16 , items = c( com_a4_5,com_a4_3,com_a4_6 )),
compare_communication(quest = 18 , items = c( com_a4_3,com_a4_4,com_a4_6 )),
compare_communication(quest = 20 , items = c( com_a4_2,com_a4_3,com_a4_6 )),
compare_communication(quest = 22 , items = c( com_a4_3,com_a4_5 ,com_a4_6 )),
compare_communication(quest = 24 , items = c( com_a4_3,com_a4_5 ,com_a4_6 )),
compare_communication(quest = 27 , items = c( com_a4_2,com_a4_4,com_a4_5 )),
compare_communication(quest = 30 , items = c( com_a4_2,com_a4_4,com_a4_6 ))
)ds %>% filter(quest == 4) %>%
select(com_a4_1,com_a4_4,com_a4_6) %>%
DataExplorer::profile_missing()
ds %>%
filter(quest == 4) %>%
rowid_to_column() %>%
filter(is.na(com_a4_4)) #its missing because asq3 com 4 was missing and 2 and 4 months are equal
left_join(
ds %>% #get ds
select(quest, ends_with("sum")) %>% #select the focus variables
pivot_longer(-quest) %>% #transpose to long formata
nest_by(quest, name) %>% #group
mutate(mean = list(map_dbl(data, ~mean(.))),
sd = list(map_dbl(data, ~sd(.))),
m_1sd = mean-sd,
m_1_half_sd = mean-1.5*sd,
m_2sd = mean-2*sd) %>%
unnest(-data) %>% #unnest
pivot_wider(id_cols = quest, names_from = name, values_from = mean:m_2sd) %>%
mutate_if(is.numeric, round, 2)
,
ds %>% count(quest)
) %>% select(quest, n, everything())ds %>%
select(quest, ends_with("sum")) %>% #get variable names
pivot_longer(-quest) %>% #tranform into the long format
nest_by(quest, name) %>% #group or nest
mutate(
questionnaire = quest,#compute questionnaire
n = map_dbl(data, ~nrow(data.frame(.))), #compute sample size
mean = map_dbl(data, ~mean(.)), #get the means
sd = map_dbl(data, ~sd(.)), #get sd
m_1sd = mean-sd, #1 below
m_1_half_sd = mean-1.5*sd, #1.5 below
m_2sd = mean-2*sd, #2 below
how_many_c = map(data, ~ ifelse(name == "c_sum", #i'll count how many participants are below
ds %>%
filter(quest == questionnaire) %>%
filter(c_sum > m_2sd & c_sum <= m_1sd) %>%
summarise(n()),
NA_integer_)), # percentage of communication
how_many_gm = map(data, ~ ifelse(name == "gm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(gm_sum > m_2sd & gm_sum <= m_1sd) %>%
summarise(n()), #percentage of gross motor
NA_integer_)),
how_many_fm = map(data, ~ ifelse(name == "fm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(fm_sum > m_2sd & fm_sum <= m_1sd) %>%
summarise(n()),#percentage of fine motor
NA_integer_)),
how_many_cg = map(data, ~ ifelse(name == "cg_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(cg_sum > m_2sd & cg_sum <= m_1sd) %>%
summarise(n()), #percentage of problem-solving
NA_integer_)),
how_many_ps = map(data, ~ ifelse(name == "ps_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(ps_sum > m_2sd & ps_sum <= m_1sd) %>%
summarise(n()), #percentage of personal social
NA_integer_))) %>%
mutate(n_com=
paste0(how_many_c[[1]],"-",
round(how_many_c[[1]]/n*100,0),"%"), #I'll get the percentages
n_gm=
paste0(how_many_gm[[1]],"-",
round(how_many_gm[[1]]/n*100,0),"%"),
n_fm=
paste0(how_many_fm[[1]],"-",
round(how_many_fm[[1]]/n*100,0),"%"),
n_cg=
paste0(how_many_cg[[1]],"-",
round(how_many_cg[[1]]/n*100,0),"%"),
n_ps=
paste0(how_many_ps[[1]],"-",
round(how_many_ps[[1]]/n*100,0),"%")
) %>%
unnest_wider(-quest) %>% #particular trickes (there is a lot of pseudo missing values)
select(quest,how_many_ps:last_col(), -how_many_ps) %>% #get some variables
pivot_longer(-quest) %>% #now transform to the long format
filter(value != "NA-NA%") %>% #filter
pivot_wider(id_cols = quest, names_from = name, values_from = value)
ds %>% filter(quest == 2 & c_sum <= 27.82241 & c_sum > 13.4012)
ds %>%
filter(quest == 4) %>%
filter(gm_sum > 35.2 & gm_sum <= 44.29)ds %>%
filter(quest == 6) %>%
filter(fm_sum > 21.6 & fm_sum <= 34.71)ds %>%
filter(quest == 27) %>%
filter(ps_sum > 22.02 & ps_sum <= 34.88)ds %>%
select(quest, ends_with("sum")) %>% #get variable names
pivot_longer(-quest) %>% #tranform into the long format
nest_by(quest, name) %>% #group or nest
mutate(
questionnaire = quest,#compute questionnaire
n = list(map_dbl(data, ~nrow(data.frame(.)))), #compute sample size
mean = list(map_dbl(data, ~mean(.))), #get the means
sd = list(map_dbl(data, ~sd(.))), #get sd
m_1sd = mean-sd, #1 below
m_1_half_sd = mean-1.5*sd, #1.5 below
m_2sd = mean-2*sd, #2 below
how_many_c = map(data, ~ ifelse(name == "c_sum", #i'll count how many participants are below
ds %>%
filter(quest == questionnaire) %>%
filter(c_sum <= m_2sd) %>%
summarise(n()),
NA_integer_)), # percentage of communication
how_many_gm = map(data, ~ ifelse(name == "gm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(gm_sum <= m_2sd) %>%
summarise(n()), #percentage of gross motor
NA_integer_)),
how_many_fm = map(data, ~ ifelse(name == "fm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(fm_sum <= m_2sd) %>%
summarise(n()),#percentage of fine motor
NA_integer_)),
how_many_cg = map(data, ~ ifelse(name == "cg_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(cg_sum <= m_2sd) %>%
summarise(n()), #percentage of problem-solving
NA_integer_)),
how_many_ps = map(data, ~ ifelse(name == "ps_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(ps_sum <= m_2sd) %>%
summarise(n()), #percentage of personal social
NA_integer_))) %>%
mutate(n_com=
paste0(how_many_c[[1]],"-",
round(how_many_c[[1]]/n*100,0),"%"), #I'll get the percentages
n_gm=
paste0(how_many_gm[[1]],"-",
round(how_many_gm[[1]]/n*100,0),"%"),
n_fm=
paste0(how_many_fm[[1]],"-",
round(how_many_fm[[1]]/n*100,0),"%"),
n_cg=
paste0(how_many_cg[[1]],"-",
round(how_many_cg[[1]]/n*100,0),"%"),
n_ps=
paste0(how_many_ps[[1]],"-",
round(how_many_ps[[1]]/n*100,0),"%")
) %>%
unnest_wider(-quest) %>% #particular trickes (there is a lot of pseudo missing values)
select(quest,how_many_ps:last_col(), -how_many_ps) %>% #get some variables
pivot_longer(-quest) %>% #now transform to the long format
filter(value != "NA-NA%") %>% #filter
pivot_wider(id_cols = quest, names_from = name, values_from = value)ds %>%
filter(quest == 2) %>%
filter(c_sum > 13.4012211 & c_sum <= 27.82241 ) %>%
count()
ds %>%
filter(quest == 4) %>%
filter(c_sum > 30.542897 & c_sum <= 40.38866) %>%
count()
ds %>%
filter(quest == 6) %>%
filter(c_sum > 28.3859573 & c_sum <= 38.04648) %>%
count()
ds %>%
filter(quest == 2) %>%
filter(gm_sum > 38.2619752 & gm_sum <= 46.0505 ) %>%
count()
ds %>%
filter(quest == 36) %>%
filter(cg_sum > 23.3887282 & cg_sum <= 35.95029 ) %>%
count()
ds %>%
filter(quest == 24) %>%
filter(cg_sum < 18.14)ds %>%
filter(quest %in% c(18,20,22)) %>%
group_by(quest) %>%
mutate(czero = ifelse(c_sum <= 0,1,0)) %>%
select(czero, everything()) %>%
summarise(n=n(), sum(czero), sum(czero)/n)
select()ds %>%
select(id, quest, ends_with("sum")) %>%
pivot_longer(-c(id, quest)) %>% #tranform into long format
nest_by(quest, name) %>% #nest (id will not be used this time)
#here is the trick. all counts are here
mutate(
n = map_dbl(data[2], ~nrow(data.frame(.))), #compute sample size
mean = map_dbl(data[2], ~mean(.)), #get the ROBUST means
sd = map_dbl(data[2], ~sd(.)), #get the ROBUST sd
one_below = mean-sd, #1 below
two_below = mean - 2 * sd,
monitor = sum(one_below >= data[[2]] & two_below < data[[2]])/n,
below = sum(two_below > data[[2]])/n) %>%
#outside of the map, add a new variable to nested data
mutate(data = list(data %>% mutate(area_monitor = ifelse(data[[2]] > two_below & data[[2]] <= one_below, paste0(name),0)))) %>%
unnest(data) %>% #I'll unnest and get again the ids. my df will be very long
pivot_wider(id_cols = c(id, quest), names_from = c(area_monitor), values_fn = length) %>% #now I'll group by ids!! my df will be sample size * quest
ungroup() %>% #end of function
#now I'll count values
mutate(how_many = select(.,ends_with("_sum")) %>% rowSums(., na.rm = T)) %>%
#now I'll present the results
tabyl(quest, how_many)%>%
as.data.frame() %>%
# adorn_totals(c("row", "col")) %>%
#adorn_percentages("row") %>%
#adorn_pct_formatting(digits = 2) %>%
#adorn_ns() #
pivot_longer(-quest) %>%
group_by(quest) %>%
mutate(prop = prop.table(value)) %>%
ungroup()%>%
ggplot(., aes(x = name, y = prop, group=quest, color = quest, label=quest)) +
stat_summary(geom = "line", fun = mean, size=1) +
#ggrepel::geom_text_repel(data = . %>% filter(name == "0"), box.padding = 0.5, max.overlaps = Inf) +
labs(x = "Domains at the monitoring zone", y = "Proportion", color = "Questionnaire") +
theme_bw() +
theme(legend.position = "bottom")ds %>%
select(id, quest, ends_with("sum")) %>%
pivot_longer(-c(id, quest)) %>% #tranform into long format
nest_by(quest, name) %>% #nest (id will not be used this time)
#here is the trick. all counts are here
mutate(
n = map_dbl(data[2], ~nrow(data.frame(.))), #compute sample size
mean = map_dbl(data[2], ~mean(.)), #get the ROBUST means
sd = map_dbl(data[2], ~sd(.)), #get the ROBUST sd
one_below = mean-sd, #1 below
two_below = mean - 2 * sd,
monitor = sum(one_below >= data[[2]] & two_below < data[[2]])/n,
below = sum(two_below > data[[2]])/n) %>%
#outside of the map, add a new variable to nested data
mutate(data = list(data %>% mutate(area_monitor = as.character(ifelse(data[[2]] <= two_below, paste0(name),0))))) %>% #attention here
unnest(data) %>% #I'll unnest and get again the ids. my df will be very long
pivot_wider(id_cols = c(id, quest), names_from = c(area_monitor), values_fn = length) %>% #now I'll group by ids!! my df will be sample size * quest
ungroup() %>% #end of function
#now I'll count values
mutate(how_many = select(.,ends_with("_sum")) %>% rowSums(., na.rm = T)) %>%
#now I'll present the results
tabyl(quest, how_many)%>%
as.data.frame() %>%
adorn_totals(c("row", "col")) %>%
adorn_percentages("row") %>%
adorn_pct_formatting(digits = 2) %>%
adorn_ns() #
pivot_longer(-quest) %>%
group_by(quest) %>%
mutate(prop = prop.table(value)) %>%
ungroup()%>%
ggplot(., aes(x = name, y = prop, group=quest, color = quest, label=quest)) +
stat_summary(geom = "line", fun = mean, size=1) +
#ggrepel::geom_text_repel(data = . %>% filter(name == "0"), box.padding = 0.5, max.overlaps = Inf) +
labs(x = "Domains below the cutoff", y = "Proportion", color = "Questionnaire") +
theme_bw() +
theme(legend.position = "bottom")ds %>%
select(quest, ends_with("sum")) %>%
pivot_longer(-quest) %>%
nest_by(quest, name) %>%
mutate(x=list(map(data, ~psych::describe(.)))) %>%
unnest_wider(x) %>%
unnest_wider(value) %>%
mutate_if(is.numeric, round, 2)
ds %>%
select(quest, ends_with("sum")) %>%
pivot_longer(-quest) %>%
nest_by(quest, name) %>%
mutate(
questionnaire = quest,
mean_r = list(map_dbl(data, ~mean(., trim=0.25))),
sd_r = list(map_dbl(data, ~chemometrics::sd_trim(., trim=0.25, const = TRUE))),
m_1sd = mean_r-sd_r,
m_1_half_sd = mean_r-1.5*sd_r,
m_2sd = mean_r-2*sd_r) %>%
unnest(-data) %>%
#unnest_wider(-quest) %>%
pivot_wider(id_cols = quest, names_from = name, values_from = mean_r:m_2sd) %>%
mutate_if(is.numeric, round, 2)ds %>%
select(quest, ends_with("sum")) %>%
pivot_longer(-quest) %>%
nest_by(quest, name) %>%
mutate(
questionnaire = quest,
n = list(map_dbl(data, ~nrow(data.frame(.)))),
mean_r = list(map_dbl(data, ~mean(., trim=0.25))),
sd_r = list(map_dbl(data, ~chemometrics::sd_trim(., trim=0.25, const = TRUE))),
m_1sd = mean_r-sd_r,
m_1_half_sd = mean_r-1.5*sd_r,
m_2sd = mean_r-2*sd_r,
how_many_c = map(data, ~ ifelse(name == "c_sum", #i'll count how many participants are below
ds %>%
filter(quest == questionnaire) %>%
filter(c_sum > m_2sd & c_sum <= m_1sd) %>%
summarise(n()),
NA_integer_)), # percentage of communication
how_many_gm = map(data, ~ ifelse(name == "gm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(gm_sum > m_2sd & gm_sum <= m_1sd) %>%
summarise(n()), #percentage of gross motor
NA_integer_)),
how_many_fm = map(data, ~ ifelse(name == "fm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(fm_sum > m_2sd & fm_sum <= m_1sd) %>%
summarise(n()),#percentage of fine motor
NA_integer_)),
how_many_cg = map(data, ~ ifelse(name == "cg_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(cg_sum > m_2sd & cg_sum <= m_1sd) %>%
summarise(n()), #percentage of problem-solving
NA_integer_)),
how_many_ps = map(data, ~ ifelse(name == "ps_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(ps_sum > m_2sd & ps_sum <= m_1sd) %>%
summarise(n()), #percentage of personal social
NA_integer_))) %>%
mutate(n_com=
paste0(how_many_c[[1]],"-",
round(how_many_c[[1]]/n*100,0),"%"), #I'll get the percentages
n_gm=
paste0(how_many_gm[[1]],"-",
round(how_many_gm[[1]]/n*100,0),"%"),
n_fm=
paste0(how_many_fm[[1]],"-",
round(how_many_fm[[1]]/n*100,0),"%"),
n_cg=
paste0(how_many_cg[[1]],"-",
round(how_many_cg[[1]]/n*100,0),"%"),
n_ps=
paste0(how_many_ps[[1]],"-",
round(how_many_ps[[1]]/n*100,0),"%")
) %>%
unnest_wider(-quest) %>% #particular trickes (there is a lot of pseudo missing values)
select(quest,how_many_ps:last_col(), -how_many_ps) %>% #get some variables
pivot_longer(-quest) %>% #now transform to the long format
filter(value != "NA-NA%") %>% #filter
pivot_wider(id_cols = quest, names_from = name, values_from = value)
ds %>%
select(quest, ends_with("sum")) %>%
pivot_longer(-quest) %>%
nest_by(quest, name) %>%
mutate(
questionnaire = quest,
n = list(map_dbl(data, ~nrow(data.frame(.)))),
mean_r = list(map_dbl(data, ~mean(., trim=0.25))),
sd_r = list(map_dbl(data, ~chemometrics::sd_trim(., trim=0.25, const = TRUE))),
m_1sd = mean_r-sd_r,
m_1_half_sd = mean_r-1.5*sd_r,
m_2sd = mean_r-2*sd_r,
how_many_c = map(data, ~ ifelse(name == "c_sum", #i'll count how many participants are below
ds %>%
filter(quest == questionnaire) %>%
filter(c_sum <= m_2sd) %>% #attention here!
summarise(n()),
NA_integer_)), # percentage of communication
how_many_gm = map(data, ~ ifelse(name == "gm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(gm_sum <= m_2sd) %>%
summarise(n()), #percentage of gross motor
NA_integer_)),
how_many_fm = map(data, ~ ifelse(name == "fm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(fm_sum <= m_2sd) %>%
summarise(n()),#percentage of fine motor
NA_integer_)),
how_many_cg = map(data, ~ ifelse(name == "cg_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(cg_sum <= m_2sd) %>%
summarise(n()), #percentage of problem-solving
NA_integer_)),
how_many_ps = map(data, ~ ifelse(name == "ps_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(ps_sum <= m_2sd) %>%
summarise(n()), #percentage of personal social
NA_integer_))) %>%
mutate(n_com=
paste0(how_many_c[[1]],"-",
round(how_many_c[[1]]/n*100,0),"%"), #I'll get the percentages
n_gm=
paste0(how_many_gm[[1]],"-",
round(how_many_gm[[1]]/n*100,0),"%"),
n_fm=
paste0(how_many_fm[[1]],"-",
round(how_many_fm[[1]]/n*100,0),"%"),
n_cg=
paste0(how_many_cg[[1]],"-",
round(how_many_cg[[1]]/n*100,0),"%"),
n_ps=
paste0(how_many_ps[[1]],"-",
round(how_many_ps[[1]]/n*100,0),"%")
) %>%
unnest_wider(-quest) %>% #particular trickes (there is a lot of pseudo missing values)
select(quest,how_many_ps:last_col(), -how_many_ps) %>% #get some variables
pivot_longer(-quest) %>% #now transform to the long format
filter(value != "NA-NA%") %>% #filter
pivot_wider(id_cols = quest, names_from = name, values_from = value)ds %>%
filter(quest == 16) %>%
filter(gm_sum > 55.49 & gm_sum <= 57.24)ds %>%
filter(quest == 2) %>%
filter(c_sum <= 49.17 & c_sum > 43.62)ds %>%
filter(quest == 18) %>%
filter(fm_sum < 45.38 & fm_sum <= 49.8)ds %>%
select(id, quest, ends_with("sum")) %>%
pivot_longer(-c(id, quest)) %>% #tranform into long format
nest_by(quest, name) %>% #nest (id will not be used this time)
#here is the trick. all counts are here
mutate(
n = map_dbl(data[2], ~nrow(data.frame(.))), #compute sample size
mean = map_dbl(data[2], ~mean(., trim = 0.25)), #get the ROBUST means
sd = map_dbl(data[2], ~chemometrics::sd_trim(., trim=0.25, const = TRUE)), #get the ROBUST sd
one_below = mean-sd, #1 below
two_below = mean - 2 * sd,
monitor = sum(one_below >= data[[2]] & two_below < data[[2]])/n,
below = sum(two_below > data[[2]])/n) %>%
#outside of the map, add a new variable to nested data
mutate(data = list(data %>% mutate(area_monitor = as.character(ifelse(data[[2]] > two_below & data[[2]] <= one_below, paste0(name),0))))) %>%
unnest(data) %>% #I'll unnest and get again the ids. my df will be very long
pivot_wider(id_cols = c(id, quest), names_from = c(area_monitor), values_fn = length) %>% #now I'll group by ids!! my df will be sample size * quest
ungroup() %>% #end of function
#now I'll count values
mutate(how_many = select(.,ends_with("_sum")) %>% rowSums(., na.rm = T)) %>%
#now I'll present the results
tabyl(quest, how_many)%>%
adorn_totals(c("row", "col")) %>%
adorn_percentages("row") %>%
adorn_pct_formatting(digits = 2) %>%
adorn_ns() #
pivot_longer(-quest) %>%
group_by(quest) %>%
mutate(prop = prop.table(value)) %>%
ungroup()%>%
ggplot(., aes(x = name, y = prop, group=quest, color = quest, label=quest)) +
stat_summary(geom = "line", fun = mean, size=1) +
#ggrepel::geom_text_repel(data = . %>% filter(name == "0"), box.padding = 0.5, max.overlaps = Inf) +
labs(x = "Domains at the monitoring zone", y = "Proportion", color = "Questionnaire") +
theme_bw() +
theme(legend.position = "bottom")ds %>%
select(id, quest, ends_with("sum")) %>%
pivot_longer(-c(id, quest)) %>% #tranform into long format
nest_by(quest, name) %>% #nest (id will not be used this time)
#here is the trick. all counts are here
mutate(
n = map_dbl(data[2], ~nrow(data.frame(.))), #compute sample size
mean = map_dbl(data[2], ~mean(., trim = 0.25)), #get the ROBUST means
sd = map_dbl(data[2], ~chemometrics::sd_trim(., trim=0.25, const = TRUE)), #get the ROBUST sd
one_below = mean-sd, #1 below
two_below = mean - 2 * sd) %>%
#outside of the map, add a new variable to nested data
mutate(data = list(data %>% mutate(area_monitor = as.character(ifelse(data[[2]] <= two_below, paste0(name),0))))) %>% #here
unnest(data) %>% #I'll unnest and get again the ids. my df will be very long
pivot_wider(id_cols = c(id, quest), names_from = c(area_monitor), values_fn = length) %>% #now I'll group by ids!! my df will be sample size * quest
ungroup() %>% #end of function
#now I'll count values
mutate(how_many = select(.,ends_with("_sum")) %>% rowSums(., na.rm = T)) %>%
#now I'll present the results
tabyl(quest, how_many)%>%
adorn_totals(c("row", "col")) %>%
adorn_percentages("row") %>%
adorn_pct_formatting(digits = 2) %>%
adorn_ns() #
pivot_longer(-quest) %>%
group_by(quest) %>%
mutate(prop = prop.table(value)) %>%
ungroup()%>%
ggplot(., aes(x = name, y = prop, group=quest, color = quest, label=quest)) +
stat_summary(geom = "line", fun = mean, size=1) +
#ggrepel::geom_text_repel(data = . %>% filter(name == "0"), box.padding = 0.5, max.overlaps = Inf) +
labs(x = "Domains at the monitoring zone", y = "Proportion", color = "Questionnaire") +
theme_bw() +
theme(legend.position = "bottom")ds %>% #get ds
select(quest, ends_with("sum")) %>% #select the focus variables
pivot_longer(-quest) %>% #transpose to long formata
nest_by(quest, name) %>% #group
mutate(mean = list(map_dbl(data, ~mean(.))),
p_05 = list(map_dbl(data, ~quantile(., prob = 0.05))),
p_10 = list(map_dbl(data, ~quantile(., prob = 0.10)))
) %>%
unnest(-data) %>% #unnest
pivot_wider(id_cols = quest, names_from = name, values_from = mean:p_10) %>%
mutate_if(is.numeric, round, 2)Double check
ds %>% filter(quest==10) %>% summarise(x=ecdf(fm_sum)(30))
ds %>% filter(quest==10) %>% summarise(x=quantile(fm_sum, probs=c(.05)))ds %>%
select(quest, ends_with("sum")) %>% #get variable names
pivot_longer(-quest) %>% #tranform into the long format
nest_by(quest, name) %>% #group or nest
mutate(
questionnaire = quest,#compute questionnaire
n = list(map_dbl(data, ~nrow(data.frame(.)))), #compute sample size
m_2sd = list(map_dbl(data, ~quantile(., prob = 0.05))), #2 below (now, 10th percentile)
how_many_c = map(data, ~ ifelse(name == "c_sum", #i'll count how many participants are below
ds %>%
filter(quest == questionnaire) %>%
filter(c_sum <= m_2sd) %>%
summarise(n()),
NA_integer_)), # percentage of communication
how_many_gm = map(data, ~ ifelse(name == "gm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(gm_sum <= m_2sd) %>%
summarise(n()), #percentage of gross motor
NA_integer_)),
how_many_fm = map(data, ~ ifelse(name == "fm_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(fm_sum <= m_2sd) %>%
summarise(n()),#percentage of fine motor
NA_integer_)),
how_many_cg = map(data, ~ ifelse(name == "cg_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(cg_sum <= m_2sd) %>%
summarise(n()), #percentage of problem-solving
NA_integer_)),
how_many_ps = map(data, ~ ifelse(name == "ps_sum",
ds %>%
filter(quest == questionnaire) %>%
filter(ps_sum <= m_2sd) %>%
summarise(n()), #percentage of personal social
NA_integer_))) %>%
mutate(n_com=
paste0(how_many_c[[1]],"-",
round(how_many_c[[1]]/n*100,0),"%"), #I'll get the percentages
n_gm=
paste0(how_many_gm[[1]],"-",
round(how_many_gm[[1]]/n*100,0),"%"),
n_fm=
paste0(how_many_fm[[1]],"-",
round(how_many_fm[[1]]/n*100,0),"%"),
n_cg=
paste0(how_many_cg[[1]],"-",
round(how_many_cg[[1]]/n*100,0),"%"),
n_ps=
paste0(how_many_ps[[1]],"-",
round(how_many_ps[[1]]/n*100,0),"%")
) %>%
unnest_wider(-quest) %>% #particular trickes (there is a lot of pseudo missing values)
select(quest,how_many_ps:last_col(), -how_many_ps) %>% #get some variables
pivot_longer(-quest) %>% #now transform to the long format
filter(value != "NA-NA%") %>% #filter
pivot_wider(id_cols = quest, names_from = name, values_from = value)ds %>%
select(id, quest, ends_with("sum")) %>%
pivot_longer(-c(id, quest)) %>% #tranform into long format
nest_by(quest, name) %>% #nest (id will not be used this time)
#here is the trick. all counts are here
mutate(mean = list(map_dbl(data[2], ~mean(.))),
two_below = list(map_dbl(data[2], ~quantile(., prob = 0.05))),
one_below = list(map_dbl(data[2], ~quantile(., prob = 0.10)))
) %>%
#outside of the map, add a new variable to nested data
mutate(data = list(data %>% mutate(area_monitor = as.character(ifelse(data[[2]] > two_below & data[[2]] <= one_below, paste0(name),0))))) %>%
unnest(data) %>% #I'll unnest and get again the ids. my df will be very long
pivot_wider(id_cols = c(id, quest), names_from = c(area_monitor), values_fn = length) %>% #now I'll group by ids!! my df will be sample size * quest
ungroup() %>% #end of function
#now I'll count values
mutate(how_many = select(.,ends_with("_sum")) %>% rowSums(., na.rm = T)) %>%
#now I'll present the results
tabyl(quest, how_many)%>%
as.data.frame() %>%
#adorn_totals(c("row", "col")) %>%
#adorn_percentages("row") %>%
#adorn_pct_formatting(digits = 2) %>%
#adorn_ns()
pivot_longer(-quest) %>%
group_by(quest) %>%
mutate(prop = prop.table(value)) %>%
ungroup()%>%
ggplot(., aes(x = name, y = prop, group=quest, color = quest, label=quest)) +
stat_summary(geom = "line", fun = mean, size=1) +
#ggrepel::geom_text_repel(data = . %>% filter(name == "0"), box.padding = 0.5, max.overlaps = Inf) +
labs(x = "Domains at the monitoring zone", y = "Proportion", color = "Questionnaire") +
theme_bw() +
theme(legend.position = "bottom")ds %>%
select(id, quest, ends_with("sum")) %>%
pivot_longer(-c(id, quest)) %>% #tranform into long format
nest_by(quest, name) %>% #nest (id will not be used this time)
#here is the trick. all counts are here
mutate(mean = list(map_dbl(data[2], ~mean(.))),
two_below = list(map_dbl(data[2], ~quantile(., prob = 0.05))),
one_below = list(map_dbl(data[2], ~quantile(., prob = 0.10)))
) %>%
#outside of the map, add a new variable to nested data
mutate(data = list(data %>% mutate(area_monitor = as.character(ifelse(data[[2]] <= two_below, paste0(name),0))))) %>%
unnest(data) %>% #I'll unnest and get again the ids. my df will be very long
pivot_wider(id_cols = c(id, quest), names_from = c(area_monitor), values_fn = length) %>% #now I'll group by ids!! my df will be sample size * quest
ungroup() %>% #end of function
#now I'll count values
mutate(how_many = select(.,ends_with("_sum")) %>% rowSums(., na.rm = T)) %>%
#now I'll present the results
tabyl(quest, how_many)%>%
as.data.frame() %>%
#adorn_totals(c("row", "col")) %>%
#adorn_percentages("row") %>%
#adorn_pct_formatting(digits = 2) %>%
#adorn_ns()
pivot_longer(-quest) %>%
group_by(quest) %>%
mutate(prop = prop.table(value)) %>%
ungroup()%>%
ggplot(., aes(x = name, y = prop, group=quest, color = quest, label=quest)) +
stat_summary(geom = "line", fun = mean, size=1) +
#ggrepel::geom_text_repel(data = . %>% filter(name == "0"), box.padding = 0.5, max.overlaps = Inf) +
labs(x = "Domains below the cutoff", y = "Proportion", color = "Questionnaire") +
theme_bw() +
theme(legend.position = "bottom")Testing to double check the differences between percentiles
#percentile 10
df_percentil_10 <- ds %>%
filter(quest == 2) %>%
select(id, quest, ends_with("sum")) %>%
pivot_longer(-c(id, quest)) %>% #tranform into long format
nest_by(quest, name) %>% #nest (id will not be used this time)
#here is the trick. all counts are here
mutate(mean = list(map_dbl(data[2], ~mean(.))),
two_below = list(map_dbl(data[2], ~quantile(., prob = 0.05))),
one_below = list(map_dbl(data[2], ~quantile(., prob = 0.10)))
) %>%
#outside of the map, add a new variable to nested data
mutate(data = list(data %>% mutate(area_monitor = as.character(ifelse(data[[2]] > two_below & data[[2]] <= one_below, paste0(name),0))))) %>%
unnest(data) %>%
unnest(id) %>%
pivot_wider(id_cols = id, names_from = area_monitor, values_fn = length) %>%
mutate(how_many = select(.,ends_with("_sum")) %>% rowSums(., na.rm = T)) %>%
filter(how_many == 2) %>%
rename_at(vars(ends_with("_sum")), ~paste(.,"_r")) %>%
mutate_all(~replace_na(., 0)) %>%
ungroup()
df_percentil_10 <- left_join(
df_percentil_10,
ds %>%
select(id, ends_with("_sum")),
by = "id") %>%
pivot_longer(-id) %>%
arrange(name) %>%
pivot_wider(id_cols = id, names_from = name, values_from = value, values_fill = 0) %>%
relocate(how_many, .before = "0") %>%
select(-"0",-how_many)
df_percentil_05 <- ds %>%
filter(quest == 2) %>%
select(id, quest, ends_with("sum")) %>%
pivot_longer(-c(id, quest)) %>% #tranform into long format
nest_by(quest, name) %>% #nest (id will not be used this time)
#here is the trick. all counts are here
mutate(mean = list(map_dbl(data[2], ~mean(.))),
two_below = list(map_dbl(data[2], ~quantile(., prob = 0.05))),
one_below = list(map_dbl(data[2], ~quantile(., prob = 0.10)))
) %>%
#outside of the map, add a new variable to nested data
mutate(data = list(data %>% mutate(area_monitor = as.character(ifelse(data[[2]] <= two_below, paste0(name),0))))) %>%
unnest(data) %>%
unnest(id) %>%
pivot_wider(id_cols = id, names_from = area_monitor, values_fn = length) %>%
mutate(how_many = select(.,ends_with("_sum")) %>% rowSums(., na.rm = T)) %>%
filter(how_many == 2) %>%
rename_at(vars(ends_with("_sum")), ~paste(.,"_r")) %>%
mutate_all(~replace_na(., 0)) %>%
ungroup()
df_percentil_05 <- left_join(
df_percentil_05,
ds %>%
select(id, ends_with("_sum")),
by = "id") %>%
pivot_longer(-id) %>%
arrange(name) %>%
pivot_wider(id_cols = id, names_from = name, values_from = value, values_fill = 0) %>%
relocate(how_many, .before = "0") %>%
select(-"0",-how_many)
rbind(df_percentil_10 %>% mutate(base="p10"),
df_percentil_05 %>% mutate(base="p05")) %>%
filter(id == "0931c359-d313-4f18-9d26-fc4bb192c0ed") %>% t()
arrange(id) %>%
select(base, everything())decided_cutoff <- function(quest, domain, type) {
quest = enquo(quest)
domain = enquo(domain)
ds2 = ds %>% mutate(quest = if_else(quest == "9", "10", as.character(quest))) #to combine 9 and 10 months questionnaires
if (type == "trad") {
#table 1
x = ds2 %>%
filter(quest == !!quest) %>%
summarise(n = n(),
mean = mean(!!domain, na.rm=T), sd = sd(!!domain,na.rm=T),
cutoff = mean-2*sd) %>%
mutate(quest = !!quest,domain = quo_name(domain)) #to identify which quest in the output list
#table2
y = ds2 %>%
filter(quest == !!quest) %>%
mutate(n = n(),
mean = mean(!!domain, na.rm=T), sd = sd(!!domain,na.rm=T),
cutoff = mean-2*sd,
class = ifelse(!!domain <= cutoff,"below","above")) %>%
tabyl(class)
j = cbind.data.frame(x,y)
j = data.frame(lapply(j, function(j) if(is.numeric(j)) round(j, 2) else j)) #round
j = j %>% select(quest, domain,everything()) #change position
#https://stackoverflow.com/questions/9063889/how-to-round-a-data-frame-in-r-that-contains-some-character-variables
return(j)
}
if (type == "10") {
#table 1
x = ds2 %>%
filter(quest == !!quest) %>%
summarise(n = n(),
mean = mean(!!domain, na.rm=T), sd = sd(!!domain,na.rm=T),
p10 = quantile(x = !!domain, prob = 0.1),
cutoff = p10) %>%
mutate(quest = !!quest,domain = quo_name(domain)) #to identify which quest in the output list
#table2
y = ds2 %>%
filter(quest == !!quest) %>%
mutate(cutoff = quantile(x = !!domain, prob = 0.1),
class = ifelse(!!domain <= cutoff,"below","above")) %>%
tabyl(class)
j = cbind.data.frame(x,y)
j = data.frame(lapply(j, function(j) if(is.numeric(j)) round(j, 2) else j)) #round
j = j %>% select(quest, domain,everything()) #change position
return(j)
}
if (type == "20p") {
#table 1
x = ds2 %>%
filter(quest == !!quest) %>%
summarise(n = n(),
mean = mean(!!domain, na.rm=T), sd = sd(!!domain,na.rm=T),
cutoff = 20) %>%
mutate(quest = !!quest,domain = quo_name(domain)) #to identify which quest in the output list
#table2
y = ds2 %>%
filter(quest == !!quest) %>%
mutate(cutoff = 20,
class = ifelse(!!domain <= cutoff,"below","above")) %>%
tabyl(class)
j = cbind.data.frame(x,y)
j = data.frame(lapply(j, function(j) if(is.numeric(j)) round(j, 2) else j)) #round
j = j %>% select(quest, domain,everything()) #change position
return(j)
}
if (type == "25p") {
#table 1
x = ds2 %>%
filter(quest == !!quest) %>%
summarise(n = n(),
mean = mean(!!domain, na.rm=T), sd = sd(!!domain,na.rm=T),
cutoff = 25) %>%
mutate(quest = !!quest,domain = quo_name(domain)) #to identify which quest in the output list
#table2
y = ds2 %>%
filter(quest == !!quest) %>%
mutate(cutoff = 25,
class = ifelse(!!domain <= cutoff,"below","above")) %>%
tabyl(class)
j = cbind.data.frame(x,y)
j = data.frame(lapply(j, function(j) if(is.numeric(j)) round(j, 2) else j)) #round
j = j %>% select(quest, domain,everything()) #change position
return(j)
}
if (type == "30p") {
#table 1
x = ds2 %>%
filter(quest == !!quest) %>%
summarise(n = n(),
mean = mean(!!domain, na.rm=T), sd = sd(!!domain,na.rm=T),
cutoff = 30) %>%
mutate(quest = !!quest,domain = quo_name(domain)) #to identify which quest in the output list
#table2
y = ds2 %>%
filter(quest == !!quest) %>%
mutate(cutoff = 30,
class = ifelse(!!domain <= cutoff,"below","above")) %>%
tabyl(class)
j = cbind.data.frame(x,y)
j = data.frame(lapply(j, function(j) if(is.numeric(j)) round(j, 2) else j)) #round
j = j %>% select(quest, domain,everything()) #change position
return(j)
}
if (type == "35p") {
#table 1
x = ds2 %>%
filter(quest == !!quest) %>%
summarise(n = n(),
mean = mean(!!domain, na.rm=T), sd = sd(!!domain,na.rm=T),
cutoff = 35) %>%
mutate(quest = !!quest,domain = quo_name(domain)) #to identify which quest in the output list
#table2
y = ds2 %>%
filter(quest == !!quest) %>%
mutate(cutoff = 35,
class = ifelse(!!domain <= cutoff,"below","above")) %>%
tabyl(class)
j = cbind.data.frame(x,y)
j = data.frame(lapply(j, function(j) if(is.numeric(j)) round(j, 2) else j)) #round
j = j %>% select(quest, domain,everything()) #change position
return(j)
}
}print(list(
decided_cutoff(type = "trad", quest = 2, domain = c_sum),
decided_cutoff(type = "trad", quest = 2, domain = gm_sum),
decided_cutoff(type = "trad", quest = 2, domain = fm_sum),
decided_cutoff(type = "trad", quest = 2, domain = cg_sum),
decided_cutoff(type = "trad", quest = 2, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 2 c_sum 1404 42.24 14.42 13.4 above 1348 0.96
2 2 c_sum 1404 42.24 14.42 13.4 below 56 0.04
[[2]]
quest domain n mean sd cutoff class n.1 percent
1 2 gm_sum 1404 53.84 7.79 38.26 above 1338 0.95
2 2 gm_sum 1404 53.84 7.79 38.26 below 66 0.05
[[3]]
quest domain n mean sd cutoff class n.1 percent
1 2 fm_sum 1404 47.36 10.08 27.2 above 1355 0.97
2 2 fm_sum 1404 47.36 10.08 27.2 below 49 0.03
[[4]]
quest domain n mean sd cutoff class n.1 percent
1 2 cg_sum 1404 45.38 10.98 23.43 above 1355 0.97
2 2 cg_sum 1404 45.38 10.98 23.43 below 49 0.03
[[5]]
quest domain n mean sd cutoff class n.1 percent
1 2 ps_sum 1404 42.98 13.24 16.51 above 1342 0.96
2 2 ps_sum 1404 42.98 13.24 16.51 below 62 0.04print(list(
decided_cutoff(type = "trad", quest = 4, domain = c_sum),
decided_cutoff(type = "trad", quest = 4, domain = gm_sum),
decided_cutoff(type = "trad", quest = 4, domain = fm_sum),
decided_cutoff(type = "trad", quest = 4, domain = cg_sum),
decided_cutoff(type = "trad", quest = 4, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 4 c_sum 1365 50.23 9.85 30.54 above 1281 0.94
2 4 c_sum 1365 50.23 9.85 30.54 below 84 0.06
[[2]]
quest domain n mean sd cutoff class n.1 percent
1 4 gm_sum 1365 53.38 9.09 35.2 above 1267 0.93
2 4 gm_sum 1365 53.38 9.09 35.2 below 98 0.07
[[3]]
quest domain n mean sd cutoff class n.1 percent
1 4 fm_sum 1365 46.14 12.9 20.34 above 1277 0.94
2 4 fm_sum 1365 46.14 12.9 20.34 below 88 0.06
[[4]]
quest domain n mean sd cutoff class n.1 percent
1 4 cg_sum 1365 48.82 10.79 27.25 above 1292 0.95
2 4 cg_sum 1365 48.82 10.79 27.25 below 73 0.05
[[5]]
quest domain n mean sd cutoff class n.1 percent
1 4 ps_sum 1365 50.28 10.8 28.68 above 1307 0.96
2 4 ps_sum 1365 50.28 10.8 28.68 below 58 0.04print(list(
decided_cutoff(type = "trad", quest = 6, domain = c_sum),
decided_cutoff(type = "trad", quest = 6, domain = gm_sum),
decided_cutoff(type = "trad", quest = 6, domain = fm_sum),
decided_cutoff(type = "trad", quest = 6, domain = cg_sum),
decided_cutoff(type = "trad", quest = 6, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 6 c_sum 1099 47.71 9.66 28.39 above 1058 0.96
2 6 c_sum 1099 47.71 9.66 28.39 below 41 0.04
[[2]]
quest domain n mean sd cutoff class n.1 percent
1 6 gm_sum 1099 42.3 12.51 17.27 above 1052 0.96
2 6 gm_sum 1099 42.3 12.51 17.27 below 47 0.04
[[3]]
quest domain n mean sd cutoff class n.1 percent
1 6 fm_sum 1099 47.56 12.85 21.85 above 1042 0.95
2 6 fm_sum 1099 47.56 12.85 21.85 below 57 0.05
[[4]]
quest domain n mean sd cutoff class n.1 percent
1 6 cg_sum 1099 45.71 11.65 22.41 above 1054 0.96
2 6 cg_sum 1099 45.71 11.65 22.41 below 45 0.04
[[5]]
quest domain n mean sd cutoff class n.1 percent
1 6 ps_sum 1099 47 12.76 21.49 above 1034 0.94
2 6 ps_sum 1099 47 12.76 21.49 below 65 0.06print(list(
decided_cutoff(type = "10", quest = 8, domain = c_sum),
decided_cutoff(type = "10", quest = 8, domain = gm_sum),
decided_cutoff(type = "10", quest = 8, domain = fm_sum),
decided_cutoff(type = "10", quest = 8, domain = cg_sum),
decided_cutoff(type = "10", quest = 8, domain = ps_sum)
))[[1]]
quest domain n mean sd p10 cutoff class n.1 percent
1 8 c_sum 1020 49.42 11.02 35 35 above 862 0.85
2 8 c_sum 1020 49.42 11.02 35 35 below 158 0.15
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 8 gm_sum 1020 48.48 12.75 30 30 above 895 0.88
2 8 gm_sum 1020 48.48 12.75 30 30 below 125 0.12
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 8 fm_sum 1020 53.28 10.15 40 40 above 893 0.88
2 8 fm_sum 1020 53.28 10.15 40 40 below 127 0.12
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 8 cg_sum 1020 49.24 11 35 35 above 869 0.85
2 8 cg_sum 1020 49.24 11 35 35 below 151 0.15
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 8 ps_sum 1020 50.84 10.41 35 35 above 896 0.88
2 8 ps_sum 1020 50.84 10.41 35 35 below 124 0.12print(list(
decided_cutoff(type = "25p", quest = 10, domain = c_sum),
decided_cutoff(type = "25p", quest = 10, domain = gm_sum),
decided_cutoff(type = "10", quest = 10, domain = fm_sum),
decided_cutoff(type = "10", quest = 10, domain = cg_sum),
decided_cutoff(type = "10", quest = 10, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 10 c_sum 937 41.86 13.9 25 above 793 0.85
2 10 c_sum 937 41.86 13.9 25 below 144 0.15
[[2]]
quest domain n mean sd cutoff class n.1 percent
1 10 gm_sum 937 45.3 14.65 25 above 797 0.85
2 10 gm_sum 937 45.3 14.65 25 below 140 0.15
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 10 fm_sum 937 51.84 10.51 40 40 above 768 0.82
2 10 fm_sum 937 51.84 10.51 40 40 below 169 0.18
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 10 cg_sum 937 47.62 10.9 30 30 above 837 0.89
2 10 cg_sum 937 47.62 10.9 30 30 below 100 0.11
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 10 ps_sum 937 48.86 11.41 35 35 above 787 0.84
2 10 ps_sum 937 48.86 11.41 35 35 below 150 0.16print(list(
decided_cutoff(type = "10", quest = 12, domain = c_sum),
decided_cutoff(type = "10", quest = 12, domain = gm_sum),
decided_cutoff(type = "10", quest = 12, domain = fm_sum),
decided_cutoff(type = "10", quest = 12, domain = cg_sum),
decided_cutoff(type = "10", quest = 12, domain = ps_sum)
))[[1]]
quest domain n mean sd p10 cutoff class n.1 percent
1 12 c_sum 1018 44.11 13.55 25 25 above 885 0.87
2 12 c_sum 1018 44.11 13.55 25 25 below 133 0.13
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 12 gm_sum 1018 48.21 13.31 30 30 above 888 0.87
2 12 gm_sum 1018 48.21 13.31 30 30 below 130 0.13
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 12 fm_sum 1018 44.95 11.2 30 30 above 874 0.86
2 12 fm_sum 1018 44.95 11.2 30 30 below 144 0.14
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 12 cg_sum 1018 41.85 13.87 23.5 23.5 above 916 0.9
2 12 cg_sum 1018 41.85 13.87 23.5 23.5 below 102 0.1
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 12 ps_sum 1018 45.49 12.39 30 30 above 865 0.85
2 12 ps_sum 1018 45.49 12.39 30 30 below 153 0.15print(list(
decided_cutoff(type = "25p", quest = 14, domain = c_sum),
decided_cutoff(type = "10", quest = 14, domain = gm_sum),
decided_cutoff(type = "10", quest = 14, domain = fm_sum),
decided_cutoff(type = "10", quest = 14, domain = cg_sum),
decided_cutoff(type = "10", quest = 14, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 14 c_sum 798 38.75 15.34 25 above 602 0.75
2 14 c_sum 798 38.75 15.34 25 below 196 0.25
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 14 gm_sum 798 51.64 13.98 30 30 above 682 0.85
2 14 gm_sum 798 51.64 13.98 30 30 below 116 0.15
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 14 fm_sum 798 45.64 12.05 30 30 above 677 0.85
2 14 fm_sum 798 45.64 12.05 30 30 below 121 0.15
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 14 cg_sum 798 43.2 14.04 20 20 above 713 0.89
2 14 cg_sum 798 43.2 14.04 20 20 below 85 0.11
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 14 ps_sum 798 41.74 13.9 20 20 above 707 0.89
2 14 ps_sum 798 41.74 13.9 20 20 below 91 0.11print(list(
decided_cutoff(type = "25p", quest = 16, domain = c_sum),
decided_cutoff(type = "10", quest = 16, domain = gm_sum),
decided_cutoff(type = "10", quest = 16, domain = fm_sum),
decided_cutoff(type = "10", quest = 16, domain = cg_sum),
decided_cutoff(type = "10", quest = 16, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 16 c_sum 915 39.12 14.96 25 above 721 0.79
2 16 c_sum 915 39.12 14.96 25 below 194 0.21
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 16 gm_sum 915 53.2 13.73 35 35 above 818 0.89
2 16 gm_sum 915 53.2 13.73 35 35 below 97 0.11
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 16 fm_sum 915 49.64 11.08 35 35 above 793 0.87
2 16 fm_sum 915 49.64 11.08 35 35 below 122 0.13
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 16 cg_sum 915 46.61 13.62 25 25 above 809 0.88
2 16 cg_sum 915 46.61 13.62 25 25 below 106 0.12
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 16 ps_sum 915 44.22 14.4 20 20 above 820 0.9
2 16 ps_sum 915 44.22 14.4 20 20 below 95 0.1print(list(
decided_cutoff(type = "25p", quest = 18, domain = c_sum),
decided_cutoff(type = "10", quest = 18, domain = gm_sum),
decided_cutoff(type = "10", quest = 18, domain = fm_sum),
decided_cutoff(type = "10", quest = 18, domain = cg_sum),
decided_cutoff(type = "10", quest = 18, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 18 c_sum 1087 33.56 17 25 above 705 0.65
2 18 c_sum 1087 33.56 17 25 below 382 0.35
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 18 gm_sum 1087 54.52 9.09 45 45 above 932 0.86
2 18 gm_sum 1087 54.52 9.09 45 45 below 155 0.14
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 18 fm_sum 1087 51.69 9.53 40 40 above 923 0.85
2 18 fm_sum 1087 51.69 9.53 40 40 below 164 0.15
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 18 cg_sum 1087 48.22 11.64 35 35 above 921 0.85
2 18 cg_sum 1087 48.22 11.64 35 35 below 166 0.15
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 18 ps_sum 1087 47.71 11.55 30 30 above 960 0.88
2 18 ps_sum 1087 47.71 11.55 30 30 below 127 0.12print(list(
decided_cutoff(type = "25p", quest = 20, domain = c_sum),
decided_cutoff(type = "10", quest = 20, domain = gm_sum),
decided_cutoff(type = "10", quest = 20, domain = fm_sum),
decided_cutoff(type = "10", quest = 20, domain = cg_sum),
decided_cutoff(type = "10", quest = 20, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 20 c_sum 756 36.12 19.87 25 above 492 0.65
2 20 c_sum 756 36.12 19.87 25 below 264 0.35
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 20 gm_sum 756 54.11 9.23 45 45 above 630 0.83
2 20 gm_sum 756 54.11 9.23 45 45 below 126 0.17
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 20 fm_sum 756 50.17 10.23 35 35 above 673 0.89
2 20 fm_sum 756 50.17 10.23 35 35 below 83 0.11
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 20 cg_sum 756 43.23 13.55 20 20 above 675 0.89
2 20 cg_sum 756 43.23 13.55 20 20 below 81 0.11
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 20 ps_sum 756 41.54 13.49 20 20 above 679 0.9
2 20 ps_sum 756 41.54 13.49 20 20 below 77 0.1print(list(
decided_cutoff(type = "25p", quest = 22, domain = c_sum),
decided_cutoff(type = "10", quest = 22, domain = gm_sum),
decided_cutoff(type = "10", quest = 22, domain = fm_sum),
decided_cutoff(type = "10", quest = 22, domain = cg_sum),
decided_cutoff(type = "10", quest = 22, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 22 c_sum 708 37.81 19.49 25 above 477 0.67
2 22 c_sum 708 37.81 19.49 25 below 231 0.33
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 22 gm_sum 708 50.13 10.9 35 35 above 618 0.87
2 22 gm_sum 708 50.13 10.9 35 35 below 90 0.13
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 22 fm_sum 708 45.68 10.52 30 30 above 618 0.87
2 22 fm_sum 708 45.68 10.52 30 30 below 90 0.13
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 22 cg_sum 708 44.22 13.51 25 25 above 616 0.87
2 22 cg_sum 708 44.22 13.51 25 25 below 92 0.13
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 22 ps_sum 708 49.4 11.9 30 30 above 634 0.9
2 22 ps_sum 708 49.4 11.9 30 30 below 74 0.1print(list(
decided_cutoff(type = "25p", quest = 24, domain = c_sum),
decided_cutoff(type = "10", quest = 24, domain = gm_sum),
decided_cutoff(type = "10", quest = 24, domain = fm_sum),
decided_cutoff(type = "10", quest = 24, domain = cg_sum),
decided_cutoff(type = "10", quest = 24, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 24 c_sum 1043 41.49 19.31 25 above 775 0.74
2 24 c_sum 1043 41.49 19.31 25 below 268 0.26
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 24 gm_sum 1043 52.64 10.02 40 40 above 905 0.87
2 24 gm_sum 1043 52.64 10.02 40 40 below 138 0.13
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 24 fm_sum 1043 41.89 11.92 25 25 above 921 0.88
2 24 fm_sum 1043 41.89 11.92 25 25 below 122 0.12
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 24 cg_sum 1043 44.88 13.37 25 25 above 918 0.88
2 24 cg_sum 1043 44.88 13.37 25 25 below 125 0.12
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 24 ps_sum 1043 48.17 11.94 30 30 above 922 0.88
2 24 ps_sum 1043 48.17 11.94 30 30 below 121 0.12print(list(
decided_cutoff(type = "25p", quest = 27, domain = c_sum),
decided_cutoff(type = "10", quest = 27, domain = gm_sum),
decided_cutoff(type = "10", quest = 27, domain = fm_sum),
decided_cutoff(type = "10", quest = 27, domain = cg_sum),
decided_cutoff(type = "10", quest = 27, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 27 c_sum 778 44.95 16.78 25 above 642 0.83
2 27 c_sum 778 44.95 16.78 25 below 136 0.17
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 27 gm_sum 778 51.25 10.37 35 35 above 692 0.89
2 27 gm_sum 778 51.25 10.37 35 35 below 86 0.11
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 27 fm_sum 778 39.81 14.34 20 20 above 687 0.88
2 27 fm_sum 778 39.81 14.34 20 20 below 91 0.12
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 27 cg_sum 778 44.16 14.03 25 25 above 668 0.86
2 27 cg_sum 778 44.16 14.03 25 25 below 110 0.14
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 27 ps_sum 778 47.74 12.86 30 30 above 672 0.86
2 27 ps_sum 778 47.74 12.86 30 30 below 106 0.14print(list(
decided_cutoff(type = "30p", quest = 30, domain = c_sum),
decided_cutoff(type = "10", quest = 30, domain = gm_sum),
decided_cutoff(type = "10", quest = 30, domain = fm_sum),
decided_cutoff(type = "10", quest = 30, domain = cg_sum),
decided_cutoff(type = "10", quest = 30, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 30 c_sum 795 48.29 15.47 30 above 676 0.85
2 30 c_sum 795 48.29 15.47 30 below 119 0.15
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 30 gm_sum 795 50.75 10.69 35 35 above 703 0.88
2 30 gm_sum 795 50.75 10.69 35 35 below 92 0.12
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 30 fm_sum 795 41.79 14.85 20 20 above 684 0.86
2 30 fm_sum 795 41.79 14.85 20 20 below 111 0.14
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 30 cg_sum 795 46.38 13.11 27 27 above 715 0.9
2 30 cg_sum 795 46.38 13.11 27 27 below 80 0.1
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 30 ps_sum 795 45.87 13.07 27 27 above 715 0.9
2 30 ps_sum 795 45.87 13.07 27 27 below 80 0.1print(list(
decided_cutoff(type = "10", quest = 33, domain = c_sum),
decided_cutoff(type = "10", quest = 33, domain = gm_sum),
decided_cutoff(type = "10", quest = 33, domain = fm_sum),
decided_cutoff(type = "10", quest = 33, domain = cg_sum),
decided_cutoff(type = "10", quest = 33, domain = ps_sum)
))[[1]]
quest domain n mean sd p10 cutoff class n.1 percent
1 33 c_sum 672 48.78 15.48 25 25 above 580 0.86
2 33 c_sum 672 48.78 15.48 25 25 below 92 0.14
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 33 gm_sum 672 52.72 9.7 40 40 above 582 0.87
2 33 gm_sum 672 52.72 9.7 40 40 below 90 0.13
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 33 fm_sum 672 38.93 16.47 15 15 above 579 0.86
2 33 fm_sum 672 38.93 16.47 15 15 below 93 0.14
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 33 cg_sum 672 47.47 12.46 30 30 above 587 0.87
2 33 cg_sum 672 47.47 12.46 30 30 below 85 0.13
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 33 ps_sum 672 48.42 12.33 30 30 above 583 0.87
2 33 ps_sum 672 48.42 12.33 30 30 below 89 0.13print(list(
decided_cutoff(type = "30p", quest = 36, domain = c_sum),
decided_cutoff(type = "10", quest = 36, domain = gm_sum),
decided_cutoff(type = "10", quest = 36, domain = fm_sum),
decided_cutoff(type = "10", quest = 36, domain = cg_sum),
decided_cutoff(type = "10", quest = 36, domain = ps_sum)
))[[1]]
quest domain n mean sd cutoff class n.1 percent
1 36 c_sum 1055 48.55 14.24 30 above 911 0.86
2 36 c_sum 1055 48.55 14.24 30 below 144 0.14
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 36 gm_sum 1055 51.42 10.51 35 35 above 940 0.89
2 36 gm_sum 1055 51.42 10.51 35 35 below 115 0.11
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 36 fm_sum 1055 39.57 17.24 15 15 above 908 0.86
2 36 fm_sum 1055 39.57 17.24 15 15 below 147 0.14
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 36 cg_sum 1055 48.51 12.56 30 30 above 922 0.87
2 36 cg_sum 1055 48.51 12.56 30 30 below 133 0.13
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 36 ps_sum 1055 50.26 11.2 35 35 above 920 0.87
2 36 ps_sum 1055 50.26 11.2 35 35 below 135 0.13print(list(
decided_cutoff(type = "10", quest = 42, domain = c_sum),
decided_cutoff(type = "10", quest = 42, domain = gm_sum),
decided_cutoff(type = "25p", quest = 42, domain = fm_sum),
decided_cutoff(type = "10", quest = 42, domain = cg_sum),
decided_cutoff(type = "10", quest = 42, domain = ps_sum)
))[[1]]
quest domain n mean sd p10 cutoff class n.1 percent
1 42 c_sum 1256 51.29 12.39 35 35 above 1094 0.87
2 42 c_sum 1256 51.29 12.39 35 35 below 162 0.13
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 42 gm_sum 1256 52.52 9.65 40 40 above 1091 0.87
2 42 gm_sum 1256 52.52 9.65 40 40 below 165 0.13
[[3]]
quest domain n mean sd cutoff class n.1 percent
1 42 fm_sum 1256 42.17 16.18 25 above 1012 0.81
2 42 fm_sum 1256 42.17 16.18 25 below 244 0.19
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 42 cg_sum 1256 49.42 12.27 30 30 above 1121 0.89
2 42 cg_sum 1256 49.42 12.27 30 30 below 135 0.11
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 42 ps_sum 1256 50.53 11.54 35 35 above 1099 0.88
2 42 ps_sum 1256 50.53 11.54 35 35 below 157 0.12print(list(
decided_cutoff(type = "10", quest = 48, domain = c_sum),
decided_cutoff(type = "10", quest = 48, domain = gm_sum),
decided_cutoff(type = "10", quest = 48, domain = fm_sum),
decided_cutoff(type = "10", quest = 48, domain = cg_sum),
decided_cutoff(type = "10", quest = 48, domain = ps_sum)
))[[1]]
quest domain n mean sd p10 cutoff class n.1 percent
1 48 c_sum 1210 51.49 12.59 35 35 above 1068 0.88
2 48 c_sum 1210 51.49 12.59 35 35 below 142 0.12
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 48 gm_sum 1210 52.09 10.71 35 35 above 1085 0.9
2 48 gm_sum 1210 52.09 10.71 35 35 below 125 0.1
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 48 fm_sum 1210 39.92 16.07 15 15 above 1067 0.88
2 48 fm_sum 1210 39.92 16.07 15 15 below 143 0.12
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 48 cg_sum 1210 47.63 12.43 30 30 above 1059 0.88
2 48 cg_sum 1210 47.63 12.43 30 30 below 151 0.12
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 48 ps_sum 1210 51.54 10.65 35 35 above 1080 0.89
2 48 ps_sum 1210 51.54 10.65 35 35 below 130 0.11print(list(
decided_cutoff(type = "10", quest = 54, domain = c_sum),
decided_cutoff(type = "35p", quest = 54, domain = gm_sum),
decided_cutoff(type = "10", quest = 54, domain = fm_sum),
decided_cutoff(type = "10", quest = 54, domain = cg_sum),
decided_cutoff(type = "10", quest = 54, domain = ps_sum)
))[[1]]
quest domain n mean sd p10 cutoff class n.1 percent
1 54 c_sum 1193 52.92 10.91 40 40 above 1031 0.86
2 54 c_sum 1193 52.92 10.91 40 40 below 162 0.14
[[2]]
quest domain n mean sd cutoff class n.1 percent
1 54 gm_sum 1193 52.7 10.6 35 above 1077 0.9
2 54 gm_sum 1193 52.7 10.6 35 below 116 0.1
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 54 fm_sum 1193 42.26 16.14 15 15 above 1066 0.89
2 54 fm_sum 1193 42.26 16.14 15 15 below 127 0.11
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 54 cg_sum 1193 48.16 11.99 30 30 above 1062 0.89
2 54 cg_sum 1193 48.16 11.99 30 30 below 131 0.11
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 54 ps_sum 1193 50.49 10.4 35 35 above 1053 0.88
2 54 ps_sum 1193 50.49 10.4 35 35 below 140 0.12print(list(
decided_cutoff(type = "10", quest = 60, domain = c_sum),
decided_cutoff(type = "10", quest = 60, domain = gm_sum),
decided_cutoff(type = "10", quest = 60, domain = fm_sum),
decided_cutoff(type = "10", quest = 60, domain = cg_sum),
decided_cutoff(type = "10", quest = 60, domain = ps_sum)
))[[1]]
quest domain n mean sd p10 cutoff class n.1 percent
1 60 c_sum 921 50.62 11.14 35 35 above 813 0.88
2 60 c_sum 921 50.62 11.14 35 35 below 108 0.12
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 60 gm_sum 921 51.06 10.76 35 35 above 809 0.88
2 60 gm_sum 921 51.06 10.76 35 35 below 112 0.12
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 60 fm_sum 921 48.26 14.94 25 25 above 811 0.88
2 60 fm_sum 921 48.26 14.94 25 25 below 110 0.12
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 60 cg_sum 921 51.05 10.8 35 35 above 824 0.89
2 60 cg_sum 921 51.05 10.8 35 35 below 97 0.11
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 60 ps_sum 921 49.95 12.24 35 35 above 777 0.84
2 60 ps_sum 921 49.95 12.24 35 35 below 144 0.16print(list(
decided_cutoff(type = "10", quest = 72, domain = c_sum),
decided_cutoff(type = "10", quest = 72, domain = gm_sum),
decided_cutoff(type = "10", quest = 72, domain = fm_sum),
decided_cutoff(type = "10", quest = 72, domain = cg_sum),
decided_cutoff(type = "10", quest = 72, domain = ps_sum)
))[[1]]
quest domain n mean sd p10 cutoff class n.1 percent
1 72 c_sum 791 50.71 12.2 35 35 above 687 0.87
2 72 c_sum 791 50.71 12.2 35 35 below 104 0.13
[[2]]
quest domain n mean sd p10 cutoff class n.1 percent
1 72 gm_sum 791 49.32 12.71 30 30 above 701 0.89
2 72 gm_sum 791 49.32 12.71 30 30 below 90 0.11
[[3]]
quest domain n mean sd p10 cutoff class n.1 percent
1 72 fm_sum 791 53.84 10.43 40 40 above 700 0.88
2 72 fm_sum 791 53.84 10.43 40 40 below 91 0.12
[[4]]
quest domain n mean sd p10 cutoff class n.1 percent
1 72 cg_sum 791 48.02 11.91 30 30 above 697 0.88
2 72 cg_sum 791 48.02 11.91 30 30 below 94 0.12
[[5]]
quest domain n mean sd p10 cutoff class n.1 percent
1 72 ps_sum 791 49.9 12.21 35 35 above 683 0.86
2 72 ps_sum 791 49.9 12.21 35 35 below 108 0.14ds_eligible %>%
select(quest,ends_with("_sum")) %>%
tableby(quest ~ ., control = tableby.control(numeric.stats=c("mean", "sd")), data = .) %>%
summary(. , digits = 2)ds_eligible %>% tabyl(gender)ds_eligible %>%
select(quest,ends_with("_sum")) %>%
group_by(quest) %>%
filter(n()>1) %>%
ungroup() %>%
pivot_longer(-quest) %>%
group_by(quest, name) %>%
mutate(
n=n(),
m=mean(value),
sd=sd(value),
monitor=m-sd,
below=m-2*sd)%>%
select(-value) %>%
distinct() %>% #remove duplicate
pivot_wider(id_cols = quest, names_from = name, values_from = n:below, names_glue = "{name}_{.value}") %>%
mutate_if(is.numeric, round, 2) %>%
.[gtools::mixedorder(.$quest), ]
#check which bases I can compare (2 groups)
bind_rows(
ds_eligible %>%
select(quest,ends_with("_sum")) %>%
group_by(quest) %>%
filter(n()>1) %>%
ungroup() %>%
mutate(base="eligible")
,
ds %>%
select(quest,ends_with("_sum")) %>%
mutate(base="original")
) %>%
group_by(quest) %>%
count(base)%>% #check each age inerval
filter(n()>1) %>% #remove if just one group
ungroup() %>%
select(quest) %>% distinct() %>%
pull(quest) -> questionnaires_to_compare
#define a df to compare
bind_rows(
ds_eligible %>%
select(quest,ends_with("_sum")) %>%
group_by(quest) %>%
filter(n()>1) %>%
ungroup() %>%
mutate(base="eligible")
,
ds %>%
select(quest,ends_with("_sum")) %>%
mutate(base="original")
) %>%
#filter those ages in which I have just one group
filter(quest %in% c(questionnaires_to_compare)) %>%
#need to have the long format to nest the all questionnaires
pivot_longer(c_sum:ps_sum) %>%
nest_by(quest,name) %>%
summarise(model = list(t.test(value ~ base, data = data, var.equal=T , alternative = "less"))) %>% #https://stackoverflow.com/questions/51074328/perform-several-t-tests-simultaneously-on-tidy-data-in-r
mutate(model = map(model, broom::tidy)) %>%
unnest(cols = c(model)) %>%
select(quest, name, estimate1, estimate2, p.value) %>%
mutate_if(is.numeric, round, 2)%>%
#presenting
pivot_wider(id_cols = quest, names_from =name, values_from = c(estimate1, estimate2, p.value), names_glue = "{name}_{.value}") %>%
select("quest", sort(colnames(.)))
#library(cutpointr)
set.seed(13)
bind_rows(
#typical
ds %>%
filter(quest == 24) %>%
select(c_sum) %>%
mutate(group = 0) %>%
sample_n(., 24)
,
#eligible
ds_eligible %>%
filter(quest == 36) %>%
mutate(group =1) %>%
select(group, c_sum)
) %>%
cutpointr(., c_sum, group,
pos_class = 1,
method = maximize_metric,
metric = youden) %>%
#plot_x()
#plot_roc() + geom_abline(slope = 1) + theme_bw()
summary(.)ds_retest_analysis ds_retest_analysis %>%
select_if(is.numeric) %>%
group_by(quest) %>%
nest() %>%
mutate(
correlations = map(data, corrr::correlate)
) %>%
unnest(correlations)# pearson correlation
retest_cor <- ds_retest_analysis %>%
split(list(.$quest)) %>%
map(~Hmisc::rcorr(as.matrix(.))$r)
# create a column with questionnaire
retest_cor <- do.call(rbind.data.frame, retest_cor)
# P vaues of pearson correlation
retest_cor_pval <- ds_retest_analysis %>%
split(list(.$quest)) %>%
map(~Hmisc::rcorr(as.matrix(.))$P)
# create a column with questionnaire
retest_cor_pval <- do.call(rbind.data.frame, retest_cor_pval)
test_rest_table <- left_join(
# R COEF
retest_cor %>%
select(-quest) %>% # remove questionnaires
rownames_to_column(var = "quest") %>% #add real questionnaires
separate(., col = "quest", into = c("quest","domain","sum", "time")) %>%#rename and separate
filter(time == "x") %>%
select(-sum, -time) %>%
select(quest, domain, ends_with("y")) %>% #select everything
pivot_longer(-c("quest","domain")) %>% #to match
mutate(name =str_extract_all(name, "\\w+(?=_)", simplify = T)) %>% #filter same domain
filter(domain == name)
,
# P values
retest_cor_pval %>%
select(-quest) %>% # remove questionnaires
rownames_to_column(var = "quest") %>% #add real questionnaires
separate(., col = "quest", into = c("quest","domain","sum", "time")) %>%#rename and separate
filter(time == "x") %>%
select(-sum, -time) %>%
select(quest, domain, ends_with("y")) %>% #select everything
pivot_longer(-c("quest","domain"))%>% #to match
mutate(name =str_extract_all(name, "\\w+(?=_)", simplify = T)) %>% #filter same domain
filter(domain == name) %>%
rename(pval=value) %>% mutate(pval = round(pval,3))
) %>%
mutate_at(vars(domain, name), ~str_replace(., "c","Communication") %>%
str_replace_all(., "gm","Gross Motor") %>%
str_replace_all(., "Communicationg","Problem Solving") %>% #gambiarra
str_replace_all(., "ps","Personal-Social") %>%
str_replace_all(., "fm","Fine Motor"))
left_join(test_rest_table,
ds_retest_analysis %>%
count(quest) %>% mutate(quest=as.character(quest)))ds_retest_analysis %>%
filter(quest == 4) %>%
{cor.test(.$c_sum.x, .$c_sum.y)}ds_retest_analysis %>%
filter(quest == 16) %>%
{cor.test(.$fm_sum.x, .$fm_sum.y)}test_rest_table %>%
filter(pval <= 0.05) %>%
tableby(domain~value, data = .) %>%
summary()test_rest_table %>%
mutate(quest= as.numeric(quest)) %>%
filter(pval <= 0.05) %>%
ggplot(aes(x = quest, y = value)) +
geom_point(aes(color = domain), alpha = 0.5, show.legend = FALSE) +
#geom_smooth(method = "lm", color = "darkgray", se = FALSE) +
facet_wrap(. ~ domain, ncol = 2) +
labs(x="", y = "r") +
ylim(0,1)+
theme_bw()ds_rater_analysis ds_rater_analysis %>%
select_if(is.numeric) %>%
group_by(quest) %>%
nest() %>%
mutate(
correlations = map(data, corrr::correlate)
) %>%
unnest(correlations)ds_rater_analysis %>%
count(quest)# pearson correlation
rater_cor <- ds_rater_analysis %>%
split(list(.$quest)) %>%
keep(~nrow(.) > 4) %>%
map(~Hmisc::rcorr(as.matrix(.))$r)
# create a column with questionnaire
rater_cor <- do.call(rbind.data.frame, rater_cor)
# P vaues of pearson correlation
rater_cor_pval <- ds_rater_analysis %>%
split(list(.$quest)) %>%
keep(~nrow(.) > 4) %>%
map(~Hmisc::rcorr(as.matrix(.))$P)
# create a column with questionnaire
rater_cor_pval <- do.call(rbind.data.frame, rater_cor_pval)
rater_table <- left_join(
# R COEF
rater_cor %>%
select(-quest) %>% # remove questionnaires
rownames_to_column(var = "quest") %>% #add real questionnaires
separate(., col = "quest", into = c("quest","domain","sum", "time")) %>%#rename and separate
filter(time == "x") %>%
select(-sum, -time) %>%
select(quest, domain, ends_with("y")) %>% #select everything
pivot_longer(-c("quest","domain")) %>% #to match
mutate(name =str_extract_all(name, "\\w+(?=_)", simplify = T)) %>% #filter same domain
filter(domain == name)
,
# P values
rater_cor_pval %>%
select(-quest) %>% # remove questionnaires
rownames_to_column(var = "quest") %>% #add real questionnaires
separate(., col = "quest", into = c("quest","domain","sum", "time")) %>%#rename and separate
filter(time == "x") %>%
select(-sum, -time) %>%
select(quest, domain, ends_with("y")) %>% #select everything
pivot_longer(-c("quest","domain"))%>% #to match
mutate(name =str_extract_all(name, "\\w+(?=_)", simplify = T)) %>% #filter same domain
filter(domain == name) %>%
rename(pval=value) %>% mutate(pval = round(pval,3))
) %>%
mutate_at(vars(domain, name), ~str_replace(., "c","Communication") %>%
str_replace_all(., "gm","Gross Motor") %>%
str_replace_all(., "Communicationg","Problem Solving") %>% #gambiarra
str_replace_all(., "ps","Personal-Social") %>%
str_replace_all(., "fm","Fine Motor"))
left_join(rater_table,
ds_rater_analysis %>%
count(quest) %>% mutate(quest=as.character(quest))) %>% write.csv(., "icc.csv")ds_rater_analysis %>%
filter(quest == 72) %>%
{cor.test(.$gm_sum.x, .$gm_sum.y)}ds_rater_analysis %>%
filter(quest == 16) %>%
{cor.test(.$fm_sum.x, .$fm_sum.y)}rater_table %>%
filter(pval <= 0.05) %>%
tableby(domain~value, data = .) %>%
summary()rater_table %>%
mutate(quest= as.numeric(quest)) %>%
filter(pval <= 0.05) %>%
ggplot(aes(x = quest, y = value)) +
geom_point(aes(color = domain), alpha = 0.5, show.legend = FALSE) +
#geom_smooth(method = "lm", color = "darkgray", se = FALSE) +
facet_wrap(. ~ domain, ncol = 2) +
labs(x="", y = "r") +
ylim(0,1)+
theme_bw()#mirt(data = ds_com_2[,-c(1:2)], model = 1, itemtype = "graded") %>%
# itemplot(., 1)
#https://groups.google.com/g/mirt-package/c/V0AX2aIXS10
plogis(3.119, location = -5)
plogis(3.119)-plogis(1.286)
help(plogis)library(mirt)
apply_irt_cfa <- function(data) {
#IRT
set.seed(123)
mod_irt <- mirt(data = data, model = 1, itemtype = "graded")
mod_coef <- coef(mod_irt, IRTpars = T, simplify = T) #get classical IRT parameterization
mod_fit <- M2(mod_irt,na.rm=TRUE)
mod_plot_trace <- plot(mod_irt, type = "trace")
#CFA
library(lavaan)
mod_cfa <- cfa(model = paste("f1=~", paste(names(data), collapse=" + ")),
data=data,
estimator = 'WLSM',
ordered=names(data))
mod_cfa_result <- summary(mod_cfa, standardized=TRUE, fit.measures = TRUE)
#Return
return(list(mod_plot_trace,mod_fit, mod_coef))
}ds %>% tabyl(summative_risk) %>% adorn_totals("row") %>%
adorn_pct_formatting(digits = 2)level_order <- c('Communication', 'Gross Motor', 'Fine Motor',"Problem Solving", "Personal-Social")
ds %>%
filter(!is.na(summative_risk)) %>% #don't use
filter(quest !=9 & quest != 72) %>%
select(quest, summative_risk, c_sum:ps_sum) %>%
mutate(summative_risk = as.numeric(as.character(summative_risk))) %>%
pivot_longer(cols = -c(quest, summative_risk))%>%
group_by(quest, name, summative_risk) %>%
nest() %>%
mutate(mean = map_dbl(data, ~mean(.x$value))) %>%
mutate(summative_risk = as.factor(summative_risk)) %>%
#plot
mutate(name = case_when(
name == "c_sum" ~ "Communication",
name == "gm_sum" ~ "Gross Motor",
name == "fm_sum" ~ "Fine Motor",
name == "cg_sum" ~ "Problem Solving",
name == "ps_sum" ~ "Personal-Social")) %>%
ggplot(.,
aes(x=factor(name, levels = level_order), y=mean, group = summative_risk, fill=summative_risk)) +
stat_summary(fun.y=mean,position=position_dodge(width=0.95),geom="bar") +
stat_summary(fun.data=mean_cl_normal,position=position_dodge(0.95),geom="errorbar") +
labs(x = "Domain", y = "Mean scores", fill = "Risk") +
theme_bw() #+ facet_wrap(~quest)Nice plot but not used
ds %>%
filter(!is.na(summative_risk)) %>%
filter(quest !=9 & quest != 72) %>%
select(quest, summative_risk, c_sum:ps_sum) %>%
mutate(summative_risk = as.numeric(as.character(summative_risk))) %>%
pivot_longer(cols = -c(quest, summative_risk))%>%
group_by(quest, name, summative_risk) %>%
nest() %>%
mutate(mean = map_dbl(data, ~mean(.x$value))) %>%
mutate(summative_risk = as.factor(summative_risk)) %>%
ggplot(.,
aes(x=quest, y=mean, group = interaction(name ,summative_risk), color=name)) +
#stat_summary(geom="line", size=1.5, aes(linetype=summative_risk)) +
stat_summary(fun.y=mean,position=position_dodge(width=0.95),geom="bar", aes(fill = summative_risk)) +
stat_summary(geom="errorbar", size=0.2, width = .2)
ds %>%
filter(!is.na(summative_risk)) %>% #don`t use missing cases on risk
filter(summative_risk %in% c(0,3)) %>% #extreme groups (no risk vs high risk)
mutate(summative_risk = if_else(summative_risk == 0,"Nonrisk","Risk")) %>%
select(quest, summative_risk, ends_with("sum"))%>% #get all variables
pivot_longer(c_sum:ps_sum) %>% #tranform to long format
mutate(name = case_when(
name == "c_sum" ~ "Communication",
name == "gm_sum" ~ "Gross Motor",
name == "fm_sum" ~ "Fine Motor",
name == "cg_sum" ~ "Problem Solving",
name == "ps_sum" ~ "Personal-Social")) %>%
#plot
ggplot(., aes(x = quest, y = value, group = summative_risk)) +
stat_summary(geom = "line", fun = mean, aes(linetype = summative_risk), size=1) +
theme_bw() +
labs(x = "", y = "Mean scores", linetype = "Risk group") +
facet_wrap(~name) +
theme(legend.position = "bottom")
#stat_summary(geom="errorbar", size=0.1, width = .2)ds %>%
filter(quest != 9) %>% #no risk here
filter(!is.na(summative_risk)) %>% #don`t use missing cases on risk
filter(summative_risk %in% c(0,3)) %>% #extreme groups (no risk vs high risk)
mutate(summative_risk = if_else(summative_risk == 0,"Nonrisk","Risk")) %>%
select(quest, summative_risk, ends_with("sum"))%>% #get all variables
pivot_longer(c_sum:ps_sum) %>%
nest_by(name,quest) %>% #group!!!
mutate(model = list(t.test(value ~ summative_risk, data = data, var.equal=T)$p.value)) %>% #compute p values
unnest(model) %>%
#split(.$quest)%>% #ungroup
filter(model <= 0.05) %>%
pivot_wider(id_cols = quest, names_from = name, values_from = model) %>% #unnest based on p values
#present
mutate_if(is.numeric, round, 2) %>%
arrange(quest)
ds %>%
filter(!is.na(summative_risk)) %>%
filter(summative_risk %in% c(0,3)) %>% #extreme groups
mutate(summative_risk = if_else(summative_risk == 0,"Nonrisk","Risk")) %>%
ggplot(., aes(x=quest, y = c_sum, linetype=summative_risk, group=summative_risk)) +
stat_summary(geom = "line", fun = mean, size=1) +
stat_summary(geom="errorbar", size=0.1, width = .2)+
theme_bw() +
labs(x = "", y = "Mean scores", title = "Communication", linetype = "Risk group") +
theme(legend.position = "bottom")ds %>%
filter(!is.na(summative_risk)) %>%
filter(summative_risk %in% c(0,3)) %>% #extreme groups
mutate(summative_risk = if_else(summative_risk == 0,"Nonrisk","Risk")) %>%
ggplot(., aes(x=quest, y = gm_sum, linetype=summative_risk, group=summative_risk)) + #attention here!!!
stat_summary(geom = "line", fun = mean, size=1) +
stat_summary(geom="errorbar", size=0.1, width = .2)+
theme_bw() +
labs(x = "", y = "Mean scores", title = "Gross Motor", linetype = "Risk group") +
theme(legend.position = "bottom")ds %>%
filter(!is.na(summative_risk)) %>%
filter(summative_risk %in% c(0,3)) %>% #extreme groups
mutate(summative_risk = if_else(summative_risk == 0,"Nonrisk","Risk")) %>%
ggplot(., aes(x=quest, y = fm_sum, linetype=summative_risk, group=summative_risk)) + ## Fine motor!!
stat_summary(geom = "line", fun = mean, size=1) +
stat_summary(geom="errorbar", size=0.1, width = .2)+
theme_bw() +
labs(x = "", y = "Mean scores", title = "Fine Motor", linetype = "Risk group") +
theme(legend.position = "bottom")ds %>%
filter(!is.na(summative_risk)) %>%
filter(summative_risk %in% c(0,3)) %>% #extreme groups
mutate(summative_risk = if_else(summative_risk == 0,"Nonrisk","Risk")) %>%
ggplot(., aes(x=quest, y = cg_sum, linetype=summative_risk, group=summative_risk)) + ## Problem Solving = Cognition
stat_summary(geom = "line", fun = mean, size=1) +
stat_summary(geom="errorbar", size=0.1, width = .2)+
theme_bw() +
labs(x = "", y = "Mean scores", title = "Problem Solving", linetype = "Risk group") +
theme(legend.position = "bottom")ds %>%
filter(!is.na(summative_risk)) %>%
select(quest, summative_risk, c_sum:ps_sum) %>%
pivot_longer(c_sum:ps_sum) %>%
group_by(quest,summative_risk,name) %>%
summarise(mean=mean(value, na.rm=T), sd = sd(value, na.rm = T)) %>%
#first level
pivot_wider(id_cols = quest, names_from = summative_risk:name, values_from = mean:sd, names_glue = "{name}_{summative_risk}_{.value}") %>%
#second level
pivot_longer(cols = -c(quest)) %>%
arrange(quest,name) %>%
separate(name, into = c("domain","constant","risk","result")) %>%
#third level
select(-constant) %>%
pivot_wider(id_cols = quest, names_from = domain:result, values_from = value) %>%
mutate_if(is.numeric, round, 2) ds %>%
filter(!is.na(summative_risk)) %>%
filter(summative_risk %in% c(0,3)) %>% #contrasting groups T TEST HERE
filter(quest !=9) %>%
select(quest, summative_risk, c_sum:ps_sum) %>%
pivot_longer(c_sum:ps_sum) %>%
group_by(quest,name) %>%
nest() %>%
mutate(p_risk = map(data, ~aov(value ~ summative_risk, data = .) %>% {summary(.)[[1]][["Pr(>F)"]][1]})) %>%
unnest_wider(p_risk) %>%
rename(p_val = "...1") %>%
pivot_wider(id_cols = quest, names_from = name, values_from = p_val) %>%
mutate_if(is.numeric, round, 3)
ds %>%
filter(quest == 2) %>%
filter(summative_risk %in% c(0,3)) %>%
t.test(c_sum ~ summative_risk, data = ., var.equal=T)
#aov(c_sum ~ summative_risk, data = .) %>% summary()library(ggpubr)
ds %>%
filter(!is.na(summative_risk)) %>%
filter(quest != "9") %>%
filter(summative_risk %in% c(0,3)) %>% #extreme groups
mutate(summative_risk = if_else(summative_risk == 0,"Nonrisk","Risk")) %>%
select(quest, summative_risk, c_sum)%>%
ggboxplot(., x = "summative_risk", y = "c_sum",
facet.by = "quest", short.panel.labs = FALSE) +
stat_compare_means(label = "p.format", method = "t.test",label.x = c("Risk"), label.y = 1.5)mod_cfa <- 'com =~ com_a4_1 + com_a4_2 + com_a4_3 + com_a4_4 + com_a4_5 + com_a4_6
gross =~ gm_a4_1 + gm_a4_2 + gm_a4_3 + gm_a4_4 + gm_a4_5 + gm_a4_6
fine =~ fm_a4_1 + fm_a4_2 + fm_a4_3 + fm_a4_4 + fm_a4_5 + fm_a4_6
cog =~ cg_a4_1 + cg_a4_2 + cg_a4_3 + cg_a4_4 + cg_a4_5 + cg_a4_6
ps =~ ps_a4_1 + ps_a4_2 + ps_a4_3 + ps_a4_4 + ps_a4_5 + ps_a4_6
'mod_cfa_asq3 <- 'com =~ com_a3_1 + com_a3_2 + com_a3_3 + com_a3_4 + com_a3_5 + com_a3_6
gross =~ gm_a3_1 + gm_a3_2 + gm_a3_3 + gm_a3_4 + gm_a3_5 + gm_a3_6
fine =~ fm_a3_1 + fm_a3_2 + fm_a3_3 + fm_a3_4 + fm_a3_5 + fm_a3_6
cog =~ cg_a3_1 + cg_a3_2 + cg_a3_3 + cg_a3_4 + cg_a3_5 + cg_a3_6
ps =~ ps_a3_1 + ps_a3_2 + ps_a3_3 + ps_a3_4 + ps_a3_5 + ps_a3_6
'Dataset
ds_2_full <- ds %>%
filter(quest == 2) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
ds_2_full %>% mutate_all(replace_na,-99)%>%
write.table(., file="ds_2_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_2 <- cfa(model = mod_cfa,
data=ds_2_full,
estimator = 'WLSM',
ordered=names(ds_2_full))CFA results
summary(mod_cfa_2, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_2, what = "est")$psi)cfa(model = mod_cfa_asq3,
data =
ds_final_merged %>%
filter(quest == 2) %>%
select(
com_a3_1, com_a3_2, com_a3_3, com_a3_4, com_a3_5, com_a3_6,
gm_a3_1, gm_a3_2, gm_a3_3, gm_a3_4, gm_a3_5, gm_a3_6,
fm_a3_1, fm_a3_2, fm_a3_3, fm_a3_4, fm_a3_5, fm_a3_6,
cg_a3_1, cg_a3_2, cg_a3_3, cg_a3_4, cg_a3_5, cg_a3_6,
ps_a3_1, ps_a3_2, ps_a3_3, ps_a3_4, ps_a3_5, ps_a3_6)
,
estimator = 'WLSM',
ordered=names(
ds_final_merged %>%
filter(quest == 2) %>%
select(
com_a3_1, com_a3_2, com_a3_3, com_a3_4, com_a3_5, com_a3_6,
gm_a3_1, gm_a3_2, gm_a3_3, gm_a3_4, gm_a3_5, gm_a3_6,
fm_a3_1, fm_a3_2, fm_a3_3, fm_a3_4, fm_a3_5, fm_a3_6,
cg_a3_1, cg_a3_2, cg_a3_3, cg_a3_4, cg_a3_5, cg_a3_6,
ps_a3_1, ps_a3_2, ps_a3_3, ps_a3_4, ps_a3_5, ps_a3_6)
)) %>% summary(., standardized=TRUE, fit.measures = TRUE)
Dataset
ds_4_full <- ds %>% #attention here to specify the correct dataset
filter(quest == 4) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_4_full, file="ds_4_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_4 <- cfa(model = mod_cfa,
data=ds_4_full,
estimator = 'WLSM',
ordered=names(ds_4_full))CFA results
summary(mod_cfa_4, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_4, what = "est")$psi)Dataset
ds_6_full <- ds %>%
filter(quest == 6) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_6_full, file="ds_6_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) CFA model
mod_cfa_6 <- cfa(model = mod_cfa,
data=ds_6_full,
estimator = 'WLSM',
ordered=names(ds_6_full))CFA results
summary(mod_cfa_6, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_6, what = "est")$psi)Dataset
ds_8_full <- ds %>%
filter(quest == 8) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_8_full, file="ds_8_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) CFA model
mod_cfa_8 <- cfa(model = mod_cfa,
data=ds_8_full,
estimator = 'WLSM',
ordered=names(ds_8_full))CFA results
summary(mod_cfa_8, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_8, what = "est")$psi)Dataset
ds_10_full <- ds %>%
filter(quest == 10) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_10_full, file="ds_10_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_10 <- cfa(model = mod_cfa,
data=ds_10_full,
estimator = 'WLSM',
ordered=names(ds_10_full))CFA results
summary(mod_cfa_10, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_10, what = "est")$psi)Dataset
ds_12_full <- ds %>%
filter(quest == 12) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_12_full, file="ds_12_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_12 <- cfa(model = mod_cfa,
data=ds_12_full,
estimator = 'WLSM',
ordered=names(ds_12_full))CFA results
summary(mod_cfa_12, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_12, what = "est")$psi)Dataset
ds_14_full <- ds %>%
filter(quest == 14) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_14_full, file="ds_14_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) CFA model
mod_cfa_14 <- cfa(model = mod_cfa,
data=ds_14_full,
estimator = 'WLSM',
ordered=names(ds_14_full))CFA results
summary(mod_cfa_14, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_14, what = "est")$psi)Dataset
ds_16_full <- ds %>%
filter(quest == 16) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_16_full, file="ds_16_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_16 <- cfa(model = mod_cfa,
data=ds_16_full,
estimator = 'WLSM',
ordered=names(ds_16_full))CFA results
summary(mod_cfa_16, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_16, what = "est")$psi)Dataset
ds_18_full <- ds %>%
filter(quest == 18) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_18_full, file="ds_18_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_18 <- cfa(model = mod_cfa,
data=ds_18_full,
estimator = 'WLSM',
ordered=names(ds_18_full))CFA results
summary(mod_cfa_18, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_18, what = "est")$psi)Dataset
ds_20_full <- ds %>%
filter(quest == 20) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_20_full, file="ds_20_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_20 <- cfa(model = mod_cfa,
data=ds_20_full,
estimator = 'WLSM',
ordered=names(ds_20_full))CFA results
summary(mod_cfa_20, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_20, what = "est")$psi)Dataset
ds_22_full <- ds %>%
filter(quest == 22) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_22_full, file="ds_22_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_22 <- cfa(model = mod_cfa,
data=ds_22_full,
estimator = 'WLSM',
ordered=names(ds_22_full))CFA results
summary(mod_cfa_22, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_22, what = "est")$psi)Dataset
ds_24_full <- ds %>%
filter(quest == 24) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_24_full, file="ds_24_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_24 <- cfa(model = mod_cfa,
data=ds_24_full,
estimator = 'WLSM',
ordered=names(ds_24_full))CFA results
summary(mod_cfa_24, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_24, what = "est")$psi)Dataset
ds_27_full <- ds %>%
filter(quest == 27) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_27_full, file="ds_27_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) CFA model
mod_cfa_27 <- cfa(model = mod_cfa,
data=ds_27_full,
estimator = 'WLSM',
ordered=names(ds_27_full))CFA results
summary(mod_cfa_27, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_27, what = "est")$psi)Dataset
ds_30_full <- ds %>%
filter(quest == 30) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_30_full, file="ds_33_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) CFA model
mod_cfa_30 <- cfa(model = mod_cfa,
data=ds_30_full,
estimator = 'WLSM',
ordered=names(ds_30_full))CFA results
summary(mod_cfa_30, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_30, what = "est")$psi)Dataset
ds_33_full <- ds %>%
filter(quest == 33) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_33_full, file="ds_33_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) CFA model
mod_cfa_33 <- cfa(model = mod_cfa,
data=ds_33_full,
estimator = 'WLSM',
ordered=names(ds_33_full))CFA results
summary(mod_cfa_33, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_33, what = "est")$psi)Dataset
ds_36_full <- ds %>%
filter(quest == 36) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_36_full, file="ds_36_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) CFA model
mod_cfa_36 <- cfa(model = mod_cfa,
data=ds_36_full,
estimator = 'WLSM',
ordered=names(ds_36_full))CFA results
summary(mod_cfa_36, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_36, what = "est")$psi)Dataset
ds_42_full <- ds %>%
filter(quest == 42) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_42_full, file="ds_42_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) CFA model
mod_cfa_42 <- cfa(model = mod_cfa,
data=ds_42_full,
estimator = 'WLSM',
ordered=names(ds_42_full))CFA results
summary(mod_cfa_42, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_42, what = "est")$psi)Dataset
ds_48_full <- ds %>%
filter(quest == 48) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_48_full, file="ds_48_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) CFA model
mod_cfa_48 <- cfa(model = mod_cfa,
data=ds_48_full,
estimator = 'WLSM',
ordered=names(ds_48_full))CFA results
summary(mod_cfa_48, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_48, what = "est")$psi)Dataset
ds_54_full <- ds %>%
filter(quest == 54) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_54_full, file="ds_54_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_54 <- cfa(model = mod_cfa,
data=ds_54_full,
estimator = 'WLSM',
ordered=names(ds_54_full))CFA results
summary(mod_cfa_54, standardized=TRUE, fit.measures = TRUE) Correlations
cov2cor(inspect(mod_cfa_54, what = "est")$psi)Dataset
ds_60_full <- ds %>%
filter(quest == 60) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_60_full, file="ds_60_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE)
CFA model
mod_cfa_60 <- cfa(model = mod_cfa,
data=ds_60_full,
estimator = 'WLSM',
ordered=names(ds_60_full))CFA results
summary(mod_cfa_60, standardized=TRUE, fit.measures = TRUE) Correlation s
cov2cor(inspect(mod_cfa_60, what = "est")$psi)ds_72_full <- ds %>%
filter(quest == 72) %>%
select(
com_a4_1, com_a4_2, com_a4_3, com_a4_4, com_a4_5, com_a4_6,
gm_a4_1, gm_a4_2, gm_a4_3, gm_a4_4, gm_a4_5, gm_a4_6,
fm_a4_1, fm_a4_2, fm_a4_3, fm_a4_4, fm_a4_5, fm_a4_6,
cg_a4_1, cg_a4_2, cg_a4_3, cg_a4_4, cg_a4_5, cg_a4_6,
ps_a4_1, ps_a4_2, ps_a4_3, ps_a4_4, ps_a4_5, ps_a4_6)
#write.table(ds_72_full, file="ds_72_full.dat", row.names=FALSE, col.names=FALSE, sep="\t", quote=FALSE) mod_cfa_72 <- cfa(model = mod_cfa,
data=ds_72_full,
estimator = 'WLSM',
ordered=names(ds_72_full))summary(mod_cfa_72, standardized=TRUE, fit.measures = TRUE) cov2cor(inspect(mod_cfa_72, what = "est")$psi)set.seed(123)
ds_1 %>%
select(quest,id,com_a4_1:com_a4_6) %>%
group_split(quest) %>%
map(. %>%
sample_n(500)) -> xfor (i in 1:length(x)) {
assign(paste0("ds_com_", unique(x[[i]][1])), as.data.frame(x[[i]]))
}apply_irt_cfa(ds_com_2[,-c(1:2)])ds_final_merged %>%
filter(quest == 2) %>%
select(
com_a3_1, com_a3_2, com_a3_3, com_a3_4, com_a3_5, com_a3_6,
gm_a3_1, gm_a3_2, gm_a3_3, gm_a3_4, gm_a3_5, gm_a3_6,
fm_a3_1, fm_a3_2, fm_a3_3, fm_a3_4, fm_a3_5, fm_a3_6,
cg_a3_1, cg_a3_2, cg_a3_3, cg_a3_4, cg_a3_5, cg_a3_6,
ps_a3_1, ps_a3_2, ps_a3_3, ps_a3_4, ps_a3_5, ps_a3_6) %>%
apply_irt_cfa(.)apply_irt_cfa(ds_com_4[,-c(1:2)])apply_irt_cfa(ds_com_6[,-c(1:2)])apply_irt_cfa(ds_com_8[,-c(1:2)])apply_irt_cfa(ds_com_10[,-c(1:2)])apply_irt_cfa(ds_com_12[,-c(1:2)])apply_irt_cfa(ds_com_14[,-c(1:2)])apply_irt_cfa(ds_com_16[,-c(1:2)])apply_irt_cfa(ds_com_18[,-c(1:2)])apply_irt_cfa(ds_com_20[,-c(1:2)])apply_irt_cfa(ds_com_22[,-c(1:2)])apply_irt_cfa(ds_com_24[,-c(1:2)])apply_irt_cfa(ds_com_27[,-c(1:2)])apply_irt_cfa(ds_com_30[,-c(1:2)])apply_irt_cfa(ds_com_33[,-c(1:2)])apply_irt_cfa(ds_com_36[,-c(1:2)])apply_irt_cfa(ds_com_42[,-c(1:2)])apply_irt_cfa(ds_com_48[,-c(1:2)])apply_irt_cfa(ds_com_54[,-c(1:2)])apply_irt_cfa(ds_com_60[,-c(1:2)])apply_irt_cfa(ds_com_72[,-c(1:2)])ds_1 %>% names