अजगर/नकाबपोश सरणियों (और/या चयनात्मक रैंकिंग)

Question

को देखते हुए एक 2D numpy सरणी

MyArray = np.array([[ 8.02, 9.54, 0.82, 7.56, 2.26, 9.47], 
      [ 2.68, 7.3 , 2.74, 3.03, 2.25, 8.84], 
      [ 2.21, 3.62, 0.55, 2.94, 5.77, 0.21], 
      [ 5.78, 5.72, 8.85, 0.24, 5.37, 9.9 ], 
      [ 9.1 , 7.21, 4.14, 9.95, 6.73, 6.08], 
      [ 1.8 , 5.14, 5.02, 6.52, 0.3 , 6.11]]) 

और एक मुखौटा सरणी पर 2 डी गिरी रैंक छानने के लिए numpy सबसे तेजी से विधि

MyMask = np.array([[ 0., 0., 1., 1., 0., 1.], 
      [ 1., 0., 0., 0., 0., 1.], 
      [ 0., 0., 0., 1., 0., 0.], 
      [ 0., 1., 1., 1., 1., 0.], 
      [ 0., 1., 0., 1., 0., 0.], 
      [ 0., 1., 0., 0., 1., 1.]]) 

मैं चलाना चाहते हैं एक 'छेद' मध्य फिल्टर जो मुखौटा तत्वों को अनदेखा करता है।

उदाहरण के लिए, एक कर्नेल

k = np.array([[ 1, 1, 1], 
       [ 1, 0, 1], 
       [ 1, 1, 1]]); 

पर चलाना शामिल है के साथ एक रैंक फिल्टर: पड़ोस के प्रत्येक तत्व के लिए कर्नेल द्वारा परिभाषित छंटाई और केवल गैर नकाबपोश तत्वों की औसत लौटने (औसत अगर सरणी एक संख्या भी है)।

अब, वर्तमान में मैं नापसंदों को मुखौटा मैप करके बोटलनेक.नमेडियन का उपयोग करके, अवांछित लूप में ऐसा कर रहा हूं। यह मुझे बिल्कुल वही दे रहा है जो मुझे चाहिए, लेकिन मैं 2 डी सरणी मैनिपुलेशन दिनचर्या पर भरोसा करने की उम्मीद कर रहा था।

scipy.signal.order_filter और scipy.ndimage.filters.rank_filter दोनों उपलब्ध (rank_filter बहुत तेजी से प्रतीत होता है) कर रहे हैं, लेकिन ऐसा लगता है कि वे तरह NaN और Inf सरणी के शीर्ष पर रैंक लौटने और परिणाम biasing से पहले। ऐसा लगता है कि इनमें से कोई भी तरीका numpy.ma arrays (मास्किंग) का समर्थन नहीं करता है, न ही वे चुनिंदा रैंकों की एक सरणी स्वीकार करते हैं (फिर मैं 0 के साथ सभी मास्क भर सकता हूं और मेरी रैंक ऑफसेट कर सकता हूं), न ही कर्नेल को बदलने के लिए एक स्पष्ट तरीका है प्रत्येक स्थान।

मुझे आश्चर्य है कि क्या मुझे संयोजन और/या पायथन सुविधा याद आई है, या यदि मुझे साइथन में एक नई दिनचर्या लागू करने की तलाश करनी चाहिए।

सीमा से निपटने की अनदेखी करते हुए उपरोक्त सवाल के आंतरिक अंक होगा

[[ 0.  0.  0.  0.  0.  0. ] 
[ 0.  3.18 3.62 2.26 2.645 0. ] 
[ 0.  2.74 3.325 2.74 2.64 0. ] 
[ 0.  3.88 3.62 4.955 6.08 0. ] 
[ 0.  5.02 5.77 5.77 6.52 0. ] 
[ 0.  0.  0.  0.  0.  0. ]] 

Answer 1

एक तरह से RAM उपयोग बलिदान करने के लिए त्याग करने के लिए अजगर लूप है। अर्थात। हम मूल सरणी को उड़ाते हैं ताकि हम एक ही समय में सभी उप-सरणी पर फ़िल्टर लागू कर सकें। जो कि Numpy broadcasting.

1000x1000 की एक सरणी के लिए वेक्टरकृत फ़ंक्शन लगभग 100x तेज करता है, मेरे परीक्षण में।

मेरे कोड में मैंने मास्किंग के लिए NaN का उपयोग किया, लेकिन कोड की कुछ अतिरिक्त पंक्तियों के साथ आप numpy.ma सरणी का भी उपयोग कर सकते हैं। और मेरे पास nanmedian फ़ंक्शन नहीं था इसलिए मैंने nanmean का उपयोग किया, प्रदर्शन हालांकि तुलनीय होना चाहिए।

import numpy as np 
from numpy.lib.stride_tricks import as_strided 

# test data 
N = 1000 
A = np.random.rand(N, N)*10 
mask = np.random.choice([True, False], size=(N, N)) 

def filter_loop(A, mask): 
    kernel = np.array([[1,1,1],[1,0,1],[1,1,1]], bool) 
    A = A.copy() 
    A[mask] = np.nan 
    N = A.shape[0] - 2 # assuming square matrix 
    out = np.empty((N, N)) 
    for i in xrange(N): 
     for j in xrange(N): 
      out[i,j] = np.nanmean(A[i:i+3, j:j+3][kernel]) 
    return out  

def filter_broadcast(A, mask): 
    A = A.copy() 
    A[mask] = np.nan 
    N = A.shape[0] - 2 
    B = as_strided(A, (N, N, 3, 3), A.strides+A.strides) 
    B = B.copy().reshape((N, N, 3*3)) 
    B[:,:,4] = np.nan 
    return np.nanmean(B, axis=2)