Validation

validation

Validation and quality assurance module for HeavyTails library.

This module provides comprehensive mathematical validation, numerical accuracy testing, and property-based testing for all distributions.

GoodnessOfFitTests

GoodnessOfFitTests(alpha_level=0.05)

Statistical tests for distribution goodness-of-fit.

Both tests answer a question that AIC and BIC cannot. Information criteria rank candidate models against each other, so the best of a bad set still ranks first. A goodness-of-fit test asks whether the winner is compatible with the data at all.

For heavy tails the Anderson-Darling test is the more informative of the two, because it weights the tails of the distribution. The Kolmogorov-Smirnov statistic is driven by the centre, which is exactly where these families agree with each other.

Examples

from heavytails import Pareto data = Pareto(alpha=2.5, xm=1.0).rvs(500, seed=42) tests = GoodnessOfFitTests() result = tests.kolmogorov_smirnov_test( ... data, "pareto", alpha=2.5, xm=1.0 ... ) result["reject"] False

Parameters:

Name	Type	Description	Default
alpha_level	float	Significance level for the reject field.	0.05

Source code in heavytails/validation.py

def __init__(self, alpha_level: float = 0.05) -> None:
    """Initialise the test suite.

    Args:
        alpha_level: Significance level for the ``reject`` field.
    """
    if not (0.0 < alpha_level < 1.0):
        raise ValueError("alpha_level must be in (0,1).")
    self.alpha_level = alpha_level

anderson_darling_test

anderson_darling_test(
    data,
    distribution,
    *,
    parameters_estimated=False,
    **params,
)

Anderson-Darling test against a named distribution.

The statistic is A^2 = -n - (1/n) sum_i (2i-1) [ln F(x_i) + ln(1 - F(x_{n+1-i}))].

Unlike the Kolmogorov-Smirnov statistic this weights the tails, which is what makes it the more useful of the two here: heavy-tailed families differ from one another in the tail and agree in the middle.

Parameters:

Name	Type	Description	Default
data	list[float]	Sample values.	required
distribution	str	Distribution name, as accepted by the fitting helpers.	required
parameters_estimated	bool	Set when the parameters came from this same sample; see :meth:kolmogorov_smirnov_test.	False
**params	Any	Distribution parameters.	{}

Returns:

Type	Description
dict[str, Any]	Dictionary with the same shape as
dict[str, Any]	meth:kolmogorov_smirnov_test.

Raises:

Type	Description
ValueError	If the sample is empty or the distribution is unknown.

Source code in heavytails/validation.py

def anderson_darling_test(
    self,
    data: list[float],
    distribution: str,
    *,
    parameters_estimated: bool = False,
    **params: Any,
) -> dict[str, Any]:
    """Anderson-Darling test against a named distribution.

    The statistic is
    ``A^2 = -n - (1/n) sum_i (2i-1) [ln F(x_i) + ln(1 - F(x_{n+1-i}))]``.

    Unlike the Kolmogorov-Smirnov statistic this weights the tails, which
    is what makes it the more useful of the two here: heavy-tailed families
    differ from one another in the tail and agree in the middle.

    Args:
        data: Sample values.
        distribution: Distribution name, as accepted by the fitting helpers.
        parameters_estimated: Set when the parameters came from this same
            sample; see :meth:`kolmogorov_smirnov_test`.
        **params: Distribution parameters.

    Returns:
        Dictionary with the same shape as
        :meth:`kolmogorov_smirnov_test`.

    Raises:
        ValueError: If the sample is empty or the distribution is unknown.
    """
    values = sorted(float(x) for x in data)
    n = len(values)
    if n == 0:
        raise ValueError("data must not be empty.")

    dist = _resolve_distribution(distribution, params)

    # Clamp away from 0 and 1: the statistic takes logs of both F and 1-F,
    # and a single saturated value would send the whole sum to infinity.
    eps = 1e-15
    cdfs = [min(max(dist.cdf(x), eps), 1.0 - eps) for x in values]

    total = 0.0
    for i in range(1, n + 1):
        total += (2 * i - 1) * (math.log(cdfs[i - 1]) + math.log(1.0 - cdfs[n - i]))
    statistic = -n - total / n

    p_value = _anderson_darling_p_value(statistic)
    result: dict[str, Any] = {
        "test": "anderson-darling",
        "statistic": statistic,
        "p_value": p_value,
        "reject": p_value < self.alpha_level,
        "alpha_level": self.alpha_level,
        "n": n,
        "distribution": distribution,
        "parameters": dict(params),
        "method": "asymptotic",
    }
    if parameters_estimated:
        result["caveat"] = (
            "Parameters were estimated from this sample, so the asymptotic "
            "null distribution does not apply and the p-value is "
            "conservative: the test rejects less often than its nominal "
            "level. Use a parametric bootstrap for a calibrated p-value."
        )
    return result

kolmogorov_smirnov_test

kolmogorov_smirnov_test(
    data,
    distribution,
    *,
    parameters_estimated=False,
    **params,
)

Kolmogorov-Smirnov test against a named distribution.

The statistic is the largest vertical distance between the empirical distribution function and the fitted one, D = max_i max(i/n - F(x_i), F(x_i) - (i-1)/n).

Parameters:

Name	Type	Description	Default
data	list[float]	Sample values.	required
distribution	str	Distribution name, as accepted by the fitting helpers.	required
parameters_estimated	bool	Set when the parameters came from this same sample. The reported p-value is then conservative, because the fitted distribution is closer to the data than the null assumes, and the result carries a caveat.	False
**params	Any	Distribution parameters.	{}

Returns:

Type	Description
dict[str, Any]	Dictionary with statistic, p_value, reject, n,
dict[str, Any]	distribution, parameters and method, plus caveat
dict[str, Any]	when the p-value should not be read at face value.

Raises:

Type	Description
ValueError	If the sample is empty or the distribution is unknown.

Source code in heavytails/validation.py

def kolmogorov_smirnov_test(
    self,
    data: list[float],
    distribution: str,
    *,
    parameters_estimated: bool = False,
    **params: Any,
) -> dict[str, Any]:
    """Kolmogorov-Smirnov test against a named distribution.

    The statistic is the largest vertical distance between the empirical
    distribution function and the fitted one,
    ``D = max_i max(i/n - F(x_i), F(x_i) - (i-1)/n)``.

    Args:
        data: Sample values.
        distribution: Distribution name, as accepted by the fitting helpers.
        parameters_estimated: Set when the parameters came from this same
            sample. The reported p-value is then conservative, because the
            fitted distribution is closer to the data than the null
            assumes, and the result carries a ``caveat``.
        **params: Distribution parameters.

    Returns:
        Dictionary with ``statistic``, ``p_value``, ``reject``, ``n``,
        ``distribution``, ``parameters`` and ``method``, plus ``caveat``
        when the p-value should not be read at face value.

    Raises:
        ValueError: If the sample is empty or the distribution is unknown.
    """
    values = sorted(float(x) for x in data)
    n = len(values)
    if n == 0:
        raise ValueError("data must not be empty.")

    dist = _resolve_distribution(distribution, params)

    d_plus = 0.0
    d_minus = 0.0
    for i, x in enumerate(values, start=1):
        cdf = dist.cdf(x)
        d_plus = max(d_plus, i / n - cdf)
        d_minus = max(d_minus, cdf - (i - 1) / n)
    statistic = max(d_plus, d_minus)

    p_value = _kolmogorov_p_value(statistic, n)
    result: dict[str, Any] = {
        "test": "kolmogorov-smirnov",
        "statistic": statistic,
        "p_value": p_value,
        "reject": p_value < self.alpha_level,
        "alpha_level": self.alpha_level,
        "n": n,
        "distribution": distribution,
        "parameters": dict(params),
        "method": "asymptotic",
    }
    if parameters_estimated:
        result["caveat"] = (
            "Parameters were estimated from this sample, so the asymptotic "
            "null distribution does not apply and the p-value is "
            "conservative: the test rejects less often than its nominal "
            "level. Use a parametric bootstrap for a calibrated p-value."
        )
    return result

NumericalValidation

NumericalValidation(tolerance=1e-10)

Comprehensive numerical validation against scipy and known results.

Validates accuracy of PDF, CDF, PPF, and sampling for all distributions against scipy implementations where available.

Parameters:

Name	Type	Description	Default
tolerance	float	Maximum allowed relative error (default: 1e-10)	1e-10

Source code in heavytails/validation.py

def __init__(self, tolerance: float = 1e-10) -> None:
    """
    Initialize numerical validation.

    Args:
        tolerance: Maximum allowed relative error (default: 1e-10)
    """
    self.tolerance = tolerance
    self.test_results: dict[str, dict[str, Any]] = {}

validate_against_scipy

validate_against_scipy(distribution, params=None)

Compare distribution against SciPy implementation.

Parameters:

Name	Type	Description	Default
distribution	str	Distribution name	required
params	dict[str, float] | None	Optional specific parameters to test (uses defaults if None)	None

Returns:

Type	Description
dict[str, Any]	Dictionary with validation results including errors and pass/fail

Examples:

>>> validator = NumericalValidation()
>>> if SCIPY_AVAILABLE:
...     result = validator.validate_against_scipy("pareto", {"alpha": 2.5, "xm": 1.0})
...     result["pass"] or result["max_error"] < 0.01
... else:
...     True  # Skip if scipy not available
True

Source code in heavytails/validation.py

def validate_against_scipy(
    self, distribution: str, params: dict[str, float] | None = None
) -> dict[str, Any]:
    """
    Compare distribution against SciPy implementation.

    Args:
        distribution: Distribution name
        params: Optional specific parameters to test (uses defaults if None)

    Returns:
        Dictionary with validation results including errors and pass/fail

    Examples:
        >>> validator = NumericalValidation()
        >>> if SCIPY_AVAILABLE:
        ...     result = validator.validate_against_scipy("pareto", {"alpha": 2.5, "xm": 1.0})
        ...     result["pass"] or result["max_error"] < 0.01
        ... else:
        ...     True  # Skip if scipy not available
        True
    """
    if not SCIPY_AVAILABLE:
        return {
            "pass": False,
            "error": "scipy not available",
            "max_error": float("inf"),
        }

    dist_lower = distribution.lower()

    # Get test parameters
    if params is None:
        params = self._get_default_params(dist_lower)

    # Test cases for evaluation
    test_points = [0.1, 0.5, 1.0, 2.0, 5.0, 10.0]
    pdf_errors = []
    cdf_errors = []

    try:
        # Create our distribution
        our_dist = self._create_heavytails_distribution(dist_lower, params)

        # Create scipy equivalent
        scipy_dist = self._create_scipy_distribution(dist_lower, params)

        if scipy_dist is None:
            return {
                "pass": False,
                "error": f"No scipy equivalent for {distribution}",
            }

        # Test PDF at multiple points
        for x in test_points:
            if x > 0:  # Most distributions require x > 0
                try:
                    our_pdf = our_dist.pdf(x)
                    scipy_pdf = scipy_dist.pdf(x)

                    if scipy_pdf > 1e-10 and math.isfinite(scipy_pdf):
                        rel_error = abs(our_pdf - scipy_pdf) / scipy_pdf
                        pdf_errors.append(rel_error)

                    # Test CDF
                    our_cdf = our_dist.cdf(x)
                    scipy_cdf = scipy_dist.cdf(x)

                    if scipy_cdf > 1e-10 and scipy_cdf < 1 - 1e-10:
                        rel_error_cdf = abs(our_cdf - scipy_cdf) / max(
                            scipy_cdf, 1 - scipy_cdf
                        )
                        cdf_errors.append(rel_error_cdf)

                except (ValueError, OverflowError, ZeroDivisionError):
                    continue

        if not pdf_errors and not cdf_errors:
            return {
                "pass": False,
                "error": "No valid comparison points",
                "max_error": float("inf"),
            }

        max_pdf_error = max(pdf_errors) if pdf_errors else 0.0
        max_cdf_error = max(cdf_errors) if cdf_errors else 0.0
        max_error = max(max_pdf_error, max_cdf_error)

        return {
            "pass": max_error < self.tolerance,
            "max_error": float(max_error),
            "max_pdf_error": float(max_pdf_error),
            "max_cdf_error": float(max_cdf_error),
            "mean_pdf_error": (
                float(sum(pdf_errors) / len(pdf_errors)) if pdf_errors else 0.0
            ),
            "mean_cdf_error": (
                float(sum(cdf_errors) / len(cdf_errors)) if cdf_errors else 0.0
            ),
            "num_pdf_tests": len(pdf_errors),
            "num_cdf_tests": len(cdf_errors),
            "distribution": distribution,
            "parameters": params,
        }

    except Exception as e:
        return {
            "pass": False,
            "error": str(e),
            "max_error": float("inf"),
        }

PropertyBasedTests

PropertyBasedTests()

Property-based testing for mathematical correctness using Hypothesis.

Tests fundamental mathematical properties that all distributions should satisfy: - PDF non-negativity - CDF monotonicity - PPF/CDF inverse relationship - Probability axioms

Source code in heavytails/validation.py

def __init__(self) -> None:
    """Initialize property-based tester."""
    self.test_results: dict[str, bool] = {}

test_cdf_monotonicity

test_cdf_monotonicity(distribution)

Test that CDF is monotonically increasing.

Parameters:

Name	Type	Description	Default
distribution	str	Distribution name to test	required

Returns:

Type	Description
dict[str, Any]	Dictionary with test results

Source code in heavytails/validation.py

def test_cdf_monotonicity(self, distribution: str) -> dict[str, Any]:
    """
    Test that CDF is monotonically increasing.

    Args:
        distribution: Distribution name to test

    Returns:
        Dictionary with test results
    """
    if not HYPOTHESIS_AVAILABLE:
        return {
            "pass": False,
            "error": "Hypothesis not available",
            "property": "cdf_monotonicity",
        }

    violations = []

    try:
        dist_lower = distribution.lower()
        test_cases = self._generate_test_cases(dist_lower, n_cases=50)

        for params, x_values in test_cases:
            try:
                # Create distribution
                try:
                    dist = create(dist_lower, **params)
                except ValueError:
                    continue

                # Test monotonicity: CDF(x1) <= CDF(x2) for x1 < x2
                sorted_x = sorted([x for x in x_values if x > 0])
                for i in range(len(sorted_x) - 1):
                    x1, x2 = sorted_x[i], sorted_x[i + 1]
                    cdf1, cdf2 = dist.cdf(x1), dist.cdf(x2)

                    if cdf1 > cdf2 + 1e-10:  # Allow small numerical error
                        violations.append(
                            {
                                "params": params,
                                "x1": x1,
                                "x2": x2,
                                "cdf1": cdf1,
                                "cdf2": cdf2,
                            }
                        )

            except Exception:
                continue

        return {
            "pass": len(violations) == 0,
            "property": "cdf_monotonicity",
            "distribution": distribution,
            "num_tests": sum(len(x_values) - 1 for _, x_values in test_cases),
            "violations": violations[:5],
            "num_violations": len(violations),
        }

    except Exception as e:
        return {
            "pass": False,
            "error": str(e),
            "property": "cdf_monotonicity",
        }

test_pdf_nonnegativity

test_pdf_nonnegativity(distribution)

Test that PDF is non-negative for all valid inputs.

Parameters:

Name	Type	Description	Default
distribution	str	Distribution name to test	required

Returns:

Type	Description
dict[str, Any]	Dictionary with test results

Examples:

>>> tester = PropertyBasedTests()
>>> result = tester.test_pdf_nonnegativity("pareto")
>>> result["property"]
'pdf_nonnegativity'

Source code in heavytails/validation.py

def test_pdf_nonnegativity(self, distribution: str) -> dict[str, Any]:
    """
    Test that PDF is non-negative for all valid inputs.

    Args:
        distribution: Distribution name to test

    Returns:
        Dictionary with test results

    Examples:
        >>> tester = PropertyBasedTests()
        >>> result = tester.test_pdf_nonnegativity("pareto")
        >>> result["property"]
        'pdf_nonnegativity'
    """
    if not HYPOTHESIS_AVAILABLE:
        return {
            "pass": False,
            "error": "Hypothesis not available",
            "property": "pdf_nonnegativity",
        }

    violations = []

    try:
        dist_lower = distribution.lower()

        # Generate test cases
        test_cases = self._generate_test_cases(dist_lower, n_cases=50)

        for params, x_values in test_cases:
            try:
                # Create distribution
                try:
                    dist = create(dist_lower, **params)
                except ValueError:
                    continue

                # Test PDF non-negativity
                for x in x_values:
                    if x > 0:  # Most distributions require x > 0
                        pdf_val = dist.pdf(x)
                        if pdf_val < 0 or math.isnan(pdf_val):
                            violations.append(
                                {
                                    "params": params,
                                    "x": x,
                                    "pdf": pdf_val,
                                }
                            )

            except Exception:
                continue

        return {
            "pass": len(violations) == 0,
            "property": "pdf_nonnegativity",
            "distribution": distribution,
            "num_tests": len(test_cases) * len(x_values) if test_cases else 0,
            "violations": violations[:5],  # Return first 5 violations
            "num_violations": len(violations),
        }

    except Exception as e:
        return {
            "pass": False,
            "error": str(e),
            "property": "pdf_nonnegativity",
        }

test_ppf_cdf_inverse

test_ppf_cdf_inverse(distribution)

Test that PPF and CDF are inverse functions: CDF(PPF(u)) ≈ u.

Parameters:

Name	Type	Description	Default
distribution	str	Distribution name to test	required

Returns:

Type	Description
dict[str, Any]	Dictionary with test results

Source code in heavytails/validation.py

def test_ppf_cdf_inverse(self, distribution: str) -> dict[str, Any]:
    """
    Test that PPF and CDF are inverse functions: CDF(PPF(u)) ≈ u.

    Args:
        distribution: Distribution name to test

    Returns:
        Dictionary with test results
    """

    violations = []

    try:
        dist_lower = distribution.lower()
        test_cases = self._generate_test_cases(dist_lower, n_cases=20)

        # Test quantiles
        u_values = [0.01, 0.05, 0.1, 0.25, 0.5, 0.75, 0.9, 0.95, 0.99]

        for params, _ in test_cases:
            try:
                # Create distribution
                try:
                    dist = create(dist_lower, **params)
                except ValueError:
                    continue

                # Test CDF(PPF(u)) ≈ u
                for u in u_values:
                    x = dist.ppf(u)
                    u_recovered = dist.cdf(x)

                    error = abs(u - u_recovered)
                    if error > 1e-6:  # Tolerance for numerical error
                        violations.append(
                            {
                                "params": params,
                                "u": u,
                                "x": x,
                                "u_recovered": u_recovered,
                                "error": error,
                            }
                        )

            except Exception:
                continue

        return {
            "pass": len(violations) == 0,
            "property": "ppf_cdf_inverse",
            "distribution": distribution,
            "num_tests": len(test_cases) * len(u_values),
            "violations": violations[:5],
            "num_violations": len(violations),
            "max_error": max((v["error"] for v in violations), default=0.0),
        }

    except Exception as e:
        return {
            "pass": False,
            "error": str(e),
            "property": "ppf_cdf_inverse",
        }

convergence_validation

convergence_validation(distribution, method='ppf')

Validate convergence of numerical algorithms.

Tests convergence properties of iterative algorithms used in the library, such as PPF computation via bisection/Newton-Raphson.

Parameters:

Name	Type	Description	Default
distribution	str	Distribution name to test	required
method	str	Method to test ("ppf", "cdf", or "pdf")	'ppf'

Returns:

Type	Description
dict[str, Any]	Dictionary with convergence diagnostics

Examples:

>>> result = convergence_validation("pareto", "ppf")
>>> "converged" in result
True

Source code in heavytails/validation.py

def convergence_validation(distribution: str, method: str = "ppf") -> dict[str, Any]:
    """
    Validate convergence of numerical algorithms.

    Tests convergence properties of iterative algorithms used in the library,
    such as PPF computation via bisection/Newton-Raphson.

    Args:
        distribution: Distribution name to test
        method: Method to test ("ppf", "cdf", or "pdf")

    Returns:
        Dictionary with convergence diagnostics

    Examples:
        >>> result = convergence_validation("pareto", "ppf")
        >>> "converged" in result
        True
    """

    try:
        dist_lower = distribution.lower()

        # A representative instance of the family, for a check that does not
        # depend on the particular parameters.
        if dist_lower not in _REPRESENTATIVE:
            return {
                "converged": False,
                "error": f"Unknown distribution: {distribution}",
            }
        dist = create(dist_lower, **_REPRESENTATIVE[dist_lower])

        if method == "ppf":
            # Test PPF convergence for various quantiles
            u_values = [0.01, 0.1, 0.5, 0.9, 0.99]
            convergence_info = []

            for u in u_values:
                try:
                    # Compute PPF
                    x = dist.ppf(u)

                    # Verify convergence: CDF(PPF(u)) should equal u
                    u_recovered = dist.cdf(x)
                    error = abs(u - u_recovered)

                    convergence_info.append(
                        {
                            "u": u,
                            "x": x,
                            "error": error,
                            "converged": error < 1e-6,
                        }
                    )

                except Exception as e:
                    convergence_info.append(
                        {
                            "u": u,
                            "error": str(e),
                            "converged": False,
                        }
                    )

            all_converged = all(
                info.get("converged", False) for info in convergence_info
            )
            max_error = max(
                (
                    info["error"]
                    for info in convergence_info
                    if isinstance(info["error"], (int, float))
                ),
                default=float("inf"),
            )

            return {
                "converged": all_converged,
                "method": method,
                "distribution": distribution,
                "convergence_info": convergence_info,
                "max_error": float(max_error),
                "num_tests": len(u_values),
            }

        else:
            return {
                "converged": False,
                "error": f"Method {method} not implemented",
            }

    except Exception as e:
        return {
            "converged": False,
            "error": str(e),
        }

parameter_stability_check

parameter_stability_check(distribution, **params)

Check parameter combinations for numerical stability with automatic fixes.

Analyzes parameters for potential numerical issues and provides specific warnings and suggested fixes.

Parameters:

Name	Type	Description	Default
distribution	str	Distribution name	required
**params	Any	Distribution parameters to check	{}

Returns:

Type	Description
dict[str, Any]	Dictionary with warnings, suggested fixes, and stability assessment

Examples:

>>> result = parameter_stability_check("pareto", alpha=1e-8, xm=1.0)
>>> len(result["warnings"]) > 0
True
>>> result["stable"]
False

Source code in heavytails/validation.py

def parameter_stability_check(distribution: str, **params: Any) -> dict[str, Any]:
    """
    Check parameter combinations for numerical stability with automatic fixes.

    Analyzes parameters for potential numerical issues and provides
    specific warnings and suggested fixes.

    Args:
        distribution: Distribution name
        **params: Distribution parameters to check

    Returns:
        Dictionary with warnings, suggested fixes, and stability assessment

    Examples:
        >>> result = parameter_stability_check("pareto", alpha=1e-8, xm=1.0)
        >>> len(result["warnings"]) > 0
        True
        >>> result["stable"]
        False
    """
    warnings_list = []
    fixes = []
    severity = "low"

    dist_lower = distribution.lower()

    if dist_lower == "pareto":
        alpha = params.get("alpha", 1.0)
        xm = params.get("xm", 1.0)

        if alpha < 1e-6:
            warnings_list.append("Alpha too small (< 1e-6), may cause overflow in PDF")
            fixes.append("Use alpha >= 1e-6")
            severity = "high"

        if alpha > 1e6:
            warnings_list.append("Alpha too large (> 1e6), may cause underflow in tail")
            fixes.append("Use alpha <= 1e6")
            severity = "medium"

        # Test numerical stability
        try:
            test_x = xm * 2.0
            pdf_val = (alpha * (xm**alpha)) / (test_x ** (alpha + 1))
            if math.isnan(pdf_val) or math.isinf(pdf_val):
                warnings_list.append("PDF computation unstable with these parameters")
                severity = "high"
        except (OverflowError, ZeroDivisionError):
            warnings_list.append("PDF computation failed with these parameters")
            severity = "high"

    elif dist_lower == "lognormal":
        mu = params.get("mu", 0.0)
        sigma = params.get("sigma", 1.0)

        if abs(mu) > 100:
            warnings_list.append("Very large |mu| (> 100) may cause numerical overflow")
            fixes.append("Use |mu| <= 100")
            severity = "medium"

        if sigma > 10:
            warnings_list.append("Very large sigma (> 10) may cause numerical issues")
            fixes.append("Use sigma <= 10")
            severity = "medium"

        if sigma < 1e-6:
            warnings_list.append(
                "Very small sigma (< 1e-6) approaches degenerate distribution"
            )
            fixes.append("Use sigma >= 1e-6")

    elif dist_lower == "cauchy":
        gamma = params.get("gamma", 1.0)

        if gamma < 1e-6:
            warnings_list.append("Very small gamma (< 1e-6) may cause numerical issues")
            fixes.append("Use gamma >= 1e-6")

        if gamma > 1e6:
            warnings_list.append("Very large gamma (> 1e6) may cause numerical issues")
            fixes.append("Use gamma <= 1e6")

    elif dist_lower == "studentt":
        nu = params.get("nu", 5.0)

        if nu < 1e-6:
            warnings_list.append("Nu too small (< 1e-6), Student-t undefined")
            fixes.append("Use nu >= 0.1")
            severity = "high"

        if nu > 1e6:
            warnings_list.append(
                "Very large nu (> 1e6): consider using normal distribution instead"
            )
            fixes.append("For nu > 30, Normal approximation often sufficient")

    elif dist_lower == "weibull":
        k = params.get("k", 1.0)
        lam = params.get("lam", 1.0)

        if k < 1e-6 or lam < 1e-6:
            warnings_list.append("Very small shape/scale parameters may cause issues")
            fixes.append("Use k, lam >= 1e-6")

        if k > 100 or lam > 1e6:
            warnings_list.append("Very large shape/scale parameters may cause overflow")
            fixes.append("Consider rescaling parameters")

    elif dist_lower == "frechet":
        alpha = params.get("alpha", 2.0)
        s = params.get("s", 1.0)

        if alpha < 1e-6 or s < 1e-6:
            warnings_list.append("Very small parameters may cause numerical issues")
            fixes.append("Use alpha, s >= 1e-6")

    # General checks for all distributions
    for param_name, param_value in params.items():
        if not math.isfinite(param_value):
            warnings_list.append(f"Parameter {param_name} is not finite")
            severity = "high"

    return {
        "warnings": warnings_list,
        "suggested_fixes": fixes,
        "stable": len(warnings_list) == 0,
        "severity": severity,
        "distribution": distribution,
        "parameters": params,
    }

ppf_edge_case_handler

ppf_edge_case_handler(distribution, u, **params)

Handle edge cases in quantile function calculation.

Problematic cases: - u very close to 0 or 1 - Parameters at boundary values - Distributions with bounded support - Numerical overflow/underflow

Should provide graceful degradation and informative errors.

Source code in heavytails/validation.py

def ppf_edge_case_handler(distribution: str, u: float, **params: Any) -> float:
    """
    Handle edge cases in quantile function calculation.

    Problematic cases:
    - u very close to 0 or 1
    - Parameters at boundary values
    - Distributions with bounded support
    - Numerical overflow/underflow

    Should provide graceful degradation and informative errors.
    """
    _ = (distribution, params)  # Reserved for future implementation
    if not (0 < u < 1):
        if u == 0:
            # TODO: Return theoretical minimum (support lower bound)
            pass
        elif u == 1:
            # TODO: Return theoretical maximum (support upper bound)
            pass
        else:
            raise ValueError(f"u must be in (0,1), got {u}")

    # TODO: Implement robust edge case handling for all distributions
    raise NotImplementedError("PPF edge case handling not fully implemented")