2010-08-09 43 views
11

ढूँढना मैं एक मैट्रिक्स जो काफी बड़ी है (50K पंक्तियां आसपास) है, और मुझे मैट्रिक्स में प्रत्येक पंक्ति के बीच सहसंबंध गुणांक प्रिंट करना चाहते हैं। मैं इस तरह अजगर कोड लिखा है:सहसंबंध मैट्रिक्स

for i in xrange(rows): # rows are the number of rows in the matrix. 
    for j in xrange(i, rows): 
     r = scipy.stats.pearsonr(data[i,:], data[j,:]) 
     print r 

कृपया ध्यान दें कि मैं pearsonr scipy मॉड्यूल (http://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.pearsonr.html) से उपलब्ध समारोह का इस्तेमाल कर रही हूँ।

मेरे सवाल यह है: वहाँ ऐसा करने का एक तेज तरीका है? क्या कोई मैट्रिक्स विभाजन तकनीक है जिसका मैं उपयोग कर सकता हूं?

धन्यवाद!

उत्तर

0

आप, 10 सेट में अजगर मल्टीप्रोसेस मॉड्यूल, हिस्सा ऊपर अपनी पंक्तियों का उपयोग कर सकते अपने परिणामों बफ़र और फिर सामान प्रिंट आउट (यह केवल यह एक मल्टीकोर मशीन पर तेजी लाने के हैं, हालांकि)

http://docs.python.org/library/multiprocessing.html

बीटीडब्ल्यू: आपको अपने स्निपेट को फ़ंक्शन में बदलना होगा और यह भी विचार करना होगा कि डेटा रीसाइक्लिंग कैसे करें। प्रत्येक उपप्रक्रिया इस प्रकार की सूची है होने ... [startcord, stopcord, शौकीन] .. अच्छी तरह से

def myfunc(thelist): 
    for i in xrange(thelist[0]:thelist[1]): 
    .... 
    thelist[2] = result 
+0

मैं यहां आपके मतलब का एक और गहन उदाहरण देखना चाहता हूं। – vgoklani

+0

मुझे लगता है कि इस बिंदु पर मेरा उत्तर इस प्रश्न से बहुत दूर है, लेकिन यदि आप मल्टीप्रोसेशन में रुचि रखते हैं, तो जांचें: http://docs.python.org/library/multiprocessing.html ... अनिवार्य रूप से पंक्तियों के माध्यम से लूपिंग के बजाय , आप एक फ़ंक्शन और थ्रेड पूल बनाते हैं और बस p.map (myfunc, xrange (पंक्तियां) करते हैं) – pyInTheSky

10

नई समाधान

काम कर सकते हैं जो Kington के जवाब को देखने के बाद, मैं इस पर गौर करने का निर्णय लिया corrcoef() कोड और निम्न कार्यान्वयन करने के लिए इसके द्वारा प्रेरित किया गया था।

ms = data.mean(axis=1)[(slice(None,None,None),None)] 
datam = data - ms 
datass = np.sqrt(scipy.stats.ss(datam,axis=1)) 
for i in xrange(rows): 
    temp = np.dot(datam[i:],datam[i].T) 
    rs = temp/(datass[i:]*datass[i]) 

प्रत्येक लूप पंक्ति पंक्ति I और पंक्तियों के बीच मैं अंतिम पंक्ति के माध्यम से पियरसन गुणांक उत्पन्न करता हूं। यह बहुत तेज़ है। यह corrcoef() अकेले उपयोग करने के रूप में कम से कम 1.5x है क्योंकि यह गुणांक और कुछ अन्य चीजों की अनावश्यक गणना नहीं करता है। यह तेज़ भी होगा और आपको 50,000 पंक्ति मैट्रिक्स के साथ मेमोरी समस्याएं नहीं देगा क्योंकि तब आप किसी अन्य सेट को उत्पन्न करने से पहले आर के प्रत्येक सेट को स्टोर कर सकते हैं या प्रक्रिया कर सकते हैं। आर के किसी भी लंबे समय तक स्टोर किए बिना, मैं अपने काफी नए लैपटॉप पर एक मिनट के भीतर यादृच्छिक रूप से जेनरेट किए गए डेटा के 50,000 x 10 सेट पर चलाने के लिए उपर्युक्त कोड प्राप्त करने में सक्षम था।

पुरानी समाधान

पहले, मैं बाहर मुद्रण स्क्रीन के लिए आर के की सिफारिश नहीं होगा। 100 पंक्तियों (10 कॉलम) के लिए, यह आपके कोड का उपयोग किए बिना प्रिंटिंग बनाम 0.301 सेकंड के साथ 19.7 9 सेकेंड का अंतर है। बस आर को स्टोर करें और अगर आप चाहें तो उन्हें बाद में इस्तेमाल करें, या उन पर कुछ प्रोसेसिंग करें जैसे आप सबसे बड़े आर की तलाश में जाते हैं।

दूसरे, आप प्रचुरता से कुछ मात्रा की गणना नहीं द्वारा कुछ बचत कर सकते हैं। पियरसन गुणांक की गणना कुछ मात्राओं का उपयोग करके की जाती है जिसे आप पंक्ति के उपयोग के हर बार गणना करने के बजाय सटीक कर सकते हैं। इसके अलावा, आप पी-मूल्य (जो भी pearsonr() द्वारा दिया जाता है उपयोग नहीं कर रहे तो खरोंच जाने वह भी नीचे दिए गए कोड का उपयोग करना:। मैं सीधे scipy से अधिक 4.8x के बारे में की रफ्तार-अप प्राप्त

r = np.zeros((rows,rows)) 
ms = data.mean(axis=1) 

datam = np.zeros_like(data) 
for i in xrange(rows): 
    datam[i] = data[i] - ms[i] 
datass = scipy.stats.ss(datam,axis=1) 
for i in xrange(rows): 
    for j in xrange(i,rows): 
     r_num = np.add.reduce(datam[i]*datam[j]) 
     r_den = np.sqrt(datass[i]*datass[j]) 
     r[i,j] = min((r_num/r_den), 1.0) 

कोड जब मैंने पी-वैल्यू स्टफ को हटा दिया है - 8.8x अगर मैं वहां पी-वैल्यू स्टफ छोड़ देता हूं (मैंने सैकड़ों पंक्तियों के साथ 10 कॉलम इस्तेमाल किए हैं)। मैंने यह भी जांच की है कि यह वही परिणाम देता है। यह नहीं है वास्तव में एक बड़ा सुधार, लेकिन यह मदद कर सकता है।

आखिरकार, आप कंप्यूटिंग (50000) * (50001)/2 = 1,250,025,000 पियरसन गुणांक (यदि मैं सही ढंग से गिन रहा हूं) की समस्या से फंस गया हूं। यह बहुत ज्यादा है। वैसे, प्रत्येक पंक्ति के पियरसन गुणांक की गणना करने की वास्तव में कोई आवश्यकता नहीं है (यह 1 के बराबर होगी), लेकिन यह आपको केवल 50,000 पियरसन गुणांक की गणना करने से बचाता है। उपर्युक्त कोड के साथ, मुझे उम्मीद है कि यदि आपके छोटे डेटासेट पर मेरे परिणामों के आधार पर आपके डेटा में 10 कॉलम हैं तो आपकी गणना करने में लगभग 4 1/4 घंटे लगेंगे।

आप ऊपर दिए गए कोड को साइथन या कुछ समान में ले कर कुछ सुधार प्राप्त कर सकते हैं। मैं उम्मीद करता हूं कि यदि आप भाग्यशाली हैं तो आप शायद सीधे Scipy पर 10x सुधार प्राप्त करेंगे। साथ ही, जैसा कि pyInTheSky द्वारा सुझाया गया है, आप कुछ मल्टीप्रोसेसिंग कर सकते हैं।

6

क्या आपने अभी numpy.corrcoef का उपयोग करने का प्रयास किया है? यह देखते हुए कि आप पी-मानों का उपयोग क्यों नहीं कर रहे हैं, यह वही करना चाहिए जो आप चाहते हैं, जितना संभव हो उतना छोटा झगड़ा। (जब तक कि मैं वास्तव में याद नहीं कर रहा हूं कि वास्तव में क्या पियरसन आर है, जो काफी संभव है।)

बस यादृच्छिक डेटा पर परिणामों की तुरंत जांच कर रहा है, यह ठीक उसी तरह की चीज़ देता है जैसे @ जस्टिन छील का कोड ऊपर और 100x तेज चलाता है ।

उदाहरण के लिए, 1000 पंक्तियों और यादृच्छिक डेटा के 10 कॉलम के साथ चीजों का परीक्षण ...:

import numpy as np 
import scipy as sp 
import scipy.stats 

def main(): 
    data = np.random.random((1000, 10)) 
    x = corrcoef_test(data) 
    y = justin_peel_test(data) 
    print 'Maximum difference between the two results:', np.abs((x-y)).max() 
    return data 

def corrcoef_test(data): 
    """Just using numpy's built-in function""" 
    return np.corrcoef(data) 

def justin_peel_test(data): 
    """Justin Peel's suggestion above""" 
    rows = data.shape[0] 

    r = np.zeros((rows,rows)) 
    ms = data.mean(axis=1) 

    datam = np.zeros_like(data) 
    for i in xrange(rows): 
     datam[i] = data[i] - ms[i] 
    datass = sp.stats.ss(datam,axis=1) 
    for i in xrange(rows): 
     for j in xrange(i,rows): 
      r_num = np.add.reduce(datam[i]*datam[j]) 
      r_den = np.sqrt(datass[i]*datass[j]) 
      r[i,j] = min((r_num/r_den), 1.0) 
      r[j,i] = r[i,j] 
    return r 

data = main() 

की एक अधिकतम निरपेक्ष अंतर पैदावार ~ 3.3e -16 दो परिणाम

और समय के बीच :

In [44]: %timeit corrcoef_test(data) 
10 loops, best of 3: 71.7 ms per loop 

In [45]: %timeit justin_peel_test(data) 
1 loops, best of 3: 6.5 s per loop 

numpy.corrcoef सिर्फ तुम क्या चाहते क्या करना चाहिए, और यह बहुत तेजी से है।

+0

आप बिल्कुल सही हैं। मैंने पहले 'कोरकोफ' के बारे में सोचा, लेकिन कुछ कारणों से मुझे याद आया कि यह धीमा है। थोड़ा भेड़िया लग रहा था कि मैंने कोशिश की बजाय मेरी बुरी याददाश्त पर भरोसा किया। यह तेज़ है क्योंकि यह पाइथन लूप को खत्म करने के लिए मैट्रिक्स गुणाओं का उपयोग करता है। मुझसे +1 –

+0

हालांकि कॉर्कोकोफ के साथ समस्या यह है कि यह आवश्यकतानुसार लगभग दोगुनी स्मृति का उपयोग करता है। यह लगभग सभी गुणांक की गणना दो बार कर रहा है। हालांकि, बड़ी समस्या स्मृति है और ओपी को स्मृति समस्याओं से बचने के लिए डेटा को तोड़ना होगा। यह अनिवार्य रूप से एक combinatorics गड़बड़ बन जाएगा। –

+0

@ जस्टिन छील - सच है, corrcoef इनपुट सरणी की एक अतिरिक्त अस्थायी प्रति बना रहा है। यह गति और उपयोग की गई स्मृति की मात्रा के बीच एक व्यापार है। यदि स्मृति मुख्य बाधा है, और 50,000 पंक्तियों के साथ, आपका समाधान बहुत बेहतर है, यह होने की संभावना है। –

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