Replicating numpy.random.randn in AFL

Extracted from numpy help https://docs.scipy.org/doc/numpy-1.15.1/reference/generated/numpy.random.randn.html

numpy.random.randn

numpy.random. randn ( d0 , d1 , ... , dn )

Return a sample (or samples) from the “standard normal” distribution.

My goal is to achieve something similar in AFL.
I would appreciate any math expert opinion , if AFL below is doing the job and if there is a simple/faster approach. Thanks in advance.

image

backgroundColor = colorWhite;
SetChartBkColor( backgroundColor );
GfxSetOverlayMode( 2 );

gl_Random_i = BarCount;
gl_Random = 0;

function NextRandom()
{
	global gl_Random_i;
	global gl_Random;
	
	gl_Random_i++;
	if(gl_Random_i >= BarCount) {
		gl_Random_i = 0;
		gl_Random = mtRandomA();
	}
	return gl_Random[gl_Random_i];
}

function MatrixRandom( rows , cols)
{
	global mxrand;
	global mxrandNorm;
	
	nextNextGaussian = 0;
	haveNextNextGaussian = false;
	
	mxrandNorm = Matrix( rows , cols );
	mxrand = Matrix( rows , cols );
	for( i = 0 ; i < rows ; i ++ ) {
		for( j = 0; j < cols ; j ++ ) {

			if (haveNextNextGaussian) {
				haveNextNextGaussian = false;
				r = nextNextGaussian;
				r1 = nextNextGaussianRand;
			} else {
				do {
					r1 = NextRandom();
					r2 = NextRandom();
					v1 = 2 * r1 - 1;   // between -1.0 and 1.0
					v2 = 2 * r2 - 1;   // between -1.0 and 1.0
					s = v1 * v1 + v2 * v2;
					r1 = v1; r2 = v2;
				} while (s >= 1 || s == 0);
				multiplier = sqrt(-2 * log(s)/s);
				nextNextGaussian = v2 * multiplier;
				nextNextGaussianRand = r2;
				haveNextNextGaussian = true;
				r = v1 * multiplier;
			}
			mxrandNorm[i][j] = r * 0.25;
			if( mxrandNorm[i][j] > 1 OR mxrandNorm[i][j] < -1 ) {
				mxrandNorm[i][j] = 2 * NextRandom() - 1;
			}
			mxrand[i][j] = r1;
		}
	}
}

function Draw(X, Y, dot_size, dd)
{
    rows = MxGetSize( X, 0 );
    cols = MxGetSize( X, 1 );
    origX = 300;
    origY = 40;

	GfxSelectFont( "Tahoma", 12 );
	GfxSetTextColor( colorBlue );
	GfxTextOut( "(0,0)", -0.1 * dd + origX, -0.12 * dd + origY);
	GfxTextOut( "(0,1)", -0.1 * dd + origX, 1 * dd + origY);

    GfxSelectSolidBrush( backgroundColor );
	GfxSelectPen( colorBlack, 1 ); 
	GfxMoveTo( -1 * dd + origX, 0 * dd + origY );
	GfxLineTo( 1 * dd + origX, 0  * dd + origY);    
	
	GfxMoveTo( 0 * dd + origX, 0 * dd + origY );
	GfxLineTo( 0 * dd + origX, 1  * dd + origY);    

	GfxSelectPen( colorRed, 1 ); 
	GfxSelectSolidBrush( colorRed );
    for( j = 0; j < cols; j++ )
    {
        GfxCircle( X[0][j] * dd + origX, abs(Y[0][j]) * dd + origY, dot_size );
    }

}

MatrixRandom( 1, 3000);

dot_size = 1;
dd = 200;
Draw(mxrandNorm, mxrand, dot_size, dd);

Sorry sir, what is this for? Does it function similar to rrg. Thanks

No , has nothing to do with rrg.

in the new 6.35 beta there is inverf() function that should aid you (I have not tried it myself)

@awilson

I'm not sure this will be relevant due to the way specific functions are implemented in AFL, but maybe it is worth to take a look at.

Did you try to use the Box-Muller scheme using the AFL trigonometric functions?

Please, see this sample article that compares the speed of execution (in Python):

1 Like

Hi @beppe , thanks for the link, it looks better ....

Is it possible to to fill a matrix with random numbers without using loop ?

image

function DrawRandomA(dot_size, dd)
{
	global	NextGaussian;
	global	Gaussian;
	global	v1;
	global	v2;

    origX = 500;
    origY = 40;

	GfxSelectFont( "Tahoma", 12 );
	GfxSetTextColor( colorBlue );
	GfxTextOut( "(0,0)", -0.05 * dd + origX, -0.05 * dd + origY);
	GfxTextOut( "(0,1)", -0.05 * dd + origX, 1 * dd + origY);

    GfxSelectSolidBrush( backgroundColor );
	GfxSelectPen( colorBlack, 1 ); 
	GfxMoveTo( -1 * dd + origX, 0 * dd + origY );
	GfxLineTo( 1 * dd + origX, 0  * dd + origY);    
	
	GfxMoveTo( 0 * dd + origX, 0 * dd + origY );
	GfxLineTo( 0 * dd + origX, 1  * dd + origY);    

	GfxSelectPen( colorRed, 1 ); 
	GfxSelectSolidBrush( colorRed );
	rows = MxGetSize(Gaussian, 0);
	cols = MxGetSize(Gaussian, 1);
    for( i = 0; i < rows; i++ ) {
		for( j = 0; j < cols; j++ )
		{
			GfxCircle( Gaussian[i][j] * dd + origX, abs(v1[i][j]) * dd + origY, dot_size );
			//GfxCircle( z1[j] * dd + origX, abs(u2[j]) * dd + origY, dot_size );
			GfxCircle( NextGaussian[i][j] * dd + origX, abs(v2[i][j]) * dd + origY, dot_size );
			//GfxCircle( z2[j] * dd + origX, abs(u2[j]) * dd + origY, dot_size );
		}
	}
}

function BoxMMatrixRandomA( rows , cols)
{
	global	NextGaussian;
	global	Gaussian;
	global	v1;
	global	v2;
	
	pi = 4 * atan( 1 ) ; //Pi
	
	r1 = Matrix( rows , cols );
	r2 = Matrix( rows , cols );
	for( i = 0 ; i < rows ; i ++ ) {
		for( j = 0; j < cols ; j ++ ) {
			do {
				rd = mtRandom();
			} while (rd == 0 OR rd == 1);
			r1[i][j] = rd;
			do {
				rd = mtRandom();
			} while (rd == 0 OR rd == 1);
			r2[i][j] = rd;
		}
	}
	v1 = r1;
	v2 = r1;
    r = sqrt(-0.1666 * log(r1));
    x = cos(2*pi*r2);
    y = sin(2*pi*r2);
	s = r * x;
	// lets make sure s only have correct numbers
	InRange = s > -1 AND s < 1;
	s = s * InRange; // an number not in range , s = zero
	EqualZero = s == 0;
	yy = Matrix(rows , cols, 0.001 );
	yy = yy * EqualZero;
	// S = 0 => s = 0.001
	Gaussian = s + yy;
	
    s = r * y;
	// lets make sure s only have correct numbers
	InRange = s > -1 AND s < 1;
	s = s * InRange; // an number not in range , s = zero
	EqualZero = s == 0;
	yy = yy * EqualZero;
	// S = 0 => s = 0.001
	NextGaussian = s + yy;
}

dot_size = 1;
dd = 300;
BoxMMatrixRandomA( 100 , 100);
DrawRandomA(dot_size, dd);

code correction...

function DrawRandomA(dot_size, dd)
{
	global	NextGaussian;
	global	Gaussian;
	global	v1;
	global	v2;

    origX = 500;
    origY = 40;

	GfxSelectFont( "Tahoma", 12 );
	GfxSetTextColor( colorBlue );
	GfxTextOut( "(0,0)", -0.05 * dd + origX, -0.05 * dd + origY);
	GfxTextOut( "(0,1)", -0.05 * dd + origX, 1 * dd + origY);

    GfxSelectSolidBrush( backgroundColor );
	GfxSelectPen( colorBlack, 1 ); 
	GfxMoveTo( -1 * dd + origX, 0 * dd + origY );
	GfxLineTo( 1 * dd + origX, 0  * dd + origY);    
	
	GfxMoveTo( 0 * dd + origX, 0 * dd + origY );
	GfxLineTo( 0 * dd + origX, 1  * dd + origY);    

	GfxSelectPen( colorRed, 1 ); 
	GfxSelectSolidBrush( colorRed );
	rows = MxGetSize(Gaussian, 0);
	cols = MxGetSize(Gaussian, 1);
    for( i = 0; i < rows; i++ ) {
		for( j = 0; j < cols; j++ )
		{
			GfxCircle( Gaussian[i][j] * dd + origX, abs(v1[i][j]) * dd + origY, dot_size );
			//GfxCircle( z1[j] * dd + origX, abs(u2[j]) * dd + origY, dot_size );
			GfxCircle( NextGaussian[i][j] * dd + origX, abs(v2[i][j]) * dd + origY, dot_size );
			//GfxCircle( z2[j] * dd + origX, abs(u2[j]) * dd + origY, dot_size );
		}
	}
}

function BoxMMatrixRandomA( rows , cols)
{
	global	NextGaussian;
	global	Gaussian;
	global	v1;
	global	v2;
	
	pi = 4 * atan( 1 ) ; //Pi
	
	r1 = Matrix( rows , cols );
	r2 = Matrix( rows , cols );
	for( i = 0 ; i < rows ; i ++ ) {
		for( j = 0; j < cols ; j ++ ) {
			do {
				rd = mtRandom();
			} while (rd == 0 OR rd == 1);
			r1[i][j] = rd;
			do {
				rd = mtRandom();
			} while (rd == 0 OR rd == 1);
			r2[i][j] = rd;
		}
	}
	v1 = r1;
	v2 = r1;
    r = sqrt(-0.1666 * log(r1));
    x = cos(2*pi*r2);
    y = sin(2*pi*r2);
	s = r * x;
	// lets make sure s only have correct numbers
	InRange = s > -1 AND s < 1;
	s = s * InRange; // an number not in range , s = zero
	EqualZero = s == 0;
	yy = Matrix(rows , cols, 0.001 );
	yy = yy * EqualZero;
	// S = 0 => s = 0.001
	Gaussian = s + yy;
	
    s = r * y;
	// lets make sure s only have correct numbers
	InRange = s > -1 AND s < 1;
	s = s * InRange; // an number not in range , s = zero
	EqualZero = s == 0;
	yy = Matrix(rows , cols, 0.001 );
	yy = yy * EqualZero;
	// S = 0 => s = 0.001
	NextGaussian = s + yy;
}

dot_size = 1;
dd = 300;
BoxMMatrixRandomA( 100 , 100);
DrawRandomA(dot_size, dd);


if you need array (not matrix) filled with normally distributed samples you can achieve it in one line of code:

_ArrStdDist = inverf(2 * mtRandomA() - 1.0);

this is how it looks on histogram:

image

1 Like

Thanks for the sugestion.
inverf works with Matrix, but when I plot it does not look as I expected

image

function inverfMatrixRandomA( rows , cols)
{
	global	inv;
	global	v1;
	
	r1 = Matrix( rows , cols );
	for( i = 0 ; i < rows ; i ++ ) {
		for( j = 0; j < cols ; j ++ ) {
			do {
				rd = mtRandom();
			} while (rd == 0 OR rd == 1);
			r1[i][j] = rd;
		}
	}
	inv = inverf(2 * r1 - 1.0);
	v1 = r1;
}