2009-12-15 12 views
5

मैं 300000 सूची (फाइबर पटरियों), जहां प्रत्येक ट्रैक (एक्स, वाई, जेड) tuples/की एक सूची है की एक सूची है निर्देशांक:घुसपैठ गिनने के लिए अधिक कुशल तरीका?

tracks= 
[[(1,2,3),(3,2,4),...] 
[(4,2,1),(5,7,3),...] 
... 
] 

मैं भी मास्क का एक समूह है, जहां प्रत्येक मुखौटा है (एक्स, वाई, जेड) tuples/की एक सूची के रूप में परिभाषित किया गया है निर्देशांक:

mask_coords_list= 
[[(1,2,3),(8,13,4),...] 
[(6,2,2),(5,7,3),...] 
... 
] 

मैं मास्क के सभी संभव जोड़े के लिए, खोजने की कोशिश कर रहा हूँ:

  1. पटरियों की संख्या है कि प्रत्येक एक दूसरे को काटना मास्क-मास्क जोड़ी (एक कॉन बनाने के लिए ectivity मैट्रिक्स)
  2. पटरियों कि प्रत्येक मुखौटा एक दूसरे को काटना, प्रत्येक (एक्स, वाई, जेड के लिए 1 जोड़ने के लिए की सबसेट) सबसेट में प्रत्येक ट्रैक के लिए समन्वय (एक "घनत्व" छवि बनाने के लिए)
इस प्रक्रिया को

def mask_tracks(tracks,masks,masks_coords_list): 
    vox_tracks_img=zeros((xdim,ydim,zdim,len(masks))) 
    for track in tracks: 
     for count,mask in enumerate(masks_coords_list): 
      if any(set(track) & set(mask)): 
       for x,y,z in track: 
        vox_tracks_img[x,y,z,count] += 1 

चौराहों को खोजने के लिए सेट का उपयोग करते हुए तेज़ी आई है काफी लेकिन दोनों भागों stil:

मैं वर्तमान में भाग 1 इसलिए की तरह कर रहा हूँ:

def mask_connectivity_matrix(tracks,masks,masks_coords_list): 
    connect_mat=zeros((len(masks),len(masks))) 
    for track in tracks: 
     cur=[] 
     for count,mask_coords in enumerate(masks_coords_list): 
      if any(set(track) & set(mask_coords)): 
       cur.append(count) 
      for x,y in list(itertools.combinations(cur,2)): 
       connect_mat[x,y] += 1 

और भाग 2 इसलिए की तरह मैं एक घंटे से अधिक समय लेता हूं जब मेरे पास 70 या अधिक मास्क की सूची होती है। क्या प्रत्येक ट्रैक के लिए पुनरावृत्ति करने से ऐसा करने का एक और अधिक प्रभावी तरीका है?

+0

सभी उत्तरों मामूली सुधार लगते हैं, लेकिन मुझे लगता है कि आपको उससे अधिक की आवश्यकता है। – McPherrinM

+0

यदि आप किसी नमूना डेटा सेट और पेस्टबिन में सही उत्तरों को कहीं भी पोस्ट कर सकते हैं, तो आपको और सहायता मिल सकती है। –

+0

क्या मुझे यह सही लगता है कि चौराहे केवल दो समन्वय tuples के रूप में परिभाषित किया गया है, और समन्वय के बीच की रेखाओं के रूप में नहीं है? – Svante

उत्तर

3

वोक्सल निर्देशांक रैखिकरण, और उन्हें दो scipy.sparse.sparse.csc matrices में डाल दें।

चलो voxels की संख्या, मास्क की संख्या और टी ट्रैक की संख्या होने दें।
एम मास्क सीएससी मैट्रिक्स, आकार (एम एक्स वी) हो, जहां 1 (i, j) का मतलब मास्क मैं voxel j ओवरलैप करता हूं।
टी को ट्रैक सीएससी मैट्रिक्स, आकार (टी एक्स वी) होना चाहिए, जहां एक पर (के, जे) का अर्थ है ट्रैक के ओवरलैप वोक्सेल जे।

Overlap = (M * T.transpose() > 0) # track T overlaps mask M 
Connected = (Overlap * Overlap.tranpose() > 0) # Connected masks 
Density[mask_idx] = numpy.take(T, nonzero(Overlap[mask_idx, :])[0], axis=0).sum(axis=0) 

मैं पिछले एक पर गलत हो सकता है, और मुझे यकीन है कि css_matrices अशून्य & ले द्वारा पर संचालित किया जा सकता नहीं हूँ। आपको प्रत्येक कॉलम को लूप में खींचने और इसे पूर्ण मैट्रिक्स में बदलने की आवश्यकता हो सकती है।


मैंने कुछ प्रयोगों को अनुकरण करने की कोशिश की जो मैंने सोचा था कि एक उचित मात्रा में डेटा था। नीचे दिए गए कोड में 2-वर्षीय मैकबुक पर लगभग 2 मिनट लगते हैं। यदि आप csr_matrices का उपयोग करते हैं, तो इसमें लगभग 4 मिनट लगते हैं। संभवतः एक ट्रेडऑफ है कि प्रत्येक ट्रैक कितनी देर तक है।

from numpy import * 
from scipy.sparse import csc_matrix 

nvox = 1000000 
ntracks = 300000 
nmask = 100 

# create about 100 entries per track 
tcoords = random.uniform(0, ntracks, ntracks * 100).astype(int) 
vcoords = random.uniform(0, nvox, ntracks * 100).astype(int) 
d = ones(ntracks * 100) 
T = csc_matrix((d, vstack((tcoords, vcoords))), shape=(ntracks, nvox), dtype=bool) 

# create around 10000 entries per mask 
mcoords = random.uniform(0, nmask, nmask * 10000).astype(int) 
vcoords = random.uniform(0, nvox, nmask * 10000).astype(int) 
d = ones(nmask * 10000) 
M = csc_matrix((d, vstack((mcoords, vcoords))), shape=(nmask, nvox), dtype=bool) 

Overlap = (M * T.transpose()).astype(bool) # mask M overlaps track T 
Connected = (Overlap * Overlap.transpose()).astype(bool) # mask M1 and M2 are connected 
Density = Overlap * T.astype(float) # number of tracks overlapping mask M summed across voxels 
+0

यदि मैट्रिस का प्रकार टाइप करने के लिए सेट है, तो मुझे नहीं लगता कि "> 0" बिट्स अब आवश्यक हैं। –

+2

वास्तव में, सच नहीं है। कम से कम स्पैस मैट्रिस के लिए, गुणा उन्हें एक बाइट को बढ़ावा देता है। (?) मुझे उम्मीद है कि इसका मतलब यह नहीं है कि रैपरराउंड मुद्दे भी हैं। –

+0

इसके लिए धन्यवाद, मुझे लगभग 10 मिनट तक औसत ट्रैक लंबाई और 500 के औसत मास्क आकार के साथ एक मिनट तक चलाया गया। – jbrown

0

आप दोनों परिणामों को एक बार में बनाने के लिए संभवतः दो कार्यों को जोड़कर शुरू कर सकते हैं। लूपिंग से पहले संयोजनों की सूची बनाने की कोई आवश्यकता नहीं है, क्योंकि यह पहले से ही जनरेटर है, और यह आपको कुछ समय बचा सकता है।

def mask_connectivity_matrix_and_tracks(tracks,masks,masks_coords_list): 
    connect_mat=zeros((len(masks),len(masks))) 
    vox_tracks_img=zeros((xdim,ydim,zdim,len(masks))) 
    for track in tracks: 
     cur=[] 
     for count,mask_coords in enumerate(masks_coords_list): 
      if any(set(track) & set(mask_coords)): 
       cur.append(count) 
       for x,y,z in track: 
        vox_tracks_img[x,y,z,count] += 1 
      for x,y in itertools.combinations(cur,2): 
       connect_mat[x,y] += 1 

इसके अलावा, यह शायद कभी नहीं के रूप में "तेज" हो जाएगा "समाप्त हो गया इससे पहले कि हम मर जाते हैं", तो सबसे अच्छा तरीका है अंत में अजगर के लिए एक सी मॉड्यूल के रूप में Cython के साथ संकलित करने के लिए है।

0

यदि आपने प्रत्येक मास्क सेट पॉइंट्स को संग्रहीत किया है: (1,2,3), (1,2,4), (1,3,1) इस तरह के शब्दकोश के रूप में: {1: [{2: set([3, 4])}, {3: set([1])}]}, आप समाप्त हो सकते हैं तेजी से मैचों की जांच करने में सक्षम ... लेकिन शायद नहीं।

0

एक छोटा सा अनुकूलन (एक ही बड़ा-ओ, sligthly छोटे गुणक) अनावश्यक आपरेशन को हटाने के द्वारा किया जा सकता है:

  1. set तो कई बार प्रत्येक ट्रैक और मुखौटा पर कॉल नहीं है: ट्रैक प्रति एक बार इसे कहते और प्रति मुखौटा एक बार, ऊपर सेट की सहायक "समांतर" सूचियों ही if someset: के रूप में स्थापित करने के लिए, तो उन
  2. if any(someset): पर काम कर रहा है शब्दार्थ लेकिन थोड़ा धीमा

, एक नाटकीय फर्क नहीं करेंगे लेकिन शायद मदद कर सकते हैं।

0

लेम अभी तक एक और वृद्धिशील सुधार है कि बनाया जा सकता, मैं जानता हूँ कि सुझाव देने के लिए है, लेकिन:

छोटे पूर्णांक के सेट अजगर के लंबे ints का उपयोग कर बिट वैक्टर के रूप में तैयार किया जा सकता। मान लीजिए कि आप प्रत्येक टुपल को एक छोटे पूर्णांक आईडी के साथ प्रतिस्थापित करते हैं, फिर प्रत्येक ट्रैक और मास्क-कॉर्ड के प्रत्येक सेट को उन छोटे आईडी के सेट में परिवर्तित करें। आप उन सेटों को लंबे समय तक स्याही के रूप में प्रस्तुत कर सकते हैं, जिससे छेड़छाड़ ऑपरेशन थोड़ा तेज हो सकता है (लेकिन असम्बद्ध रूप से तेज़ नहीं)।

1

ठीक है, मुझे लगता है कि आखिर में ऐसा कुछ है जो जटिलता को कम करेगा। यह कोड वास्तव में आपके पास जो मिला है उसकी तुलना में उड़ना चाहिए।

ऐसा लगता है कि आपको पता होना चाहिए कि कौन से ट्रैक मास्क, incidence matrix के साथ मेल खाते हैं।

import numpy 
from collections import defaultdict 

def by_point(sets): 
    d = defaultdict(list) 
    for i, s in enumerate(sets): 
     for pt in s: 
      d[pt].append(i) 
    return d 

def calc(xdim, ydim, zdim, mask_coords_list, tracks): 
    masks_by_point = by_point(mask_coords_list) 
    tracks_by_point = by_point(tracks) 

    a = numpy.zeros((len(mask_coords_list), len(tracks)), dtype=int) 
    for pt, maskids in masks_by_point.iteritems(): 
     for trackid in tracks_by_point.get(pt,()): 
      a[maskids, trackid] = 1 
    m = numpy.matrix(a) 

adjacency matrix आप देख रहे हैं m * m.T है।

जो कोड आपने अभी तक ऊपरी त्रिकोण की गणना की है। आप उस आधे हिस्से को पकड़ने के लिए triu का उपयोग कर सकते हैं।

am = m * m.T # calculate adjacency matrix 
    am = numpy.triu(am, 1) # keep only upper triangle 
    am = am.A # convert matrix back to array 

वोक्सेल गणना घटना मैट्रिक्स का भी उपयोग कर सकती है।

vox_tracks_img = numpy.zeros((xdim, ydim, zdim, len(mask_coords_list)), dtype=int) 
    for trackid, track in enumerate(tracks): 
     for x, y, z in track: 
      vox_tracks_img[x, y, z, :] += a[:,trackid] 
    return am, vox_tracks_img 

मेरे लिए यह डेटा सेट के लिए एक सेकंड के नीचे चलता है जिसमें सैकड़ों मास्क और ट्रैक होते हैं।

यदि आपके पास मास्क में दिखाई देने वाले कई बिंदु हैं लेकिन किसी भी ट्रैक पर नहीं हैं, तो लूप में प्रवेश करने से पहले masks_by_point से उन बिंदुओं के लिए प्रविष्टियों को हटाने के लिए उपयुक्त हो सकता है।

संबंधित मुद्दे