# Linear Predictive Coding and Block Size

I'm experimenting with linear predictive coding as a low-complexity audio compression method, and I'm observing something very odd. I'm compressing a 20 MB mono 16-bit PCM song in different block sizes using 2-sample LPC (increasing the number of samples appears to benefit little compression-wise). The optimal block size appears to be 32k samples, which yields the lowest standard deviation of output samples; going lower than that increases the standard deviation and by 64 samples it's 4x as large as for 32k samples. I don't understand why this would be. Shouldn't a smaller sample size yield a tighter fit to the local data and give higher predictive accuracy?

Here are the functions I wrote to compute the coefficients. Note that it's a prototype that uses floats for convenience:

#define ADD_INDEX(x, add, size) ((x) + (add) - ((x) + (add) >= (size) ? (size) : 0))
#define SUB_INDEX(x, sub, size) ((x) - (sub) + ((x) < (sub) ? (size) : 0))

void ComputeCovars(const short samples[], size_t numSamples, short prevSamples[], size_t prevIdx, size_t numCovars, float covars[])
{
float mean = 0.0f;
for (size_t i = 0; i < numCovars; i++)
mean += prevSamples[i];
for (size_t i = 0; i < numSamples; i++)
mean += samples[i];
mean /= numSamples + numCovars;

for (size_t i = 0; i < numSamples; i++)
{
float sample = (float)samples[i] - mean;

covars[0] += sample * sample;
for (size_t j = 0; j < numCovars; j++)
covars[j + 1] += sample * ((float)prevSamples[SUB_INDEX(prevIdx, j + 1, numCovars)] - mean);

prevSamples[prevIdx] = samples[i];
}
}

void ComputeLinearCoeffs(const float covars[], size_t numCovars, float coeffs[])
{
for (size_t i = 0; i < numCovars; i++)
coeffs[i] = 0.0f;

switch (numCovars)
{
case 1:
coeffs[0] = covars[1] / covars[0];
break;

case 2:
{
float a = covars[0], b = covars[1], c = covars[1], d = covars[0],
y = covars[1], z = covars[2],
det = (a * d - b * c),
A = d, B = -b, C = -c, D = a;

coeffs[0] = (A * y + B * z) / det;
coeffs[1] = (C * y + D * z) / det;

break;
}

case 3:
{
float a = covars[0], b = covars[1], c = covars[2],
d = covars[1], e = covars[0], f = covars[1],
g = covars[2], h = covars[1], i = covars[0],
x = covars[1], y = covars[2], z = covars[3],
det = a * (e * i - f * h) - b * (i * d - f * g) + c * (d * h - e * g),
A = e * i - f * h, B = f * g - d * i, C = d * h - e * g,
D = c * h - b * i, E = a * i - c * g, F = b * g - a * h,
G = b * f - c * e, H = c * d - a * f, I = a * e - b * d;

coeffs[0] = (A * x + D * y + G * z) / det;
coeffs[1] = (B * x + E * y + H * z) / det;
coeffs[2] = (C * x + F * y + I * z) / det;

break;
}
}
}


And the code that gets executed on each block:

fill_n(autovars, nPredSamples + 1, 0.0f);
ComputeCovars(pnSrcBuffer, nSamples, workBuff, iSample, nPredSamples, autovars);
ComputeLinearCoeffs(autovars, numCoeffs, coeffs);

for (int i = 0; i < nSamples; i++)
{
float pred = 0;
for (int j = 0; j < numCoeffs; j++)
pred += coeffs[j] * samples[SUB_INDEX(iSample, j + 1, numCoeffs)];

samples[iSample] = pnSrcBuffer[i];
pnSrcBuffer[i] -= short(pred + (pred >= 0.0f ? 0.5f : -0.5f));

var += pnSrcBuffer[i] * pnSrcBuffer[i];
}


Is there some rational reason for this behavior, or is there a stupid bug in my code, or what?

EDIT: I've rewritten part of the method to prevent knees at the buffer borders, but the same general problem remains.

Making some progress on this. So far I've found that the method for computing the coefficients (from http://nptel.ac.in/courses/IIT-MADRAS/Principles_of_Communication1/Pdfs/1_9.pdf) appears to have issues. It will create coefficients with a sum of almost 1, but the coefficients can sometimes be quite large (e.g. 15 and -14 in one case). When this is applied to samples with sizable values (e.g. > 9000), this can lead to absurd prediction values such as 190,000.

Adding this simple function dramatically reduces the error, bringing the standard deviation for 64 block size down to 20% more than the uncorrected version with 32k block size:

void SanitizeCoeffs(float coeffs[], size_t numCoeffs)
{
switch (numCoeffs)
{
case 2:
{
if (abs(coeffs[0]) + abs(coeffs[1]) <= 3.0f)
return;

coeffs[0] = 2.0f;
coeffs[1] = -1.0f;

break;
}
}
}


Still leaves unanswered questions, but I think I can live with this, so probably won't investigate further.

I think there is a relation between LP order (which you supposed 2) and frame length. Selecting low order means just low frequency is predicted that's why block size with large length gives you better results. The low frequency better predicted in larger block size. I hope it helps.