old_contents: import numpy as np
  from scipy import percentile
  #from scipy.stats import norm
  
  def _select_sigma(X, axis=0):
      """
      Returns the smaller of std(X, ddof=1) or normalized IQR(X)
  
      References
      ----------
      Silverman (1986) p.47
      """
  #    normalize = norm.ppf(.75) - norm.ppf(.25)
      normalize = 1.349
      IQR = np.subtract.reduce(percentile(X, [75,25],
                               axis=axis), axis=axis)/normalize
      return np.minimum(np.std(X, axis=axis, ddof=1), IQR)
  
  
  ## Univariate Rule of Thumb Bandwidths ##
  def bw_scott(x, axis=0):
      """
      Scott's Rule of Thumb
  
      Parameter
      ---------
      x : array-like
          Array for which to get the bandwidth
  
      Returns
      -------
      bw : float
          The estimate of the bandwidth
  
      Notes
      -----
      Returns 1.059 * A * n ** (-1/5.)
  
      A = min(std(x, ddof=1), IQR/1.349)
      IQR = np.subtract.reduce(scipy.percentile(x, [75,25]))
  
      References
      ---------- ::
  
      Scott, D.W. (1992) `Multivariate Density Estimation: Theory, Practice, and
          Visualization.`
      """
      A = _select_sigma(x, axis)
      n = len(x)
      return 1.059 * A * n ** -.2
  
  def bw_silverman(x, axis=0):
      """
      Silverman's Rule of Thumb
  
      Parameter
      ---------
      x : array-like
          Array for which to get the bandwidth
  
      Returns
      -------
      bw : float
          The estimate of the bandwidth
  
      Notes
      -----
      Returns .9 * A * n ** (-1/5.)
  
      A = min(std(x, ddof=1), IQR/1.349)
      IQR = np.subtract.reduce(scipy.percentile(x, [75,25]))
  
      References
      ---------- ::
  
      Silverman, B.W. (1986) `Density Estimation.`
      """
      A = _select_sigma(x, axis)
      n = len(x)
      return .9 * A * n ** -.2
  
  ## Plug-In Methods ##
  
  ## Least Squares Cross-Validation ##
  
  ## Helper Functions ##
  
  bandwidth_funcs = dict(scott=bw_scott,silverman=bw_silverman)
  
  def select_bandwidth(X, bw, kernel):
      """
      Selects bandwidth
      """
      bw = bw.lower()
      if bw not in ["scott","silverman"]:
          raise ValueError("Bandwidth %s not understood" % bw)
  #TODO: uncomment checks when we have non-rule of thumb bandwidths for diff. kernels
  #    if kernel == "gauss":
      return bandwidth_funcs[bw](X)
  #    else:
  #        raise ValueError("Only Gaussian Kernels are currently supported")
  
  
new_contents: import numpy as np
  from scipy.stats import scoreatpercentile as sap
  #from scipy.stats import norm
  
  def _select_sigma(X):
      """
      Returns the smaller of std(X, ddof=1) or normalized IQR(X) over axis 0.
  
      References
      ----------
      Silverman (1986) p.47
      """
  #    normalize = norm.ppf(.75) - norm.ppf(.25)
      normalize = 1.349
  #    IQR = np.subtract.reduce(percentile(X, [75,25],
  #                             axis=axis), axis=axis)/normalize
      IQR = (sap(X, 75) - sap(X, 25))/normalize
      return np.minimum(np.std(X, axis=0, ddof=1), IQR)
  
  
  ## Univariate Rule of Thumb Bandwidths ##
  def bw_scott(x):
      """
      Scott's Rule of Thumb
  
      Parameter
      ---------
      x : array-like
          Array for which to get the bandwidth
  
      Returns
      -------
      bw : float
          The estimate of the bandwidth
  
      Notes
      -----
      Returns 1.059 * A * n ** (-1/5.)
  
      A = min(std(x, ddof=1), IQR/1.349)
      IQR = np.subtract.reduce(np.percentile(x, [75,25]))
  
      References
      ---------- ::
  
      Scott, D.W. (1992) `Multivariate Density Estimation: Theory, Practice, and
          Visualization.`
      """
      A = _select_sigma(x)
      n = len(x)
      return 1.059 * A * n ** -.2
  
  def bw_silverman(x):
      """f
      Silverman's Rule of Thumb
  
      Parameter
      ---------
      x : array-like
          Array for which to get the bandwidth
  
      Returns
      -------
      bw : float
          The estimate of the bandwidth
  
      Notes
      -----
      Returns .9 * A * n ** (-1/5.)
  
      A = min(std(x, ddof=1), IQR/1.349)
      IQR = np.subtract.reduce(np.percentile(x, [75,25]))
  
      References
      ---------- ::
  
      Silverman, B.W. (1986) `Density Estimation.`
      """
      A = _select_sigma(x)
      n = len(x)
      return .9 * A * n ** -.2
  
  ## Plug-In Methods ##
  
  ## Least Squares Cross-Validation ##
  
  ## Helper Functions ##
  
  bandwidth_funcs = dict(scott=bw_scott,silverman=bw_silverman)
  
  def select_bandwidth(X, bw, kernel):
      """
      Selects bandwidth
      """
      bw = bw.lower()
      if bw not in ["scott","silverman"]:
          raise ValueError("Bandwidth %s not understood" % bw)
  #TODO: uncomment checks when we have non-rule of thumb bandwidths for diff. kernels
  #    if kernel == "gauss":
      return bandwidth_funcs[bw](X)
  #    else:
  #        raise ValueError("Only Gaussian Kernels are currently supported")
  
  
subject: Replace numpy.percentile with stats.scoreatpercentile