Message Boards Message Boards

GROUPS:

[✓]  Simplest method for computing cross correlation?

Posted 11 months ago
986 Views
|
2 Replies
|
1 Total Likes
|

I have two lists, x and y, of equal length, Nx, and I wish to construct the cross-correlation of the two. To compute the autocorrelation of just x (or y) I simply used

Rxx = CorrelationFunction[x, {Nx - 1}]

For the cross-correlation, I first attempted to construct a function using the following:

Rxy[n_] = Sum[x[[m]]*y[[m + n]], {m, 1, Nx - n}],

where n denotes the prescribed lag. This resulted in an error message "The expression m cannot be used as a part specification."

I'm not sure as to what this error means. Also, is there a simple way to use the ListConvolve function to create the desired cross-correlation? (I'm a bit confused by the description at http://reference.wolfram.com/language/ref/ListConvolve.html)

2 Replies

Hi,

you can use something like

ListLinePlot[CorrelationFunction[TemporalData[{Transpose[{x, y}]}, Automatic], {-200, 200}]["Values"][[;; , 1, 2]], PlotRange -> {All, {-1, 1}}]

I guess. It does something slightly different from your low level implementation.

You can get yours to work (fix standard deviations and normalisations etc) like this:

ClearAll["Global`*"]
Nx = 10000;
x = RandomVariate[NormalDistribution[], Nx];
y = RandomVariate[NormalDistribution[], Nx];
Rxy[n_] := 1/(Nx-n) Sum[x[[m]]*y[[m + n]], {m, 1, Nx - n}]
ListLinePlot[Table[Rxy[z], {z, 0, 200, 1}], PlotRange -> {All, {-1, 1}}]

enter image description here

It is somewhat more elegant to do it like this:

xcorr[m_] := Mean[Drop[x, -m]*Drop[y, m]]
ListLinePlot[ Table[xcorr[m], {m, 1, 200}], PlotRange -> {All, {-1, 1}}]

which gives the same image as above. Don't forget that for the correlation you need to divide the expected value by the standard deviations of the respective time series and before that subtract the mean. You can do that my hand, but there are functions for that:

xcorr[m_] := Mean[Drop[Standardize[x], -m]*Drop[Standardize[y], m]]
ListLinePlot[Table[xcorr[m], {m, 1, 200}], PlotRange -> {All, {-1, 1}}]

This doesn't make much of a difference for the time series I have chosen above, but is a bit more general. There are more things to say, but that's the basic idea. Also this is numerically probably not the most efficient way.

Cheers,

Marco

Wonderful! The Standardize and Drop commands will help me in many other situations involving convolutions.

Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract