1

我正在尝试使用以下功能：将 SNR 降低在此处输入图像描述到我的数据中。C1、C2和h是我需要从leastsq的方法中获取的参数。C1 和 C2 很简单，但问题是我的 h(t) 实际上是：基于 B 样条的函数。我想要获得的是该函数内的系数 hj （在我的例子中，有 35 个不同的 hj's ）。该函数是不同基 B-Spline 的总和，每个 B-Spline 的权重不同，系数的数量等于 B-Spline 的节点数。因为我想获得 C1、C2 和 h1..35，所以我执行以下操作：

funcLine = lambda tpl, eix_x: (tpl[0]*np.sin((4*math.pi*np.sum(bsplines_evaluades * np.transpose([tpl[2],tpl[3],tpl[4],tpl[5],tpl[6],tpl[7],tpl[8],tpl[9],tpl[10],tpl[11],tpl[12],tpl[13],tpl[14],tpl[15],tpl[16],tpl[17],tpl[18],tpl[19],tpl[20],tpl[21],tpl[22],tpl[23],tpl[24],tpl[25],tpl[26],tpl[27],tpl[28],tpl[29],tpl[30],tpl[31],tpl[32],tpl[33],tpl[34],tpl[35],tpl[36],tpl[37]]) , axis=0))*eix_x/lambda1) + tpl[1]*np.cos((4*math.pi*np.sum(bsplines_evaluades * np.transpose([tpl[2],tpl[3],tpl[4],tpl[5],tpl[6],tpl[7],tpl[8],tpl[9],tpl[10],tpl[11],tpl[12],tpl[13],tpl[14],tpl[15],tpl[16],tpl[17],tpl[18],tpl[19],tpl[20],tpl[21],tpl[22],tpl[23],tpl[24],tpl[25],tpl[26],tpl[27],tpl[28],tpl[29],tpl[30],tpl[31],tpl[32],tpl[33],tpl[34],tpl[35],tpl[36],tpl[37]]) , axis=0))*eix_x/lambda1))*np.exp(-4*np.power(k, 2)*lambda_big*np.power(eix_x, 2))
func = funcLine
ErrorFunc = lambda tpl, eix_x, ydata: np.power(func(tpl, eix_x) - ydata,2)
tplFinal1, success = leastsq(ErrorFunc, [2, -2, 8.2*np.ones(35)], args=(eix_x, ydata))

tpl(0)=C1, tpl(1)=C2 和 tpl(2..35)=我的系数。bsplines_evalades 是一个矩阵 [35,86000]，其中每一行是每个基本 b 样条的时间函数，所以我用其单独的系数对每一行进行加权，86000 是 eix_x 的长度。ydata(eix_x) 是我想要近似的函数。λ1= 0.1903 ; lambda_big = 2; k=2*pi/λ1。输出是相同的初始参数，不是逻辑。谁能帮我？我也尝试过curvefit，但它不起作用。数据在：http://www.filedropper.com/data_5>http://www.filedropper.com/download_button.png width=127 height=145 border=0/>
http://www.filedropper.com >在线备份存储

编辑 现在的代码是：

lambda1 = 0.1903
k = 2 * math.pi / lambda1
lambda_big = 2
def funcLine(tpl, eix_x):
    C1, C2, h = tpl[0], tpl(1), tpl[2:]
    hsum = np.sum(bsplines_evaluades * h, axis=1)  # weight each
    theta = 4 * np.pi * np.array(hsum) * np.array(eix_x) / lambda1
    return (C1*np.sin(theta)+C2*np.cos(theta))*np.exp(-4*lambda_big*(k*eix_x)**2)  # lambda_big = 2
if len(eix_x) != 0:
    ErrorFunc = lambda tpl, eix_x, ydata: funcLine(tpl, eix_x) - ydata
    param_values = 7.5 * np.ones(37)
    param_values[0] = 2
    param_values(1) = -2
    tplFinal2, success = leastsq(ErrorFunc, param_values, args=(eix_x, ydata))

问题是输出参数不会相对于初始参数发生变化。数据（x_axis,ydata,bsplines_evaluates）：gist.github.com/hect1995/dcd36a4237fe57791d996bd70e7a9fc7 gist.github.com/hect1995/39ae4768ebb32c27f1ddea97e24d96af gist.github.com/hect1995/bdddf02debfc716

4

2 回答 2

0

这里有几件事可能是错误的：我不确定您tpl是否正确索引了数组（如果它有 37 个条目，则索引应该是0:36）。并且您的 errorFunc 可能应该返回残差而不是平方残差。

最后，我认为您的 h-sum 可能不正确：您想在 $N$ 轴上求和，而不是在 $x$ 轴上求和，对吧？

您可以按如下方式整理您的代码，看看它是否有帮助（没有一些数据很难为自己测试）：

def funcLine(tpl, eix_x):
    C1, C2, h = tpl[0], tpl[1], tpl[2:]
    hsum = np.sum(bsplines_evaluades * h, axis=1)
    theta = 4 * np.pi * hsum * eix_x / lambda1
    return (C1 * np.sin(theta) + C2 * np.cos(theta)) * np.exp(-4  *lambda_big *
                 (k * eix_x)**2)

errorFunc = lambda tpl, eix_x, ydata: funcLine(tpl, eix_x) - ydata
tplFinal2, success = leastsq(errorFunc, [2, -2, 8.2*np.ones(35)],
                             args=(eix_x, ydata))
于 2017-04-10T09:35:03.420 回答
0

提供一个足够完整可运行的可读示例总是有帮助的（对您自己和我们都是如此）。您带有lambdas 和很长的行的示例绝对不可读，这使您很容易错过简单的错误。使用 Python 的要点之一是使代码更易于阅读。

将样条系数作为拟合变量很好，但您希望变量的 np.ndarray 恰好是一维的。所以你的参数数组应该是

param_values = 8.2 * np.ones(37)
param_values[0] = 2
param_values[1] = -2.
result_params, success = leastsq(errorFunc, param_values, ....)

使用 curve_fit() 也应该没问题。除此之外，很难提供太多帮助，因为您既没有提供完整的可运行程序（留下许多未定义的术语），也没有提供运行代码的输出或错误消息。

于 2017-04-11T01:05:25.400 回答