主页 / coding / 问题 / 79352276
Accepted
IzaeDA
Asked: 2025-01-13 21:06:00 +0800 CST

对二维数组的下三角应用函数

  • 772

我有一个数组：

U = np.array([3, 5, 7, 9, 11])

我想要得到如下结果：

result = np.array([
    [  np.nan,      np.nan,      np.nan,      np.nan,    np.nan],
    [U[0] - U[1],   np.nan,      np.nan,      np.nan,    np.nan],
    [U[0] - U[2], U[1] - U[2],   np.nan,      np.nan,    np.nan],
    [U[0] - U[3], U[1] - U[3], U[2] - U[3],   np.nan,    np.nan],
    [U[0] - U[4], U[1] - U[4], U[2] - U[4], U[3] - U[4], np.nan]
])

我可以用来np.tril_indices(4, k=-1)获取没有对角线的下三角的索引，但下一步是什么？

python

2 个回答

  1. Best Answer

    一种简单的方法是计算整个差异并选择所需的元素，但这样做会做更多不必要的工作：

    np.where(np.arange(U.size)[:, None] > np.arange(U.size), U[:, None] - U, np.nan)

    这是一个实际上比简单的面具更有用的时候np.where，虽然也可以用面具来做：

    result = np.full((U.size, U.size), np.nan)
index = np.arange(U.size)
mask = index[:, None] > index
result[mask] = [U[:, None] - U][mask]

    更有效的方法可能是更直接地使用索引来索引源：

    result = np.full((U.size, U.size), np.nan)
r, c = np.tril_indices(U.size, k=-1)
result[r, c] = U[c] - U[r]
    • 2

  2. 虽然不是最有效的，但只有几行代码：

    • 通过广播计算完整矩阵u。
    • 将上三角和对角线设为nan。
    import numpy as np

u = np.array([3, 5, 7, 9, 11], dtype=float)

result = u - u[:, np.newaxis]
result[np.triu_indices(len(u), k=0)] = np.nan
print(result)
# [[nan nan nan nan nan]
#  [-2. nan nan nan nan]
#  [-4. -2. nan nan nan]
#  [-6. -4. -2. nan nan]
#  [-8. -6. -4. -2. nan]]

    以下是使用接受的答案作为基线的计时结果- 使用where和使用掩码的版本似乎对于较小的尺寸具有优势; 对于较大的尺寸，差异似乎越来越可以忽略不计：

    计时结果

    用于计时的代码（请注意，我交换了减法的顺序baseline_where()和baseline_mask()以符合OP的问题，同时我还修复了中的拼写错误（(…)[mask]而不是[…][mask]）baseline_mask()：

    import matplotlib.pyplot as plt
import numpy as np
from timeit import Timer

def proposed(U):
    result = U.astype(float) - U[:, np.newaxis]
    result[np.triu_indices(len(U), k=0)] = np.nan
    return result

def baseline_where(U):
    return np.where(np.arange(U.size)[:, None] > np.arange(U.size),
                    U - U[:, None], np.nan)

def baseline_mask(U):
    result = np.full((U.size, U.size), np.nan)
    index = np.arange(U.size)
    mask = index[:, None] > index
    result[mask] = (U - U[:, None])[mask]
    return result

def baseline_direct_indexing(U):
    result = np.full((U.size, U.size), np.nan)
    r, c = np.tril_indices(U.size, k=-1)
    result[r, c] = U[c] - U[r]
    return result

rand = np.random.default_rng(seed=42)
n = 100  # number of timing runs
sizes = np.round(np.logspace(0, 4, num=9)).astype(int)
timings = {}
for size in sizes:
    u = rand.integers(100, size=size)
    for fct in proposed, baseline_where, baseline_mask, baseline_direct_indexing:
        fct_name = fct.__name__
        fct_timings = timings.get(fct_name, [])
        fct_timings.append(Timer(lambda: fct(u)).timeit(n))
        timings[fct_name] = fct_timings

plt.loglog()
plt.xlabel("len(U)")
plt.ylabel(f"calculation time ({n=} runs) (s)")
for fct_name, fct_timings in timings.items():
    plt.plot(sizes, fct_timings, "x", label=fct_name)
plt.legend()
    • 2

相关问题