Эффективно применять переход к матрице состояний с учетом матрицы переходов - PullRequest
Новая
       4

Эффективно применять переход к матрице состояний с учетом матрицы переходов

3 голосов
/ 02 ноября 2019

Я хочу применить изменение состояния к большой категориальной матрице (M) с k категориями, где я знаю вероятности перехода каждой категории к другой категории в k (T)

По сути, я хочучтобы иметь возможность эффективно взять каждый элемент в M, смоделировать изменение состояния с учетом вероятностей в T и заменить элемент вычисленным изменением.

Я испробовал несколько решений:

  • Грубая сила, вложенная для циклов с индексацией (слишком много, слишком длинная)
  • Вложенная с помощью Numba вложенная для циклов (~ 500 мсэто слишком долго для моих целей)
  • предварительно рассчитанные тиражи для каждой категории и замены (~ 400 мс) 
import numpy as np


def categorical_transition(mat, t_mat, k=4):

    transformed_mat = mat.copy()
    cat_counts = np.bincount(mat.reshape(-1,))

    for i in range(k):
        rand_vec = np.random.multinomial(1, t_mat[i], cat_counts[i])

        choice = np.where(rand_vec)[1]

        transformed_mat[mat == i] = choice

    return transformed_mat


# load data
mat = np.random.choice(4, (16000, 256))
t_mat = np.random.random((4, 4))

# normalize transition matrix
for i in range(t_mat.shape[0]):
    t_mat[i] = t_mat[i] / t_mat[i].sum()

transformed_mat = categorical_transition(mat, t_mat)

Этот метод работает, но он медленный, и я был бы признателен за любыесоветы по более эффективным способам его реализации

Ответы [ 2 ]

1 голос
/ 04 ноября 2019

Всегда предоставляйте все реализации, которые вы пробовали до сих пор

Я попробовал, например, простую реализацию, описанную здесь . Она должна быть примерно в 20-80 раз быстрее, чемВаше решение, в зависимости от того, сколько ядер у вас есть.

Реализация 

@nb.njit(parallel=True)  
def categorical_transition_nb(mat_in, t_mat):
    mat=np.reshape(mat_in,-1)
    transformed_mat = np.empty_like(mat)
    for i in nb.prange(mat.shape[0]):
        rand_number=np.random.rand()
        probabilities=t_mat[mat[i],:]
        if rand_number<probabilities[0]:
            transformed_mat[i]=0
        else:
            for j in range(1,probabilities.shape[0]):
                if rand_number>=probabilities[j-1] and rand_number<probabilities[j]:
                    transformed_mat[i]=j

    return transformed_mat.reshape(mat_in.shape)

Время 

import numpy as np
import numba as nb

# load data
mat = np.random.choice(4, (16_000,256))
t_mat = np.random.random((4, 4))

# normalize transition matrix
for i in range(t_mat.shape[0]):
    t_mat[i] = t_mat[i] / t_mat[i].sum()

t_mat_2=np.cumsum(t_mat,axis=1)
%timeit transformed_mat_2 = categorical_transition_nb(mat, t_mat_2)
21.7 ms ± 1.85 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
0 голосов
/ 02 ноября 2019

Вот метод, который дает мне 5-кратное ускорение на примере 

import numpy as np

def categorical_transition(mat, t_mat, k=4):

    transformed_mat = mat.copy()
    cat_counts = np.bincount(mat.reshape(-1,))

    for i in range(k):
        rand_vec = np.random.multinomial(1, t_mat[i], cat_counts[i])

        choice = np.where(rand_vec)[1]

        transformed_mat[mat == i] = choice

    return transformed_mat

def pp(mat,t_mat):
    ps = t_mat.cumsum(1)
    ps /= ps[:,-1:]
    return (np.random.random(mat.shape+(1,))<ps[mat]).argmax(-1)


# load data
mat = np.random.choice(4, (16000, 256))
t_mat = np.random.random((4, 4))

# normalize transition matrix
for i in range(t_mat.shape[0]):
    t_mat[i] = t_mat[i] / t_mat[i].sum()

transformed_mat = categorical_transition(mat, t_mat)
transformed_mat_pp = pp(mat, t_mat)

# check correctness
from pprint import pprint
np.set_printoptions(3)

cnts = np.bincount(mat.ravel())
pprint([[np.bincount(tm[mat==i])/cnts[i] for tm in (transformed_mat,transformed_mat_pp)] + [t_mat[i]] for i in range(4)])

from timeit import timeit

print('OP',timeit(lambda:categorical_transition(mat, t_mat),number=10)*100,'ms')
print('pp',timeit(lambda:pp(mat, t_mat),number=10)*100,'ms')

Пример выполнения: 

[[array([0.186, 0.1  , 0.078, 0.637]),
  array([0.186, 0.099, 0.078, 0.637]),
  array([0.186, 0.099, 0.078, 0.637])],
 [array([0.303, 0.517, 0.088, 0.092]),
  array([0.303, 0.517, 0.089, 0.092]),
  array([0.303, 0.517, 0.088, 0.092])],
 [array([0.319, 0.27 , 0.329, 0.082]),
  array([0.319, 0.271, 0.328, 0.083]),
  array([0.318, 0.27 , 0.329, 0.082])],
 [array([0.408, 0.106, 0.264, 0.222]),
  array([0.409, 0.107, 0.263, 0.221]),
  array([0.408, 0.107, 0.264, 0.221])]]
OP 872.7993675973266 ms
pp 170.54899749346077 ms
...