化工数值计算-6.ppt

回归分析和曲线拟合生产过程和科学实验中 常用的变量大体可分两类 一类为确定性变量 另一类为随机变量 确定性变量是指两个或多个变量之间有确定的关系 即其中某个变量的每个值 都与一变量的一个或几个完全确定的值相对应 即它们之间存在着 函数关系 例如 理想气体的压力P与摩尔体积V间 存在着确定的函数关系 但在实际问题中 由于变量之间的关系比较复杂 或由于生产或实验过程中不可避免地存在着误差 使变量之间的关系具有不确定性 也就是说 某个变量对应的 不是一个或几个确定的值 而是整个集合的值 这时 变量x和y间的关系 就称为相关关系 例如 流体在圆形直管中做湍流时的情形 通过量纲分析可知 努塞尔特数Nu 普兰特 数Pr和雷诺数Re之间存在着如下相关关系 这种关系的不确定性 表现为式中a和b的数值 在每次测量中不尽相同 不确定的原因 首先是影响该过程的因素甚多 有些因素至今尚未弄清 其次是受到实验过程中的偶然因素影响 这 种不确定性关系并不说明上述三个量纲为1的数群之间无规律可循 相反 通过大量试验 人们发现 a和b的数值总是围绕着某一定值波动 而且随着试验次数的增多 a b的数值趋于稳定 a b的稳定值 可作为a和b的最佳估计值 在一定条件下 a 0 023 b 0 8 由此可见 通过大量试验 是可以找到隐藏在随机性后面的统计规律性的 回归分析和曲线拟合是一种处理变量相关关系的数理统计方法 用它可以寻找隐藏在随机性后面的统计规律性 函数与相关是两种不同类型的变量关系 它们之间并无严格界限 一方面 相关的变量之间 并无确定的关系 但在一定的条件下 从一定的统计意 义上看 它们之间又可能存在着某种确定的函数关系 另一方面 由于实际测定的数据中 总存在着误差 即使是确定性变量 也会出现某些非确定性结果 6 1一元线性回归 一元线性回归处理的是两个变量之间的线性关系 所用的数学模型为一元线性代数模型 其模型方程式是对这种模型参数的估计 就是根据原始数据点 x1 y1 x2 y2 xi yi xn yn 确定式 6 1 中a b的估计值 在实际体系中 自变量x与因变量y之间服从线性关系的情况虽然不多 但在不少情况下 x y之间存在着某种函数组合关系 例如f1 x y f2 x y 设两个函数之间服从线性关系 f1与f2是不含待定系数的已知函数 若把f1 x y 与f2 x y 分别视为自变量与因变量 则仍可以借用线性模型去估计其参数值 这种方法称为化直法 它在化学化工的实际问题中是常见的 例如单分子基元反应AB的动力学方程式为 对上式积分得式中 cA t是不呈线性关系的函数 若对方程两边取对数 上式可化为lncA t的线性函数 又例如 按照阿仑尼乌斯定律 反应速率常数k与温度T之间不呈线性关系 但lnk与1 T则呈线性关系 这些都是属于可化为线性关系的例子 一元线性代数模型中的待定参数a和b 称为 估计值 之所以称为 估计 值 是因为a b的值是从实验值中通过数理统计方法确定的 图6 1一元线性回归 6 1 1方法概述设有一组实验数据 x1 y1 x2 y2 xn yn 自变量x与因变量y存在着式 6 1 的关系 当x取值为xi时 y的测定值为yi 计算值为yi 并有由于参数a b为未知值 故yi 也是未知值 若将全部实验 数据标绘在x y图中 见图6 1 由于各种因素的影响 它们不会全部落在一条直线上 即n个yi不会与n个yi 完全重合 它们将随机地分布在与xi呈线性关系的yi 的周围 以 i表示它们之间的差值 则有这里 i就是误差 它反映了xi使yi偏离直线的各种影响因素的总和 现在 要寻找一条最靠近各个数据点的直线 这条直线称为回归直线 由于回归直线是一切直线中最接近各数据点 xi yi 的 用它代表x与y之间的线性关系 比任何其他直线更为可靠 究竟如何确定回归曲线中的参数a和b呢 目前最常用的方法就是最小二乘法 即残差平方和最小法 式 6 3 中的误差 i又称为残差 表示第i个数据与回归直线的偏离程度 则残差平方和 Q表示全部数据与回归直线的总偏离程度 显然Q是a和b的函数 不用残差和 i的原因是 i有正有负 相加时可能彼此抵消 从而不能反映总的偏离程度 而用残差的平方和不会发生 这种现象 由多元函数的极值理论可知 要使Q值最小 a b必须满足下列条件 即得式 6 6 称为一元线性回归的正规方程组 通过求解该方程组 可得 式 6 7 中等号右侧的量全部取自原始数据 因此 就可以确定回归系数a和b 完成参数估计 为了简化a和b的表达式 定义 式中 分别为xi和yi的平均值 xi与之差 xi 称为xi的离差 全部xi的离差平方和 称为x的 离差平方和 记为Lxx yi与之差 yi 称为yi的离差 全部yi的离差平方和 称为y的离差平方和 记为Lyy 同理再令Lxy为全部xi的离差与yi的离差乘积的总和 将以上关系式代入式 6 7 得 由式 6 12 第二式可以看出 回归直线是通过点 的 从力学观点看 相当于n个实验点 xi yi 的重心 回归直线是通过重心的 应当指出 残差 i只用yi y i表示时 表明yi有测量误差 而xi无测量误差 或表示与yi相比 xi的误差很小 因此 测量误差使实验点偏离回归直线 都表现为yi偏离y i 如果xi的误差与yi的误差相比 不可忽略 则两者都必须考虑 这种情况比较复杂 此处不予介绍 求回归方程的计算过程中 不需要事先假定两个变量之间必须有相关关系 即使是一组杂乱无章的数据 也可以用最小二乘法绘制一条直线 以表示x与y的关系 显然 这种情况下 绘制的直线并无实际意义 为了判断两个变量间线性关系的优劣程度 引入一个新的指标R 称为简单相关系数 它的定义为R值不同时 数据点的分布情况如下 1 R 0 图6 2R 0的数据点分布 此时Lxy 0 b 0 即回归直线平行于x轴 y的变化与x无关 表示数据点的分布是无规则的 如图6 2所示 但亦有当R 0时 x与y确实存在明显相关性的情况 这种情形 不能应用线性回归方法 只能用化直线法或曲线拟合法处理 2 0 R 1 绝大多数属于这种情况 它表示x与y之间存在着一定的线性关系 如图6 3所示 R 0时 b 0 数据点的y值随着x增加而增加 这种情况称为x与y正相关 R 0时 b 0 数据点的y值随着x增加而减小 这种情况称为x与y负相关 R的绝对值越小 数据点沿回归直线越分散 图6 30 1的数据点分布 1的数据点分布 3 R 1x与y完全相关 全部数据点均落在回归直线上 若x与y为非线性相关 但经变量变换后 用回归直线的方法处理 所求得的回归系数仅对变换后的变量是最佳的 而对原变量来说则并非最佳 但通常还能令人满意 此时应注意原变量 的残差平方和并非最小 由以上讨论可知 相关系数R的绝对值在0与1之间 而且越接近于1 其线性关系越密切 那么 R 与1接近到什么程度 才能说明x与y之间存在线性相关关系呢 要回答这个问题 就要对相关系数进行显著性检验 由于篇幅所限 有关相关系数的显著性检验和回归方程的方差分析等问题将不在此讨论 如有需要 可参考有关数理统计方面的书籍 6 1 2程序框图图6 4是一元线性回归的通用计算程序框图 程序框图中的主要变量 N数据点数X Y一维数组 用于存放原始数据中的x和y值 XXLx离差平方和LxxYYLy离差平方和LyyXYLx离差与y离差乘积总和LxyA回归直线截距aB回归直线斜率b R简单相关系数6 1 3计算实例 6 2多元线性回归一元线性代数模型中 只有一个自变量 若有多个因素影响体系的性质时 必须考虑因变量y与多个自变量xl x2 xn之间的关系 例如 化学反应速率要受到温度 压力和浓度的影响 在 气相反应动力学中 反应动力学方程可表示为式中 r为反应速率 pA pB pC分别为反应物A B C组分的分压 a b c是方程式中的待定指数 k为反应速率常数 若将上式取对数得到再令y lnr d lnk x1 lnpA x2 lnpB x3 lnpC 则得 可见该式具有多元线性方程式的特征 a b c d为系数 x1 x2 x3为自变量 多元线性方程的普遍式为它是一个含有m 1个系数的m元线性方程式 下 面介绍多元线性回归的最小二乘法 6 2 1方法概述设x取值为xi1 xi2 xim时 实验测定的y值为yi i 1 2 n 由于测定值yi存在着误差 所以会偏离线性关系 现在要寻找一组aj的估计值以构成回归方程 确定aj的原则 仍然是使yi的实验值与回归方程计算值的残差平方和最小 即使 最小 式中i表示实验点序号 i 1 2 n j表示自变量标号 j 1 2 m 自变量xij为第j个自变量的第i次测定值 此外注意仅n m 1才能求出上式中的m 1个回归系数 同样由多元函数的极值理论可知 要使Q值最小 a0和aj必须满 足下列条件 式 6 15 经整理可得 式 6 16 称为多元线性回归模型的正规方程组 它是一个m 1元的线性代数方程组 由于xij和yi已知 故可求得m 1个待定 系数a0 a1 am 实际计算时 一般作如下处理 先将式 6 16 的第一式写成然后将式 6 17 代入方程组 6 16 的第2至第m 1式 重新组成一个m元线性方程组 其中有a1 a2 am等m个待定系数 通过求解此m元线性方程组 获得系数a1 a2 am 再代回式 6 17 求得a0 为简化计算 用表示第j个x的平均值 表示y的平均值 则用Ljk表示第j个x离差与第k个x离差乘积之和 则 用Lyy表示y离差的平方和 则 用Ljy表示第j个x离差与y离差乘积之和 则将式 6 17 分别代入式 6 16 的第2至m 1式 经简化整理可得如下m元线性方程组 可用主元素消去法求解此式 然后将求得的a1 a2 am代入式 6 17 求出a0 从而完成对多元线性回归模型的参数估计 多元线性回归的计算中 常用复相关系数衡量数据点之间的线性优劣 复相关系数定义如下 式中 U称为回归平方和 应当指出 并非所有曲线都可以按这种方法处理 例如抛物线就不能通过变量变换把它化为直线 但是如果令x1 x x2 x2 则上式就化成一个包含两个自变量的线性方程 从而将抛物线按二元线性回归计算 对于含多变量的任意多项式也可以通过类似的变换 把它们转化成多元线性回归计算 6 2 2程序框图图6 6是多元线性回归的通用计算程序框图 图6 6 a 多元线性回归的通用计算程序框图 1 图6 6 b 多元线性回归的通用计算程序框图 2 程序框图中的主要变量 N数据点数M多元线性模型元数X二维数组 用于存放原始数据的x值Y一维数组 用于存放原始数据的y值 YP值YYLLyy值XP一维数组 用于存放值A二维数组 用于存放m元线性方程组的系数LjkB一维数组 用于存放m元线性方程组的常数项Ljy C一维数组 用于存放多元线性模型的系数aj j 0 1 M R复相关系数R0U回归平方和Q残差平方和子程序XYF为列主元消去法求解线性方程组的程序 可参见 图5 2和图5 3 6 2 3计算实例 6 3剔除可疑数据及其计算程序6 3 1剔除可疑数据的方法在线性回归计算中 假定每个测定数据与回归结果之间的误差均在随机误差允许的范围之内 然而 由于测量误差或过失误差等多种原因 在一组实验值中 误差往往会超出随机误 差的允许范围 这些数据 称为可疑数据 为保证回归结果的可靠性 必须剔除这些可疑的数据 剔除可疑数据 应当有一个科学的标准 这个标准就是统计判据 属于统计判据的剔除准则有多种 以一元线性回归为例 其代数模型为y a bx 若自变量x无测量误差 则y的标准偏差为 式中 n为原始数据点数 m为回归模型中自变量的个数 对一元线性回归m 1 i为残差 即 i yi a b是按最小二乘法求出的最佳估计值 根据数理统计分析 合理的数据 其残差不应超出 的k倍 若取k 3 便是常 用的3 准则 据此 可以把残差绝对值超过3 的个别数据 xi yi 判为可疑数据而加以剔除 必须指出 3 准则是以数据点数n 为前提的 当n为有限值时 3 判据并不十分可靠 下面介绍一种广泛采用的判据 即所谓肖维奈特准则 按肖维奈特准则 若n次等精度测量中 有某个测量值yi 其残差的绝对值超出k 就可以认为是可疑数据而予以剔除 表6 1列出了肖维奈特准则中与n相对应的k值 表6 1肖维奈特准则的n和k值 使用这个准则时 可根据回归结果 对全部实验值进行逐级检查 把属于可疑数据的实验值选出 若发现不止一个可疑数据 则应把其中残差绝对值最大者剔除 然后重新计算 值 根据新的 值 再次用肖维奈特准则进行检查 每次只剔除一个可疑数据 其余数据重新进行回归 直至回归所用的数据中不再含有可疑数据为止 6 3 2剔除可疑数据的计算程序框图 图6 7是具有剔除可疑数据功能的一元线性回归通用计算程序框图 整个计算过程分为输入原始数据 一元线性回归计算 确定肖维奈特准则的k值 确定残差绝对值最大的数据点 剔除最可疑数据点 即残差绝对值最大的数据点 图6 7具有剔除可疑数据功能的一元线性回归通用计算程序框图 程序框图中的主要变量 N原始数据点数或剔除可疑数据后的合格数据点数N1可疑数据点数X一维数组 用于存放原始数据及合格数据中的x值Y一维数组 用于存放原始数据或合格数据中的y值 X1一维数组 用于存放可疑数据点的x值Y1一维数组 用于存放可疑数据点的y值A回归直线截距B回归直线斜率R简单相关系数 SD标准偏差 ER平均相对误差DALTA绝对值最大的残差ID残差绝对值最大的数据点序号U肖维奈特准则的k值 子程序LINEAR1A为一元线性回归计算子程序 比例6 1中的子程序LINEAR1增加了标准偏差和平均相对误差的计算 子程序RULES为肖维奈特准则中k值的计算程序 采用类似的方法 可以编写能剔除可疑数据的多元线性回归计算程序框图 6 3 3计算实例 6 4多项式拟合在化学化工的实验或科研中 经常需要从一组测定数据 例如从n对 xi yi 数据 去求自变量x和因变量y的近似函数关系式y p x 从图形上看 这是由给定的n个点 xi yi i 1 2 n 作曲线 拟合 在曲线拟合中 多项式拟合问题占特殊的地位 任何函数在一个比较小的范围内 可以用多项式任意逼近 因此 在比较复杂的实际问题中 可以不问y与各因素的确切关系 而用多项式拟合进行分析和计算 下面以多项式拟合为例 说明曲线拟合的方法和计算程序 6 4 1方法概述设用下列m次多项式 拟合一组数据 xi yi i 1 2 n 即曲线y f x 上已给定n个点 用多项式求作该曲线的近似图形 这一问题与前述的插值问题有类似之处 但插值问题要求近似曲线y p x 严格地通过所给的n个点 这一要求将会使近似曲线y p x 保留数据的全部测试点的测量误差 如果个别数据的误差很大 那么插值的效果显然是不够理想的 鉴于这种情况 考虑放弃严格通过所有结点 xi yi 这一要求 而采用别的方法去构造近似曲线 以尽可能反映所给数据的总趋势 曲线拟合的常 用方法仍然是最小二乘法 即残差平方和最小法 若以 i代表结点处的残差 则残差的平方和为 由于xi与yi为已知值 故Q是aj j 0 1 2 m