docs: def initialize_means(data, clusters, k):
      """
      Initializes the M matrix given the data and a set of cluster labels.
      Cluster centers are set to the mean of each cluster.
  
      Args:
          data (array): genes x cells
          clusters (array): 1d array of ints (0...k-1)
          k (int): number of clusters
      """
      init_w = np.zeros((data.shape[0], k))
      if sparse.issparse(data):
          for i in range(k):
              if data[:,clusters==i].shape[1]==0:
                  point = np.random.randint(0, data.shape[1])
                  init_w[:,i] = data[:,point].toarray().flatten()
              else:
                  # memory usage might be a problem here?
                  init_w[:,i] = np.array(data[:,clusters==i].mean(1)).flatten() + eps
      else:
          for i in range(k):
              if data[:,clusters==i].shape[1]==0:
                  point = np.random.randint(0, data.shape[1])
                  init_w[:,i] = data[:,point].flatten()
              else:
                  init_w[:,i] = data[:,clusters==i].mean(1) + eps
      return init_w