xcov

Cross-covariance

Syntax

  • [c,lags] = xcov(___)

Description

c = xcov(x,y) returns the cross-covariance of two discrete-time sequences, x and y. Cross-covariance measures the similarity between x and shifted (lagged) copies of y as a function of the lag. If x and y have different lengths, the function appends zeros at the end of the shorter vector so it has the same length as the other.

example

c = xcov(x) returns the autocovariance sequence of x. If x is a matrix, then c is a matrix whose columns contain the autocovariance and cross-covariance sequences for all combinations of the columns of x.

example

c = xcov(___,maxlag) limits the lag range from –maxlag to maxlag. This syntax accepts one or two input sequences. maxlag defaults to N – 1.

example

c = xcov(___,scaleopt) additionally specifies a normalization option for the cross-covariance or autocovariance. Any option other than 'none' (the default) requires x and y to have the same length.

[c,lags] = xcov(___) also outputs a vector with the lags at which the covariances are computed.

Examples

expand all

Cross-Covariance of Two Shifted Signals

Create $s$, a two-channel 150-sample signal consisting of $s_1$, a uniform random sequence, and $s_2$, a copy of $s_1$ shifted circularly by 50 samples. Reset the random number generator for reproducible results. Plot the sequences.

rng default
shft = 50;

s1 = rand(150,1);
s2 = circshift(s1,[shft 0]);
s = [s1 s2];

subplot(2,1,1)
plot(s1)
title('s_1')
subplot(2,1,2)
plot(s2)
title('s_2')
hold on
plot([shft shft],[0 1])

Compute biased estimates of the autocovariance and mutual cross-covariance sequences. The output array is organized as ${\texttt c}=\pmatrix{c_{s_1s_1}&c_{s_1s_2}&c_{s_2s_1}&c_{s_2s_2}\cr}$. Plot the result. The maxima at $\mp50$ and $\pm100$ are a result of the circular shift.

[c,lg] = xcov(s,'biased');

figure
plot(lg,c)
legend('c_{s_1s_1}','c_{s_1s_2}','c_{s_2s_1}','c_{s_2s_2}')

Change the normalization so that the cross-covariance and autocovariance sequences are unity at zero lag. Plot each sequence in its own subplot.

[c,lg] = xcov(s,'coeff');

for a = 1:2
    for b = 1:2
        nm = 2*(a-1)+b;
        subplot(2,2,nm)
        plot(lg,c(:,nm))
        title(sprintf('c_{s_%ds_%d}',a,b))
        axis([-150 150 -0.2 1])
    end
end

Autocovariance of White Gaussian Noise

Display the estimated autocovariance of white Gaussian noise, $c_{ww}(m)$, for $-10 \le m \le 10$. Reset the random number generator for reproducible results. Normalize the sequence so that it is unity at zero lag.

rng default

ww = randn(1000,1);
[cov_ww,lags] = xcov(ww,10,'coeff');

stem(lags,cov_ww)

GPU Acceleration for Autocovariance Sequence Estimation

This example requires Parallel Computing Toolbox™ software and a CUDA-enabled NVIDIA GPU with compute capability 1.3 or above. See GPU System Requirements for details.

Create a signal consisting of a 10 Hz sine wave in additive noise, sampled at 1 kHz. Use gpuArray to create a gpuArray object stored on your computer's GPU.

t = 0:0.001:10-0.001;
x = cos(2*pi*10*t) + randn(size(t));
X = gpuArray(x);

Compute the autocovariance sequence to lag 200.

[xc,lags] = xcov(X,200);

The output, xc, is a gpuArray object.

Use gather to transfer the data from the GPU to the MATLAB® workspace as a double-precision vector.

xc = gather(xc);

Input Arguments

expand all

x — Input arrayvector | matrix | gpuArray object

Input array, specified as a vector, a matrix, or a gpuArray object.

See GPU Computing and GPU System Requirements for details on using xcov with gpuArray objects.

Example: sin(2*pi*(0:9)/10) + randn([1 10])/10 specifies a noisy sinusoid as a row vector.

Example: sin(2*pi*[0.1;0.3]*(0:39))' + randn([40 2])/10 specifies a two-channel noisy sinusoid.

Example: gpuArray(sin(2*pi*(0:9)/10) + randn([1 10])/10) specifies a noisy sinusoid as a gpuArray object.

Data Types: single | double
Complex Number Support: Yes

y — Input arrayvector | gpuArray object

Input array, specified as a vector or a gpuArray object.

Data Types: single | double
Complex Number Support: Yes

maxlag — Maximum laginteger scalar

Maximum lag, specified as an integer scalar. If you specify maxlag, the returned cross-covariance sequence ranges from –maxlag to maxlag. If you do not specify maxlag, the lag range equals 2N – 1, where N is the greater of the lengths of x and y.

Data Types: single | double

scaleopt — Normalization option'none' (default) | 'biased' | 'unbiased' | 'coeff'

Normalization option, specified as one of these strings:

  • 'none' — Raw, unscaled cross-covariance. This is the only allowed option when x and y have different lengths.

  • 'biased' — Biased estimate of the cross-covariance.

  • 'unbiased' — Unbiased estimate of the cross-covariance.

  • 'coeff' — Normalizes the sequence so that the autocovariances at zero lag equal 1.

Data Types: char

Output Arguments

expand all

c — Cross-covariance or autocovariance sequencevector | matrix | gpuArray object

Cross-covariance or autocovariance sequence, returned as a vector, a matrix, or a gpuArray object.

If x is an M × N signal matrix representing N channels in its columns, then xcov(x) returns a (2M – 1) × N2 matrix with the autocovariances and mutual cross-covariances of the channels of x. If you specify maxlag, then c has size (2 × maxlag – 1) × N2.

For example, if S is a three-channel signal, S=(x1x2x3), then the result of C = xcov(S) is organized as

c=(cx1x1cx1x2cx1x3cx2x1cx2x2cx2x3cx3x1cx3x2cx3x3).

lags — Lag indicesvector

Lag indices, returned as a vector.

More About

expand all

Cross-Covariance and Autocovariance

xcov computes the mean of its inputs, subtracts the mean, and then calls xcorr. xcov does not check for any errors other than the correct number of input arguments. Instead, it relies on the error checking in xcorr.

The result of xcov can be interpreted as an estimate of the covariance between two random sequences or as the deterministic covariance between two deterministic signals.

The true cross-covariance sequence of two jointly stationary random processes, xn and yn, is the cross-correlation of mean-removed sequences,

ϕxy(m)=E{(xn+mμx)(ynμy))},

where μx and μy are the mean values of the two stationary random processes, the asterisk denotes complex conjugation, and E is the expected value operator. xcov can only estimate the sequence because, in practice, only a finite segment of one realization of the infinite-length random process is available.

By default, xcov computes raw covariances with no normalization:

cxy(m)={n=0Nm1(xn+m1Ni=0N1xi)(yn1Ni=0N1yi),m0,cyx(m),m<0.

The output vector, c, has elements given by

c(m)=cxy(mN),m=1,,2N1.

The covariance function requires normalization to estimate the function properly. You can control the normalization of the correlation by using the input argument scaleopt.

References

[1] Orfanidis, Sophocles J. Optimum Signal Processing: An Introduction. 2nd Edition. New York: McGraw-Hill, 1996.

[2] Larsen, Jan. "Correlation Functions and Power Spectra." November, 2009. http://www2.imm.dtu.dk/pubdb/views/edoc_download.php/4932/pdf/imm4932.pdf

See Also

| | | |

Was this topic helpful?