function ydot = algaesys(t,y,theta,xdata) % ode system function for MCMC algae example A = y(1); Z = y(2); P = y(3); % control variables are assumed to be saved % at each time unit interval QpV = xdata(ceil(t),2); T = xdata(ceil(t),3); Pin = xdata(ceil(t),4); mumax = theta(1); rhoa = theta(2); rhoz = theta(3); k = theta(4); alpha = theta(5); th = theta(6); mu = mumax*th^(T-20)*P/(k+P); dotA = (mu - rhoa - QpV -alpha*Z)*A; dotZ = alpha*Z*A - rhoz*Z; dotP = -QpV*(P-Pin) + (rhoa-mu)*A + rhoz*Z; ydot=[dotA;dotZ;dotP];