備忘録とか日常とか

学んだこととかを書きます。

【曲紹介とかしてみる】Tombo in 7/4 (A Bossa Electrica)

今年は暖冬と聞いてすっかり油断していましたが、ここ数日で急激に寒くなりましたね。
寒すぎて外に出るのが億劫です・・・

こんな寒い日こそ、熱い曲を聴いてあったまりましょう!
今回はサンバの曲をご紹介します。季節外れとか言わない。

A Bossa Electricaで、「Tombo in 7/4」です。


[Count the Beats] A Bossa Electrica - Tombo in 7/4


タイトルにもある通り、7/4拍子です。ベースラインもいいのですが、ドラムとパーカッションが刻むビートが秀逸です。7拍子とか聴いてると気持ち悪さが先行しそうなものなのに、この曲に関しては延々とドラムソロを聴いていたくなります。。
サビの部分は誰もが聴いたことがあるでしょう。Samba de janeiroという曲で、世界中でダンスミュージックとして古くから親しまれています。日本では野球とかサッカーの応援歌で有名でしょう。
あろうことか、Samba de janeiro のサビ部分をそのまま丸パクリ引用している曲です。どうしてこうなった。


Tombo in 7/4のオリジナルはAirto Moreira(1941-)という方が作曲しました。サンバの本場ブラジルのパーカッショニストだそうです。この曲は色んなアレンジがされていますので、他のが知りたければyoutubeで調べてみるといいと思います。
で、上のアレンジをしたのがA Bossa Electricaというスウェーデンの6人グループ。調べてみましたが情報が少ないのであまりよくわかりませんでした・・・音楽的なジャンルとしては北欧ブラジリアンフュージョン?というらしいです。スウェーデン人が北欧ブラジリアンとはこれいかに。



あともう一つ、jazz界でも屈指の腕を持つピアニストであるMichel Camilo(1954-)もこの曲をアレンジしています。


Michel Camilo - Tombo7/4

さすがの超絶技巧です。。こっちのほうが気持ち悪さ3割り増しですねー
セカンドリフが最高にきもかっこいい(褒め言葉です)
ドラムはDave Weckl 、ベースはAnthony Jacksonが担当してます。個人技が突出しているので、フュージョン好きな方や何か楽器をされてる方はこのアレンジも気に入ると思います。興味ない人がBGMとして聴くと意味不明になりそうですが・・・


聴いてても暖かくならないかもしれないですけど、演奏者は汗だくになること請け合いです。
興味があれば是非。。

matlabとnumpy reshapeの違い

matlabは特に行列計算が容易に書けるが、ライセンスが有料なので気軽には使えなかった。そこで無料で導入できるpythonに、matlabのような使い勝手の良い数値計算ライブラリ(numpyとかmatplotlibとか)が開発されたわけである。機械学習を利用する上で便利な言語としてpythonが使われだしたのも、これらのライブラリが開発された頃からであるとどっかの本で読んだ気がしなくもない。

matlabを使っている時に、行列の構造を変えるreshapeという関数でつまずいたのでメモしておく。
pythonにもnumpy.ndarrayの中に全く同じ名前のメソッドがあるが、若干仕様が違う。

例えば以下のようなベクトルがあるとする

% in matlab
>> a = [0 1 2 3 4 5 6 7 8 9 10 11]

a =

     0     1     2     3     4     5     6     7     8     9    10    11

これをreshapeして(3,4)の行列にするのだが、pythonだと以下のようになる

# in python
>> a.reshape((3,4))
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])

これに対し、matlabで普通にreshape関数を使うと

% in matlab
>> reshape(a,[3,4])

ans =

     0     3     6     9
     1     4     7    10
     2     5     8    11

こうなってしまう。
列から埋めるか、行から埋めるかという違いがあるらしい。

pythonmatlabと同じような結果にするには、orderという引数を与える。

# in python
>> a.reshape((3,4),order="F")
array([[ 0,  3,  6,  9],
       [ 1,  4,  7, 10],
       [ 2,  5,  8, 11]])

こんな感じ。デフォルトではorder="C"となっている。詳しくはドキュメントを参照

matlabでの結果をpythonorder="C"の場合にする方法は、調べたところよくわからない。。
なので以下のようにする

% in matlab
>> reshape(a,[4,3]).'

ans =

     0     1     2     3
     4     5     6     7
     8     9    10    11

普通にreshapeしてから.'で転置をとるだけである。


pythonに慣れてたので不意打ちを食らった

Machine Learning: An Algorithmic Perspective 読んでいく(3)

完全に忘れてたけど前回の続き。
読んでいるのはこれ

例によってほぼ自分用です。

続きを読む

theano CNNで抽出した特徴をファイルに出力

タイトルの通り。
CNNで特徴抽出して、今まではそのままMLPで識別していた。だが抽出した特徴を使って別の手法を試す必要が出てきたので、特徴をファイルに出力する方法をメモしておく。

CNNは畳み込みとプーリングの繰り返しにより、画像から自動的に特徴を学習し、抽出する。
オーソドックスなCNNだと、抽出した特徴はそのまま多層パーセプトロン(MLP)に入力されて識別が行われる。なので、MLPに入力される前のデータをそのまま取り出せばそれが特徴となっているわけである。
厳密に言うと別にMLP直前じゃなくても、任意の層の出力を取り出せば特徴が得られるのだが、識別器の違いによる性能の比較をするならば同じ箇所の特徴を使わなければならない。
ここの論文によると、CNNのどの階層から特徴を取るかによって性能がかなり変わってくるらしい。CNNをハナから特徴抽出器として使う場合は注意する必要がある。

参考元は以下
Pythonの基礎 ファイル(CSV)に書き込む編 - Pythonの学習の過程とか


早速やってみたいところだが、とりあえず学習済みネットワークがないと特徴抽出できない。
データセットは適当に用意。以前の記事参照
で、学習。各layerはここ参照
DeepLearningTutorialに準拠してるとこ多し

# train_set_x, train_set_yにデータとラベルを格納してある
# conv-pool-conv-pool-conv-pool-MLP-softmax
# image_size: 150x150  color

def gradient_updates_momentum(cost, params, learning_rate, momentum):
    assert momentum < 1 and momentum >= 0
    updates = []
    for param in params:
        param_update = theano.shared(param.get_value()*0., broadcastable=param.broadcastable)
        updates.append((param, param - learning_rate*param_update))
        updates.append((param_update, momentum*param_update + (1. - momentum)*T.grad(cost, param)))
    return updates

def training(learning_rate=0.001, n_epochs=300, momentum=0.9,
             nkerns=[16, 32, 48], batch_size=50, L1_reg=0.00, L2_reg=0.00):
    
    index = T.lscalar()
    
    x = T.matrix("x")
    y = T.ivector("y")
    
    print "building the model..."
    layer0_input = x.reshape((batch_size, 3, 150, 150))
    
    layerC0 = ConvLayer(
        input=layer0_input,
        image_shape=(batch_size, 3, 150, 150),
        filter_shape=(nkerns[0], 3, 11, 11),
        padding=0,
        st=3
        )
    
    layerP0 = PoolLayer(
        input=layerC0.output,
        poolsize=3,
        st=2
        )    
    
    layerC1 = ConvLayer(
        input=layerN_out,
        image_shape=(batch_size, nkerns[0], 24, 24),
        filter_shape=(nkerns[1], nkerns[0], 5, 5),
        padding=2
        )
        
    layerP1 = PoolLayer(
        input=layerC1.output,
        poolsize=3,
        st=2
        )

    
    layerC2 = ConvLayer(
        input=layerP1.output,
        image_shape=(batch_size, nkerns[1], 12, 12),
        filter_shape=(nkerns[2], nkerns[1], 5, 5),
        padding=2
        )
    
    layerP2 = PoolLayer(
        input=layerC2.output,
        poolsize=3,
        st=2
        )
        
        
    # hidden layer 全結合層を定義
    layer3_input = layerP2.output.flatten(2)
    layer3 = HiddenLayer(
        rng=numpy.random.RandomState(23455),
        input=layer3_input,
        n_in=nkerns[2] * 6 * 6,
        n_out=500,
        activation=T.tanh,
        )
        
    # logistic regression 層を定義
    layer4 = LogisticRegression(input=layer3.output, n_in=500, n_out=3)
    
    # 誤差関数NLLの数式を定義
    L1 = ( abs(layerC0.W).sum() + abs(layerC1.W).sum() + abs(layerC2.W).sum() 
            + abs(layer3.W).sum() + abs(layer4.W).sum() )
    L2_sqr = ( abs(layerC0.W).sum() + abs(layerC1.W).sum() + abs(layerC2.W).sum() 
            + abs(layer3.W ** 2).sum()  +  abs(layer4.W ** 2).sum() )

    cost = ( layer4.negative_log_likelihood(y)
        +  L1 * L1_reg
        +  L2_sqr * L2_reg
        )
    params = (layer4.params + layer3.params + layerC2.params 
            + layerC1.params + layerC0.params)

    updates = gradient_updates_momentum(cost, params, learning_rate, momentum)
    train_model = theano.function(
        [index],
        cost,
        updates=updates,
        givens={
            x: train_set_x[index * batch_size: (index + 1) * batch_size],
            y: train_set_y[index * batch_size: (index + 1) * batch_size]
        }
    )

    print 'training ...'
    epoch = 0
    while (epoch < n_epochs):
        epoch = epoch + 1
        for minibatch_index in xrange(n_train_batches):

            iter = (epoch - 1) * n_train_batches + minibatch_index

            if iter % 100 == 0:
                print 'training @ iter = ', iter
            cost_ij = train_model(minibatch_index)
            
    # 学習したネットワークのモデルを保存
    cPickle.dump(layerC0, open("layerC0.pkl","wb"))
    cPickle.dump(layerC1, open("layerC1.pkl", "wb"))
    cPickle.dump(layerC2, open("layerC2.pkl", "wb"))
    cPickle.dump(layer3, open("layer3.pkl", "wb"))
    cPickle.dump(layer4, open("layer4.pkl", "wb"))

パラメータはかなり適当。

その後、モデルを使って特徴を抽出する。
csvwriterを使うと簡単にできる。

# 要csv, PIL, pylearn2

def feedforward_conv_pool(input, Clayer, pad, conv_st, pool_size, pool_st):
    input_shuffled = input.dimshuffle(1,2,3,0)
    filters_shuffled = Clayer.W.dimshuffle(1,2,3,0)
    conv_op = FilterActs(stride=conv_st, pad=pad, partial_sum=1)   #"full" convolution
    contiguous_input = gpu_contiguous(input_shuffled)
    contiguous_filters = gpu_contiguous(filters_shuffled)
    conv_out_shuffled = conv_op(contiguous_input, contiguous_filters)  
    conv_out = conv_out_shuffled.dimshuffle(3,0,1,2)
    
    pool_input = relu(conv_out + Clayer.b.dimshuffle('x', 0, 'x', 'x'))
    
    input_shuffled = pool_input.dimshuffle(1,2,3,0)
    pool_op = MaxPool(ds=pool_size, stride=pool_st)
    pooled_out_shuffled = pool_op(input_shuffled)
        
    pooled_out = pooled_out_shuffled.dimshuffle(3,0,1,2)
    
    return pooled_out

def feedforward_hidden(input, layer, activation=T.tanh, drop_rate=0):
    hidden_input = input.flatten(2)
    lin_output = T.dot(hidden_input, layer.W) + layer.b
    if drop_rate>0:
        output = lin_output * drop_rate
        return output
    else:
        output = activation(lin_output)
        return output

def feedforward_logistic(input, layer):
    output = T.nnet.softmax(T.dot(input,layer.params[0]) + layer.params[1])
    return output

def feature():
    layerC0 = cPickle.load(open("layerC0.pkl", "rb"))
    layerC1 = cPickle.load(open("layerC1.pkl", "rb"))
    layerC2 = cPickle.load(open("layerC2.pkl", "rb"))
    layer3 = cPickle.load(open("layer3.pkl", "rb"))
    layer4  = cPickle.load(open("layer4.pkl", "rb"))
    
    # 特徴抽出する画像を用意
    img = numpy.array(Image.open(...), dtype=theano.config.floatX)
    img = img.transpose(2,0,1)
    img = img.flatten()

    # 拡張子は.csvとか
    f = open("feature.csv","w")   
    dataWriter = csv.writer(f)
    x = T.matrix("x")
        
    layer0_input = x.reshape((1,3,150,150))
    layer0_out = feedforward_conv_pool(layer0_input, layerC0,
                                       pad=0, conv_st=3, pool_size=3, pool_st=2)
    layer1_out = feedforward_conv_pool(layer0_out, layerC1,
                                       pad=2, conv_st=1, pool_size=3, pool_st=2)
    layer2_out = feedforward_conv_pool(layer1_out, layerC2,
                                       pad=2, conv_st=1, pool_size=3, pool_st=2)
    layer3_out = feedforward_hidden(layer2_out, layer3)
    layer4_out = feedforward_logistic(layer3_out, layer4)
    # theano.functionで途中のlayerの出力を取り出す
    test_model = theano.function(
                    inputs=[layer0_input],
                    outputs=layer2_out,
                    )
    
    ft = test_model(img.reshape((1,3,150,150)))  # typeはnumpy.ndarray
    ft = ft.flatten()
    dataWriter.writerow(ft)
    f.close()

theano.functionoutputsに途中の層を指定することで、入力データをネットワークの最後まで通さなくても特徴を取り出せる。

ftにpool layerの出力がndarrayで入る。
この時点では[1, channel, row, column]の形になっているので、一次元配列になおしてからwriterow()でファイルに書き込む。
二次元配列の場合はwriterows()を使うとループを使わず一発で書き込める。
実際はこっちのが使うかも。

出力したファイルはcsv形式なので、

[特徴1,特徴2,特徴3,   ....]

のようにカンマ(,)で区切られた数字の羅列になる。
自分はこれをmatlabで使いたかったので、csvread('feature.txt')ですぐ読み込めた。

追記
numpyのsavetxt()loadtxt()でもできた。てかこっちのほうが簡単。

import numpy
a = numpy.zeros((3,3))
numpy.savetxt("test.txt",a)

読み込むときは

import numpy
a = numpy.loadtxt("test.txt")

でおけ
ただしこっちの方法だと全て小数でっぽい?桁数がすごいことになる
多分やり方はあると思うけど必要にかられてから調べよう。。



やっとこさ少しだけtheanoに慣れてきた気がしないでもない。

theanoでLocal Response Normalization(LRN)を使う

Deeplearning Tutorialでtheanoによる実装、アルゴリズムを勉強中。
前回のLCNに引き続いて、LRNの正規化についても試す。

今回はpylearn2内のコードがそのまま流用できるので、新しくコードを書いたりする必要はない。

参考元は以下
pylearn2/normalize.py at master · lisa-lab/pylearn2 · GitHub
MIRU2014 tutorial deeplearning


pylearn2内のコードではLocalResponseNormalizationではなくCrossChannelNormalizationという名前になっている。LCNとごっちゃになるのを避けてるのかもしれない。

LRNが最初に世に知られたのが、Imagenetの論文である。オーバーラッププーリングや活性化関数のReLUなどと並べて説明されており、このCNNの性能の要になっているように思える。

LRNは端的に述べると、「同一位置(ピクセル)において複数の特徴マップ間で正規化する」ということだそうだ。元の論文にも書いてあるが、LRNは"brightness normalization"であり、LCNのように輝度の平均を減算して0にしないことがミソらしい。
式で書くと以下のようになる(論文の式そのままです)


{\displaystyle
b^i_{x,y}=a^i_{x,y}/ \left( k+\alpha \sum^{min(N-1,i+\frac{n}{2})}_{j=max(0,i-\frac{n}{2})} (a^j_{x,y})^2 \right)^\beta
}


k, n, {\alpha}, {\beta}がパラメータである。{a^i_{x,y}}はi番目の特徴マップの(x,y)のピクセルを、Nは特徴マップの総数を表す。
summationの部分は、「i番目の特徴マップに対して、n近傍の特徴マップの二乗和をとる」という意味である。
下みたいなイメージ
f:id:may46onez:20160109161610p:plain


端っこの特徴マップに関しては片側だけ和をとる。

この処理をCNN内の1層として組み込む。Imagenetの論文では畳みこみ層の直後にLRNを配置しているが、Caffeで用意されているモデルではプーリングの直後に配置されていた。岡谷先生の[asin:4061529021:title]でもそのモデルが紹介されていたので、どっちでもいいのかもしれない。


実際に使うのはすごく簡単。
pylearn2からnormalize.pyのCrossChannelNormalizationを引っ張ってくる。
Deeplearning Tutorialでは画像データの扱いは[batch_size, channel, row, column](bc01)となっていたので、それにあったものをインポートする。

from pylearn2.expr.normalize import CrossChannelNormalizationBC01 as Normalize
# Normalizeは引数にalpha, k, beta, nが指定できる。
# デフォルトでは alpha=1e-4, k=2, beta=0.75, n=5
Normlayer = Normalize()

名前は適当に決めた。
あとはconvolutionかpoolingの出力をNormlayerに入れるだけである。
以前の記事に載せてあるConvLayerを使って書くと、

layerC0 = ConvLayer(...)
inputP0 = Normlayer(layerC0.output)
layerP0 = PoolLayer(inputP0, ...)

これでおけ。


theanoじゃなくてpylearn2をそのまま使いたいところだけど、pylearn2でCNN使ったら11epochで止まっちゃって学習できなかった。
同じ症状の人も結構いるっぽいけど直し方がよくわからないし・・・ドキュメント全然更新されてないし・・・
そのうちchainerとかtensorflowとかも触ってみようかな。