16 Pre-Class Assignment: Linear Dynamical Systems

1. Linear Dynamical Systems

A linear dynamical system is a simple model of how a system changes with time. These systems can be represented by the following “dynamics” or “update equation”:

\[x_{(t+1)} = A_tx_t\]

Where \(t\) is an integer representing th progress of time and \(A_t\) are an \(n \times n\) matrix called the dynamics matrices. Often the above matrix does not change with \(t\). In this case the system is called “time-invariant”.

We have seen a few “time-invarient” examples in class.

DO THIS: Review Chapter 9 in the Boyd and Vandenberghe text and become familiar with the contents and the basic terminology.


2. Markov Models

This is not the first time we have used Dynamical Linear Systems.

DO THIS: Review markov models in 10–Eigenproblems_pre-class-assignment.ipynb. See how this is a special type of linear dynamical systems that work with state probabilities.

Example

The dynamics of infection and the spread of an epidemic can be modeled as a linear dynamical system.

We count the fraction of the population in the following four groups:

  • Susceptible: the individuals can be infected next day

  • Infected: the infected individuals

  • Recovered (and immune): recovered individuals from the disease and will not be infected again

  • Decreased: the individuals died from the disease

We denote the fractions of these four groups in \(x(t)\). For example \(x(t)=(0.8,0.1,0.05,0.05)\) means that at day \(t\), 80% of the population are susceptible, 10% are infected, 5% are recovered and immuned, and 5% died.

We choose a simple model here. After each day,

  • 5% of the susceptible individuals will get infected

  • 3% of infected inviduals will die

  • 10% of infected inviduals will recover and immuned to the disease

  • 4% of infected inviduals will recover but not immuned to the disease

  • 83% of the infected inviduals will remain

Do this: Write the dynamics matrix for the above markov linear dynamical system. Come to class ready to discuss the matrix. (hint the columns of the matrix should add to 1).

# Put your matrix here

Do this: Review how we solved for the long term steady state of the markov system. See if you can find these probabilities for your dyamics matrix.

# Put your matrix here

3. Ordinary Differential Equations

Ordinary Differential Equations (ODEs) are yet another for of linear dynamical systems and are a scientific model used in a wide range of problems of the basic form:

$\(\dot{x} = A x\)$

These are equations such that the is the instantaneous rate of change in \(x\) (i.e. \(\dot{x}\) is the derivative of \(x\)) is dependent on \(x\). Many systems can be modeled with these types of equations.

Here is a quick video that introduces the concepts of Differential Equations. The following is a good review of general ODEs.

from IPython.display import YouTubeVideo
YouTubeVideo("8QeCQn7uxnE",width=640,height=360, cc_load_policy=True)

Now consider an ODE as a system of linear equations:

\[\dot{x_t} = A x_t\]

Based on the current \(x\) vector at time \(t\) and the matrix \(A\), we can calculate the derivative at \(\dot{x}\) at time \(t\). Once we know the derivative, we can increment the time to by some small amount \(dt\) and calculate a new value of \(x\) as follows:

\[x_{t+1} = x_t + \dot{x_t}dt\]

Then we can do the exact sequence of calculations again for \(t+2\). The following function has the transition matrix (\(A\)), the starting state vector (\(x_0\)) and a number of time steps (\(N\)) and uses the above equations to calculate each state and return all of the \(x\) statues:

The following code generates a trajectory of points starting from x_0, applying the matrix \(A\) to get \(x_1\) and then applying \(A\) again to see how the system progresses from the start state.

%matplotlib inline
import matplotlib.pylab as plt
import numpy as np
import sympy as sym
sym.init_printing()
---------------------------------------------------------------------------
ModuleNotFoundError                       Traceback (most recent call last)
<ipython-input-4-97c04e533e83> in <module>
----> 1 get_ipython().run_line_magic('matplotlib', 'inline')
      2 import matplotlib.pylab as plt
      3 import numpy as np
      4 import sympy as sym
      5 sym.init_printing()

~/REPOS/MTH314_Textbook/MakeTextbook/envs/lib/python3.9/site-packages/IPython/core/interactiveshell.py in run_line_magic(self, magic_name, line, _stack_depth)
   2342                 kwargs['local_ns'] = self.get_local_scope(stack_depth)
   2343             with self.builtin_trap:
-> 2344                 result = fn(*args, **kwargs)
   2345             return result
   2346 

~/REPOS/MTH314_Textbook/MakeTextbook/envs/lib/python3.9/site-packages/decorator.py in fun(*args, **kw)
    230             if not kwsyntax:
    231                 args, kw = fix(args, kw, sig)
--> 232             return caller(func, *(extras + args), **kw)
    233     fun.__name__ = func.__name__
    234     fun.__doc__ = func.__doc__

~/REPOS/MTH314_Textbook/MakeTextbook/envs/lib/python3.9/site-packages/IPython/core/magic.py in <lambda>(f, *a, **k)
    185     # but it's overkill for just that one bit of state.
    186     def magic_deco(arg):
--> 187         call = lambda f, *a, **k: f(*a, **k)
    188 
    189         if callable(arg):

~/REPOS/MTH314_Textbook/MakeTextbook/envs/lib/python3.9/site-packages/IPython/core/magics/pylab.py in matplotlib(self, line)
     97             print("Available matplotlib backends: %s" % backends_list)
     98         else:
---> 99             gui, backend = self.shell.enable_matplotlib(args.gui.lower() if isinstance(args.gui, str) else args.gui)
    100             self._show_matplotlib_backend(args.gui, backend)
    101 

~/REPOS/MTH314_Textbook/MakeTextbook/envs/lib/python3.9/site-packages/IPython/core/interactiveshell.py in enable_matplotlib(self, gui)
   3511         """
   3512         from IPython.core import pylabtools as pt
-> 3513         gui, backend = pt.find_gui_and_backend(gui, self.pylab_gui_select)
   3514 
   3515         if gui != 'inline':

~/REPOS/MTH314_Textbook/MakeTextbook/envs/lib/python3.9/site-packages/IPython/core/pylabtools.py in find_gui_and_backend(gui, gui_select)
    278     """
    279 
--> 280     import matplotlib
    281 
    282     if gui and gui != 'auto':

ModuleNotFoundError: No module named 'matplotlib'
def traj(A, x, n):
    dt = 0.01
    x_all = np.matrix(np.zeros((len(x),n)))   # Store all points on the trajectory
    for i in range(n):  
        x_dot = A*x         # First we transform x into the derrivative
        x = x + x_dot*dt    # Then we estimate x based on the previous value and a small increment of time.
        x_all[:,i] = x[:,0] 
    return x_all

For example the following code uses the matrix \(A= \begin{bmatrix}1 & 1 \\ 1 & -2\end{bmatrix}\) and the starting point (0,0) over 50 timesteps to get a graph:

A = np.matrix([[1,1],[1,-2]])
x0 = np.matrix([[1],[1]])

x_all = traj(A, x0, 50)
plt.scatter(np.asarray(x_all[0,:]),np.asarray(x_all[1,:]))

plt.scatter(list(x0[0,:]),list(x0[1,:])) #Plot the start point as a refernce

Do this: Let $\(A= \begin{bmatrix}2 & 3 \\ 4 & -2\end{bmatrix}\)$

Write a loop over the points \((1.5,1)\), \((-1.5,-1)\), \((-1,2)\) and plot the results of the traj function:

A = np.matrix([[2,3],[4,-2]])
x0 = np.matrix([[1.5, -1.5, -1, 1, 2],[1, -1, 2, -2, -2]])
# Put your code here

Do this: Let $\(A= \begin{bmatrix}6 & -1 \\ 1 & 4\end{bmatrix}\)$

Write a loop over the points \((1.5,1)\), \((-1.5,-1)\), \((-1,2)\), \((1,-2)\) and \((2,-2)\) and plot the results of the traj function:

# Put your code here

Do this: Let $\(A= \begin{bmatrix}5 & 2 \\ -4 & 1\end{bmatrix}\)$

Write a loop over the points \((1.5,1)\), \((-1.5,-1)\), \((-1,2)\), \((1,-2)\) and \((2,-2)\) and plot the results of the traj function:

# Put your code here

4. Assignment wrap up

Assignment-Specific QUESTION: Where you able to get the ODE code working in the above example. If not, where did you get stuck?

Put your answer to the above question here

QUESTION: Summarize what you did in this assignment.

Put your answer to the above question here

QUESTION: What questions do you have, if any, about any of the topics discussed in this assignment after working through the jupyter notebook?

Put your answer to the above question here

QUESTION: How well do you feel this assignment helped you to achieve a better understanding of the above mentioned topic(s)?

Put your answer to the above question here

QUESTION: What was the most challenging part of this assignment for you?

Put your answer to the above question here

QUESTION: What was the least challenging part of this assignment for you?

Put your answer to the above question here

QUESTION: What kind of additional questions or support, if any, do you feel you need to have a better understanding of the content in this assignment?

Put your answer to the above question here

QUESTION: Do you have any further questions or comments about this material, or anything else that’s going on in class?

Put your answer to the above question here

QUESTION: Approximately how long did this pre-class assignment take?

Put your answer to the above question here


Written by Dr. Dirk Colbry, Michigan State University Creative Commons License
This work is licensed under a Creative Commons Attribution-NonCommercial 4.0 International License.