evaluation.R

# this is old code! I do not advocate doing any of the below :-)

rm(list=ls())
library(psych)
count.NAS=function(x) length(which(is.na(x))) 


ASPECT=read.csv("ASPECT.csv")
HONOS=read.csv("HONOS.csv")
ESSENCES=read.csv("ESSENCES.csv")
WOAQ=read.csv("WOAQ.csv")
teammatrix=data.frame(matrix(c("AH", "ATU", "Horizon", "HS", "LC", "LY", "OH", "SC", "SG", "Alexander House", "A/T Unit", "Horizon", "Hillside", "Lindsay C.", "Leylands", "Old Hall", "Stone Cross", "Sycamore Grove"), ncol=2))

ESSENCES2=ESSENCES

ESSENCES$PC=apply(as.data.frame(matrix(c(ESSENCES$item2, ESSENCES$item5, ESSENCES$item8, ESSENCES$item11, ESSENCES$item14), ncol=5)), 1, mean, na.rm=TRUE)
ESSENCES$ES=apply(as.data.frame(matrix(c(ESSENCES$item3, ESSENCES$item6, ESSENCES$item9, ESSENCES$item12, ESSENCES$item15), ncol=5)), 1, mean, na.rm=TRUE)
ESSENCES$TH=apply(as.data.frame(matrix(c(ESSENCES$item4, ESSENCES$item7, ESSENCES$item10, ESSENCES$item13, ESSENCES$item16), ncol=5)), 1, mean, na.rm=TRUE)
ESSENCES$Total17=apply(ESSENCES[8:24], 1, mean, na.rm=TRUE)


HONOS$CogComm=apply(as.data.frame(matrix(c(HONOS$X4, HONOS$X5, HONOS$X6, HONOS$X7), ncol=4)), 1, mean, na.rm=TRUE)
HONOS$Behav=apply(as.data.frame(matrix(c(HONOS$X1, HONOS$X3a, HONOS$X17, HONOS$X2, HONOS$X9), ncol=5)), 1, mean, na.rm=TRUE)
HONOS$LossFunc=apply(as.data.frame(matrix(c(HONOS$X16, HONOS$X12, HONOS$X18, HONOS$X14, HONOS$X15), ncol=5)), 1, mean, na.rm=TRUE)
HONOS$IntDys=apply(as.data.frame(matrix(c(HONOS$X8, HONOS$X3e, HONOS$X10, HONOS$X3c), ncol=4)), 1, mean, na.rm=TRUE)
HONOS$Total=apply(HONOS[which(names(HONOS)==c("X1")):which(names(HONOS)==c("X18"))], 1, mean, na.rm=TRUE)

WOAQ$Total=apply(WOAQ[which(names(WOAQ)==c("work1")):which(names(WOAQ)==c("work28"))], 1, mean, na.rm=TRUE)
WOAQ$Manag=apply(as.data.frame(matrix(c(WOAQ$work3, WOAQ$work5, WOAQ$work7, WOAQ$work11, WOAQ$work16, WOAQ$work17, WOAQ$work22, WOAQ$work26, WOAQ$work27), ncol=9)), 1, mean, na.rm=TRUE)
WOAQ$Colleagues=apply(as.data.frame(matrix(c(WOAQ$work10, WOAQ$work28), ncol=2)), 1, mean, na.rm=TRUE)
WOAQ$Valued=apply(as.data.frame(matrix(c(WOAQ$work12, WOAQ$work13, WOAQ$work14, WOAQ$work21, WOAQ$work23, WOAQ$work24, WOAQ$work25), ncol=7)), 1, mean, na.rm=TRUE)
WOAQ$Environment=apply(as.data.frame(matrix(c(WOAQ$work1, WOAQ$work2, WOAQ$work4, WOAQ$work9, WOAQ$work18, WOAQ$work20), ncol=6)), 1, mean, na.rm=TRUE)
WOAQ$Workload=apply(as.data.frame(matrix(c(WOAQ$work6, WOAQ$work8, WOAQ$work15, WOAQ$work19), ncol=4)), 1, mean, na.rm=TRUE)

#### ASPECT

# draw averages for all

with(ASPECT, barplot(tapply(Average, Unit, mean, na.rm=TRUE), ylim=c(0,6), col=c(rep("blue", 0), rep("red", 1), rep("blue", 7)), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(ASPECT, tapply(Average,list(Unit), mean, na.rm=TRUE)))),2], ylab=c("Average ASPECT item score"), xlab=c("Unit")))
with(ASPECT, barplot(tapply(Average, Unit, mean, na.rm=TRUE), ylim=c(0,6), col=c(rep("blue", 1), rep("red", 1), rep("blue", 6)), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(ASPECT, tapply(Average,list(Unit), mean, na.rm=TRUE)))),2], ylab=c("Average ASPECT item score"), xlab=c("Unit")))
with(ASPECT, barplot(tapply(Average, Unit, mean, na.rm=TRUE), ylim=c(0,6), col=c(rep("blue", 2), rep("red", 1), rep("blue", 5)), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(ASPECT, tapply(Average,list(Unit), mean, na.rm=TRUE)))),2], ylab=c("Average ASPECT item score"), xlab=c("Unit")))
with(ASPECT, barplot(tapply(Average, Unit, mean, na.rm=TRUE), ylim=c(0,6), col=c(rep("blue", 3), rep("red", 1), rep("blue", 4)), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(ASPECT, tapply(Average,list(Unit), mean, na.rm=TRUE)))),2], ylab=c("Average ASPECT item score"), xlab=c("Unit")))
with(ASPECT, barplot(tapply(Average, Unit, mean, na.rm=TRUE), ylim=c(0,6), col=c(rep("blue", 4), rep("red", 1), rep("blue", 3)), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(ASPECT, tapply(Average,list(Unit), mean, na.rm=TRUE)))),2], ylab=c("Average ASPECT item score"), xlab=c("Unit")))
with(ASPECT, barplot(tapply(Average, Unit, mean, na.rm=TRUE), ylim=c(0,6), col=c(rep("blue", 5), rep("red", 1), rep("blue", 2)), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(ASPECT, tapply(Average,list(Unit), mean, na.rm=TRUE)))),2], ylab=c("Average ASPECT item score"), xlab=c("Unit")))
with(ASPECT, barplot(tapply(Average, Unit, mean, na.rm=TRUE), ylim=c(0,6), col=c(rep("blue", 6), rep("red", 1), rep("blue", 1)), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(ASPECT, tapply(Average,list(Unit), mean, na.rm=TRUE)))),2], ylab=c("Average ASPECT item score"), xlab=c("Unit")))
with(ASPECT, barplot(tapply(Average, Unit, mean, na.rm=TRUE), ylim=c(0,6), col=c(rep("blue", 7), rep("red", 1), rep("blue", 0)), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(ASPECT, tapply(Average,list(Unit), mean, na.rm=TRUE)))),2], ylab=c("Average ASPECT item score"), xlab=c("Unit")))

# unit summaries across different scales

barplot(c(mean(ASPECT$Privacy[ASPECT$Unit=="AH"], na.rm=TRUE), mean(ASPECT$Views[ASPECT$Unit=="AH"], na.rm=TRUE),mean(ASPECT$Nature[ASPECT$Unit=="AH"], na.rm=TRUE),mean(ASPECT$Comfort[ASPECT$Unit=="AH"], na.rm=TRUE),mean(ASPECT$Place[ASPECT$Unit=="AH"], na.rm=TRUE),mean(ASPECT$Interior[ASPECT$Unit=="AH"], na.rm=TRUE),mean(ASPECT$Staff[ASPECT$Unit=="AH"], na.rm=TRUE)), names.arg=c("Privacy", "Views", "Nature", "Comfort", "Place", "Interior", "Staff"), col=terrain.colors(7), main=c("Alexander House"), ylim=c(0,6), ylab=c("Average ASPECT item score"), xlab=c("ASPECT score"))

barplot(c(mean(ASPECT$Privacy[ASPECT$Unit=="ATU"], na.rm=TRUE), mean(ASPECT$Views[ASPECT$Unit=="ATU"], na.rm=TRUE),mean(ASPECT$Nature[ASPECT$Unit=="ATU"], na.rm=TRUE),mean(ASPECT$Comfort[ASPECT$Unit=="ATU"], na.rm=TRUE),mean(ASPECT$Place[ASPECT$Unit=="ATU"], na.rm=TRUE),mean(ASPECT$Interior[ASPECT$Unit=="ATU"], na.rm=TRUE),mean(ASPECT$Staff[ASPECT$Unit=="ATU"], na.rm=TRUE)), names.arg=c("Privacy", "Views", "Nature", "Comfort", "Place", "Interior", "Staff"), col=terrain.colors(7), ylim=c(0,6), main=c("Assessment/ Treatment Unit"), ylab=c("Average ASPECT item score"), xlab=c("ASPECT score"))

barplot(c(mean(ASPECT$Privacy[ASPECT$Unit=="HS"], na.rm=TRUE), mean(ASPECT$Views[ASPECT$Unit=="HS"], na.rm=TRUE),mean(ASPECT$Nature[ASPECT$Unit=="HS"], na.rm=TRUE),mean(ASPECT$Comfort[ASPECT$Unit=="HS"], na.rm=TRUE),mean(ASPECT$Place[ASPECT$Unit=="HS"], na.rm=TRUE),mean(ASPECT$Interior[ASPECT$Unit=="HS"], na.rm=TRUE),mean(ASPECT$Staff[ASPECT$Unit=="HS"], na.rm=TRUE)), names.arg=c("Privacy", "Views", "Nature", "Comfort", "Place", "Interior", "Staff"), col=terrain.colors(7), ylim=c(0,6), main=c("Hillside House"), ylab=c("Average ASPECT item score"), xlab=c("ASPECT score"))

barplot(c(mean(ASPECT$Privacy[ASPECT$Unit=="LC"], na.rm=TRUE), mean(ASPECT$Views[ASPECT$Unit=="LC"], na.rm=TRUE),mean(ASPECT$Nature[ASPECT$Unit=="LC"], na.rm=TRUE),mean(ASPECT$Comfort[ASPECT$Unit=="LC"], na.rm=TRUE),mean(ASPECT$Place[ASPECT$Unit=="LC"], na.rm=TRUE),mean(ASPECT$Interior[ASPECT$Unit=="LC"], na.rm=TRUE),mean(ASPECT$Staff[ASPECT$Unit=="LC"], na.rm=TRUE)), names.arg=c("Privacy", "Views", "Nature", "Comfort", "Place", "Interior", "Staff"), col=terrain.colors(7), ylim=c(0,6), main=c("Lindsay Close"), ylab=c("Average ASPECT item score"), xlab=c("ASPECT score"))

barplot(c(mean(ASPECT$Privacy[ASPECT$Unit=="LY"], na.rm=TRUE), mean(ASPECT$Views[ASPECT$Unit=="LY"], na.rm=TRUE),mean(ASPECT$Nature[ASPECT$Unit=="LY"], na.rm=TRUE),mean(ASPECT$Comfort[ASPECT$Unit=="LY"], na.rm=TRUE),mean(ASPECT$Place[ASPECT$Unit=="LY"], na.rm=TRUE),mean(ASPECT$Interior[ASPECT$Unit=="LY"], na.rm=TRUE),mean(ASPECT$Staff[ASPECT$Unit=="LY"], na.rm=TRUE)), names.arg=c("Privacy", "Views", "Nature", "Comfort", "Place", "Interior", "Staff"), col=terrain.colors(7), ylim=c(0,6), main=c("Leylands"), ylab=c("Average ASPECT item score"), xlab=c("ASPECT score"))

barplot(c(mean(ASPECT$Privacy[ASPECT$Unit=="SG"], na.rm=TRUE), mean(ASPECT$Views[ASPECT$Unit=="SG"], na.rm=TRUE),mean(ASPECT$Nature[ASPECT$Unit=="SG"], na.rm=TRUE),mean(ASPECT$Comfort[ASPECT$Unit=="SG"], na.rm=TRUE),mean(ASPECT$Place[ASPECT$Unit=="SG"], na.rm=TRUE),mean(ASPECT$Interior[ASPECT$Unit=="SG"], na.rm=TRUE),mean(ASPECT$Staff[ASPECT$Unit=="SG"], na.rm=TRUE)), names.arg=c("Privacy", "Views", "Nature", "Comfort", "Place", "Interior", "Staff"), col=terrain.colors(7), ylim=c(0,6), main=c("Sycamore Grove"), ylab=c("Average ASPECT item score"), xlab=c("ASPECT score"))

barplot(c(mean(ASPECT$Privacy[ASPECT$Unit=="Horizon"], na.rm=TRUE), mean(ASPECT$Views[ASPECT$Unit=="Horizon"], na.rm=TRUE),mean(ASPECT$Nature[ASPECT$Unit=="Horizon"], na.rm=TRUE),mean(ASPECT$Comfort[ASPECT$Unit=="Horizon"], na.rm=TRUE),mean(ASPECT$Place[ASPECT$Unit=="Horizon"], na.rm=TRUE),mean(ASPECT$Interior[ASPECT$Unit=="Horizon"], na.rm=TRUE),mean(ASPECT$Staff[ASPECT$Unit=="Horizon"], na.rm=TRUE)), names.arg=c("Privacy", "Views", "Nature", "Comfort", "Place", "Interior", "Staff"), col=terrain.colors(7), ylim=c(0,6), main=c("Horizon"), ylab=c("Average ASPECT item score"), xlab=c("ASPECT score"))

#### ESSENCES

barplot(with(ESSENCES, tapply(Total17,list(Time,Unit), mean, na.rm=TRUE)), beside=TRUE, col=rainbow(2), ylim=c(0,4), names.arg=teammatrix[c(1,2,4:9),2], xlab=c("Unit"), ylab=c("Average item score"), main=c("EssenCES average item score"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barx=barplot(c(tapply(ESSENCES$TH[ESSENCES$Unit=="AH"], ESSENCES$Time[ESSENCES$Unit=="AH"], mean, na.rm=TRUE), tapply(ESSENCES$ES[ESSENCES$Unit=="AH"], ESSENCES$Time[ESSENCES$Unit=="AH"], mean, na.rm=TRUE), tapply(ESSENCES$PC[ESSENCES$Unit=="AH"], ESSENCES$Time[ESSENCES$Unit=="AH"], mean, na.rm=TRUE)), main=c("Alexander House"), col=rainbow(2), xaxt="n", ylim=c(0,4), xlab=c("Unit"), ylab=c("Average item score"))
axis(1, at=c(1.3, 3.7, 6.1), labels=c("TH", "ES", "PC"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barx=barplot(c(tapply(ESSENCES$TH[ESSENCES$Unit=="ATU"], ESSENCES$Time[ESSENCES$Unit=="ATU"], mean, na.rm=TRUE), tapply(ESSENCES$ES[ESSENCES$Unit=="ATU"], ESSENCES$Time[ESSENCES$Unit=="ATU"], mean, na.rm=TRUE), tapply(ESSENCES$PC[ESSENCES$Unit=="ATU"], ESSENCES$Time[ESSENCES$Unit=="ATU"], mean, na.rm=TRUE)), main=c("A/T Unit"), col=rainbow(2), xaxt="n", ylim=c(0,4), xlab=c("Unit"), ylab=c("Average item score"))
axis(1, at=c(1.3, 3.7, 6.1), labels=c("TH", "ES", "PC"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barx=barplot(c(tapply(ESSENCES$TH[ESSENCES$Unit=="HS"], ESSENCES$Time[ESSENCES$Unit=="HS"], mean, na.rm=TRUE), tapply(ESSENCES$ES[ESSENCES$Unit=="HS"], ESSENCES$Time[ESSENCES$Unit=="HS"], mean, na.rm=TRUE), tapply(ESSENCES$PC[ESSENCES$Unit=="HS"], ESSENCES$Time[ESSENCES$Unit=="HS"], mean, na.rm=TRUE)), main=c("Hillside House"), col=rainbow(1), names.arg=c("TH", "ES", "PC"), ylim=c(0,4), xlab=c("Unit"), ylab=c("Average item score"))

barx=barplot(c(tapply(ESSENCES$TH[ESSENCES$Unit=="LC"], ESSENCES$Time[ESSENCES$Unit=="LC"], mean, na.rm=TRUE), tapply(ESSENCES$ES[ESSENCES$Unit=="LC"], ESSENCES$Time[ESSENCES$Unit=="LC"], mean, na.rm=TRUE), tapply(ESSENCES$PC[ESSENCES$Unit=="LC"], ESSENCES$Time[ESSENCES$Unit=="LC"], mean, na.rm=TRUE)), main=c("Lindsay Close"), col=rainbow(1), names.arg=c("TH", "ES", "PC"), ylim=c(0,4), xlab=c("Unit"), ylab=c("Average item score"))

barx=barplot(c(tapply(ESSENCES$TH[ESSENCES$Unit=="LY"], ESSENCES$Time[ESSENCES$Unit=="LY"], mean, na.rm=TRUE), tapply(ESSENCES$ES[ESSENCES$Unit=="LY"], ESSENCES$Time[ESSENCES$Unit=="LY"], mean, na.rm=TRUE), tapply(ESSENCES$PC[ESSENCES$Unit=="LY"], ESSENCES$Time[ESSENCES$Unit=="LY"], mean, na.rm=TRUE)), main=c("Leylands"), col=rainbow(2), xaxt="n", ylim=c(0,4), xlab=c("Unit"), ylab=c("Average item score"))
axis(1, at=c(1.3, 3.7, 6.1), labels=c("TH", "ES", "PC"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barx=barplot(c(tapply(ESSENCES$TH[ESSENCES$Unit=="OH"], ESSENCES$Time[ESSENCES$Unit=="OH"], mean, na.rm=TRUE), tapply(ESSENCES$ES[ESSENCES$Unit=="OH"], ESSENCES$Time[ESSENCES$Unit=="OH"], mean, na.rm=TRUE), tapply(ESSENCES$PC[ESSENCES$Unit=="OH"], ESSENCES$Time[ESSENCES$Unit=="OH"], mean, na.rm=TRUE)), main=c("Old Hall"), col=rainbow(2), xaxt="n", ylim=c(0,4), xlab=c("Unit"), ylab=c("Average item score"))
axis(1, at=c(1.3, 3.7, 6.1), labels=c("TH", "ES", "PC"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barx=barplot(c(tapply(ESSENCES$TH[ESSENCES$Unit=="SC"], ESSENCES$Time[ESSENCES$Unit=="SC"], mean, na.rm=TRUE), tapply(ESSENCES$ES[ESSENCES$Unit=="SC"], ESSENCES$Time[ESSENCES$Unit=="SC"], mean, na.rm=TRUE), tapply(ESSENCES$PC[ESSENCES$Unit=="SC"], ESSENCES$Time[ESSENCES$Unit=="SC"], mean, na.rm=TRUE)), main=c("Stone Cross Lane"), col=rainbow(1), names.arg=c("TH", "ES", "PC"), ylim=c(0,4), xlab=c("Unit"), ylab=c("Average item score"))

barx=barplot(c(tapply(ESSENCES$TH[ESSENCES$Unit=="SG"], ESSENCES$Time[ESSENCES$Unit=="SG"], mean, na.rm=TRUE), tapply(ESSENCES$ES[ESSENCES$Unit=="SG"], ESSENCES$Time[ESSENCES$Unit=="SG"], mean, na.rm=TRUE), tapply(ESSENCES$PC[ESSENCES$Unit=="SG"], ESSENCES$Time[ESSENCES$Unit=="SG"], mean, na.rm=TRUE)), main=c("Sycamore Grove"), col=rainbow(1), names.arg=c("TH", "ES", "PC"), ylim=c(0,4), xlab=c("Unit"), ylab=c("Average item score"))

#### HONOS

barplot(with(HONOS, tapply(Total,list(Measure, ID), mean, na.rm=TRUE))[,1:5], beside=TRUE, col=rainbow(2), ylim=c(0,2), ylab=c("Average item score"), xlab=c("Service user code"), main=c("Alexander House"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barplot(with(HONOS, tapply(Total,list(Measure, ID), mean, na.rm=TRUE))[,6:10], beside=TRUE, col=rainbow(2), ylim=c(0,2), ylab=c("Average item score"), xlab=c("Service user code"), main=c("Assessment/ Treatment Unit"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barplot(with(HONOS, tapply(Total,list(Measure, ID), mean, na.rm=TRUE))[,11:17], beside=TRUE, col=rainbow(2), ylim=c(0,2), ylab=c("Average item score"), xlab=c("Service user code"), main=c("Hillside House"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barplot(with(HONOS, tapply(Total,list(Measure, ID), mean, na.rm=TRUE))[,18:23], beside=TRUE, col=rainbow(2), ylim=c(0,3), ylab=c("Average item score"), xlab=c("Service user code"), main=c("Leylands"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barplot(with(HONOS, tapply(Total,list(Measure, ID), mean, na.rm=TRUE))[,24:27], beside=TRUE, col=rainbow(2), ylim=c(0,2), ylab=c("Average item score"), xlab=c("Service user code"), main=c("Old Hall"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

barplot(with(HONOS, tapply(Total,list(Measure, ID), mean, na.rm=TRUE))[,28:40], beside=TRUE, col=rainbow(2), ylim=c(0,2), ylab=c("Average item score"), xlab=c("Service user code"), main=c("Sycamore Grove"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

#### WOAQ

# WOAQ$Total

barplot(with(WOAQ, tapply(Total, list(Time, Unit), mean, na.rm=TRUE)), beside=TRUE, col=rainbow(2), ylim=c(0, 5), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(WOAQ, tapply(Total,list(Unit), mean, na.rm=TRUE)))),2], main=c("WOAQ Total"), xlab=c("Unit"), ylab=c("Average item score"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# WOAQ$Manag

barplot(with(WOAQ, tapply(Manag, list(Time, Unit), mean, na.rm=TRUE)), beside=TRUE, col=rainbow(2), ylim=c(0, 5), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(WOAQ, tapply(Total,list(Unit), mean, na.rm=TRUE)))),2], main=c("WOAQ Manager"), xlab=c("Unit"), ylab=c("Average item score"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# WOAQ$Colleagues

barplot(with(WOAQ, tapply(Colleagues, list(Time, Unit), mean, na.rm=TRUE)), beside=TRUE, col=rainbow(2), ylim=c(0, 5), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(WOAQ, tapply(Total,list(Unit), mean, na.rm=TRUE)))),2], main=c("WOAQ Colleagues"), xlab=c("Unit"), ylab=c("Average item score"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# WOAQ$Valued

barplot(with(WOAQ, tapply(Valued, list(Time, Unit), mean, na.rm=TRUE)), beside=TRUE, col=rainbow(2), ylim=c(0, 5), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(WOAQ, tapply(Total,list(Unit), mean, na.rm=TRUE)))),2], main=c("WOAQ Valued"), xlab=c("Unit"), ylab=c("Average item score"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# WOAQ$Environment

barplot(with(WOAQ, tapply(Environment, list(Time, Unit), mean, na.rm=TRUE)), beside=TRUE, col=rainbow(2), ylim=c(0, 5), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(WOAQ, tapply(Total,list(Unit), mean, na.rm=TRUE)))),2], main=c("WOAQ Environment"), xlab=c("Unit"), ylab=c("Average item score"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# WOAQ$Workload

barplot(with(WOAQ, tapply(Workload, list(Time, Unit), mean, na.rm=TRUE)), beside=TRUE, col=rainbow(2), ylim=c(0, 5), names.arg=teammatrix[which(teammatrix[,1] %in% names(with(WOAQ, tapply(Total,list(Unit), mean, na.rm=TRUE)))),2], main=c("WOAQ Workload"), xlab=c("Unit"), ylab=c("Average item score"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# WOAQ by unit

# WOAQ$Total
# WOAQ$Manag
# WOAQ$Colleagues
# WOAQ$Valued
# WOAQ$Environment
# WOAQ$Workload

barx=barplot(c(tapply(WOAQ$Total[WOAQ$Unit=="AH"], WOAQ$Time[WOAQ$Unit=="AH"], mean, na.rm=TRUE), tapply(WOAQ$Manag[WOAQ$Unit=="AH"], WOAQ$Time[WOAQ$Unit=="AH"], mean, na.rm=TRUE), tapply(WOAQ$Colleagues[WOAQ$Unit=="AH"], WOAQ$Time[WOAQ$Unit=="AH"], mean, na.rm=TRUE), tapply(WOAQ$Valued[WOAQ$Unit=="AH"], WOAQ$Time[WOAQ$Unit=="AH"], mean, na.rm=TRUE), tapply(WOAQ$Environment[WOAQ$Unit=="AH"], WOAQ$Time[WOAQ$Unit=="AH"], mean, na.rm=TRUE), tapply(WOAQ$Workload[WOAQ$Unit=="AH"], WOAQ$Time[WOAQ$Unit=="AH"], mean, na.rm=TRUE)), main=c("Alexander House"), col=rainbow(2), xaxt="n", ylim=c(0,5), ylab=c("Average item score"))
axis(1, at=c(1.3, 3.7, 6.1, 8.5, 10.9, 13.3), labels=c("Total", "Manager", "Colleagues", "Valued", "Environment", "Workload"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# ATU

barx=barplot(c(tapply(WOAQ$Total[WOAQ$Unit=="ATU"], WOAQ$Time[WOAQ$Unit=="ATU"], mean, na.rm=TRUE), tapply(WOAQ$Manag[WOAQ$Unit=="ATU"], WOAQ$Time[WOAQ$Unit=="ATU"], mean, na.rm=TRUE), tapply(WOAQ$Colleagues[WOAQ$Unit=="ATU"], WOAQ$Time[WOAQ$Unit=="ATU"], mean, na.rm=TRUE), tapply(WOAQ$Valued[WOAQ$Unit=="ATU"], WOAQ$Time[WOAQ$Unit=="ATU"], mean, na.rm=TRUE), tapply(WOAQ$Environment[WOAQ$Unit=="ATU"], WOAQ$Time[WOAQ$Unit=="ATU"], mean, na.rm=TRUE), tapply(WOAQ$Workload[WOAQ$Unit=="ATU"], WOAQ$Time[WOAQ$Unit=="ATU"], mean, na.rm=TRUE)), main=c("Assessment and treatment unit"), col=rainbow(2), xaxt="n", ylim=c(0,5), ylab=c("Average item score"))
axis(1, at=c(1.3, 3.7, 6.1, 8.5, 10.9, 13.3), labels=c("Total", "Manager", "Colleagues", "Valued", "Environment", "Workload"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# HS

barx=barplot(c(tapply(WOAQ$Total[WOAQ$Unit=="HS"], WOAQ$Time[WOAQ$Unit=="HS"], mean, na.rm=TRUE), tapply(WOAQ$Manag[WOAQ$Unit=="HS"], WOAQ$Time[WOAQ$Unit=="HS"], mean, na.rm=TRUE), tapply(WOAQ$Colleagues[WOAQ$Unit=="HS"], WOAQ$Time[WOAQ$Unit=="HS"], mean, na.rm=TRUE), tapply(WOAQ$Valued[WOAQ$Unit=="HS"], WOAQ$Time[WOAQ$Unit=="HS"], mean, na.rm=TRUE), tapply(WOAQ$Environment[WOAQ$Unit=="HS"], WOAQ$Time[WOAQ$Unit=="HS"], mean, na.rm=TRUE), tapply(WOAQ$Workload[WOAQ$Unit=="HS"], WOAQ$Time[WOAQ$Unit=="HS"], mean, na.rm=TRUE)), main=c("Hillside House"), col=rainbow(1), xaxt="n", ylim=c(0,5), ylab=c("Average item score"))
axis(1, at=barx, labels=c("Total", "Manager", "Colleagues", "Valued", "Environment", "Workload"))

#LY

barx=barplot(c(tapply(WOAQ$Total[WOAQ$Unit=="LY"], WOAQ$Time[WOAQ$Unit=="LY"], mean, na.rm=TRUE), tapply(WOAQ$Manag[WOAQ$Unit=="LY"], WOAQ$Time[WOAQ$Unit=="LY"], mean, na.rm=TRUE), tapply(WOAQ$Colleagues[WOAQ$Unit=="LY"], WOAQ$Time[WOAQ$Unit=="LY"], mean, na.rm=TRUE), tapply(WOAQ$Valued[WOAQ$Unit=="LY"], WOAQ$Time[WOAQ$Unit=="LY"], mean, na.rm=TRUE), tapply(WOAQ$Environment[WOAQ$Unit=="LY"], WOAQ$Time[WOAQ$Unit=="LY"], mean, na.rm=TRUE), tapply(WOAQ$Workload[WOAQ$Unit=="LY"], WOAQ$Time[WOAQ$Unit=="LY"], mean, na.rm=TRUE)), main=c("Leylands"), col=rainbow(2), xaxt="n", ylim=c(0,5), ylab=c("Average item score"))
axis(1, at=c(1.3, 3.7, 6.1, 8.5, 10.9, 13.3), labels=c("Total", "Manager", "Colleagues", "Valued", "Environment", "Workload"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# OH

barx=barplot(c(tapply(WOAQ$Total[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Manag[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Colleagues[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Valued[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Environment[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Workload[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE)), main=c("Old Hall"), col=rainbow(2), xaxt="n", ylim=c(0,5), ylab=c("Average item score"))
axis(1, at=c(1.3, 3.7, 6.1, 8.5, 10.9, 13.3), labels=c("Total", "Manager", "Colleagues", "Valued", "Environment", "Workload"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

# SC

barx=barplot(c(tapply(WOAQ$Total[WOAQ$Unit=="SC"], WOAQ$Time[WOAQ$Unit=="SC"], mean, na.rm=TRUE), tapply(WOAQ$Manag[WOAQ$Unit=="SC"], WOAQ$Time[WOAQ$Unit=="SC"], mean, na.rm=TRUE), tapply(WOAQ$Colleagues[WOAQ$Unit=="SC"], WOAQ$Time[WOAQ$Unit=="SC"], mean, na.rm=TRUE), tapply(WOAQ$Valued[WOAQ$Unit=="SC"], WOAQ$Time[WOAQ$Unit=="SC"], mean, na.rm=TRUE), tapply(WOAQ$Environment[WOAQ$Unit=="SC"], WOAQ$Time[WOAQ$Unit=="SC"], mean, na.rm=TRUE), tapply(WOAQ$Workload[WOAQ$Unit=="SC"], WOAQ$Time[WOAQ$Unit=="SC"], mean, na.rm=TRUE)), main=c("Stone Cross"), col=rainbow(1), xaxt="n", ylim=c(0,5), ylab=c("Average item score"))
axis(1, at=barx, labels=c("Total", "Manager", "Colleagues", "Valued", "Environment", "Workload"))

# SG

barx=barplot(c(tapply(WOAQ$Total[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Manag[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Colleagues[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Valued[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Environment[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE), tapply(WOAQ$Workload[WOAQ$Unit=="OH"], WOAQ$Time[WOAQ$Unit=="OH"], mean, na.rm=TRUE)), main=c("Sycamore Grove"), col=rainbow(2), xaxt="n", ylim=c(0,5), ylab=c("Average item score"))
axis(1, at=c(1.3, 3.7, 6.1, 8.5, 10.9, 13.3), labels=c("Total", "Manager", "Colleagues", "Valued", "Environment", "Workload"))
legend("topright", c("Time 1", "Time 2"), fill=rainbow(2))

#### Decriptives

ASPECT$Total=ASPECT$Average*39

HONOS$CogComm2=apply(as.data.frame(matrix(c(HONOS$X4, HONOS$X5, HONOS$X6, HONOS$X7), ncol=4)), 1, sum, na.rm=TRUE)/(4-apply(as.data.frame(matrix(c(HONOS$X4, HONOS$X5, HONOS$X6, HONOS$X7), ncol=4)), 1, count.NAS))*4

HONOS$Behav2=apply(as.data.frame(matrix(c(HONOS$X1, HONOS$X3a, HONOS$X17, HONOS$X2, HONOS$X9), ncol=5)), 1, sum, na.rm=TRUE)/(5-apply(as.data.frame(matrix(c(HONOS$X1, HONOS$X3a, HONOS$X17, HONOS$X2, HONOS$X9), ncol=5)), 1, count.NAS))*5
HONOS$LossFunc2=apply(as.data.frame(matrix(c(HONOS$X16, HONOS$X12, HONOS$X18, HONOS$X14, HONOS$X15), ncol=5)), 1, sum, na.rm=TRUE)/(5-apply(as.data.frame(matrix(c(HONOS$X16, HONOS$X12, HONOS$X18, HONOS$X14, HONOS$X15), ncol=5)), 1, count.NAS))*5
HONOS$IntDys2=apply(as.data.frame(matrix(c(HONOS$X8, HONOS$X3e, HONOS$X10, HONOS$X3c), ncol=4)), 1, sum, na.rm=TRUE)/(4-apply(as.data.frame(matrix(c(HONOS$X8, HONOS$X3e, HONOS$X10, HONOS$X3c), ncol=4)), 1, count.NAS))*4
HONOS$Total2=apply(HONOS[which(names(HONOS)==c("X1")):which(names(HONOS)==c("X18"))], 1, sum, na.rm=TRUE)/(22-apply(HONOS[which(names(HONOS)==c("X1")):which(names(HONOS)==c("X18"))], 1, count.NAS))*22

ESSENCES$PC2=apply(as.data.frame(matrix(c(ESSENCES$item2, ESSENCES$item5, ESSENCES$item8, ESSENCES$item11, ESSENCES$item14), ncol=5)), 1, sum, na.rm=TRUE)/(5-apply(as.data.frame(matrix(c(ESSENCES$item2, ESSENCES$item5, ESSENCES$item8, ESSENCES$item11, ESSENCES$item14), ncol=5)), 1, count.NAS))*5

ESSENCES$ES2=apply(as.data.frame(matrix(c(ESSENCES$item3, ESSENCES$item6, ESSENCES$item9, ESSENCES$item12, ESSENCES$item15), ncol=5)), 1, sum, na.rm=TRUE)/(5-apply(as.data.frame(matrix(c(ESSENCES$item3, ESSENCES$item6, ESSENCES$item9, ESSENCES$item12, ESSENCES$item15), ncol=5)), 1, count.NAS))*5

ESSENCES$TH2=apply(as.data.frame(matrix(c(ESSENCES$item4, ESSENCES$item7, ESSENCES$item10, ESSENCES$item13, ESSENCES$item16), ncol=5)), 1, sum, na.rm=TRUE)/(5-apply(as.data.frame(matrix(c(ESSENCES$item4, ESSENCES$item7, ESSENCES$item10, ESSENCES$item13, ESSENCES$item16), ncol=5)), 1, count.NAS))*5
ESSENCES$Total17_2=apply(ESSENCES[8:24], 1, sum, na.rm=TRUE)/(17-apply(ESSENCES[8:24], 1, count.NAS))*17

WOAQ$Total2=apply(WOAQ[which(names(WOAQ)==c("work1")):which(names(WOAQ)==c("work28"))], 1, sum, na.rm=TRUE)/(28-apply(WOAQ[which(names(WOAQ)==c("work1")):which(names(WOAQ)==c("work28"))], 1, count.NAS))*28

WOAQ$Manag2=apply(as.data.frame(matrix(c(WOAQ$work3, WOAQ$work5, WOAQ$work7, WOAQ$work11, WOAQ$work16, WOAQ$work17, WOAQ$work22, WOAQ$work26, WOAQ$work27), ncol=9)), 1, sum, na.rm=TRUE)/(9-apply(as.data.frame(matrix(c(WOAQ$work3, WOAQ$work5, WOAQ$work7, WOAQ$work11, WOAQ$work16, WOAQ$work17, WOAQ$work22, WOAQ$work26, WOAQ$work27), ncol=9)), 1, count.NAS))*9

WOAQ$Colleagues2=apply(as.data.frame(matrix(c(WOAQ$work10, WOAQ$work28), ncol=2)), 1, sum, na.rm=TRUE)/(2-apply(as.data.frame(matrix(c(WOAQ$work10, WOAQ$work28), ncol=2)), 1, count.NAS))*2

WOAQ$Valued2=apply(as.data.frame(matrix(c(WOAQ$work12, WOAQ$work13, WOAQ$work14, WOAQ$work21, WOAQ$work23, WOAQ$work24, WOAQ$work25), ncol=7)), 1, sum, na.rm=TRUE)/(7-apply(as.data.frame(matrix(c(WOAQ$work12, WOAQ$work13, WOAQ$work14, WOAQ$work21, WOAQ$work23, WOAQ$work24, WOAQ$work25), ncol=7)), 1, count.NAS))*7

WOAQ$Environment2=apply(as.data.frame(matrix(c(WOAQ$work1, WOAQ$work2, WOAQ$work4, WOAQ$work9, WOAQ$work18, WOAQ$work20), ncol=6)), 1, sum, na.rm=TRUE)/(6-apply(as.data.frame(matrix(c(WOAQ$work1, WOAQ$work2, WOAQ$work4, WOAQ$work9, WOAQ$work18, WOAQ$work20), ncol=6)), 1, count.NAS))*6

WOAQ$Workload2=apply(as.data.frame(matrix(c(WOAQ$work6, WOAQ$work8, WOAQ$work15, WOAQ$work19), ncol=4)), 1, sum, na.rm=TRUE)/(4-apply(as.data.frame(matrix(c(WOAQ$work6, WOAQ$work8, WOAQ$work15, WOAQ$work19), ncol=4)), 1, count.NAS))*4

describe.by(ASPECT, ASPECT$Unit)

describe.by(ESSENCES, ESSENCES$Unit)

ESSENCES2$ES-ESSENCES$ES2

describe.by(HONOS, HONOS$Unit)

# ESSENCES
with(subset(ESSENCES, Time==1), tapply(Total, Unit, mean, na.rm=TRUE))
with(subset(ESSENCES, Time==1), tapply(Total, Unit, sd, na.rm=TRUE))

with(subset(ESSENCES, Time==1), tapply(PC, Unit, mean, na.rm=TRUE))
with(subset(ESSENCES, Time==1), tapply(PC, Unit, sd, na.rm=TRUE))

with(subset(ESSENCES, Time==1), tapply(TH, Unit, mean, na.rm=TRUE))
with(subset(ESSENCES, Time==1), tapply(TH, Unit, sd, na.rm=TRUE))

with(subset(ESSENCES, Time==1), tapply(ES, Unit, mean, na.rm=TRUE))
with(subset(ESSENCES, Time==1), tapply(ES, Unit, sd, na.rm=TRUE))

# WOAQ
with(subset(WOAQ, Time==1), tapply(Total2, Unit, mean, na.rm=TRUE))
with(subset(WOAQ, Time==1), tapply(Total2, Unit, sd, na.rm=TRUE))

with(subset(WOAQ, Time==1), tapply(Manag2, Unit, mean, na.rm=TRUE))
with(subset(WOAQ, Time==1), tapply(Manag2, Unit, sd, na.rm=TRUE))

with(subset(WOAQ, Time==1), tapply(Colleagues2, Unit, mean, na.rm=TRUE))
with(subset(WOAQ, Time==1), tapply(Colleagues2, Unit, sd, na.rm=TRUE))

with(subset(WOAQ, Time==1), tapply(Valued2, Unit, mean, na.rm=TRUE))
with(subset(WOAQ, Time==1), tapply(Valued2, Unit, sd, na.rm=TRUE))

with(subset(WOAQ, Time==1), tapply(Environment2, Unit, mean, na.rm=TRUE))
with(subset(WOAQ, Time==1), tapply(Environment2, Unit, sd, na.rm=TRUE))

with(subset(WOAQ, Time==1), tapply(Workload2, Unit, mean, na.rm=TRUE))
with(subset(WOAQ, Time==1), tapply(Workload2, Unit, sd, na.rm=TRUE))

# HONOS
with(subset(HONOS, Measure==1), tapply(Total2, Unit, mean, na.rm=TRUE))
with(subset(HONOS, Measure==1), tapply(Total2, Unit, sd, na.rm=TRUE))

with(subset(HONOS, Measure==2), tapply(Total2, Unit, mean, na.rm=TRUE))
with(subset(HONOS, Measure==2), tapply(Total2, Unit, sd, na.rm=TRUE))

with(subset(HONOS, Measure==1), tapply(CogComm2, Unit, mean, na.rm=TRUE))
with(subset(HONOS, Measure==1), tapply(CogComm2, Unit, sd, na.rm=TRUE))

with(subset(HONOS, Measure==2), tapply(CogComm2, Unit, mean, na.rm=TRUE))
with(subset(HONOS, Measure==2), tapply(CogComm2, Unit, sd, na.rm=TRUE))

with(subset(HONOS, Measure==2), tapply(Behav2, Unit, mean, na.rm=TRUE))
with(subset(HONOS, Measure==2), tapply(Behav2, Unit, sd, na.rm=TRUE))

with(subset(HONOS, Measure==2), tapply(Behav2, Unit, mean, na.rm=TRUE))
with(subset(HONOS, Measure==2), tapply(Behav2, Unit, sd, na.rm=TRUE))

### mean differences in unit variables by ANOVA

# ASPECT

ASPECTmodel=aov(Average~Unit, data=ASPECT)
summary(ASPECTmodel)
TukeyHSD(ASPECTmodel)

# ESSENCES

ESSENCESmodel=aov(Total~Unit, data=ESSENCES)
summary(ESSENCESmodel)
TukeyHSD(ESSENCESmodel)

ESSENCESmodel=aov(PC~Unit, data=ESSENCES)
summary(ESSENCESmodel)

ESSENCESmodel=aov(TH~Unit, data=ESSENCES)
summary(ESSENCESmodel)

ESSENCESmodel=aov(ES~Unit, data=ESSENCES)
summary(ESSENCESmodel)

# HONOS time 1

HONOSmodel=aov(Total~Unit, data=subset(HONOS, Measure==1))
summary(HONOSmodel)
TukeyHSD(HONOSmodel)

# HONOS time 2

HONOSmodel=aov(Total~Unit, data=subset(HONOS, Measure==2))
summary(HONOSmodel)
TukeyHSD(HONOSmodel)

# WOAQ

WOAQmodel=aov(Total~Unit, data=WOAQ)
summary(WOAQmodel)
TukeyHSD(WOAQmodel)

WOAQmodel=aov(Manag~Unit, data=WOAQ)
summary(WOAQmodel)

WOAQmodel=aov(Colleagues~Unit, data=WOAQ)
summary(WOAQmodel)

WOAQmodel=aov(Environment~Unit, data=WOAQ)
summary(WOAQmodel)

WOAQmodel=aov(Workload~Unit, data=WOAQ)
summary(WOAQmodel)

### ICC
library(psych)

# HONOS

alpha(HONOS[7:28], na.rm=TRUE)

# WOAQ

alpha(WOAQ[7:34], na.rm=TRUE)