Python:寻找两个高斯曲线的交点

12 投票
2 回答
14092 浏览
提问于 2025-04-17 23:25

我有两个高斯曲线的图:

x = np.linspace(-5,9,10000)
plot1=plt.plot(x,mlab.normpdf(x,2.5,1))
plot2=plt.plot(x,mlab.normpdf(x,5,1))

我想找出这两条曲线交叉的点。有没有什么方法可以做到这一点?特别是我想知道它们相遇时的x坐标值。

2 个回答

2

这里有一个完全基于numpy的解决方案,这个方法也可以用在除了高斯曲线以外的其他曲线。

def get_intersection_locations(y1,y2,test=False,x=None): 
    """
    return indices of the intersection point/s.
    """
    idxs=np.argwhere(np.diff(np.sign(y1 - y2))).flatten()
    if test:
        x=range(len(y1)) if x is None else x
        plt.figure(figsize=[2.5,2.5])
        ax=plt.subplot()
        ax.plot(x,y1,color='r',label='line1',alpha=0.5)
        ax.plot(x,y2,color='b',label='line2',alpha=0.5)
        _=[ax.axvline(x[i],color='k') for i in idxs]
        _=[ax.text(x[i],ax.get_ylim()[1],f"{x[i]:1.1f}",ha='center',va='bottom') for i in idxs]
        ax.legend(bbox_to_anchor=[1,1])
        ax.set(xlabel='x',ylabel='density')
    return idxs
# single intersection
x = np.arange(-10, 10, 0.001)
y1=sc.stats.norm.pdf(x,-2,2)
y2=sc.stats.norm.pdf(x,2,3)
get_intersection_locations(y1=y1,y2=y2,x=x,test=True) # returns indice/s array([10173])

在这里输入图片描述

# double intersection
x = np.arange(-10, 10, 0.001)
y1=sc.stats.norm.pdf(x,-2,1)
y2=sc.stats.norm.pdf(x,2,3)
get_intersection_locations(y1=y1,y2=y2,x=x,test=True)

在这里输入图片描述

这个内容是基于对一个类似问题的回答

22

你想找到x的值,使得两个高斯函数的高度相同,也就是它们相交的地方。

你可以通过把这两个高斯函数相等来求解x。最后你会得到一个二次方程,这个方程的系数和高斯函数的均值及方差有关。下面是最终的结果:

import numpy as np

def solve(m1,m2,std1,std2):
  a = 1/(2*std1**2) - 1/(2*std2**2)
  b = m2/(std2**2) - m1/(std1**2)
  c = m1**2 /(2*std1**2) - m2**2 / (2*std2**2) - np.log(std2/std1)
  return np.roots([a,b,c])

m1 = 2.5
std1 = 1.0
m2 = 5.0
std2 = 1.0

result = solve(m1,m2,std1,std2)

输出结果是:

array([ 3.75])

你可以画出找到的交点:

x = np.linspace(-5,9,10000)
plot1=plt.plot(x,mlab.normpdf(x,m1,std1))
plot2=plt.plot(x,mlab.normpdf(x,m2,std2))
plot3=plt.plot(result,mlab.normpdf(result,m1,std1),'o')

绘制的图像将是:

enter image description here

如果你的高斯函数有多个交点,代码也会找到所有的交点(比如m1=2.5,std1=3.0,m2=5.0,std2=1.0):

enter image description here

撰写回答