Вычисление автокорреляции с использованием БПФ ⇐ C#

[Anonymous]

Это был мой первоначальный класс для вычисления автокорреляции.
Он написан на C# 4.7.0.

Код: Выделить всё

using System;
using System.Collections.Generic;
using System.Linq;

namespace MyNameSpace
{
public sealed class TimeSeriesAnalysis
{
private readonly IEnumerable _source;
private int _count;

public TimeSeriesAnalysis(IEnumerable vectorList)
{
_source = vectorList;
_count = vectorList.Count();
}

public double AutoCorr(int lag)
{
int n = _count;

if (n == 0 || lag >= n || lag < 0)
throw new ArgumentException($"Lag must be between 0 and {n - 1}, and the list cannot be empty.");

double sum = 0;

int targetIndex = n - lag;

for (int i = 0; i < targetIndex; i++)
{
var vecA = _source.Skip(i).Take(1).First();

var vecB = _source.Skip(i + lag).Take(1).First();

sum += vecA.Dot(vecB);
}

return sum / (n - lag);
}

public IEnumerable GetCResults(int maxLag, double c0)
{
for (int lag = 0; lag  (double)lag);

return new PairReturn(values, GetCResults(maxLag, c0));
}

}
}

Этот класс работает хорошо.

Код: Выделить всё

using System;
using System.Collections.Generic;
using Microsoft.VisualStudio.TestTools.UnitTesting;
using MyNameSpace;

namespace MyNameSpace____unit_test
{
[TestClass]
public class TimeSeriesAnalysisVec3AutoCorrelationMemoryLess____unit_test
{

[TestInitialize]
public void Initialize()
{
}

[TestMethod]
public void GetAutoCorrelationTestMethod()
{
IList vect3List = new List()
{
new Vector3(1,1,1),
new Vector3(2,2,2),
new Vector3(3,3,3),
new Vector3(4,4,4),
new Vector3(5,5,5)
};

var timeSeriesAnalysis = new TimeSeriesAnalysis(vect3List);

double c0 = timeSeriesAnalysis.AutoCorr(0);//33
double c1 = timeSeriesAnalysis.AutoCorr(1);//30
double c2 = timeSeriesAnalysis.AutoCorr(2);//26
double c3 = timeSeriesAnalysis.AutoCorr(3);//21
double c4 = timeSeriesAnalysis.AutoCorr(4);//15

Assert.AreEqual(33, c0);
Assert.AreEqual(30, c1);
Assert.AreEqual(26, c2);
Assert.AreEqual(21, c3);
Assert.AreEqual(15, c4);
}

[TestMethod]
public void GetAutoCorrelationPointsTestMethod()
{
List vect3List = new List()
{
new Vector3(1,1,1),
new Vector3(2,2,2),
new Vector3(3,3,3),
new Vector3(4,4,4),
new Vector3(5,5,5)
};

var timeSeriesAnalysis = new TimeSeriesAnalysis(vect3List);

var autocorr = timeSeriesAnalysis.GetAutoCorrelationPoints(10);

IList tausList = new List(autocorr.XData);
IList  autocorrList = new List(autocorr.YData);

Assert.AreEqual(5, tausList.Count);
Assert.AreEqual(5, autocorrList.Count);

Assert.AreEqual(1.0, autocorrList[0]);
Assert.AreEqual(0.90909090909090906, autocorrList[1]);
Assert.AreEqual(0.78787878787878785, autocorrList[2]);
Assert.AreEqual(0.63636363636363635, autocorrList[3]);
Assert.AreEqual(0.45454545454545453, autocorrList[4]);
}
}
}

Затем мне пришлось переписать этот класс, используя БПФ (быстрое преобразование Фурье).

Код: Выделить всё

using System;
using System.Collections.Generic;
using System.Linq;
using MathNet.Numerics;
using MathNet.Numerics.IntegralTransforms;

namespace MyNameSpace
{
public sealed class TimeSeriesAnalysisFFT
{
private readonly List _source;
private int _count;

public TimeSeriesAnalysisFFT(IEnumerable vectorList)
{
_source = vectorList.ToList();
_count = _source.Count;
}

public double AutoCorr(int lag)
{
int n = _count;

if (n == 0 || lag >= n || lag < 0)
throw new ArgumentException($"Lag must be between 0 and {n - 1}, and the list cannot be empty.");

// Using FFT to compute autocorrelation
var realValues = _source.Select(v => v.X).ToArray();

// Remove mean to center the data
double mean = realValues.Average();
for (int i = 0; i < n; i++)
{
realValues[i] -= mean;
}

var complexValues = new Complex32[n];
for (int i = 0; i < n; i++)
{
complexValues[i] = new Complex32((float)realValues[i], 0);
}

// Forward FFT
Fourier.Forward(complexValues, FourierOptions.Matlab);

// Compute power spectrum (magnitude squared)
for (int i = 0; i < n; i++)
{
complexValues[i] = complexValues[i] * Complex32.Conjugate(complexValues[i]);
}

// Inverse FFT to get the autocorrelation
Fourier.Inverse(complexValues, FourierOptions.Matlab);

// Extract the real part and normalize by number of samples
var autocorr = complexValues.Select(c => c.Real / n).ToArray();

// Return the normalized value for the given lag
double variance = realValues.Sum(x => x * x) / n;

if (variance == 0)
return 0;

return autocorr[lag] / variance;
}

public IEnumerable GetCResults(int maxLag, double c0)
{
for (int lag = 0; lag   (double)lag);

return new PairReturn(values, GetCResults(maxLag, c0));
}
}
}

Однако следующий модульный тест завершился неудачно:

Код: Выделить всё

using System;
using System.Collections.Generic;
using Microsoft.VisualStudio.TestTools.UnitTesting;
using MyNameSpace;

namespace MyNameSpace____unit_test
{
[TestClass]
public class TimeSeriesAnalysisFFT____unit_test
{

[TestInitialize]
public void Initialize()
{
}

[TestMethod]
public void GetAutoCorrelationTestMethod()
{
IList vect3List = new List()
{
new Vector3(1,1,1),
new Vector3(2,2,2),
new Vector3(3,3,3),
new Vector3(4,4,4),
new Vector3(5,5,5)
};

var timeSeriesAnalysis = new TimeSeriesAnalysisFFT(vect3List);

double c0 = timeSeriesAnalysis.AutoCorr(0);//33
double c1 = timeSeriesAnalysis.AutoCorr(1);//30
double c2 = timeSeriesAnalysis.AutoCorr(2);//26
double c3 = timeSeriesAnalysis.AutoCorr(3);//21
double c4 = timeSeriesAnalysis.AutoCorr(4);//15

Assert.AreEqual(33, c0);
Assert.AreEqual(30, c1);
Assert.AreEqual(26, c2);
Assert.AreEqual(21, c3);
Assert.AreEqual(15, c4);
}

[TestMethod]
public void GetAutoCorrelationPointsTestMethod()
{
List vect3List = new List()
{
new Vector3(1,1,1),
new Vector3(2,2,2),
new Vector3(3,3,3),
new Vector3(4,4,4),
new Vector3(5,5,5)
};

var timeSeriesAnalysis = new TimeSeriesAnalysisFFT(vect3List);

var autocorr = timeSeriesAnalysis.GetAutoCorrelationPoints(10);

IList tausList = new List(autocorr.XData);
IList autocorrList = new List(autocorr.YData);

Assert.AreEqual(5, tausList.Count);
Assert.AreEqual(5, autocorrList.Count);

Assert.AreEqual(1.0, autocorrList[0]);
Assert.AreEqual(0.90909090909090906, autocorrList[1]);
Assert.AreEqual(0.78787878787878785, autocorrList[2]);
Assert.AreEqual(0.63636363636363635, autocorrList[3]);
Assert.AreEqual(0.45454545454545453, autocorrList[4]);
}
}
}

Что я делаю не так?

Подробнее здесь: https://stackoverflow.com/questions/786 ... -using-fft