8

これは別のフォーラムへの質問かもしれません。もしそうなら、私に知らせてください。ウェーブレットタグをフォローしているのは14人だけであることに気づきました。

ここでは、pywt(pyWaveletsパッケージ)のウェーブレット分解を複数の次元に拡張するエレガントな方法を紹介します。pywtがインストールされている場合、これは箱から出して実行する必要があります。テスト1は、3Dアレイの分解と再構成を示しています。必要なのは、次元数を増やすことだけです。コードは、4、6、または18次元のデータで分解/再構成する際に機能します。

ここでpywt.wavedec関数とpywt.waverec関数を置き換えました。また、fn_decでは、新しいwavedec関数が古い関数と同じように機能することを示しています。

ただし、1つの落とし穴があります。それは、ウェーブレット係数をデータと同じ形状の配列として表します。結果として、ウェーブレットに関する知識が限られているため、ハールウェーブレットにしか使用できませんでした。DB4のような他の例では、この厳密な境界の端で係数をブリードします(配列のリストとしての係数の現在の表現[CA、CD1 ... CDN]では問題ありません。もう1つの問題は、これを2でしか処理していないことです。データの^Nエッジ直方体。

理論的には、「出血」が起こらないようにすることは可能だと思います。この種のウェーブレット分解と再構成のアルゴリズムは、William Press、Saul A teukolsky、William T. Vetterling、Brian P. Flannery(Second Edition)による「Cの数値レシピ」で説明されています。このアルゴリズムは、他の形式のエッジ拡張(zpdなど)ではなく、エッジでの反射を想定していますが、この方法は、他の形式の拡張で機能するのに十分一般的です。

この作業を他のウェーブレットに拡張する方法について何か提案はありますか?

注:このクエリはhttp://groups.google.com/group/pywaveletsにも投稿されています

ありがとう、阿城

import pywt
import sys
import numpy as np

def waveFn(wavelet):
    if not isinstance(wavelet, pywt.Wavelet):
        return pywt.Wavelet(wavelet)
    else:
        return wavelet

# given a single dimensional array ... returns the coefficients.
def wavedec(data, wavelet, mode='sym'):
    wavelet = waveFn(wavelet)

    dLen = len(data)
    coeffs = np.zeros_like(data)
    level = pywt.dwt_max_level(dLen, wavelet.dec_len)

    a = data    
    end_idx = dLen
    for idx in xrange(level):
        a, d = pywt.dwt(a, wavelet, mode)
        begin_idx = end_idx/2
        coeffs[begin_idx:end_idx] = d
        end_idx = begin_idx

    coeffs[:end_idx] = a
    return coeffs

def waverec(data, wavelet, mode='sym'):
    wavelet = waveFn(wavelet)

    dLen = len(data)
    level = pywt.dwt_max_level(dLen, wavelet.dec_len)

    end_idx = 1
    a = data[:end_idx] # approximation ... also the original data 
    d = data[end_idx:end_idx*2]    
    for idx in xrange(level):
        a = pywt.idwt(a, d, wavelet, mode)
        end_idx *= 2
        d = data[end_idx:end_idx*2]
    return a

def fn_dec(arr):
    return np.array(map(lambda row: reduce(lambda x,y : np.hstack((x,y)), pywt.wavedec(row, 'haar', 'zpd')), arr))
    # return np.array(map(lambda row: row*2, arr))

if __name__ == '__main__':
    test  = 1
    np.random.seed(10)
    wavelet = waveFn('haar')
    if test==0:
        # SIngle dimensional test.
        a = np.random.randn(1,8)
        print "original values A"
        print a
        print "decomposition of A by method in pywt"
        print fn_dec(a)
        print " decomposition of A by my method"
        coeffs =  wavedec(a[0], 'haar', 'zpd')
        print coeffs
        print "recomposition of A by my method"
        print waverec(coeffs, 'haar', 'zpd')
        sys.exit()
    if test==1:
        a = np.random.randn(4,4,4)
        # 2 D test
        print "original value of A"
        print a

        # decompose the signal into wavelet coefficients.
        dimensions = a.shape
        for dim in dimensions:
            a = np.rollaxis(a, 0, a.ndim)
            ndim = a.shape
            #a = fn_dec(a.reshape(-1, dim))
            a = np.array(map(lambda row: wavedec(row, wavelet), a.reshape(-1, dim)))
            a = a.reshape(ndim)
        print " decomposition of signal into coefficients"
        print a

        # re-composition of the coefficients into original signal
        for dim in dimensions:
            a = np.rollaxis(a, 0, a.ndim)
            ndim = a.shape
            a = np.array(map(lambda row: waverec(row, wavelet), a.reshape(-1, dim)))
            a = a.reshape(ndim)
        print "recomposition of coefficients to signal"
        print a
4

1 に答える 1

5

まず、シングルレベル多次元変換ソース)をすでに実装している関数を紹介します。n次元係数配列の辞書を返します。係数は、各ディメンションに適用される変換のタイプ(近似/詳細)を説明するキーによってアドレス指定されます。

たとえば、2Dの場合、結果は近似係数と詳細係数の配列を持つ辞書になります。

>>> pywt.dwtn([[1,2,3,4],[3,4,5,6],[5,6,7,8],[7,8,9,10]], 'db1')
{'aa': [[5.0, 9.0], [13.0, 17.0]],
 'ad': [[-1.0, -1.0], [-1.0, -1.0]],
 'da': [[-2.0, -2.0], [-2.0, -2.0]],
 'dd': [[0.0, 0.0], [0.0, -0.0]]}

ここaaで、は両方の次元(LL)に適用された近似変換を使用しdaた係数配列であり、は最初の次元に適用された詳細変換と2番目の次元(HL)に適用された近似変換を使用した係数配列です(dwt2出力と比較してください)。

これに基づいて、マルチレベルのケースに拡張するのはかなり簡単なはずです。

これが分解部分についての私の見解です:https ://gist.github.com/934166 。

また、質問で言及した1つの問題についても取り上げたいと思います。

ただし、1つの落とし穴があります。それは、ウェーブレット係数をデータと同じ形状の配列として表します。

入力データと同じ形状/サイズの配列として結果を表すというアプローチは、私の意見では有害です。とにかく、出力配列の係数にアクセスして逆変換を実行できるようにするには、仮定を行うか、インデックスを使用して2次データ構造を維持する必要があるため、全体を理解して操作するのが不必要に複雑になります(wavedec / waverecについてはMatlabのドキュメントを参照してください)。 )。

また、紙の上ではうまく機能しますが、前述の問題のために、実際のアプリケーションに常に適合するとは限りません。ほとんどの場合、入力データサイズは2 ^ nではなく、ウェーブレットフィルターを使用した畳み込み信号のデシメート結果は大きくなります。その「ストレージスペース」は、データの損失や不完全な再構築につながる可能性があります。

これらの問題を回避するには、Pythonのリスト、辞書、タプル(利用可能な場合)など、より自然なデータ構造を使用して結果データ階層を表すことをお勧めします。

于 2011-04-21T12:19:10.017 に答える