如何在Python中高效地计算两个高斯分布的热图?

How to efficiently compute the heat map of two Gaussian distribution in Python?(如何在Python中高效地计算两个高斯分布的热图?)
本文介绍了如何在Python中高效地计算两个高斯分布的热图?的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我正在尝试生成一个热图,其中像素值由两个独立的2D高斯分布控制。设它们分别是Kernel1(muX1,muY1,sigmaX1,sigmaY1)和Kernel2(muX2,muY2,sigmaX2,sigmaY2)。更具体地说,每个内核的长度是其标准差的三倍。第一个内核具有sigmaX1=sigmaY1,第二个内核具有sigmaX2<;sigmaY2。两个核函数的协方差矩阵是对角的(X和Y是独立的)。Kernel1通常完全在Kernel2中。

我尝试了以下两种方法,结果都不满意。有人能给我一些建议吗?

方法1:

迭代地图上的所有像素值对(i,j)并计算i(i,j)的值,公式为i(i,j)=P(i,j|Kernel1,Kernel2)=1-(1-P(i,j|Kernel1))*(1-P(i,j|Kernel2))。然后我得到了下面的结果,就流畅度而言,这是很好的。但在我的电脑上运行需要10秒钟,太慢了。

编码:

def genDensityBox(self, height, width, muY1, muX1, muY2, muX2, sigmaK1, sigmaY2, sigmaX2):
    densityBox = np.zeros((height, width))
    for y in range(height):
        for x in range(width):
            densityBox[y, x] += 1. - (1. - multivariateNormal(y, x, muY1, muX1, sigmaK1, sigmaK1)) * (1. - multivariateNormal(y, x, muY2, muX2, sigmaY2, sigmaX2))
    return densityBox

def multivariateNormal(y, x, muY, muX, sigmaY, sigmaX):
    return norm.pdf(y, loc=muY, scale=sigmaY) * norm.pdf(x, loc=muX, scale=sigmaX)

方法2:

分别生成两个内核对应的两个图像,然后使用某个Alpha值将它们混合在一起。每幅图像由两个一维高斯滤波的外积生成。然后我得到了下面的结果,这是非常粗糙的。但这种方法的优点是,由于在两个向量之间使用外积,所以速度非常快。

由于第一种方法速度慢,第二种方法粗糙,所以我正在努力寻找一种新的方法,同时获得良好的光滑度和较低的时间复杂度。有谁能帮帮我吗?

谢谢!

对于第二种方法,可以很容易地生成二维高斯映射,如here:

def gkern(self, sigmaY, sigmaX, yKernelLen, xKernelLen, nsigma=3):
    """Returns a 2D Gaussian kernel array."""
    yInterval = (2*nsigma+1.)/(yKernelLen)
    yRow = np.linspace(-nsigma-yInterval/2.,nsigma+yInterval/2.,yKernelLen + 1)
    kernelY = np.diff(st.norm.cdf(yRow, 0, sigmaY))
    xInterval = (2*nsigma+1.)/(xKernelLen)
    xRow = np.linspace(-nsigma-xInterval/2.,nsigma+xInterval/2.,xKernelLen + 1)
    kernelX = np.diff(st.norm.cdf(xRow, 0, sigmaX))    
    kernelRaw = np.sqrt(np.outer(kernelY, kernelX))
    kernel = kernelRaw / (kernelRaw.sum())
    return kernel

推荐答案

您的方法很好,只是您不应该循环norm.pdf,而是只推入您想要评估内核的所有值,然后将输出重塑为所需的图像形状。

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import multivariate_normal

# create 2 kernels
m1 = (-1,-1)
s1 = np.eye(2)
k1 = multivariate_normal(mean=m1, cov=s1)

m2 = (1,1)
s2 = np.eye(2)
k2 = multivariate_normal(mean=m2, cov=s2)

# create a grid of (x,y) coordinates at which to evaluate the kernels
xlim = (-3, 3)
ylim = (-3, 3)
xres = 100
yres = 100

x = np.linspace(xlim[0], xlim[1], xres)
y = np.linspace(ylim[0], ylim[1], yres)
xx, yy = np.meshgrid(x,y)

# evaluate kernels at grid points
xxyy = np.c_[xx.ravel(), yy.ravel()]
zz = k1.pdf(xxyy) + k2.pdf(xxyy)

# reshape and plot image
img = zz.reshape((xres,yres))
plt.imshow(img); plt.show()

此方法不应花费太长时间:

In [26]: %timeit zz = k1.pdf(xxyy) + k2.pdf(xxyy)
1000 loops, best of 3: 1.16 ms per loop

这篇关于如何在Python中高效地计算两个高斯分布的热图?的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持编程学习网!

本站部分内容来源互联网,如果有图片或者内容侵犯您的权益请联系我们删除!

相关文档推荐

Leetcode 234: Palindrome LinkedList(Leetcode 234:回文链接列表)
How do I read an Excel file directly from Dropbox#39;s API using pandas.read_excel()?(如何使用PANDAS.READ_EXCEL()直接从Dropbox的API读取Excel文件?)
subprocess.Popen tries to write to nonexistent pipe(子进程。打开尝试写入不存在的管道)
I want to realize Popen-code from Windows to Linux:(我想实现从Windows到Linux的POpen-code:)
Reading stdout from a subprocess in real time(实时读取子进程中的标准输出)
How to call type safely on a random file in Python?(如何在Python中安全地调用随机文件上的类型?)