Repository: scipy/scipy
Pull request: BUG: stats.tukeylambda: Fix the support() method.
Closes gh-21369.


Problem statement:
BUG: stats.tukeylambda: `support(lam)` is incorrect when `lam > 0`
### Describe your issue.

When λ > 0, the support of the Tukey lambda distribution is (-1/λ, 1/λ).  Currently, `tukeylambda.support(lam)` returns (-inf, inf) even when `lam` > 0.

### Reproducing Code Example

```python
In [1]: from scipy.stats import tukeylambda

In [2]: tukeylambda.support(0.5)  # Should be (-2, 2)
Out[2]: (np.float64(-inf), np.float64(inf))
```


### Error message

There is an error related to this problem.  The `pdf` method will generate a spurious warning when evaluated at the bounds of the support when `lam` > 0, e.g.

```
In [3]: tukeylambda.pdf(2, 0.5)
/home/warren/repos/git/forks/scipy/build-install/lib/python3.12/site-packages/scipy/stats/_continuous_distns.py:10492: RuntimeWarning: divide by zero encountered in power
  Px = Fx**(lam-1.0) + (np.asarray(1-Fx))**(lam-1.0)
Out[3]: np.float64(0.0)
```

<details>

<summary>
SciPy/NumPy/Python version and system information
</summary>

```shell
1.15.0.dev0+1377.235e778 2.0.0 sys.version_info(major=3, minor=12, micro=4, releaselevel='final', serial=0)
Build Dependencies:
  blas:
    detection method: pkgconfig
    found: true
    include directory: /usr/include/x86_64-linux-gnu/openblas-pthread/
    lib directory: /usr/lib/x86_64-linux-gnu/openblas-pthread/
    name: openblas
    openblas configuration: USE_64BITINT= DYNAMIC_ARCH=1 DYNAMIC_OLDER=1 NO_CBLAS=
      NO_LAPACK= NO_LAPACKE=1 NO_AFFINITY=1 USE_OPENMP=0 generic MAX_THREADS=64
    pc file directory: /usr/lib/x86_64-linux-gnu/pkgconfig
    version: 0.3.20
  lapack:
    detection method: pkgconfig
    found: true
    include directory: /usr/include/x86_64-linux-gnu/openblas-pthread/
    lib directory: /usr/lib/x86_64-linux-gnu/openblas-pthread/
    name: openblas
    openblas configuration: USE_64BITINT= DYNAMIC_ARCH=1 DYNAMIC_OLDER=1 NO_CBLAS=
      NO_LAPACK= NO_LAPACKE=1 NO_AFFINITY=1 USE_OPENMP=0 generic MAX_THREADS=64
    pc file directory: /usr/lib/x86_64-linux-gnu/pkgconfig
    version: 0.3.20
  pybind11:
    detection method: config-tool
    include directory: unknown
    name: pybind11
    version: 2.13.1
Compilers:
  c:
    commands: cc
    linker: ld.bfd
    name: gcc
    version: 11.4.0
  c++:
    commands: c++
    linker: ld.bfd
    name: gcc
    version: 11.4.0
  cython:
    commands: cython
    linker: cython
    name: cython
    version: 3.0.10
  fortran:
    commands: gfortran
    linker: ld.bfd
    name: gcc
    version: 11.4.0
  pythran:
    include directory: ../../../../../py3.12.4/lib/python3.12/site-packages/pythran
    version: 0.16.1
Machine Information:
  build:
    cpu: x86_64
    endian: little
    family: x86_64
    system: linux
  cross-compiled: false
  host:
    cpu: x86_64
    endian: little
    family: x86_64
    system: linux
Python Information:
  path: /home/warren/py3.12.4/bin/python3
  version: '3.12'
```

</details>

Hints from issue discussion:
I'll submit a PR to fix this soon.

Commit message: BUG: stats.tukeylambda: Fix the support() method.

Closes gh-21369.

Patch to review:
diff --git a/scipy/stats/_continuous_distns.py b/scipy/stats/_continuous_distns.py
index f8a478bab47b..84fd0aa1191c 100644
--- a/scipy/stats/_continuous_distns.py
+++ b/scipy/stats/_continuous_distns.py
@@ -10481,12 +10481,20 @@ class tukeylambda_gen(rv_continuous):
     %(example)s
 
     """
+    _support_mask = rv_continuous._open_support_mask
+
     def _argcheck(self, lam):
         return np.isfinite(lam)
 
     def _shape_info(self):
         return [_ShapeInfo("lam", False, (-np.inf, np.inf), (False, False))]
 
+    def _get_support(self, lam):
+        b = _lazywhere(lam > 0, (lam,),
+                       f=lambda lam: 1/lam,
+                       f2=lambda lam: np.inf)
+        return -b, b
+
     def _pdf(self, x, lam):
         Fx = np.asarray(sc.tklmbda(x, lam))
         Px = Fx**(lam-1.0) + (np.asarray(1-Fx))**(lam-1.0)
diff --git a/scipy/stats/tests/test_distributions.py b/scipy/stats/tests/test_distributions.py
index 61235b81b27c..d029345cdf5f 100644
--- a/scipy/stats/tests/test_distributions.py
+++ b/scipy/stats/tests/test_distributions.py
@@ -7870,6 +7870,22 @@ def test_regression_tukey_lambda():
     assert_((p[2] == 0.0).any())
 
 
+def test_tukeylambda_support():
+    lam = np.array([-1.75, -0.5, 0.0, 0.25, 0.5, 2.0])
+    a, b = stats.tukeylambda.support(lam)
+    expected_b = np.array([np.inf, np.inf, np.inf, 4, 2, 0.5])
+    assert_equal(b, expected_b)
+    assert_equal(a, -expected_b)
+
+
+def test_tukeylambda_pdf_support_boundary():
+    # Verify that tukeylambda.pdf() doesn't generate a
+    # warning when evaluated at the bounds of the support.
+    # For lam=0.5, the support is (-2, 2).
+    p = stats.tukeylambda.pdf([-2.0, 2.0], 0.5)
+    assert_equal(p, [0.0, 0.0])
+
+
 @pytest.mark.skipif(DOCSTRINGS_STRIPPED, reason="docstrings stripped")
 def test_regression_ticket_1421():
     assert_('pdf(x, mu, loc=0, scale=1)' not in stats.poisson.__doc__)


Code review comments:
diff_hunk: @@ -10481,12 +10481,20 @@ class tukeylambda_gen(rv_continuous):
     %(example)s
 
     """
+    _support_mask = rv_continuous._open_support_mask
+
     def _argcheck(self, lam):
         return np.isfinite(lam)
 
     def _shape_info(self):
         return [_ShapeInfo("lam", False, (-np.inf, np.inf), (False, False))]
 
+    def _get_support(self, lam):
+        b = _lazywhere(lam > 0, (lam,),
+                       f=lambda lam: 1/lam,
+                       f2=lambda lam: np.inf)
path: scipy/stats/_continuous_distns.py
text: @user1:
Just to show that I'm actually thinking about this - not that there's any real motivation to change it - it would be *slightly* simpler as:
```suggestion
        b = _lazywhere(lam <= 0, (lam,), np.inf, f2=1/lam)
```
right?
And `_argcheck` should already be taking care of when `lam` is invalid, but I suppose only this way would do the right thing if `lam = np.nan` (because `cond` would fail and it would have to calculate `1/np.nan` rather than returning `np.inf`.).
Again, zero effect in practice, so it's fine as-is.

@author:
That's not quite right, but you are right that it can be simplified: replace `f2=lambda lam: np.inf` with `fillvalue=np.inf`.  I'll push that change.


@user1:
Oops, should have checked the signature. OK, I'll merge after that.

diff_hunk: @@ -7870,6 +7870,22 @@ def test_regression_tukey_lambda():
     assert_((p[2] == 0.0).any())
 
 
+def test_tukeylambda_support():
path: scipy/stats/tests/test_distributions.py
text: @user1:
I'd be willing to review a commit that brought all the scattered `tukeylambda` tests into one class if you're interested.

@author:
I almost did that originally, so I'll do that now.  I'm seeing some quirks (possible bugs) in `test_regression_tukey_lambda` that I'll clean up while I'm at it.

@author:
The bug is the test `assert_(~np.isnan(p).all())`.  The intent was almost certainly to check that there are no nans in `p`, but all it takes for that test to succeed is at least one non-nan value in `p`.  I'll change it to `assert np.isfinite(p).all()`.