Skip to content
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
86 changes: 86 additions & 0 deletions CLAUDE.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,86 @@
# CLAUDE.md

This file provides guidance to Claude Code (claude.ai/code) when working with code in this repository.

## Project Overview

SMLMFrameConnection.jl performs frame-connection on 2D localization microscopy data, combining repeated localizations of blinking fluorophores into higher-precision localizations. Part of the JuliaSMLM ecosystem.

**Reference:** Schodt & Lidke (2021), "Spatiotemporal Clustering of Repeated Super-Resolution Localizations via Linear Assignment Problem"

## Commands

```bash
# Run tests
julia --project -e 'using Pkg; Pkg.test()'

# Run a specific test file
julia --project test/test_frameconnect.jl

# Interactive development
julia --project
using SMLMFrameConnection, SMLMData
```

## Architecture

### Algorithm Pipeline

The main entry point `frameconnect()` executes four stages:

1. **Preclustering** (`precluster.jl`): Groups spatiotemporally nearby localizations using KD-tree nearest-neighbor search. Distance threshold: `nsigmadev × mean(σ_x, σ_y)`. Temporal threshold: `maxframegap` frames.

2. **Parameter Estimation** (`estimateparams.jl`, `estimatedensities.jl`): Estimates photophysics (k_on, k_off, k_bleach, p_miss) and emitter densities from precluster statistics.

3. **LAP Solving** (`create_costmatrix.jl`, `solveLAP.jl`, `linkclusters.jl`): Constructs cost matrices for each precluster and solves linear assignment problems via Hungarian algorithm to determine which localizations belong to the same physical emitter.

4. **Combination** (`combinelocalizations.jl`): Merges connected localizations using MLE weighted mean: `x = Σ(x/σ²) / Σ(1/σ²)`, uncertainty: `σ = √(1/Σ(1/σ²))`.

### Key Types

- **Input**: `BasicSMLD{T, Emitter2DFit{T}}` from SMLMData.jl
- **Output**: `(combined::BasicSMLD, info::ConnectInfo)` tuple
- **ConnectInfo{T}** (`connectinfo.jl`): Secondary output with connected SMLD, statistics, and photophysics
- **ParamStruct** (`structdefinitions.jl`): Internal algorithm parameters

### File Organization

```
src/
├── frameconnect.jl # Main entry point
├── connectinfo.jl # ConnectInfo output struct
├── structdefinitions.jl # ParamStruct, helper conversions
├── precluster.jl # Stage 1: spatiotemporal clustering
├── organizeclusters.jl # Prepare cluster data structures
├── estimateparams.jl # Stage 2: photophysics estimation
├── estimatedensities.jl # Stage 2: density estimation
├── create_costmatrix.jl # Stage 3: LAP cost matrix construction
├── solveLAP.jl # Stage 3: Hungarian algorithm wrapper
├── linkclusters.jl # Stage 3: apply LAP solution
├── combinelocalizations.jl # Stage 4: MLE combination with covariance
├── defineidealFC.jl # Ground-truth validation utility
└── compress_connectID.jl # Normalize track_id to 1:N
```

### Dependencies

- **SMLMData.jl** (v0.6+): Core types (BasicSMLD, Emitter2DFit with σ_xy)
- **Hungarian.jl**: LAP solver
- **NearestNeighbors.jl**: KD-tree for spatial queries
- **Optim.jl**: Parameter optimization

## Testing

Tests use helper fixtures in `test/test_helpers.jl`:
- `make_test_camera()`, `make_emitter()`, `make_test_smld()`
- `make_blinking_molecule()`: Creates emitter sequence for one molecule
- `make_two_molecules_smld()`, `make_single_molecule_smld()`

## API

Five exports:
- `frameconnect(smld)` → `(combined::BasicSMLD, info::ConnectInfo)` tuple
- `ConnectInfo{T}` → struct with connected SMLD, statistics (n_input, n_tracks, n_combined), photophysics (k_on, k_off, k_bleach, p_miss), timing
- `combinelocalizations(smld)` → combines emitters by `track_id` using MLE with full covariance
- `defineidealFC(smld)` → ground-truth frame connection for simulated data
- `ParamStruct` → internal algorithm parameters
4 changes: 2 additions & 2 deletions Project.toml
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "SMLMFrameConnection"
uuid = "1517b1a9-81e8-461f-b994-92eb29599690"
authors = ["klidke@unm.edu"]
version = "0.2.0"
version = "0.3.0"

[deps]
Hungarian = "e91730f6-4275-51fb-a7a0-7064cfbd3b39"
Expand All @@ -15,7 +15,7 @@ StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
Hungarian = "0.6, 0.7"
NearestNeighbors = "0.4"
Optim = "1, 2"
SMLMData = "0.5"
SMLMData = "0.6"
StatsBase = "0.33, 0.34"
julia = "1.6"

Expand Down
75 changes: 38 additions & 37 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -24,9 +24,10 @@ Pkg.add("SMLMFrameConnection")
using SMLMData, SMLMFrameConnection

# Run frame connection on your data
smld_connected, smld_preclustered, smld_combined, params = frameconnect(smld)
(combined, info) = frameconnect(smld)

# smld_combined is the main output - higher precision localizations
# combined is the main output - higher precision localizations
# info contains connected SMLD, statistics, and estimated photophysics
```

## Input Requirements
Expand All @@ -52,41 +53,53 @@ using SMLMData, SMLMFrameConnection
cam = IdealCamera(1:512, 1:512, 0.1)

# Create emitters representing the same molecule blinking across 3 frames
# Constructor: Emitter2DFit{T}(x, y, photons, bg, σ_x, σ_y, σ_photons, σ_bg; frame, dataset, track_id)
# Constructor: Emitter2DFit{T}(x, y, photons, bg, σ_x, σ_y, σ_xy, σ_photons, σ_bg, frame, dataset, track_id, id)
emitters = [
Emitter2DFit{Float64}(
5.0, 5.0, # x, y position (μm)
1000.0, 10.0, # photons, background
0.02, 0.02, # σ_x, σ_y uncertainties (μm)
50.0, 2.0; # σ_photons, σ_bg
frame=1
),
Emitter2DFit{Float64}(5.01, 5.01, 1200.0, 12.0, 0.02, 0.02, 60.0, 2.0; frame=2),
Emitter2DFit{Float64}(5.02, 4.99, 1100.0, 11.0, 0.02, 0.02, 55.0, 2.0; frame=3),
Emitter2DFit{Float64}(5.0, 5.0, 1000.0, 10.0, 0.02, 0.02, 0.0, 50.0, 2.0, 1, 1, 0, 1),
Emitter2DFit{Float64}(5.01, 5.01, 1200.0, 12.0, 0.02, 0.02, 0.0, 60.0, 2.0, 2, 1, 0, 2),
Emitter2DFit{Float64}(5.02, 4.99, 1100.0, 11.0, 0.02, 0.02, 0.0, 55.0, 2.0, 3, 1, 0, 3),
]

# Create SMLD: BasicSMLD(emitters, camera, n_frames, n_datasets)
smld = BasicSMLD(emitters, cam, 3, 1)

# Run frame connection
_, _, smld_combined, _ = frameconnect(smld)
(combined, info) = frameconnect(smld)

# Result: localizations connected based on spatial/temporal proximity
println("Combined $(info.n_input) localizations into $(info.n_combined)")
println("Formed $(info.n_tracks) tracks from $(info.n_preclusters) preclusters")

# Combined uncertainty: σ_combined ≈ σ_individual / √n_connected
```

## Outputs Explained

```julia
smld_connected, smld_preclustered, smld_combined, params = frameconnect(smld)
(combined, info) = frameconnect(smld)
```

| Output | Description | When to use |
|--------|-------------|-------------|
| `smld_combined` | **Main output.** Combined high-precision localizations | Standard analysis |
| `smld_connected` | Original localizations with `track_id` populated | When you need per-frame data with connection labels |
| `smld_preclustered` | Intermediate preclustering result | Debugging, algorithm tuning |
| `params` | Estimated photophysics + input parameters | Inspecting estimated k_on, k_off, density |
| Output | Type | Description |
|--------|------|-------------|
| `combined` | `BasicSMLD` | **Main output.** Combined high-precision localizations |
| `info` | `ConnectInfo` | Metadata: connected SMLD, statistics, estimated photophysics |

### ConnectInfo Fields

| Field | Type | Description |
|-------|------|-------------|
| `info.connected` | `BasicSMLD` | Original localizations with `track_id` populated |
| `info.n_input` | `Int` | Number of input localizations |
| `info.n_tracks` | `Int` | Number of tracks formed |
| `info.n_combined` | `Int` | Number of output localizations |
| `info.n_preclusters` | `Int` | Number of preclusters |
| `info.k_on` | `Float64` | Estimated on rate (1/frame) |
| `info.k_off` | `Float64` | Estimated off rate (1/frame) |
| `info.k_bleach` | `Float64` | Estimated bleach rate (1/frame) |
| `info.p_miss` | `Float64` | Estimated miss probability |
| `info.initialdensity` | `Vector{Float64}` | Density estimate per cluster (emitters/μm²) |
| `info.elapsed_ns` | `UInt64` | Processing time in nanoseconds |
| `info.algorithm` | `Symbol` | Algorithm used (`:lap`) |

## Parameters

Expand All @@ -104,36 +117,24 @@ frameconnect(smld;
- `maxframegap`: Set based on expected blinking duration. For dSTORM with long dark states, increase to 10-20.
- Defaults work well for standard dSTORM/PALM data with typical blinking kinetics.

## Estimated Parameters (ParamStruct)

The algorithm estimates fluorophore photophysics from your data:

| Field | Description |
|-------|-------------|
| `k_on` | Rate of transitioning from dark to visible state (1/frame) |
| `k_off` | Rate of transitioning from visible to dark state (1/frame) |
| `k_bleach` | Photobleaching rate (1/frame) |
| `p_miss` | Probability of missing a localization when fluorophore is on |
| `initialdensity` | Estimated emitter density per cluster (emitters/μm²) |

## Combination Method

Connected localizations are combined using **maximum likelihood estimation (MLE) weighted mean**:
- Position: inverse-variance weighted average → `x_combined = Σ(x/σ²) / Σ(1/σ²)`
- Uncertainty: `σ_combined = √(1/Σ(1/σ²))` ≈ `σ_individual / √n`
Connected localizations are combined using **maximum likelihood estimation (MLE) weighted mean** with full covariance propagation:
- Position: inverse-variance weighted average using 2x2 covariance matrix
- Uncertainty: properly propagated including σ_xy correlation
- Photons: summed across connected localizations

## Utility Functions

### combinelocalizations
```julia
smld_combined = combinelocalizations(smld)
combined = combinelocalizations(smld)
```
Combines emitters that share the same `track_id`. Use when you have pre-labeled data.

### defineidealFC
```julia
smld_connected, smld_combined = defineidealFC(smld; maxframegap=5)
(connected, combined) = defineidealFC(smld; maxframegap=5)
```
For **simulated data** where `track_id` already contains ground-truth emitter IDs. Useful for validating frame-connection performance against known truth.

Expand Down
71 changes: 56 additions & 15 deletions api_overview.md
Original file line number Diff line number Diff line change
Expand Up @@ -2,16 +2,20 @@

Frame-connection for 2D SMLM data: combines repeated localizations of blinking fluorophores into higher-precision localizations.

## Exports

4 exports: `frameconnect`, `combinelocalizations`, `defineidealFC`, `ConnectInfo`, `ParamStruct`

## Main Functions

### frameconnect
```julia
frameconnect(smld::BasicSMLD{T,Emitter2DFit{T}};
(combined, info) = frameconnect(smld::BasicSMLD{T,Emitter2DFit{T}};
nnearestclusters::Int=2,
nsigmadev::Float64=5.0,
maxframegap::Int=5,
nmaxnn::Int=2
) where T -> (smld_connected, smld_preclustered, smld_combined, params)
) where T
```
Main entry point. Connects repeated localizations and combines them.

Expand All @@ -21,28 +25,44 @@ Main entry point. Connects repeated localizations and combines them.
- `maxframegap`: Maximum frame gap for temporal adjacency in preclusters
- `nmaxnn`: Maximum nearest-neighbors inspected for precluster membership

**Returns:**
- `smld_connected`: Input with `track_id` populated (localizations uncombined)
- `smld_preclustered`: Intermediate preclustering result (for debugging)
- `smld_combined`: **Main output** - combined high-precision localizations
- `params::ParamStruct`: Algorithm parameters (input + estimated photophysics)
**Returns tuple:**
- `combined::BasicSMLD`: **Main output** - combined high-precision localizations
- `info::ConnectInfo`: Metadata including connected SMLD, statistics, and estimated photophysics

### combinelocalizations
```julia
combinelocalizations(smld::BasicSMLD{T,Emitter2DFit{T}}) where T -> BasicSMLD
combined = combinelocalizations(smld::BasicSMLD{T,Emitter2DFit{T}}) where T -> BasicSMLD
```
Combines emitters sharing the same `track_id` using MLE weighted mean. Use when `track_id` is already populated.
Combines emitters sharing the same `track_id` using MLE weighted mean with full covariance propagation. Use when `track_id` is already populated.

### defineidealFC
```julia
defineidealFC(smld::BasicSMLD{T,Emitter2DFit{T}};
(connected, combined) = defineidealFC(smld::BasicSMLD{T,Emitter2DFit{T}};
maxframegap::Int=5
) where T -> (smld_connected, smld_combined)
) where T
```
For simulated data where `track_id` indicates ground-truth emitter ID. Useful for validation/benchmarking.

## Types

### ConnectInfo{T}
```julia
struct ConnectInfo{T}
connected::BasicSMLD{T} # Input with track_id assigned (uncombined)
n_input::Int # Number of input localizations
n_tracks::Int # Number of tracks formed
n_combined::Int # Number of output localizations
k_on::Float64 # Estimated on rate (1/frame)
k_off::Float64 # Estimated off rate (1/frame)
k_bleach::Float64 # Estimated bleach rate (1/frame)
p_miss::Float64 # Probability of missed detection
initialdensity::Vector{Float64} # Density estimate per cluster (emitters/μm²)
elapsed_ns::UInt64 # Processing time in nanoseconds
algorithm::Symbol # Algorithm used (:lap)
n_preclusters::Int # Number of preclusters formed
end
```

### ParamStruct
```julia
mutable struct ParamStruct
Expand Down Expand Up @@ -76,20 +96,41 @@ Input `BasicSMLD` must contain `Emitter2DFit` emitters.
**Optional fields:**
- `photons`, `σ_photons`: Photon count (summed in output)
- `bg`, `σ_bg`: Background
- `σ_xy`: Position covariance (microns², propagated in output)
- `dataset`: Dataset identifier (default: 1)
- `track_id`: Set to 0 for input (populated by algorithm)

## Output

Output `smld_combined` contains `Emitter2DFit` emitters with:
- Combined position via MLE weighted mean: `x = Σ(x/σ²) / Σ(1/σ²)`
- Reduced uncertainties: `σ = √(1/Σ(1/σ²))`
Output `combined` contains `Emitter2DFit` emitters with:
- Combined position via MLE weighted mean with full 2x2 covariance
- Properly propagated uncertainties including σ_xy
- Summed photons
- `track_id` indicating connection group

## Example

```julia
using SMLMData, SMLMFrameConnection

# Load or create SMLD
smld = ...

# Run frame connection
(combined, info) = frameconnect(smld; maxframegap=10)

# Access results
println("Combined $(info.n_input) → $(info.n_combined) localizations")
println("Estimated k_on=$(info.k_on), k_off=$(info.k_off)")
println("Processing time: $(info.elapsed_ns / 1e9) seconds")

# Access connected (uncombined) localizations
connected_smld = info.connected
```

## Dependencies

- SMLMData.jl (v0.5+): BasicSMLD, Emitter2DFit types
- SMLMData.jl (v0.6+): BasicSMLD, Emitter2DFit types
- Hungarian.jl (v0.6-0.7): Linear assignment problem solver
- NearestNeighbors.jl: Spatial clustering
- Optim.jl: Parameter estimation
Expand Down
Loading
Loading