具有可变大小矩阵输入的 Cython 函数答案

【问题标题】：Cython function with variable sized matrix input具有可变大小矩阵输入的 Cython 函数
【发布时间】：2014-04-01 07:30:19
【问题描述】：

我正在尝试将本机 python 函数的一部分转换为 cython 以缩短计算时间。我想为占用时间的循环组件编写一个 cython 函数（正如 ipython lprun 好心告诉我的那样）。然而，这个函数接受可变大小的矩阵..我看不出如何轻松地将它带到静态类型的 cython。

for index1 in range(0,num_products):
    for index2 in range(0,num_products):
        cond_prob = (data[index1] * data[index2]).sum() / max(col_sums[index1], col_sums[index2])
        prox[index1][index2] = cond_prob

这个问题是 num_products 每年都在变化，因此矩阵（数据）大小是可变的。

这里最好的策略是什么？

我应该写两个 C 函数吗？一个使用 memalloc 创建一个特定维度的矩阵，然后一个对创建的矩阵进行循环？
在这种情况下是否有一些花哨的 cython/numpy 魔法可以提供帮助？我可以编写一个 C 函数，在内存中接收可变大小的 Numpy 数组并传递大小吗？

【问题讨论】：

您好，在 cython 中，您也可以执行 num_products = data.shape[0] 或类似的操作，以便您的循环长度适合您。但是，我不太确定您将要改进循环的方式是否会有很大帮助，因为您在 cython 循环中使用了 numpy 函数。您正在尝试做的事情在 numpy 中似乎也非常简单......上面的表达式是您的原始 python 代码吗？因为如果是这样，您可以通过省略 for 循环来优化它，从而使其更快。
几个问题：数组data和prox的shape是什么？它们都是numpy数组吗？ num_products 通常有多大？
感谢 cmets。好点re：矢量化而不是循环遍历numpy数组！数据 = (170 x 800)，因此 prox 为 800 x 800

标签： python numpy cython

【解决方案1】：

Cython 代码（战略上）是静态类型的，但这并不意味着数组必须具有固定大小。在直接 C 中，将多维数组传递给函数可能有点尴尬，但在 Cython 中，您应该能够执行以下操作：

请注意，我从您的 follow-up question. 中获取了函数和变量名称

import numpy as np
cimport numpy as np
cimport cython

@cython.boundscheck(False)
@cython.cdivision(True)
def cooccurance_probability_cy(double[:,:] X):
    cdef int P, i, j, k
    P = X.shape[0]
    cdef double item
    cdef double [:] CS = np.sum(X, axis=1)
    cdef double [:,:] D = np.empty((P, P), dtype=np.float)

    for i in range(P):
        for j in range(P):
            item = 0
            for k in range(P):
                item += X[i,k] * X[j,k]
            D[i,j] = item / max(CS[i], CS[j])
    return D

另一方面，如果你使用正确的函数和一些广播，只使用 Numpy 也应该很快解决这个问题。事实上，由于计算复杂度以矩阵乘法为主，我发现下面的代码比上面的 Cython 代码要快得多（np.inner 使用了高度优化的 BLAS 例程）：

def new(X):
    CS = np.sum(X, axis=1, keepdims=True)
    D = np.inner(X,X) / np.maximum(CS, CS.T)
    return D

【讨论】：

很好的解释：Cython ...我已经用你的 Cython 函数更新了后续问题 (stackoverflow.com/questions/22853837/…)
@sanguineturtle，与此同时，我编辑了我的答案：事实证明，禁用 Cython 函数的边界检查会显着提高性能。我还发现普通的 Numpy 更快，你一定要检查一下:)

【解决方案2】：

您是否尝试过摆脱 numpy 中的 for 循环？

例如，对于等式的第一部分，您可以尝试：

(data[ np.newaxis,:] * data[:,np.newaxis]).sum(2)

如果内存是个问题，您也可以使用 np.einsum() 函数。对于第二部分，如果你还没有尝试过的话，你可能还可以制作一个 numpy 表达式（有点困难）。

【讨论】：

好主意！ ...我明天会试试这个并更新这篇文章