65.9K
CodeProject 正在变化。 阅读更多。
Home

C#.NET 中的时间序列分析:第二部分

starIconstarIconstarIconstarIcon
emptyStarIcon
starIcon

4.60/5 (4投票s)

2012年5月21日

CPOL

3分钟阅读

viewsIcon

29400

在 R 和 C#.NET 中使用开源库进行格兰杰因果关系分析。

引言

在本系列的第 2 部分中,我们将研究如何使用 C#.NET 和 R 统计语言执行线性格兰杰因果关系检验。 通过同时使用两者,我们可以使用 Pfaff (2008) 的工作在 R 中原型化统计算法,然后转换为 C# 和/或 C++ 代码进行优化。 首先,我们希望在 C# 中拥有多个统计视角,因此我们可以处理与 R 的通信,并建立在第一部分的基础上。 其次,我们希望能够将我们的工作代码从开源社区迁移到高性能 CUDA 库。 这对于神经科学和商品价格研究中所需的大数据和实时预测是必要的。 这是本系列的主题之一。

背景

再次,我参考了许多关于多变量时间序列分析的好书,这些书都以格兰杰因果关系检验为特色。 当然,我的书——不是那本……等等……由 Springer Science 出版,今年将发布 “R 和 Model Maker 3 中的神经信息学和计算神经科学研究” 将更详细地介绍这项工作——无耻的宣传。 无论如何,您可以阅读 Library 1 的 PDF 文件

Library 1:CRAN 上提供的 R 包

1. tSeries - 诊断测试和 ARMA 建模

2. vars - 因果关系测试和多变量时间序列建模(良好的指导)

Library 2:R 到 .NET 的接口

R.NET 是一个到 C#.NET 的良好接口,可在 codeplex 上获得:http://rdotnet.codeplex.com/。 这比 StatConnector 好得多,StatConnector 是之前一篇文章的主题,并且占据了我一两年的人生。

使用代码

在图 1 的列表中,您必须用您自己的数据集和变量替换几行代码。 R 中的主要方法是causality(),它计算 VAR(p) 模型的格兰杰 (1969) 和瞬时因果关系的 F 和 Wald 检验统计量。 REngine 的酷之处在于您可以在方法中进行字符串处理,并且可以根据决策变量构建字符串。 这里有很多 AI 的可能性,以及我的神经对角线思维的东西(http://neuraldiagonalthinking.wordpress.com/) 与交叉频率耦合相结合,用于使用 HTML5 和 ASP.NET 预测脑电波状态——太令人兴奋了。 再次,这是本系列中的另一篇文章。

图1。

使用 Brainwaver R 包中的大脑数据在 R 中进行格兰杰因果关系检验

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using ABMath.ModelFramework.Models;
using ABMath.ModelFramework.Data;
using MathNet.Numerics.Distributions;
using MathNet.Numerics.RandomSources;
using MathNet.Numerics;
using MathNet.Numerics.LinearAlgebra;
using RDotNet;

namespace TCW_ModelMaker3
{
    public class R_GrangerCausality
    {

        public R_GrangerCausality()
        {
        }

        public void runR_CausalityTest()
        {
            REngine engine;
            REngine.SetDllDirectory(@"C:\Program Files\R\R-2.13.0\bin\i386");
            engine = REngine.CreateInstance("RDotNet", new[] { "-q" });
            engine.EagerEvaluate("load('C:/2011/The Cromwell Workshop/Research/Data/fmri/brain.RData')");
            engine.EagerEvaluate("library(tseries)");
            engine.EagerEvaluate("install.packages('vars',,'http://cran.ms.unimelb.edu.au')");
            engine.EagerEvaluate("library(vars)");
            engine.EagerEvaluate("x<-brain$V2");
            engine.EagerEvaluate("y<-brain$V3");
            engine.EagerEvaluate("X1<-data.frame(x,y)");
            engine.EagerEvaluate("var_Model<-VAR(X1,p=1)");
            engine.EagerEvaluate("causTest1<-causality(var_Model, cause='x')");
            engine.EagerEvaluate("causStat<-causTest1$Granger$statistic[1,1]");
            engine.EagerEvaluate("causValue<-causTest1$Granger$p.value[1,1]");
 
            string statistic = engine.GetSymbol("causStat").AsNumeric().First().ToString();
            string pvalue = engine.GetSymbol("causValue").AsNumeric().First().ToString();
            string hypothesis = "HO: X Does Not Granger Cause Y";
            string decision;
            double signif = 0.05;
            if (Convert.ToDouble(pvalue) < signif)
            {
                decision = "Reject";
 
            }
            else
            {
 
                decision = "Accept";
            }
            Console.WriteLine("Hypothesis:" + hypothesis);
            Console.WriteLine("Statistic:" + statistic + " and " + 
               "P-Value:" + pvalue + " and  " + decision);
 
            Console.ReadKey();
        }
	}
}

在 C#.NET 中,我不会像上面那样重复。 这是你的任务。 这里的过程是将 VAR(p) 模型拟合到数据,然后估计变量 X 和 Y 的 AR(p) 模型,并获得方差并计算检验统计量。 其他方法用于数据生成。

图2。

C# 中的格兰杰因果关系类。

public class GrangerCausality
{

    public double[] sigma2 { get; set; }
    public Matrix s2 { set; get; }
    public double[,] resids { set; get; }
    public double covariance { set; get; }
    public double[,] _data { set; get; }
    StatTest st;
 
    public GrangerCausality(int order, int dim, int obs)
    {
        List<TimeSeries> zt = generateTimeSeries(obs);
        VARModel model1 = new VARModel(order, dim);
        MVTimeSeries m1 = new MVTimeSeries(zt, false);

        resids = model1.resid_Data;
        model1.SetInput(0, m1, null);
        model1.FitByMethodOfMoments();
        model1.ComputeResidualsAndOutputs();
        s2 = model1.Sigma;
        sigma2 = model1.sigMa2;
        covariance = model1.coVariance;
        Console.WriteLine(model1.Description);
        //
        //Fit AR Model for X
        //
        ARMAModel model2 = new ARMAModel(order, 0);
        TimeSeries ts1 = new TimeSeries();
        ts1 = zt[0];
        model2.SetInput(0, ts1, null);
        model2.FitByMLE(10, 10, 10, null);
        model2.ComputeResidualsAndOutputs();
        Console.WriteLine(model2.Description);
        var sigmaVAR_X = sigma2[0];
        var sigma_AR_X = model2.Sigma;
        //
        //Fit AR Model for Y
        //
        ARMAModel model3 = new ARMAModel(order, 0);
        TimeSeries ts2 = new TimeSeries();
        ts2 = zt[1];
        model3.SetInput(0, ts2, null);
        model3.FitByMLE(10, 10, 10, null);
        model3.ComputeResidualsAndOutputs();
        Console.WriteLine(model3.Description);
        var sigmaVAR_Y = sigma2[1];
        var sigma_AR_Y = model3.Sigma; 
        //
        //  Compute Statistics
        //
        st = new StatTest();
        var GC_1 = obs*Math.Log(sigma_AR_X / sigmaVAR_X);
        var GC_2 = obs*Math.Log(sigma_AR_Y / sigmaVAR_Y);
        var GC_3 = obs*Math.Log((sigmaVAR_X * sigmaVAR_Y) / ((sigmaVAR_X * sigmaVAR_Y) - model1.coVariance));
        var GC_TL = GC_1 + GC_2 + GC_3;
        st._FXToY = GC_1;
        st._FYToX = GC_2;
        st._FYX = GC_3;
        st._FLinear = GC_TL;
        Console.WriteLine("GC_XToY:"+st.computeHypothTest(0.10,st._FXToY,1));
        Console.WriteLine("GC_YToX:"+st.computeHypothTest(0.10, st._FYToX, 1));
        Console.WriteLine("GC_YX:"+st.computeHypothTest(0.10, st._FYX, 1));
        Console.WriteLine("GC_FLinear:"+st.computeHypothTest(0.10, st._FLinear, 3));
        Console.ReadKey();
    }
    public List<TimeSeries> generateTimeSeries(int obs)
    {
        //Generate Data
        TimeSeries ts1 = new TimeSeries();
        TimeSeries ts2 = new TimeSeries();
        List<TimeSeries> _lstTime = new List<TimeSeries>();

        var current = new DateTime(2001, 1, 1);
        _data = generateData(obs, 2);
        for (int t = 0; t < obs; ++t)
        {

            ts1.Add(current, _data[t, 0], false);
            ts2.Add(current, _data[t, 1], false);
            current = new DateTime(current.AddDays(1).Ticks);
        }
        _lstTime.Add(ts1);
        _lstTime.Add(ts2);
        return _lstTime;

    }
    public double[,] generateData(int obs, int dim)
    {
        var normalSim = new StandardDistribution();

        _data = new double[obs, dim];

        for (int t = 0; t < obs; ++t)
        {
            double s1 = normalSim.NextDouble();
            double s2 = normalSim.NextDouble();
            _data[t, dim - 2] = s1;
            _data[t, dim - 1] = s2;
        }
        return _data;
    }
}

然后,我们可以使用我们在第一部分中的开源数学库,并将其包装到一个类中以获得检验统计量。 我仍然需要写一篇关于这个的文章……使用基于代理的模型设计指定、估计和预测线性和非线性时间序列模型的方法。 哦,好吧,我正在接近它——宁愿阅读它也不愿做它。 我不会在这里抱怨太多。

图3。

StatTest

public class StatTest
{
    public double _FTest { set; get; }
    public double _pvalue { set; get; }
    public double _FXToY { set; get; }
    public double _FYToX { set; get; }
    public double _FYX { set; get; }
    public double _FLinear { set; get; }
    public string decision { set; get; }

    public StatTest()
    {
    }
    public string computeHypothTest(double pvalue, double _Fstat, int dgFreedom)
    {
        if (_Fstat > 0)
        {
            ChiSquareDistribution chi = new ChiSquareDistribution();
            chi.SetDistributionParameters(dgFreedom);

            double x = chi.ProbabilityDensity(_Fstat);
            double pval = 1 - chi.CumulativeDistribution(_Fstat);
            if (pval > pvalue)
            {
                decision = "Accept H0: No Granger Causality.";
            }
            else
            {
                decision = "Reject H0: There is Granger Causality.";
            }
        }
        else
        {

            pvalue = 0;
            decision = "Accept H0: No Granger Causality.";
        }
        return decision;
    }
}

当然,这可以扩展到多个多变量关系,并将上述内容放入人工智能框架中,如前所述。 换句话说,对于给定的实时任务,有多少个布罗德曼区域是连接的,并且在什么频率上连接? 代码位于 http://modelmaker3.codeplex.com/ 可供下载。

参考文献

Granger, C.W.J. (1969), "Investigating causal relations be econometric models and cross-spectral methods. Econometrica 37:424-438.

Pfaff, B. (2008), "VAR, SVAR and SVEC Models: Implementation Within R Package vars." Journal of Statistical Software 27(4) at http://www.jstatsoft.org/v27/i04/.

历史

  • 2012 年 5 月 21 日 代码已添加到 CodePlex。
© . All rights reserved.