Factor of 2?

[JonHowMIT]

Hi - is is possible that this equation should be modified as noted below?

// return mean + variance * cos( 2 * M_PI * R1) * sqrt(-log(R2)); return mean + variance * cos( 2 * M_PI * R1) * sqrt(-2*log(R2));

This appears to be the Box-Muller Transformation, and thus I think there needs to be a 2 inside the sqrt. thx

[drf5n]

Yes. The Box-Muller transformation does mean that:

https://github.com/ivanseidel/Gaussian/blob/e891a2a01f4bf2f15a4e7196b08660135697f58e/Gaussian.h#L130

should read

         return mean + variance * cos( 2 * M_PI * R1) * sqrt(-2*log(R2));

This example of summing samples of Gaussian(0,1) produces approximately mean(0), var=0.5:

// https://wokwi.com/projects/360282019322936321
// for https://github.com/ivanseidel/Gaussian/issues/6

#include "Gaussian.h"


const int sample_size = 10000;
const float var_set = 1.0;
const float mean_set = 0;

void setup() {
  Serial.begin(9600);
  Serial.println("\n\nStarting Random Gaussian Distribution...");
  Serial.print("Gaussian(");
  Serial.print(mean_set);
  Serial.print(",");
  Serial.print(var_set);
  Serial.println(")");
  delay(20);
  randomSeed(20230326);
}


void loop() {
  Gaussian g1 = Gaussian(mean_set, var_set);
  float sum = 0, sumsq = 0, mean = 0, variance = 0;

  // Zero the array
  for (int i = 0; i < sample_size; i++) {
    double x = g1.random();
    sum += x;
    sumsq += x*x;
  }
  mean = sum/sample_size;
  variance = (sumsq-sample_size*mean*mean)/(sample_size-1);

  Serial.print("n:");
  Serial.print(sample_size);
  
  Serial.print(" sum:");
  Serial.print(sum);
  
  Serial.print(" ssq:");
  Serial.print(sumsq);
  
  Serial.print(" mean:");
  Serial.print(mean);
  
  Serial.print(" var:");
  Serial.print(variance);
  Serial.println();
}

produces this output:

Starting Random Gaussian Distribution...
Gaussian(0.00,1.00)
n:10000 sum:-21.59 ssq:5016.92 mean:-0.00 var:0.50
n:10000 sum:-67.49 ssq:4974.67 mean:-0.01 var:0.50
n:10000 sum:54.96 ssq:5090.27 mean:0.01 var:0.51
n:10000 sum:55.73 ssq:4858.84 mean:0.01 var:0.49
n:10000 sum:68.58 ssq:5038.47 mean:0.01 var:0.50
n:10000 sum:-31.55 ssq:4904.44 mean:-0.00 var:0.49
n:10000 sum:56.46 ssq:5057.28 mean:0.01 var:0.51
n:10000 sum:-19.36 ssq:4876.92 mean:-0.00 var:0.49
n:10000 sum:32.28 ssq:4994.59 mean:0.00 var:0.50

See this for a simulation https://wokwi.com/projects/360282019322936321

For Gaussian(0,10) the results are mean=0, var = 50:

Starting Random Gaussian Distribution...
Gaussian(0.00,10.00)
n:10000 sum:-215.91 ssq:501691.40 mean:-0.02 var:50.17
n:10000 sum:-674.85 ssq:497466.50 mean:-0.07 var:49.75
n:10000 sum:549.60 ssq:509027.06 mean:0.05 var:50.90
n:10000 sum:557.26 ssq:485883.06 mean:0.06 var:48.59
n:10000 sum:685.77 ssq:503847.90 mean:0.07 var:50.39
n:10000 sum:-315.49 ssq:490444.37 mean:-0.03 var:49.05
n:10000 sum:564.61 ssq:505727.43 mean:0.06 var:50.57

The input spread parameter appears to be sigma rather than variance.

Gaussian( mu, sigma) versus Gaussian(mu, sigma^2)

[drf5n]

Looking at it further, the line

https://github.com/ivanseidel/Gaussian/blob/master/Gaussian.h#L130

should be:

return mean + sqrt(variance) * cos( 2 * M_PI * R1) * sqrt(-2*log(R2));