Repository: rmarkello/abagen

Issue:
Memory issues with `abagen.get_expression_data`
## The issue

Using the fantastic [`memory_profiler`](https://pypi.org/project/memory-profiler/) package and running a few tests with `abagen` reveals that, in the best case, the `get_expression_data()` functions uses ~4GiB RAM and, in the worst case, uses up to ~7GiB RAM.

This isn't *great*, and it would be nice to find ways to optimize this down a bit. Unfortunately, I think at a minimum we're going to be dealing with ~2GiB of RAM given that's about how much loading and storing the raw microarray + pacall data uses, but beyond that we can definitely try and get things down.

Currently, the worst of the memory usage is happening when `lr_mirror` is set to `True`. The step where we duplicate the samples has, for a brief period of time, both the original microarray data AND the duplicated microarray data in memory, leading to huge jumps until the original is removed.

## Proposed solution

I tried to resolve this issue by adding an optional keyword argument `inplace` to the `samples.mirror_samples` command, but that didn't seem to help a ton (at least not with the current implementation). I'm going to have to dig into the internals of things a bit more and try and figure out if there are stages where we can offload holding things in memory. 

The flipside is that offloading things from memory might increase runtime (due to the need to load them back in whenever they're needed), so perhaps this might warrant a `low_memory` parameter or something?

Discussion:

Patch:
diff --git a/.circleci/config.yml b/.circleci/config.yml
index 448b188..2287221 100644
--- a/.circleci/config.yml
+++ b/.circleci/config.yml
@@ -117,9 +117,9 @@ jobs:
             apt-get install -yqq make
             source activate abagen_py37
             pip install -e .
-            # python -c "import abagen; abagen.fetch_microarray(donors='all', n_proc=6)"
+            python -c "import abagen; abagen.fetch_microarray(donors='all', n_proc=6)"
             cd docs/
-            make html  # doctest
+            make html doctest
       - store_artifacts:
           path: /tmp/src/abagen/docs/_build
   unittest_py35:
diff --git a/.travis.yml b/.travis.yml
deleted file mode 100644
index 478fb53..0000000
--- a/.travis.yml
+++ /dev/null
@@ -1,31 +0,0 @@
-language: python
-sudo: false
-dist: xenial
-notifications:
-  email: change
-
-branches:
-  only:
-    - master
-
-python:
-  - 3.5
-  - 3.6
-  - 3.7
-
-env:
-  matrix:
-    - CHECK_TYPE=docdoctest
-
-before_install:
-  - python -m pip install --upgrade pip
-
-install:
-  - python setup.py install;
-
-script:
-  - |
-    cd docs;
-    pip install -r ./requirements.txt;
-    make html;
-    travis_wait 40 make doctest;
diff --git a/abagen/allen.py b/abagen/allen.py
index 29e89f4..1291ea0 100644
--- a/abagen/allen.py
+++ b/abagen/allen.py
@@ -338,44 +338,24 @@ def get_expression_data(atlas, atlas_info=None, *, exact=True,
                  .format(len(all_labels)))
         centroids = utils.get_centroids(atlas, labels=all_labels)
 
-    # are we using corrected MNI coordinates? update the annotation "files"
-    if corrected_mni:
-        for n, annot in enumerate(annotation):
-            annotation[n] = samples.update_mni_coords(annot)
+    # update the annotation "files"
+    # this handles updating the MNI coordinates, dropping mistmatched samples
+    # (where MNI coordinates don't match ontology), and mirror samples across
+    # hemisphere, as applicable.
+    for n, (annot, ont) in enumerate(zip(annotation, ontology)):
+        annotation[n] = samples.update_samples(annot, ont, lr_mirror=lr_mirror,
+                                               corrected_mni=corrected_mni)
 
     # get dataframe of probe information (reannotated or otherwise)
     probe_info = io.read_probes(probe_info, copy=True)
     if reannotated:
-        lgr.info('Reannotating microarray probes with information from '
-                 'Arnatkevic̆iūtė et al., 2019, NeuroImage')
+        lgr.info('Reannotating probes with information from Arnatkevic̆iūtė '
+                 'et al., 2019, NeuroImage')
         probe_info = probes.reannotate_probes(probe_info)
 
-    # when probes are reannotated there are a few that are removed due to their
-    # not being matched to a gene. in this case, it's good to subset microarray
-    # and pacall data to save a bit on memory
-    for n in range(len(microarray)):
-        microarray[n] = io.read_microarray(microarray[n]).loc[probe_info.index]
-        pacall[n] = io.read_pacall(pacall[n]).loc[probe_info.index]
-
-    # if we're mirroring samples across hemispheres we need to modify the
-    # microarray, pacall, and annotation dataframes. the other data (ontology
-    # and probes) are redundant across subjects and have no sample-specific
-    # information.
-    # if this is done with uncorrected MNI coords (i.e., corrected_mni=False)
-    # then mirroring will change the number of probes / genes in the generated
-    # microarray expression dataframes. with corrected MNI coords the number of
-    # probes and genes is consistent regardless of whether samples are mirrored
-    if lr_mirror:
-        lgr.info('Left/right mirroring requested; mirroring samples across '
-                 'hemispheres')
-        microarray, pacall, annotation = samples.mirror_samples(microarray,
-                                                                pacall,
-                                                                annotation,
-                                                                ontology,
-                                                                inplace=True)
-
     # intensity-based filtering of probes
-    probe_info = probes.filter_probes(pacall, probe_info, ibf_threshold)
+    probe_info = probes.filter_probes(pacall, annotation, probe_info,
+                                      threshold=ibf_threshold)
     lgr.info('{} probes survive intensity-based filtering with threshold of {}'
              .format(len(probe_info), ibf_threshold))
 
@@ -395,17 +375,13 @@ def get_expression_data(atlas, atlas_info=None, *, exact=True,
                                    dtype=int),
                           index=np.append([0], all_labels))
     for subj in range(len(microarray)):
-        # get rid of samples whose coordinates don't match ontological profile
-        annotation[subj] = samples.drop_mismatch_samples(annotation[subj],
-                                                         ontology[subj])
+        if lr_mirror:  # reset index (duplicates will cause issues if we don't)
+            annotation[subj] = annotation[subj].reset_index(drop=True)
+            microarray[subj] = microarray[subj].reset_index(drop=True)
 
         # assign samples to regions
         labels = samples.label_samples(annotation[subj], atlas,
-                                       atlas_info=atlas_info,
-                                       tolerance=tolerance)
-
-        # subset representative probes + samples from microarray data
-        microarray[subj] = microarray[subj].loc[annotation[subj].index]
+                                       atlas_info, tolerance=tolerance)
 
         if exact:  # remove all samples not assigned a label before norming
             nz = np.asarray(labels != 0).squeeze()
@@ -474,6 +450,8 @@ def get_expression_data(atlas, atlas_info=None, *, exact=True,
         expression = pd.concat(expression).groupby('label').aggregate(metric)
         # some genes may have been poorly normalized; remove these
         drop = expression.dropna(axis=0, how='all').isna().any(axis=0)
+        lgr.info('Dropping {} gene from concatenated expression data due to '
+                 'poor normalization'.format(drop.sum()))
         expression = expression.drop(drop[drop].index, axis=1)
 
     # drop the "zero" label from the counts dataframe (this is background)
diff --git a/abagen/probes.py b/abagen/probes.py
index 620204f..0f61bc7 100644
--- a/abagen/probes.py
+++ b/abagen/probes.py
@@ -58,7 +58,7 @@ def reannotate_probes(probes):
     return reannotated
 
 
-def filter_probes(pacall, probes, threshold=0.5):
+def filter_probes(pacall, samples, probes, threshold=0.5):
     """
     Performs intensity based filtering (IBF) of expression probes
 
@@ -72,6 +72,8 @@ def filter_probes(pacall, probes, threshold=0.5):
         List of filepaths to PAcall files from Allen Brain Institute (i.e.,
         as obtained by calling :func:`abagen.fetch_microarray` and accessing
         the `pacall` attribute on the resulting object).
+    samples : list of str or pandas.DataFrame
+
     probes : str or pandas.DataFrame
         Dataframe containing information on microarray probes that should be
         considered in filtering (probes not in this dataframe will be ignored)
@@ -89,16 +91,17 @@ def filter_probes(pacall, probes, threshold=0.5):
 
     threshold = np.clip(threshold, 0.0, 1.0)
 
-    probes = io.read_probes(probes, copy=True)
-    signal, samples = [], 0
-    for fname in pacall:
-        data = io.read_pacall(fname).loc[probes.index]
-        samples += data.shape[-1]
+    probes = io.read_probes(probes)
+    signal, n_samp = np.zeros(len(probes), dtype=int), 0
+    for pa, annot in zip(pacall, samples):
+        data = io.read_pacall(pa).loc[probes.index,
+                                      io.read_annotation(annot).index]
+        n_samp += data.shape[-1]
         # sum binary expression indicator across samples for current subject
-        signal.append(data.sum(axis=1).values)
+        signal += np.asarray(data.sum(axis=1))
 
     # calculate proportion of signal to noise for given probe across samples
-    keep = (np.sum(signal, axis=0) / samples) >= threshold
+    keep = (signal / n_samp) >= threshold
 
     return probes[keep]
 
@@ -187,7 +190,7 @@ def _max_loading(df):
     """
 
     if len(df) == 1:
-        return df.index[0]
+        return np.squeeze(df.index)
 
     data = np.asarray(df)
     data = data - data.mean(axis=0, keepdims=True)
@@ -389,11 +392,12 @@ def _collapse(expression, probes, *args, method='variance', **kwargs):
     # determine aggregation function based on provided method; also reduce
     # probe expression if required (i.e., max_variance, max_intensity)
     if method == 'max_variance':
-        probe_exp = pd.DataFrame(probe_exp.std(axis=1), columns=[method])
+        probe_exp = probe_exp.std(axis=1)
+        probe_exp.name = method
         agg = functools.partial(_max_idx, column=method)
     elif method == 'max_intensity':
-        probe_exp = pd.DataFrame(probe_exp.mean(axis=1), columns=[method])
-        # probe_exp.name = method
+        probe_exp = probe_exp.mean(axis=1)
+        probe_exp.name = method
         agg = functools.partial(_max_idx, column=method)
     elif method == 'pc_loading':
         agg = _max_loading
@@ -403,7 +407,8 @@ def _collapse(expression, probes, *args, method='variance', **kwargs):
         raise ValueError('Provided method {} invalid. Please check inputs '
                          'and try again.'.format(method))
 
-    info = pd.merge(probes[['gene_symbol']], probe_exp, on='probe_id')
+    info = pd.merge(probes['gene_symbol'], probe_exp, on='probe_id')
+
     return _groupby_and_apply(expression, probes, info, agg)
 
 
@@ -608,8 +613,9 @@ def collapse_probes(microarray, annotation, probes, method='diff_stability',
     # read in microarray data for all subjects; this can be quite slow...
     probes = io.read_probes(probes)
     exp = [
-        io.read_microarray(m, copy=not inplace).loc[probes.index]
-        for m in microarray
+        (io.read_microarray(m, copy=not inplace)
+           .loc[probes.index, io.read_annotation(a, copy=not inplace).index])
+        for m, a in zip(microarray, annotation)
     ]
 
     return [e.T for e in collfunc(exp, probes, annotation)]
diff --git a/abagen/samples.py b/abagen/samples.py
index a3b19d0..ce7764b 100644
--- a/abagen/samples.py
+++ b/abagen/samples.py
@@ -141,8 +141,6 @@ def drop_mismatch_samples(annotation, ontology):
     hemisphere designation in `ontology`; samples who do not match (e.g., L
     hemisphere designation in `ontology` but X coordinate < 0) are removed.
 
-    Optionally updates MNI coordinates in `annotation` (see `corrected`).
-
     Parameters
     ----------
     annotation : str
@@ -153,10 +151,6 @@ def drop_mismatch_samples(annotation, ontology):
         Ontology file from Allen Brain Institute. Optimally obtained by
         calling `abagen.fetch_microarray()` and accessing the `ontology`
         attribute on the resulting object
-    corrected : bool, optional
-        Whether to use the "corrected" MNI coordinates shipped with the
-        `alleninf` package instead of the coordinates provided with the AHBA
-        data. Default: True
 
     Returns
     -------
@@ -165,25 +159,168 @@ def drop_mismatch_samples(annotation, ontology):
     """
 
     # read in data files
-    annot = io.read_annotation(annotation)
-    ont = io.read_ontology(ontology).set_index('id')
+    annotation = io.read_annotation(annotation)
+    ontology = io.read_ontology(ontology).set_index('id')
+    sid = np.asarray(annotation['structure_id'])
 
     # get hemisphere and structure path
-    hemisphere = np.asarray(ont.loc[annot['structure_id'], 'hemisphere'])
-    structure = np.asarray(ont.loc[annot['structure_id'], 'structure_id_path']
-                              .apply(_get_struct))
+    hemisphere = np.asarray(ontology.loc[sid, 'hemisphere'])
+    structure = np.asarray(ontology.loc[sid, 'structure_id_path']
+                                   .apply(_get_struct))
 
     # add hemisphere + brain "structure" designation to annotation data and
     # only keep samples with consistent hemisphere + MNI coordinate designation
-    annot = annot.assign(hemisphere=hemisphere, structure=structure) \
-                 .query('(hemisphere == "L" & mni_x < 0) '
-                        '| (hemisphere == "R" & mni_x > 0) '
-                        '| (hemisphere.isna() & mni_x == 0)',
-                        engine='python')
+    annot = annotation.assign(hemisphere=hemisphere, structure=structure) \
+                      .query('(hemisphere == "L" & mni_x < 0) '
+                             '| (hemisphere == "R" & mni_x > 0) '
+                             '| (hemisphere.isna() & mni_x == 0)',
+                             engine='python')
 
     return annot
 
 
+def _mirror_ontology(annotation, ontology):
+    """
+    Assumes all hemisphere assignments of structures in `annotation` are wrong
+
+    Parameters
+    ----------
+    annotation,ontology : str
+        Filepath to {annotation,ontology} file
+
+    Returns
+    -------
+    annotation : pandas.DataFrame
+        Loaded annotation with updated structure information
+    """
+
+    HEMI_SWAP = dict(L='R', R='L')
+
+    annotation = io.read_annotation(annotation, copy=True)
+    ontology = io.read_ontology(ontology)
+
+    # structure IDs are specific to structure + hemisphere, so we can use this
+    # to grab the hemisphere designation for each sample and flip that L<->R
+    hemi = ontology.set_index('id') \
+                   .loc[annotation['structure_id'], 'hemisphere'] \
+                   .replace(HEMI_SWAP)
+
+    # structure acronyms are distinct to structure but not to hemisphere, so we
+    # can use this to grab all variations of a given structure and, with the
+    # flipped hemisphere designations we just generated, find the ID/name of
+    # the relevant structure in the new hemisphere
+    acr = ontology.set_index(['acronym', 'hemisphere']) \
+                  .loc[zip(annotation['structure_acronym'], hemi)] \
+                  .reset_index()
+
+    # update the original annotation with the new values
+    annotation['structure_acronym'] = np.asarray(acr['acronym'])
+    annotation['structure_id'] = np.asarray(acr['id'])
+    annotation['structure_name'] = np.asarray(acr['name'])
+
+    return annotation
+
+
+def mirror_samples(annotation, ontology):
+    """
+    Mirrors all tissue samples across hemispheres (L->R and R->L)
+
+    Methodology follows [SA1]_ in that samples are mirrored bi-directionally
+    rather than uni-directionally (R->L) as in [SA2]_.
+
+    Parameters
+    ----------
+    annotation : str or pandas.DataFrame
+        Filepath to annotation file from Allen Brain Institute (i.e., as
+        obtained by calling :func:`abagen.fetch_microarray` and accessing the
+        `annotation` attribute on the resulting object).
+    ontology : str or pandas.DataFrame
+        Filepath to ontology file from Allen Brain Institute (i.e., as
+        obtained by calling :func:`abagen.fetch_microarray` and accessing the
+        `ontology` attribute on the resulting object).
+
+    Returns
+    -------
+    annotation : pandas.DataFrame
+        Loaded input `annotation` data with all samples duplicated across
+        hemispheres (where structure IDs and MNI coordinates are updated,
+        accordingly)
+
+    References
+    ----------
+    .. [SA1] Gryglewski, G., Seiger, R., James, G. M., Godbersen, G. M.,
+       Komorowski, A., Unterholzner, J., ... & Kasper, S. (2018). Spatial
+       analysis and high resolution mapping of the human whole-brain
+       transcriptome for integrative analysis in neuroimaging. NeuroImage, 176,
+       259-267.
+
+    .. [SA2] Romero-Garcia, R., Whitaker, K. J., Váša, F., Seidlitz, J., Shinn,
+       M., Fonagy, P., ... & Vértes, P. E. (2018). Structural covariance
+       networks are coupled to expression of genes enriched in supragranular
+       layers of the human cortex. NeuroImage, 171, 256-267.
+    """
+
+    annotation = io.read_annotation(annotation)
+    ontology = io.read_ontology(ontology)
+
+    # take all lh and rh samples and flip x-coordinate
+    # also update ontology information (structure_id/acronym/name) as this is
+    # used when dropping mismatched samples later in the workflow
+    lh = _mirror_ontology(annotation[annotation['mni_x'] < 0], ontology)
+    rh = _mirror_ontology(annotation[annotation['mni_x'] > 0], ontology)
+
+    # update mni-x coordinate (simply flip the sign)
+    for df in [lh, rh]:
+        df['mni_x'] *= -1
+
+    # grow microarray and pacall based on duplicated samples
+    annotation = pd.concat([annotation, lh, rh])
+
+    return annotation
+
+
+def update_samples(annotation, ontology, lr_mirror=False, corrected_mni=True):
+    """
+    Updates samples in provided `annotation` file
+
+    Sample MNI coordinates are updated (if `corrected_mni`; see
+    :func:`~.update_mni_coords`), mismatched samples are dropped (where MNI
+    x-coordinate does not match hemisphere designation of structure in
+    `ontology`; see :func:`~.drop_mismatch_samples`), and samples are mirrored
+    across hemispheres (if `lr_mirror`; see :func:`mirror_samples`)
+
+    Parameters
+    ----------
+    annotation : str or pandas.DataFrame
+        List of filepaths to annotation files from Allen Brain Institute (i.e.,
+        as obtained by calling :func:`abagen.fetch_microarray` and accessing
+        the `annotation` attribute on the resulting object).
+    ontology : str or pandas.DataFrame
+        List of filepaths to ontology files from Allen Brain Institute (i.e.,
+        as obtained by calling :func:`abagen.fetch_microarray` and accessing
+        the `ontology` attribute on the resulting object).
+    lr_mirror : bool, optional
+    corrected_mni : bool, optional
+
+    Returns
+    -------
+    annotation : pandas.DataFrame
+    """
+
+    # update MNI coordinates first; this will affect all later stages
+    if corrected_mni:
+        annotation = update_mni_coords(annotation)
+
+    # drop samples where (new?) MNI coords don't match hemisphere in ontology
+    annotation = drop_mismatch_samples(annotation, ontology)
+
+    # optionally duplicate samples across hemisphere
+    if lr_mirror:
+        annotation = mirror_samples(annotation, ontology)
+
+    return annotation
+
+
 def _assign_sample(sample, atlas, sample_info=None, atlas_info=None,
                    tolerance=2):
     """
@@ -362,152 +499,3 @@ def label_samples(annotation, atlas, atlas_info=None, tolerance=2):
 
     return pd.DataFrame(labelled, dtype=int,
                         columns=['label'], index=annotation.index)
-
-
-def mirror_samples(microarray, pacall, annotation, ontology, inplace=False):
-    """
-    Mirrors all tissue samples across hemispheres (L->R and R->L)
-
-    Methodology follows [SA1]_ in that samples are mirrored bi-directionally
-    rather than uni-directionally (R->L) as in [SA2]_.
-
-    Parameters
-    ----------
-    microarray : list of str or pandas.DataFrame
-        List of filepaths to microarray expression files from Allen Brain
-        Institute (i.e., as obtained by calling :func:`abagen.fetch_microarray`
-        and accessing the `microarray` attribute on the resulting object).
-    pacall : list of str or pandas.DataFrame
-        List of filepaths to pacall files from Allen Brain Institute (i.e., as
-        obtained by calling :func:`abagen.fetch_microarray` and accessing the
-        `pacall` attribute on the resulting object).
-    annotation : list of str or pandas.DataFrame
-        List of filepaths to annotation files from Allen Brain Institute (i.e.,
-        as obtained by calling :func:`abagen.fetch_microarray` and accessing
-        the `annotation` attribute on the resulting object).
-    ontology : list of str or pandas.DataFrame
-        List of filepaths to ontology files from Allen Brain Institute (i.e.,
-        as obtained by calling :func:`abagen.fetch_microarray` and accessing
-        the `ontology` attribute on the resulting object).
-    inplace : bool, optional
-        Whether to conserve memory by editing dataframes in-place instead of
-        returning edited copies. Default: False
-
-    Returns
-    -------
-    microarray,pacall,annotation : list of pandas.DataFrame
-        Loaded input data with all samples duplicated across hemispheres
-
-    References
-    ----------
-    .. [SA1] Gryglewski, G., Seiger, R., James, G. M., Godbersen, G. M.,
-       Komorowski, A., Unterholzner, J., ... & Kasper, S. (2018). Spatial
-       analysis and high resolution mapping of the human whole-brain
-       transcriptome for integrative analysis in neuroimaging. NeuroImage, 176,
-       259-267.
-
-    .. [SA2] Romero-Garcia, R., Whitaker, K. J., Váša, F., Seidlitz, J., Shinn,
-       M., Fonagy, P., ... & Vértes, P. E. (2018). Structural covariance
-       networks are coupled to expression of genes enriched in supragranular
-       layers of the human cortex. NeuroImage, 171, 256-267.
-    """
-
-    # FIXME: seems like there should be a better way to do this, but both str
-    # and dataframes _are_ iterable just...not the way we want them to be here
-    inp = [microarray, pacall, annotation, ontology]
-    for n, i in enumerate(inp):
-        if isinstance(i, (str, pd.DataFrame)):
-            inp[n] = [i]
-
-    # flipped will be a list of len-3 tuples: (microarray, pacall, annotation)
-    flipped = [_mirror_samples(*i, inplace=inplace) for i in zip(*inp)]
-
-    # unpack so that even if user assigns output to 1 variable they get a tuple
-    microarray, pacall, annotation = [list(f) for f in zip(*flipped)]
-
-    return microarray, pacall, annotation
-
-
-def _mirror_samples(microarray, pacall, annotation, ontology, inplace=False):
-    """
-    Mirrors tissue samples across hemispheres for single donor
-
-    microarray,pacall,annotation,ontology : str or pandas.DataFrame
-        Filepath to {microarray,pacall,annotation,ontology} file
-    inplace : bool, optional
-        Whether to conserve memory by editing dataframes in-place instead of
-        returning edited copies. Default: False
-
-    Returns
-    -------
-    microarray,pacall,annotation : pandas.DataFrame
-        Loaded input data with all samples duplicated across hemispheres
-    """
-
-    microarray = io.read_microarray(microarray, copy=not inplace)
-    pacall = io.read_pacall(pacall, copy=not inplace)
-    annotation = io.read_annotation(annotation, copy=not inplace)
-    ontology = io.read_ontology(ontology)
-
-    # take all lh and rh samples and flip x-coordinate
-    # also update ontology information (structure_id/acronym/name) as this is
-    # used when dropping mismatched samples later in the workflow
-    lh = _mirror_ontology(annotation[annotation['mni_x'] < 0], ontology)
-    rh = _mirror_ontology(annotation[annotation['mni_x'] > 0], ontology)
-    for df in [lh, rh]:
-        df['mni_x'] *= -1
-
-    # grow microarray and pacall based on duplicated samples
-    annotation = pd.concat([annotation, lh, rh])
-    microarray = microarray.loc[:, annotation.index]
-    pacall = pacall.loc[:, annotation.index]
-
-    # now reset index / columns of all dataframes to be consecutive
-    sids = pd.Series(range(len(annotation)), name='sample_id')
-    microarray.columns = sids
-    pacall.columns = sids
-    annotation.index = sids
-
-    return microarray, pacall, annotation
-
-
-def _mirror_ontology(annotation, ontology):
-    """
-    Assumes all hemisphere assignments of structures in `annotation` are wrong
-
-    Parameters
-    ----------
-    annotation,ontology : str
-        Filepath to {annotation,ontology} file
-
-    Returns
-    -------
-    annotation : pandas.DataFrame
-        Loaded annotation with updated structure information
-    """
-
-    HEMI_SWAP = dict(L='R', R='L')
-
-    annotation = io.read_annotation(annotation, copy=True)
-    ontology = io.read_ontology(ontology)
-
-    # structure IDs are specific to structure + hemisphere, so we can use this
-    # to grab the hemisphere designation for each sample and flip that L<->R
-    hemi = ontology.set_index('id') \
-                   .loc[annotation['structure_id'], 'hemisphere'] \
-                   .replace(HEMI_SWAP)
-
-    # structure acronyms are distinct to structure but not to hemisphere, so we
-    # can use this to grab all variations of a given structure and, with the
-    # flipped hemisphere designations we just generated, find the ID/name of
-    # the relevant structure in the new hemisphere
-    acr = ontology.set_index(['acronym', 'hemisphere']) \
-                  .loc[zip(annotation['structure_acronym'], hemi)] \
-                  .reset_index()
-
-    # update the original annotation with the new values
-    annotation['structure_acronym'] = acr['acronym'].values
-    annotation['structure_id'] = acr['id'].values
-    annotation['structure_name'] = acr['name'].values
-
-    return annotation
diff --git a/docs/user_guide/expression.rst b/docs/user_guide/expression.rst
index 8df7e3b..9e70e16 100644
--- a/docs/user_guide/expression.rst
+++ b/docs/user_guide/expression.rst
@@ -82,15 +82,15 @@ expression data normalized and aggregated across donors:
     >>> print(expression)
     gene_symbol      A1BG  A1BG-AS1       A2M  ...       ZYX     ZZEF1      ZZZ3
     label                                      ...
-    1            0.584804  0.680186  0.407087  ...  0.628665  0.371389  0.488942
-    2            0.518708  0.647558  0.418932  ...  0.553781  0.290940  0.564853
+    1            0.584886  0.680115  0.456676  ...  0.628942  0.372817  0.486182
+    2            0.520235  0.648150  0.390960  ...  0.554345  0.294109  0.563199
     3                 NaN       NaN       NaN  ...       NaN       NaN       NaN
     ...               ...       ...       ...  ...       ...       ...       ...
-    81           0.495329  0.341752  0.508363  ...  0.266736  0.573451  0.490982
-    82           0.756134  0.391646  0.421420  ...  0.292576  0.379677  0.674880
-    83           0.428424  0.222963  0.660437  ...  0.411786  0.661593  0.439605
+    81           0.496183  0.342802  0.426420  ...  0.266740  0.574540  0.492474
+    82           0.756227  0.392104  0.453205  ...  0.291817  0.380244  0.673324
+    83           0.429359  0.224082  0.652792  ...  0.411622  0.661333  0.440884
     <BLANKLINE>
-    [83 rows x 15655 columns]
+    [83 rows x 15661 columns]
 
 
 By default the data are normalized using a scaled robust sigmoid function such
@@ -159,15 +159,15 @@ will return a dense matrix (at the expense of some anatomical precision):
     >>> print(exp_exact)
     gene_symbol      A1BG  A1BG-AS1       A2M  ...       ZYX     ZZEF1      ZZZ3
     label                                      ...
-    1            0.578008  0.688770  0.406520  ...  0.556295  0.374005  0.497330
-    2            0.510950  0.655648  0.420663  ...  0.481134  0.296428  0.569186
-    3            0.664877  0.818036  0.319912  ...  0.538502  0.147872  0.624574
+    1            0.578278  0.689091  0.460607  ...  0.556139  0.374723  0.494891
+    2            0.512639  0.656650  0.397684  ...  0.481831  0.298676  0.568115
+    3            0.670337  0.819654  0.262418  ...  0.540717  0.153573  0.627640
     ...               ...       ...       ...  ...       ...       ...       ...
-    81           0.486826  0.354870  0.517001  ...  0.251989  0.566561  0.498105
-    82           0.749263  0.412805  0.440311  ...  0.270990  0.376689  0.667931
-    83           0.416738  0.239668  0.661322  ...  0.389969  0.654765  0.448254
+    81           0.487488  0.355369  0.437644  ...  0.252065  0.567333  0.499693
+    82           0.749459  0.412666  0.468690  ...  0.270398  0.376335  0.667313
+    83           0.417314  0.240221  0.648598  ...  0.389826  0.654373  0.449491
     <BLANKLINE>
-    [83 rows x 15655 columns]
+    [83 rows x 15661 columns]
 
 
 .. _usage_expression_lrmirror:
@@ -197,15 +197,16 @@ but it will dramatically increase the likelihood that this will happen:
     >>> print(exp_mirror)
     gene_symbol      A1BG  A1BG-AS1       A2M  ...       ZYX     ZZEF1      ZZZ3
     label                                      ...
-    1            0.633077  0.680281  0.447366  ...  0.606646  0.378124  0.543360
-    2            0.578483  0.609816  0.478465  ...  0.597044  0.294332  0.540773
-    3            0.828558  0.643389  0.431069  ...  0.584006  0.436933  0.801877
+    1            0.626998  0.685507  0.451278  ...  0.623884  0.392830  0.547517
+    2            0.580470  0.623133  0.486455  ...  0.614684  0.301301  0.541552
+    3            0.819186  0.690949  0.486005  ...  0.643147  0.512646  0.833094
     ...               ...       ...       ...  ...       ...       ...       ...
-    81           0.478099  0.335455  0.427816  ...  0.262703  0.577574  0.492169
-    82           0.762843  0.409607  0.423332  ...  0.249643  0.389806  0.686732
-    83           0.417837  0.219333  0.676730  ...  0.404555  0.675958  0.458144
+    81           0.479263  0.349313  0.420825  ...  0.273995  0.575442  0.498685
+    82           0.766295  0.425379  0.430502  ...  0.268913  0.391012  0.687172
+    83           0.413256  0.220973  0.672984  ...  0.417636  0.669784  0.454275
     <BLANKLINE>
-    [83 rows x 15655 columns]
+    [83 rows x 15661 columns]
+
 
 
 Note that since this effectively duplicates the number of tissue samples the


Test patch:
diff --git a/abagen/tests/test_probes.py b/abagen/tests/test_probes.py
index a1bf33b..06756be 100644
--- a/abagen/tests/test_probes.py
+++ b/abagen/tests/test_probes.py
@@ -44,13 +44,15 @@ def test_filter_probes(testfiles, threshold, expected_length):
     # set up a few useful variables
     pacall = testfiles['pacall']
     probe_file = testfiles['probes'][0]
+    samples = testfiles['annotation']
     probe_df = abagen.io.read_probes(probe_file)
 
     # should work with either a filename _or_ a dataframe
-    filtered = probes.filter_probes(pacall, probe_file, threshold=threshold)
+    filtered = probes.filter_probes(pacall, samples, probe_file,
+                                    threshold=threshold)
     pd.testing.assert_frame_equal(
         filtered,
-        probes.filter_probes(pacall, probe_df, threshold=threshold)
+        probes.filter_probes(pacall, samples, probe_df, threshold=threshold)
     )
 
     # provided threshold returns expected output
diff --git a/abagen/tests/test_samples.py b/abagen/tests/test_samples.py
index d23508a..9b368aa 100644
--- a/abagen/tests/test_samples.py
+++ b/abagen/tests/test_samples.py
@@ -69,25 +69,6 @@ def annotation(mm_annotation):
     return out
 
 
-@pytest.fixture(scope='module')
-def microarray():
-    """ Fake microarray dataframe
-    """
-    data = np.arange(9).reshape(3, 3)
-    cols = pd.Series(range(3), name='sample_id')
-    ind = pd.Series([1058685, 1058684, 1058683], name='probe_id')
-    return pd.DataFrame(data, columns=cols, index=ind)
-
-
-@pytest.fixture(scope='module')
-def pacall(microarray):
-    """ Fake PACall dataframe
-    """
-    data = [[0, 0, 1], [0, 1, 1], [1, 1, 1]]
-    return pd.DataFrame(data, columns=microarray.columns,
-                        index=microarray.index)
-
-
 # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
 # test all the functions on our generated fake data so we know what to expect #
 # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
@@ -172,8 +153,7 @@ def test_label_samples():
     assert False
 
 
-def test_mirror_samples(microarray, pacall, annotation, ontology):
-    sids = pd.Series(range(5), name='sample_id')
+def test_mirror_samples(annotation, ontology):
     # we're changing quite a bit of stuff in the annotation dataframe
     aexp = pd.DataFrame(dict(mni_x=[-10, 30, 0, 10, -30],
                              structure_acronym=['S', 'Cl', 'CC', 'S', 'Cl'],
@@ -183,23 +163,11 @@ def test_mirror_samples(microarray, pacall, annotation, ontology):
                                              'central canal',
                                              'subiculum, right',
                                              'claustrum, left']),
-                        index=sids)
-    # and far less stuff in the microarray/pacall dataframes
-    mexp = pd.DataFrame(microarray.loc[:, [0, 1, 2, 0, 1]].values,
-                        columns=sids, index=microarray.index)
-    pexp = pd.DataFrame(pacall.loc[:, [0, 1, 2, 0, 1]].values,
-                        columns=sids, index=pacall.index)
+                        index=[0, 1, 2, 0, 1])
 
     # but let's confirm all the outputs are as-expected
-    m, p, a = samples.mirror_samples(microarray, pacall, annotation, ontology)
-    pd.testing.assert_frame_equal(a[0], aexp, check_like=True)
-    pd.testing.assert_frame_equal(m[0], mexp)
-    pd.testing.assert_frame_equal(p[0], pexp)
-
-    m, p, a = samples._mirror_samples(microarray, pacall, annotation, ontology)
+    a = samples.mirror_samples(annotation, ontology)
     pd.testing.assert_frame_equal(a, aexp, check_like=True)
-    pd.testing.assert_frame_equal(m, mexp)
-    pd.testing.assert_frame_equal(p, pexp)
 
 
 def test__mirror_ontology(annotation, ontology):
@@ -241,21 +209,12 @@ def test_drop_mismatch_samples_real(testfiles):
 
 
 def test_mirror_samples_real(testfiles):
-    testfiles = {k: v for k, v in testfiles.items() if k != 'probes'}
-
-    original = [363, 470]  # number of original probes
-    for m, p, a, o in zip(*samples.mirror_samples(**testfiles), original):
-        # all the outputs should have same number of samples
-        assert m.shape[-1] == p.shape[-1] == len(a)
-        # which should be more than the original # of samples but less than or
+    orig = [363, 470]
+    for an, on, o in zip(testfiles['annotation'], testfiles['ontology'], orig):
+        out = samples.mirror_samples(an, on)
+        # there should be more than the original # of samples but less than or
         # equal to 2x that number (we can't MORE than duplicate)
-        assert len(a) > o and len(a) <= o * 2
-
-    # we can also supply single filenames, too, and that should work w/o issue
-    samples.mirror_samples(**{k: v[0] for k, v in testfiles.items()})
-
-    # check that the sub-function works too
-    samples._mirror_samples(**{k: v[0] for k, v in testfiles.items()})
+        assert len(out) > o and len(out) <= o * 2
 
 
 def test__mirror_ontology_real(testfiles):