library (foreign)
library (ggplot2)
library (nlme)
library (lattice)
library (MASS)
library (Matrix)
library (lme4)
library (multilevel)
library (psy)
library(ordinal)
library(pspearman)
library(readxl)
library(car)


#REGION CRONBACH ALPHA ALLIANCE

data_alliLSHM <- read_excel("PATH/Alliance dataLSHM.xlsx")
View(data_alliLSHM)
cronbach(data_alliLSHM)

data_alliHSHM <- read_excel("PATH/Alliance dataHSHM.xlsx")
View(data_alliHSHM)
cronbach(data_alliHSHM)

data_alliLSMM <- read_excel("PATH/Alliance dataLSMM.xlsx")
View(data_alliLSMM)
cronbach(data_alliLSMM)

data_alliMSMM <- read_excel("PATH/Alliance dataMSMM.xlsx")
View(data_alliMSMM)
cronbach(data_alliMSMM)

data_alliHSMM <- read_excel("PATH/Alliance dataHSMM.xlsx")
View(data_alliHSMM)
cronbach(data_alliHSMM)

data_alliLSLM <- read_excel("PATH/Alliance dataLSLM.xlsx")
View(data_alliLSLM)
cronbach(data_alliLSLM)

data_alliMSLM <- read_excel("PATH/Alliance dataMSLM.xlsx")
View(data_alliMSLM)
cronbach(data_alliMSLM)

data_alliMSHM <- read_excel("PATH/Alliance dataMSHM.xlsx")
View(data_alliMSHM)
cronbach(data_alliMSHM)

data_alliHSLM <- read_excel("PATH/Alliance dataHSLM.xlsx")
View(data_alliHSLM)
cronbach(data_alliHSLM)



#Region Data analysis & Descriptive outcomes

data <- read_excel("PATH/Datatable.xlsx")

#Specify levels of MatchSys as 0 or 1 
data$MatchSys = factor(data$MatchSys,levels=c(0,1))
#Specify levels for Route as Motiv, Accept or Refer
data$Route = factor(data$Route,levels=c("motivate","accept","refer"))

#Specify levels for SubSev as Low, Med, High
data$SubSev = factor(data$SubSev,levels=c("Low","Med","High"))
data$SubSev <- ordered(data$SubSev)
#Specify levels of SubMoti as Low, Med, High
data$SubMoti = factor(data$SubMoti,levels=c("Low","Med","High"))
data$SubMoti <- ordered(data$SubMoti)

#Specify levels for ScenSev as Low, Med, High
data$ScenSev = factor(data$ScenSev,levels=c("Low","Med","High"))
data$ScenSev <- ordered(data$ScenSev)
#Specify levels of ScenMoti as Low, Med, High
data$ScenMoti = factor(data$ScenMoti,levels=c("Low","Med","High"))
data$ScenMoti <- ordered(data$ScenMoti)


data$ScenSit = factor(data$ScenSit,levels=c("LsLm","MsLm","HsLm","LsMm","MsMm","HsMm","LsHm","MsHm","HsHm"))
data$SubSit = factor(data$SubSit,levels=c("LsLm","MsLm","HsLm","LsMm","MsMm","HsMm","LsHm","MsHm","HsHm"))

data$WhyTxt = factor(data$WhyTxt,levels=c("Mon","Stig","Time"))

data$ppCode = factor(data$ppCode)

#Display data
View(data)


#EXPLORATION Graphs
### data exploration
stem(data$FBH)
hist(data$FBH)

ggplot(data, aes(x=Help, fill=MatchSys)) + geom_density(alpha = 0.3)
ggplot(data, aes(x=Help, fill=Route)) + geom_density(alpha = 0.3)

ggplot(data, aes(x=Return, fill=MatchSys)) + geom_density(alpha = 0.3)
ggplot(data, aes(x=Return, fill=Route)) + geom_density(alpha = 0.3)

ggplot(data, aes(x=FBH, fill=MatchSys)) + geom_density(alpha = 0.3)
ggplot(data, aes(x=FBH, fill=Route)) + geom_density(alpha = 0.3)
ggplot(data, aes(x=FBH, fill=Gen)) + geom_density(alpha = 0.3)

#MEAN, 0 Modell revealing deviation from 0 & moti + severity added to check if these have an overall effect.
mean(data$Help)
sd(data$Help)
Model_null_h <- lme(fixed = Help~1, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_null_h)
anova(Model_null_h)
Model_scen_h <- lme(fixed = Help~ScenSit, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_h,Model_scen_h)
anova(Model_scen_h)


mean(data$Return)
sd(data$Return)
Model_null_r <- lme(fixed = Return~1, random = ~1 | ppCode, data = data, method = "ML")

anova(Model_null_r)
Model_scen_r <- lme(fixed = Return~ScenSit, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_r,Model_scen_r)
anova(Model_scen_r)

mean(data$FBH)
sd(data$FBH)
Model_null_a <- lme(fixed = FBH~1, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_a)
Model_scen_a <- lme(fixed = FBH~ScenSit, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_scen_a)
anova(Model_null_a,Model_scen_a)


#REGION COMPARE ScenJECTIVE - SCENE SITUATION
fmodelSnull <- clmm2(ScenSev~1,random = ppCode,data=data,Hess=TRUE,nAGQ = 10)
summary(fmodelSnull)

fmodelS <- clmm2(ScenSev ~ SubSev, random = ppCode, data=data, Hess=TRUE,nAGQ = 10)
summary(fmodelS)
anova(fmodelSnull,fmodelS)

( McF.pR2 <- 1 - fmodelS$logLik/fmodelSnull$logLik )


fmodelMnull <- clmm2(ScenMoti~1,random = ppCode,data=data,Hess=TRUE,nAGQ = 10)
summary(fmodelMnull)
fmodelM <- clmm2(SubMoti ~ ScenMoti , random = ppCode,data=data,Hess=TRUE,nAGQ = 10)
summary(fmodelM)
anova(fmodelMnull,fmodelM)

( McF.pR2 <- 1 - fmodelM$logLik/fmodelMnull$logLik )






#Covariate checks
Model_age_h <- lme(fixed = Help~Age, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_h,Model_age_h)
Model_age_r <- lme(fixed = Return~Age, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_r,Model_age_r)
Model_age_a <- lme(fixed = FBH~Age, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_a,Model_age_a)

Model_gen_h <- lme(fixed = Help~Gen, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_h,Model_gen_h)
Model_gen_r <- lme(fixed = Return~Gen, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_r,Model_gen_r)
Model_gen_a <- lme(fixed = FBH~Gen, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_gen_a)
anova(Model_null_a,Model_gen_a)
anova(Model_gen_a)

Model_isi_h <- lme(fixed = Help~ISITot, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_h,Model_isi_h)
Model_isi_r <- lme(fixed = Return~ISITot, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_r,Model_isi_r)
Model_isi_a <- lme(fixed = FBH~ISITot, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_a,Model_isi_a)

Model_why_h <- lme(fixed = Help~WhyTxt, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_h,Model_why_h)
Model_why_r <- lme(fixed = Return~WhyTxt, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_r,Model_why_r)
Model_why_a <- lme(fixed = FBH~WhyTxt, random = ~1 | ppCode, data = data, method = "ML")
anova(Model_null_a,Model_why_a)

#Best Route generally
Model_route_h <- lme(fixed = Help~Route, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_route_h)
anova(Model_null_h,Model_route_h)
anova(Model_route_h)
Model_route_r <- lme(fixed = Return~Route, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_route_r)
anova(Model_null_r,Model_route_r)
anova(Model_route_r)
Model_route_a <- lme(fixed = FBH~Route, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_route_a)
anova(Model_null_a,Model_route_a)
anova(Model_route_a)


ggplot(data, aes(x=Route, y=FBH)) + geom_density(alpha = 0.3)
qplot(x=Route, y=FBH,data=data, geom="bar")


#HYPOTHESIS1

data2  <-subset(data, (Route == "accept"| Route == "refer"))
View(data2)


ModelH1_0 <- lme(fixed = Help~1, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_0)

ModelH1_1 <- lme(fixed = Help~ScenSit, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_1)
anova(ModelH1_0,ModelH1_1)

ModelH1_2 <- lme(fixed = Help~Route, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_2)
anova(ModelH1_0,ModelH1_2)

ModelH1_3 <- lme(fixed = Help~Route+ScenSit, random = ~1 | ppCode, data = data, method = "ML")
ModelH1_4 <- lme(fixed = Help~Route+ScenSit+Route:ScenSit, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_4)
anova(ModelH1_3,ModelH1_4)

ModelH1_3 <- lme(fixed = Help~Route+ScenSit, random = ~1 | ppCode, data = data2, method = "ML")
ModelH1_4 <- lme(fixed = Help~Route+ScenSit+Route:ScenSit, random = ~1 | ppCode, data = data2, method = "ML")
summary(ModelH1_4)
anova(ModelH1_3,ModelH1_4)


ModelH1_0 <- lme(fixed = Return~1, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_0)

ModelH1_1 <- lme(fixed = Return~ScenSit, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_1)
anova(ModelH1_0,ModelH1_1)

ModelH1_2 <- lme(fixed = Return~Route, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_2)
anova(ModelH1_0,ModelH1_2)

ModelH1_3 <- lme(fixed = Return~Route+ScenSit, random = ~1 | ppCode, data = data, method = "ML")
ModelH1_4 <- lme(fixed = Return~Route+ScenSit+Route:ScenSit, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_4)
anova(ModelH1_3,ModelH1_4)

ModelH1_3 <- lme(fixed = Return~Route+ScenSit, random = ~1 | ppCode, data = data2, method = "ML")
ModelH1_4 <- lme(fixed = Return~Route+ScenSit+Route:ScenSit, random = ~1 | ppCode, data = data2, method = "ML")
summary(ModelH1_4)
anova(ModelH1_3,ModelH1_4)



ModelH1_0 <- lme(fixed = FBH~1, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_0)

ModelH1_1 <- lme(fixed = FBH~ScenSit, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_1)
anova(ModelH1_0,ModelH1_1)

ModelH1_2 <- lme(fixed = FBH~Route, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_2)
anova(ModelH1_0,ModelH1_2)

ModelH1_3 <- lme(fixed = FBH~Route+ScenSit, random = ~1 | ppCode, data = data, method = "ML")
ModelH1_4 <- lme(fixed = FBH~Route+ScenSit+Route:ScenSit, random = ~1 | ppCode, data = data, method = "ML")
summary(ModelH1_4)
anova(ModelH1_3,ModelH1_4)

ModelH1_3 <- lme(fixed = FBH~Route+ScenSit, random = ~1 | ppCode, data = data2, method = "ML")
ModelH1_4 <- lme(fixed = FBH~Route+ScenSit+Route:ScenSit, random = ~1 | ppCode, data = data2, method = "ML")
summary(ModelH1_4)
anova(ModelH1_3,ModelH1_4)


#HYPOTHESIS2
dataH2 <- subset(data, (Systemroute == "motivate"))
#View(dataH2)

ModelH2_0 <- lme(fixed = Help~1, random = ~1 | ppCode, data = dataH2, method = "ML")
summary(ModelH2_0)

ModelH2_1 <- lme(fixed = Help~Route, random = ~1 | ppCode, data = dataH2, method = "ML")
summary(ModelH2_1)
anova(ModelH2_0,ModelH2_1)

dataH2B <- subset(dataH2, (Route == "motivate" | Route == "refer"))

ModelH2_0B <- lme(fixed = Help~1, random = ~1 | ppCode, data = dataH2B, method = "ML")
summary(ModelH2_0B)

ModelH2_1B <- lme(fixed = Help~Route, random = ~1 | ppCode, data = dataH2B, method = "ML")
summary(ModelH2_1B)
anova(ModelH2_0B,ModelH2_1B)

dataH2C <- subset(dataH2, (Route == "motivate" | Route == "accept"))

ModelH2_0C <- lme(fixed = Help~1, random = ~1 | ppCode, data = dataH2C, method = "ML")
summary(ModelH2_0C)

ModelH2_1C <- lme(fixed = Help~Route, random = ~1 | ppCode, data = dataH2C, method = "ML")
summary(ModelH2_1C)
anova(ModelH2_0C,ModelH2_1C)


#HYPOTHESIS3
dataH3 <- subset(data, (Systemroute == "refer" | Route == "refer"))

ModelH3_0 <- lme(fixed = Help~1, random = ~1 | ppCode, data = dataH3, method = "ML")
summary(ModelH3_0)
anova(ModelH3_0)

#HYPOTHESIS4
dataH4 <- subset(data, (Systemroute == "accept"))

ModelH4_0R <- lme(fixed = Return~1, random = ~1 | ppCode, data = dataH4, method = "ML")
summary(ModelH4_0R)

ModelH4_1R <- lme(fixed = Return~Route, random = ~1 | ppCode, data = dataH4, method = "ML")
summary(ModelH4_1R)
anova(ModelH4_0R, ModelH4_1R)

ModelH4_0A <- lme(fixed = FBH~1, random = ~1 | ppCode, data = dataH4, method = "ML")
summary(ModelH4_0A)

ModelH4_1A <- lme(fixed = FBH~Route, random = ~1 | ppCode, data = dataH4, method = "ML")
summary(ModelH4_1A)
anova(ModelH4_0A, ModelH4_1A)



dataH4B <- subset(dataH4, (Route == "accept" | Route == "motivate"))

View(dataH4B)

ModelH4_0RB <- lme(fixed = Return~1, random = ~1 | ppCode, data = dataH4B, method = "ML")
summary(ModelH4_0RB)

ModelH4_1RB <- lme(fixed = Return~Route, random = ~1 | ppCode, data = dataH4B, method = "ML")
summary(ModelH4_1RB)
anova(ModelH4_0RB, ModelH4_1RB)

ModelH4_0AB <- lme(fixed = FBH~1, random = ~1 | ppCode, data = dataH4B, method = "ML")
summary(ModelH4_0AB)

ModelH4_1AB <- lme(fixed = FBH~Route, random = ~1 | ppCode, data = dataH4B, method = "ML")
summary(ModelH4_1AB)
anova(ModelH4_0AB, ModelH4_1AB)



dataH4C <- subset(dataH4, (Route == "accept" | Route == "refer"))


ModelH4_0RC <- lme(fixed = Return~1, random = ~1 | ppCode, data = dataH4C, method = "ML")
summary(ModelH4_0RC)

ModelH4_1RC <- lme(fixed = Return~Route, random = ~1 | ppCode, data = dataH4C, method = "ML")
summary(ModelH4_1RC)
anova(ModelH4_0RC, ModelH4_1RC)

ModelH4_0AC <- lme(fixed = FBH~1, random = ~1 | ppCode, data = dataH4C, method = "ML")
summary(ModelH4_0AC)

ModelH4_1AC <- lme(fixed = FBH~Route, random = ~1 | ppCode, data = dataH4C, method = "ML")
summary(ModelH4_1AC)
anova(ModelH4_0AC, ModelH4_1AC)


#PRIMARY OUTCOMES
#Comparison Motivation & System

datacon <- read_excel("PATH/Data_conditionsSM.xlsx")
#View(datacon)

#Specify levels of MatchSys as 0 or 1 
datacon$MatchSys = factor(datacon$MatchSys,levels=c(0,1))
#Specify levels for Route as Motiv, Accept or Refer
datacon$Route = factor(datacon$Route,levels=c("motivate","accept","refer"))

#Specify levels for SubSev as Low, Med, High
datacon$SubSev = factor(datacon$SubSev,levels=c("Low","Med","High"))
datacon$SubSev <- ordered(datacon$SubSev)
#Specify levels of SubMoti as Low, Med, High
datacon$SubMoti = factor(datacon$SubMoti,levels=c("Low","Med","High"))
datacon$SubMoti <- ordered(datacon$SubMoti)

#Specify levels for ScenSev as Low, Med, High
datacon$ScenSev = factor(datacon$ScenSev,levels=c("Low","Med","High"))
datacon$ScenSev <- ordered(datacon$ScenSev)
#Specify levels of ScenMoti as Low, Med, High
datacon$ScenMoti = factor(datacon$ScenMoti,levels=c("Low","Med","High"))
datacon$ScenMoti <- ordered(datacon$ScenMoti)

datacon$Condition = factor(datacon$Condition,levels=c("motivate","system"))
datacon$ppCode = factor(datacon$ppCode)

#HELP
Modelc_h <- lme(fixed = Help~1, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Modelc_h)

ggplot(datacon, aes(x=Help, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on Help
Modelc_hCon <- lme(fixed = Help~Condition, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Modelc_hCon)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_h,Modelc_hCon)
anova(Modelc_hCon)

#Return
Modelc_r <- lme(fixed = Return~1, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Modelc_r)

ggplot(datacon, aes(x=Return, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on Return
Modelc_rCon <- lme(fixed = Return~Condition, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Modelc_rCon)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_r,Modelc_rCon)
anova(Modelc_rCon)

#ALLI
Modelc_a <- lme(fixed = FBH~1, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Modelc_a)

ggplot(datacon, aes(x=FBH, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on FBH
Modelc_aCon <- lme(fixed = FBH~Condition, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Modelc_aCon)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_a,Modelc_aCon)
anova(Modelc_aCon)

#Model effect of Condition on FBH
Modelc_aConGen <- lme(fixed = FBH~Condition+Gen, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Modelc_aConGen)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_aCon,Modelc_aConGen)



#Comparison per Accept & Rever route
#Accept
dataconA <- read_excel("PATH/Data_conditionsSMaccept.xlsx")
#View(datacon)

#Specify levels of MatchSys as 0 or 1 
dataconA$MatchSys = factor(dataconA$MatchSys,levels=c(0,1))
#Specify levels for Route as Motiv, Accept or Refer
dataconA$Route = factor(dataconA$Route,levels=c("motivate","accept","refer"))

#Specify levels for SubSev as Low, Med, High
dataconA$SubSev = factor(dataconA$SubSev,levels=c("Low","Med","High"))
dataconA$SubSev <- ordered(dataconA$SubSev)
#Specify levels of SubMoti as Low, Med, High
dataconA$SubMoti = factor(dataconA$SubMoti,levels=c("Low","Med","High"))
dataconA$SubMoti <- ordered(dataconA$SubMoti)

#Specify levels for ScenSev as Low, Med, High
dataconA$ScenSev = factor(dataconA$ScenSev,levels=c("Low","Med","High"))
dataconA$ScenSev <- ordered(dataconA$ScenSev)
#Specify levels of ScenMoti as Low, Med, High
dataconA$ScenMoti = factor(dataconA$ScenMoti,levels=c("Low","Med","High"))
dataconA$ScenMoti <- ordered(dataconA$ScenMoti)

dataconA$Condition = factor(dataconA$Condition,levels=c("motivate","system"))
dataconA$ppCode = factor(dataconA$ppCode)

#HELP
Modelc_hA <- lme(fixed = Help~1, random = ~1 | ppCode, data = dataconA, method = "ML")
summary(Modelc_hA)

ggplot(dataconA, aes(x=Help, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on Help
Modelc_hConA <- lme(fixed = Help~Condition, random = ~1 | ppCode, data = dataconA, method = "ML")
summary(Modelc_hConA)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_hA,Modelc_hConA)
anova(Modelc_hConA)

#Return
Modelc_rA <- lme(fixed = Return~1, random = ~1 | ppCode, data = dataconA, method = "ML")
summary(Modelc_rA)

ggplot(dataconA, aes(x=Return, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on Return
Modelc_rConA <- lme(fixed = Return~Condition, random = ~1 | ppCode, data = dataconA, method = "ML")
summary(Modelc_rConA)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_rA,Modelc_rConA)
anova(Modelc_rConA)

#ALLI
Modelc_aA <- lme(fixed = FBH~1, random = ~1 | ppCode, data = dataconA, method = "ML")
summary(Modelc_aA)

ggplot(dataconA, aes(x=FBH, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on FBH
Modelc_aConA <- lme(fixed = FBH~Condition, random = ~1 | ppCode, data = dataconA, method = "ML")
summary(Modelc_aConA)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_aA,Modelc_aConA)
anova(Modelc_aConA)


#REFER
dataconR <- read_excel("PATH/Data_conditionsSMrefer.xlsx")
#View(datacon)

#Specify levels of MatchSys as 0 or 1 
dataconR$MatchSys = factor(dataconR$MatchSys,levels=c(0,1))
#Specify levels for Route as Motiv, Accept or Refer
dataconR$Route = factor(dataconR$Route,levels=c("motivate","accept","refer"))

#Specify levels for SubSev as Low, Med, High
dataconR$SubSev = factor(dataconR$SubSev,levels=c("Low","Med","High"))
dataconR$SubSev <- ordered(dataconR$SubSev)
#Specify levels of SubMoti as Low, Med, High
dataconR$SubMoti = factor(dataconR$SubMoti,levels=c("Low","Med","High"))
dataconR$SubMoti <- ordered(dataconR$SubMoti)

#Specify levels for ScenSev as Low, Med, High
dataconR$ScenSev = factor(dataconR$ScenSev,levels=c("Low","Med","High"))
dataconR$ScenSev <- ordered(dataconR$ScenSev)
#Specify levels of ScenMoti as Low, Med, High
dataconR$ScenMoti = factor(dataconR$ScenMoti,levels=c("Low","Med","High"))
dataconR$ScenMoti <- ordered(dataconR$ScenMoti)

dataconR$Condition = factor(dataconR$Condition,levels=c("motivate","system"))
dataconR$ppCode = factor(dataconR$ppCode)

#HELP
Modelc_hR <- lme(fixed = Help~1, random = ~1 | ppCode, data = dataconR, method = "ML")
summary(Modelc_hR)

ggplot(dataconR, aes(x=Help, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on Help
Modelc_hConR <- lme(fixed = Help~Condition, random = ~1 | ppCode, data = dataconR, method = "ML")
summary(Modelc_hConR)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_hR,Modelc_hConR)
anova(Modelc_hConR)

#Return
Modelc_rR <- lme(fixed = Return~1, random = ~1 | ppCode, data = dataconR, method = "ML")
summary(Modelc_rR)

ggplot(dataconR, aes(x=Return, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on Return
Modelc_rConR <- lme(fixed = Return~Condition, random = ~1 | ppCode, data = dataconR, method = "ML")
summary(Modelc_rConR)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_rR,Modelc_rConR)
anova(Modelc_rConR)

#ALLI
Modelc_aR <- lme(fixed = FBH~1, random = ~1 | ppCode, data = dataconR, method = "ML")
summary(Modelc_aR)

ggplot(dataconR, aes(x=FBH, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on FBH
Modelc_aConR <- lme(fixed = FBH~Condition, random = ~1 | ppCode, data = dataconR, method = "ML")
summary(Modelc_aConR)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Modelc_aR,Modelc_aConR)
anova(Modelc_aConR)










#REGION WHY NOT
ggplot(data, aes(x=Help, fill=WhyTxt)) + geom_density(alpha = 0.3)
ggplot(data, aes(x=FBH, fill=WhyTxt)) + geom_density(alpha = 0.3)
ggplot(data, aes(x=Return, fill=WhyTxt)) + geom_density(alpha = 0.3)

Model_null_h <- lme(fixed = Help~1, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_null_h)
Model_W_h<- lme(fixed = Help~WhyTxt, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_W_h)
anova(Model_null_h,Model_W_h)

Model_null_a <- lme(fixed = FBH~1, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_null_a)
Model_W_a<- lme(fixed = FBH~WhyTxt, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_W_a)
anova(Model_null_a,Model_W_a)

Model_null_r <- lme(fixed = Return~1, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_null_r)
Model_W_r<- lme(fixed = Return~WhyTxt, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_W_r)
anova(Model_null_r,Model_W_r)




#REGION ISI + age Stats
mean(data$ISITot)
sd(data$ISITot)

mean(data$Age)
sd(data$Age)

#Correlations ISI&Gender + ISI  & Age
datapp <- read_excel("PATH/ppData.xlsx")
View(datapp)
datapp$Gen = factor(datapp$Gen,levels=c("F","M"))
spearman.test(datapp$Gen,datapp$ISITot)
spearman.test(datapp$Age,datapp$ISITot)

#Correlation Return & Allia & HELP
spearman.test(data$FBH,data$Return)
spearman.test(data$Help,data$Return)
spearman.test(data$FBH,data$Help)





#REGION ALLIANCE

#Null model for Agent Attitude
Model_null_aa <- lme(fixed = FBH~1, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_null_aa)

Model_ISI_aa<- lme(fixed = FBH~Age, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_ISI_aa)

anova(Model_null_aa,Model_ISI_aa)

#Model effect of MatchSys on Agent Attitude
Model_MatchSys_aa <- lme(fixed = FBH~MatchSys, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_MatchSys_aa)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_null_aa,Model_MatchSys_aa)

#Model effect of Route on Agent Attitude
Model_route_aa <- lme(fixed = FBH~Route, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_route_aa)
#Vergelijk of model route_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_null_aa,Model_route_aa)

#Model effect of Route & MatchSys on Agent Attitude
Model_routeMatchSys_aa<- lme(fixed = FBH~Route+MatchSys, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_routeMatchSys_aa)
#Vergelijk of model routeMatchSys_aa beter 'MatchSys' heeft dan model route. Lagere logLik is betere MatchSys
anova(Model_route_aa,Model_routeMatchSys_aa)

#Model effect of Route & MatchSys & ScenSev on Agent Attitude
Model_routeMatchSyssev_aa<- lme(fixed = FBH~Route+MatchSys+ScenSev, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_routeMatchSyssev_aa)
#Vergelijk of model routeMatchSyssev_aa beter 'MatchSys' heeft dan model routeMatchSys_aa. Lagere logLik is betere MatchSys
anova(Model_routeMatchSys_aa,Model_routeMatchSyssev_aa)

#Model effect of Route & MatchSys & ScenSev on Agent Attitude
Model_routeMatchSevMoti_aa<- lme(fixed = FBH~Route+MatchSys+ScenSev+ScenMoti, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_routeMatchSevMoti_aa)
#Vergelijk of model routeMatchSyssev_aa beter 'MatchSys' heeft dan model routeMatchSys_aa. Lagere logLik is betere MatchSys
anova(Model_routeMatchSys_aa,Model_routeMatchSevMoti_aa)


#InterHelp effects lme (doesn't work...)
#Model_interHelp1 <- lme(fixed = FBH~Route:ScenSev:ScenMoti, random = ~1 | ppCode, data = data, method = "ML")
#Model_interHelp2 <- lme(fixed = FBH~Route:ScenMoti, random = ~1 | ppCode, data = data, method = "ML")

#InterHelp effects lmer (werkt maar met waarschuwing) - 
Model_interHelp3 <- lmer(FBH~Route:ScenSev:ScenMoti +(1|ppCode), data = data)
summary(Model_interHelp3)

#Full model Route, MatchSys & ScenSev + interaction effects on Agent Attitude
Model_full_aa <- lme(fixed = FBH~Route+ScenSev+ScenMoti+Route:ScenSev+Route:ScenMoti+ScenSev:ScenMoti+Route:ScenSev:ScenMoti, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_full_aa)



#REGION Help

#Null model for Help
Model_null_a <- lme(fixed = Help~1, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_null_a)

Model_ISI_h<- lme(fixed = Help~Age, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_ISI_h)

anova(Model_null_a,Model_ISI_h)

#Model effect of MatchSys on Help
Model_MatchSys_a <- lme(fixed = Help~MatchSys, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_MatchSys_a)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_null_a,Model_MatchSys_a)

#Model effect of Route on Help
Model_route_a <- lme(fixed = Help~Route, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_route_a)
#Vergelijk of model route_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_null_a,Model_route_a)

#Model effect of Route & MatchSys on Help
Model_routeMatchSys_a<- lme(fixed = Help~Route+MatchSys, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_routeMatchSys_a)
#Vergelijk of model routeMatchSys_a beter 'MatchSys' heeft dan model route. Lagere logLik is betere MatchSys
anova(Model_route_a,Model_routeMatchSys_a)

#Model effect of Route & MatchSys & ScenSevon Help
Model_routeMatchSyssev_a<- lme(fixed = Help~Route+MatchSys+ScenSev, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_routeMatchSyssev_a)
#Vergelijk of model routeMatchSyssev_a beter 'MatchSys' heeft dan model routeMatchSys_a. Lagere logLik is betere MatchSys
anova(Model_routeMatchSys_a,Model_routeMatchSyssev_a)

#Model effect of Route & MatchSys & ScenSevon Help
Model_routeMatchSevMoti_a<- lme(fixed = Help~Route+MatchSys+ScenSev+ScenMoti, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_routeMatchSevMoti_a)
#Vergelijk of model routeMatchSyssev_a beter 'MatchSys' heeft dan model routeMatchSys_a. Lagere logLik is betere MatchSys
anova(Model_routeMatchSys_a,Model_routeMatchSevMoti_a)
 
#InterHelp effects lme (doesn't work...)
#Model_interHelp1 <- lme(fixed = Help~Route:ScenSev:ScenMoti, random = ~1 | ppCode, data = data, method = "ML")
#Model_interHelp2 <- lme(fixed = Help~Route:ScenMoti, random = ~1 | ppCode, data = data, method = "ML")

#InterHelp effects lmer (werkt maar met waarschuwing) - 
Model_interHelp3 <- lmer(Help~Route:ScenSev:ScenMoti +(1|ppCode), data = data)
summary(Model_interHelp3)

#Full model Route, MatchSys & ScenSev + interHelp effects on Help
Model_full_a <- lme(fixed = Help~Route+ScenSev+ScenMoti+Route:ScenSev+Route:ScenMoti+ScenSev:ScenMoti+Route:ScenSev:ScenMoti, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_full_a)


#REGION 
#Null model for Return
Model_null_r <- lme(fixed = Return~1, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_null_r)

Model_ISI_r<- lme(fixed = Return~Age, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_ISI_r)

anova(Model_null_r,Model_ISI_r)

#Model effect of MatchSys on Help
Model_MatchSys_a <- lme(fixed = Help~MatchSys, random = ~1 | ppCode, data = data, method = "ML")
summary(Model_MatchSys_a)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_null_a,Model_MatchSys_a)





#Conditions SYSTEM VS MOTIVATION
datacon <- read_excel("PATH/Data_conditionsSM.xlsx")
#View(datacon)

#Specify levels of MatchSys as 0 or 1 
datacon$MatchSys = factor(datacon$MatchSys,levels=c(0,1))
#Specify levels for Route as Motiv, Accept or Refer
datacon$Route = factor(datacon$Route,levels=c("motivate","accept","refer"))

#Specify levels for SubSev as Low, Med, High
datacon$SubSev = factor(datacon$SubSev,levels=c("Low","Med","High"))
datacon$SubSev <- ordered(datacon$SubSev)
#Specify levels of SubMoti as Low, Med, High
datacon$SubMoti = factor(datacon$SubMoti,levels=c("Low","Med","High"))
datacon$SubMoti <- ordered(datacon$SubMoti)

#Specify levels for ScenSev as Low, Med, High
datacon$ScenSev = factor(datacon$ScenSev,levels=c("Low","Med","High"))
datacon$ScenSev <- ordered(datacon$ScenSev)
#Specify levels of ScenMoti as Low, Med, High
datacon$ScenMoti = factor(datacon$ScenMoti,levels=c("Low","Med","High"))
datacon$ScenMoti <- ordered(datacon$ScenMoti)

datacon$Condition = factor(datacon$Condition,levels=c("motivate","system"))
datacon$ppCode = factor(datacon$ppCode)


#HELP

#Null
Model_h <- lme(fixed = Help~1, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_h)

ggplot(datacon, aes(x=Help, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on Help
Model_hCon <- lme(fixed = Help~Condition, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_hCon)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_h,Model_hCon)

#Model effect of Condition on Help
Model_hConISI <- lme(fixed = Help~Condition+Gen, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_hConISI)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_hCon,Model_hConISI)


ggplot(datacon, aes(x=Help, fill=Systemroute)) + geom_density(alpha = 0.3)

#Model effect of Condition+route on Help
Model_hConSys <- lme(fixed = Help~Condition+Systemroute, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_hConSys)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_hCon,Model_hConSys)

#Model effect of Condition+route on Help
Model_hConSysInt <- lme(fixed = Help~Condition+Systemroute+Condition:Systemroute, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_hConSysInt)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_hConSys,Model_hConSysInt)

#Model effect of Condition+route on Help
Model_hConSysInta <- lmer(Help~ScenMoti:ScenSev:Route+ (1|ppCode), data = data)
summary(Model_hConSysInta)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_hConSysInta)



anova(Model_hConSysInt)


#RETURN

#Null
Model_r <- lme(fixed = Return~1, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_r)

ggplot(datacon, aes(x=Return, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on Help
Model_rCon <- lme(fixed = Return~Condition, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_rCon)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_r,Model_rCon)

#Model effect of Condition on Help
Model_rConISI <- lme(fixed = Return~Condition+Gen, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_rConISI)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_rCon,Model_rConISI)

ggplot(datacon, aes(x=Return, fill=Systemroute)) + geom_density(alpha = 0.3)

#Model effect of Condition+route on Help
Model_rConSys <- lme(fixed = Return~Condition+Systemroute, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_rConSys)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_rCon,Model_rConSys)


#ALLIANCE

#Null
Model_a <- lme(fixed = FBH~1, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_a)

ggplot(datacon, aes(x=FBH, fill=Condition)) + geom_density(alpha = 0.3)

#Model effect of Condition on Help
Model_aCon <- lme(fixed = FBH~Condition, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_aCon)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_a,Model_aCon)

#Model effect of Condition on Help
Model_aConISI <- lme(fixed = FBH~Condition+Gen, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_aConISI)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_aCon,Model_aConISI)

ggplot(datacon, aes(x=FBH, fill=Systemroute)) + geom_density(alpha = 0.3)

#Model effect of Condition+route on Help
Model_aConSys <- lme(fixed = FBH~Condition+Systemroute, random = ~1 | ppCode, data = datacon, method = "ML")
summary(Model_aConSys)
#Vergelijk of model MatchSys_aa beter 'MatchSys' heeft dan model test. Lagere logLik is betere MatchSys
anova(Model_aCon,Model_aConSys)