4.) refactoring allen smFish with new spot finding

notebooks/osmFISH.ipynb

@@ -161,15 +161,22 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "from starfish.spots import DetectSpots\n",
+    "from starfish.spots import DecodeSpots, FindSpots\n",
+    "from starfish.types import TraceBuildingStrategies\n",
     "\n",
-    "lmp = DetectSpots.LocalMaxPeakFinder(\n",
+    "\n",
+    "lmp = FindSpots.LocalMaxPeakFinder(\n",
     "    min_distance=6,\n",
     "    stringency=0,\n",
     "    min_obj_area=6,\n",
     "    max_obj_area=600,\n",
+    "    is_volume=True\n",
     ")\n",
-    "spot_intensities = lmp.run(mp)"
+    "spots = lmp.run(mp)\n",
+    "\n",
+    "decoder = DecodeSpots.PerRoundMaxChannel(codebook=experiment.codebook,\n",
+    "                                         trace_building_strategy=TraceBuildingStrategies.SEQUENTIAL)\n",
+    "decoded_intensities = decoder.run(spots=spots)"
    ]
   },
   {
@@ -240,11 +247,11 @@
    "outputs": [],
    "source": [
     "benchmark_spot_count = len(benchmark_peaks)\n",
-    "starfish_spot_count = len(spot_intensities)\n",
+    "starfish_spot_count = len(decoded_intensities)\n",
     "\n",
     "plt.figure(figsize=(10,10))\n",
     "plt.plot(benchmark_peaks.x, -benchmark_peaks.y, \"o\")\n",
-    "plt.plot(spot_intensities[Axes.X.value], -spot_intensities[Axes.Y.value], \"x\")\n",
+    "plt.plot(decoded_intensities[Axes.X.value], -decoded_intensities[Axes.Y.value], \"x\")\n",
     "\n",
     "plt.legend([\"Benchmark: {} spots\".format(benchmark_spot_count),\n",
     "            \"Starfish: {} spots\".format(starfish_spot_count)])\n",

notebooks/py/osmFISH.py

@@ -102,15 +102,22 @@
 # EPY: END markdown
 
 # EPY: START code
-from starfish.spots import DetectSpots
+from starfish.spots import DecodeSpots, FindSpots
+from starfish.types import TraceBuildingStrategies
 
-lmp = DetectSpots.LocalMaxPeakFinder(
+
+lmp = FindSpots.LocalMaxPeakFinder(
     min_distance=6,
     stringency=0,
     min_obj_area=6,
     max_obj_area=600,
+    is_volume=True
 )
-spot_intensities = lmp.run(mp)
+spots = lmp.run(mp)
+
+decoder = DecodeSpots.PerRoundMaxChannel(codebook=experiment.codebook,
+                                         trace_building_strategy=TraceBuildingStrategies.SEQUENTIAL)
+decoded_intensities = decoder.run(spots=spots)
 # EPY: END code
 
 # EPY: START markdown
@@ -162,11 +169,11 @@ def get_benchmark_peaks(loaded_results, redo_flag=False):
 
 # EPY: START code
 benchmark_spot_count = len(benchmark_peaks)
-starfish_spot_count = len(spot_intensities)
+starfish_spot_count = len(decoded_intensities)
 
 plt.figure(figsize=(10,10))
 plt.plot(benchmark_peaks.x, -benchmark_peaks.y, "o")
-plt.plot(spot_intensities[Axes.X.value], -spot_intensities[Axes.Y.value], "x")
+plt.plot(decoded_intensities[Axes.X.value], -decoded_intensities[Axes.Y.value], "x")
 
 plt.legend(["Benchmark: {} spots".format(benchmark_spot_count),
             "Starfish: {} spots".format(starfish_spot_count)])

notebooks/py/smFISH.py

@@ -30,7 +30,8 @@
 
 import starfish
 import starfish.data
-from starfish import FieldOfView, IntensityTable
+from starfish import FieldOfView, DecodedIntensityTable
+from starfish.types import TraceBuildingStrategies
 
 # equivalent to %gui qt
 ipython = get_ipython()
@@ -72,7 +73,7 @@
 # EPY: END markdown
 
 # EPY: START code
-tlmpf = starfish.spots.DetectSpots.TrackpyLocalMaxPeakFinder(
+tlmpf = starfish.spots.FindSpots.TrackpyLocalMaxPeakFinder(
     spot_diameter=5,  # must be odd integer
     min_mass=0.02,
     max_size=2,  # this is max radius
@@ -96,6 +97,7 @@
 import sys
 print = partial(print, file=sys.stderr)
 
+
 def processing_pipeline(
     experiment: starfish.Experiment,
     fov_name: str,
@@ -123,6 +125,11 @@ def processing_pipeline(
     print("Loading images...")
     images = enumerate(experiment[fov_name].get_images(FieldOfView.PRIMARY_IMAGES))
 
+    decoder = starfish.spots.DecodeSpots.PerRoundMaxChannel(
+        codebook=codebook,
+        trace_building_strategy=TraceBuildingStrategies.SEQUENTIAL
+    )
+
     for image_number, primary_image in images:
         print(f"Filtering image {image_number}...")
         filter_kwargs = dict(
@@ -140,15 +147,14 @@ def processing_pipeline(
         clip2.run(primary_image, **filter_kwargs)
 
         print("Calling spots...")
-        spot_attributes = tlmpf.run(primary_image)
-        all_intensities.append(spot_attributes)
+        spots = tlmpf.run(primary_image)
+        print("Decoding spots...")
+        decoded_intensities = decoder.run(spots)
+        all_intensities.append(decoded_intensities)
 
-    spot_attributes = IntensityTable.concatenate_intensity_tables(all_intensities)
-
-    print("Decoding spots...")
-    decoded = codebook.decode_per_round_max(spot_attributes)
+    decoded = DecodedIntensityTable.concatenate_intensity_tables(all_intensities)
     decoded = decoded[decoded["total_intensity"] > .025]
-    
+
     print("Processing complete.")
 
     return primary_image, decoded

notebooks/smFISH.ipynb

@@ -40,7 +40,8 @@
     "\n",
     "import starfish\n",
     "import starfish.data\n",
-    "from starfish import FieldOfView, IntensityTable\n",
+    "from starfish import FieldOfView, DecodedIntensityTable\n",
+    "from starfish.types import TraceBuildingStrategies\n",
     "\n",
     "# equivalent to %gui qt\n",
     "ipython = get_ipython()\n",
@@ -98,7 +99,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "tlmpf = starfish.spots.DetectSpots.TrackpyLocalMaxPeakFinder(\n",
+    "tlmpf = starfish.spots.FindSpots.TrackpyLocalMaxPeakFinder(\n",
     "    spot_diameter=5,  # must be odd integer\n",
     "    min_mass=0.02,\n",
     "    max_size=2,  # this is max radius\n",
@@ -130,6 +131,7 @@
     "import sys\n",
     "print = partial(print, file=sys.stderr)\n",
     "\n",
+    "\n",
     "def processing_pipeline(\n",
     "    experiment: starfish.Experiment,\n",
     "    fov_name: str,\n",
@@ -157,6 +159,11 @@
     "    print(\"Loading images...\")\n",
     "    images = enumerate(experiment[fov_name].get_images(FieldOfView.PRIMARY_IMAGES))\n",
     "\n",
+    "    decoder = starfish.spots.DecodeSpots.PerRoundMaxChannel(\n",
+    "        codebook=codebook,\n",
+    "        trace_building_strategy=TraceBuildingStrategies.SEQUENTIAL\n",
+    "    )\n",
+    "\n",
     "    for image_number, primary_image in images:\n",
     "        print(f\"Filtering image {image_number}...\")\n",
     "        filter_kwargs = dict(\n",
@@ -174,15 +181,14 @@
     "        clip2.run(primary_image, **filter_kwargs)\n",
     "\n",
     "        print(\"Calling spots...\")\n",
-    "        spot_attributes = tlmpf.run(primary_image)\n",
-    "        all_intensities.append(spot_attributes)\n",
+    "        spots = tlmpf.run(primary_image)\n",
+    "        print(\"Decoding spots...\")\n",
+    "        decoded_intensities = decoder.run(spots)\n",
+    "        all_intensities.append(decoded_intensities)\n",
     "\n",
-    "    spot_attributes = IntensityTable.concatenate_intensity_tables(all_intensities)\n",
-    "\n",
-    "    print(\"Decoding spots...\")\n",
-    "    decoded = codebook.decode_per_round_max(spot_attributes)\n",
+    "    decoded = DecodedIntensityTable.concatenate_intensity_tables(all_intensities)\n",
     "    decoded = decoded[decoded[\"total_intensity\"] > .025]\n",
-    "    \n",
+    "\n",
     "    print(\"Processing complete.\")\n",
     "\n",
     "    return primary_image, decoded"

starfish/core/spots/DecodeSpots/__init__.py

@@ -1,4 +1,5 @@
 from ._base import DecodeSpotsAlgorithm
+from .metric_decoder import MetricDistance
 from .per_round_max_channel_decoder import PerRoundMaxChannel
 
 # autodoc's automodule directive only captures the modules explicitly listed in __all__.

starfish/core/spots/DecodeSpots/metric_decoder.py

@@ -0,0 +1,92 @@
+from starfish.core.codebook.codebook import Codebook
+from starfish.core.intensity_table.decoded_intensity_table import DecodedIntensityTable
+from starfish.core.intensity_table.intensity_table_coordinates import \
+    transfer_physical_coords_to_intensity_table
+from starfish.core.spots.DecodeSpots.trace_builders import TRACE_BUILDERS
+from starfish.core.types import Number, SpotFindingResults, TraceBuildingStrategies
+from ._base import DecodeSpotsAlgorithm
+
+
+class MetricDistance(DecodeSpotsAlgorithm):
+    """
+    Normalizes both the magnitudes of the codes and the spot intensities, then decodes spots by
+    assigning each spot to the closest code, measured by the provided metric.
+
+    Codes greater than max_distance from the nearest code, or dimmer than min_intensity, are
+    discarded.
+
+    Parameters
+    ----------
+    codebook : Codebook
+        codebook containing targets the experiment was designed to quantify
+    max_distance : Number
+        spots greater than this distance from their nearest target are not decoded
+    min_intensity : Number
+        spots dimmer than this intensity are not decoded
+    metric : str
+        the metric to use to measure distance. Can be any metric that satisfies the triangle
+        inequality that is implemented by :py:mod:`scipy.spatial.distance` (default "euclidean")
+    norm_order : int
+        the norm to use to normalize the magnitudes of spots and codes (default 2, L2 norm)
+    trace_building_strategy: TraceBuildingStrategies
+        Defines the strategy for building spot traces to decode across rounds and chs of spot
+        finding results.
+    anchor_round : Optional[int]
+        Only applicable if trace_building_strategy is TraceBuildingStrategies.NEAREST_NEIGHBORS.
+        The imaging round against which other rounds will be checked for spots in the same
+        approximate pixel location.
+    search_radius : Optional[int]
+        Only applicable if trace_building_strategy is TraceBuildingStrategies.NEAREST_NEIGHBORS.
+        Number of pixels over which to search for spots in other rounds and channels.
+    """
+
+    def __init__(
+        self,
+        codebook: Codebook,
+        max_distance: Number,
+        min_intensity: Number,
+        norm_order: int = 2,
+        metric: str = "euclidean",
+        trace_building_strategy: TraceBuildingStrategies = TraceBuildingStrategies.EXACT_MATCH,
+        anchor_round: int = 1,
+        search_radius: int = 3,
+    ) -> None:
+        self.codebook = codebook
+        self.max_distance = max_distance
+        self.min_intensity = min_intensity
+        self.norm_order = norm_order
+        self.metric = metric
+        self.trace_builder = TRACE_BUILDERS[trace_building_strategy]
+        self.anchor_round = anchor_round
+        self.search_radius = search_radius
+
+    def run(
+        self,
+        spots: SpotFindingResults,
+        *args
+    ) -> DecodedIntensityTable:
+        """Decode spots by selecting the max-valued channel in each sequencing round
+
+        Parameters
+        ----------
+        spots : SpotFindingResults
+            A Dict of tile indices and their corresponding measured spots
+
+        Returns
+        -------
+        DecodedIntensityTable :
+            IntensityTable decoded and appended with Features.TARGET and Features.QUALITY values.
+
+        """
+
+        intensities = self.trace_builder(spot_results=spots,
+                                         anchor_round=self.anchor_round,
+                                         search_radius=self.search_radius)
+        transfer_physical_coords_to_intensity_table(intensity_table=intensities, spots=spots)
+        return self.codebook.decode_metric(
+            intensities,
+            max_distance=self.max_distance,
+            min_intensity=self.min_intensity,
+            norm_order=self.norm_order,
+            metric=self.metric,
+        )

starfish/core/spots/DecodeSpots/trace_builders.py

@@ -1,7 +1,10 @@
 from typing import Callable, Mapping
 
+import pandas as pd
+
 from starfish.core.intensity_table.intensity_table import IntensityTable
-from starfish.core.types import Axes, Features, SpotFindingResults, TraceBuildingStrategies
+from starfish.core.types import Axes, Features, SpotAttributes, SpotFindingResults, \
+    TraceBuildingStrategies
 from .util import _build_intensity_table, _match_spots, _merge_spots_by_round
 
 
@@ -30,6 +33,40 @@ def build_spot_traces_exact_match(spot_results: SpotFindingResults, **kwargs) ->
     return intensity_table
 
 
+def build_traces_sequential(spot_results: SpotFindingResults, **kwargs) -> IntensityTable:
+    """
+    Build spot traces  without merging across channels and imaging rounds. Used for sequential
+    methods like smFIsh.
+
+    Parameters
+    ----------
+    spot_results: SpotFindingResults
+        Spots found across rounds/channels of an ImageStack
+
+    Returns
+    -------
+    IntensityTable :
+        concatenated input SpotAttributes, converted to an IntensityTable object
+
+    """
+
+    all_spots = pd.concat([sa.data for sa in spot_results.values()], sort=True)
+
+    intensity_table = IntensityTable.zeros(
+        spot_attributes=SpotAttributes(all_spots),
+        ch_labels=spot_results.ch_labels,
+        round_labels=spot_results.round_labels,
+    )
+
+    i = 0
+    for (r, c), attrs in spot_results.items():
+        for _, row in attrs.data.iterrows():
+            selector = dict(features=i, c=c, r=r)
+            intensity_table.loc[selector] = row[Features.INTENSITY]
+            i += 1
+    return intensity_table
+
+
 def build_traces_nearest_neighbors(spot_results: SpotFindingResults, anchor_round: int=1,
                                    search_radius: int=3):
     """
@@ -65,5 +102,6 @@ def build_traces_nearest_neighbors(spot_results: SpotFindingResults, anchor_roun
 
 TRACE_BUILDERS: Mapping[TraceBuildingStrategies, Callable] = {
     TraceBuildingStrategies.EXACT_MATCH: build_spot_traces_exact_match,
-    TraceBuildingStrategies.NEAREST_NEIGHBOR: build_traces_nearest_neighbors
+    TraceBuildingStrategies.NEAREST_NEIGHBOR: build_traces_nearest_neighbors,
+    TraceBuildingStrategies.SEQUENTIAL: build_traces_sequential,
 }

starfish/core/spots/DetectSpots/trackpy_local_max_peak_finder.py

@@ -163,6 +163,12 @@ def run(
         n_processes : Optional[int] = None,
             Number of processes to devote to spot finding.
         """
+        DeprecationWarning("Starfish is embarking on a SpotFinding data structures refactor"
+                           "(See https://github.com/spacetx/starfish/issues/1514) This version of "
+                           "TrackpyLocalMaxPeakFinder will soon be deleted. To find and decode your"
+                           "spots please instead use FindSpots.TrackpyLocalMaxPeakFinder then "
+                           "DecodeSpots.PerRoundMaxChannel with the parameter "
+                           "TraceBuildingStrategies.SEQUENTIAL. See example in smFISH.py")
         intensity_table = detect_spots(
             data_stack=primary_image,
             spot_finding_method=self.image_to_spots,

starfish/core/spots/FindSpots/__init__.py

@@ -1,5 +1,7 @@
 from ._base import FindSpotsAlgorithm
 from .blob import BlobDetector
+from .local_max_peak_finder import LocalMaxPeakFinder
+from .trackpy_local_max_peak_finder import TrackpyLocalMaxPeakFinder
 
 # autodoc's automodule directive only captures the modules explicitly listed in __all__.
 all_filters = {