tsai
/
mathnet-numerics
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity
Browse Source

Implementation of Hybrid Monte Carlo method.

pull/37/head
manyue 14 years ago
parent
7447ddbab7
commit
483cd7638f
2 changed files with 0 additions and 232 deletions
Whitespace
Show all changes Ignore whitespace when comparing lines Ignore changes in amount of whitespace Ignore changes in whitespace at EOL
Split View
Diff Options
Show Stats Download Patch File Download Diff File
  1. 91
      src/Numerics/Statistics/MCMC/MCMCDiagonistics.cs
  2. 141
      src/UnitTests/StatisticsTests/MCMCTests/MCMCDiagonisticsTest.cs

91
src/Numerics/Statistics/MCMC/MCMCDiagonistics.cs
View File

@ -1,91 +0,0 @@
// <copyright file="MCMCDiagonistics.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
// Copyright (c) 2009-2010 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Numerics;
namespace MathNet.Numerics.Statistics.Mcmc.Diagonistics
{
/// <summary>
/// Provides utilities to analysis the convergence of a set of samples from
/// a <seealso cref="McmcSampler{T}"/>.
/// </summary>
static public class MCMCDiagonistics
{
/// <summary>
/// Computes the auto correlations of a series evaluated by a function f.
/// </summary>
/// <param name="Series">The series for computing the auto correlation.</param>
/// <param name="lag">The lag in the series</param>
/// <param name="f">The function used to evaluate the series.</param>
/// <returns>The auto correlation.</returns>
/// <exception cref="ArgumentOutOfRangeException">Throws if lag is zero or if lag is
/// greater than or equal to the length of Series.</exception>
static public double ACF<T>(IEnumerable<T> Series, int lag, Func<T,double> f)
{
if (lag < 0)
throw new ArgumentOutOfRangeException("Lag must be positive");
int Length = Series.Count();
if (lag >= Length)
throw new ArgumentOutOfRangeException("Lag must be smaller than the sample size");
var TransformedSeries = from data in Series
select f(data);
var FirstSeries = TransformedSeries.Take(Length-lag);
var SecondSeries = TransformedSeries.Skip(lag);
return Correlation.Pearson(FirstSeries, SecondSeries);
}
/// <summary>
/// Computes the effective size of the sample when evaluated by a function f.
/// </summary>
/// <param name="Series">The samples.</param>
/// <param name="f">The function use for evaluating the series.</param>
/// <returns>The effective size when auto correlation is taken into account.</returns>
static public double EffectiveSize<T>(IEnumerable<T> Series, Func<T,double> f)
{
int Length = Series.Count();
double rho = ACF(Series, 1, f);
return ((1 - rho) / (1 + rho)) * Length;
}
}
}

141
src/UnitTests/StatisticsTests/MCMCTests/MCMCDiagonisticsTest.cs
View File

@ -1,141 +0,0 @@
// <copyright file="MCMCDiagonistics.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
// Copyright (c) 2009-2010 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using NUnit.Framework;
using MathNet.Numerics.Statistics;
using MathNet.Numerics.Distributions;
using MathNet.Numerics.Statistics.Mcmc.Diagonistics;
using MathNet.Numerics.Random;
namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
{
/// <summary>
/// MCMCDiagonistics testing.
/// </summary>
[TestFixture]
public class MCMCDiagonisticsTest
{
/// <summary>
/// For generation of a random series to test the methods.
/// </summary>
private System.Random rnd = new System.Random();
/// <summary>
/// Distribution to sample the entries of the random series from.
/// </summary>
private Normal dis = new Normal(0, 1);
/// <summary>
/// Testing the ACF function using a randomly generated series with a range
/// of lags.
/// </summary>
/// <param name="startlag">Minimum value of lag in the test.</param>
/// <param name="endlag">Maximum value of lag in the test.</param>
[TestCase(0, 10)]
[TestCase(11, 20)]
[TestCase(21, 30)]
[TestCase(31, 40)]
public void TestACF(int startlag, int endlag)
{
for (int lag = startlag; lag < endlag; lag++)
{
int Length = 10000;
double[] firstSeries = new double[Length - lag];
double[] secondSeries = new double[Length - lag];
double[] Series = new double[Length];
for (int i = 0; i < Length; i++)
{ Series[i] = RandomSeries(); }
double[] TransformedSeries = new double[Length];
for (int i = 0; i < Length; i++)
{ TransformedSeries[i] = Series[i] * Series[i]; }
Array.Copy(TransformedSeries, firstSeries, Length - lag);
Array.Copy(TransformedSeries, lag, secondSeries, 0, Length - lag);
double result = MCMCDiagonistics.ACF(Series, lag, x=>x*x);
double correlation = Correlation.Pearson(firstSeries, secondSeries);
Assert.AreEqual(result, correlation, 10e-13);
}
}
/// <summary>
/// Set lag to be greater than the length of the series throws a
/// <c>ArgumentOutOfRangeException</c>.
/// </summary>
[Test]
public void LagOutOfRange()
{
int Length = 10;
double[] Series = new double[Length];
Assert.Throws<ArgumentOutOfRangeException>(() => MCMCDiagonistics.ACF(Series, 11, x=>x));
}
/// <summary>
/// Set lag to be negative throws a <c>ArgumentOutOfRangeException</c>.
/// </summary>
[Test]
public void LagNegative()
{
Assert.Throws<ArgumentOutOfRangeException>(() => MCMCDiagonistics.ACF(new double[10], -1, x=>x));
}
/// <summary>
/// Generating a random number used for the entry of the series.
/// </summary>
/// <returns>A random number.</returns>
private double RandomSeries()
{ return rnd.NextDouble() + rnd.NextDouble() * (dis.Sample()); }
/// <summary>
/// Testing the effective size using a random series.
/// </summary>
[Test]
public void EffectiveSizeTest()
{
int Length = 10;
double[] Series = new double[Length];
for (int i = 0; i < Length; i++)
{ Series[i] = RandomSeries(); }
double rho = MCMCDiagonistics.ACF(Series, 1,x=>x*x);
double ESS = (1 - rho) / (1 + rho) * Length;
Assert.AreEqual(ESS, MCMCDiagonistics.EffectiveSize(Series,x=>x*x), 10e-13);
}
}
}
Write Preview
Loading…
Cancel
Save