feat: add stats/base/dists/bradford/skewness

[vivekmaurya001]

Progresses #168.

This pull request:

• add stats/base/dists/bradford/skewness

Related Issues

Does this pull request have any related issues?

This pull request:

• resolves #{{TODO: add issue number}}

Questions

Any questions for reviewers of this pull request?

No.

Other

Any other information relevant to this pull request? This may include screenshots, references, and/or implementation notes.

No.

Checklist

Please ensure the following tasks are completed before submitting this pull request.

• Read, understood, and followed the contributing guidelines

---

@stdlib-js/reviewers

[vivekmaurya001]

/stdlib update-copyright-years

[vivekmaurya001]

currently there are 10 test cases which fail at tolerance factor of 50, but these cases persists even I increase the factor to 1e4

[vivekmaurya001]

@anandkaranubc , @kgryte , Any suggestion on this

[anandkaranubc]

@vivekmaurya001, that's the main reason I didn’t submit my PRs on skewness and kurtosis — they lead to large divergent results due to ("My Hypothesis"):

1. A lot of calculations being performed, which causes floating-point precision errors.
2. The shape parameter being very small, leading to differences in the ln calculations compared to the reference (SciPy) implementation.
3. The expected values being very small.

These are my local results:

• For small c: tol = 930226.0 * EPS * abs(expected[i])
• For large c: tol = 600.0 * EPS * abs(expected[i])

There are actually two ways to implement this:

k = ln(1.0 + c);
out = SQRT2;
out *= ( (12.0*c*c) - (9.0*k*c*(c+2.0)) + (2.0*k*k*(c*(c+3.0)+3.0)) );
out /= sqrt(c * (c*(k-2.0) + (2.0*k)));
out /= ( (3.0*c*(k-2.0)) + (6.0*k) );
return out;

or

var g1 = SQRT2 * (12*c*c - 9*c*k*(c+2) + 2*k*k*(c*(c+3)+3));
g1 /= sqrt(c*(c*(k-2)+2*k)) * (3*c*(k-2) + 6*k);
return g1;

The second one gives better results than the first. Here, better results mean less divergence from SciPy.

That said, I think we might need to check if there have been any changes in the upstream implementations of ln or log1p, and update ours if needed.

cc: @kgryte

[vivekmaurya001]

Yess, Also the scipy docs are using log(1+c) instead of log1p

[anandkaranubc]

Yess, Also the scipy docs are using log(1+c) instead of log1p

Would you mind checking the upstream implementation of ln to see if there's any difference? If there is, create an issue and we can sync up on that.

FreeBSD implementation
Our implementation

[vivekmaurya001]

Sure I will try it first

[vivekmaurya001]

@anandkaranubc , I found some slight differences in lg2, lg3 etc values in free bsd and our, I tried correcting it and running but still not improved, apart from this all looked same

[anandkaranubc]

@anandkaranubc , I found some slight differences in lg2, lg3 etc values in free bsd and our, I tried correcting it and running but still not improved, apart from this all looked same

Could you create a separate issue for this so we can tackle it there and get some feedback? It might also be worthwhile to check how SciPy has implemented ln.

Edit: Are you sure about that?

The values seem the same to me:

stdlib/lib/node_modules/@stdlib/math/base/special/ln/scripts/evalpoly.js

Lines 40 to 51 in 14a166d

	// Polynomial coefficients ordered in ascending degree...
	var P = [
	3.999999999940941908e-01, // 3FD99999 9997FA04
	2.222219843214978396e-01, // 3FCC71C5 1D8E78AF
	1.531383769920937332e-01 // 3FC39A09 D078C69F
	];
	var Q = [
	6.666666666666735130e-01, // 3FE55555 55555593
	2.857142874366239149e-01, // 3FD24924 94229359
	1.818357216161805012e-01, // 3FC74664 96CB03DE
	1.479819860511658591e-01 // 3FC2F112 DF3E5244
	];

[vivekmaurya001]

Oops yes, I realised some last digits were not matching

[vivekmaurya001]

Added ln logic instead of log1p and TODO comment in test file