生物信息学 第八章 数学模型 毛理凯

本课目录 概述 差分方程 微分方程 应用 E-Cell

一、概述

数学模型的例子(米氏方程) 酶促反应机制 根据稳态/定态(steady state)假设 和反应动力学 推导出米氏方程

为什么要使用数学模型？ 通常利用数学模型来作为所关心的系统工作原理的假设 通过模拟(simulation)的结果可以证明假设是否正确  理解生命现象的机制 正确的模型可以进一步预测生命系统的其他未知特性  预言试验结果,指导实验设计,减少实验成本 善于在短时间内完成复杂的实验,甚至某些当前实验条件尚无法达到的

定义、构成元素 数学模型(mathematical model)是用数学语言来描述一个系统的抽象模型 这个数学语言通常是包含一些方程 例如一个群体增长模型 这个数学语言通常是包含一些方程 这些方程(equation)用来建立一些变量之间的关系 这些变量(variable)通常代表了系统的某些属性(property) 如某群体的大小

构成元素关系 系统 关系/规律 属性 数学模型 方程 变量

参数 模型还包括参数(parameter) 参数通常是常数,用于描述系统的某个相对不变的属性 参数在模型中相对于变量为从属地位 如某群体的生殖率(以群体大小为变量) 参数在模型中相对于变量为从属地位 一个属性是变量还是参数没有明显界限,由具体问题的性质决定 如果以生殖率为研究对象(变量),那么生殖率就不是参数,而是变量

数学模型的分类(1) 静态的(static)和动态的(dynamic) 区别在于是否考虑时间 动态模型常由差分方程或微分方程来表示 确定性的(deterministic)和随机性的(stochastic) 看是否唯一参数决定唯一结果 注意: 确定性模型可能产生貌似随机的结果,如混沌(chaos)

数学模型的分类(2) (时间)离散的(discrete)和连续的(continuous) 如差分方程(离散)和微分方程(连续) 线性(linear)和非线性的(nonlinear) y=ax+b (线性) y=ax2+bx+c (非线性) 对于方程组来说,只有全部方程都是线性的,该模型才是线性模型

数学模型的分类(3) 集总/中(lumped)参数和分布(distributed)参数模型 看参数是(集总)否(分布)均一分布 分布参数模型常用偏微分方程表示

一个离散模型的具体例子 生命游戏(life game) 演示… 属于细胞自动机(cellular automaton)的一种 给定某初始条件和繁衍条件 根据这些条件,观察群体的演化 定态,周期解,混沌… 演示…

二、差分方程 (difference equation)

例: 逻辑斯蒂映射(logistic map) 方程 Xn+1=rXn(1-Xn) Xn是变量,范围[0,1],代表某群体中第n代的个体数(已归一化) r是参数,表示增长率 如果知道前一项Xn,我们就可以推出后一项Xn+1 所以差分方程也叫递归(recursion)

解差分方程 要解这个差分方程,或者说进行模拟(run a simulation), 需要知道参数值(parameter values)、(变量)初值(initial values) 令 r = 1.0 X0 = 0.5 这样可以通过迭代(iteration)来求解差分方程

不同参数的效果(1) 周期一 周期一 周期二

不同参数的效果(2) 混沌(Chaos) 周期四…

迭代 对于本例(参数r=1.0) 用Excel操作、三维演示… X0=0.5 X1=0.25 X2=0.1875 X3=0.152344

换个方式演示迭代过程 用笔和尺

混沌的初值敏感性(sensitivity to initial conditions)

分岔图(bifurcation diagram) 就是横轴为参数、纵轴为变量的图,显示整个系统随参数的变化

丰富多彩的分岔图 – 前分岔、后分岔 后分岔(r<0) 前分岔(r>0)

前分岔局部放大 程序、动画演示… 丰富多彩的分岔图 – 自相似

丰富多彩的分岔图 – 三维 前后分岔、r为复数

三、微分方程 (differential equation)

(微分基础)微分/导数就是速度 从导数的定义开始 导数表示在x的某一点的切线的斜率,也就是变化率 变化率就是速度 Δx0

两种主要的微分方程 常微分方程(ordinary differential equation) u是x的函数(都是变量) 该方程的解为u(x)=c c为任意常数

两种主要的微分方程 偏微分方程(partial differential equation) u是x,y的函数 该方程暗示u独立于x 所以该方程的解为u(x,y)=f(y) f是y的任意函数

(生态学例子)群体增长模型(1) 方程 x是变量,代表某群体的个体数,即该群体大小,对时间t求导 m是参数,表示增长率 求导表示上变量对下变量变化的速度,所以这里的求导代表某群体大小的变化速度

群体增长模型(2) 这样上述方程就表示某群体的增长速度跟现有的群体大小成正比(这意味着指数增长！) 该方程其实就是著名的马尔萨斯人口方程,m是马尔萨斯参数(Malthusian parameter)

群体增长模型(3) 该方程的(解析)解(analytic solution)是 m=1, x0=1

(混沌例子)Lorenz奇怪吸引子 微分方程也可以产生混沌！而且更漂亮！ 例如Lorenz奇怪吸引子(strange attrator)

微分方程的数值解 这个方程不易得出解析解 需转化成差分方程并借助计算机求得数值解(numerical solution) 欧拉折线法(Euler method) dy/dx=f(x,y) (yn+1-yn)/h=f(xn,yn) yn+1=yn+h f(xn,yn) 转化成了差分方程 用Excel也可以解（演示…）!

用软件Euler解Lorenz方程 Euler 免费Matlab克隆 几乎可做常见的任何数学操作,甚至可以符号运算! ~2M! Homepage 演示…

(例子)Logistic映射的微分形式(单物种增长) [差分] Xn+1 = rXn(1-Xn) [微分] dX/dt = rX(1-X/K) X : 群体大小(变量) t : 时间 r : 增值率(参数) K : 群体大小极限(参数) 该方程比Malthus模型更接近现实,考虑了资源限制

单物种增长模型的解 变量初值 参数值 (变化) Euler演示解的演化、解受参数的影响 X0=1 r=1 (1…10) 参数值 (变化) r=1 (1…10) K=10000 (1000…10000) Euler演示解的演化、解受参数的影响 不再指数增长(资源限制K起作用了!) 还不如差分方程的解丰富 只有定态解(steady states, fixed points, equilibria)

定态解及其稳定性 令方程右边rX(1-X/K)=0,即可得定态解 求这些定态解的稳定性(stability) X1=0, X2=K 对方程右边求导 [rX(1-X/K)]’=r-2rX/K 将定态解代入 r-2rX1/K=r >0  X1不稳定  不可见 r-2rX2/K=-r <0  X2稳定  可见

丰富多彩的混沌 分形学

Dynamics Solver 免费数学运算、作图软件 特别擅长于非线性动力学、混沌、分形 ~7M 软件自带混沌示例 bifurcation.ds (Logistic) circle.ds, Crutchfield.ds, tent.ds (不同的分岔图) Henon4.ds (初值敏感) Henon1.ds, baker.ds, Lozi.ds, Julia.ds, Mandelbrot.ds, Newton.ds, von Koch.ds, snowflake.ds, tree.ds (自相似,丰富的细节,分形)

四、应用

应用广泛(仅生命科学方面的部分列举) 生态学 酶动力学(生化) 神经系统 细胞代谢系统 信号转导系统 传染病 群体遗传学 捕食-被捕食模型 米氏方程 神经系统 细胞代谢系统 信号转导系统 传染病 群体遗传学

群体遗传学 – 模拟突变 研究对象/假设 代与代不重叠,随机交配,群体无限大 1个位点,2个等位基因(A1, A2),pn和qn=1-pn是它们在第n代时的基因频率 A1变异为A2的突变率是u,A2变异为A1的突变率是v 设一代中一个等位基因只能变异一次 u pn A1 A2 v qn

突变方程及其解 这样下一代的A1为 这个差分方程的解为 通常u, v很小(10-6或10-5的量级) pn+1=(1-u)pn+v(1-pn) 这个差分方程的解为 这里p0是开始时(第0代)A1的频率 通常u, v很小(10-6或10-5的量级) 当n∞, pnv/(u+v), qnu/(u+v) 达到平衡(实际很难达到)

predator-prey模型 Malthus和Logistic模型是单物种模型 predator-prey模型是一类双物种模型

Lotka-Volterra模型 Lotka-Volterra模型是最早的predator-prey模型 [美]生物物理学家Alfred Lotka (1925) [意]数学家Vito Volterra (1926) 基于一阶非线性常微分方程 被捕食者 捕食者 Euler数值解演示…

定态解 求定态解 -αx-βxy=0 -δxy-γy=0 得 {x=y=0} (定态解1) {x=α/β, y=γ/δ} (定态解2)

定态解的稳定性 用偏导数线性化方程右端 得Jacobian matrix 该矩阵的本征值(eigenvalue)是 λ1=α>0, λ2 =-γ<0 (定态解1) 该定态解是鞍点(saddle point,不稳定) λ1=i√αγ>0, λ2 =-i√αγ<0 (定态解2) 该定态解是焦点(focus, 稳定周期)

五、E-Cell

E-Cell简介 功能: 在分子水平上全细胞模拟 免费/Gnu General Public License (GPL)、开源 跨平台(Linux, Windows, Mac) 程序架构: 前端/界面python,核心C++ 支持各类数学模型,参数估计,分析,便于自动化 E-Cell 3D (for Mac) 演示…

考试 不定项选择题 30 (15) 是非题 名词解释题 20 (5) 综合分析题 20 (2)

完