Repository: rmarkello/abagen

Issue:
CLI is out-of-sync with main workflow
## The issue

Due to some recent updates the CLI has become out-of-sync with the main `abagen.get_expression_data()` workflow and needs to be updated. It would also be great to find a way to make it so that they could be kept in-sync more readily / without manual intervention.

## Proposed solution

Other than simply manually adding / updating all the new or changed parameters and help text in the CLI, I don't have much of a solution for this 🤷‍♂️ 

Discussion:

Patch:
diff --git a/abagen/allen.py b/abagen/allen.py
index 5fb7497..22c7d5f 100644
--- a/abagen/allen.py
+++ b/abagen/allen.py
@@ -42,7 +42,7 @@ def get_expression_data(atlas,
                         return_donors=False,
                         donors='all',
                         data_dir=None,
-                        verbose=1,
+                        verbose=0,
                         n_proc=1):
     """
     Assigns microarray expression data to ROIs defined in `atlas`
@@ -54,10 +54,13 @@ def get_expression_data(atlas,
 
         1. Intensity-based filtering of microarray probes to remove probes that
            do not exceed a certain level of background noise (specified via the
-           `ibf_threshold` parameter), and
+           `ibf_threshold` parameter),
         2. Selection of a single, representative probe (or collapsing across
            probes) for each gene, specified via the `probe_selection`
-           parameter.
+           parameter (and influenced by the `donor_probes` parameter), and
+        3. Optional mirroring of the tissue samples across the left/right
+           hemisphere boundary, as specified via the `lr_mirror` parameter
+           (turned off by default).
 
     Tissue samples are then matched to parcels in the defined `atlas` for each
     donor. If `atlas_info` is provided then this matching is constrained by
@@ -78,29 +81,33 @@ def get_expression_data(atlas,
            closest parcel, as determined by the parcel centroid.
 
     If at any step a sample can be assigned to a parcel the matching process is
-    terminated. If multiple sample are assigned to the same parcel they are
-    aggregated with the metric specified via the `metric` parameter. More
-    control over the sample matching can be obtained by setting the `exact`
-    parameter; see the parameter description for more information.
+    terminated. When the provided atlas is not volumetric (i.e., surface-based)
+    the samples are simply matched to the nearest vertex, and `tolerance` is
+    used as a standard deviation threshold. More control over the sample
+    matching can be obtained by setting the `exact` parameter; see the
+    parameter description for more information.
 
     Once all samples have been matched to parcels for all supplied donors, the
     microarray expression data are optionally normalized via the provided
-    `sample_norm` and `gene_norm` functions before being combined within
-    parcels and across donors via the supplied `agg_metric`.
+    `sample_norm` and `gene_norm` functions (which are influenced by the
+    `norm_matched` and `norm_structures` parameters) before being aggregated
+    across donors via the supplied `region_agg` and `agg_metric` parameters.
 
     Parameters
     ----------
     atlas : niimg-like object or dict
-        A parcellation image in MNI space, where each parcel is identified by a
-        unique integer ID. Alternatively, a dictionary where keys are donor IDs
-        and values are parcellation images in the native space of each donor.
+        A parcellation image in MNI space or a tuple of GIFTI images in
+        fsaverage5 space, where each parcel is identified by a unique integer
+        ID. Alternatively, a dictionary where keys are donor IDs and values are
+        parcellation images (or surfaces) in the native space of each donor.
     atlas_info : os.PathLike or pandas.DataFrame, optional
         Filepath to or pre-loaded dataframe containing information about
         `atlas`. Must have at least columns 'id', 'hemisphere', and 'structure'
-        containing information mapping atlas IDs to hemisphere (i.e, "L", "R")
-        and broad structural class (i.e., "cortex", "subcortex/brainstem",
+        containing information mapping atlas IDs to hemisphere (i.e, "L", "R",
+        "B") and broad structural class (i.e., "cortex", "subcortex/brainstem",
         "cerebellum"). If provided, this will constrain matching of tissue
-        samples to regions in `atlas`. Default: None
+        samples to regions in `atlas`. If `atlas` is a tuple of GIFTI images
+        with valid label tables this will be intuited. Default: None
     ibf_threshold : [0, 1] float, optional
         Threshold for intensity-based filtering. This number specifies the
         ratio of samples, across all supplied donors, for which a probe must
@@ -123,7 +130,7 @@ def get_expression_data(atlas,
         increase spatial coverage. Using 'bidirectional' will mirror samples
         across both hemispheres, 'leftright' will mirror samples in the left
         hemisphere to the right, and 'rightleft' will mirror the right to the
-        left. Default: 'bidirectional'
+        left. Default: None
     exact : bool, optional
         Whether to use exact matching of donor tissue samples to parcels in
         `atlas`. If True, this function will ONLY match tissue samples to
@@ -136,8 +143,10 @@ def get_expression_data(atlas,
         Euclidean distance to the parcel centroid). Default: True
     tolerance : int, optional
         Distance (in mm) that a sample must be from a parcel for it to be
-        matched to that parcel. This is only considered if the sample is not
-        directly within a parcel. Default: 2
+        matched to that parcel. If `atlas` is a tuple of surface files then
+        this measure is a standard deviation threshold (i.e., samples greater
+        than `tolerance` SDs away from the mean matched distance are ignored).
+        Default: 2
     sample_norm : {'rs', 'srs', 'minmax', 'center', 'zscore', None}, optional
         Method by which to normalize microarray expression values for each
         sample. Expression values are normalized separately for each sample and
diff --git a/abagen/cli/run.py b/abagen/cli/run.py
index 8b359f5..b67b88c 100644
--- a/abagen/cli/run.py
+++ b/abagen/cli/run.py
@@ -4,16 +4,25 @@ import argparse
 import logging
 import os
 from pathlib import Path
+from typing import Iterable
 import sys
 
 lgr = logging.getLogger('abagen')
 
 
+def isiterable(val):
+    """ Helper function to check whether value is iterable (but not string)
+    """
+    return isinstance(val, Iterable) and not isinstance(val, str)
+
+
 def _resolve_path(path):
     """ Helper function for get_parser() to resolve paths
     """
 
     if path is not None:
+        if isiterable(path):
+            return [_resolve_path(p) for p in path]
         try:
             return str(Path(path).expanduser().resolve())
         except FileNotFoundError:
@@ -34,7 +43,14 @@ class CheckExists(argparse.Action):
     """
     def __call__(self, parser, namespace, values, option_string=None):
         values = self.type(values)
-        if not os.path.exists(_resolve_path(values)):
+        missing = False
+        if isiterable(values):
+            missing = any(not os.path.exists(val) for val in values)
+            if len(values) == 1:
+                values = values[0]
+        else:
+            missing = not os.path.exists(values)
+        if missing:
             parser.error('Provided value for {} does not exist: {}'
                          .format(option_string, values))
         setattr(namespace, self.dest, values)
@@ -52,7 +68,7 @@ def get_parser():
     parser = argparse.ArgumentParser(
         formatter_class=argparse.RawDescriptionHelpFormatter,
         description="""
-Assigns microarray expression data to ROIs defined in the specific `atlas`
+Assigns microarray expression data to ROIs defined in the specified `atlas`
 
 This command aims to provide a workflow for generating pre-processed microarray
 expression data from the Allen Human Brain Atlas for arbitrary atlas
@@ -61,12 +77,16 @@ including:
 
     1. Intensity-based filtering of microarray probes to remove probes that do
        not exceed a certain level of background noise (specified via the
-       `ibf_threshold` parameter), and
+       `--ibf_threshold` parameter),
     2. Selection of a single, representative probe (or collapsing across
-       probes) for each gene, specified via the `probe_selection` parameter.
+       probes) for each gene, specified via the `--probe_selection` parameter
+       (and influenced by the `--donor_probes` parameter), and
+    3. Optional mirroring of the tissue samples across the left/right
+       hemisphere boundary, as specified via the `--lr_mirror` parameter
+       (turned off by default).
 
 Tissue samples are then matched to parcels in the defined `atlas` for each
-donor. If `atlas_info` is provided then this matching is constrained by both
+donor. If `--atlas_info` is provided then this matching is constrained by both
 hemisphere and tissue class designation (e.g., cortical samples from the left
 hemisphere are only matched to ROIs in the left cortex, subcortical samples
 from the right hemisphere are only matched to ROIs in the left subcortex); see
@@ -78,26 +98,29 @@ process:
     1. Determine if the sample falls directly within a parcel,
     2. Check to see if there are nearby parcels by slowly expanding the search
        space to include nearby voxels, up to a specified distance (specified
-       via the `tolerance` parameter),
+       via the `--tolerance` parameter),
     3. If there are multiple nearby parcels, the sample is assigned to the
        closest parcel, as determined by the parcel centroid.
 
 If at any step a sample can be assigned to a parcel the matching process is
-terminated. If multiple sample are assigned to the same parcel they are
-aggregated with the metric specified via the `agg_metric` parameter. More
-control over the sample matching can be obtained by setting the `inexact`
-parameter; see the parameter description for more information.
+terminated. When the provided atlas is not volumetric (i.e., is surface-based)
+the samples are simply matched to the nearest vertex, and `--tolerance` is used
+as a standard deviation threshold. More control over the sample matching can be
+obtained by setting the `--inexact` parameter.
 
 Once all samples have been matched to parcels for all supplied donors, the
 microarray expression data are optionally normalized via the provided
-`sample_norm` and `gene_norm` functions before being combined within parcels
-and across donors via the supplied `agg_metric`.
+`--sample_norm` and `--gene_norm` functions (which are influenced by the
+`--norm_matched` and `--norm_structures` parameters) before being aggregated
+across donors via the supplied `--region_agg` and `--agg_metric` parameters.
 """
     )
 
     parser.add_argument('atlas', action=CheckExists, type=_resolve_path,
-                        help='An image in MNI space, where each parcel in the '
-                             'image is identified by a unique integer ID.')
+                        nargs='+',
+                        help='A NIFTI image in MNI152 space or two GIFTI '
+                             'images in fsaverage5 space, where each parcel '
+                             'is identified by a unique integer ID.')
 
     # because I like consistency in capitalization and punctuation...
     for act in parser._actions:
@@ -106,25 +129,27 @@ and across donors via the supplied `agg_metric`.
             break
 
     parser.add_argument('--version', action='version', version=verstr,
-                        help='Show program\'s version number and exit.')
-    parser.add_argument('-v', '--verbose', action='count', default=1,
+                        help='Show program version and exit.')
+    parser.add_argument('-v', '--verbose', action='count', default=0,
                         help='Increase verbosity of status messages to '
                              'display during workflow.')
-    parser.add_argument('-q', '--quiet', action='store_true', default=False,
-                        help='Suppress all status messages during workflow.')
     parser.add_argument('--debug', action='store_true', help=argparse.SUPPRESS)
 
     a_data = parser.add_argument_group('Options to specify information about '
                                        'the atlas used')
     a_data.add_argument('--atlas_info', '--atlas-info', action=CheckExists,
                         type=_resolve_path, default=None, metavar='PATH',
-                        help='Filepath to CSV files containing information '
+                        help='Filepath to CSV file containing information '
                              'about `atlas`. The CSV file must have at least '
                              'columns ["id", "hemisphere", "structure"] which '
                              'contain information mapping the atlas IDs to '
-                             'hemispheres (i.e, "L", "R") and broad '
-                             'structural groups (i.e., "cortex", "subcortex", '
-                             '"cerebellum", "brainstem").')
+                             'hemispheres (i.e, "L", "R", or "B") and broad '
+                             'structural groups (i.e., "cortex", "subcortex/'
+                             'brainstem", "cerebellum"). If provided, this '
+                             'will constrain matching of tissue samples to '
+                             'regions in `atlas`. If the supplied `atlas` is '
+                             'a pair of GIFTI files with valid label tables '
+                             'this information will be intuited.')
 
     g_data = parser.add_argument_group('Options to specify which AHBA data to '
                                        'use during processing')
@@ -146,8 +171,39 @@ and across donors via the supplied `agg_metric`.
                              'exist at one of those locations then it will be '
                              'downloaded to the first of these location that '
                              'exists and for which write access is enabled.')
+    g_data.add_argument('--n_proc', '--n-proc', action='store', type=int,
+                        default=1,
+                        help='Number of processors to use to download AHBA '
+                             'data. Can paralellize up to six times if all '
+                             'donors are requested. Default: 1')
 
     w_data = parser.add_argument_group('Options to specify processing options')
+    w_data.add_argument('--ibf_threshold', '--ibf-threshold', action='store',
+                        default=0.5, metavar='THRESHOLD',
+                        help='Threshold for intensity-based filtering of '
+                             'probes. This number should specify the ratio of '
+                             'samples, across all supplied donors, for which '
+                             'a probe must have signal above background noise '
+                             'in order to be retained. Default: 0.5')
+    w_data.add_argument('--probe_selection', '--probe-selection',
+                        action='store', default='diff_stability',
+                        metavar='METHOD', choices=sorted(SELECTION_METHODS),
+                        help='Selection method for subsetting (or collapsing '
+                             'across) probes that index the same gene. Must '
+                             'be one of {"average", "mean", "max_intensity", '
+                             '"max_variance", "pc_loading", "corr_variance", '
+                             '"corr_intensity", "diff_stability", "rnaseq"}. '
+                             'Default: "diff_stability"')
+    w_data.add_argument('--lr_mirror', '--lr-mirror', metavar='METHOD',
+                        type=_resolve_none, default=None, choices=(
+                            None, 'bidirectional', 'leftright', 'rightleft'),
+                        help='Whether to mirror microarray expression samples '
+                             'across hemispheres to increase spatial coverage.'
+                             ' Using "bidirectional" will mirror samples '
+                             'across both hemispheres, "leftright" will '
+                             'mirror samples in the left hemisphere to the '
+                             'right, and "rightleft" will mirror the right to '
+                             'the left. Default: None')
     w_data.add_argument('--inexact', dest='exact', action='store_false',
                         default=True,
                         help='Whether to use inexact matching of donor tissue '
@@ -167,16 +223,45 @@ and across donors via the supplied `agg_metric`.
     w_data.add_argument('--tol', '--tolerance', dest='tolerance',
                         action='store', type=float, default=2,
                         help='Distance (in mm) that a sample can be from a '
-                             'parcel for it to be matched to that parcel. '
-                             'This is only considered if the sample is not '
-                             'directly within a parcel. Default: 2')
-    w_data.add_argument('--ibf_threshold', '--ibf-threshold', action='store',
-                        default=0.5, metavar='THRESHOLD',
-                        help='Threshold for intensity-based filtering of '
-                             'probes. This number should specify the ratio of '
-                             'samples, across all supplied donors, for which '
-                             'a probe must have signal above background noise '
-                             'in order to be retained. Default: 0.5')
+                             'parcel for it to be matched to that parcel. If '
+                             '`atlas` is GIFTI files then this measure is a '
+                             'standard deviation threshold (i.e., samples '
+                             'greater than `tolerance` SDs away from the mean '
+                             'matched distance are ignored). Default: 2')
+    w_data.add_argument('--sample_norm', '--sample-norm', action='store',
+                        default='srs', metavar='METHOD', type=_resolve_none,
+                        choices=sorted(NORMALIZATION_METHODS) + ['None', None],
+                        help='Method by which to normalize microarray '
+                             'expression values for each sample prior to '
+                             'collapsing into regions in `atlas`. Expression '
+                             'values are normalized separately for each '
+                             'sample and donor across genes. If None is '
+                             'specified then no normalization is performed. '
+                             'Default: "srs"')
+    w_data.add_argument('--gene_norm', '--gene-norm', action='store',
+                        default='srs', metavar='METHOD', type=_resolve_none,
+                        choices=sorted(NORMALIZATION_METHODS) + ['None', None],
+                        help='Method by which to normalize microarray '
+                             'expression values for each donor prior to '
+                             'collapsing across donors. Expression values are '
+                             'normalized separately for each gene for each '
+                             'donor across all expression samples. If None is '
+                             'specified then no normalization is performed. '
+                             'Default: "srs"')
+    w_data.add_argument('--norm_all', '--norm-all', dest='norm_matched',
+                        action='store_false', default=True,
+                        help='Whether to perform gene normalization '
+                             '(`gene_norm`) across all available samples '
+                             'instead of only across samples that were '
+                             'matched to regions in `atlas`. If `atlas` is '
+                             'very small (i.e., only a few regions of '
+                             'interest) using `--norm_all` is suggested.')
+    w_data.add_argument('--norm_structures', '--norm-structures',
+                        action='store_true', default=False,
+                        help='Whether to perform gene normalization '
+                             '(`gene_norm`) within structural classes (i.e., '
+                             '"cortex", "subcortex/brainstem", "cerebellum") '
+                             'instead of across all available samples.')
     w_data.add_argument('--region_agg', '--region-agg', action='store',
                         default='donors', metavar='METHOD',
                         choices=['donors', 'samples'],
@@ -200,42 +285,6 @@ and across donors via the supplied `agg_metric`.
                              'data into a single "group" expression '
                              'dataframe. Must be one of {"mean", "median"}. '
                              'Default: "mean"')
-    w_data.add_argument('--probe_selection', '--probe-selection',
-                        action='store', default='diff_stability',
-                        metavar='METHOD', choices=sorted(SELECTION_METHODS),
-                        help='Selection method for subsetting (or collapsing '
-                             'across) probes that index the same gene. Must '
-                             'be one of {"average", "mean", "max_intensity", '
-                             '"max_variance", "pc_loading", "corr_variance", '
-                             '"corr_intensity", "diff_stability", "rnaseq"}. '
-                             'Default: "diff_stability"')
-    w_data.add_argument('--lr_mirror', '--lr-mirror', action='store_true',
-                        help='Whether to mirror microarray expression samples '
-                             'across hemispheres to increase spatial '
-                             'coverage. This will duplicate samples across '
-                             'both hemispheres (i.e., L->R and R->L), '
-                             'approximately doubling the number of available '
-                             'samples. Default: False (i.e., no mirroring)')
-    w_data.add_argument('--gene_norm', '--gene-norm', action='store',
-                        default='srs', metavar='METHOD', type=_resolve_none,
-                        choices=sorted(NORMALIZATION_METHODS) + ['None', None],
-                        help='Method by which to normalize microarray '
-                             'expression values for each donor prior to '
-                             'collapsing across donors. Expression values are '
-                             'normalized separately for each gene for each '
-                             'donor across all expression samples. If None is '
-                             'specified then no normalization is performed. '
-                             'Default: "srs"')
-    w_data.add_argument('--sample_norm', '--sample-norm', action='store',
-                        default='srs', metavar='METHOD', type=_resolve_none,
-                        choices=sorted(NORMALIZATION_METHODS) + ['None', None],
-                        help='Method by which to normalize microarray '
-                             'expression values for each sample prior to '
-                             'collapsing into regions in `atlas`. Expression '
-                             'values are normalized separately for each '
-                             'sample and donor across genes. If None is '
-                             'specified then no normalization is performed. '
-                             'Default: "srs"')
 
     p_data = parser.add_argument_group('Options to modify the AHBA data used')
     p_data.add_argument('--no-reannotated', '--no_reannotated',
@@ -299,10 +348,6 @@ def main(args=None):
 
     opts = get_parser().parse_args(args)
 
-    # quiet overrides any verbosity setting
-    if opts.quiet:
-        opts.verbose = 0
-
     # debugging is fun
     if opts.debug:
         print(opts)
@@ -311,22 +356,25 @@ def main(args=None):
     # run the workflow
     expression = get_expression_data(atlas=opts.atlas,
                                      atlas_info=opts.atlas_info,
-                                     exact=opts.exact,
-                                     tolerance=opts.tolerance,
-                                     region_agg=opts.region_agg,
-                                     agg_metric=opts.agg_metric,
                                      ibf_threshold=opts.ibf_threshold,
                                      probe_selection=opts.probe_selection,
                                      lr_mirror=opts.lr_mirror,
-                                     gene_norm=opts.gene_norm,
+                                     exact=opts.exact,
+                                     tolerance=opts.tolerance,
                                      sample_norm=opts.sample_norm,
+                                     gene_norm=opts.gene_norm,
+                                     norm_matched=opts.norm_matched,
+                                     norm_structures=opts.norm_structures,
+                                     region_agg=opts.region_agg,
+                                     agg_metric=opts.agg_metric,
                                      corrected_mni=opts.corrected_mni,
                                      reannotated=opts.reannotated,
                                      return_counts=opts.save_counts,
                                      return_donors=opts.save_donors,
                                      donors=opts.donors,
                                      data_dir=opts.data_dir,
-                                     verbose=opts.verbose)
+                                     verbose=opts.verbose,
+                                     n_proc=opts.n_proc)
 
     output_path = os.path.dirname(opts.output_file)
     fname_pref = os.path.splitext(os.path.basename(opts.output_file))[0]


Test patch:
diff --git a/abagen/tests/cli/test_run.py b/abagen/tests/cli/test_run.py
index bf00d13..1ef7dc8 100644
--- a/abagen/tests/cli/test_run.py
+++ b/abagen/tests/cli/test_run.py
@@ -40,20 +40,15 @@ def test_run_get_parser(capsys, atlas, datadir):
     assert 'abagen {}'.format(version) == capsys.readouterr().out.strip()
 
     # arguments with invalid choices (probe_selection) raise errors
-    with pytest.raises(SystemExit):
-        parser.parse_args(['--probe-selection', 'notamethod', atlas['image']])
-    assert "invalid choice: 'notamethod'" in capsys.readouterr().err
-    with pytest.raises(SystemExit):
-        parser.parse_args(['--gene-norm', 'notamethod', atlas['image']])
-    assert "invalid choice: 'notamethod'" in capsys.readouterr().err
-    with pytest.raises(SystemExit):
-        parser.parse_args(['--sample-norm', 'notamethod', atlas['image']])
-    assert "invalid choice: 'notamethod'" in capsys.readouterr().err
+    for arg in ('--probe-selection', '--gene-norm',
+                '--sample-norm', '--lr-mirror'):
+        with pytest.raises(SystemExit):
+            parser.parse_args([arg, 'notamethod', atlas['image']])
+        assert "invalid choice: 'notamethod'" in capsys.readouterr().err
 
     # just test every option
     args = parser.parse_args([
         '-v', '-v', '-v',
-        '--quiet',
         '--debug',
         '--atlas-info', atlas['info'],
         '--donors', '12876', '15496',
@@ -64,7 +59,9 @@ def test_run_get_parser(capsys, atlas, datadir):
         '--region-agg', 'donors',
         '--agg-metric', 'median',
         '--probe-selection', 'average',
-        '--lr_mirror',
+        '--norm-all',
+        '--norm-structures',
+        '--lr_mirror', 'bidirectional',
         '--gene-norm', 'srs',
         '--sample-norm', 'srs',
         '--no-reannotated', '--no-corrected-mni',