code: def _compute_individual_regions(
      domain, stats, threshold=3.0, smin=5, method="gauss_mixture"
  ):
      """Compute the individual regions that are real activation candidates
  
      Parameters
      ----------
      domain : StructuredDomain instance,
            generic descriptor of the space domain
      stats: an array of shape (n_voxels, n_subjects)
             the multi-subject statistical maps
      threshold: float, optional
             first level threshold
      smin: int, optional
            minimal size of the regions to validate them
      method: {'gauss_mixture', 'emp_null', 'gam_gauss', 'prior'}, optional,
             'gauss_mixture' A Gaussian Mixture Model is used
             'emp_null' a null mode is fitted to test
             'gam_gauss' a Gamma-Gaussian mixture is used
             'prior' a hard-coded function is used
  
      Returns
      -------
      hrois: list of nipy.labs.spatial_models.hroi.HierrachicalROI instances
             that represent individual ROIs
             let nr be the number of terminal regions across subjects
      prior_h0: array of shape (nr),
                the mixture-based prior probability
                that the terminal regions are false positives
      subjects: array of shape (nr),
                the subject index associated with the terminal regions
      coords: array of shape (nr, coord.shape[1]),
              the coordinates of the of the terminal regions
  
      Notes
      -----
      FIXME: should allow for subject specific domains
      """
      hrois = []
      coords = []
      prior_h0 = []
      subjects = []
      n_subjects = stats.shape[1]
      nvox = stats.shape[0]
  
      for subject in range(n_subjects):
          # description in terms of blobs
          stats_ = np.reshape(stats[:, subject], (nvox, 1))
          hroi = HROI_as_discrete_domain_blobs(
              domain, stats_, threshold=threshold, smin=smin
          )
  
          if hroi is not None and hroi.k > 0:
              # get the leave regions (i.e. the local maxima)
              leaves = [
                  hroi.select_id(id) for id in hroi.get_leaves_id()
              ]
  
              # get the region mean statistical value
              mean_val = hroi.representative_feature(
                  "signal", "weighted mean"
              )
              mean_val = mean_val[leaves]
  
              # get the regions position
              mean_pos = np.asarray(
                  [np.mean(coord, 0) for coord in hroi.get_coord()]
              )
              hroi.set_roi_feature("position", mean_pos)
              coords.append(mean_pos[leaves])
  
              # compute the prior proba of being null
              learning_set = np.squeeze(stats_[stats_ != 0])
              prior_h0.append(
                  _stat_to_proba(mean_val, learning_set, method)
              )
              subjects.append(
                  subject * np.ones(mean_val.size).astype(np.int_)
              )
          else:
              subjects.append([])
              prior_h0.append([])
              coords.append(np.empty((0, domain.dim)))
          hrois.append(hroi)
  
      prior_h0 = np.concatenate(prior_h0)
      subjects = np.concatenate(subjects)
      coords = np.concatenate(coords)
      return hrois, prior_h0, subjects, coords
  
file_path: nipy/nipy/labs/spatial_models/bayesian_structural_analysis.py
imports: import numpy
  import scipy.stats
  from nipy.algorithms.graph import wgraph_from_coo_matrix
  from ...algorithms.statistics.empirical_pvalue import NormalEmpiricalNull, gamma_gaussian_fit, three_classes_GMM_fit
  from .hroi import HROI_as_discrete_domain_blobs
  from .structural_bfls import build_landmarks
label_desc: Use image node_id 2 for calling a global function with example usage: _compute_individual_regions(domain, stats, threshold, smin, method) and returns: hrois,  prior_h0,  subjects,  coords