#set working directory to the folder to which files created in R should be saved (e.g., "C:/PCmatch") setwd("C:/PCmatch") #install packages install.packages(c("psych","data.table","car","ggplot2","gridExtra","flextable","officer","tidyverse","dplyr","numform")) #load packages library(psych) library(data.table) library(car) library(ggplot2) library(gridExtra) library(flextable) library(officer) library(tidyverse) library(tibble) library(dplyr) library(numform) #####################################Analyses with full OOS Dataset (not publicly available)##################################### # ##install additional packages #install.packages(c("lavaan","semTools","multilevel")) # ##load additional packages #library(lavaan) #library(semTools) #library(multilevel) # ##read data #dat <- as.data.frame(fread("04_OOS_MASTER DATASET_2016_03_03.csv", header = T, sep = ',')) # ##select variables #dat <- subset(dat, select = c(record_id, cnt_now, age, female, langu, orelig, ose, oext01:oopn10)) # ##recode items to compute oagy_ac and oagy_bf #dat$oneu02 <- car::recode(dat$oneu02r, "1=5; 2=4; 3=3; 4=2; 5=1;") #dat$oneu07 <- car::recode(dat$oneu07r, "1=5; 2=4; 3=3; 4=2; 5=1;") # ##compute personality traits #dat$oagr <- rowMeans(dat[c("oagr01r","oagr02","oagr03r","oagr04","oagr05","oagr06r","oagr07","oagr08r","oagr09")], na.rm=TRUE) #dat$ocns <- rowMeans(dat[c("ocns01","ocns02r","ocns03","ocns04r","ocns05r","ocns06","ocns07","ocns08","ocns09r")], na.rm=TRUE) #dat$oext <- rowMeans(dat[c("oext01","oext02r","oext03","oext04","oext05r","oext06","oext07r","oext08")], na.rm=TRUE) #dat$oopn <- rowMeans(dat[c("oopn01","oopn02","oopn03","oopn04","oopn05","oopn06","oopn07r","oopn08","oopn09r","oopn10")], na.rm=TRUE) #dat$oneu <- rowMeans(dat[c("oneu01","oneu02r","oneu03","oneu04","oneu05r","oneu06","oneu07r","oneu08")], na.rm=TRUE) #dat$oagy_ts <- rowMeans(dat[c("oext01","oext03","oext04","oext05r","oext06","oext07r","oext08","oopn01")], na.rm=TRUE) #dat$ocom_ts <- rowMeans(dat[c("oagr01r","oagr02","oagr03r","oagr04","oagr07","oagr08r","oagr09","ocns07")], na.rm=TRUE) #dat$oagy_ac <- rowMeans(dat[c("oext03","oext06","oext07r","oopn04","oopn07r","oneu07")], na.rm=TRUE) #dat$ocom_ac <- rowMeans(dat[c("oagr01r","oagr03r","oagr05","oagr07","oagr08r","oagr09","ocns01")], na.rm=TRUE) #dat$oagy_bf <- rowMeans(dat[c("ocns06","oext02r","oext05r","oopn01","oopn03","oopn05","oopn10","oneu02")], na.rm=TRUE) #dat$ocom_bf <- rowMeans(dat[c("oagr02","oagr03r","oagr04","oagr06r","oagr07","oagr08r","oagr09","ocns01")], na.rm=TRUE) # ##N = 9,328,610 # ##select data #dat[is.na(dat)] <- NA #dat <- subset(dat, cnt_now != -999 & is.na(orelig)==F & is.na(ose)==F & is.na(oagr)==F & is.na(ocns)==F & is.na(oext)==F & is.na(oopn)==F & is.na(oneu)==F # & is.na(oagy_ts)==F & is.na(ocom_ts)==F & is.na(oagy_ac)==F & is.na(ocom_ac)==F & is.na(oagy_bf)==F & is.na(ocom_bf)==F) #dat <- transform(dat, ocnt_frq = ave(seq(nrow(dat)), cnt_now, FUN=length)) #dat <- subset(dat, ocnt_frq >= 300) # ##N = 850,877 # ##measurement invariance #dat2 <- setNames(data.frame(matrix(ncol = 13, nrow = 6)), c("model","fit_index","oagy_ts","ocom_ts","oagy_ac","ocom_ac","oagy_bf","ocom_bf","oagr","ocns","oext","oopn","oneu")) #dat2$model[1:3] <- "config" #dat2$model[4:6] <- "metric" #dat2$fit_index[c(1,4)] <- "cfi" #dat2$fit_index[c(2,5)] <- "rmsea" #dat2$fit_index[c(3,6)] <- "srmr" # ##BIG TWO #mod_oagy_ts <- 'oagy =~ oext01 + oext03 + oext04 + oext05r + oext06 + oext07r + oext08 + oopn01 # meth =~ 1*oext05r + 1*oext07r # oagy ~~ 0*meth' #config_oagy_ts <- cfa(mod_oagy_ts, data = dat, missing = "ML", group = "cnt_now") #dat2$oagy_ts[1] <- fitMeasures(config_oagy_ts, "cfi") #dat2$oagy_ts[2] <- fitMeasures(config_oagy_ts, "rmsea") #dat2$oagy_ts[3] <- fitMeasures(config_oagy_ts, "srmr") #metric_oagy_ts <- cfa(mod_oagy_ts, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$oagy_ts[4] <- fitMeasures(metric_oagy_ts, "cfi") #dat2$oagy_ts[5] <- fitMeasures(metric_oagy_ts, "rmsea") #dat2$oagy_ts[6] <- fitMeasures(metric_oagy_ts, "srmr") #rm(mod_oagy_ts,config_oagy_ts,metric_oagy_ts) # #mod_ocom_ts <- 'ocom =~ oagr01r + oagr02 + oagr03r + oagr04 + oagr07 + oagr08r + oagr09 + ocns07 # meth =~ 1*oagr01r + 1*oagr03r + 1*oagr08r # ocom ~~ 0*meth' #config_ocom_ts <- cfa(mod_ocom_ts, data = dat, missing = "ML", group = "cnt_now") #dat2$ocom_ts[1] <- fitMeasures(config_ocom_ts, "cfi") #dat2$ocom_ts[2] <- fitMeasures(config_ocom_ts, "rmsea") #dat2$ocom_ts[3] <- fitMeasures(config_ocom_ts, "srmr") #metric_ocom_ts <- cfa(mod_ocom_ts, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$ocom_ts[4] <- fitMeasures(metric_ocom_ts, "cfi") #dat2$ocom_ts[5] <- fitMeasures(metric_ocom_ts, "rmsea") #dat2$ocom_ts[6] <- fitMeasures(metric_ocom_ts, "srmr") #rm(mod_ocom_ts,config_ocom_ts,metric_ocom_ts) # #mod_oagy_ac <- 'oagy =~ oext03 + oext06 + oext07r + oopn04 + oopn07r + oneu07 # meth =~ 1*oext07r + 1*oopn07r # oagy ~~ 0*meth' #config_oagy_ac <- cfa(mod_oagy_ac, data = dat, missing = "ML", group = "cnt_now") #dat2$oagy_ac[1] <- fitMeasures(config_oagy_ac, "cfi") #dat2$oagy_ac[2] <- fitMeasures(config_oagy_ac, "rmsea") #dat2$oagy_ac[3] <- fitMeasures(config_oagy_ac, "srmr") #metric_oagy_ac <- cfa(mod_oagy_ac, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$oagy_ac[4] <- fitMeasures(metric_oagy_ac, "cfi") #dat2$oagy_ac[5] <- fitMeasures(metric_oagy_ac, "rmsea") #dat2$oagy_ac[6] <- fitMeasures(metric_oagy_ac, "srmr") #rm(mod_oagy_ac,config_oagy_ac,metric_oagy_ac) # #mod_ocom_ac <- 'ocom =~ oagr01r + oagr03r + oagr05 + oagr07 + oagr08r + oagr09 + ocns01 # meth =~ 1*oagr01r + 1*oagr03r + 1*oagr08r # ocom ~~ 0*meth' #config_ocom_ac <- cfa(mod_ocom_ac, data = dat, missing = "ML", group = "cnt_now") #dat2$ocom_ac[1] <- fitMeasures(config_ocom_ac, "cfi") #dat2$ocom_ac[2] <- fitMeasures(config_ocom_ac, "rmsea") #dat2$ocom_ac[3] <- fitMeasures(config_ocom_ac, "srmr") #metric_ocom_ac <- cfa(mod_ocom_ac, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$ocom_ac[4] <- fitMeasures(metric_ocom_ac, "cfi") #dat2$ocom_ac[5] <- fitMeasures(metric_ocom_ac, "rmsea") #dat2$ocom_ac[6] <- fitMeasures(metric_ocom_ac, "srmr") #rm(mod_ocom_ac,config_ocom_ac,metric_ocom_ac) # #mod_oagy_bf <- 'oagy =~ ocns06 + oext02r + oext05r + oopn01 + oopn03 + oopn05 + oopn10 + oneu02 # meth =~ 1*oext02r + 1*oext05r # oagy ~~ 0*meth' #config_oagy_bf <- cfa(mod_oagy_bf, data = dat, missing = "ML", group = "cnt_now") #dat2$oagy_bf[1] <- fitMeasures(config_oagy_bf, "cfi") #dat2$oagy_bf[2] <- fitMeasures(config_oagy_bf, "rmsea") #dat2$oagy_bf[3] <- fitMeasures(config_oagy_bf, "srmr") #metric_oagy_bf <- cfa(mod_oagy_bf, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$oagy_bf[4] <- fitMeasures(metric_oagy_bf, "cfi") #dat2$oagy_bf[5] <- fitMeasures(metric_oagy_bf, "rmsea") #dat2$oagy_bf[6] <- fitMeasures(metric_oagy_bf, "srmr") #rm(mod_oagy_bf,config_oagy_bf,metric_oagy_bf) # #mod_ocom_bf <- 'ocom =~ oagr02 + oagr03r + oagr04 + oagr06r + oagr07 + oagr08r + oagr09 + ocns01 # meth =~ 1*oagr03r + 1*oagr06r + 1*oagr08r # ocom ~~ 0*meth' #config_ocom_bf <- cfa(mod_ocom_bf, data = dat, missing = "ML", group = "cnt_now") #dat2$ocom_bf[1] <- fitMeasures(config_ocom_bf, "cfi") #dat2$ocom_bf[2] <- fitMeasures(config_ocom_bf, "rmsea") #dat2$ocom_bf[3] <- fitMeasures(config_ocom_bf, "srmr") #metric_ocom_bf <- cfa(mod_ocom_bf, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$ocom_bf[4] <- fitMeasures(metric_ocom_bf, "cfi") #dat2$ocom_bf[5] <- fitMeasures(metric_ocom_bf, "rmsea") #dat2$ocom_bf[6] <- fitMeasures(metric_ocom_bf, "srmr") #rm(mod_ocom_bf,config_ocom_bf,metric_ocom_bf) # ##BIG FIVE #mod_oagr <- 'oagr.lat =~ oagr01r + oagr02 + oagr03r + oagr04 + oagr05 + oagr06r + oagr07 + oagr08r + oagr09 # meth =~ 1*oagr01r + 1*oagr03r + 1*oagr06r + 1*oagr08r # oagr.lat ~~ 0*meth' #config_oagr <- cfa(mod_oagr, data = dat, missing = "ML", group = "cnt_now") #dat2$oagr[1] <- fitMeasures(config_oagr, "cfi") #dat2$oagr[2] <- fitMeasures(config_oagr, "rmsea") #dat2$oagr[3] <- fitMeasures(config_oagr, "srmr") #metric_oagr <- cfa(mod_oagr, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$oagr[4] <- fitMeasures(metric_oagr, "cfi") #dat2$oagr[5] <- fitMeasures(metric_oagr, "rmsea") #dat2$oagr[6] <- fitMeasures(metric_oagr, "srmr") #rm(mod_oagr,config_oagr,metric_oagr) # #mod_ocns <- 'ocns.lat =~ ocns01 + ocns02r + ocns03 + ocns04r + ocns05r + ocns06 + ocns07 + ocns08 + ocns09r # meth =~ 1*ocns02r + 1*ocns04r + 1*ocns05r + 1*ocns09r # ocns.lat ~~ 0*meth' #config_ocns <- cfa(mod_ocns, data = dat, missing = "ML", group = "cnt_now") #dat2$ocns[1] <- fitMeasures(config_ocns, "cfi") #dat2$ocns[2] <- fitMeasures(config_ocns, "rmsea") #dat2$ocns[3] <- fitMeasures(config_ocns, "srmr") #metric_ocns <- cfa(mod_ocns, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$ocns[4] <- fitMeasures(metric_ocns, "cfi") #dat2$ocns[5] <- fitMeasures(metric_ocns, "rmsea") #dat2$ocns[6] <- fitMeasures(metric_ocns, "srmr") #rm(mod_ocns,config_ocns,metric_ocns) # #mod_oext <- 'oext.lat =~ oext01 + oext02r + oext03 + oext04 + oext05r + oext06 + oext07r + oext08 # meth =~ 1*oext02r + 1*oext05r + 1*oext07r # oext.lat ~~ 0*meth' #config_oext <- cfa(mod_oext, data = dat, missing = "ML", group = "cnt_now") #dat2$oext[1] <- fitMeasures(config_oext, "cfi") #dat2$oext[2] <- fitMeasures(config_oext, "rmsea") #dat2$oext[3] <- fitMeasures(config_oext, "srmr") #metric_oext <- cfa(mod_oext, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$oext[4] <- fitMeasures(metric_oext, "cfi") #dat2$oext[5] <- fitMeasures(metric_oext, "rmsea") #dat2$oext[6] <- fitMeasures(metric_oext, "srmr") #rm(mod_oext,config_oext,metric_oext) # #mod_oopn <- 'oopn.lat =~ oopn01 + oopn02 + oopn03 + oopn04 + oopn05 + oopn06 + oopn07r + oopn08 + oopn09r + oopn10 # meth =~ 1*oopn07r + 1*oopn09r # oopn.lat ~~ 0*meth' #config_oopn <- cfa(mod_oopn, data = dat, missing = "ML", group = "cnt_now") #dat2$oopn[1] <- fitMeasures(config_oopn, "cfi") #dat2$oopn[2] <- fitMeasures(config_oopn, "rmsea") #dat2$oopn[3] <- fitMeasures(config_oopn, "srmr") #metric_oopn <- cfa(mod_oopn, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$oopn[4] <- fitMeasures(metric_oopn, "cfi") #dat2$oopn[5] <- fitMeasures(metric_oopn, "rmsea") #dat2$oopn[6] <- fitMeasures(metric_oopn, "srmr") #rm(mod_oopn,config_oopn,metric_oopn) # #mod_oneu <- 'oneu.lat =~ oneu01 + oneu02r + oneu03 + oneu04 + oneu05r + oneu06 + oneu07r + oneu08 # meth =~ 1*oneu02r + 1*oneu05r + 1*oneu07r # oneu.lat ~~ 0*meth' #config_oneu <- cfa(mod_oneu, data = dat, missing = "ML", group = "cnt_now") #dat2$oneu[1] <- fitMeasures(config_oneu, "cfi") #dat2$oneu[2] <- fitMeasures(config_oneu, "rmsea") #dat2$oneu[3] <- fitMeasures(config_oneu, "srmr") #metric_oneu <- cfa(mod_oneu, data = dat, missing = "ML", group = "cnt_now", group.equal = c("loadings")) #dat2$oneu[4] <- fitMeasures(metric_oneu, "cfi") #dat2$oneu[5] <- fitMeasures(metric_oneu, "rmsea") #dat2$oneu[6] <- fitMeasures(metric_oneu, "srmr") #rm(mod_oneu,config_oneu,metric_oneu) # ##select variables #dat <- subset(dat, select = c(record_id, cnt_now, age, female, langu, ose, orelig, oagy_ts, ocom_ts, oagy_ac, ocom_ac, oagy_bf, ocom_bf, oagr, ocns, oext, oopn, oneu)) # ##save data #fwrite(dat, file = "PCmatch_csv-data_SOM-U1.csv") #fwrite(dat2, file = "PCmatch_mi_SOM-U1.csv") # #rm(dat2) #####################################Analyses with publicly available dataset##################################### #read data dat <- as.data.frame(fread("https://madata.bib.uni-mannheim.de/330/5/PCmatch_csv-data_SOM-U1.csv", header = T, sep = ",")) #descripitives dat$female <- car::recode(dat$female, "-999=NA; 0=0; 1=1") 100*prop.table(table(dat$female)) dat$age <- as.numeric(dat$age) dat$age2 <- ifelse(dat$age > 150, NA, ifelse(dat$age < 0, NA, dat$age)) round(mean(dat$age2, na.rm=T),digits=2) round(sd(dat$age2, na.rm=T),digits=2) 100*prop.table(table(dat$langu)) #Measurement invariance Table dat1 <- fread("https://madata.bib.uni-mannheim.de/330/10/PCmatch_mi_SOM-U1.csv") dat1 <- dat1[, c(2,1,3:10,12,11,13)] dat1 <- add_row(dat1) dat1 <- add_row(dat1) dat1 <- dat1[c(1,4,7,10,2,5,8,11,3,6,9), ] dat1[,3:13] <- round(dat1[,3:13], digits=2) dat1[1,1] <- "CFI" dat1[2,1] <- "" dat1[3,1] <- "" dat1[3,2] <- intToUtf8(916) dat1[3,3:13] <- dat1[1,3:13] - dat1[2,3:13] dat1[5,1] <- "RMSEA" dat1[6,1] <- "" dat1[7,1] <- "" dat1[7,2] <- intToUtf8(916) dat1[7,3:13] <- dat1[5,3:13] - dat1[6,3:13] dat1[9,1] <- "SRMR" dat1[10,1] <- "" dat1[11,1] <- "" dat1[11,2] <- intToUtf8(916) dat1[11,3:13] <- dat1[9,3:13] - dat1[10,3:13] dat1[,3:13] <- mutate_all(data.frame(dat1[,3:13]), list(~f_num(unlist(.), digits=2))) dat1$blank1 <- NA dat1$blank2 <- NA dat1$blank3 <- NA dat1 <- dat1[, c(1:4,14,5:6,15,7:8,16,9:13)] col_keys <- c("fit_index","model","oagy_ts","ocom_ts","blank1","oagy_ac","ocom_ac","blank2","oagy_bf","ocom_bf","blank3","oagr","ocns","oopn","oext","oneu") head1 <- c("","","Big Two","Big Two","Big Two","Big Two","Big Two","Big Two","Big Two","Big Two","","Big Five","Big Five","Big Five","Big Five","Big Five") head2 <- c("","","TS","TS","","AC","AC","","BF","BF","","","","","","") head3 <- c("","","agy","com","","agy","com","","agy","com","","agr","cns","opn","ext","neu") head <- data.frame(col_keys,head1,head2,head3, stringsAsFactors = FALSE) rm(col_keys,head1,head2,head3) tbl <- flextable(dat1) tbl <- set_header_df(tbl, mapping=head, key="col_keys") tbl <- colformat_lgl(tbl, j =~ blank1 + blank2 + blank3, na_str="") tbl <- hline_top(tbl, j=1:16, border=fp_border(width=2), part="header") tbl <- merge_at(tbl, i=1, j=3:10, part="header") tbl <- merge_at(tbl, i=1, j=12:16, part="header") tbl <- merge_at(tbl, i=2, j=3:4, part="header") tbl <- merge_at(tbl, i=2, j=6:7, part="header") tbl <- merge_at(tbl, i=2, j=9:10, part="header") tbl <- hline(tbl, i=1, j=3:10, border=fp_border(width=1), part="header") tbl <- hline(tbl, i=1, j=12:16, border=fp_border(width=1), part="header") tbl <- hline(tbl, i=2, j=3:4, border=fp_border(width=1), part="header") tbl <- hline(tbl, i=2, j=6:7, border=fp_border(width=1), part="header") tbl <- hline(tbl, i=2, j=9:10, border=fp_border(width=1), part="header") tbl <- hline(tbl, i=3, j=3:10, border=fp_border(width=1), part="header") tbl <- hline(tbl, i=3, j=12:16, border=fp_border(width=1), part="header") tbl <- font(tbl, fontname="Times", part="all") tbl <- fontsize(tbl, size=9, part="all") tbl <- bold(tbl, part="header") tbl <- align(tbl, align="center", part="all") tbl <- align(tbl, j = c("fit_index","model"), align="left", part="body") tbl <- width(tbl, j =~ fit_index, width=1) tbl <- width(tbl, j =~ model, width=.7) tbl <- width(tbl, j =~ oagy_ts + ocom_ts + oagy_ac + ocom_ac + oagy_bf + ocom_bf + oagr + ocns + oopn + oext + oneu, width=.5) tbl <- width(tbl, j =~ blank1 + blank2, width=.1) tbl <- width(tbl, j =~ blank3, width=.3) tbl <- height(tbl, i=c(4,8), height=.1, part="body") tbl doc <- read_docx() doc <- body_add_flextable(doc, value = tbl) print(doc, target = "s1_mi_table.docx") rm(dat1,head,tbl,doc) #prepare data for mixed-effects models ##z-standardize self-esteem dat$zose <- scale(dat$ose) ##group-mean center and z-standardize religiosity and personality traits dat$zorelgrp <- scale(dat$orelig - ave(dat$orelig, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zoagytsgrp <- scale(dat$oagy_ts - ave(dat$oagy_ts, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zocomtsgrp <- scale(dat$ocom_ts - ave(dat$ocom_ts, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zoagyacgrp <- scale(dat$oagy_ac - ave(dat$oagy_ac, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zocomacgrp <- scale(dat$ocom_ac - ave(dat$ocom_ac, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zoagybfgrp <- scale(dat$oagy_bf - ave(dat$oagy_bf, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zocombfgrp <- scale(dat$ocom_bf - ave(dat$ocom_bf, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zoagrgrp <- scale(dat$oagr - ave(dat$oagr, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zocnsgrp <- scale(dat$ocns - ave(dat$ocns, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zoextgrp <- scale(dat$oext - ave(dat$oext, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zoopngrp <- scale(dat$oopn - ave(dat$oopn, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) dat$zoneugrp <- scale(dat$oneu - ave(dat$oneu, dat$cnt_now, FUN = function(x) mean(x, na.rm=T))) ##compute country-level religiosiy (z-standardized) dat1 <- aggregate(orelig ~ cnt_now, dat, mean) dat1$zocrel <- scale(dat1$orelig) dat1$orelig <- NULL dat <- merge(dat, dat1, by = "cnt_now", all = TRUE) ###correlation between computed country-level religiosiy and ###country-level religiosity index based on Gallup World Poll data (Joshanloo & Gebauer, 2020) dat2 <- as.data.frame(fread("https://madata.bib.uni-mannheim.de/330/17/crel_gallup.csv", header = T, sep = ",")) dat3 <- merge(dat1, dat2, by = "cnt_now") dat3 <- subset(dat3, select = c(zocrel,crel_gallup)) print(corr.test(dat3), short=F) rm(dat1,dat2,dat3) ##re-center personality traits and country-level religiosity for simple slope analyses dat$zoagytsgrpverylow <- dat$zoagytsgrp + 2 dat$zoagytsgrplow <- dat$zoagytsgrp + 1 dat$zoagytsgrphigh <- dat$zoagytsgrp - 1 dat$zoagytsgrpveryhigh <- dat$zoagytsgrp - 2 dat$zocomtsgrpverylow <- dat$zocomtsgrp + 2 dat$zocomtsgrplow <- dat$zocomtsgrp + 1 dat$zocomtsgrphigh <- dat$zocomtsgrp - 1 dat$zocomtsgrpveryhigh <- dat$zocomtsgrp - 2 dat$zoagyacgrpverylow <- dat$zoagyacgrp + 2 dat$zoagyacgrplow <- dat$zoagyacgrp + 1 dat$zoagyacgrphigh <- dat$zoagyacgrp - 1 dat$zoagyacgrpveryhigh <- dat$zoagyacgrp - 2 dat$zocomacgrpverylow <- dat$zocomacgrp + 2 dat$zocomacgrplow <- dat$zocomacgrp + 1 dat$zocomacgrphigh <- dat$zocomacgrp - 1 dat$zocomacgrpveryhigh <- dat$zocomacgrp - 2 dat$zoagybfgrpverylow <- dat$zoagybfgrp + 2 dat$zoagybfgrplow <- dat$zoagybfgrp + 1 dat$zoagybfgrphigh <- dat$zoagybfgrp - 1 dat$zoagybfgrpveryhigh <- dat$zoagybfgrp - 2 dat$zocombfgrpverylow <- dat$zocombfgrp + 2 dat$zocombfgrplow <- dat$zocombfgrp + 1 dat$zocombfgrphigh <- dat$zocombfgrp - 1 dat$zocombfgrpveryhigh <- dat$zocombfgrp - 2 dat$zoagrgrpverylow <- dat$zoagrgrp + 2 dat$zoagrgrplow <- dat$zoagrgrp + 1 dat$zoagrgrphigh <- dat$zoagrgrp - 1 dat$zoagrgrpveryhigh <- dat$zoagrgrp - 2 dat$zocnsgrpverylow <- dat$zocnsgrp + 2 dat$zocnsgrplow <- dat$zocnsgrp + 1 dat$zocnsgrphigh <- dat$zocnsgrp - 1 dat$zocnsgrpveryhigh <- dat$zocnsgrp - 2 dat$zoextgrpverylow <- dat$zoextgrp + 2 dat$zoextgrplow <- dat$zoextgrp + 1 dat$zoextgrphigh <- dat$zoextgrp - 1 dat$zoextgrpveryhigh <- dat$zoextgrp - 2 dat$zoopngrpverylow <- dat$zoopngrp + 2 dat$zoopngrplow <- dat$zoopngrp + 1 dat$zoopngrphigh <- dat$zoopngrp - 1 dat$zoopngrpveryhigh <- dat$zoopngrp - 2 dat$zoneugrpverylow <- dat$zoneugrp + 2 dat$zoneugrplow <- dat$zoneugrp + 1 dat$zoneugrphigh <- dat$zoneugrp - 1 dat$zoneugrpveryhigh <- dat$zoneugrp - 2 dat$zocrelmax <- dat$zocrel - abs(max(dat$zocrel)) dat$zocrelmin <- dat$zocrel + abs(min(dat$zocrel)) ##select variables dat <- subset(dat, select = c(record_id, age, female, langu, cnt_now, zose:zocrelmin)) #####Before you continue with this R-Syntax, download and install Julia (https://julialang.org/downloads/)##### #install JuliaCall install.packages("JuliaCall") #load JuliaCall library(JuliaCall) #start Julia julia <- julia_setup() #set working directory to the Julia folder indicated in julia_setup (e.g., "C:/AppData/Local/JULIA-~1.1/bin") setwd("C:/AppData/Local/JULIA-~1.1/bin") #install Julia packages julia_install_package("MixedModels") julia_install_package("DataFrames") julia_install_package("StatsBase") #load Julia packages julia_library("MixedModels") julia_library("DataFrames") julia_library("StatsBase") #conduct mixed-effects models in Julia dat$cnt_now <- as.factor(dat$cnt_now) julia_assign("dat", dat) ##BIG TWO TS julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrel * zocomtsgrp + zorelgrp * zocrel * zoagytsgrp + (1 + zorelgrp + zocomtsgrp + zoagytsgrp | cnt_now)), dat)") BigTwoTS <- julia_eval("BigTwoTS = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomtsgrphigh + zorelgrp * zocrelmax * zoagytsgrplow + (1 + zorelgrp + zocomtsgrphigh + zoagytsgrplow | cnt_now)), dat)") BigTwoTS_cmax_asshigh <- julia_eval("BigTwoTS_cmax_asshigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomtsgrpveryhigh + zorelgrp * zocrelmax * zoagytsgrpverylow + (1 + zorelgrp + zocomtsgrpveryhigh + zoagytsgrpverylow | cnt_now)), dat)") BigTwoTS_cmax_assveryhigh <- julia_eval("BigTwoTS_cmax_assveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomtsgrp + zorelgrp * zocrelmax * zoagytsgrp + (1 + zorelgrp + zocomtsgrp + zoagytsgrp | cnt_now)), dat)") BigTwoTS_cmax_average <- julia_eval("BigTwoTS_cmax_average = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomtsgrplow + zorelgrp * zocrelmax * zoagytsgrphigh + (1 + zorelgrp + zocomtsgrplow + zoagytsgrphigh | cnt_now)), dat)") BigTwoTS_cmax_conthigh <- julia_eval("BigTwoTS_cmax_conthigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomtsgrpverylow + zorelgrp * zocrelmax * zoagytsgrpveryhigh + (1 + zorelgrp + zocomtsgrpverylow + zoagytsgrpveryhigh | cnt_now)), dat)") BigTwoTS_cmax_contveryhigh <- julia_eval("BigTwoTS_cmax_contveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomtsgrphigh + zorelgrp * zocrelmin * zoagytsgrplow + (1 + zorelgrp + zocomtsgrphigh + zoagytsgrplow | cnt_now)), dat)") BigTwoTS_cmin_asshigh <- julia_eval("BigTwoTS_cmin_asshigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomtsgrpveryhigh + zorelgrp * zocrelmin * zoagytsgrpverylow + (1 + zorelgrp + zocomtsgrpveryhigh + zoagytsgrpverylow | cnt_now)), dat)") BigTwoTS_cmin_assveryhigh <- julia_eval("BigTwoTS_cmin_assveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomtsgrp + zorelgrp * zocrelmin * zoagytsgrp + (1 + zorelgrp + zocomtsgrp + zoagytsgrp | cnt_now)), dat)") BigTwoTS_cmin_average <- julia_eval("BigTwoTS_cmin_average = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomtsgrplow + zorelgrp * zocrelmin * zoagytsgrphigh + (1 + zorelgrp + zocomtsgrplow + zoagytsgrphigh | cnt_now)), dat)") BigTwoTS_cmin_conthigh <- julia_eval("BigTwoTS_cmin_conthigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomtsgrpverylow + zorelgrp * zocrelmin * zoagytsgrpveryhigh + (1 + zorelgrp + zocomtsgrpverylow + zoagytsgrpveryhigh | cnt_now)), dat)") BigTwoTS_cmin_contveryhigh <- julia_eval("BigTwoTS_cmin_contveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") ##BIG TWO AC julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrel * zocomacgrp + zorelgrp * zocrel * zoagyacgrp + (1 + zorelgrp + zocomacgrp + zoagyacgrp | cnt_now)), dat)") BigTwoAC <- julia_eval("BigTwoAC = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomacgrphigh + zorelgrp * zocrelmax * zoagyacgrplow + (1 + zorelgrp + zocomacgrphigh + zoagyacgrplow | cnt_now)), dat)") BigTwoAC_cmax_asshigh <- julia_eval("BigTwoAC_cmax_asshigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomacgrpveryhigh + zorelgrp * zocrelmax * zoagyacgrpverylow + (1 + zorelgrp + zocomacgrpveryhigh + zoagyacgrpverylow | cnt_now)), dat)") BigTwoAC_cmax_assveryhigh <- julia_eval("BigTwoAC_cmax_assveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomacgrp + zorelgrp * zocrelmax * zoagyacgrp + (1 + zorelgrp + zocomacgrp + zoagyacgrp | cnt_now)), dat)") BigTwoAC_cmax_average <- julia_eval("BigTwoAC_cmax_average = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomacgrplow + zorelgrp * zocrelmax * zoagyacgrphigh + (1 + zorelgrp + zocomacgrplow + zoagyacgrphigh | cnt_now)), dat)") BigTwoAC_cmax_conthigh <- julia_eval("BigTwoAC_cmax_conthigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocomacgrpverylow + zorelgrp * zocrelmax * zoagyacgrpveryhigh + (1 + zorelgrp + zocomacgrpverylow + zoagyacgrpveryhigh | cnt_now)), dat)") BigTwoAC_cmax_contveryhigh <- julia_eval("BigTwoAC_cmax_contveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomacgrphigh + zorelgrp * zocrelmin * zoagyacgrplow + (1 + zorelgrp + zocomacgrphigh + zoagyacgrplow | cnt_now)), dat)") BigTwoAC_cmin_asshigh <- julia_eval("BigTwoAC_cmin_asshigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomacgrpveryhigh + zorelgrp * zocrelmin * zoagyacgrpverylow + (1 + zorelgrp + zocomacgrpveryhigh + zoagyacgrpverylow | cnt_now)), dat)") BigTwoAC_cmin_assveryhigh <- julia_eval("BigTwoAC_cmin_assveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomacgrp + zorelgrp * zocrelmin * zoagyacgrp + (1 + zorelgrp + zocomacgrp + zoagyacgrp | cnt_now)), dat)") BigTwoAC_cmin_average <- julia_eval("BigTwoAC_cmin_average = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomacgrplow + zorelgrp * zocrelmin * zoagyacgrphigh + (1 + zorelgrp + zocomacgrplow + zoagyacgrphigh | cnt_now)), dat)") BigTwoAC_cmin_conthigh <- julia_eval("BigTwoAC_cmin_conthigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocomacgrpverylow + zorelgrp * zocrelmin * zoagyacgrpveryhigh + (1 + zorelgrp + zocomacgrpverylow + zoagyacgrpveryhigh | cnt_now)), dat)") BigTwoAC_cmin_contveryhigh <- julia_eval("BigTwoAC_cmin_contveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") ##BIG TWO BF julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrel * zocombfgrp + zorelgrp * zocrel * zoagybfgrp + (1 + zorelgrp + zocombfgrp + zoagybfgrp | cnt_now)), dat)") BigTwoBF <- julia_eval("BigTwoBF = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocombfgrphigh + zorelgrp * zocrelmax * zoagybfgrplow + (1 + zorelgrp + zocombfgrphigh + zoagybfgrplow | cnt_now)), dat)") BigTwoBF_cmax_asshigh <- julia_eval("BigTwoBF_cmax_asshigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocombfgrpveryhigh + zorelgrp * zocrelmax * zoagybfgrpverylow + (1 + zorelgrp + zocombfgrpveryhigh + zoagybfgrpverylow | cnt_now)), dat)") BigTwoBF_cmax_assveryhigh <- julia_eval("BigTwoBF_cmax_assveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocombfgrp + zorelgrp * zocrelmax * zoagybfgrp + (1 + zorelgrp + zocombfgrp + zoagybfgrp | cnt_now)), dat)") BigTwoBF_cmax_average <- julia_eval("BigTwoBF_cmax_average = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocombfgrplow + zorelgrp * zocrelmax * zoagybfgrphigh + (1 + zorelgrp + zocombfgrplow + zoagybfgrphigh | cnt_now)), dat)") BigTwoBF_cmax_conthigh <- julia_eval("BigTwoBF_cmax_conthigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zocombfgrpverylow + zorelgrp * zocrelmax * zoagybfgrpveryhigh + (1 + zorelgrp + zocombfgrpverylow + zoagybfgrpveryhigh | cnt_now)), dat)") BigTwoBF_cmax_contveryhigh <- julia_eval("BigTwoBF_cmax_contveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocombfgrphigh + zorelgrp * zocrelmin * zoagybfgrplow + (1 + zorelgrp + zocombfgrphigh + zoagybfgrplow | cnt_now)), dat)") BigTwoBF_cmin_asshigh <- julia_eval("BigTwoBF_cmin_asshigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocombfgrpveryhigh + zorelgrp * zocrelmin * zoagybfgrpverylow + (1 + zorelgrp + zocombfgrpveryhigh + zoagybfgrpverylow | cnt_now)), dat)") BigTwoBF_cmin_assveryhigh <- julia_eval("BigTwoBF_cmin_assveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocombfgrp + zorelgrp * zocrelmin * zoagybfgrp + (1 + zorelgrp + zocombfgrp + zoagybfgrp | cnt_now)), dat)") BigTwoBF_cmin_average <- julia_eval("BigTwoBF_cmin_average = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocombfgrplow + zorelgrp * zocrelmin * zoagybfgrphigh + (1 + zorelgrp + zocombfgrplow + zoagybfgrphigh | cnt_now)), dat)") BigTwoBF_cmin_conthigh <- julia_eval("BigTwoBF_cmin_conthigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zocombfgrpverylow + zorelgrp * zocrelmin * zoagybfgrpveryhigh + (1 + zorelgrp + zocombfgrpverylow + zoagybfgrpveryhigh | cnt_now)), dat)") BigTwoBF_cmin_contveryhigh <- julia_eval("BigTwoBF_cmin_contveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") ##BIG FIVE julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrel * zoagrgrp + zorelgrp * zocrel * zocnsgrp + zorelgrp * zocrel * zoopngrp + zorelgrp * zocrel * zoextgrp + zorelgrp * zocrel * zoneugrp + (1 + zorelgrp + zoagrgrp + zocnsgrp + zoopngrp + zoextgrp + zoneugrp | cnt_now)), dat)") BigFive <- julia_eval("BigFive = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zoagrgrplow + zorelgrp * zocrelmax * zocnsgrphigh + zorelgrp * zocrelmax * zoopngrphigh + zorelgrp * zocrelmax * zoextgrphigh + zorelgrp * zocrelmax * zoneugrplow + (1 + zorelgrp + zoagrgrplow + zocnsgrphigh + zoopngrphigh + zoextgrphigh + zoneugrplow | cnt_now)), dat)") BigFive_cmax_conthigh <- julia_eval("BigFive_cmax_conthigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zoagrgrpverylow + zorelgrp * zocrelmax * zocnsgrpveryhigh + zorelgrp * zocrelmax * zoopngrpveryhigh + zorelgrp * zocrelmax * zoextgrpveryhigh + zorelgrp * zocrelmax * zoneugrpverylow + (1 + zorelgrp + zoagrgrpverylow + zocnsgrpveryhigh + zoopngrpveryhigh + zoextgrpveryhigh + zoneugrpverylow | cnt_now)), dat)") BigFive_cmax_contveryhigh <- julia_eval("BigFive_cmax_contveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zoagrgrp + zorelgrp * zocrelmax * zocnsgrp + zorelgrp * zocrelmax * zoopngrp + zorelgrp * zocrelmax * zoextgrp + zorelgrp * zocrelmax * zoneugrp + (1 + zorelgrp + zoagrgrp + zocnsgrp + zoopngrp + zoextgrp + zoneugrp | cnt_now)), dat)") BigFive_cmax_average <- julia_eval("BigFive_cmax_average = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zoagrgrphigh + zorelgrp * zocrelmax * zocnsgrplow + zorelgrp * zocrelmax * zoopngrplow + zorelgrp * zocrelmax * zoextgrplow + zorelgrp * zocrelmax * zoneugrphigh + (1 + zorelgrp + zoagrgrphigh + zocnsgrplow + zoopngrplow + zoextgrplow + zoneugrphigh | cnt_now)), dat)") BigFive_cmax_asshigh <- julia_eval("BigFive_cmax_asshigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmax * zoagrgrpveryhigh + zorelgrp * zocrelmax * zocnsgrpverylow + zorelgrp * zocrelmax * zoopngrpverylow + zorelgrp * zocrelmax * zoextgrpverylow + zorelgrp * zocrelmax * zoneugrpveryhigh + (1 + zorelgrp + zoagrgrpveryhigh + zocnsgrpverylow + zoopngrpverylow + zoextgrpverylow + zoneugrpveryhigh | cnt_now)), dat)") BigFive_cmax_assveryhigh <- julia_eval("BigFive_cmax_assveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zoagrgrplow + zorelgrp * zocrelmin * zocnsgrphigh + zorelgrp * zocrelmin * zoopngrphigh + zorelgrp * zocrelmin * zoextgrphigh + zorelgrp * zocrelmin * zoneugrplow + (1 + zorelgrp + zoagrgrplow + zocnsgrphigh + zoopngrphigh + zoextgrphigh + zoneugrplow | cnt_now)), dat)") BigFive_cmin_conthigh <- julia_eval("BigFive_cmin_conthigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zoagrgrpverylow + zorelgrp * zocrelmin * zocnsgrpveryhigh + zorelgrp * zocrelmin * zoopngrpveryhigh + zorelgrp * zocrelmin * zoextgrpveryhigh + zorelgrp * zocrelmin * zoneugrpverylow + (1 + zorelgrp + zoagrgrpverylow + zocnsgrpveryhigh + zoopngrpveryhigh + zoextgrpveryhigh + zoneugrpverylow | cnt_now)), dat)") BigFive_cmin_contveryhigh <- julia_eval("BigFive_cmin_contveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zoagrgrp + zorelgrp * zocrelmin * zocnsgrp + zorelgrp * zocrelmin * zoopngrp + zorelgrp * zocrelmin * zoextgrp + zorelgrp * zocrelmin * zoneugrp + (1 + zorelgrp + zoagrgrp + zocnsgrp + zoopngrp + zoextgrp + zoneugrp | cnt_now)), dat)") BigFive_cmin_average <- julia_eval("BigFive_cmin_average = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zoagrgrphigh + zorelgrp * zocrelmin * zocnsgrplow + zorelgrp * zocrelmin * zoopngrplow + zorelgrp * zocrelmin * zoextgrplow + zorelgrp * zocrelmin * zoneugrphigh + (1 + zorelgrp + zoagrgrphigh + zocnsgrplow + zoopngrplow + zoextgrplow + zoneugrphigh | cnt_now)), dat)") BigFive_cmin_asshigh <- julia_eval("BigFive_cmin_asshigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") julia_eval("fm = fit(LinearMixedModel, @formula(zose ~ 1 + zorelgrp * zocrelmin * zoagrgrpveryhigh + zorelgrp * zocrelmin * zocnsgrpverylow + zorelgrp * zocrelmin * zoopngrpverylow + zorelgrp * zocrelmin * zoextgrpverylow + zorelgrp * zocrelmin * zoneugrpveryhigh + (1 + zorelgrp + zoagrgrpveryhigh + zocnsgrpverylow + zoopngrpverylow + zoextgrpverylow + zoneugrpveryhigh | cnt_now)), dat)") BigFive_cmin_assveryhigh <- julia_eval("BigFive_cmin_assveryhigh = rename!(DataFrame([(coeftable(fm).rownms) (round.(fixef(fm, false),digits=5)) (round.(((fixef(fm, false))-1.96*(stderror(fm))),digits=5)) (round.(((fixef(fm, false))+1.96*(stderror(fm))),digits=5))]), :x1 => :predictor, :x2 => :zPe, :x3 => :lowerCI, :x4 => :upperCI)") #set working directory to the folder to which files created in R should be saved (e.g., "C:/PCmatch") setwd("C:/PCmatch") #Table SOM-U1 BigTwoTS[, c(2:4)] <- sapply(BigTwoTS[, c(2:4)], as.numeric) BigTwoAC[, c(2:4)] <- sapply(BigTwoAC[, c(2:4)], as.numeric) BigTwoBF[, c(2:4)] <- sapply(BigTwoBF[, c(2:4)], as.numeric) BigFive[, c(2:4)] <- sapply(BigFive[, c(2:4)], as.numeric) dat1 <- data.frame(BigTwoTS,BigTwoAC,BigTwoBF,BigFive) dat1[,5] <- NULL dat1[,8] <- NULL dat1[,1] <- c("(intercept)","rel", "cnt rel", "com", "agy", "rel x cnt rel", "rel x com", "cnt rel x com", "rel x agy", "cnt rel x agy", "rel x cnt rel x com", "rel x cnt rel x agy", "", "", "", "", "", "", "", "", "", "", "", "") dat1[,11] <- c("(intercept)","rel", "cnt rel", "agr", "cns", "opn", "ext", "neu", "rel x cnt rel", "rel x agr", "cnt rel x agr", "rel x cns", "cnt rel x cns", "rel x opn", "cnt rel x opn", "rel x ext", "cnt rel x ext", "rel x neu", "cnt rel x neu", "rel x cnt rel x agr", "rel x cnt rel x cns", "rel x cnt rel x opn", "rel x cnt rel x ext", "rel x cnt rel x neu") trunc0e <- function(x){ ifelse(abs(x)< 0.0005, formatC(x, format = "e", digits = 0), f_num(unlist(x), digits=3)) } dat1 <- dat1 %>% mutate_if(is.numeric, trunc0e) dat1$ci <- paste0("[", dat1$lowerCI, ", ", dat1$upperCI, "]") dat1$ci.1 <- paste0("[", dat1$lowerCI.1, ", ", dat1$upperCI.1, "]") dat1$ci.2 <- paste0("[", dat1$lowerCI.2, ", ", dat1$upperCI.2, "]") dat1$ci.3 <- paste0("[", dat1$lowerCI.3, ", ", dat1$upperCI.3, "]") dat1$blank <- NA dat1$blank.1 <- NA dat1$blank.2 <- NA dat1 <- dat1[c("predictor","zPe","ci","blank","zPe.1","ci.1","blank.1","zPe.2","ci.2","blank.2","predictor.3","zPe.3","ci.3")] dat1[13:24,1:10] <- NA col_keys <- c("predictor","zPe","ci","blank","zPe.1","ci.1","blank.1","zPe.2","ci.2","blank.2","predictor.3","zPe.3","ci.3") head1 <- c("Big Two","Big Two","Big Two","","Big Two","Big Two","","Big Two","Big Two","","Big Five","Big Five","Big Five") head2 <- c("","TS","TS","","AC","AC","","BF","BF","","","","") head3 <- c(""," zPE"," 95% CI",""," zPE"," 95% CI",""," zPE"," 95% CI","",""," zPE"," 95% CI") head <- data.frame(col_keys,head1,head2,head3, stringsAsFactors = FALSE) rm(col_keys,head1,head2,head3) tbl <- flextable(dat1) tbl <- set_header_df(tbl, mapping=head, key="col_keys") tbl <- hline_top(tbl, j=1:13, border=fp_border(width=2), part="header") tbl <- merge_at(tbl, i=1, j=1:9, part="header") tbl <- merge_at(tbl, i=1, j=11:13, part="header") tbl <- merge_at(tbl, i=2, j=2:3, part="header") tbl <- merge_at(tbl, i=2, j=5:6, part="header") tbl <- merge_at(tbl, i=2, j=8:9, part="header") tbl <- hline(tbl, i=1, j=1:3, border=fp_border(width=1.2), part="header") tbl <- hline(tbl, i=1, j=5:6, border=fp_border(width=1.2), part="header") tbl <- hline(tbl, i=1, j=8:9, border=fp_border(width=1.2), part="header") tbl <- hline(tbl, i=1, j=11:13, border=fp_border(width=1.2), part="header") tbl <- hline(tbl, i=2, j=2:3, border=fp_border(width=1.2), part="header") tbl <- hline(tbl, i=2, j=5:6, border=fp_border(width=1.2), part="header") tbl <- hline(tbl, i=2, j=8:9, border=fp_border(width=1.2), part="header") tbl <- hline(tbl, i=3, j=1:9, border=fp_border(width=1.2), part="header") tbl <- hline(tbl, i=3, j=11:13, border=fp_border(width=1.2), part="header") tbl <- font(tbl, fontname="Times", part="all") tbl <- fontsize(tbl, size=11, part="all") tbl <- bold(tbl, part="header") tbl <- italic(tbl, i=3, j=c("zPe","zPe.1","zPe.2","zPe.3"), part="header") tbl <- align(tbl, align="center", part="all") tbl <- align(tbl, j = c("ci","ci.1","ci.2","ci.3"), align="left", part="header") tbl <- align(tbl, j = c("predictor","predictor.3","ci","ci.1","ci.2","ci.3"), align="left", part="body") tbl <- align(tbl, j = c("zPe","zPe.1","zPe.2","zPe.3"), align="right", part="body") tbl <- width(tbl, j =~ predictor + predictor.3, width=1.3) tbl <- width(tbl, j =~ zPe + zPe.1 + zPe.2 + zPe.3, width=.67) tbl <- width(tbl, j =~ ci + ci.1 + ci.2 + ci.3, width=1.1) tbl <- width(tbl, j =~ blank + blank.1 + blank.2, width=.1) tbl <- colformat_lgl(tbl, j =~ blank + blank.1 + blank.2, na_str="") tbl <- height_all(tbl, height=.1, part="all") tbl doc <- read_docx() doc <- body_add_flextable(doc, value = tbl) print(doc, target = "Table_SOM-U1.docx") rm(dat1,head,tbl,doc) #simple slopes for the religiosity x country-level religiosity interactions ##mean association between religiosity and self-esteem in least religious countries BigTwoTS_cmin_average[, c(2:4)] <- sapply(BigTwoTS_cmin_average[, c(2:4)], as.numeric) BigTwoAC_cmin_average[, c(2:4)] <- sapply(BigTwoAC_cmin_average[, c(2:4)], as.numeric) BigTwoBF_cmin_average[, c(2:4)] <- sapply(BigTwoBF_cmin_average[, c(2:4)], as.numeric) BigFive_cmin_average[, c(2:4)] <- sapply(BigFive_cmin_average[, c(2:4)], as.numeric) cmin_zpe <- (BigTwoTS_cmin_average[2,2] + BigTwoAC_cmin_average[2,2] + BigTwoBF_cmin_average[2,2] + BigFive_cmin_average[2,2]) / 4 cmin_lowerCI <- (BigTwoTS_cmin_average[2,3] + BigTwoAC_cmin_average[2,3] + BigTwoBF_cmin_average[2,3] + BigFive_cmin_average[2,3]) / 4 cmin_upperCI <- (BigTwoTS_cmin_average[2,4] + BigTwoAC_cmin_average[2,4] + BigTwoBF_cmin_average[2,4] + BigFive_cmin_average[2,4]) / 4 paste0("least religious countries: zPE = ", round(cmin_zpe, digits=2), " 95% CI [", round(cmin_lowerCI, digits=4), ", ", round(cmin_upperCI, digits=2), "]") ##mean association between religiosity and self-esteem in most religious countries BigTwoTS_cmax_average[, c(2:4)] <- sapply(BigTwoTS_cmax_average[, c(2:4)], as.numeric) BigTwoAC_cmax_average[, c(2:4)] <- sapply(BigTwoAC_cmax_average[, c(2:4)], as.numeric) BigTwoBF_cmax_average[, c(2:4)] <- sapply(BigTwoBF_cmax_average[, c(2:4)], as.numeric) BigFive_cmax_average[, c(2:4)] <- sapply(BigFive_cmax_average[, c(2:4)], as.numeric) cmax_zpe <- (BigTwoTS_cmax_average[2,2] + BigTwoAC_cmax_average[2,2] + BigTwoBF_cmax_average[2,2] + BigFive_cmax_average[2,2]) / 4 cmax_lowerCI <- (BigTwoTS_cmax_average[2,3] + BigTwoAC_cmax_average[2,3] + BigTwoBF_cmax_average[2,3] + BigFive_cmax_average[2,3]) / 4 cmax_upperCI <- (BigTwoTS_cmax_average[2,4] + BigTwoAC_cmax_average[2,4] + BigTwoBF_cmax_average[2,4] + BigFive_cmax_average[2,4]) / 4 paste0("most religious countries: zPE = ", round(cmax_zpe, digits=2), " 95% CI [", round(cmax_lowerCI, digits=2), ", ", round(cmax_upperCI, digits=2), "]") #mean size of the person-culture match effect power_culture <- round(cmax_zpe, digits=2) - round(cmin_zpe, digits=2) power_lowerCI <- round(cmax_lowerCI, digits=2) - round(cmin_upperCI, digits=2) power_upperCI <- round(cmax_upperCI, digits=2) - round(cmin_lowerCI, digits=3) paste0("power of culture: delta_zPE = ", round(power_culture, digits=2), " 95% CI [", round(power_lowerCI, digits=2), ", ", round(power_upperCI, digits=2), "]") rm(cmin_zpe,cmin_lowerCI,cmin_upperCI,cmax_zpe,cmax_lowerCI,cmax_upperCI,power_culture,power_lowerCI,power_upperCI) #Figure SOM-U1 BigTwoTS_cmax_assveryhigh[, c(2:4)] <- sapply(BigTwoTS_cmax_assveryhigh[, c(2:4)], as.numeric) BigTwoTS_cmax_asshigh[, c(2:4)] <- sapply(BigTwoTS_cmax_asshigh[, c(2:4)], as.numeric) BigTwoTS_cmax_average[, c(2:4)] <- sapply(BigTwoTS_cmax_average[, c(2:4)], as.numeric) BigTwoTS_cmax_conthigh[, c(2:4)] <- sapply(BigTwoTS_cmax_conthigh[, c(2:4)], as.numeric) BigTwoTS_cmax_contveryhigh[, c(2:4)] <- sapply(BigTwoTS_cmax_contveryhigh[, c(2:4)], as.numeric) BigTwoTS_cmin_assveryhigh[, c(2:4)] <- sapply(BigTwoTS_cmin_assveryhigh[, c(2:4)], as.numeric) BigTwoTS_cmin_asshigh[, c(2:4)] <- sapply(BigTwoTS_cmin_asshigh[, c(2:4)], as.numeric) BigTwoTS_cmin_average[, c(2:4)] <- sapply(BigTwoTS_cmin_average[, c(2:4)], as.numeric) BigTwoTS_cmin_conthigh[, c(2:4)] <- sapply(BigTwoTS_cmin_conthigh[, c(2:4)], as.numeric) BigTwoTS_cmin_contveryhigh[, c(2:4)] <- sapply(BigTwoTS_cmin_contveryhigh[, c(2:4)], as.numeric) BigTwoAC_cmax_assveryhigh[, c(2:4)] <- sapply(BigTwoAC_cmax_assveryhigh[, c(2:4)], as.numeric) BigTwoAC_cmax_asshigh[, c(2:4)] <- sapply(BigTwoAC_cmax_asshigh[, c(2:4)], as.numeric) BigTwoAC_cmax_average[, c(2:4)] <- sapply(BigTwoAC_cmax_average[, c(2:4)], as.numeric) BigTwoAC_cmax_conthigh[, c(2:4)] <- sapply(BigTwoAC_cmax_conthigh[, c(2:4)], as.numeric) BigTwoAC_cmax_contveryhigh[, c(2:4)] <- sapply(BigTwoAC_cmax_contveryhigh[, c(2:4)], as.numeric) BigTwoAC_cmin_assveryhigh[, c(2:4)] <- sapply(BigTwoAC_cmin_assveryhigh[, c(2:4)], as.numeric) BigTwoAC_cmin_asshigh[, c(2:4)] <- sapply(BigTwoAC_cmin_asshigh[, c(2:4)], as.numeric) BigTwoAC_cmin_average[, c(2:4)] <- sapply(BigTwoAC_cmin_average[, c(2:4)], as.numeric) BigTwoAC_cmin_conthigh[, c(2:4)] <- sapply(BigTwoAC_cmin_conthigh[, c(2:4)], as.numeric) BigTwoAC_cmin_contveryhigh[, c(2:4)] <- sapply(BigTwoAC_cmin_contveryhigh[, c(2:4)], as.numeric) BigTwoBF_cmax_assveryhigh[, c(2:4)] <- sapply(BigTwoBF_cmax_assveryhigh[, c(2:4)], as.numeric) BigTwoBF_cmax_asshigh[, c(2:4)] <- sapply(BigTwoBF_cmax_asshigh[, c(2:4)], as.numeric) BigTwoBF_cmax_average[, c(2:4)] <- sapply(BigTwoBF_cmax_average[, c(2:4)], as.numeric) BigTwoBF_cmax_conthigh[, c(2:4)] <- sapply(BigTwoBF_cmax_conthigh[, c(2:4)], as.numeric) BigTwoBF_cmax_contveryhigh[, c(2:4)] <- sapply(BigTwoBF_cmax_contveryhigh[, c(2:4)], as.numeric) BigTwoBF_cmin_assveryhigh[, c(2:4)] <- sapply(BigTwoBF_cmin_assveryhigh[, c(2:4)], as.numeric) BigTwoBF_cmin_asshigh[, c(2:4)] <- sapply(BigTwoBF_cmin_asshigh[, c(2:4)], as.numeric) BigTwoBF_cmin_average[, c(2:4)] <- sapply(BigTwoBF_cmin_average[, c(2:4)], as.numeric) BigTwoBF_cmin_conthigh[, c(2:4)] <- sapply(BigTwoBF_cmin_conthigh[, c(2:4)], as.numeric) BigTwoBF_cmin_contveryhigh[, c(2:4)] <- sapply(BigTwoBF_cmin_contveryhigh[, c(2:4)], as.numeric) BFIBigTwo <- c("TS/least religious countries","AC/least religious countries","BF/least religious countries","TS/most religious countries","AC/most religious countries","BF/most religious countries", "TS/least religious countries","AC/least religious countries","BF/least religious countries","TS/most religious countries","AC/most religious countries","BF/most religious countries", "TS/least religious countries","AC/least religious countries","BF/least religious countries","TS/most religious countries","AC/most religious countries","BF/most religious countries", "TS/least religious countries","AC/least religious countries","BF/least religious countries","TS/most religious countries","AC/most religious countries","BF/most religious countries", "TS/least religious countries","AC/least religious countries","BF/least religious countries","TS/most religious countries","AC/most religious countries","BF/most religious countries") ass_cont <- c("strong\ncontrasters","strong\ncontrasters","strong\ncontrasters","strong\ncontrasters","strong\ncontrasters","strong\ncontrasters", "contrasters","contrasters","contrasters","contrasters","contrasters","contrasters", "average","average","average","average","average","average", "assimilators","assimilators","assimilators","assimilators","assimilators","assimilators", "strong\nassimilators","strong\nassimilators","strong\nassimilators","strong\nassimilators","strong\nassimilators","strong\nassimilators") zPE <- c(BigTwoTS_cmin_contveryhigh$zPe[2],BigTwoAC_cmin_contveryhigh$zPe[2],BigTwoBF_cmin_contveryhigh$zPe[2],BigTwoTS_cmax_contveryhigh$zPe[2],BigTwoAC_cmax_contveryhigh$zPe[2],BigTwoBF_cmax_contveryhigh$zPe[2], BigTwoTS_cmin_conthigh$zPe[2],BigTwoAC_cmin_conthigh$zPe[2],BigTwoBF_cmin_conthigh$zPe[2],BigTwoTS_cmax_conthigh$zPe[2],BigTwoAC_cmax_conthigh$zPe[2],BigTwoBF_cmax_conthigh$zPe[2], BigTwoTS_cmin_average$zPe[2],BigTwoAC_cmin_average$zPe[2],BigTwoBF_cmin_average$zPe[2],BigTwoTS_cmax_average$zPe[2],BigTwoAC_cmax_average$zPe[2],BigTwoBF_cmax_average$zPe[2], BigTwoTS_cmin_asshigh$zPe[2],BigTwoAC_cmin_asshigh$zPe[2],BigTwoBF_cmin_asshigh$zPe[2],BigTwoTS_cmax_asshigh$zPe[2],BigTwoAC_cmax_asshigh$zPe[2],BigTwoBF_cmax_asshigh$zPe[2], BigTwoTS_cmin_assveryhigh$zPe[2],BigTwoAC_cmin_assveryhigh$zPe[2],BigTwoBF_cmin_assveryhigh$zPe[2],BigTwoTS_cmax_assveryhigh$zPe[2],BigTwoAC_cmax_assveryhigh$zPe[2],BigTwoBF_cmax_assveryhigh$zPe[2]) lowerCI <- c(BigTwoTS_cmin_contveryhigh$lowerCI[2],BigTwoAC_cmin_contveryhigh$lowerCI[2],BigTwoBF_cmin_contveryhigh$lowerCI[2],BigTwoTS_cmax_contveryhigh$lowerCI[2],BigTwoAC_cmax_contveryhigh$lowerCI[2],BigTwoBF_cmax_contveryhigh$lowerCI[2], BigTwoTS_cmin_conthigh$lowerCI[2],BigTwoAC_cmin_conthigh$lowerCI[2],BigTwoBF_cmin_conthigh$lowerCI[2],BigTwoTS_cmax_conthigh$lowerCI[2],BigTwoAC_cmax_conthigh$lowerCI[2],BigTwoBF_cmax_conthigh$lowerCI[2], BigTwoTS_cmin_average$lowerCI[2],BigTwoAC_cmin_average$lowerCI[2],BigTwoBF_cmin_average$lowerCI[2],BigTwoTS_cmax_average$lowerCI[2],BigTwoAC_cmax_average$lowerCI[2],BigTwoBF_cmax_average$lowerCI[2], BigTwoTS_cmin_asshigh$lowerCI[2],BigTwoAC_cmin_asshigh$lowerCI[2],BigTwoBF_cmin_asshigh$lowerCI[2],BigTwoTS_cmax_asshigh$lowerCI[2],BigTwoAC_cmax_asshigh$lowerCI[2],BigTwoBF_cmax_asshigh$lowerCI[2], BigTwoTS_cmin_assveryhigh$lowerCI[2],BigTwoAC_cmin_assveryhigh$lowerCI[2],BigTwoBF_cmin_assveryhigh$lowerCI[2],BigTwoTS_cmax_assveryhigh$lowerCI[2],BigTwoAC_cmax_assveryhigh$lowerCI[2],BigTwoBF_cmax_assveryhigh$lowerCI[2]) upperCI <- c(BigTwoTS_cmin_contveryhigh$upperCI[2],BigTwoAC_cmin_contveryhigh$upperCI[2],BigTwoBF_cmin_contveryhigh$upperCI[2],BigTwoTS_cmax_contveryhigh$upperCI[2],BigTwoAC_cmax_contveryhigh$upperCI[2],BigTwoBF_cmax_contveryhigh$upperCI[2], BigTwoTS_cmin_conthigh$upperCI[2],BigTwoAC_cmin_conthigh$upperCI[2],BigTwoBF_cmin_conthigh$upperCI[2],BigTwoTS_cmax_conthigh$upperCI[2],BigTwoAC_cmax_conthigh$upperCI[2],BigTwoBF_cmax_conthigh$upperCI[2], BigTwoTS_cmin_average$upperCI[2],BigTwoAC_cmin_average$upperCI[2],BigTwoBF_cmin_average$upperCI[2],BigTwoTS_cmax_average$upperCI[2],BigTwoAC_cmax_average$upperCI[2],BigTwoBF_cmax_average$upperCI[2], BigTwoTS_cmin_asshigh$upperCI[2],BigTwoAC_cmin_asshigh$upperCI[2],BigTwoBF_cmin_asshigh$upperCI[2],BigTwoTS_cmax_asshigh$upperCI[2],BigTwoAC_cmax_asshigh$upperCI[2],BigTwoBF_cmax_asshigh$upperCI[2], BigTwoTS_cmin_assveryhigh$upperCI[2],BigTwoAC_cmin_assveryhigh$upperCI[2],BigTwoBF_cmin_assveryhigh$upperCI[2],BigTwoTS_cmax_assveryhigh$upperCI[2],BigTwoAC_cmax_assveryhigh$upperCI[2],BigTwoBF_cmax_assveryhigh$upperCI[2]) dat1 <- data.frame(BFIBigTwo,ass_cont,zPE,lowerCI,upperCI) rm(BFIBigTwo,ass_cont,zPE,lowerCI,upperCI) dat1$BFIBigTwo <- factor(dat1$BFIBigTwo, levels=unique(as.character(dat1$BFIBigTwo))) dat1$ass_cont <- factor(dat1$ass_cont, levels=unique(as.character(dat1$ass_cont))) BigTwoPlot <- ggplot(data=dat1, aes(fill=BFIBigTwo, x=ass_cont, y=zPE, group=BFIBigTwo, color=BFIBigTwo)) + geom_bar(stat="identity", position=position_dodge(), width=0.9, color="black") + geom_errorbar(aes(ymin=lowerCI, ymax=upperCI), position=position_dodge(width=0.9), width=0.3, size=0.5, color="black") + ylim(-0.075, 0.33) + scale_fill_manual(values=c("#f4a460", "#cd853f", "#8b4513", "#AAEEFF", "#2E9AFE", "#2E64FE")) + theme(axis.title.x = element_blank(), legend.position=c(0.17, 0.8), legend.title=element_blank(), axis.text = element_text(color="black", size=12), legend.text = element_text(size=12), axis.title.y = element_text(size=14), panel.background = element_blank(), panel.grid = element_blank(), axis.line = element_line(colour = "black")) + geom_hline(yintercept=0) + ggtitle("Figure SOM-U1a. Big Two") + ylab("religiosity-esteem association (zPE)") BigFive_cmax_assveryhigh[, c(2:4)] <- sapply(BigFive_cmax_assveryhigh[, c(2:4)], as.numeric) BigFive_cmax_asshigh[, c(2:4)] <- sapply(BigFive_cmax_asshigh[, c(2:4)], as.numeric) BigFive_cmax_average[, c(2:4)] <- sapply(BigFive_cmax_average[, c(2:4)], as.numeric) BigFive_cmax_conthigh[, c(2:4)] <- sapply(BigFive_cmax_conthigh[, c(2:4)], as.numeric) BigFive_cmax_contveryhigh[, c(2:4)] <- sapply(BigFive_cmax_contveryhigh[, c(2:4)], as.numeric) BigFive_cmin_assveryhigh[, c(2:4)] <- sapply(BigFive_cmin_assveryhigh[, c(2:4)], as.numeric) BigFive_cmin_asshigh[, c(2:4)] <- sapply(BigFive_cmin_asshigh[, c(2:4)], as.numeric) BigFive_cmin_average[, c(2:4)] <- sapply(BigFive_cmin_average[, c(2:4)], as.numeric) BigFive_cmin_conthigh[, c(2:4)] <- sapply(BigFive_cmin_conthigh[, c(2:4)], as.numeric) BigFive_cmin_contveryhigh[, c(2:4)] <- sapply(BigFive_cmin_contveryhigh[, c(2:4)], as.numeric) cnt <- c("least religious countries","most religious countries","least religious countries","most religious countries","least religious countries","most religious countries", "least religious countries","most religious countries","least religious countries","most religious countries") ass_cont <- c("strong\ncontrasters","strong\ncontrasters","contrasters","contrasters","average","average","assimilators","assimilators","strong\nassimilators","strong\nassimilators") zPE <- c(BigFive_cmin_contveryhigh$zPe[2],BigFive_cmax_contveryhigh$zPe[2],BigFive_cmin_conthigh$zPe[2],BigFive_cmax_conthigh$zPe[2],BigFive_cmin_average$zPe[2],BigFive_cmax_average$zPe[2], BigFive_cmin_asshigh$zPe[2],BigFive_cmax_asshigh$zPe[2],BigFive_cmin_assveryhigh$zPe[2],BigFive_cmax_assveryhigh$zPe[2]) lowerCI <- c(BigFive_cmin_contveryhigh$lowerCI[2],BigFive_cmax_contveryhigh$lowerCI[2],BigFive_cmin_conthigh$lowerCI[2],BigFive_cmax_conthigh$lowerCI[2],BigFive_cmin_average$lowerCI[2],BigFive_cmax_average$lowerCI[2], BigFive_cmin_asshigh$lowerCI[2],BigFive_cmax_asshigh$lowerCI[2],BigFive_cmin_assveryhigh$lowerCI[2],BigFive_cmax_assveryhigh$lowerCI[2]) upperCI <- c(BigFive_cmin_contveryhigh$upperCI[2],BigFive_cmax_contveryhigh$upperCI[2],BigFive_cmin_conthigh$upperCI[2],BigFive_cmax_conthigh$upperCI[2],BigFive_cmin_average$upperCI[2],BigFive_cmax_average$upperCI[2], BigFive_cmin_asshigh$upperCI[2],BigFive_cmax_asshigh$upperCI[2],BigFive_cmin_assveryhigh$upperCI[2],BigFive_cmax_assveryhigh$upperCI[2]) dat1 <- data.frame(cnt,ass_cont,zPE,lowerCI,upperCI) rm(cnt,ass_cont,zPE,lowerCI,upperCI) dat1$cnt <- factor(dat1$cnt, levels=unique(as.character(dat1$cnt))) dat1$ass_cont <- factor(dat1$ass_cont, levels=unique(as.character(dat1$ass_cont))) BigFivePlot <- ggplot(data=dat1, aes(fill=cnt, x=ass_cont, y=zPE, color=cnt, group=cnt)) + geom_bar(stat="identity", position=position_dodge(), width=0.9, color="black") + geom_errorbar(aes(ymin=lowerCI, ymax=upperCI), position=position_dodge(width=0.9), width=0.2, size=0.5, color="black") + ylim(-0.075, 0.33) + scale_fill_manual(values=c("#cd853f", "#2E9AFE")) + theme(axis.title.x = element_blank(), legend.position=c(0.15, 0.8), legend.title=element_blank(), axis.text = element_text(color="black", size=12), legend.text = element_text(size=12), axis.title.y = element_text(size=14), panel.background = element_blank(), panel.grid = element_blank(), axis.line = element_line(colour = "black")) + geom_hline(yintercept=0) + ggtitle("Figure SOM-U1b. Big Five") + ylab("religiosity-esteem association (zPE)") figure <- grid.arrange(BigTwoPlot,BigFivePlot) ggsave("Figure_SOM-U1.jpg", plot=figure, device="jpg", width=210, height=210, units="mm", dpi=720) ggsave("Figure_SOM-U1.pdf", plot=figure, device="pdf", width=210, height=210, units="mm", dpi=720) rm(dat1,BigTwoPlot,BigFivePlot,figure) #size of the person-culture match effect in each model (power of culture) delta_zPE_b2ts <- BigTwoTS_cmax_average[2,2] - BigTwoTS_cmin_average[2,2] lowerCI_b2ts <- BigTwoTS_cmax_average[2,3] - BigTwoTS_cmin_average[2,4] upperCI_b2ts <- BigTwoTS_cmax_average[2,4] - BigTwoTS_cmin_average[2,3] delta_zPE_b2ac <- BigTwoAC_cmax_average[2,2] - BigTwoAC_cmin_average[2,2] lowerCI_b2ac <- BigTwoAC_cmax_average[2,3] - BigTwoAC_cmin_average[2,4] upperCI_b2ac <- BigTwoAC_cmax_average[2,4] - BigTwoAC_cmin_average[2,3] delta_zPE_b2bf <- BigTwoBF_cmax_average[2,2] - BigTwoBF_cmin_average[2,2] lowerCI_b2bf <- BigTwoBF_cmax_average[2,3] - BigTwoBF_cmin_average[2,4] upperCI_b2bf <- BigTwoBF_cmax_average[2,4] - BigTwoBF_cmin_average[2,3] delta_zPE_b5 <- BigFive_cmax_average[2,2] - BigFive_cmin_average[2,2] lowerCI_b5 <- BigFive_cmax_average[2,3] - BigFive_cmin_average[2,4] upperCI_b5 <- BigFive_cmax_average[2,4] - BigFive_cmin_average[2,3] paste0("Big Two TS Model: delta_zPE = ", round(delta_zPE_b2ts, digits=2), " 95% CI [", round(lowerCI_b2ts, digits=2), ", ", round(upperCI_b2ts, digits=2), "]") paste0("Big Two AC Model: delta_zPE = ", round(delta_zPE_b2ac, digits=2), " 95% CI [", round(lowerCI_b2ac, digits=2), ", ", round(upperCI_b2ac, digits=2), "]") paste0("Big Two BF Model: delta_zPE = ", round(delta_zPE_b2bf, digits=2), " 95% CI [", round(lowerCI_b2bf, digits=2), ", ", round(upperCI_b2bf, digits=2), "]") paste0("Big Five Model: delta_zPE = ", round(delta_zPE_b5, digits=2), " 95% CI [", round(lowerCI_b5, digits=2), ", ", round(upperCI_b5, digits=2), "]") rm(delta_zPE_b2ts,lowerCI_b2ts,upperCI_b2ts,delta_zPE_b2ac,lowerCI_b2ac,upperCI_b2ac,delta_zPE_b2bf,lowerCI_b2bf,upperCI_b2bf,delta_zPE_b5,lowerCI_b5,upperCI_b5) #extent to which personality alters the person-culture match effect in each model (power of personality) ##BIG TWO TS delta_zPE_ass <- BigTwoTS_cmax_assveryhigh$zPe[2] - BigTwoTS_cmin_assveryhigh$zPe[2] lowerCI_ass <- BigTwoTS_cmax_assveryhigh$lowerCI[2] - BigTwoTS_cmin_assveryhigh$upperCI[2] upperCI_ass <- BigTwoTS_cmax_assveryhigh$upperCI[2] - BigTwoTS_cmin_assveryhigh$lowerCI[2] delta_zPE_cont <- BigTwoTS_cmax_contveryhigh$zPe[2] - BigTwoTS_cmin_contveryhigh$zPe[2] lowerCI_cont <- BigTwoTS_cmax_contveryhigh$lowerCI[2] - BigTwoTS_cmin_contveryhigh$upperCI[2] upperCI_cont <- BigTwoTS_cmax_contveryhigh$upperCI[2] - BigTwoTS_cmin_contveryhigh$lowerCI[2] paste0("Big Two TS Model: strong assimilators: delta_zPE = ", round(delta_zPE_ass, digits=2), " 95% CI [", round(lowerCI_ass, digits=2), ", ", round(upperCI_ass, digits=2), "]") paste0("Big Two TS Model: strong contrasters: delta_zPE = ", round(delta_zPE_cont, digits=2), " 95% CI [", round(lowerCI_cont, digits=2), ", ", round(upperCI_cont, digits=2), "]") power_personality <- round(delta_zPE_ass, digits=2) - round(delta_zPE_cont, digits=2) power_lowerCI <- round(lowerCI_ass, digits=2) - round(upperCI_cont, digits=2) power_upperCI <- round(upperCI_ass, digits=2) - round(lowerCI_cont, digits=2) paste0("Power of Big Two TS: delta_zPE = ", power_personality, " 95% CI [", power_lowerCI, ", ", power_upperCI, "]") ##BIG TWO AC delta_zPE_ass <- BigTwoAC_cmax_assveryhigh$zPe[2] - BigTwoAC_cmin_assveryhigh$zPe[2] lowerCI_ass <- BigTwoAC_cmax_assveryhigh$lowerCI[2] - BigTwoAC_cmin_assveryhigh$upperCI[2] upperCI_ass <- BigTwoAC_cmax_assveryhigh$upperCI[2] - BigTwoAC_cmin_assveryhigh$lowerCI[2] delta_zPE_cont <- BigTwoAC_cmax_contveryhigh$zPe[2] - BigTwoAC_cmin_contveryhigh$zPe[2] lowerCI_cont <- BigTwoAC_cmax_contveryhigh$lowerCI[2] - BigTwoAC_cmin_contveryhigh$upperCI[2] upperCI_cont <- BigTwoAC_cmax_contveryhigh$upperCI[2] - BigTwoAC_cmin_contveryhigh$lowerCI[2] paste0("Big Two AC Model: strong assimilators: delta_zPE = ", round(delta_zPE_ass, digits=2), " 95% CI [", round(lowerCI_ass, digits=2), ", ", round(upperCI_ass, digits=2), "]") paste0("Big Two AC Model: strong contrasters: delta_zPE = ", round(delta_zPE_cont, digits=2), " 95% CI [", round(lowerCI_cont, digits=2), ", ", round(upperCI_cont, digits=2), "]") power_personality <- round(delta_zPE_ass, digits=2) - round(delta_zPE_cont, digits=2) power_lowerCI <- round(lowerCI_ass, digits=2) - round(upperCI_cont, digits=2) power_upperCI <- round(upperCI_ass, digits=2) - round(lowerCI_cont, digits=2) paste0("Power of Big Two AC: delta_zPE = ", power_personality, " 95% CI [", power_lowerCI, ", ", power_upperCI, "]") ##BIG TWO BF delta_zPE_ass <- BigTwoBF_cmax_assveryhigh$zPe[2] - BigTwoBF_cmin_assveryhigh$zPe[2] lowerCI_ass <- BigTwoBF_cmax_assveryhigh$lowerCI[2] - BigTwoBF_cmin_assveryhigh$upperCI[2] upperCI_ass <- BigTwoBF_cmax_assveryhigh$upperCI[2] - BigTwoBF_cmin_assveryhigh$lowerCI[2] delta_zPE_cont <- BigTwoBF_cmax_contveryhigh$zPe[2] - BigTwoBF_cmin_contveryhigh$zPe[2] lowerCI_cont <- BigTwoBF_cmax_contveryhigh$lowerCI[2] - BigTwoBF_cmin_contveryhigh$upperCI[2] upperCI_cont <- BigTwoBF_cmax_contveryhigh$upperCI[2] - BigTwoBF_cmin_contveryhigh$lowerCI[2] paste0("Big Two BF Model: strong assimilators: delta_zPE = ", round(delta_zPE_ass, digits=2), " 95% CI [", round(lowerCI_ass, digits=2), ", ", round(upperCI_ass, digits=2), "]") paste0("Big Two BF Model: strong contrasters: delta_zPE = ", round(delta_zPE_cont, digits=2), " 95% CI [", round(lowerCI_cont, digits=2), ", ", round(upperCI_cont, digits=2), "]") power_personality <- round(delta_zPE_ass, digits=2) - round(delta_zPE_cont, digits=2) power_lowerCI <- round(lowerCI_ass, digits=2) - round(upperCI_cont, digits=2) power_upperCI <- round(upperCI_ass, digits=2) - round(lowerCI_cont, digits=2) paste0("Power of Big Two BF: delta_zPE = ", power_personality, " 95% CI [", power_lowerCI, ", ", power_upperCI, "]") ##BIG FIVE delta_zPE_ass <- BigFive_cmax_assveryhigh$zPe[2] - BigFive_cmin_assveryhigh$zPe[2] lowerCI_ass <- BigFive_cmax_assveryhigh$lowerCI[2] - BigFive_cmin_assveryhigh$upperCI[2] upperCI_ass <- BigFive_cmax_assveryhigh$upperCI[2] - BigFive_cmin_assveryhigh$lowerCI[2] delta_zPE_cont <- BigFive_cmax_contveryhigh$zPe[2] - BigFive_cmin_contveryhigh$zPe[2] lowerCI_cont <- BigFive_cmax_contveryhigh$lowerCI[2] - BigFive_cmin_contveryhigh$upperCI[2] upperCI_cont <- BigFive_cmax_contveryhigh$upperCI[2] - BigFive_cmin_contveryhigh$lowerCI[2] paste0("Big Five Model: strong assimilators: delta_zPE = ", round(delta_zPE_ass, digits=2), " 95% CI [", round(lowerCI_ass, digits=2), ", ", round(upperCI_ass, digits=2), "]") paste0("Big Five Model: strong contrasters: delta_zPE = ", round(delta_zPE_cont, digits=2), " 95% CI [", round(lowerCI_cont, digits=2), ", ", round(upperCI_cont, digits=2), "]") power_personality <- round(delta_zPE_ass, digits=2) - round(delta_zPE_cont, digits=2) power_lowerCI <- round(lowerCI_ass, digits=2) - round(upperCI_cont, digits=2) power_upperCI <- round(upperCI_ass, digits=2) - round(lowerCI_cont, digits=2) paste0("Power of Big Five: delta_zPE = ", power_personality, " 95% CI [", power_lowerCI, ", ", power_upperCI, "]")