code: def h_hat(ac):
      
  
      # check inputs
      ac = asarray_ndim(ac, 2)
      assert ac.shape[1] == 2, 'only biallelic variants supported'
  
      # compute allele number
      an = ac.sum(axis=1)
  
      # compute estimator
      x = (ac[:, 0] * ac[:, 1]) / (an * (an - 1))
  
      return x
docs: Unbiased estimator for h, where 2*h is the heterozygosity
      of the population.
  
      Parameters
      ----------
      ac : array_like, int, shape (n_variants, 2)
          Allele counts array for a single population.
  
      Returns
      -------
      h_hat : ndarray, float, shape (n_variants,)
  
      Notes
      -----
      Used in Patterson (2012) for calculation of various statistics.