jax.scipy.stats.beta.pdf handles x=0, x=1 incorrectly
```
a = b = 1.
x = np.array([0., 1.])
scipy.stats.beta.pdf(x, a, b)
# array([1., 1.])
jax.scipy.stats.beta.pdf(x, a, b)
# DeviceArray([nan, nan], dtype=float32) 
```

Same issue and resolution as #7256 basically--use xlogy/xlog1py here:
https://github.com/google/jax/blob/f6f1debf700b55ad8c0ecc1daa0ed1f677c2bf00/jax/_src/scipy/stats/beta.py#L30-L31


diff --git a/jax/_src/scipy/stats/beta.py b/jax/_src/scipy/stats/beta.py
index 64a9e59cb..eb3135de3 100644
--- a/jax/_src/scipy/stats/beta.py
+++ b/jax/_src/scipy/stats/beta.py
@@ -18,7 +18,7 @@ from jax import lax
 from jax._src.numpy.util import _wraps
 from jax._src.numpy.lax_numpy import (_promote_args_inexact, _constant_like,
                                       where, inf, logical_or)
-from jax.scipy.special import betaln
+from jax.scipy.special import betaln, xlogy, xlog1py
 
 
 @_wraps(osp_stats.beta.logpdf, update_doc=False)
@@ -27,8 +27,8 @@ def logpdf(x, a, b, loc=0, scale=1):
   one = _constant_like(x, 1)
   shape_term = lax.neg(betaln(a, b))
   y = lax.div(lax.sub(x, loc), scale)
-  log_linear_term = lax.add(lax.mul(lax.sub(a, one), lax.log(y)),
-                            lax.mul(lax.sub(b, one), lax.log1p(lax.neg(y))))
+  log_linear_term = lax.add(xlogy(lax.sub(a, one), y),
+                            xlog1py(lax.sub(b, one), lax.neg(y)))
   log_probs = lax.sub(lax.add(shape_term, log_linear_term), lax.log(scale))
   return where(logical_or(lax.gt(x, lax.add(loc, scale)),
                           lax.lt(x, loc)), -inf, log_probs)
diff --git a/tests/scipy_stats_test.py b/tests/scipy_stats_test.py
index 548a5e373..154f29b03 100644
--- a/tests/scipy_stats_test.py
+++ b/tests/scipy_stats_test.py
@@ -148,6 +148,13 @@ class LaxBackedScipyStatsTests(jtu.JaxTestCase):
     self._CompileAndCheck(lax_fun, args_maker,
                           rtol={np.float32: 2e-3, np.float64: 1e-4})
 
+  def testBetaLogPdfZero(self):
+    # Regression test for https://github.com/google/jax/issues/7645
+    a = b = 1.
+    x = np.array([0., 1.])
+    self.assertAllClose(
+      osp_stats.beta.pdf(x, a, b), lsp_stats.beta.pdf(x, a, b), atol=1E-6)
+
   @genNamedParametersNArgs(3)
   def testCauchyLogPdf(self, shapes, dtypes):
     rng = jtu.rand_default(self.rng())
['tests/scipy_stats_test.py', 'jax/_src/scipy/stats/beta.py']