2010-07-13 17 views
7

के साथ एक पेटी वितरण को फ़िट करना मेरे पास एक डेटा सेट है जो मुझे पता है कि एक पेटी वितरण है। क्या कोई मुझे सिसि में इस डेटा सेट को फिट करने के लिए इंगित कर सकता है? मुझे नीचे कोड चलाने के लिए मिला लेकिन मुझे नहीं पता कि मुझे क्या लौटाया जा रहा है (ए, बी, सी)। इसके अलावा, ए, बी, सी प्राप्त करने के बाद, मैं उनका उपयोग कर भिन्नता की गणना कैसे करूं?(python) Scipy

import scipy.stats as ss 
import scipy as sp 

a,b,c=ss.pareto.fit(data) 

उत्तर

1

बहुत सावधान फिटिंग कानून बनें !! कई रिपोर्ट किए गए बिजली कानून वास्तव में एक बिजली कानून द्वारा बुरी तरह फिट किए जाते हैं। सभी विवरणों के लिए Clauset et al. देखें (यदि आपके पास जर्नल तक पहुंच नहीं है तो arxiv पर भी)। उनके पास आलेख में companion website है जो अब पाइथन कार्यान्वयन से लिंक है। पता नहीं है कि यह Scipy का उपयोग करता है क्योंकि मैंने अपने आर कार्यान्वयन का उपयोग किया था जब मैंने इसे आखिरी बार इस्तेमाल किया था।

+1

पायथन कार्यान्वयन (http://code.google.com/p/agpy/wiki/PowerLaw) में दो संस्करण शामिल हैं; एक numpy पर निर्भर करता है, कोई नहीं करता है। (मैने यह लिखा) – keflavich

3

यहां एक त्वरित लिखित संस्करण है, जो कि संदर्भ पृष्ठ से कुछ संकेत लेता है। यह वर्तमान में scipy और figuresmodels में प्रगति पर काम कर रहा है और कुछ निश्चित या जमे हुए पैरामीटर के साथ एमएलई की आवश्यकता है, जो केवल ट्रंक संस्करणों में उपलब्ध है। पैरामीटर अनुमानों या अन्य परिणाम आंकड़ों पर कोई मानक त्रुटियां अभी तक उपलब्ध नहीं हैं।

'''estimating pareto with 3 parameters (shape, loc, scale) with nested 
minimization, MLE inside minimizing Kolmogorov-Smirnov statistic 

running some examples looks good 
Author: josef-pktd 
''' 

import numpy as np 
from scipy import stats, optimize 
#the following adds my frozen fit method to the distributions 
#scipy trunk also has a fit method with some parameters fixed. 
import scikits.statsmodels.sandbox.stats.distributions_patch 

true = (0.5, 10, 1.) # try different values 
shape, loc, scale = true 
rvs = stats.pareto.rvs(shape, loc=loc, scale=scale, size=1000) 

rvsmin = rvs.min() #for starting value to fmin 


def pareto_ks(loc, rvs): 
    est = stats.pareto.fit_fr(rvs, 1., frozen=[np.nan, loc, np.nan]) 
    args = (est[0], loc, est[1]) 
    return stats.kstest(rvs,'pareto',args)[0] 

locest = optimize.fmin(pareto_ks, rvsmin*0.7, (rvs,)) 
est = stats.pareto.fit_fr(rvs, 1., frozen=[np.nan, locest, np.nan]) 
args = (est[0], locest[0], est[1]) 
print 'estimate' 
print args 
print 'kstest' 
print stats.kstest(rvs,'pareto',args) 
print 'estimation error', args - np.array(true)