2012-12-15 14 views
5

क्या कोई अजगर पैकेज है जो मुझे छवि की स्केवनेस और कुर्टोसिस को क्लैक्लुएट करने का एक तरीका प्रदान करेगा? कोई भी उदाहरण महान होगा।पाइथन में छवि skewness और kurtosis

बहुत बहुत धन्यवाद।

+0

छवि से क्या आपका मतलब है कि आपके पास अंक का एक सेट है? इन वर्णनकर्ताओं की गणना करने के लिए आपको क्या समस्या है? – mmgp

उत्तर

7

मुझे लगता है कि आपके पास एक ऐसी छवि है जो किसी प्रकार की चोटी दिखाती है, और आप एक्स और वाई दोनों दिशाओं में उस चोटी के स्केवनेस और कुर्टोसिस (और शायद मानक विचलन और केंद्र) प्राप्त करने में रूचि रखते हैं।

मैं इसके बारे में भी सोच रहा था। उत्सुकता से, मुझे यह किसी भी पायथन छवि विश्लेषण पैकेज में लागू नहीं हुआ। ओपनसीवी में moments function है, और हमें इनसे स्केवनेस प्राप्त करने में सक्षम होना चाहिए, लेकिन क्षण केवल तीसरे क्रम में जाते हैं, और हमें कुर्टोसिस प्राप्त करने के लिए चौथे आदेश की आवश्यकता होती है।

चीजों को आसान और तेज़ बनाने के लिए, मुझे लगता है कि एक्स और वाई दिशाओं में छवि का प्रक्षेपण करना और इन अनुमानों से आंकड़े ढूंढना गणितीय रूप से पूर्ण छवि का उपयोग करके आंकड़ों को ढूंढने के बराबर है। निम्नलिखित कोड में, मैं दोनों विधियों का उपयोग करता हूं और दिखाता हूं कि वे इस चिकनी उदाहरण के लिए समान हैं। एक असली, शोर छवि का उपयोग करके, मैंने पाया कि दोनों विधियां भी एक ही परिणाम प्रदान करती हैं, लेकिन केवल तभी जब आप छवि डेटा को फ्लोट 64 पर मैन्युअल रूप से डालें (यह फ्लोट 32 के रूप में आयात किया गया है, और "संख्यात्मक सामग्री" के परिणामस्वरूप परिणाम थोड़ा अलग हो गए हैं।

नीचे, एक उदाहरण के उत्पादन की

import numpy as np 
import matplotlib.pyplot as plt 
import time 

plt.figure(figsize=(10,10)) 

ax1 = plt.subplot(221) 
ax2 = plt.subplot(222) 
ax4 = plt.subplot(224) 

#Make some sample data as a sum of two elliptical gaussians: 
x = range(200) 
y = range(200) 

X,Y = np.meshgrid(x,y) 

def twoD_gaussian(X,Y,A=1,xo=100,yo=100,sx=20,sy=10): 
    return A*np.exp(-(X-xo)**2/(2.*sx**2)-(Y-yo)**2/(2.*sy**2)) 

Z = twoD_gaussian(X,Y) + twoD_gaussian(X,Y,A=0.4,yo=75) 

ax2.imshow(Z) #plot it 

#calculate projections along the x and y axes for the plots 
yp = np.sum(Z,axis=1) 
xp = np.sum(Z,axis=0) 

ax1.plot(yp,np.linspace(0,len(yp),len(yp))) 
ax4.plot(np.linspace(0,len(xp),len(xp)),xp) 

#Here is the business: 
def image_statistics(Z): 
    #Input: Z, a 2D array, hopefully containing some sort of peak 
    #Output: cx,cy,sx,sy,skx,sky,kx,ky 
    #cx and cy are the coordinates of the centroid 
    #sx and sy are the stardard deviation in the x and y directions 
    #skx and sky are the skewness in the x and y directions 
    #kx and ky are the Kurtosis in the x and y directions 
    #Note: this is not the excess kurtosis. For a normal distribution 
    #you expect the kurtosis will be 3.0. Just subtract 3 to get the 
    #excess kurtosis. 
    import numpy as np 

    h,w = np.shape(Z) 

    x = range(w) 
    y = range(h) 


    #calculate projections along the x and y axes 
    yp = np.sum(Z,axis=1) 
    xp = np.sum(Z,axis=0) 

    #centroid 
    cx = np.sum(x*xp)/np.sum(xp) 
    cy = np.sum(y*yp)/np.sum(yp) 

    #standard deviation 
    x2 = (x-cx)**2 
    y2 = (y-cy)**2 

    sx = np.sqrt(np.sum(x2*xp)/np.sum(xp)) 
    sy = np.sqrt(np.sum(y2*yp)/np.sum(yp)) 

    #skewness 
    x3 = (x-cx)**3 
    y3 = (y-cy)**3 

    skx = np.sum(xp*x3)/(np.sum(xp) * sx**3) 
    sky = np.sum(yp*y3)/(np.sum(yp) * sy**3) 

    #Kurtosis 
    x4 = (x-cx)**4 
    y4 = (y-cy)**4 
    kx = np.sum(xp*x4)/(np.sum(xp) * sx**4) 
    ky = np.sum(yp*y4)/(np.sum(yp) * sy**4) 


    return cx,cy,sx,sy,skx,sky,kx,ky 

#We can check that the result is the same if we use the full 2D data array 
def image_statistics_2D(Z): 
    h,w = np.shape(Z) 

    x = range(w) 
    y = range(h) 

    X,Y = np.meshgrid(x,y) 

    #Centroid (mean) 
    cx = np.sum(Z*X)/np.sum(Z) 
    cy = np.sum(Z*Y)/np.sum(Z) 

    ###Standard deviation 
    x2 = (range(w) - cx)**2 
    y2 = (range(h) - cy)**2 

    X2,Y2 = np.meshgrid(x2,y2) 

    #Find the variance 
    vx = np.sum(Z*X2)/np.sum(Z) 
    vy = np.sum(Z*Y2)/np.sum(Z) 

    #SD is the sqrt of the variance 
    sx,sy = np.sqrt(vx),np.sqrt(vy) 

    ###Skewness 
    x3 = (range(w) - cx)**3 
    y3 = (range(h) - cy)**3 

    X3,Y3 = np.meshgrid(x3,y3) 

    #Find the thid central moment 
    m3x = np.sum(Z*X3)/np.sum(Z) 
    m3y = np.sum(Z*Y3)/np.sum(Z) 

    #Skewness is the third central moment divided by SD cubed 
    skx = m3x/sx**3 
    sky = m3y/sy**3 

    ###Kurtosis 
    x4 = (range(w) - cx)**4 
    y4 = (range(h) - cy)**4 

    X4,Y4 = np.meshgrid(x4,y4) 

    #Find the fourth central moment 
    m4x = np.sum(Z*X4)/np.sum(Z) 
    m4y = np.sum(Z*Y4)/np.sum(Z) 

    #Kurtosis is the fourth central moment divided by SD to the fourth power 
    kx = m4x/sx**4 
    ky = m4y/sy**4 

    return cx,cy,sx,sy,skx,sky,kx,ky 


#Calculate the image statistics using the projection method 
stats_pr = image_statistics(Z) 

#Confirm that they are the same by using a 2D calculation 
stats_2d = image_statistics_2D(Z) 

names = ('Centroid x','Centroid y','StdDev x','StdDev y','Skewness x','Skewness y','Kurtosis x','Kurtosis y') 

print 'Statistis\t1D\t2D' 
for name,i1,i2 in zip(names, stats_2d,stats_pr): 
    print '%s \t%.2f \t%.2f'%(name, i1,i2) 

plt.show() 

स्क्रीन शॉट है। आप में कटौती और "image_statistics()" अपने खुद के कोड में कार्य करते हैं। उम्मीद है कि यह किसी के लिए काम करता पेस्ट करने के लिए सक्षम होना चाहिए! :) बस मज़े के लिए:

screen shot of output

एक और बात: आप छवियों के साथ वास्तव में क्या कर रहे हैं, इस पर निर्भर करते हुए, आप अपने छवि विश्लेषण के लिए ImageJ का उपयोग करने पर विचार कर सकते हैं - लेकिन सावधान रहें! moments plugin आपको स्क्व्यूनेस, कुर्टोसिस इत्यादि की गणना करने देगा। छविजे में विश्लेषण में "स्काईनेस" और "कुर्टोसिस" है >> मापन मेनू सेट करें, लेकिन मुझे लगता है कि यह वास्तव में तीव्रता हिस्टोग्राम की स्काईनेस और कुर्टोसिस पाता है (मैं था एक मिनट के लिए मूर्ख)।