From 78508b911c7923c74d982f6e8a0145ebfbbd5ce0 Mon Sep 17 00:00:00 2001 From: Abhist17 Date: Mon, 28 Sep 2026 04:12:51 +0530 Subject: [PATCH 1/2] fix: avoid overflow in `stats/base/dists/poisson/entropy` for large `lambda` The series `sum_k lambda^k ln(k!)/k!` overflowed once `lambda^k` or `k!` reached `Infinity`, so the function returned `Infinity` for `lambda >= ~80` and `NaN` for `lambda >= ~745`. Evaluate `-sum p_k ln(p_k)` with weights relative to the mode instead, and use the asymptotic expansion for `lambda >= 5000`. Regenerate the fixtures in extended precision and add fixtures for medium and large `lambda`. --- .../base/dists/poisson/entropy/lib/main.js | 101 +++++++++++++----- .../base/dists/poisson/entropy/manifest.json | 22 ++-- .../base/dists/poisson/entropy/src/main.c | 90 ++++++++++++---- .../entropy/test/fixtures/julia/REQUIRE | 2 +- .../entropy/test/fixtures/julia/data.json | 2 +- .../entropy/test/fixtures/julia/large.json | 1 + .../entropy/test/fixtures/julia/medium.json | 1 + .../entropy/test/fixtures/julia/runner.jl | 44 +++++++- .../base/dists/poisson/entropy/test/test.js | 42 +++++++- .../dists/poisson/entropy/test/test.native.js | 42 +++++++- 10 files changed, 284 insertions(+), 63 deletions(-) create mode 100644 lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/large.json create mode 100644 lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/medium.json diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/lib/main.js b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/lib/main.js index 80006e45bf83..5cebcadbffa5 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/lib/main.js +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/lib/main.js @@ -20,41 +20,44 @@ // MODULES // -var sumSeries = require( '@stdlib/math/base/tools/sum-series' ); var isnan = require( '@stdlib/math/base/assert/is-nan' ); -var factorialln = require( '@stdlib/math/base/special/factorialln' ); -var factorial = require( '@stdlib/math/base/special/factorial' ); +var floor = require( '@stdlib/math/base/special/floor' ); +var log1p = require( '@stdlib/math/base/special/log1p' ); var exp = require( '@stdlib/math/base/special/exp' ); var ln = require( '@stdlib/math/base/special/ln' ); +var LN_TWO_PI = require( '@stdlib/constants/float64/ln-two-pi' ); +var EPS = require( '@stdlib/constants/float64/eps' ); + + +// VARIABLES // + +// Value of `lambda` above which the asymptotic expansion is accurate to within double-precision (relative error < 1e-16): +var THRESHOLD = 5000.0; + +// Relative weight below which the probability mass of a term no longer affects the sum: +var TOL = EPS * 1.0e-3; // FUNCTIONS // /** -* Returns a function to retrieve elements of the series \\( \sum_{k=0}^{\infty} \frac{ \lambda^k \log(k!) }{ k! } \\). +* Returns `ln(a/b)` for positive `a` and `b`. +* +* ## Notes +* +* - When `a` and `b` are within a factor of two of each other, `a-b` is exact (Sterbenz lemma), and `log1p` avoids the cancellation which `ln(a/b)` suffers when `a/b` is close to one. * * @private -* @param {NonNegativeNumber} lambda - mean parameter -* @returns {Function} function to retrieve series elements +* @param {PositiveNumber} a - numerator +* @param {PositiveNumber} b - denominator +* @returns {number} natural logarithm of the ratio */ -function seriesClosure( lambda ) { - var lk; - var k; - k = 1; - lk = lambda; - return seriesElement; - - /** - * Returns the current series element. - * - * @private - * @returns {number} series element - */ - function seriesElement() { - k += 1; - lk *= lambda; - return lk * factorialln( k ) / factorial( k ); +function lnratio( a, b ) { + var r = a / b; + if ( r > 0.5 && r < 2.0 ) { + return log1p( ( a-b ) / b ); } + return ln( r ); } @@ -83,18 +86,58 @@ function seriesClosure( lambda ) { * // returns NaN */ function entropy( lambda ) { - var gen; - var out; + var sum; + var wm1; + var lw; + var l2; + var w; + var m; + var k; if ( isnan( lambda ) || lambda < 0.0 ) { return NaN; } if ( lambda === 0.0 ) { return 0.0; } - gen = seriesClosure( lambda ); - out = lambda * ( 1.0-ln(lambda) ); - out += exp( -lambda ) * sumSeries( gen ); - return out; + if ( lambda >= THRESHOLD ) { + // Asymptotic expansion: H = ln(2πeλ)/2 - 1/(12λ) - 1/(24λ²) - 19/(360λ³) + O(1/λ⁴) + l2 = lambda * lambda; + return ( 0.5 * ( LN_TWO_PI + 1.0 + ln( lambda ) ) ) - ( 1.0 / ( 12.0*lambda ) ) - ( 1.0 / ( 24.0*l2 ) ) - ( 19.0 / ( 360.0*l2*lambda ) ); // eslint-disable-line max-len + } + /* + * Evaluate H = -Σ p_k ln(p_k) using weights relative to the mode, w_k = p_k/p_m. As Σ p_k = 1, we have p_m = 1/W, where W = Σ w_k, and thus + * + * H = ln(W) - Σ w_k ln(w_k) / W + * + * Computing `ln(w_k)` by recurrence avoids evaluating `λ^k/k!`, which overflows for large `k`, and avoids the cancellation between `k*ln(λ)`, `λ`, and `ln(k!)`, each of which is much larger than the entropy. + */ + m = floor( lambda ); + wm1 = 0.0; // W - 1 + sum = 0.0; // Σ w_k ln(w_k) + + // Terms above the mode, using w_{k+1} = w_k * λ/(k+1): + lw = 0.0; + for ( k = m; ; k++ ) { + lw += lnratio( lambda, k+1 ); + w = exp( lw ); + wm1 += w; + sum += w * lw; + if ( w < TOL ) { + break; + } + } + // Terms below the mode, using w_{k-1} = w_k * k/λ: + lw = 0.0; + for ( k = m; k > 0; k-- ) { + lw += lnratio( k, lambda ); + w = exp( lw ); + wm1 += w; + sum += w * lw; + if ( w < TOL ) { + break; + } + } + return log1p( wm1 ) - ( sum / ( 1.0+wm1 ) ); } diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/manifest.json b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/manifest.json index bc14a938f4e2..7f3d9702b45c 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/manifest.json +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/manifest.json @@ -40,10 +40,12 @@ "dependencies": [ "@stdlib/math/base/napi/unary", "@stdlib/math/base/assert/is-nan", - "@stdlib/math/base/special/factorialln", - "@stdlib/math/base/special/factorial", + "@stdlib/math/base/special/floor", + "@stdlib/math/base/special/log1p", "@stdlib/math/base/special/exp", - "@stdlib/math/base/special/ln" + "@stdlib/math/base/special/ln", + "@stdlib/constants/float64/ln-two-pi", + "@stdlib/constants/float64/eps" ] }, { @@ -59,10 +61,12 @@ "libpath": [], "dependencies": [ "@stdlib/math/base/assert/is-nan", - "@stdlib/math/base/special/factorialln", - "@stdlib/math/base/special/factorial", + "@stdlib/math/base/special/floor", + "@stdlib/math/base/special/log1p", "@stdlib/math/base/special/exp", "@stdlib/math/base/special/ln", + "@stdlib/constants/float64/ln-two-pi", + "@stdlib/constants/float64/eps", "@stdlib/constants/float64/eps" ] }, @@ -79,10 +83,12 @@ "libpath": [], "dependencies": [ "@stdlib/math/base/assert/is-nan", - "@stdlib/math/base/special/factorialln", - "@stdlib/math/base/special/factorial", + "@stdlib/math/base/special/floor", + "@stdlib/math/base/special/log1p", "@stdlib/math/base/special/exp", - "@stdlib/math/base/special/ln" + "@stdlib/math/base/special/ln", + "@stdlib/constants/float64/ln-two-pi", + "@stdlib/constants/float64/eps" ] } ] diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/src/main.c b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/src/main.c index 7bbe6179292e..c02abc9a439a 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/src/main.c +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/src/main.c @@ -18,10 +18,37 @@ #include "stdlib/stats/base/dists/poisson/entropy.h" #include "stdlib/math/base/assert/is_nan.h" -#include "stdlib/math/base/special/factorialln.h" -#include "stdlib/math/base/special/factorial.h" +#include "stdlib/math/base/special/floor.h" +#include "stdlib/math/base/special/log1p.h" #include "stdlib/math/base/special/exp.h" #include "stdlib/math/base/special/ln.h" +#include "stdlib/constants/float64/ln_two_pi.h" +#include "stdlib/constants/float64/eps.h" + +// Value of `lambda` above which the asymptotic expansion is accurate to within double-precision (relative error < 1e-16): +static const double THRESHOLD = 5000.0; + +// Relative weight below which the probability mass of a term no longer affects the sum: +static const double TOL = STDLIB_CONSTANT_FLOAT64_EPS * 1.0e-3; + +/** +* Returns `ln(a/b)` for positive `a` and `b`. +* +* ## Notes +* +* - When `a` and `b` are within a factor of two of each other, `a-b` is exact (Sterbenz lemma), and `log1p` avoids the cancellation which `ln(a/b)` suffers when `a/b` is close to one. +* +* @param a numerator +* @param b denominator +* @return natural logarithm of the ratio +*/ +static double lnratio( const double a, const double b ) { + double r = a / b; + if ( r > 0.5 && r < 2.0 ) { + return stdlib_base_log1p( ( a-b ) / b ); + } + return stdlib_base_ln( r ); +} /** * Returns the entropy of a Poisson distribution. @@ -34,11 +61,13 @@ * // returns ~2.508 */ double stdlib_base_dists_poisson_entropy( const double lambda ) { - double term; double sum; - double out; - double lk; - int k; + double wm1; + double lw; + double l2; + double w; + double m; + double k; if ( stdlib_base_is_nan( lambda ) || lambda < 0.0 ) { return 0.0 / 0.0; // NaN @@ -46,22 +75,43 @@ double stdlib_base_dists_poisson_entropy( const double lambda ) { if ( lambda == 0.0 ) { return 0.0; } - out = lambda * ( 1.0 - stdlib_base_ln( lambda ) ); - - // Compute: sum_{k=2}^{inf} lambda^k * ln(k!) / k! - // k starts at 1, lk = lambda^1, first iteration makes k=2, lk=lambda^2 - sum = 0.0; - lk = lambda; // lambda^1 - for ( k = 2; k < 1000; k++ ) { - lk *= lambda; - term = lk * stdlib_base_factorialln( k ) / stdlib_base_factorial( k ); - sum += term; + if ( lambda >= THRESHOLD ) { + // Asymptotic expansion: H = ln(2πeλ)/2 - 1/(12λ) - 1/(24λ²) - 19/(360λ³) + O(1/λ⁴) + l2 = lambda * lambda; + return ( 0.5 * ( STDLIB_CONSTANT_FLOAT64_LN_TWO_PI + 1.0 + stdlib_base_ln( lambda ) ) ) - ( 1.0 / ( 12.0*lambda ) ) - ( 1.0 / ( 24.0*l2 ) ) - ( 19.0 / ( 360.0*l2*lambda ) ); + } + /* + * Evaluate H = -Σ p_k ln(p_k) using weights relative to the mode, w_k = p_k/p_m. As Σ p_k = 1, we have p_m = 1/W, where W = Σ w_k, and thus + * + * H = ln(W) - Σ w_k ln(w_k) / W + * + * Computing `ln(w_k)` by recurrence avoids evaluating `λ^k/k!`, which overflows for large `k`, and avoids the cancellation between `k*ln(λ)`, `λ`, and `ln(k!)`, each of which is much larger than the entropy. + */ + m = stdlib_base_floor( lambda ); + wm1 = 0.0; // W - 1 + sum = 0.0; // Σ w_k ln(w_k) - // Convergence: check relative size of term - if ( term < 1.0e-16 * sum ) { + // Terms above the mode, using w_{k+1} = w_k * λ/(k+1): + lw = 0.0; + for ( k = m; ; k++ ) { + lw += lnratio( lambda, k+1.0 ); + w = stdlib_base_exp( lw ); + wm1 += w; + sum += w * lw; + if ( w < TOL ) { + break; + } + } + // Terms below the mode, using w_{k-1} = w_k * k/λ: + lw = 0.0; + for ( k = m; k > 0.0; k-- ) { + lw += lnratio( k, lambda ); + w = stdlib_base_exp( lw ); + wm1 += w; + sum += w * lw; + if ( w < TOL ) { break; } } - out += stdlib_base_exp( -lambda ) * sum; - return out; + return stdlib_base_log1p( wm1 ) - ( sum / ( 1.0+wm1 ) ); } diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/REQUIRE b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/REQUIRE index 98be20b58ed3..b10f12e99e14 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/REQUIRE +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/REQUIRE @@ -1,3 +1,3 @@ -Distributions 0.23.8 julia 1.5 JSON 0.21 +SpecialFunctions 2.4.0 diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/data.json b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/data.json index 47cb39a16fec..49c396081b81 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/data.json +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/data.json @@ -1 +1 @@ -{"expected":[1.0575457765396141,0.6470436116370292,0.03327965566239778,2.799401973338398,2.6196306581496795,2.2366184227234824,2.8531659700120997,2.312487216716022,2.6321705389324563,1.9419612209248767,2.477223457676331,2.6152660326049784,2.170658516362966,2.1254772650134997,2.724881399670359,2.5867874797490256,2.627480112888172,1.3091810799871255,2.458617686155188,2.71276132498215,2.720758084634763,2.6920618659542335,2.8688189360816523,2.6119570283323803,2.140534512162709,2.0998067114179992,2.7721366011626074,2.6440719325061792,2.6571057932324322,2.2352489072657225,2.8485810002839855,2.4043086031422973,2.4672093357702565,1.7892080304755866,2.686425423697422,1.6365753856511867,2.2486062853834694,2.8688927739242303,2.792090530874738,2.7159828379667275,2.7803923253324925,0.4196166468632322,2.396648461847926,2.638405958997282,2.71833651447216,2.891984802921556,2.2814204893411114,2.8439983406727762,0.7986289762786074,2.467947913560158,2.493445617382333,2.7498116632181,2.406641656639583,2.7080565453849346,2.3526484777117123,1.8420709374345652,2.393528478308421,2.227934743126822,2.746988597651086,2.7034492415884834,1.1874911077878751,2.6326734534073566,2.524674936529836,0.8670440970120721,2.9095700383438015,2.896062637108905,0.9648408493623966,2.4588935622718946,2.665543441999141,2.445634463562767,2.62140945937783,2.849122485500402,2.116865944372655,2.8462104520713893,2.526060719547587,2.307677851342901,2.7845115235707105,2.036334759358266,2.7896610873334033,2.0874246913515777,1.5511918320768934,2.797630105287407,2.1348295179386643,2.8666401229788647,1.9861054120667785,2.7289212685181106,2.4207167800189957,2.357569338236999,0.08590818323979812,2.80800773661624,2.849397386563126,1.2335961034408194,2.240216232298917,2.9101922123820003,2.4267410516461014,2.7250456625479593,2.694275766249433,2.9060655744151163,2.773468151531276,2.8069421103710006,2.209037776392863,2.8187061934885342,1.9242640155098103,2.744022544525162,1.8336808653159506,1.6938416500111417,2.621025426721408,2.905197808127781,1.5194825530243603,2.846269409891214,2.8390087644888666,2.607940710579424,2.5431018333482953,2.506511929334952,2.8452017849160924,2.283909847514825,2.559323603189455,2.5083333906025427,2.8890261094468173,2.52121390708621,2.516356916570272,2.753632165299038,2.8849088551850244,2.274766257909241,2.745521839633291,2.8262550958635586,2.4461006563922787,2.802362120627773,2.842886998120214,2.485085926309033,2.3831749317068116,2.381507487715136,1.6359915489637777,2.4138867474008716,2.9088367626540617,2.8704779179115363,2.7209634840188563,2.6689759064187903,1.8904754849744638,2.8867457825989433,2.5934789824358084,2.911824860135546,2.715620068356653,2.008557199584022,2.737987499212796,2.580920344821939,2.5793499798224264,2.2282633512261683,2.240150400981382,2.7588664170022525,0.4015271733819643,1.8968287258023246,2.8005656552625666,2.797999222777502,2.690839682986862,2.7905472986081605,2.6581554810003,2.811900188773226,2.5920645189278897,2.8150392794827077,2.540851461400349,0.3503234864782779,2.6023283881789254,2.434891592731123,1.353780018879244,2.7916480832678943,2.632410025921814,2.1741401523839454,2.6573856512030645,2.8224749750097686,2.159748633875626,2.418079408904135,2.797287090601433,2.4705617320750424,2.8397003856891345,2.820744224205775,1.5857656576403216,1.7025898350220694,2.5854787164590185,2.6030840927620336,2.895578369113082,2.170160892053644,2.7307656334806225,2.336469755999895,2.2697483693976634,2.658518814380912,1.9433274743815578,0.9746945533344427,2.7036332266689804,1.3917184607382058,2.465760133462675,1.5830926450630547,2.818863455121658,2.775269109960931,2.773929500651839,1.6865808706952246,2.7183508402758925,2.5269562322010035,2.0752802151308485,0.2666288178619782,2.7057109390702188,2.2856620140245765,2.8918995400370022,2.8103556033356156,2.1950800594316378,2.771951457501004,2.8189573751229666,2.4666514431552766,2.1830318989330486,2.84142067687673,2.305679710830278,1.9076970796921384,2.5063516281276534,2.299715374427814,2.3304841306826933,2.147216642907506,2.1564831136240086,2.5900409963496394,2.2829817225790032,2.887555419229045,2.6104564623877415,2.312677297081154,2.677677112874772,2.6307446133119186,2.593219749368462,2.680409732239081,2.318829014513735,2.778871628614972,1.4671790785989076,1.5979225419430487,2.622803427786213,1.5669745923086968,2.0505729057965785,2.147383587517248,2.857749080389155,2.3788827465976183,2.7480985122586894,2.8830910580214777,2.265151513541934,2.6902882542425353,2.4332436211135224,2.8834385884308418,2.547503387211604,2.5844816940010595,2.5040254197460694,2.600307605615118,2.903070032811179,1.854903951335478,1.6935803213778595,2.1393157503017415,2.1400276533638345,1.4047836624925794,2.708394725322801,1.46476760109723,0.8564718087615242,2.7523637634826663,1.9106349599380676,2.615974630567985,2.3101729369873327,2.199358970388438,2.4226267742791086,2.6287411143661252,1.5766524660268195,1.9696150989932235,2.802342284096266,2.620269776606788,2.4033443334335196,2.8657970187924775,2.7038534858300736,2.5747036066266435,2.8513755271852332,1.397589441963741,2.6426764771549536,2.4079721223205492,1.1433525675664593,2.624281140369572,1.3395339629681207,2.5455254782973835,2.6366642339596886,2.8265479807814913,2.7223880348228064,2.8596602497624843,2.685044685422014,1.4526484483474456,2.4859783203842927,2.9075481008009376,2.8616169701768897,2.666119837399627,1.5956047949058538,1.3121557127245467,2.6567515480750217,1.667875743837941,2.7831467841858775,2.5794166716197715,2.804513291040319,2.7135302480939174,2.750583924491771,2.7670797390390725,2.4328818951044084,1.5679943288585467,2.653773705309632,2.9037177175325795,2.5079549824980543,2.4422883757679337,2.8722843168925394,1.741110770824861,2.606306254849727,2.3345727445441726,2.422109626186412,1.729754479701369,1.3141970482576348,2.134053087943086,1.4056958802918293,2.8907398963425237,1.3825274631847126,2.8889744601361613,2.7944301060951275,2.442343904587487,2.8055860338752687,1.9490267506034533,2.2522578417054513,2.84952456804416,2.755159451514295,2.5110583450163038,2.7088303209598976,2.903844205331154,2.561110028617037,2.7902815875913483,2.5563564408206894,2.6140271091337137,1.0134402976242562,2.3478781074157276,2.524258127891148,2.9058084493359075,2.18634188833494,1.8616787182281556,2.446023074287556,1.817730589237057,1.9379485083344492,1.7837013347824806,2.491315311876134,2.793755033194607,2.7112970465201975,2.774568663604935,2.823989373648768,2.114551302564699,2.7884171822715302,0.3301465477636736,2.369034123423952,2.870265485854091,2.6960304282783163,2.837878210025025,2.603323662751464,1.9044667967857374,2.3054636966082835,2.7011295073694903,2.26679654669905,1.5706029126070673,2.766268888280571,2.8163688775043667,2.870738376913202,1.9004834747133534,1.4538168232733746,2.4779708155249445,2.7784152598874527,1.9038067312291447,2.853897994884157,2.441131537727241,2.4633912488836707,1.8143954530664448,2.7816525380206087,2.615052406103107,2.4538067884277606,2.0769254572681213,2.7824641899328135,2.105846728927757,0.5081762581527878,2.1092756428275816,2.9122565691286013,2.3100783400747105,1.3183005992917165,2.010235813810286,1.6711750456995391,2.609246637202089,0.98050550443473,2.6779189020412275,1.912200413888608,2.57273206827478,2.2069888938585254,2.7050703799501505,2.8027972474729523,2.90781612526753,2.485935207670268,2.289122208149081,2.764661867926975,2.891103020081381,2.259191860715658,2.305838372635206,2.8099489774022643,2.8043815387104054,2.147591527705324,1.1806493889592922,2.7361011209503836,2.8851829378708445,2.885620704942667,2.3522853818106153,2.801831761803964,2.7757490640770257,2.641494490102964,2.6171026454713857,2.739860722533436,2.809862746140549,2.8998903084023198,2.13375835150578,2.8269373978405916,2.48756101736536,2.3952428546131,2.3178378411918628,2.604006601254927,2.7926066862111227,2.0811064232718897,1.518464898349316,2.880937948806128,1.9932571071774992,1.5813903113186716,2.433052863826779,1.750059167905037,2.4491553959149233,1.5122525186906746,2.6585805776380056,2.879117001436086,2.6809450708382236,2.365729364615496,2.780913760722516,2.5473125381794866,2.3875595034076547,2.7292150290876194,2.5765315770929895,1.8558621210382582,2.79061006072957,2.228223932625453,1.939925020079063,2.720068421111005,2.721244704227555,2.6745599306916077,2.686747105860377,2.5397558102358886,0.04763364287661594,2.660693580535632,2.804319654172115,2.4252861239653365,2.192898844975382,2.648132573354733,2.667144802700605,2.8608889676056393,2.837773023025253,2.248991088667581,2.4522619201455917,2.195369243329398,2.347999461792199,2.023621847239233,2.536462368147756,2.5615382455052345,2.7094725824256933,2.6783799727765327,2.4191307794296817,2.515167953876853,1.8524663071441283,1.82982568503696,2.8864049724891956,2.5305270031179052,1.7886861239214624,2.2787678983721014,2.1421669822563256,2.870553305774564,0.7692039922860429,2.4408624040210807,2.2186776009173306,2.6089446515530206,2.888640148759528,2.8313786673198926,2.114151481435882,2.800359204689439,2.82376091055211,2.0505384909793434,2.1123681059708974,2.837083045369681,2.0917279461708773,2.8187141034757133,2.378342834147147,2.6404217895148143,1.6674501110038182,2.6913659589388246,2.5471334047729517,2.8495430814328913,2.321885895096801,2.608840351593317,2.886623712618892,2.0537208395635833,2.0004530541107437],"lambda":[0.6417250092665938,0.26943920029941726,0.005337722116496302,15.985930662238328,11.212261523308712,5.318171219583796,17.780697487727966,6.154810517198985,11.492398294300292,3.058120203586596,8.479408066270079,11.116396801579649,4.688495962872374,4.303490710214795,13.79700979664031,10.51106223259234,11.386790932248175,1.0075931571814678,8.176637400509549,13.470869959535943,13.685165884987072,12.931845503499527,18.340567281465194,11.044274484135093,4.427878504259404,4.100042166375224,15.146940289721451,11.764860548171967,12.070796320874546,5.304220473332095,17.61999725809719,7.354970765491027,8.315042987863063,2.321334296819404,12.78889872740255,1.7749894910044262,5.44195753971163,18.34325006358828,15.75643674942186,13.556787084152578,15.396164076587882,0.13888418038666028,7.246135815780983,11.63433737795113,13.619910458691638,19.202012584993582,5.7966077864179955,17.46084461313594,0.38280646985544564,8.327052879372312,8.752784105578314,14.493262928476351,7.388453596542446,13.346385888249523,6.652506707865142,2.551415916853048,7.202286220751333,5.230367391487412,14.412664397541711,13.225612160518349,0.8131332831855431,11.503780693087538,9.304801501941839,0.44232172987561214,19.88312272205833,19.357827356387762,0.5379341482629307,8.181044541813002,12.273149714001722,7.971969149593585,11.251572251147767,17.63889927085453,4.234051539526789,17.53748751974131,9.33010775931152,6.097874093533653,15.522062991553597,3.6404705060077847,15.680920940561297,4.005655852221928,1.5308457660176922,15.930005993060092,4.380299473662133,18.26158168652377,3.316661993854666,13.907489372414577,7.5938233283392975,6.716310979302231,0.016887773945302165,16.260389227502003,17.648503280960952,0.8825022785367098,5.355006173795243,19.90766244257973,7.683516751898689,13.801484539321205,12.988435417480932,19.745469235690713,15.186859092426182,16.226146843915284,5.044565084572503,16.608246462828156,2.9607848143233495,14.328472829628357,2.5132786096729465,1.9616426142115273,11.243073434688565,19.71153267088431,1.449309575259572,17.53953487367855,17.28921170370052,10.957375777455013,9.647122355637823,8.979545033687879,17.502498197590835,5.82448733619632,9.959130040225247,9.011630813582459,19.089753605527264,9.24190493759243,9.15437325525903,14.603066134444855,18.934638223482523,5.7227703979372535,14.370967915598914,16.858217438512174,7.979226150948264,16.079804022051384,17.422467841390343,8.61079446781336,7.058732266999095,7.035891562505148,1.7731873115981855,7.493439765259495,19.854240134278477,18.400939483397977,13.690715468557837,12.356451313369572,2.7843107070013406,19.00368563066564,10.650206031238604,19.972202607805258,13.547084377261447,3.4573165352377444,14.158703462054433,10.390586588524972,10.358580395427177,5.233661943361678,5.354329762212093,14.75487114797021,0.13025609306086317,2.8165852147743475,16.02276747366854,15.94163993249435,12.900712370123495,15.708425028019416,12.095784282296332,16.38609005856924,10.620637782187622,16.488177945266436,9.604635198609976,0.10711975000242902,10.837110832828412,7.806610824825118,1.088800417547211,15.742656530990526,11.497817142319517,4.719661917802576,12.077453242823047,16.73257215211963,4.59225563774456,7.554897268848211,15.919202470462697,8.369699177494248,17.31290010397768,16.675361179081527,1.625169315232431,1.991968350364437,10.48406519876487,10.853225925809324,19.33925679444947,4.684059496182735,13.958225671987714,6.4471439141002485,5.667745222152396,12.104445730268116,3.0657802957624547,0.5483217026530474,13.230413732411442,1.1626897448705842,8.291529479887675,1.6176654335989005,16.613415529864714,15.241019915851863,15.200714802363148,1.9368515375551931,13.62029557469237,9.346498488282645,3.915368658572218,0.07327680684219295,13.28476002969495,5.84419518450658,19.198768225164482,16.336092793265784,4.911840434865327,15.141398253103091,16.61650338296516,8.305982996553007,4.800267023823812,17.371963401909163,6.074381185571949,2.872759574894248,8.976726875053757,6.004818420792453,6.37285044572065,4.484315311858107,4.563862892859709,10.578483107383185,5.814076366564542,19.03419935935938,11.011725755769861,6.157072210552861,12.57020239578766,11.460187788922305,10.644780691921806,12.638103276930055,6.230740469162059,15.349947314756847,1.32434608790581,1.6597598405193592,11.282476410644705,1.5731773897527335,3.7384479558774464,4.485735185082715,17.94281287960805,7.000092928368651,14.444298136894966,18.866559282936276,5.617825916299024,12.886690540435048,7.781559190156631,18.879555624122318,9.73078038708493,10.463546198344233,8.935933105383459,10.794138163569475,19.62856868413109,2.610990279071048,1.960744421989551,4.41766760567075,4.423628911524697,1.1892291050644888,13.355294749632346,1.318853075109483,0.43275023981676153,14.566518449319773,2.88815614444895,11.131903509878246,6.127343415120334,4.952128898382959,7.622142991948073,11.415085563638868,1.599734366158856,3.2173654956314257,16.079173105344186,11.226369637257744,7.341177917257085,18.231110050881007,13.236164298210502,10.264468898084637,17.71776809576083,1.174544025713229,11.732576758444232,7.407618322206928,0.7511428520518804,11.315331556857933,1.062220283880042,9.693095921584316,11.594511244752903,16.867992207989445,13.729267563288223,18.010855943308886,12.754126932771763,1.2915885809856675,8.625838417757116,19.803584312548242,18.080790882488124,12.287098128216552,1.6531059606691034,1.0128294781320202,12.06237531616333,1.8745259573955542,15.480236145592464,10.359937618101181,16.14837257363756,13.491326666531279,14.515390449806977,14.996305262614689,7.776071518464311,1.575953889770183,11.99182270382909,19.65378514403396,9.00495534397599,7.920081345922636,18.46690468424649,2.131682505298458,10.922211756966608,6.42350093373429,7.614464440239219,2.089417808120042,1.0164373781149205,4.373866831921442,1.191103787269041,19.1546972411064,1.144362462538453,19.087799819461743,15.829507030039554,7.920939586105824,16.182676723289752,3.09796330511654,5.480265855962156,17.652948311884874,14.647196625620854,9.05985153957169,13.366778779637048,19.65871353825716,9.994114485659523,15.700173409071635,9.901298874973229,11.089337481624465,0.5906157827265757,6.591255664805171,9.297203771445153,19.735407455838384,4.8306474051909865,2.6430588440775793,7.978017995492737,2.4424977457089136,3.0357442516806366,2.2987114713714796,8.71637417918846,15.808387828211297,13.432000755688529,15.219931954813184,16.782794293538117,4.215594075111864,15.642397820566654,0.09851503157033203,6.867444604567239,18.393197650284442,13.033465211977528,17.25056013352567,10.858339755423149,2.8559352268781213,6.071847106874948,13.165224136992574,5.635636630270331,1.5830796548793113,14.972292976822082,16.531612179597506,18.410436093437795,2.8353385692884148,1.2941925208043425,8.491807800204985,15.336104853716247,2.8525108100702568,17.80649125871384,7.902223310143737,8.253241131412633,2.4279771488398527,15.434570989491556,11.111726192261049,8.10017452294729,3.9274694403367327,15.459358665835747,4.1469571273597605,0.18469164254584225,4.17384919984205,19.989303762029582,6.126223421435588,1.0237259763137363,3.4680929519766757,1.8853607360027924,10.985554647505896,0.5545162451555186,12.576195466102117,2.896397454614381,10.224799868108313,5.024845345282705,13.267980942244577,16.0936498813171,19.81410929217474,8.625111005078022,5.883320456315531,14.924818178917722,19.168486363966828,5.553796244899276,6.076243164437396,16.32295624580077,16.14416447627468,4.487504404124412,0.803250750832345,14.106057919871592,18.944924487668903,18.96136549704393,6.647823831030739,16.062944151601066,15.255486715732228,11.705302053346362,11.156634057897739,14.211178982282968,16.320171813252685,19.50524466535413,4.371427637850802,16.88099752682124,8.652585871094356,7.226346691317218,6.21880892248702,10.872931153897106,15.772528055990186,3.958403950799405,1.4467677531752798,18.786242163307357,3.3607705761544127,1.6129053997505505,7.778664752571749,2.165642408842272,8.026946499703662,1.4313491322513405,12.10591872188012,18.71858518630713,12.651449102474484,6.823520769119269,15.41204387035721,9.727137681627443,7.119159829081179,13.915557840340341,10.301389564955926,2.6154996004803133,15.710374737824374,5.233266625096009,3.0467431354910524,13.666548897415623,13.698317319510572,12.493198149586707,12.797013622144103,9.584018554215877,0.008204298669394205,12.156421629798757,16.142188300940198,7.661755124957486,4.8914376383955505,11.859318285836643,12.311941469003965,18.054739425450524,17.24696842880276,5.4459811247600465,8.075777300318032,4.914552230581397,6.592806541988159,3.555359651207204,9.522318634512299,10.00251910371274,13.383729623373334,12.587631771930035,7.5703901252998,9.133075686747908,2.5995569360434656,2.495965997791343,18.99085589567262,9.412148565780791,2.3191797272994252,5.767054440371462,4.4415953176849765,18.403687695476002,0.3588982714092914,7.898074663233641,5.138455610900321,10.979031940862441,19.075158397131343,17.03004130588626,4.212414634339883,16.016225937252557,16.77520799809891,3.738207689817812,4.198264500191318,17.22342725728618,4.038188842430808,16.608506419072807,6.992752563170286,11.680605005001583,1.873133069050672,12.914109057937498,9.723719858808618,17.653595454707105,6.267689522868869,10.976780041988254,18.999089322625654,3.7604966438698018,3.405804007846518]} +{"expected":[2.4629889243949132,0.6490466839180405,2.8486843584085615,2.8315757076578763,2.432513441863522,2.156440673495002,2.88228238537559,2.732582636409852,2.702131175619679,2.2390875951429994,2.2577622719574744,2.7372014946653698,2.6873441896799894,2.829643495198632,2.8318288774340883,2.6780887407063143,2.254859268116933,2.630181949485469,2.5484096042251188,2.7407362787101337,2.555421330008701,2.8664059592233495,1.1752638696297892,2.452268346231772,2.8998897482938033,2.476405592784444,1.8401701390742453,2.7760683590951776,2.4913123929584367,2.847898893534935,2.7655634012316703,2.1893441725624707,2.6653063883386006,1.966634165139636,2.401054829585707,1.6143469074426506,2.911205132062997,2.874167477814133,2.7553173777642757,2.460091646544711,2.682127662946939,2.627495131098585,2.4282688003235147,2.495278498520831,2.59201374422612,2.4552420591509914,2.8101751685298555,2.0209971507902007,2.5948203066831717,2.6212053578327748,2.3405599703627082,2.7818365748316864,2.781571758330083,2.296785261125609,2.81915480868648,2.6852317280098585,1.749253578527715,2.826211268386401,1.9395160675653293,2.4239578061765013,1.835098615399897,2.854123616965943,1.9632509919990133,1.9337238243209283,2.395037307781092,2.7523094210569314,2.474046884682675,2.761420645609898,2.6568106849347055,1.2786336101817755,2.7170355920442244,2.4194698541375628,2.1068644478890057,2.801204306759726,2.890029854812679,2.825344822788038,1.8193252442775603,2.553741579493865,2.605172202081624,2.36818638628223,2.4394615090936034,2.27969776432806,2.7679416087909714,2.6484115536679456,2.3009493854915486,2.774742349191659,2.6747368076087805,2.907324656701854,2.850961953733042,2.731061150266681,2.169935220474274,2.687635527223334,0.03510614846263468,1.564354480132895,2.8428851838280984,1.9743963059749186,1.8197924108446952,2.8487435060677155,2.2172393231793155,2.8475422058251496,2.888758755806742,2.4436276186491273,2.8776226923604176,2.1282082146834216,2.103816185884134,2.4230555559253366,2.7726489227129587,2.602774841204006,1.3165052114761484,2.724371512196043,2.7263903486253094,2.3939488888821137,2.225808067978727,2.718653733461336,2.666302263453715,2.4895536580716264,2.2228333469562536,2.6996593098022608,2.7085301091301854,2.7393099577824302,2.635345160616009,2.755811447760469,2.787270402441578,2.769183569185367,2.442899911260982,2.7432409207504347,2.700119449181164,2.429163178871473,1.8268318192810842,2.2627215154619744,2.403522458830795,1.8676265653728565,2.498354074141832,2.3266504324435426,2.7347402262886322,2.5255919442276684,2.329787449734913,2.69020351893427,2.843298975421366,2.159832604603069,2.6014553625826173,2.7075688258310295,2.537979555657497,2.739621246091708,2.5580580115354805,0.11058331365389691,2.852929529339114,1.9738047480124166,1.9646587697096742,2.8629040055934207,2.343687427622616,2.186789535349518,2.562713906557213,2.2086961133869925,2.720100652444552,2.243179086156389,2.582842653995932,2.7836751454404953,1.8295125675015496,2.850416432746892,1.3884984472153077,2.2013850490026057,2.4613869214366657,2.5826394596232567,2.8043622617813257,2.4492956422679417,2.675719087682463,2.750038482461368,2.2068804989628257,2.6978384412813865,2.8658324408371247,2.8266305352325753,1.9992592212102023,2.615018851979965,2.6341012339655614,2.381886897393737,2.7925994488656953,2.3624592497382753,2.627098895578477,2.669026511583448,2.4197588379752233,2.755056293316249,2.6768450651003586,2.218564537317703,2.7559557766770784,2.7429172201545375,2.8349656691073197,2.4802219264168923,2.5571084779854547,2.792710864478139,1.8870076225434695,2.3467880445917437,0.8788114393142497,2.167368525382205,2.814562644979234,2.078557727632965,2.1335833435841565,1.6349871860772718,1.9527631311000282,2.4248714022001208,2.6909405960565715,2.687264294789066,2.873437465853311,2.045303838588533,2.5308215707245814,2.0886508800230668,2.6415026096596352,2.7422459878857763,2.8119509523934982,2.619582188058162,2.9005012606618976,2.007943385006302,2.570111211314345,1.9566217512488009,1.642117199400533,2.810997638399134,2.6296825879639414,2.3103051058170156,1.911310535713477,1.0890846361023283,1.3925791699029315,2.0324556639249893,1.7378409869823674,2.799093256373254,2.6963539753341648,1.7669128291062464,2.743595932464178,2.1917196230196403,2.2774279736567253,2.4546297676326096,2.491125838248946,2.1238370003631504,2.810263813251357,2.784551245413457,2.820578618114591,2.292232526172072,2.730658956766529,2.9005305705936197,2.1772149277075226,2.3493498132821986,2.8578402530021783,2.118805384371396,2.8167633140708705,2.8859655275805,1.5630009048575004,2.49373519786252,2.313913660466024,2.6463364676054635,2.2354313527834964,0.5539501286374111,1.9319003440244462,2.764838594329933,2.412775703761378,1.9358198850412998,2.7701296187681286,2.297610981066386,2.4148515961242634,2.384727712341081,2.657188841588983,1.3944307923559627,2.881424294248621,1.4330408798365721,1.3741872238577764,1.3610830742421178,2.833115682525379,2.871685815106433,2.7843136224751395,1.7957312845612314,1.9323284123509876,2.7161409686088023,2.799518016573457,2.6669363195285656,2.08000664192582,2.7986908692121704,2.409933694852188,2.5475675519132315,1.9906583553147559,2.667548484222183,2.534295050040528,2.561133882962611,1.846273751228519,2.7970376380187285,2.6579320318693513,2.615979602239625,2.8660672318487856,2.8627465771606837,2.6920728934904647,2.6456230443018467,2.423705015332769,2.3942230813994705,2.8153197982721734,1.0995317361727568,2.7170889123714828,2.792530772053622,2.1958774968734422,2.0454269712253565,2.8367336462200408,2.877687583934946,2.758154427852524,2.771755086876001,0.7757344465423851,2.227610764365306,2.6269574110657405,2.907522826500943,2.826113150255876,2.605908266846303,2.7516106719486184,2.489385188349811,2.095459961275355,2.4326617824698076,2.819685482488171,2.0537462778495112,2.513561553655822,2.7866667705614323,2.5953059930563014,2.8466168075527727,2.3016596287135944,2.8187748877840666,1.8397740245098109,2.743983105135291,1.630440180766161,2.6087160604852273,2.8251446687528565,1.4514633166764335,2.6550400291738416,2.8526121012025607,1.9074498685428731,1.887585733841246,2.7593105631089516,2.8893348311335725,2.643296906293974,2.8029364655091613,2.3621246017332735,2.213880842385037,2.099782668699738,1.55274584730499,2.674422469146063,2.0295571069034963,2.001136187488677,2.846143626405479,2.8827902994230037,2.3702132412118555,2.462823929891343,2.061016475599032,2.7513685915824992,2.796925225811483,2.844844829548351,2.0806567623200163,2.5668819630642234,2.723274346380575,2.8932026168558553,1.6693717503890033,2.1032240805833173,2.740810701583887,2.856899928225367,1.9500055599652864,2.7520518579194495,2.632858134184459,2.647740844133197,2.127120086442942,2.2165080025143715,1.9542412708259127,2.819670888850559,2.8368300754749485,2.300299071355497,2.6499955840264127,2.4193574770787527,2.727646456486939,1.0949244457359408,1.9067786995725966,2.836061821167255,0.6647586185924594,2.600619466662586,2.0785404183731937,2.8408019853055797,0.7023977448743433,1.2599901759012289,2.6654020598657917,1.8051629009398769,2.099864621159009,2.772618968684151,2.8244581015904964,2.596469064351948,2.5769257362169267,0.9302487998079609,1.3906369499308062,2.8757330583684095,2.474444388708199,2.712630578774029,2.1799028300333227,2.2881058069602447,1.8970348860515225,2.6226036405727355,2.0331384608908247,0.4320616825825225,2.8934015273524625,0.20984258193149333,2.5868768177227364,2.834934010739965,1.6517099504516926,2.015715787791685,1.6118981530550984,2.570938206703597,2.7516407689507405,2.744965778661612,1.9963954978061174,2.3006578976331853,1.9230058811590616,2.685448275374594,2.690684546476746,2.800239429142434,2.866866677367626,2.8797963239299054,2.889800970578229,1.5784242742981807,1.7653203393276364,2.692470876397463,2.8099412688685366,2.900115164323219,2.8150556257865618,2.1118740591953125,2.726352351302078,2.33715334495968,2.7493572012191914,0.9608632607028332,2.6724505499776896,2.910868033175069,2.72113795846131,2.701653021633275,2.564293020581775,2.795606704270818,2.649873514247681,2.3963949539878824,2.6809774377115865,2.6350273131249526,2.054274166440598,2.2480036928929956,2.631935438436224,2.731235161871219,2.8268017144418147,2.4053483321590976,2.459327802340719,2.147672757693727,2.767855318216674,2.5755560815312606,2.0328911378402736,2.7811060079684737,1.4653017084246984,2.8483372894901886,2.777920732267622,1.9195072092710947,2.700446403619882,2.4034812656278906,2.8961701460356024,2.8350106105599093,2.602474140600425,2.7769095895336933,2.6488280083711984,2.90118289626075,2.656867244762428,2.547570253276452,2.2753358693106662,2.502753365381319,0.13949752700180174,2.905668161295761,2.5437680271087815,2.431114095032338,2.34768043136307,1.520585000223515,2.6469778988315467,1.9773983775941721,2.654674308962112,2.7286787139234487,2.9009873801968262,2.2477183660601616,2.6762958242523163,2.632420747841483,2.630062924895635,2.718752792197358,2.8247945701183217,2.7292469316596395,2.274061550635082,2.87162697620915,2.0523457339765034,2.0601357011290427,2.6647791535827965,2.732526126469152,2.7765135879734335,1.9299046944008558,2.5475500521383982,2.145852622229778,2.8162701911884316,2.5761068128343343,2.893476416748988,2.136553900296061,2.619973629671793,2.7640769468306403,0.9237661638955924,2.3288055567515387,2.398712180204382],"lambda":[8.246756435708562,0.2707869976718702,17.623603672769445,17.036684487308058,7.770485892644034,4.5634951310927345,18.83635290151152,14.008392804709187,13.191265461749417,5.3434220208226435,5.5385517000570665,14.13674327558631,12.812089782772873,16.971653242918133,17.045223909294315,12.58040688176026,5.507730976621654,11.4475029658937,9.748096250529766,14.23577393316894,9.883144391649248,18.25311333750156,0.7955436726273546,8.075878626991397,19.50522301087118,8.46585983453365,2.542720170413939,15.265118622718717,8.716324397536042,17.596215601746444,14.951432582709003,4.858381198544635,12.2674178326237,3.199767666729254,7.308536098179124,1.7077264606687304,19.947679320881,18.53592736985725,14.651767560421241,8.200212399850951,12.680981520836285,11.387127493752061,7.706436241183711,8.784235513173602,10.619577934681836,8.122908923285376,16.3302623070252,3.5380606293244043,10.678322903355213,11.247054577401016,6.498430672993257,15.44018790470836,15.432106197623206,5.970950504457518,16.622996394066558,12.758831687926246,2.1625610726916844,16.856755228830462,3.0444637935478047,7.641942777721487,2.5196786310476837,17.814448943015144,3.1799236506646,3.0123814527841453,7.2234575517303075,14.564954700707318,8.426911401548772,14.82952896264967,12.0637806835862,0.9552341622535132,13.58498418423057,7.575393679919666,4.154919200964966,16.043020732891268,19.127763326140823,16.827874446347227,2.4494743758006465,9.850618660478785,10.89788110252573,6.856149325332863,7.876516346959919,5.777396292806269,15.02187117490007,11.865836115118668,6.019142036673893,15.225158256789243,12.497554715187674,19.794814260415517,17.703264096413708,13.966372417113558,4.682049049938288,12.819452501603246,0.005688660447991101,1.56606672802321,17.422405260190235,3.24581428475913,2.4515223896179528,17.625667813094136,5.12432991384566,17.58379260223925,19.07964235524777,7.940807102465309,18.66324826072353,4.325768333232938,4.13112060403134,7.628515551170823,15.162286683609121,10.846628311756392,1.020531148534698,13.783129143164643,13.838170989580282,7.208178799251619,5.209098924810394,13.628440741646248,12.291516082258344,8.686382458115196,5.179503104165059,13.127096476512055,13.358862931478209,14.195729710952712,11.564441766914339,14.666076999644652,15.60696749142794,15.058789493864573,7.929538378199679,14.30636950027586,13.139017527167391,7.719886554257638,2.4826120343994607,5.591624879450148,7.343723794554384,2.671572059790992,8.837271739829845,6.325736501634445,14.068200977497254,9.32153941646809,6.36426147985492,12.884537247822623,17.436684246660658,4.592988234099229,10.81852427514522,13.333548240814672,9.550691257427534,14.204459435770094,9.934421070078741,0.023168727358360286,17.772374288005274,3.2422794634011454,3.18816444897712,18.126940199970317,6.537932317442028,4.834771826499407,10.025630956495213,5.041270991843218,13.667418387364158,5.385543061351116,10.429903290697226,15.496415929597504,2.494565653604035,17.684151054089828,1.1562370705485359,4.97132835126348,8.22098717208135,10.425740222126336,16.14354887475055,8.029144436047579,12.521776860172748,14.499758408729136,5.023804371594611,13.080029952121038,18.23238924700755,16.870748423421897,3.398287136227889,11.110992765501742,11.536158033956402,7.04108198538578,15.772302314516164,6.780344627186363,11.378251156007623,12.357683725875566,7.579660764538724,14.644211648822656,12.549601132272656,5.137343687633362,14.670259779156446,14.297225780664473,17.151387638059198,8.52926961588889,9.91592383778852,15.775777863279018,2.766866545293043,6.577341967894865,0.45314276158500144,4.659248001348983,16.472635755500402,3.9395150314752625,4.369979992618713,1.7700916081398543,3.1192643633782136,7.655563711475701,12.903280082331616,12.810071419499485,18.50913968115955,3.701857773130337,9.417585339924715,4.014896857369267,11.705489194338831,14.278284046212752,16.387735869938897,11.211192315198836,19.528879286475863,3.4533851338215427,10.172308172015159,3.1414320373758176,1.792193257126522,16.35685629863847,11.43625721510501,6.128908605811631,2.8917095208030363,0.6803397420141377,1.164420396979453,3.614266451875421,2.1194186726427433,15.976172495689358,13.041785672790862,2.231219774229376,14.3164044945578,4.880445016755273,5.752186193251689,8.113202326214992,8.713143343380715,4.290169856749904,16.33312646814439,15.523282107673387,16.66989735430217,5.918724331823859,13.955286014374193,19.53001386442567,4.747370271962838,6.6100893393495745,17.946052976275205,4.249584331992651,16.544519559926222,18.974325980170935,1.562406372721843,8.757745479503726,6.171804462400767,11.817442888862642,5.306076738467002,0.2108012295851247,3.002358880769671,14.930031581911066,7.477240438083712,3.0239479047629514,15.086973262573997,5.980474328080363,7.507536950223006,7.080071440428841,12.072771385091098,1.1681519106037408,18.80435399483967,1.2486537341071235,1.1279711748088639,1.1026678899770048,17.088694818640494,18.445022467919756,15.515990588245472,2.348450346503528,3.004708402184826,13.561018683224532,15.989600206728012,12.306884074129247,3.9502411877424493,15.96346257631857,7.435967643004053,9.732005401881034,3.3446683965162016,12.321740000053529,9.48193754099941,9.9945824854046,2.570759638469813,15.911350418305252,12.090460642660535,11.132012385989679,18.240870546428805,18.121288895740214,12.93212675913532,11.80085150818982,7.638178303277876,7.212024638953238,16.49733208896843,0.6935305397065639,13.586413904675451,15.77016036520633,4.9193222015124665,3.702708368768093,17.21151858919765,18.665647876122918,14.7341282853177,15.135522395836858,0.3641204808666343,5.227121405181851,11.375083373424946,19.802592114153622,16.8534821929773,10.913666783219453,14.544862761517532,8.683519856907589,4.066635492541188,7.7727341879742795,16.64046142061235,3.7606753954121164,9.104381266600122,15.588350630328216,10.68852252784823,17.551603420872617,6.02740232163073,16.610504190550984,2.54091215172346,14.327356704945542,1.7561482550087781,10.974097136751249,16.821209933081086,1.2889525978625938,12.021773693232738,17.761206353123036,2.8714681569704226,2.7697661715647737,14.767825738866364,19.101436116729765,11.74691921375497,16.09808238535895,6.7759422291270965,5.091505806805232,4.099856571284222,1.5349610062310948,12.489813469564703,3.5948209196902337,3.4101135187514653,17.535167241067583,18.855319325815117,6.88318735378912,8.244098557680669,3.812144170941585,14.537908497860892,15.907813285462058,17.490132800415033,3.9550642769853073,10.10800949132689,13.753309033587541,19.2484123977703,1.879430480359714,4.126515027428043,14.237866509437415,17.91266404479767,3.103527926838119,14.557545420061576,11.507963429229779,11.850172444658435,4.316877087754547,5.117163250121739,3.1277357833161568,16.63998088075048,17.214804382344255,6.011589145361961,11.902913019080954,7.573735010004816,13.872530137810035,0.6876881698388182,2.8679652402433153,17.188643974794488,0.2814905421462921,10.8007586499095,3.939387084509247,17.350698252032974,0.30810325771643665,0.9245022842334727,12.269730811901614,2.388273619042074,4.100489231570624,15.161388989026829,16.798369633382674,10.712987923500922,10.309368392182039,0.502597921522796,1.1605186264210143,18.593508250428506,8.43346230849512,13.467394672236795,4.771734893334538,5.871799302984271,2.8176394140065053,11.278041844639386,3.618863868526281,0.1449605436051149,19.256001815433528,0.053064099279027666,10.512907670604331,17.150312838031724,1.82239395850196,3.5035308710163626,1.7004832748294008,10.188841830416065,14.54572759627042,14.355192387995054,3.380330153053097,6.015755403400824,2.95399643182805,12.764280789840864,12.896765965536874,16.01243200081835,18.2697786263729,18.74379646846746,19.119089242324804,1.6046467387490604,2.224931425556491,12.942281410632361,16.322707315442305,19.51393983215101,16.488711231686636,4.19435360303625,13.837132969543632,6.455685902410952,14.480257256272377,0.5337821760219374,12.441361834265868,19.93435272941721,13.695431299536036,13.178827920952951,10.056759991377206,15.866384391134755,11.900051584316046,7.24256261185559,12.652256457382293,11.55720801668159,3.7643868943120373,5.4356630146642555,11.487081227365527,13.971171787896932,16.87646497267294,7.369873015500355,8.187986510911466,4.488195737634031,15.019309523210413,10.281670033519406,3.6171978509791214,15.417902751130255,1.3200677369008118,17.611496552197888,15.321119234764032,2.9352092399327145,13.147494769526526,7.34313495489058,19.361952517479292,17.152913515132443,10.840217044605623,15.290524815031535,11.875572544064134,19.555282407144198,12.065124964471147,9.732056979021309,5.729052087099986,8.9137060375459,0.031137156414311473,19.729919907667025,9.659736968122308,7.749309947446632,6.58873022991334,1.4520683271488743,11.832380277047488,3.263818547433883,12.013115905005634,13.900830910488548,19.54770541866093,5.432685249134915,12.536020884251569,11.498059806928271,11.444821479966327,13.631105623865059,16.809559142844886,13.916434366427843,5.715008896811256,18.44287263490946,3.7508477437912235,3.8058680856286586,12.254679218734953,14.006829824019544,15.278559752357047,2.9914321850485948,9.731671285369002,4.4727323314907945,16.52838440233083,10.292798263652191,19.258859998531065,4.39462229480077,11.219830144814097,14.90757624001563,0.4961661579745913,6.352176814179254,7.275291993760115]} diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/large.json b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/large.json new file mode 100644 index 000000000000..056800c8425c --- /dev/null +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/large.json @@ -0,0 +1 @@ +{"expected":[7.636170468247301,7.972537561030521,6.7610278260932,6.917374059516255,6.727728616304167,6.750677706556752,7.244541238871671,6.509536989419537,7.821177896107505,7.269680611911745,7.970309910988387,6.524486872450528,7.21123753063748,6.453430815347789,6.710260269147566,7.591407693115048,7.719570138214441,6.60409923648532,7.275309925886296,6.569444099784174,7.220662662609307,6.195700902330057,6.549083061504953,7.650034426380673,6.943061039654225,5.758624091480383,6.314255620964442,6.854093883650863,7.092309409089726,6.267512362927276,7.234898393348585,7.645304883181197,6.77246469767916,6.690957199656208,6.370033457060966,6.414296508435013,7.346668967419823,6.369298474087861,6.498332839194182,7.529607127859304,7.78693387327993,7.0160844635493245,6.928044277357777,7.139826340202468,6.062510702123153,6.886249261145232,7.1787683992812275,7.110396468179512,6.486337927437986,6.644103590157126,7.433638895070484,6.570855743395127,7.434170413341826,7.797867116769563,7.102138962481776,7.11346739600779,7.962627037349753,6.861612790136011,7.27952065202966,7.193231082968424,7.239885488716427,6.893343446896546,7.485664259242629,7.957547318248099,7.470112622815086,7.475431581084728,7.787995523408879,6.753984879486759,7.6407055038951635,6.089444968267598,6.024594986004413,6.359695233109921,7.13377710402349,7.710607191114309,6.798541224414055,7.564456723270363,7.13716031513703,7.397293948193246,7.121284650541771,6.727903509452129,6.663396008687282,6.975561358374909,7.545200793938018,5.746414776130449,6.33421494448371,7.499186674598618,7.764468119811527,5.97851871881829,7.0678479033698505,6.733128250255487,6.848576538876062,6.072626808270496,6.7933007812668835,7.575245446153155,7.819236694491278,6.661459830661217,6.581854920708985,6.638681032060693,6.631801879734862,7.259650282374132],"lambda":[251315.53326565857,492474.71298286394,43659.771794731394,59687.2581321705,40846.819797937875,42765.2973878073,114829.86335154864,26402.21600383054,363843.52720370906,120750.96191720419,490285.47158018046,27203.5512726292,107430.51508539419,23599.755402424704,39444.41313241093,229794.11735901757,296933.5844250128,31898.89445166201,122118.13198356815,29762.864363011704,109474.81275048848,14094.452562317174,28575.211472947598,258381.49418485272,62833.752321240754,5880.514053612631,17865.799193184277,52591.714641078666,84689.230158319,16271.296811890477,112636.51143423203,255948.96576958572,44669.93950410649,37950.64294699449,19974.229406768616,21823.08400408184,140851.5885728032,19944.889784149582,25817.170369697844,203076.49506414964,339758.77880560927,72714.35888680838,60974.69504154627,93132.39887869543,10798.470393349922,56085.04214822674,100675.8660019734,87808.86722918529,25205.19494818494,34555.946838829004,167611.1040471175,29847.01173730285,167789.37533616516,347269.9280158415,86370.61460762519,88349.83516636364,482809.45586904045,53388.55283910989,123150.88513890035,103630.47313789162,113765.58927913052,56886.46721094308,185990.77378900204,477929.216950017,180294.89658037858,182223.0946912322,340480.95478128234,43049.098297129065,253605.35026226818,11396.11168158752,10009.89674047982,19565.477074651004,92012.4297604441,291658.21154573106,47061.42413460863,217735.68359812972,92637.13470783619,155859.75871650188,89741.99335667417,40861.109896531,35915.33318731494,67053.63737486127,209509.6997880244,5738.662975894696,18593.396811852715,191089.49451916435,324830.7921852627,9128.70049245127,80645.70336010086,41290.324199144525,52014.57337985769,11019.169049149645,46570.7562549823,222484.91377096222,362433.67913489684,35776.52580414858,30510.86832457014,34183.21022702833,33716.12998381677,118352.75659560131]} diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/medium.json b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/medium.json new file mode 100644 index 000000000000..8d9345c8ddc7 --- /dev/null +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/medium.json @@ -0,0 +1 @@ +{"expected":[5.297932442873073,5.640605404492273,5.543771616854158,5.370819361586786,5.360382244826263,4.709628668970456,5.547902675556583,4.96695917727993,5.152682686632906,5.069992559903182,5.041338700302777,5.256902537127248,5.639958789397951,5.338899287998377,5.4862229661241555,4.548564108910101,5.610494915920906,5.37585177373591,3.7489247539088346,5.166621821454751,5.030803802220767,4.993260275532558,5.5075636285470075,5.225770815285207,5.420541970034563,4.683389894380692,5.152705617969719,5.543571583244232,5.244271784831994,5.159547869887345,4.889244274426916,5.371282405880172,4.889083757430596,3.7816562933309044,5.142966577155905,5.071688578701656,5.653474411603716,5.125991108549598,5.631783909958582,4.808676922482843,5.5997664125837145,5.02109432558016,5.332623074160816,5.66659573479357,5.274960414794218,4.998539279177382,5.476844399618951,5.236174497701167,5.271251924405752,5.175219041445591,5.673278655738039,5.466184153517774,5.437179214188343,5.586299185844098,5.573965212126386,5.13588300406301,5.308319809752897,5.239345701577919,5.360993766001567,5.1247520477297215,5.075356030764413,4.9222144882271115,5.187978875717892,5.666673749797058,5.399694066829916,4.877864115763167,5.310157404917509,3.4238260347132203,5.160331877328516,5.432343059633697,4.874065041740192,5.213556633802033,5.080970240937916,5.026646317175687,3.1539891438054304,5.176099024264579,4.966183062144444,5.413886688412867,5.31918007643341,5.519869365397794,5.586723471483822,4.531426386923642,4.956912637177755,4.401164646198089,5.589007587367674,5.6703810229628235,4.571150515579051,5.512747350922097,5.510488983067623,5.415089191659981,5.506833438082915,5.510722376065468,4.872584962579213,5.630662041700444,5.424258008636119,5.562628288509823,5.312529061804464,4.628025960242309,5.478764468722338,3.567045650989211,5.641238606725956,5.592866169797883,5.465217163527386,5.324507176232086,4.980492704353633,4.753950029216426,5.144882113717389,5.311383841107571,5.546066107512675,5.5354790820814825,4.314628937644934,5.370263659399351,5.548313556155101,4.866284330022998,5.433959764955782,5.607781661487631,5.380984581212988,5.588730519484812,5.439346207792007,5.40208328393351,5.051882724673507,4.640192991997508,5.251944165504522,5.631434531961809,5.194877276550785,5.294777570561486,4.968224505716788,5.551737621413624,4.8626523217802085,5.152514885302616,5.595587695054463,5.337449815148776,4.719106410104614,5.65520062374582,5.672557476799,4.685693743745963,5.195795740652329,5.397539989450434,4.944885376781741,5.384882685616704,5.458572436801836,5.469706625628674,4.946764117131257,5.6395733243617885,5.357634730639558,5.505011335686712,5.368636506110417,5.6354945163302785,5.342199159442905,5.37804300081126,5.272149560583296,5.557627154935549,5.4875781275601225,5.521964351366361,5.353588851229649,5.1285785273587186,5.198294523798764,4.82766341936424,5.220871005848583,5.414449781540382,5.37865959644923,4.196093594364468,5.5225095869472,5.065862667992173,5.566958085440327,5.509200831922085,5.666442571835966,5.272540721703526,4.784143653221107,5.189008281909669,5.557123422532671,5.151187008416531,5.428942543057168,5.306174284504015,5.017511282351347,5.237290424987225,5.613510218927985,4.534659380187328,5.6563413613799804,5.272177374212716,5.48744702197241,5.565463286012549,5.518992246618705,5.6493592572407465,5.1434801047470895,5.192239489303117,5.411423728163354,5.421727135386185,4.46121382447911,5.001444319778764,4.555181925515497,3.383135491854309,5.322623889696441,5.248664375550273,5.292012143221108,4.6275358207816835,5.663894485423631,5.639751146834365,4.22813099738241,5.418053825909536],"lambda":[2340.3576738539855,4644.17782611079,3826.514932854369,2707.614378217217,2651.6842635317357,721.7013189241054,3858.2596297881296,1207.3454544131355,1750.3718882915923,1483.590483849412,1400.9688920413707,2155.9900775913134,4638.175932473216,2540.1721539454124,3410.512122396125,522.9940747986507,4372.765320880077,2735.001957605865,105.80050248987388,1799.8510716640246,1371.7630379504308,1272.5455039231038,3559.221157717916,2025.8544765970753,2990.6995771303173,684.8133200902868,1750.452159219515,3824.9844425257734,2102.2128429489476,1774.5685809615331,1033.5686628797375,2710.1228759609353,1033.236959001079,112.94689438268276,1716.689837152039,1488.6308566262217,4765.25693458066,1659.3904572052531,4562.962138959306,879.7751194432647,4279.941861182236,1345.3850197502084,2508.4881908094103,4891.960763165716,2235.2723177020403,1286.050455049055,3347.1400173514467,2068.445264899669,2218.7559091752523,1831.0632044246024,4957.7826115484,3276.536227270597,3091.8820906231276,4166.207024821573,4064.6965362242377,1692.543053000416,2389.4830340735134,2081.604818205917,2654.929165703917,1655.2837895452553,1499.588736441516,1104.008519243578,1878.3883681788293,4892.724089377631,2868.570721685812,1010.3137570270702,2398.2803823295167,55.302417979297836,1777.3530519230903,3062.122219150958,1002.6675960539138,1976.9697589807865,1516.5197530844046,1360.4055400584864,32.30890689548534,1834.288356530179,1205.4730886826324,2951.1576468045496,2441.9479793388814,3647.9015338715617,4169.743707341614,505.3775738849921,1183.3315808177542,389.50654504628943,4188.834874722779,4929.13500124589,547.1532078577239,3596.31139137181,3580.1051833850406,2958.2633403576606,3554.027375846935,3581.7766386075514,999.7044228910427,4552.7359049877205,3013.00824544156,3973.5749499276876,2409.682408189455,613.0502710831641,3360.0175673942313,73.5889390639579,4650.062748110987,4221.284566340912,3270.2059175321233,2468.1022834001,1240.4664467136774,788.5804863596213,1723.2785755363611,2404.1698692120503,3844.114324147936,3763.578030449759,327.6317582470885,2704.606980569967,3861.431363875705,987.1880113352477,3072.038802232697,4349.101639343529,2763.2213198810546,4186.514426857288,3105.3106240514117,2882.3100014907895,1430.822749030851,628.1471534645827,2134.7169852186735,4559.774971952711,1904.4814057939307,2325.638153647375,1210.4042782032143,3887.9645559249525,980.0442506535879,1749.7846133403496,4244.322966821464,2532.819478773259,735.5088029384757,4781.736478658021,4950.637109598563,687.9752436726112,1907.9827099494562,2856.239774201832,1155.2106918358738,2784.846859590006,3227.036368891014,3299.6995641835597,1159.5589110743347,4634.601729584859,2637.1540579013067,3541.0999491634757,2695.820208015054,4596.949572640294,2556.990975640223,2747.0135485547603,2222.74245839733,3934.029958614029,3419.7677968428707,3663.217504777043,2615.902290394913,1667.9989249620323,1917.5410109648608,913.8187586221379,2006.1004577925635,2954.482884299023,2750.4030256107108,258.5156107338719,3667.2141347472298,1471.3881881170273,4008.1324626644914,3570.894050109291,4890.462509359505,2224.481909138184,837.6575174322203,1882.259257474847,3930.0687255688326,1745.1442236834937,3041.368411896223,2379.2523230830043,1335.7795274709035,2073.066516858599,4399.214412296707,508.65484299725046,4792.657965675557,2222.866097631534,3418.8712567682533,3996.1681459763295,3641.5081498288832,4726.199717757124,1718.4537067274423,1894.4615055156344,2936.657043805421,2997.7965367448296,439.1892892416472,1293.5432895607858,529.960020761461,50.993348173617086,2458.824108208546,2120.761054825679,2312.8118121324837,612.4497688054169,4865.60411187112,4636.250236102088,275.61108160281594,2975.8547901818274]} diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/runner.jl b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/runner.jl index ace6801dd0ad..234be06a5d97 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/runner.jl +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/fixtures/julia/runner.jl @@ -16,9 +16,43 @@ # See the License for the specific language governing permissions and # limitations under the License. -import Distributions: entropy, Poisson +import SpecialFunctions: loggamma import JSON +# Evaluate the entropy in extended precision, as the double-precision series used by `Distributions.entropy( Poisson( lambda ) )` is only accurate to tens of ULPs (and switches to an asymptotic expansion for `lambda >= 50`): +setprecision( BigFloat, 256 ); + +""" + poisson_entropy( lambda ) + +Compute the entropy of a Poisson distribution as `-sum( p_k * ln(p_k) )` in extended precision. + +# Arguments + +* `lambda`: mean parameter + +# Examples + +``` julia +julia> poisson_entropy( 9.0 ) +2.507673896080289 +``` +""" +function poisson_entropy( lambda ) + l = BigFloat( lambda ); + ll = log( l ); + + # Terms more than 40 standard deviations from the mean are smaller than `exp(-800)`: + lo = max( 0, floor( Int, lambda - 40.0*sqrt( lambda ) ) ); + hi = ceil( Int, lambda + 40.0*sqrt( lambda ) + 40.0 ); + h = zero( BigFloat ); + for k in lo:hi + lp = k*ll - l - loggamma( BigFloat( k+1 ) ); + h -= exp( lp ) * lp; + end + return Float64( h ); +end + """ gen( lambda, name ) @@ -39,7 +73,7 @@ julia> gen( lambda, \"data.json\" ); function gen( lambda, name ) z = Array{Float64}( undef, length(lambda) ); for i in eachindex(lambda) - z[ i ] = entropy( Poisson( lambda[ i ] ) ); + z[ i ] = poisson_entropy( lambda[ i ] ); end # Store data to be written to file as a collection: @@ -68,3 +102,9 @@ dir = dirname( file ); lambda = rand( 500 ) .* 20.0; gen( lambda, "data.json" ); +lambda = ( rand( 200 ) .* 4980.0 ) .+ 20.0; +gen( lambda, "medium.json" ); + +lambda = exp10.( ( rand( 100 ) .* 2.0 ) .+ 3.7 ); +gen( lambda, "large.json" ); + diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/test.js b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/test.js index 3b28ee1adbc6..d10e10000381 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/test.js +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/test.js @@ -30,6 +30,8 @@ var entropy = require( './../lib' ); // FIXTURES // var data = require( './fixtures/julia/data.json' ); +var medium = require( './fixtures/julia/medium.json' ); +var large = require( './fixtures/julia/large.json' ); // TESTS // @@ -80,7 +82,45 @@ tape( 'the function returns the entropy of a Poisson distribution', function tes if ( y === expected[i] ) { t.strictEqual( y, expected[i], 'lambda: '+lambda[i]+', y: '+y+', expected: '+expected[i] ); } else { - t.ok( isAlmostSameValue( y, expected[i], 20 ), 'within tolerance. lambda: '+lambda[i]+'. y: '+y+'. E: '+expected[ i ]+'.' ); + t.ok( isAlmostSameValue( y, expected[i], 4 ), 'within tolerance. lambda: '+lambda[i]+'. y: '+y+'. E: '+expected[ i ]+'.' ); + } + } + t.end(); +}); + +tape( 'the function returns the entropy of a Poisson distribution (medium `lambda`)', function test( t ) { + var expected; + var lambda; + var i; + var y; + + expected = medium.expected; + lambda = medium.lambda; + for ( i = 0; i < expected.length; i++ ) { + y = entropy( lambda[i] ); + if ( y === expected[i] ) { + t.strictEqual( y, expected[i], 'lambda: '+lambda[i]+', y: '+y+', expected: '+expected[i] ); + } else { + t.ok( isAlmostSameValue( y, expected[i], 4 ), 'within tolerance. lambda: '+lambda[i]+'. y: '+y+'. E: '+expected[ i ]+'.' ); + } + } + t.end(); +}); + +tape( 'the function returns the entropy of a Poisson distribution (large `lambda`)', function test( t ) { + var expected; + var lambda; + var i; + var y; + + expected = large.expected; + lambda = large.lambda; + for ( i = 0; i < expected.length; i++ ) { + y = entropy( lambda[i] ); + if ( y === expected[i] ) { + t.strictEqual( y, expected[i], 'lambda: '+lambda[i]+', y: '+y+', expected: '+expected[i] ); + } else { + t.ok( isAlmostSameValue( y, expected[i], 4 ), 'within tolerance. lambda: '+lambda[i]+'. y: '+y+'. E: '+expected[ i ]+'.' ); } } t.end(); diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/test.native.js b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/test.native.js index 602c2c02e2e7..17e7b6c8c219 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/test.native.js +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/test/test.native.js @@ -39,6 +39,8 @@ var opts = { // FIXTURES // var data = require( './fixtures/julia/data.json' ); +var medium = require( './fixtures/julia/medium.json' ); +var large = require( './fixtures/julia/large.json' ); // TESTS // @@ -89,7 +91,45 @@ tape( 'the function returns the entropy of a Poisson distribution', opts, functi if ( y === expected[i] ) { t.strictEqual( y, expected[i], 'lambda: '+lambda[i]+', y: '+y+', expected: '+expected[i] ); } else { - t.ok( isAlmostSameValue( y, expected[i], 20 ), 'within tolerance. lambda: '+lambda[i]+'. y: '+y+'. E: '+expected[ i ]+'.' ); + t.ok( isAlmostSameValue( y, expected[i], 4 ), 'within tolerance. lambda: '+lambda[i]+'. y: '+y+'. E: '+expected[ i ]+'.' ); + } + } + t.end(); +}); + +tape( 'the function returns the entropy of a Poisson distribution (medium `lambda`)', opts, function test( t ) { + var expected; + var lambda; + var y; + var i; + + expected = medium.expected; + lambda = medium.lambda; + for ( i = 0; i < expected.length; i++ ) { + y = entropy( lambda[i] ); + if ( y === expected[i] ) { + t.strictEqual( y, expected[i], 'lambda: '+lambda[i]+', y: '+y+', expected: '+expected[i] ); + } else { + t.ok( isAlmostSameValue( y, expected[i], 4 ), 'within tolerance. lambda: '+lambda[i]+'. y: '+y+'. E: '+expected[ i ]+'.' ); + } + } + t.end(); +}); + +tape( 'the function returns the entropy of a Poisson distribution (large `lambda`)', opts, function test( t ) { + var expected; + var lambda; + var y; + var i; + + expected = large.expected; + lambda = large.lambda; + for ( i = 0; i < expected.length; i++ ) { + y = entropy( lambda[i] ); + if ( y === expected[i] ) { + t.strictEqual( y, expected[i], 'lambda: '+lambda[i]+', y: '+y+', expected: '+expected[i] ); + } else { + t.ok( isAlmostSameValue( y, expected[i], 4 ), 'within tolerance. lambda: '+lambda[i]+'. y: '+y+'. E: '+expected[ i ]+'.' ); } } t.end(); From 48c6c29ee103724c6e608fd836482d89151db571 Mon Sep 17 00:00:00 2001 From: Abhist17 Date: Mon, 28 Sep 2026 09:33:42 +0530 Subject: [PATCH 2/2] docs: document the method in the JSDoc and C doc comments Move the explanation of the evaluation method out of the function body, as requested in review of #15611. --- .../base/dists/poisson/entropy/lib/main.js | 26 +++++++++++++------ .../base/dists/poisson/entropy/src/main.c | 26 +++++++++++++------ 2 files changed, 36 insertions(+), 16 deletions(-) diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/lib/main.js b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/lib/main.js index 5cebcadbffa5..241099a91aa5 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/lib/main.js +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/lib/main.js @@ -66,6 +66,24 @@ function lnratio( a, b ) { /** * Returns the entropy of a Poisson distribution. * +* ## Method +* +* For \\(\lambda < 5000\\), the function evaluates the definition \\(H = -\sum_k p_k \ln p_k\\) using weights relative to the mode \\(m = \lfloor \lambda \rfloor\\), \\(w_k = p_k / p_m\\). As the probabilities sum to one, \\(p_m = 1/W\\), where \\(W = \sum_k w_k\\), and thus +* +* ```tex +* H = \ln W - \frac{1}{W} \sum_k w_k \ln w_k +* ``` +* +* The logarithms \\(\ln w_k\\) are accumulated using the recurrences \\(w_{k+1} = w_k \lambda/(k+1)\\) and \\(w_{k-1} = w_k k/\lambda\\). This avoids evaluating \\(\lambda^k/k!\\), which overflows for large \\(k\\), and avoids the cancellation between \\(k \ln \lambda\\), \\(\lambda\\), and \\(\ln k!\\), each of which is much larger than the entropy. +* +* For \\(\lambda \geq 5000\\), the function uses the asymptotic expansion +* +* ```tex +* H = \frac{1}{2} \ln(2 \pi e \lambda) - \frac{1}{12\lambda} - \frac{1}{24\lambda^2} - \frac{19}{360\lambda^3} + O\left(\lambda^{-4}\right) +* ``` +* +* whose truncation error is less than \\(10^{-16}\\) relative to the result. +* * @param {NonNegativeNumber} lambda - mean parameter * @returns {PositiveNumber} entropy * @@ -100,17 +118,9 @@ function entropy( lambda ) { return 0.0; } if ( lambda >= THRESHOLD ) { - // Asymptotic expansion: H = ln(2πeλ)/2 - 1/(12λ) - 1/(24λ²) - 19/(360λ³) + O(1/λ⁴) l2 = lambda * lambda; return ( 0.5 * ( LN_TWO_PI + 1.0 + ln( lambda ) ) ) - ( 1.0 / ( 12.0*lambda ) ) - ( 1.0 / ( 24.0*l2 ) ) - ( 19.0 / ( 360.0*l2*lambda ) ); // eslint-disable-line max-len } - /* - * Evaluate H = -Σ p_k ln(p_k) using weights relative to the mode, w_k = p_k/p_m. As Σ p_k = 1, we have p_m = 1/W, where W = Σ w_k, and thus - * - * H = ln(W) - Σ w_k ln(w_k) / W - * - * Computing `ln(w_k)` by recurrence avoids evaluating `λ^k/k!`, which overflows for large `k`, and avoids the cancellation between `k*ln(λ)`, `λ`, and `ln(k!)`, each of which is much larger than the entropy. - */ m = floor( lambda ); wm1 = 0.0; // W - 1 sum = 0.0; // Σ w_k ln(w_k) diff --git a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/src/main.c b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/src/main.c index c02abc9a439a..677c646aca57 100644 --- a/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/src/main.c +++ b/lib/node_modules/@stdlib/stats/base/dists/poisson/entropy/src/main.c @@ -53,6 +53,24 @@ static double lnratio( const double a, const double b ) { /** * Returns the entropy of a Poisson distribution. * +* ## Method +* +* For `\lambda < 5000`, the function evaluates the definition `H = -\sum_k p_k \ln p_k` using weights relative to the mode `m = \lfloor \lambda \rfloor`, `w_k = p_k / p_m`. As the probabilities sum to one, `p_m = 1/W`, where `W = \sum_k w_k`, and thus +* +* ```tex +* H = \ln W - \frac{1}{W} \sum_k w_k \ln w_k +* ``` +* +* The logarithms `\ln w_k` are accumulated using the recurrences `w_{k+1} = w_k \lambda/(k+1)` and `w_{k-1} = w_k k/\lambda`. This avoids evaluating `\lambda^k/k!`, which overflows for large `k`, and avoids the cancellation between `k \ln \lambda`, `\lambda`, and `\ln k!`, each of which is much larger than the entropy. +* +* For `\lambda \geq 5000`, the function uses the asymptotic expansion +* +* ```tex +* H = \frac{1}{2} \ln(2 \pi e \lambda) - \frac{1}{12\lambda} - \frac{1}{24\lambda^2} - \frac{19}{360\lambda^3} + O\left(\lambda^{-4}\right) +* ``` +* +* whose truncation error is less than `10^{-16}` relative to the result. +* * @param lambda mean parameter * @return entropy * @@ -76,17 +94,9 @@ double stdlib_base_dists_poisson_entropy( const double lambda ) { return 0.0; } if ( lambda >= THRESHOLD ) { - // Asymptotic expansion: H = ln(2πeλ)/2 - 1/(12λ) - 1/(24λ²) - 19/(360λ³) + O(1/λ⁴) l2 = lambda * lambda; return ( 0.5 * ( STDLIB_CONSTANT_FLOAT64_LN_TWO_PI + 1.0 + stdlib_base_ln( lambda ) ) ) - ( 1.0 / ( 12.0*lambda ) ) - ( 1.0 / ( 24.0*l2 ) ) - ( 19.0 / ( 360.0*l2*lambda ) ); } - /* - * Evaluate H = -Σ p_k ln(p_k) using weights relative to the mode, w_k = p_k/p_m. As Σ p_k = 1, we have p_m = 1/W, where W = Σ w_k, and thus - * - * H = ln(W) - Σ w_k ln(w_k) / W - * - * Computing `ln(w_k)` by recurrence avoids evaluating `λ^k/k!`, which overflows for large `k`, and avoids the cancellation between `k*ln(λ)`, `λ`, and `ln(k!)`, each of which is much larger than the entropy. - */ m = stdlib_base_floor( lambda ); wm1 = 0.0; // W - 1 sum = 0.0; // Σ w_k ln(w_k)