# Model evaluation: Class work

setwd("C:/Users/Tim/Desktop/NANP Nutrition Workshop files/Exercise files")

data = read.csv(file = "Lesson 4 exercise data file.csv")
names(data) = c("Obs", "Pred")

#Always plot your data first to get a feel for the data
plot(data$Pred, data$Obs, type='p')
abline (0, 1)

obs=data$Obs
pred=data$Pred

#You need a funciton to write the equation.
#MSEP

if(length(obs) != length(pred)) stop("The length of obs and pred are not equal or incorrect names")
MSEP = sum((obs-pred)^2)/length(pred)
RMSEP = sqrt(MSEP)
STDEVP <- sd(obs)*sqrt((length(obs)-1)^2/length(obs)^2)
RSR = RMSEP/STDEVP
CCC = (2*cov(obs,pred))/
  (sd(obs)^2+sd(pred)^2 + (mean(obs)-mean(pred))^2)
SEobs = sqrt(sum((obs-mean(obs))^2)/length(obs))
SEpred = sqrt(sum((pred-mean(pred))^2)/length(pred))
r = sum((pred - mean(pred))*(obs - mean(obs)))/(length(obs)*SEobs*SEpred) 
v = sd(obs)/sd(pred)
u = (mean(obs)-mean(pred))/(sqrt(sd(obs)*sd(pred)))
Cb = 2/(v + 1/v + u^2)
Results = matrix(c(MSEP, RMSEP, RSR, CCC, r, Cb ), byrow = TRUE, ncol = 6)
colnames(Results) = c("MSEP", "RMSEP", "RSR",  "CCC", "r", "Cb")

Results