MFCC coefficients

Question

I have read this, this, this, this and this as a reference for computing the MFCC for a given wave file. Although, I am sure the values look wrong.

In short I followed the procedure in link 5.

Frame the signal into short frames.
For each frame calculate the periodogram estimate of the power spectrum.
Apply the mel filterbank to the power spectra, sum the energy in each filter.
Take the logarithm of all filterbank energies.
Take the DCT of the log filterbank energies.
Keep DCT coefficients 2-13, discard the rest.

For step 1&2: I compute the STFT (Complex STFT matrix (numFrames x NFFT)) and then abs(magnitude) square yields periodogram estimate of the power spectrum.

for step 3: I used the numpy.fft.fftfreq which is where I am a little skeptical. The frequencies (Frequency axis values in Hz (NFFT) )to get the MEL scale were the ones which I got from the numpy.fft.fftfreq.

Are these frequencies right?

I somehow feel the MFCC values are incorrect because they are in a cycle.

Here are the first five columns of the 12 rows (since I consider the 12 coefficients)

Row 1:

-121.14120041896194,36.31415193982116,-14.643619564107524,-57.269690625660424,-234.4386822674871,32.84089116534659,

Row2:

0.0,0.0,0.0,0.0,0.0,0.0

Row3:

-121.22760208014765,36.066028260522074,-14.639033011310872,-57.2746737936947,-234.79568545601416,32.93127531287391

Row4:

0.0,0.0,0.0,0.0,0.0,0.0

Row5:

9.318848744263936e-14,-1.8994827565125353e-13,-4.1272764363905036e-14,-1.146561606635031e-14,-1.2768410873828004e-13,

Row6:

0.0,0.0,0.0,0.0,0.0,0.0

Row7:

-121.4016310115022,35.5644565194515,-14.62969158843071,-57.28437605058768,-235.51590947244588,33.11399608942877

Row8:

3.511817254704773e-14,2.9096677134632476e-14,-5.78167494726597e-14,-1.169310279972752e-14,-2.271276447305103e-14

Row9:

-122.02396968016332,33.75128464855077,-14.595168239618843,-57.315447372114335,-238.1039561000439,33.77468522727183,

Row10:

0.0,0.0,0.0,0.0,0.0,0.0

Row11:

-122.48161817552142,32.39877413428054,-14.568682122506447,-57.334732903026804,-240.01937917255765,34.26769359167513

Row12:

3.6307708978933654e-14,2.0687378665389866e-13,3.2738059255324807e-14,-3.1731562727398686e-14,3.216220044843834e-13

Since, I cannot post the entire code, I will post the important parts:

part 1 of the code is with the wave file:

self.samples, sampleRate= loadWAV(fileName, mono=True, startSec=None, endSec=None)

part2:

self.X, self.freqHz, timeSec = stft(self.samples)

where X: (2D ndarray) Complex STFT matrix (numFrames x NFFT) timeVec: (ndarray) Time axis values in seconds (numFrames) freqHz: (ndarray) Frequency axis values in Hz (NFFT)

part 3:

self.mfccFeature2 = run(self.X, self.freqHz)

where run leads to this :

 def run(self, mag, freq, nFilters=26, nCoeff=12):
     mag = mag**2

     # create filterbank for given frequency axis
     filterBank = melFilterbankLinFreq(freq, nFilters=26)

     """change 1 : since there is a mismatch of operands, the * function is changed to dot """
     # filter original magnitude spectrogram (original)
     magFilt = np.tile(mag, (nFilters, 1)) * filterBank.T

     **New version where I think the problem exists:**
     magFilt_1 = np.tile(mag, (nFilters, 1))
     magFilt = np.dot(magFilt_1, filterBank.T)

     # accumulate energy values over filterbanks
     fbEnergyLog = np.log(np.sum(magFilt, axis=1)+.0000001)

     # perform DCT and keep only first nCoeff coefficients
     # DCT can work only on real values and the values are complex 128 instead of float
     return dct(fbEnergyLog.real)[1:nCoeff+1]

Answer 1

Your code should look like this:

 filterEnergies = numpy.dot(power, filters) // power is 1-d array, filters are 2-d arrays. Result is 1-d array of energies

 logFilterEnergies = np.log(magFilt) // take a lot from them

You can find complete numpy-based implementation here:

https://github.com/cmusphinx/sphinxtrain/blob/master/python/cmusphinx/mfcc.py

Answer 2

I don't think that those frequencies look right. Here is a running example I also use from another Stackoverflow thread. You could it to calculate the filter bank again and then calculate the dot-product on this.

Source: How to get MFCC from FFT

import numpy
from scipy.fftpack import dct
from scipy.io import wavfile

sampleRate, signal = wavfile.read("file.wav")
numCoefficients = 13 # choose the sive of mfcc array
minHz = 0
maxHz = 22.000  

complexSpectrum = numpy.fft(signal)
powerSpectrum = abs(complexSpectrum) ** 2
filteredSpectrum = numpy.dot(powerSpectrum, melFilterBank())
logSpectrum = numpy.log(filteredSpectrum)
dctSpectrum = dct(logSpectrum, type=2)  # MFCC :)

def melFilterBank(blockSize):
    numBands = int(numCoefficients)
    maxMel = int(freqToMel(maxHz))
    minMel = int(freqToMel(minHz))

    # Create a matrix for triangular filters, one row per filter
    filterMatrix = numpy.zeros((numBands, blockSize))

    melRange = numpy.array(xrange(numBands + 2))

    melCenterFilters = melRange * (maxMel - minMel) / (numBands + 1) + minMel

    # each array index represent the center of each triangular filter
    aux = numpy.log(1 + 1000.0 / 700.0) / 1000.0
    aux = (numpy.exp(melCenterFilters * aux) - 1) / 22050
    aux = 0.5 + 700 * blockSize * aux
    aux = numpy.floor(aux)  # Arredonda pra baixo
    centerIndex = numpy.array(aux, int)  # Get int values

    for i in xrange(numBands):
        start, centre, end = centerIndex[i:i + 3]
        k1 = numpy.float32(centre - start)
        k2 = numpy.float32(end - centre)
        up = (numpy.array(xrange(start, centre)) - start) / k1
        down = (end - numpy.array(xrange(centre, end))) / k2

        filterMatrix[i][start:centre] = up
        filterMatrix[i][centre:end] = down

    return filterMatrix.transpose()

def freqToMel(freq):
    return 1127.01048 * math.log(1 + freq / 700.0)

def melToFreq(mel):
    return 700 * (math.exp(freq / 1127.01048 - 1))

Having every second row/ coefficient zeroed would suppose that your input signal (which seems to be a wav-file here) is already scaled based on your mel-filter bank. By using real signals it seems quite unrealistic to me. Furthermore row 1 and 3 for example look quite similar which of course can happen but (if he is using real-world audio files) might not hold true