Typo in ci.confidence_interval

[ethanclem]

Hello, I noticed a discrepancy in the the confidence intervals determined using allantools vs. those determined with matlab and checked against the NIST reference noted. It looks like there is a typo in ci.confidence_interval where the chi2_h and the chi2_l got swapped between var_l and var_h. The code snippet is below where the issue is in the last three lines of code

def confidence_interval(dev, edf, ci=one_sigma_ci):
    """ returns confidence interval (dev_min, dev_max) 
        for a given deviation dev, equivalent degrees of freedom edf,
        and degree of confidence ci.
        
    Parameters
    ----------
    dev: float
        Mean value (e.g. adev) around which we produce the confidence interval
    edf: float
        Equivalent degrees of freedon
    ci: float, defaults to scipy.special.erf(1/math.sqrt(2))
        for 1-sigma standard error set
        ci = scipy.special.erf(1/math.sqrt(2))
            = 0.68268949213708585
        
    Returns
    -------
    (dev_min, dev_max): (float, float)
        Confidence interval
    """
    ci_l = min(np.abs(ci), np.abs((ci-1))) / 2
    ci_h = 1 - ci_l

    # function from scipy, works OK, but scipy is large and slow to build
    chi2_l = scipy.stats.chi2.ppf(ci_l, edf)
    chi2_h = scipy.stats.chi2.ppf(ci_h, edf)

    variance = dev*dev
    var_l = float(edf) * variance / chi2_h  # NIST SP1065 eqn (45)
    var_h = float(edf) * variance / chi2_l
    return (np.sqrt(var_l), np.sqrt(var_h))

In the output, this error is causing the confidence intervals to be swapped where the upper confidence interval is smaller than the lower confidence interval, where I would expect that due to the chi squared distribution being asymmetric when the mean is near zero, the lower confidence interval should be the smaller one.

[aewallin]

the definition of ci_l with min() might be the problem here?
however the example code with one_sigma_ci seems to work: https://github.com/aewallin/allantools/blob/master/examples/ci_demo.py

could you please provide a code example that shows this problem?

[ethanclem]

My apologies, in the process of making the example code I realized it was an error in my understanding and the original was correct. Thank you for your discussion, I will close the comment!