Southern-Summer School on Mathematical Biology

User Tools

Site Tools

Trace:
2016:groups:g4:start

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision		Previous revision
2016:groups:g4:start [2016/01/03 19:33] – [Introduction] mendes2016:groups:g4:start [2024/01/09 18:45] (current) – external edit 127.0.0.1
Line 25: Line 25:
    *The authors make no conclusions about which plant species is being mimicked. How would the dynamics of this system change if there were more than one plant species serving as model?    *The authors make no conclusions about which plant species is being mimicked. How would the dynamics of this system change if there were more than one plant species serving as model?
    *Is it possible for orchid mantis to maintain higher population densities than their flower model? If so, what mechanisms would allow this surprising outcome?    *Is it possible for orchid mantis to maintain higher population densities than their flower model? If so, what mechanisms would allow this surprising outcome?
 +
 +===== Proposed Model =====
 +
 +Our model is in this {{:2016:groups:g4:model.pdf|}} pdf file. 
 +
 +===== For Testing the Model =====
 +
 +If you are in windows download you can just download unzip and run this file:
 +   *{{:2016:groups:g4:group4_model.rar|}}
 +
 +For testing our model during the presentation you can copy and paste this python source code and run it for an interactive model simulator.
 +If it doesn't work for you follow this instructions:
 +
 +If you already have python:
 +
 +   *Run in a console (this will download the packages needed and automatically install them):
 +<code bash>
 +pip install PyQt4 pyqtgraph
 +</code>
 +
 +Otherwise follow the instructions on the tutorials to install python. And then you should be able to run the line above. If the line above doesn't work, maybe you dont have pip.
 +
 +To install pip (python package manager):
 +  *In Ubuntu/Linux run:
 +<code bash>
 +  sudo apt-get install python-pip
 +</code>
 +  *In another operating system you should find how to install pip and then follow the instructions above.
 +
 +<file python model.py>
 +
 +#!/bin/python
 +
 +from scipy.integrate import odeint
 +import numpy as np
 +from PyQt4.QtCore import *
 +from PyQt4.QtGui import *
 +from pyqtgraph.widgets.PlotWidget import PlotWidget
 +from pyqtgraph.graphicsItems.InfiniteLine import InfiniteLine
 +from pyqtgraph.graphicsItems.ScatterPlotItem import ScatterPlotItem
 +from pyqtgraph import mkPen, mkBrush, setConfigOption
 +from pyqtgraph.parametertree import ParameterTree, Parameter
 +import pyqtgraph as pg
 +from pylab import plot, show, xlabel, ylabel, quiver
 +
 +t_init, t_step, t_end = 0, 0.1, 10000
 +t = t_init
 +t_values = np.arange(t_init, t_end, t_step)
 +
 +
 +class MyWidget(QWidget):
 +    def __init__(self):
 +        super(MyWidget, self).__init__()
 +        self.parameter_tree = ParameterTree()
 +        self.parameter_tree.setMaximumWidth(300)
 +
 +        self.k1, self.k2, self.k3, self.k4 = 0.4, 0.6, 0.1, 0.005
 +
 +        self.alpha, self.beta = 0.01, 0.1
 +        self.Kf, self.Kp = 0, 0
 +        self.d1, self.d2, self.d3 = 0.3, 0.2, 0.044
 +
 +        self.F, self.P, self.M = 100, 50, 2
 +
 +        self.parameter = Parameter(name='Parameters',
 +                                   children=[
 +                                       { 'name':'Initial Conditions', 'type':'group', 'children':
 +                                           [{ 'name':'F', 'type':'float', 'value':self.F },
 +                                            { 'name':'P', 'type':'float', 'value':self.P },
 +                                            { 'name':'M', 'type':'float', 'value':self.M }] },
 +                                       { 'name':'K', 'type':'group', 'children':
 +                                           [{ 'name':'k1', 'type':'float', 'value':self.k1 },
 +                                            { 'name':'k2', 'type':'float', 'value':self.k2 },
 +                                            { 'name':'k3', 'type':'float', 'value':self.k3 },
 +                                            { 'name':'k4', 'type':'float', 'value':self.k4 }] },
 +                                       { 'name':'D', 'type':'group', 'children':
 +                                           [{ 'name':'d1', 'type':'float', 'value':self.d1 },
 +                                            { 'name':'d2', 'type':'float', 'value':self.d2 },
 +                                            { 'name':'d3', 'type':'float', 'value':self.d3 }] },
 +                                       { 'name':'Other paramters', 'type':'group', 'children':
 +                                           [{ 'name':'alpha', 'type':'float', 'value':self.alpha },
 +                                            { 'name':'beta', 'type':'float', 'value':self.beta },
 +                                            { 'name':'Kf', 'type':'float', 'value':self.Kf },
 +                                            { 'name':'Kp', 'type':'float', 'value':self.Kp },
 +                                            ] }]
 +                                   )
 +        self.parameter_tree.setParameters(self.parameter)
 +        self.plot_widget = PlotWidget()
 +        self.plt = self.plot_widget.getPlotItem()
 +
 +        l = self.plt.addLegend(size=(100, 100), offset=(500, 30))
 +        self.plt.showGrid(True, True)
 +
 +        res = self.solve()
 +        self.f_plt = self.plt.plot(t_values, res[0], name='Flowers', antialias=True)
 +        self.f_plt.setPen(mkPen('#EE02FF', width=3))
 +        self.p_plt = self.plt.plot(t_values, res[1], name='Pollinators', antialias=True)
 +        self.p_plt.setPen(mkPen('#FFD800', width=3))
 +        self.m_plt = self.plt.plot(t_values, res[2], name='Mantis', antialias=True)
 +        self.m_plt.setPen(mkPen('#FF0000', width=3))
 +        print('Flowers {0}\nPollinators {1}\nMantis {2}'.format(res[0][-1], res[1][-1], res[2][-1]))
 +
 +        self.fp_plt = PlotWidget()
 +        fp_plt = self.fp_plt.getPlotItem()
 +        fp_plt.showGrid(True, True)
 +        l = fp_plt.addLegend(size=(100, 100), offset=(500, 30))
 +
 +        self.fp_f_plt = fp_plt.plot([], [], name='Flowers', antialias=True)  # , symbol='+')
 +        self.fp_f_plt.setPen(mkPen('#EE02FF', width=3))
 +        self.fp_p_plt = fp_plt.plot([], [], name='Pollinators', antialias=True)  # , symbol='s')
 +        self.fp_p_plt.setPen(mkPen('#FFD800', width=3))
 +        self.fp_m_plt = fp_plt.plot([], [], name='Mantis', antialias=True)  # , symbol='t')
 +        self.fp_m_plt.setPen(mkPen('#FF0000', width=3))
 +
 +        wid = QWidget()
 +        l2 = QVBoxLayout()
 +        self.min_value = QDoubleSpinBox()
 +        self.max_value = QDoubleSpinBox()
 +        self.parameter_name = QLineEdit('d3')
 +        run_button = QPushButton('Run')
 +        l3 = QHBoxLayout()
 +        l3.addWidget(QLabel('Parameter: '))
 +        l3.addWidget(self.parameter_name)
 +        l3.addWidget(QLabel('Min value: '))
 +        l3.addWidget(self.min_value)
 +        l3.addWidget(QLabel('Max value: '))
 +        l3.addWidget(self.max_value)
 +        l3.addWidget(run_button)
 +        l2.addLayout(l3)
 +        l2.addWidget(self.fp_plt)
 +        wid.setLayout(l2)
 +
 +        self.fpa_plt = PlotWidget()
 +        fpa_plt = self.fpa_plt.getPlotItem()
 +        fpa_plt.showGrid(True, True)
 +        yaxis = fpa_plt.getAxis('left')
 +        yaxis.setTicks([[(1, '0 0 0'), (2, 'CO'), (3, 'NM -'), (4, ('NM +'))]])
 +        line1 = InfiniteLine(QPointF(0, 1), angle=0, pen=mkPen('#f00', width=2))
 +        line2 = InfiniteLine(QPointF(0, 2), angle=0, pen=mkPen('#f00', width=2))
 +        line3 = InfiniteLine(QPointF(0, 3), angle=0, pen=mkPen('#f00', width=2))
 +        line4 = InfiniteLine(QPointF(0, 4), angle=0, pen=mkPen('#f00', width=2))
 +        fpa_plt.addItem(line1)
 +        fpa_plt.addItem(line2)
 +        fpa_plt.addItem(line3)
 +        fpa_plt.addItem(line4)
 +        self.sc1 = ScatterPlotItem(antialias=True, brush=mkBrush('#00f'), size=20)
 +        self.sc2 = ScatterPlotItem(antialias=True, brush=mkBrush('#ff0'), size=20)
 +        self.index = 0
 +        fpa_plt.addItem(self.sc1)
 +        fpa_plt.addItem(self.sc2)
 +        wid2 = QWidget()
 +        l4 = QVBoxLayout()
 +        l5 = QHBoxLayout()
 +        self.min_fxa = QDoubleSpinBox()
 +        self.max_fxa = QDoubleSpinBox()
 +        self.param_fxa = QLineEdit('d3')
 +        run_button_2 = QPushButton('Run')
 +        l5.addWidget(QLabel('Parameter: '))
 +        l5.addWidget(self.param_fxa)
 +        l5.addWidget(QLabel('Min value: '))
 +        l5.addWidget(self.min_fxa)
 +        l5.addWidget(QLabel('Max value: '))
 +        l5.addWidget(self.max_fxa)
 +        l5.addWidget(run_button_2)
 +        l4.addLayout(l5)
 +        l4.addWidget(self.fpa_plt)
 +        wid2.setLayout(l4)
 +
 +        lay = QHBoxLayout()
 +        self.tab = QTabWidget()
 +
 +        l6 = QVBoxLayout()
 +        l6.addWidget(self.parameter_tree)
 +        phase_button = QPushButton('Plot phase diagram')
 +        l6.addWidget(phase_button)
 +        lay.addLayout(l6)
 +        lay.addWidget(self.tab)
 +        self.tab.addTab(self.plot_widget, 'Simulation')
 +        self.tab.addTab(wid, 'Fixed Points')
 +        self.tab.addTab(wid2, 'Fixed Points Analysis')
 +
 +        self.setLayout(lay)
 +        self.parameter.sigTreeStateChanged.connect(self.upd)
 +        self.update_fp()
 +
 +        run_button.clicked.connect(self.run_a_lot)
 +        run_button_2.clicked.connect(self.run_a_lot_of_points)
 +        phase_button.clicked.connect(self.plot_phase)
 +
 +    def run_a_lot(self):
 +        interval = self.min_value.value(), self.max_value.value()
 +        values = np.linspace(interval[0], interval[1], 200)
 +        parameter_name = self.parameter_name.text()
 +        old_val = eval(str('self.' + parameter_name))
 +        f, p, m = [], [], []
 +        self.fp_plt.getPlotItem().setLabels(bottom=str(parameter_name), left='Fixed Points')
 +        for i in values:
 +            exec (str('self.' + parameter_name + '=' + str(i)))
 +            res = self.solve()
 +            f.append(res[0][-1])
 +            p.append(res[1][-1])
 +            m.append(res[2][-1])
 +            l = len(f)
 +            print('Finished simulation for {0} = {1}'.format(parameter_name, eval(str('self.' + parameter_name))))
 +
 +        self.fp_f_plt.setData(values, f)
 +        self.fp_p_plt.setData(values, p)
 +        self.fp_m_plt.setData(values, m)
 +        self.fp_plt.update()
 +        exec (str('self.' + parameter_name + '=' + str(old_val)))
 +
 +    def run_a_lot_of_points(self):
 +        interval = self.min_fxa.value(), self.max_fxa.value()
 +        values = np.linspace(interval[0], interval[1], 100)
 +        parameter_name = self.param_fxa.text()
 +        old_val = eval(str('self.' + parameter_name))
 +        s1x, s1y, s2x, s2y = [], [], [], []
 +        self.sc1.clear()
 +        self.sc2.clear()
 +        self.index = 0
 +
 +        for i in values:
 +            exec (str('self.' + parameter_name + '=' + str(i)))
 +            self.update_fp()
 +            res = self.solve()
 +
 +            print('Finished simulation for {0} = {1}'.format(parameter_name, eval(str('self.' + parameter_name))))
 +
 +            f, p, m = res[0][-1], res[1][-1], res[2][-1]
 +            eps = 5
 +            if np.abs(f - self.fp0[0]) < eps and np.abs(p - self.fp0[1]) < eps and np.abs(m - self.fp0[2]) < eps:
 +                val = 1
 +            elif np.abs(f - self.fp1[0]) < eps and np.abs(p - self.fp1[1]) < eps and np.abs(m - self.fp1[2]) < eps:
 +                val = 2
 +            elif np.abs(f - self.fp2[0]) < eps and np.abs(p - self.fp2[1]) < eps and np.abs(m - self.fp2[2]) < eps:
 +                val = 3
 +            elif np.abs(f - self.fp3[0]) < eps and np.abs(p - self.fp3[1]) < eps and np.abs(m - self.fp3[2]) < eps:
 +                val = 4
 +            else:
 +                val = 5
 +            if ((self.k4 / self.d3) - (self.k1 / self.d1)) < self.beta \
 +                    and (self.k4 / self.d3) > self.beta \
 +                    and (self.alpha * self.k1 * self.d2) / (self.d1 * self.k2) \
 +                            < ((self.k4 / self.d3) - self.beta) * \
 +                                    (self.beta + (self.k1 / self.d1) - (self.k4 / self.d3)):
 +                s2x.append(self.index)
 +                s2y.append(val)
 +            else:
 +                s1x.append(self.index)
 +                s1y.append(val)
 +            show = (4 * self.alpha * self.d1 * self.d2) / (self.k1 * self.k2)
 +            self.fpa_plt.setXRange(0, self.index, 0.5, True)
 +            self.index += 1
 +
 +        self.sc1.setData(s1x, s1y)
 +        self.sc2.setData(s2x, s2y)
 +
 +    def upd(self, *args):
 +        exec ('self.' + args[1][0][0].name() + '=' + str(args[1][0][2]))
 +        self.update_fp()
 +        res = self.solve()
 +        print('Flowers {0}\nPollinators {1}\nMantis {2}'.format(res[0][-1], res[1][-1], res[2][-1]))
 +
 +        if self.tab.currentIndex() == 0:
 +            self.f_plt.setData(t_values, res[0])
 +            self.p_plt.setData(t_values, res[1])
 +            self.m_plt.setData(t_values, res[2])
 +            self.plot_widget.update()
 +
 +    def update_fp(self):
 +        P0 = F0 = M0 = 0
 +
 +        P1 = self.d3 / (self.k4 - self.d3 * self.beta)
 +        F1 = (self.k1 * P1 - self.d1) / (self.alpha * P1 * self.d1)
 +        M1 = ((self.k2 * F1) / (1 + self.alpha * P1 * F1) - self.d2) * ((1 + self.beta * P1) / self.k3)
 +
 +        M2 = 0
 +        P2 = (self.k1 * self.k2
 +              - np.sqrt(self.k1 * self.k2 * (-4 * self.alpha * self.d1 * self.d2 + self.k1 * self.k2))) \
 +             / (2 * self.alpha * self.d2 * self.k1)
 +        F2 = (self.k1 * self.k2
 +              - np.sqrt(self.k1 * self.k2 * (-4 * self.alpha * self.d1 * self.d2 + self.k1 * self.k2))) \
 +             / (2 * self.alpha * self.d1 * self.k2)
 +        M3 = 0
 +        P3 = (self.k1 * self.k2
 +              + np.sqrt(self.k1 * self.k2 * (-4 * self.alpha * self.d1 * self.d2 + self.k1 * self.k2))) \
 +             / (2 * self.alpha * self.d2 * self.k1)
 +        F3 = (self.k1 * self.k2
 +              + np.sqrt(self.k1 * self.k2 * (-4 * self.alpha * self.d1 * self.d2 + self.k1 * self.k2))) \
 +             / (2 * self.alpha * self.d1 * self.k2)
 +        self.fp0, self.fp1, self.fp2, self.fp3 = (F0, P0, M0), (F1, P1, M1), (F2, P2, M2), (F3, P3, M3)
 +        print('F0 = {0}\nP0 = {1}\nM0 = {2}\n'
 +              .format(F0, P0, M0))
 +        print('F1 = {0}\nP1 = {1}\nM1 = {2}\n'
 +              .format(F1, P1, M1))
 +        print('F2 = {0}\nP2 = {1}\nM2 = {2}\n'
 +              .format(F2, P2, M2))
 +        print('F3 = {0}\nP3 = {1}\nM3 = {2}\n'
 +              .format(F3, P3, M3))
 +        print()
 +
 +    def solve(self):
 +        res = odeint(self.func(), (self.F, self.P, self.M), t_values).T
 +        return res
 +
 +    def func(self):
 +        flower = lambda t, f, p, m:(self.k1 * p * f) / (1 + self.alpha * f * p) \
 +                                   - self.d1 * f \
 +                                   - self.Kf * f * f
 +
 +        polli = lambda t, f, p, m:(self.k2 * p * f) / (1 + self.alpha * f * p) \
 +                                  - (self.k3 * m * p) / (1 + self.beta * p) \
 +                                  - self.d2 * p \
 +                                  - self.Kp * p * p
 +
 +        mantis = lambda t, f, p, m:((self.k4 * m * p)) / (1 + self.beta * p) \
 +                                   - self.d3 * m
 +
 +        return lambda y, t:(flower(t, *y), polli(t, *y), mantis(t, *y))
 +
 +    def plot_phase(self):
 +        # Equations without the mantis (simpler model for the phase space)
 +        flower = lambda t, f, p:(self.k1 * p * f) / (1 + self.alpha * f * p) \
 +                                - self.d1 * f
 +
 +        polli = lambda t, f, p:(self.k2 * p * f) / (1 + self.alpha * f * p) \
 +                               - self.d2 * p
 +
 +        f = lambda y, t:np.array([flower(t, *y), polli(t, *y)])
 +
 +        # res = odeint(f, (self.F, self.P), t_values).T
 +        # plot(res[0], res[1])
 +        res = odeint(f, (20, 20), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (150, 20), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (150, 350), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (200, 20), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (150, 20), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (150, 0), t_values).T
 +        plot(res[0], res[1])
 +
 +        xlabel('Flowers')
 +        ylabel('Pollinators')
 +        res = odeint(f, (0, 50), t_values).T
 +        plot(res[0], res[1])
 +        R, C = np.meshgrid(np.arange(-30, 200, 30), np.arange(-30, 330, 10))
 +        dy = f(np.array([R, C]), 0)
 +        plot([self.fp0[0]], [self.fp0[1]], 'o')
 +        plot([self.fp2[0]], [self.fp2[1]], 'o')
 +        plot([self.fp3[0]], [self.fp3[1]], 'o')
 +        quiver(R, C, dy[0, :], dy[1, :], scale_units='xy', angles='xy')
 +        show()
 +
 +        xlabel('Flowers')
 +        ylabel('Pollinators')
 +        res = odeint(f, (0.5, 0.3), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (0.3, 0.8), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (.1, 1), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (0.3, 0.7), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (0.2, 0.75), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (0.2, 0.7), t_values).T
 +        plot(res[0], res[1])
 +        res = odeint(f, (0.334170, 0.75188), t_values).T
 +        plot(res[0], res[1])
 +
 +        plot([self.fp0[0]], [self.fp0[1]], 'o')
 +        plot([self.fp2[0]], [self.fp2[1]], 'o')
 +        R, C = np.meshgrid(np.arange(-1, 2, .2), np.arange(-1, 2, .2))
 +        dy = f(np.array([R, C]), 0)
 +        quiver(R, C, dy[0, :], dy[1, :], scale_units='xy', angles='xy')
 +        show()
 +
 +
 +app = QApplication([])
 +
 +wid = MyWidget()
 +wid.show()
 +
 +app.exec_()
 +
 +
 +</file>
  
 ===== References ===== ===== References =====
2016/groups/g4/start.1451849627.txt.gz · Last modified: 2024/01/09 18:45 (external edit)