Hi. We are now in our fourth lecture on mathematical modeling of the cell cycle, and this our last lecture in this series. What we're going to talk about today in this lecture Are some of the prac, practical aspects of implementing this mathematical model we've been talking about. So we've been discussing the 1993 Novak and Tyson model of the cell cycle. We're going to talk a little bit about the model equations and, and the notation that we use. And then I'm going to discuss some of the practical considerations of implementing this model in MATLAB. And that's going to force us to review a concept from one of the early lectures on introduction to MATLAB the idea of MATLAB scripts versus MATLAB functions. Because that's critical to how we have implemented this Novak and Tyson cell cycle model. and then that's going to allow us to discuss how you can use MATLAB's built-in ordinary differential equation solvers. And we're gong discuss that in the context of, of scripts versus functions. The Novak and Tyson cell cycle model as we've discussed, consists of seven ordinary differential equations. For the seven unique molecular species in this model. By unique in this case I mean is there are some species that can be either the phospholated form or the unphospholated form. For instance we want one of these examples, but because the total amount of We one is constant, you only have to keep track of either the phospholated or the unphosfolated. And you can always get the whatever the, the other state is by taking the total amount and subtracting that. So, we have, only have one equation for the phosphorylated we want because if we want to calculate unphosphorylated we want, it's just this term here. The total we want minus the phosphorylated form. So once Novak and Tyson made those simplifications based on conservation of mass, they come up, came up with seven ODEs for seven molecular species. And we reviewed structure of a couple of these equations in, in one of our previous lectures. But I just have all the equations here for completeness, and this set of equations came from this very nice review article by Sible and Tyson that we've encountered previously. Although we've discussed some of the equations in the Novak-Tyson 1993 cell cycle model. But we haven't spent a lot of time discussing all of the parameters of the model. All of the different numbers that have to be defined in order to implement this model and, and to run it. And I'm just going to summarize these here in, in a few categories. There are a bunch of maximal rates. These are usually kcats of enzymes. They can also be rates of of synthesis or degradation. And these are, are denoted with, with small ks here. So you've got little k1, little k3, little k a, b, c etcetera. And these are the values that are used in the implementation, and that we're going to provide you with. There are also a, a number of Michaelis constants, and these are denoted with a capital K. Rather than lowercase k, over here. Then something else that needs to be defined in the model are a number of total protein concentrations, total amount of CDK cdc25. These are the species that don't get synthesized or degraded. We want intermediate. And so. And you see that cycling, for instance, is not out here as a parameter of total protein concentration. That's because cycling gets synthesized, cycling gets degraded. So that's a variable that changes as a function of time. These are total amounts of the proteins that, that don't change as a function of time. And therefore, whatever the, the mac, the total amount is has to be defined. At the beginning before the simulation runs. The Novak-Tyson model also can aide several parameters that I've grouped down here that I call weighting parameters. What I mean by that is that some of the rate constants are weighted averages that are calculated using these parameters down here. For instance Kwee, which is one of the rate constants that involves Wee1, depends on the phosphorylated amount of phosphorylated Wee1. It also depends on the amount of non-phosphorylated Wee1. And these are grouped together, they're weighted together based on these two terms Vwee. Prime, V Wee prime [INAUDIBLE]. These, this is a notation used by Tyson. The, the terminology that I used in the MATLAB implementation is vwee_1 and vwee_2. And these are the two that are used to, to weight the phosphorylated wee1, the non-phosphorylated wee1. You come up with this overall rate constant here. So you can see that there's a, a number or parameters that are necessary to simulate the Novak-Tyson cell cycle model. Some are k cats, some are Michaelis constants, some are total protein concentrations, and a series of them are these weighting parameters. Now let's talk, more specifically. About some of the things we need to consider when implementing the Novac Tyson model in the MATLAB programming environment. We encountered this slide before when we discussed Oyler's method in the context of dynamical systems, and we said this is the general way you structure a script that implements a mathematical model that solves a system of ODEs. we, first we defined constants, then we set the time step and the simulation time, then we defined the initial conditions for our system. We had a for loop that implemented Euler's method. And at each time step in the for loop, we saved the output into some vector that kept all the output, or some array that kept all the output. We computed the derivative of all the variables with respect to time, and then we computed all the variables at the next time step. And when our "for" loop was over, we plot it and we output the results. And when we encountered this before, this was a specific example that we discussed. Very simple ordinary differential equation, dx dt equals a minus bx for these values of little a, little b, and little c. And we wrote this MATLAB script that we called euler.m. This solves this differential equation numerically. But, one of the things we discussed when we talked about Euler's method before. Is that, this is really that the simplest and most straight forward way to solve ODE systems. In the over 200 years since Euler died there have been many improvements made of, of ways to solve ODE systems numerically. And we talked a little bit about the Runge-Kutta method or the variable time-step methods. And what I've discussed last time is I said these algorithms are available in MATLAB as the built-in ODE solvers. So we dis-, I think it was worth talking about Euler's method. So that, conceptually, you can understand how ODE systems are, are solved. But, in practice Euler's method is not usually used. Usually you'll use one of these more sophisticated algorithms. And that's what we're going to do to to solve the Novak and Tyson model. But we're not going to write Sophisticated algorithm ourselves. We're going to exploit the fact that MATLAB has already done this, and these ode solvers that are, that are provided with the software, things like ode 23. Ode 15s. Now, when we use MATLAB's built-in solvers to solve to integrate ode systems, we're going to have to make a couple of changes to the structure of our program. We're still going to find, find the constants at the beginning. We're going to set the initial conditions in set three and we're going to apply and output the results in step five. But you can see steps two and four are changed a little bit. For one thing we're not going to explicitly set a time step like we did when we used Oiler's method. MATLAB is going to compute the steps automatically for us. And then at step four, rather than having a four loop. That goes through, calculate the derivatives, integrates moves to the next time step. We're going to do all of this using MATLAB's solvers. And so step four is going to be altered a little bit when we do this for the Novak and Tyson model. But in order to do number four, we're going to have to write a function that is going to compute the set of derivatives in our In our model. Given the current collection of state variables. And that term function that I just mentioned should sound familiar to you, because that's something we also discussed in the introduction to MATLAB lectures at the very beginning. This is the slide we showed back then to illustrate the schematic relationship. Between, between scripts and functions. Remember that if you call a function, sample function, A, B, C. What this is going is it's taking the variables, A, B, and C, that are defined in your script. Sending them to your function, sample function .m. Sample function.m is going to perform a lot of calculations, some of which are going to involve the vari the variables, A, B, and C that you sent. And then it's going to return results. And a couple of things we discussed when we talked about scripts versus functions is that the variables that are called capital A, capital B and capital C in the script are called x, y and z in the function and that's okay. It's okay for the variables to have different names in the script and in the function. And similarly, sample function is going to calculate something called little a and little a is the thing that is going to get returned to whatever whatever the function was called. And when, when it comes back to this script it's going to be called results and, again, it's a different name but this is okay. So, these are things that we discussed in the beginning when we talked about the differences between MATLAB's scripts and MATLAB functions. And another lesson we learned when we talked about scripts versus function is, is that after the function is called the variables defined within the function are gone. So, we called sample function. We sent it these variables: capital A, capital B, capital C. It referred to these variables as little x, little y, and little z, but then by the time the function is over little x, little y, and little z didn't exist anymore. Why did we just go back and review the difference between MATLAB scripts and MATLAB functions? Well, this distinction is critical to understand in order to use MATLAB's built-in solvers to integrate systems of, of ODEs. For one thing the ode solvers are themselves functions. So ode23 or ode15s is a MATLAB function, and in order to know how to use a MATLAB function you need to know what variables should get sent to the, to the MATLAB function. For instance this is for instance a way that you could. Integrate the Novak-Tyson and cell cycle model, using Matlab's solver ode23. So let's review what these, what these variables represent here. You call ode23 and then you have a parenthesis and you say, at dydt_novak comma. Brack, you know, square bracket, zero comma t last, end square bracket, and then a comma and then state variable underscore i. The three things that are getting sent to ODE two three. The first is what's called a function handle. And the way that you know this is a function handle, rather than just say, a regular old variable, It's because of this at symbol in front of it. So dydt_novak is not just a variable, it's equal to the number four or something like that. This at symbol in front of it, tells MATLAB that that you should look fora file, called dydt_novak. So that at symbol here is is a special code. That MatLab understands that this is not just a regular old variable, this is a file I should look for. The second variable that is getting sent here, you can probably deduce this looking at the variable names I have choose. This is saying the time over which to perform this integration start at time 0, end at time t last. Meaning it's defined that way, because it's like the last time. And then this final one, is the set of initial conditions. So if you have a bunch of variables in your model, that are sometimes called state variables of that model, then statevar_i is the variable name that represents initial conditions for the state variable. Little i, in this case Meaning initial. And then what are the variables that get computed by ode23, and get sent back to the script, he call to ode23? Well, you can see that there's two variables that are going to get sent back. One is time, this is going to be a vector with all the time points in it. And then statevars is going to be a two dimensional array. Where each row is going to represent a different time point, and each column and statevars is going to represent a different variable with seven new variables in the case of the the Novak and Tyson cell cycle model. And then the second reason for going over scripts versus functions again, is that in order to MATLAB ODE solver, you need to create a function. In this function will compute the derivatives of your variables, so in the MATLAB code we will provide you with to do your homework assignment, there's going to be two files. And one of them is going to be a function called dydt_novac and you know its a function because the first line of that function will look like this, function deriv equals dydt_novac T comma. Stavar, meaning that if you send dydt novak time point and a vector containing all the state variables in it, then dydt_novak is going to compute seven derivatives and send them back to whatever script called it. So it's worth reviewing scripts versus functions to understand how all this. Model is implemented with a script called novak.m and then a function called dydt_novak.m. Now that we've discussed scripts versus functions again, let's go back to our simple example, of the var the one variable ordinary differential equation that we previously solved using Euler.m. Knowing what we know now about how MATLAB's ode solvers work. This is how we would divide it, into a script and a function to solve it using ODE through three in MATLAB. Ode.m would set our would define our, our constants here. It would set the last time. We want to integrate this for up to, up to four minutes. Set the initial condition x0 equals c, and then this line, instead of a four loop we have this single line here that would perform the integration, call Ode23 dxdt is the name of the function the computes the derivative, integrate it at time zero up to time t last tlast in this case equal to 4 minutes. And then send it x0 as my as my initial condition. What about dxdt, what's that going to contain? Well, this is going to be a function where we compute the derivative. If we know t and we know x. And the derivative is just computed in this one line here a minus b times x. And, there's one more thing in this in this arrangement that you may have noticed which is that we like to define our, our parameters in the very beginning. Whatever our, our, whatever is controlling the behavior of our system, whatever our total concentrations are, or our rate constants etcetera. It's good to define them at the very beginning. That way, they're all in one place if you want to go change them. So, I prefer to structure things this way, where at the very beginning you say, okay, a, b, and c are the parameters that control the system. I'm going to define them at the beginning, but if they're defined in the script here A and B have to get used in the in the function over here that computes the derivative. So the way I've set this up, is to have, is to define A and B as global in the script. And then, in the function, I also defined A and B as global. And we talked about global variables versus local variables before, in the introduction to MATLAB lectures. Global variables can be used in both a script. And the function. But they have to be defined in both in order for MATLAB to understand them as, as variables that are truly global. And this general structure like we've shown here with ode.m and dxdt.m has been followed in in the two files novak.m and dydt_novak.m. To solve the Novak and Tyson cell cycle model. To summarize this lecture, one way to use MATLAB's built-in ODE solvers is to write a script that sets up the model and then a function that computes the ordinary differential equations. I want to emphasize there are other ways to solve this, and we're not going to go into all the different ways to do this. There are things called embedded functions. You can have the whole thing in one file. With one function that sets up the model and a second function that just computes the, the derivatives. But this is a structure that I like and this is the structure that you're going to be provided with. Where the script, Novac.m, sets up the parameters, tells you how long to integrate it for, and then actually tells MATLAB to solve it and does the post processing. The plotting the computations etcetera. And then there's a function and the purpose of the function is to compute the derivatives and send them back. So this is one way to use MATLAB's ODE solvers to solve a system of ODEs. And then the second thing that goes along with that is that. Parameters can be defined in the script and then used in the function and so, in that case, it can be very helpful to define these as global variables, and that's also the pattern that I've used in solving the Novak and Tyson cell cycle model. So this concludes our lectures on the cell cycle model. Particularly the the cell cycle models designed by Novak at all. And in your homework assignment you're going to use this, the MATLAB script and the MATLAB function that we provide you with. And you're going to reproduce some of the plots that we talked about in, in the lectures when we discussed what the Novak and Tyson cell cycle model actually was able to show. [BLANK_AUDIO]