# CorrelationFunction for vectors

Mathematica 9 has a new CorrelationFunction. Sadly the joy does not last long, as I can't get it to work with vectors. I would like to make a velocity correlation function. That is, given a list of n vectors I would like to calculate:

For example, this is almost instantaneous:

testdata = RandomVariate[BinormalDistribution[{10, 20}, 0.4], {10000}];
AbsoluteTiming[
x = CorrelationFunction[testdata[[All, 1]], {0, 1000}];
y = CorrelationFunction[testdata[[All, 2]], {0, 1000}];
]
ListPlot[{x, y}, PlotRange -> {0, 1}]


(the test data here is not the best as it is by definition uncorrelated)

Here I would like the autocorrelation of 2D vectors evenly spaced in time. But this does not work:

r = CorrelationFunction[testdata, {0, 100}];


How can I make an autocorrelation for vectors? Or do I have to go the manual way?

ClearAll[AbsoluteAutocorrelationT, AutocorrelationT]

AbsoluteAutocorrelationT[l_, 0] := Total[(#.#) & /@ l];
AbsoluteAutocorrelationT[l_, n_] :=
Total@Table[l[[i]].l[[i + n]], {i, Length[l] - n}];
AutocorrelationT[l_, range_] := Block[{t0, ti},
t0 = AbsoluteAutocorrelationT[l, 0];
ti = AbsoluteAutocorrelationT[l, #] & /@ (Range @@ range);
(ti/t0)];

AbsoluteTiming[r = AutocorrelationT[testdata, {1, 1000}];]


Which is 200 times slower. I know I can compile this (but than I have to worry about overflows in the sum, etc.)

-
Easy to guess why this isn't supported directly, but it does look like something of an oversight not to at least mention it in the documentation. One obvious way around it would be to interpolate and resample the data. – Oleksandr R. Feb 25 '13 at 17:45
@OleksandrR. and what is the reason (must be a bit slow today:) – Ajasja Feb 25 '13 at 17:48
Correlation is probably being done using an FFT, which requires evenly spaced data. The good news is an FFT is fast, so you can happily upsample your data significantly, especially if you're able to arrive at a power of two length without absurdly excessive padding. – Oleksandr R. Feb 25 '13 at 17:51
No, the 1-d case is assumed to be evenly spaced to start with. Your 2-d input is assumed to be of the form {{t1, x1}, ..., {tn, xn}}. If I misunderstood your question and you want the n-dimensional autocorrelation for values sampled on a regular grid, rather than that for 1-d but unevenly spaced data, you can get it yourself using ListCorrelate. – Oleksandr R. Feb 25 '13 at 17:59
Given your code, isn't this what you are trying to do? absCorr[lists_, n_] := Length@lists AbsoluteCorrelationFunction[#, n] & /@ Transpose@lists // Total – Rojo Feb 25 '13 at 19:56

When each $v_i = (v_{i,1}, v_{i,2}, \ldots, v_{i,k})$ is a $k$-vector we may interchange the summations to obtain

$$\sum_{i=1}^{n-t} v_i\cdot v_{i+t} = \sum_{j=1}^k\sum_{i=1}^{n-t} v_{i,j} v_{i+t,j},$$

whence the desired correlation is the sum of the correlations of the components. Thus:

kernel = PDF[NormalDistribution[0, 1]] /@ Range[-5, 5, .005]; (* Induces autocorrelation *)
testdata = ListConvolve[kernel, #] & /@
Transpose[RandomVariate[BinormalDistribution[{10, 20}, 0.4], {10^4}]];
c = Sum[ListCorrelate[v, v, {1, 1}], {v, testdata}]; // AbsoluteTiming


$\{0.0030002, \text{Null}\}$

ListLinePlot[c[[;; Length[c]/2]]/c[[1]]]


The timings scale as expected ($k=2$ shown):

-

The following is (only) six times faster:

l = RandomVariate[BinormalDistribution[{10, 20}, 0.4], {10000}];

c[l_, 0] := ListCorrelate[l, l, {1, -1}, {}, Dot, Plus, 1];
c[l_, s_List] := Module[{m = c[l, 0][[1]]},
ListCorrelate[l[[;;Length@l - #]], l[[1+#;;]], {1, -1}, {}, Dot, Plus, 1] / m & /@ Range @@ s];

ListLinePlot@Table[ AbsoluteTiming[AutocorrelationT[l, {1, n}];][[1]]/
AbsoluteTiming[c[l, {1, n}];][[1]], {n, 30}]


-