В numpy вычисляем матрицу, в которой каждая ячейка содержит произведение всех других записей в этой строке

Ответ 1

Если вы готовы терпеть один цикл:

B = np.empty_like(A)
for col in range(A.shape[1]):
    B[:,col] = np.prod(np.delete(A, col, 1), 1)

Это вычисляет то, что вам нужно, по одному столбцу за раз. Это не так эффективно, как теоретически возможно, потому что np.delete() создает копию; если вы много заботитесь о распределении памяти, используйте вместо этого маску:

B = np.empty_like(A)
mask = np.ones(A.shape[1], dtype=bool)
for col in range(A.shape[1]):
    mask[col] = False
    B[:,col] = np.prod(A[:,mask], 1)
    mask[col] = True

Ответ 2

Вариант вашего решения с использованием repeat, использует [:,None].

np.prod(A,axis=1)[:,None]/A

Мой первый удар при обработке 0s:

In [21]: B
array([[ 0.2,  0.4,  0.6],
       [ 0. ,  0.5,  0.5],
       [ 0.6,  0.4,  0.2]])

In [22]: np.prod(B,axis=1)[:,None]/(B+np.where(B==0,1,0))
array([[ 0.24,  0.12,  0.08],
       [ 0.  ,  0.  ,  0.  ],
       [ 0.08,  0.12,  0.24]])

Но, как отметил комментарий; ячейка [0,1] должна быть 0,25.

Это исправляет эту проблему, но теперь имеет проблемы, когда в строке есть несколько 0.

In [30]: I=B==0
In [31]: B1=B+np.where(I,1,0)
In [32]: B2=np.prod(B1,axis=1)[:,None]/B1
In [33]: B3=np.prod(B,axis=1)[:,None]/B1
In [34]: np.where(I,B2,B3)
Out[34]: 
array([[ 0.24,  0.12,  0.08],
       [ 0.25,  0.  ,  0.  ],
       [ 0.08,  0.12,  0.24]])

In [55]: C
array([[ 0.2,  0.4,  0.6],
       [ 0. ,  0.5,  0. ],
       [ 0.6,  0.4,  0.2]])
In [64]: np.where(I,sum1[:,None],sum[:,None])/C1
array([[ 0.24,  0.12,  0.08],
       [ 0.5 ,  0.  ,  0.5 ],
       [ 0.08,  0.12,  0.24]])

Блаз Братанический epsilon подход - лучшее нетеративное решение (пока):

In [74]: np.prod(C+eps,axis=1)[:,None]/(C+eps)

Различное решение, повторяющееся по столбцам:

def paulj(A):
    P = np.ones_like(A)
    for i in range(1,A.shape[1]):
        P *= np.roll(A, i, axis=1)
    return P

In [130]: paulj(A)
array([[ 0.24,  0.12,  0.08],
       [ 0.25,  0.25,  0.25],
       [ 0.08,  0.12,  0.24]])
In [131]: paulj(B)
array([[ 0.24,  0.12,  0.08],
       [ 0.25,  0.  ,  0.  ],
       [ 0.08,  0.12,  0.24]])
In [132]: paulj(C)
array([[ 0.24,  0.12,  0.08],
       [ 0.  ,  0.  ,  0.  ],
       [ 0.08,  0.12,  0.24]])

Я пробовал некоторые тайминги на большой матрице

In [13]: A=np.random.randint(0,100,(1000,1000))*0.01

In [14]: timeit paulj(A)
1 loops, best of 3: 23.2 s per loop

In [15]: timeit blaz(A)
10 loops, best of 3: 80.7 ms per loop

In [16]: timeit zwinck1(A)
1 loops, best of 3: 15.3 s per loop

In [17]: timeit zwinck2(A)
1 loops, best of 3: 65.3 s per loop

Эпсилонное приближение, вероятно, является лучшей скоростью, которую мы можем ожидать, но имеет некоторые проблемы округления. Необходимость итерации по многим колонкам вредит скорости. Я не уверен, почему подход np.prod(A[:,mask], 1) самый медленный.

eeclo fooobar.com/questions/493340/... предложил использовать as_strided. Вот что я думаю, что он имеет в виду (адаптированный из перекрывающегося блочного вопроса, fooobar.com/questions/64269/...)

def strided(A):
    h,w = A.shape
    A2 = np.hstack([A,A])
    x,y = A2.strides
    strides = (y,x,y)
    shape = (w, h, w-1)
    blocks = np.lib.stride_tricks.as_strided(A2[:,1:], shape=shape, strides=strides)
    P = blocks.prod(2).T # faster to prod on last dim
    # alt: shape = (w-1, h, w), and P=blocks.prod(0)
    return P

Сроки для массива (1000,1000) значительно улучшаются по сравнению с итерациями столбцов, хотя все еще намного медленнее, чем подход epsilon.

In [153]: timeit strided(A)
1 loops, best of 3: 2.51 s per loop

Другой подход к индексированию, хотя и относительно прямолинейный, медленнее и быстрее производит ошибки памяти.

def foo(A):
    h,w = A.shape
    I = (np.arange(w)[:,None]+np.arange(1,w))
    I1 = np.array(I)%w
    P = A[:,I1].prod(2)
    return P

Ответ 3

Im на ходу, поэтому у меня нет времени для разработки этого решения; но какой id делает создание непрерывного кругового обзора по последней оси, путем объединения массива к себе вдоль последней оси, а затем используйте np.lib.index_tricks.as_strided, чтобы выбрать соответствующие элементы, чтобы взять np.prod поверх, Нет циклов питона, нет числового приближения.

edit: здесь вы идете:

import numpy as np

A = np.array([[0.2, 0.4, 0.6],
              [0.5, 0.5, 0.5],
              [0.5, 0.0, 0.5],
              [0.6, 0.4, 0.2]])

B = np.concatenate((A,A),axis=1)
C = np.lib.index_tricks.as_strided(
        B,
        A.shape  +A.shape[1:],
        B.strides+B.strides[1:])
D = np.prod(C[...,1:], axis=-1)

print D

Примечание: этот метод не идеален, так как это O (n ^ 3). См. Мое другое опубликованное решение, которое является O (n ^ 2)

Ответ 4

Если вы готовы терпеть небольшую ошибку, вы можете использовать решение, которое вы впервые предложили.

A += 1e-10
np.around(np.repeat(np.prod(A, 1, keepdims = True), 3, axis = 1) / A, 9)

Ответ 5

Вот метод O (n ^ 2) без петонов питона или численное приближение:

def double_cumprod(A):
    B = np.empty((A.shape[0],A.shape[1]+1),A.dtype)
    B[:,0] = 1
    B[:,1:] = A
    L = np.cumprod(B, axis=1)
    B[:,1:] = A[:,::-1]
    R = np.cumprod(B, axis=1)[:,::-1]
    return L[:,:-1] * R[:,1:]

Примечание: он, по-видимому, примерно в два раза медленнее, чем метод численной аппроксимации, который соответствует ожиданиям.