matlab - 使用 FFT 在 Matlab 中自相关多个时间序列

Question

Matlab 函数 autocorr(文档 here)使用快速傅立叶变换 (FFT) 算法计算单个时间序列的样本自相关:

nFFT = 2^(nextpow2(length(y))+1);
F = fft(y-mean(y),nFFT);
F = F.*conj(F);
acf = ifft(F);
acf = acf(1:(numLags+1)); % Retain non-negative lags
acf = acf./acf(1); % Normalize
acf = real(acf);

假设我对同一个随机过程(时间序列)有多个实现，并且我想计算样本自相关。一种天真的方法是为每个时间序列调用 autocorr 并对每个滞后的相关性进行平均。然而，最好同时考虑由滞后 k 分隔的所有点对。我通过构建一个 (n-iLag)*nSamples 来实现这一点，所有对的矩阵由滞后大小 iLag 分隔，然后计算样本相关性。

% input: samplesMat (nSamples (times series) by n (time points) matrix)

sizeMat = size(samplesMat);
nSamples = sizeMat(1); 
n = sizeMat(2);

ACF = zeros(n,1);

muHat = mean(mean(samplesMat)); % sample mean

for iLag = 0:(n-1)
    index = 1;
    nPairs = (n-iLag)*nSamples;
    pairs = zeros(nPairs,2);
    for iSample = 1:nSamples
        for ix = 1:(n-iLag)
            pairs(index,1) = samplesMat(iSample,ix);
            pairs(index,2) = samplesMat(iSample,ix+iLag);
            index = index + 1;
        end
    end

    X2 = pairs(:,1)-muHat;
    Y2 = pairs(:,2)-muHat;
    ACF(iLag+1) = sum(X2.*Y2)/n*nSamples % calculate covariance
end
ACF = ACF/ACF(1); % divide by variance

对于单个时间序列 nSamples = 1，此代码提供与 autocorr 相同的输出，尽管速度明显较慢。有没有办法利用FFT算法计算多个时间序列的ACF？

Answer 1

在傅立叶域中取算术平均得到相同的结果

function [ acf ] = corrFFT(samplesMat,denomType,mu)
% Adaptation of inbuilt function autocorr to allow for multiple
% realisations, known mean and variable denominator.
%
% Input:
% samplesMat - nx by nSamples matrix of realisations
% denomType (optional) - 'cnst' or 'var'
% mu (optional) - if true mean is known

if nargin < 2
    denomType = 'cnst';
end
if nargin < 3
    mu = mean(mean(samplesMat));
end

[nx,~] = size(samplesMat);

nFFT = 2^(nextpow2(nx)+1);
F = fft(samplesMat-mu,nFFT);
F = F.*conj(F);
acf = ifft(mean(F,2));
acf = acf(1:nx);
if strcmp(denomType,'var')
    acf = (nx./[nx:-1:1]').*acf;
end

acf = acf./acf(1);
acf = real(acf);
end