Fix failing doctests on Python 3.13

- Fix discrimination doctest coordinate order by enforcing consistent ordering
- Suppress NumPy scalar conversion warnings in multipletests
- Update pearson_r_eff_p_value doctest to reflect behavior change from #437

🤖 Generated with [Claude Code](https://claude.com/claude-code)

Co-Authored-By: Claude <noreply@anthropic.com>

Author: aaronspring

xskillscore/core/deterministic.py

@@ -465,7 +465,7 @@ def pearson_r_eff_p_value(
     <xarray.DataArray (x: 3, y: 3)> Size: 72B
     array([[0.82544245,        nan, 0.25734167],
            [0.78902959, 0.57503354, 0.8059353 ],
-           [0.79242625, 0.66792245, 1.        ]])
+           [0.79242625, 0.66792245,        nan]])
     Dimensions without coordinates: x, y
     """
     _fail_if_dim_empty(dim)

xskillscore/core/probabilistic.py

@@ -1006,8 +1006,8 @@ def discrimination(
     array([[0.00437637, 0.15536105, 0.66739606, 0.12472648, 0.04814004],
            [0.00451467, 0.16704289, 0.66365688, 0.1241535 , 0.04063205]])
     Coordinates:
-      * forecast_probability  (forecast_probability) float64 40B 0.1 0.3 0.5 0.7 0.9
       * event                 (event) bool 2B True False
+      * forecast_probability  (forecast_probability) float64 40B 0.1 0.3 0.5 0.7 0.9
 
     References
     ----------
@@ -1032,9 +1032,20 @@ def discrimination(
         dim=dim,
     ) / (np.logical_not(observations)).sum(dim=dim)
 
-    return xr.concat([hist_event, hist_no_event], dim="event").assign_coords(
+    result = xr.concat([hist_event, hist_no_event], dim="event").assign_coords(
         {"event": [True, False]}
     )
+    # Ensure consistent dimension and coordinate order across versions
+    result = result.transpose("event", FORECAST_PROBABILITY_DIM, ...)
+    # Reconstruct to ensure coordinate order
+    return xr.DataArray(
+        result.values,
+        dims=result.dims,
+        coords={
+            "event": result.coords["event"],
+            FORECAST_PROBABILITY_DIM: result.coords[FORECAST_PROBABILITY_DIM],
+        },
+    )
 
 
 def reliability(

xskillscore/core/stattests.py

@@ -1,5 +1,6 @@
 from __future__ import annotations
 
+import warnings
 from typing import Literal, Mapping, Optional, Tuple, Union
 
 import xarray as xr
@@ -159,24 +160,32 @@ def multipletests(
             f"Expected `return_results` from {allowed_return_results}, found {return_results}"
         )
 
-    ret = xr.apply_ufunc(
-        statsmodels_multipletests,
-        p.stack(s=p.dims),
-        input_core_dims=[[]],
-        vectorize=True,
-        output_core_dims=[[]] * 4,
-        output_dtypes=[bool, float, float, float],
-        kwargs=dict(method=method, alpha=alpha, **multipletests_kwargs),
-        dask="parallelized",
-        keep_attrs=keep_attrs,
-    )
+    # Suppress NumPy scalar conversion deprecation warning from internal numpy operations
+    with warnings.catch_warnings():
+        warnings.filterwarnings(
+            "ignore",
+            message="Conversion of an array with ndim > 0 to a scalar is deprecated",
+            category=DeprecationWarning,
+        )
+        ret = xr.apply_ufunc(
+            statsmodels_multipletests,
+            p.stack(s=p.dims),
+            input_core_dims=[[]],
+            vectorize=True,
+            output_core_dims=[[]] * 4,
+            output_dtypes=[bool, float, float, float],
+            kwargs=dict(method=method, alpha=alpha, **multipletests_kwargs),
+            dask="parallelized",
+            keep_attrs=keep_attrs,
+        )
 
     ret = tuple(r.unstack("s").transpose(*p.dims, ...) for r in ret)
 
     def _add_kwargs_as_coords(r: XArray):
-        r.coords["multipletests_method"] = method
-        r.coords["multipletests_alpha"] = alpha
-        return r
+        return r.assign_coords(
+            multipletests_method=method,
+            multipletests_alpha=alpha
+        )
 
     ret = tuple(_add_kwargs_as_coords(r) for r in ret)