Блог

Главная - Вопросы по программированию - Быстрое (разреженное) построение матрицы с помощью определенных векторных элементов

Быстрое (разреженное) построение матрицы с помощью определенных векторных элементов

#python #numpy #scipy

#python #numpy #scipy

Вопрос:

У меня есть два numpy-вектора a и b, и я хочу построить (идеально разреженную) матрицу C с помощью чего-то вроде:

C[i,j] = (a[j]>b[j])*(a[i] b[j])

Таким образом, только те элементы, в которых соблюдается определенное соотношение между j-индексами, будут вычисляться с использованием некоторых других элементов. Опять же, это легко сделать с помощью цикла i amp; j, но я хотел бы знать, есть ли более эффективный способ «numpy / scipy»? Проще говоря, я не мог понять, как справиться с различиями в используемых индексах.

Заранее спасибо!

Ответ №1:

Я должен был настаивать на [MCVE], но я сделаю один

Для массивов одинакового размера:

 In [232]: a = np.random.randint(0,10,10)
In [233]: b = np.random.randint(0,10,10)

a[j]>b[j] часть вашего вычисления:

 In [234]: a>b
Out[234]: 
array([ True, False, False, False, False,  True, False, False, False,
       False])

часть суммирования, a[i] b[j] :

 In [235]: a[:,None] b
Out[235]: 
array([[ 9,  9, 10, 14, 13, 11, 11, 12,  9, 14],
       [ 4,  4,  5,  9,  8,  6,  6,  7,  4,  9],
       [ 7,  7,  8, 12, 11,  9,  9, 10,  7, 12],
       [10, 10, 11, 15, 14, 12, 12, 13, 10, 15],
       [12, 12, 13, 17, 16, 14, 14, 15, 12, 17],
       [12, 12, 13, 17, 16, 14, 14, 15, 12, 17],
       [ 6,  6,  7, 11, 10,  8,  8,  9,  6, 11],
       [11, 11, 12, 16, 15, 13, 13, 14, 11, 16],
       [ 7,  7,  8, 12, 11,  9,  9, 10,  7, 12],
       [12, 12, 13, 17, 16, 14, 14, 15, 12, 17]])

Установка желаемых столбцов в 0:

 In [236]: _[:,np.nonzero(a>b)] = 0
In [237]: _
Out[237]: 
array([[ 0,  9, 10, 14, 13,  0, 11, 12,  9, 14],
       [ 0,  4,  5,  9,  8,  0,  6,  7,  4,  9],
       [ 0,  7,  8, 12, 11,  0,  9, 10,  7, 12],
       [ 0, 10, 11, 15, 14,  0, 12, 13, 10, 15],
       [ 0, 12, 13, 17, 16,  0, 14, 15, 12, 17],
       [ 0, 12, 13, 17, 16,  0, 14, 15, 12, 17],
       [ 0,  6,  7, 11, 10,  0,  8,  9,  6, 11],
       [ 0, 11, 12, 16, 15,  0, 13, 14, 11, 16],
       [ 0,  7,  8, 12, 11,  0,  9, 10,  7, 12],
       [ 0, 12, 13, 17, 16,  0, 14, 15, 12, 17]])

Упс, я переключил это, я должен был установить для других столбцов значение 0.

Но нам не нужно делать это отдельными шагами:

 In [238]: (a>b)*(a[:,None] b)
Out[238]: 
array([[ 9,  0,  0,  0,  0, 11,  0,  0,  0,  0],
       [ 4,  0,  0,  0,  0,  6,  0,  0,  0,  0],
       [ 7,  0,  0,  0,  0,  9,  0,  0,  0,  0],
       [10,  0,  0,  0,  0, 12,  0,  0,  0,  0],
       [12,  0,  0,  0,  0, 14,  0,  0,  0,  0],
       [12,  0,  0,  0,  0, 14,  0,  0,  0,  0],
       [ 6,  0,  0,  0,  0,  8,  0,  0,  0,  0],
       [11,  0,  0,  0,  0, 13,  0,  0,  0,  0],
       [ 7,  0,  0,  0,  0,  9,  0,  0,  0,  0],
       [12,  0,  0,  0,  0, 14,  0,  0,  0,  0]])

Это вычисление сильно зависит от broadcasting . Если вы не знакомы с этим, это не будет иметь особого смысла.

Массивы:

 In [239]: a
Out[239]: array([5, 0, 3, 6, 8, 8, 2, 7, 3, 8])
In [240]: b
Out[240]: array([4, 4, 5, 9, 8, 6, 6, 7, 4, 9])

Тот факт, что многие (или только некоторые) столбцы будут равны 0, не имеет большого значения при выполнении numpy вычислений всего массива. Обычно проще и быстрее передавать целые массивы.

Вы можете выполнять итерации по столбцам и просто выполнять добавление для выбранных. Нет необходимости выполнять итерации по строкам. Но обычно это будет медленнее, чем [238]:

 In [245]: c =  np.zeros((10,10),int)
     ...: for j in range(10):
     ...:     if a[j]>b[j]:
     ...:         c[:,j] = a   b[j]

или перемещение сравнения за пределы цикла:

 In [249]: c =  np.zeros((10,10),int)
     ...: m = a>b
     ...: for j,v in enumerate(m):
     ...:     if v:
     ...:         c[:,j] = a   b[j]

Еще лучше:

 In [253]: c =  np.zeros((10,10),int)
     ...: idx = np.nonzero(a>b)[0]
     ...: print(idx)
     ...: c[:,idx] = a[:,None]   b[idx]

Для некоторой степени разреженности это последнее может быть быстрее, но не здесь:

 In [256]: %%timeit
     ...: c =  np.zeros((10,10),int)
     ...: idx = np.nonzero(a>b)[0]
     ...: c[:,idx] = a[:,None]   b[idx]
17.7 µs ± 772 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

In [258]: timeit (a>b)*(a[:,None] b)
10.7 µs ± 139 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

Ответ №2:

Нет ничего плохого в использовании циклов. Есть одно свойство, которое вы можете использовать, подумав: Условие a[j]>b[j] сообщает нам, может ли весь столбец быть равен нулю, поэтому, если мы используем scipy.sparse.coo_matrix((data,(i,j))) , имеет смысл создать внешний цикл над столбцами:

 data=[]
i=[]
j=[]
for jx = ...:
    if a[j]>b[j]:
        for ix = ...: # this whole loop can be skipped if the condition is not satisfied
            data.append(a[i] a[j])
            i.append(ix)
            j.append(jx)

Комментарии:

1. Спасибо за предложение. Я должен был подчеркнуть, что этот процесс будет выполняться очень часто, и, насколько я думал, добавление является одним из довольно медленных процессов, не так ли? И это также будет выполняться однопоточно, верно?

2. В этом случае должно быть легко распараллелить цикл jx . Добавление происходит не медленно, поскольку это операция O (1), но если это проблема, вы также можете рассмотреть возможность предварительного выделения i,j как np.array s и в зависимости от того, как a и b распределены, вы могли бы ускорить его, векторизовав второй цикл, возможно, даже сделав уступку разреженности, которая для каждого j , где a[j]>b[j] вы простовведите полный столбец. Но это то, что можете решить только вы, поскольку мы недостаточно знаем о вашей проблеме и ваших инструментах!

3. В зависимости от вашего варианта использования может даже не потребоваться реализовать создание этих матриц, но, возможно, реализовать только метод, который вычисляет результат умножения матрицы. Но, как я уже сказал, это не то, о чем мы можем судить со стороны. В качестве альтернативы вы можете попытаться векторизовать все, что может быть даже быстрее, но это опять же зависит от того, насколько велики a и b , какое у них распределение чисел и т. Д. То, что я предоставил здесь, было именно тем, что я бы счел наиболее «numpy / scipy»-способом.

4. ОП сказал, что у него были numpy массивы. Обычно мы стараемся избегать циклов (python) при работе с массивами.

5. @hpaulj В целом это верно, однако в данном случае это звучит так, как будто мы имеем дело с разреженными данными. Допустим, условие выполняется только в 0,01% случаев a[j]>b[j] , и a и b достаточно велики. Тогда нам, возможно, действительно придется работать с циклами, поскольку вся матрица может даже не поместиться в память. И память, вероятно, является проблемой, если OP уже хочет работать с разреженными структурами данных. Но это всего лишь предположение, поскольку мы не знаем точной ситуации, в которой находится OP.