#install packages
install.packages(c("psych","car","flextable","officer","data.table","dplyr","numform","lavaan","ggplot2","ggpubr","effectsize","multilevel","MplusAutomation","lme4","r2mlm","simr"))


#load packages
library(psych)
library(car)
library(flextable)
library(officer)
library(data.table)
library(dplyr)
library(numform)
library(lavaan)
library(ggplot2)
library(ggpubr)
library(effectsize)
library(multilevel)
library(MplusAutomation)
library(lme4)
library(r2mlm)
library(simr)


#save directory to the folder to which files created in R should be saved (e.g., "C:/NarcSub/")
files_wd <- "C:/NarcSub/"


#read data
dat2a <- as.data.frame(fread("https://madata.bib.uni-mannheim.de/427/5/NarcSub_Study2a_Data.csv", header = T, sep = ","))
dat2b <- as.data.frame(fread("https://madata.bib.uni-mannheim.de/427/7/NarcSub_Study2b_Data.csv", header = T, sep = ","))


#apply data-exclusion criteria for Study 2a
dat2a <- subset(dat2a, subset = exclude == 0)


#sample demographics
##Study 2a
###sex
dat2a_1 <- table(dat2a$sex)
length(which(is.na(dat2a$sex))) #no missing values
###years with supervisor
####check for unrealistic values
dat2a$check_yearswithsupervisor <- dat2a$work_experience - dat2a$years_with_supervisor
table(dat2a$check_yearswithsupervisor) #6 values are negative, but all values are still realistic when considering that professional training can be seen as either work experience or education
dat2a[is.na(dat2a["years_with_supervisor"]),"id"] #no missing values
dat2a$years_with_supervisor[dat2a$id==77] #participant with missing value for work experience has realistic value
###print sample demographics
paste0(dat2a_1[2], " women, ", dat2a_1[1], " men; age: ", min(dat2a$age,na.rm=TRUE), "-", max(dat2a$age,na.rm=TRUE)," years, M = " , (round(mean(dat2a$age,na.rm=TRUE),2)), ", SD = ", (round(sd(dat2a$age,na.rm=TRUE),2)))
paste0("Participants had been working with their supervisor between ", min(dat2a$years_with_supervisor), " and ", max(dat2a$years_with_supervisor), " years (M = ", (round(mean(dat2a$years_with_supervisor),2)), ", SD = ", (round(sd(dat2a$years_with_supervisor),2)), ")")
paste0("Most participants (", round(length(which(dat2a$contact_with_supervisor >= 3))/length(dat2a$id)*100,0), "%) reported seeing their supervisor at least once a week.")
rm(dat2a_1)
##Study 2b
###sex
dat2b_1 <- table(dat2b$sex)
length(which(is.na(dat2b$sex))) #one missing value
###years with supervisor
####check for unrealistic values
dat2b$check_yearswithsupervisor <- dat2b$work_experience - dat2b$years_with_supervisor
table(dat2b$check_yearswithsupervisor) #21 values are negative, but 19 values are still realistic when considering that professional training can be seen as either work experience or education
dat2b[is.na(dat2b["years_with_supervisor"]),"id"] #no missing values
dat2b$years_with_supervisor[dat2b$id==49] #participant with missing value for work experience has realistic value
dat2b$years_with_supervisor[dat2b$id==70] #participant with missing value for work experience has realistic value
dat2b$years_with_supervisor[dat2b$id==120] #participant with missing value for work experience has realistic value
####exclude participants with unrealistic values
dat2b$id[dat2b$check_yearswithsupervisor == -20]
dat2b$id[dat2b$check_yearswithsupervisor == -32]
dat2b_2 <- subset(dat2b, subset = c(id != 164 & id != 60))
###team size
dat2b_3 <- dat2b %>% count(team)
###print sample demographics
paste0(dat2b_1[2], " women, ", dat2b_1[1], " men, ", dat2b_1[3]+length(which(is.na(dat2b$sex)))," undisclosed; age: ", min(dat2b$age,na.rm=TRUE), "-", max(dat2b$age,na.rm=TRUE)," years, M = " , (round(mean(dat2b$age,na.rm=TRUE),2)), ", SD = ", (round(sd(dat2b$age,na.rm=TRUE),2)))
paste0(length(dat2b_3$team)," teams with ", min(dat2b_3$n), " to ", max(dat2b_3$n)," team members (M = ", (round(mean(dat2b_3$n),2)), ", SD = ", (round(sd(dat2b_3$n),2)), ")")
paste0("Except for ", length(which(dat2b$years_with_supervisor==0)) , " participants (", round(length(which(dat2b$years_with_supervisor==0))/nrow(dat2b),2)*100, "%), they had been working with their supervisor for at least one year (1-", max(dat2b_2$years_with_supervisor), " years, M = ", (round(mean(dat2b_2$years_with_supervisor),2)), ", SD = ", (round(sd(dat2b_2$years_with_supervisor),2)), ")")
paste0("Most participants (", round(length(which(dat2b$contact_with_supervisor >= 4))/length(dat2b$id)*100,0), "%) reported seeing their supervisor at least once a week.")
rm(dat2b_1,dat2b_2,dat2b_3)


#recode variables
##Study 2a
dat2a$description_aut_supervisor_r <- car::recode(dat2a$description_aut_supervisor, "1=7; 2=6; 3=5; 4=4; 5=3; 6=2; 7=1")
dat2a$lbdq2r_supervisor <- car::recode(dat2a$lbdq2_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq4r_supervisor <- car::recode(dat2a$lbdq4_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq6r_supervisor <- car::recode(dat2a$lbdq6_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq8r_supervisor <- car::recode(dat2a$lbdq8_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq10r_supervisor <- car::recode(dat2a$lbdq10_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq12r_supervisor <- car::recode(dat2a$lbdq12_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq14r_supervisor <- car::recode(dat2a$lbdq14_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq16r_supervisor <- car::recode(dat2a$lbdq16_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq18r_supervisor <- car::recode(dat2a$lbdq18_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq20r_supervisor <- car::recode(dat2a$lbdq20_supervisor, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$description_aut_self_r <- car::recode(dat2a$description_aut_self, "1=7; 2=6; 3=5; 4=4; 5=3; 6=2; 7=1")
dat2a$lbdq2r_self <- car::recode(dat2a$lbdq2_self, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq4r_self <- car::recode(dat2a$lbdq4_self, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq6r_self <- car::recode(dat2a$lbdq6_self, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq8r_self <- car::recode(dat2a$lbdq8_self, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq10r_self <- car::recode(dat2a$lbdq10_self, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq12r_self <- car::recode(dat2a$lbdq12_self, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq14r_self <- car::recode(dat2a$lbdq14_self, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq16r_self <- car::recode(dat2a$lbdq16_self, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq18r_self <- car::recode(dat2a$lbdq18_self, "1=5; 2=4; 3=3; 4=2; 5=1")
dat2a$lbdq20r_self <- car::recode(dat2a$lbdq20_self, "1=5; 2=4; 3=3; 4=2; 5=1")
##Study 2b
dat2b$description_aut_supervisor_r <- car::recode(dat2b$description_aut_supervisor, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq2r_supervisor <- car::recode(dat2b$lbdq2_supervisor, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq4r_supervisor <- car::recode(dat2b$lbdq4_supervisor, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq6r_supervisor <- car::recode(dat2b$lbdq6_supervisor, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq8r_supervisor <- car::recode(dat2b$lbdq8_supervisor, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq10r_supervisor <- car::recode(dat2b$lbdq10_supervisor, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$description_aut_self_r <- car::recode(dat2b$description_aut_self, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq2r_self <- car::recode(dat2b$lbdq2_self, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq4r_self <- car::recode(dat2b$lbdq4_self, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq6r_self <- car::recode(dat2b$lbdq6_self, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq8r_self <- car::recode(dat2b$lbdq8_self, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")
dat2b$lbdq10r_self <- car::recode(dat2b$lbdq10_self, "1=6; 2=5; 3=4; 4=3; 5=2; 6=1")


#compute variables 
##Study 2a
###predictors
dat2a$narq_adm <- rowMeans(dat2a[c("narq1", "narq2", "narq3", "narq5", "narq7", "narq8", "narq15", "narq16", "narq18")], na.rm=TRUE)
dat2a$Znarq_adm <- scale(dat2a$narq_adm)
dat2a$narq_riv <- rowMeans(dat2a[c("narq4", "narq6", "narq9", "narq10", "narq11", "narq12", "narq13", "narq14", "narq17")], na.rm=TRUE)
dat2a$Znarq_riv <- scale(dat2a$narq_riv)
###moderator
dat2a$descriptions_supervisor <- rowMeans(dat2a[c("description_dem_supervisor","description_aut_supervisor_r")], na.rm = TRUE) 
dat2a$Zdescriptions_supervisor <- scale(dat2a$descriptions_supervisor)
dat2a$lbdq_supervisor <- rowMeans(dat2a[c("lbdq1_supervisor","lbdq3_supervisor","lbdq5_supervisor","lbdq7_supervisor","lbdq9_supervisor","lbdq11_supervisor","lbdq13_supervisor", 
                                          "lbdq15_supervisor","lbdq17_supervisor","lbdq19_supervisor","lbdq2r_supervisor","lbdq4r_supervisor","lbdq6r_supervisor","lbdq8r_supervisor",
                                          "lbdq10r_supervisor","lbdq12r_supervisor","lbdq14r_supervisor","lbdq16r_supervisor","lbdq18r_supervisor","lbdq20r_supervisor")], na.rm = TRUE)
dat2a$Zlbdq_supervisor <- scale(dat2a$lbdq_supervisor)
dat2a$style_index_supervisor <- rowMeans(dat2a[c("Zdescriptions_supervisor","Zlbdq_supervisor")], na.rm = TRUE)
dat2a$Zstyle_index_supervisor <- scale(dat2a$style_index_supervisor)
###mediator
dat2a$descriptions_self <- rowMeans(dat2a[c("description_dem_self","description_aut_self_r")], na.rm=TRUE)
dat2a$Zdescriptions_self <- scale(dat2a$descriptions_self)
dat2a$lbdq_self <- rowMeans(dat2a[c("lbdq1_self","lbdq3_self","lbdq5_self","lbdq7_self","lbdq9_self","lbdq11_self","lbdq13_self", 
                                    "lbdq15_self","lbdq17_self","lbdq19_self","lbdq2r_self","lbdq4r_self","lbdq6r_self","lbdq8r_self",
                                    "lbdq10r_self","lbdq12r_self","lbdq14r_self","lbdq16r_self","lbdq18r_self","lbdq20r_self")], na.rm=TRUE)
dat2a$Zlbdq_self <- scale(dat2a$lbdq_self)
dat2a$style_index_self <- rowMeans(dat2a[c("Zdescriptions_self","Zlbdq_self")], na.rm = TRUE)
dat2a$Zstyle_index_self <- scale(dat2a$style_index_self)
###interactions
dat2a$ZADMxSLS <- scale(dat2a$Znarq_adm*dat2a$Zstyle_index_supervisor)
dat2a$ZRIVxSLS <- scale(dat2a$Znarq_riv*dat2a$Zstyle_index_supervisor)
dat2a$ZSHLSxSLS <- scale(dat2a$Zstyle_index_self*dat2a$Zstyle_index_supervisor)
###dependent variable
dat2a$zvalence <- scale(rowMeans(dat2a[c("valence1","valence2","valence3")], na.rm=TRUE))
dat2a$zsuccess <- scale(rowMeans(dat2a[c("success1","success2","success3")], na.rm=TRUE))
dat2a$boss_desirability3r <- 8 - dat2a$boss_desirability3
dat2a$zboss_desirability <- scale(rowMeans(dat2a[c("boss_desirability1","boss_desirability2","boss_desirability3r")], na.rm=TRUE))
dat2a$zliking <- scale(rowMeans(dat2a[c("liking1","liking2","liking3")], na.rm=TRUE))
dat2a$zreputation <- scale(rowMeans(dat2a[c("reputation1","reputation2")], na.rm=TRUE))
dat2a$zcompetence <- scale(rowMeans(dat2a[c("competence1","competence2","competence3")], na.rm=TRUE))
dat2a$zwarmth <- scale(rowMeans(dat2a[c("warmth1","warmth2","warmth3")], na.rm=TRUE))
dat2a$zworking_conditions <- scale(rowMeans(dat2a[c("working_conditions1","working_conditions2","working_conditions3")], na.rm=TRUE))
dat2a$attraction <- rowMeans(dat2a[c("zvalence","zsuccess","zboss_desirability","zliking","zreputation","zcompetence","zwarmth","zworking_conditions")], na.rm=TRUE)
dat2a$Zattraction <- scale(dat2a$attraction)
##Study 2b
###predictors
dat2b$narq_adm <- rowMeans(dat2b[c("narq1", "narq2", "narq3", "narq5", "narq7", "narq8", "narq15", "narq16", "narq18")], na.rm=TRUE)
dat2b$Znarq_adm <- scale(dat2b$narq_adm)
dat2b$narq_riv <- rowMeans(dat2b[c("narq4", "narq6", "narq9", "narq10", "narq11", "narq12", "narq13", "narq14", "narq17")], na.rm=TRUE)
dat2b$Znarq_riv <- scale(dat2b$narq_riv)
###moderator
dat2b$descriptions_supervisor <- rowMeans(dat2b[c("description_dem_supervisor","description_aut_supervisor_r")], na.rm = TRUE)
dat2b$lbdq_supervisor <- rowMeans(dat2b[c("lbdq1_supervisor","lbdq3_supervisor","lbdq5_supervisor","lbdq7_supervisor","lbdq9_supervisor",
                                          "lbdq2r_supervisor","lbdq4r_supervisor","lbdq6r_supervisor","lbdq8r_supervisor","lbdq10r_supervisor")], na.rm = TRUE)
dat2b_team <- aggregate(cbind(descriptions_supervisor,lbdq_supervisor) ~ team, dat2b, mean)
colnames(dat2b_team) <- c("team","descriptions_supervisor_team","lbdq_supervisor_team")
dat2b_team$Zdescriptions_supervisor_team <- scale(dat2b_team$descriptions_supervisor_team)
dat2b_team$Zlbdq_supervisor_team <- scale(dat2b_team$lbdq_supervisor_team)
dat2b_team$style_index_supervisor_team <- rowMeans(dat2b_team[c("Zdescriptions_supervisor_team","Zlbdq_supervisor_team")], na.rm=TRUE)
dat2b_team$Zstyle_index_supervisor_team <- scale(dat2b_team$style_index_supervisor_team)
dat2b <- merge(dat2b,dat2b_team, by = "team")
###mediator
dat2b$descriptions_self <- rowMeans(dat2b[c("description_dem_self","description_aut_self_r")], na.rm=TRUE) 
dat2b$Zdescriptions_self <- scale(dat2b$descriptions_self)
dat2b$lbdq_self <- rowMeans(dat2b[c("lbdq1_self","lbdq3_self","lbdq5_self","lbdq7_self","lbdq9_self","lbdq2r_self",
                                    "lbdq4r_self","lbdq6r_self","lbdq8r_self","lbdq10r_self")], na.rm=TRUE) 
dat2b$Zlbdq_self <- scale(dat2b$lbdq_self)
dat2b$style_index_self <- rowMeans(dat2b[c("Zdescriptions_self","Zlbdq_self")], na.rm = TRUE)
dat2b$Zstyle_index_self <- scale(dat2b$style_index_self)
###interactions
dat2b$ZADMxSLS <- scale(dat2b$Znarq_adm*dat2b$Zstyle_index_supervisor_team)
dat2b$ZRIVxSLS <- scale(dat2b$Znarq_riv*dat2b$Zstyle_index_supervisor_team)
dat2b$ZSHLSxSLS <- scale(dat2b$Zstyle_index_self*dat2b$Zstyle_index_supervisor_team)
###dependent variable
dat2b$zvalence <- scale(rowMeans(dat2b[c("valence1","valence2")], na.rm=TRUE))
dat2b$zsuccess <- scale(rowMeans(dat2b[c("success1","success2")], na.rm=TRUE))
dat2b$boss_desirability2r <- 8 - dat2b$boss_desirability2
dat2b$zboss_desirability <- scale(rowMeans(dat2b[c("boss_desirability1","boss_desirability2r")], na.rm=TRUE))
dat2b$zliking <- scale(rowMeans(dat2b[c("liking1","liking2")], na.rm=TRUE))
dat2b$zcompetence <- scale(rowMeans(dat2b[c("competence1","competence2")], na.rm=TRUE))
dat2b$zwarmth <- scale(rowMeans(dat2b[c("warmth1","warmth2")], na.rm=TRUE))
dat2b$attraction <- rowMeans(dat2b[c("zvalence","zsuccess","zboss_desirability","zliking","zcompetence","zwarmth")], na.rm=TRUE)
dat2b$Zattraction <- scale(dat2b$attraction)


#internal consistencies of variables
##Study 2a
dat2a.narq_adm <- subset(dat2a, select = c(narq1,narq2,narq3,narq5,narq7,narq8,narq15,narq16,narq18))
dat2a.narq_riv <- subset(dat2a, select = c(narq4,narq6,narq9,narq10,narq11,narq12,narq13,narq14,narq17))
dat2a.style_index_supervisor <- subset(dat2a, select = c(Zdescriptions_supervisor,Zlbdq_supervisor))
dat2a.style_index_self <- subset(dat2a, select = c(Zdescriptions_self,Zlbdq_self))
dat2a.attraction <- subset(dat2a, select = c(zvalence, zsuccess, zboss_desirability, zliking, zreputation, zcompetence, zwarmth, zworking_conditions))
round(psych::alpha(dat2a.narq_adm, na.rm=T)$total[1], 2)
round(psych::alpha(dat2a.narq_riv, na.rm=T)$total[1], 2)
round(psych::alpha(dat2a.style_index_supervisor, na.rm=T)$total[1], 2)
round(psych::alpha(dat2a.style_index_self, na.rm=T)$total[1], 2)
round(psych::alpha(dat2a.attraction, na.rm=T)$total[1], 2)
##Study 2b
dat2b.narq_adm <- subset(dat2b, select = c(narq1,narq2,narq3,narq5,narq7,narq8,narq15,narq16,narq18))
dat2b.narq_riv <- subset(dat2b, select = c(narq4,narq6,narq9,narq10,narq11,narq12,narq13,narq14,narq17))
dat2b.style_index_supervisor_team <- subset(dat2b_team, select = c(Zdescriptions_supervisor_team,Zlbdq_supervisor_team))
dat2b.style_index_self <- subset(dat2b, select = c(Zdescriptions_self,Zlbdq_self))
dat2b.attraction <- subset(dat2b, select = c(zvalence, zsuccess, zboss_desirability, zliking, zcompetence, zwarmth))
round(psych::alpha(dat2b.narq_adm, na.rm=T)$total[1], 2)
round(psych::alpha(dat2b.narq_riv, na.rm=T)$total[1], 2)
round(ICC2(aov(lbdq_supervisor ~ as.factor(team),data=dat2b)),2)
round(ICC2(aov(descriptions_supervisor ~ as.factor(dat2b$team),data = dat2b)),2)
round(psych::alpha(dat2b.style_index_supervisor_team, na.rm=T)$total[1], 2)
round(psych::alpha(dat2b.style_index_self, na.rm=T)$total[1], 2)
round(psych::alpha(dat2b.attraction, na.rm=T)$total[1], 2)
rm(dat2a.narq_adm,dat2a.narq_riv,dat2a.style_index_supervisor,dat2a.style_index_self,
   dat2b_team,dat2b.narq_adm,dat2b.narq_riv,dat2b.style_index_supervisor_team,dat2b.style_index_self)


#Footnote 10
round(ICC1(aov(lbdq_supervisor ~ as.factor(team),data=dat2b)),2)
round(ICC1(aov(descriptions_supervisor ~ as.factor(dat2b$team),data = dat2b)),2)
rwg_lbdq <- as.data.frame(rwg(dat2b$lbdq_supervisor, dat2b$team))
round(mean(rwg_lbdq$rwg),2)
rwg_descriptions <- as.data.frame(rwg(dat2b$descriptions_supervisor, dat2b$team))
round(mean(rwg_descriptions$rwg),2)
rm(rwg_lbdq,rwg_descriptions)


#principal component analysis for attraction to supervisor
##Study2a
scree(dat2a.attraction)
dat2a.model <- principal(dat2a.attraction, nfactors = 1)
prop.table(dat2a.model$values)[1]*100 #factor explains 71.17% of the total variance
dat2a.model #all measures manifest loadings of at least .68
##Study2b
scree(dat2b.attraction)
dat2b.model <- principal(dat2b.attraction, nfactors = 1)
prop.table(dat2b.model$values)[1]*100 #factor explains 75.64% of the total variance
dat2b.model #all measures manifest loadings of at least .80
rm(dat2a.attraction,dat2b.attraction,dat2a.model,dat2b.model)


#Table 3: Descriptive Statistics and Zero-Order Correlations in Studies 2a-2b
dat_table3 <- data.frame(matrix(nrow=5, ncol=5))
dat_table3[,1] <- c("ADM","RIV","SHLS","SLS","Attraction")
dat_table3[,2] <- c(mean(dat2a$narq_adm, na.rm=T),mean(dat2a$narq_riv, na.rm=T),mean(dat2a$style_index_self, na.rm=T), mean(dat2a$style_index_supervisor, na.rm=T), mean(dat2a$attraction, na.rm=T))
dat_table3[,3] <- c(sd(dat2a$narq_adm, na.rm=T),sd(dat2a$narq_riv, na.rm=T),sd(dat2a$style_index_self, na.rm=T), sd(dat2a$style_index_supervisor, na.rm=T), sd(dat2a$attraction, na.rm=T))
dat_table3[,4] <- c(mean(dat2b$narq_adm, na.rm=T),mean(dat2b$narq_riv, na.rm=T),mean(dat2b$style_index_self, na.rm=T), mean(dat2b$style_index_supervisor_team, na.rm=T), mean(dat2b$attraction, na.rm=T))
dat_table3[,5] <- c(sd(dat2b$narq_adm, na.rm=T),sd(dat2b$narq_riv, na.rm=T),sd(dat2b$style_index_self, na.rm=T), sd(dat2b$style_index_supervisor_team, na.rm=T), sd(dat2b$attraction, na.rm=T))
dat2a_1 <- subset(dat2a, select = c(narq_adm,narq_riv,style_index_self,style_index_supervisor,attraction))
cor_dat2a_1 <- corr.test(dat2a_1)
dat_cor_2a <- as.data.frame(matrix(unlist(cor_dat2a_1), nrow=250, ncol = 1, byrow=F), stringsAsFactors=FALSE)
dat_cor_2a_r <- as.data.frame(matrix(dat_cor_2a[1:25,], nrow=5))
dat_cor_2a_p <- as.data.frame(matrix(dat_cor_2a[52:76,], nrow=5))
dat2b_1 <- subset(dat2b, select = c(narq_adm,narq_riv,style_index_self,style_index_supervisor_team,attraction))
cor_dat2b_1 <- corr.test(dat2b_1)
dat_cor_2b <- as.data.frame(matrix(unlist(cor_dat2b_1), nrow=250, ncol = 1, byrow=F), stringsAsFactors=FALSE)
dat_cor_2b_r <- as.data.frame(matrix(dat_cor_2b[1:25,], nrow=5))
dat_cor_2b_p <- as.data.frame(matrix(dat_cor_2b[52:76,], nrow=5))
dat_cor_r <- cbind(dat_cor_2a_r,dat_cor_2b_r)
dat_cor_p <- cbind(dat_cor_2a_p,dat_cor_2b_p)

dat_table3[,2:5] <- sapply(dat_table3[,2:5], as.numeric)
format_msd <- function(x){
  formatC(x, format = "f", digits = 2)
}
dat_table3[,2:5] <- dat_table3[,2:5] %>% mutate_if(is.numeric, format_msd)
dat_table3[3,2] <- "0.00"

dat_cor_r[,1:10] <- sapply(dat_cor_r[,1:10], as.numeric)
format_r <- function(x){
  ifelse(abs(x) < 0.005,  f_num(unlist(x), digits=3), f_num(unlist(x), digits=2))
}
dat_cor_r[,1:10] <- dat_cor_r[,1:10] %>% mutate_if(is.numeric, format_r)

dat_cor_p[,1:10] <- sapply(dat_cor_p[,1:10], as.numeric)
format_p <- function(x){
  ifelse(abs(x) < 0.05 & abs(x) >= 0.01,  "*", ifelse(abs(x) < 0.01 & abs(x) >= 0.001,  "**", ifelse(abs(x) < 0.001,  "***", "")))
}
dat_cor_p[,1:10] <- dat_cor_p[,1:10] %>% mutate_if(is.numeric, format_p)

dat_table3$r1 <- paste0(dat_cor_r[,1],dat_cor_p[,1])
dat_table3$r2 <- paste0(dat_cor_r[,2],dat_cor_p[,2])
dat_table3$r3 <- paste0(dat_cor_r[,3],dat_cor_p[,3])
dat_table3$r4 <- paste0(dat_cor_r[,4],dat_cor_p[,4])
dat_table3$r5 <- paste0(dat_cor_r[,5],dat_cor_p[,5])
dat_table3[1,6] <- NA
dat_table3[1,7] <- paste0(dat_cor_r[2,6],dat_cor_p[2,6])
dat_table3[1,8] <- paste0(dat_cor_r[3,6],dat_cor_p[3,6])
dat_table3[1,9] <- paste0(dat_cor_r[4,6],dat_cor_p[4,6])
dat_table3[1,10] <- paste0(dat_cor_r[5,6],dat_cor_p[5,6])
dat_table3[2,7] <- NA
dat_table3[2,8] <- paste0(dat_cor_r[3,7],dat_cor_p[3,7])
dat_table3[2,9] <- paste0(dat_cor_r[4,7],dat_cor_p[4,7])
dat_table3[2,10] <- paste0(dat_cor_r[5,7],dat_cor_p[5,7])
dat_table3[3,8] <- NA
dat_table3[3,9] <- paste0(dat_cor_r[4,8],dat_cor_p[4,8])
dat_table3[3,10] <- paste0(dat_cor_r[5,8],dat_cor_p[5,8])
dat_table3[4,9] <- NA
dat_table3[4,10] <- paste0(dat_cor_r[5,9],dat_cor_p[5,9])
dat_table3[5,10] <- NA

dat_table3$blank1 <- NA
dat_table3$blank2 <- NA
dat_table3 <- dat_table3[,c(1:3,11,4:5,12,6:10)]

col_keys <- c("X1","X2","X3","blank1","X4","X5","blank2","r1","r2","r3","r4","r5")
head1 <- c("","Study 2a","Study 2a","","Study 2b","Study 2b","","","","","","")
head2 <- c("Variable","M","SD","","M","SD","","1","2","3","4","5")
head <- data.frame(col_keys,head1,head2, stringsAsFactors = FALSE)
rm(col_keys,head1,head2)

tbl <- flextable(dat_table3)
tbl <- set_header_df(tbl, mapping=head, key="col_keys")
tbl <- hline_top(tbl, j=1:12, border=fp_border(width=2), part="header")
tbl <- merge_at(tbl, i=1, j=2:3, part="header")
tbl <- merge_at(tbl, i=1, j=5:6, part="header")
tbl <- hline(tbl, i=1, j=c(2:3,5:6), border=fp_border(width=1.2), part="header")
tbl <- hline(tbl, i=2, j=1:12, border=fp_border(width=1.2), part="header")
tbl <- flextable::font(tbl, fontname="Times", part="all")
tbl <- fontsize(tbl, size=12, part="all")
tbl <- align(tbl, align="center", part="all")
tbl <- align(tbl, j = c("X1"), align="left", part="body")
tbl <- italic(tbl, i=2, j=c("X2","X3","X4","X5"), part="header")
tbl <- width(tbl, j =~ X1, width=1)
tbl <- width(tbl, j =~ X2 + X3 + X4 + X5, width=.55)
tbl <- width(tbl, j =~ blank1 + blank2, width=.1)
tbl

setwd(files_wd)
doc <- read_docx()
doc <- body_add_flextable(doc, value = tbl)
print(doc, target = "Table_3.docx")

rm(cor_dat2a_1,cor_dat2b_1, dat_cor_2a,dat_cor_2a_p,dat_cor_2a_r,dat_cor_2b,dat_cor_2b_p,dat_cor_2b_r,
   dat_cor_p,dat_cor_r,dat_table3, dat2a_1, dat2b_1,doc,head,tbl, format_msd, format_p, format_r)


#Note of Table 3
round(mean(dat2a$descriptions_self),2)
round(sd(dat2a$descriptions_self),2)
round(mean(dat2a$lbdq_self),2)
round(sd(dat2a$lbdq_self),2)
round(mean(dat2a$descriptions_supervisor),2)
round(sd(dat2a$descriptions_supervisor),2)
round(mean(dat2a$lbdq_supervisor),2)
round(sd(dat2a$lbdq_supervisor),2)
round(mean(dat2b$descriptions_self),2)
round(sd(dat2b$descriptions_self),2)
round(mean(dat2b$lbdq_self),2)
round(sd(dat2b$lbdq_self),2)
round(mean(dat2b$descriptions_supervisor_team),2)
round(sd(dat2b$descriptions_supervisor_team),2)
round(mean(dat2b$lbdq_supervisor_team),2)
round(sd(dat2b$lbdq_supervisor_team),2)
t.test(dat2a$descriptions_self, mu = 4, alternative = "two.sided")
cohens_d(dat2a$descriptions_self, mu=4)
t.test(dat2a$lbdq_self, mu = 3, alternative = "two.sided")
cohens_d(dat2a$lbdq_self, mu=3)
t.test(dat2a$descriptions_supervisor, mu = 4, alternative = "two.sided")
cohens_d(dat2a$descriptions_supervisor, mu=4)
t.test(dat2a$lbdq_supervisor, mu = 3, alternative = "two.sided")
cohens_d(dat2a$lbdq_supervisor, mu=3)
t.test(dat2b$descriptions_self, mu = 3.5, alternative = "two.sided")
cohens_d(dat2b$descriptions_self, mu=3.5)
t.test(dat2b$lbdq_self, mu = 3.5, alternative = "two.sided")
cohens_d(dat2b$lbdq_self, mu=3.5)
t.test(dat2b$descriptions_supervisor_team, mu = 3.5, alternative = "two.sided")
cohens_d(dat2b$descriptions_supervisor_team, mu=3.5)
t.test(dat2b$lbdq_supervisor_team, mu = 3.5, alternative = "two.sided")
cohens_d(dat2b$lbdq_supervisor_team, mu=3.5)


#Table 4: Unique Effects of Narcissistic Admiration and Rivalry on Attraction to Supervisor Moderated by Supervisor's Leadership Style in Studies 2a-2b

##Study 2a
###regression model
moderation_model2a <- lm(Zattraction ~ Znarq_adm + Znarq_riv + Zstyle_index_supervisor + ZADMxSLS + ZRIVxSLS, data=dat2a)
dat2a_mod <- data.frame(summary(moderation_model2a)$coef,confint(moderation_model2a))
moderation_model2a_rsquared_p <- pf(summary(moderation_model2a)$fstatistic[1],summary(moderation_model2a)$fstatistic[2],summary(moderation_model2a)$fstatistic[3],lower.tail=FALSE)
moderation_model2a_rsquared_p <- ifelse(moderation_model2a_rsquared_p < 0.001,  "< .001", paste0("= ", substr(round(moderation_model2a_rsquared_p, 3),2,5)))
moderation_model2a_rsquared <- paste0("R2 = ",substr(round(summary(moderation_model2a)$r.squared,2),2,4),", F(",summary(moderation_model2a)$fstatistic[2],", ",summary(moderation_model2a)$fstatistic[3],") = ",round(summary(moderation_model2a)$fstatistic[1],2),", p ",moderation_model2a_rsquared_p)

###effect size of interaction between narcissistic rivalry and supervisor’s leadership style 
moderation_model2a_no_RIVxSLS <- lm(Zattraction ~ Znarq_adm + Znarq_riv + Zstyle_index_supervisor + ZADMxSLS, data=dat2a)
cohens_f_squared_RIVxSLS <- cohens_f_squared(moderation_model2a, model2 = moderation_model2a_no_RIVxSLS)
fstatistic_RIVxSLS <- anova(moderation_model2a,moderation_model2a_no_RIVxSLS)
fstatistic_RIVxSLS_p <- ifelse(fstatistic_RIVxSLS[2,6] < 0.001,  "< .001", paste0("= ", substr(round(fstatistic_RIVxSLS[2,6], 3),2,5)))
paste0("f2 = ",substr(round(cohens_f_squared_RIVxSLS[1,1],2),2,4),", ΔR2 = ",substr(round(cohens_f_squared_RIVxSLS[1,5],2),2,4),", F(5, 150) = ",round(fstatistic_RIVxSLS[2,5],2),", p ",fstatistic_RIVxSLS_p)

###simple slopes
dat2a$Zsupervisor_dem <- dat2a$Zstyle_index_supervisor - abs(max(dat2a$Zstyle_index_supervisor))
dat2a$Zsupervisor_aut <- dat2a$Zstyle_index_supervisor + abs(min(dat2a$Zstyle_index_supervisor))
dat2a$ZADMxSLS_dem <- scale(dat2a$Znarq_adm*dat2a$Zsupervisor_dem)
dat2a$ZRIVxSLS_dem <- scale(dat2a$Znarq_riv*dat2a$Zsupervisor_dem)
dat2a$ZADMxSLS_aut <- scale(dat2a$Znarq_adm*dat2a$Zsupervisor_aut)
dat2a$ZRIVxSLS_aut <- scale(dat2a$Znarq_riv*dat2a$Zsupervisor_aut)
model_dem <- lm(Zattraction ~ Znarq_adm + Znarq_riv + Zsupervisor_dem + ZADMxSLS_dem + ZRIVxSLS_dem, data=dat2a)
model_aut <- lm(Zattraction ~ Znarq_adm + Znarq_riv + Zsupervisor_aut + ZADMxSLS_aut + ZRIVxSLS_aut, data=dat2a)
dat2a_mod_dem <- data.frame(summary(model_dem)$coef,confint(model_dem))
dat2a_mod_aut <- data.frame(summary(model_aut)$coef,confint(model_aut))
paste0("relation between rivalry and atraction to democratic supervisors: beta = ", round(dat2a_mod_dem[3,1], digits = 2), ", 95% CI [", round(dat2a_mod_dem[3,5], digits = 2), ", ", round(dat2a_mod_dem[3,6], digits = 2), "], t(150) = ", round(dat2a_mod_dem[3,3], digits = 2), ", p = ", round(dat2a_mod_dem[3,4], digits = 3))
paste0("relation between rivalry and atraction to autocratic supervisors: beta = ", round(dat2a_mod_aut[3,1], digits = 2), ", 95% CI [", round(dat2a_mod_aut[3,5], digits = 2), ", ", round(dat2a_mod_aut[3,6], digits = 2), "], t(150) = ", round(dat2a_mod_aut[3,3], digits = 2), ", p = ", round(dat2a_mod_aut[3,4], digits = 3))


##Study2b
##################################################################################################################################

##### ACCESS TO MPLUS #####

#select variables for Mplus analyses
dat2b_mplus <- subset(dat2b, select = c(id,team,Znarq_adm,Znarq_riv,Zstyle_index_supervisor_team,Zstyle_index_self,Zattraction))

#save dataset for Mplus analyses to folder "C:/NarcSub/" (or change path of dataset in FILE command of Mplus script)
setwd("C:/NarcSub/")
fwrite(dat2b_mplus, file = "NarcSub_Mplus-Data.csv", col.names = FALSE)
rm(dat2b_mplus)

#download the folder "NarcSub_Mplus-Scripts" via https://madata.bib.uni-mannheim.de/427/15/NarcSub_Mplus-Scripts.zip and set working directory to that folder (e.g., "C:/NarcSub/NarcSub_Mplus-Scripts/")
setwd("C:/NarcSub/NarcSub_Mplus-Scripts/")

#run moderation analysis in Mplus
runModels("MainText/narq_index_moderation.inp")


##### NO ACCESS TO MPLUS #####

#download the folder "NarcSub_Mplus-Outputs" via https://madata.bib.uni-mannheim.de/427/16/NarcSub_Mplus-Outputs.zip and set working directory to that folder (e.g., "C:/NarcSub/NarcSub_Mplus-Outputs/")
setwd("C:/NarcSub/NarcSub_Mplus-Outputs/")

##################################################################################################################################

###extract coefficients from Mplus-Output
output <- readModels("MainText/narq_index_moderation.out", what = "parameters")
output_est_ci <- sapply(output, "[", "ci.unstandardized")
est_ci <- sapply(output_est_ci, "[", c("paramHeader","param","est","low2.5","up2.5"))
output_t_p <- sapply(output, "[", "unstandardized")
t_p <- sapply(output_t_p, "[", c("paramHeader","param","est_se","pval"))
dat2b_est_ci <- data.frame(matrix(unlist(est_ci), nrow=12, ncol = 5, byrow=F), stringsAsFactors=FALSE)
dat2b_t_p <- data.frame(matrix(unlist(t_p), nrow=12, ncol = 4, byrow=F), stringsAsFactors=FALSE)
dat2b_mod <- cbind(dat2b_est_ci,dat2b_t_p[,c(3:4)])
dat2b_mod_simpleslopes <- dat2b_mod[c(11:12),]
dat2b_mod <- dat2b_mod[c(5,7,6,4,2,3),]
dat2b_mod_simpleslopes[,3:7] <- sapply(dat2b_mod_simpleslopes[,3:7], as.numeric)
paste0("relation between rivalry and atraction to democratic supervisors: zPE = ", round(dat2b_mod_simpleslopes[2,3], digits = 2), ", 95% CI [", round(dat2b_mod_simpleslopes[2,4], digits = 2), ", ", round(dat2b_mod_simpleslopes[2,5], digits = 2), "], z = ", round(dat2b_mod_simpleslopes[2,6], digits = 2), ", p = ", round(dat2b_mod_simpleslopes[2,7], digits = 3))
paste0("relation between rivalry and atraction to autocratic supervisors: zPE = ", round(dat2b_mod_simpleslopes[1,3], digits = 2), ", 95% CI [", round(dat2b_mod_simpleslopes[1,4], digits = 2), ", ", round(dat2b_mod_simpleslopes[1,5], digits = 2), "], z = ", round(dat2b_mod_simpleslopes[1,6], digits = 2), ", p = ", round(dat2b_mod_simpleslopes[1,7], digits = 3))

###estimate R-squared total
dat2b$Zattraction <- as.numeric(dat2b$Zattraction)
dat2b$Znarq_adm <- as.numeric(dat2b$Znarq_adm)
dat2b$Znarq_riv <- as.numeric(dat2b$Znarq_riv)
dat2b$Zstyle_index_supervisor_team <- as.numeric(dat2b$Zstyle_index_supervisor_team)
dat2b$ZADMxSLS <- as.numeric(dat2b$ZADMxSLS)
dat2b$ZRIVxSLS <- as.numeric(dat2b$ZRIVxSLS)
dat2b$team <- as.factor(dat2b$team)
moderation_model2b <- lmer(Zattraction ~ Znarq_adm + Znarq_riv + Zstyle_index_supervisor_team + ZADMxSLS + ZRIVxSLS + (1 + Znarq_adm + Znarq_riv || team), data=dat2b, verbose=2, REML=FALSE)
moderation_model2b_rsquared <- paste0("R2 = ",substr(round(r2mlm(moderation_model2b)$R2s[7,1],2),2,4))


##create table
dat2a_mod$Predictor <- c("(Intercept)","ADM","RIV","SLS","ADM x SLS","RIV x SLS")
dat2a_mod$blank <- NA
dat_table4 <- cbind(dat2a_mod[,c(7,1,5:6,3:4,8)],dat2b_mod[,c(3:7)])

colnames(dat_table4) <- c("Predictor","beta1","ci_lower1","ci_upper1","t","p1","blank","beta2","ci_lower2","ci_upper2","z","p2")

dat_table4[,2:12] <- sapply(dat_table4[,2:12], as.numeric)

format_beta <- function(x){
  ifelse(abs(x) < 0.005,  f_num(unlist(x), digits=3), f_num(unlist(x), digits=2))
}
dat_table4[,c(2:4,8:10)] <- dat_table4[,c(2:4,8:10)] %>% mutate_if(is.numeric, format_beta)

format_t <- function(x){
  formatC(x, format = "f", digits = 2)
}
dat_table4[,c(5,11)] <- dat_table4[,c(5,11)] %>% mutate_if(is.numeric, format_t)

format_p <- function(x){
  ifelse(abs(x) < 0.001,  "< .001", f_num(unlist(x), digits=3))
}
dat_table4[,c(6,12)] <- dat_table4[,c(6,12)] %>% mutate_if(is.numeric, format_p)

dat_table4$ci1 <- paste0("[", dat_table4$ci_lower1, ", ", dat_table4$ci_upper1, "]")
dat_table4$ci2 <- paste0("[", dat_table4$ci_lower2, ", ", dat_table4$ci_upper2, "]")

dat_table4 <- dat_table4[,c("Predictor","beta1","ci1","t","p1","blank","beta2","ci2","z","p2")]

dat_table4[1,c(2:5)] <- NA
dat_table4[7,2] <- moderation_model2a_rsquared
dat_table4[7,7] <- moderation_model2b_rsquared

col_keys <- c("Predictor","beta1","ci1","t","p1","blank","beta2","ci2","z","p2")
head1 <- c("","Study 2a","Study 2a","Study 2a","Study 2a","","Study 2b","Study 2b","Study 2b","Study 2b")
head2 <- c("Predictor","\u03B2","95% CI","t","p","","\u03B2","95% CI","z","p")
head <- data.frame(col_keys,head1,head2, stringsAsFactors = FALSE)
rm(col_keys,head1,head2)

tbl <- flextable(dat_table4)
tbl <- set_header_df(tbl, mapping=head, key="col_keys")
tbl <- hline_top(tbl, j=1:10, border=fp_border(width=2), part="header")
tbl <- merge_at(tbl, i=1, j=2:5, part="header")
tbl <- merge_at(tbl, i=1, j=7:10, part="header")
tbl <- merge_at(tbl, i=7, j=2:5, part="body")
tbl <- merge_at(tbl, i=7, j=7:10, part="body")
tbl <- hline(tbl, i=1, j=c(2:5,7:10), border=fp_border(width=1.2), part="header")
tbl <- hline(tbl, i=2, j=1:10, border=fp_border(width=1.2), part="header")
tbl <- flextable::font(tbl, fontname="Times", part="all")
tbl <- fontsize(tbl, size=12, part="all")
tbl <- fontsize(tbl, i=7, size=11, part="body")
tbl <- align(tbl, align="center", part="all")
tbl <- align(tbl, j = c("Predictor"), align="left", part="body")
tbl <- italic(tbl, i=2, j=c("t","p1","z","p2"), part="header")
tbl <- width(tbl, j =~ Predictor, width=1)
tbl <- width(tbl, j =~ beta1 + beta2, width=.45)
tbl <- width(tbl, j =~ ci1 + ci2, width=.95)
tbl <- width(tbl, j =~ t + z, width=.6)
tbl <- width(tbl, j =~ p1 + p2, width=.6)
tbl <- width(tbl, j =~ blank, width=.3)
tbl <- height_all(tbl, height=.2, part="all")
tbl

setwd(files_wd)
doc <- read_docx()
doc <- body_add_flextable(doc, value = tbl)
print(doc, target = "Table_4.docx")

rm(moderation_model2a,dat2a_mod,moderation_model2a_rsquared,moderation_model2a_rsquared_p,moderation_model2a_no_RIVxSLS,cohens_f_squared_RIVxSLS,
   fstatistic_RIVxSLS,fstatistic_RIVxSLS_p,model_aut,model_dem,dat2a_mod_aut,dat2a_mod_dem,output,output_est_ci,est_ci,output_t_p,t_p,
   dat2b_t_p,dat2b_mod,dat2b_mod_simpleslopes,moderation_model2b_rsquared,format_beta,format_p,format_t,head,tbl,doc)


#Note of Table 4
var_rand_intercept <- dat2b_est_ci[8,3]
var_rand_slope_adm <- dat2b_est_ci[10,3] #less than 0.001
var_rand_slope_riv <- dat2b_est_ci[9,3]
paste0("In Study 2b, the variance of the random intercept was ", var_rand_intercept, " and the variance of the random slopes was less than 0.001 for narcissistic admiration and ", var_rand_slope_riv, " for narcissistic rivalry.")
rm(dat2b_est_ci,var_rand_intercept,var_rand_slope_adm,var_rand_slope_riv)


#Figure 1: Attraction to Supervisor as a Function of Subordinate's Narcissistic Rivalry and Supervisor's Leadership Style in Studies 2a-2b
dat_figure1 <- dat_table4[-c(7), ]
dat_figure1[,c(2,7)] <- sapply(dat_figure1[,c(2,7)], as.numeric)

dat2a$y_plot_dem <- dat_figure1[3,2] * dat2a$Znarq_riv + dat_figure1[4,2] * max(dat2a$Zstyle_index_supervisor) + dat_figure1[6,2] * dat2a$Znarq_riv * max(dat2a$Zstyle_index_supervisor)
dat2a$y_plot_aut <- dat_figure1[3,2] * dat2a$Znarq_riv + dat_figure1[4,2] * min(dat2a$Zstyle_index_supervisor) + dat_figure1[6,2] * dat2a$Znarq_riv * min(dat2a$Zstyle_index_supervisor)
dat2a_plot_dem <- subset(dat2a, select = c(Znarq_riv,y_plot_dem))
dat2a_plot_aut <- subset(dat2a, select = c(Znarq_riv,y_plot_aut))
dat2a_plot_dem$supervisor <- "democratic"
dat2a_plot_aut$supervisor <- "autocratic"
names(dat2a_plot_dem)[2] <- "y_plot"
names(dat2a_plot_aut)[2] <- "y_plot"
dat2a_plot <- rbind(dat2a_plot_dem,dat2a_plot_aut)

dat2b$y_plot_dem <- dat_figure1[1,7] + dat_figure1[3,7] * dat2b$Znarq_riv + dat_figure1[4,7] * max(dat2b$Zstyle_index_supervisor_team) + dat_figure1[6,7] * dat2b$Znarq_riv * max(dat2b$Zstyle_index_supervisor_team)
dat2b$y_plot_aut <- dat_figure1[1,7] + dat_figure1[3,7] * dat2b$Znarq_riv + dat_figure1[4,7] * min(dat2b$Zstyle_index_supervisor_team) + dat_figure1[6,7] * dat2b$Znarq_riv * min(dat2b$Zstyle_index_supervisor_team)
dat2b_plot_dem <- subset(dat2b, select = c(Znarq_riv,y_plot_dem))
dat2b_plot_aut <- subset(dat2b, select = c(Znarq_riv,y_plot_aut))
dat2b_plot_dem$supervisor <- "democratic"
dat2b_plot_aut$supervisor <- "autocratic"
names(dat2b_plot_dem)[2] <- "y_plot"
names(dat2b_plot_aut)[2] <- "y_plot"
dat2b_plot <- rbind(dat2b_plot_dem,dat2b_plot_aut)

max(dat2a_plot$y_plot) #maximum value on y-axis in Study 2a: 1.794562 
min(dat2a_plot$y_plot) #minimum value on y-axis in Study 2a: -2.438056
max(dat2b_plot$y_plot) #maximum value on y-axis in Study 2b: 1.370678
min(dat2b_plot$y_plot) #minimum value on y-axis in Study 2b: -2.034231

legend_title <- "Supervisor's Leadership Style"
legend_labels <- c("Autocratic     ","Democratic")

Figure1A <- ggplot(data=dat2a_plot, aes(x=Znarq_riv, y=y_plot, color=factor(supervisor))) +
  geom_line(aes(linetype = factor(supervisor)), linewidth = 1) +
  scale_linetype_manual(name=legend_title, labels=legend_labels, values=c("solid","solid"), 
                        guide = guide_legend(direction = "horizontal", title.position = "top")) +
  scale_color_manual(name=legend_title, labels=legend_labels, values=c("black","grey"), 
                     guide = guide_legend(direction = "horizontal", title.position = "top")) +
  theme(legend.position = "top", legend.title = element_text(size=11), legend.title.align=0.5, 
        legend.text = element_text(size=11), legend.key = element_blank(), legend.key.width = unit(1,"cm"), legend.box.margin=margin(10,10,5,10),
        axis.text = element_text(colour = "black", size=9), axis.title.x = element_text(size=11), axis.title.y = element_text(size=11, margin=margin(0,5,0,0)),
        axis.ticks = element_blank(), panel.background = element_blank(), axis.line = element_line(colour = "black"),
        panel.grid = element_blank(), plot.tag = element_text(color="black", size=11)) +
  ylim(-2.5, 2) +
  labs(tag="A") +
  ylab("Attraction to Supervisor") +
  xlab("Subordinate's Rivalry")

Figure1B <- ggplot(data=dat2b_plot, aes(x=Znarq_riv, y=y_plot, color=factor(supervisor))) +
  geom_line(aes(linetype = factor(supervisor)), linewidth = 1) +
  scale_linetype_manual(name=legend_title, labels=legend_labels, values=c("solid","solid"), 
                        guide = guide_legend(direction = "horizontal", title.position = "top")) +
  scale_color_manual(name=legend_title, labels=legend_labels, values=c("black","grey"), 
                     guide = guide_legend(direction = "horizontal", title.position = "top")) +
  theme(legend.position = "top", legend.title = element_text(size=11), legend.title.align=0.5, 
        legend.text = element_text(size=11), legend.key = element_blank(), legend.key.width = unit(1,"cm"), legend.box.margin=margin(10,10,5,10),
        axis.text = element_text(colour = "black", size=9), axis.title.x = element_text(size=11), axis.title.y = element_text(size=11, margin=margin(0,5,0,0)),
        axis.ticks = element_blank(), panel.background = element_blank(), axis.line = element_line(colour = "black"),
        panel.grid = element_blank(), plot.tag = element_text(color="black", size=11)) +
  ylim(-2.5, 2) +
  labs(tag="B") +
  ylab("Attraction to Supervisor") +
  xlab("Subordinate's Rivalry")

Figure1 <- ggarrange(Figure1A,Figure1B, ncol=2, common.legend = TRUE, legend = "top")
ggsave("Figure_1.jpg", plot=Figure1, device="jpg", path="C:/NarcSub/", width=160, height=90, units="mm", dpi=720)

rm(dat_table4,dat_figure1,dat2a_plot,dat2a_plot_aut,dat2a_plot_dem,dat2b_plot,dat2b_plot_aut,dat2b_plot_dem,
   Figure1,Figure1A,Figure1B,legend_labels,legend_title)


#Footnote 8
max(dat2a$Zstyle_index_supervisor) #leadership style value of the most democratic supervisors in Study 2a: 1.826888
dat2a$descriptions_supervisor[dat2a$Zstyle_index_supervisor==max(dat2a$Zstyle_index_supervisor)] #leadership descriptions value of the most democratic supervisors in Study 2a: 7
dat2a$lbdq_supervisor[dat2a$Zstyle_index_supervisor==max(dat2a$Zstyle_index_supervisor)] #leadership behaviors value of the most democratic supervisors in Study 2a: 4.75
min(dat2a$Zstyle_index_supervisor) #leadership style value of the most autocratic supervisors in Study 2a: -2.404079
dat2a$descriptions_supervisor[dat2a$Zstyle_index_supervisor==min(dat2a$Zstyle_index_supervisor)] #leadership descriptions value of the most autocratic supervisors in Study 2a: 1
dat2a$lbdq_supervisor[dat2a$Zstyle_index_supervisor==min(dat2a$Zstyle_index_supervisor)] #leadership behaviors value of the most autocratic supervisors in Study 2a: 1.20
boxplot(dat2a$Zstyle_index_supervisor, ylab = "Zstyle_index_supervisor") #no leadership style values were outliers in Study 2a


#Footnote 12
max(dat2b$Zstyle_index_supervisor_team) #leadership style value of the most democratic supervisors in Study 2b: 1.541836
dat2b$descriptions_supervisor_team[dat2b$Zstyle_index_supervisor_team==max(dat2b$Zstyle_index_supervisor_team)] #leadership descriptions value of the most democratic supervisors in Study 2b: 5.75
dat2b$lbdq_supervisor_team[dat2b$Zstyle_index_supervisor_team==max(dat2b$Zstyle_index_supervisor_team)] #leadership behaviors value of the most democratic supervisors in Study 2b: 5.5
min(dat2b$Zstyle_index_supervisor_team) #leadership style value of the most autocratic supervisors in Study 2b: -2.2365
dat2b$descriptions_supervisor_team[dat2b$Zstyle_index_supervisor_team==min(dat2b$Zstyle_index_supervisor_team)] #leadership descriptions value of the most autocratic supervisors in Study 2b: 2.5
dat2b$lbdq_supervisor_team[dat2b$Zstyle_index_supervisor_team==min(dat2b$Zstyle_index_supervisor_team)] #leadership behaviors value of the most autocratic supervisors in Study 2b: 1.63
boxplot(dat2b$Zstyle_index_supervisor_team, ylab = "Zstyle_index_supervisor_team") #no leadership style values were extreme outliers in Study 2b


#Table 5: Results of the Moderated Mediation Model in Studies 2a-2b

##Study 2a
###moderated mediation model
set.seed(123) #setting the seed for reproducibility
mediation_model2a <- '
Zattraction ~ Znarq_adm + c*Znarq_riv + Zstyle_index_supervisor + b*Zstyle_index_self + ZADMxSLS + h*ZRIVxSLS + g*ZSHLSxSLS
Zstyle_index_self ~ a1*Znarq_adm + a2*Znarq_riv
## index of moderated mediaion
index_adm := a1 * g
index_riv := a2 * g
## conditional effects of RIV
indirect_aut := a2*(b-2.404079*g)
indirect_dem := a2*(b+1.826888*g)
direct_aut := c-2.404079*h
direct_dem := c+1.826888*h
## conditional effects of SHLS
SHLS_aut := b-2.404079*g
SHLS_dem := b+1.826888*g'
results <- sem(mediation_model2a, data = dat2a, meanstructure = TRUE,
               se = "bootstrap", bootstrap = 10000)
dat2a_modmed <- data.frame(parameterEstimates(results, ci = TRUE, level = 0.95,
                                              boot.ci.type = "bca.simple")) #bias-corrected percentile method
dat2a_modmed_index_simpleslopes <- dat2a_modmed[c(41:48),]
dat2a_modmed <- dat2a_modmed[c(8:9,33,1:2,4,3,5:7),]

mediation_model2a_path_a <- lm(Zstyle_index_self ~ Znarq_adm + Znarq_riv, data=dat2a)
mediation_model2a_path_a_rsquared_p <- pf(summary(mediation_model2a_path_a)$fstatistic[1],summary(mediation_model2a_path_a)$fstatistic[2],summary(mediation_model2a_path_a)$fstatistic[3],lower.tail=FALSE)
mediation_model2a_path_a_rsquared_p <- ifelse(mediation_model2a_path_a_rsquared_p < 0.001,  "< .001", paste0("= ", substr(round(mediation_model2a_path_a_rsquared_p, 3),2,5)))
mediation_model2a_path_a_rsquared <- paste0("R2 = ",substr(round(summary(mediation_model2a_path_a)$r.squared,2),2,4),", F(",summary(mediation_model2a_path_a)$fstatistic[2],", ",summary(mediation_model2a_path_a)$fstatistic[3],") = ",round(summary(mediation_model2a_path_a)$fstatistic[1],2),", p ",mediation_model2a_path_a_rsquared_p)

mediation_model2a_path_b <- lm(Zattraction ~ Znarq_adm + Znarq_riv + Zstyle_index_supervisor + Zstyle_index_self + ZADMxSLS + ZRIVxSLS + ZSHLSxSLS, data=dat2a)
mediation_model2a_path_b_rsquared_p <- pf(summary(mediation_model2a_path_b)$fstatistic[1],summary(mediation_model2a_path_b)$fstatistic[2],summary(mediation_model2a_path_b)$fstatistic[3],lower.tail=FALSE)
mediation_model2a_path_b_rsquared_p <- ifelse(mediation_model2a_path_b_rsquared_p < 0.001,  "< .001", paste0("= ", substr(round(mediation_model2a_path_b_rsquared_p, 3),2,5)))
mediation_model2a_path_b_rsquared <- paste0("R2 = ",substr(round(summary(mediation_model2a_path_b)$r.squared,2),2,4),", F(",summary(mediation_model2a_path_b)$fstatistic[2],", ",summary(mediation_model2a_path_b)$fstatistic[3],") = ",round(summary(mediation_model2a_path_b)$fstatistic[1],2),", p ",mediation_model2a_path_b_rsquared_p)

###effect size of relation between narcissistic rivalry and subordinate’s hypothetical leadership style 
mediation_model2a_path_a_no_riv <- lm(Zstyle_index_self ~ Znarq_adm, data=dat2a)
cohens_f_squared_path_a_riv <- cohens_f_squared(mediation_model2a_path_a, model2 = mediation_model2a_path_a_no_riv)
fstatistic_path_a_riv <- anova(mediation_model2a_path_a,mediation_model2a_path_a_no_riv)
fstatistic_path_a_riv_p <- ifelse(fstatistic_path_a_riv[2,6] < 0.001,  "< .001", paste0("= ", substr(round(fstatistic_path_a_riv[2,6], 3),2,5)))
paste0("f2 = ",substr(round(cohens_f_squared_path_a_riv[1,1],2),2,4),", ΔR2 = ",substr(round(cohens_f_squared_path_a_riv[1,5],2),2,4),", F(1, 153) = ",round(fstatistic_path_a_riv[2,5],2),", p ",fstatistic_path_a_riv_p)

###effect size of interaction between subordinate’s hypothetical leadership style  and supervisor’s leadership style 
mediation_model2a_path_b_no_SHLSxSLS <- lm(Zattraction ~ Znarq_adm + Znarq_riv + Zstyle_index_supervisor + Zstyle_index_self + ZADMxSLS + ZRIVxSLS, data=dat2a)
cohens_f_squared_SHLSxSLS <- cohens_f_squared(mediation_model2a_path_b, model2 = mediation_model2a_path_b_no_SHLSxSLS)
fstatistic_SHLSxSLS <- anova(mediation_model2a_path_b,mediation_model2a_path_b_no_SHLSxSLS)
fstatistic_SHLSxSLS_p <- ifelse(fstatistic_SHLSxSLS[2,6] < 0.001,  "< .001", paste0("= ", substr(round(fstatistic_SHLSxSLS[2,6], 3),2,5)))
paste0("f2 = ",substr(round(cohens_f_squared_SHLSxSLS[1,1],2),2,4),", ΔR2 = ",substr(round(cohens_f_squared_SHLSxSLS[1,5],2),2,4),", F(1, 148) = ",round(fstatistic_SHLSxSLS[2,5],2),", p ",fstatistic_SHLSxSLS_p)


##Study2b
##################################################################################################################################

##### ACCESS TO MPLUS #####

#select variables for analyses in Mplus
dat2b_mplus <- subset(dat2b, select = c(id,team,Znarq_adm,Znarq_riv,Zstyle_index_supervisor_team,Zstyle_index_self,Zattraction))

#save dataset for analyses in Mplus to folder "C:/NarcSub/" (or change path of dataset in FILE command of Mplus script)
setwd("C:/NarcSub/")
fwrite(dat2b_mplus, file = "NarcSub_Mplus-Data.csv", col.names = FALSE)

#download the folder "NarcSub_Mplus-Scripts" via https://madata.bib.uni-mannheim.de/427/15/NarcSub_Mplus-Scripts.zip and set working directory to that folder (e.g., "C:/NarcSub/NarcSub_Mplus-Scripts/")
setwd("C:/NarcSub/NarcSub_Mplus-Scripts/")

#run mediation analysis in Mplus
runModels("MainText/narq_index_mediation.inp")


##### NO ACCESS TO MPLUS #####

#download the folder "NarcSub_Mplus-Outputs" via https://madata.bib.uni-mannheim.de/427/16/NarcSub_Mplus-Outputs.zip and set working directory to that folder (e.g., "C:/NarcSub/NarcSub_Mplus-Outputs/")
setwd("C:/NarcSub/NarcSub_Mplus-Outputs/")

##################################################################################################################################

###extract coefficients from Mplus-Output
output <- readModels("MainText/narq_index_mediation.out", what = "parameters")
output_est_ci <- sapply(output, "[", "ci.unstandardized")
est_ci <- sapply(output_est_ci, "[", c("paramHeader","param","est","low2.5","up2.5"))
output_t_p <- sapply(output, "[", "unstandardized")
t_p <- sapply(output_t_p, "[", c("paramHeader","param","est_se","pval"))
dat2b_est_ci <- data.frame(matrix(unlist(est_ci), nrow=25, ncol = 5, byrow=F), stringsAsFactors=FALSE)
dat2b_t_p <- data.frame(matrix(unlist(t_p), nrow=25, ncol = 4, byrow=F), stringsAsFactors=FALSE)
dat2b_modmed <- cbind(dat2b_est_ci,dat2b_t_p[,c(3:4)])
dat2b_modmed_index_simpleslopes <- dat2b_modmed[c(18:25),]
dat2b_modmed <- dat2b_modmed[c(2,1,10,12,11,13,9,6,8,7),]

###estimate R-squared total
dat2b$Zattraction <- as.numeric(dat2b$Zattraction)
dat2b$Znarq_adm <- as.numeric(dat2b$Znarq_adm)
dat2b$Znarq_riv <- as.numeric(dat2b$Znarq_riv)
dat2b$Zstyle_index_self <- as.numeric(dat2b$Zstyle_index_self)
dat2b$Zstyle_index_supervisor_team <- as.numeric(dat2b$Zstyle_index_supervisor_team)
dat2b$ZADMxSLS <- as.numeric(dat2b$ZADMxSLS)
dat2b$ZRIVxSLS <- as.numeric(dat2b$ZRIVxSLS)
dat2b$ZSHLSxSLS <- as.numeric(dat2b$ZSHLSxSLS)
dat2b$team <- as.factor(dat2b$team)
mediation_model2b_path_a <- lmer(Zstyle_index_self ~ Znarq_adm + Znarq_riv + (1 + Znarq_adm + Znarq_riv || team), data=dat2b, verbose=2, REML=FALSE)
mediation_model2b_path_a_rsquared <- paste0("R2 = ",substr(round(r2mlm(mediation_model2b_path_a)$R2s[7,1],2),2,4))
mediation_model2b_path_b <- lmer(Zattraction ~ Znarq_adm + Znarq_riv + Zstyle_index_self + Zstyle_index_supervisor_team + ZADMxSLS + ZRIVxSLS + ZSHLSxSLS + (1 + Znarq_adm + Znarq_riv + Zstyle_index_self || team), data=dat2b, verbose=2, REML=FALSE)
mediation_model2b_path_b_rsquared <- paste0("R2 = ",substr(round(r2mlm(mediation_model2b_path_b)$R2s[7,1],2),2,4))


##create table
dat2a_modmed$Predictor <- c("ADM","RIV","(Intercept)","ADM","RIV","SHLS","SLS","ADM x SLS","RIV x SLS","SHLS x SLS")
dat2a_modmed$blank <- NA
dat_table5 <- cbind(dat2a_modmed[,c(11,5,9:10,7:8,12)],dat2b_modmed[,c(3:7)])

colnames(dat_table5) <- c("Predictor","beta1","ci_lower1","ci_upper1","z1","p1","blank","beta2","ci_lower2","ci_upper2","z2","p2")

dat_table5[,2:12] <- sapply(dat_table5[,2:12], as.numeric)

format_beta <- function(x){
  ifelse(abs(x) < 0.005,  f_num(unlist(x), digits=3), f_num(unlist(x), digits=2))
}
dat_table5[,c(2:4,8:10)] <- dat_table5[,c(2:4,8:10)] %>% mutate_if(is.numeric, format_beta)

format_t <- function(x){
  formatC(x, format = "f", digits = 2)
}
dat_table5[,c(5,11)] <- dat_table5[,c(5,11)] %>% mutate_if(is.numeric, format_t)

format_p <- function(x){
  ifelse(abs(x) < 0.001,  "< .001", f_num(unlist(x), digits=3))
}
dat_table5[,c(6,12)] <- dat_table5[,c(6,12)] %>% mutate_if(is.numeric, format_p)

dat_table5$ci1 <- paste0("[", dat_table5$ci_lower1, ", ", dat_table5$ci_upper1, "]")
dat_table5$ci2 <- paste0("[", dat_table5$ci_lower2, ", ", dat_table5$ci_upper2, "]")

dat_table5 <- dat_table5[,c("Predictor","beta1","ci1","z1","p1","blank","beta2","ci2","z2","p2")]

dat_table5[3,c(2:5)] <- NA
dat_table5[11,2] <- "Subordinate's Hypothetical Leadership Style"
dat_table5[12,2] <- mediation_model2a_path_a_rsquared
dat_table5[12,7] <- mediation_model2b_path_a_rsquared
dat_table5[13,2] <- "Attraction to Supervisor"
dat_table5[14,2] <- mediation_model2a_path_b_rsquared
dat_table5[14,7] <- mediation_model2b_path_b_rsquared
dat_table5 <- dat_table5[c(11,1:2,12,13,3:10,14),]

col_keys <- c("Predictor","beta1","ci1","z1","p1","blank","beta2","ci2","z2","p2")
head1 <- c("","Study 2a","Study 2a","Study 2a","Study 2a","","Study 2b","Study 2b","Study 2b","Study 2b")
head2 <- c("Predictor","\u03B2","95% CI","z","p","","\u03B2","95% CI","z","p")
head <- data.frame(col_keys,head1,head2, stringsAsFactors = FALSE)
rm(col_keys,head1,head2)

tbl <- flextable(dat_table5)
tbl <- set_header_df(tbl, mapping=head, key="col_keys")
tbl <- hline_top(tbl, j=1:10, border=fp_border(width=2), part="header")
tbl <- merge_at(tbl, i=1, j=2:5, part="header")
tbl <- merge_at(tbl, i=1, j=7:10, part="header")
tbl <- merge_at(tbl, i=1, j=2:10, part="body")
tbl <- merge_at(tbl, i=4, j=2:5, part="body")
tbl <- merge_at(tbl, i=4, j=7:10, part="body")
tbl <- merge_at(tbl, i=5, j=2:10, part="body")
tbl <- merge_at(tbl, i=14, j=2:5, part="body")
tbl <- merge_at(tbl, i=14, j=7:10, part="body")
tbl <- hline(tbl, i=1, j=c(2:5,7:10), border=fp_border(width=1.2), part="header")
tbl <- hline(tbl, i=2, j=1:10, border=fp_border(width=1.2), part="header")
tbl <- flextable::font(tbl, fontname="Times", part="all")
tbl <- fontsize(tbl, size=12, part="all")
tbl <- fontsize(tbl, i=c(4,14), size=11, part="body")
tbl <- align(tbl, align="center", part="all")
tbl <- align(tbl, j = c("Predictor"), align="left", part="body")
tbl <- italic(tbl, i=2, j=c("z1","p1","z2","p2"), part="header")
tbl <- width(tbl, j =~ Predictor, width=1.1)
tbl <- width(tbl, j =~ beta1 + beta2, width=.5)
tbl <- width(tbl, j =~ ci1 + ci2, width=.95)
tbl <- width(tbl, j =~ z1 + z2, width=.6)
tbl <- width(tbl, j =~ p1 + p2, width=.6)
tbl <- width(tbl, j =~ blank, width=.3)
tbl <- height_all(tbl, height=.2, part="all")
tbl

setwd(files_wd)
doc <- read_docx()
doc <- body_add_flextable(doc, value = tbl)
print(doc, target = "Table_5.docx")

rm(mediation_model2a,dat2a_modmed,results,mediation_model2a_path_a,mediation_model2a_path_a_no_riv,mediation_model2a_path_b_no_SHLSxSLS,mediation_model2a_path_a_rsquared_p,mediation_model2a_path_a_rsquared,mediation_model2a_path_b,mediation_model2a_path_b_rsquared_p,mediation_model2a_path_b_rsquared,
   dat2b_est_ci,est_ci,dat2b_t_p,dat2b_modmed,mediation_model2b_path_a_rsquared,mediation_model2b_path_b_rsquared,cohens_f_squared_path_a_riv,fstatistic_path_a_riv,fstatistic_path_a_riv_p ,cohens_f_squared_SHLSxSLS,fstatistic_SHLSxSLS,fstatistic_SHLSxSLS_p,
   output,output_est_ci,output_t_p,t_p,head,tbl,format_beta,format_p,format_t,doc)


#Figure 2: Moderated Mediation Model in Studies 2a-2b
##path coefficients in Figure 2a
paste0("effect of rivalry on SHLS: beta = ", dat_table5$beta1[3], ", 95% CI [", dat_table5$ci1[3], "]")
paste0("direct effect of rivalry on attraction to democratic supervisors: beta = ", round(dat2a_modmed_index_simpleslopes[6,5], digits = 2), ", 95% CI [", round(dat2a_modmed_index_simpleslopes[6,9], digits = 2), ", ", round(dat2a_modmed_index_simpleslopes[6,10], digits = 2), "]")
paste0("direct effect of rivalry on attraction to autocratic supervisors: beta = ", round(dat2a_modmed_index_simpleslopes[5,5], digits = 2), ", 95% CI [", round(dat2a_modmed_index_simpleslopes[5,9], digits = 2), ", ", round(dat2a_modmed_index_simpleslopes[5,10], digits = 2), "]")
paste0("relation between SHLS and attraction to democratic supervisors: beta = ", round(dat2a_modmed_index_simpleslopes[8,5], digits = 2), ", 95% CI [", round(dat2a_modmed_index_simpleslopes[8,9], digits = 2), ", ", round(dat2a_modmed_index_simpleslopes[8,10], digits = 2), "]")
paste0("relation between SHLS and attraction to autocratic supervisors: beta = ", round(dat2a_modmed_index_simpleslopes[7,5], digits = 2), ", 95% CI [", round(dat2a_modmed_index_simpleslopes[7,9], digits = 2), ", ", round(dat2a_modmed_index_simpleslopes[7,10], digits = 2), "]")
##path coefficients in Figure 2b
dat2b_modmed_index_simpleslopes[,3:7] <- sapply(dat2b_modmed_index_simpleslopes[,3:7], as.numeric)
paste0("effect of rivalry on SHLS: zPE = ", dat_table5$beta2[3], ", 95% CI [", dat_table5$ci2[3], "]")
paste0("direct effect of rivalry on attraction to democratic supervisors: zPE = ", round(dat2b_modmed_index_simpleslopes[6,3], digits = 2), ", 95% CI [", round(dat2b_modmed_index_simpleslopes[6,4], digits = 2), ", ", round(dat2b_modmed_index_simpleslopes[6,5], digits = 2), "]")
paste0("direct effect of rivalry on attraction to autocratic supervisors: zPE = ", round(dat2b_modmed_index_simpleslopes[5,3], digits = 2), ", 95% CI [", round(dat2b_modmed_index_simpleslopes[5,4], digits = 2), ", ", round(dat2b_modmed_index_simpleslopes[5,5], digits = 2), "]")
paste0("relation between SHLS and attraction to democratic supervisors: zPE = ", round(dat2b_modmed_index_simpleslopes[8,3], digits = 2), ", 95% CI [", round(dat2b_modmed_index_simpleslopes[8,4], digits = 2), ", ", round(dat2b_modmed_index_simpleslopes[8,5], digits = 2), "]")
paste0("relation between SHLS and attraction to autocratic supervisors: zPE = ", round(dat2b_modmed_index_simpleslopes[7,3], digits = 2), ", 95% CI [", round(dat2b_modmed_index_simpleslopes[7,4], digits = 2), ", ", round(dat2b_modmed_index_simpleslopes[7,5], digits = 2), "]")


#indirect effects
##Study 2a
paste0("Study 2a: index of moderated mediation for rivalry: estimate = ", round(dat2a_modmed_index_simpleslopes[2,5], digits = 2), ", 95% CI [", round(dat2a_modmed_index_simpleslopes[2,9], digits = 2), ", ", round(dat2a_modmed_index_simpleslopes[2,10], digits = 2), "], z = ", round(dat2a_modmed_index_simpleslopes[2,7], digits = 2), ", p = ", round(dat2a_modmed_index_simpleslopes[2,8], digits = 3))
paste0("Study 2a: indirect effect of rivalry for democratic supervisors: estimate = ", round(dat2a_modmed_index_simpleslopes[4,5], digits = 2), ", 95% CI [", round(dat2a_modmed_index_simpleslopes[4,9], digits = 2), ", ", round(dat2a_modmed_index_simpleslopes[4,10], digits = 2), "], z = ", round(dat2a_modmed_index_simpleslopes[4,7], digits = 2), ", p = ", round(dat2a_modmed_index_simpleslopes[4,8], digits = 3))
paste0("Study 2a: indirect effect of rivalry for autocratic supervisors: estimate = ", round(dat2a_modmed_index_simpleslopes[3,5], digits = 2), ", 95% CI [", round(dat2a_modmed_index_simpleslopes[3,9], digits = 2), ", ", round(dat2a_modmed_index_simpleslopes[3,10], digits = 2), "], z = ", round(dat2a_modmed_index_simpleslopes[3,7], digits = 2), ", p = ", round(dat2a_modmed_index_simpleslopes[3,8], digits = 3))
paste0("Study 2a: index of moderated mediation for admiration: estimate = ", round(dat2a_modmed_index_simpleslopes[1,5], digits = 2), ", 95% CI [", round(dat2a_modmed_index_simpleslopes[1,9], digits = 2), ", ", round(dat2a_modmed_index_simpleslopes[1,10], digits = 2), "], z = ", round(dat2a_modmed_index_simpleslopes[1,7], digits = 2), ", p = ", round(dat2a_modmed_index_simpleslopes[1,8], digits = 3))
##Study 2b
dat2b_modmed_index_simpleslopes[,3:7] <- sapply(dat2b_modmed_index_simpleslopes[,3:7], as.numeric)
paste0("Study 2b: index of moderated mediation for rivalry: estimate = ", round(dat2b_modmed_index_simpleslopes[1,3], digits = 2), ", 95% CI [", round(dat2b_modmed_index_simpleslopes[1,4], digits = 2), ", ", round(dat2b_modmed_index_simpleslopes[1,5], digits = 2), "], z = ", round(dat2b_modmed_index_simpleslopes[1,6], digits = 2), ", p = ", round(dat2b_modmed_index_simpleslopes[1,7], digits = 3))
paste0("Study 2b: indirect effect of rivalry for democratic supervisors: estimate = ", round(dat2b_modmed_index_simpleslopes[4,3], digits = 2), ", 95% CI [", round(dat2b_modmed_index_simpleslopes[4,4], digits = 2), ", ", round(dat2b_modmed_index_simpleslopes[4,5], digits = 2), "], z = ", round(dat2b_modmed_index_simpleslopes[4,6], digits = 2), ", p = ", round(dat2b_modmed_index_simpleslopes[4,7], digits = 3))
paste0("Study 2b: indirect effect of rivalry for autocratic supervisors: estimate = ", round(dat2b_modmed_index_simpleslopes[3,3], digits = 2), ", 95% CI [", round(dat2b_modmed_index_simpleslopes[3,4], digits = 2), ", ", round(dat2b_modmed_index_simpleslopes[3,5], digits = 2), "], z = ", round(dat2b_modmed_index_simpleslopes[3,6], digits = 2), ", p = ", round(dat2b_modmed_index_simpleslopes[3,7], digits = 3))
paste0("Study 2b: index of moderated mediation for admiration: estimate = ", round(dat2b_modmed_index_simpleslopes[2,3], digits = 2), ", 95% CI [", round(dat2b_modmed_index_simpleslopes[2,4], digits = 2), ", ", round(dat2b_modmed_index_simpleslopes[2,5], digits = 2), "], z = ", round(dat2b_modmed_index_simpleslopes[2,6], digits = 2), ", p = ", round(dat2b_modmed_index_simpleslopes[2,7], digits = 3))
rm(dat2a_modmed_index_simpleslopes,dat2b_modmed_index_simpleslopes)


#observed power in Study 2b
set.seed(456) #setting the seed for reproducibility
sim.ef_ZRIVxSLS <- powerSim(moderation_model2b, fixed("ZRIVxSLS"))
print(sim.ef_ZRIVxSLS)
set.seed(456) #setting the seed for reproducibility
sim.ef_Znarq_riv <- powerSim(mediation_model2b_path_a, fixed("Znarq_riv"))
print(sim.ef_Znarq_riv)
set.seed(456) #setting the seed for reproducibility
sim.ef_ZSHLSxSLS <- powerSim(mediation_model2b_path_b, fixed("ZSHLSxSLS"))
print(sim.ef_ZSHLSxSLS)


#clear environment
rm(list = ls())