Digital Signal Processing

Representing Linear Combinations

The dotProduct performs the math of a linear combination. A linear combination has the following form:

(a1*v1)+(a2*v2)...

In the above example a1 and a2 are random variables that change. v1 and v2 are constant values.

When computing the dot product the elements of two vectors are multiplied together and the results are added. If the first vector contains random variables and the second vector contains constant values then the dot product is performing a linear combination.

This scenario comes up again and again in machine learning. For example both linear and logistic regression solve for a vector of constant weights. In order to perform a prediction, a dot product is calculated between a random observation vector and the constant weight vector. That dot product is a linear combination because one of the vectors holds constant weights.

Lets look at simple example of how a linear combination can be used to find the mean of a vector of numbers.

In the example below two arrays are set to variables a and b and then operated on by the dotProduct function. The output of the dotProduct function is set to variable c.

The mean function is then used to compute the mean of the first array which is set to the variable d.

Both the dot product and the mean are included in the output.

When we look at the output of this expression we see that the dot product and the mean of the first array are both 30.

The dotProduct function calculated the mean of the first array.

let(echo="c, d",
    a=array(10, 20, 30, 40, 50),
    b=array(.2, .2, .2, .2, .2),
    c=dotProduct(a, b),
    d=mean(a))

When this expression is sent to the /stream handler it responds with:

{
  "result-set": {
    "docs": [
      {
        "c": 30,
        "d": 30
      },
      {
        "EOF": true,
        "RESPONSE_TIME": 0
      }
    ]
  }
}

To get a better understanding of how the dot product calculated the mean we can perform the steps of the calculation using vector math and look at the output of each step.

In the example below the ebeMultiply function performs an element-by-element multiplication of two arrays. This is the first step of the dot product calculation. The result of the element-by-element multiplication is assigned to variable c.

In the next step the add function adds all the elements of the array in variable c.

Notice that multiplying each element of the first array by .2 and then adding the results is equivalent to the formula for computing the mean of the first array. The formula for computing the mean of an array is to add all the elements and divide by the number of elements.

The output includes the output of both the ebeMultiply function and the add function.

let(echo="c, d",
    a=array(10, 20, 30, 40, 50),
    b=array(.2, .2, .2, .2, .2),
    c=ebeMultiply(a, b),
    d=add(c))

When this expression is sent to the /stream handler it responds with:

{
  "result-set": {
    "docs": [
      {
        "c": [
          2,
          4,
          6,
          8,
          10
        ],
        "d": 30
      },
      {
        "EOF": true,
        "RESPONSE_TIME": 0
      }
    ]
  }
}

In the example above two arrays were combined in a way that produced the mean of the first. In the second array each value was set to .2. Another way of looking at this is that each value in the second array is applying the same weight to the values in the first array. By varying the weights in the second array we can produce a different result. For example if the first array represents a time series, the weights in the second array can be set to add more weight to a particular element in the first array.

The example below creates a weighted average with the weight decreasing from right to left. Notice that the weighted mean of 36.666 is larger than the previous mean which was 30. This is because more weight was given to last element in the array.

let(echo="c, d",
    a=array(10, 20, 30, 40, 50),
    b=array(.066666666666666,.133333333333333,.2, .266666666666666, .33333333333333),
    c=ebeMultiply(a, b),
    d=add(c))

When this expression is sent to the /stream handler it responds with:

{
  "result-set": {
    "docs": [
      {
        "c": [
          0.66666666666666,
          2.66666666666666,
          6,
          10.66666666666664,
          16.6666666666665
        ],
        "d": 36.66666666666646
      },
      {
        "EOF": true,
        "RESPONSE_TIME": 0
      }
    ]
  }
}

Representing Correlation

Often when we think of correlation, we are thinking of Pearson correlation in the field of statistics. But the definition of correlation is actually more general: a mutual relationship or connection between two or more things. In the field of digital signal processing the dot product is used to represent correlation. The examples below demonstrates how the dot product can be used to represent correlation.

In the example below the dot product is computed for two vectors. Notice that the vectors have different values that fluctuate together. The output of the dot product is 190, which is hard to reason about because it’s not scaled.

let(echo="c, d",
    a=array(10, 20, 30, 20, 10),
    b=array(1, 2, 3, 2, 1),
    c=dotProduct(a, b))

When this expression is sent to the /stream handler it responds with:

{
  "result-set": {
    "docs": [
      {
        "c": 190
      },
      {
        "EOF": true,
        "RESPONSE_TIME": 0
      }
    ]
  }
}

One approach to scaling the dot product is to first scale the vectors so that both vectors have a magnitude of 1. Vectors with a magnitude of 1, also called unit vectors, are used when comparing only the angle between vectors rather than the magnitude. The unitize function can be used to unitize the vectors before calculating the dot product.

Notice in the example below the dot product result, set to variable e, is effectively 1. When applied to unit vectors the dot product will be scaled between 1 and -1. Also notice in the example cosineSimilarity is calculated on the unscaled vectors and the answer is also effectively 1. This is because cosine similarity is a scaled dot product.

let(echo="e, f",
    a=array(10, 20, 30, 20, 10),
    b=array(1, 2, 3, 2, 1),
    c=unitize(a),
    d=unitize(b),
    e=dotProduct(c, d),
    f=cosineSimilarity(a, b))

When this expression is sent to the /stream handler it responds with:

{
  "result-set": {
    "docs": [
      {
        "e": 0.9999999999999998,
        "f": 0.9999999999999999
      },
      {
        "EOF": true,
        "RESPONSE_TIME": 0
      }
    ]
  }
}

If we transpose the first two numbers in the first array, so that the vectors are not perfectly correlated, we see that the cosine similarity drops. This illustrates how the dot product represents correlation.

let(echo="c, d",
    a=array(20, 10, 30, 20, 10),
    b=array(1, 2, 3, 2, 1),
    c=cosineSimilarity(a, b))

When this expression is sent to the /stream handler it responds with:

{
  "result-set": {
    "docs": [
      {
        "c": 0.9473684210526314
      },
      {
        "EOF": true,
        "RESPONSE_TIME": 0
      }
    ]
  }
}

Moving Average Function

Before looking at an example of convolution its useful to review the movingAvg function. The moving average function computes a moving average by sliding a window across a vector and computing the average of the window at each shift. If that sounds similar to convolution, that’s because the movingAvg function is syntactic sugar for convolution.

Below is an example of a moving average with a window size of 5. Notice that original vector has 13 elements but the result of the moving average has only 9 elements. This is because the movingAvg function only begins generating results when it has a full window. In this case because the window size is 5 so the moving average starts generating results from the 4^th index of the original array.

let(a=array(1, 2, 3, 4, 5, 6, 7, 6, 5, 4, 3, 2, 1),
    b=movingAvg(a, 5))

When this expression is sent to the /stream handler it responds with:

{
  "result-set": {
    "docs": [
      {
        "b": [
          3,
          4,
          5,
          5.6,
          5.8,
          5.6,
          5,
          4,
          3
        ]
      },
      {
        "EOF": true,
        "RESPONSE_TIME": 0
      }
    ]
  }
}

Convolutional Smoothing

The moving average can also be computed using convolution. In the example below the conv function is used to compute the moving average of the first array by applying the second array as the filter.

Looking at the result, we see that it is not exactly the same as the result of the movingAvg function. That is because the conv pads zeros to the front and back of the first vector so that the window size is always full.

let(a=array(1, 2, 3, 4, 5, 6, 7, 6, 5, 4, 3, 2, 1),
    b=array(.2, .2, .2, .2, .2),
    c=conv(a, b))

When this expression is sent to the /stream handler it responds with:

{
  "result-set": {
    "docs": [
      {
        "c": [
          0.2,
          0.6000000000000001,
          1.2,
          2.0000000000000004,
          3.0000000000000004,
          4,
          5,
          5.6000000000000005,
          5.800000000000001,
          5.6000000000000005,
          5.000000000000001,
          4,
          3,
          2,
          1.2000000000000002,
          0.6000000000000001,
          0.2
        ]
      },
      {
        "EOF": true,
        "RESPONSE_TIME": 0
      }
    ]
  }
}

We achieve the same result as the movingAvg function by using the copyOfRange function to copy a range of the result that drops the first and last 4 values of the convolution result. In the example below the precision function is also also used to remove floating point errors from the convolution result. When this is added the output is exactly the same as the movingAvg function.

let(a=array(1, 2, 3, 4, 5, 6, 7, 6, 5, 4, 3, 2, 1),
    b=array(.2, .2, .2, .2, .2),
    c=conv(a, b),
    d=copyOfRange(c, 4, 13),
    e=precision(d, 2))

When this expression is sent to the /stream handler it responds with:

{
  "result-set": {
    "docs": [
      {
        "e": [
          3,
          4,
          5,
          5.6,
          5.8,
          5.6,
          5,
          4,
          3
        ]
      },
      {
        "EOF": true,
        "RESPONSE_TIME": 0
      }
    ]
  }
}

Sine Wave Interpolation, Extrapolation

The oscillate function returns a function which can be used by the predict function to interpolate or extrapolate a sine wave. The example below extrapolates the sine wave to an x-axis sequence of 0-256.

let(a=oscillate(1, 0.28, 1.57),
    b=predict(a, sequence(256, 0, 1)))

The extrapolated sine wave is plotted below:

Complex Result

The fft function performs the discrete Fourier Transform on a vector of real data. The result of the fft function is returned as complex numbers. A complex number has two parts, real and imaginary. The imaginary part of the complex number is ignored in the examples below, but there are many tutorials on the FFT and that include complex numbers available online.

But before diving into the examples it is important to understand how the fft function formats the complex numbers in the result.

The fft function returns a matrix with two rows. The first row in the matrix is the real part of the complex result. The second row in the matrix is the imaginary part of the complex result.

The rowAt function can be used to access the rows so they can be processed as vectors. This approach was taken because all of the vector math functions operate on vectors of real numbers. Rather then introducing a complex number abstraction into the expression language, the fft result is represented as two vectors of real numbers.

Fast Fourier Transform Examples

In the first example the fft function is called on the sine wave used in the autocorrelation example.

The results of the fft function is a matrix. The rowAt function is used to return the first row of the matrix which is a vector containing the real values of the fft response.

The plot of the real values of the fft response is shown below. Notice there are two peaks on opposite sides of the plot. The plot is actually showing a mirrored response. The right side of the plot is an exact mirror of the left side. This is expected when the fft is run on real rather than complex data.

Also notice that the fft has accumulated significant power in a single peak. This is the power associated with the specific frequency of the sine wave. The vast majority of frequencies in the plot have close to 0 power associated with them. This fft shows a clear signal with very low levels of noise.

let(a=oscillate(1, 0.28, 1.57),
    b=predict(a, sequence(256, 0, 1)),
    c=fft(b),
    d=rowAt(c, 0))

In the second example the fft function is called on a vector of random data similar to one used in the autocorrelation example. The plot of the real values of the fft response is shown below.

Notice that in is this response there is no clear peak. Instead all frequencies have accumulated a random level of power. This fft shows no clear sign of signal and appears to be noise.

let(a=sample(uniformDistribution(-1.5, 1.5), 256),
    b=fft(a),
    c=rowAt(b, 0))

In the third example the fft function is called on the same signal hidden within noise that was used for the autocorrelation example. The plot of the real values of the fft response is shown below.

Notice that there are two clear mirrored peaks, at the same locations as the fft of the pure signal. But there is also now considerable noise on the frequencies. The fft has found the signal and but also shows that there is considerable noise along with the signal.

let(a=oscillate(1, 0.28, 1.57),
    b=predict(a, sequence(256, 0, 1)),
    c=sample(uniformDistribution(-1.5, 1.5), 256),
    d=ebeAdd(b, c),
    e=fft(d),
    f=rowAt(e, 0))