Originally reported by @stefanv on mailing list, for an issue that was detected in skimage. My local testing suggests that the issue is reproducible for NumPy 1.26.4 and NumPy 2.0.0b1 though, so the reason may be due to a code change we've made directly.

I confirmed locally that the test below passes at v1.12.0 tag, and fails on main. I think Stefan's argument was that the docstring still says float for the statistic.

--- a/scipy/stats/tests/test_stats.py
+++ b/scipy/stats/tests/test_stats.py
@@ -8495,3 +8495,8 @@ class TestExpectile:
         for alpha in np.r_[0, alpha_seq, 1 - alpha_seq[:-1:-1], 1]:
             e_list.append(stats.expectile(x, alpha=alpha))
         assert np.all(np.diff(e_list) > 0)
+
+
+def test_pearsonr_stat_ret_type():
+    res = stats.pearsonr([1, 2, 3], [4, 5, 6])
+    assert isinstance(res.statistic, float)

Maybe one of the stats folks can comment re: what we might want to do here for 1.13.0? I'll see if I can bisect this for you folks.