Data Analysis of Psychology Meditation Study
Charles_Li
2022-08-14
Tables(except descriptive statistics containing confidence intervals) and figures (except interactive figures) are in APA format. Some errors are caused by Rstudio when displaying latex code. The APA tables and figures are currently in html format.
We only use non-parametric tests
AEQ scores are essentially in ordinal scale
Chi-squared test of independence (nominal and ordinal (except AEQ score))
- Cramer’s V (effect size): examine the association between two categorical variables when at least one of the variables has more than two values (i.e. more than a \(2 \times 2\) contingency table).
Kruskal Wallis test and post hoc (pairwise Wilcox test) (AEQ score and all other numeric data)
- Spearman’s rank or Kendall’s tau (correlation coefficient)
- Eta-squared(\(\eta^2\)) (effect size)
The effect sizes of all hypothesis tests except the post-hoc test of Chi-squared test are included.
Chi-squared – Cramer’s V or Cramer’s \(\phi\) (\(\phi_c\))
\(\phi_c \in [0.1,0.3]\): weak association
\(\phi_c \in [0.4,0.5]\): medium association
\(\phi_c>0.5\): strong association
There is a disagreement between the criteria of Cramer’s V https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6107969/table/tbl2/?report=objectonly
Kruskal Wallis Test by Rank – Eta-squared (\(\eta^2\))
- \(\eta^2 = 0.01\) indicates a small effect;
- \(\eta^2 = 0.06\) indicates a medium effect;
- \(\eta^2 = 0.14\) indicates a large effect.
Pairwise Wilcoxon Rank Sum Test – (Cohen) \(r = z/\sqrt{n}\); z statistic; n: sample size number of pairs
- 0.1 (small effect),
- 0.3 (moderate effect) and
- 0.5 and above (large effect).
- P adjustment method “BH” is used.
\(\phi_c\), \(\eta^2\), and \(r\) of pairwise Wilcoxon are incomparable.
The effect size of Wilcoxon test is not explained well in the original source, but this method of calculation is provided by both R and SPSS:
Rosenthal, R. (1991). Applied Social Research Methods: Meta-analytic procedures for social research Thousand Oaks, CA: SAGE Publications Ltd doi: 10.4135/9781412984997
Spearman Correlation – \(\rho\) (Coefficient)
Interpretation of the Pearson’s and Spearman’s correlation coefficients.
Akoglu H. (2018). User’s guide to correlation coefficients. Turkish journal of emergency medicine, 18(3), 91–93. https://doi.org/10.1016/j.tjem.2018.08.001
You can view significant data by limiting p value and effect size (e.g. p <0.05, effect.size \(\in [0.4,\infty]\)). All of the interpretation reports generated are in APA format.
Example_APA_format_Kruskal
# function that converts group notations (P,S,C) to more readable texts
group.legend = function(data){
new = data %>%
mutate(group = case_when(
group == "control" ~ "Control",
group == "single" ~ "20 minuted 1X/day",
TRUE ~ "10 minutes 2X/day"
),
group = factor(group, levels = c("Control","20 minuted 1X/day","10 minutes 2X/day"))) %>%
return(new)
}
# remove leading 0 from tables, strings, including negative sign
rm_0 <- function(id){
return(
ifelse(id == 0,
"0",ifelse(abs(id) < 1,
ifelse(id > 0,sub("^0+", "", id),
paste0("-",sub("^-0+", "", id))),
id))
)
}
p.format <- function(p, sig = 0.05, table = TRUE){
# APA format of p values:
# 1. 3 decimal places and no leading zeros for values less than 1
# 2. use "p < .001"
if(table){
str <- ifelse(
is.nan(p) | is.na(p),"NA",ifelse(p < 0.001,"< .001",rm_0(round(p, digits = 3)))
)
}
else{
str <- ifelse(
is.nan(p) | is.na(p),"NA",ifelse(p < 0.001,"p < .001",sprintf("p = %s",rm_0(round(p, digits = 3))))
)
}
return(str)
}
num.format <- function(n){
# APA format of numeric values:
# 1. 2 decimal places and no leading zeros for values less than 1
str <- ifelse(
is.nan(n) | is.na(n),"NA",ifelse(abs(n) < 1,rm_0(round(n, digits = 2)), round(n, digits = 2))
)
return(str)
}
# apply change of group names
i = 1
for(data in datasets){
datasets[[i]] <- group.legend(data)
i = i+1
}
# APA formatting variables (no CI)
## figure names
test.type <- c("chisq","kruskal", "wilcox")
var.name <- c("demographic", "productivity", "grade_hours","AEQ", "confounding")
var.type <- c("cat","num")
table.names <- c("chisq.demographic","chisq.demographic.post", "demographic.cat", # 1.1 Group differences in all demographic variables, categorical
"kruskal.demographic", # 1.2 Group differences in all demographic variables, numeric
"productivity", # 2.1 productivity
"kruskal.grade_hours","wilcox.grade", "grade_hours", # 2.2 exam score and total hours of studying
"kruskal.AEQ", "wilcox.AEQ", "AEQ", # 3.0 AEQs
"spearman", # 5.1 correlation
"chisq.confounding", "chisq.confounding.post","confounding.cat" ,# 5.2 confounding, categorical
"kruskal.confounding", "wilcox.confounding","confounding.num" # 5.3 confounding, categorical
) # object names
APA.table.names <- paste("APA", table.names,sep = ".")
# Change Table number or orders here, maybe use a excel file later
names(APA.table.names) <- c("Table 1", "Table 2","Table 1 and 2",
"Table 3",
"Table 4",
"Table 5", "Table 6", "Table 5 and 6",
"Table 7", "Table 8", "Table 7 and 8",
"Table 9",
"Table 10", "Table 11", "Table 10 and 11",
"Table 12", "Table 13", "Table 12 and 13"
) # store table numbers
# Change Table names here, maybe use a excel file later
APA.table.titles <- c("Chi-squared Test Results of Meditation Group vs Demographic Factors (Categorical)", "Frequency Table and Post-hoc Chi-squared Test Result of Meditation Group vs Demographic Factors (Categorical)",
"Kruskal-Wallis Test Results of Meditation Group vs Demographic Factors (Numeric)",
"Chi-square Test Results of Meditation Group vs Level of Productivity on Day 6 and 7",
"Kruskal-Wallis Test Results of Meditation Group vs Exam Score and Total Hours of Studying", "Post-hoc Wilcoxon Test Results of Meditation Group vs Exam Score",
"Kruskal-Wallis Test Results of Meditation Group vs AEQ Score", "Post-hoc Wilcoxon Test Results of Meditation Group vs AEQ Score",
"Spearman Correlation between Baseline and Outcome Variables",
"Chi-squared Test Results of Nominal (Except Meditation group) vs Ordinal Variables", "Frequency Table and Post-hoc Chi-squared Test Result of Nominal (Except Meditation group) vs Ordinal Variables",
"Kruskal-Wallis Test Results of Ordinal vs Numeric Variables", "Post-hoc Wilcoxon Test Results of Ordinal vs Numeric Variables"
) # store table titles
names(APA.table.titles) = paste("Table", c(1:13)) # store table numbers
# APA.table.notes # add manually
# APA figures
APA.figure.names <- paste("APA",c("AEQ.fig", "mood.fig", "spearman.fig"),sep = ".") # object names
names(APA.figure.names) <- paste("Figure", c(1:3))# store figure numbers
APA.figure.titles <- c("AEQ Scores","Mood Scores", "Spearman Correlation")
names(APA.figure.titles) <- paste("Figure", c(1:3))# store figure numbers
cat("APA table names:\n")## APA table names:print(APA.table.names)## Table 1 Table 2
## "APA.chisq.demographic" "APA.chisq.demographic.post"
## Table 1 and 2 Table 3
## "APA.demographic.cat" "APA.kruskal.demographic"
## Table 4 Table 5
## "APA.productivity" "APA.kruskal.grade_hours"
## Table 6 Table 5 and 6
## "APA.wilcox.grade" "APA.grade_hours"
## Table 7 Table 8
## "APA.kruskal.AEQ" "APA.wilcox.AEQ"
## Table 7 and 8 Table 9
## "APA.AEQ" "APA.spearman"
## Table 10 Table 11
## "APA.chisq.confounding" "APA.chisq.confounding.post"
## Table 10 and 11 Table 12
## "APA.confounding.cat" "APA.kruskal.confounding"
## Table 13 Table 12 and 13
## "APA.wilcox.confounding" "APA.confounding.num"########################################################################
# prevent rendering of latex code
latex.escape <- function(str, report = FALSE,escape = escape.format){ # string vector or a string
if(escape == TRUE){
if(report == TRUE){
str <- gsub(pattern="\\", replacement="\\\\", str,fixed = T) #
}
else{
str <- gsub(pattern="\\", replacement="\\\\", str, fixed = T)
}
str <- str_replace_all(str, pattern=fixed("$"), replacement="\\$")
# str <- knitr:::escape_latex(str)
# gsub(pattern="$", replacement="\\$", str) # add escape char \ in front of $
# }
}
return(str)
}
# gsub(pattern="\\$", replacement="\\\\$", c("", "$df$", "$\\chi^2$", "$\\phi_c$", "$p$", "$N$"))
latex.escape(c("", "$df$", "$\\chi^2$", "$\\phi_c$", "$p$", "$N$"))## [1] "" "$df$" "$\\chi^2$" "$\\phi_c$" "$p$" "$N$"latex.escape("For the first Chi-squared test, only variable relationships which $p \\le .05$ are included. For the post-hoc Chi-squared test, only variable relationships which $p_{post-hoc} \\le .05$ are included. Some relationships shows significance in the first test while showing no significance in the post-hoc test. Please compare the number of incidences (N) and the proportion in the post-hoc test to analyze the direction of the significance.", report = TRUE)## [1] "For the first Chi-squared test, only variable relationships which $p \\le .05$ are included. For the post-hoc Chi-squared test, only variable relationships which $p_{post-hoc} \\le .05$ are included. Some relationships shows significance in the first test while showing no significance in the post-hoc test. Please compare the number of incidences (N) and the proportion in the post-hoc test to analyze the direction of the significance."# report.demographic.cat
#
#
# cat(report.demographic.cat)
########################################################################
APA.tbl.caption <- function(name, format = output.format){ # enter object name
number <- names(which(APA.table.names == name)) # enter object name
title <- APA.table.titles[which(names(APA.table.titles) == number)]
if(format == "html"){
caption <- sprintf("<b>%s</b><br><i>%s</i>", number, title)
}
else if(format == "latex"){
caption <- sprintf("\\textbf{%s}\\\\\\emph{%s}", number, title) # need double lines, 4 in total to create a blank line
}
else{
print("Error!!!")
}
return(caption)
}
APA.fig.title <- function(name){ # enter object name
number <- names(which(APA.figure.names == name)) # title
title <- APA.figure.titles[which(names(APA.figure.titles) == number)] # subtitle
# caption <- sprintf("<b>%s</b><br><i>%s</i>", number, title)
return(c(number, title))
}
preserve <- append(preserve,c(APA.table.names, APA.table.titles, APA.figure.names, APA.figure.titles,
"APA.table.names","APA.table.titles", "APA.figure.names", "APA.figure.titles", "APA.tbl.caption", "APA.fig.title"))# Modify factor Levels
# Define factor levels of 3; this step is useful for future correlation analysis and data display
levels = as.character(c(1:3))
labels = c("Barely concerned", "Somewhat concerned", "Very Concerned")
# The variable "data" is used for 1.1
# Do not remove it until 1.2 Numeric Data
datasets.raw$pre_study <- datasets.raw$pre_study %>%
rename(improve.score = exam.score) %>%
left_join(select(datasets$other, emory.id, group)) %>%
mutate(
# Separate Race, religious affiliation, family affiliation, and types of meditation practiced using ';'
race = strsplit(race,";"),
religion = strsplit(religion,";"),
religion.family = strsplit(religion.family,";"),
meditation.type = strsplit(meditation.type,";"),
major = strsplit(major,";"),
worry = factor(worry, levels = as.character(c(1:5))),
concern.1 = factor(concern.1, levels=levels),
concern.2 = factor(concern.2, levels=levels),
concern.3 = factor(concern.3, levels=levels),
concern.4 = factor(concern.4, levels=levels),
concern.5 = factor(concern.5, levels=levels),
concern.6 = factor(concern.6, levels=levels),
religious.level = factor(religious.level , levels=as.character(c(1:5))),
meditation.freq = factor(meditation.freq, levels=as.character(c(1:3))),
meditation.times = factor(meditation.times, levels=as.character(c(1:4))),
relax = factor(relax, levels=as.character(c(1:5))),
focus = factor(focus, levels=as.character(c(1:5))),
improve.score = factor(improve.score, levels=as.character(c(1:5))),
)
concern.level = c("Barely concerned", "Somewhat concerned", "Very Concerned")
levels(datasets.raw$pre_study$concern.1) <- concern.level
levels(datasets.raw$pre_study$concern.2) <- concern.level
levels(datasets.raw$pre_study$concern.3) <- concern.level
levels(datasets.raw$pre_study$concern.4) <- concern.level
levels(datasets.raw$pre_study$concern.5) <- concern.level
levels(datasets.raw$pre_study$concern.6) <- concern.level
freq <- c("Never", "Weekly", "Monthly") # delete other because there is 0 response
times <- c("Never tried it", "Tried it 1-10 times", "Tried it 10-100 times", "Tried it 100-1000 times")
levels(datasets.raw$pre_study$meditation.freq) <- freq
levels(datasets.raw$pre_study$meditation.times) <- times
likert <- c("Stongly disagree", "Disagree", "Neutral", "Agree", "Strongly agree")
levels(datasets.raw$pre_study$relax) <- likert
levels(datasets.raw$pre_study$focus) <- likert
levels(datasets.raw$pre_study$improve.score) <- likert
demographic <- datasets.raw$pre_study %>%
select(1:35,72)
# remove unnecessary data
preserve <- append(preserve, c("group.legend", "demographic", "concern.level","freq", "times", "likert", "ordinal", "nominal", "ordinal.label", "nominal.label"))
rm(list=ls()[! ls() %in% preserve])1.0 Group differences in all demographic variables
1.1 Chi-Squared Test on Categorical Demographic Variables
Ordinal variables:
- Meditation frequency, meditation times, Difficulty, worry, and concerns (6 questions), religious level,
- Attitude (Likert Scale): meditation helps me relax, focus, improve exam score
Nested variables
- race, major, religion, religion.family, meditation.type (with multiple values stored as vectors inside a cell)
You can interact with the following tables.
meditation.habit <- c("meditation.freq","meditation.times")
meditation.habit.label = c("Meditation frequency","Meditation times")
demographic.cat.1 <- c("gender", "grade","ethnicity",meditation.habit)
demographic.cat.1.label = c("Gender", "Grade","Ethnicity",meditation.habit.label)
concerns = paste0("concern.",c(1:6))
concerns.label = paste("Concern:",
c("Disappointing your parents",
"Embarrassing yourself amongst friends",
"Embarrassing yourself amongst professors",
"Not getting the job or into the grad school that you hope to",
"Feeling unintelligent",
"Feeling lazy"))
attitude = c("relax","focus","improve.score")
attitude.label = c(paste("Meditation helps me", c("relax", "focus","improve my exam score")))
demographic.cat.1 = c(demographic.cat.1, "difficulty", "worry", concerns, "religious.level",
attitude)
demographic.cat.1.label = c(demographic.cat.1.label, "Difficulty of the selected course", "How worried are you about this exam?",
concerns.label,
"Religious level", attitude.label
)
demographic.cat.2 = c("race", "major", "religion", "religion.family", "meditation.type")
demographic.cat.2.label = c("Race", "Major", "Religion affiliation", "Religion affiliation of family members", "Meditation type")
baseline.vars <- list(demographic.cat.1, demographic.cat.2)
baseline.vars.label <- list(demographic.cat.1.label, demographic.cat.2.label)
ordinal <- c("difficulty", "worry", concerns, attitude,meditation.habit, "religious.level", "productivity6", "productivity7", "mood.min.1_6", "mood.min.7_21")
nominal <- c("race", "ethnicity", "gender", "grade","ethnicity","major","religion", "religion.family", "meditation.type","course")
ordinal.label <- c("Difficulty", "Worry", concerns.label, attitude.label,meditation.habit.label, "Religious level", "Productivity (Day 6)", "Productivity (Day 7)", "Min mood score (Day 1-6)", "Min mood score (Day 7-21)")
nominal.label <- c("Race", "Ethnicity", "Gender", "Grade","Ethnicity","Major","Religious affiliation", "Religious affiliation of family members", "Meditation type","Course")
preserve <- append(preserve, c("concerns","concerns.label","demographic.cat.1", "demographic.cat.2", "demographic.cat.1.label", "demographic.cat.2.label", "baseline.vars", "baseline.vars.label"))# function to calculate CI
CI.categorical = function(data, x="group", y){ # colname is a string # y is dependent variable, x is the grouping variable
new <- data %>%
group_by_(x,y) %>% # use the 'standard evaluation' counterpart of `group_by`
# ---*** Warning! some subjects left the question `meditation.freq` in blank ***---
filter(!is.na(.data[[y]])) %>% # To refer to column names that are stored as strings, use the `.data` pronoun
summarise(count = n()) %>%
ungroup() %>%
mutate(sum = sum(count),
proportion = count / sum) %>%
group_by_(x,y) %>% # use the 'standard evaluation' counterpart of `group_by`
mutate(lower = prop.test(count, sum)$conf.int[1], # CI.formula(sum, proportion)[1],
upper = prop.test(count, sum)$conf.int[2], # CI.formula(sum,proportion)[2],
) %>%
mutate(independent = x,dependent = y) %>%
rename_(.dots=setNames(x,"independent.value")) %>%
rename_(.dots=setNames(y,"dependent.value")) %>%
# rename_(.dots=setNames("var1","independent")) %>%
# rename_(.dots=setNames("var2","dependent")) %>%
select(independent,independent.value, dependent, dependent.value, proportion, lower, upper) %>% # remove original count of frequences of categorical value and the sum of all instences
arrange_('dependent.value') %>%
ungroup()
return(new)
}
# function for Chi-squared Post-hoc Testing and contingency table
post.chi <- function(tbl,x = "group",y, method = "bonferroni"){ # Post-hoc test of Chi-squared test
post.df <- chisq.posthoc.test(tbl, method = method)
post.df <- post.df %>%
rename(setNames("Dimension",x), value.type = Value) %>%
pivot_longer(cols = names(post.df)[-c(1, 2)], names_to = y, values_to = "value") %>%
pivot_wider(names_from = "value.type", values_from = "value") %>%
rename(p.value = `p values`, stdres = Residuals)
contingent.df <- as.data.frame(tbl) %>%
# rename(n = Freq) %>%
rename(setNames("Var1",x), setNames("Var2",y), n = Freq) %>%
mutate(proportion = n/sum(n))
return(post.df %>%
full_join(contingent.df, by = c(x,y)) %>%
mutate(relationship = paste(.data[[x]],.data[[y]], sep = " and "),
stdres = round(stdres, 2),
p.value = round(p.value, 3),
proportion = paste0(round(proportion,3) * 100,"%")) %>%
# rename(setNames("relationship",paste(x,y, sep = " and "))) %>%
# select(-.data[[x]],-.data[[y]]) %>%
select(7, 1:6)
)
}
# post.chi(tbl, y = "ethnicity")
# function for Chi-squared Hypothesis Testing
chi <- function(data, x="group", y, bonf = TRUE){# pass the variable name instead of the variable x, y
# data <- data %>%
# na.omit(cols = y)
tbl = table(data[[x]], data[[y]])
result <- chisq.test(data[[x]], data[[y]])
n = sum(tbl)
q = min(ncol(tbl),nrow(tbl))
effect = num.format(sqrt((result$statistic)/(n * (q-1)))) # parameter is the degree of freedom, calculates Cramer's V
chi.df <- data.frame("Chi-squared", x,y, num.format(result$statistic), result$parameter, n, p.format(result$p.value), effect)
colnames(chi.df) <- c("test.type","independent", "dependent","statistic", "df", "n", "p.value", "effect.size")
rownames(chi.df) = NULL
# Post-hoc test
if(!is.na(result$p.value) & result$p.value <= 0.05){
# define levels of factor
if(y %in% concerns){levels = concern.level
}else if(y %in% attitude){levels = likert
}else if(y == "meditation.freq"){levels = freq
}else if(y == "meditation.times"){levels = times
}else if(y %in% ordinal){levels = 1:5
}else{levels= NA_character_
}
report <- sprintf("A chi-square test of independence with continuity correction revealed a significant association between %s and %s, $\\chi^2(%d, N = %d) = %s$, $%s$, $\\phi_c = %s$.\n", x, y, result$parameter, n, num.format(result$statistic), p.format(result$p.value, table = F), effect)
report <- latex.escape(report)
post.df <- post.chi(tbl,x,y) %>%
mutate(stdres = num.format(stdres),
report = case_when(
p.value <= 0.05 ~ sprintf('Post-hoc chi-square test of independence with Bonferroni adjustment revealed a significant difference in the proportions of subjects who answered "%s" between "%s" and other %s counterparts, $Std\\ Residual(%d, N = %d, proportion = %s) = %s$, $%s$.\n',.[[3]],.[[2]],x, 1, n, proportion, stdres, p.format(p.value, table = F)),
# , $\\phi_c = %.2f$.\n don't know how to calculate the effect size
TRUE ~ NA_character_
),
p.value= p.format(p.value))
if(!TRUE %in% is.na(levels)){
post.df <- arrange(post.df, factor(.data[[y]], levels = levels))
}
# concatenate report from post-hoc test
for(i in post.df$report){
if(!is.na(i)) report <- paste0(report, latex.escape(i))# cat(i)
}
post.df <- post.df %>%
select(-report, -.data[[x]],-.data[[y]])
# colnames(chisq.demographic.post) <- c("relationship","stdres", "p.value", "n", "proportion")
rownames(chisq.demographic.post) = NULL
return(list(chi.df,post.df, report))
}
else{
return(list(chi.df))
}
}# data for first test
CI.demographic.cat = data.frame(matrix(ncol = 7, nrow = 0))
chisq.demographic = data.frame(matrix(ncol = 8, nrow = 0))
chisq.demographic.post <- data.frame(matrix(ncol = 6, nrow = 0))
report.demographic.cat <- ""
# analyze
i = 1
for(var in demographic.cat.1){
var.label = demographic.cat.1.label[i]
# CI.categorical
result = CI.categorical(demographic, y=var)
CI.demographic.cat = rbind(CI.demographic.cat, result)
# Hypothesis Testing
result = chi(demographic, y=var)
chisq.demographic = rbind(chisq.demographic, result[[1]])
if(length(result) > 1){
chisq.demographic.post <- rbind(chisq.demographic.post, result[[2]] %>%
mutate(relation.group = paste("Group",var.label,sep = " and ")) %>%
# rename(relationship = .data[[1]]) %>%
select(relation.group, 1:5)
)
report.demographic.cat <- paste(report.demographic.cat, result[[3]])
}
i = i + 1
gc(verbose = FALSE) # free some RAM
}
# Analyze nested variables
i = 1
for(var in demographic.cat.2){
var.label = demographic.cat.2.label[i]
# CI
tmp = demographic %>% select(group, var) %>%
unnest(var) # unnest the variable
result = CI.categorical(tmp, y=var)
CI.demographic.cat = rbind(CI.demographic.cat, result)
# Hypothesis testing
result = chi(tmp, y=var)
chisq.demographic = rbind(chisq.demographic, result[[1]])
if(length(result) > 1){
chisq.demographic.post = rbind(chisq.demographic.post, result[[2]] %>%
mutate(relation.group = paste("Group",var.label,sep = " and ")) %>%
select(relation.group, 1:5)
)
report.demographic.cat <- paste(report.demographic.cat, result[[3]])
}
i = i + 1
gc(verbose = FALSE)
}
## calculate group indentation in kable() display
APA.group.index <- function(data){
index.df <<- data %>%
rename(setNames(names(data)[[1]],"group"),setNames(names(data)[[2]],"relationship")) %>%
arrange(group, relationship) %>%
group_by(group) %>%
summarise(n = n()) %>%
mutate(group = case_when(is.na(group) ~ " ",
TRUE ~ group))
index <- index.df$n
names(index) <- index.df$group
return(index)
}
# Format dataframes into APA format
concerns.label <- c("Disappointing your parents",
"Embarrassing yourself amongst friends",
"Embarrassing yourself amongst professors",
"Not getting the job or into the grad school that you hope to",
"Feeling unintelligent",
"Feeling lazy")
## Chi-squared Table, APA format, theme color not changed (still gray)
chisq.tbl <- chisq.demographic %>% select(-independent, -test.type) %>%
mutate(dependent.group = case_when(startsWith(dependent, "concern") ~ "Concern",
startsWith(dependent, "meditation") ~ "Meditation habit",
startsWith(dependent, "religi") ~ "Religious affiliation",
dependent == "focus" ~ "Attitudes toward meditation",
dependent == "relax" ~ "Attitudes toward meditation",
dependent == "improve.score" ~ "Attitudes toward meditation",
),
dependent = case_when(
startsWith(dependent, "concern") & endsWith(dependent, "1") ~ concerns.label[1],
startsWith(dependent, "concern") & endsWith(dependent, "2") ~ concerns.label[2],
startsWith(dependent, "concern") & endsWith(dependent, "3") ~ concerns.label[3],
startsWith(dependent, "concern") & endsWith(dependent, "4") ~ concerns.label[4],
startsWith(dependent, "concern") & endsWith(dependent, "5") ~ concerns.label[5],
startsWith(dependent, "concern") & endsWith(dependent, "6") ~ concerns.label[6],
dependent == "religion.family" ~ "religious affliation of family members",
dependent == "religion" ~ "religious affliation",
dependent == "religious.level" ~ "level of religiosity",
TRUE ~ dependent
)
# statistic = as.numeric(statistic),
# effect.size = as.numeric(effect.size),
# p.value = as.numeric(p.value)
) %>%
select(dependent.group, dependent, df, statistic, effect.size, p.value, n) %>%
arrange(dependent.group,dependent)
index = APA.group.index(chisq.tbl)
APA.chisq.demographic <- chisq.tbl %>%
select(-dependent.group) %>%
# format table
kable(
# digits = c(NA, 0, 2, 2, 3, 0),
format = output.format, # stored before, at the beginning
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "$df$", "$\\chi^2$", "$\\phi_c$", "$p$", "$N$")),
align = c("l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.chisq.demographic")) %>%
# use html code to format the title
pack_rows(index = index) %>%
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
column_spec(column = 1, width = "2in")
## Frequency Table with proportion and post-hoc chisq table
index = APA.group.index(chisq.demographic.post)
APA.chisq.demographic.post <- chisq.demographic.post %>% select(-1) %>%
kable(
format = output.format,
longtable = TRUE,
col.names = latex.escape(c("", "Standard Residual", "$p$", "$N$", "proportion")),
caption = APA.tbl.caption("APA.chisq.demographic.post")) %>%
pack_rows(index = index) %>%
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
# kable_styling(full_width = TRUE) %>%
footnote(
general_title = "Note. ",
general = report.demographic.cat,
threeparttable = TRUE,
footnote_as_chunk = TRUE,
escape = escape.format
) %>%
column_spec(column = 1, width = "2in")
# Display
datatable(
CI.demographic.cat, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)APA.demographic.cat<- structure(paste0(APA.chisq.demographic,APA.chisq.demographic.post), format = format, class = ifelse(output.format == "latex", "knitr_kable", "kableExtra")) # "kableExtra"
APA.demographic.cat| \(df\) | \(\chi^2\) | \(\phi_c\) | \(p\) | \(N\) | |
|---|---|---|---|---|---|
| Attitudes toward meditation | |||||
| focus | 6 | 6.52 | .23 | .368 | 61 |
| improve.score | 8 | 4.63 | .19 | .797 | 61 |
| relax | 6 | 5.35 | .21 | .5 | 61 |
| Concern | |||||
| Disappointing your parents | 4 | .75 | .08 | .945 | 60 |
| Embarrassing yourself amongst friends | 4 | 3.09 | .16 | .544 | 60 |
| Embarrassing yourself amongst professors | 4 | .37 | .06 | .985 | 60 |
| Feeling lazy | 4 | 1.02 | .09 | .907 | 60 |
| Feeling unintelligent | 4 | 8.17 | .26 | .086 | 60 |
| Not getting the job or into the grad school that you hope to | 4 | 4.08 | .19 | .395 | 59 |
| Meditation habit | |||||
| meditation.freq | 4 | 11.33 | .31 | .023 | 58 |
| meditation.times | 6 | 9.93 | .29 | .128 | 61 |
| meditation.type | 8 | 3.6 | .14 | .892 | 89 |
| Religious affiliation | |||||
| level of religiosity | 8 | 5.9 | .22 | .659 | 61 |
| religious affliation | 14 | 8.93 | .27 | .836 | 61 |
| religious affliation of family members | 14 | 12.77 | .31 | .545 | 66 |
| difficulty | 4 | 3.38 | .17 | .496 | 60 |
| ethnicity | 6 | 6.98 | .24 | .323 | 61 |
| gender | 2 | .65 | .1 | .722 | 61 |
| grade | 2 | .09 | .04 | .957 | 61 |
| major | 44 | 46.39 | .56 | .374 | 73 |
| race | 6 | 8.94 | .26 | .177 | 65 |
| worry | 4 | 4.52 | .19 | .34 | 60 |
| Standard Residual | \(p\) | \(N\) | proportion | |
|---|---|---|---|---|
| Group and Meditation frequency | ||||
| Control and Never | .51 | 1 | 13 | 22.4% |
| 20 minuted 1X/day and Never | -.23 | 1 | 13 | 22.4% |
| 10 minutes 2X/day and Never | -.29 | 1 | 11 | 19% |
| Control and Weekly | -2.13 | .301 | 0 | 0% |
| 20 minuted 1X/day and Weekly | -.71 | 1 | 2 | 3.4% |
| 10 minutes 2X/day and Weekly | 2.89 | .034 | 6 | 10.3% |
| Control and Monthly | 1.17 | 1 | 6 | 10.3% |
| 20 minuted 1X/day and Monthly | .85 | 1 | 6 | 10.3% |
| 10 minutes 2X/day and Monthly | -2.07 | .35 | 1 | 1.7% |
| Note. A chi-square test of independence with continuity correction revealed a significant association between group and meditation.freq, \(\chi^2(4, N = 58) = 11.33\), \(p = .023\), \(\phi_c = .31\). Post-hoc chi-square test of independence with Bonferroni adjustment revealed a significant difference in the proportions of subjects who answered “Weekly” between “10 minutes 2X/day” and other group counterparts, \(Std\ Residual(1, N = 6, proportion = 10.3%) = 2.89\), \(p = .034\). | ||||
# remove unnecessary data
preserve <- append(preserve,c("APA.group.index", "CI.categorical","p.format","num.format","post.chi", "chi", "chisq.demographic","chisq.demographic.post","report.demographic.cat", "prop.tables.demographic","CI.demographic.cat", "baseline.vars"))
# rm(list=ls()[! ls() %in% preserve])- Most feared class and Reason to meditate are not included currently.
1.2 Analysis of the differences of numeric demographic variables across groups
- Discrete variables: credits, age
- Continuous variables: expected grade for the course, and estimated GPA
# Define functions
CI.numeric <- function(data, x = "group", y){
new = data %>%
group_by_(x) %>%
filter(!is.na(.data[[y]])) %>%
summarise(
mean = t.test(.data[[y]], conf.level = 0.95)$estimate,
lower = t.test(.data[[y]], conf.level = 0.95)$conf.int[1],
lower = case_when(lower < 0 ~ 0, TRUE ~ lower),
upper = t.test(.data[[y]], conf.level = 0.95)$conf.int[2]) %>%
ungroup() %>%
mutate(variable = y) %>%
select_(x, "variable", "mean", "lower", "upper")
return(new)
}
pairwise.wilcox.test2 <- function(x, g, p.adjust.method = p.adjust.methods, paired = FALSE,
...) # Modified according to the stats package as the original function does not provide test statistic
{
p.adjust.method <- match.arg(p.adjust.method)
DNAME <- paste(deparse1(substitute(x)), "and", deparse1(substitute(g)))
g <- factor(g) #
METHOD <- NULL
# unmodified
compare.levels <- function(i, j) {
xi <- x[as.integer(g) == i]
xj <- x[as.integer(g) == j]
if (is.null(METHOD)) {
wt <- wilcox.test(xi, xj, paired = paired, ...)
METHOD <<- wt$method
wt$p.value
}
else wilcox.test(xi, xj, paired = paired, ...)$p.value
}
PVAL <- pairwise.table(compare.levels, levels(g), p.adjust.method)
# modified
PVAL <- as.data.frame(as.table(PVAL)) %>%
filter(!is.na(Freq)) %>%
rename(p.value = Freq)
# The following 2 functions are used to display test statistic
# modified according to original function
compare.levels.stat <- function(i, j) {
xi <- x[as.integer(g) == i]
xj <- x[as.integer(g) == j]
wilcox.test(xi, xj, paired = paired, ...)$statistic
}
pairwise.table.stat <- function(compare.levels, level.names)
{
ix <- setNames(seq_along(level.names), level.names)
stat <- outer(ix[-1L], ix[-length(ix)], function(ivec, jvec) vapply(seq_along(ivec),
function(k) {
i <- ivec[k]
j <- jvec[k]
if (i > j)
compare.levels.stat(i, j)
else NA_real_
}, 0.05))
# il.tri <- lower.tri(pp, TRUE)
# pp[il.tri] <- p.adjust(pp[il.tri], p.adjust.method)
stat
}
STAT <- pairwise.table.stat(compare.levels, levels(g))
STAT <- as.data.frame(as.table(STAT)) %>%
filter(!is.na(Freq)) %>%
rename(statistic = Freq)
ans <- list(method = METHOD, data.name = DNAME, statistic = STAT, p.value = PVAL,
p.adjust.method = p.adjust.method)
class(ans) <- "pairwise.htest2" # self-modified, avoid repetition of class name
ans
}
kruskal_wilcox <- function(data, x = "group", y, sig = 0.05){ # y is the dependent variable (numeric), x is the grouping variable (categorical)
# Kruskal-Wallis rank sum test
data = data %>% select(x,y) %>% na.omit()
kruskal = kruskal_test(reformulate(x,y), data = data) %>% # y~x
rename(test.type = method,
dependent = .y.,
p.value = p) %>%
mutate(independent = x,
effect.size = kruskal_effsize(data, reformulate(x,y))$effsize) %>% # incorporate effect size
select(test.type, independent, dependent, statistic, df, n, p.value, effect.size)
# Multiple pairwise Wilcox unpaired test between groups
if(kruskal$p.value <= sig){
report <- kruskal.report(kruskal)
wilcox = pairwise.wilcox.test2(data[[y]],data[[x]], p.adjust.method = "BH") # y,x
effect = wilcox_effsize(data, reformulate(x,y))
med <- data %>%
group_by(.data[[x]]) %>%
summarise(mdn = median(.data[[y]], na.rm = T))
wilcox = effect %>%
full_join(wilcox$statistic, by = c("group2"= "Var1","group1" = "Var2")) %>%
full_join(wilcox$p.value, by = c("group2"= "Var1","group1" = "Var2")) %>%
left_join(med, by = c("group1" = x)) %>% rename(mdn1 = mdn) %>%
left_join(med, by = c("group2" = x)) %>% rename(mdn2 = mdn) %>%
rename(dependent = .y.,
independent1 = group1,
independent2 = group2,
effect.size = effsize) %>%
mutate(independent = x,
test.type = "pairwise Wilcox") %>%
select(test.type, independent, independent1, independent2, dependent, statistic, n1, n2, mdn1, mdn2, p.value, effect.size)
report <- paste(report, wilcox.report(wilcox))
kruskal <- kruskal %>%
mutate(statistic = num.format(statistic),
effect.size = num.format(effect.size),
p.value = p.format(p.value),
)
wilcox <- wilcox %>%
mutate(statistic = num.format(statistic),
effect.size = num.format(effect.size),
p.value = p.format(p.value),
mdn1 = num.format(mdn1),
mdn2 = num.format(mdn2),
)
return(list(kruskal, wilcox, report))
}
else{
kruskal <- kruskal %>%
mutate(statistic = num.format(statistic),
effect.size = num.format(effect.size),
p.value = p.format(p.value),
)
return(list(kruskal))
}
}
kruskal.report <- function(data){
result <- NULL
# x,y,df, n, statistic, p, effsize
for (i in 1:nrow(data)) {
row = data[i,]
# print(row)
str <- sprintf("A Kruskal-Wallis test showed that the mean rank of %s was significantly different between the levels of %s, $H(%d,N=%d)=%s$, $%s$, $\\eta^2=%s$.", row$dependent,row$independent, row$df, row$n, num.format(row$statistic), p.format(row$p.value, table = F), num.format(row$effect.size))
str <- latex.escape(str)
# print(str)
result <- append(result, str)
}
return(result)
# $\phi_c$, $\eta^2$, and $r$
}
wilcox.report <- function(data){
result <- c()
# x,y,df, n, statistic, p, effsize
for (i in 1:nrow(data)) {
row = data[i,]
# print(row)
if(row$p.value < 0.05){
str.size <- ifelse(row$n1 == row$n2, sprintf("n_1 = n_2 = %d", row$n1), sprintf("n_1 = %d, n_2 = %d", row$n1, row$n2))
str <- sprintf("Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of %s between %s: %s ($Mdn_1=%.2f$) and %s: %s ($Mdn_2=%.2f$), $W(%s)=%s$, $%s$, $r = %s$.", row$dependent,row$independent,row$independent1,row$mdn1, row$independent,row$independent2, row$mdn2, str.size, num.format(row$statistic), p.format(row$p.value, table = F), num.format(row$effect.size))
str <- latex.escape(str)
result <- append(result, str)
}
}
return(result)
}
# Define variable names and labels
demographic.num = c("credits","age",
"GPA.previous", "grade.predict")
demographic.num.labels = c("Credits of the current semester","Age", "Expected grade for the most challenging course", "Estimated GPA of the previous semester")
# Define datasets
CI.demographic.num = data.frame(matrix(ncol = 5, nrow = 0))
kruskal.demographic = data.frame(matrix(ncol = 8, nrow = 0))
wilcox.demographic = data.frame(matrix(ncol = 12, nrow = 0))
report.demographic.num = ""
for(var in demographic.num){
# CI
result = CI.numeric(demographic, y = var)
CI.demographic.num = rbind(CI.demographic.num, result)
# Hypothesis Testing
result = kruskal_wilcox(demographic %>% select(group, var), y = var)
kruskal.demographic = rbind(kruskal.demographic, result[[1]])
if(length(result) > 1){
wilcox.demographic = rbind(wilcox.demographic, result[[2]])
report.demographic.num = paste(report.demographic.num, result[[3]])
}
gc(verbose = FALSE) # free some RAM
}
baseline.vars <- append(baseline.vars, list(demographic.num))
# CI_func <- CI.numeric
# remove unnecessary data
preserve <- append(preserve,c("pairwise.wilcox.test2","kruskal_wilcox", "kruskal.report", "wilcox.report", "CI.numeric", "CI.demographic.num", "kruskal.demographic", "wilcox.demographic"))
rm(list=ls()[! ls() %in% preserve])
# Display
# datatable(
# kruskal.demographic, filter = 'top', extensions = 'Buttons', options = list(
# dom = 'Bfrtip',
# buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
# )
# )
APA.kruskal.demographic <- kruskal.demographic %>% select(-independent, -test.type) %>%
mutate(dependent = case_when(
dependent == "grade.predict" ~ "Predicted grade of the chosen challenging course",
dependent == "GPA.previous" ~ "GPA of previous semester",
TRUE ~ dependent
)
) %>%
select(dependent, df, statistic, effect.size, p.value, n) %>%
arrange(dependent) %>%
# select(-dependent.group) %>%
# format table
kable(
# digits = c(NA, 0, 2, 2, 3, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "df", "$H$", "$\\eta^2$", "p", "N")),
align = c("l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.kruskal.demographic")) %>%
# use html code to format the title
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
footnote(
general_title = "Note.",
general = latex.escape('No significant differences were found in Kruskal-Wallis Test.', report = TRUE),
threeparttable = TRUE,
footnote_as_chunk = TRUE,
escape = escape.format
) %>%
column_spec(column = 1, width = "2in")
APA.kruskal.demographic| df | \(H\) | \(\eta^2\) | p | N | |
|---|---|---|---|---|---|
| age | 2 | .65 | -.02 | .722 | 61 |
| credits | 2 | 4.18 | .04 | .124 | 60 |
| GPA of previous semester | 2 | 1.19 | -.02 | .552 | 57 |
| Predicted grade of the chosen challenging course | 2 | .31 | -.03 | .856 | 60 |
| Note. No significant differences were found in Kruskal-Wallis Test. |
2.0 Group differences in self-reported study productivity on day 6 and 7, exam grade, and hours studying
2.1 Group differences in self-reported study productivity on day 6 and 7
# modified mood.6 and mood.7
productivity67 <- datasets.raw$daily %>%
filter(day == 6) %>% select(emory.id, productivity, day) %>%
bind_rows(select(datasets.raw$daily, emory.id, productivity, day) %>% filter(day == 7) ) %>%
left_join(datasets$other %>% select(emory.id, group), by = "emory.id")
# Create Productivity Dataset separately for day 6 and 7
productivity6 <- productivity67 %>%
filter(day == 6) %>%
select(emory.id,group, productivity) %>%
rename(productivity6 = productivity)
productivity7 <- productivity67 %>%
filter(day == 7) %>%
select(emory.id,group, productivity) %>%
rename(productivity7 = productivity)
productivity6$duplicate <- duplicated(productivity6$emory.id)
productivity6$productivity6[which(productivity6$duplicate == TRUE)-1] <- 3.5 # vary between 3 and 4
productivity6 <- productivity6[-which(productivity6$duplicate == TRUE),]
productivity6$duplicate = NULL
datasets.raw$pre_study <- datasets.raw$pre_study %>%
left_join(select(productivity6,-group),by="emory.id") %>%
left_join(select(productivity7,-group),by="emory.id")
productivity67 <- list(productivity6,productivity7)
# Define Variables Names
productivity = c("productivity6", "productivity7")
productivity.label = c("productivity day 6", "productivity day 7")
baseline.vars <- append(baseline.vars, list(productivity))
# Define data sets
CI.productivity <- data.frame(matrix(ncol = 7, nrow = 0))
chisq.productivity = data.frame(matrix(ncol = 8, nrow = 0))
chisq.productivity.post <- data.frame(matrix(ncol = 6, nrow = 0))
report.productivity <- ""
i = 1
for(var in demographic.cat.1){
var.label = demographic.cat.1.label[i]
# CI.categorical
result = CI.categorical(demographic, y=var)
CI.demographic.cat = rbind(CI.demographic.cat, result)
# Hypothesis Testing
result = chi(demographic, y=var)
chisq.demographic = rbind(chisq.demographic, result[[1]])
if(length(result) > 1){
chisq.demographic.post <- rbind(chisq.demographic.post, result[[2]] %>%
mutate(relation.group = paste("Group",var.label,sep = " and ")) %>%
# rename(relationship = .data[[1]]) %>%
select(relation.group, 1:5)
)
report.demographic.cat <- paste(report.demographic.cat, result[[3]])
}
i = i + 1
gc(verbose = FALSE) # free some RAM
}
# analyze
i = 1
for(var in productivity){
# CI.categorical
index.df = which(productivity == var)
result = CI.categorical(productivity67[[index.df]], y=var)
CI.productivity = rbind(CI.productivity, result)
# Hypothesis Testing
result = chi(productivity67[[index.df]], x="group", y=var)
# chisq.productivity[nrow(chisq.productivity) + 1,] = result
chisq.productivity = rbind(chisq.productivity, result[[1]])
if(length(result) > 1){
chisq.productivity.post <- rbind(chisq.productivity.post, result[[2]] %>%
mutate(relation.group = paste("Group",var.label,sep = " and ")) %>%
# rename(relationship = .data[[1]]) %>%
select(relation.group, 1:5)
)
report.productivity <- paste(report.productivity, result[[3]])
}
i = i + 1
gc(verbose = FALSE) # free some RAM
}
# Display
datatable(
CI.productivity, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)datatable(
chisq.productivity, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)# Chi-squared Table, APA format, theme color not changed (still gray)
APA.productivity <- chisq.productivity %>% select(-independent, -test.type) %>%
mutate(
dependent = case_when(
dependent == "productivity6" ~ "Productivity (Day 6)",
dependent == "productivity7" ~ "Productivity (Day 7)",
TRUE ~ dependent
) # ,
# statistic = as.numeric(statistic),
# effect.size = as.numeric(effect.size),
# p.value = as.numeric(p.value)
) %>%
select(dependent, df, statistic, effect.size, p.value, n) %>%
arrange(dependent) %>%
# select(-dependent.group) %>%
# format table
kable(
digits = c(NA, 0, 2, 2, 3, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "df", "$\\chi^2$", "$\\phi_c$", "p", "N")),
align = c("l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.productivity")) %>%
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
# kable_styling(full_width = TRUE) %>%
footnote(
general_title = "Note.",
general = latex.escape('A chi-square test of independence with continuity correction revealed no significant association between meditation group and level of productivity on Day 6 and 7.', report = TRUE),
threeparttable = TRUE,
footnote_as_chunk = TRUE,
escape = escape.format
) %>%
column_spec(column = 1, width = "2in")
APA.productivity| df | \(\chi^2\) | \(\phi_c\) | p | N | |
|---|---|---|---|---|---|
| Productivity (Day 6) | 10 | 5.29 | .24 | .871 | 46 |
| Productivity (Day 7) | 8 | 8.69 | .3 | .369 | 48 |
| Note. A chi-square test of independence with continuity correction revealed no significant association between meditation group and level of productivity on Day 6 and 7. |
# remove unnecessary data
preserve <- append(preserve,c("chisq.productivity", "CI.productivity", "productivity6", "productivity7"))
rm(list=ls()[! ls() %in% preserve])2.2 Group differences in exam grade and hours studying
# Define variable names and labels
vars = c("exam.score", "study.time.total")
labels = c("Exam score", "Hours studying")
baseline.vars <- append(baseline.vars, list(vars))
# Define data sets
CI.grade_hours = data.frame(matrix(ncol = 5, nrow = 0))
kruskal.grade_hours = data.frame(matrix(ncol = 8, nrow = 0))
wilcox.grade_hours = data.frame(matrix(ncol = 12, nrow = 0))
report.grade_hours = ""
# analyze
for(var in vars){
# CI
result = CI.numeric(datasets$other, y = var)
CI.grade_hours = rbind(CI.grade_hours, result)
# hypothesis testing
result = kruskal_wilcox(datasets$other, y = var)
kruskal.grade_hours = rbind(kruskal.grade_hours, result[[1]])
if(length(result) > 1){
wilcox.grade_hours = rbind(wilcox.grade_hours, result[[2]])
report.grade_hours = paste(report.grade_hours, result[[3]])
}
gc(verbose = FALSE) # free some RAM
}
# Display
datatable(
CI.grade_hours, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)# datatable(
# kruskal.grade_hours, filter = 'top', extensions = 'Buttons', options = list(
# dom = 'Bfrtip',
# buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
# )
# )
# datatable(
# wilcox.grade_hours, filter = 'top', extensions = 'Buttons', options = list(
# dom = 'Bfrtip',
# buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
# )
# )
APA.kruskal.grade_hours <- kruskal.grade_hours %>% select(-independent, -test.type) %>%
mutate(dependent = case_when(
dependent == "exam.score" ~ "exam score",
# dependent == "GPA.previous" ~ "GPA of previous semester",
TRUE ~ "total hours of studying"
)
) %>%
select(dependent, df, statistic, effect.size, p.value, n) %>%
arrange(dependent) %>%
# select(-dependent.group) %>%
# format table
kable(
# digits = c(NA, 0, 2, 2, 3, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "df", "$H$", "$\\eta^2$", "p", "N")),
align = c("l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.kruskal.grade_hours")) %>%
# use html code to format the title
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
column_spec(column = 1, width = "2in")
APA.tbl.caption("APA.wilcox.grade")## [1] "<b>Table 6</b><br><i>Post-hoc Wilcoxon Test Results of Meditation Group vs Exam Score</i>"APA.wilcox.grade <- wilcox.grade_hours %>% mutate(pair = sprintf("%s - %s", independent1, independent2),
n = sprintf("%s - %s", n1, n2),
mdn = sprintf("%s - %s", mdn1, mdn2)
) %>%
select(pair, statistic, effect.size, p.value, n, mdn) %>%
arrange(pair) %>%
# select(-dependent.group) %>%
# format table
kable(
# digits = c(NA, 2, 2, 3, 0, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "$W$", "$r$", "p", "N", "Median")),
align = c("l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.wilcox.grade")) %>%
# use html code to format the title
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
footnote(
general_title = "Note.",
general = report.grade_hours,
threeparttable = TRUE,
footnote_as_chunk = TRUE,
escape = escape.format
) %>%
column_spec(column = 1, width = "2in")
APA.grade_hours <- structure(paste0(APA.kruskal.grade_hours,APA.wilcox.grade), format = output.format, class = ifelse(output.format == "latex", "knitr_kable", "kableExtra"))
APA.grade_hours| df | \(H\) | \(\eta^2\) | p | N | |
|---|---|---|---|---|---|
| exam score | 2 | 6.39 | .11 | .041 | 44 |
| total hours of studying | 2 | .72 | -.02 | .697 | 61 |
| \(W\) | \(r\) | p | N | Median | |
|---|---|---|---|---|---|
| 20 minuted 1X/day - 10 minutes 2X/day | 90.5 | .12 | .539 | 15 - 14 | 93 - 92 |
| Control - 10 minutes 2X/day | 141.0 | .29 | .181 | 15 - 14 | 81 - 92 |
| Control - 20 minuted 1X/day | 174.0 | .47 | .034 | 15 - 15 | 81 - 93 |
| Note. A Kruskal-Wallis test showed that the mean rank of exam.score was significantly different between the levels of group, \(H(2,N=44)=6.39\), \(p = .041\), \(\eta^2=.11\). Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of exam.score between group: Control (\(Mdn_1=81.00\)) and group: 20 minuted 1X/day (\(Mdn_2=93.00\)), \(W(n_1 = n_2 = 15)=174\), \(p = .034\), \(r = .47\). |
# Report
# cat(kruskal.report(kruskal.grade_hours))
# cat(wilcox.report(wilcox.grade_hours))
# remove unnecessary data
preserve <- append(preserve,c("CI.grade_hours", "kruskal.grade_hours", "wilcox.grade_hours"))
rm(list=ls()[! ls() %in% preserve])3.0 AEQs
State Measures are AEQs on day 6, 7, and 14. Trait Measure is the difference between AEQ score on Day 21 and before the study began (day 21 - pre study).
CI.AEQ <- function(data){
new = data %>%
group_by(subscale, group) %>%
summarise(
avg = t.test(score, conf.level = 0.95)$estimate,
L = t.test(score, conf.level = 0.95)$conf.int[1],
U = t.test(score, conf.level = 0.95)$conf.int[2]) %>%
ungroup()
return(new)
}
plot.aeq <- function(data, title = "AEQ", x = "Group", y = "Score"){ # plot a dotplot
subscales = unique(data$subscale)
p = data %>%
# filter(subscale == subscales[i]) %>%
ggplot(mapping = aes(x=group, y=score)) +
geom_dotplot(binaxis='y', stackdir='center', alpha = 0.6) +
facet_grid(rows = vars(subscale)) +
stat_summary(fun.data = "mean_cl_normal",
geom = "errorbar",size = 0.6,
width = .3, color = "black") +
stat_summary(fun = "mean", geom = "point", color = "red", size = 2) +
# geom_errorbar(aes(ymin = , ymax = , fill = group), width = 0.1, size = 0.8) +
labs(# title=title,
y =y, x = x # , caption = "Error bars indicate 95% CI"
)
return(p)
}
hypo.aeq <- function(data){
results = vector(mode = "list", length = 2)
kruskal = data.frame(matrix(ncol = 8, nrow = 0))
wilcox = data.frame(matrix(ncol = 12, nrow = 0))
report = ""
subscales = unique(data$subscale) # get the subscale names by removing duplicates
for(var in subscales){
new = data %>%
filter(subscale == var) %>%
mutate_(.dots=setNames(list(quote(score)),var)) %>% # standard evaluation https://stackoverflow.com/questions/37259581/r-understanding-standard-evaluation-in-mutate
select_('group', var)
result = kruskal_wilcox(new, y = var)
if(length(result) > 1){
kruskal = rbind(kruskal, result[[1]])
wilcox = rbind(wilcox, result[[2]])
report = paste(report, result[[3]])
}
else{
kruskal = rbind(kruskal, result[[1]])
}
results = list(kruskal, wilcox, report)
}
return(results)
}# CI in tables
measure <- c(paste("Total State AEQ, Day", c(6,7,14)), "Trait AEQ change (Day 21 - Day 0)")
AEQ_6.CI <- datasets$AEQ_6 %>% CI.AEQ() %>%
mutate(AEQ.type = "state", measure = measure[1])
AEQ_7.CI <- datasets$AEQ_7 %>% CI.AEQ() %>%
mutate(AEQ.type = "state", measure = measure[2])
AEQ_14.CI <- datasets$AEQ_14 %>% CI.AEQ() %>%
mutate(AEQ.type = "state", measure = measure[3])
difference.CI <- datasets$difference %>%
group_by(subscale, group) %>%
summarise(
avg = t.test(score, conf.level = 0.95)$estimate,
L = t.test(score, conf.level = 0.95)$conf.int[1],
U = t.test(score, conf.level = 0.95)$conf.int[2]) %>%
ungroup() %>%
mutate(AEQ.type = "trait", measure = measure[4])
CI.AEQ <- rbind(AEQ_6.CI, AEQ_7.CI, AEQ_14.CI, difference.CI)
# CI in graphs
# titles = paste("Total State AEQ Score on Day", c(6, 7, 14))
i = 1
AEQ.graphs = list()
for(data in datasets[c(1,2,3)]){
graph = data %>% plot.aeq(title = measure[i])
AEQ.graphs = append(AEQ.graphs, list(graph))
names(AEQ.graphs)[length(AEQ.graphs)] = c(6, 7, 14)[i]
i = i+1
gc(verbose = FALSE) # free some RAM
}
difference.graph <- datasets$difference %>%
plot.aeq(title = measure[4])
AEQ.graphs = append(AEQ.graphs, list(difference.graph))
names(AEQ.graphs)[length(AEQ.graphs)] = "difference"
AEQ.6714 <- ggarrange(AEQ.graphs[[1]] + labs(subtitle = "Day 6") +
theme(plot.subtitle = element_text(size = 12, hjust = 0.2,face= "italic"),
plot.title.position = "plot"),
AEQ.graphs[[2]] + labs(subtitle = "Day 7") +
theme(plot.subtitle = element_text(size = 12, hjust = 0.2,face= "italic"),
plot.title.position = "plot"),
AEQ.graphs[[3]] + labs(subtitle = "Day 14") +
theme(plot.subtitle = element_text(size = 12, hjust = 0.2,face= "italic"),
plot.title.position = "plot") ,
# labels = c(paste("Day", c(6, 7, 14))),
ncol = 1, nrow = 3, # hjust = c(-0.9, -0.9,-0.9),
vjust = 0,
heights = c(0.8,0.8,0.8))
APA.AEQ.fig <- ggarrange(AEQ.6714,AEQ.graphs[[4]] + labs(subtitle = "Day 21 - pre-study") +
theme(plot.subtitle = element_text(size = 12, hjust = 0.2,face= "italic"),
plot.title.position = "plot"
),
# labels = c(paste("Day", c(6, 7, 14)), "21-0"),
ncol = 2, hjust = 0,
vjust = 0) +
labs(title = APA.fig.title("APA.AEQ.fig")[1],
subtitle = APA.fig.title("APA.AEQ.fig")[2],
caption = expression(paste(italic("Note. "), "Error bars indicate 95% CI."))) +
theme(# text = element_text(size = 12),
plot.title = element_text(size = 12, face = "bold",
hjust = 0,
margin = margin(b = 15)),
plot.subtitle = element_text(size = 12, face = "italic",
hjust = 0,
margin = margin(b = 15)
),
plot.caption = element_text(size = 12, hjust = 0,face= "italic"),
# axis.text = element_text(size = 8, color = "black"),
plot.title.position = "plot",
plot.caption.position = "plot",
)
APA.AEQ.fig
# hypothesis testing
kruskal.AEQ = data.frame(matrix(ncol = 9, nrow = 0))
wilcox.AEQ = data.frame(matrix(ncol = 13, nrow = 0))
report.AEQ = ""
i = 1
for(data in datasets[1:4]){
result = hypo.aeq(data)
result[[1]] = result[[1]] %>%
mutate(measure = measure[i])
kruskal.AEQ = rbind(kruskal.AEQ, result[[1]])
if(length(result[[2]]) != 0){
result[[2]] = result[[2]] %>%
mutate(measure = measure[i])
wilcox.AEQ = rbind(wilcox.AEQ, result[[2]])
report.AEQ = paste(report.AEQ, result[[3]])
}
i = i+1
gc(verbose = FALSE) # free some RAM
}
# Display
datatable(
CI.AEQ, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)# datatable(
# kruskal.AEQ, filter = 'top', extensions = 'Buttons', options = list(
# dom = 'Bfrtip',
# buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
# )
# )
# datatable(
# wilcox.AEQ, filter = 'top', extensions = 'Buttons', options = list(
# dom = 'Bfrtip',
# buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
# )
# )
APA.kruskal.AEQ <- kruskal.AEQ %>% select(-test.type) %>%
select(measure, dependent, df, statistic, effect.size, p.value, n) %>%
arrange(dependent) %>%
# select(-dependent.group) %>%
# format table
kable(
# digits = c(NA, 0, 2, 2, 3, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "", "$df$", "$H$", "$ \\eta^2$", "$p$", "$N$")),
align = c("l","l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.kruskal.AEQ")) %>%
# use html code to format the title
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
column_spec(column = 1, width = "2in")
APA.wilcox.AEQ <- wilcox.AEQ %>% mutate(pair = sprintf("%s - %s", independent1, independent2),
n = sprintf("%s - %s", n1, n2),
mdn = sprintf("%s - %s", mdn1, mdn2)
) %>%
select(measure, dependent, pair, statistic, effect.size, p.value, n, mdn) %>%
arrange(dependent,pair) %>%
# select(-dependent.group) %>%
# format table
kable(
# digits = c(NA, 2, 2, 3, 0, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "", "", "$W$", "$r$", "$p$", "$N$", "Median")),
align = c("l", "l", "l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.wilcox.AEQ")) %>%
# use html code to format the title
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
footnote(
general_title = "Note.",
general = report.AEQ,
threeparttable = TRUE,
footnote_as_chunk = TRUE,
escape = escape.format
) %>%
column_spec(column = 1, width = "2in")
APA.AEQ <- structure(paste0(APA.kruskal.AEQ,APA.wilcox.AEQ), format = output.format, class = ifelse(output.format == "latex", "knitr_kable", "kableExtra"))
APA.AEQ| \(df\) | \(H\) | $ ^2$ | \(p\) | \(N\) | ||
|---|---|---|---|---|---|---|
| Total State AEQ, Day 7 | anger | 2 | .73 | -.03 | .696 | 47 |
| Total State AEQ, Day 14 | anger | 2 | 1.68 | -.01 | .431 | 27 |
| Trait AEQ change (Day 21 - Day 0) | anger | 2 | .88 | -.04 | .644 | 31 |
| Total State AEQ, Day 6 | anxiety | 2 | 6.4 | .11 | .041 | 42 |
| Trait AEQ change (Day 21 - Day 0) | anxiety | 2 | 2.94 | .03 | .23 | 31 |
| Total State AEQ, Day 6 | hope | 2 | 2.45 | .01 | .294 | 42 |
| Trait AEQ change (Day 21 - Day 0) | hope | 2 | 1.5 | -.02 | .473 | 30 |
| Total State AEQ, Day 6 | hopeless | 2 | 3.41 | .04 | .182 | 42 |
| Trait AEQ change (Day 21 - Day 0) | hopeless | 2 | 1.62 | -.01 | .444 | 30 |
| Total State AEQ, Day 6 | joy | 2 | .62 | -.04 | .734 | 42 |
| Trait AEQ change (Day 21 - Day 0) | joy | 2 | .94 | -.04 | .624 | 31 |
| Total State AEQ, Day 7 | pride | 2 | 7.04 | .11 | .03 | 48 |
| Total State AEQ, Day 14 | pride | 2 | 2 | 0 | .368 | 27 |
| Trait AEQ change (Day 21 - Day 0) | pride | 2 | .83 | -.04 | .66 | 31 |
| Total State AEQ, Day 7 | relief | 2 | 1.06 | -.02 | .589 | 48 |
| Total State AEQ, Day 14 | relief | 2 | 2.56 | .02 | .278 | 27 |
| Trait AEQ change (Day 21 - Day 0) | relief | 2 | .26 | -.06 | .878 | 31 |
| Total State AEQ, Day 7 | shame | 2 | 3.29 | .03 | .193 | 48 |
| Total State AEQ, Day 14 | shame | 2 | 1.44 | -.02 | .486 | 27 |
| Trait AEQ change (Day 21 - Day 0) | shame | 2 | 1.36 | -.02 | .507 | 31 |
| \(W\) | \(r\) | \(p\) | \(N\) | Median | |||
|---|---|---|---|---|---|---|---|
| Total State AEQ, Day 6 | anxiety | 20 minuted 1X/day - 10 minutes 2X/day | 145.0 | .43 | .071 | 16 - 12 | 13 - 15.5 |
| Total State AEQ, Day 6 | anxiety | Control - 10 minutes 2X/day | 86.0 | .02 | .938 | 14 - 12 | 16 - 15.5 |
| Total State AEQ, Day 6 | anxiety | Control - 20 minuted 1X/day | 64.0 | .37 | .071 | 14 - 16 | 16 - 13 |
| Total State AEQ, Day 7 | pride | 20 minuted 1X/day - 10 minutes 2X/day | 86.0 | .28 | .179 | 17 - 15 | 15 - 13 |
| Total State AEQ, Day 7 | pride | Control - 10 minutes 2X/day | 142.5 | .16 | .382 | 16 - 15 | 12.5 - 13 |
| Total State AEQ, Day 7 | pride | Control - 20 minuted 1X/day | 209.0 | .46 | .026 | 16 - 17 | 12.5 - 15 |
| Note. A Kruskal-Wallis test showed that the mean rank of anxiety was significantly different between the levels of group, \(H(2,N=42)=6.4\), \(p = .041\), \(\eta^2=.11\). A Kruskal-Wallis test showed that the mean rank of pride was significantly different between the levels of group, \(H(2,N=48)=7.04\), \(p = .03\), \(\eta^2=.11\). Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of pride between group: Control (\(Mdn_1=12.50\)) and group: 20 minuted 1X/day (\(Mdn_2=15.00\)), \(W(n_1 = 16, n_2 = 17)=209\), \(p = .026\), \(r = .46\). |
# remove unnecessary data
preserve <- append(preserve,c("CI.AEQ", "kruskal.AEQ", "wilcox.AEQ", "figure.AEQ", "report.AEQ"))
rm(list=ls()[! ls() %in% preserve])
# Report
cat(report.AEQ)## A Kruskal-Wallis test showed that the mean rank of anxiety was significantly different between the levels of group, $H(2,N=42)=6.4$, $p = .041$, $\eta^2=.11$. A Kruskal-Wallis test showed that the mean rank of pride was significantly different between the levels of group, $H(2,N=48)=7.04$, $p = .03$, $\eta^2=.11$. Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of pride between group: Control ($Mdn_1=12.50$) and group: 20 minuted 1X/day ($Mdn_2=15.00$), $W(n_1 = 16, n_2 = 17)=209$, $p = .026$, $r = .46$.A Kruskal-Wallis test showed that the mean rank of anxiety was significantly different between the levels of group, \(H(2, N = 42) = 6.40\), \(p = .041\), $^2 = 0.11 $.
A Kruskal-Wallis test showed that the mean rank of pride was significantly different between the levels of group, \(H(2, N = 48) = 7.04\), \(p = .03\), $^2 = 0.11 $.
- Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of pride between group: Control (\(Mdn1 = 12.50\)) and group: 20 minuted 1X/day (\(Mdn2 = 15.00\)), \(W(n1 = 16, n2 = 17) = 209.00\), \(p = .026\), \(r = 0.46\).
4.0 Daily Mood Score
95% CI for daily mood
group by duration (1~6 or 7~21) and group: 2 * 3 = 6 CI in total
CI.mood <- datasets$mood %>%
group_by(mood.type, group) %>%
summarise(
avg.mean = t.test(mood.avg, conf.level = 0.95)$estimate,
min.mean = t.test(mood.min, conf.level = 0.95)$estimate,
avg.L = t.test(mood.avg, conf.level = 0.95)$conf.int[1],
avg.U = t.test(mood.avg, conf.level = 0.95)$conf.int[2],
min.L = t.test(mood.min, conf.level = 0.95)$conf.int[1],
min.U = t.test(mood.min, conf.level = 0.95)$conf.int[2]) %>%
ungroup()
avg.graph <- datasets$mood %>%
ggplot(mapping = aes(x=mood.type, y=mood.avg)) +
geom_dotplot(binaxis='y', stackdir='center') +
stat_summary(fun.data = "mean_cl_normal",
geom = "errorbar",size = 1,
width = .3, color = "red") +
stat_summary(fun = "mean", geom = "point", color = "red") +
# geom_errorbar(aes(ymin = , ymax = , fill = group), width = 0.1, size = 0.8) +
facet_wrap(~group) +
labs(# title="4.1 Average Daily Mood across Groups",
x ="Day", y = "Mood (Likert Scale)" #, caption = "Error bars indicate 95% CI"
)
min.graph <- datasets$mood %>%
ggplot(mapping = aes(x=mood.type, y=mood.min)) +
geom_dotplot(binaxis='y', stackdir='center') +
stat_summary(fun.data = "mean_cl_normal",
geom = "errorbar",size = 1,
width = .3, color = "red") +
stat_summary(fun = "mean", geom = "point", color = "red") +
# geom_errorbar(aes(ymin = , ymax = , fill = group), width = 0.1, size = 0.8) +
facet_wrap(~group) +
labs( # title="4.2 Minimum Daily Mood across Groups",
x ="Day", y = "Mood (Likert Scale)" # , caption = "Error bars indicate 95% CI"
)
# combine average mood and min mood score into one graph
APA.mood.fig <- ggarrange(avg.graph,min.graph,
labels = c("Average","Minimum"),
ncol = 1, nrow = 2, hjust = c(-1, -0.9), vjust = 0.2,font.label = list(size = 10, face = "italic")) +
labs(title = APA.fig.title("APA.mood.fig")[1],
subtitle = APA.fig.title("APA.mood.fig")[2],
caption = expression(paste(italic("Note. "), "Error bars indicate 95% CI "))
) +
theme(# text = element_text(size = 12),
plot.title = element_text(size = 12, face = "bold",
hjust = 0,
margin = margin(b = 15)),
plot.subtitle = element_text(size = 12, face = "italic",
hjust = 0,
margin = margin(b = 15)
),
plot.caption = element_text(size = 12, hjust = 0),
# axis.text = element_text(size = 8, color = "black"),
plot.title.position = "plot",
plot.caption.position = "plot",
)
APA.mood.fig
# use plotly to display the gradual change of mood score
p <- datasets.raw$daily %>%
left_join(select(datasets$other, emory.id, group), by = "emory.id") %>%
filter(day > 0) %>%
ggplot(aes(group, mood, frame = day)) +
geom_jitter(position=position_jitter(0.2), cex=1.2) +
stat_summary(fun.data = "mean_cl_normal",
geom = "errorbar",size = 1,
width = .3, color = "red") +
stat_summary(fun = "mean", geom = "point", color = "red") +
labs(title = "Mood Score", x ="Group", y = "Score (Likert Scale)")
figure.mood.ani <- ggplotly(p) %>%
animation_opts(
1000, easing = "elastic", redraw = FALSE
) %>%
animation_slider(
currentvalue = list(prefix = "Day ", font = list(color="black"))
) %>%
layout(annotations = list(x = 1, y = -0.13, text = "Error bars indicate 95% CI",
showarrow = F, xref='paper', yref='paper',
xanchor='right', yanchor='auto', xshift=0, yshift=0,
font=list(size=12, color="gray"))
)
figure.mood.ani# remove unnecessary data
preserve <- append(preserve,c("CI.mood", "figure.mood", "figure.mood.ani"))
rm(list=ls()[! ls() %in% preserve])No significant differences in the daily, mean, and minimum mood scores were found across groups.
5.0 Who did meditation help? and Potential Confounding Variables
Since there are many variables (thousands of combinations), insignificant data are excluded from all graphs and tables in this section.
5.1 Spearman Correlation (Numeric vs Numeric)
# continuous vs continuous
# continuous <- c("GPA.previous", "grade.predict", "study.time.total")
# mean and min, treat as continuous
# continuous.2 <- c("mood.avg", "mood.min", "productivity.avg")
# discrete <- c("age", "credits", "exam.score")
# list variable names
# move AEQ subscales from values in `subscale` to column names
pivot_wider_AEQ <- function(data, name){
AEQ.name = paste0("AEQ.",name, ".")
new = select(data, emory.id, subscale, score) %>%
pivot_wider(names_prefix = AEQ.name,
names_from = "subscale",
values_from = "score",
names_sep = ".")
return(new)
}
# prepare data for spearman correlation
spearman.ori <- datasets.raw$pre_study %>%
select(1:35) %>%
# get AEQ scores
left_join(pivot_wider_AEQ(datasets$AEQ_6,6), by = "emory.id") %>%
left_join(pivot_wider_AEQ(datasets$AEQ_7,7), by = "emory.id") %>%
left_join(pivot_wider_AEQ(datasets$AEQ_14,14), by = "emory.id") %>%
left_join(pivot_wider_AEQ(datasets$difference,"difference"), by = "emory.id") %>%
# get mood avg and min scores
left_join(select(datasets$mood, emory.id, mood.type, mood.avg, mood.min) %>%
pivot_wider(names_from = "mood.type",
values_from = c("mood.avg","mood.min"),
names_sep = "."
),
by = "emory.id") %>%
# get other:
left_join(select(datasets$other,-start.date,-group),by = "emory.id") %>%
select(-(1:6)) %>%
select(-meditation.reason, -exam.date) %>%
select(!one_of(nominal)) %>%
mutate_if(is.factor, as.numeric) # the factor can affect the construction of correlation matrix
# function of correlation and contructring correlation matrix plot
# adapted from
# https://towardsdatascience.com/how-to-create-a-correlation-matrix-with-too-many-variables-309cc0c0a57
corr_simple <- function(data, sig=0.05,type = "spearman"){
#run a correlation and drop the insignificant ones
corr <- rcorr(as.matrix(data), type = type)
#prepare to drop duplicates and correlations of 1
corr$P[lower.tri(corr$P,diag=TRUE)] <- NA
#drop perfect correlations
# corr$P[corr$P == 1] <- NA
corr$P <- as.data.frame(as.table(corr$P)) %>% #turn into a 3-column table, important
na.omit() %>% #remove the NA values from above
arrange(Var1,Var2) %>% # arrange the varaibles in alphabetical order
filter(abs(Freq) <= sig) %>% #select significant values
# exclude relationships between AEQ, mood, and concern
filter(!(substr(Var1,1,4) == "AEQ." & substr(Var2,1,4) == "AEQ.")) %>%
filter(!(substr(Var1,1,5) == "mood." & substr(Var2,1,5) == "mood.")) %>%
filter(!(substr(Var1,1,8) == "concern." & substr(Var2,1,8) == "concern."))
#sort by lowest P value
corr$P <- corr$P[order(abs(corr$P$Freq)),]
# remove values from matrix P and n and convert p and n to dataframe
# corr$P <- corr$r %>%
# select(Var1,Var2)%>%
# inner_join(as.data.frame(as.table(corr$P)), by = c("Var1","Var2"))
# corr$n <- corr$r %>%
# select(Var1,Var2)%>%
# inner_join(as.data.frame(as.table(corr$n)) , by = c("Var1","Var2"))
# concatenate r, P, and n into one dataframe
result <- corr$P %>%
rename(p = Freq) %>%
inner_join(as.data.frame(as.table(corr$r)) %>% rename(r = Freq), by = c("Var1","Var2")) %>%
inner_join(as.data.frame(as.table(corr$n)) %>% rename(n = Freq), by = c("Var1","Var2")) %>% # number of pairwise cases
mutate(df = n-2)
# plot the correlation coefficient, create as a global variable
spearman.plt <<- result %>%
select(Var1,Var2,r) %>%
reshape2::acast(Var1~Var2, value.var="r") %>% # shape correlation coefficient into matrix
ggcorrplot(title = "Figure 3",
outline.col = "white", lab = F,
tl.cex = 10)
return(result)
}
spearman <- corr_simple(spearman.ori)
m <- list(
l = 0,
r = 0,
b = 100,
t = 0,
pad = 50
)# These two lists corresponds to values in columns Var1 and Var2. A for loop is used to discard these relationships.
useless1 <- c("AEQ.14.pride","difficulty","worry","AEQ.7.shame", "AEQ.14.pride","AEQ.14.pride",
"difficulty", "meditation.freq",
# not sure
"grade.predict", "AEQ.7.pride", "AEQ.difference.shame")
useless2 <- c("exam.score", "grade.predict","grade.predict","exam.score","mood.min.7~21","mood.min.1~6",
"worry", "meditation.times",
# not sure
"worry", "exam.score", "exam.score")
# remove meaningless pairs
spearman$useless = FALSE
for(i in 1:length(useless1)){
for(j in 1:nrow(spearman)) {
if(spearman[j,]$Var1 == useless1[i] & spearman[j,]$Var2 == useless2[i]){
spearman[j,]$useless = TRUE
}
}
}
spearman <- spearman[-which(spearman$useless == TRUE),] %>%
filter(abs(r) > 0.35) # filter out weak associations
# reconstruct plot
APA.spearman.fig <- spearman %>%
select(Var1,Var2,r) %>%
reshape2::acast(Var1~Var2, value.var="r") %>% # shape correlation coefficient into matrix
ggcorrplot(
outline.col = "white", lab = F,
# type = "upper",
tl.cex = 10) +
# scale_x_discrete(guide = guide_axis(n.dodge=3)) +
labs(title = APA.fig.title("APA.spearman.fig")[1],
subtitle = APA.fig.title("APA.spearman.fig")[2],
caption = expression(paste(italic("Note. "), "<Add Text Here>"))) +
theme(# text = element_text(size = 12),
plot.title = element_text(size = 12, face = "bold",
# hjust = 0.5,
margin = margin(b = 15)),
plot.subtitle = element_text(size = 12, face = "italic",
# hjust = 0.5,
margin = margin(b = 15)),
plot.caption = element_text(size = 12, hjust = 0),
# axis.text = element_text(size = 8, color = "black"),
axis.text.y = element_text(hjust = 1),
axis.text.x = element_text(angle = 30, vjust = 1, hjust = 1),
plot.title.position = "plot",
plot.caption.position = "plot",
) # + facet_wrap(~ type, ncol = 1)
ggsave(file="APA.spearman.svg", plot=APA.spearman.fig, width=10, height=10)
include_graphics("APA.spearman.svg")
# animated plt
ggplotly(APA.spearman.fig) %>% # layout(autosize = F, width = 800, height = 700, margin = m) %>%
layout(yaxis= list(showticklabels = FALSE),xaxis= list(showticklabels = FALSE))Spearman Correlation Coefficient in Dataframe and report in APA format (same as the plot)
# datatable(
# spearman, filter = 'top', extensions = 'Buttons', options = list(
# dom = 'Bfrtip',
# buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
# )
# )
spearman.report <- function(x,y,r,df,p){
str <- sprintf("Spearman’s rank correlation was computed to assess the relationship between %s and %s. There was a %s correlation between the two variables, $\\rho(%d) = %s$, $%s$.\n", x,y,ifelse(r > 0, "positive","negative"),df,num.format(r),p.format(p, table = F))
return(latex.escape(str))
}
spm.rep <- spearman %>% mutate(report = spearman.report(Var1,Var2,r,df,p))
APA.spearman <- spearman %>%
mutate(p = p.format(p),
r = num.format(r),
relationship = paste(Var1, Var2, sep = " and ")
) %>%
arrange(Var1) %>%
select(relationship, df, r, p, n) %>%
# select(-dependent.group) %>%
# format table
kable(
# digits = c(NA, 0, 2, 2, 3, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "$df$", "$\\rho$", "$p$", "$N$")),
align = c("l", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.spearman")) %>%
# use html code to format the title
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
footnote(
general_title = "Note.",
general = latex.escape("Correlations which $p\\le .05$ and $|\\rho|>0.35$ are included.", report = TRUE),
# paste("Correlations which $p < .05$ and $|\\rho| > 0.35$ are included.",paste(spm.rep$report,collapse = ""))
threeparttable = TRUE,
footnote_as_chunk = TRUE,
escape = escape.format
) %>%
column_spec(column = 1, width = "2in")
APA.spearman| \(df\) | \(\rho\) | \(p\) | \(N\) | |
|---|---|---|---|---|
| age and AEQ.6.joy | 40 | .41 | .007 | 42 |
| age and AEQ.7.pride | 46 | .35 | .014 | 48 |
| age and AEQ.6.hope | 40 | .36 | .02 | 42 |
| age and AEQ.6.hopeless | 40 | -.35 | .022 | 42 |
| age and AEQ.14.pride | 25 | .4 | .037 | 27 |
| GPA.previous and religious.level | 55 | -.38 | .004 | 57 |
| GPA.previous and AEQ.14.shame | 22 | -.46 | .023 | 24 |
| difficulty and mood.min.7~21 | 45 | -.36 | .014 | 47 |
| grade.predict and meditation.freq | 55 | .38 | .004 | 57 |
| grade.predict and AEQ.14.pride | 24 | .46 | .017 | 26 |
| worry and AEQ.difference.hopeless | 28 | -.39 | .033 | 30 |
| concern.1 and AEQ.14.anger | 24 | .52 | .006 | 26 |
| concern.2 and AEQ.14.anger | 24 | .5 | .009 | 26 |
| concern.5 and AEQ.7.shame | 45 | .35 | .015 | 47 |
| meditation.freq and relax | 56 | .56 | < .001 | 58 |
| meditation.freq and focus | 56 | .35 | .007 | 58 |
| meditation.freq and AEQ.difference.pride | 29 | -.39 | .029 | 31 |
| meditation.times and relax | 59 | .44 | < .001 | 61 |
| meditation.times and AEQ.6.hopeless | 40 | -.41 | .006 | 42 |
| meditation.times and AEQ.6.anxiety | 40 | -.38 | .013 | 42 |
| meditation.times and AEQ.difference.joy | 29 | -.39 | .031 | 31 |
| relax and focus | 59 | .46 | < .001 | 61 |
| relax and improve.score | 59 | .37 | .004 | 61 |
| focus and improve.score | 59 | .56 | < .001 | 61 |
| AEQ.6.hope and mood.avg.7~21 | 31 | .55 | < .001 | 33 |
| AEQ.6.hope and exam.score | 36 | .45 | .005 | 38 |
| AEQ.6.hopeless and exam.score | 36 | -.43 | .007 | 38 |
| AEQ.6.joy and study.time.total | 40 | .46 | .002 | 42 |
| AEQ.6.joy and productivity.avg | 40 | .43 | .004 | 42 |
| AEQ.6.joy and mood.avg.7~21 | 31 | .37 | .033 | 33 |
| AEQ.difference.anger and mood.min.7~21 | 29 | .36 | .049 | 31 |
| AEQ.difference.joy and study.time.total | 29 | .37 | .04 | 31 |
| study.time.total and productivity.avg | 59 | .35 | .005 | 61 |
| Note. Correlations which \(p\le .05\) and \(|\rho|>0.35\) are included. |
All results shown here are filtered by \(p <= 0.05\) and \(|\rho| > 0.35\) . \(\rho\) is correlation coefficient of spearman correlation.
Attitudes toward the effectiveness of meditation are interrelated. There are significant positive associations among the three types of attitudes:
- “Meditation can help me focus” and “Meditation can help me relax” (+)
- “Meditation can help me focus” and “Meditation can help me improve score” (+)
- “Meditation can help me improve score” and “Meditation can help me relax” (+)
Possible effect of more experience related to meditation:
- meditation frequency and the attitude “Meditation can help me relax” (+)
- meditation times and the attitude “Meditation can help me relax” (+)
- meditation frequency and the attitude “Meditation can help me focus”, (+)
which may indicate a possible common reason to meditate for those people who believe that meditation can help them relax or mediatation can help those who meditate more frequently and more times.
- meditation times and AEQ.6.anxiety (-)
- meditation times and AEQ.difference.joy (-) (Exception)
- meditation frequency and AEQ.difference.pride (-) (Exception)
- meditation times and AEQ.6.hopeless (-)
indicates that people with more experience in mediation may benefit more by having more positive emotions.
predicted grade and meditation frequency. People who meditate more frequently tend to predict a higher score on exam.
Moderate positive relationship related to accuracy of prediction (exam score)
- between AEQ.6.hope and average mood from day 7 to 21 (+)
- between AEQ.6.joy and average mood from day 7 to 21 (+)
- between AEQ.6.hope and exam.score (+)
- between difficulty (estimated by subjects) and minimum mood score from day 7 to 21 (-)
- concern.5 (same order as the document in google folder) and AEQ.7.shame (+)
- concern.2 and AEQ.14.anger (+)
- concern.1 and AEQ.14.anger (+)
- AEQ.6.hopeless and exam.score (-)
- predicted grade and AEQ.14.pride (+) may suggest that students on average made a relatively accurate assumption for their exam score the day before exam.
Moderate positive relationship related to study time and productivity
- AEQ.6.joy and total study time (+)
- AEQ.difference.joy and total study time (+)
- AEQ.6.joy and productivity (+)
- productivity and total study time (+)
Relationship between age and many AEQ subscales:
- AEQ.6.joy (+)
- AEQ.14.pride (+)
- AEQ.6.hope (+)
- AEQ.7.pride (+)
- AEQ.6.hopeless (-)
worry and AEQ.difference.hopeless (day 21-0) (-)
Examples:
- level 1 vs (day 21 - 0) (5-1)
- level 2 vs (2-1)
- level 3 vs (3-5)
Students who scored a higher level of worriness tend to have less increment in the level of trait hopelessness.
previous GPA and religious level (-)
- Perhaps waste more time on praying?
- Believe that supernatural forces can affect score?
- pessimistic determinism?
previous GPA and AEQ.14.shame (-)
AEQ.difference.anger and minimum mood score from day 7 to 21
5.2 Association between Nominal variables and Ordinal variables
The relationships between the variable group and other
numeric variables have been analyzed before
confounding.ori <- datasets.raw$pre_study %>%
select(1:35,73:74) %>%
# get AEQ scores
left_join(pivot_wider_AEQ(datasets$AEQ_6,6), by = "emory.id") %>%
left_join(pivot_wider_AEQ(datasets$AEQ_7,7), by = "emory.id") %>%
left_join(pivot_wider_AEQ(datasets$AEQ_14,14), by = "emory.id") %>%
left_join(pivot_wider_AEQ(datasets$difference,"difference"), by = "emory.id") %>%
# get mood avg and min scores
left_join(select(datasets$mood, emory.id, mood.type, mood.avg, mood.min) %>%
mutate(mood.type = case_when(mood.type == "1~6" ~ "1_6", TRUE ~ "7_21")) %>%
pivot_wider(names_from = "mood.type",
values_from = c("mood.avg","mood.min"),
names_sep = "."
),
by = "emory.id") %>%
# get other:
left_join(select(datasets$other,-start.date,-group),by = "emory.id") %>%
# left_join(select(productivity6,-group, -duplicate), by = "emory.id") %>%
# left_join(select(productivity7,-group), by = "emory.id") %>%
select(-(1:6)) %>%
select(-meditation.reason, -exam.date)
numeric <- names(
confounding.ori %>%
select(!one_of(nominal)) %>%
select(!one_of(ordinal))
)
##############################
nominal <- nominal[nominal != "course"] # only delete this variable temporarily
nominal.label <- nominal.label[nominal.label != "Course"]
##############################
col.nested = c("race", "major", "religion", "religion.family", "meditation.type")CI.confounding.cat = data.frame(matrix(ncol = 7, nrow = 0))
chisq.confounding = data.frame(matrix(ncol = 8, nrow = 0))
# colnames(chisq.confounding) = c("test.type","independent", "dependent","statistic", "df", "n", "p.value", "effect.size")
chisq.confounding.post <- data.frame(matrix(ncol = 6, nrow = 0))
# colnames(chisq.confounding.post) = c("independent", "independent.value", "dependent", "dependent.value", "std_residual", "p.adj", "report")
report.confounding.cat <- ""
i = 1
for(y in ordinal){ # dependent variable
var.label.ordinal <- ordinal.label[i]
j = 1
for(x in nominal){ # grouping variable
var.label.nominal <- nominal.label[j]
data <- confounding.ori %>%
select(x,y) %>%
mutate_if(is.factor, as.integer) %>% # can affect sd()
na.omit() # apply to both variables
if(x %in% col.nested){ # these are all part of nominal variables
data = data %>%
unnest(x)
}
# CI.categorical
result = CI.categorical(data, x=x, y=y)
CI.confounding.cat = rbind(CI.confounding.cat, result)
# Hypothesis Testing
result = chi(data, x=x, y=y)
chisq.confounding = rbind(chisq.confounding, result[[1]])
if(length(result) > 1){
chisq.confounding.post = rbind(chisq.confounding.post, result[[2]] %>%
mutate(relation.group = paste(nominal.label[j],ordinal.label[i],sep = " and ")) %>%
select(relation.group, 1:5)
)
report.confounding.cat <- paste(report.confounding.cat, result[[3]])
}
j = j + 1
gc(verbose = FALSE) # free some RAM
}
i = i + 1
gc(verbose = FALSE) # free some RAM
}
datatable(
CI.confounding.cat, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)# datatable(
# chisq.confounding , filter = 'top', extensions = 'Buttons', options = list(
# dom = 'Bfrtip',
# buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
# )
# )
# chisq.confounding.post <- chisq.confounding.post[-c(351:700),]Filter out insignificant data same as the previous two tables
chisq.sig.confounding <- chisq.confounding %>%
mutate(sig = case_when(p.value == "< .001" ~ TRUE,
as.numeric(p.value) <= 0.05 ~ TRUE,
TRUE ~ FALSE
)) %>%
filter(sig == TRUE) %>%
select(-sig, -test.type)
sig.group <- chisq.confounding.post %>%
mutate(# sig.group = FALSE,
sig = case_when(p.value == "< .001" ~ TRUE,
as.numeric(p.value) <= 0.05 ~ TRUE,
TRUE ~ FALSE
)) %>%
filter(sig == TRUE)
# unique(sig.group$relation.group)
# select(-sig)
chisq.sig.confounding.post <- chisq.confounding.post %>%
filter(relation.group %in% sig.group$relation.group)
discard <- chisq.sig.confounding %>% select(independent, dependent)
CI.sig.confounding.cat <- CI.confounding.cat %>%
left_join(chisq.sig.confounding %>% select(independent, dependent, p.value), by = c("independent", "dependent")) %>%
filter(!is.na(p.value))
datatable(
CI.sig.confounding.cat, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)## Chi-squared Table, APA format, theme color not changed (still gray)
chisq.tbl <- chisq.sig.confounding %>% # select(-test.type) %>%
mutate(dependent.group = case_when(startsWith(dependent, "concern") ~ "Concern",
startsWith(dependent, "meditation") ~ "Meditation habit",
startsWith(dependent, "religi") ~ "Religious affiliation",
dependent == "focus" ~ "Attitudes toward meditation",
dependent == "relax" ~ "Attitudes toward meditation",
dependent == "improve.score" ~ "Attitudes toward meditation",
startsWith(dependent, "productivity") ~ "Productivity",
startsWith(dependent, "mood") ~ "Mood score",
TRUE ~ " "
),
dependent = case_when(
startsWith(dependent, "concern") & endsWith(dependent, "1") ~ concerns.label[1],
startsWith(dependent, "concern") & endsWith(dependent, "2") ~ concerns.label[2],
startsWith(dependent, "concern") & endsWith(dependent, "3") ~ concerns.label[3],
startsWith(dependent, "concern") & endsWith(dependent, "4") ~ concerns.label[4],
startsWith(dependent, "concern") & endsWith(dependent, "5") ~ concerns.label[5],
startsWith(dependent, "concern") & endsWith(dependent, "6") ~ concerns.label[6],
# dependent == "religion.family" ~ "religious affliation of family members",
# dependent == "religion" ~ "religious affliation",
dependent == "religious.level" ~ "Level of religiosity",
dependent == "productivity6" ~ "Productivity (Day 6)",
dependent == "productivity7" ~ "Productivity (Day 7)",
dependent == "meditation.times" ~ "Meditation times",
dependent == "meditation.freq" ~ "Meditation frequency",
startsWith(dependent, "mood") & endsWith(dependent, "6") ~ "Mood score (Day 1-6)",
startsWith(dependent, "mood") & endsWith(dependent, "7") ~ "Mood score (Day 7-21)",
TRUE ~ str_to_title(dependent)
),
independent.group = case_when(
startsWith(independent, "meditation") ~ "Meditation habit",
startsWith(independent, "religi") ~ "Religious affiliation",
TRUE ~ " "
),
independent = case_when(
independent == "meditation.type" ~ "Meditation type",
independent == "religion.family" ~ "Religious affliation of family members",
independent == "religion" ~ "Religious affliation",
# dependent == "religious.level" ~ "level of religiosity",
TRUE ~ str_to_title(independent)
),
relationship.group = paste(independent.group, dependent.group, sep = " and "),
relationship = paste(independent, dependent, sep = " and ")
# statistic = as.numeric(statistic),
# effect.size = as.numeric(effect.size),
# p.value = as.numeric(p.value)
) %>%
select(dependent.group, relationship, df, statistic, effect.size, p.value, n, relationship.group) %>%
arrange(dependent.group,relationship)
index = APA.group.index(chisq.tbl)
APA.chisq.confounding <- chisq.tbl %>%
select(-relationship.group, -dependent.group) %>%
# format table
kable(
# digits = c(NA, 0, 2, 2, 3, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "$df$", "$\\chi^2$", "$\\phi_c$", "$p$", "$N$")),
align = c("l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.chisq.confounding")) %>%
# use html code to format the title
pack_rows(index = index) %>%
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
column_spec(column = 1, width = "2in")
## Frequency Table with proportion and post-hoc chisq table
index = APA.group.index(chisq.sig.confounding.post)
APA.chisq.confounding.post<- chisq.sig.confounding.post %>%
arrange(relation.group,relationship) %>%
select(-1) %>%
kable(
format = output.format,
longtable = TRUE,
col.names = latex.escape(c("", "Standard Residual", "$p$", "$N$", "proportion")),
caption = APA.tbl.caption("APA.chisq.confounding.post")) %>%
pack_rows(index = index) %>%
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
# kable_styling(full_width = TRUE) %>%
footnote(
general_title = "Note. ",
general = latex.escape("For the first Chi-squared test, only variable relationships which $p \\le .05$ are included. For the post-hoc Chi-squared test, only variable relationships which $p_{post-hoc} \\le .05$ are included. Some relationships shows significance in the first test while showing no significance in the post-hoc test. Please compare the number of incidences (N) and the proportion in the post-hoc test to analyze the direction of the significance.", report = TRUE),
threeparttable = TRUE,
footnote_as_chunk = TRUE,
escape = escape.format
) %>%
column_spec(column = 1, width = "2in")
APA.confounding.cat <- structure(paste0(APA.chisq.confounding,APA.chisq.confounding.post), format = output.format, class = ifelse(output.format == "latex", "knitr_kable", "kableExtra"))
APA.confounding.cat| \(df\) | \(\chi^2\) | \(\phi_c\) | \(p\) | \(N\) | |
|---|---|---|---|---|---|
| Attitudes toward meditation | |||||
| Meditation type and Relax | 12 | 31.37 | .34 | .002 | 89 |
| Concern | |||||
| Gender and Not getting the job or into the grad school that you hope to | 2 | 6.26 | .33 | .044 | 59 |
| Grade and Embarrassing yourself amongst friends | 2 | 8.58 | .38 | .014 | 60 |
| Meditation habit | |||||
| Gender and Meditation times | 3 | 10.03 | .41 | .018 | 61 |
| Major and Meditation frequency | 44 | 62.79 | .66 | .033 | 71 |
| Major and Meditation times | 66 | 92.41 | .65 | .018 | 73 |
| Meditation type and Meditation frequency | 8 | 18.3 | .33 | .019 | 85 |
| Meditation type and Meditation times | 12 | 63.49 | .49 | < .001 | 89 |
| Mood score | |||||
| Race and Mood.min.7_21 | 9 | 32.05 | .45 | < .001 | 52 |
| Religious affliation of family members and Mood score (Day 1-6) | 14 | 27.23 | .51 | .018 | 53 |
| Productivity | |||||
| Gender and Productivity (Day 7) | 4 | 10.13 | .46 | .038 | 48 |
| Race and Productivity (Day 6) | 15 | 27.29 | .43 | .026 | 50 |
| Religious affiliation | |||||
| Religious affliation and Level of religiosity | 28 | 74.58 | .55 | < .001 | 61 |
| Standard Residual | \(p\) | \(N\) | proportion | |
|---|---|---|---|---|
| Gender and Meditation times | ||||
| Female and 1 | 1.4 | 1 | 7 | 11.5% |
| Female and 2 | 1.63 | .829 | 24 | 39.3% |
| Female and 3 | -3.09 | .016 | 6 | 9.8% |
| Female and 4 | .41 | 1 | 3 | 4.9% |
| Male and 1 | -1.4 | 1 | 1 | 1.6% |
| Male and 2 | -1.63 | .829 | 8 | 13.1% |
| Male and 3 | 3.09 | .016 | 11 | 18% |
| Male and 4 | -.41 | 1 | 1 | 1.6% |
| Gender and Productivity (Day 7) | ||||
| Female and 1 | -3.01 | .026 | 1 | 2.1% |
| Female and 2 | .08 | 1 | 2 | 4.2% |
| Female and 3 | .43 | 1 | 11 | 22.9% |
| Female and 4 | 1.89 | .591 | 16 | 33.3% |
| Female and 5 | -.44 | 1 | 1 | 2.1% |
| Male and 1 | 3.01 | .026 | 6 | 12.5% |
| Male and 2 | -.08 | 1 | 1 | 2.1% |
| Male and 3 | -.43 | 1 | 5 | 10.4% |
| Male and 4 | -1.89 | .591 | 4 | 8.3% |
| Male and 5 | .44 | 1 | 1 | 2.1% |
| Grade and Concern: Embarrassing yourself amongst friends | ||||
| Freshman and 1 | -2.49 | .076 | 4 | 6.7% |
| Freshman and 2 | 2.83 | .028 | 23 | 38.3% |
| Freshman and 3 | -.81 | 1 | 4 | 6.7% |
| Sophomore and 1 | 2.49 | .076 | 12 | 20% |
| Sophomore and 2 | -2.83 | .028 | 11 | 18.3% |
| Sophomore and 3 | .81 | 1 | 6 | 10% |
| Major and Meditation times | ||||
| Biology and 1 | -.38 | 1 | 0 | 0% |
| Biology and 2 | .99 | 1 | 1 | 1.4% |
| Biology and 3 | -.66 | 1 | 0 | 0% |
| Biology and 4 | -.27 | 1 | 0 | 0% |
| Anthropology and 1 | -.38 | 1 | 0 | 0% |
| Anthropology and 2 | .99 | 1 | 1 | 1.4% |
| Anthropology and 3 | -.66 | 1 | 0 | 0% |
| Anthropology and 4 | -.27 | 1 | 0 | 0% |
| Biology and 1 | .54 | 1 | 1 | 1.4% |
| Biology and 2 | 1.36 | 1 | 4 | 5.5% |
| Biology and 3 | -1.52 | 1 | 0 | 0% |
| Biology and 4 | -.63 | 1 | 0 | 0% |
| Business and 1 | -.66 | 1 | 0 | 0% |
| Business and 2 | -.61 | 1 | 1 | 1.4% |
| Business and 3 | .12 | 1 | 1 | 1.4% |
| Business and 4 | 1.85 | 1 | 1 | 1.4% |
| Chemistry and 1 | 2.92 | .319 | 2 | 2.7% |
| Chemistry and 2 | -.61 | 1 | 1 | 1.4% |
| Chemistry and 3 | -1.16 | 1 | 0 | 0% |
| Chemistry and 4 | -.48 | 1 | 0 | 0% |
| Computer Science and 1 | -.38 | 1 | 0 | 0% |
| Computer Science and 2 | -1.02 | 1 | 0 | 0% |
| Computer Science and 3 | 1.53 | 1 | 1 | 1.4% |
| Computer Science and 4 | -.27 | 1 | 0 | 0% |
| Economics and 1 | -.77 | 1 | 0 | 0% |
| Economics and 2 | -1.06 | 1 | 1 | 1.4% |
| Economics and 3 | 2.01 | 1 | 3 | 4.1% |
| Economics and 4 | -.56 | 1 | 0 | 0% |
| Finance and 1 | 1.13 | 1 | 1 | 1.4% |
| Finance and 2 | -.61 | 1 | 1 | 1.4% |
| Finance and 3 | .12 | 1 | 1 | 1.4% |
| Finance and 4 | -.48 | 1 | 0 | 0% |
| Human Biology and 1 | -.38 | 1 | 0 | 0% |
| Human Biology and 2 | .99 | 1 | 1 | 1.4% |
| Human Biology and 3 | -.66 | 1 | 0 | 0% |
| Human Biology and 4 | -.27 | 1 | 0 | 0% |
| Human Health and 1 | 1.64 | 1 | 1 | 1.4% |
| Human Health and 2 | -.02 | 1 | 1 | 1.4% |
| Human Health and 3 | -.94 | 1 | 0 | 0% |
| Human Health and 4 | -.39 | 1 | 0 | 0% |
| International Studies and 1 | -.38 | 1 | 0 | 0% |
| International Studies and 2 | .99 | 1 | 1 | 1.4% |
| International Studies and 3 | -.66 | 1 | 0 | 0% |
| International Studies and 4 | -.27 | 1 | 0 | 0% |
| Medial Anthropology and 1 | -.38 | 1 | 0 | 0% |
| Medial Anthropology and 2 | -1.02 | 1 | 0 | 0% |
| Medial Anthropology and 3 | -.66 | 1 | 0 | 0% |
| Medial Anthropology and 4 | 3.71 | .019 | 1 | 1.4% |
| Music and 1 | -.38 | 1 | 0 | 0% |
| Music and 2 | -1.02 | 1 | 0 | 0% |
| Music and 3 | 1.53 | 1 | 1 | 1.4% |
| Music and 4 | -.27 | 1 | 0 | 0% |
| Neuroscience and 1 | -.54 | 1 | 0 | 0% |
| Neuroscience and 2 | -.02 | 1 | 1 | 1.4% |
| Neuroscience and 3 | .62 | 1 | 1 | 1.4% |
| Neuroscience and 4 | -.39 | 1 | 0 | 0% |
| Neuroscience and Behavioral Biology and 1 | 1.01 | 1 | 3 | 4.1% |
| Neuroscience and Behavioral Biology and 2 | 1.39 | 1 | 10 | 13.7% |
| Neuroscience and Behavioral Biology and 3 | -1.59 | 1 | 2 | 2.7% |
| Neuroscience and Behavioral Biology and 4 | -1.18 | 1 | 0 | 0% |
| Nursing and 1 | -.66 | 1 | 0 | 0% |
| Nursing and 2 | 1.74 | 1 | 3 | 4.1% |
| Nursing and 3 | -1.16 | 1 | 0 | 0% |
| Nursing and 4 | -.48 | 1 | 0 | 0% |
| Physical Therapy and 1 | -.38 | 1 | 0 | 0% |
| Physical Therapy and 2 | -1.02 | 1 | 0 | 0% |
| Physical Therapy and 3 | 1.53 | 1 | 1 | 1.4% |
| Physical Therapy and 4 | -.27 | 1 | 0 | 0% |
| Political Science and 1 | -.66 | 1 | 0 | 0% |
| Political Science and 2 | -.61 | 1 | 1 | 1.4% |
| Political Science and 3 | 1.41 | 1 | 2 | 2.7% |
| Political Science and 4 | -.48 | 1 | 0 | 0% |
| Pre-6(Biology Major (Area 6: Pre-Med, Pre-Dental, Pre-Vet)) and 1 | -.38 | 1 | 0 | 0% |
| Pre-6(Biology Major (Area 6: Pre-Med, Pre-Dental, Pre-Vet)) and 2 | -1.02 | 1 | 0 | 0% |
| Pre-6(Biology Major (Area 6: Pre-Med, Pre-Dental, Pre-Vet)) and 3 | -.66 | 1 | 0 | 0% |
| Pre-6(Biology Major (Area 6: Pre-Med, Pre-Dental, Pre-Vet)) and 4 | 3.71 | .019 | 1 | 1.4% |
| Psychology and 1 | -.75 | 1 | 1 | 1.4% |
| Psychology and 2 | .23 | 1 | 8 | 11% |
| Psychology and 3 | .93 | 1 | 6 | 8.2% |
| Psychology and 4 | -1.18 | 1 | 0 | 0% |
| Quantitative Science and 1 | -.77 | 1 | 0 | 0% |
| Quantitative Science and 2 | -1.06 | 1 | 1 | 1.4% |
| Quantitative Science and 3 | .89 | 1 | 2 | 2.7% |
| Quantitative Science and 4 | 1.48 | 1 | 1 | 1.4% |
| Spanish and 1 | -.38 | 1 | 0 | 0% |
| Spanish and 2 | -1.02 | 1 | 0 | 0% |
| Spanish and 3 | -.66 | 1 | 0 | 0% |
| Spanish and 4 | 3.71 | .019 | 1 | 1.4% |
| WGS and 1 | -.38 | 1 | 0 | 0% |
| WGS and 2 | -1.02 | 1 | 0 | 0% |
| WGS and 3 | 1.53 | 1 | 1 | 1.4% |
| WGS and 4 | -.27 | 1 | 0 | 0% |
| Meditation type and Meditation frequency | ||||
| Focused attention and 1 | -.4 | 1 | 11 | 12.9% |
| Focused attention and 2 | .49 | 1 | 4 | 4.7% |
| Focused attention and 3 | .03 | 1 | 5 | 5.9% |
| Loving-kindness and 1 | -1.7 | 1 | 2 | 2.4% |
| Loving-kindness and 2 | .9 | 1 | 2 | 2.4% |
| Loving-kindness and 3 | 1.16 | 1 | 3 | 3.5% |
| Mindfulness and 1 | -2.31 | .314 | 16 | 18.8% |
| Mindfulness and 2 | 1.23 | 1 | 8 | 9.4% |
| Mindfulness and 3 | 1.58 | 1 | 12 | 14.1% |
| Never tried meditation and 1 | 2.31 | .313 | 7 | 8.2% |
| Never tried meditation and 2 | -1.23 | 1 | 0 | 0% |
| Never tried meditation and 3 | -1.58 | 1 | 0 | 0% |
| Not sure and 1 | 2.99 | .041 | 14 | 16.5% |
| Not sure and 2 | -1.9 | .871 | 0 | 0% |
| Not sure and 3 | -1.79 | 1 | 1 | 1.2% |
| Meditation type and Meditation helps me relax | ||||
| Focused attention and 2 | -.98 | 1 | 0 | 0% |
| Focused attention and 3 | -1.85 | 1 | 2 | 2.2% |
| Focused attention and 4 | 1.61 | 1 | 15 | 16.9% |
| Focused attention and 5 | .48 | 1 | 4 | 4.5% |
| Loving-kindness and 2 | -.51 | 1 | 0 | 0% |
| Loving-kindness and 3 | .25 | 1 | 2 | 2.2% |
| Loving-kindness and 4 | -.74 | 1 | 3 | 3.4% |
| Loving-kindness and 5 | .97 | 1 | 2 | 2.2% |
| Mindfulness and 2 | -.33 | 1 | 1 | 1.1% |
| Mindfulness and 3 | -2.18 | .582 | 5 | 5.6% |
| Mindfulness and 4 | 1.15 | 1 | 24 | 27% |
| Mindfulness and 5 | 1.19 | 1 | 8 | 9% |
| Never tried meditation and 2 | -.51 | 1 | 0 | 0% |
| Never tried meditation and 3 | 3.9 | .002 | 6 | 6.7% |
| Never tried meditation and 4 | -2.33 | .399 | 1 | 1.1% |
| Never tried meditation and 5 | -1.19 | 1 | 0 | 0% |
| Not sure and 2 | 2.23 | .509 | 2 | 2.2% |
| Not sure and 3 | 1.95 | 1 | 7 | 7.9% |
| Not sure and 4 | -1.11 | 1 | 7 | 7.9% |
| Not sure and 5 | -1.91 | 1 | 0 | 0% |
| Meditation type and Meditation times | ||||
| Focused attention and 1 | -1.65 | 1 | 0 | 0% |
| Focused attention and 2 | .19 | 1 | 11 | 12.4% |
| Focused attention and 3 | .21 | 1 | 7 | 7.9% |
| Focused attention and 4 | .97 | 1 | 3 | 3.4% |
| Loving-kindness and 1 | -.87 | 1 | 0 | 0% |
| Loving-kindness and 2 | -1.21 | 1 | 2 | 2.2% |
| Loving-kindness and 3 | 1.52 | 1 | 4 | 4.5% |
| Loving-kindness and 4 | .51 | 1 | 1 | 1.1% |
| Mindfulness and 1 | -2.56 | .21 | 0 | 0% |
| Mindfulness and 2 | .34 | 1 | 20 | 22.5% |
| Mindfulness and 3 | .94 | 1 | 14 | 15.7% |
| Mindfulness and 4 | .44 | 1 | 4 | 4.5% |
| Never tried meditation and 1 | 7.39 | < .001 | 6 | 6.7% |
| Never tried meditation and 2 | -2 | .91 | 1 | 1.1% |
| Never tried meditation and 3 | -1.87 | 1 | 0 | 0% |
| Never tried meditation and 4 | -.87 | 1 | 0 | 0% |
| Not sure and 1 | .54 | 1 | 2 | 2.2% |
| Not sure and 2 | 1.61 | 1 | 11 | 12.4% |
| Not sure and 3 | -1.21 | 1 | 3 | 3.4% |
| Not sure and 4 | -1.39 | 1 | 0 | 0% |
| Race and Min mood score (Day 7-21) | ||||
| Asian and 1 | -1.45 | 1 | 5 | 9.6% |
| Asian and 2 | 1.09 | 1 | 17 | 32.7% |
| Asian and 3 | 1.46 | 1 | 7 | 13.5% |
| Asian and 4 | -2 | .722 | 0 | 0% |
| Black or African American and 1 | -1.03 | 1 | 0 | 0% |
| Black or African American and 2 | -.66 | 1 | 1 | 1.9% |
| Black or African American and 3 | -.82 | 1 | 0 | 0% |
| Black or African American and 4 | 4.66 | < .001 | 2 | 3.8% |
| Other and 1 | .34 | 1 | 1 | 1.9% |
| Other and 2 | -.66 | 1 | 1 | 1.9% |
| Other and 3 | -.82 | 1 | 0 | 0% |
| Other and 4 | 2.11 | .559 | 1 | 1.9% |
| White and 1 | 1.88 | .967 | 7 | 13.5% |
| White and 2 | -.49 | 1 | 8 | 15.4% |
| White and 3 | -.74 | 1 | 2 | 3.8% |
| White and 4 | -1.24 | 1 | 0 | 0% |
| Race and Productivity (Day 6) | ||||
| Asian and 1 | -1.14 | 1 | 0 | 0% |
| Asian and 2 | -.17 | 1 | 1 | 2% |
| Asian and 3 | .89 | 1 | 11 | 22% |
| Asian and 3.5 | .38 | 1 | 2 | 4% |
| Asian and 4 | -1.14 | 1 | 12 | 24% |
| Asian and 5 | 1.28 | 1 | 2 | 4% |
| Black or African American and 1 | 4 | .002 | 1 | 2% |
| Black or African American and 2 | 2.67 | .18 | 1 | 2% |
| Black or African American and 3 | -1.28 | 1 | 0 | 0% |
| Black or African American and 3.5 | -.45 | 1 | 0 | 0% |
| Black or African American and 4 | -.6 | 1 | 1 | 2% |
| Black or African American and 5 | -.36 | 1 | 0 | 0% |
| Other and 1 | -.26 | 1 | 0 | 0% |
| Other and 2 | -.36 | 1 | 0 | 0% |
| Other and 3 | -.03 | 1 | 1 | 2% |
| Other and 3.5 | -.45 | 1 | 0 | 0% |
| Other and 4 | .6 | 1 | 2 | 4% |
| Other and 5 | -.36 | 1 | 0 | 0% |
| White and 1 | -.69 | 1 | 0 | 0% |
| White and 2 | -.99 | 1 | 0 | 0% |
| White and 3 | -.28 | 1 | 5 | 10% |
| White and 3.5 | .05 | 1 | 1 | 2% |
| White and 4 | 1.21 | 1 | 10 | 20% |
| White and 5 | -.99 | 1 | 0 | 0% |
| Religious affiliation and Religious level | ||||
| Buddhism and 1 | -2.15 | 1 | 0 | 0% |
| Buddhism and 2 | .23 | 1 | 1 | 1.6% |
| Buddhism and 3 | 2.85 | .173 | 4 | 6.6% |
| Buddhism and 4 | -.62 | 1 | 0 | 0% |
| Buddhism and 5 | -.53 | 1 | 0 | 0% |
| Christianity and 1 | -3.37 | .03 | 0 | 0% |
| Christianity and 2 | -.72 | 1 | 1 | 1.6% |
| Christianity and 3 | 1.6 | 1 | 5 | 8.2% |
| Christianity and 4 | 3.07 | .087 | 3 | 4.9% |
| Christianity and 5 | 2.25 | .986 | 2 | 3.3% |
| Hinduism and 1 | -2.38 | .7 | 0 | 0% |
| Hinduism and 2 | 1.18 | 1 | 2 | 3.3% |
| Hinduism and 3 | 1.39 | 1 | 3 | 4.9% |
| Hinduism and 4 | 1.05 | 1 | 1 | 1.6% |
| Hinduism and 5 | -.59 | 1 | 0 | 0% |
| Islam and 1 | .12 | 1 | 1 | 1.6% |
| Islam and 2 | -.64 | 1 | 0 | 0% |
| Islam and 3 | -.86 | 1 | 0 | 0% |
| Islam and 4 | -.38 | 1 | 0 | 0% |
| Islam and 5 | 3 | .109 | 1 | 1.6% |
| Judaism and 1 | -.45 | 1 | 1 | 1.6% |
| Judaism and 2 | .81 | 1 | 1 | 1.6% |
| Judaism and 3 | .29 | 1 | 1 | 1.6% |
| Judaism and 4 | -.47 | 1 | 0 | 0% |
| Judaism and 5 | -.4 | 1 | 0 | 0% |
| None and 1 | 6.01 | < .001 | 25 | 41% |
| None and 2 | -1.21 | 1 | 3 | 4.9% |
| None and 3 | -3.85 | .005 | 1 | 1.6% |
| None and 4 | -1.97 | 1 | 0 | 0% |
| None and 5 | -1.69 | 1 | 0 | 0% |
| Other and 1 | -1.32 | 1 | 0 | 0% |
| Other and 2 | -.64 | 1 | 0 | 0% |
| Other and 3 | 2.41 | .635 | 2 | 3.3% |
| Other and 4 | -.38 | 1 | 0 | 0% |
| Other and 5 | -.33 | 1 | 0 | 0% |
| Prefer not to specify and 1 | -.45 | 1 | 1 | 1.6% |
| Prefer not to specify and 2 | 2.41 | .634 | 2 | 3.3% |
| Prefer not to specify and 3 | -1.06 | 1 | 0 | 0% |
| Prefer not to specify and 4 | -.47 | 1 | 0 | 0% |
| Prefer not to specify and 5 | -.4 | 1 | 0 | 0% |
| Religious affiliation of family members and Min mood score (Day 1-6) | ||||
| Buddhism and 1 | 3.67 | .006 | 4 | 7.5% |
| Buddhism and 2 | -2.39 | .4 | 0 | 0% |
| Buddhism and 3 | -.52 | 1 | 1 | 1.9% |
| Christianity and 1 | -1.04 | 1 | 2 | 3.8% |
| Christianity and 2 | .48 | 1 | 10 | 18.9% |
| Christianity and 3 | .36 | 1 | 6 | 11.3% |
| Hinduism and 1 | -1 | 1 | 0 | 0% |
| Hinduism and 2 | -1.08 | 1 | 1 | 1.9% |
| Hinduism and 3 | 2.03 | 1 | 3 | 5.7% |
| Islam and 1 | 1.15 | 1 | 1 | 1.9% |
| Islam and 2 | -.03 | 1 | 1 | 1.9% |
| Islam and 3 | -.95 | 1 | 0 | 0% |
| Judaism and 1 | -.86 | 1 | 0 | 0% |
| Judaism and 2 | .56 | 1 | 2 | 3.8% |
| Judaism and 3 | .12 | 1 | 1 | 1.9% |
| None and 1 | -.29 | 1 | 3 | 5.7% |
| None and 2 | 1.64 | 1 | 12 | 22.6% |
| None and 3 | -1.54 | 1 | 3 | 5.7% |
| Other and 1 | -.7 | 1 | 0 | 0% |
| Other and 2 | -1.47 | 1 | 0 | 0% |
| Other and 3 | 2.19 | .681 | 2 | 3.8% |
| Prefer not to specify and 1 | -.49 | 1 | 0 | 0% |
| Prefer not to specify and 2 | .99 | 1 | 1 | 1.9% |
| Prefer not to specify and 3 | -.66 | 1 | 0 | 0% |
| Note. For the first Chi-squared test, only variable relationships which \(p \le .05\) are included. For the post-hoc Chi-squared test, only variable relationships which \(p_{post-hoc} \le .05\) are included. Some relationships shows significance in the first test while showing no significance in the post-hoc test. Please compare the number of incidences (N) and the proportion in the post-hoc test to analyze the direction of the significance. | ||||
5.3 Association between Nominal variables and Numeric variables (including AEQ)
kruskal.confounding = data.frame(matrix(ncol = 8, nrow = 0))
# colnames(kruskal.confounding) = c("test.type","independent","dependent","statistic","df","n","p.value","effect.size")
wilcox.confounding = data.frame(matrix(ncol = 12, nrow = 0))
# colnames(wilcox.confounding) = c("test.type","independent","independent1","independent2","dependent","statistic","n1","n2","mdn1", "mdn2","p.value","effect.size")
CI.confounding.num = data.frame(matrix(ncol = 6, nrow = 0))
report.confounding.num <- ""
# define new CI
CI.confounding <- function(data, x = "group", y){
new = data %>%
mutate_if(is.factor, as.integer) %>% # can affect sd()
group_by_(x) %>%
filter(!is.na(.data[[y]])) %>%
summarise(
mean = mean(.data[[y]]),
sd = sd(.data[[y]]),
n = n(),
error = qnorm(0.975) - sd/sqrt(n),
lower = mean - error,
lower = case_when(lower < 0 ~ 0, TRUE ~ lower),
upper = mean + error) %>%
ungroup() %>%
mutate(independent = x, dependent = y) %>%
rename(independent.value = x) %>% # need to include a column of grouping variable name since we are dealing with multiple grouping variables instead of one
select(independent,independent.value,dependent, "mean", "lower", "upper")
return(new)
}
# beware of your RAM
for(y in numeric){ # dependent variable
for(x in nominal){ # grouping variable
# subset data set
data = confounding.ori %>%
select(x,y) %>%
mutate_if(is.factor, as.integer) %>% # can affect sd()
na.omit() # apply to both variables
# unnest variables
if(x %in% col.nested){
data = data %>%
unnest(x)
}
# CI
tmp = CI.confounding(data, x, y)
CI.confounding.num = rbind(CI.confounding.num, tmp)
# Hypothesis Testing
tmp = kruskal_wilcox(data, x, y)
kruskal.confounding = rbind(kruskal.confounding, tmp[[1]])
if(length(tmp) > 1){
wilcox.confounding = rbind(wilcox.confounding, tmp[[2]])
report.confounding.num = paste(report.confounding.num, tmp[[3]])
}
gc(verbose = FALSE) # free some RAM
}
gc(verbose = FALSE) # free some RAM
}
# Display
datatable(
CI.confounding.num, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)# datatable(
# kruskal.confounding %>% filter(p.value <= 0.05), filter = 'top', extensions = 'Buttons', options = list(
# dom = 'Bfrtip',
# buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
# )
# )
#
# datatable(
# wilcox.confounding %>% filter(p.value <= 0.05), filter = 'top', extensions = 'Buttons', options = list(
# dom = 'Bfrtip',
# buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
# )
# )
# Format tables into APA format
kruskal.sig.confounding <- kruskal.confounding %>%
# exclude insignificant data
mutate(sig = case_when(p.value == "< .001" ~ TRUE,
as.numeric(p.value) <= 0.05 ~ TRUE,
TRUE ~ FALSE
)) %>%
filter(sig == TRUE) %>%
select(-sig)
APA.kruskal.confounding <- kruskal.sig.confounding %>%
select(-test.type) %>%
mutate(dependent = case_when(dependent == "age" ~ "Age", TRUE ~ dependent),
relationship = paste(str_to_title(independent), dependent, sep = " vs. ")
) %>%
select(relationship,df, statistic, effect.size, p.value, n) %>%
arrange(relationship) %>%
# format table
kable(
# digits = c(NA, 0, 2, 2, 3, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "$df$", "$H$", "$\\eta^2$", "$p$", "$N$")),
align = c("l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.kruskal.confounding")) %>%
# use html code to format the title
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
column_spec(column = 1, width = "2in")
# Wilcox
wilcox.tbl <- wilcox.confounding %>% mutate(
dependent = case_when(dependent == "age" ~ "Age", TRUE ~ dependent),
group = sprintf("%s vs %s", str_to_title(independent), dependent),
pair = sprintf("%s - %s", independent1, independent2),
n = sprintf("%s - %s", n1, n2),
mdn = sprintf("%s - %s", mdn1, mdn2)
) %>%
select(group,pair, statistic, effect.size, p.value, n, mdn) %>%
arrange(group)
index = APA.group.index(wilcox.tbl)
# wilcox.example <- wilcox.tbl
APA.wilcox.confounding <- wilcox.tbl %>%
select(-group) %>%
# select(-dependent.group) %>%
# format table
kable(
# digits = c(NA, 2, 2, 3, 0, 0),
format = output.format,
booktabs = TRUE,
escape = escape.format,
longtable = TRUE,
col.names = latex.escape(c("", "$W$", "$r$", "$p$", "$N$", "Median")),
align = c("l", "c", "c", "c", "c", "c"), # left, centered, right
caption = APA.tbl.caption("APA.wilcox.confounding")) %>%
# use html code to format the title
kable_classic(full_width = TRUE, html_font = "Cambria") %>%
footnote(
general_title = "Note.",
general = latex.escape("For Kruskal-Wallis Test, only results which $p \\le .05$ are included. All of the results of the post-hoc test are included. The significant difference in AEQ pride score on Day 14 between Asian and White students can be due to non-representative or biased selection on Asian students (International Asian students can bias the result). ", report = TRUE),
threeparttable = TRUE,
footnote_as_chunk = TRUE,
escape = escape.format
) %>%
pack_rows(index = index) %>%
column_spec(column = 1, width = "2in")
APA.confounding.num <- structure(paste0(APA.kruskal.confounding,APA.wilcox.confounding), format = output.format, class = ifelse(output.format == "latex", "knitr_kable", "kableExtra"))# concatenate the html code and stored still as a kableExtra object, easier for display
APA.confounding.num| \(df\) | \(H\) | \(\eta^2\) | \(p\) | \(N\) | |
|---|---|---|---|---|---|
| Gender vs. AEQ.6.anxiety | 1 | 4.23 | .08 | .04 | 42 |
| Grade vs. AEQ.7.anger | 1 | 4.08 | .07 | .043 | 47 |
| Grade vs. AEQ.7.pride | 1 | 3.88 | .06 | .049 | 48 |
| Grade vs. AEQ.difference.anger | 1 | 4.13 | .11 | .042 | 31 |
| Grade vs. Age | 1 | 30.13 | .49 | < .001 | 61 |
| Grade vs. GPA.previous | 1 | 9.78 | .16 | .002 | 57 |
| Race vs. AEQ.14.pride | 3 | 10.13 | .31 | .018 | 27 |
| \(W\) | \(r\) | \(p\) | \(N\) | Median | |
|---|---|---|---|---|---|
| Gender vs AEQ.6.anxiety | |||||
| Female - Male | 124.5 | .32 | .041 | 27 - 15 | 15 - 13 |
| Grade vs AEQ.7.anger | |||||
| Freshman - Sophomore | 182 | .29 | .044 | 24 - 23 | 7 - 6 |
| Grade vs AEQ.7.pride | |||||
| Freshman - Sophomore | 382.5 | .28 | .05 | 25 - 23 | 13 - 15 |
| Grade vs AEQ.difference.anger | |||||
| Freshman - Sophomore | 167.5 | .36 | .044 | 18 - 13 | -3.5 - -1 |
| Grade vs Age | |||||
| Freshman - Sophomore | 808 | .7 | < .001 | 32 - 29 | 18 - 19 |
| Grade vs GPA.previous | |||||
| Freshman - Sophomore | 211.5 | .41 | .002 | 28 - 29 | 3.92 - 3.7 |
| Race vs AEQ.14.pride | |||||
| Asian - Black or African American | 25 | .07 | .974 | 15 - 3 | 13 - 15 |
| Asian - Other | .5 | .38 | .312 | 15 - 1 | 13 - 8 |
| Asian - White | 18 | .57 | .043 | 15 - 8 | 13 - 8 |
| Black or African American - Other | 0 | .71 | .519 | 3 - 1 | 15 - 8 |
| Black or African American - White | 1.5 | .65 | .12 | 3 - 8 | 15 - 8 |
| Other - White | 4 | 0 | 1 | 1 - 8 | 8 - 8 |
| Note. For Kruskal-Wallis Test, only results which \(p \le .05\) are included. All of the results of the post-hoc test are included. The significant difference in AEQ pride score on Day 14 between Asian and White students can be due to non-representative or biased selection on Asian students (International Asian students can bias the result). | |||||
# Interpretation, APA format
cat(report.confounding.num)## A Kruskal-Wallis test showed that the mean rank of age was significantly different between the levels of grade, $H(1,N=61)=30.13$, $p < .001$, $\eta^2=.49$. Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of age between grade: Freshman ($Mdn_1=18.00$) and grade: Sophomore ($Mdn_2=19.00$), $W(n_1 = 32, n_2 = 29)=808$, $p < .001$, $r = .7$. A Kruskal-Wallis test showed that the mean rank of GPA.previous was significantly different between the levels of grade, $H(1,N=57)=9.78$, $p = .002$, $\eta^2=.16$. Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of GPA.previous between grade: Freshman ($Mdn_1=3.92$) and grade: Sophomore ($Mdn_2=3.70$), $W(n_1 = 28, n_2 = 29)=211.5$, $p = .002$, $r = .41$. A Kruskal-Wallis test showed that the mean rank of AEQ.6.anxiety was significantly different between the levels of gender, $H(1,N=42)=4.23$, $p = .04$, $\eta^2=.08$. Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of AEQ.6.anxiety between gender: Female ($Mdn_1=15.00$) and gender: Male ($Mdn_2=13.00$), $W(n_1 = 27, n_2 = 15)=124.5$, $p = .041$, $r = .32$. A Kruskal-Wallis test showed that the mean rank of AEQ.7.anger was significantly different between the levels of grade, $H(1,N=47)=4.08$, $p = .043$, $\eta^2=.07$. Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of AEQ.7.anger between grade: Freshman ($Mdn_1=7.00$) and grade: Sophomore ($Mdn_2=6.00$), $W(n_1 = 24, n_2 = 23)=182$, $p = .044$, $r = .29$. A Kruskal-Wallis test showed that the mean rank of AEQ.7.pride was significantly different between the levels of grade, $H(1,N=48)=3.88$, $p = .049$, $\eta^2=.06$. Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of AEQ.7.pride between grade: Freshman ($Mdn_1=13.00$) and grade: Sophomore ($Mdn_2=15.00$), $W(n_1 = 25, n_2 = 23)=382.5$, $p = .05$, $r = .28$. A Kruskal-Wallis test showed that the mean rank of AEQ.14.pride was significantly different between the levels of race, $H(3,N=27)=10.13$, $p = .018$, $\eta^2=.31$. Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of AEQ.14.pride between race: Asian ($Mdn_1=13.00$) and race: White ($Mdn_2=8.00$), $W(n_1 = 15, n_2 = 8)=18$, $p = .043$, $r = .57$. A Kruskal-Wallis test showed that the mean rank of AEQ.difference.anger was significantly different between the levels of grade, $H(1,N=31)=4.13$, $p = .042$, $\eta^2=.11$. Post-hoc pairwise Wilcoxon rank sum test with Bonferroni adjustment revealed a significant difference in the mean rank of AEQ.difference.anger between grade: Freshman ($Mdn_1=-3.50$) and grade: Sophomore ($Mdn_2=-1.00$), $W(n_1 = 18, n_2 = 13)=167.5$, $p = .044$, $r = .36$.# remove unnecessary data
preserve <- append(preserve,c("pivot_wider_AEQ","spearman", "spearman.plt", "CI.confounding.num", "kruskal.confounding", "wilcox.confounding", "CI.confounding.cat", "chisq.confounding"))
rm(list=ls()[! ls() %in% preserve])Summary and Export of APA Tables and Figures
You can copy tables and figures here into word document.
html_format <- function(tbl.name){ # adjust the Note. section of tables
pattern <- '<tfoot><tr><td style=\"padding: 0; \" colspan=\"100%\">'
replacement <- '<tfoot><tr><td style=\"padding: 0; \" colspan=\"10\">' # The largest column number is 8. 0 (as flexible and largest) is not supported by chrome and Rstudio
tbl <- get(tbl.name) %>%
kable_classic(full_width = T, html_font = 'sans-serif')
return(sub(pattern, replacement, tbl))
# structure(str, format = format, class = "kableExtra")
}
if(escape.format == TRUE){
tbls <- lapply(APA.table.names[names(APA.table.names)%in% names(APA.table.titles)]
, html_format) # revise all of the tables (except the tables that are already connected), store as a html string
names(tbls) <- APA.table.names[names(APA.table.names)%in% names(APA.table.titles)]
list2env(tbls, envir = .GlobalEnv) # revise all of the tables (except the tables that are already connected), revise the individual tables
structure(paste(tbls, collapse = ""), format = "html", class = "kableExtra") # print all tables together
}
# Save figures into svg files
for(fig in APA.figure.names[2:3]){
ggsave(file=paste0(fig, ".svg"), plot=get(fig), width=10, height=10)
}
fig <- "APA.AEQ.fig"
ggsave(file=paste0(fig, ".svg"), plot=get(fig), width=10, height=14)
include_graphics(paste0("APA.AEQ.fig", ".svg"))
include_graphics(paste0(APA.figure.names[2], ".svg"))
include_graphics(paste0(APA.figure.names[3], ".svg"))
# store tables and figures into Rdata file
# if(output.format == "latex"){
save("APA.table.names", "APA.figure.names","APA.table.titles", "APA.figure.titles",list = c(APA.table.names,APA.figure.names) ,file = "APA.RData")
#}6.0 Summary of All Data (without daily mood score)
combined_data <- datasets.raw$pre_study %>% # without daily surveys
select(1:35) %>%
# get AEQ scores
left_join(pivot_wider_AEQ(datasets$AEQ_6,6), by = "emory.id") %>%
left_join(pivot_wider_AEQ(datasets$AEQ_7,7), by = "emory.id") %>%
left_join(pivot_wider_AEQ(datasets$AEQ_14,14), by = "emory.id") %>%
left_join(pivot_wider_AEQ(datasets$difference,"difference"), by = "emory.id") %>%
# get mood avg and min scores
left_join(select(datasets$mood, emory.id, mood.type, mood.avg, mood.min) %>%
pivot_wider(names_from = "mood.type",
values_from = c("mood.avg","mood.min"),
names_sep = "."
),
by = "emory.id") %>%
# get other:
left_join(select(datasets$other,emory.id, start.date,group, exam.score, study.time.total),by = "emory.id") %>%
select(-(1:3),-(5:6)) %>%
left_join(select(productivity6,-group), by = "emory.id") %>%
left_join(select(productivity7,-group), by = "emory.id")
# mutate(group = factor(group, levels = c("Control","20 minuted 1X/day","10 minutes 2X/day")))
# select(meditation.reason, exam.date) %>%
# select(!one_of(nominal))
datatable(
combined_data, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)7.0 Experimental Exploration of Non-linear Regression accoring to 2C analysis document (Not encouraged for formal use)
Non-parametric methods of regression are non-linear regression models. Here, generalized additive model is used.
The model generally follows the format of y ~ x1 + x2.
This is a simple analysis of two independent variables and one dependent
variable, namely, AEQ.scores ~ group + baseline.data.
Generalized (Approximate) Cross Validation (GCV or GACV) method is
used.
bs(define smooth terms):use
bs = "fs"for all categorical variables (coerce them into factors, including nominal variables)use default setting
bs = "tp"for exam score and hours studying totalFrom Rdocumentation: “These are low rank isotropic smoothers of any number of covariates. By isotropic is meant that rotation of the covariate co-ordinate system will not change the result of smoothing. By low rank is meant that they have far fewer coefficients than there are data to smooth. They are reduced rank versions of the thin plate splines and use the thin plate spline penalty. They are the default smooth for s terms because there is a defined sense in which they are the optimal smoother of any given basis dimension/rank (Wood, 2003).”
Further model checks should be runned.
gam.check(model)plot(b,pages=1,residuals=TRUE) # show partial residualsplot(b,pages=1,seWithMean=TRUE) # 'with intercept' CIs
More detailed information on this website.
Only significant regression models are included in the output. The data set includes all models.
# AEQs <- c("AEQ_6","AEQ_7","AEQ_14","difference")
AEQs <- names(combined_data)[31:50]
var.types <- c("NA","demographic.cat.1","demographic.cat.2","demographic.num", "productivity", "exam.score_study.time")
regression_data <- combined_data
baseline.vars = append("NA", baseline.vars)
names(baseline.vars) <- var.types
# models <- NULL # 5 dimension
# baseline.vars = baseline.vars[-c(1:11)]
models.df <- data.frame(matrix(ncol = 6, nrow = 0))
colnames(models.df) = c("baseline.var.type","baseline.var","AEQ.subscale","coefficient","smooth.sig")
i = 1
report <- c()
for(baseline.var.type in baseline.vars){
for(baseline.var in baseline.var.type){
# print(baseline.var)
if(var.types[i] == "demographic.cat.2"){
regression_unnest <- regression_data %>%
unnest(baseline.var) %>% na.omit(cols=baseline.var)
}
else regression_unnest <- regression_data %>% na.omit(cols=baseline.var)
if(var.types[i] != "NA" & var.types[i] != "exam.score_study.time"){
regression_unnest[baseline.var] = lapply(regression_unnest[baseline.var], factor)
}
# model.var <- list() # grouped by a pack of AEQ subscales
for(AEQ in AEQs){ # reformulate(c("`A var`", "B"), "Y") AEQ ~ group + baseline.var
if(var.types[i] == "NA"){
formula <- as.formula(paste(AEQ,"s(group, bs = 'fs')", sep = " ~ "))
title <- sprintf("Response variable: %s\nPredictor Variable: group",AEQ)
}
else if(var.types[i] == "exam.score_study.time"){ # numeric
formula <- as.formula(
paste(AEQ,
paste("s(group, bs = 'fs')",paste0("s(",baseline.var,")"), sep = " + "), # use default setting
sep = " ~ "))
title <- sprintf("Response variable: %s\nPredictor Variables: group and %s\n",AEQ, baseline.var)
}
else{ # categorical: nominal and ordinal
formula <- as.formula(
paste(AEQ,
paste("s(group, bs = 'fs')",paste0("s(",baseline.var,", bs = 'fs')"), sep = " + "),
sep = " ~ "))
title <- sprintf("Response variable: %s\nPredictor Variables: group and %s\n",AEQ, baseline.var)
}
model <- mgcv::gam(formula, data = regression_unnest %>% na.omit(cols="AEQ"))
model.sum <- summary(model)
coefficient <- round(model$coefficients,2) # estimate, effect, slope of each value of smoothed variable
# p.coeff.smooth <- round(model$p.coeff,3) # p.value of each value of smoothed variable
sig.smooth.table <- model.sum$s.table # Approximate significance of smooth terms; combined effect of a variable (use F statistics, related to ANOVA)
# print(model.sum$p.table)
isSig <- ifelse(model.sum$s.table[,4] <= 0.05,T,F)
if(TRUE %in% isSig){
report <- append(report, paste0(strrep("*",times=60),"\n",title,"Approximate significance of smooth terms:\n"))
report <- append(report, paste0(capture.output(print(sig.smooth.table)),"\n"))
report <- append(report, "Coefficients:\n")
report <- append(report, paste0(capture.output(print(coefficient)),"\n"))
report <- append(report, "\n\n")
}
models.df[nrow(models.df) + 1,] = list(var.types[i], baseline.var, AEQ, list(coefficient), list(sig.smooth.table))
}
# names(model.group)[length(model.group)] <- baseline.var
}
i = i + 1
# names(models)[length(models)] <- baseline.var.type
}
report.regression <- report
cat(report.regression,sep = "")## ************************************************************
## Response variable: AEQ.6.anxiety
## Predictor Variables: group and gender
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.963646e-11 3 4.080340e-12 0.536745919
## s(gender) 8.924445e-01 1 8.297529e+00 0.008117999
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(gender).1 s(gender).2
## 15.67 0.00 0.00 0.00 3.59 -3.59
##
##
## ************************************************************
## Response variable: AEQ.6.hopeless
## Predictor Variables: group and gender
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.537428e-11 3 1.575858e-12 0.72813180
## s(gender) 9.206092e-01 1 1.159592e+01 0.00291672
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(gender).1 s(gender).2
## 11.23 0.00 0.00 0.00 3.08 -3.08
##
##
## ************************************************************
## Response variable: AEQ.7.shame
## Predictor Variables: group and ethnicity
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.005383 2 1.275182 0.18833951
## s(ethnicity) 1.760722 2 6.137339 0.01274067
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(ethnicity).1
## 9.80 -0.45 -0.66 1.11 4.83
## s(ethnicity).2 s(ethnicity).3
## -1.67 -3.15
##
##
## ************************************************************
## Response variable: AEQ.6.anxiety
## Predictor Variables: group and meditation.times
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.7103534 2 0.7138021 0.28588164
## s(meditation.times) 1.7697671 2 5.4975622 0.01753378
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 18.03 -0.68 -0.51
## s(group).3 s(meditation.times).1 s(meditation.times).2
## 1.20 7.53 -2.61
## s(meditation.times).3
## -4.92
##
##
## ************************************************************
## Response variable: AEQ.14.pride
## Predictor Variables: group and meditation.times
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.125754e-10 3 3.219509e-11 0.43984921
## s(meditation.times) 1.233176e+00 2 2.478647e+00 0.04073988
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 12.99 0.00 0.00
## s(group).3 s(meditation.times).1 s(meditation.times).2
## 0.00 1.14 -2.08
## s(meditation.times).3
## 0.94
##
##
## ************************************************************
## Response variable: AEQ.14.anger
## Predictor Variables: group and concern.2
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.411756e-11 3 1.972596e-12 0.66472807
## s(concern.2) 1.655418e+00 2 3.657819e+00 0.03255441
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(concern.2).1
## 7.11 0.00 0.00 0.00 -2.05
## s(concern.2).2 s(concern.2).3
## 0.01 2.04
##
##
## ************************************************************
## Response variable: AEQ.6.hope
## Predictor Variables: group and concern.4
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 3.276738e-11 3 4.792431e-12 0.64941233
## s(concern.4) 1.514079e+00 2 3.218125e+00 0.03613567
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(concern.4).1
## 11.03 0.00 0.00 0.00 0.75
## s(concern.4).2 s(concern.4).3
## 2.94 -3.69
##
##
## ************************************************************
## Response variable: AEQ.14.anger
## Predictor Variables: group and concern.6
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.1649865 2 0.1127946 0.35435630
## s(concern.6) 1.7128591 2 4.1515520 0.02871288
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(concern.6).1
## 8.50 -0.10 0.11 -0.01 3.07
## s(concern.6).2 s(concern.6).3
## -0.33 -2.74
##
##
## ************************************************************
## Response variable: AEQ.6.hope
## Predictor Variables: group and religious.level
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.873650e-10 3 5.559132e-11 0.42878503
## s(religious.level) 1.156671e+00 3 1.355876e+00 0.04686367
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 11.73 0.00 0.00
## s(group).3 s(religious.level).1 s(religious.level).2
## 0.00 -2.73 0.92
## s(religious.level).3 s(religious.level).4
## 1.97 -0.16
##
##
## ************************************************************
## Response variable: AEQ.6.joy
## Predictor Variables: group and religious.level
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.7090335 2 0.4981796 0.29359816
## s(religious.level) 1.6078022 3 2.3868580 0.02545027
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 5.92 -0.36 0.03
## s(group).3 s(religious.level).1 s(religious.level).2
## 0.33 -1.15 -0.70
## s(religious.level).3 s(religious.level).4
## 1.34 0.52
##
##
## ************************************************************
## Response variable: AEQ.7.pride
## Predictor Variables: group and religious.level
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.775846 2 7.956171 0.01848977
## s(religious.level) 2.534293 3 7.525949 0.01462344
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 11.32 -0.53 2.55
## s(group).3 s(religious.level).1 s(religious.level).2
## -2.02 -1.11 -3.74
## s(religious.level).3 s(religious.level).4
## 2.38 2.46
##
##
## ************************************************************
## Response variable: AEQ.14.pride
## Predictor Variables: group and religious.level
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 3.541847e-11 3 8.767240e-12 0.49166972
## s(religious.level) 2.062603e+00 3 2.368111e+00 0.04729275
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 13.92 0.00 0.00
## s(group).3 s(religious.level).1 s(religious.level).2
## 0.00 -3.04 3.23
## s(religious.level).3 s(religious.level).4
## -0.76 0.57
##
##
## ************************************************************
## Response variable: AEQ.difference.joy
## Predictor Variables: group and religious.level
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.518268e-09 2 8.706893e-10 0.34431256
## s(religious.level) 2.196038e+00 3 3.393932e+00 0.01943028
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 1.55 0.00 0.00
## s(group).3 s(religious.level).1 s(religious.level).2
## 0.00 -1.72 -0.91
## s(religious.level).3 s(religious.level).4
## 0.60 2.03
##
##
## ************************************************************
## Response variable: AEQ.difference.pride
## Predictor Variables: group and religious.level
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 4.238554e-11 3 1.019828e-11 0.49939220
## s(religious.level) 1.318848e+00 3 1.636651e+00 0.03936183
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 0.98 0.00 0.00
## s(group).3 s(religious.level).1 s(religious.level).2
## 0.00 -1.65 0.53
## s(religious.level).3 s(religious.level).4
## 0.60 0.53
##
##
## ************************************************************
## Response variable: AEQ.difference.relief
## Predictor Variables: group and religious.level
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.1212395 2 0.07899692 0.31850355
## s(religious.level) 2.5952509 3 3.21635525 0.04057952
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 0.91 -0.03 0.05
## s(group).3 s(religious.level).1 s(religious.level).2
## -0.02 -1.07 3.16
## s(religious.level).3 s(religious.level).4
## -0.61 -1.47
##
##
## ************************************************************
## Response variable: AEQ.7.relief
## Predictor Variables: group and relax
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 2.517715e-11 3 3.146659e-12 0.69179466
## s(relax) 2.654209e+00 3 3.834007e+00 0.02537426
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(relax).1 s(relax).2
## 12.91 0.00 0.00 0.00 -7.18 2.99
## s(relax).3 s(relax).4
## 0.09 4.10
##
##
## ************************************************************
## Response variable: AEQ.6.hopeless
## Predictor Variables: group and improve.score
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.441772 2 2.850944 0.08649617
## s(improve.score) 1.555476 2 5.992362 0.01461130
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 12.59 -0.03 -2.17 2.20
## s(improve.score).1 s(improve.score).2 s(improve.score).3
## -2.74 -2.10 4.84
##
##
## ************************************************************
## Response variable: AEQ.difference.joy
## Predictor Variables: group and improve.score
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.0749976 2 0.0444193 0.33873609
## s(improve.score) 1.6180180 2 3.7230578 0.02996132
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 0.71 0.04 -0.03 -0.01
## s(improve.score).1 s(improve.score).2 s(improve.score).3
## -1.75 -0.46 2.21
##
##
## ************************************************************
## Response variable: AEQ.difference.joy
## Predictor Variables: group and race
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.929364 2 26.00661 0.01423183
## s(race) 1.694602 3 67.10021 0.02271690
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(race).1 s(race).2
## 0.17 2.03 -1.66 -0.38 1.36 0.47
## s(race).3 s(race).4
## -0.33 -1.50
##
##
## ************************************************************
## Response variable: AEQ.7.pride
## Predictor Variables: group and major
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.536758e+00 2 2.838733e+00 0.04536525
## s(major) 4.079446e-10 9 4.678239e-11 0.41872698
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(major).1 s(major).2
## 11.47 0.03 1.82 -1.84 0.00 0.00
## s(major).3 s(major).4 s(major).5 s(major).6 s(major).7 s(major).8
## 0.00 0.00 0.00 0.00 0.00 0.00
## s(major).9 s(major).10
## 0.00 0.00
##
##
## ************************************************************
## Response variable: AEQ.6.joy
## Predictor Variables: group and religion
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.6216313 2 0.4690953 0.29604160
## s(religion) 2.8145536 4 3.3254006 0.01670994
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(religion).1
## 5.70 -0.31 0.07 0.23 2.44
## s(religion).2 s(religion).3 s(religion).4 s(religion).5
## -2.28 0.37 -0.89 0.37
##
##
## ************************************************************
## Response variable: AEQ.7.pride
## Predictor Variables: group and religion
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.554670 2 3.457742 0.045157181
## s(religion) 1.809944 4 2.818659 0.008763601
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(religion).1
## 12.29 -0.06 1.65 -1.59 1.43
## s(religion).2 s(religion).3 s(religion).4 s(religion).5
## 0.69 -0.08 -2.29 0.26
##
##
## ************************************************************
## Response variable: AEQ.difference.anxiety
## Predictor Variables: group and religion
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.1250956 2 0.08206911 0.32255156
## s(religion) 3.3844217 4 3.63593183 0.02266672
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(religion).1
## -3.40 0.09 -0.07 -0.03 3.11
## s(religion).2 s(religion).3 s(religion).4 s(religion).5
## -6.70 -1.33 1.54 3.38
##
##
## ************************************************************
## Response variable: AEQ.difference.joy
## Predictor Variables: group and religion
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.107924 2 1.656441 0.105493556
## s(religion) 2.079555 4 3.294438 0.006621289
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(religion).1
## 1.47 0.65 -0.43 -0.22 1.02
## s(religion).2 s(religion).3 s(religion).4 s(religion).5
## 0.39 -0.45 -1.71 0.75
##
##
## ************************************************************
## Response variable: AEQ.difference.anxiety
## Predictor Variables: group and religion.family
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 2.656787e-11 3 6.852469e-12 0.47140680
## s(religion.family) 3.473816e+00 4 4.089946e+00 0.01295566
## Coefficients:
## (Intercept) s(group).1 s(group).2
## -3.94 0.00 0.00
## s(group).3 s(religion.family).1 s(religion.family).2
## 0.00 3.51 -6.27
## s(religion.family).3 s(religion.family).4 s(religion.family).5
## -0.83 1.45 2.13
##
##
## ************************************************************
## Response variable: AEQ.difference.shame
## Predictor Variables: group and religion.family
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 2.333945e-11 3 4.244912e-12 0.58225320
## s(religion.family) 3.345417e+00 4 2.830841e+00 0.03927289
## Coefficients:
## (Intercept) s(group).1 s(group).2
## -3.88 0.00 0.00
## s(group).3 s(religion.family).1 s(religion.family).2
## 0.00 3.27 -7.37
## s(religion.family).3 s(religion.family).4 s(religion.family).5
## -1.54 1.10 4.53
##
##
## ************************************************************
## Response variable: AEQ.6.anxiety
## Predictor Variables: group and meditation.type
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.147921 2 1.534574 0.14022336
## s(meditation.type) 2.992005 4 2.537193 0.04021984
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 17.55 -0.90 -1.04
## s(group).3 s(meditation.type).1 s(meditation.type).2
## 1.94 0.29 -2.14
## s(meditation.type).3 s(meditation.type).4 s(meditation.type).5
## -2.99 6.04 -1.20
##
##
## ************************************************************
## Response variable: AEQ.7.shame
## Predictor Variables: group and meditation.type
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.709489e+00 2 4.155214e+00 0.01905394
## s(meditation.type) 3.977563e-11 5 5.215105e-12 0.60223393
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 9.42 -1.03 -2.19
## s(group).3 s(meditation.type).1 s(meditation.type).2
## 3.22 0.00 0.00
## s(meditation.type).3 s(meditation.type).4 s(meditation.type).5
## 0.00 0.00 0.00
##
##
## ************************************************************
## Response variable: AEQ.14.anger
## Predictor Variables: group and meditation.type
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.9945681 2 1.403846 0.11275876
## s(meditation.type) 2.7329805 4 2.582344 0.02595243
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 7.60 -0.73 0.27
## s(group).3 s(meditation.type).1 s(meditation.type).2
## 0.46 0.70 -0.57
## s(meditation.type).3 s(meditation.type).4 s(meditation.type).5
## -0.95 -1.38 2.21
##
##
## ************************************************************
## Response variable: AEQ.14.shame
## Predictor Variables: group and meditation.type
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.625554 2 3.955241 0.06164407
## s(meditation.type) 3.298442 4 4.229381 0.01739870
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 8.22 -1.02 -0.98
## s(group).3 s(meditation.type).1 s(meditation.type).2
## 2.01 1.46 -0.17
## s(meditation.type).3 s(meditation.type).4 s(meditation.type).5
## -1.23 -3.50 3.44
##
##
## ************************************************************
## Response variable: AEQ.difference.hope
## Predictor Variables: group and meditation.type
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.570956e+00 2 5.132338e+00 0.006582657
## s(meditation.type) 1.008619e-10 5 1.690587e-11 0.492909946
## Coefficients:
## (Intercept) s(group).1 s(group).2
## -0.07 2.15 -1.18
## s(group).3 s(meditation.type).1 s(meditation.type).2
## -0.97 0.00 0.00
## s(meditation.type).3 s(meditation.type).4 s(meditation.type).5
## 0.00 0.00 0.00
##
##
## ************************************************************
## Response variable: AEQ.difference.joy
## Predictor Variables: group and meditation.type
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.355606e+00 2 2.828541e+00 0.02992754
## s(meditation.type) 8.357097e-11 5 9.980704e-12 0.64139394
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 0.55 1.11 -0.66
## s(group).3 s(meditation.type).1 s(meditation.type).2
## -0.45 0.00 0.00
## s(meditation.type).3 s(meditation.type).4 s(meditation.type).5
## 0.00 0.00 0.00
##
##
## ************************************************************
## Response variable: AEQ.6.anxiety
## Predictor Variables: group and credits
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.617952e-11 3 3.252343e-12 0.55178696
## s(credits) 5.145533e+00 6 3.817766e+00 0.01588944
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(credits).1 s(credits).2
## 17.97 0.00 0.00 0.00 -1.57 15.17
## s(credits).3 s(credits).4 s(credits).5 s(credits).6 s(credits).7
## 0.03 -5.14 -0.89 -3.64 -3.96
##
##
## ************************************************************
## Response variable: AEQ.6.hope
## Predictor Variables: group and credits
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.7607605 2 1.024785 0.40204462
## s(credits) 5.1511350 6 3.302533 0.04885178
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(credits).1 s(credits).2
## 11.76 -0.40 1.10 -0.70 -0.30 13.73
## s(credits).3 s(credits).4 s(credits).5 s(credits).6 s(credits).7
## -5.16 -0.34 -3.39 -3.03 -1.51
##
##
## ************************************************************
## Response variable: AEQ.6.hopeless
## Predictor Variables: group and credits
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.229261 2 1.696196 0.15505283
## s(credits) 2.236800 6 1.416187 0.03488754
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(credits).1 s(credits).2
## 10.87 0.43 -1.83 1.40 -1.13 -0.07
## s(credits).3 s(credits).4 s(credits).5 s(credits).6 s(credits).7
## 0.18 -1.42 -0.79 3.49 -0.27
##
##
## ************************************************************
## Response variable: AEQ.6.joy
## Predictor Variables: group and credits
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.607523e-11 3 5.773268e-14 0.98863051
## s(credits) 5.081982e+00 6 3.050695e+00 0.03182184
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(credits).1 s(credits).2
## 5.96 0.00 0.00 0.00 0.03 4.73
## s(credits).3 s(credits).4 s(credits).5 s(credits).6 s(credits).7
## -3.10 -0.66 -0.57 -1.24 0.81
##
##
## ************************************************************
## Response variable: AEQ.14.pride
## Predictor Variables: group and age
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.220416 2 1.475823 0.174047799
## s(age) 1.444133 2 6.381231 0.007180118
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(age).1 s(age).2
## 11.61 0.12 1.02 -1.14 -2.06 2.42
## s(age).3
## -0.36
##
##
## ************************************************************
## Response variable: AEQ.14.relief
## Predictor Variables: group and age
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.748432e-11 3 1.894520e-12 0.72515335
## s(age) 1.103376e+00 2 2.112201e+00 0.04666311
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(age).1 s(age).2
## 14.33 0.00 0.00 0.00 -1.54 2.02
## s(age).3
## -0.48
##
##
## ************************************************************
## Response variable: AEQ.7.pride
## Predictor Variables: group and GPA.previous
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 2.00002 2 2.190591e+02 0.01553894
## s(GPA.previous) 10.46619 12 8.420023e+11 0.02479328
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 11.10 -1.37 3.43 -2.06
## s(GPA.previous).1 s(GPA.previous).2 s(GPA.previous).3 s(GPA.previous).4
## -0.04 3.75 -1.52 2.69
## s(GPA.previous).5 s(GPA.previous).6 s(GPA.previous).7 s(GPA.previous).8
## 2.91 -4.01 1.29 -0.04
## s(GPA.previous).9 s(GPA.previous).10 s(GPA.previous).11 s(GPA.previous).12
## 3.75 -2.22 -1.06 -0.64
## s(GPA.previous).13
## -4.86
##
##
## ************************************************************
## Response variable: AEQ.7.shame
## Predictor Variables: group and GPA.previous
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.7393867 2 2.069981 0.22943482
## s(GPA.previous) 9.1638228 12 3.970002 0.03952446
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 8.49 0.15 -0.72 0.58
## s(GPA.previous).1 s(GPA.previous).2 s(GPA.previous).3 s(GPA.previous).4
## 6.65 -1.96 3.99 -0.60
## s(GPA.previous).5 s(GPA.previous).6 s(GPA.previous).7 s(GPA.previous).8
## 1.41 3.49 -2.80 -1.54
## s(GPA.previous).9 s(GPA.previous).10 s(GPA.previous).11 s(GPA.previous).12
## -1.96 0.18 -3.58 -3.45
## s(GPA.previous).13
## 0.18
##
##
## ************************************************************
## Response variable: AEQ.difference.anger
## Predictor Variables: group and GPA.previous
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 2.879492e-12 3 1.942965e-14 0.979456033
## s(GPA.previous) 1.139887e+01 12 1.787596e+01 0.002247874
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## -2.36 0.00 0.00 0.00
## s(GPA.previous).1 s(GPA.previous).2 s(GPA.previous).3 s(GPA.previous).4
## -10.01 -0.60 4.11 2.78
## s(GPA.previous).5 s(GPA.previous).6 s(GPA.previous).7 s(GPA.previous).8
## -5.95 -4.50 -6.25 4.11
## s(GPA.previous).9 s(GPA.previous).10 s(GPA.previous).11 s(GPA.previous).12
## 0.34 5.05 0.84 3.16
## s(GPA.previous).13
## 6.93
##
##
## ************************************************************
## Response variable: AEQ.difference.pride
## Predictor Variables: group and GPA.previous
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 2.185669e-11 3 5.553739e-12 0.49026982
## s(GPA.previous) 7.284372e+00 12 1.893067e+00 0.04810839
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 0.17 0.00 0.00 0.00
## s(GPA.previous).1 s(GPA.previous).2 s(GPA.previous).3 s(GPA.previous).4
## -0.09 1.58 -1.21 3.11
## s(GPA.previous).5 s(GPA.previous).6 s(GPA.previous).7 s(GPA.previous).8
## -0.48 -2.27 -0.09 -0.09
## s(GPA.previous).9 s(GPA.previous).10 s(GPA.previous).11 s(GPA.previous).12
## 2.14 -1.77 -1.20 -1.21
## s(GPA.previous).13
## 1.58
##
##
## ************************************************************
## Response variable: AEQ.6.anxiety
## Predictor Variables: group and grade.predict
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.430244 2 75.48078 0.090016847
## s(grade.predict) 8.548368 9 29.91289 0.004968139
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 16.29 1.39 -2.28 0.88
## s(grade.predict).1 s(grade.predict).2 s(grade.predict).3 s(grade.predict).4
## -1.88 0.32 5.66 -2.32
## s(grade.predict).5 s(grade.predict).6 s(grade.predict).7 s(grade.predict).8
## 1.60 18.74 -8.21 -0.96
## s(grade.predict).9 s(grade.predict).10
## -7.26 -5.69
##
##
## ************************************************************
## Response variable: AEQ.6.hope
## Predictor Variables: group and grade.predict
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.3064873 2 0.2416603 0.3962420
## s(grade.predict) 6.8992217 9 2.5714508 0.0498696
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 11.36 -0.33 0.34 -0.01
## s(grade.predict).1 s(grade.predict).2 s(grade.predict).3 s(grade.predict).4
## -2.65 -5.21 -5.48 -0.33
## s(grade.predict).5 s(grade.predict).6 s(grade.predict).7 s(grade.predict).8
## -0.72 11.86 -0.02 -0.50
## s(grade.predict).9 s(grade.predict).10
## -0.02 3.06
##
##
## ************************************************************
## Response variable: AEQ.6.joy
## Predictor Variables: group and grade.predict
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 7.464011e-11 3 2.162902e-11 0.43593760
## s(grade.predict) 7.154198e+00 9 2.842487e+00 0.03579758
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 5.75 0.00 0.00 0.00
## s(grade.predict).1 s(grade.predict).2 s(grade.predict).3 s(grade.predict).4
## -0.22 -1.87 -2.78 -0.69
## s(grade.predict).5 s(grade.predict).6 s(grade.predict).7 s(grade.predict).8
## 0.21 4.62 0.92 0.92
## s(grade.predict).9 s(grade.predict).10
## -1.30 0.18
##
##
## ************************************************************
## Response variable: AEQ.difference.hopeless
## Predictor Variables: group and grade.predict
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.3729402 2 0.4882811 0.33759603
## s(grade.predict) 7.5653496 9 3.9959974 0.02444592
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## -0.90 0.32 -0.19 -0.13
## s(grade.predict).1 s(grade.predict).2 s(grade.predict).3 s(grade.predict).4
## 1.60 -0.42 -7.95 -0.15
## s(grade.predict).5 s(grade.predict).6 s(grade.predict).7 s(grade.predict).8
## 3.42 -3.99 -1.14 3.29
## s(grade.predict).9 s(grade.predict).10
## 3.67 1.68
##
##
## ************************************************************
## Response variable: AEQ.6.anxiety
## Predictor Variables: group and productivity6
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 8.774577e-10 3 3.207978e-10 0.3615621055
## s(productivity6) 4.785494e+00 5 1.035719e+01 0.0005892683
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 17.06 0.00 0.00 0.00
## s(productivity6).1 s(productivity6).2 s(productivity6).3 s(productivity6).4
## -7.53 -0.05 0.44 18.64
## s(productivity6).5 s(productivity6).6
## -3.03 -8.46
##
##
## ************************************************************
## Response variable: AEQ.6.hope
## Predictor Variables: group and productivity6
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.300778 2 7.172321 0.137565173
## s(productivity6) 4.653697 5 8.370785 0.006635131
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 12.90 -0.01 1.66 -1.65
## s(productivity6).1 s(productivity6).2 s(productivity6).3 s(productivity6).4
## -1.70 -6.24 -5.33 15.15
## s(productivity6).5 s(productivity6).6
## -3.19 1.31
##
##
## ************************************************************
## Response variable: AEQ.6.joy
## Predictor Variables: group and productivity6
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.968938e-11 3 2.702484e-12 0.66371958
## s(productivity6) 4.384781e+00 5 4.173470e+00 0.01259658
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 5.96 0.00 0.00 0.00
## s(productivity6).1 s(productivity6).2 s(productivity6).3 s(productivity6).4
## 0.85 -2.22 -1.62 4.93
## s(productivity6).5 s(productivity6).6
## -0.33 -1.60
##
##
## ************************************************************
## Response variable: AEQ.14.anger
## Predictor Variables: group and productivity6
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 0.6524268 2 0.5417918 0.25320991
## s(productivity6) 2.1429294 5 1.5496847 0.03465215
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 6.43 -0.49 0.27 0.22
## s(productivity6).1 s(productivity6).2 s(productivity6).3 s(productivity6).4
## -0.27 -1.04 0.11 -0.19
## s(productivity6).5 s(productivity6).6
## 2.18 -0.79
##
##
## ************************************************************
## Response variable: AEQ.14.pride
## Predictor Variables: group and productivity6
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.967035 2 53.58624 0.002265159
## s(productivity6) 4.359018 5 85.43736 0.004438901
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 13.71 1.43 2.26 -3.69
## s(productivity6).1 s(productivity6).2 s(productivity6).3 s(productivity6).4
## -0.11 0.39 -5.53 4.99
## s(productivity6).5 s(productivity6).6
## -3.09 3.36
##
##
## ************************************************************
## Response variable: AEQ.6.hopeless
## Predictor Variables: group and productivity7
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 9.758472e-12 3 1.105974e-12 0.71331338
## s(productivity7) 2.502778e+00 3 2.872017e+00 0.04916358
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 12.30 0.00 0.00 0.00
## s(productivity7).1 s(productivity7).2 s(productivity7).3 s(productivity7).4
## -4.88 4.97 -1.19 1.10
##
##
## ************************************************************
## Response variable: AEQ.7.anger
## Predictor Variables: group and productivity7
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 5.245249e-11 3 1.384312e-11 0.46802089
## s(productivity7) 2.793599e+00 3 4.790408e+00 0.01461879
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3
## 8.41 0.00 0.00 0.00
## s(productivity7).1 s(productivity7).2 s(productivity7).3 s(productivity7).4
## -0.15 3.85 -1.53 -2.17
##
##
## ************************************************************
## Response variable: AEQ.7.anger
## Predictor Variables: group and exam.score
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 2.629128e-11 3.00000 2.534469e-12 0.73127260
## s(exam.score) 2.076956e+00 2.58593 4.043894e+00 0.03767037
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(exam.score).1
## 7.29 0.00 0.00 0.00 -0.41
## s(exam.score).2 s(exam.score).3 s(exam.score).4 s(exam.score).5 s(exam.score).6
## 0.45 -0.13 0.52 -0.20 0.50
## s(exam.score).7 s(exam.score).8 s(exam.score).9
## 0.00 -3.02 0.95
##
##
## ************************************************************
## Response variable: AEQ.difference.anger
## Predictor Variables: group and exam.score
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.232233 2 1.827157 0.084401806
## s(exam.score) 1.000000 1 10.773325 0.005454788
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(exam.score).1
## -3.07 2.01 -0.45 -1.56 0.00
## s(exam.score).2 s(exam.score).3 s(exam.score).4 s(exam.score).5 s(exam.score).6
## 0.00 0.00 0.00 0.00 0.00
## s(exam.score).7 s(exam.score).8 s(exam.score).9
## 0.00 0.00 3.37
##
##
## ************************************************************
## Response variable: AEQ.difference.relief
## Predictor Variables: group and exam.score
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.308655 2.000000 3.462043 0.04608096
## s(exam.score) 8.427342 8.897769 4.165225 0.05169583
## Coefficients:
## (Intercept) s(group).1 s(group).2 s(group).3 s(exam.score).1
## 0.21 -0.86 0.99 -0.13 -7.33
## s(exam.score).2 s(exam.score).3 s(exam.score).4 s(exam.score).5 s(exam.score).6
## 4.97 1.09 2.24 2.11 0.05
## s(exam.score).7 s(exam.score).8 s(exam.score).9
## -1.41 -4.74 11.03
##
##
## ************************************************************
## Response variable: AEQ.7.pride
## Predictor Variables: group and study.time.total
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.700853 2.000000 4.727285 0.01715818
## s(study.time.total) 1.276414 1.494543 5.235040 0.02018309
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 10.98 0.39 2.05
## s(group).3 s(study.time.total).1 s(study.time.total).2
## -2.43 -0.03 0.21
## s(study.time.total).3 s(study.time.total).4 s(study.time.total).5
## 0.06 0.15 0.06
## s(study.time.total).6 s(study.time.total).7 s(study.time.total).8
## -0.17 -0.12 0.77
## s(study.time.total).9
## 1.72
##
##
## ************************************************************
## Response variable: AEQ.7.shame
## Predictor Variables: group and study.time.total
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.690449 2 3.970114 0.02918300
## s(study.time.total) 1.000000 1 4.388741 0.05633225
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 8.93 -1.13 -2.17
## s(group).3 s(study.time.total).1 s(study.time.total).2
## 3.30 0.00 0.00
## s(study.time.total).3 s(study.time.total).4 s(study.time.total).5
## 0.00 0.00 0.00
## s(study.time.total).6 s(study.time.total).7 s(study.time.total).8
## 0.00 0.00 0.00
## s(study.time.total).9
## -1.75
##
##
## ************************************************************
## Response variable: AEQ.14.anger
## Predictor Variables: group and study.time.total
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.893573 2.000000 5.157426 0.04035853
## s(study.time.total) 8.035068 8.737403 4.698256 0.03564614
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 7.66 -2.26 -0.24
## s(group).3 s(study.time.total).1 s(study.time.total).2
## 2.50 11.81 -42.79
## s(study.time.total).3 s(study.time.total).4 s(study.time.total).5
## 1.50 -26.48 -3.28
## s(study.time.total).6 s(study.time.total).7 s(study.time.total).8
## 24.32 14.47 -61.53
## s(study.time.total).9
## 20.60
##
##
## ************************************************************
## Response variable: AEQ.14.shame
## Predictor Variables: group and study.time.total
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.979829 2.000000 4.001508 0.04857260
## s(study.time.total) 6.214019 7.304922 2.748287 0.08876767
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 8.09 -2.51 -1.51
## s(group).3 s(study.time.total).1 s(study.time.total).2
## 4.03 12.89 -36.87
## s(study.time.total).3 s(study.time.total).4 s(study.time.total).5
## 6.73 -22.62 -4.20
## s(study.time.total).6 s(study.time.total).7 s(study.time.total).8
## 15.71 12.09 -51.19
## s(study.time.total).9
## 23.61
##
##
## ************************************************************
## Response variable: AEQ.difference.hope
## Predictor Variables: group and study.time.total
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.486005 2 3.011461 0.04093086
## s(study.time.total) 1.000000 1 3.736281 0.07373060
## Coefficients:
## (Intercept) s(group).1 s(group).2
## -0.46 1.85 -0.41
## s(group).3 s(study.time.total).1 s(study.time.total).2
## -1.44 0.00 0.00
## s(study.time.total).3 s(study.time.total).4 s(study.time.total).5
## 0.00 0.00 0.00
## s(study.time.total).6 s(study.time.total).7 s(study.time.total).8
## 0.00 0.00 0.00
## s(study.time.total).9
## 1.26
##
##
## ************************************************************
## Response variable: AEQ.difference.joy
## Predictor Variables: group and study.time.total
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(group) 1.380504 2 2.477424 0.055260845
## s(study.time.total) 1.000000 1 10.330455 0.006246406
## Coefficients:
## (Intercept) s(group).1 s(group).2
## 0.51 1.03 -0.35
## s(group).3 s(study.time.total).1 s(study.time.total).2
## -0.68 0.00 0.00
## s(study.time.total).3 s(study.time.total).4 s(study.time.total).5
## 0.00 0.00 0.00
## s(study.time.total).6 s(study.time.total).7 s(study.time.total).8
## 0.00 0.00 0.00
## s(study.time.total).9
## 1.30report <- NULL
datatable(
models.df, filter = 'top', extensions = 'Buttons', options = list(
dom = 'Bfrtip',
buttons = c('copy', 'csv', 'excel', 'pdf', 'print')
)
)