Time-frequency ridges from wavelet synchrosqueezing

`[___]= wsstridge(___,`

returns the time-frequency ridge with additional options specified by one or more
`Name,Value`

)`Name,Value`

pair arguments.

The function uses a penalized forward-backward greedy algorithm to extract the maximum-energy
ridges from a time-frequency matrix. The algorithm finds the maximum time-frequency ridge by
minimizing –ln *A* at each time point, where *A* is the absolute value of the
matrix. Minimizing –ln *A* is equivalent to maximizing the value of *A*. The algorithm
optionally constrains jumps in frequency with a penalty that is proportional to the distance
between frequency bins.

The following example illustrates the time-frequency ridge algorithm using a penalty that is
two times the distance between frequency bins. Specifically, the distance between the elements
`(j,k)`

and `(m,n)`

is defined as
`(j-m)`

. The time-frequency matrix has three
frequency bins and three time steps. The matrix columns correspond to time steps, and the matrix
rows correspond to frequency bins. The values in the second row represent a sine wave.^{2}

Suppose you have the matrix:

1 4 4 2 2 2 5 5 4

Update the value for the (1,2) element as follows.

Leave the values at the first time point unaltered. Begin the algorithm with the (1,2) element of the matrix, which presents the first frequency bin at the second time point. The bin value is 4. Penalize the values in the first column based on their distance from the (1,2) element. Applying the penalty to the first column produces

original value + penalty × distance 1 + 2 × 0 = 1 2 + 2 × 1 = 4 5 + 2 × 4 = 13

The minimum value of the first column is 1, which is in bin 1.1 4 4 2 13 5

Add the minimum value in column 1 to the current bin value, 4. The updated value for (1,2) becomes 5, which came from bin 1.

Update the values for the remaining elements in column 2 as follows.

Recompute the original column 1 values with the penalty factor using the same process as in Step 2a. Obtain the remaining second column values using the same process as in Step 2b. For example, when updating the (2,2) element, which has bin value 2, applying the penalty to the column yields

Add the minimum value, 2, to the current bin value. The updated value for (2,2) becomes 4. After updating the (3,2) element, the matrix isoriginal value + penalty × distance 1 + 2 × 1 = 3 2 + 2 × 0 = 2 5 + 2 × 1 = 7

Only the second column has been updated. The subscripts indicate the index of the bin in the previous column from which a value came.1 5

_{(1)}4 2 4_{(2)}2 5 9_{(2)}4Repeat Step 2 for the third column. But now the penalty is applied to the updated second column. For example, when updating the (1,3) element, the penalty is

The minimum value, 5, which is in the first bin, is added to the (1,3) bin value. After updating all the values in the third column, the final matrix is5 + 2 × 0 = 5 4 + 2 × 1 = 6 9 + 2 × 4 = 17

1 5

_{(1)}9_{(1)}2 4_{(2)}6_{(2)}5 9_{(2)}10_{(2)}Starting at the last column of the matrix, find the minimum value. Walk back in time through the matrix by going from the current bin to the origin of that bin at the previous time point. Keep track of the bin indices, which form the path composing the ridge. The algorithm smooths the transition by using the origin bin instead of the bin with the minimum value. For this example, the ridge indices are

`2`

,`2`

,`2`

, which matches the energy path of the sine wave in row 2 of the matrix shown in Step 1.

If you are extracting multiple ridges, the algorithm removes the first ridge from the time-frequency matrix and repeats the process.

[1] Daubechies, I., J. Lu, and H.-T. Wu. "Synchrosqueezed wavelet transforms: an
empirical mode decomposition-like tool." *Applied and Computational Harmonic
Analysis*. Vol. 30, Number 2, 2011, pp. 243–261.

[2] Thakur, G., E. Brevdo, N. S. Fučkar, and H.-T. Wu. "The Synchrosqueezing algorithm
for time-varying spectral analysis: Robustness properties and new paleoclimate applications."
*Signal Processing*. Vol. 93, Number 4, 2013, pp. 1079–1094.