Added updated Duet tests, transformer NMF and demo for it.
ADD
demos/demo_transformer_nmf.py
+import time
+import numpy as np
+import nussl
+
+import nimfa
+
+def simple_example():
+ """
+ A simple example of NMF in nussl
+ """
+ # Make two simple matrices
+ n = 4
+ a = np.arange(n ** 2).reshape((n, n))
+ b = 2 * a + 3
+
+ # Mix them together
+ mixture = np.dot(b, a)
+
+ # Set up NU NMF
+ n_bases = 2
+ n_iter = 100
+ nussl_nmf = nussl.TransformerNMF(mixture, n_bases)
+ nussl_nmf.should_use_epsilon = False
+ nussl_nmf.max_num_iterations = n_iter
+ nussl_nmf.distance_measure = nussl.DistanceType.DIVERGENCE
+
+ # run NU NMF
+ nussl_nmf.run()
+
+ signals = nussl_nmf.recombine_calculated_matrices()
+ # for sig in signals:
+ # print sig
+
+ V = mixture
+ nimfa_nmf = nimfa.Nmf(V, max_iter=n_iter, rank=n_bases, update='divergence')
+ nmf_fit = nimfa_nmf()
+ W = nmf_fit.basis()
+ H = nmf_fit.coef()
+
+ nimfa_est = np.dot(W, H)
+ nussl_est = np.dot(nussl_nmf.templates, nussl_nmf.activation_matrix)
+ error = np.abs(nussl_est - nimfa_est) / nimfa_est * 100
+
+ print('nimfa estimate:\n%s\n\n' % nimfa_est)
+ print('nussl estimate:\n%s\n\n' % nussl_est)
+ print('Errors:\n%s' % error)
+
+ sm = nmf_fit.summary()
+ # print( np.dot(W, H))
+
+ x = 1
+ y = 2
+ z = x + y
+
+if __name__ == '__main__':
+ simple_example()

MODIFY
nussl/Nmf.py
nussl/Nmf.py
@@ -26,7 +26,6 @@ class NMF(separation_base.SeparationBase):
def __str__(self):
return "Nmf"
- # TODO: Change this so that NMF accepts an AudioSignal object and not a raw stft
def __init__(self, input_audio_signal, num_templates,
activation_matrix=None, templates=None, distance_measure=None,
should_update_template=None, should_update_activation=None):

MODIFY
nussl/__init__.py
nussl/__init__.py
@@ -9,6 +9,7 @@ from .audio_signal import AudioSignal
from separation_base import SeparationBase
from Duet import Duet
from Nmf import NMF, DistanceType
+from transformer_nmf import TransformerNMF, DistanceType
from repet import Repet
from repet_sim import RepetSim
from ft2d import FT2D

ADD
nussl/transformer_nmf.py
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+import random
+import math
+import numpy as np
+
+import constants
+
+class TransformerNMF(object):
+ def __init__(self, input_matrix, num_templates,
+ activation_matrix=None, templates=None, distance_measure=None,
+ should_update_template=None, should_update_activation=None):
+
+ if num_templates <= 0:
+ raise Exception('Need more than 0 bases!')
+
+ self.input_matrix = input_matrix
+ self.num_templates = num_templates
+
+ if self.input_matrix.size <= 0:
+ raise Exception('Input matrix size must be > 0!')
+
+ if np.iscomplexobj(input_matrix):
+ raise Exception('Input matrix must be real')
+
+ if np.min(input_matrix) < 0.0:
+ raise Exception('Input matrix must be nonnegative!')
+
+ self.activation_matrix = None
+ self.templates = None
+
+ if activation_matrix is None and templates is None:
+ self.templates = np.zeros((self.input_matrix.shape[0], num_templates)) # W, in literature
+ self.activation_matrix = np.zeros((num_templates, self.input_matrix.shape[1])) # H, in literature
+ self.randomize_input_matrices()
+ elif activation_matrix is not None and templates is not None:
+ self.templates = templates
+ self.activation_matrix = activation_matrix
+ else:
+ raise Exception('Must provide both activation matrix and template vectors or nothing at all!')
+
+ self.distance_measure = distance_measure if distance_measure is not None else DistanceType.DEFAULT
+ self.should_update_template = True if should_update_template is None else should_update_template
+ self.should_update_activation = True if should_update_activation is None else should_update_activation
+
+ self.should_use_epsilon = False # Replace this with something more general
+ self.epsilon_euclidean_type = True
+ self.stopping_epsilon = 1e10
+ self.max_num_iterations = 20
+
+ self.reconstruction_error = []
+
+ def run(self):
+ """
+ This runs the NMF separation algorithm. This function assumes that all
+ parameters have been set prior to running.
+
+ No inputs. N must be set prior to calling this function.
+
+ Returns an activation matrix (in a 2d numpy array)
+ and a set of template vectors (also 2d numpy array).
+ """
+
+ if self.input_matrix is None or self.input_matrix.size == 0:
+ raise Exception('Cannot do NMF with an empty input_matrix!')
+
+ if self.num_templates is None or self.num_templates == 0:
+ raise Exception('Cannot do NMF with no bases!')
+
+ if self.should_use_epsilon:
+ print('Warning: User is expected to have set stopping_epsilon prior to using'
+ ' this function. Expect this to take a long time if you have not set'
+ ' a suitable epsilon.')
+
+ should_stop = False
+ num_iterations = 0
+ while not should_stop:
+
+ self.update()
+
+ self.reconstruction_error.append(self._euclidean_distance() if self.epsilon_euclidean_type else
+ self._divergence())
+
+ # Stopping conditions
+ num_iterations += 1
+ if self.should_use_epsilon:
+ if self.epsilon_euclidean_type:
+ should_stop = self._euclidean_distance() <= self.stopping_epsilon
+ else:
+ should_stop = self._divergence() <= self.stopping_epsilon
+ else:
+ should_stop = num_iterations >= self.max_num_iterations
+
+ return self.activation_matrix, self.templates
+
+ def update(self):
+ """
+ Computes a single update using the update function specified.
+ :return: nothing
+ """
+ # update activation matrix
+ if self.should_update_activation:
+ if self.distance_measure == DistanceType.EUCLIDEAN:
+ self.activation_matrix = self._update_activation_euclidean()
+
+ elif self.distance_measure == DistanceType.DIVERGENCE:
+ self.activation_matrix = self._update_activation_divergent()
+
+ # update template vectors
+ if self.should_update_template:
+ if self.distance_measure == DistanceType.EUCLIDEAN:
+ self.templates = self._update_template_euclidean()
+
+ elif self.distance_measure == DistanceType.DIVERGENCE:
+ self.templates = self._update_template_divergence()
+
+ def _update_activation_euclidean(self):
+ # make a new matrix to store results
+ activation_copy = np.empty_like(self.activation_matrix)
+
+ # store in memory so we don't have to do n*m calculations.
+ template_T = self.templates.T
+ temp_T_matrix = np.dot(template_T, self.input_matrix)
+ temp_T_act = np.dot(np.dot(template_T, self.templates), self.activation_matrix)
+
+ # Eq. 4, H update from [1]
+ for indices, val in np.ndenumerate(self.activation_matrix):
+ result = temp_T_matrix[indices]
+ result /= temp_T_act[indices]
+ result *= self.activation_matrix[indices]
+ activation_copy[indices] = result
+
+ return activation_copy
+
+ def _update_template_euclidean(self):
+ # make a new matrix to store results
+ template_copy = np.empty_like(self.templates)
+
+ # store in memory so we don't have to do n*m calculations.
+ activation_T = self.activation_matrix.T
+ input_matrix_act_T = np.dot(self.input_matrix, activation_T)
+ temp_act = np.dot(np.dot(self.templates, self.activation_matrix), activation_T)
+
+ # Eq. 4, W update from [1]
+ for indices, val in np.ndenumerate(self.templates):
+ result = input_matrix_act_T[indices]
+ result /= temp_act[indices]
+ result *= self.templates[indices]
+ template_copy[indices] = result
+
+ return template_copy
+
+ def _update_activation_divergent(self):
+ # make a new matrix to store results
+ activation_copy = np.empty_like(self.activation_matrix)
+
+ dot = np.dot(self.templates, self.activation_matrix)
+
+ # Eq. 5, H update from [1]
+ for indices, val in np.ndenumerate(self.activation_matrix):
+ (a, mu) = indices
+ result = sum((self.templates[i][a] * self.input_matrix[i][mu]) / dot[i][mu]
+ for i in range(self.templates.shape[0]))
+ result /= sum(self.templates[k][a] for k in range(self.templates.shape[0]))
+ result *= self.activation_matrix[indices]
+ activation_copy[indices] = result
+
+ return activation_copy
+
+ def _update_template_divergence(self):
+ # make a new matrix to store results
+ template_copy = np.empty_like(self.templates)
+
+ dot = np.dot(self.templates, self.activation_matrix)
+
+ # Eq. 5, W update from [1]
+ for indices, val in np.ndenumerate(self.templates):
+ (i, a) = indices
+ result = sum((self.activation_matrix[a][mu] * self.input_matrix[i][mu]) / dot[i][mu]
+ for mu in range(self.activation_matrix.shape[1]))
+ result /= sum(self.activation_matrix[a][nu] for nu in range(self.activation_matrix.shape[1]))
+ result *= self.templates[indices]
+ template_copy[indices] = result
+
+ return template_copy
+
+ def _euclidean_distance(self):
+ try:
+ mixture = np.dot(self.templates, self.activation_matrix)
+ except:
+ print(self.activation_matrix.shape, self.templates.shape)
+ return
+
+ if mixture.shape != self.input_matrix.shape:
+ raise Exception('Something went wrong! Recombining the activation matrix '
+ 'and template vectors is not the same size as the input matrix!')
+
+ return sum((self.input_matrix[index] - val) ** 2 for index, val in np.ndenumerate(mixture))
+
+ def _divergence(self):
+ mixture = np.dot(self.activation_matrix, self.templates)
+
+ if mixture.shape != self.input_matrix.shape:
+ raise Exception('Something went wrong! Recombining the activation matrix '
+ 'and template vectors is not the same size as the input matrix!')
+
+ return sum(
+ (self.input_matrix[index] * math.log(self.input_matrix[index] / val, 10) + self.input_matrix[index] - val)
+ for index, val in np.ndenumerate(mixture))
+
+ def recombine_calculated_matrices(self):
+ new_matrices = []
+ for n in range(self.num_templates):
+ matrix = np.empty_like(self.activation_matrix)
+ matrix[n,] = self.activation_matrix[n,]
+
+ new_matrix = np.dot(self.templates, matrix)
+ new_matrices.append(new_matrix)
+ return new_matrices
+
+ def randomize_input_matrices(self, shouldNormalize=False):
+ self._randomize_matrix(self.activation_matrix, shouldNormalize)
+ self._randomize_matrix(self.templates, shouldNormalize)
+
+ @staticmethod
+ def _randomize_matrix(M, shouldNormalize=False):
+ for i, row in enumerate(M):
+ for j, col in enumerate(row):
+ M[i][j] = random.random()
+
+ if not shouldNormalize:
+ M[i][j] *= constants.DEFAULT_MAX_VAL
+ return M
+
+ def plot(self, outputFile, **kwargs):
+ raise NotImplementedError('Sorry, you cannot do this yet.')
+
+class DistanceType:
+ EUCLIDEAN = 'euclidean'
+ DIVERGENCE = 'divergence'
+ DEFAULT = EUCLIDEAN
+
+ def __init__(self):
+ pass
\ No newline at end of file

MODIFY
tests/test_duet.py
tests/test_duet.py
@@ -10,13 +10,24 @@ import scipy.io
import os
import warnings
-
class DuetUnitTests(unittest.TestCase):
+
def setUp(self):
- # Use the same audio file that the benchmarks were run with
path = os.path.join('..', 'Input', 'dev1_female3_inst_mix.wav')
- signal = nussl.AudioSignal(path)
- self.duet = nussl.Duet(signal, 3)
+ self.signal = nussl.AudioSignal(path)
+ #Call benchmarks
+ self.duet = nussl.Duet(self.signal, 3)
+ self.benchmark_dict = self.load_benchmarks()
+
+ def load_benchmarks(self):
+ benchmark_dict = {}
+ directory = 'duet_reference/duet_benchmarks'
+ for filename in os.listdir(directory):
+ key = os.path.splitext(filename)[0]
+ file_path = os.path.join('duet_reference', 'duet_benchmarks', filename)
+ value = np.load(file_path)
+ benchmark_dict[key] = value
+ return benchmark_dict
# def test_duet_final_outputs(self):
# #Test final outputs
@@ -40,78 +51,75 @@ class DuetUnitTests(unittest.TestCase):
def test_compute_spectrogram_wmat(self):
# Load Duet values to benchmark against
- stft0_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_stft_ch0.npy')
- stft1_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_stft_ch1.npy')
- wmat_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_wmat.npy')
-
- benchmark_stft_ch0 = np.load(stft0_path)
- benchmark_stft_ch1 = np.load(stft1_path)
- benchmark_wmat = np.load(wmat_path)
+ duet = nussl.Duet(self.signal, 3)
+ duet_sft0, duet_sft1, duet_wmat = duet._compute_spectrogram(self.duet.sample_rate)
+ assert np.all(np.array(self.benchmark_dict['benchmark_stft_ch0']) == duet_sft0)
+ assert np.all(np.array(self.benchmark_dict['benchmark_stft_ch1']) == duet_sft1)
+ assert np.all(np.array(self.benchmark_dict['benchmark_wmat']) == duet_wmat)
- duet_sft0, duet_sft1, duet_wmat = self.duet._compute_spectrogram(self.duet.sample_rate)
- assert np.all(benchmark_stft_ch0 == duet_sft0)
- assert np.all(benchmark_stft_ch1 == duet_sft1)
- assert np.all(benchmark_wmat == duet_wmat)
def test_compute_atn_delay(self):
+ #Use the same stfts for comparing the two functions' outputs
+ duet = nussl.Duet(self.signal, 3)
+ duet.stft_ch0 = self.benchmark_dict['benchmark_stft_ch0']
+ duet.stft_ch1 = self.benchmark_dict['benchmark_stft_ch1']
+ duet.frequency_matrix = self.benchmark_dict['benchmark_wmat']
+
+ symmetric_atn, delay = duet._compute_atn_delay()
+
+ assert np.all(np.array(self.benchmark_dict['benchmark_sym_atn'])== symmetric_atn)
+ assert np.all(np.array(self.benchmark_dict['benchmark_delay'])== delay)
+
+ def test_make_histogram(self):
#Use the same stfts for comparing this function's outputs
- stft0_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_stft_ch0.npy')
- stft1_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_stft_ch1.npy')
- wmat_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_wmat.npy')
- self.duet.stft_ch0 = np.load(stft0_path)
- self.duet.stft_ch1 = np.load(stft1_path)
- self.duet.frequency_matrix = np.load(wmat_path)
+ duet = nussl.Duet(self.signal, 3)
+ duet.stft_ch0 = self.benchmark_dict['benchmark_stft_ch0']
+ duet.stft_ch1 = self.benchmark_dict['benchmark_stft_ch1']
+ duet.frequency_matrix = self.benchmark_dict['benchmark_wmat']
+ duet.symmetric_atn = self.benchmark_dict['benchmark_sym_atn']
+ duet.delay = self.benchmark_dict['benchmark_delay']
#Load benchmarks
- sym_atn_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_sym_atn.npy')
- delay_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_delay.npy')
- benchmark_sym_atn = np.load(sym_atn_path)
- benchmark_delay = np.load(delay_path)
-
- symmetric_atn, delay = self.duet._compute_atn_delay()
-
- assert np.all(benchmark_sym_atn == symmetric_atn)
- assert np.all(benchmark_delay == delay)
-
- # def test_make_histogram(self):
- # #Use the same stfts for comparing this function's outputs
- # stft0_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_stft_ch0.npy')
- # stft1_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_stft_ch1.npy')
- # wmat_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_wmat.npy')
- # sym_atn_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_sym_atn.npy')
- # delay_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_delay.npy')
- # self.duet.stft_ch0 = np.load(stft0_path)
- # self.duet.stft_ch1 = np.load(stft1_path)
- # self.duet.frequency_matrix = np.load(wmat_path)
- # self.symmetric_atn = np.load(sym_atn_path)
- # self.delay = np.load(delay_path)
- #
- # #Load benchmarks
- # hist_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_hist.npy')
- # atn_bins_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_atn_bins.npy')
- # delay_bins_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_delay_bins.npy')
- # benchmark_hist = np.load(hist_path)
- # benchmark_atn_bins = np.load(atn_bins_path)
- # benchmark_delay_bins = np.load(delay_bins_path)
- #
- # hist, atn_bins, delay_bins = self.duet.make_histogram(self.duet.p, self.duet.q)
- #
- # assert np.all(benchmark_hist == hist)
- # assert np.all(benchmark_atn_bins == atn_bins)
- # assert np.all(benchmark_delay_bins == delay_bins)
- #
- # def test_peak_indices(self):
- # hist_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_hist.npy')
- # peak_indices_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_peak_indices.npy')
- # benchmark_hist = np.load(hist_path)
- # benchmark_peak_indices = np.load(peak_indices_path)
- #
- # duet_peak_indices = nussl.utils.find_peak_indices(benchmark_hist, self.duet.num_sources, threshold=self.duet.peak_threshold,
- # min_dist=[self.duet.attenuation_min_distance, self.duet.delay_min_distance])
- #
- # assert np.all(benchmark_peak_indices == duet_peak_indices)
+ hist_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_hist.npy')
+ atn_bins_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_atn_bins.npy')
+ delay_bins_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_delay_bins.npy')
+ benchmark_hist = np.load(hist_path)
+ benchmark_atn_bins = np.load(atn_bins_path)
+ benchmark_delay_bins = np.load(delay_bins_path)
+
+ hist, atn_bins, delay_bins = duet.make_histogram(duet.p, duet.q)
+
+ assert np.all(benchmark_hist == hist)
+ assert np.all(benchmark_atn_bins == atn_bins)
+ assert np.all(benchmark_delay_bins == delay_bins)
+
+ def test_peak_indices(self):
+ duet = nussl.Duet(self.signal, 3)
+ duet.stft_ch0 = self.benchmark_dict['benchmark_stft_ch0']
+ duet.stft_ch1 = self.benchmark_dict['benchmark_stft_ch1']
+ duet.frequency_matrix = self.benchmark_dict['benchmark_wmat']
+ duet.symmetric_atn = self.benchmark_dict['benchmark_sym_atn']
+ duet.delay = self.benchmark_dict['benchmark_delay']
+
+ hist_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_hist.npy')
+ peak_indices_path = os.path.join('duet_reference', 'duet_benchmarks', 'benchmark_peak_indices.npy')
+ benchmark_hist = np.load(hist_path)
+ benchmark_peak_indices = np.load(peak_indices_path)
+
+ duet_peak_indices = nussl.utils.find_peak_indices(benchmark_hist, self.duet.num_sources, threshold=self.duet.peak_threshold,
+ min_dist=[self.duet.attenuation_min_distance, self.duet.delay_min_distance])
+
+ # original DUET returned peaks as [[28 21 44] [24 25 24]] instead of [(28, 24), (21,25), (44,24)]
+ # Convert current peaks being tested to original format
+ atn_indices = [x[0] for x in duet_peak_indices]
+ delay_indices = [x[1] for x in duet_peak_indices]
+ duet_peak_indices = np.row_stack((atn_indices, delay_indices))
+
+ assert np.all(benchmark_peak_indices == duet_peak_indices)
+
+ def test_convert_peaks(self):
+ pass
- # def test_convert_peaks
# def test_compute_masks
@@ -119,6 +127,7 @@ class DuetUnitTests(unittest.TestCase):
# test_smooth_matrix
+ #TODO: fix this function before writing test for it
if __name__ == '__main__':
unittest.main()

Added updated Duet tests, transformer NMF and demo for it.