2017-01-05 4 views
2

Рассмотрим матрицу M1, задающую значения для всех комбинаций x,y. Рассмотрим раздел f(x)->X и раздел g(y)->Y. Кроме того, рассмотрите операцию p(A) на наборе A цифр, то есть max(A) или sum(A).numpy: матрица сжимающего блока

Отображения f,g могут быть использованы для создания из M1 блок-матрицы M2 где все x, которые отображаются в одной и той же X являются смежными, и то же самое для всех y.

Эта матрица M2 имеет блок для каждой комбинации «наборы» X,Y.

Теперь я хотел бы сконденсировать эту матрицу M2 в другую матрицу M3 путем применения p на каждом блоке отдельно. M3 имеет одно значение для каждой комбинации X,Y.

В идеале, я хотел бы пропустить преобразование M1 в M2 с использованием f и g на лету.

Каким будет наиболее эффективный способ выполнения такой операции, и можно ли было бы установить для него numpy или scipy?

Особый случай: На самом деле, в моем случае x и y идентичны, и есть только одна функция f применяется к обоим из них. Меня интересует только часть M2, которая находится под диагональю.

+0

Do '' f' и G' работают только со скалярными входами? В идеале для использования 'numpy' вы хотите записать их таким образом, чтобы работать с массивом (может быть 1d) значений, возвращая массив соответствующего размера. В противном случае вы застряли с повторением, так или иначе, над элементами 'M1'. Что вы надеетесь получить, пропустив «M2»? – hpaulj

ответ

3

Самый простой способ сделать это, хотя, возможно, и не самый эффективный (особенно если ваша матрица огромна) состоит в том, чтобы преобразовать вашу матрицу в одномерный массив, а затем иметь соответствующие массивы для индексы групп групп X и Y. Затем вы можете группировать индексы группы разделов и, наконец, реструктурировать матрицу обратно в ее исходную форму.

Например, если матрица

>>> M1 = np.arange(25).reshape((5,5)) 
>>> M1 
array([[ 0, 1, 2, 3, 4], 
     [ 5, 6, 7, 8, 9], 
     [10, 11, 12, 13, 14], 
     [15, 16, 17, 18, 19], 
     [20, 21, 22, 23, 24]]) 

и ваши перегородки

>>> def f(x): 
...  return np.array([1,1,1,2,2])[x] 
>>> def g(y): 
...  return np.array([3,4,4,4,5])[y] 

С этого момента, существует несколько способов реализации на перепрофилирование и последующую группировку. Вы можете сделать это с помощью Pandas, например, построив DataFrame и используя его метод stack() для «стека» всех строк друг над другом в одном столбце, индексированных их исходными индексами строк и столбцов.

>>> st = pd.DataFrame(M1).stack().to_frame('M1') 
>>> st 
    M1 
0 0 0 
    1 1 
    2 2 
    3 3 
    4 4 
1 0 5 
... 
4 3 23 
    4 24 

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

>>> st['X'] = f(st.index.get_level_values(0)) 
>>> st['Y'] = g(st.index.get_level_values(1)) 

Затем вы можете группировать эти индексы и применять свою функцию агрегации по выбору.

>>> stp = st.groupby(['X', 'Y']).agg(p) 

Вы должны определить p (или найти существующее определение) таким образом, что он принимает одномерный массив Numpy и возвращает один номер. Если вы хотите использовать что-то вроде sum(), вы можете просто использовать st.groupby(...).sum(), потому что у Pandas есть встроенная поддержка этого и несколько других стандартных функций, но agg является общим и работает для любой функции уменьшения p, которую вы можете предоставить.

И наконец, метод unstack() преобразует DataFrame обратно в надлежащую двумерную «матричную форму», а затем, если вы хотите, вы можете использовать метод as_matrix(), чтобы превратить его обратно в чистый массив Numpy.

>>> M3 = stp.unstack().as_matrix() 
>>> M3 
array([[ 15, 63, 27], 
     [ 35, 117, 43]]) 

Если вы не хотите вводить Pandas, есть и другие библиотеки, которые делают то же самое. Например, вы можете посмотреть на numpy-groupies, но по мере того, как я пишу это, он кажется совместимым только с Python 2. Однако я не нашел ни одной библиотеки, которая бы выполняла настоящую двумерную группировку, что вам может понадобиться, если вы работаете с очень больших матриц, достаточно больших, чтобы иметь дополнительные 2 или 3 копии из них исчерпала доступную память.

+0

Я бы описал 'vectorize' как простой способ, а не быстрый способ (что подразумевает улучшение скорости). – hpaulj

+0

Ваш 'f (x)' может быть записан как 'np.array ([1,1,1,2,2]) [x]' и будет работать намного быстрее, чем версия с «векторизованным». – hpaulj

+0

Это первый раз, когда я использовал '@ vectorize', используемый в качестве декоратора. Он работает, хотя он не позволяет использовать такие параметры, как 'otypes'. Часто, когда у плакатов есть проблемы с 'vectorize' (кроме ожиданий скорости), это потому, что им нужен один из дополнительных параметров. – hpaulj

1

Да M1 be numpy n x m массив. Вы можете начать с определения, какие разделы у вас есть. Конструктор set удаляет повторяющиеся записи, но упорядочивает их произвольно. Сортировать их только, чтобы иметь четко определенный порядок:

xs = sorted(set(f(i) for i in range(n))) 
ys = sorted(set(g(i) for i in range(m))) 

Чтобы построить матрицу блока для каждого X,Y вы можете использовать Numpy булевой индексацию вместе с сеткой-строительным помощником ix_ выбрать только строку и столбцы, принадлежащие X и Y, соответственно. И, наконец, применить p к выбранной подматрице:

from numpy import zeros, arange, ix_ 

ii, jj = arange(n), arange(m) 
M3 = zeros((len(xs), len(ys))) 

for k, X in enumerate(xs): 
    for l, Y in enumerate(ys): 
     M3[k,l] = p(M1[ix_(f(ii) == X, g(jj) == Y)]) 

Перегородка f и g должна применять поэлементно в Numpy массивов для этой работы. Как уже упоминалось в другом ответе, для достижения этого можно использовать декоратор numpy.vectorize.

Для примера:

from __future__ import division 
n = m = 5 
M1 = np.arange(25).reshape(5,5) 
f = lambda x: x // 3  # f(ii) = [0, 0, 0, 1, 1] 
g = lambda x: (x+2) // 3 # g(jj) = [0, 1, 1, 1, 2] 
p = numpy.sum 

M3 = [[ 15., 63., 27.], 
     [ 35., 117., 43.]] 
+0

Проблема в том, что у нас все еще есть две петли python. Я надеялся сделать это на C, используя numpy ... –

+1

Я не думаю, что это можно сделать без циклов. Альтернативой было бы иметь 4D-массив, индексированный с помощью '(X, Y, x, y)' и использовать 'apply_over_axes (p, a, (2,3))', но это просто сдвигает петли в конструкцию массива (и менее эффективен и явно строит «M2»). В любом случае, петли проходят только над _partitions_. Пока каждый содержит больше, чем просто несколько элементов, основная часть вычисления должна идти в оценку 'p'. –