数值分析第3章5-7节.ppt

1 3 5曲线拟合的最小二乘法 3 5 1最小二乘法及其计算 在函数的最佳平方逼近中如果只在一组离散点集上给定 这就是科学实验中经常见到的实验数据的曲线拟合 2 记误差 则的各分量分别为个数据点上的误差 3 设是上线性无关函数族 在中找一函数 使误差平方和 5 1 这里 5 2 4 这个问题称为最小二乘逼近 几何上称为曲线拟合的最小二乘法 用最小二乘求拟合曲线时 首先要确定的形式 确定的形式问题不仅是数学问题 还与问题的实际背景有关 通常要用问题的运动规律及给定的数据进行数据描图 确定的形式 然后通过实际计算选出较好的结果 5 为了使问题的提法更有一般性 通常在最小二乘法中考虑加权平方和 5 3 这里是上的权函数 它表示不同点处的数据比重不同 就是次多项式 若是次多项式 的一般表达式为 5 2 表示的线性形式 6 这样 最小二乘问题就转化为求多元函数 5 4 的极小点问题 由求多元函数极值的必要条件 有 使 5 3 取得最小 7 若记 5 5 上式可改写为 5 6 这个方程称为法方程 可写成矩阵形式 8 其中 5 7 要使法方程 5 6 有惟一解 就要求矩阵非奇异 9 显然在任意个点上满足哈尔条件 函数的最小二乘解为 定义10 方程 5 6 存在惟一的解 从而得到 于是 10 这样得到的 对任何形如 5 2 的 都有 故确是所求最小二乘解 11 一般可取 但这样做当时 通常对的简单情形都可通过求法方程 5 6 得到 给定的离散数据 例如 求解法方程 5 6 将出现系数矩阵为病态的问题 若两边取对数得 12 例7 这样就变成了形如 5 2 的线性模型 此时 若令 则 已知一组实验数据如下 求它的拟合曲线 13 解 从图中看到各点在一条直线附近 故可选择线性函数作拟合曲线 将所给数据在坐标纸上标出 见图3 4 图3 4 14 令 这里 故 15 解得 由 5 6 得方程组 于是所求拟合曲线为 16 关于多项式拟合 Matlab中有现成的程序 其中输入参数为要拟合的数据 为拟合多项式的次数 输出参数为拟合多项式的系数 利用下面的程序 可在Matlab中完成上例的多项式拟合 17 x 11233345 f 444 566688 5 aa poly x f 1 y polyval aa x plot x f r x y k xlabel x ylabel y gtext y s1 x 18 结果如下 19 例8 设数据由表3 1给出 用最小二乘法确定及 解 表中第4行为 通过描点可以看出数学模型为 它不是线性形式 用给定数据描图可确定拟合曲线方程为 两边取对数得 20 若令 先将转化为 为确定 根据最小二乘法 取 则得 数据表见表3 1 得 21 故有法方程 解得 于是得最小二乘拟合曲线为 22 利用下面的程序 可在Matlab中完成曲线拟合 x 1 001 251 501 752 00 y 5 105 796 537 458 46 y1 log y aa poly x y1 1 a aa 1 b exp aa 2 y2 b exp a x plot x y r x y2 k xlabel x ylabel y gtext y a exp bx 23 结果如下 24 3 5 2用正交多项式做最小二乘拟合 如果是关于点集 5 8 用最小二乘法得到的法方程组 5 6 其系数矩阵是病态的 25 5 9 则方程 5 6 的解为 且平方误差为 26 接下来根据给定节点及权函数 构造带权正交的多项式 注意 用递推公式表示 即 5 10 这里是首项系数为1的次多项式 27 5 11 下面用归纳法证明这样给出的是正交的 28 假定对及 要证对均成立 由 5 10 有 由 5 10 第二式及 5 11 中的表达式 有 均成立 5 12 29 而 于是由 5 12 当时 另外 是首项系数为1的次多项式 它可由 由归纳法假定 当时 的线性组合表示 由归纳法假定又有 30 由假定有 再考虑 5 13 利用 5 11 中表达式及以上结果 得 31 至此已证明了由 5 10 及 5 11 确定的多项式 组成一个关于点集的正交系 用正交多项式的线性组合作最小二乘曲线拟合 只要根据公式 5 10 及 5 11 逐步求的同时 相应计算出系数 最后 由和的表达式 5 11 有 32 并逐步把累加到中去 最后就可得到所求的 用这种方法编程序不用解方程组 只用递推公式 并且当逼近次数增加一次时 只要把程序中循环数加1 其余不用改变 这里可事先给定或在计算过程中根据误差确定 拟合曲线 33 3 6最佳平方三角逼近与快速傅里叶变换 当是周期函数时 显然用三角多项式逼近比用代数多项式更合适 本节主要讨论用三角多项式做最小平方逼近及快速傅里叶变换 简称FFT算法 34 3 6 1最佳平方三角逼近与三角插值 设是以 为周期的平方可积函数 用三角多项式 6 1 作为最佳平方逼近函数 由于三角函数族 在上是正交函数族 于是在上的最小平方三角逼近多项式的系数是 35 称为傅里叶系数 函数按傅里叶系数展开得到的级数 6 3 就称为傅里叶级数 6 2 36 只要在上分段连续 则级数 6 3 一致收敛到 对于最佳平方逼近多项式 6 1 有 由此可以得到相应于 4 11 的贝塞尔不等式 因为右边不依赖于 左边单调有界 所以级数 37 当只在给定的离散点集 上已知时 则可类似得到离散点集正交性与相应的离散傅里叶系数 下面只给出奇数个点的情形 收敛 并有 38 可以证明对任何成立 令 39 这表明函数族在点集 上正交 若令 其中 40 当时 于是 6 4 就是三角插值多项式 系数仍由 6 4 表示 41 由于 所以函数族在区间上是正交的 42 当时 个复向量具有如下正交性 6 5 43 事实上 令 于是 即 若 若 则有 则 从而 44 于是 若 这就证明了 6 5 成立 即是正交的 则 于是 因此 在个点上的最小二乘傅里叶逼近为 45 6 6 其中 6 7 于是由 6 6 得 即 6 8 46 而 6 8 是由求的过程 称为反变换 47 3 6 2快速傅氏变换 FFT 不论是按 6 7 式由求 由求 6 9 其中 正变换 或 反变换 还是由 6 4 计算傅里叶逼近系数 都可归结为计算 是已知复数序列 或是按 6 8 48 当较大且处理数据很多时 就是用高速的电子计算机 很多实际问题仍然无法计算 直到20世纪60年代中期产生了FFT算法 大大提高了运算速度 从而使傅氏变换得以广泛应用 FFT算法的基本思想就是尽量减少乘法次数 49 用 6 9 计算全部 表面看要做个乘法 特别当时 只有个不同的值 因此可把同一个对应的相加后再乘 这就能大量减少乘法次数 50 设正整数除以后得商及余数 则 称为的同余数 以表示 由于 因此计算时可用的同余数代替 从而推出FFT算法 以为例 说明FFT的计算方法 由于则 6 9 的和是 6 10 故有 51 将用二进制表示为 其中只能取0或1 例如 根据表示法 有 公式 6 10 可表示为 52 6 11 若引入记号 6 12 53 则 6 11 变成 它说明利用同余数可把计算分为步 用公式 6 12 计算 每计算一个只用2次复数乘法 计算一个用 次复数乘法 计算全部共用次复数乘法 若注意公式 6 12 还可进一步简化为 54 将这表达式中二进制表示还原为十进制表示 55 6 13 同样 6 12 中的也可简化为 即 即得 56 把二进制表示还原为十进制表示 得 6 14 同理 6 12 中可简化为 即 57 表示为十进制 有 6 15 58 根据公式 6 13 6 14 6 15 由 逐次计算到 见表3 2 略 上面推导的的计算公式可类似地推广到的情形 根据公式 6 13 6 14 6 15 一般情况的FFT计算公式如下 59 6 16 其中 从出发 由到算到 一组占用个复数单元 计算时需给出两组单元 括号内的数代表它的位置 在计算机中代表存放数的地址 即为所求 60 这个计算公式除了具有不倒地址的优点外 计算只有两重循环 计算过程中只要按地址号存放则最后得到的 就是所求离散频谱的次序 外循环由计算到 内循环由计算到 由计算到 更重要的是整个计算过程省计算量 由公式看到算一个共做次复数乘法 而最后一步计算时 由于 61 当时比值是它比一般FFT的计算量 次乘法 也快一倍 我们称 6 16 的计算公式为改进的FFT算法 62 3 7有理逼近 3 7 1有理逼近与连分式 有理函数逼近是指用形如 的函数逼近 与前面讨论一样 如果最小就可得到最佳有理一致逼近 7 1 63 如果最小则可得到最佳有理平方逼近函数 本节主要讨论利用函数的泰勒展开获得有理逼近函数的方法 对函数用泰勒展开得 7 2 取部分和 64 7 3 65 7 4 7 3 右端为的无穷连分式的前5项 最后式子 若取 7 3 的前2 4 6 8项 则可分别得到的以下有理逼近 是它的紧凑形式 66 若用同样多项的泰勒展开部分和逼近 并计算处的值及 计算结果见表3 2 的准确值为 从表3 1可以看出 67 但它们的计算量是相当的 这说明用有理逼近比多项式逼近好得多 由此看出的精度比高出近10万倍 例9 用辗转相除法将它化为连分式并写成紧凑形式 解 给出有理函数 用辗转相除可逐步得到 68 本例中用连分式计算的值只需3次除法 1次乘法和7次加法 69 若直接用多项式计算的秦九韶算法则需6次乘法和1次除法及7次加法 可见将化成连分式可节省计算乘除法次数 对一般的有理函数 7 1 可转化为一个连分式 它的乘除法运算只需次 而直接用有理函数 7 1 计算乘除法次数为次 70 3 7 2帕德逼近 利用函数的泰勒展开可以得到它的有理逼近 设在的泰勒展开为 7 5 它的部分和记作 7 6 71 定义11 设 其中无公因式 且满足条件 7 8 则称为函数在处的阶帕德逼近 记作 简称的帕德逼近 如果有理函数 7 7 72 根据定义 若令 则满足条件 7 8 等价于 即 由于应用莱布尼茨求导公式得 73 这里是由 7 6 得到的 上式两端除 并由可得 7 9 及 7 10 注意当时 故 7 10 可写成 74 7 11 其中时 若记 7 12 75 则方程组 7 11 的矩阵形式为 定理10 76 根据定理10 求的帕