numpy Python -在给定函数值的情况下,查找2D高斯函数的x和y值

ix0qys7i  于 2022-12-04  发布在  Python
关注(0)|答案(2)|浏览(176)

我有一个2D高斯函数f(x,y),我知道函数的峰值g₀出现时的值x₀y₀,但我想找到f(xₑ, yₑ) = g₀ / e¹出现时的值xₑyₑ,我知道这个问题有多个解,但至少有一个解就足够了。
到目前为止我已经

def f(x, y, g0,x0,y0,sigma_x,sigma_y,offset):
    return offset + g0* np.exp(-(((x-x0)**(2)/(2*sigma_x**(2))) + ((y-y0)**(2)/(2*sigma_y**(2)))))
  • 作为参数的所有变量都是已知的 *,因为它们是从曲线拟合中提取的。

我知道在x中取导数并设置f() = 0,类似地在y中,给出了(x,y)的可解线性系统,但手动实现这似乎有点大材小用,必须有一些库或工具可以做我试图实现的事情?

uurv41yg

uurv41yg1#

您可以使用scipy.optimize.fsolve函数来查找满足方程f(x,y)= g/ e¹的(x,y)值。此函数使用数值求根算法来查找非线性方程组的根。
下面是一个例子,说明如何使用scipy.optimize.fsolve来查找满足等式f(x,y)= g/ e¹的(x,y)值:

from scipy.optimize import fsolve

# Define the function f(x, y)
def f(x, y, g0, x0, y0, sigma_x, sigma_y, offset):
    return offset + g0 * np.exp(-(((x-x0)**(2)/(2*sigma_x**(2))) + ((y-y0)**(2)/(2*sigma_y**(2)))))

# Define the function that we want to find the roots of
def g(xy, g0, x0, y0, sigma_x, sigma_y, offset):
    x, y = xy
    return f(x, y, g0, x0, y0, sigma_x, sigma_y, offset) - g0 / np.exp(1)

# Define the initial guess for the root
x0_guess = x0
y0_guess = y0

# Find the roots of the function g(x, y)
x, y = fsolve(g, (x0_guess, y0_guess), args=(g0, x0, y0, sigma_x, sigma_y, offset))

# Print the result
print(f"x = {x}, y = {y}")

在此示例中,fsolve将使用初始猜测值(x0_guess,y0_guess)作为起点,并迭代尝试查找函数g(x,y)的根。如果函数g(x,y)有多个根,则fsolve将仅返回其中一个根(最接近初始猜测值的根)。

nwwlzxa7

nwwlzxa72#

有无穷多种可能性(在g0的特殊情况下可能是1个平凡的或没有)。一个解可以用直接方法在常数时间内解析计算。不需要近似或迭代方法来求给定函数的根。它只是纯粹的数学。
高斯核具有有趣的对称性,其中之一是当峰值被平移到(0,0)时对旋转的不变性,另一个是2D高斯曲面的1D截面是高斯曲线。
让我们暂时忽略offset:它实际上并不改变问题(它只是Z轴平移),并且为分辨率添加了额外的无用项。
该问题的几何解是椭圆,因此解(xe, ye)遵循二次曲线表达式(xe-x0)² / a² + (ye-y0)² / b² = 1。如果sigma_x = sigma_y,则解更简单:这是一个表达式为(xe-x0)² + (ye-y0)² = r的圆。请注意,abr取决于搜索值和内核参数(例如sigma_x)。更改sigma_xsigma_y就像拉伸空间,所以解也是一样的,改变x0y0就像平移空间一样,解也是一样的。
实际上,我们可以解决x0=0y0=0sigma_x=1sigma_y=1的简单情况。然后我们可以应用平移,然后使用变换矩阵进行线性变换。基本乘法4x 4矩阵可以做到这一点。解决简单情况要容易得多,因为需要考虑的参数较少。实际上,g0offset也可以部分地从f中去掉,因为它在表达式的两边,我们只需要解线性方程offset + g0 * h(xe,ye) = g0 / e,所以h(x,y) = 1 / e - offset / g0,其中h(xe, ye) = exp(-(xe² + ye²)/2)。假设我们暂时忘记了平移和线性变换,问题可以很容易地解决:
x1米25英寸
exp(-(xe² + ye²)/2) = 1 / e - offset / g0
-(xe² + ye²)/2 = ln(1 / e - offset / g0)
x1米28英寸1x
就是这样!我们得到了半径为r的圆的表达式,-2*ln(1 / e - offset / g0)!注意,表达式中的ln基本上是自然对数。
现在,我们可以尝试找到4x 4矩阵系数,或者实际上尝试直接求解完整表达式,这最终并不那么困难。
x1米32英寸
exp(-((x-x0)²/(2*sigma_x²) + (y-y0)²/(2*sigma_y²))) = 1 / e - offset / g0
x1米34英寸
x1米35英寸
x1米36英寸1x
就是这样!我们得到了圆锥曲线表达式,其中r = -2 * ln(1 / e - offset / g0)是一个常数,a = sigma_xb = sigma_y是上面表达式中的未知参数。它可以用a = sigma_x/sqrt(r)b = sigma_y/sqrt(r)进行归一化,所以右手边是1,完全符合上面的表达式,但这只是一些数学细节。
你可以很容易地找到椭圆的一个点,因为你知道椭圆的中心(x0, y0),并且在直线y=y0和上面的二次曲线表达式的交点上至少有一个点。
(x-x0)²/sigma_x² + (y0-y0)²/sigma_y² = -2 * ln(1 / e - offset / g0)
x1米45英寸
(x-x0)² = -2 * ln(1 / e - offset / g0) * sigma_x²
x = sqrt(-2 * ln(1 / e - offset / g0) * sigma_x²) + x0
注意,这里有两个解(-sqrt(...) + x0),但你只需要其中一个。我希望我没有在计算中犯任何错误(至少细节应该足够容易地找到它),而且在你的情况下,解不是一个复数。这个解的好处是它计算起来非常非常快

最终解决方案是:
一米四十九分一秒

相关问题