library(cusp)
attach(S4_Cusp_data)
names(S4_Cusp_data)
# Run seven models
fit1 <- cusp(y~uSB+cFV+uME+cAN,alpha~f+b+m,beta~f+b+m,data=Cusp_data)
fit2 <- cusp(y~uSB+cFV+uME+cAN,alpha~f,beta~b+m,data=Cusp_data)
fit3 <- cusp(y~uSB+cFV+uME+cAN,alpha~b+m,beta~f,data=Cusp_data)
fit4 <- cusp(y~uSB+cFV+uME+cAN,alpha~f+b,beta~m,data=Cusp_data)
fit5 <- cusp(y~uSB+cFV+uME+cAN,alpha~m,beta~f+b,data=Cusp_data)
fit6 <- cusp(y~uSB+cFV+uME+cAN,alpha~f,beta~m,data=Cusp_data)
fit7 <- cusp(y~uSB+cFV+uME+cAN,alpha~m,beta~f,data=Cusp_data)
# Print out the AICc for the seven models
A=numeric(7)
A[1] <- summary(fit1)[["r2cusp.aicc"]]
A[2] <- summary(fit2)[["r2cusp.aicc"]]
A[3] <- summary(fit3)[["r2cusp.aicc"]]
A[4] <- summary(fit4)[["r2cusp.aicc"]]
A[5] <- summary(fit5)[["r2cusp.aicc"]]
A[6] <- summary(fit6)[["r2cusp.aicc"]]
A[7] <- summary(fit7)[["r2cusp.aicc"]]
A
# Model 6 has smallest AICc, so get coefficients and plot
summary(fit6)
plot(fit6)
cusp3d(fit6, B=6, Y=3, Yfloor= -5, w=.06,surf.alpha=.4, surf.gamma=2, n.surface=50, surf.chroma=90)

# Get alpha and beta values for data
predict(fit6)
# Get predicted position of Y
fitted(fit6)
#Get the boundaries of the folds
B <- seq(0,4,by=.02)
cusp.bifset(B)