Generate Random variable dataset

Been reading and trying some things out.
Might be useful for people.
Code and comments should be explanation.

pi    = 3.14159265;
TWOPI = 6.28318530718;

N     = 5000;//Param("Dataset Size", 5000, 20, 5000);
hist  = 50;//Param("No of bins", 50, 10, 200);

function mzRandGauss( N ) {
/* - generate gauss variables f(x|0,1)~N(0,1), from uniform mtRandom~[0,1]
   converted from Python 3.6 code (Box-Muller transform)
	z = self-gauss_next
	self.gauss_next = none
	if z is None:
		x2pi = random() * TWOPI
		g2rad = _sqrt(-2.0 * _log(1.0 - random()))
		z = _cos(x2pi) * g2rad
		self.gauss_next = _sin(x2pi) * g2rad
	return mu + z*sigma
*/
	x2pi  = g2rad = Matrix( N, 1 );
	x2pi  = MxSetBlock( x2pi, 0, n-1, 0, 0, mtRandomA()*TWOPI );
	g2rad = MxSetBlock( g2rad, 0, n-1, 0, 0, sqrt(-2.0*log(1.0-mtRandomA())) );
	return (g2rad * cos(x2pi));
}

// generate rand.gauss numbers
ym    = MxGetBlock( mzRandGauss(N), 0, N-1, 0, 0 );
ym    = MxSortrows( ym, True, 0 );                    // order the variables
zm    = MxSetBlock( ym, 0, N-1, 0, 0, 1 );            // initalise "Vol" or frequency value

nS    = MxGetSize( ym, 0 );

all
I’m sorry but the mtRandomA() has been used incorrectly.
I’ve found that mtRandomA has created barcount-1 array and applied 15 groups (5000/331) to the 5000-element matrix.

If there are solutions, I would appreciate.
In the meantime, back to the drawing board.

Hello,
re-worked code and lessons learnt,

  1. always explicitly declare matrix
    if you don’t, then mult matrix by scalar results in 1D vector
    eg. x2pi = matrix(N,1); // z not declared
    z = 2.0 * x2pi // z is 1D vector

    then,
    z[j][0] = … will have error since using 2 subscripts

then,
z[j] = … will have error, array subscript out of range, only use 0…barcount-1
This has been my experience.
corrected code as separate post for clarity.

corrected code

pi    = 3.14159265;
TWOPI = 6.28318530718;

N     = 5000;//Param("Dataset Size", 5000, 20, 5000);
hist  = 50;//Param("No of bins", 50, 10, 200);

function mzRandGauss( N ) {
/* - generate gauss variables f(x|0,1)~N(0,1), from uniform mtRandom~[0,1]
   converted from Python 3.6 code (Box-Muller transform)
	z = self-gauss_next
	self.gauss_next = none
	if z is None:
		x2pi = random() * TWOPI
		g2rad = _sqrt(-2.0 * _log(1.0 - random()))
		z = _cos(x2pi) * g2rad
		self.gauss_next = _sin(x2pi) * g2rad
	return mu + z*sigma
*/
	x2pi  = g2rad = z = Matrix( N, 1 );
	for(j=0; j<N; j++) { x2pi[j][0] = mtRandom() * TWOPI; }
	z = cos( x2pi);

	for(j=0; j<N; j++) { g2rad[j][0] = sqrt( -2.0*log( 1.0-mtRandom() ) ); }
	return g2rad * z;
}

// generate rand.gauss numbers
ym = zm = yt = Matrix( N, 1 );
ym    = MxGetBlock( mzRandGauss(N), 0, N-1, 0, 0 );
yt    = MxSortrows( ym, True, 0 );                    // order the variables

zm    = MxSetBlock( yt, 0, N-1, 0, 0, 1.0 );       // initalise "Vol" or frequency value

nS    = MxGetSize( yt, 0 );

1 Like

Generate mtRandom numbers in a matrix, which produces a matrix of Gaussian random numbers, then plot the distribution.
Hope it’s useful or put in helper code folder.

pi    = 3.14159265;
TWOPI = 6.28318530718;

// --- parameters ---
N     = 5000;//Param("Dataset Size", 5000, 20, 5000);
hist  = 50;//Param("No of bins", 50, 10, 200);

// --- functions ---
function mzRand( yx, M ) {
	for( j=0; j<M; j++ ) { yx[j][0]  = mtRandom(); }
}

function mzRandGauss( yx, M ) {
// - generate gauss variables f(x|0,1)~N(0,1), from uniform mtRandom~[0,1]
	rand  = Matrix( M, 1, 0 );
	mzRand( &rand, M );
	yx  = cos ( rand * TWOPI );
	
	mzRand( &rand, M );
	yx  = yx * sqrt( -2 * log( 1-rand ) );
	return 0;
}

xt    = Matrix( N, 1, 0 );                            // use xt as working matrix
ok    = mzRandGauss( &xt, N );
// end of generating distribution 1D-vector, unordered
/****************************************************/
// --- distribution ---

// combine x-variable with frequency data
yt    = yOrd = Matrix( N, 2, 0 );                     // xy-data
for( i=0; i<N; i++ ) { yt[i][0] = xt[i][0]; }
yt    = MxSetBlock( yt, 0, N-1, 1, 1, 1 );
// sort the matrix, keeping x & y-axis values together
yOrd  = MxSortRows( yt, True, 0 );                   // ordered data

// Min, Max values; x-axis
xMin  = yOrd[0][0];
xMax  = yOrd[N-1][0];
rangX = (xMax - xMin) / hist;

/****************************************************/
// arrange ordered x-y dataset (yOrd) into bins
// ie. the distribution
mx    = Matrix( hist+1, 2, 0 );
yMax  = 0;
for( i=0; i<N; i++ ) {
	bin = floor( ( yOrd[i][0] - xMin ) / rangX );
	if( mx[bin][0] ==0 ) {
		mx[bin][0] = xMin+i*rangX + rangX/2;          // mid-point value for each bin
	}
	mx[bin][1] = mx[bin][1] + yOrd[i][1];             // frequency for each bin
	yMax  = Max( yMax, mx[bin][1] );                  // running max frequency / yScale
}

// scale y-data ~ [0,1]
for( i=0; i<=bin; i++ ) {
	mx[i][1] = mx[i][1] / yMax;
}

/****************************************************/
// Mean
EX    = MxSum( xt ) / N;

// Standard deviation:- VAR(X) = E(X^2)- (EX)^2
xSD   = sqrt( ( MxSum( xt * xt ) - EX * EX ) / N );

// percentile or quartile(1/4) or quintile(1/5):- pick from ordered data
q25tile   = yOrd[ ceil(0.25*N) ][0];
p50cent   = yOrd[ ceil(0.5*N) ][0];
q75tile   = yOrd[ ceil(0.75*N) ][0];
p5cent    = yOrd[ ceil(0.05*N) ][0];
p95cent   = yOrd[ ceil(0.95*N) ][0];

/****************************************************/
// --- Plot distribution ---
bi    = BarIndex();
selB  = SelectedValue( bi );
Nbar  = Max( Min(selB, N), 1 );

fvb   = FirstVisibleValue( bi ); 
lvb   = LastVisibleValue( bi );

GfxSetCoordsMode( 3 );                                // mode=3: x-bar, y-pixel
pH     = Status("pxChartHeight");                     // returns chart pixel height (bottom-top)
scaleY = 0.95;                                        // space at top of chart

bins   = MxGetSize( mx, 0 );
width  = (lvb-fvb) / bins;                            // no. of bars per bin *visible bars*
for( i=0; i<bins; i++ ) { 
	price = mx[i][0];                                 // price level 
	relvolume = scaleY * mx[i][1];                    // relative 0..1

	X1 = fvb+i*width;             X2 = X1 + width;
	Y1 = ph * ( 1 - relvolume);   Y2 = ph;
	
	GfxGradientRect( X1, Y1, X2, Y2, ColorRGB(70,255,255), ColorRGB(70,20,255) );
}


/****************************************************/
// Legend
_N( Title =   StrFormat("Mean[%g]        = %2.5f", N, EX )
			+ StrFormat("\nSDev                   = %2.5f", xSD )
			+ StrFormat("\nMin                       = %2.5f", xMin )
			+ StrFormat("\nMax                      = %2.5f", xMax )
			+ StrFormat("\nquartile: 25%%     = %2.5f (%g)", q25tile, ceil(0.25*N) )
			+ StrFormat("\nquartile: 50%%     = %2.5f (%g)", p50cent, ceil(0.5*N) )
			+ StrFormat("\nquartile: 75%%     = %2.5f (%g)", q75tile, ceil(0.75*N) )
			+ StrFormat("\npercentile: 5%%   = %2.5f (%g)", p5cent, ceil(0.05*N) )
			+ StrFormat("\npercentile: 95%% = %2.5f (%g)", p95cent, ceil(0.95*N) ) );



3 Likes