11

私が持っているデータに基づいて分布を作成し、その分布からランダムに描画しようとしています。ここに私が持っているものがあります:

from scipy import stats
import numpy

def getDistribution(data):
    kernel = stats.gaussian_kde(data)
    class rv(stats.rv_continuous):
        def _cdf(self, x):
            return kernel.integrate_box_1d(-numpy.Inf, x)
    return rv()

if __name__ == "__main__":
    # pretend this is real data
    data = numpy.concatenate((numpy.random.normal(2,5,100), numpy.random.normal(25,5,100)))
    d = getDistribution(data)

    print d.rvs(size=100) # this usually fails

これは私がやりたいことをやっていると思いますが、やろうとするとエラーが頻繁に発生し (以下を参照) d.rvs()、d.rvs(100)まったく機能しません。私は何か間違ったことをしていますか？これを行うためのより簡単またはより良い方法はありますか? scipy のバグである場合、回避する方法はありますか?

最後に、カスタム ディストリビューションの作成に関するドキュメントは他にありますか? 私が見つけた最高のものは scipy.stats.rv_continuous のドキュメントです。これは非常に質素で、有用な例が含まれていません。

トレースバック:

トレースバック (最新の呼び出しが最後): ファイル "testDistributions.py"、19 行目、出力 d.rvs(size=100) 内 ファイル "/usr/local/lib/python2.6/dist-packages/scipy-0.10.0 -py2.6-linux-x86_64.egg/scipy/stats/distributions.py"、696 行目、rvs vals = self._rvs(*args) ファイル "/usr/local/lib/python2.6/dist-packages /scipy-0.10.0-py2.6-linux-x86_64.egg/scipy/stats/distributions.py"、1193 行目、_rvs Y = self._ppf(U,*args) ファイル"/usr/local/lib /python2.6/dist-packages/scipy-0.10.0-py2.6-linux-x86_64.egg/scipy/stats/distributions.py"、1212 行目、_ppf で self.vecfunc(q,*args) ファイルを返す「/usr/local/lib/python2.6/dist-packages/numpy-1.6.1-py2.6-linux-x86_64.egg/numpy/lib/function_base.py」、1862行目、呼び出し中 theout = self.thefunc(*newargs) ファイル "/usr/local/lib/python2.6/dist-packages/scipy-0.10.0-py2.6-linux-x86_64.egg/scipy/stats/distributions.py" 、1158 行目、_ppf_single_call で return optimize.brentq(self._ppf_to_solve, self.xa, self.xb, args=(q,)+args, xtol=self.xtol) ファイル "/usr/local/lib/python2.6 /dist-packages/scipy-0.10.0-py2.6-linux-x86_64.egg/scipy/optimize/zeros.py"、366 行目、brentq r = _zeros._brentq(f,a,b,xtol,maxiter 内) ,args,full_output,disp) ValueError: f(a) と f(b) には異なる符号が必要です

編集

好奇心旺盛な人のために、以下の回答のアドバイスに従って、動作するコードを次に示します。

from scipy import stats
import numpy

def getDistribution(data):
    kernel = stats.gaussian_kde(data)
    class rv(stats.rv_continuous):
        def _rvs(self, *x, **y):
            # don't ask me why it's using self._size 
            # nor why I have to cast to int
            return kernel.resample(int(self._size)) 
        def _cdf(self, x):
            return kernel.integrate_box_1d(-numpy.Inf, x)
        def _pdf(self, x):
            return kernel.evaluate(x)
    return rv(name='kdedist', xa=-200, xb=200)
4

1 に答える 1

7

具体的には、トレースバックに:

rvs は、cdf の逆数 ppf を使用して乱数を作成します。ppf を指定していないため、根探索アルゴリズムによって計算されますbrentq。brentq関数がゼロの値を検索する必要がある場所の下限と上限を使用します (cdf(x)=q、q が分位数となるような x を見つけます)。

制限のデフォルト と は、このxa例xbでは小さすぎます。xa以下は、scipy 0.9.0 でxb機能します。関数インスタンスの作成時に設定できます。

def getDistribution(data):
    kernel = stats.gaussian_kde(data)
    class rv(stats.rv_continuous):
        def _cdf(self, x):
            return kernel.integrate_box_1d(-numpy.Inf, x)
    return rv(name='kdedist', xa=-200, xb=200)

現在、これを改善するための scipy のプル リクエストがあるため、次のリリースxaでは、例外xbを回避するために自動的に展開されます。f(a) and f(b) must have different signs

これに関するドキュメントはあまりありません。最も簡単なのは、いくつかの例に従うことです (メーリング リストで質問してください)。

編集：追記

pdf : gaussian_kde によっても指定された密度関数_pdfがあるため、いくつかの計算をより効率的にするメソッドを追加します。

edit2: 追加

rvs : 乱数の生成に関心がある場合は、gaussian_kde に resample メソッドがあります。ランダム サンプルは、データからサンプリングし、ガウス ノイズを追加することによって生成できます。したがって、これは ppf メソッドを使用する一般的な rvs よりも高速になります。gaussian_kde の resample メソッドを呼び出すだけの ._rvs メソッドを作成します。

ppfの事前計算: ppf を事前計算する一般的な方法を知りません。ただし、私が考えた方法 (ただし、これまで試したことはありません) は、多くのポイントで ppf を事前計算し、線形補間を使用して ppf 関数を近似することです。

edit3:_rvsコメントで Srivatsan の質問に答えようとしています

_rvspublic メソッドによって呼び出されるディストリビューション固有のメソッドrvsです。rvsいくつかの引数チェックを行い、場所とスケールを追加しself._size、要求された確率変数の配列のサイズである属性を設定してから、分布固有._rvsのメソッドまたは対応する汎用メソッドを呼び出すジェネリック メソッドです。の追加の引数._rvsは形状パラメーターですが、この場合は何もないため、*x冗長**yで未使用です。

sizeメソッドのor shape が.rvs多変量の場合にどの程度うまく機能するかはわかりません。これらの分布は一変量分布用に設計されており、多変量の場合には完全に機能しないか、何らかの形の変更が必要になる場合があります。

于 2012-05-21T04:10:40.723 に答える