Как эффективно рассчитать пользовательские парные метри c в numpy или python для строковых данных? - PullRequest
Новая
       50

Как эффективно рассчитать пользовательские парные метри c в numpy или python для строковых данных?

0 голосов
/ 07 марта 2020

Я пытаюсь вычислить попарно взвешенный индекс Жакара в python списках строк. Взвешенная Jaccard включает в себя смежные ** пары значений в данном списке и взвешивает их так, чтобы пара имела вес 0,5 как для пересечения, так и для объединения пар между двумя массивами.

Взвешенная Jaccard рассчитывается следующим образом: (number_intersect + 0.5 * (num_pairs_intersect)) / (number_union + 0.5 * (num_pairs_union))

** = для моего сценария смежные пары могут включать элементы вверх на три позиции («частичная синхронность»), но в приведенном ниже примере рассматриваются только непосредственно смежные элементы как пары

Пример: 

list1 = ['a', 'b', 'c']  #pairs = ab, bc
list2 = ['a', 'b', 'c', 'd', 'e'] #pairs = ab, bc, cd, de

weight_jaccard(list1, list2) = (3+(0.5*2))/(5+(0.5*4)) = 0.57

У меня есть следующий сценарий python для расчета Взвешенный индекс Жакара. Проблема, с которой мне нужна помощь, заключается в том, что для ее выполнения потребуется ~ 270 дней, поскольку у меня более 1,3 миллиона списков / пар для выполнения парных вычислений, просто проходя по парам. Я попытался реализовать в numpy векторизованную функцию, которая вычисляет по паре массивов, сгенерированных из np.meshgrid, но выполнение на N = 1000 уже занимает больше времени, чем мой оригинальный скрипт (то есть ~ 17 se c против ~ 34 se c). 

from sys import argv
import numpy as np

script, input = argv

with open(input, "r") as f:
        dict = [i.strip().split("\t") for i in f]

#number of CGCs
dict_init = len(dict) #initial length of dict

#function to calculate weighted jaccard
def jaccard(cgcA, cgcB):
        A = cgcA[1:] #first element is the list tag/accession
        B = cgcB[1:] #first element is the list tag/accession
        inter = len(list(set(A) & set(B)))
        pairsA = makepairs(A)
        pairsB = makepairs(B)
        inter_pair = len(list(set(pairsA) & set(pairsB)))*0.5
        union = len(list(set(A) | set(B)))
        union_pair = len(list(set(pairsA) | set(pairsB)))*0.5
        jaccard_in = (inter+inter_pair)/(union+union_pair)
        return jaccard_in

#function to create pairs
def makepairs(cgc):
        distance = 4 ###max number of elements away that a item can be part of a pair
        pairs = []
        cgclen = len(cgc)
        a = 0
        while a < cgclen - 1: #create adjacent pairs
                pair=str(cgc[a]+"-"+cgc[a+1])
                pairs.append(pair)
                a = a + 1
        if distance != 1: #create nonadjacent pairs
                if distance == 2:
                        a = 0
                        while a < cgclen-3:
                                pair=str(cgc[a]+"-"+cgc[a+2])
                                pairs.append(pair)
                                a = a + 1
                else:
                        a = 0
                        for r in range(2, distance+1):
                                while a < cgclen-1-r:
                                        pair = str(cgc[a]+"-"+cgc[a+r])
                                        pairs.append(pair)
                                        a = a + 1
                                a = 0
        return pairs

#calculate weighted jaccard index for pairs of CGCs
for y in range(0, dict_init-1):
        A=dict[0]
        dict_len=len(dict)
        x=1
        while x < dict_len:
                B=dict[x]
                rhino = jaccard(A, B)
                if rhino > 0.0:
                        mongoose=[dict[0][0], dict[x][0], str(rhino)]
                        mongoose2 = "\t".join(mongoose)
                        print(mongoose2)
                x=x+1

        dict.pop(0) #remove first element in list to prevent redundant calculations


#calculate weighted jaccard index for pairs of CGCs with NUMPY
vfunc = np.vectorize(jaccard)
arr = np.array(dict)
xs, ys = np.meshgrid(arr, arr)
z = vfunc(xs, ys)
print(z)

Любое понимание или предложения будут оценены.

...