Losowe zmienne Gaussa

Czy w standardowej bibliotece platformy .NET znajduje się klasa, która daje mi możliwość tworzenia zmiennych losowych zgodnych z rozkładem Gaussa?

Sugestia Jarretta dotycząca użycia transformacji Boxa-Mullera jest dobra dla szybkiego i brudnego rozwiązania. Prosta implementacja:

Random rand = new Random(); //reuse this if you are generating many
double u1 = 1.0-rand.NextDouble(); //uniform(0,1] random doubles
double u2 = 1.0-rand.NextDouble();
double randStdNormal = Math.Sqrt(-2.0 * Math.Log(u1)) *
             Math.Sin(2.0 * Math.PI * u2); //random normal(0,1)
double randNormal =
             mean + stdDev * randStdNormal; //random normal(mean,stdDev^2)

Wydaje się, że to pytanie przeniosło się do Google w przypadku generacji Gaussa .NET, więc pomyślałem, że opublikuję odpowiedź.

Zrobiłem kilka metod rozszerzających dla klasy .NET Random , w tym implementację transformacji Boxa-Mullera. Ponieważ są to rozszerzenia, o ile projekt jest uwzględniony (lub odwołujesz się do skompilowanej biblioteki DLL), nadal możesz to zrobić

var r = new Random();
var x = r.NextGaussian();

Mam nadzieję, że nikomu nie przeszkadza bezwstydna wtyczka.

Przykładowy histogram wyników (dołączona jest aplikacja demonstracyjna do rysowania):

Math.NET zapewnia tę funkcjonalność. Oto jak:

double mean = 100;
double stdDev = 10;

MathNet.Numerics.Distributions.Normal normalDist = new Normal(mean, stdDev);
double randomGaussianValue=   normalDist.Sample();

Dokumentację można znaleźć tutaj: http://numerics.mathdotnet.com/api/MathNet.Numerics.Distributions/Normal.htm

Utworzyłem prośbę o taką funkcję na Microsoft Connect. Jeśli jest to coś, czego szukasz, zagłosuj na to i zwiększ jego widoczność.

https://connect.microsoft.com/VisualStudio/feedback/details/634346/guassian-normal-distribution-random-numbers

Ta funkcja jest zawarta w Java SDK. Jego implementacja jest dostępna jako część dokumentacji i można ją łatwo przenieść na C # lub inne języki .NET.

Jeśli szukasz czystej prędkości, algorytm Zigorat jest ogólnie uznawany za najszybsze podejście.

Nie jestem jednak ekspertem w tym temacie - natknąłem się na taką potrzebę podczas wdrażania filtra cząstek stałych do mojej biblioteki symulowanej robotycznej piłki nożnej RoboCup 3D i byłem zaskoczony, że nie zostało to uwzględnione we frameworku.

W międzyczasie, oto opakowanie, Randomktóre zapewnia wydajną implementację metody biegunowej Box Muller:

public sealed class GaussianRandom
{
    private bool _hasDeviate;
    private double _storedDeviate;
    private readonly Random _random;

    public GaussianRandom(Random random = null)
    {
        _random = random ?? new Random();
    }

    /// <summary>
    /// Obtains normally (Gaussian) distributed random numbers, using the Box-Muller
    /// transformation.  This transformation takes two uniformly distributed deviates
    /// within the unit circle, and transforms them into two independently
    /// distributed normal deviates.
    /// </summary>
    /// <param name="mu">The mean of the distribution.  Default is zero.</param>
    /// <param name="sigma">The standard deviation of the distribution.  Default is one.</param>
    /// <returns></returns>
    public double NextGaussian(double mu = 0, double sigma = 1)
    {
        if (sigma <= 0)
            throw new ArgumentOutOfRangeException("sigma", "Must be greater than zero.");

        if (_hasDeviate)
        {
            _hasDeviate = false;
            return _storedDeviate*sigma + mu;
        }

        double v1, v2, rSquared;
        do
        {
            // two random values between -1.0 and 1.0
            v1 = 2*_random.NextDouble() - 1;
            v2 = 2*_random.NextDouble() - 1;
            rSquared = v1*v1 + v2*v2;
            // ensure within the unit circle
        } while (rSquared >= 1 || rSquared == 0);

        // calculate polar tranformation for each deviate
        var polar = Math.Sqrt(-2*Math.Log(rSquared)/rSquared);
        // store first deviate
        _storedDeviate = v2*polar;
        _hasDeviate = true;
        // return second deviate
        return v1*polar*sigma + mu;
    }
}

Math.NET Iridium twierdzi również, że implementuje „niejednorodne generatory losowe (normalne, poissona, dwumianowe, ...)”.

Oto kolejne szybkie i brudne rozwiązanie do generowania zmiennych losowych o rozkładzie normalnym . Rysuje losowy punkt (x, y) i sprawdza, czy ten punkt leży pod krzywą twojej funkcji gęstości prawdopodobieństwa, w przeciwnym razie powtórz.

Bonus: Możesz wygenerować zmienne losowe dla dowolnego innego rozkładu (np. Rozkładu wykładniczego lub rozkładu poissona ), zastępując funkcję gęstości.

    static Random _rand = new Random();

    public static double Draw()
    {
        while (true)
        {
            // Get random values from interval [0,1]
            var x = _rand.NextDouble(); 
            var y = _rand.NextDouble(); 

            // Is the point (x,y) under the curve of the density function?
            if (y < f(x)) 
                return x;
        }
    }

    // Normal (or gauss) distribution function
    public static double f(double x, double μ = 0.5, double σ = 0.5)
    {
        return 1d / Math.Sqrt(2 * σ * σ * Math.PI) * Math.Exp(-((x - μ) * (x - μ)) / (2 * σ * σ));
    }

Ważne: Wybierz przedział y oraz parametry σ i μ tak, aby krzywa funkcji nie była odcięta w jej punktach maksimum / minimum (np. Przy x = średnia). Pomyśl o przedziałach x i y jak o obwiedni, w której krzywa musi pasować.

Chciałbym rozwinąć odpowiedź @ yoyoyoyosef, czyniąc ją jeszcze szybszą i pisząc klasę opakowującą. Poniesione narzuty mogą nie oznaczać dwa razy szybciej, ale myślę, że powinno być prawie dwa razy szybsze. Nie jest to jednak bezpieczne dla wątków.

public class Gaussian
{
     private bool _available;
     private double _nextGauss;
     private Random _rng;

     public Gaussian()
     {
         _rng = new Random();
     }

     public double RandomGauss()
     {
        if (_available)
        {
            _available = false;
            return _nextGauss;
        }

        double u1 = _rng.NextDouble();
        double u2 = _rng.NextDouble();
        double temp1 = Math.Sqrt(-2.0*Math.Log(u1));
        double temp2 = 2.0*Math.PI*u2;

        _nextGauss = temp1 * Math.Sin(temp2);
        _available = true;
        return temp1*Math.Cos(temp2);
     }

    public double RandomGauss(double mu, double sigma)
    {
        return mu + sigma*RandomGauss();
    }

    public double RandomGauss(double sigma)
    {
        return sigma*RandomGauss();
    }
}

Rozwijając odpowiedzi @Noakes i @ Hameer, zaimplementowałem również klasę 'Gaussian', ale aby uprościć przestrzeń pamięci, uczyniłem ją dzieckiem klasy Random, aby można było również wywołać podstawowe Next (), NextDouble () , etc z klasy Gaussa bez konieczności tworzenia dodatkowego obiektu Random do obsługi tego. Wyeliminowałem również właściwości klas globalnych _available i _nextgauss, ponieważ nie uważałem ich za konieczne, ponieważ ta klasa jest oparta na instancji, powinna być bezpieczna dla wątków, jeśli nadasz każdemu wątkowi własny obiekt Gaussa. Usunąłem również wszystkie zmienne przydzielone w czasie wykonywania z funkcji i uczyniłem je właściwościami klas, co zmniejszy liczbę wywołań menedżera pamięci, ponieważ te 4 dublety teoretycznie nigdy nie powinny zostać cofnięte, dopóki obiekt nie zostanie zniszczony.

public class Gaussian : Random
{

    private double u1;
    private double u2;
    private double temp1;
    private double temp2;

    public Gaussian(int seed):base(seed)
    {
    }

    public Gaussian() : base()
    {
    }

    /// <summary>
    /// Obtains normally (Gaussian) distrubuted random numbers, using the Box-Muller
    /// transformation.  This transformation takes two uniformly distributed deviates
    /// within the unit circle, and transforms them into two independently distributed normal deviates.
    /// </summary>
    /// <param name="mu">The mean of the distribution.  Default is zero</param>
    /// <param name="sigma">The standard deviation of the distribution.  Default is one.</param>
    /// <returns></returns>

    public double RandomGauss(double mu = 0, double sigma = 1)
    {
        if (sigma <= 0)
            throw new ArgumentOutOfRangeException("sigma", "Must be greater than zero.");

        u1 = base.NextDouble();
        u2 = base.NextDouble();
        temp1 = Math.Sqrt(-2 * Math.Log(u1));
        temp2 = 2 * Math.PI * u2;

        return mu + sigma*(temp1 * Math.Cos(temp2));
    }
}

Rozwijając odpowiedź Drew Noakesa, jeśli potrzebujesz lepszej wydajności niż Box-Muller (około 50-75% szybciej), Colin Green udostępnił implementację algorytmu Ziggurat w C #, którą możesz znaleźć tutaj:

http://heliosphan.org/zigguratalgorithm/zigguratalgorithm.html

Ziggurat używa tabeli przeglądowej do obsługi wartości, które znajdują się wystarczająco daleko od krzywej, którą szybko akceptuje lub odrzuca. W około 2,5% przypadków musi wykonać dalsze obliczenia, aby określić, po której stronie krzywej znajduje się liczba.

Możesz wypróbować Infer.NET. Jednak nie jest to jeszcze licencja komercyjna. Tutaj jest link

Jest to probabilistyczny framework dla .NET, który opracowałem w moim badaniu Microsoft. Mają typy .NET dla dystrybucji Bernoulliego, Beta, Gamma, Gaussa, Poissona i prawdopodobnie kilka innych, które pominąłem.

Może osiągnąć to, czego chcesz. Dzięki.

To moja prosta realizacja inspirowana Box Mullerem. Możesz zwiększyć rozdzielczość, aby dopasować ją do swoich potrzeb. Chociaż dla mnie to działa świetnie, jest to przybliżenie o ograniczonym zakresie, więc pamiętaj, że ogony są zamknięte i skończone, ale z pewnością możesz je rozszerzyć w razie potrzeby.

    //
    // by Dan
    // islandTraderFX
    // copyright 2015
    // Siesta Key, FL
    //    
// 0.0  3231 ********************************
// 0.1  1981 *******************
// 0.2  1411 **************
// 0.3  1048 **********
// 0.4  810 ********
// 0.5  573 *****
// 0.6  464 ****
// 0.7  262 **
// 0.8  161 *
// 0.9  59 
//Total: 10000

double g()
{
   double res = 1000000;
   return random.Next(0, (int)(res * random.NextDouble()) + 1) / res;
}

public static class RandomProvider
{
   public static int seed = Environment.TickCount;

   private static ThreadLocal<Random> randomWrapper = new ThreadLocal<Random>(() =>
       new Random(Interlocked.Increment(ref seed))
   );

   public static Random GetThreadRandom()
   {
       return randomWrapper.Value;
   }
} 

Myślę, że nie ma. I naprawdę mam nadzieję, że tak nie jest, ponieważ framework jest już wystarczająco rozdęty, bez tak wyspecjalizowanej funkcjonalności wypełniającej go jeszcze bardziej.

Spójrz na http://www.extremeoptimization.com/Statistics/UsersGuide/ContinuousDistributions/NormalDistribution.aspx i http://www.vbforums.com/showthread.php?t=488959 dla rozwiązań .NET innej firmy.