code: def bayesfactor_pearson(
      r, n, alternative="two-sided", method="ly", kappa=1.0
  ):
      """
      Bayes Factor of a Pearson correlation.
  
      Parameters
      ----------
      r : float
          Pearson correlation coefficient.
      n : int
          Sample size.
      alternative : string
          Defines the alternative hypothesis, or tail of the correlation. Must be one of
          "two-sided" (default), "greater" or "less". Both "greater" and "less" return a one-sided
          p-value. "greater" tests against the alternative hypothesis that the correlation is
          positive (greater than zero), "less" tests against the hypothesis that the correlation is
          negative.
      method : str
          Method to compute the Bayes Factor. Can be "ly" (default) or
          "wetzels". The former has an exact analytical solution, while the
          latter requires integral solving (and is therefore slower). "wetzels"
          was the default in Pingouin <= 0.2.5. See Notes for details.
      kappa : float
          Kappa factor. This is sometimes called the *rscale* parameter, and
          is only used when ``method`` is "ly".
  
      Returns
      -------
      bf : float
          Bayes Factor (BF10).
          The Bayes Factor quantifies the evidence in favour of the alternative
          hypothesis.
  
      See also
      --------
      corr : (Robust) correlation between two variables
      pairwise_corr : Pairwise correlation between columns of a pandas DataFrame
      bayesfactor_ttest : Bayes Factor of a T-test
      bayesfactor_binom : Bayes Factor of a binomial test
  
      Notes
      -----
      To compute the Bayes Factor directly from the raw data, use the
      :py:func:`pingouin.corr` function.
  
      The two-sided **Wetzels Bayes Factor** (also called *JZS Bayes Factor*)
      is calculated using the equation 13 and associated R code of [1]_:
  
      .. math::
  
          \\text{BF}_{10}(n, r) = \\frac{\\sqrt{n/2}}{\\gamma(1/2)}*
          \\int_{0}^{\\infty}e((n-2)/2)*
          log(1+g)+(-(n-1)/2)log(1+(1-r^2)*g)+(-3/2)log(g)-n/2g
  
      where :math:`n` is the sample size, :math:`r` is the Pearson correlation
      coefficient and :math:`g` is is an auxiliary variable that is integrated
      out numerically. Since the Wetzels Bayes Factor requires solving an
      integral, it is slower than the analytical solution described below.
  
      The two-sided **Ly Bayes Factor** (also called *Jeffreys
      exact Bayes Factor*) is calculated using equation 25 of [2]_:
  
      .. math::
  
          \\text{BF}_{10;k}(n, r) = \\frac{2^{\\frac{k-2}{k}}\\sqrt{\\pi}}
          {\\beta(\\frac{1}{k}, \\frac{1}{k})} \\cdot
          \\frac{\\Gamma(\\frac{2+k(n-1)}{2k})}{\\Gamma(\\frac{2+nk}{2k})}
          \\cdot 2F_1(\\frac{n-1}{2}, \\frac{n-1}{2}, \\frac{2+nk}{2k}, r^2)
  
      The one-sided version is described in eq. 27 and 28 of Ly et al, 2016.
      Please take note that the one-sided test requires the
      `mpmath <http://mpmath.org/>`_ package.
  
      Results have been validated against JASP and the BayesFactor R package.
  
      References
      ----------
      .. [1] Ly, A., Verhagen, J. & Wagenmakers, E.-J. Harold Jeffreys’s default
         Bayes factor hypothesis tests: Explanation, extension, and
         application in psychology. J. Math. Psychol. 72, 19–32 (2016).
  
      .. [2] Wetzels, R. & Wagenmakers, E.-J. A default Bayesian hypothesis test
         for correlations and partial correlations. Psychon. Bull. Rev. 19,
         1057–1064 (2012).
  
      Examples
      --------
      Bayes Factor of a Pearson correlation
  
      >>> from pingouin import bayesfactor_pearson
      >>> r, n = 0.6, 20
      >>> bf = bayesfactor_pearson(r, n)
      >>> print("Bayes Factor: %.3f" % bf)
      Bayes Factor: 10.634
  
      Compare to Wetzels method:
  
      >>> bf = bayesfactor_pearson(r, n, method='wetzels')
      >>> print("Bayes Factor: %.3f" % bf)
      Bayes Factor: 8.221
  
      One-sided test
  
      >>> bf10pos = bayesfactor_pearson(r, n, alternative='greater')
      >>> bf10neg = bayesfactor_pearson(r, n, alternative='less')
      >>> print("BF-pos: %.3f, BF-neg: %.3f" % (bf10pos, bf10neg))
      BF-pos: 21.185, BF-neg: 0.082
      """
      from scipy.special import gamma, betaln, hyp2f1
  
      assert method.lower() in [
          "ly",
          "wetzels",
      ], "Method not recognized."
      assert alternative in [
          "two-sided",
          "greater",
          "less",
      ], "Alternative must be one of 'two-sided' (default), 'greater' or 'less'."
  
      # Wrong input
      if not np.isfinite(r) or n < 2:
          return np.nan
      assert -1 <= r <= 1, "r must be between -1 and 1."
  
      if alternative != "two-sided" and method.lower() == "wetzels":
          warnings.warn(
              "One-sided Bayes Factor are not supported by the "
              "Wetzels's method. Switching to method='ly'."
          )
          method = "ly"
  
      if method.lower() == "wetzels":
          # Wetzels & Wagenmakers, 2012. Integral solving
  
          def fun(g, r, n):
              return exp(
                  ((n - 2) / 2) * log(1 + g)
                  + (-(n - 1) / 2) * log(1 + (1 - r**2) * g)
                  + (-3 / 2) * log(g)
                  + -n / (2 * g)
              )
  
          integr = quad(fun, 0, np.inf, args=(r, n))[0]
          bf10 = np.sqrt(n / 2) / gamma(1 / 2) * integr
  
      else:
          # Ly et al, 2016. Analytical solution.
          k = kappa
          lbeta = betaln(1 / k, 1 / k)
          log_hyperterm = log(
              hyp2f1(
                  ((n - 1) / 2),
                  ((n - 1) / 2),
                  ((n + 2 / k) / 2),
                  r**2,
              )
          )
          bf10 = exp(
              (1 - 2 / k) * log(2)
              + 0.5 * log(pi)
              - lbeta
              + lgamma((n + 2 / k - 1) / 2)
              - lgamma((n + 2 / k) / 2)
              + log_hyperterm
          )
  
          if alternative != "two-sided":
              # Directional test.
              # We need mpmath for the generalized hypergeometric function
              from .utils import _is_mpmath_installed
  
              _is_mpmath_installed(raise_error=True)
              from mpmath import hyp3f2
  
              hyper_term = float(
                  hyp3f2(
                      1,
                      n / 2,
                      n / 2,
                      3 / 2,
                      (2 + k * (n + 1)) / (2 * k),
                      r**2,
                  )
              )
              log_term = (
                  2 * (lgamma(n / 2) - lgamma((n - 1) / 2)) - lbeta
              )
              C = (
                  2 ** ((3 * k - 2) / k)
                  * k
                  * r
                  / (2 + (n - 1) * k)
                  * exp(log_term)
                  * hyper_term
              )
  
              bf10neg = bf10 - C
              bf10pos = 2 * bf10 - bf10neg
              if alternative == "greater":
                  # We expect the correlation to be positive
                  bf10 = bf10pos
              else:
                  # We expect the correlation to be negative
                  bf10 = bf10neg
  
      return bf10
  
imports: import warnings
  import numpy
  from scipy.integrate import quad
  from math import pi, exp, log, lgamma
label_desc: Use image node_id 3 for calling a global function with example usage: bayesfactor_pearson(r, n, alternative, method, kappa) and returns: bf10, np, exp