任取一个周期（25 个点）的数据，取模型中常数项初值为 \(c=1.0\)，初始步长 \(h=0.1\)，进行下列迭代：
[*]在 \(\) 范围内均匀取 11 个点 \(c_i\)。
[*]在每个点 \(c_i\) 上用 DFT 计算 \(y-c_i\) 的常数项。
[*]以误差绝对值最小的 \(c_i\) 作为新的迭代参数 \(c\)，同时步长 \(h\) 缩小到原来的 1/10 。
重复上述迭代过程，直到误差绝对值最小的 \(c_i\) 满足指定的数度要求。
【Mathematica 代码】
加载数据
ClearAll["Global`*"]
SetDirectory["C:\\"];
tbl = Import["Data.txt", "Table"];
n = Table][], {i, Length}];
y = SetPrecision][], {i, Length}], 100];
常数分量计算函数
calConstant :=
(*=======================================================*)
(*计算 y 中的常数分量。*)
(*=======================================================*)
(*【输入参数】*)
(*y: {Subscript,Subscript,...,Subscript} 原始数据。*)
(*【输出数据】*)
(*常数分量。*)
Module[{ferr, optimal, step, offset, err},
(*-------------------------------------------------------*)
(*初始化。*)
(*-------------------------------------------------------*)
ferr := 1/5 Fourier[];
optimal = SetPrecision;
step = SetPrecision;
(*-------------------------------------------------------*)
(*迭代。*)
(*-------------------------------------------------------*)
While[step > 10^-48,
offset = Subdivide;
err = Abs]], {i, Length}]];
optimal = Part][][]];
step /= 10;
]; (*end while*)
(*-------------------------------------------------------*)
(*返回结果。*)
(*-------------------------------------------------------*)
optimal
] (*end Module *)
计算并保存结果：
result = Table[{n[], calConstant]]}, {i, 25, Length, 25}];
Export["Constant.txt", result, "TSV"]
【结果】
小数点后稳定的数字比原数据多了 1~2 位。
这是小数点后稳定的 33 位数字：0.998953960449742112505647176797204
如果数据中的常数项信息都已基本被提取了出来，
则意味着这个数据并没有提供预测未来某个时间步常数项小数点后 40 位数字所需的信息。
如果假设前提不成立，还有什么方法来提取更多的信息？ 本帖最后由 Jack315 于 2024-4-12 05:59 编辑
原数据减去 0.99895396044974211 再乘以 1e+18 作图看到了与 8# 形状相同的趋势。
因此模型应该修改为：
\(\begin{matrix}y(t)=\frac{p_1t^2+p_2t+p_3}{t^2+q_1t+q_2}+\exp(\beta_0+\beta_1t)\cos(\omega_ct+\varphi)&t=1,2,3,\cdots\end{matrix}\)
据此有：\(y(\infty)=p_1\)
建模数据取自中心左右各四分之一的部分，模型参数：
\(p_1\): 0.998953960449742112505647162231022030432
\(p_2\): 61.92624846613714685707263143501205487642
\(p_3\): 7240.718073573050381480974427123042069732
\(q_1\): 61.99109360180837709884047211271530449650
\(q_2\): 7248.300082131097003177689864784403599322
\(\beta_0\): -51.01143750080362404210707554419013802146
\(\beta_1\): -0.01004441927896885604364847479566826014080
\(\omega_c\): 0.2509257588303796389235689874214851411490
\(\varphi_c\): 0.1300513307914659036485461634070844812088
残差图：
离真相应该不太远了。
估计模型再加上 \(f(t)\cos[\omega(t)+\varphi]\) 这一项可能会有更理想的结果。
【讨论】
如果 \(f(x)\) 是一个衰减指数或幂函数，则 \(y(\infty)=p_1=\) 0.998953960449742112505647162231022030432...
需要经过非常长的时间才能得到这个极限，否则真实的有效数字只有小数点后 26~27 位的样子。
如果 \(f(x)\) 也是一个有理分式函数，则 \(y(\infty)\) 很有可能根本就不存在。
【更新】
非常可能：
\(f_(x)\) 是一个衰减的幂函数。
\(\omega(t)\) 是一个自然对数函数。
\(y(\infty)=p_1\) 由于我要的是振幅衰减到无限小后，整个数列的极限值，所以我是优先考虑尾部的数据的，因为尾部的数据振幅最小，它们最接近极限值。
我把尾部的100个数据喂给文心大模型4.0后，他返回给我的是一段python代码：
https://yiyan.baidu.com/share/FUanbmPf8z
我以为他是随便写的伪代码，但神奇的是，我把y值填进去之后，这段python代码竟然可以运行出结果：
现在的大语言模型都这么厉害了吗？他竟然能够根据我的问题描述，编写出一段能解决这个问题的程序代码，
而且代码的运行结果的第1个数值和我在3楼里给出的预测值39295561几乎一样！ 本帖最后由 Jack315 于 2024-4-12 08:03 编辑
这是残差的数据：
喂点给无所不能的 AI 不知道会出个啥:lol
会不会是类似这样的函数：\(x^{-0.5}\cos(5\ln x)\) 接 12#...
模型简记为：\(f_1(t)+f_2(t)\cdot f_3(t)\)
继续寻找 \(f_1(t)\) 的模型感到有点难度了。
理论上只有完全确定准确的模型形式、
或者三个函数分量用傅立叶变换后在频率上能完全分离，
才可能得到 \(y(\infty)\) 比较准确的估计。
以我自己目前的结果，或许只能得到小数点后约 34 位有效数字的估计：
\(y(\infty)\approx\) 0.9989 5396 0449 7421 1250 5647 1767 9720 39 指数衰减：
https://bbs.emath.ac.cn/data/attachment/forum/202404/02/224223awhjdzdtsdooco9w.png
从图形上看这是一条直线，但查看细节，会有所发现。
首选用傅立叶变换取出幅度值：
amp = Log]][], {i, 25, Length}]]];然后用下列函数建立多项式模型：
buildAmpPolyModel :=
(*=======================================================*)
(*生成并返回多项式模型参数。*)
(*=======================================================*)
(*【输入参数】*)
(*data: {Subscript,Subscript,...} 数据。*)
(*order: 多项式阶数。*)
(*【输出数据】*)
(*建立的模型。*)
Module[{fun, var},
(*-------------------------------------------------------*)
(*初始化。*)
(*-------------------------------------------------------*)
fun = Plus @@ Table" <> ToString] t^i, {i, 0, order}];
var = Table" <> ToString], {i, 0, order}];
(*-------------------------------------------------------*)
(*拟合并返回模型。*)
(*-------------------------------------------------------*)
NonlinearModelFit
] (*end Module *)
mdlAmp = buildAmpPolyModel模型残差如下图所示：
在不知道准确模型的情况下，无法知道这个残差
究竟是数据生成过程本身的特性，还是因为耦合了振荡所造成。
取最后 25 个数据点拟合下列模型：
\(y(t)=\alpha+\exp(\sum_{i=0}^n\beta_it^i)\cos(\omega t+\varphi)\)
当 \(n=1\) 时，得到 \(\alpha\) 95% 置信区间为：
置信下限：0.9989 5396 0449 7421 1250 5647 1767 9720 3929 5351
置信上限：0.9989 5396 0449 7421 1250 5647 1767 9720 3929 5712
小数点后 37 位有效数字。
增加 \(n\) 收效甚微。这是当 \(n=4\) 时的残差图：
不知道还有什么好的办法能获得更理想的结果。 所以你具体在算什么？
也许有更好的算法。而且我有空闲算力，如果问题有趣的话也可以重算一遍。 用时间序列的 AI 来做预测，这个或许还可以再玩一下:D Ickiverar 发表于 2024-4-15 04:43
所以你具体在算什么？
也许有更好的算法。而且我有空闲算力，如果问题有趣的话也可以重算一遍。 ...
我算的是下面这个问题的答案在n趋于无穷大时的极限值：
https://bbs.emath.ac.cn/thread-1819-1-1.html
这个问题早在1963年就已经有人给出0.593这3位数的近似结果了，我（This paper）在2011年发表了9位小数的计算结果，但在2015年被13位小数的计算结果赶超了，如下表所示：
https://bbs.emath.ac.cn/data/attachment/forum/201709/28/223158l4jkz0hdzkvunh1j.png
我正在计算的是极限值的第14位小数及以后的数据，企图在2024年赶超2015年的计算结果。
页:
1
[2]