200字范文,内容丰富有趣,生活中的好帮手!
200字范文 > 2月18 多项式拟合原理 全过程 loss函数 正则化 迭代与梯度下降 plotly

2月18 多项式拟合原理 全过程 loss函数 正则化 迭代与梯度下降 plotly

时间:2018-12-11 02:48:41

相关推荐

2月18 多项式拟合原理 全过程 loss函数 正则化 迭代与梯度下降 plotly

资料库:暂未更新

文章目录

1.0多项式的定义:1.1问题等价:多项式拟合=关于多项式系数 W 的线性函数的求解1.2问题实质:通过误差函数来对拟合进行评估,并得出最优的多项式系数2.0误差函数的定义2.1误差函数也称为损失函数lost或者代价cost函数2.2误差函数2.2.1 **第一类:** 适用于回归问题(Regression)的误差函数,这种误差函数的目标是量化推测值和真实值的逻辑距离,理论上我们可以使用任何距离计算公式作为误差函数。实践中为常用的是以下两种距离:2.2.2 **第二类:** 适用于分类问题(classification)的误差函数,分类问题的目标是推测出正确的类型,一般使用概率描述推测属于某种类型的可能性,因此误差函数就需要能够计算两个概率分布之间的"距离".2.2.2.1 最常用的此类方法是交叉熵损失函数「Cross Entropy Loss」2.2.2.2 其他的常见loss2.2.3 **第三类:** 多任务问题 此类问题是指使用一个模型同时学习多个指标,比如使用深度学习解决计算机视觉中的目标定位(object localization)问题,模型需要同时学习对象类型和对象位置两个指标,因此误差函数需要具备同时衡量类型误差和位置误差的能力,常见的做法是先分别设计类型误差函数 Lclass 和位置误差函数 Lposition,在按一定比例(由超参数 α 和 β 确定)合并这两个误差从而形成一个新的复合误差函数 L 来综合的反应总误差。2.2.4 交叉熵函数与最大似然函数的联系和区别2.2.5 损失函数和激活函数的组合2.2.1如何设计误差函数任务目标替代误差函数(Surrogate loss function)拓展:(该部分不是本文章主题内容):基于神经网络的误差函数设计2.3基于梯度下降的优化问题2.3.1可导性2.3.2非凸性2.3.4 泛化(Generalization)2.3.5 正则化3.0 最小二乘法—多项式拟合非线性函数3.1 函数的近似表示—高次多项式3.2误差函数—最小二乘法MSE3.3引出案例函数曲线3.4目标函数3.5 优化目标函数3.5.1 有两种方法可以求解3.5.2 优化目标函数—梯度下降法3.5.3 优化目标函数—求解线性方程组3.5.3.1 求解线性方程组,大概有如下算法附录:python科学计算中的常用方法-曲线拟合plotly库附录测试代码未完成项目:三种梯度下降方式的理解和优劣,梯度下降(Gradient Descent)小结,迭代

1.0多项式的定义:

其中 M是多项式的阶数,w0,w1,w2…wM是多项式的系数,记做W,可以看到虽然多项式函数y(x, W) 是关于 x 非线性函数,但是却是关于多项式系数 W 的线性函数。有了多项式以后我们还需要一个误差函数来对我们拟合出的多项式进行评估。误差函数很多类型。

1.1问题等价:多项式拟合=关于多项式系数 W 的线性函数的求解

1.2问题实质:通过误差函数来对拟合进行评估,并得出最优的多项式系数

2.0误差函数的定义

2.1误差函数也称为损失函数lost或者代价cost函数

参考:/developerworks/cn/cognitive/library/cc-lo-talking-about-the-loss-function/index.html

基于梯度下降算法的机器学习问题本质上都可以抽象成误差函数的优化过程。

误差函数一种量化模型拟合程度的工具,机器学习(基于梯度下降算法)的基本思想是设计一个由参数 θ 决定的模型,使得输入 x 经过模型计算后能够得到预测值,模型的训练过程是用训练数据 xi 输入模型计算得到预测值并计算预测值和真实值 yi 的差距 L,通过调整模型的参数 θ 来减小差距 L 直到这个差距不再减小为止,这时模型达到最优状态。误差函数就是被设计出来的一个关于 θ 的用来衡量预测值和真实值之间的差距的函数,这样通过求的最小值,就可以获得一个确定的 θ, 从而够得到最佳预测结果的模型。所以机器学习实际上是一个通过优化误差函数(求最小值)来确定 θ 的过程。

2.2误差函数

损失函数用来评价模型的预测值和真实值不一样的程度,损失函数越好,通常模型的性能越好。不同的模型用的损失函数一般也不一样。可以采用例如:MSE(最小均方误差)、LMS(最小均方)和LS(最小二乘)

损失函数分为经验风险损失函数结构风险损失函数。经验风险损失函数指预测结果和实际结果的差别,结构风险损失函数是指经验风险损失函数加上正则项

根据不同类型机器学习任务可以将误差函数主要可以以下三类:

2.2.1第一类:适用于回归问题(Regression)的误差函数,这种误差函数的目标是量化推测值和真实值的逻辑距离,理论上我们可以使用任何距离计算公式作为误差函数。实践中为常用的是以下两种距离:

MSE(Mean Squared Error,L2) 均方误差MAE(Mean Absolute Error, L1)

例如:在小批量梯度下降(Mini-batch gradient descent)算法中,通常求平均值作为批量的误差来作为反向传递

(back propagation) 的输入,即 batch_size=m 时,

2.2.2第二类:适用于分类问题(classification)的误差函数,分类问题的目标是推测出正确的类型,一般使用概率描述推测属于某种类型的可能性,因此误差函数就需要能够计算两个概率分布之间的"距离".

2.2.2.1 最常用的此类方法是交叉熵损失函数「Cross Entropy Loss」

代码:

"""实现常见的损失函数"""import numpy as npdef cross_entropy(Y, P):"""Cross-Entropy loss function.以向量化的方式实现交叉熵函数Y and P are lists of labels and estimationsreturns the float corresponding to their cross-entropy."""Y = np.float_(Y)P = np.float_(P)return -np.sum(Y * np.log(P) + (1 - Y) * np.log(1 - P)) / len(Y)

参考:/red_stone1/article/details/80735068

参考:2月7日 SVM&线性回归&逻辑回归 /djfjkj52/article/details/104210888

从上面两种图,可以帮助我们对交叉熵损失函数有更直观的理解。无论真实样本标签 y 是 0 还是 1,L 都表征了预测输出与 y 的差距。

另外,重点提一点的是,从图形中我们可以发现:预测输出与 y 差得越多,L 的值越大,也就是说对当前模型的 “ 惩罚 ” 越大,而且是非线性增大,是一种类似指数增长的级别。这是由 log 函数本身的特性所决定的。这样的好处是模型会倾向于让预测输出更接近真实样本标签 y。

交叉熵损失函数的其它形式

什么?交叉熵损失函数还有其它形式?没错!我刚才介绍的是一个典型的形式。接下来我将从另一个角度推导新的交叉熵损失函数

第一种形式在实际应用中更加常见,例如神经网络等复杂模型;第二种多用于简单的逻辑回归模型。

2.2.2.2 其他的常见loss

参考:/p/58883095

2.2.3第三类:多任务问题 此类问题是指使用一个模型同时学习多个指标,比如使用深度学习解决计算机视觉中的目标定位(object localization)问题,模型需要同时学习对象类型和对象位置两个指标,因此误差函数需要具备同时衡量类型误差和位置误差的能力,常见的做法是先分别设计类型误差函数 Lclass 和位置误差函数 Lposition,在按一定比例(由超参数 α 和 β 确定)合并这两个误差从而形成一个新的复合误差函数 L 来综合的反应总误差。

2.2.4 交叉熵函数与最大似然函数的联系和区别

参考:/wydbyxr/article/details/83212703

对数似然函数值/最大近然估计/log likelihood:在参数估计中有一类方法叫做“最大似然估计”,因为涉及到的估计函数往往是是指数型族,取对数后不影响它的单调性但会让计算过程变得简单,所以就采用了似然函数的对数,称“对数似然函数”。

根据涉及的模型不同,对数函数会不尽相同,但是原理是一样的,都是从因变量的密度函数的到来,并涉及到对随机干扰项分布的假设。

最大似然估计法的基本思想:极大似然原理的直观想法是:一个随机试验如有若干个可能的结果A,B,C,…。若在一次试验中,结果A出现,则一般认为试验条件对A出现有利,也即A出现的概率很大。

最大似然估计法的思想很简单:在已经得到试验结果的情况下,我们应该寻找使这个结果出现的可能性最大的那个X作为真X的估计。求X的极大似然估计就归结为求L(X)的最大值点,而由于对数函数是单调增函数,所以对L(X)取log。

对log(L(X))关于X求导数,并命其等于零,得到的方程组称为似然方程组。解方程组log(L(X)),又能验证它是一个极大值点,则它必是L(X)的最大值点,即为所求的最大似然估计。

理解对数似然估计函数值

对数似然估计函数值一般取负值,实际值(不是绝对值)越大越好。

1)第一,基本推理。对于似然函数,如果是离散分布,最后得到的数值直接就是概率,取值区间为0-1,对数化之后的值就是负数了;如果是连续变量,因为概率密度函数的取值区间并不局限于0-1,所以最后得到的似然函数值不是概率而只是概率密度函数值,这样对数化之后的正负就不确定了。

2)第二,Eviews的计算公式解释。公式值的大小关键取之于残差平方和(以及样本容量),只有当残差平方和与样本容量的比之很小时,括号内的值才可能为负,从而公式值为正,这时说明参数拟合效度很高;反之公式值为负,但其绝对值越小表示残差平方和越小,因而参数拟合效度越高。

1.交叉熵函数与最大似然函数的联系和区别?

区别:交叉熵函数使用来描述模型预测值和真实值的差距大小,越大代表越不相近;似然函数的本质就是衡量在某个参数下,整体的估计和真实的情况一样的概率,越大代表越相近。

联系:交叉熵函数可以由最大似然函数在伯努利分布的条件下推导出来,或者说最小化交叉熵函数的本质就是对数似然函数的最大化

推导过程:/p/58883095

2.2.5 损失函数和激活函数的组合

参考:/vivian_ll/article/details/78498432

1.均方差损失函数+Sigmoid激活函数(不推荐)

对于Sigmoid,当z的取值越来越大后,函数曲线变得越来越平缓,意味着此时的导数σ′(z)也越来越小。同样的,当z的取值越来越小时,也有这个问题。仅仅在z取值为0附近时,导数σ′(z)的取值较大。

在反向传播算法中,每一层向前递推都要乘以σ′(z),得到梯度变化值。Sigmoid的这个曲线意味着在大多数时候,我们的梯度变化值很小,导致我们的W,b更新到极值的速度较慢,也就是我们的算法收敛速度较慢。

使用交叉熵损失函数+Sigmoid激活函数(提高收敛速度

在用sigmoid作为激活函数的时候,为什么要用交叉熵损失函数,而不用均方误差损失函数?请参考:/p/58883095

2.2.1如何设计误差函数

误差函数不仅仅是误差的衡量工具,更重要的是梯度下降算法会为了不断缩小误差而沿着误差函数的导数方向调整模型参数,因此可以说误差函数它决定了模型学习的方向。使用不同的误差函数我们可以在完全相同的模型架构(model architecture)的基础上学习出不同的模型参数,来达到不同的目的。比如,在多层卷积神经网络架构上使用 cross entropy loss 可以判断图像对象的类型(是猫还是狗),而同样的网络架构配合复合误差函数则可以用来进行目标定位。从某种程度上说模型设计决定了能力,误差函数设计决定了方向。 那么如何设计(选择)误差函数呢?可以从以下几个方面考虑。

任务目标

决定误差函数设计的最重要的因素就是任务目标,有时任务目标和误差函数的关系很直接,比如预测价格。 有时他们的关系就不那么明显,需要一些的专业知识才能和误差函数建立联系。另外对目标的理解程度也很关键,有时一些细节会对误差函数的设计起到关键的作用,比如 MSE 和 MAE 是相似的误差函数,如何选择取决于任务目标,如果需要避免较大误差则应选择 MSE。无论如何,任务都应该是设计误差函数最先考虑的东西。

替代误差函数(Surrogate loss function)

有些情况下根据问题目标得到的误差函数很难使用梯度下降,例如图 5-1 中的误差函数在所有可导处导数都是 0(flat plateau),意味着梯度下降无法工作。这时可以使用一个如图 5-2 所示的近似函数进行模拟,通过求替代误差函数的最小值来实现优化原来的误差函数的目的。当然这只是一个理想化的例子,并且建立在了解误差函数的形态的基础上,但是实际中这种信息通常不容易得到,因此替代误差函数在实际中如何应用是一个比较复杂的问题。

拓展:(该部分不是本文章主题内容):基于神经网络的误差函数设计

下面通过介绍人脸识别领域常用的误差函数来展示如何结合实际的场景来设计误差函数。

模型设计

传统的面部识别技术的基本思想是使用人工设计的一系列指标,如双眼的距离,鼻尖到嘴的距离等作为标识不同个体的特征。使用相同的思路,基于机器学习的方法通过神经网络学习类似的特征来分离不同的个体。面部图像的特征提取可使用典型的多层卷积神经网络(convolutional neural network)来实现,经过多层卷积运算得到的大量图像特征再通过一层或多层的全连接运算汇聚成最终的代表样本的特征向量。

误差函数设计

人脸识别的基本目标是判断输入图像的类型(不同个体属于的不同类型,每个个体有多个不同的样本),这里误差函数的作用就是引导模型提取能够区别不同类型的特征向量。如何设计误差函数来引导模型选取合适的关键特征呢?

Contrastive loss对比误差

最直接的方法就是根据不同样本特征向量的"距离"做出判断,当两个样本之间的特征距离足够近的时候,就认为是同一个体,反之就是不同的个体。这就是对比误差函数(contrastive loss)。

Triplet loss

为了解决 contrastive loss 的局限性问题,Triplet loss 被提出,如下式所示:

Large margin softmax loss(L-Softmax)

Contrastive loss 和 triplet loss 都是以欧式距离衡量样本间的差异,Large Margin softmax loss(L-Softmax) 则提出了另外一种思路,那就是将欧式距离转化为极坐标下的角度差,同类样本角度差较小而异类样本角度差较大。

传统的 softmax loss 虽然能够分离不同类型,但是由于它能放大输入(特征向量)的差异,使较小的特征向量的差异就能产生较大的最终差异而使得训练误差迅速减小,这样就难以学习到类型区分度大的特征。为此 L-Softmax 引入了角度差距超参数 m,使用 ‖W1‖ ‖x‖ cos(mθ1) > ‖W2‖ ‖x‖ cos(θ2) 强制产生类间间隔(decision margin)来分离不同类型(如图 9),引导模型学习类间距离大,类内距离小的特征。

2.3基于梯度下降的优化问题

数学意义上的优化问题一般有两类解法,

一个是解析方法(analytical optimization),适用于存在解析解(closed-form solution);另一种是迭代优化(iterative optimization)方法用于不存在解析解的情况。

梯度下降(Gradient Descent)属于迭代优化的一种典型方法。

现实中的机器学习模型中由于使用非线性激活函数(nolinear activation)所以无法直接计算导数为零时误差函数的解,因此只能使用基于梯度下降的优化方法,如图 1 所示对于一个可导的凸函数,从任意一点出发,沿着倒数下降的方向前进直到倒数为零的点,就是函数的最小值。

2.3.1可导性

虽然从数学原理上梯度下降要求误差函数连续可导,但在实践中误差函数允许存在不可导的点,这是因为随机梯度下降算法(Stochastic Gradient Descent,SGD)使用小批量数据的误差均值进行求导,这样使得落在不可导的点上的概率很低。因此可以对 Hinge loss 这样的不连续的函数使用梯度下降来进行优化。事实上即使在某组数据真的发生小概率事件导致求导失效,由于小批量梯度下降算法使用了大量的分组,绝大多数可求导的分组仍然可以保证在整个数据集上有效运行。【需要修改,我没读懂】

2.3.2非凸性

对于所有的凸函数,使用梯度下降都可以找到最小值,但是在实际的机器学习任务中,由于模型参数的数量很大(如 VGG16 ),这时的误差函数是凸函数的概率非常低,函数表面会呈现出非常复杂的形态。图 2 展示了某个二元模型的误差函数的形态,可见其中有很多区域导数为零,但显然他们并不都是最小值,甚至不是局部最小值。数学上将导数为零的点成为关键点,误差函数表面主要有以下几类关键点:

图 2. 误差函数表面上的关键点

水平区域(flat plateau),是指在一个其中所有的点的导数都为 0 的区域,在三维空间中的水平区域就是水平面, 如何避免陷入水平区域一直是一个讨论较多的话题,目前常见的方法有:

适当的参数初始化使用替代误差函数(surrogate function)调整优化器(如增加动量 momentum)

鞍点(saddle point),是指一阶导数为零二阶导数符号相反的点,如图 2 中所示。由于在高维度误差函数的所有导数为 0 的点中,只有在所有维度具有相同的凹凸性(即二阶导数都大于 0 或小于 0)的时候误差函数才会处于极值(局部最小值或者全局最小值),任何一个维度的凹凸性不同于其他的维度都会使误差函数处于鞍点,考虑到实际的误差函数通常会有万或十万(甚至百万)级的维度数量,因此鞍点是非常普遍的。

大量的试验表明,使用随机梯度下降 SGD 选择具有动量或可变学习率(adaptive learning rate)的优化器(如 Adagrad, Rmsprop)就可以有效的脱离鞍点,如图 3 所示。图三

局部最小值(local minimum), 通过大量实践发现在高维度的优化问题中,局部最小值和全局最小值通常没有太大的区别,甚至在有些情况下比全局最小值有更好的归纳能力(泛化能力)。

总结来说,使用随机梯度下降算法 SGD,结合合适的优化器(例如 momentum) 以及合适的参数初始化可以有效找到非凸函数的最优解(接近最优解)。

2.3.4 泛化(Generalization)

机器学习的最终目的是提高对未知的输入进行预测的准确率,也叫泛化能力。误差函数虽然可以引导梯度下降算法进行模型的优化,但是一个常见的问题是模型虽然达到了很高的训练准确率,但是泛化能力并没有提高甚至反而降低,这就是过拟合(overfitting)。图 4 中,左中右三个小图分别表示了欠拟合,正常和过拟合的情况。右侧的过拟合模型具有三个模型中最复杂的函数形态,可以看到随着训练准确率的提升,这个复杂函数对陌生数据的预测准确率反而下降了。这种问题源自于模型为了提高训练准确率学习了训练数据中的噪声从而导致模型和真实规律产生偏差。

正则化(regularization)是解决过拟合问题的常见方法之一,它把原误差函数和参数的模(norm)相加形成新的目标函数,再使用梯度下降对这个目标函数求最小值。这样做的目的在于通过降低参数维度来增加模型的泛化能力,由于参数的模变成了目标函数的一部分,因此梯度下降也会尽力降低参数的模,而参数的模和参数的维度正相关,因此梯度下降会降低参数的维度,从而避免生成过度复杂的函数,最终增强泛化能力。

2.3.5 正则化

参考:/p/53056358

目标函数

其中M是多项式的阶数

误差函数

拟合数据的目的是最小化误差函数,因为误差函数是多项式系数W的二次函数,因此它关于系数 W的导数是线性函数,所以误差函数的最小值有一个唯一解,我们记做W*,下面我们通过选择M阶数=1,3,9,不同的阶来对数据进行拟合。阶数过低或者过高,导致过拟合和欠拟合,都不能代表我们的目标函数,并且对新数据不具备很好的泛化能力。

第一:值得注意的是当我们一个模型的复杂度给定了之后,数据规模增加能够有效的减轻模型的过拟合问题,所以这也是一个防止模型过拟合的方式。

第二:除了通过增加数据量的方式来减轻过拟合的影响,还可以通过正则化的方式来实现,这种技术就是在我们定义为误差函数增加一个惩罚项,是的多项式系数被有效的控制。均方误差函数修改后的形式为:

目标函数是最终需要优化的函数,其中包括经验损失和结构损失。

obj=loss+Ω

经验损失(loss)就是传说中的损失函数或者代价函数。结构损失(Ω)就是正则项之类的来控制模型复杂程度的函数。

其中,这种二次正则项的应用也叫作山脊回归(Ridge Regression),下面我们同样展示当 多项式的阶M-9时,取λ=0.001和λ=1时过拟合被压制的情况。如下图:

当λ 过大的时候会过度抑制模型,所以根据模型复杂度来选择合适的λ对于模型最终的拟合结果有着重要的影响。

在多项式拟合曲线的过程中,**如果当我们模型复杂度被限制了,我们可以通过增加训练数据的方式防止模型的过拟合以及减轻被噪声数据的影响。同样当我们数据有限,想通过分析数据来选择模型的复杂度的话也可以通过正则化的方式来抑制模型的过拟合。**当然也可以同时使用两种方式,合理的调整才能达到最佳拟合。

fig = plt.figure(figsize=(12,4))for i,lamb in enumerate([0.001,1]):X_train,y_train = create_data(10)plt.subplot(1,2,i+1)poly = PolynomialFeatures(9)X_train_ploy = poly.fit_transform(X_train)X_test_ploy = poly.fit_transform(X_test)lr = Ridge(alpha=(lamb/2))lr.fit(X_train_ploy,y_train)y_pred = lr.predict(X_test_ploy)plt.scatter(X_train,y_train,facecolor="none", edgecolor="b", s=50, label="training data")plt.plot(X_test,y_pred,c="r",label="fitting")plt.plot(X_test,y_test,c="g",label="$\sin(2\pi x)$")plt.title("$\lambda$={}".format(lamb))plt.legend()plt.show()

3.0 最小二乘法—多项式拟合非线性函数

基于均方误差最小化来进行模型求解的方法称为“最小二乘法”——机器学习。均方误差(MSE)是一种加权最小二乘,它的权值是概率。

参考:/p/af0a4f71c05a

流程:

1、函数的近似表示—高次多项式

2、误差函数—最小二乘法

3、引出案例函数曲线

4、目标函数

5、优化目标函数

6、优化目标函数—梯度下降法

7、优化目标函数—求解线性方程组

8、python编程实现拟合曲线函数

9、结果分析

3.1 函数的近似表示—高次多项式

为了研究一些复杂的函数,我们希望对函数自变量进行有限的加、减、乘法三种算数运算,便可以求出原函数值,因此我们通常使用高次多项式来近似表示函数

可以看到,我们可以用n阶多项式的系数线性组合来近似表示f(x)

特别说明,由泰勒展开式可知,当pn(x)的各阶导数和f(x)的各阶导数都相等,则f(x)和pn(x)的误差只是(x-x0)^n的高阶无穷小

3.2误差函数—最小二乘法MSE

我们用f(x)的真实函数值减去多项式函数的结果的平方,来表示f(x)和多项式函数的误差关系,即

我们令x0=0,则最小二乘法表示的误差为

3.3引出案例函数曲线

下面我们用一个案例来介绍最小二乘法拟合非线性函数的算法步骤

案例:求一个多项式来拟合下列函数的函数值

其中我们加入了噪点数据noise,使得函数在定义域内随机的震荡

我们画出f(x)在[-1, 1]之间的图像,即

案例问题即为:已知上述N个数据点,求其函数f(x)的表达式?

生成带有噪点的待拟合的数据集合

3.4目标函数

显然,f(x)也就是机器学习要学习的目标

机器学习的目标为:

(1) 学习一个f(x)多项式,可以拟合真实数据的变化趋势

(2)f(x)的目标:使每一个真实数值到f(x)的拟合数值的距离之和最小

翻译为数学语言,即

计算最小二乘法当前的误差:

3.5 优化目标函数

为了使目标函数最优,我们对每个系数求其偏导数为0,即

至此,我们需要根据上述k的等式,求出所有的系数ak,有两种方法可以求解

1)梯度下降法

(2)求解线性方程

3.5.1 有两种方法可以求解

1)梯度下降法

(2)求解线性方程

3.5.2 优化目标函数—梯度下降法

采用梯度下降法,可以避免我们用数学方法直接一步来求解上述k个方程的最优解,而是采用迭代逼近最优解的思路,其具体步骤为:

(1)初始化k个系数值,开始迭代

(2)每次迭代,分别求出各个系数ak对应的梯度值

(3)用梯度值和学习率来更新每一个系数ak

(4)保证每次更新完所有系数,对应的损失函数值都在减少(说明梯度一直在下降)

翻译为数学语言为:

计算ak对应的梯度值

更新ak

迭代解法:最小二乘法+梯度下降法

3.5.3 优化目标函数—求解线性方程组

求解线性方程组让我们直接一步求出偏导数最优解,其具体步骤为:

(1)将最优化偏导数的方程组写为矩阵乘法形式:XA=Y

(2)利用数学消元算法来求解XA=Y

我们对k个偏导数=0的等式进行化简处理,即

下面我们将其写为矩阵相乘的形式,将所有只含有x的系数写为第一个矩阵X,即

将只含有待求解的多项式系数ak写为第二个矩阵A,即

最后将含有y的系数写为第三个矩阵Y,即

此时,我们就可以用矩阵的乘法来表示最初的k个偏导数为0的等式,即

注意:其中X是kk的矩阵,A是k1的矩阵,Y是k*1的矩阵,符合矩阵的乘法规则

从矩阵表达式可以看出,X和Y矩阵是已知的,A矩阵只包含待求解的f(x)的所有多项式系数,则问题即转化为线性方程组:XA=Y

求解线性方程组,大概有如下算法

(1)高斯消元法(全主元、列主消元)

(2)LU三角分解法

(3)雅克比迭代法

(4)逐次超松弛(SOR)迭代法

(5)高斯-赛德尔迭代法

3.5.3.1 求解线性方程组,大概有如下算法

(1)高斯消元法(全主元、列主消元)

(2)LU三角分解法

(3)雅克比迭代法

(4)逐次超松弛(SOR)迭代法

(5)高斯-赛德尔迭代法

数学解法:最小二乘法+求解线性方程组(我们采用高斯消元法来求解A)

备注:高斯消元法是线性代数中的一种算法,可用于求解线性方程组的解,以及求出矩阵的秩,和求出可逆矩阵的逆矩阵。

高斯消元法的原理是:通过用初等行变换将增广矩阵化为行阶梯阵,然后通过回带求解线性方程组的解。

参考:/u010182633/article/details/45225179

比较上述两种解法,以及多项式拟合的思想,我们可以总结出:

(1)利用函数的近似多项式表示,和最小二乘法来表示目标函数,我们可以高效的拟合非线性函数

(2)梯度下降法可以一步步迭代逼近目标函数的极值,而不用直接求出偏导数最优解

(3)求解方程组,使用数学消元的算法,直接一步优化完成了目标函数的偏导数最优解

(4)从结果上比较,求解方程组的效率和拟合结果要比梯度下降法好

参考:/p/af0a4f71c05a

附录:python科学计算中的常用方法-曲线拟合

最小二乘法是一种数学优化技术,它通过最小化误差的平方和寻找数据的最佳函数匹配。优化是找到最小值或等式的数值解的问题。而线性回归就是要求样本回归函数尽可能好地拟合目标函数值,也就是说,这条直线应该尽可能的处于样本数据的中心位置。因此,选择最佳拟合曲线的标准可以确定为:使总的拟合误差(即总残差)达到最小

在机器学习或者是科学计算项目中,有时可以涉及到一些趋势预测的工作,或者是曲线拟合的工作,现在想整理一下这些基础的方法。

方法会基于Python,以下是github地址

/PeterJaq/curve-fitting

然后就是定义调用函数了,为了方便后续的使用,我们定义函数Fitting。

def Fitting(self, model="line"):info = []if model is "line":A1, B1 = optimize.curve_fit(self.f_1, self.x_0, self.y_0)[0]info = [A1, B1]y_1 = A1 * self.x_1 + B1elif model is "square":A1, B1, C1 = optimize.curve_fit(self.f_2, self.x_0, self.y_0)[0]info = [A1, B1, C1]y_1 = A1 * self.x_1 * self.x_1 + B1*self.x_1 + C1elif model is "cube":A1, B1, C1, D1 = optimize.curve_fit(self.f_3, self.x_0, self.y_0)[0]info = [A1, B1, C1, D1]y_1 = A1 * self.x_1 * self.x_1 * self.x_1 + B1 * self.x_1 * self.x_1 + C1* self.x_1 + D1elif model is "gauss":A1, B1, C1, sigma = optimize.curve_fit(self.f_gauss, self.x_0, self.y_0)[0]info = [A1, B1, C1, sigma]y_1 = A1 * np.exp(-(self.x_1 - B1) ** 2 / (2 * sigma ** 2)) + C1elif model is "ln":A1, B1 = optimize.curve_fit(self.f_ln, self.x_0, self.y_0)[0]info = [A1, B1]y_1 = A1 * np.log(self.x_1) + B1return y_1, info

参考:/p/72241280

参考:/p/53056358

参考:/p/34115676

参考:/p/af0a4f71c05a

参考:/p/50347696

参考:/p/53977691

plotly库

bug1:plotly.offline.iplot(data)和plotly.offline.plot(data)。如果我cmd状态下运行下面的代码,它只是在浏览器中打开。如果我将最后一行更改为pyo.iplot而不是plot,则它不会执行任何操作。plotly更适合jupyter 。

import plotly.offline as pyo from plotly.graph_objs import Figure, Data pyo.offline.init_notebook_mode() x = list(range(10)) y = [22,23,19,26,24,23,29,30,33,29] trace = {'type':'scatter', 'x' : x, 'y' : y, 'name':'simple trace', 'mode':'lines'} layout = {'title':'basic line', 'xaxis':{'title':'x'}, 'yaxis':{'title':'y'}} data = Data([trace]) fig = Figure(data=data, layout=layout) pyo.plot(fig)

附录测试代码

'''最小二乘法拟合函数曲线f(x)1、拟合多项式为:y = a0 + a1*x + a2*x^2 + ... + ak*x^k2、求每个点到曲线的距离之和:Loss = ∑(yi - (a0 + a1*x + a2*x^2 + ... + ak*x^k))^23、最优化Loss函数,即求Loss对a0,a1,...ak的偏导数为03.1、数学解法——求解线性方程组:整理最优化的偏导数矩阵为:X:含有xi的系数矩阵,A:含有ai的系数矩阵,Y:含有yi的系数矩阵求解:XA=Y中的A矩阵3.2、迭代解法——梯度下降法:计算每个系数矩阵A[k]的梯度,并迭代更新A[k]的梯度A[k] = A[k] - (learn_rate * gradient)'''import numpy as npimport matplotlib.pyplot as pltplt.rcParams['font.sans-serif'] = ['SimHei']plt.rcParams['axes.unicode_minus'] = False'''高斯列主消元算法'''# 得到增广矩阵def get_augmented_matrix(matrix, b):row, col = np.shape(matrix)matrix = np.insert(matrix, col, values=b, axis=1)return matrix# 取出增广矩阵的系数矩阵(第一列到倒数第二列)def get_matrix(a_matrix):return a_matrix[:, :a_matrix.shape[1] - 1]# 选列主元,在第k行后的矩阵里,找出最大值和其对应的行号和列号def get_pos_j_max(matrix, k):max_v = np.max(matrix[k:, :])pos = np.argwhere(matrix == max_v)i, _ = pos[0]return i, max_v# 矩阵的第k行后,行交换def exchange_row(matrix, r1, r2, k):matrix[[r1, r2], k:] = matrix[[r2, r1], k:]return matrix# 消元计算(初等变化)def elimination(matrix, k):row, col = np.shape(matrix)for i in range(k + 1, row):m_ik = matrix[i][k] / matrix[k][k]matrix[i] = -m_ik * matrix[k] + matrix[i]return matrix# 回代求解def backToSolve(a_matrix):matrix = a_matrix[:, :a_matrix.shape[1] - 1] # 得到系数矩阵b_matrix = a_matrix[:, -1] # 得到值矩阵row, col = np.shape(matrix)x = [None] * col # 待求解空间X# 先计算上三角矩阵对应的最后一个分量x[-1] = b_matrix[col - 1] / matrix[col - 1][col - 1]# 从倒数第二行开始回代x分量for _ in range(col - 1, 0, -1):i = _ - 1sij = 0xidx = len(x) - 1for j in range(col - 1, i, -1):sij += matrix[i][j] * x[xidx]xidx -= 1x[xidx] = (b_matrix[i] - sij) / matrix[i][i]return x# 求解非齐次线性方程组:ax=bdef solve_NLQ(a, b):a_matrix = get_augmented_matrix(a, b)for k in range(len(a_matrix) - 1):# 选列主元max_i, max_v = get_pos_j_max(get_matrix(a_matrix), k=k)# 如果A[ik][k]=0,则矩阵奇异退出if a_matrix[max_i][k] == 0:print('矩阵A奇异')return None, []if max_i != k:a_matrix = exchange_row(a_matrix, k, max_i, k=k)# 消元计算a_matrix = elimination(a_matrix, k=k)# 回代求解X = backToSolve(a_matrix)return a_matrix, X'''最小二乘法多项式拟合曲线'''# 生成带有噪点的待拟合的数据集合def init_fx_data():# 待拟合曲线f(x) = sin2x * [(x^2 - 1)^3 + 0.5]xs = np.arange(-1, 1, 0.01) # 200个点ys = [((x ** 2 - 1) ** 3 + 0.5) * np.sin(x * 2) for x in xs]ys1 = []for i in range(len(ys)):z = np.random.randint(low=-10, high=10) / 100 # 加入噪点ys1.append(ys[i] + z)return xs, ys1# 计算最小二乘法当前的误差def last_square_current_loss(xs, ys, A):error = 0.0for i in range(len(xs)):y1 = 0.0for k in range(len(A)):y1 += A[k] * xs[i] ** kerror += (ys[i] - y1) ** 2return error# 迭代解法:最小二乘法+梯度下降法def last_square_fit_curve_Gradient(xs, ys, order, iternum=1000, learn_rate=0.001):A = [0.0] * (order + 1)for r in range(iternum + 1):for k in range(len(A)):gradient = 0.0for i in range(len(xs)):y1 = 0.0for j in range(len(A)):y1 += A[j] * xs[i]**jgradient += -2 * (ys[i] - y1) * xs[i]**k # 计算A[k]的梯度A[k] = A[k] - (learn_rate * gradient) # 更新A[k]的梯度# 检查误差变化if r % 100 == 0:error = last_square_current_loss(xs=xs, ys=ys, A=A)print('最小二乘法+梯度下降法:第{}次迭代,误差下降为:{}'.format(r, error))return A# 数学解法:最小二乘法+求解线性方程组def last_square_fit_curve_Gauss(xs, ys, order):X, Y = [], []# 求解偏导数矩阵里,含有xi的系数矩阵Xfor i in range(0, order + 1):X_line = []for j in range(0, order + 1):sum_xi = 0.0for xi in xs:sum_xi += xi ** (j + i)X_line.append(sum_xi)X.append(X_line)# 求解偏导数矩阵里,含有yi的系数矩阵Yfor i in range(0, order + 1):Y_line = 0.0for j in range(0, order + 1):sum_xi_yi = 0.0for k in range(len(xs)):sum_xi_yi += (xs[k] ** i * ys[k])Y_line = sum_xi_yiY.append(Y_line)a_matrix, A = solve_NLQ(np.array(X), Y) # 高斯消元:求解XA=Y的A# A = np.linalg.solve(np.array(X), np.array(Y)) # numpy API 求解XA=Y的Aerror = last_square_current_loss(xs=xs, ys=ys, A=A)print('最小二乘法+求解线性方程组,误差下降为:{}'.format(error))return A# 可视化多项式曲线拟合结果def draw_fit_curve(xs, ys, A, order):fig = plt.figure()ax = fig.add_subplot(111)fit_xs, fit_ys = np.arange(min(xs) * 0.8, max(xs) * 0.8, 0.01), []for i in range(0, len(fit_xs)):y = 0.0for k in range(0, order + 1):y += (A[k] * fit_xs[i] ** k)fit_ys.append(y)ax.plot(fit_xs, fit_ys, color='g', linestyle='-', marker='', label='多项式拟合曲线')ax.plot(xs, ys, color='m', linestyle='', marker='.', label='曲线真实数据')#matplotlib.pyplot.title(label, fontdict=None, loc='center', pad=None, **kwargs)plt.title('最小二乘法拟合多项式N={}的函数曲线f(x)'.format(order))plt.legend()plt.show()if __name__ == '__main__':order = 10 # 拟合的多项式项数xs, ys = init_fx_data() # 曲线数据# 数学解法:最小二乘法+求解线性方程组#A = last_square_fit_curve_Gauss(xs=xs, ys=ys, order=order)# 迭代解法:最小二乘法+梯度下降A = last_square_fit_curve_Gradient(xs=xs, ys=ys, order=order, iternum=10000, learn_rate=0.001)draw_fit_curve(xs=xs, ys=ys, A=A, order=order) # 可视化多项式曲线拟合结果

未完成项目:三种梯度下降方式的理解和优劣,梯度下降(Gradient Descent)小结,迭代

之后更新

本内容不代表本网观点和政治立场,如有侵犯你的权益请联系我们处理。
网友评论
网友评论仅供其表达个人看法,并不表明网站立场。