Эффективный способ зацикливания на 2D массив

Ответ 1

Как правило, мы избегаем использования двойного цикла for с numpy; это медленно и с умной индексацией (массив [:, iN]...) мы можем сделать многое в одном цикле.
Но для вашей проблемы со сверткой, это может быть самый простой ~~ (и единственный?) ~~ способ сделать то, что вы хотите. (редактировать: это не так. См. ниже @Masoud answer).

Bpmap_data=np.array(a[row-N:row+N,column-N:column+N].ravel()) создает новый массив в каждом цикле.
Вычисление медианы и стандартного ввода без создания нового массива (т.е. С использованием непосредственного представления) будет намного быстрее.

На самом деле, это в 40 раз быстрее (на моем экземпляре Google Colab)

Чтобы получить алгоритм в 400 раз быстрее, посмотрите ответ @Masoud, в котором используется фильтр scipy для 2D-массива.

import numpy as np
from astropy.stats import sigma_clipped_stats


N=10
a=np.random.random((80,40))


def f():
  """Your original code"""
  BPmap=[]
  for row in range(N,a.shape[0]-N):
      BPmap_row=[]
      for column in range(N,a.shape[1]-N):
          Bpmap_data=np.array(a[row-N:row+N,column-N:column+N].ravel()) 
          mean, median, std = sigma_clipped_stats(Bpmap_data, sigma=3,iters=5)
          BPmap_Nsigma=float(a[row][column]-median)/std                 
          BPmap_row.append(BPmap_Nsigma)
      BPmap.append(BPmap_row)
  return BPmap

def f2():
  """this little guy is improving a lot your work"""
  BPmap=[]
  for row in range(N,a.shape[0]-N):
      BPmap_row=[]
      for column in range(N,a.shape[1]-N):
          # the next 3 lines do not need any more memory
          view = a[row-N:row+N,column-N:column+N]
          std_without_outliers = view[view - view.mean() < 3*view.std()].std()
          median = np.median(view)
          # back to your workflow
          BPmap_Nsigma=float(a[row][column]-median)/std_without_outliers                 
          BPmap_row.append(BPmap_Nsigma)
      BPmap.append(BPmap_row)
  return BPmap

%time _ = f()
%time _ = f2()

f() == f2()

>>>CPU times: user 39.7 s, sys: 14.2 ms, total: 39.7 s
Wall time: 39.7 s
CPU times: user 969 ms, sys: 2.99 ms, total: 972 ms
Wall time: 976 ms
True

редактировать
Фактически, sigma_clipped_stats(a[row-N:row+N,column-N:column+N]) действительно замедляет цикл. Я подозреваю, что sigma_clipped_stats создает копию своего аргумента.

я беру стандартное после избавления от выбросов от 3 сигма

Я показываю здесь, как сделать это с помощью чистого NumPy; это действительно быстрее, чем функция, которую вы использовали раньше.

В конце концов, f() = f2(), так зачем больше использовать эту функцию astropy?

Ответ 2

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

1. Как описано здесь, избегайте использования for-циклов.
2. Вы на самом деле квадраты каждого числа 10 * 10 раз.

Вместо цикла for вы можете использовать библиотекарей Scipy.ndimage и opencv для выполнения свертки. В то время как эти библиотеки используются для обработки изображений, они настолько эффективны для обработки любого 2D-массива. Вот код, который выполняет ту же задачу, что и при использовании инструментов Scipy.ndimage, но в 1000 раз быстрее (23 мс против 27 с для массива 200X400). Я использовал приведенный здесь алгоритм для расчета стандартного отклонения:

import numpy as np
from scipy.ndimage.filters import uniform_filter, median_filter

a=np.random.random((200,400))
c1 = uniform_filter(a, size = (10,10))
c2 = uniform_filter(a*a, size = (10,10))
std = ((c2 - c1*c1)**.5)
med = median_filter(a, size=(10, 10))

BPmap = (a - med)/std

Ответ 3

Есть два способа сделать что-то быстрее: на самом деле ищите nop операторы /nop и исправляйте их ИЛИ бросайте на это деньги (больше вычислительной мощности).

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

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

Очередь может поддерживаться внутри python или во внешнем хранилище данных в памяти, например, redis.

Оригинальный однострочный ответ:

По моему мнению, единственный способ сделать это быстрее - это разделить работу между несколькими процессами в многоядерной системе.