NOTE: Run all the code included on the “Base Code” page prior to running the following.
Load eigenworms and worm data of choice.
# working directory
setwd("") # enter appropriate working directory
# load eigenworms [Broekmans et al. 2016]
ew = read.csv("eigenworms.csv", header=F, sep=",")
# load coefficients from worm of choice
data = read.csv("12 Foraging Worms/1.txt", header=F, sep="")[10001:12000,1:5]
These are the optimal parameters for the data produced by Broekmans et al. 2016.
E = 5 # embedding dimension
tp = 1
theta = 2 # linearity
tau = 1
# library and prediction selection
lib = c(1,1000)
pred = c(1001,2000)
# convert frames to seconds (i.e. indicate frames per second as fps)
fps = 16 # foraging worms
Run the functions created on the “Base Code” page.
# run the embedding function for separate lib and pred
matrix1 <- make_embed(data[lib[1]:lib[2],], E, tau, tp)
matrix2 <- make_embed(data[pred[1]:pred[2],], E, tau, tp)
matrix = rbind(matrix1, matrix2)
# run the prediction function
new_pred <- make_pred(E, matrix, theta, lib, pred)
observations_total = new_pred$obs
predictions_total = new_pred$pred
# run the EDM error function for sections
errors <- edm_error_section(pred, observations_total, ew, predictions_total)
# run the cp error function
cp_mean <- cp_error(pred, data, ew)
Create an RMS error vs time plot, broken into worm body angle segments.
# RMS error by body segment for all 5 segments
par(mai=c(0.9,0.9,0.15,0.15))
plot(c(1:length(errors[,1]))/fps, errors[,1], col="purple", type="l", xlab="time (s)",
ylab="prediction error (rad)", cex.lab=1.5, cex.axis=1.5)
lines(c(1:length(errors[,2]))/fps, errors[,2], col="blue")
lines(c(1:length(errors[,3]))/fps, errors[,3], col="green")
lines(c(1:length(errors[,4]))/fps, errors[,4], col="orange")
lines(c(1:length(errors[,5]))/fps, errors[,5], col="red")
abline(h=cp_mean, lty="dashed", col="grey40") # average cp error
# add a plot legend
legend(53,0.12, title="worm body angles", lty=1, cex=1, bty="n",
legend=c("1-20", "21-40", "41-60", "61-80", "81-100"),
col=c("purple", "blue", "green", "orange", "red"))