How can I replicate this plot using these differential equations and parameters?

How can I replicate the plot above given these parameters?

[Torsten]

How can I replicate the plot above given these parameters?

By using ODE45 to solve the 5 differential equations.

[faraz khan]

I'm kind of confused on how to start. How would I enter the differential equations into MATLAB?

[Torsten]

Take a look at the example "Solve Nonstiff Equation" under

https://de.mathworks.com/help/matlab/ref/ode45.html

It should not be a problem to enter five instead of two differential equations.

[faraz khan]

This doesn't have any information regarding entering ODE's that look similar to mine. I still don't understand how to enter them. It would be very helpful if you could use an example

[Sam Chak]

Aha.. I see. You are looking for an example that has the same number of ODEs and lots of rational terms.

The second example "Solve Non-stiff ODEs" is relevant. Since you want 5, this example has 5.

https://www.mathworks.com/help/matlab/math/solve-stiff-transistor-dae.html

Answer 1

[William Rose]

@faraz khan,

@Torsten is exactly correct, as usual. But since you are still stuck, here's some additional help.

Define a state vector x with 5 elements. The five elements correspond to the 5 variables:

, , , , ,

Then write the five differential equations in terms of the elements of the state vector. For example, the first diff. eq, for dm/dt, becomes

and the last diff.eq., for pN(t), becomes

,

and you can work out the three in between. Use the differential equations you have written to create a function that returns a column vector of length 5, called dxdt. You must define the various constants inside the function, or make them global. (I recommend defining them inside the function, because global variables can cause slow execution.)

function dxdt = circOsc(t,x)
vs=0.76; vm=0.65; KI=1; Km1=0.5;
k1=1.9; k2=1.3; n=4;
%Add definitions of additional constants which will be used below.
ks=0; V1=0; V2=0;  %fix these and add the others

dxdt = zeros(5,1);  %allocate column vector dxdt
dxdt(1) = vs/(1+(x(5)/KI)^n)-vm*x(1)/(Km1+x(1));
dxdt(2) = 0; %fill this in with the right equation
dxdt(3) = 0; %fill this in with the right equation
dxdt(4) = 0; %fill this in with the right equation
dxdt(5) = k1*x(4)-k2*x(5);
end

You can save this function as a file, circOsc.m, or you can include circOsc() as code after the body of your main script, if you prefer to have everything in a single file. The main script calls ode45(). You must pass to ode45() the initial condition (vector x0), the time span (0 to 70 hours, to reproduce Fig.7.19A), and the name of the function that computes the derivative of the state vector (circOsc, in the example above). All of this is illustrated clearly in the many examples associated with the ode45() help. The van der Pol oscillator example in the ode45 help is probably the most similar to your model.

%Goldbeter.m
%Simulate the Goldbeter circadian oscillator
tspan=[0,70];
x0=[2.5; 0.5; 0.5; 0.5; 0.7];   %estimated from Fig. 7.19A, assuming p0=p1=p2 at t=0, which may not be true
[t,x]=ode45(@circOsc,tspan,x0);
%% Display results graphically
figure;
plot(t,x(:,1),'-r',t,x(:,2),'-g',t,x(:,3),'-b',t,x(:,4),'-c',t,x(:,5),'-m');
legend('m','p_0','p_1','p_2','p_N'); grid on
xlabel('Time (h)'); ylabel('Concentration (\muM)'); 

Try it. Check my equations and the constants, since I might have made mistakes. Check everything and add the missing bits.

[faraz khan]

it worked, thanks so much. I never thought to use circOsc.

[William Rose]

@faraz khan, great! Good luck.