TTC analysis
07 June, 2022
library(tidyverse)
library(zoo)
library(knitr)
library(tibble)
library(plyr)
library(dplyr)
library(data.table)
library(Hmisc)
#library(grt)
library(reshape) #reshape data from long to wide
library(MASS) #backward elimination of non-significant predictors in linear models
#library(ggcorrplot)
#library(car) #Anova
library(emmeans) #multiple comparisons
library(Rmisc)
library(lubridate) #time stamp conversion
library(RcppRoll)
library(ltm)
library(kableExtra)Qualtrics Data
setwd('/home/tzhu/Desktop/new_TTC_data')
df <- read.csv("new_ATQ.csv", header=T)Cleaning - Extract Demographic and ATQ items
df[13:74]<-df[13:74] %>% mutate_all(~ifelse(is.na(.x), mean(.x, na.rm = TRUE), .x)) #fill na with col average
cols.num <- c(13:74) #find ATQ columns
df[cols.num] <- sapply(df[cols.num],as.numeric) #convert to numerics
colnames(df)[which(names(df) == "PA.3_1")] <- "EX.3_1"
colnames(df)[which(names(df) == "EX.40_1")] <- "EC.40_1"
#sapply(df, class)Reverse coding
indx <- grepl('R', colnames(df)) #find reverse coded items
df[,indx]<-abs(df[,indx]-8) #reverse coding
# df<-df[,c(1:12,15,26,35,53,69,21,29,32,40,49,59,66,17,36,41,45,54)] #when looking at subfactor within each trait factorGet average factor scale scores
NA_indx <- grepl('NA', colnames(df)) #find Neg_Affect items
NA_sum <- rowSums(df[, NA_indx]) #total score of NA items for each participant
# df$NA_score <- NA_sum/25 #append average NA scores to df
df$NA_score <- NA_sum/25
cronbach.alpha(df[,NA_indx])##
## Cronbach's alpha for the 'df[, NA_indx]' data-set
##
## Items: 26
## Sample units: 353
## alpha: 0.804EC_indx <- grepl('EC', colnames(df)) #find EC items
EC_sum <- rowSums(df[, EC_indx]) #total score of EC items for each participant
# df$EC_score <- EC_sum/19 #append EC scores to df
df$EC_score <- EC_sum/19
cronbach.alpha(df[,EC_indx])##
## Cronbach's alpha for the 'df[, EC_indx]' data-set
##
## Items: 17
## Sample units: 353
## alpha: 0.815EX_indx <- grepl('EX', colnames(df)) #find EX items
EX_sum <- rowSums(df[, EX_indx]) #total score of EX items for each participant
# df$EX_score <- EX_sum/17 #append EX scores to df
df$EX_score <- EX_sum/17
cronbach.alpha(df[,EX_indx])##
## Cronbach's alpha for the 'df[, EX_indx]' data-set
##
## Items: 19
## Sample units: 353
## alpha: 0.749df_temp <-df[,c(1:12, 75:77)]
# df_temp <-df[,c(1:12, 30:32)]
#generate columns that take temperament factors as low/medium/high levels for potential analyses
df_temp$EC_level <- as.factor(ifelse(df_temp$EC_score<=median(df_temp$EC_score),"low","high"))
df_temp$NA_level <- as.factor(ifelse(df_temp$NA_score<=median(df_temp$NA_score),"low","high"))
df_temp$EX_level <- as.factor(ifelse(df_temp$EX_score<=median(df_temp$EX_score),"low","high"))
df_temp$UniqueID<-toupper(df_temp$UniqueID)
# df_temp<-df_temp[!df_temp$UniqueID == "8803SA",] #8803SA only exists in temperament but not cat learn
getwd()## [1] "/home/tzhu/Desktop/new_TTC_data"Pavlovia Data
Get CR task data
because of changes in experimental design, data columns has been changed, thus outputs with same columns are put in a folder (batch)
# CR task from CR first condition
setwd("/home/tzhu/Desktop/new_TTC_data/CR_first/first_batch")
csv <- list.files(full.names = TRUE, pattern = '*.csv')
df_Pav <- as_tibble(rbind.fill(lapply (csv, fread)))
cols1 <- c(1, 8:11, 14, 17:18, 21, 23:25)
df_cat1<-df_Pav[,cols1]
df_cat1$UniqueID <- toupper(df_cat1$UniqueID)
names(df_cat1)[names(df_cat1) == 'CR_loop.thisN'] <- 'trial_num'
df_cat1$trial_num<-df_cat1$trial_num + 1
#remove empty rows
df_cat1<-df_cat1[!is.na(df_cat1$key_resp_25.corr),]
#add column to indicate category type
df_cat1$cat_type<-'CR'
names(df_cat1)[names(df_cat1) == 'CR_cond_(DO_NOT_CHANGE)'] <- 'CR_rule_cond'
names(df_cat1)[names(df_cat1) == 'key_resp_25.corr'] <- 'corr'
names(df_cat1)[names(df_cat1) == 'key_resp_25.keys'] <- 'key_resp'
names(df_cat1)[names(df_cat1) == 'key_resp_25.rt'] <- 'rt'
names(df_cat1)[names(df_cat1) == 'corr_count1'] <- 'corr_count'
df_cat1$cat_order <- 'CR_first'
df_cat1<-df_cat1[, c(1,9,13,10,14,8,2,3,4,5,6,7)]
df_cat1$UniqueID[df_cat1$UniqueID==' 2503ND'] <- '2503ND'setwd("/home/tzhu/Desktop/new_TTC_data/CR_first/second_batch")
csv_second <- list.files(full.names = TRUE, pattern = '*.csv')
df_Pav_second <- as_tibble(rbind.fill(lapply (csv_second, fread)))
cols1_second <- c(1,6:9, 12, 15, 16, 19:21)
df_cat1_second<-df_Pav_second[,cols1_second]
df_cat1_second$UniqueID <- toupper(df_cat1_second$UniqueID)
names(df_cat1_second)[names(df_cat1_second) == 'CR_loop.thisN'] <- 'trial_num'
df_cat1_second$trial_num<-df_cat1_second$trial_num + 1
#remove empty rows
df_cat1_second<-df_cat1_second[!is.na(df_cat1_second$key_resp_25.corr),]
#add column to indicate category type
df_cat1_second$cat_type <- 'CR'
names(df_cat1_second)[names(df_cat1_second) == 'CR_cond_(DO_NOT_CHANGE)'] <- 'CR_rule_cond'
names(df_cat1_second)[names(df_cat1_second) == 'key_resp_25.corr'] <- 'corr'
names(df_cat1_second)[names(df_cat1_second) == 'key_resp_25.keys'] <- 'key_resp'
names(df_cat1_second)[names(df_cat1_second) == 'key_resp_25.rt'] <- 'rt'
names(df_cat1_second)[names(df_cat1_second) == 'corr_count1'] <- 'corr_count'
df_cat1_second<-df_cat1_second[, c(1,9,12,11,10,7,2:6,8)]setwd("/home/tzhu/Desktop/new_TTC_data/CR_first/third_batch")
csv_third <- list.files(full.names = TRUE, pattern = '*.csv')
df_Pav_third <- as_tibble(rbind.fill(lapply (csv_third, fread)))
cols1_third <- c(1,6:9, 12, 15,16, 19, 21,22)
df_cat1_third<-df_Pav_third[,cols1_third]
df_cat1_third$UniqueID <- toupper(df_cat1_third$UniqueID)
names(df_cat1_third)[names(df_cat1_third) == 'CR_loop.thisN'] <- 'trial_num'
df_cat1_third$trial_num<-df_cat1_third$trial_num + 1
#remove empty rows
df_cat1_third<-df_cat1_third[!is.na(df_cat1_third$key_resp_25.corr),]
#add column to indicate category type
df_cat1_third$cat_type <- 'CR'
names(df_cat1_third)[names(df_cat1_third) == 'CR_cond_(DO_NOT_CHANGE)'] <- 'CR_rule_cond'
names(df_cat1_third)[names(df_cat1_third) == 'key_resp_25.corr'] <- 'corr'
names(df_cat1_third)[names(df_cat1_third) == 'key_resp_25.keys'] <- 'key_resp'
names(df_cat1_third)[names(df_cat1_third) == 'key_resp_25.rt'] <- 'rt'
names(df_cat1_third)[names(df_cat1_third) == 'corr_count1'] <- 'corr_count'
df_cat1_third<-df_cat1_third[, c(1,9,12,11,10,7,2:6,8)]# CR task from II first condition
setwd("/home/tzhu/Desktop/new_TTC_data/II_first")
csv_fourth <- list.files(full.names = TRUE, pattern = '*.csv')
df_Pav_fourth <- as_tibble(rbind.fill(lapply (csv_fourth, fread)))
cols1_fourth <- c(1,15,16, 19, 21, 22,135:138, 141)
df_cat1_fourth<-df_Pav_fourth[,cols1_fourth]
df_cat1_fourth$UniqueID <- toupper(df_cat1_fourth$UniqueID)
names(df_cat1_fourth)[names(df_cat1_fourth) == 'CR_loop.thisN'] <- 'trial_num'
df_cat1_fourth$trial_num<-df_cat1_fourth$trial_num + 1
#remove empty rows
df_cat1_fourth<-df_cat1_fourth[!is.na(df_cat1_fourth$key_resp_26.corr),]
#add column to indicate category type
df_cat1_fourth$cat_type <- 'CR'
names(df_cat1_fourth)[names(df_cat1_fourth) == 'CR_cond_(DO_NOT_CHANGE)'] <- 'CR_rule_cond'
names(df_cat1_fourth)[names(df_cat1_fourth) == 'key_resp_26.corr'] <- 'corr'
names(df_cat1_fourth)[names(df_cat1_fourth) == 'key_resp_26.keys'] <- 'key_resp'
names(df_cat1_fourth)[names(df_cat1_fourth) == 'key_resp_26.rt'] <- 'rt'
df_cat1_fourth$cat_order <- 'II_first'
df_cat1_fourth<-df_cat1_fourth[, c(1,4,12,6,5,2,7:11,3)]cat_df_cr1 <- bind_rows(df_cat1, df_cat1_second, df_cat1_third, df_cat1_fourth)Get II data
# II task from II first condition
setwd("/home/tzhu/Desktop/new_TTC_data/II_first")
csv_ii <- list.files(full.names = TRUE, pattern = '*.csv')
df_Pav_ii <- as_tibble(rbind.fill(lapply (csv_ii, fread)))
cols1_ii <- c(1,6:9, 12, 15,16, 19, 21,23)
df_cat1_ii<-df_Pav_ii[,cols1_ii]
df_cat1_ii$UniqueID <- toupper(df_cat1_ii$UniqueID)
names(df_cat1_ii)[names(df_cat1_ii) == 'II_loop.thisN'] <- 'trial_num'
df_cat1_ii$trial_num<-df_cat1_ii$trial_num + 1
#remove empty rows
df_cat1_ii<-df_cat1_ii[!is.na(df_cat1_ii$key_resp_25.corr),]
#add column to indicate category type
df_cat1_ii$cat_type <- 'II'
names(df_cat1_ii)[names(df_cat1_ii) == 'II_cond_(DO_NOT_CHANGE)'] <- 'II_rule_cond'
names(df_cat1_ii)[names(df_cat1_ii) == 'key_resp_25.corr'] <- 'corr'
names(df_cat1_ii)[names(df_cat1_ii) == 'key_resp_25.keys'] <- 'key_resp'
names(df_cat1_ii)[names(df_cat1_ii) == 'key_resp_25.rt'] <- 'rt'
names(df_cat1_ii)[names(df_cat1_ii) == 'corr_count1'] <- 'corr_count'
df_cat1_ii$cat_order <- 'II_first'
df_cat1_ii<-df_cat1_ii[, c(1,9,12,11,10,7,2:6,8)]# II task from CR first condition first batch
setwd("/home/tzhu/Desktop/new_TTC_data/CR_first/first_batch")
csv_ii1 <- list.files(full.names = TRUE, pattern = '*.csv')
df_Pav_ii1 <- as_tibble(rbind.fill(lapply (csv_ii1, fread)))
cols1_ii1 <- c(1, 17:18, 21, 24, 135:138, 141)
df_cat1_ii1<-df_Pav_ii1[,cols1_ii1]
df_cat1_ii1$UniqueID <- toupper(df_cat1_ii1$UniqueID)
names(df_cat1_ii1)[names(df_cat1_ii1) == 'II_loop.thisN'] <- 'trial_num'
df_cat1_ii1$trial_num<-df_cat1_ii1$trial_num + 1
#remove empty rows
df_cat1_ii1<-df_cat1_ii1[!is.na(df_cat1_ii1$key_resp_26.corr),]
#add column to indicate category type
df_cat1_ii1$cat_type<-'II'
names(df_cat1_ii1)[names(df_cat1_ii1) == 'II_cond_(DO_NOT_CHANGE)'] <- 'II_rule_cond'
names(df_cat1_ii1)[names(df_cat1_ii1) == 'key_resp_26.corr'] <- 'corr'
names(df_cat1_ii1)[names(df_cat1_ii1) == 'key_resp_26.keys'] <- 'key_resp'
names(df_cat1_ii1)[names(df_cat1_ii1) == 'key_resp_26.rt'] <- 'rt'
df_cat1_ii1$cat_order <- 'CR_first'
df_cat1_ii1<-df_cat1_ii1[, c(1,4,11,5,12,3,6:10,2)]
df_cat1_ii1$UniqueID[df_cat1_ii1$UniqueID==' 2503ND'] <- '2503ND'# II task from CR first condition second batch
setwd("/home/tzhu/Desktop/new_TTC_data/CR_first/second_batch")
csv_ii2 <- list.files(full.names = TRUE, pattern = '*.csv')
df_Pav_ii2 <- as_tibble(rbind.fill(lapply (csv_ii2, fread)))
cols1_ii2 <- c(1,15, 16,19, 20, 22,133:136, 139)
df_cat1_ii2<-df_Pav_ii2[,cols1_ii2]
df_cat1_ii2$UniqueID <- toupper(df_cat1_ii2$UniqueID)
names(df_cat1_ii2)[names(df_cat1_ii2) == 'II_loop.thisN'] <- 'trial_num'
df_cat1_ii2$trial_num<-df_cat1_ii2$trial_num + 1
#remove empty rows
df_cat1_ii2<-df_cat1_ii2[!is.na(df_cat1_ii2$key_resp_26.corr),]
#add column to indicate category type
df_cat1_ii2$cat_type <- 'II'
names(df_cat1_ii2)[names(df_cat1_ii2) == 'II_cond_(DO_NOT_CHANGE)'] <- 'II_rule_cond'
names(df_cat1_ii2)[names(df_cat1_ii2) == 'key_resp_26.corr'] <- 'corr'
names(df_cat1_ii2)[names(df_cat1_ii2) == 'key_resp_26.keys'] <- 'key_resp'
names(df_cat1_ii2)[names(df_cat1_ii2) == 'key_resp_26.rt'] <- 'rt'
df_cat1_ii2<-df_cat1_ii2[, c(1,4,12,6,5,2,7:11,3)]# II task from CR first condition third batch
setwd("/home/tzhu/Desktop/new_TTC_data/CR_first/third_batch")
csv_ii3 <- list.files(full.names = TRUE, pattern = '*.csv')
df_Pav_ii3 <- as_tibble(rbind.fill(lapply (csv_ii3, fread)))
cols1_ii3 <- c(1,15, 16,19, 21, 23,135:138, 141)
df_cat1_ii3<-df_Pav_third[,cols1_ii3]
df_cat1_ii3$UniqueID <- toupper(df_cat1_ii3$UniqueID)
names(df_cat1_ii3)[names(df_cat1_ii3) == 'II_loop.thisN'] <- 'trial_num'
df_cat1_ii3$trial_num<-df_cat1_ii3$trial_num + 1
#remove empty rows
df_cat1_ii3<-df_cat1_ii3[!is.na(df_cat1_ii3$key_resp_26.corr),]
#add column to indicate category type
df_cat1_ii3$cat_type <- 'II'
names(df_cat1_ii3)[names(df_cat1_ii3) == 'II_cond_(DO_NOT_CHANGE)'] <- 'II_rule_cond'
names(df_cat1_ii3)[names(df_cat1_ii3) == 'key_resp_26.corr'] <- 'corr'
names(df_cat1_ii3)[names(df_cat1_ii3) == 'key_resp_26.keys'] <- 'key_resp'
names(df_cat1_ii3)[names(df_cat1_ii3) == 'key_resp_26.rt'] <- 'rt'
df_cat1_ii3<-df_cat1_ii3[, c(1,4,12,6,5,2,7:11,3)]cat_df_ii1 <- bind_rows(df_cat1_ii,df_cat1_ii1, df_cat1_ii2, df_cat1_ii3)One dataframe for Category learning
cat_df <- bind_rows(cat_df_cr1,cat_df_ii1)
cat_df<-cat_df[,c(1:3,5:12)]
cat_df$UniqueID[cat_df$UniqueID==' 2503ND'] <- '2503ND'perf_df <- cat_df %>%
group_by(UniqueID, cat_type) %>%
dplyr::summarize(average = mean(corr))## `summarise()` has grouped output by 'UniqueID'. You can override using the
## `.groups` argument.cat_df<- merge(cat_df, perf_df, by = c("UniqueID", "cat_type"))RSpan task
RSpan from CR first
#first batch
df_Pav_rspan1 <- df_Pav %>%
mutate_at(vars(colnames(.)),
.funs = funs(ifelse(.=="", NA, as.character(.)))) #fill empty cells with NA
# which(colnames(df_Pav)=="sentence_num_total")
df_rspan1<-df_Pav_rspan1[, c(1, 106:116, 124, 125)] #subset relevant cols
df_rspan1<-df_rspan1[, c(1,2, 4, 6, 12:14)]
df_rspan1<-df_rspan1[!is.na(df_rspan1$exp_Sentence) | !is.na(df_rspan1$exp_typedWord),] #keep rows with either judgment or recall responses
cols4 <- c(3, 4, 6,7)
df_rspan1[cols4] <- sapply(df_rspan1[cols4],as.numeric)#second batch
df_Pav_rspan2 <- df_Pav_second %>%
mutate_at(vars(colnames(.)),
.funs = funs(ifelse(.=="", NA, as.character(.)))) #fill empty cells with NA
# which(colnames(df_Pav)=="sentence_num_total")
df_rspan2<-df_Pav_rspan2[, c(1, 104:114, 122, 123)] #subset relevant cols
df_rspan2<-df_rspan2[, c(1,2, 4, 6, 12:14)]
df_rspan2<-df_rspan2[!is.na(df_rspan2$exp_Sentence) | !is.na(df_rspan2$exp_typedWord),] #keep rows with either judgment or recall responses
cols4 <- c(3, 4, 6,7)
df_rspan2[cols4] <- sapply(df_rspan2[cols4],as.numeric)#third batch
df_Pav_rspan3 <- df_Pav_third %>%
mutate_at(vars(colnames(.)),
.funs = funs(ifelse(.=="", NA, as.character(.)))) #fill empty cells with NA
# which(colnames(df_Pav)=="sentence_num_total")
df_rspan3<-df_Pav_rspan3[, c(1, 106:116, 124, 125)] #subset relevant cols
df_rspan3<-df_rspan3[, c(1,2, 4, 6, 12:14)]
df_rspan3<-df_rspan3[!is.na(df_rspan3$exp_Sentence) | !is.na(df_rspan3$exp_typedWord),] #keep rows with either judgment or recall responses
cols4 <- c(3, 4, 6,7)
df_rspan3[cols4] <- sapply(df_rspan3[cols4],as.numeric)RSpan from II first
df_Pav_rspan4 <- df_Pav_ii %>%
mutate_at(vars(colnames(.)),
.funs = funs(ifelse(.=="", NA, as.character(.)))) #fill empty cells with NA
## Create RSpan DataFrame for II first
# which(colnames(df_Pav)=="sentence_num_total")
df_rspan4<-df_Pav_rspan4[, c(1, 106:116, 124, 125)] #subset relevant cols
df_rspan4<-df_rspan4[, c(1,2, 4,6, 12:14)]
df_rspan4<-df_rspan4[!is.na(df_rspan4$exp_Sentence) | !is.na(df_rspan4$exp_typedWord),] #keep rows with either judgment or recall responses
cols4 <- c(3, 4, 6, 7)
df_rspan4[cols4] <- sapply(df_rspan4[cols4],as.numeric)One dataframe for RSpan
rspan_df <- bind_rows(df_rspan1, df_rspan2, df_rspan3, df_rspan4)
rspan_df$UniqueID <- toupper(rspan_df$UniqueID)
rspan_df$UniqueID[rspan_df$UniqueID==' 2503ND'] <- '2503ND'RSpan judgment part
#keep only judgment related rows
Proportion_Correct <- rspan_df %>%
filter(!is.na(exp_isCorrect))
#late responses on judgment are filled with rt of 4.5
Proportion_Correct$key_resp_10.rt[is.na(Proportion_Correct$key_resp_10.rt)] <- 4.5
Proportion_Correct<-Proportion_Correct[,1:4]RSpan recall part
#keep only recalled_word related rows
recall_num <- rspan_df %>%
filter(!is.na(numberWords_recalled))
recall_num$number_answered<-sapply(strsplit(recall_num$exp_typedWord, " "), length)
recall_num$intrusion_error<-recall_num$number_answered - recall_num$numberWords_recalled #intrusion are the incorrect words typed
recall_num<-recall_num[,c(1,5:8)]RSpan data ready for later merging
#Proportion_Correct dataframe
df1<-aggregate(Proportion_Correct[, 3], list(Proportion_Correct$UniqueID), mean) #average for each Ss performance
#recalled number dataframe
df2<-aggregate(recall_num[, c(3:5)], list(recall_num$UniqueID), sum)
df_rspan_perf <- merge(df1,df2,by="Group.1")
names(df_rspan_perf)[1] <- "UniqueID"
#head(df_rspan_perf)Merge temperament df with Rspan df
df_tem_rspan1<-merge(df_temp, df_rspan_perf, by = "UniqueID")
# df_tem_rspan<-df_tem_rspan1[df_tem_rspan1$numberWords_recalled > 15,]
#head(df_tem_rspan)Merge category learning with tem_rspan
a <- as.data.frame(table(cat_df$UniqueID, cat_df$cat_type)) #check whether same ID did study more than once
cat_tem_rspan<-merge(cat_df, df_tem_rspan1, by.x="UniqueID",by.y="UniqueID")full_df <- cat_tem_rspan
#append temperament profile
full_df$profile1 <- ifelse(full_df$NA_level == "high" & full_df$EX_level == "low" , "highNA, lowEX",
ifelse(full_df$NA_level == "low" & full_df$EX_level == "high" ,"lowNA, highEX",
ifelse(full_df$NA_level == "low" & full_df$EX_level == "low" ,"lowNA, lowEX", 'highNA, highEX')))cat_df1<-as.data.table(cat_df)
df_10correct<-cat_df1[cat_df1[, .I[corr_count==10], by=UniqueID]$V1]
df_10correct_cr<-df_10correct[df_10correct$cat_type=='CR',]
df_10correct_ii<-df_10correct[df_10correct$cat_type=='II',]
first_10correct_cr<-df_10correct_cr[ , .SD[which.min(trial_num)], by = UniqueID]
first_10correct_ii<-df_10correct_ii[ , .SD[which.min(trial_num)], by = UniqueID]
cr_max<-cat_df1[cat_df1[cat_df1$cat_type=='CR', .I[corr_count==max(corr_count)], by=UniqueID]$V1]
cr_max<-unique(cr_max[,c(1,2,9)])
ii_max<-cat_df1[cat_df1[cat_df1$cat_type=='II', .I[corr_count==max(corr_count)], by=UniqueID]$V1]
ii_max<-unique(ii_max[,c(1,2,9)])
max_df <- bind_rows(cr_max,ii_max)
names(max_df)[names(max_df) == 'corr_count'] <- 'max_corr_count'
full_df1 <-merge(full_df, max_df, by=c('UniqueID','cat_type'))
first_10correct<-bind_rows(first_10correct_cr,first_10correct_ii)
first_10correct<-first_10correct[,c(1,2,10)]
names(first_10correct)[names(first_10correct) == 'trial_num'] <- 'first_TTC'
full_df2<-merge(full_df1, first_10correct, by = c('UniqueID', 'cat_type'), all.x = TRUE)
full_df3<-full_df2[order(full_df2[,1],full_df2[,2],full_df2[,8]),]d <-full_df3[,c(1:4, 12:36)] #remove columns with cat stimuli information
d<-distinct(d) #remove duplicate rows
d<-d[d$numberWords_recalled>0,]
d[is.na(d)]<-200
d<-d[,c(1:5,9:11,13,17:29)]d_wide<-reshape(d, idvar = c("UniqueID", "cat_order","Ethnicity","Age","Gender","Country_origin", "NA_score",'EC_score','EX_score','NA_level','EC_level','EX_level', 'exp_isCorrect','numberWords_recalled','sentence_num_total','number_answered','profile1'), timevar = c('cat_type'), dir = "wide")
cr_better<-d_wide[d_wide$first_TTC.CR<d_wide$first_TTC.II,]
cr_better$group <- 'cr_better'
ii_better<-d_wide[d_wide$first_TTC.II<d_wide$first_TTC.CR,]
ii_better$group<-'ii_better'
same<-d_wide[d_wide$first_TTC.II==200 & d_wide$first_TTC.CR==200,]
same$group<-'both_nonlearner'
d_wide<-rbind(cr_better,ii_better, same)
d_wide$cat_perf_CR<-ifelse(d_wide$average.CR<=median(d_wide$average.CR),"low","high")
d_wide$cat_perf_II<-ifelse(d_wide$average.II<=median(d_wide$average.II),"low","high")
d_wide$cat_perf_group<-ifelse(d_wide$cat_perf_CR=='high' & d_wide$cat_perf_II=='high', 'hh',
ifelse(d_wide$cat_perf_CR=='high' & d_wide$cat_perf_II=='low', 'hl',
ifelse(d_wide$cat_perf_CR=='low' & d_wide$cat_perf_II=='high', 'lh','ll')))
d_wide$Ethnicity[d_wide$Ethnicity=='South Asian'] <- 'Asian or Pacific Islander'
d_wide$intrusion<-d_wide$number_answered-d_wide$numberWords_recalled
duration<-read.csv('duration.csv')
setnames(duration, old = c('Q102','Duration..in.seconds.'), new = c('UniqueID','duration_seconds'))
duration$UniqueID<-toupper(duration$UniqueID)
d_wide<-merge(d_wide,duration,by=c('UniqueID'))Ethnic distribution
library(RColorBrewer)
mytable <- table(d_wide$Ethnicity)
lbls <- paste(names(mytable), " - ", mytable, sep="")
coul <- brewer.pal(9, "Set3")
pie(mytable, labels = lbls,col=coul, cex = 1.2,
main="")
summary(d_wide[,c("NA_score","EX_score","EC_score","numberWords_recalled")]) #because different subjects may be removed in each task, this summary reflects the complete 300 people data unlike in Experiment 1, post-cleaning descriptives were reported.## NA_score EX_score EC_score numberWords_recalled
## Min. :2.560 Min. :2.588 Min. :1.842 Min. : 7.00
## 1st Qu.:4.110 1st Qu.:4.471 1st Qu.:3.368 1st Qu.:41.00
## Median :4.520 Median :5.059 Median :3.789 Median :47.00
## Mean :4.523 Mean :5.003 Mean :3.850 Mean :45.44
## 3rd Qu.:5.000 3rd Qu.:5.529 3rd Qu.:4.368 3rd Qu.:52.25
## Max. :6.520 Max. :7.235 Max. :5.842 Max. :60.00x<-d_wide[,c(7:9, 14)]
sapply(x, function(x) c( "Stand dev" = sd(x),
"Mean"= mean(x,na.rm=TRUE),
"Median" = median(x),
"Minimum" = min(x),
"Maximun" = max(x),
"Upper Quantile" = quantile(x,1),
"LowerQuartile" = quantile(x,0)
)
)## NA_score EC_score EX_score numberWords_recalled
## Stand dev 0.712845 0.7494327 0.7978208 9.50863
## Mean 4.523176 3.8495482 5.0030125 45.43667
## Median 4.520000 3.7894737 5.0588235 47.00000
## Minimum 2.560000 1.8421053 2.5882353 7.00000
## Maximun 6.520000 5.8421053 7.2352941 60.00000
## Upper Quantile.100% 6.520000 5.8421053 7.2352941 60.00000
## LowerQuartile.0% 2.560000 1.8421053 2.5882353 7.00000Compare temperament scores across largest ethnicity groups
race_data<-d_wide[d_wide$Ethnicity=='White or Caucasian'|d_wide$Ethnicity=='Asian or Pacific Islander',]
table(race_data$Ethnicity)##
## Asian or Pacific Islander White or Caucasian
## 131 132t.test(race_data$NA_score~race_data$Ethnicity)##
## Welch Two Sample t-test
##
## data: race_data$NA_score by race_data$Ethnicity
## t = 0.37557, df = 260.88, p-value = 0.7075
## alternative hypothesis: true difference in means between group Asian or Pacific Islander and group White or Caucasian is not equal to 0
## 95 percent confidence interval:
## -0.1382115 0.2033608
## sample estimates:
## mean in group Asian or Pacific Islander mean in group White or Caucasian
## 4.513581 4.481006t.test(race_data$EX_score~race_data$Ethnicity)##
## Welch Two Sample t-test
##
## data: race_data$EX_score by race_data$Ethnicity
## t = -4.8605, df = 260.58, p-value = 2.025e-06
## alternative hypothesis: true difference in means between group Asian or Pacific Islander and group White or Caucasian is not equal to 0
## 95 percent confidence interval:
## -0.6648945 -0.2814884
## sample estimates:
## mean in group Asian or Pacific Islander mean in group White or Caucasian
## 4.791503 5.264695t.test(race_data$EC_score~race_data$Ethnicity)##
## Welch Two Sample t-test
##
## data: race_data$EC_score by race_data$Ethnicity
## t = 0.17344, df = 260.64, p-value = 0.8624
## alternative hypothesis: true difference in means between group Asian or Pacific Islander and group White or Caucasian is not equal to 0
## 95 percent confidence interval:
## -0.1685590 0.2011209
## sample estimates:
## mean in group Asian or Pacific Islander mean in group White or Caucasian
## 3.863570 3.847289Compare temperament scores across genders
gender_data<-d_wide[d_wide$Gender=='Female'|d_wide$Gender=='Male',]
table(gender_data$Gender)##
## Female Male
## 233 62t.test(gender_data$NA_score~gender_data$Gender)##
## Welch Two Sample t-test
##
## data: gender_data$NA_score by gender_data$Gender
## t = 5.9422, df = 102.53, p-value = 3.901e-08
## alternative hypothesis: true difference in means between group Female and group Male is not equal to 0
## 95 percent confidence interval:
## 0.3656005 0.7319323
## sample estimates:
## mean in group Female mean in group Male
## 4.633411 4.084644Analysis with ALL participants using 3 temperament factors as predictors
CR task
Overall performance
setnames(d_wide, old = c('NA_score','EC_score','EX_score','numberWords_recalled'), new = c('negative_affect','effortful_control','extraversion','Rspan_score'))
summary(d_wide[,c("negative_affect","extraversion","effortful_control","Rspan_score")])## negative_affect extraversion effortful_control Rspan_score
## Min. :2.560 Min. :2.588 Min. :1.842 Min. : 7.00
## 1st Qu.:4.110 1st Qu.:4.471 1st Qu.:3.368 1st Qu.:41.00
## Median :4.520 Median :5.059 Median :3.789 Median :47.00
## Mean :4.523 Mean :5.003 Mean :3.850 Mean :45.44
## 3rd Qu.:5.000 3rd Qu.:5.529 3rd Qu.:4.368 3rd Qu.:52.25
## Max. :6.520 Max. :7.235 Max. :5.842 Max. :60.00avg_wide<-read.csv('avg_wide.csv')#avg_wide is created in strategy_analysis.Rmd file, it has 4 block performance for CR and II
DATA<-merge(d_wide,avg_wide,by='UniqueID')
# find and remove bad data
dat<-DATA[DATA$exp_isCorrect>0.80 & DATA$number_answered<100 & DATA$intrusion<60&(DATA$performance_block_1_CR<DATA$performance_block_4_CR|DATA$max.CR>0.599)&DATA$duration_seconds>300,]
dat$group2 <- ifelse(dat$first_TTC.CR==200,"cr_nlearner","cr_learner")descriptives after general data cleaning
dat_general_removal<-DATA[DATA$exp_isCorrect>0.80 & DATA$number_answered<100 & DATA$intrusion<60&DATA$duration_seconds>300,]
x<-dat_general_removal[,c(7:9, 14)]
sapply(x, function(x) c( "Stand dev" = sd(x),
"Mean"= mean(x,na.rm=TRUE),
"Median" = median(x),
"Minimum" = min(x),
"Maximun" = max(x),
"Upper Quantile" = quantile(x,1),
"LowerQuartile" = quantile(x,0)
)
)## negative_affect effortful_control extraversion Rspan_score
## Stand dev 0.7241304 0.780895 0.8296661 7.694052
## Mean 4.5246679 3.890690 5.0145589 46.972222
## Median 4.5200000 3.894737 5.0882353 48.000000
## Minimum 2.5600000 2.000000 2.5882353 16.000000
## Maximun 6.5200000 5.842105 7.2352941 60.000000
## Upper Quantile.100% 6.5200000 5.842105 7.2352941 60.000000
## LowerQuartile.0% 2.5600000 2.000000 2.5882353 16.000000hist(dat$average.CR, col='yellow', xlab='Average CR Task Performance', main="", labels = T, ylim=c(0,50))
# df_for_tables<- dat[dat$Rspan_score>20,]
df_for_tables<- dat
df.corr<- rcorr(as.matrix(df_for_tables[,c(7:9, 14,19)]))
df.corr$r<-round(df.corr$r,3)
df.corr$p<-round(df.corr$P, 3)
df.corr$r[upper.tri(df.corr$r)] <- "-"
df.corr$p[upper.tri(df.corr$r)] <- "-"
df.corr$p## negative_affect effortful_control extraversion Rspan_score
## negative_affect NA "-" "-" "-"
## effortful_control "0" NA "-" "-"
## extraversion "0" "0.655" NA "-"
## Rspan_score "0.617" "0.359" "0.161" NA
## average.CR "0.775" "0.86" "0.67" "0.042"
## average.CR
## negative_affect "-"
## effortful_control "-"
## extraversion "-"
## Rspan_score "-"
## average.CR NAhist(dat$first_TTC.CR, col='yellow', xlab='Trials to Criterion on the CR Task', main="", labels=T, ylim=c(0,120))
m_cr<-lm(average.CR ~ negative_affect+extraversion+effortful_control+negative_affect*effortful_control+extraversion*effortful_control+Rspan_score+cond.CR+cat_order, dat)
summary(m_cr)##
## Call:
## lm(formula = average.CR ~ negative_affect + extraversion + effortful_control +
## negative_affect * effortful_control + extraversion * effortful_control +
## Rspan_score + cond.CR + cat_order, data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.32641 -0.13074 0.02095 0.12984 0.29065
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.748893 0.546473 3.200 0.00157 **
## negative_affect -0.125631 0.065629 -1.914 0.05688 .
## extraversion -0.115573 0.067388 -1.715 0.08774 .
## effortful_control -0.284526 0.128273 -2.218 0.02756 *
## Rspan_score 0.002599 0.001403 1.852 0.06531 .
## cond.CRCR2 -0.011628 0.035951 -0.323 0.74667
## cond.CRCR3 0.034159 0.036080 0.947 0.34479
## cond.CRCR4 0.009133 0.037326 0.245 0.80694
## cond.CRCR5 -0.058765 0.037179 -1.581 0.11540
## cond.CRCR6 0.065757 0.037280 1.764 0.07913 .
## cat_orderII_first -0.028928 0.021559 -1.342 0.18105
## negative_affect:effortful_control 0.031604 0.015981 1.978 0.04922 *
## extraversion:effortful_control 0.029791 0.016742 1.779 0.07654 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1579 on 221 degrees of freedom
## Multiple R-squared: 0.09485, Adjusted R-squared: 0.0457
## F-statistic: 1.93 on 12 and 221 DF, p-value: 0.03213m_cr_stan<- lm(scale(average.CR)~scale(negative_affect)+scale(extraversion)+scale(effortful_control)+scale(negative_affect*effortful_control)+scale(extraversion*effortful_control)+scale(Rspan_score)+ cond.CR + cat_order, dat)
summary(m_cr_stan)##
## Call:
## lm(formula = scale(average.CR) ~ scale(negative_affect) + scale(extraversion) +
## scale(effortful_control) + scale(negative_affect * effortful_control) +
## scale(extraversion * effortful_control) + scale(Rspan_score) +
## cond.CR + cat_order, data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.0192 -0.8088 0.1296 0.8032 1.7980
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.04033 0.16804 0.240 0.8105
## scale(negative_affect) -0.55185 0.28828 -1.914 0.0569
## scale(extraversion) -0.59810 0.34874 -1.715 0.0877
## scale(effortful_control) -1.36096 0.61356 -2.218 0.0276
## scale(negative_affect * effortful_control) 0.65641 0.33192 1.978 0.0492
## scale(extraversion * effortful_control) 0.96093 0.54000 1.779 0.0765
## scale(Rspan_score) 0.12368 0.06677 1.852 0.0653
## cond.CRCR2 -0.07193 0.22240 -0.323 0.7467
## cond.CRCR3 0.21132 0.22320 0.947 0.3448
## cond.CRCR4 0.05650 0.23091 0.245 0.8069
## cond.CRCR5 -0.36353 0.23000 -1.581 0.1154
## cond.CRCR6 0.40679 0.23062 1.764 0.0791
## cat_orderII_first -0.17895 0.13337 -1.342 0.1810
##
## (Intercept)
## scale(negative_affect) .
## scale(extraversion) .
## scale(effortful_control) *
## scale(negative_affect * effortful_control) *
## scale(extraversion * effortful_control) .
## scale(Rspan_score) .
## cond.CRCR2
## cond.CRCR3
## cond.CRCR4
## cond.CRCR5
## cond.CRCR6 .
## cat_orderII_first
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9769 on 221 degrees of freedom
## Multiple R-squared: 0.09485, Adjusted R-squared: 0.0457
## F-statistic: 1.93 on 12 and 221 DF, p-value: 0.03213Regression stats table in R
# create table depicting main stats to export later
# use tidy() in broom package & DescTools::StdCoef()
library(broom)
library(DescTools)##
## Attaching package: 'DescTools'## The following objects are masked from 'package:Hmisc':
##
## %nin%, Label, Mean, Quantile## The following object is masked from 'package:data.table':
##
## %like%m_cr_coef <- cbind(tidy(m_cr)[,1:3],
tidy(m_cr_stan)[c(1:4,7,8:13,5,6),2:3],
tidy(m_cr)[,4:5])
# remove rownames
rownames(m_cr_coef) <- c(1:nrow(m_cr_coef))
# rename cols
colnames(m_cr_coef) <- c("term", "beta", "beta_SE", "std_beta", "std_beta_SE", "t_val", "p_val")
# finally, create df
mod1_t <-
rbind(m_cr_coef,
tibble("term" = c("f_stat", "R2/Adj.R2", "N"),
"beta" = c(summary(m_cr)$fstatistic[2], summary(m_cr)$r.squared, nrow(m_cr$model)),
"beta_SE" = c(summary(m_cr)$fstatistic[3], summary(m_cr)$adj.r.squared, NA),
"std_beta" = c(summary(m_cr)$fstatistic[1], rep(NA, 2)),
"std_beta_SE" = c(pf(summary(m_cr)$fstatistic[1], summary(m_cr)$fstatistic[2], summary(m_cr)$fstatistic[3], lower.tail = F), rep(NA, 2)),
"t_val" = rep(NA, 3),
"p_val" = rep(NA, 3)
)
)
# remove NA from kable tables
options(knitr.kable.NA = '')
# show kable table
kable(mod1_t, caption = "Multiple Regression Model Predicting CR overall average performance", align = rep('lcccc'), digits = 3, row.names = T, col.names = c("", "Beta", "Beta's SE", "Std. Beta", "Std. Beta's SE", "t-value", "p-value")) %>%
kable_paper(full_width = F) %>%
# add_header_above(c(" ", "mturk" = 2, "Mturk" = 2)) %>%
# pack_rows("Criterion(s)", 1, 6) %>%
# # column_spec(1, width = "10em") %>%
# pack_rows("Predictor(s)", 7, 19) %>%
footnote(general = "...") %>%
column_spec(c(4), border_left = T) %>%
kable_classic()| Beta | Beta’s SE | Std. Beta | Std. Beta’s SE | t-value | p-value | ||
|---|---|---|---|---|---|---|---|
| 1 | (Intercept) | 1.749 | 0.546 | 0.040 | 0.168 | 3.200 | 0.002 |
| 2 | negative_affect | -0.126 | 0.066 | -0.552 | 0.288 | -1.914 | 0.057 |
| 3 | extraversion | -0.116 | 0.067 | -0.598 | 0.349 | -1.715 | 0.088 |
| 4 | effortful_control | -0.285 | 0.128 | -1.361 | 0.614 | -2.218 | 0.028 |
| 5 | Rspan_score | 0.003 | 0.001 | 0.124 | 0.067 | 1.852 | 0.065 |
| 6 | cond.CRCR2 | -0.012 | 0.036 | -0.072 | 0.222 | -0.323 | 0.747 |
| 7 | cond.CRCR3 | 0.034 | 0.036 | 0.211 | 0.223 | 0.947 | 0.345 |
| 8 | cond.CRCR4 | 0.009 | 0.037 | 0.056 | 0.231 | 0.245 | 0.807 |
| 9 | cond.CRCR5 | -0.059 | 0.037 | -0.364 | 0.230 | -1.581 | 0.115 |
| 10 | cond.CRCR6 | 0.066 | 0.037 | 0.407 | 0.231 | 1.764 | 0.079 |
| 11 | cat_orderII_first | -0.029 | 0.022 | -0.179 | 0.133 | -1.342 | 0.181 |
| 12 | negative_affect:effortful_control | 0.032 | 0.016 | 0.656 | 0.332 | 1.978 | 0.049 |
| 13 | extraversion:effortful_control | 0.030 | 0.017 | 0.961 | 0.540 | 1.779 | 0.077 |
| 14 | f_stat | 12.000 | 221.000 | 1.930 | 0.032 | ||
| 15 | R2/Adj.R2 | 0.095 | 0.046 | ||||
| 16 | N | 234.000 | |||||
| Note: | |||||||
| … |
write.csv(mod1_t,"TTC_CR_reg.csv", row.names = T)# reset NA vals in kable tables to "NA"
options(knitr.kable.NA = NA)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
# plot lm
# plot standardized betas to compare cooefficients
ttc_cr_reg_plot <-
sjPlot::plot_model(model = m_cr_stan, # lm output
title = "", # to remove title
# type = "std", # type "std" = standardized betas
#sort.est = T, # sort based on betas (highest on top)
# wrap.title = 55, # char length of each title line
#axis.lim = c(-.25,.25)
# grid.breaks = .25, # x-axis breaks
#se = T # use se instead of CI
colors = reg_colours,
show.values = T, # show betas
show.p = T, # show * for sig p
p.threshold = c(.049, .009, .0009), # state up to 3 cut-offs for pvals, the number is inclusive**
value.offset = 0.4, # offset to beta labels
value.size = 3, # size of beta labels
digits = 2, # num decimals in beta labels
dot.size = 2.5, line.size = 1,
vline.color = "#F98B99",
# group.terms = c(1, 2, 3, 3, 3, 3, 3, 4, 5, 6, 7, 7, 7),
axis.labels = c("Cat Order : II first", "CR Condition 6", "CR Condition 5","CR Condition 4","CR Condition 3","CR Condition 2", "Working Memory", "Extraversion x Effortful Control", "Negative Affect x Effortful Control", "Effortful Control", "Extraversion", "Negative Affect"),
axis.title = "",
width = 0.2, transform = NULL
) +
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(color = "black"),
axis.text = element_text(colour = "black", size = 11),
axis.ticks.y = element_blank(),
)
# show plot
ttc_cr_reg_plot
Scatterplot of correlational relationship between each predictor and DV
library(tidyr)
predictor_corr_data<-gather(dat, predictor, Score, c(7:9,14), factor_key = TRUE)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
predictor_names <- list(
'negative_affect'="Negative Affect",
'effortful_control'="Effortful Control",
'extraversion'="Extraversion",
'Rspan_score'="Working Memory"
)
predictor_labeller <- function(variable,value){
return(predictor_names[value])
}
predictor_scatter <-
ggplot(predictor_corr_data, aes(y = average.CR, x = Score, colour = factor(predictor))) +
geom_point(shape = 1) +
scale_color_manual(values = reg_colours[c(2,2,2,2)]) +
facet_wrap(predictor ~ ., labeller=predictor_labeller, scales = "free_x")+
# facet_wrap(~ predictor, ncol = 2, scales = "free_x") +
geom_smooth(method = "lm") +
xlab("Predictor Value")+
ylab("Overall Average Performance") +
# geom_line(intercept = 83.731, slopr = .5)+
# add theme
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(size=20),
axis.text.x = element_text(angle=0, hjust=2)) +
theme_classic() +
theme(
legend.position = "none",
text = element_text(color = "black"),
axis.text = element_text(colour = "black"),
axis.ticks.y = element_blank(),
strip.text.x = element_text(size=12),
axis.title = element_text(size=12)
)## Warning: The labeller API has been updated. Labellers taking `variable` and
## `value` arguments are now deprecated. See labellers documentation.# show fig
predictor_scatter ## `geom_smooth()` using formula 'y ~ x'
step_mcr<-step(m_cr)
step_mcr_stan<-step(m_cr_stan)
summary(step_mcr)
summary(step_mcr_stan)# create table depicting main stats to export later
# use tidy() in broom package & DescTools::StdCoef()
library(broom)
library(DescTools)
m_cr_red_coef <- cbind(tidy(step_mcr)[c(1:4, 11:12,5, 6:10),1:3],
tidy(step_mcr_stan)[,2:3],
tidy(step_mcr)[c(1:4, 11:12,5, 6:10),4:5])
# remove rownames
rownames(m_cr_red_coef) <- c(1:nrow(m_cr_red_coef))
# rename cols
colnames(m_cr_red_coef) <- c("term", "beta", "beta_SE", "std_beta", "std_beta_SE", "t_val", "p_val")
# finally, create df
reduced_TTC_CR_mod1_t <-
rbind(m_cr_red_coef,
tibble("term" = c("f_stat", "R2/Adj.R2", "N"),
"beta" = c(summary(step_mcr)$fstatistic[2], summary(step_mcr)$r.squared, nrow(step_mcr$model)),
"beta_SE" = c(summary(step_mcr)$fstatistic[3], summary(step_mcr)$adj.r.squared, NA),
"std_beta" = c(summary(step_mcr)$fstatistic[1], rep(NA, 2)),
"std_beta_SE" = c(pf(summary(step_mcr)$fstatistic[1], summary(step_mcr)$fstatistic[2], summary(step_mcr)$fstatistic[3], lower.tail = F), rep(NA, 2)),
"t_val" = rep(NA, 3),
"p_val" = rep(NA, 3)
)
)
# remove NA from kable tables
options(knitr.kable.NA = '')
# show kable table
kable(reduced_TTC_CR_mod1_t, caption = "Multiple Regression Model Predicting CR overall average performance", align = rep('lcccc'), digits = 3, row.names = T, col.names = c("", "Beta", "Beta's SE", "Std. Beta", "Std. Beta's SE", "t-value", "p-value")) %>%
kable_paper(full_width = F) %>%
# add_header_above(c(" ", "mturk" = 2, "Mturk" = 2)) %>%
# pack_rows("Criterion(s)", 1, 6) %>%
# # column_spec(1, width = "10em") %>%
# pack_rows("Predictor(s)", 7, 19) %>%
footnote(general = "...") %>%
column_spec(c(4), border_left = T) %>%
kable_classic()| Beta | Beta’s SE | Std. Beta | Std. Beta’s SE | t-value | p-value | ||
|---|---|---|---|---|---|---|---|
| 1 | (Intercept) | 1.756 | 0.547 | -0.022 | 0.162 | 3.208 | 0.002 |
| 2 | negative_affect | -0.127 | 0.066 | -0.556 | 0.289 | -1.926 | 0.055 |
| 3 | extraversion | -0.119 | 0.067 | -0.618 | 0.349 | -1.771 | 0.078 |
| 4 | effortful_control | -0.288 | 0.128 | -1.375 | 0.615 | -2.238 | 0.026 |
| 5 | negative_affect:effortful_control | 0.031 | 0.016 | 0.646 | 0.332 | 1.943 | 0.053 |
| 6 | extraversion:effortful_control | 0.030 | 0.017 | 0.984 | 0.541 | 1.819 | 0.070 |
| 7 | Rspan_score | 0.003 | 0.001 | 0.137 | 0.066 | 2.072 | 0.039 |
| 8 | cond.CRCR2 | -0.012 | 0.036 | -0.071 | 0.223 | -0.321 | 0.749 |
| 9 | cond.CRCR3 | 0.029 | 0.036 | 0.177 | 0.222 | 0.798 | 0.426 |
| 10 | cond.CRCR4 | 0.004 | 0.037 | 0.027 | 0.230 | 0.116 | 0.908 |
| 11 | cond.CRCR5 | -0.063 | 0.037 | -0.391 | 0.229 | -1.704 | 0.090 |
| 12 | cond.CRCR6 | 0.058 | 0.037 | 0.361 | 0.229 | 1.581 | 0.115 |
| 13 | f_stat | 11.000 | 222.000 | 1.935 | 0.036 | ||
| 14 | R2/Adj.R2 | 0.087 | 0.042 | ||||
| 15 | N | 234.000 | |||||
| Note: | |||||||
| … |
write.csv(reduced_TTC_CR_mod1_t,"reduced_TTC_CR_mod1_t.csv", row.names = T)# reset NA vals in kable tables to "NA"
options(knitr.kable.NA = NA)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
# plot lm
# plot standardized betas to compare cooefficients
ttc_cr_reg_reduced_plot <-
sjPlot::plot_model(model = step_mcr_stan, # lm output
title = "", # to remove title
# type = "std", # type "std" = standardized betas
#sort.est = T, # sort based on betas (highest on top)
# wrap.title = 55, # char length of each title line
#axis.lim = c(-.25,.25)
# grid.breaks = .25, # x-axis breaks
#se = T # use se instead of CI
colors = reg_colours,
show.values = T, # show betas
show.p = T, # show * for sig p
p.threshold = c(.049, .009, .0009), # state up to 3 cut-offs for pvals, the number is inclusive**
value.offset = 0.4, # offset to beta labels
value.size = 3, # size of beta labels
digits = 2, # num decimals in beta labels
dot.size = 2.5, line.size = 1,
vline.color = "#F98B99",
# group.terms = c(1, 2, 3, 3, 3, 3, 3, 4, 5, 6, 7, 7, 7),
axis.labels = c("CR Condition 6", "CR Condition 5","CR Condition 4","CR Condition 3","CR Condition 2", "Working Memory", "Extraversion x Effortful Control", "Negative Affect x Effortful Control", "Effortful Control", "Extraversion", "Negative Affect"),
axis.title = "",
width = 0.2, transform = NULL
) +
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(color = "black"),
axis.text = element_text(colour = "black", size = 11),
axis.ticks.y = element_blank(),
)
# show plot
ttc_cr_reg_reduced_plot
anova_cond<-aov(average.CR~cond.CR, dat)
summary(anova_cond)## Df Sum Sq Mean Sq F value Pr(>F)
## cond.CR 5 0.237 0.04734 1.844 0.105
## Residuals 228 5.852 0.02567anova_cond<-aov(first_TTC.CR~cond.CR, dat[dat$first_TTC.CR<200,])
summary(anova_cond)## Df Sum Sq Mean Sq F value Pr(>F)
## cond.CR 5 24476 4895 2.338 0.0442 *
## Residuals 161 337118 2094
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tuk<-TukeyHSD(anova_cond)Tukey plot comparing every 2 CR conditions
tuk_plot <- function (x, xlab, ylab, ylabels = NULL, ...) {
for (i in seq_along(x)) {
xi <- x[[i]][, -4L, drop = FALSE]
yvals <- nrow(xi):1L
dev.hold()
on.exit(dev.flush())
plot(c(xi[, "lwr"], xi[, "upr"]), rep.int(yvals, 2L),
type = "n", axes = FALSE, xlab = "", ylab = "", main = NULL,
...)
axis(1, ...)
# change for custom axis labels
if (is.null(ylabels)) ylabels <- dimnames(xi)[[1L]]
axis(2, at = nrow(xi):1, labels = ylabels,
srt = 0, ...)
abline(h = yvals, lty = 1, lwd = 1, col = "lightgray")
abline(v = 0, lty = 2, lwd = 1, ...)
segments(xi[, "lwr"], yvals, xi[, "upr"], yvals, ...)
segments(as.vector(xi), rep.int(yvals - 0.1, 3L), as.vector(xi),
rep.int(yvals + 0.1, 3L), ...)
title(main ="",
# change for custom axis titles
xlab = xlab, ylab = ylab)
box()
dev.flush()
on.exit()
}
}
with(par(mai=c(0.7,0.8,0.5,0.5), cex.lab=1.2, cex.axis=0.8),{tuk_plot(tuk, "", "condition",c("2-1","3-1","4-1","5-1","6-1","3-2","4-2","5-2","6-2","4-3","5-3","6-3","5-4","6-4","6-5"),las=2)})
### TTC (learners)
m_cr1<-lm(first_TTC.CR ~ negative_affect+extraversion+effortful_control+negative_affect*effortful_control+extraversion*effortful_control+Rspan_score + cond.CR + cat_order, dat[dat$first_TTC.CR<200,])
m_cr1_stan<- lm(scale(first_TTC.CR)~scale(negative_affect)+scale(extraversion)+scale(effortful_control)+scale(negative_affect*effortful_control)+scale(extraversion*effortful_control)+scale(Rspan_score)+ cond.CR + cat_order, dat[dat$first_TTC.CR<200,])
summary(m_cr1)##
## Call:
## lm(formula = first_TTC.CR ~ negative_affect + extraversion +
## effortful_control + negative_affect * effortful_control +
## extraversion * effortful_control + Rspan_score + cond.CR +
## cat_order, data = dat[dat$first_TTC.CR < 200, ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -70.475 -32.906 -9.668 27.791 138.266
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 93.2786 205.5752 0.454 0.651
## negative_affect -19.7379 25.5799 -0.772 0.442
## extraversion 22.4527 24.7298 0.908 0.365
## effortful_control 14.9193 51.2095 0.291 0.771
## Rspan_score -0.5948 0.5235 -1.136 0.258
## cond.CRCR2 14.4275 11.5980 1.244 0.215
## cond.CRCR3 -18.2986 12.0870 -1.514 0.132
## cond.CRCR4 -3.6443 12.7055 -0.287 0.775
## cond.CRCR5 10.7225 13.8462 0.774 0.440
## cond.CRCR6 -5.1941 12.0720 -0.430 0.668
## cat_orderII_first -6.4547 7.3828 -0.874 0.383
## negative_affect:effortful_control 3.8709 6.4964 0.596 0.552
## extraversion:effortful_control -6.9734 6.1676 -1.131 0.260
##
## Residual standard error: 45.89 on 154 degrees of freedom
## Multiple R-squared: 0.103, Adjusted R-squared: 0.03311
## F-statistic: 1.474 on 12 and 154 DF, p-value: 0.1395summary(m_cr1_stan)##
## Call:
## lm(formula = scale(first_TTC.CR) ~ scale(negative_affect) + scale(extraversion) +
## scale(effortful_control) + scale(negative_affect * effortful_control) +
## scale(extraversion * effortful_control) + scale(Rspan_score) +
## cond.CR + cat_order, data = dat[dat$first_TTC.CR < 200, ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.5100 -0.7050 -0.2072 0.5955 2.9625
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.07421 0.18873 0.393 0.695
## scale(negative_affect) -0.28882 0.37430 -0.772 0.442
## scale(extraversion) 0.39868 0.43912 0.908 0.365
## scale(effortful_control) 0.24598 0.84429 0.291 0.771
## scale(negative_affect * effortful_control) 0.29282 0.49144 0.596 0.552
## scale(extraversion * effortful_control) -0.77407 0.68463 -1.131 0.260
## scale(Rspan_score) -0.09125 0.08030 -1.136 0.258
## cond.CRCR2 0.30912 0.24850 1.244 0.215
## cond.CRCR3 -0.39207 0.25898 -1.514 0.132
## cond.CRCR4 -0.07808 0.27223 -0.287 0.775
## cond.CRCR5 0.22974 0.29667 0.774 0.440
## cond.CRCR6 -0.11129 0.25866 -0.430 0.668
## cat_orderII_first -0.13830 0.15818 -0.874 0.383
##
## Residual standard error: 0.9833 on 154 degrees of freedom
## Multiple R-squared: 0.103, Adjusted R-squared: 0.03311
## F-statistic: 1.474 on 12 and 154 DF, p-value: 0.1395Regression stats table in R
library(broom)
library(DescTools)
m_cr_ttc_coef <- cbind(tidy(m_cr1)[c(1:4, 12:13,5, 6:11),1:3],
tidy(m_cr1_stan)[,2:3],
tidy(m_cr1)[c(1:4, 12:13,5, 6:11),4:5])
# remove rownames
rownames(m_cr_ttc_coef) <- c(1:nrow(m_cr_ttc_coef))
# rename cols
colnames(m_cr_ttc_coef) <- c("term", "beta", "beta_SE", "std_beta", "std_beta_SE", "t_val", "p_val")
# finally, create df
learner_TTC_CR_mod1_t <-
rbind(m_cr_ttc_coef,
tibble("term" = c("f_stat", "R2/Adj.R2", "N"),
"beta" = c(summary(m_cr1)$fstatistic[2], summary(m_cr1)$r.squared, nrow(m_cr1$model)),
"beta_SE" = c(summary(m_cr1)$fstatistic[3], summary(m_cr1)$adj.r.squared, NA),
"std_beta" = c(summary(m_cr1)$fstatistic[1], rep(NA, 2)),
"std_beta_SE" = c(pf(summary(m_cr1)$fstatistic[1], summary(m_cr1)$fstatistic[2], summary(m_cr1)$fstatistic[3], lower.tail = F), rep(NA, 2)),
"t_val" = rep(NA, 3),
"p_val" = rep(NA, 3)
)
)
# remove NA from kable tables
options(knitr.kable.NA = '')
# show kable table
kable(learner_TTC_CR_mod1_t, caption = "Multiple Regression Model Predicting TTC in the CR task", align = rep('lcccc'), digits = 3, row.names = T, col.names = c("", "Beta", "Beta's SE", "Std. Beta", "Std. Beta's SE", "t-value", "p-value")) %>%
kable_paper(full_width = F) %>%
# add_header_above(c(" ", "mturk" = 2, "Mturk" = 2)) %>%
# pack_rows("Criterion(s)", 1, 6) %>%
# # column_spec(1, width = "10em") %>%
# pack_rows("Predictor(s)", 7, 19) %>%
footnote(general = "...") %>%
column_spec(c(4), border_left = T) %>%
kable_classic()| Beta | Beta’s SE | Std. Beta | Std. Beta’s SE | t-value | p-value | ||
|---|---|---|---|---|---|---|---|
| 1 | (Intercept) | 93.279 | 205.575 | 0.074 | 0.189 | 0.454 | 0.651 |
| 2 | negative_affect | -19.738 | 25.580 | -0.289 | 0.374 | -0.772 | 0.442 |
| 3 | extraversion | 22.453 | 24.730 | 0.399 | 0.439 | 0.908 | 0.365 |
| 4 | effortful_control | 14.919 | 51.209 | 0.246 | 0.844 | 0.291 | 0.771 |
| 5 | negative_affect:effortful_control | 3.871 | 6.496 | 0.293 | 0.491 | 0.596 | 0.552 |
| 6 | extraversion:effortful_control | -6.973 | 6.168 | -0.774 | 0.685 | -1.131 | 0.260 |
| 7 | Rspan_score | -0.595 | 0.523 | -0.091 | 0.080 | -1.136 | 0.258 |
| 8 | cond.CRCR2 | 14.427 | 11.598 | 0.309 | 0.248 | 1.244 | 0.215 |
| 9 | cond.CRCR3 | -18.299 | 12.087 | -0.392 | 0.259 | -1.514 | 0.132 |
| 10 | cond.CRCR4 | -3.644 | 12.706 | -0.078 | 0.272 | -0.287 | 0.775 |
| 11 | cond.CRCR5 | 10.723 | 13.846 | 0.230 | 0.297 | 0.774 | 0.440 |
| 12 | cond.CRCR6 | -5.194 | 12.072 | -0.111 | 0.259 | -0.430 | 0.668 |
| 13 | cat_orderII_first | -6.455 | 7.383 | -0.138 | 0.158 | -0.874 | 0.383 |
| 14 | f_stat | 12.000 | 154.000 | 1.474 | 0.140 | ||
| 15 | R2/Adj.R2 | 0.103 | 0.033 | ||||
| 16 | N | 167.000 | |||||
| Note: | |||||||
| … |
write.csv(learner_TTC_CR_mod1_t,"learner_TTC_CR_mod1_t.csv", row.names = T)# reset NA vals in kable tables to "NA"
options(knitr.kable.NA = NA)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
# plot lm
# plot standardized betas to compare cooefficients
ttc_cr_learner_plot <-
sjPlot::plot_model(model = m_cr1_stan, # lm output
title = "", # to remove title
# type = "std", # type "std" = standardized betas
#sort.est = T, # sort based on betas (highest on top)
# wrap.title = 55, # char length of each title line
#axis.lim = c(-.25,.25)
# grid.breaks = .25, # x-axis breaks
#se = T # use se instead of CI
colors = reg_colours,
show.values = T, # show betas
show.p = T, # show * for sig p
p.threshold = c(.049, .009, .0009), # state up to 3 cut-offs for pvals, the number is inclusive**
value.offset = 0.4, # offset to beta labels
value.size = 3, # size of beta labels
digits = 3, # num decimals in beta labels
dot.size = 2.5, line.size = 1,
vline.color = "#F98B99",
# group.terms = c(1, 2, 3, 3, 3, 3, 3, 4, 5, 6, 7, 7, 7),
axis.labels = c("Cat Order : II first", "CR Condition 6", "CR Condition 5","CR Condition 4","CR Condition 3","CR Condition 2", "Working Memory", "Extraversion x Effortful Control", "Negative Affect x Effortful Control", "Effortful Control", "Extraversion", "Negative Affect"),
axis.title = "",
width = 0.2, transform = NULL
) +
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(color = "black"),
axis.text = element_text(colour = "black", size = 11),
axis.ticks.y = element_blank(),
)
# show plot
ttc_cr_learner_plot
Scatterplot of correlational relationship between each predictor and DV
library(tidyr)
predictor_corr_data<-gather(dat[dat$first_TTC.CR<200,], predictor, Score, c(7:9,14), factor_key = TRUE)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
predictor_names <- list(
'negative_affect'="Negative Affect",
'effortful_control'="Effortful Control",
'extraversion'="Extraversion",
'Rspan_score'="Working Memory"
)
predictor_labeller <- function(variable,value){
return(predictor_names[value])
}
predictor_scatter <-
ggplot(predictor_corr_data, aes(y = first_TTC.CR, x = Score, colour = factor(predictor))) +
geom_point(shape = 1) +
scale_color_manual(values = reg_colours[c(2,2,2,2)]) +
facet_wrap(~ predictor, labeller=predictor_labeller, ncol = 2, scales = "free_x") +
geom_smooth(method = "lm") +
xlab("Predictor Value")+
ylab("Trials to Criterion") +
# geom_line(intercept = 83.731, slopr = .5)+
# add theme
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(size=2),
axis.text.x = element_text(angle=0, hjust=2)) +
theme_classic() +
theme(
legend.position = "none",
text = element_text(color = "black"),
axis.text = element_text(colour = "black"),
axis.ticks.y = element_blank(),
strip.text.x = element_text(size=12),
axis.title = element_text(size=15)
)## Warning: The labeller API has been updated. Labellers taking `variable` and
## `value` arguments are now deprecated. See labellers documentation.# show fig
predictor_scatter ## `geom_smooth()` using formula 'y ~ x'
step_mcr1<-step(m_cr1)
step_mcr1_stan<- lm(scale(first_TTC.CR)~scale(extraversion)+scale(effortful_control)+scale(extraversion*effortful_control), dat[dat$first_TTC.CR<200,])summary(step_mcr1)##
## Call:
## lm(formula = first_TTC.CR ~ extraversion + effortful_control +
## extraversion:effortful_control, data = dat[dat$first_TTC.CR <
## 200, ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -68.33 -34.48 -13.22 27.24 123.30
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -119.875 110.422 -1.086 0.2793
## extraversion 39.825 21.886 1.820 0.0707 .
## effortful_control 55.534 27.751 2.001 0.0470 *
## extraversion:effortful_control -11.365 5.481 -2.073 0.0397 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 46.32 on 163 degrees of freedom
## Multiple R-squared: 0.03285, Adjusted R-squared: 0.01505
## F-statistic: 1.846 on 3 and 163 DF, p-value: 0.1409summary(step_mcr1_stan)##
## Call:
## lm(formula = scale(first_TTC.CR) ~ scale(extraversion) + scale(effortful_control) +
## scale(extraversion * effortful_control), data = dat[dat$first_TTC.CR <
## 200, ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.4641 -0.7387 -0.2832 0.5836 2.6418
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -9.366e-17 7.680e-02 0.000 1.0000
## scale(extraversion) 7.071e-01 3.886e-01 1.820 0.0707
## scale(effortful_control) 9.156e-01 4.575e-01 2.001 0.0470
## scale(extraversion * effortful_control) -1.262e+00 6.084e-01 -2.073 0.0397
##
## (Intercept)
## scale(extraversion) .
## scale(effortful_control) *
## scale(extraversion * effortful_control) *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9924 on 163 degrees of freedom
## Multiple R-squared: 0.03285, Adjusted R-squared: 0.01505
## F-statistic: 1.846 on 3 and 163 DF, p-value: 0.1409Regression stats table in R
library(broom)
library(DescTools)
m_cr_ttc_red_coef <- cbind(tidy(step_mcr1)[,1:3],
tidy(step_mcr1_stan)[,2:3],
tidy(step_mcr1)[,4:5])
# remove rownames
rownames(m_cr_ttc_red_coef) <- c(1:nrow(m_cr_ttc_red_coef))
# rename cols
colnames(m_cr_ttc_red_coef) <- c("term", "beta", "beta_SE", "std_beta", "std_beta_SE", "t_val", "p_val")
# finally, create df
learner_TTC_CR_red_mod1_t <-
rbind(m_cr_ttc_red_coef,
tibble("term" = c("f_stat", "R2/Adj.R2", "N"),
"beta" = c(summary(step_mcr1)$fstatistic[2], summary(step_mcr1)$r.squared, nrow(step_mcr1$model)),
"beta_SE" = c(summary(step_mcr1)$fstatistic[3], summary(step_mcr1)$adj.r.squared, NA),
"std_beta" = c(summary(step_mcr1)$fstatistic[1], rep(NA, 2)),
"std_beta_SE" = c(pf(summary(step_mcr1)$fstatistic[1], summary(step_mcr1)$fstatistic[2], summary(step_mcr1)$fstatistic[3], lower.tail = F), rep(NA, 2)),
"t_val" = rep(NA, 3),
"p_val" = rep(NA, 3)
)
)
# remove NA from kable tables
options(knitr.kable.NA = '')
# show kable table
kable(learner_TTC_CR_red_mod1_t, caption = "Multiple Regression Model Predicting TTC in the CR task", align = rep('lcccc'), digits = 3, row.names = T, col.names = c("", "Beta", "Beta's SE", "Std. Beta", "Std. Beta's SE", "t-value", "p-value")) %>%
kable_paper(full_width = F) %>%
# add_header_above(c(" ", "mturk" = 2, "Mturk" = 2)) %>%
# pack_rows("Criterion(s)", 1, 6) %>%
# # column_spec(1, width = "10em") %>%
# pack_rows("Predictor(s)", 7, 19) %>%
footnote(general = "...") %>%
column_spec(c(4), border_left = T) %>%
kable_classic()| Beta | Beta’s SE | Std. Beta | Std. Beta’s SE | t-value | p-value | ||
|---|---|---|---|---|---|---|---|
| 1 | (Intercept) | -119.875 | 110.422 | 0.000 | 0.077 | -1.086 | 0.279 |
| 2 | extraversion | 39.825 | 21.886 | 0.707 | 0.389 | 1.820 | 0.071 |
| 3 | effortful_control | 55.534 | 27.751 | 0.916 | 0.458 | 2.001 | 0.047 |
| 4 | extraversion:effortful_control | -11.365 | 5.481 | -1.262 | 0.608 | -2.073 | 0.040 |
| 5 | f_stat | 3.000 | 163.000 | 1.846 | 0.141 | ||
| 6 | R2/Adj.R2 | 0.033 | 0.015 | ||||
| 7 | N | 167.000 | |||||
| Note: | |||||||
| … |
write.csv(learner_TTC_CR_red_mod1_t,"learner_TTC_CR_red_mod1_t.csv", row.names = T)# reset NA vals in kable tables to "NA"
options(knitr.kable.NA = NA)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
# plot lm
# plot standardized betas to compare cooefficients
ttc_cr_learner_red_plot <-
sjPlot::plot_model(model = step_mcr1_stan, # lm output
title = "", # to remove title
# type = "std", # type "std" = standardized betas
#sort.est = T, # sort based on betas (highest on top)
# wrap.title = 55, # char length of each title line
#axis.lim = c(-.25,.25)
# grid.breaks = .25, # x-axis breaks
#se = T # use se instead of CI
colors = reg_colours,
show.values = T, # show betas
show.p = T, # show * for sig p
p.threshold = c(.049, .009, .0009), # state up to 3 cut-offs for pvals, the number is inclusive**
value.offset = 0.4, # offset to beta labels
value.size = 3, # size of beta labels
digits = 3, # num decimals in beta labels
dot.size = 2.5, line.size = 1,
vline.color = "#F98B99",
# group.terms = c(1, 2, 3, 3, 3, 3, 3, 4, 5, 6, 7, 7, 7),
axis.labels = c("Extraversion x Effortful Control","Effortful Control","Extraversion"),
axis.title = "",
width = 0.2, transform = NULL
) +
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(color = "black"),
axis.text = element_text(colour = "black", size = 12),
axis.ticks.y = element_blank(),
)
# show plot
ttc_cr_learner_red_plot
t.test(dat$extraversion~dat$group2)##
## Welch Two Sample t-test
##
## data: dat$extraversion by dat$group2
## t = -0.82717, df = 118.11, p-value = 0.4098
## alternative hypothesis: true difference in means between group cr_learner and group cr_nlearner is not equal to 0
## 95 percent confidence interval:
## -0.3449043 0.1416612
## sample estimates:
## mean in group cr_learner mean in group cr_nlearner
## 4.971024 5.072646t.test(dat$negative_affect~dat$group2)##
## Welch Two Sample t-test
##
## data: dat$negative_affect by dat$group2
## t = 0.010737, df = 108.86, p-value = 0.9915
## alternative hypothesis: true difference in means between group cr_learner and group cr_nlearner is not equal to 0
## 95 percent confidence interval:
## -0.2145995 0.2169373
## sample estimates:
## mean in group cr_learner mean in group cr_nlearner
## 4.539875 4.538706t.test(dat$effortful_control~dat$group2)##
## Welch Two Sample t-test
##
## data: dat$effortful_control by dat$group2
## t = -0.50125, df = 119.38, p-value = 0.6171
## alternative hypothesis: true difference in means between group cr_learner and group cr_nlearner is not equal to 0
## 95 percent confidence interval:
## -0.2805634 0.1672093
## sample estimates:
## mean in group cr_learner mean in group cr_nlearner
## 3.861330 3.918007t.test(dat$Rspan_score~dat$group2)##
## Welch Two Sample t-test
##
## data: dat$Rspan_score by dat$group2
## t = 1.5193, df = 102.79, p-value = 0.1318
## alternative hypothesis: true difference in means between group cr_learner and group cr_nlearner is not equal to 0
## 95 percent confidence interval:
## -0.5614353 4.2377245
## sample estimates:
## mean in group cr_learner mean in group cr_nlearner
## 47.77844 45.94030II task
Overall performance
m_ii_learner<-d_wide[!is.na(d_wide$first_TTC.II),]
# find and remove bad data
dat1<-DATA[DATA$exp_isCorrect>0.80 & DATA$number_answered<100 & DATA$intrusion<60& (DATA$performance_block_1_II<DATA$performance_block_4_II|DATA$max.II>0.599)&DATA$duration_seconds>300,]
dat1$group2 <- ifelse(dat1$first_TTC.II==200,"ii_nlearner","ii_learner")Check correlation matrix among predictors
ii.corr<- rcorr(as.matrix(dat1[,c(7:9, 14)]))
ii.r<-ii.corr$r
ii.p<-ii.corr$P
ii.corr$r## negative_affect effortful_control extraversion Rspan_score
## negative_affect 1.0000000 -0.47376482 -0.36693610 0.0149425
## effortful_control -0.4737648 1.00000000 0.09243824 -0.1189793
## extraversion -0.3669361 0.09243824 1.00000000 0.1161389
## Rspan_score 0.0149425 -0.11897932 0.11613891 1.0000000ii.p## negative_affect effortful_control extraversion Rspan_score
## negative_affect NA 6.217249e-14 1.516122e-08 0.82400244
## effortful_control 6.217249e-14 NA 1.679866e-01 0.07555600
## extraversion 1.516122e-08 1.679866e-01 NA 0.08285476
## Rspan_score 8.240024e-01 7.555600e-02 8.285476e-02 NAhist(dat1$average.II, col='yellow', xlab='Average II Task Performance', main="", labels=T, ylim=c(0,50))
hist(dat1$first_TTC.II, col='yellow', xlab='Trails to Criterion', main="", labels = T, ylim=c(0,120), xlim=c(0,200))
m_ii<-lm(average.II ~ negative_affect+extraversion+effortful_control+negative_affect*effortful_control+extraversion*effortful_control+Rspan_score+cond.II+cat_order, dat1)
m_ii_stan<- lm(scale(average.II)~scale(negative_affect)+scale(extraversion)+scale(effortful_control)+scale(negative_affect*effortful_control)+scale(extraversion*effortful_control)+scale(Rspan_score)+ cond.II + cat_order, dat1)
summary(m_ii)##
## Call:
## lm(formula = average.II ~ negative_affect + extraversion + effortful_control +
## negative_affect * effortful_control + extraversion * effortful_control +
## Rspan_score + cond.II + cat_order, data = dat1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.20713 -0.07389 -0.01456 0.05036 0.31874
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.5680389 0.4053035 3.869 0.000145 ***
## negative_affect -0.1077590 0.0470301 -2.291 0.022925 *
## extraversion -0.0966423 0.0490720 -1.969 0.050204 .
## effortful_control -0.2191160 0.0946734 -2.314 0.021596 *
## Rspan_score 0.0005461 0.0010359 0.527 0.598611
## cond.IIII2 0.0038787 0.0202534 0.192 0.848311
## cond.IIII3 0.0239685 0.0216430 1.107 0.269349
## cond.IIII4 -0.0050570 0.0210230 -0.241 0.810140
## cat_orderII_first -0.0253353 0.0150793 -1.680 0.094395 .
## negative_affect:effortful_control 0.0270300 0.0112759 2.397 0.017387 *
## extraversion:effortful_control 0.0205774 0.0120214 1.712 0.088404 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1085 on 213 degrees of freedom
## Multiple R-squared: 0.0612, Adjusted R-squared: 0.01713
## F-statistic: 1.389 on 10 and 213 DF, p-value: 0.187summary(m_ii_stan)##
## Call:
## lm(formula = scale(average.II) ~ scale(negative_affect) + scale(extraversion) +
## scale(effortful_control) + scale(negative_affect * effortful_control) +
## scale(extraversion * effortful_control) + scale(Rspan_score) +
## cond.II + cat_order, data = dat1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.8927 -0.6752 -0.1331 0.4602 2.9125
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.06474 0.14542 0.445 0.6566
## scale(negative_affect) -0.71903 0.31381 -2.291 0.0229
## scale(extraversion) -0.74514 0.37836 -1.969 0.0502
## scale(effortful_control) -1.55460 0.67170 -2.314 0.0216
## scale(negative_affect * effortful_control) 0.82908 0.34586 2.397 0.0174
## scale(extraversion * effortful_control) 1.02069 0.59630 1.712 0.0884
## scale(Rspan_score) 0.03627 0.06880 0.527 0.5986
## cond.IIII2 0.03544 0.18507 0.192 0.8483
## cond.IIII3 0.21901 0.19776 1.107 0.2693
## cond.IIII4 -0.04621 0.19210 -0.241 0.8101
## cat_orderII_first -0.23150 0.13779 -1.680 0.0944
##
## (Intercept)
## scale(negative_affect) *
## scale(extraversion) .
## scale(effortful_control) *
## scale(negative_affect * effortful_control) *
## scale(extraversion * effortful_control) .
## scale(Rspan_score)
## cond.IIII2
## cond.IIII3
## cond.IIII4
## cat_orderII_first .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9914 on 213 degrees of freedom
## Multiple R-squared: 0.0612, Adjusted R-squared: 0.01713
## F-statistic: 1.389 on 10 and 213 DF, p-value: 0.187# create table depicting main stats to export later
# use tidy() in broom package & DescTools::StdCoef()
library(broom)
library(DescTools)
m_ii_reg_coef <- cbind(tidy(m_ii)[c(1:4, 10:11,5, 6:9),1:3],
tidy(m_ii_stan)[,2:3],
tidy(m_ii)[c(1:4, 10:11,5, 6:9),4:5])
# remove rownames
rownames(m_ii_reg_coef) <- c(1:nrow(m_ii_reg_coef))
# rename cols
colnames(m_ii_reg_coef) <- c("term", "beta", "beta_SE", "std_beta", "std_beta_SE", "t_val", "p_val")
# finally, create df
reg_TTC_II_mod1_t <-
rbind(m_ii_reg_coef,
tibble("term" = c("f_stat", "R2/Adj.R2", "N"),
"beta" = c(summary(m_ii)$fstatistic[2], summary(m_ii)$r.squared, nrow(m_ii$model)),
"beta_SE" = c(summary(m_ii)$fstatistic[3], summary(m_ii)$adj.r.squared, NA),
"std_beta" = c(summary(m_ii)$fstatistic[1], rep(NA, 2)),
"std_beta_SE" = c(pf(summary(m_ii)$fstatistic[1], summary(m_ii)$fstatistic[2], summary(m_ii)$fstatistic[3], lower.tail = F), rep(NA, 2)),
"t_val" = rep(NA, 3),
"p_val" = rep(NA, 3)
)
)
# remove NA from kable tables
options(knitr.kable.NA = '')
# show kable table
kable(reg_TTC_II_mod1_t, caption = "Multiple Regression Model Predicting II overall average performance", align = rep('lcccc'), digits = 3, row.names = T, col.names = c("", "Beta", "Beta's SE", "Std. Beta", "Std. Beta's SE", "t-value", "p-value")) %>%
kable_paper(full_width = F) %>%
# add_header_above(c(" ", "mturk" = 2, "Mturk" = 2)) %>%
# pack_rows("Criterion(s)", 1, 6) %>%
# # column_spec(1, width = "10em") %>%
# pack_rows("Predictor(s)", 7, 19) %>%
footnote(general = "...") %>%
column_spec(c(4), border_left = T) %>%
kable_classic()| Beta | Beta’s SE | Std. Beta | Std. Beta’s SE | t-value | p-value | ||
|---|---|---|---|---|---|---|---|
| 1 | (Intercept) | 1.568 | 0.405 | 0.065 | 0.145 | 3.869 | 0.000 |
| 2 | negative_affect | -0.108 | 0.047 | -0.719 | 0.314 | -2.291 | 0.023 |
| 3 | extraversion | -0.097 | 0.049 | -0.745 | 0.378 | -1.969 | 0.050 |
| 4 | effortful_control | -0.219 | 0.095 | -1.555 | 0.672 | -2.314 | 0.022 |
| 5 | negative_affect:effortful_control | 0.027 | 0.011 | 0.829 | 0.346 | 2.397 | 0.017 |
| 6 | extraversion:effortful_control | 0.021 | 0.012 | 1.021 | 0.596 | 1.712 | 0.088 |
| 7 | Rspan_score | 0.001 | 0.001 | 0.036 | 0.069 | 0.527 | 0.599 |
| 8 | cond.IIII2 | 0.004 | 0.020 | 0.035 | 0.185 | 0.192 | 0.848 |
| 9 | cond.IIII3 | 0.024 | 0.022 | 0.219 | 0.198 | 1.107 | 0.269 |
| 10 | cond.IIII4 | -0.005 | 0.021 | -0.046 | 0.192 | -0.241 | 0.810 |
| 11 | cat_orderII_first | -0.025 | 0.015 | -0.232 | 0.138 | -1.680 | 0.094 |
| 12 | f_stat | 10.000 | 213.000 | 1.389 | 0.187 | ||
| 13 | R2/Adj.R2 | 0.061 | 0.017 | ||||
| 14 | N | 224.000 | |||||
| Note: | |||||||
| … |
write.csv(reg_TTC_II_mod1_t,"reg_TTC_II_mod1_t.csv", row.names = T)Scatterplot of correlational relationship between each predictor and DV
library(tidyr)
predictor_corr_data<-gather(dat1, predictor, Score, c(7:9,14), factor_key = TRUE)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
predictor_corr_data<-predictor_corr_data[,c(19,41,42)]
predictor_names <- list(
'negative_affect'="Negative Affect",
'effortful_control'="Effortful Control",
'extraversion'="Extraversion",
'Rspan_score'="Working Memory"
)
predictor_labeller <- function(variable,value){
return(predictor_names[value])
}
predictor_scatter <-
ggplot(predictor_corr_data, aes(y = average.II, x = Score, colour = factor(predictor))) +
geom_point(shape = 1) +
scale_color_manual(values = reg_colours[c(2,2,2,2)]) +
facet_wrap(predictor ~ ., labeller=predictor_labeller, scales = "free_x")+
# facet_wrap(~ predictor, ncol = 2, scales = "free_x") +
geom_smooth(method = "lm") +
xlab("Predictor Value")+
ylab("Overall Average Performance") +
# geom_line(intercept = 83.731, slopr = .5)+
# add theme
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(size=20),
axis.text.x = element_text(angle=0, hjust=2)) +
theme_classic() +
theme(
legend.position = "none",
text = element_text(color = "black"),
axis.text = element_text(colour = "black"),
axis.ticks.y = element_blank(),
strip.text.x = element_text(size=12),
axis.title = element_text(size=12)
)## Warning: The labeller API has been updated. Labellers taking `variable` and
## `value` arguments are now deprecated. See labellers documentation.# show fig
predictor_scatter## `geom_smooth()` using formula 'y ~ x'
# reset NA vals in kable tables to "NA"
options(knitr.kable.NA = NA)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
# plot lm
# plot standardized betas to compare cooefficients
reg_ii_plot <-
sjPlot::plot_model(model = m_ii_stan, # lm output
title = "", # to remove title
# type = "std", # type "std" = standardized betas
#sort.est = T, # sort based on betas (highest on top)
# wrap.title = 55, # char length of each title line
#axis.lim = c(-.25,.25)
# grid.breaks = .25, # x-axis breaks
#se = T # use se instead of CI
colors = reg_colours,
show.values = T, # show betas
show.p = T, # show * for sig p
p.threshold = c(.049, .009, .0009), # state up to 3 cut-offs for pvals, the number is inclusive**
value.offset = 0.4, # offset to beta labels
value.size = 3, # size of beta labels
digits = 3, # num decimals in beta labels
dot.size = 2.5, line.size = 1,
vline.color = "#F98B99",
# group.terms = c(1, 2, 3, 3, 3, 3, 3, 4, 5, 6, 7, 7, 7),
axis.labels = c("Cat Order: II First", "II Condition 4", "II Condition 3","II Condition 2","Working Memory", "Extraversion x Effortful Control", "Negative Affect x Effortful Control", "Effortful Control", "Extraversion", "Negative Affect"),
axis.title = "",
width = 0.2, transform = NULL
) +
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(color = "black"),
axis.text = element_text(colour = "black", size = 11),
axis.ticks.y = element_blank(),
)
# show plot
reg_ii_plot
step_mii<-step(m_ii)
step_mii_stan<-step(m_ii_stan)summary(step_mii)##
## Call:
## lm(formula = average.II ~ negative_affect + extraversion + effortful_control +
## cat_order + negative_affect:effortful_control + extraversion:effortful_control,
## data = dat1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.22445 -0.07445 -0.01492 0.05527 0.33401
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.60190 0.39416 4.064 6.74e-05 ***
## negative_affect -0.10549 0.04639 -2.274 0.0239 *
## extraversion -0.10030 0.04880 -2.055 0.0411 *
## effortful_control -0.22339 0.09370 -2.384 0.0180 *
## cat_orderII_first -0.02554 0.01457 -1.753 0.0810 .
## negative_affect:effortful_control 0.02670 0.01112 2.402 0.0171 *
## extraversion:effortful_control 0.02200 0.01193 1.844 0.0665 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1081 on 217 degrees of freedom
## Multiple R-squared: 0.05068, Adjusted R-squared: 0.02443
## F-statistic: 1.931 on 6 and 217 DF, p-value: 0.07706summary(step_mii_stan)##
## Call:
## lm(formula = scale(average.II) ~ scale(negative_affect) + scale(extraversion) +
## scale(effortful_control) + scale(negative_affect * effortful_control) +
## scale(extraversion * effortful_control) + cat_order, data = dat1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.0509 -0.6803 -0.1364 0.5051 3.0520
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.11459 0.09289 1.234 0.2187
## scale(negative_affect) -0.70388 0.30954 -2.274 0.0239
## scale(extraversion) -0.77334 0.37630 -2.055 0.0411
## scale(effortful_control) -1.58491 0.66480 -2.384 0.0180
## scale(negative_affect * effortful_control) 0.81899 0.34095 2.402 0.0171
## scale(extraversion * effortful_control) 1.09115 0.59172 1.844 0.0665
## cat_orderII_first -0.23334 0.13311 -1.753 0.0810
##
## (Intercept)
## scale(negative_affect) *
## scale(extraversion) *
## scale(effortful_control) *
## scale(negative_affect * effortful_control) *
## scale(extraversion * effortful_control) .
## cat_orderII_first .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9877 on 217 degrees of freedom
## Multiple R-squared: 0.05068, Adjusted R-squared: 0.02443
## F-statistic: 1.931 on 6 and 217 DF, p-value: 0.07706# create table depicting main stats to export later
# use tidy() in broom package & DescTools::StdCoef()
library(broom)
library(DescTools)
red_m_ii_reg_coef <- cbind(tidy(step_mii)[c(1:4, 6,7,5),1:3],
tidy(step_mii_stan)[,2:3],
tidy(step_mii)[c(1:4,6,7,5),4:5])
# remove rownames
rownames(red_m_ii_reg_coef) <- c(1:nrow(red_m_ii_reg_coef))
# rename cols
colnames(red_m_ii_reg_coef) <- c("term", "beta", "beta_SE", "std_beta", "std_beta_SE", "t_val", "p_val")
# finally, create df
red_reg_TTC_II_mod1_t <-
rbind(red_m_ii_reg_coef,
tibble("term" = c("f_stat", "R2/Adj.R2", "N"),
"beta" = c(summary(step_mii)$fstatistic[2], summary(step_mii)$r.squared, nrow(step_mii$model)),
"beta_SE" = c(summary(step_mii)$fstatistic[3], summary(step_mii)$adj.r.squared, NA),
"std_beta" = c(summary(step_mii)$fstatistic[1], rep(NA, 2)),
"std_beta_SE" = c(pf(summary(step_mii)$fstatistic[1], summary(step_mii)$fstatistic[2], summary(step_mii)$fstatistic[3], lower.tail = F), rep(NA, 2)),
"t_val" = rep(NA, 3),
"p_val" = rep(NA, 3)
)
)
# remove NA from kable tables
options(knitr.kable.NA = '')
# show kable table
kable(red_reg_TTC_II_mod1_t, caption = "Multiple Regression Model Predicting II overall average performance", align = rep('lcccc'), digits = 3, row.names = T, col.names = c("", "Beta", "Beta's SE", "Std. Beta", "Std. Beta's SE", "t-value", "p-value")) %>%
kable_paper(full_width = F) %>%
# add_header_above(c(" ", "mturk" = 2, "Mturk" = 2)) %>%
# pack_rows("Criterion(s)", 1, 6) %>%
# # column_spec(1, width = "10em") %>%
# pack_rows("Predictor(s)", 7, 19) %>%
footnote(general = "...") %>%
column_spec(c(4), border_left = T) %>%
kable_classic()| Beta | Beta’s SE | Std. Beta | Std. Beta’s SE | t-value | p-value | ||
|---|---|---|---|---|---|---|---|
| 1 | (Intercept) | 1.602 | 0.394 | 0.115 | 0.093 | 4.064 | 0.000 |
| 2 | negative_affect | -0.105 | 0.046 | -0.704 | 0.310 | -2.274 | 0.024 |
| 3 | extraversion | -0.100 | 0.049 | -0.773 | 0.376 | -2.055 | 0.041 |
| 4 | effortful_control | -0.223 | 0.094 | -1.585 | 0.665 | -2.384 | 0.018 |
| 5 | negative_affect:effortful_control | 0.027 | 0.011 | 0.819 | 0.341 | 2.402 | 0.017 |
| 6 | extraversion:effortful_control | 0.022 | 0.012 | 1.091 | 0.592 | 1.844 | 0.067 |
| 7 | cat_orderII_first | -0.026 | 0.015 | -0.233 | 0.133 | -1.753 | 0.081 |
| 8 | f_stat | 6.000 | 217.000 | 1.931 | 0.077 | ||
| 9 | R2/Adj.R2 | 0.051 | 0.024 | ||||
| 10 | N | 224.000 | |||||
| Note: | |||||||
| … |
write.csv(red_reg_TTC_II_mod1_t,"red_reg_TTC_II_mod1_t.csv", row.names = T)# reset NA vals in kable tables to "NA"
options(knitr.kable.NA = NA)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
# plot lm
# plot standardized betas to compare cooefficients
red_reg_ii_plot <-
sjPlot::plot_model(model = step_mii_stan, # lm output
title = "", # to remove title
# type = "std", # type "std" = standardized betas
#sort.est = T, # sort based on betas (highest on top)
# wrap.title = 55, # char length of each title line
#axis.lim = c(-.25,.25)
# grid.breaks = .25, # x-axis breaks
#se = T # use se instead of CI
colors = reg_colours,
show.values = T, # show betas
show.p = T, # show * for sig p
p.threshold = c(.049, .009, .0009), # state up to 3 cut-offs for pvals, the number is inclusive**
value.offset = 0.4, # offset to beta labels
value.size = 3, # size of beta labels
digits = 3, # num decimals in beta labels
dot.size = 2.5, line.size = 1,
vline.color = "#F98B99",
# group.terms = c(1, 2, 3, 3, 3, 3, 3, 4, 5, 6, 7, 7, 7),
axis.labels = c("Cat Order: II First","Extraversion x Effortful Control", "Negative Affect x Effortful Control", "Effortful Control", "Extraversion", "Negative Affect"),
axis.title = "",
width = 0.2, transform = NULL
) +
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(color = "black"),
axis.text = element_text(colour = "black", size = 12),
axis.ticks.y = element_blank(),
)
# show plot
red_reg_ii_plot
TTC (learners)
Scatterplot of correlational relationship between each predictor and DV
library(tidyr)
predictor_corr_data<-gather(dat1[dat1$first_TTC.II<200,], predictor, Score, c(7:9,14), factor_key = TRUE)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
predictor_names <- list(
'negative_affect'="Negative Affect",
'effortful_control'="Effortful Control",
'extraversion'="Extraversion",
'Rspan_score'="Working Memory"
)
predictor_labeller <- function(variable,value){
return(predictor_names[value])
}
predictor_scatter <-
ggplot(predictor_corr_data, aes(y = first_TTC.II, x = Score, colour = factor(predictor))) +
geom_point(shape = 1) +
scale_color_manual(values = reg_colours[c(2,2,2,2)]) +
facet_wrap(~ predictor, labeller=predictor_labeller, ncol = 2, scales = "free_x") +
geom_smooth(method = "lm") +
xlab("Predictor Value")+
ylab("Trials to Criterion") +
# geom_line(intercept = 83.731, slopr = .5)+
# add theme
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(size=2),
axis.text.x = element_text(angle=0, hjust=2)) +
theme_classic() +
theme(
legend.position = "none",
text = element_text(color = "black"),
axis.text = element_text(colour = "black"),
axis.ticks.y = element_blank(),
strip.text.x = element_text(size=12),
axis.title = element_text(size=15)
)## Warning: The labeller API has been updated. Labellers taking `variable` and
## `value` arguments are now deprecated. See labellers documentation.# show fig
predictor_scatter ## `geom_smooth()` using formula 'y ~ x'
m_ii1<-lm(first_TTC.II ~ negative_affect+extraversion+effortful_control+negative_affect*effortful_control+extraversion*effortful_control+Rspan_score+cond.II+cat_order, dat1[dat1$first_TTC.II<200,])
m_ii1_stan<- lm(scale(first_TTC.II)~scale(negative_affect)+scale(extraversion)+scale(effortful_control)+scale(negative_affect*effortful_control)+scale(extraversion*effortful_control)+scale(Rspan_score)+ cond.II + cat_order, dat1[dat1$first_TTC.II<200,])
summary(m_ii1)##
## Call:
## lm(formula = first_TTC.II ~ negative_affect + extraversion +
## effortful_control + negative_affect * effortful_control +
## extraversion * effortful_control + Rspan_score + cond.II +
## cat_order, data = dat1[dat1$first_TTC.II < 200, ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -79.309 -41.293 -4.436 37.242 116.307
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -262.5186 296.3043 -0.886 0.3775
## negative_affect 27.8263 34.7790 0.800 0.4254
## extraversion 59.1527 34.8931 1.695 0.0928 .
## effortful_control 117.5066 73.3542 1.602 0.1120
## Rspan_score -0.4322 0.6923 -0.624 0.5337
## cond.IIII2 -15.5093 12.5444 -1.236 0.2189
## cond.IIII3 -38.0886 13.3909 -2.844 0.0053 **
## cond.IIII4 -27.3858 14.0118 -1.954 0.0532 .
## cat_orderII_first 4.7749 9.4736 0.504 0.6152
## negative_affect:effortful_control -9.4683 8.5824 -1.103 0.2723
## extraversion:effortful_control -15.8562 8.7278 -1.817 0.0720 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 50.07 on 111 degrees of freedom
## Multiple R-squared: 0.113, Adjusted R-squared: 0.03313
## F-statistic: 1.415 on 10 and 111 DF, p-value: 0.183summary(m_ii1_stan)##
## Call:
## lm(formula = scale(first_TTC.II) ~ scale(negative_affect) + scale(extraversion) +
## scale(effortful_control) + scale(negative_affect * effortful_control) +
## scale(extraversion * effortful_control) + scale(Rspan_score) +
## cond.II + cat_order, data = dat1[dat1$first_TTC.II < 200,
## ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.55744 -0.81089 -0.08711 0.73135 2.28400
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.34824 0.20962 1.661 0.0995
## scale(negative_affect) 0.38749 0.48431 0.800 0.4254
## scale(extraversion) 0.98068 0.57848 1.695 0.0928
## scale(effortful_control) 1.75332 1.09452 1.602 0.1120
## scale(negative_affect * effortful_control) -0.68888 0.62443 -1.103 0.2723
## scale(extraversion * effortful_control) -1.68222 0.92595 -1.817 0.0720
## scale(Rspan_score) -0.05937 0.09509 -0.624 0.5337
## cond.IIII2 -0.30457 0.24634 -1.236 0.2189
## cond.IIII3 -0.74797 0.26297 -2.844 0.0053
## cond.IIII4 -0.53779 0.27516 -1.954 0.0532
## cat_orderII_first 0.09377 0.18604 0.504 0.6152
##
## (Intercept) .
## scale(negative_affect)
## scale(extraversion) .
## scale(effortful_control)
## scale(negative_affect * effortful_control)
## scale(extraversion * effortful_control) .
## scale(Rspan_score)
## cond.IIII2
## cond.IIII3 **
## cond.IIII4 .
## cat_orderII_first
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9833 on 111 degrees of freedom
## Multiple R-squared: 0.113, Adjusted R-squared: 0.03313
## F-statistic: 1.415 on 10 and 111 DF, p-value: 0.183# create table depicting main stats to export later
# use tidy() in broom package & DescTools::StdCoef()
library(broom)
library(DescTools)
TTC_ii_reg_coef <- cbind(tidy(m_ii1)[c(1:4, 10:11,5, 6:9),1:3],
tidy(m_ii1_stan)[,2:3],
tidy(m_ii1)[c(1:4, 10:11,5, 6:9),4:5])
# remove rownames
rownames(TTC_ii_reg_coef) <- c(1:nrow(TTC_ii_reg_coef))
# rename cols
colnames(TTC_ii_reg_coef) <- c("term", "beta", "beta_SE", "std_beta", "std_beta_SE", "t_val", "p_val")
# finally, create df
TTC_II_mod1_t <-
rbind(TTC_ii_reg_coef,
tibble("term" = c("f_stat", "R2/Adj.R2", "N"),
"beta" = c(summary(m_ii1)$fstatistic[2], summary(m_ii1)$r.squared, nrow(m_ii$model)),
"beta_SE" = c(summary(m_ii1)$fstatistic[3], summary(m_ii1)$adj.r.squared, NA),
"std_beta" = c(summary(m_ii1)$fstatistic[1], rep(NA, 2)),
"std_beta_SE" = c(pf(summary(m_ii1)$fstatistic[1], summary(m_ii1)$fstatistic[2], summary(m_ii1)$fstatistic[3], lower.tail = F), rep(NA, 2)),
"t_val" = rep(NA, 3),
"p_val" = rep(NA, 3)
)
)
# remove NA from kable tables
options(knitr.kable.NA = '')
# show kable table
kable(TTC_II_mod1_t, caption = "Multiple Regression Model Predicting II overall average performance", align = rep('lcccc'), digits = 3, row.names = T, col.names = c("", "Beta", "Beta's SE", "Std. Beta", "Std. Beta's SE", "t-value", "p-value")) %>%
kable_paper(full_width = F) %>%
# add_header_above(c(" ", "mturk" = 2, "Mturk" = 2)) %>%
# pack_rows("Criterion(s)", 1, 6) %>%
# # column_spec(1, width = "10em") %>%
# pack_rows("Predictor(s)", 7, 19) %>%
footnote(general = "...") %>%
column_spec(c(4), border_left = T) %>%
kable_classic()| Beta | Beta’s SE | Std. Beta | Std. Beta’s SE | t-value | p-value | ||
|---|---|---|---|---|---|---|---|
| 1 | (Intercept) | -262.519 | 296.304 | 0.348 | 0.210 | -0.886 | 0.378 |
| 2 | negative_affect | 27.826 | 34.779 | 0.387 | 0.484 | 0.800 | 0.425 |
| 3 | extraversion | 59.153 | 34.893 | 0.981 | 0.578 | 1.695 | 0.093 |
| 4 | effortful_control | 117.507 | 73.354 | 1.753 | 1.095 | 1.602 | 0.112 |
| 5 | negative_affect:effortful_control | -9.468 | 8.582 | -0.689 | 0.624 | -1.103 | 0.272 |
| 6 | extraversion:effortful_control | -15.856 | 8.728 | -1.682 | 0.926 | -1.817 | 0.072 |
| 7 | Rspan_score | -0.432 | 0.692 | -0.059 | 0.095 | -0.624 | 0.534 |
| 8 | cond.IIII2 | -15.509 | 12.544 | -0.305 | 0.246 | -1.236 | 0.219 |
| 9 | cond.IIII3 | -38.089 | 13.391 | -0.748 | 0.263 | -2.844 | 0.005 |
| 10 | cond.IIII4 | -27.386 | 14.012 | -0.538 | 0.275 | -1.954 | 0.053 |
| 11 | cat_orderII_first | 4.775 | 9.474 | 0.094 | 0.186 | 0.504 | 0.615 |
| 12 | f_stat | 10.000 | 111.000 | 1.415 | 0.183 | ||
| 13 | R2/Adj.R2 | 0.113 | 0.033 | ||||
| 14 | N | 224.000 | |||||
| Note: | |||||||
| … |
write.csv(TTC_II_mod1_t,"TTC_II_mod1_t.csv", row.names = T)# reset NA vals in kable tables to "NA"
options(knitr.kable.NA = NA)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
# plot lm
# plot standardized betas to compare cooefficients
TTC_ii_plot <-
sjPlot::plot_model(model = m_ii1_stan, # lm output
title = "", # to remove title
# type = "std", # type "std" = standardized betas
#sort.est = T, # sort based on betas (highest on top)
# wrap.title = 55, # char length of each title line
#axis.lim = c(-.25,.25)
# grid.breaks = .25, # x-axis breaks
#se = T # use se instead of CI
colors = reg_colours,
show.values = T, # show betas
show.p = T, # show * for sig p
p.threshold = c(.049, .009, .0009), # state up to 3 cut-offs for pvals, the number is inclusive**
value.offset = 0.4, # offset to beta labels
value.size = 3, # size of beta labels
digits = 2, # num decimals in beta labels
dot.size = 2.5, line.size = 1,
vline.color = "#F98B99",
# group.terms = c(1, 2, 3, 3, 3, 3, 3, 4, 5, 6, 7, 7, 7),
axis.labels = c("Cat Order: II First", "II Condition 4", "II Condition 3","II Condition 2","Working Memory", "Extraversion x Effortful Control", "Negative Affect x Effortful Control", "Effortful Control", "Extraversion", "Negative Affect"),
axis.title = "",
width = 0.2, transform = NULL
) +
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(color = "black"),
axis.text = element_text(colour = "black", size = 11),
axis.ticks.y = element_blank(),
)
# show plot
TTC_ii_plot
step_mii1<-step(m_ii1)
step_mii1_stan<-lm(scale(first_TTC.II)~scale(extraversion)+scale(effortful_control)+scale(extraversion*effortful_control)+cond.II, dat1[dat1$first_TTC.II<200,])summary(step_mii1_stan)##
## Call:
## lm(formula = scale(first_TTC.II) ~ scale(extraversion) + scale(effortful_control) +
## scale(extraversion * effortful_control) + cond.II, data = dat1[dat1$first_TTC.II <
## 200, ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.4270 -0.8591 -0.1678 0.6824 2.1583
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.3875 0.1857 2.087 0.03907 *
## scale(extraversion) 0.6572 0.4641 1.416 0.15946
## scale(effortful_control) 0.7487 0.5271 1.420 0.15821
## scale(extraversion * effortful_control) -1.0920 0.7456 -1.465 0.14577
## cond.IIII2 -0.3400 0.2444 -1.391 0.16686
## cond.IIII3 -0.7105 0.2590 -2.743 0.00706 **
## cond.IIII4 -0.4933 0.2703 -1.825 0.07059 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9808 on 115 degrees of freedom
## Multiple R-squared: 0.08566, Adjusted R-squared: 0.03796
## F-statistic: 1.796 on 6 and 115 DF, p-value: 0.106summary(step_mii1)##
## Call:
## lm(formula = first_TTC.II ~ extraversion + effortful_control +
## cond.II + extraversion:effortful_control, data = dat1[dat1$first_TTC.II <
## 200, ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -72.665 -43.746 -8.543 34.747 109.907
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -69.804 139.764 -0.499 0.61842
## extraversion 39.639 27.992 1.416 0.15946
## effortful_control 50.179 35.329 1.420 0.15821
## cond.IIII2 -17.314 12.445 -1.391 0.16686
## cond.IIII3 -36.178 13.189 -2.743 0.00706 **
## cond.IIII4 -25.118 13.763 -1.825 0.07059 .
## extraversion:effortful_control -10.293 7.028 -1.465 0.14577
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 49.95 on 115 degrees of freedom
## Multiple R-squared: 0.08566, Adjusted R-squared: 0.03796
## F-statistic: 1.796 on 6 and 115 DF, p-value: 0.106# create table depicting main stats to export later
# use tidy() in broom package & DescTools::StdCoef()
library(broom)
library(DescTools)
red_TTC_ii_reg_coef <- cbind(tidy(step_mii1)[c(1:3, 7, 4:6),1:3],
tidy(step_mii1_stan)[,2:3],
tidy(step_mii1)[c(1:3, 7, 4:6),4:5])
# remove rownames
rownames(red_TTC_ii_reg_coef) <- c(1:nrow(red_TTC_ii_reg_coef))
# rename cols
colnames(red_TTC_ii_reg_coef) <- c("term", "beta", "beta_SE", "std_beta", "std_beta_SE", "t_val", "p_val")
# finally, create df
red_TTC_II_mod1_t <-
rbind(red_TTC_ii_reg_coef,
tibble("term" = c("f_stat", "R2/Adj.R2", "N"),
"beta" = c(summary(step_mii1)$fstatistic[2], summary(step_mii1)$r.squared, nrow(m_ii$model)),
"beta_SE" = c(summary(step_mii1)$fstatistic[3], summary(step_mii1)$adj.r.squared, NA),
"std_beta" = c(summary(step_mii1)$fstatistic[1], rep(NA, 2)),
"std_beta_SE" = c(pf(summary(step_mii1)$fstatistic[1], summary(step_mii1)$fstatistic[2], summary(step_mii1)$fstatistic[3], lower.tail = F), rep(NA, 2)),
"t_val" = rep(NA, 3),
"p_val" = rep(NA, 3)
)
)
# remove NA from kable tables
options(knitr.kable.NA = '')
# show kable table
kable(red_TTC_II_mod1_t, caption = "Multiple Regression Model Predicting II overall average performance", align = rep('lcccc'), digits = 3, row.names = T, col.names = c("", "Beta", "Beta's SE", "Std. Beta", "Std. Beta's SE", "t-value", "p-value")) %>%
kable_paper(full_width = F) %>%
# add_header_above(c(" ", "mturk" = 2, "Mturk" = 2)) %>%
# pack_rows("Criterion(s)", 1, 6) %>%
# # column_spec(1, width = "10em") %>%
# pack_rows("Predictor(s)", 7, 19) %>%
footnote(general = "...") %>%
column_spec(c(4), border_left = T) %>%
kable_classic()| Beta | Beta’s SE | Std. Beta | Std. Beta’s SE | t-value | p-value | ||
|---|---|---|---|---|---|---|---|
| 1 | (Intercept) | -69.804 | 139.764 | 0.388 | 0.186 | -0.499 | 0.618 |
| 2 | extraversion | 39.639 | 27.992 | 0.657 | 0.464 | 1.416 | 0.159 |
| 3 | effortful_control | 50.179 | 35.329 | 0.749 | 0.527 | 1.420 | 0.158 |
| 4 | extraversion:effortful_control | -10.293 | 7.028 | -1.092 | 0.746 | -1.465 | 0.146 |
| 5 | cond.IIII2 | -17.314 | 12.445 | -0.340 | 0.244 | -1.391 | 0.167 |
| 6 | cond.IIII3 | -36.178 | 13.189 | -0.710 | 0.259 | -2.743 | 0.007 |
| 7 | cond.IIII4 | -25.118 | 13.763 | -0.493 | 0.270 | -1.825 | 0.071 |
| 8 | f_stat | 6.000 | 115.000 | 1.796 | 0.106 | ||
| 9 | R2/Adj.R2 | 0.086 | 0.038 | ||||
| 10 | N | 224.000 | |||||
| Note: | |||||||
| … |
write.csv(red_TTC_II_mod1_t,"red_TTC_II_mod1_t.csv", row.names = T)anova_condII<-aov(average.II~cond.II, dat1)
summary(anova_condII)## Df Sum Sq Mean Sq F value Pr(>F)
## cond.II 3 0.020 0.00668 0.554 0.646
## Residuals 220 2.651 0.01205anova_cond_ii<-aov(average.II~cond.II, dat1)
summary(anova_cond_ii)## Df Sum Sq Mean Sq F value Pr(>F)
## cond.II 3 0.020 0.00668 0.554 0.646
## Residuals 220 2.651 0.01205anova_condII<-aov(first_TTC.II~cond.II, dat1[dat1$first_TTC.II<200,])
summary(anova_condII)## Df Sum Sq Mean Sq F value Pr(>F)
## cond.II 3 21448 7149 2.886 0.0386 *
## Residuals 118 292318 2477
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tuk<-TukeyHSD(anova_condII, conf.level=.95)Tukey plot comparing every 2 II conditions
tuk_plot <- function (x, xlab, ylab, ylabels = NULL, ...) {
for (i in seq_along(x)) {
xi <- x[[i]][, -4L, drop = FALSE]
yvals <- nrow(xi):1L
dev.hold()
on.exit(dev.flush())
plot(c(xi[, "lwr"], xi[, "upr"]), rep.int(yvals, 2L),
type = "n", axes = FALSE, xlab = "", ylab = "", main = NULL,
...)
axis(1, ...)
# change for custom axis labels
if (is.null(ylabels)) ylabels <- dimnames(xi)[[1L]]
axis(2, at = nrow(xi):1, labels = ylabels,
srt = 0, ...)
abline(h = yvals, lty = 1, lwd = 1, col = "lightgray")
abline(v = 0, lty = 2, lwd = 1, ...)
segments(xi[, "lwr"], yvals, xi[, "upr"], yvals, ...)
segments(as.vector(xi), rep.int(yvals - 0.1, 3L), as.vector(xi),
rep.int(yvals + 0.1, 3L), ...)
title(main ="",
# change for custom axis titles
xlab = xlab, ylab = ylab)
box()
dev.flush()
on.exit()
}
}
with(par(mai=c(0.4,0.8,0.1,0.5),cex.lab=1.2, cex.axis=0.8),{tuk_plot(tuk, "", "Condition",c("2-1","3-1","4-1","3-2","4-2","4-3"),las=2)})
# reset NA vals in kable tables to "NA"
options(knitr.kable.NA = NA)
reg_colours <- c("#414535", "#08605F", "#3F826D", "#858585", "#DDA448", "#F0803C", "#98410B" ) # green > orange
# plot lm
# plot standardized betas to compare cooefficients
TTC_red_reg_ii_plot <-
sjPlot::plot_model(model = step_mii1_stan, # lm output
title = "", # to remove title
# type = "std", # type "std" = standardized betas
#sort.est = T, # sort based on betas (highest on top)
# wrap.title = 55, # char length of each title line
#axis.lim = c(-.25,.25)
# grid.breaks = .25, # x-axis breaks
#se = T # use se instead of CI
colors = reg_colours,
show.values = T, # show betas
show.p = T, # show * for sig p
p.threshold = c(.049, .009, .0009), # state up to 3 cut-offs for pvals, the number is inclusive**
value.offset = 0.4, # offset to beta labels
value.size = 3, # size of beta labels
digits = 2, # num decimals in beta labels
dot.size = 2.5, line.size = 1,
vline.color = "#F98B99",
# group.terms = c(1, 2, 3, 3, 3, 3, 3, 4, 5, 6, 7, 7, 7),
axis.labels = c("II Condition 4","II Condition3","II Condition 2", "Extraversion x Effortful Control", "Effortful Control", "Extraversion"),
axis.title = "",
width = 0.2, transform = NULL
) +
# font_size(offset.x = 10) +
# font_size(base.theme = theme_classic()) +
theme(text = element_text(color = "black"),
axis.text = element_text(colour = "black", size = 11),
axis.ticks.y = element_blank(),
)
# show plot
TTC_red_reg_ii_plot
t.test(extraversion~group2, dat1)##
## Welch Two Sample t-test
##
## data: extraversion by group2
## t = -0.2607, df = 214.76, p-value = 0.7946
## alternative hypothesis: true difference in means between group ii_learner and group ii_nlearner is not equal to 0
## 95 percent confidence interval:
## -0.2532671 0.1940975
## sample estimates:
## mean in group ii_learner mean in group ii_nlearner
## 4.983331 5.012915t.test(effortful_control~group2, dat1)##
## Welch Two Sample t-test
##
## data: effortful_control by group2
## t = -0.91897, df = 211.14, p-value = 0.3592
## alternative hypothesis: true difference in means between group ii_learner and group ii_nlearner is not equal to 0
## 95 percent confidence interval:
## -0.3024621 0.1101225
## sample estimates:
## mean in group ii_learner mean in group ii_nlearner
## 3.854616 3.950786t.test(negative_affect~group2, dat1)##
## Welch Two Sample t-test
##
## data: negative_affect by group2
## t = 0.88808, df = 209.65, p-value = 0.3755
## alternative hypothesis: true difference in means between group ii_learner and group ii_nlearner is not equal to 0
## 95 percent confidence interval:
## -0.1067825 0.2818677
## sample estimates:
## mean in group ii_learner mean in group ii_nlearner
## 4.557996 4.470453t.test(Rspan_score~group2, dat1)##
## Welch Two Sample t-test
##
## data: Rspan_score by group2
## t = 1.9274, df = 209.18, p-value = 0.05529
## alternative hypothesis: true difference in means between group ii_learner and group ii_nlearner is not equal to 0
## 95 percent confidence interval:
## -0.04290022 3.80149874
## sample estimates:
## mean in group ii_learner mean in group ii_nlearner
## 48.06557 46.18627t.test(extraversion~group, dat1[dat1$group=="cr_better"|dat1$group=="ii_better",])##
## Welch Two Sample t-test
##
## data: extraversion by group
## t = -0.079294, df = 100.06, p-value = 0.937
## alternative hypothesis: true difference in means between group cr_better and group ii_better is not equal to 0
## 95 percent confidence interval:
## -0.2827597 0.2610258
## sample estimates:
## mean in group cr_better mean in group ii_better
## 4.983477 4.994344t.test(effortful_control~group, dat1[dat1$group=="cr_better"|dat1$group=="ii_better",])##
## Welch Two Sample t-test
##
## data: effortful_control by group
## t = -0.090496, df = 90.281, p-value = 0.9281
## alternative hypothesis: true difference in means between group cr_better and group ii_better is not equal to 0
## 95 percent confidence interval:
## -0.2642527 0.2412264
## sample estimates:
## mean in group cr_better mean in group ii_better
## 3.843750 3.855263t.test(negative_affect~group, dat1[dat1$group=="cr_better"|dat1$group=="ii_better",])##
## Welch Two Sample t-test
##
## data: negative_affect by group
## t = 0.20989, df = 100.2, p-value = 0.8342
## alternative hypothesis: true difference in means between group cr_better and group ii_better is not equal to 0
## 95 percent confidence interval:
## -0.1950678 0.2412260
## sample estimates:
## mean in group cr_better mean in group ii_better
## 4.561597 4.538518t.test(Rspan_score~group, dat1[dat1$group=="cr_better"|dat1$group=="ii_better",])##
## Welch Two Sample t-test
##
## data: Rspan_score by group
## t = 0.22279, df = 102.29, p-value = 0.8241
## alternative hypothesis: true difference in means between group cr_better and group ii_better is not equal to 0
## 95 percent confidence interval:
## -1.970924 2.469722
## sample estimates:
## mean in group cr_better mean in group ii_better
## 47.71094 47.46154