Hello again. We are now at the third part of our lecture on Introduction to Dynamical Systems. In the last lecture in part two of the Dynamical Systems Lecture Series, we discussed how you would numerically solve systems of ODE's specifically we discussed using Euler's method. Now in this lecture we're going to talk about not just how you solve them numerically but how you analyze these ODE systems to get, to understand general properties of their system, and we're going to focus on this concept of stability of ODE systems. So first we'll discuss what we mean in general when we talk about stability. Then I'm going to show a couple of one-dimensional examples so you can understand how, how stability works in that context. And then we're going to introduce phase-plane techniques for analyzing two-dimensional systems. Although most of what we explained about these phase-plane techniques are going to be covered in more depth in the, in the fourth part of this series of lectures. And then I'm going to briefly introduce an example that we're going to deal with in in part four which is a mathematical model of glycolytic oscillations in yeast. The focus of this lecture as mentioned is going to be an analyzing stability of OED systems. So to introduce what we mean when we talk about stability, let's go back once again to this generic three component repressive network that we've seen before. Where you have three protein species, a, b, and c, that could either be dephosphorylated or they could be phosphorylated. And a, b, and c regulate either the dephosphorylation or the phosphorylation reactions of the other species. And as we have discussed, this scheme implies a set of differential equations, which we've gone through and looked at as follows over here. What we haven't discussed as much are the model parameters in the system. And when we look at this we can see several constants in here. That are that are model parameters that can be varied. And that can govern the behavior of the system. One is we have Kcat and km's for the phosphorylation reactions. And those are the ones that are over here for, labeled as kinases. And Kcats are the little k's in the numerator, and then the km's are the big k's for the phosphorylation reactions in the denominator. Similarly we have Kcats and Kms for the dephosphorylation reactions, those are the analogous parameters over here on these left, these leftmost terms. And then the other things we have, the other parameters we have our total amounts of A B and C that is A total B total and C total here. And as we have discussed, these equations are solved using standard numerical techniques. And I'll just give a couple examples of what you can see when you vary the parameters in this case. There are two two general categories of solutions to this system. With one set of parameters, we can see that A increases to a steady level, B goes to a different steady level that's much lower than the steady level of A, and then C goes to a level that's very, very low. But if you continue to run this for longer and longer periods of time, it would continue to stay at this very steady level. Once A gets to this level here, right around 0.9, it stays there, and it will stay there forever. But with a different parameter set, we see this sort of behavior. We can see the sustained oscillations able to go up and then go down and similarly B and C will also increase and decrease as a function of time and these oscillations will go, go on forever. So what we conclude from this in general is that the parameter values can greatly influence the system behavior. We have qualitatively different behavior over here on the left. Where A, B, and C go to steady levels. Compared to on the right, where A, B, and C continue to oscillate. How do we understand these different behaviors? Well, that's where the tools of dynamical systems come in. And that's how we can, we can understand and categorize the different types of behaviors that we observe. [SOUND] To illustrate how the tools of dynamical systems can be used to analyze stability, we will consider a one-dimensional example. This is a standard mathematical model of an isolated cardiac myocyte, where different ionic currents are responsible for changing the, the cells transmembrane potential, the details of these ionic currents are not important here. What is important here is that the differential equation describing voltage across the cell membrane. DV, dt is a negative of the ionic current. That means the sum of, of the current through all of these different channels, and pumps and transporters, divided by the capacitance. And we want to do a, a relatively simple experiment with this isolated cardiac myocyte model. We begin with the myocyte at rest, which is approximately minus 85 millivolts. Then we instantaneously change the voltage. When we instantaneously change the voltage we calculate the instantaneous ionic current, I ion. And if we, once we calculate the instantaneous ionic current, then we know what the instantaneous change is in the in the derivative. So, we can fit, we can vary voltage to whatever we want it to be instantaneously from rest. And then we want to calculate what's the resulting derivative going to be, what's the dV/dt going to be. And in that case we get a curve that looks like this. For voltages that are much lower than minus 85 millivolts, dV/dt is positive. For voltages between minus 85 millivolts and around plus 58 millivolts dV/dt is, is negative. And then for voltages above minus 58 millivolts dV/dt is positive again. So what happens when this crosses zero, right. There are two, there are two voltages here for which dV/dt instantaneously will be equal to zero. And if dV/dt is equal to zero, that means that if you change voltage and you put it at exactly that voltage right there where it's not changing, then it's going to stay there. And because these are voltages where the derivative of the voltage is zero. These are, are what are known as fixed points. So the fixed points in this case you can represent in this one dimensional example by plotting with derivative on the y axis, tying the variable on the x axis. And then seeing where the curve crosses zeri, right and so where they, where, where this is zeri is the same where derivative is zero. That means that if the voltage is at this level, and the derivative at that, in that case is is zero, then it's going to stay there. But then what we want to do is we want to say, well, what's going to happen to the voltage when we move it away from one of these fixed points. And this is where it gets more complicated and, and, and more interesting. Right, if we start at minus 85 and we go negative. This is a dV/dt is positive, that means that voltage is going to move to the right. Voltage is going to go up because the derivative is positive. So when we change to minus 85 millivolts we have a positive dV/dt. If we were to change something like minus 70, in this case the derivative is negative. And so therefore we can draw the arrow this way. And then if we were to go from minus 85 all the way to like minus 55, we would have a positive dV/dt, and so we would draw the arrow this way. So that's why it's actually useful to, to draw this kind of plot, where you have the derivative on the the y axis, and the variable on the x axis. Is because knowing if the derivative is positive or negative in particular ranges will tell you, well which way is the system going to evolve. Which way are things going to change when we when we move the variable, in this case, voltage to that particular location. So here we have an arrow going to the right, indicating that the volt, the derivative is positive. Here we have an arrow moving to the left, indicating that the derivative is negative. And here we have an arrow moving to the right again, indicating that the derivative is positive. Why is it useful to plot the the arrows this way to show how the system evolves when we change the when we change the variable, in this case the voltage. Well this is a way that we can tell whether our fix points are stable fix points or unstable fix points. What do we mean by that? Well, this fix point where the derivative crosses zero has an arrow pointing towards it. And an arrow pointing on the left, and an arrow pointing towards it on the right. That means that we start at minus 85 and we move away from it, it's going to go back to minus 85, because that's the way that it's pointed, right? Similarly, if we go from minus 85 to minus 75, the arrow is pointing negative, it's pointing to the left, it's pointing back toward minus 85. And so, the voltage as we go back is going to go back to minus 85. So, we can now delineate this fixed point as a stable fixed point because deviations away from that fixed point will always push the system back towards it. This fixed point at minus 58, in contrast, is an unstable fixed point. What happens if you move negative to minus 58, well the derivative in that case is going to be negative, and that means that it's voltage is going to continue moving negative. Similarly, if you're at minus 58 and you get nudged a little bit positive, the derivative is going to be positive and that means it's going to continue to move away from that. So, this fixed point has arrows pointing away from it. That's how we identify it as an unstable fixed point. And this fixed point has arrows pointing towards it. That's how we identify it as a stable fixed point. And we can perform this experiment numerically by changing voltage and then allowing the system to evolve, then integrating the equations numerically using Euler's Method. And we see exactly what we can see over here on the left in schematic format. We're starting at minus 85 millivolts. When we hyper-polarize the cell, when we move the voltage down to minus 95 millivolts, that's the black curve, what do we see? It comes back to minus 85. When we go up to minus 75, which is a Wrenn curve, we see, we instantaneously change it and then it goes back. Same thing with minus 65 which is the green curve. it goes up to that instantaneously and it goes back. And we can identify mi, minus 75 and minus 65. On this, on this plot, and see that yes, the system will go back to minus 85, when we're in that regime. But if we start at minus 85, and we move the cell all the way to minus 55, what do we see? The voltage continues to go up. And then it reaches a peak, and then it continues to evolve, and that's because there are multiple differential equations in this model. There's somewhere around the order of eight differential equations in this model. But we can see that small deviations from minus 85, the system will come back to minus 85 millivolts. But a large deviation, if you cross this unstable fixed point, then you're going to move away from that unstable fixed point and move away from minus 85 millivolts. And so by doing this numerical equation here, we can confirm what we see graphically over here on the left. That between mi, between minus 85 and minus 58, or something that's negative to minus 85 is going to move back to the resting state, to minus 85. So therefore that's a stable fixed point. Small deviations away from minus 58 are going to cause action potentials or return to the states and therefore this one is an unstable fixed point. What we're going to do next is we're going to learn to analyze these stable fixed points and unstable fixed points more rigorously. The example I just showed of of a cardiac action potential, we were plotting everything as the derivative of voltage, minus voltage but it is in fact a model that has eight differential equations in that case. And, what's happening with eight differential equations you can't, you cannot visualize everything graphically. So now what we are going to do is we're going to look at a two variable model where you can analyze things graphically. And that's going to help us develop the tools we need of dynamical systems, to be able to analyze these models. And the two dimensional example we're going to use describes the biological phenomenon of glycolytic oscillations in yeast. And here's an, an example from a, from a pretty famous paper showing what happens. There's an experimental procedure that you have to undertake in order to get to, in order to observe these oscillations, where you condition the cells, you starve them, and then when you add glucose, if you look at oxygen consumption in this case, what do you see? You see it goes up and down, up and down, up and down. And these oscillations will continue for a very long time. And this is a phenomenon that's been, that's been known to occur for quite a while. And in fact, many mathematical models of this process have been developed and have been published over the years. And, we're going to analyze one that was published in, in the year 2000 by Bier et al, and the citation for this for this model that we're going to consider is given right here. And this is the, the example system that we're going to use to introduce concepts of stability and to show you how to analyze stability of ODE systems. Now we want to ask what's the structure of the Bier et al model and it's summarized in this diagram here. The Bier et al model has two variables so the two variables that are computed that are evolved with respect to time are glucose inside the cell and ATP inside the cell. Now we'll go through each step in the process. Model simulates transport of glucose from the outside of the cell to the inside of the cell through this rate called, called VN. Then once glucose is inside the cell it can get converted into ATP and this step here represents glycolysis. in, in real life glycolysis repre, is occurs through the action of several different enzymes, there are many different steps in this process. But to keep things simple Bier et al represented lumped all of these steps into, into one, and used the single rate constant k1. The most important enzyme in this process, the rate limiting step, is an enzyme that some of you may have heard of called phosphofructokinase. So, this this rate k1, again it represents the action of several steps in the process, several different enzymes. But it can be generally considered primarily the action of phosphofructokinase activity within the cells. Once ATP is produced, it can get consumed. And again, this is another place where, where Bier et al made an approximation. There are several different ATPases that act, that act within cells, each one of which has, has different activity to keep things simple Bier et al lumped all the ATPases together. And then there's one more aspect of this model that's really important and that's this dashed arrow, right, I drew down here. Glycolysis is an interesting reaction in that it produces ATP from glucose, but it also requires ATP to be initiated. So when ATP goes up there are a couple of actions. One is that ATP goes up, then ATP gets consumed, but the other thing that happens is that this reaction here, this production of ATP occurs more quickly. And so we draw this as a regulatory arrow as ATP goes up. Then you get more conversion of glucose into ATP. So, we can look at the equations for the differential equations for ATP and the differential equations for glucose and we can understand the four terms here as the action as either the steps that increase glucose or decrease glucose. The steps that increase ATP or decrease ATP. We'll look at glucose, dG dt first this ODE for glucose, initially. Right, if you have more transport of glucose into the cell, then glucose concentration's going to go up. And then this term here represents conversion of glucose into AT, ATP through glycolysis, through the action of several different enzymes, but phosphofructokinese being the most important one. Now if we look at the things that can either produce or consume ATP, this positive term here represents conversion of glucose into ATP, and again it depends on both glucose and ATP. And then this term here represents consumption of ATP through the action of ATPases. And there are four four primary parameters in the Bier et al model and the default values of three of them VN, K one and KP are given here. And there is a fourth model, a fourth model parameter that we haven't really discussed yet, and that is the the KM for the action of the, of the ATPases. And let's look at what happens in the, in the Bier et al model when we change this fourth parameter, the Km for the ATPases and the kls constants for the ATPases. When Km is equal to 13, we see these sustained oscillations of both glucose and ATP. Glucose here is plotted in black. ATP is plotted in red. And, and you see that glucose goes up and goes down, goes up and goes down, and this continues indefinitely. Similarly, ATP exhibits oscillations, and they continue indefinitely. What happens if we set Km equal to 20 instead of 13. Well at the very beginning of this simulation, we see oscillations, but what do we notice about these oscillations? As time goes on the oscillations keep getting smaller, and smaller, and smaller. And this is what we would refer to as, as damped oscillations and with Km equal to 20. And then finally what we notice is the oscillations get so small you can't see them anymore. So eventually glucose settles in at a constant value and ATP settles in at a constant value. So just like what we saw with the three component repressive network, we can see qualitatively different behavior as we change one of the parameters. We see sustained oscillations over here, with Km equals 13, and we see damped oscillations, and eventual settling at steady levels of glucose and ATP over here. So what we want to do next is we want to ask, how can we understand the qualitatively different behavior that we see in these two cases. [SOUND] What we want to introduce now are what we call phase plane techniques for analyzing two dimensional systems of ODEs. And what we mean by phase plane is that we're not going to plot glucose and ATP versus time. Instead of plotting glucose versus time, and ATP versus time, we want to plot glucose versus ATP. So, glucose will be on the y axis in this case, and ATP will be on the x axis. So let's take the time courses that we simulated on the last slide, for Km equals 13 on the left, and for Km equals 20 on the right. And plot them in the phase plane where ATP is on the x axis and glucose is on the y axis. And what we see with ATP what we see in the phase plane with Km equals 13 is this loop. And this loop makes sense intuitively because we saw that glucose would oscillate, and ATP would oscillate and those would continue indefinitely. And when you have two variables that are, that are both oscillating, what you see when you plot one versus the other is something like this that looks like a loop. What about for case of Km equals 20 where you saw these damped oscillations, where you saw these small fluctuations in the beginning but eventually they settle down to a steady level. Well, what this trajectory looks like in the phase plan is a something that looks like this, it starts here and then it goes like this and it spirals. And it goes down to some level and it eventually, it keeps spiraling and then eventually it gets to some steady level where ATP stays constant and glucose stays constant forever. So the terminology we introduced to describe this is on the left. We say that glucose and ATP oscillate indef, indefinitely in something that we call a stable limit cycle. In contrast on the right, we see glucose and ATP that converge to a stable fixed point and in the next lecture we will discuss these, this terminology in, in more detail. Why do we think it's useful to plot these in the phase plane? When you get similar information when you plot glucose versus time and ATP versus time. Well, one of the reasons why the the plotting things in the phase plane is useful is that what you care about is which direction the system is moving. And the direction is determined by the derivative of ATP with respect to time and the derivative of glucose with respect to time. And so at any given location in your phase plane and any combination of ATP and glucose, you can calculate the root of ATP with respect to time and the derivative of glucose with respect to time. And this defines a vector in the phase plane, and this vector tells you how the system is moving, how ATP is changing with respect to time and how glucose is changing with respect to time. So we can understand this directionality of our of our two, of our two dimensional derivative vector, like this. This is our differential equation for glucose, this is our differential equation for ATP. And what we saw was with Km equals 13 we, we saw the stable limit cycle. One question becomes well, which way is it moving around this stable limit cycle? Is it going clockwise? Like I've drawn it here with the arrows, or is it, would it be going counterclockwise? Well, we can understand that by looking at these differential equations here. What happens if we had very large values of glucose and, and ATP? [SOUND] Glucose, the changing glucose with respect to time would be negative. Because Vn in this case is a constant. And glucose and ATP are multiplied together, in this case, on this term, that has a negative in front of it. So eventually, this term would become greater than this term Vn. And glucose would be, would would be decreasing with respect to time. We'd have a negative derivative. So that's how we know that this arrow for glucose when we are up here at the top right corner of this phase plane, should be pointed down. What about ATP? Well, when glucose and ATP both get large this term here, this first term keeps getting bigger and bigger, this positive term gets bigger and bigger. This term here, which has a negative in front of it, doesn't continue to get bigger. As ATP grows this term is eventually going to saturate. And so this term is going to go up for low levels of ATP but then it'll get to a point where it cannot continue to go up any more. So clearly this term is going to eventually overcome this negative term here. So the change in ATP with respect to time for large values of glucose and ATP is going to be positive. So, just by looking at these differential equations qualitatively we can draw this arrow up here for large levels of glucose and ATP, and we know that it's pointing to the right because the derivative of ATP is positive. And we know that it's pointing down because the derivative of glucose in this case is negative. And so now that we know that this arrow is pointing this way, we can then deduce how things are going through the rest of the system. And we could have made the same sort of argument if we had picked very low levels of glucose and ATP, for instance. In the next lecture, we'll we'll talk about how we can deduce when these arrows are going to switch directions. Why does it go from pointing down and to the right over here, to pointing down and to the left over here? So, we'll talk about that in the next lecture. But for now, we just want to introduce this concept of plotting, of looking at systems and how they involve in the phase plane. Which means you have one variable on the y axis, and another variable here on the x axis. To summarize this lecture, a fixed point of a dynamical system is a set of variables where all derivatives are equal to zero. A fixed point can be stable. Which means that small perturbations away from that fixed point will cause the system to evolve, and return back to that fixed point. Conversely, as we saw in the one dimensional example, we can also have unstable fixed points. Which means that if you have a small perturbation away from an unstable fixed point, it will continue to move away from that fixed point. Finally, with two variable systems it's helpful to plot one variable versus the other and that's what we refer to as analyzing the system in the phase plane.