| 324 | |
| 325 | @derived_from(scipy.stats) |
| 326 | def kurtosistest(a, axis=0, nan_policy="propagate"): |
| 327 | if nan_policy != "propagate": |
| 328 | raise NotImplementedError( |
| 329 | "`nan_policy` other than 'propagate' have not been implemented." |
| 330 | ) |
| 331 | |
| 332 | n = float(a.shape[axis]) |
| 333 | b2 = kurtosis(a, axis, fisher=False) |
| 334 | |
| 335 | E = 3.0 * (n - 1) / (n + 1) |
| 336 | varb2 = ( |
| 337 | 24.0 * n * (n - 2) * (n - 3) / ((n + 1) * (n + 1.0) * (n + 3) * (n + 5)) |
| 338 | ) # [1]_ Eq. 1 |
| 339 | x = (b2 - E) / np.sqrt(varb2) # [1]_ Eq. 4 |
| 340 | # [1]_ Eq. 2: |
| 341 | sqrtbeta1 = ( |
| 342 | 6.0 |
| 343 | * (n * n - 5 * n + 2) |
| 344 | / ((n + 7) * (n + 9)) |
| 345 | * np.sqrt((6.0 * (n + 3) * (n + 5)) / (n * (n - 2) * (n - 3))) |
| 346 | ) |
| 347 | # [1]_ Eq. 3: |
| 348 | A = 6.0 + 8.0 / sqrtbeta1 * (2.0 / sqrtbeta1 + np.sqrt(1 + 4.0 / (sqrtbeta1**2))) |
| 349 | term1 = 1 - 2 / (9.0 * A) |
| 350 | denom = 1 + x * np.sqrt(2 / (A - 4.0)) |
| 351 | denom = np.where(denom < 0, 99, denom) |
| 352 | term2 = np.where(denom < 0, term1, np.power((1 - 2.0 / A) / denom, 1 / 3.0)) |
| 353 | Z = (term1 - term2) / np.sqrt(2 / (9.0 * A)) # [1]_ Eq. 5 |
| 354 | Z = np.where(denom == 99, 0, Z) |
| 355 | if Z.ndim == 0: |
| 356 | Z = Z[()] |
| 357 | |
| 358 | # zprob uses upper tail, so Z needs to be positive |
| 359 | return delayed(KurtosistestResult, nout=2)(Z, 2 * distributions.norm.sf(np.abs(Z))) |
| 360 | |
| 361 | |
| 362 | @derived_from(scipy.stats) |