0

10 多天以来，我一直在努力让这个拟合操作正常工作。我使用 Levenberg-Marquardt 算法编写了一个用于拟合的 C++ 类，并用简单的多项式对其进行了测试，拟合成功；但是对于我的实验所需的真正的实验功能，它不起作用。

我在 Mathematica 上安装了相同的功能，在那里执行它是可以的（但速度很慢，每次拟合 12 秒，这就是我使用 C++ 的原因）。当我适应时，我的程序表现得很疯狂并返回错误的结果。您能否检查我使用的功能并告诉我是否有问题？

以下是我在课堂上使用的拟合过程。再一次，我想提一下它适用于多项式，所以我可能在这里犯了一个非常严重的错误，它只与一些特殊函数冲突......这是我的假设。

我使用的初始值非常接近正确值。

如果您需要任何其他信息，请告诉我。非常感谢任何帮助。

4

2 回答 2

3

好吧，快速浏览一下，你应该有

   parameters -= SamArray<T>::MatrixProduct(SamArray<T>::Inverse(JTJ + combinationMatrix),J); //add the new deltas to the parameters and test them

   instead of 

   parameters += SamArray<T>::MatrixProduct(SamArray<T>::Inverse(JTJ + combinationMatrix),J); //add the new deltas to the parameters and test them

也就是说，当然，假设J是梯度而不是负梯度。
任何牛顿法或准牛顿法（如 Levenberg-Marquardt）都是这种情况。

于 2012-05-29T13:07:46.993 回答
2

与其讨论代码没有给出预期结果的原因，我将讨论速度问题。你的代码中有

parameters += SamArray<T>::MatrixProduct(SamArray<T>::Inverse(JTJ + combinationMatrix),J);

它首先计算一个矩阵逆然后乘以J，我认为它是另一个矩阵。A X = I这是低效的，因为它与先求解，然后计算基本相同A^-1 B，其中A和X是矩阵，并且I是单位矩阵。相反，最好A X = B直接解决 , 。有几种分解方案可以让您做到这一点，例如QR和LU。在 LU 分解中，矩阵A被分解为上U、 和下、L三角形部分，然后分别有效地反转，如下所示

A X = L U X = B
U X = L^-1 B
X = U^-1 L^-1 B

这是作为算法本身的一部分执行的。如果您想要直接反转，那么您将替换B为I，但如您所见，没有必要。查看Levenberg-Marquardt 算法，可能是这样的情况Transpose[J].J + DiagonalMatrix[Diagonal@J]，使用 Mathematica 表示法，可能是正定的，然后Cholesky 分解可用，而且速度明显更快。但是，我没有仔细研究过算法，所以在这种情况下，Cholesky 分解可能不可用。

\begin{edit} 在第 1.5 章第 1.5 章。在图3中，Stewart 讨论了两种算法的相对速度：invert-then-multiply v. LU 分解。两者具有相同的渐近复杂度，O(n^3)但它们的系数不同。具体来说，在求解方程 时A X == B，invert-then-multiply 的运算计数为，5 n^3 /6 + l n^2而 LU 的n^3 /3 + l n^2wheren是行数，B并且l是列数。两者的比值总是大于 1，即使在 时l == n，LU 分解也比 invert-then-multiply 快 30%。 \结束{编辑}

其次，您要为每组参数计算两次卡方：一次是在最初确定它们时，另一次是在您将它们与下一次迭代进行比较时。该值应与先前接受的参数一起缓存。

最后，我不会低估 Mathematica 以合理的速度处理拟合的能力。

于 2012-05-29T15:16:40.610 回答