Gamma Process Non-negative Matrix Factorization (GaP-NMF)

This notebook demonstrates how GaP-NMF works on real data. The code used in the demo is a Julia translation of MATALB code for GaP-NMF as in:

Bayesian Nonparametric Matrix Factorization for Recorded Music by Matthew D. Hoffman et al. in ICML 2010

The original MATLAB code: http://www.cs.princeton.edu/~mdhoffma/code/gapnmfmatlab.tar



In [2]:

    
using PyPlot









    



INFO: Loading help data...



In [3]:

    
# Bayesian NMF package
using BNMF



In [4]:

    
# Read test file
using WAV
x, fs = wavread("arayuru.wav")
x = x[:] # monoral
println(length(x)/fs, " sec.")









    



3.79375 sec.



In [5]:

    
# Short-Time Fourier Transform
include("stft.jl")









    Out[5]:





istft (generic function with 1 method)



In [6]:

    
# Compute spectrogram
framelen = 2048
X = stft(x, framelen=framelen, hopsize=256, window=hanning(framelen))
X = abs(X[:,1:framelen/2])

figure(figsize=(16, 6), dpi=80, facecolor="w", edgecolor="k")
imshow(log(X)', origin="lower", aspect="auto")
colorbar()









    












    Out[6]:





PyObject <matplotlib.colorbar.Colorbar instance at 0x7f942299ab00>



In [7]:

    
p = GaPNMF(X, K=100, a=0.1, b=0.1, alpha=1.0, smoothness=100)
typeof(p)









    Out[7]:





GaPNMF (constructor with 3 methods)



In [8]:

    
tic()
fit!(p, epochs=100, verbose=true)
toc()









    



#1 bound: -258355.66569987283,
                    improvement: NaN,
                    activek: 100
#2 bound: -166964.96002387482,
                    improvement: 0.35373989352401103,
                    activek: 100
#3 bound: -129991.8965370445,
                    improvement: 0.22144205276091128,
                    activek: 100
#4 bound: -114130.39875258185,
                    improvement: 0.12201912739954912,
                    activek: 100
#5 bound: -106144.52246008086,
                    improvement: 0.0699715096046691,
                    activek: 100
#6 bound: -100518.69946800527,
                    improvement: 0.05300153848439386,
                    activek: 100
#7 bound: -94373.87015499652,
                    improvement: 0.0611312058903491,
                    activek: 100
#8 bound: -85407.36706304806,
                    improvement: 0.09501044173797438,
                    activek: 100
#9 bound: -71379.20170161093,
                    improvement: 0.16425006230529837,
                    activek: 100
#10 bound: -51067.294076852355,
                    improvement: 0.284563390182888,
                    activek: 99
#11 bound: -24669.156508926593,
                    improvement: 0.5169284577365425,
                    activek: 94
#12 bound: 4102.394646412344,
                    improvement: 1.166296510581052,
                    activek: 86
#13 bound: 37172.20978047181,
                    improvement: 8.061100402171185,
                    activek: 80
#14 bound: 71905.2228691249,
                    improvement: 0.9343811759853963,
                    activek: 70
#15 bound: 104135.88179713141,
                    improvement: 0.4482380784309606,
                    activek: 59
#16 bound: 133533.96505598474,
                    improvement: 0.28230503023083,
                    activek: 46
#17 bound: 157475.44680815807,
                    improvement: 0.17929132668333303,
                    activek: 39
#18 bound: 176234.83052766905,
                    improvement: 0.11912576912618193,
                    activek: 33
#19 bound: 191371.97749562288,
                    improvement: 0.08589191434310302,
                    activek: 31
#20 bound: 203320.95952385763,
                    improvement: 0.0624385146906268,
                    activek: 28
#21 bound: 212527.40919362585,
                    improvement: 0.045280376855038086,
                    activek: 27
#22 bound: 221271.98254512594,
                    improvement: 0.041145626273236285,
                    activek: 24
#23 bound: 227047.40465476713,
                    improvement: 0.026101009459990537,
                    activek: 22
#24 bound: 231900.9935011162,
                    improvement: 0.021376984483610847,
                    activek: 22
#25 bound: 236935.8543801122,
                    improvement: 0.02171125187081945,
                    activek: 22
#26 bound: 242323.12157731538,
                    improvement: 0.022737239204668788,
                    activek: 19
#27 bound: 245692.53493397307,
                    improvement: 0.013904630044073822,
                    activek: 19
#28 bound: 249330.1082934242,
                    improvement: 0.014805388207780368,
                    activek: 19
#29 bound: 253188.99310633534,
                    improvement: 0.015477010936721,
                    activek: 19
#30 bound: 257098.1937642234,
                    improvement: 0.015439852301344915,
                    activek: 18
#31 bound: 260811.65862233413,
                    improvement: 0.014443760976073729,
                    activek: 18
#32 bound: 264488.21050513384,
                    improvement: 0.014096577975923657,
                    activek: 16
#33 bound: 267381.2119991927,
                    improvement: 0.010938111337868948,
                    activek: 15
#34 bound: 269830.27669571573,
                    improvement: 0.009159449454999127,
                    activek: 15
#35 bound: 272294.75787379453,
                    improvement: 0.009133449397370501,
                    activek: 15
#36 bound: 274581.35622993053,
                    improvement: 0.008397511483477829,
                    activek: 15
#37 bound: 276529.68533708015,
                    improvement: 0.00709563509300361,
                    activek: 14
#38 bound: 278124.43559486896,
                    improvement: 0.005767012882703217,
                    activek: 14
#39 bound: 279586.42900052277,
                    improvement: 0.005256616170840243,
                    activek: 14
#40 bound: 280897.39361979894,
                    improvement: 0.004688942249316833,
                    activek: 13
#41 bound: 281907.5343082436,
                    improvement: 0.0035961198337493686,
                    activek: 13
#42 bound: 282847.9552359242,
                    improvement: 0.003335919807848574,
                    activek: 13
#43 bound: 283730.3996598744,
                    improvement: 0.0031198543514805036,
                    activek: 13
#44 bound: 284563.1679619092,
                    improvement: 0.002935069005764316,
                    activek: 13
#45 bound: 285351.8363099131,
                    improvement: 0.0027715053696248074,
                    activek: 13
#46 bound: 286101.56184193055,
                    improvement: 0.002627372375495129,
                    activek: 13
#47 bound: 286818.95426385384,
                    improvement: 0.0025074746789381106,
                    activek: 13
#48 bound: 287506.9574985963,
                    improvement: 0.002398736989011982,
                    activek: 13
#49 bound: 288154.78159329866,
                    improvement: 0.002253246670406354,
                    activek: 13
#50 bound: 288743.0419019079,
                    improvement: 0.0020414733545511578,
                    activek: 12
#51 bound: 289252.8625989089,
                    improvement: 0.0017656553510099952,
                    activek: 12
#52 bound: 289725.232744448,
                    improvement: 0.0016330699073983674,
                    activek: 12
#53 bound: 290162.90757765586,
                    improvement: 0.0015106548679310501,
                    activek: 12
#54 bound: 290569.38496781344,
                    improvement: 0.0014008592399040452,
                    activek: 12
#55 bound: 290948.291487965,
                    improvement: 0.0013040139111472632,
                    activek: 12
#56 bound: 291303.1780167969,
                    improvement: 0.0012197580780314925,
                    activek: 12
#57 bound: 291636.2081829912,
                    improvement: 0.0011432424749417106,
                    activek: 12
#58 bound: 291948.6619722085,
                    improvement: 0.001071382017905087,
                    activek: 12
#59 bound: 292242.2464596381,
                    improvement: 0.0010056031270922305,
                    activek: 12
#60 bound: 292518.53416279465,
                    improvement: 0.0009454064444946287,
                    activek: 12
#61 bound: 292778.7895722327,
                    improvement: 0.0008897057076499888,
                    activek: 12
#62 bound: 293025.04199046985,
                    improvement: 0.0008410869468957971,
                    activek: 12
#63 bound: 293259.715969421,
                    improvement: 0.0008008666336401717,
                    activek: 12
#64 bound: 293484.2688095808,
                    improvement: 0.0007657132157326826,
                    activek: 12
#65 bound: 293699.2632813068,
                    improvement: 0.000732558758934513,
                    activek: 12
#66 bound: 293905.0015859551,
                    improvement: 0.0007005067099921809,
                    activek: 12
#67 bound: 294101.767618701,
                    improvement: 0.0006694885479461332,
                    activek: 12
#68 bound: 294290.2079882837,
                    improvement: 0.0006407318497556445,
                    activek: 12
#69 bound: 294471.29677259346,
                    improvement: 0.0006153408417753195,
                    activek: 12
#70 bound: 294645.7946388076,
                    improvement: 0.0005925802213208346,
                    activek: 12
#71 bound: 294813.6967608542,
                    improvement: 0.0005698439451764352,
                    activek: 12
#72 bound: 294974.629941982,
                    improvement: 0.0005458809509054502,
                    activek: 12
#73 bound: 295128.54492309236,
                    improvement: 0.0005217905727711023,
                    activek: 12
#74 bound: 295275.94797967165,
                    improvement: 0.000499453743512694,
                    activek: 12
#75 bound: 295417.893278301,
                    improvement: 0.00048072082944974404,
                    activek: 12
#76 bound: 295555.37588700384,
                    improvement: 0.000465383484991994,
                    activek: 12
#77 bound: 295688.6802731279,
                    improvement: 0.00045103015204513733,
                    activek: 12
#78 bound: 295817.46751446003,
                    improvement: 0.0004355501239112557,
                    activek: 12
#79 bound: 295941.1726987199,
                    improvement: 0.0004181807967570829,
                    activek: 12
#80 bound: 296059.4741081589,
                    improvement: 0.00039974636972681316,
                    activek: 12
#81 bound: 296172.4748717353,
                    improvement: 0.000381682646423729,
                    activek: 12
#82 bound: 296280.6788987094,
                    improvement: 0.000365341266169173,
                    activek: 12
#83 bound: 296384.9379170988,
                    improvement: 0.00035189273487872333,
                    activek: 12
#84 bound: 296485.85837681143,
                    improvement: 0.0003405046842861189,
                    activek: 12
#85 bound: 296583.4631691834,
                    improvement: 0.00032920555775016393,
                    activek: 12
#86 bound: 296677.33816635475,
                    improvement: 0.00031652134670027217,
                    activek: 12
#87 bound: 296766.91172191023,
                    improvement: 0.00030192247277498294,
                    activek: 12
#88 bound: 296852.27510867815,
                    improvement: 0.0002876445567082202,
                    activek: 12
#89 bound: 296934.38152826263,
                    improvement: 0.0002765901644325338,
                    activek: 12
#90 bound: 297014.2363602156,
                    improvement: 0.0002689309050099011,
                    activek: 12
#91 bound: 297092.3936934883,
                    improvement: 0.00026314339080350937,
                    activek: 12
#92 bound: 297169.241321607,
                    improvement: 0.0002586657543241824,
                    activek: 12
#93 bound: 297245.1230941944,
                    improvement: 0.00025534867690174215,
                    activek: 12
#94 bound: 297320.00587926933,
                    improvement: 0.00025192267006919736,
                    activek: 12
#95 bound: 297393.48991517304,
                    improvement: 0.0002471546967934055,
                    activek: 12
#96 bound: 297465.2228586288,
                    improvement: 0.00024120549335564673,
                    activek: 12
#97 bound: 297535.04699156235,
                    improvement: 0.0002347304073483458,
                    activek: 12
#98 bound: 297602.8387909525,
                    improvement: 0.000227844753334511,
                    activek: 12
#99 bound: 297668.60083128273,
                    improvement: 0.00022097249003866942,
                    activek: 12
#100 bound: 297732.4308089611,
                    improvement: 0.00021443302215992348,
                    activek: 12
iteration finished
elapsed time: 65.55415797 seconds






    Out[8]:





65.55415797



In [9]:

    
# Get active K-components
good = goodk(p)

goodspec = p.Eh[good,:]









    Out[9]:





12x1024 Array{Float64,2}:
  0.00277722    0.00726524   0.0211036  …  0.000520452  0.000491621
  0.00280999    0.00493925   0.0108074     0.00037936   0.000340876
  0.00271666    0.00594996   0.0126983     0.000312902  0.0166263  
  0.00580218    0.0131258    0.0379031     0.172967     0.125151   
 11.8068       34.5922      84.7486        0.000627676  0.000588681
  0.00194749    0.0042438    3.31672    …  0.249485     0.274738   
  3.31153       3.91104      4.46299       0.0971082    0.0279081  
  0.0036766     0.00851538   0.0226062     0.000484815  0.000428812
  0.00377926    1.49767      4.74214       0.000353943  0.000278075
  0.000794334   0.456527     1.37461       0.867147     0.964443   
  0.00158975    0.00809176   0.0213092  …  0.446736     0.363154   
  0.464875      0.698002     0.84818       0.0435258    0.0458247



In [10]:

    
# Visualize K-active spectrum
figure(figsize=(16, 6), dpi=80, facecolor="w", edgecolor="k")

for k=1:size(goodspec,1)
    plot(goodspec[k,:][:], label=string("#", k))
end
legend()









    












    Out[10]:





PyObject <matplotlib.legend.Legend object at 0x7f94225e8510>



In [11]:

    
# Reconstruction by goodk
Y = xbar(p, good)
figure(figsize=(16, 6), dpi=80, facecolor="w", edgecolor="k")
imshow(log(Y)', origin="lower", aspect="auto", interpolation="nearest")
colorbar()









    












    Out[11]:





PyObject <matplotlib.colorbar.Colorbar instance at 0x7f941c1d69e0>



In [12]:

    
# E[templates]
figure(figsize=(16, 12), dpi=80, facecolor="w", edgecolor="k")
for k=1:length(good)
    subplot(length(good), 1, k)
    plot(p.Eh[good[k],:][:])
    xlim(0,1024)
end



In [13]:

    
# E[activations]
figure(figsize=(16, 12), dpi=80, facecolor="w", edgecolor="k")
for k=1:length(good)
    subplot(length(good), 1, k)
    plot(p.Ew[:,good[k]][:])
    xlim(0, 229)
end



In [14]:

    
# Only length(good) components are active
figure(figsize=(16, 6), dpi=80, facecolor="w", edgecolor="k")
bar(good, p.Et[good])









    












    Out[14]:





(PyObject <matplotlib.patches.Rectangle object at 0x7f94190fb050>,PyObject <matplotlib.patches.Rectangle object at 0x7f94190fb690>,PyObject <matplotlib.patches.Rectangle object at 0x7f94190fbd10>,PyObject <matplotlib.patches.Rectangle object at 0x7f94190873d0>,PyObject <matplotlib.patches.Rectangle object at 0x7f9419087a50>,PyObject <matplotlib.patches.Rectangle object at 0x7f9419087f50>,PyObject <matplotlib.patches.Rectangle object at 0x7f9419093790>,PyObject <matplotlib.patches.Rectangle object at 0x7f9419093e10>,PyObject <matplotlib.patches.Rectangle object at 0x7f94190a24d0>,PyObject <matplotlib.patches.Rectangle object at 0x7f94190a2b50>,PyObject <matplotlib.patches.Rectangle object at 0x7f94190ae210>,PyObject <matplotlib.patches.Rectangle object at 0x7f94190ae890>)



In [15]: