将上三角矩阵转换为对称矩阵的快速方法 [英] Fast way to convert upper triangular matrix into symmetric matrix

问题描述

我有一个np.float64值的上三角矩阵，如下所示:

I have an upper-triangular matrix of np.float64 values, like this:

array([[ 1.,  2.,  3.,  4.],
       [ 0.,  5.,  6.,  7.],
       [ 0.,  0.,  8.,  9.],
       [ 0.,  0.,  0., 10.]])

我想将其转换为相应的对称矩阵，如下所示:

I would like to convert this into the corresponding symmetric matrix, like this:

array([[ 1.,  2.,  3.,  4.],
       [ 2.,  5.,  6.,  7.],
       [ 3.,  6.,  8.,  9.],
       [ 4.,  7.,  9., 10.]])

转换可以就地完成，也可以作为新矩阵进行.我希望它尽快.我该如何快速做到这一点?

The conversion can be done in place, or as a new matrix. I would like it to be as fast as possible. How can I do this quickly?

推荐答案

这是迄今为止我发现的最快的例程，它不使用Cython或Numba之类的JIT.我在计算机上花费约1.6μs的时间来处理4x4阵列(整个100K 4x4阵列列表上的平均时间):

This is the fastest routine I've found so far that doesn't use Cython or a JIT like Numba. I takes about 1.6 μs on my machine to process a 4x4 array (average time over a list of 100K 4x4 arrays):

inds_cache = {}

def upper_triangular_to_symmetric(ut):
    n = ut.shape[0]
    try:
        inds = inds_cache[n]
    except KeyError:
        inds = np.tri(n, k=-1, dtype=np.bool)
        inds_cache[n] = inds
    ut[inds] = ut.T[inds]

以下是我尝试过的其他一些事情，这些事情并没有那么快:

Here are some other things I've tried that are not as fast:

上面的代码，但是没有缓存.每个4x4阵列大约需要8.3μs:

The above code, but without the cache. Takes about 8.3 μs per 4x4 array:

def upper_triangular_to_symmetric(ut):
    n = ut.shape[0]
    inds = np.tri(n, k=-1, dtype=np.bool)
    ut[inds] = ut.T[inds]

一个普通的Python嵌套循环.每个4x4阵列大约需要2.5μs:

A plain Python nested loop. Takes about 2.5 μs per 4x4 array:

def upper_triangular_to_symmetric(ut):
    n = ut.shape[0]
    for r in range(1, n):
        for c in range(r):
            ut[r, c] = ut[c, r]

使用np.triu进行浮点加法.每个4x4阵列大约需要11.9μs:

Floating point addition using np.triu. Takes about 11.9 μs per 4x4 array:

def upper_triangular_to_symmetric(ut):
    ut += np.triu(ut, k=1).T

Numba版本的Python嵌套循环.这是我发现的最快的结果(每个4x4阵列约0.4μs)，并且最终在生产中使用，至少直到我开始遇到Numba的问题并不得不恢复为纯Python版本:

Numba version of Python nested loop. This was the fastest thing I found (about 0.4 μs per 4x4 array), and was what I ended up using in production, at least until I started running into issues with Numba and had to revert back to a pure Python version:

import numba

@numba.njit()
def upper_triangular_to_symmetric(ut):
    n = ut.shape[0]
    for r in range(1, n):
        for c in range(r):
            ut[r, c] = ut[c, r]

Cython版本的Python嵌套循环.我是Cython的新手，因此可能无法完全优化.由于Cython会增加运营开销，因此我想听听Cython和纯Numpy的答案.每个4x4阵列大约需要0.6μs:

Cython version of Python nested loop. I'm new to Cython so this may not be fully optimized. Since Cython adds operational overhead, I'm interested in hearing both Cython and pure-Numpy answers. Takes about 0.6 μs per 4x4 array:

cimport numpy as np
cimport cython

@cython.boundscheck(False)
@cython.wraparound(False)
def upper_triangular_to_symmetric(np.ndarray[np.float64_t, ndim=2] ut):
    cdef int n, r, c
    n = ut.shape[0]
    for r in range(1, n):
        for c in range(r):
            ut[r, c] = ut[c, r]

这篇关于将上三角矩阵转换为对称矩阵的快速方法的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！