# Dimensionality reduction of multiple signals using Fourier transform

I have $N$ recorded signals, $x$, each of which have been sampled 672 times across a time period of a week (= 15 min intervals). I will denote $x_{ij}$ as the $j$th sample for the $i$th recorded signal. The signals are a density representation of the activity across the week, meaning that

$\sum_{j=1}^{672}(x_{ij}) = 1$ for $i = 1,2,...,N$.

Also, the mean has been subtracted from each sampling point, meaning that

$\sum_{i=1}^{N}(x_{ij}) / N = 0$ for $j = 1,2,...,672$.

As I know that a lot of the information in the signals have a periodicity of a day, it would be obvious to represent the data using a Fourier transform as a way of dimensionality reduction. However, I've only ever seen the Fourier transform used for a single signal at a time. I want to apply Fourier transform as a dimensionality reduction technique across all of my signals and I therefore have the following questions:

1. How do I decompose each recorded signal into a linear combination of the same basic sinusoidals (= fixed amplitude, period and phase)?
2. How do I quantify the amount of variance captured across the population of recorded signals with the top $k$ basic sinusoidals?

Both general advice as well as advice for packages in Python would be great!

• Let me get this right, you have 672 variables, each with its own set of 30,000 time sampled points of the respective variable's value as it changes over time? It's a bit unclear, especially when you're referring to a Fourier transform for a single observation at a time. – Falimond Oct 10 '15 at 1:44
• Please see my edit - I've tried to be more precise in my description of the problem. – pir Oct 10 '15 at 10:14
• FT is not dimension reducing in and of itself. It is a univariate procedure that decomposes any periodicity in a time series signal, xj (not xij), into frequencies, amplitude (intensity), cosines, etc. You can rerank the Nyquist frequencies (N/2) by amplitude and select the "top k basic signals" using some heuristic criteria based on a PCA of the full Nyquist set of frequencies. This would address question #2. Does your question #1 then become one of clustering the signals based on the FT output? That is certainly possible. Of course, if your data is aperiodic, then FT is useless. – Mike Hunter Oct 10 '15 at 10:54
• @DJohnson, thanks. I am simply surprised that FT hasn't been expanded in a way to be used as a dimensionality reduction technique for multiple time series, $x_{ij}$. If I apply the FT to each $x_j$ individually then I won't be able to compare the coefficients and the procedure won't make sense at all as all of the signals will have different representations. – pir Oct 10 '15 at 11:03
• A small example.. Say I have $N$ recordings of a certain type of time series signal. The underlying true signal is simply 3 sinusoidal functions and a lot of random noise. What would be the standard way of reducing the dimensionality of all of these $N$ recordings at once? – pir Oct 10 '15 at 11:05

You can think of this as a sparse signal approximation problem. First, let's just consider the case for a single signal.

Say we want to approximate a signal $x$. We have a signal dictionary $D$, which is a matrix where each column is a basis function. In your case, the basis functions would be sinusoids. We can approximate $x$ as a linear combination of basis functions: $\hat{x} = D w$, where $w$ is a weight vector. We can additionally impose the constraint that $w$ be sparse. That is, it should contain few non-zero entries. We can measure the sparsity using the $l_0$ norm, $\|w\|_0$, which just counts the number of non-zero entries. Given the signal dictionary, we want to find the weights that give a good approximation of the signal:

$$\underset{w}{\min} \|x - Dw\|_2^2 \quad s.t. \|w\|_0 \le k$$

This tries to find a good least-squares approximation to $x$, subject to the constraint that we use no more than $k$ basis functions. There's a big literature on solving this problem. It turns out that solving it is NP-hard, but people have come up with many good approximations. One solution is to relax the $l_0$ norm to the $l_1$ norm. This is called basis pursuit denoising in the signal processing literature and LASSO in the statistics literature, and has proven very successful. Another class of solutions is to keep the $l_0$ norm, but use a heuristic search strategy that finds a locally optimal solution instead of the global minimum (or that finds the global minimum under certain conditions). One such method is called Orthogonal Matching Pursuit.

Your case has an interesting twist, in that there are multiple signals to approximate, and you want to impose the condition that the sparsity pattern is the same for all signals (i.e. they all use the same set of basis functions). Also, the sparsity pattern is not known ahead of time.

Writing out the problem, let's replace the signal vector with a signal matrix $X$, where each column is a signal. Replace the weight vector with a weight matrix $W$, where column $i$ gives the weights for signal $i$. Replace the $l_2$ norm with the Frobenius norm, which is just the analog for matrices.

$$\underset{W}{\min} \|X - DW\|_F^2 \quad$$

Subject to the constraint that $W$ contains no more than $k$ nonzero rows.

This finds a least squares approximation to the entire set of signals, using a small set of basis functions (no more than $k$ of them) that are shared across all signals. Of course, the caveat about NP-completeness still applies, so we'll have to settle for approximation algorithms.

One method that solves this problem is called Simultaneous Orthogonal Matching Pursuit:

Gilbert and Strauss (2005). Simultaneous sparse approximation via greedy pursuit http://authors.library.caltech.edu/9043/1/TROicassp05.pdf

There are probably others too.