boblsturm@0: function [ Y, cost ] = nmf_beta( V, W, varargin )
boblsturm@0: 
boblsturm@0: if nargin > 2
boblsturm@0:     nmf_params = varargin{1};
boblsturm@0:     iterations = nmf_params.Iterations;
boblsturm@0:     lambda = nmf_params.Lambda;
boblsturm@0:     beta = nmf_params.Beta % 1: KL Divergence; 2: Euclidean
boblsturm@0: end
boblsturm@0: 
boblsturm@0: cost=0;
boblsturm@0: K=size(W, 2);
boblsturm@0: M=size(V, 2);
boblsturm@0: 
boblsturm@0: H=random('unif',0, 1, K, M);
boblsturm@0: 
boblsturm@0: V = V+1E-6;
boblsturm@0: W = W+1E-6;
boblsturm@0: 
boblsturm@0: for l=1:L-1    
boblsturm@0:     recon = W*H;
boblsturm@0:     num = H.*(W'*(((recon).^(beta-2)).*V));
boblsturm@0:     den = W'*((recon).^(beta-1));
boblsturm@0:     H = num./den;    
boblsturm@0: end
boblsturm@0: 
boblsturm@0: fprintf('Iterations: %i/%i\n', l, L);
boblsturm@0: fprintf('Convergence Criteria: %i\n', convergence*100);
boblsturm@0: fprintf('Repitition: %i\n', r);
boblsturm@0: fprintf('Polyphony: %i\n', p);
boblsturm@0: fprintf('Continuity: %i\n', c);
boblsturm@0: 
boblsturm@0: Y=H;
boblsturm@0: Y = Y./max(max(Y)); %Normalize activations
boblsturm@0: 
boblsturm@0: end