Вычисление переменного взаимодействия (точечный продукт векторов в матрице)

Drey спросил: 12 мая 2018 в 05:02 в: python

Если я умножу вектор x (1, n) на себя, то есть np.dot(x.T, x), я получу матрицу в квадратичной форме.

Если I имеют матрицу Xmat (k, n), как я могу эффективно вычислить продукт роллинга и выбрать только верхние треугольные элементы?

Итак, атм. У меня есть следующее решение:

def compute_interaction(x):
    xx = np.reshape(x, (1, x.size))
    return np.concatenate((x, np.dot(xx.T, xx)[np.triu_indices(xx.size)]))

Затем compute_interaction(np.asarray([2,5])) yield array([ 2, 5, 4, 10, 25]).

И когда у меня есть матрица Я использую

np.apply_along_axis(compute_interaction, axis=1, arr = np.asarray([[2,5], [3,4], [8,9]]))

, который дает то, что я хочу:

array([[ 2,  5,  4, 10, 25],
       [ 3,  4,  9, 12, 16],
       [ 8,  9, 64, 72, 81]])

Есть ли другой способ, чем вычислить это, используя apply_along_axis? Возможно, используя np.einsum?

Divakar ответил: 12 мая 2018 в 05:15

Каковы типичные значения n и k в вашем фактическом прецеденте?

Drey ответил: 12 мая 2018 в 05:18

Скажем, k будет в тысячах или миллионах и n до 100.

2 ответа

Есть решение

Divakar ответил: 12 мая 2018 в 05:45

Подход № 1

Одно из решений с np.triu_indices будет -

r,c = np.triu_indices(arr.shape[1])
out = np.concatenate((arr,arr[:,r]*arr[:,c]),axis=1)

Подход # 2

Быстрее с slicing -

def pairwise_col_mult(a):
    n = a.shape[1]
    N = n*(n+1)//2
    idx = n + np.concatenate(( [0], np.arange(n,0,-1).cumsum() ))
    start, stop = idx[:-1], idx[1:]
    out = np.empty((a.shape[0],n+N),dtype=a.dtype)
    out[:,:n] = a
    for j,i in enumerate(range(n)):
        out[:,start[j]:stop[j]] = a[:,[i]] * a[:,i:]
    return out

Сроки -

In [254]: arr = np.random.randint(0,9,(10000,100))In [255]: %%timeit
     ...: r,c = np.triu_indices(arr.shape[1])
     ...: out = np.concatenate((arr,arr[:,r]*arr[:,c]),axis=1)
1 loop, best of 3: 577 ms per loopIn [256]: %timeit pairwise_col_mult(arr)
1 loop, best of 3: 233 ms per loop

Drey ответил: 12 мая 2018 в 05:32

Это очень лаконично, спасибо! Еще разобраться, почему это работает :-)

Drey ответил: 12 мая 2018 в 04:56

# 1 уже в 100 раз быстрее: D, то решение, которое я использую

hpaulj ответил: 12 мая 2018 в 09:40

In [165]: arr = np.asarray([[2,5], [3,4], [8,9]])
In [166]: arr
Out[166]: 
array([[2, 5],
       [3, 4],
       [8, 9]])
In [167]: compute_interaction(arr[0])
Out[167]: array([ 2,  5,  4, 10, 25])

Для чего стоит apply_along_axis:

In [168]: np.array([compute_interaction(row) for row in arr])
Out[168]: 
array([[ 2,  5,  4, 10, 25],
       [ 3,  4,  9, 12, 16],
       [ 8,  9, 64, 72, 81]])

apply... - это просто удобный инструмент, чтобы сделать итерацию по нескольким осям более четкой (но не быстрее).