docs: def calculate_gap(data, clustering, km, B=50, **kwargs):
      """
      See: https://datasciencelab.wordpress.com/2013/12/27/finding-the-k-in-k-means-clustering/
  
      https://web.stanford.edu/~hastie/Papers/gap.pdf
  
      Returns two results: the gap score, and s_k.
      """
      k = len(set(clustering))
      Wk = km.inertia_
      mins, maxes = calculate_bounding_box(data)
      Wk_est = []
      for i in range(B):
          data_sample = (maxes-mins)*np.random.random(data.shape) + mins
          km_b = KMeans(k)
          km_b.fit_predict(data_sample)
          Wk_est.append(km_b.inertia_)
      Wk_est = np.log(np.array(Wk_est))
      Wk_mean = np.mean(Wk_est)
      Wk_std = np.std(Wk_est)
      gap = Wk_mean - np.log(Wk)
      sk = np.sqrt(1 + 1.0/B)*Wk_std
      return gap, sk