3 生物反应过程系统辨识与状态估计v3.1 概述v3.2 系统辩识的基本方法v3.3 参数估计的基本方法v3.4 生物反应器的状态观测 3.1 概述v模型：系统的输入向量，输出向量和状态向量之间关系的数学表达式。v输入向量：施加于系统或仪表上的独立变量；v输出向量：系统或仪表送出的变量。v输入和输出指的是信号流，而非物料流或能量流。v状态向量：确定动态系统状态所需的一组最小数目的变量。例：生物反应器内发酵液体积变化v从物料流看，Fin是流入物料的流量，Fout为流出液流量。v就信息流而言，输入变量： Fout、Fin和时间t （影响发酵液体积V变化的原因），输出变量为： dV=（Fin - Fout）dt。v在规定Fin和Fout变化的情况下，只需知道初始体积V(t0)，体积变化过程就完全确定了。因此该过程的状态变量只有一个，即V(t)。图3.1生物反应器内发酵液体积变化v建立系统数学模型的方法：理论分析和实验两种。v理论分析方法：根据系统内部已知的规律，用物理或化学的基本定律直接推导出系统的数学模型。v但对于复杂的生物反应系统而言，其内部规律往往很不清楚，因而单靠理论分析无法建立其数学模型，必须用实验方法（即系统辩识的方法）才能建立。v系统辩识（扎德，Zadch定义）：“系统辩识是在输入输出的基础上，从一类系统中确定一个与所测系统等价的系统。”根据测试结果，由输入、输出间的数值关系来确定适宜的数学模型。v参数估计：由输入、输出间的数值关系来决定数学模型中待定的未知系数。v扎德（LAZadeh）1921年2月生于苏联巴库，1942年毕业于伊朗德黑兰大学电机工程系，获学士学位。1944年获美国麻省理工学院（MIT）电机工程系硕士学位，1949年获美国哥伦比亚大学博士学位，随后在哥伦比亚、普林斯顿等著名大学工作。从1959年起，在加里福尼亚大学伯克莱分校电机工程、计算机科学系任教授至今。 v扎德在20世纪50年代从事工程控制论的研究，在非线性滤波器的设计方面取得了一系列重要成果，已被该领域视为经典并广泛引用。3.2 系统辩识的基本方法(1)依据在线或离线的反应曲线、脉冲或阶跃响应曲线等确定模型结构与参数；v图解法（粗略估计）v曲线拟合（更精确）v典型例子：Monod方程(2)依据一段时期内各个采样时刻的输入和输出数据进行推算。v最小二乘参数估计（典型的手段）。v典型例子：描述对数生长期中生物质浓度随时间变化。步骤：v明确建立数学模型的目的（不同目的的辩识所需模型的精度不同，辩识的方法也就不同）。v掌握过程的事先资料（以便预知系统是否线性，有否时变，能否自衡，过渡时间长短及纯滞后存在情况等。这些资料包括系统内在的物理、化学及生物特性、测试数据和操作记录等历史资料）。v进行实验设计（选择和确定输入信号的幅度和频谱、采样时间、测试时间、开环或闭环辩识、离线或在线辩识、信号发生数据贮存、计算装置和信号滤波等）。3.3 参数估计的基本方法v建立输入输出模型：Y=f(X,K ) (3.1) Y：系统的m维输出向量，X为系统的r维输入向量，K为需要估计的n维未知参数。v观测被控过程，取得m组数据(X1i，X2i，Xri；Yi)，i=1，2，m (3.2)v参数估计：根据数据(3.2)估计模型(3.1)中的未知参数。vY与X的函数关系不同，参数估计的方法也不完全相同。3.3.1 线性模型的参数估计线性模型的参数估计1. 一维线性模型的参数估计v根据已知数据估计式(3.3)中未知参数K0和K1 ：Y=K0+ K1X (3.3)v以离差平方和最小为最佳估计准则，即目标函数： (3.4) v由此求得的K0和K1是最小二乘估计量。求得K0和K1的最小二乘估计为： (3.5) 2多维线性模型的参数估计v多维线性模型： (3.6)v实验观测数据： (3.7)参数集K0，K1， ，KN估计：v整理后得方程组： (3.8)v代入实测数据联立求解该方程组便可得到K0KN的最小二乘估计值。v衡量参数估计结果好坏的方法：(1) 以数据的拟合误差为标准（即可以针对参加拟合计算的数据,也可以用新取得的或不参予拟合的数据来检验模型）；(2) 用数理统计的方法对估计的结果进行分析。 3.3.2 非线性模型的参数估计非线性模型的参数估计v生物反应过程特性：多为非线性v参数估计方法：最小二乘法v一般的非线性模型形式：Y=f(X,K) (3.1)观测数据：(X1i，X2i，Xri；Yi)，i=1，2，m (3.2)性能指标为： (3.9)v多元函数求极值：非常复杂、难于求解。v单纯形法：一种最常用的直接搜索方法。 单纯形法单纯形法v是一种通过比较性能指标J值的大小，直接寻找J最小值的方法。v正规单纯形：Spendley等人v改正单纯形：Nelder和Mead，操作简单，便于在计算机上实现，广为应用。v初始单纯形构成和新的单纯形选取：关系到是否收敛及收敛速度快慢。什么是单纯形？v单纯形是空间的基元v一维单纯形v二维单纯形v三维单纯形v四维单纯形v 初始单纯形构成初始单纯形构成v定义：n维参数，用n+1个顶点构成的凸包（当性能指标J(x)为二元函数时，单纯形为三角形）。v若J(x)为n维函数，则取n+1个点x(0)，x(1) ，x(n)构成单纯形，并且保证n个向量x(1)-x(0)，x(2)-x(0)，x(n)-x(0)为线性独立。v正单纯形：单纯形各顶点之间的距离相等。若选定x(0)点，则x(0) ，x(1) ，x(n)由下式确定： (3.10)其中ei为第i个坐标单位向量。x右上角标i表示第i个点的序号。 单纯形新顶点取法单纯形新顶点取法v最差点的反射点。v以图3.2为例，点H、G、L分别表示性能指标最差、次差和最好点。先找出G、L两点的重心F，则H关于F的反射点就是新的顶点。v设X(H)、X(G)和X(L)分别为H、G、L三个点的坐标，则G、L两点的重心F为：(3.11)v反射点R的坐标为： (3.12)v对于n维变量，去除H点后重心F点坐标为： (3.13)v反射点R的坐标为： (3.14)v式中X(i)为第i个坐标向量，它有n个分量。 单纯形加速单纯形加速为了更快地搜索到极小值，还可以采取：v单纯形新点扩张v单纯形压缩v单纯形收缩新点扩张新点扩张v如果J(X(R)1）。v若J(X(E)J(X(R)，新的顶点仍取为R，见图3.2(b)。新点压缩新点压缩v若J(X(R)J(X(G)，可试探将反射点压缩到S点见图3.2(c)。S点坐标为： (3.16)式中为压缩系数，0.5J(X(H)，即反射点比最坏点还要差，则应该压缩到S点： (3.17) 为压缩系数，00.5。单纯形缩边单纯形缩边v如果HR方向所有各点的性能指标都比H点的差，这表明沿着HF方向反射是不对的。为有效地搜索，可把原单纯形HGL缩成HGL后从(2)重新开始，H 和G 点的坐标为： (3.18) (3.19)单纯形缩边如图3.2（d）所示。 单纯形寻优程序框图单纯形寻优程序框图v说明：程序框中第一步设置了最多次数M，以使搜索不到极小值时能够自动停机。图3.3 单纯形法寻优程序框图另一类方法：梯度法另一类方法：梯度法v使用函数的梯度(一阶导数)或它的Hessen矩阵(二阶导数)来构造算法的。v由于导数反映了函数值变化规律，因而在导数值可求得的情况下，充分利用函数的梯度信息，一般能获得加速收敛的效果。v常见的算法：共轭梯度法，拟牛顿法(变尺法)，Marquardt 法等等。 3.3.3 微分方程模型的参数估算微分方程模型的参数估算v微分方程模型是描述生物反应过程最常见的数学模型形式。v估计微分方程模型参数的方法解析解法数值解法 (1) 解析解法解析解法v能否采用解析解法取决于微分方程模型的积分式。v对于微生物分批培养而言，描述该系统的数学模型通常由如下两个微分方程组成： (3.20) (3.21)忽略ms,(3.20)/(3.21)得： (3.22)考虑ms的存在,作者把Pirt所提出的关系式： (3.23)改写成： (3.24)v把(3.20)代入(3.24)后整理得： (3.25)v以边界条件： 积分式(3.25)便可得到X与S的函数关系式： (3.26)v把式(3.26)代入式(3.21)得积分式： (3.27)v假定BS0，则式(3.27)中的对数项可写成： (3.28)v取ln (B+S)/S0展开式的第一项得： (3.29)v经分析上述假定及简化在较大范围内给X值的计算所带来的偏差很小，可忽略不计。据此，式(3.27)可简写成： (3.30)v式中M、N、P、C、D为A、B、E、F等的复合参数。v以边界条件：积分式(3.30)便可得到S-t关系式(3.31)v同理可得服从Monod方程、用于描述单罐连续培养微生物的动力学模型的解析解： (3.32) (3.33) v发酵动力学模型： (3.34)=(3.20) (3.35) (3.36)v式中P为产物浓度，YG和YP分别最大生物质生长得率和最大产物生成得率。这是一组含有七个参数(umax、Ks、YG、YP和ms)的多元非线性微分方程组。经研究，可导出如下所示的解析方程组： (3.37) (3.38) (3.39)式中：a1a4和b1b9均为所述七个参数的复合常数。(2) 数值解法数值解法v随着计算技术的发展，用数值解法求取用微分方程形式表示的动力学参数日益受到重视。这是由于该方法即不受参数数目多少的限制，又不受动力学方程的形式和实验方法的约束。v常见的方法是龙格库特法。设用于描述生物反应过程的微分方程为： (3.40)：n维状态向量(常表示X、S、P等)K：P维参数向量f：n维向量函数初始状态给定，应用龙格库特法可求得： (3.41)无噪声输出为： (3.42)若在 ti 时刻实际观测值为Yi, 按“离差平方和最小”原则求K的估计值，其性能指标： (3.43) 可用单纯形搜索法或其它搜索法求解（图3.4）。图图3.4 3.4 数值解法求取微分方程动力学参数框图数值解法求取微分方程动力学参数框图3.4 生物反应器的状态观测生物反应器的状态观测由于生物反应系统及其测量装置常常受到无法预计的随机噪声干扰随机噪声干扰随机噪声干扰随机噪声干扰，加上检测检测检测检测生物反应过程变量的种种局限性局限性局限性局限性，往往无法观测到系统的全部状态变量。因此，工程上要实现闭环最优控制，首先必须根据生物反应中所获得的不完全且夹杂着随机噪声的过程状态变量检测信息进行全系统的状态观测状态观测状态观测状态观测。（略）