Skip to content

Instantly share code, notes, and snippets.

@zmughal
Last active June 9, 2021 04:05
Show Gist options
  • Save zmughal/fd79961a166d653a7316aef2f010a1a2 to your computer and use it in GitHub Desktop.
Save zmughal/fd79961a166d653a7316aef2f010a1a2 to your computer and use it in GitHub Desktop.

Running with SciPy:

#!/usr/bin/python3
import scipy.integrate as integ
import numpy as np
def normalProbabilityDensity(x):
    constant = 1.0 / np.sqrt(2*np.pi)
    return(constant * np.exp((-x**2) / 2.0) )
neg_percentile, _ = integ.quad(normalProbabilityDensity, np.NINF, -1.00)
pos_percentile, _ = integ.quad(normalProbabilityDensity, np.NINF, 1.00)
print('Neg: ', neg_percentile)
print('Pos: ', pos_percentile)
# Neg:  0.15865525393145707
# Pos:  0.8413447460685435

Running in R:

#!/usr/bin/env Rscript -
pnorm(c(-1,1))
# [1] 0.1586553 0.8413447

Running with PDL:

#!/usr/bin/env perl

use PDL;
use PDL::Constants qw(PI);
use PDL::GSL::INTEG;
use PDL::GSL::CDF;
use feature qw(state say);

sub norm_pdf {
  my ($x) = @_;
  state $constant = 1 / sqrt(2*PI);
  return $constant * exp(- $x**2 / 2);
}

my $xvals = [-1, 1];
my ($res) = gslinteg_qagil( \&norm_pdf, $xvals , 1e-15, 0, 1000);

my $direct = gsl_cdf_gaussian_P( $xvals, 1);

say $res->string('%.18f');
say $direct->string('%.18f');
say all(abs( $res - $direct ) < 1e-15);

# [0.158655253931457074 0.841344746068543148]
# [0.158655253931457046 0.841344746068542926]
# 1
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment