Julian/JuliaMSI
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Improved documentation, added precautions to te json generation and loading specially for windows

Browse Source
This commit is contained in:
Pixelguy14 2026-03-03 14:36:39 -06:00
parent df94d50954
commit 3da0dcaae9
15 changed files with 1500 additions and 85 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

217
JuliaMSI_MainGuide.md Normal file
View File

@ -0,0 +1,217 @@
# JuliaMSI Platform User Guide
## Core Workflow and Interactive Analysis
---
### Table of Contents
1. [Introduction](#introduction)
2. [Loading a Dataset](#loading-a-dataset)
3. [Converting mzML to imzML](#converting-mzml-to-imzml)
4. [Batch Processing](#batch-processing)
5. [Spectral Visualization](#spectral-visualization)
6. [Imaging Slice Generation](#imaging-slice-generation)
7. [Topography and Surface Plots](#topography-and-surface-plots)
8. [Small Features and Tools](#small-features-and-tools)
* [Metadata Viewing](#metadata-viewing)
* [Optical Image Overlay](#optical-image-overlay)
* [Comparison of Datasets](#comparison-of-datasets)
* [Registry and Safe Processing](#registry-and-safe-processing)
9. [Mask Generation and Usage (mask.jl)](#mask-generation-and-usage-maskjl)
10. [Preprocessing](#preprocessing)
11. [Practical Tips for Analysts](#practical-tips-for-analysts)
* [Julia Version Compatibility](#julia-version-compatibility)
* [First-Time Startup](#first-time-startup)
* [Accessing the Web Interface](#accessing-the-web-interface)
* [Memory Management](#memory-management)
* [Registry and File Locations](#registry-and-file-locations)
* [Mask Handling](#mask-handling)
12. [Troubleshooting Common Issues](#troubleshooting-common-issues)
13. [Citation and Contact](#citation-and-contact)
---
## Introduction
This document serves as an entry point for analysts using the **JuliaMSI** platform. It focuses on the user-facing modules: `app.jl` (the main application entry), `julia_imzml_visual.jl` (visualization and heavy procedures), and `mask.jl` (mask generation and filtering). The goal is to explain why each step requires specific inputs and what those inputs mean scientifically.
**JuliaMSI** is a user-friendly Graphical User Interface (GUI) developed in Julia for MSI data analyses. Our framework allows for rapid **imzML** MSI data analysis and **mzML** spectra plot loading, significantly improving analytical workflows. High-resolution instruments produce massive datasets, requiring efficient memory management and parallel processing—reasons why JuliaMSI is written in Julia, which offers high performance.
### Core Objectives:
* **Simplification**: Reducing command-line dependency through intuitive visualizations and easy parameter input.
* **Accessibility**: Breaking the proprietary software barrier across fields like plant sciences, biomedical research, and environmental studies.
* **Interactive Exploration**: Using the `Plotly.jl` package for dynamic visualizations. Users can explore multiple interpretations of images, select specific coordinates for spectra analysis, or pick an *m/z* value directly from the displayed spectrum.
### Key Capabilities:
* **Smart data handling**: Processes large files in manageable sections without overwhelming memory.
* **Batch processing**: Analyzes multiple samples simultaneously for high-throughput studies.
* **Interactive exploration**: Web interface enables real-time data visualization and analysis.
---
## Loading a Dataset
**Action**: Click the **Search** button (`@onbutton btnSearch`)
Load an **imzML** or **mzML** file into the workspace, extract metadata, and prepare the dataset for visualization and processing.
The platform checks whether a pre-processed version of the same file exists in the registry (**fast load**). If not, it performs a **full load** using `OpenMSIData`. The GUI handles multiple file uploads for MSI slice generation; however, for spectra and preprocessing, it will only use the last selected file.
---
## Converting mzML to imzML
**Action**: Click **Convert/Process** (`@onbutton convert_process`)
This tool converts an **mzML** file (from a non-imaging mass spectrometry run) into an **imzML** file suitable for imaging analysis. This is necessary when spectra were acquired in a raster pattern but saved as a continuous **mzML** file without spatial coordinates.
### Inputs Required:
* **mzML file**: The raw mass spectrometry data file.
* **Sync file**: A text file that maps each spectrum in the **mzML** to a specific pixel coordinate (*x, y*). The sync file typically contains one line per spectrum with the format `x y` or similar.
> **Why these inputs?**: **mzML** files lack spatial metadata; they contain a list of spectra in acquisition order. Without coordinate information, the data cannot be reconstructed as an image. **Sync files** provide the missing spatial mapping, often generated by instrument software from a motor stage log. The conversion process embeds these coordinates into the **imzML** structure.
---
## Batch Processing
**Action**: Click **Main Process** (`@onbutton mainProcess`)
This is the core image generation routine. It takes a list of comma-separated *m/z* values and produces ion images (and optionally **TrIQ Quantized** images) for each file loaded.
### Required Inputs:
* **m/z values**: A comma-separated list of *m/z* values to visualize spatial distribution. These are the core of MSI: each image maps the abundance of a particular ion across the sample. Without a target, there is no image.
* **Tolerance**: The mass window around each target ion. Necessary because no mass spectrometer measures masses with infinite precision. The signal from a given ion is spread over a small *m/z* range; integrating over that range (by summing intensities of all data points within ±tolerance) captures the true signal.
* **Color level**: Number of discrete color levels in the output BMP image. A low number (e.g., 16) gives a coarse view; a high number (256) shows subtle intensity variations but can be noisy.
* **TrIQ**: A toggle to enable **Threshold Intensity Quantization**. This produces images with enhanced contrast and reduced influence of extreme outliers, making spatial patterns more visible. It is a contrast normalization and quantization technique.
* **TrIQ Threshold**: Probability threshold for outlier removal in TrIQ.
* **Median filter**: Replaces each pixel’s intensity with the median of its neighbors, effectively removing isolated "hot pixels" while preserving sharp boundaries.
* **Mask**: Toggle to focus on biologically relevant regions (e.g., tissue vs. background), avoiding wasted computational resources on irrelevant pixels.
> **Note**: The mask path is determined inside the JSON file, which can be manually altered, or using the specific section in `mask.jl`.
---
## Imaging Slice Generation (Part of Main Process)
In MSI, each pixel contains a full mass spectrum. To create an ion image (or "slice"), we must decide which data points in that spectrum belong to the ion of interest.
* **Tolerance**: Defines the *m/z* window. Too narrow may miss part of the peak; too wide includes noise or interfering peaks. Optimal tolerance depends on instrument resolving power and calibration stability.
* **Colorbar Generation**: Each image is accompanied by a colorbar PNG showing the mapping between intensity and colors. The bounds are determined during quantization (TrIQ or min-max scaling) and passed to `generate_colorbar_image`. This ensures consistent color interpretation with labels indicating the intensity range and scaling factors (e.g., "Intensity ×10³").
---
## Spectral Visualization
**Functions**: `createMeanPlot`, `createSumPlot`, `createXYPlot`, `createNSpectrumPlot`
Displays spectra from the loaded dataset, optionally filtered by a mask. This generates a Plotly layout plot for quality control, identification, and *m/z* selection.
### Inputs & Features:
* **Dataset**: The last **imzML** / **mzML** file loaded in the GUI.
* **Mask**: Toggle for mask filtering. The spectrum is computed only from pixels inside the mask.
* **User click**: If a user clicks on an already generated **slice** (2D image), two vertical lines appear and values are filled (one for *x*-axis, one for *y*-axis) to generate spectrum plot indicating the position of the click in the image. This allows the creation of a spectrum plot correlating to specific spatial selection.
* **NSpectrumPlot**: Allows the user to select specific spectra based on the **ID** (the index of the spectra in the dataset).
### Spectral Types:
* **Mean spectrum**: The average of all pixel spectra. Gives an overview of ions present in the entire sample and their relative abundances. Used to identify peaks for further imaging.
* **Sum spectrum**: Total ion current summed across all pixels. Emphasizes ions that are abundant across many pixels.
* **Per-pixel spectra (XY and N)**: Reveal the distribution of ions in a specific location, useful for comparing tissue regions (e.g., tumor vs. normal).
> **Why these inputs?**: Spectra are the raw data of MSI. Visualizing them is the first step: checking for unexpected peaks, assessing signal-to-noise ratio, and verifying mass calibration. Applying a mask (e.g., only tumor region) allows researchers to extract region-specific molecular signatures.
---
## Topography and Surface Plots
**Functions**: `image3dPlot`, `imageCPlot`, `triq3dPlot`, `triqCPlot`
Convert a 2D ion image into a 3D surface plot where the height (*z*-axis) represents ion intensity.
* **Inputs**: Current ion image and optional Mask.
* **Why these inputs?**: 3D surfaces and 2D Contour plots make intensity variations more tangible, especially for features with large dynamic range. They are useful for understanding the topography of ion distributions (e.g., a lipid enriched in a specific tissue layer). The 3D plot is interactive, allowing for rotation and zoom.
---
## Small Features and Tools
### Metadata Viewing (`showMetadataBtn`)
Display a summary of the dataset’s metadata: instrument settings, acquisition parameters, image dimensions, number of pixels, *m/z* range, etc.
### Optical Image Overlay (`btnOptical`)
Overlay a histological or photographic image onto the MSI image to correlate molecular distributions with tissue anatomy.
* **Inputs**: Current ion image, Transparency slider (`imgTrans`), and Optical image file (PNG, BMP, JPG, JPEG).
* **Output**: The main image tab updates to show the MSI image (or TrIQ image) blended with the optical background.
### Comparison of Datasets
Functionality to compare two datasets (e.g., control vs. treated). This involves side-by-side visualization or ratio image generation to highlight differences in molecular abundance. It allows visualization of slices with or without TrIQ filtering.
### Registry and Safe Processing
Batch processing is orchestrated by `process_file_safely`. The registry (`registry.json`) stores critical info: source path, image dimensions, summary statistics, and mask path. This file can be edited manually following its generated structure.
---
## Mask Generation and Usage (`mask.jl`)
Create binary masks that define **Regions of Interest (ROI)** in the sample. Masks filter spectra, preprocessing, and images to include only pixels inside the ROI.
### User-Defined Inputs:
* **Drawing tools**: Draw freehand shapes on the ion image or optical overlay.
* **Otsu Scale**: Automatically generates a mask based on the Otsu thresholding algorithm with the slider percentage.
* **Noise Size**: Removes noise from the mask based on the slider percentage.
* **Hole Size**: Removes holes in the mask based on the slider percentage.
* **Smoothing**: Smooths the mask based on the slider percentage.
* **Thresholding**: Set intensity thresholds to automatically segment regions.
* **Import**: Load a pre-existing mask image (e.g., from a segmentation tool).
### Outputs:
* A binary mask matrix (same dimensions as the MSI image) saved in the public/css/masks folder as a `.png` file with the name of the dataset.
* The mask is registered in the dataset’s metadata entry for automatic application in subsequent analyses.
> **Why these inputs?**: Often only a portion of the imaged area contains tissue; the rest could be background. Masking removes background pixels, improving statistical power and avoiding artifacts. Manual delineation based on optical images is standard practice in histology-directed MSI. It can also work for a selection of an area of interest (e.g. only the tumor region) in samples with multiple areas of interest.
---
## Preprocessing
To understand the deep structure of the preprocessing pipeline, please refer to the dedicated guide: **[JuliaMSI Preprocessing Guide](JuliaMSI_Preprocessing.md)**
This document describes the rationale, methods (TIC/PQN, SNIP, Wavelets), and implementation choices in detail.
---
## Practical Tips for Analysts
### Julia Version Compatibility
* The platform is tested and compatible with **Julia 1.11**.
* **Recommendation**: Use Julia 1.11 to ensure all dependencies work correctly. Use **Juliaup** to manage multiple versions easily. Downgrading your system Julia is unnecessary; you can keep multiple versions and launch the platform with the specific executable path.
### First-Time Startup
* The script will detect a missing `Manifest.toml` and automatically instantiate packages. This may take a few minutes. Subsequent startups will be much faster.
* **Caution**: Do not close the terminal/console while the server is running; it hosts the web interface.
### Accessing the Web Interface
* The platform opens your browser to `http://127.0.0.1:1481`.
* If port 1481 is in use, change this segment in the startup script: `up(host="127.0.0.1", port=1481)`.
### Memory Management
* MSI datasets can be large (gigabytes). JuliaMSI includes automatic **garbage collection** and, on Linux, calls `malloc_trim` to return freed memory to the OS.
* **OOM Troubleshooting**: Process fewer *m/z* values at once, use a mask to exclude background, or increase system swap space.
### Registry and File Locations
* `registry.json` stores paths to your datasets. If you move original **imzML** files, update the paths in the JSON or delete the registry file to let the platform rebuild it.
* Processed images and colorbars are saved in the `public/` folder, organized by dataset name.
### Mask Handling
* Ensure manually created masks have the **exact same dimensions** (width × height) as the MSI image.
* Masks are automatically applied during batch processing if the "Mask" toggle is enabled. If you move the mask file, you must update its path in the registry.
---
## Troubleshooting Common Issues
* **"Error loading active file"**: Check that the file is a valid **imzML** or **mzML** and that the corresponding `.ibd` file is in the same directory.
* **Images appear all black**: This can happen if the tolerance is too narrow (no signal captured) or if the mask is incorrectly applied. Verify *m/z* and tolerance.
* **Slow performance**: Reduce the number of *m/z* values or disable the median filter. Ensure you have enough RAM.
* **Browser interface not responding**: Refresh the page. Check the console for errors; you may need to restart the Julia process.
---
## Citation and Contact
Any revision to the code, bug, or error should be reported to: **julian.sierrag@icloud.com**
*Note: This guide is based on collective knowledge of these algorithms; please perform your own research and cite properly. To cite the JuliaMSI platform, please use:*
> José Julián Sierra-Álvarez, Martín Orlando Camargo-Escalante, Carlos Daniel Sierra-Álvarez, Carmelo Hernández-Caricio, Juan Francisco Moreno-Luna, Isabel Buendía-Corona, Robert Winkler,
> **JuliaMSI: A high-performance graphical platform for mass spectrometry imaging data analysis**,
> *Analytica Chimica Acta*, Volume 1377, 2025, 344613, ISSN 0003-2670,
> [https://doi.org/10.1016/j.aca.2025.344613](https://doi.org/10.1016/j.aca.2025.344613)

354
JuliaMSI_Preprocessing.md Normal file
View File

@ -0,0 +1,354 @@
# A Comprehensive Guide to Mass Spectrometry Imaging Preprocessing
## Rationale, Methods, and Implementation
---
### Table of Contents
1. [Introduction](#introduction)
2. [Profile vs. Centroid Mode Data](#profile-vs-centroid-mode-data)
3. [General Principles of MSI Preprocessing](#general-principles-of-msi-preprocessing)
4. [Step-by-Step Breakdown of Preprocessing Steps](#step-by-step-breakdown-of-preprocessing-steps)
* [Validation and Quality Control](#validation-and-quality-control)
* [Intensity Transformation (Stabilization)](#intensity-transformation-stabilization)
* [Smoothing](#smoothing)
* [Baseline Correction](#baseline-correction)
* [Normalization](#normalization)
* [Peak Detection](#peak-detection)
* [Peak Selection (Filtering)](#peak-selection-filtering)
* [Calibration](#calibration)
* [Peak Alignment](#peak-alignment)
* [Peak Binning and Feature Matrix Generation](#peak-binning-and-feature-matrix-generation)
5. [Example Walkthrough: From Raw Spectra to Feature Matrix](#example-walkthrough-from-raw-spectra-to-feature-matrix)
* [Phase 1 – Instrument & Data Characteristics](#phase-1--instrument--data-characteristics)
* [Phase 2 – Noise & Signal Quality](#phase-2--noise--signal-quality)
* [Phase 3 – Mass Accuracy Analysis](#phase-3--mass-accuracy-analysis)
* [Phase 4 – Spatial Region Analysis](#phase-4--spatial-region-analysis)
* [Phase 5 – Peak Characteristics](#phase-5--peak-characteristics)
* [Phase 6 – Preprocessing Recommendations](#phase-6--preprocessing-recommendations)
6. [Extra Tips for the Analyst](#extra-tips-for-the-analyst)
* [Isotopic and Charge State Matters](#isotopic-and-charge-state-matters)
* [Always Visualize the Average Spectrum First](#always-visualize-the-average-spectrum-first)
* [Preprocessing Is Not a One-Shot Process](#preprocessing-is-not-a-one-shot-process)
---
## Introduction
Mass Spectrometry Imaging (MSI) generates thousands of mass spectra, each acquired from a distinct pixel location. In the literature, before any interpretation can be drawn from these data, a series of preprocessing steps must be applied. These steps correct for instrumental artifacts, reduce noise, and align data across pixels to produce a clean, comparable feature matrix suitable for statistical analysis and image reconstruction.
Inside the **JuliaMSI** code (`Preprocessing.jl` and `PreprocessingPipeline.jl`), we have a preprocessing workflow inspired by established R packages like **MALDIquant** and **Cardinal**.
This document explains the rationale behind each step, the required inputs, the scientific meaning of parameters, and the context behind key algorithms such as **TIC normalization**, **PQN**, **wavelet-based peak detection**, and **LOWESS alignment**.
---
## Profile vs. Centroid Mode Data
Mass spectrometry data can be acquired in two primary modes: **profile** and **centroid**. Understanding the difference is crucial because preprocessing steps must be adapted accordingly.
* **Profile mode**:
The instrument records intensity at every digitized point along the *m/z* axis, producing a continuous waveform. Peaks appear as Gaussian-like shapes. This mode retains full peak shape information, which is valuable for accurate mass measurement, resolution assessment, and advanced peak detection. However, file sizes are large and noise is present between peaks.
* **Centroid mode**:
The instrument performs real-time peak detection and reduces each peak to a single centroid (*m/z*, intensity) pair, discarding the underlying profile. This dramatically reduces file size but loses peak shape details. Centroid data are essentially a list of detected peaks with no information about peak width or baseline.
In **profile mode**, algorithms must operate on the full *m/z* and intensity vectors, performing tasks like smoothing, baseline correction, and peak picking on continuous data. In **centroid mode**, these steps are often unnecessary or simplified: smoothing and baseline correction are not applicable because the data are already peak lists; peak picking is replaced by filtering based on intensity thresholds.
> [!IMPORTANT]
> The JuliaMSI code internally distinguishes this difference in some steps by acquiring the mode from the metadata on the files, but it's up to the analyst to decide skipping smoothing/baseline steps for centroid data, for example.
---
## General Principles of MSI Preprocessing
All preprocessing steps share a common goal: **to end up with an aligned feature matrix where rows correspond to pixels (spectra) and columns correspond to *m/z* bins (features).**
The steps are applied sequentially, often in a specific order to avoid compounding errors. For example, **baseline correction should precede normalization** because baseline offsets can distort total ion current calculations. Similarly, **peak alignment should follow peak detection** because alignment requires a set of common peaks across spectra.
---
## Step-by-Step Breakdown of Preprocessing Steps
### Validation and Quality Control
*(This step doesn't require analyst input)*
* **Functions in code**: `validate_spectrum`, `qc_is_empty`, `qc_is_regular`
* **Purpose**: Ensure that each input spectrum is well-formed and contains usable data before any processing.
* **Inputs required**:
* `mz::AbstractVector{<:Real}`: The *m/z* values.
* `intensity::AbstractVector{<:Real}`: The corresponding intensities.
> **Why these inputs?**
> Every subsequent algorithm depends on the integrity of these vectors. Checking for non-emptiness, equal length, finite values, non-negativity, and monotonic *m/z* order prevents cryptic errors later.
---
### Intensity Transformation (Stabilization)
* **Function**: `transform_intensity_core`
* **Purpose**: Apply a variance-stabilizing transformation to make the data more homoscedastic (constant variance across intensity range).
* **Inputs required**:
* `intensity::AbstractVector{<:Real}`
* `method::Symbol` (e.g., `:sqrt`, `:log1p`)
> **Why these inputs?**
> The transformation operates on the intensity vector only; *m/z* values remain unchanged. The choice of method depends on the data characteristics (e.g., presence of zeros).
>
> Stabilization uses functions like **square root**, **log**, **log2**, **log10** and **log1p** helps the analyst stabilize variance for count data, or handle quadratical variance. The analyst must know which one of these methods to use according to the data they are visualizing (e.g., by the mean spectrum plot, or a plot from a specific coordinate).
>
> In MSI in general, **square root stabilization** is common because ion counting statistics approximate Poisson.
>
> **Example**: A spectrum with intense peaks and low baseline noise: after square root, the dynamic range is compressed, making low-intensity features more comparable to high-intensity ones.
---
### Smoothing
* **Function**: `smooth_spectrum_core`, `moving_average_smooth`
* **Purpose**: Reduce high-frequency noise while preserving peak shapes.
* **Inputs required**:
* `y::AbstractVector{<:Real}`: Intensity vector.
* `method::Symbol` (`:savitzky_golay`, `:moving_average`)
* `window::Int`: Size of the smoothing window.
* `order::Int`: Polynomial order for Savitzky–Golay.
> **Why these inputs?**
> Smoothing operates on the intensity vector only. The window size must be an odd integer for symmetry; the polynomial order must be less than the window. These constraints ensure the filter behaves properly (e.g., Savitzky–Golay requires the window bigger than the order parameter).
>
> **Savitzky–Golay** fits a low-degree polynomial to the local window, preserving moments of the peak (height, width) better than a simple moving average, which can distort peak shapes.
>
> **Example**: A noisy peak with 5% random noise: after Savitzky–Golay, the peak becomes smoother while its centroid and height remain accurate.
---
### Baseline Correction
* **Function**: `apply_baseline_correction_core`, `_snip_baseline_impl`, `convex_hull_baseline`, `median_baseline`
* **Purpose**: Estimate and subtract the background signal (baseline) from each spectrum.
* **Inputs required**:
* `y::AbstractVector{<:Real}`: Intensity vector.
* `method::Symbol` (`:snip`, `:convex_hull`, `:median`)
* `iterations::Int` (for SNIP)
* `window::Int` (for median)
> **Why these inputs?**
> Baseline correction requires only the intensity vector. The method-specific parameters control the aggressiveness and shape of the baseline estimate.
>
> * **SNIP (Sensitive Nonlinear Iterative Peak clipping)**: Iteratively replaces each point with the minimum of itself and the average of its neighbors. After many iterations, only the baseline remains. SNIP is the default in **MALDIquant** due to its effectiveness.
> * **Convex hull**: Constructs the upper convex hull of the spectrum (or lower hull for baseline) and interpolates.
> * **Moving median**: Estimates baseline as the local median. Simple but can be biased by dense peaks.
>
> **Example**: In a spectrum with a rising baseline after SNIP correction, the baseline is removed, revealing peaks sitting on a flat background.
---
### Normalization
* **Function**: `apply_normalization_core`, `tic_normalize`, `pqn_normalize`, `median_normalize`, `rms_normalize`
* **Purpose**: Correct for differences in total ion abundance between spectra, making intensities comparable across pixels.
* **Inputs required**:
* `y::AbstractVector{<:Real}` (per-spectrum) or a matrix `M` (for PQN).
* `method::Symbol`
> **Why these inputs?**
> Normalization operates on intensities. **PQN** requires a matrix of spectra (rows = *m/z* bins, columns = spectra) because it computes a reference spectrum across all samples.
>
> * **TIC (Total Ion Current)**: Divides each intensity by the sum of all intensities in the spectrum. A simple algorithm but sensitive to a few very intense peaks (e.g., matrix peaks).
> * **Median normalization**: Divides by the median intensity. Robust to outliers but can be unstable if many intensities are zero.
> * **RMS normalization**: Divides by root-mean-square, preserving the shape of the spectrum while scaling to unit energy.
> * **PQN (Probabilistic Quotient Normalization)**: Calculates a reference spectrum and scales each spectrum by the median of the quotients (intensity/reference). This method is preferred when some peaks vary drastically because it downweights their influence.
>
> **Example**: Two spectra from different tissue regions: one has a huge lipid peak, the other has lower overall signal. TIC would suppress the lipid peak in the first spectrum, possibly losing biological information. PQN would scale both spectra more equitably.
---
### Peak Detection
* **Functions**: `detect_peaks_profile_core`, `detect_peaks_wavelet_core`, `detect_peaks_centroid_core`
* **Purpose**: Identify peaks (signals of interest) in each spectrum.
* **Inputs required**:
* *m/z*, intensity vectors.
* Method-specific parameters: `snr_threshold`, `half_window`, `min_peak_prominence`, `merge_peaks_tolerance`, etc.
> **Why these inputs?**
> Peak detection uses both *m/z* and intensity to locate local maxima and evaluate peak quality. Parameters control sensitivity and filtering. A peak is a local maximum that rises significantly above the local noise. Detection algorithms vary by data mode:
>
> * **Profile mode (local maxima)**: Scans for points that are highest in a sliding window (`half_window`). Then filters by **SNR** (signal-to-noise ratio) and prominence (height relative to surrounding minima). Close peaks are merged if within `merge_peaks_tolerance`.
> * **Wavelet-based**: Uses **Continuous Wavelet Transform (CWT)** to identify peaks at multiple scales, making it robust to varying peak widths. Local maxima in the CWT magnitude indicate peak positions. Wavelet is often used for complex spectra with overlapping peaks. CWT-based peak detection is implemented in **MALDIquant** as the default method for profile data.
> * **Centroid mode**: Simply filters the existing peak list by SNR, as peaks are already identified.
>
> **Example**: A spectrum with a small peak on the shoulder of a large peak: wavelet detection may resolve both, while simple local maxima might miss the shoulder.
---
### Peak Selection (Filtering)
* **Function**: `apply_peak_selection`
* **Purpose**: Remove peaks that do not meet quality criteria, such as minimum SNR, acceptable width, or good Gaussian shape.
* **Inputs required**:
* peaks list (from peak detection).
* Thresholds: `min_snr`, `min_fwhm_ppm`, `max_fwhm_ppm`, `min_shape_r2`.
> **Why these inputs?**
> These thresholds reflect the analyst's confidence in peak quality. They are typically derived from instrument specifications or empirical data.
>
> Not all local maxima are true peaks. Noise spikes, baseline artifacts, or electronic interference can produce spurious peaks. Filtering by **SNR** removes noise. **FWHM** (full width at half maximum) in **ppm** ensures that peaks have physically reasonable widths (e.g., isotopic peaks have predictable spacing). **Shape R²** from a Gaussian fit assesses how well the peak resembles a theoretical peak; poor shape may indicate saturation or interference.
>
> **Example**: A peak with SNR < 3 is likely noise and is discarded. A peak with FWHM > 100 *ppm* might be background rather than a specific ion.
---
### Calibration
*(This step requires that you have a list of Internal Standards inside JuliaMSI preprocessing tab filled in)*
* **Function**: `calibrate_spectra_core`, `find_calibration_peaks_core`
* **Purpose**: Correct systematic mass errors by aligning the *m/z* axis to known reference masses (internal standards).
* **Inputs required**:
* spectra (list of `MutableSpectrum`).
* `internal_standards`: dictionary of theoretical *m/z* values.
* `ppm_tolerance`: matching tolerance.
* `fit_order`: polynomial order for calibration curve.
> **Why these inputs?**
> Calibration requires a set of known masses that are present in the spectra. The algorithm matches detected peaks to these references, then fits a transformation (e.g., linear) to map measured to theoretical *m/z*.
>
> Mass spectrometers drift; this step is made to correct this drift, ensuring that the same *m/z* value corresponds to the same ion across the entire imaging dataset. **Linear or quadratic calibration** is standard. The code currently uses linear interpolation between matched peaks, which is sufficient for most MSI applications where drift is small.
>
> **Example**: A reference peak expected at **760.585 Da** (a lipid) is found at **760.590 Da** in a spectrum. Calibration shifts the entire axis so that the peak aligns with the theoretical value.
---
### Peak Alignment
* **Function**: `align_peaks_lowess_core`
* **Purpose**: Correct residual *m/z* shifts between spectra after calibration, ensuring that the same analyte peak appears at the same *m/z* across all pixels.
* **Inputs required**:
* Reference spectrum's peak *m/z* list.
* Target spectrum's peak *m/z* list.
* Parameters: method (`:lowess`, `:linear`, `:ransac`), `tolerance`, `tolerance_unit`, `max_shift_ppm`.
> **Why these inputs?**
> Alignment operates on peak lists, not full spectra, to be efficient. It finds common peaks between a reference spectrum (e.g., the one with most peaks) and each target, then fits a warping function.
>
> Even after calibration, small nonlinear shifts may persist due to sample topography or matrix effects. Alignment brings all spectra onto a common *m/z* scale, which is essential for constructing a feature matrix.
>
> **LOWESS (Locally Weighted Scatterplot Smoothing)**: Fits a smooth curve through the matched peak pairs, allowing for nonlinear corrections. It is robust to outliers and does not assume a global polynomial shape. It has been used in **MALDIquant**. **RANSAC (Random Sample Consensus)** is a robust alternative. The code defaults to LOWESS because it is well-suited for smooth, continuous drift.
>
> **Example**: Peaks in a target spectrum are systematically shifted by **+0.01 Da** at low *m/z* and **+0.005 Da** at high *m/z*. LOWESS fits a curve that corrects each peak individually, producing aligned *m/z* values.
---
### Peak Binning and Feature Matrix Generation
* **Function**: `bin_peaks_core`
* **Purpose**: Group peaks from all spectra into common *m/z* bins, producing a final feature matrix.
* **Inputs required**:
* spectra (with detected peaks).
* `params::PeakBinning`: method (`:adaptive` or `:uniform`), `tolerance`, `frequency_threshold`, etc.
> **Why these inputs?**
> Binning requires the list of peaks from each spectrum. The method determines how bins are defined.
>
> * **Adaptive binning**: Clusters peaks based on their *m/z* proximity (using tolerance). Bins are created where peaks are dense, and bins without enough peaks are discarded.
> * **Uniform binning**: Divides the total *m/z* range into a fixed number of equally spaced bins.
>
> After alignment, peaks from different spectra still have slightly different *m/z* values due to residual noise. Binning assigns them to discrete bins, creating a matrix suitable for multivariate analysis. The **frequency threshold** removes bins that appear in too few spectra (likely noise).
>
> **Example**: Peaks around *m/z* **760.58** from 100 spectra are grouped into one bin. The intensity for each spectrum in that bin is the maximum intensity of peaks assigned to it.
---
## Example Walkthrough: From Raw Spectra to Feature Matrix
### with Intelligent Parameter Recommendations
To illustrate the entire pipeline, consider a mass spectrometry imaging dataset of a mouse brain section acquired in **profile mode**. The analyst has a list of internal standards for lipids (e.g., known phospholipids at *m/z* **760.585**, **798.540**, **834.528**) to aid in mass calibration.
Before running the full preprocessing, the **JuliaMSI** code executes the intelligent pre-analysis pipeline (`main_precalculation` which is a wrapper for `run_preprocessing_analysis`). This function samples a representative subset of spectra (e.g., 100 pixels) and performs a battery of diagnostic tests to suggest optimal parameters.
#### The pre-analysis yields a report:
* **Phase 1 – Instrument & Data Characteristics**:
* Acquisition mode: **profile** (continuous waveforms)
* *m/z* axis type: **regular** (constant step of ~0.01 Da)
* Dynamic range: ~4 orders of magnitude
* **Phase 2 – Noise & Signal Quality**:
* Estimated noise level (median absolute deviation): 45 counts
* Suggested SNR threshold: **3.0**
* TIC coefficient of variation: **35%** → indicates moderate heterogeneity, so robust normalization (e.g., **PQN**) may be beneficial.
* **Phase 3 – Mass Accuracy Analysis** (using the provided lipid standards):
* Mean PPM error: **12.3 ppm**
* Standard deviation: **8.1 ppm**
* Suggested bin tolerance (mean + 3σ): **36.6 ppm** (capped between 10 and 100 ppm).
* **Phase 4 – Spatial Region Analysis**:
* *(if region masks were supplied) – skipped in this example.*
* **Phase 5 – Peak Characteristics**:
* Mean FWHM: **22.5 ppm** (converted from Δm at *m/z* 500)
* Mean Gaussian R²: **0.85** (peaks are well-shaped)
* Mean peaks per spectrum: **187**
* **Phase 6 – Preprocessing Recommendations**:
* The system synthesizes all information into a dictionary of suggested parameters for each step. Parameters that cannot be determined automatically (e.g., frequency threshold for binning, minimum matched peaks) are left as `nothing`; the analyst must supply them based on experimental goals.
> [!NOTE]
> If you have internal standards, remember to re-run the preprocessing recommendations after you have loaded them. This is only for precalculation and is not a required step for analysts that know about their data.
#### Execution of the Pipeline:
The analyst can now execute the pipeline, possibly overriding a few parameters (e.g., increasing SNR threshold to 4 for more stringent peak picking).
1. **Validation**: All spectra pass QC; one empty spectrum is flagged and removed.
2. **Intensity transformation (sqrt)**: Stabilizes variance, making low-intensity features more comparable.
3. **Baseline correction (SNIP, 100 iterations)**: Removes background using the recommended window of 200 points.
4. **Smoothing (Savitzky–Golay, window 9, order 2)**: Reduces noise while preserving peak shape.
5. **Peak detection (profile mode)**: Detects peaks with SNR ≥ 3, prominence ≥ 0.002, and uses the recommended half-window of 5 points. Peaks closer than 18 *ppm* are merged.
6. **Peak selection**: Filters peaks: keeps only those with SNR ≥ 3, FWHM between 11 and 55 *ppm*, and shape R² ≥ 0.68.
7. **Calibration**: Using the internal standards, each spectrum’s *m/z* axis is linearly corrected based on matched peaks within 36.6 *ppm* tolerance.
8. **Peak alignment (LOWESS)**: Aligns all spectra to a reference spectrum (the one with most peaks) with a span of 0.6, using the same 36.6 *ppm* tolerance. This corrects residual nonlinear drift.
9. **Normalization (RMS)**: Divides each spectrum by its root-mean-square intensity, chosen because TIC variation was high.
10. **Peak binning (adaptive)**: Groups peaks within 36.6 *ppm* into bins. Bins containing fewer than 3 peaks or appearing in <20% of spectra are discarded. The final bin centers are intensity-weighted averages.
11. **Feature matrix**: A **999 × 425 matrix** is produced, where rows are pixels and columns are *m/z* bins. This matrix is ready for PCA, segmentation, or statistical testing.
#### How the Analyst Uses the Recommendations:
* The suggested values are data-driven, but the analyst can override them based on domain knowledge. Always verify; experimenting is encouraged as data integrity is kept.
* Parameters like `frequency_threshold` in peak selection and binning require the analyst’s judgment. A threshold of **0.2** means that a bin must contain a peak in at least 20% of all pixels to be retained. This reduces noise but may discard rare but real signals.
* This is an iterative refinement; the analyst may rerun the analysis with adjusted parameters (e.g., a lower bin tolerance).
* JuliaMSI aims for a tailored preprocessing workflow that respects both the instrument’s characteristics and the biological question.
---
## Extra Tips for the Analyst
The preprocessing pipeline is powerful, but its success ultimately depends on the analyst’s informed choices. Here are practical tips to guide you through the process, especially when working with internal standards and interpreting intermediate results.
> [!WARNING]
> This is a really resource-intensive pipeline dependent on the sample, so run your precautions.
### Isotopic and Charge State Matters
When providing a dictionary of reference masses for calibration, ensure that the theoretical *m/z* values correspond to the same charge state and isotopic composition as the peaks you expect to see in your data.
* **Example**: If your instrument detects lipids mainly as **[M+H]⁺** ions, the reference mass must be the monoisotopic mass plus a proton. Using the neutral mass will introduce a systematic offset of **1.0078 Da**, leading to poor calibration. (Make the same adjustment for negative ion mode).
* If isotopic peaks are prominent, you may want to include multiple isotopes (e.g., **M**, **M+1**, **M+2**) as separate references, but be aware that their relative intensities vary. The calibration algorithm will match the closest detected peak. This is common for MSI preprocessing.
* The pre-analysis suggests a bin tolerance based on observed mass error. A tolerance that is too large may match the wrong peak (e.g., an isobaric interference); a tolerance that is too small may miss valid matches.
### Always Visualize the Average Spectrum First
Before running any preprocessing, generate a sum or average spectrum. This single plot tells you:
1. **Baseline shape**: If it's flat, rising, or undulating, it can guide the baseline correction method (**SNIP** works well for most shapes; **convex hull** is faster but may fail on curvy baselines).
2. **Noise level**: By zooming in.
3. ***m/z* Scatter**: How scattered are the spectra *m/z*.
4. **Peak density & Intensity range**: Influence the peak detection thresholds (SNR, prominence) and binning strategy.
* If the average spectrum shows a **rising baseline** from low to high *m/z*, SNIP with ~100 iterations will handle it well. If extremely wavy, increase iterations or try a different method.
* If running from the console (`run_preprocessing.jl`), the output of `run_preprocessing_analysis` provides:
* **TIC coefficient of variation (CV)**: A CV > 30% suggests substantial pixel-to-pixel variation, favoring **RMS** or **PQM** over TIC.
* **Mean PPM error**: If >20 *ppm*, calibration is essential. If <5 *ppm* and drift is linear, you might skip alignment.
* **FWHM and Gaussian R²**: Indicate peak quality. If R² < 0.6, peaks may be poorly shaped; consider increasing the `min_peak_shape_r2` threshold.
* **Acquisition mode**: For **centroid data**, never apply smoothing or baseline correction because those steps assume continuous waveforms.
* The pipeline will set `:method = :centroid` in peak picking and use permissive filters.
* Use a low SNR threshold (e.g., 2.0).
* Disable shape-based filtering (`min_peak_shape_r2 = 0.0`).
* Centroid data do not retain peak width information; set `min_fwhm_ppm = 0.0` and `max_fwhm_ppm = 500.0` to avoid discarding real peaks.
* Adaptive binning still works; tolerance should be based on expected mass accuracy.
### Preprocessing Is Not a One-Shot Process
The recommended parameters are a starting point, not a final answer.
* **Reproducibility**: Keep a record of all parameters used. The `params` dictionary can be saved as a **JSON**. This is essential for publications.
* **Biological Context**: If looking for rare metabolites appearing in few pixels, lower the `frequency_threshold` in binning. If interested only in high-abundance lipids, increase the SNR threshold.
* **Modularity**: The pipeline is modular; you can always rerun parts of it with adjusted parameters without starting from scratch.
---
*Any revision to the code, bug, error, report to julian.sierrag@icloud.com, if you we're to cite this file, remember i provide no citation and its based on internet knowledge of the web about these algorithms, so do your own research and cite properly. to cite JuliaMSI platform, please cite:*
José Julián Sierra-Álvarez, Martín Orlando Camargo-Escalante, Carlos Daniel Sierra-Álvarez, Carmelo Hernández-Caricio, Juan Francisco Moreno-Luna, Isabel Buendía-Corona, Robert Winkler,
JuliaMSI: A high-performance graphical platform for mass spectrometry imaging data analysis,
Analytica Chimica Acta,
Volume 1377,
2025,
344613,
ISSN 0003-2670,
https://doi.org/10.1016/j.aca.2025.344613

8
README.md
View File

@ -3,7 +3,7 @@
https://codeberg.org/LabABI/JuliaMSI https://codeberg.org/LabABI/JuliaMSI
## Local Installation ## Local Installation
1. Make sure you have Julia installed and at least Julia version 1.6, If not, follow the next installation guide which contains the latest files and instructions for [juliaup](https://github.com/JuliaLang/juliaup) (recommended) alternatively, you can follow the official Julia installation instructions at [julia](https://julialang.org/downloads/platform/). 1. Make sure you have Julia installed and at least Julia version 1.6 but not higher than 1.11, If not, follow the next installation guide which contains the latest files and instructions for [juliaup](https://github.com/JuliaLang/juliaup) (recommended) alternatively, you can follow the official Julia installation instructions at [julia](https://julialang.org/downloads/platform/).
2. Download and decompress the following repository (You can download the repository from [download](https://codeberg.org/LabABI/JuliaMSI/archive/main.zip) or locate the option manually at the top right of the page, in the box with [•••] (More operations). For a more advanced method using Git, use this link: https://codeberg.org/LabABI/JuliaMSI.git). 2. Download and decompress the following repository (You can download the repository from [download](https://codeberg.org/LabABI/JuliaMSI/archive/main.zip) or locate the option manually at the top right of the page, in the box with [•••] (More operations). For a more advanced method using Git, use this link: https://codeberg.org/LabABI/JuliaMSI.git).
3. Ensure that you have decompressed the ZIP file and know the working directory where it is located (the folder on your computer where you saved the JuliaMSI files). 3. Ensure that you have decompressed the ZIP file and know the working directory where it is located (the folder on your computer where you saved the JuliaMSI files).
4. Julia Version 1.11 is required. You might need to downgrade your whole Julia installation to support it. 4. Julia Version 1.11 is required. You might need to downgrade your whole Julia installation to support it.
@ -42,3 +42,9 @@ JuliaMSI is published under the terms of the MIT License.
Robert Winkler. (2023). mzML mass spectrometry and imzML mass spectrometry imaging test data [Data set]. Robert Winkler. (2023). mzML mass spectrometry and imzML mass spectrometry imaging test data [Data set].
Zenodo. <https://doi.org/10.5281/zenodo.10084132> Zenodo. <https://doi.org/10.5281/zenodo.10084132>
## User Guides
For a general guide on how to use the JuliaMSI platform, please refer to the [JuliaMSI Main Guide](JuliaMSI_MainGuide.md).
For a more detailed guide on the preprocessing pipeline, please refer to the [JuliaMSI Preprocessing Guide](JuliaMSI_Preprocessing.md).

45
app.jl
View File

@ -706,7 +706,16 @@ end
is_processing = false is_processing = false
end end
end end
#=
@onbutton btnSearch
This reactive handler job is to open a file dialog and prompt the user to select a file
its the entry point to any process done later in the interface, it creates a new dataset entry in the json
and updates the list of available datasets in the interface. preprocesses the file and saves
metadata to the json file. if the file is already in the json file, it will skip the preprocessing
and load the metadata from the json file. creates the suggestions to the preprocessing pipeline, and makes
sure the file can be processed by later steps like mainProcess
=#
@onbutton btnSearch begin @onbutton btnSearch begin
is_processing = true is_processing = true
push!(__model__) push!(__model__)
@ -1219,6 +1228,9 @@ end
is_processing = false is_processing = false
end end
#=
This reactive handler job is to run the full preprocessing pipeline on the selected dataset.
=#
@onbutton run_full_pipeline begin @onbutton run_full_pipeline begin
is_processing = true is_processing = true
push!(__model__) push!(__model__)
@ -1915,7 +1927,11 @@ end
println("DEBUG: Memory cleanup complete.") println("DEBUG: Memory cleanup complete.")
end end
end end
#=
@onbutton recalculate_suggestions_btn
This reactive handler job is to recalculate the suggestions for the preprocessing pipeline.
=#
@onbutton recalculate_suggestions_btn begin @onbutton recalculate_suggestions_btn begin
is_processing = true is_processing = true
push!(__model__) push!(__model__)
@ -2089,6 +2105,11 @@ end
end end
end end
#=
@onbutton addReferencePeak
This reactive handler job is to add a new reference peak to the list of reference peaks.
=#
@onbutton addReferencePeak begin @onbutton addReferencePeak begin
is_processing = true is_processing = true
new_list = deepcopy(reference_peaks_list) new_list = deepcopy(reference_peaks_list)
@ -2097,6 +2118,11 @@ end
is_processing = false is_processing = false
end end
#=
@onbutton remove_peak_trigger
This reactive handler job is to remove a reference peak from the list of reference peaks.
=#
@onbutton remove_peak_trigger begin @onbutton remove_peak_trigger begin
is_processing = true is_processing = true
if action_index > -1 if action_index > -1
@ -2157,7 +2183,11 @@ end
is_processing = false is_processing = false
end end
# This new handler correctly adds the file from full_route to the batch list. #=
@onbutton btnAddBatch
This reactive handler job is to add a new file to the batch list.
=#
@onbutton btnAddBatch begin @onbutton btnAddBatch begin
is_processing = true is_processing = true
push!(__model__) push!(__model__)
@ -2179,6 +2209,11 @@ end
is_processing = false is_processing = false
end end
#=
@onbutton clear_batch_btn
This reactive handler job is to clear the batch list in the UI
=#
@onbutton clear_batch_btn begin @onbutton clear_batch_btn begin
is_processing = true is_processing = true
push!(__model__) push!(__model__)
@ -2200,7 +2235,11 @@ end
msg = "File automatically added to batch" msg = "File automatically added to batch"
end end
end end
#=
@onbutton showMetadataBtn
This reactive handler job is to show the metadata of the selected dataset.
=#
@onbutton showMetadataBtn begin @onbutton showMetadataBtn begin
if !isempty(available_folders) if !isempty(available_folders)
if !isempty(selected_folder_main) if !isempty(selected_folder_main)
@ -2214,7 +2253,11 @@ end
warning_msg = true warning_msg = true
end end
end end
#=
@onchange selected_folder_metadata
This reactive handler job is to update the metadata table when the selected folder changes.
=#
@onchange selected_folder_metadata begin @onchange selected_folder_metadata begin
if !isempty(selected_folder_metadata) if !isempty(selected_folder_metadata)
registry = load_registry(registry_path) registry = load_registry(registry_path)

112
julia_imzML_visual.jl
View File

@ -906,7 +906,7 @@ end
""" """
load_registry(registry_path) load_registry(registry_path)
Loads the dataset registry from a JSON file. Loads the dataset registry from a JSON file, with retries for robustness on Windows.
# Arguments # Arguments
- `registry_path`: Path to the `registry.json` file. - `registry_path`: Path to the `registry.json` file.
@ -916,17 +916,24 @@ Loads the dataset registry from a JSON file.
""" """
function load_registry(registry_path) function load_registry(registry_path)
return lock(REGISTRY_LOCK) do return lock(REGISTRY_LOCK) do
if isfile(registry_path) clean_path = abspath(registry_path)
if isfile(clean_path)
for attempt in 1:5
try try
JSON.parsefile(registry_path, dicttype=Dict{String,Any}) return JSON.parsefile(clean_path, dicttype=Dict{String,Any})
catch e catch e
@error "Failed to parse registry.json: $e" if attempt == 5
Dict{String,Any}() @error "Final attempt failed to parse registry.json: $(sprint(showerror, e))"
end return Dict{String,Any}()
else else
Dict{String,Any}() @warn "Attempt $attempt to read registry.json failed, retrying in 0.2s... Error: $(sprint(showerror, e))"
sleep(0.2)
end end
end end
end
end
return Dict{String,Any}()
end
end end
""" """
@ -987,6 +994,47 @@ function extract_metadata(msi_data::MSIData, source_path::String)
) )
end end
"""
save_registry(registry_path, registry_data)
Saves the dataset registry to a JSON file, ensuring thread-safe access and robustness on Windows.
# Arguments
- `registry_path`: Path to the `registry.json` file.
- `registry_data`: The dictionary containing the registry data to save.
"""
function save_registry(registry_path, registry_data)
lock(REGISTRY_LOCK) do
# Ensure the path is clean and valid for the OS
clean_path = abspath(registry_path)
dir = dirname(clean_path)
if !isdir(dir)
mkpath(dir)
end
# On Windows, opening files for writing can occasionally fail if another process
# (like an anti-virus or file indexer) holds a transient lock.
# Retry with a small delay to improve robustness.
for attempt in 1:5
try
open(clean_path, "w") do f
JSON.print(f, registry_data, 4)
end
return true # Success
catch e
if attempt == 5
@error "Final attempt failed to write to registry.json: $(sprint(showerror, e))"
return false
else
@warn "Attempt $attempt to write to registry.json failed (Error: $(sprint(showerror, e))), retrying in 0.2s..."
sleep(0.2)
end
end
end
return false
end
end
""" """
update_registry(registry_path, dataset_name, source_path, metadata=nothing, is_imzML=false) update_registry(registry_path, dataset_name, source_path, metadata=nothing, is_imzML=false)
@ -1001,16 +1049,7 @@ Adds or updates an entry in the dataset registry JSON file.
""" """
function update_registry(registry_path, dataset_name, source_path, metadata=nothing, is_imzML=false) function update_registry(registry_path, dataset_name, source_path, metadata=nothing, is_imzML=false)
lock(REGISTRY_LOCK) do lock(REGISTRY_LOCK) do
registry = if isfile(registry_path) registry = load_registry(registry_path)
try
JSON.parsefile(registry_path, dicttype=Dict{String,Any})
catch e
@error "Failed to parse registry.json while updating: $e"
Dict{String,Any}() # Start with empty if parsing fails
end
else
Dict{String,Any}()
end
# Get existing entry if it exists, otherwise create new one # Get existing entry if it exists, otherwise create new one
existing_entry = get(registry, dataset_name, Dict{String,Any}()) existing_entry = get(registry, dataset_name, Dict{String,Any}())
@ -1030,34 +1069,8 @@ function update_registry(registry_path, dataset_name, source_path, metadata=noth
registry[dataset_name] = entry registry[dataset_name] = entry
try # Use the robust save function
open(registry_path, "w") do f save_registry(registry_path, registry)
JSON.print(f, registry, 4)
end
catch e
@error "Failed to write to registry.json: $e"
end
end
end
"""
save_registry(registry_path, registry_data)
Saves the dataset registry to a JSON file, ensuring thread-safe access.
# Arguments
- `registry_path`: Path to the `registry.json` file.
- `registry_data`: The dictionary containing the registry data to save.
"""
function save_registry(registry_path, registry_data)
lock(REGISTRY_LOCK) do
try
open(registry_path, "w") do f
JSON.print(f, registry_data, 4)
end
catch e
@error "Failed to write to registry.json: $e"
end
end end
end end
@ -1208,13 +1221,8 @@ function update_registry_mask_fields(registry_path, dataset_name, has_mask, mask
) )
end end
try # Use the robust save function
open(registry_path, "w") do f save_registry(registry_path, registry)
JSON.print(f, registry, 4)
end
catch e
@error "Failed to write to registry.json: $e"
end
end end
end end

12
mask.jl
View File

@ -799,11 +799,7 @@ end
end end
# Save the updated registry # Save the updated registry
lock(REGISTRY_LOCK) do save_registry(reg_path, registry)
open(reg_path, "w") do f
JSON.print(f, registry, 4)
end
end
@info "Registry updated with mask: $(final_mask_name)" @info "Registry updated with mask: $(final_mask_name)"
@ -876,11 +872,7 @@ end
if !isempty(new_folders) || !isempty(removed_folders) if !isempty(new_folders) || !isempty(removed_folders)
println("Registry changed, saving...") println("Registry changed, saving...")
lock(REGISTRY_LOCK) do save_registry(reg_path, registry)
open(reg_path, "w") do f
JSON.print(f, registry, 4)
end
end
end end
all_folders = sort(collect(keys(registry)), lt=natural) all_folders = sort(collect(keys(registry)), lt=natural)

51
src/BloomFilters.jl
View File

@ -24,6 +24,15 @@ end
StreamingBloomFilter StreamingBloomFilter
A memory-efficient Bloom filter that doesn't store all values at once. A memory-efficient Bloom filter that doesn't store all values at once.
# Fields
- `bits::BitVector`: The underlying bit array
- `size::Int`: Number of bits in the filter
- `hash_count::Int`: Number of hash functions to use
- `seed::UInt64`: Random seed for hash functions
- `count::Int`: Number of elements added (for monitoring)
- `current_chunk::Vector{Int}`: Current chunk of data being processed
- `chunk_size::Int`: Size of each chunk
""" """
mutable struct StreamingBloomFilter mutable struct StreamingBloomFilter
bits::BitVector bits::BitVector
@ -51,6 +60,9 @@ Creates a Bloom filter optimized for the expected number of elements and desired
- `seed::Union{UInt32,UInt64}`: Keyword description - `seed::Union{UInt32,UInt64}`: Keyword description
(**Default**: `0x12345678`) (**Default**: `0x12345678`)
# Returns
- A `BloomFilter{T}`.
""" """
function BloomFilter{T}(expected_elements::Int, false_positive_rate::Float64=0.01; seed::Union{UInt32,UInt64}=0x12345678) where T function BloomFilter{T}(expected_elements::Int, false_positive_rate::Float64=0.01; seed::Union{UInt32,UInt64}=0x12345678) where T
# Convert seed to UInt64 for consistency # Convert seed to UInt64 for consistency
@ -73,6 +85,9 @@ Calculates the optimal number of bits for a Bloom filter
- `n::Int`: Expected number of elements - `n::Int`: Expected number of elements
- `p::Float64`: Desired false positive rate - `p::Float64`: Desired false positive rate
# Returns
- An `Int` to determine the optimal number of bits for a Bloom filter.
""" """
function optimal_bit_size(n::Int, p::Float64)::Int function optimal_bit_size(n::Int, p::Float64)::Int
if p <= 0.0 || p >= 1.0 if p <= 0.0 || p >= 1.0
@ -91,7 +106,10 @@ Calculates the optimal number of hash functions for a Bloom filter.
# Arguments # Arguments
- `n::Int`: Expected number of elements - `n::Int`: Expected number of elements
- `p::Float64`: Desired false positive rate - `m::Int`: Number of bits in the filter
# Returns
- An `Int` to determine the optimal number of hash functions for a Bloom filter.
""" """
function optimal_hash_count(n::Int, m::Int)::Int function optimal_hash_count(n::Int, m::Int)::Int
if n <= 0 || m <= 0 if n <= 0 || m <= 0
@ -106,6 +124,16 @@ end
hash_functions(item::T, count::Int, size::Int, seed::UInt64) -> Vector{Int} hash_functions(item::T, count::Int, size::Int, seed::UInt64) -> Vector{Int}
Generates multiple hash values for an item using double hashing technique. Generates multiple hash values for an item using double hashing technique.
# Arguments
- `item::T`: Argument description
- `count::Int`: Argument description
- `size::Int`: Argument description
- `seed::UInt64`: Argument description
# Returns
- A `Vector{Int}` of hash values.
""" """
function hash_functions(item::T, count::Int, size::Int, seed::UInt64) where T function hash_functions(item::T, count::Int, size::Int, seed::UInt64) where T
# Use Julia's built-in hash with different seeds # Use Julia's built-in hash with different seeds
@ -129,7 +157,7 @@ end
Adds an element to the Bloom filter. Adds an element to the Bloom filter.
# Arguments # Arguments:
- `bf::BloomFilter{T}`: Argument description - `bf::BloomFilter{T}`: Argument description
- `item::T`: Argument description - `item::T`: Argument description
@ -150,6 +178,11 @@ end
Checks whether an element is possibly in the Bloom filter. Checks whether an element is possibly in the Bloom filter.
Returns `true` if the element might be in the set, `false` if it's definitely not. Returns `true` if the element might be in the set, `false` if it's definitely not.
# Arguments:
- `item::T`: Argument description
- `bf::BloomFilter{T}`: Argument description
""" """
function Base.in(item::T, bf::BloomFilter{T})::Bool where T function Base.in(item::T, bf::BloomFilter{T})::Bool where T
hashes = hash_functions(item, bf.hash_count, bf.size, bf.seed) hashes = hash_functions(item, bf.hash_count, bf.size, bf.seed)
@ -166,6 +199,11 @@ end
contains(bf::BloomFilter{T}, item::T) -> Bool contains(bf::BloomFilter{T}, item::T) -> Bool
Alias for `in()` for compatibility with your existing code. Alias for `in()` for compatibility with your existing code.
# Arguments:
- `bf::BloomFilter{T}`: Argument description
- `item::T`: Argument description
""" """
contains(bf::BloomFilter{T}, item::T) where T = item in bf contains(bf::BloomFilter{T}, item::T) where T = item in bf
@ -173,6 +211,11 @@ contains(bf::BloomFilter{T}, item::T) where T = item in bf
add!(bf::BloomFilter{T}, item::T) add!(bf::BloomFilter{T}, item::T)
Alias for `push!()` for compatibility with your existing code. Alias for `push!()` for compatibility with your existing code.
# Arguments:
- `bf::BloomFilter{T}`: Argument description
- `item::T`: Argument description
""" """
add!(bf::BloomFilter{T}, item::T) where T = push!(bf, item) add!(bf::BloomFilter{T}, item::T) where T = push!(bf, item)
@ -180,6 +223,10 @@ add!(bf::BloomFilter{T}, item::T) where T = push!(bf, item)
false_positive_rate(bf::BloomFilter) -> Float64 false_positive_rate(bf::BloomFilter) -> Float64
Estimates the current false positive rate of the Bloom filter. Estimates the current false positive rate of the Bloom filter.
# Arguments:
- `bf::BloomFilter`: Argument description
""" """
function false_positive_rate(bf::BloomFilter)::Float64 function false_positive_rate(bf::BloomFilter)::Float64
if bf.count == 0 if bf.count == 0

72
src/Common.jl
View File

@ -8,7 +8,7 @@ using Base.Threads
A thread-safe pool of `Vector{UInt8}` buffers, categorized by their size. A thread-safe pool of `Vector{UInt8}` buffers, categorized by their size.
This helps reduce memory allocations and garbage collection overhead by reusing buffers. This helps reduce memory allocations and garbage collection overhead by reusing buffers.
# Fields # Arguments:
- `lock`: A `ReentrantLock` to ensure thread-safe access to the pool. - `lock`: A `ReentrantLock` to ensure thread-safe access to the pool.
- `buffers`: A dictionary mapping buffer size (in bytes) to a list of available buffers of that size. - `buffers`: A dictionary mapping buffer size (in bytes) to a list of available buffers of that size.
""" """
@ -27,11 +27,11 @@ end
Retrieves a `Vector{UInt8}` buffer of at least `size` bytes from the pool. Retrieves a `Vector{UInt8}` buffer of at least `size` bytes from the pool.
If no suitable buffer is available, a new one is allocated. If no suitable buffer is available, a new one is allocated.
# Arguments # Arguments:
- `pool`: The `BufferPool` to retrieve the buffer from. - `pool`: The `BufferPool` to retrieve the buffer from.
- `size`: The minimum desired size of the buffer in bytes. - `size`: The minimum desired size of the buffer in bytes.
# Returns # Returns:
- A `Vector{UInt8}` buffer. - A `Vector{UInt8}` buffer.
""" """
function get_buffer!(pool::BufferPool, size::Int)::Vector{UInt8} function get_buffer!(pool::BufferPool, size::Int)::Vector{UInt8}
@ -58,7 +58,7 @@ end
Returns a `Vector{UInt8}` buffer to the pool for future reuse. Returns a `Vector{UInt8}` buffer to the pool for future reuse.
# Arguments # Arguments:
- `pool`: The `BufferPool` to return the buffer to. - `pool`: The `BufferPool` to return the buffer to.
- `buffer`: The `Vector{UInt8}` buffer to release. - `buffer`: The `Vector{UInt8}` buffer to release.
""" """
@ -78,6 +78,11 @@ end
SimpleBufferPool SimpleBufferPool
A dramatically simplified buffer pool that avoids complex locking. A dramatically simplified buffer pool that avoids complex locking.
# Arguments:
- `buffers::Dict{Int, Vector{Vector{UInt8}}}`: A dictionary mapping buffer size (in bytes) to a list of available buffers of that size.
- `max_pool_size::Int`: The maximum number of buffers to keep in the pool.
- `lock::ReentrantLock`: A lock to ensure thread-safe access to the pool.
""" """
mutable struct SimpleBufferPool mutable struct SimpleBufferPool
buffers::Dict{Int, Vector{Vector{UInt8}}} buffers::Dict{Int, Vector{Vector{UInt8}}}
@ -85,8 +90,26 @@ mutable struct SimpleBufferPool
lock::ReentrantLock # Add lock for thread safety lock::ReentrantLock # Add lock for thread safety
end end
"""
SimpleBufferPool()
Creates a new `SimpleBufferPool` with a default maximum pool size of 50.
"""
SimpleBufferPool() = SimpleBufferPool(Dict{Int, Vector{Vector{UInt8}}}(), 50, ReentrantLock()) SimpleBufferPool() = SimpleBufferPool(Dict{Int, Vector{Vector{UInt8}}}(), 50, ReentrantLock())
"""
get_buffer!(pool::SimpleBufferPool, size::Int) -> Vector{UInt8}
Retrieves a `Vector{UInt8}` buffer of at least `size` bytes from the pool.
If no suitable buffer is available, a new one is allocated.
# Arguments:
- `pool`: The `SimpleBufferPool` to retrieve the buffer from.
- `size`: The minimum desired size of the buffer in bytes.
# Returns:
- A `Vector{UInt8}` buffer.
"""
function get_buffer!(pool::SimpleBufferPool, size::Int)::Vector{UInt8} function get_buffer!(pool::SimpleBufferPool, size::Int)::Vector{UInt8}
lock(pool.lock) do lock(pool.lock) do
# Check for existing buffers of exact size first # Check for existing buffers of exact size first
@ -99,6 +122,15 @@ function get_buffer!(pool::SimpleBufferPool, size::Int)::Vector{UInt8}
return Vector{UInt8}(undef, size) return Vector{UInt8}(undef, size)
end end
"""
release_buffer!(pool::SimpleBufferPool, buffer::Vector{UInt8})
Returns a `Vector{UInt8}` buffer to the pool for future reuse.
# Arguments:
- `pool`: The `SimpleBufferPool` to return the buffer to.
- `buffer`: The `Vector{UInt8}` buffer to release.
"""
function release_buffer!(pool::SimpleBufferPool, buffer::Vector{UInt8}) function release_buffer!(pool::SimpleBufferPool, buffer::Vector{UInt8})
size = length(buffer) size = length(buffer)
@ -132,15 +164,26 @@ struct InvalidSpectrumError <: MSIError
end end
# --- Shared Constants --- # --- Shared Constants ---
const DEFAULT_CACHE_SIZE = 100 const DEFAULT_CACHE_SIZE = 100 # Default cache size for spectra
const DEFAULT_NUM_BINS = 2000 const DEFAULT_NUM_BINS = 2000 # Default number of bins for spectra
const DEFAULT_BLOOM_FILTER_SIZE = 10000 # Default size for empty spectra const DEFAULT_BLOOM_FILTER_SIZE = 10000 # Default size for empty spectra
const DEFAULT_FALSE_POSITIVE_RATE = 0.01 const DEFAULT_FALSE_POSITIVE_RATE = 0.01 # Default false positive rate for bloom filters
""" """
validate_spectrum_data(mz, intensity, id) validate_spectrum_data(mz, intensity, id)
Checks a spectrum for common data integrity issues. Checks a spectrum for common data integrity issues.
# Arguments:
- `mz`: The m/z values of the spectrum.
- `intensity`: The intensity values of the spectrum.
- `id`: The ID of the spectrum.
# Returns:
- `true` if the spectrum is valid.
# Throws:
- `InvalidSpectrumError` if the spectrum is invalid.
""" """
function validate_spectrum_data(mz::AbstractVector, intensity::AbstractVector, id) function validate_spectrum_data(mz::AbstractVector, intensity::AbstractVector, id)
if length(mz) != length(intensity) if length(mz) != length(intensity)
@ -163,6 +206,12 @@ end
Memory-optimized version that processes data in chunks to reduce temporary allocations Memory-optimized version that processes data in chunks to reduce temporary allocations
and avoid holding large intermediate arrays. and avoid holding large intermediate arrays.
# Arguments:
- `mz`: The m/z values of the spectrum.
# Returns:
- A `BloomFilter{Int}` for the spectrum.
""" """
function create_bloom_filter_for_spectrum(mz::AbstractVector{<:Real})::BloomFilter{Int} function create_bloom_filter_for_spectrum(mz::AbstractVector{<:Real})::BloomFilter{Int}
if isempty(mz) if isempty(mz)
@ -193,6 +242,9 @@ end
empty_bloom_filter() -> BloomFilter{Float64} empty_bloom_filter() -> BloomFilter{Float64}
Creates an empty Bloom filter for spectra with no data. Creates an empty Bloom filter for spectra with no data.
# Returns:
- An empty `BloomFilter{Float64}`.
""" """
function empty_bloom_filter()::BloomFilter{Int} function empty_bloom_filter()::BloomFilter{Int}
return BloomFilter{Int}(size=DEFAULT_BLOOM_FILTER_SIZE, hash_count=3) return BloomFilter{Int}(size=DEFAULT_BLOOM_FILTER_SIZE, hash_count=3)
@ -202,6 +254,12 @@ end
AtomicFlag AtomicFlag
A thread-safe boolean flag using atomic operations. A thread-safe boolean flag using atomic operations.
# Fields:
- `value::Base.Threads.Atomic{Int}`: The atomic value of the flag.
# Arguments:
- `value::Base.Threads.Atomic{Int}`: The atomic value of the flag.
""" """
struct AtomicFlag struct AtomicFlag
value::Base.Threads.Atomic{Int} value::Base.Threads.Atomic{Int}

346
src/MSIData.jl
View File

@ -14,6 +14,11 @@ const FILE_HANDLE_LOCK = ReentrantLock()
ThreadSafeFileHandle ThreadSafeFileHandle
Wraps a file handle with thread-safe operations. Wraps a file handle with thread-safe operations.
# Fields
- `path::String`: The path to the file.
- `handle::IO`: The file handle.
- `lock::ReentrantLock`: The lock for thread-safe operations.
""" """
mutable struct ThreadSafeFileHandle mutable struct ThreadSafeFileHandle
path::String path::String
@ -21,6 +26,20 @@ mutable struct ThreadSafeFileHandle
lock::ReentrantLock lock::ReentrantLock
end end
"""
ThreadSafeFileHandle(path::String, mode::String="r")
Wraps a file handle with thread-safe operations.
# Arguments
- `path::String`: The path to the file.
- `mode::String`: The mode to open the file in.
# Returns
- A `ThreadSafeFileHandle`.
"""
function ThreadSafeFileHandle(path::String, mode::String="r") function ThreadSafeFileHandle(path::String, mode::String="r")
handle = lock(FILE_HANDLE_LOCK) do handle = lock(FILE_HANDLE_LOCK) do
open(path, mode) open(path, mode)
@ -285,12 +304,37 @@ mutable struct MSIData
end end
end end
# Thread-safe state management """
set_global_mz_range!(data::MSIData, min_mz::Float64, max_mz::Float64)
Sets the global m/z range for the MSIData object.
# Arguments
- `data::MSIData`: The MSIData object.
- `min_mz::Float64`: The minimum m/z value.
- `max_mz::Float64`: The maximum m/z value.
# Returns
- `nothing`
"""
function set_global_mz_range!(data::MSIData, min_mz::Float64, max_mz::Float64) function set_global_mz_range!(data::MSIData, min_mz::Float64, max_mz::Float64)
Threads.atomic_xchg!(data.global_min_mz, min_mz) Threads.atomic_xchg!(data.global_min_mz, min_mz)
Threads.atomic_xchg!(data.global_max_mz, max_mz) Threads.atomic_xchg!(data.global_max_mz, max_mz)
end end
"""
get_global_mz_range(data::MSIData)
Gets the global m/z range for the MSIData object.
# Arguments
- `data::MSIData`: The MSIData object.
# Returns
- `(min_mz, max_mz)`: The global m/z range.
"""
function get_global_mz_range(data::MSIData) function get_global_mz_range(data::MSIData)
# Atomic read using atomic_add! with zero # Atomic read using atomic_add! with zero
min_mz = Threads.atomic_add!(data.global_min_mz, 0.0) min_mz = Threads.atomic_add!(data.global_min_mz, 0.0)
@ -298,6 +342,20 @@ function get_global_mz_range(data::MSIData)
return (min_mz, max_mz) return (min_mz, max_mz)
end end
"""
set_analytics_data!(data::MSIData, stats_df::DataFrame, bloom_filters::Vector)
Sets the analytics data for the MSIData object.
# Arguments
- `data::MSIData`: The MSIData object.
- `stats_df::DataFrame`: The statistics DataFrame.
- `bloom_filters::Vector`: The bloom filters.
# Returns
- `nothing`
"""
function set_analytics_data!(data::MSIData, stats_df::DataFrame, bloom_filters::Vector) function set_analytics_data!(data::MSIData, stats_df::DataFrame, bloom_filters::Vector)
lock(data.cache_lock) do lock(data.cache_lock) do
data.spectrum_stats_df = stats_df data.spectrum_stats_df = stats_df
@ -305,12 +363,36 @@ function set_analytics_data!(data::MSIData, stats_df::DataFrame, bloom_filters::
end end
end end
"""
get_bloom_filters(data::MSIData)
Gets the bloom filters for the MSIData object.
# Arguments
- `data::MSIData`: The MSIData object.
# Returns
- `bloom_filters::Vector`: The bloom filters.
"""
function get_bloom_filters(data::MSIData) function get_bloom_filters(data::MSIData)
lock(data.cache_lock) do lock(data.cache_lock) do
return data.bloom_filters return data.bloom_filters
end end
end end
"""
get_spectrum_stats(data::MSIData)
Gets the spectrum statistics for the MSIData object.
# Arguments
- `data::MSIData`: The MSIData object.
# Returns
- `stats_df::DataFrame`: The statistics DataFrame.
"""
function get_spectrum_stats(data::MSIData) function get_spectrum_stats(data::MSIData)
lock(data.cache_lock) do lock(data.cache_lock) do
return data.spectrum_stats_df return data.spectrum_stats_df
@ -322,6 +404,13 @@ end
Explicitly closes the file handles associated with the `MSIData` object and clears the spectrum cache. Explicitly closes the file handles associated with the `MSIData` object and clears the spectrum cache.
It is good practice to call this method when you are finished with an `MSIData` object to release resources immediately. It is good practice to call this method when you are finished with an `MSIData` object to release resources immediately.
# Arguments
- `data::MSIData`: The MSIData object.
# Returns
- `nothing`
""" """
function Base.close(data::MSIData) function Base.close(data::MSIData)
if data.source isa ImzMLSource && isopen(data.source.ibd_handle) if data.source isa ImzMLSource && isopen(data.source.ibd_handle)
@ -343,6 +432,14 @@ end
Pre-loads a specified list of spectra into the cache. This is useful for warming up Pre-loads a specified list of spectra into the cache. This is useful for warming up
the cache with frequently accessed spectra to improve performance for subsequent the cache with frequently accessed spectra to improve performance for subsequent
interactive analysis. interactive analysis.
# Arguments
- `data::MSIData`: The MSIData object.
- `indices::AbstractVector{Int}`: The indices of the spectra to warm the cache with.
# Returns
- `nothing`
""" """
function warm_cache(data::MSIData, indices::AbstractVector{Int}) function warm_cache(data::MSIData, indices::AbstractVector{Int})
println("Warming cache with $(length(indices)) spectra...") println("Warming cache with $(length(indices)) spectra...")
@ -359,8 +456,12 @@ end
Calculates and returns a summary of the memory currently used by the MSIData object, Calculates and returns a summary of the memory currently used by the MSIData object,
including the spectrum cache and metadata. including the spectrum cache and metadata.
# Arguments
- `data::MSIData`: The MSIData object.
# Returns # Returns
- A `NamedTuple` with memory usage in bytes for different components. - `(total_mb, cache_mb, metadata_mb, coordinate_map_mb)`: A named tuple with memory usage in MB for different components.
""" """
function get_memory_usage(data::MSIData) function get_memory_usage(data::MSIData)
cache_bytes = 0 cache_bytes = 0
@ -488,6 +589,21 @@ function read_binary_vector(data::MSIData, io::IO, asset::SpectrumAsset)
return out_array return out_array
end end
"""
read_binary_vector(data::MSIData, ts_handle::ThreadSafeFileHandle, asset::SpectrumAsset)
Reads a single spectrum's m/z and intensity arrays directly from the `.ibd`
binary file for an `.imzML` dataset. It handles reading the raw binary data
and converting it from little-endian to the host's native byte order.
# Arguments
- `data`: The `MSIData` object.
- `ts_handle`: The `ThreadSafeFileHandle` containing the file handle and data formats.
- `asset`: The `SpectrumAsset` containing metadata for the array.
# Returns
- A `Vector` of the appropriate type containing the decoded data.
"""
function read_binary_vector(data::MSIData, ts_handle::ThreadSafeFileHandle, asset::SpectrumAsset) function read_binary_vector(data::MSIData, ts_handle::ThreadSafeFileHandle, asset::SpectrumAsset)
lock(ts_handle.lock) do lock(ts_handle.lock) do
return read_binary_vector(data, ts_handle.handle, asset) return read_binary_vector(data, ts_handle.handle, asset)
@ -679,6 +795,22 @@ single spectrum. Your logic should be inside the function passed to this helper.
For bulk processing of all spectra, consider using the `_iterate_spectra_fast` For bulk processing of all spectra, consider using the `_iterate_spectra_fast`
function, which is optimized for sequential access and avoids cache overhead. function, which is optimized for sequential access and avoids cache overhead.
# Arguments
- `f::Function`: The function to apply to the spectrum.
- `data::MSIData`: The MSIData object.
- `index::Int`: The index of the spectrum.
# Returns
- `result::Any`: The result of applying the function to the spectrum.
# Example
```julia
process_spectrum(data, 1) do mz, intensity
# Process the spectrum here
end
``` ```
""" """
function process_spectrum(f::Function, data::MSIData, index::Int) function process_spectrum(f::Function, data::MSIData, index::Int)
@ -696,6 +828,24 @@ end
A coordinate-based version of the "function barrier" helper. See `process_spectrum` A coordinate-based version of the "function barrier" helper. See `process_spectrum`
for details on the performance pattern. for details on the performance pattern.
# Arguments
- `f::Function`: The function to apply to the spectrum.
- `data::MSIData`: The MSIData object.
- `x::Int`: The x-coordinate of the spectrum.
- `y::Int`: The y-coordinate of the spectrum.
# Returns
- `result::Any`: The result of applying the function to the spectrum.
# Example
```julia
process_spectrum(data, 1, 1) do mz, intensity
# Process the spectrum here
end
```
""" """
function process_spectrum(f::Function, data::MSIData, x::Int, y::Int) function process_spectrum(f::Function, data::MSIData, x::Int, y::Int)
mz, intensity = GetSpectrum(data, x, y) mz, intensity = GetSpectrum(data, x, y)
@ -721,6 +871,13 @@ The computed statistics include:
Subsequent calls to functions like `get_total_spectrum` will be much faster Subsequent calls to functions like `get_total_spectrum` will be much faster
as they can use this cached data. This function modifies the `MSIData` object in-place as they can use this cached data. This function modifies the `MSIData` object in-place
and is idempotent. and is idempotent.
# Arguments
- `msi_data::MSIData`: The MSIData object.
# Returns
- `nothing`
""" """
function precompute_analytics(msi_data::MSIData) function precompute_analytics(msi_data::MSIData)
# Double-checked locking pattern # Double-checked locking pattern
@ -861,6 +1018,21 @@ It uses a fast, two-pass approach:
This function is highly optimized for `.imzML` by leveraging direct binary This function is highly optimized for `.imzML` by leveraging direct binary
reading and optimized binning logic. It is called by `get_total_spectrum`. reading and optimized binning logic. It is called by `get_total_spectrum`.
# Arguments
- `msi_data::MSIData`: The MSIData object.
- `num_bins::Int`: The number of bins to use for the total spectrum.
- `masked_indices::Union{Set{Int}, Nothing}`: The indices of the spectra to mask.
# Returns
- `(mz_range, total_spectrum)`: A tuple containing the m/z range and the total spectrum.
# Example
```julia
get_total_spectrum_imzml(msi_data)
```
""" """
function get_total_spectrum_imzml(msi_data::MSIData; num_bins::Int=2000, masked_indices::Union{Set{Int}, Nothing}=nothing) function get_total_spectrum_imzml(msi_data::MSIData; num_bins::Int=2000, masked_indices::Union{Set{Int}, Nothing}=nothing)
println("Calculating total spectrum for imzML (2-pass method)...") println("Calculating total spectrum for imzML (2-pass method)...")
@ -963,6 +1135,21 @@ It uses a two-pass approach analogous to the `imzML` implementation:
2. The second pass sums intensities into a pre-defined number of bins. 2. The second pass sums intensities into a pre-defined number of bins.
This function is called by `get_total_spectrum`. This function is called by `get_total_spectrum`.
# Arguments
- `msi_data::MSIData`: The MSIData object.
- `num_bins::Int`: The number of bins to use for the total spectrum.
- `masked_indices::Union{Set{Int}, Nothing}`: The indices of the spectra to mask.
# Returns
- `(mz_range, total_spectrum)`: A tuple containing the m/z range and the total spectrum.
# Example
```julia
get_total_spectrum_mzml(msi_data)
```
""" """
function get_total_spectrum_mzml(msi_data::MSIData; num_bins::Int=2000, masked_indices::Union{Set{Int}, Nothing}=nothing) function get_total_spectrum_mzml(msi_data::MSIData; num_bins::Int=2000, masked_indices::Union{Set{Int}, Nothing}=nothing)
println("Calculating total spectrum for mzML (2-pass method)...") println("Calculating total spectrum for mzML (2-pass method)...")
@ -1061,6 +1248,21 @@ This function dispatches to a specialized implementation based on the file type
(.imzML or .mzML) for optimal performance. (.imzML or .mzML) for optimal performance.
Returns a tuple containing two vectors: the binned m/z axis and the summed intensities. Returns a tuple containing two vectors: the binned m/z axis and the summed intensities.
# Arguments
- `msi_data::MSIData`: The MSIData object.
- `num_bins::Int`: The number of bins to use for the total spectrum.
- `mask_path::Union{String, Nothing}`: The path to the mask file.
# Returns
- `(mz_range, total_spectrum)`: A tuple containing the m/z range and the total spectrum.
# Example
```julia
get_total_spectrum(msi_data)
```
""" """
function get_total_spectrum(msi_data::MSIData; num_bins::Int=2000, mask_path::Union{String, Nothing}=nothing)::Tuple{Vector{Float64}, Vector{Float64}, Int} function get_total_spectrum(msi_data::MSIData; num_bins::Int=2000, mask_path::Union{String, Nothing}=nothing)::Tuple{Vector{Float64}, Vector{Float64}, Int}
local masked_indices::Union{Set{Int}, Nothing} = nothing local masked_indices::Union{Set{Int}, Nothing} = nothing
@ -1084,6 +1286,21 @@ Calculates the average of all spectra in the dataset by binning.
This is effectively the total spectrum divided by the number of spectra. This is effectively the total spectrum divided by the number of spectra.
Returns a tuple containing two vectors: the binned m/z axis and the averaged intensities. Returns a tuple containing two vectors: the binned m/z axis and the averaged intensities.
# Arguments
- `msi_data::MSIData`: The MSIData object.
- `num_bins::Int`: The number of bins to use for the total spectrum.
- `mask_path::Union{String, Nothing}`: The path to the mask file.
# Returns
- `(mz_range, average_spectrum)`: A tuple containing the m/z range and the average spectrum.
# Example
```julia
get_average_spectrum(msi_data)
```
""" """
function get_average_spectrum(msi_data::MSIData; num_bins::Int=2000, mask_path::Union{String, Nothing}=nothing)::Tuple{Vector{Float64}, Vector{Float64}} function get_average_spectrum(msi_data::MSIData; num_bins::Int=2000, mask_path::Union{String, Nothing}=nothing)::Tuple{Vector{Float64}, Vector{Float64}}
mz_bins, intensity_sum, num_spectra_processed = get_total_spectrum(msi_data, num_bins=num_bins, mask_path=mask_path) mz_bins, intensity_sum, num_spectra_processed = get_total_spectrum(msi_data, num_bins=num_bins, mask_path=mask_path)
@ -1105,6 +1322,19 @@ end
A stateful iterator for sequentially accessing spectra from an `MSIData` object. A stateful iterator for sequentially accessing spectra from an `MSIData` object.
This is created by the `IterateSpectra` function. This is created by the `IterateSpectra` function.
# Arguments
- `data::MSIData`: The MSIData object.
# Returns
- `MSIDataIterator`: An iterator that yields each spectrum.
# Example
```julia
IterateSpectra(data)
```
""" """
struct MSIDataIterator struct MSIDataIterator
data::MSIData data::MSIData
@ -1118,6 +1348,19 @@ This iterator is useful for processing all spectra in a loop and benefits from
the caching implemented in `GetSpectrum`. the caching implemented in `GetSpectrum`.
This function is the core of the "Event-driven" access pattern. This function is the core of the "Event-driven" access pattern.
# Arguments
- `data::MSIData`: The MSIData object.
# Returns
- `MSIDataIterator`: An iterator that yields each spectrum.
# Example
```julia
IterateSpectra(data)
```
""" """
function IterateSpectra(data::MSIData) function IterateSpectra(data::MSIData)
return MSIDataIterator(data) return MSIDataIterator(data)
@ -1251,6 +1494,17 @@ is significantly faster for bulk processing than reading spectra one by one.
- `f`: A function to execute for each spectrum, with the signature `f(index, mz_view, int_view)`. - `f`: A function to execute for each spectrum, with the signature `f(index, mz_view, int_view)`.
- `data`: The `MSIData` object. - `data`: The `MSIData` object.
- `source`: The `ImzMLSource`. - `source`: The `ImzMLSource`.
# Returns
- `nothing`
# Example
```julia
_iterate_uncompressed_fast(data, 1) do mz, intensity
# Process the spectrum here
end
```
""" """
function _iterate_uncompressed_fast(f::Function, data::MSIData, source::ImzMLSource, indices_to_iterate::Union{AbstractVector{Int}, Nothing}) function _iterate_uncompressed_fast(f::Function, data::MSIData, source::ImzMLSource, indices_to_iterate::Union{AbstractVector{Int}, Nothing})
# Optimized path for uncompressed data using buffer reuse # Optimized path for uncompressed data using buffer reuse
@ -1305,6 +1559,17 @@ and will be slower and allocate more memory than `_iterate_uncompressed_fast`.
- `f`: A function to execute for each spectrum, with the signature `f(index, mz_array, int_array)`. - `f`: A function to execute for each spectrum, with the signature `f(index, mz_array, int_array)`.
- `data`: The `MSIData` object. - `data`: The `MSIData` object.
- `source`: The `ImzMLSource`. - `source`: The `ImzMLSource`.
# Returns
- `nothing`
# Example
```julia
_iterate_compressed_fast(data, 1) do mz, intensity
# Process the spectrum here
end
```
""" """
function _iterate_compressed_fast(f::Function, data::MSIData, source::ImzMLSource, indices_to_iterate::Union{AbstractVector{Int}, Nothing}) function _iterate_compressed_fast(f::Function, data::MSIData, source::ImzMLSource, indices_to_iterate::Union{AbstractVector{Int}, Nothing})
# Path for datasets containing at least one compressed spectrum. # Path for datasets containing at least one compressed spectrum.
@ -1342,6 +1607,23 @@ This function acts as a dispatcher:
- If any spectrum is compressed, it uses a slower path that decompresses each spectrum individually. - If any spectrum is compressed, it uses a slower path that decompresses each spectrum individually.
- If all spectra are uncompressed, it uses a highly optimized path that reuses pre-allocated - If all spectra are uncompressed, it uses a highly optimized path that reuses pre-allocated
buffers to minimize memory allocations and overhead. buffers to minimize memory allocations and overhead.
# Arguments
- `f::Function`: A function to execute for each spectrum, with the signature `f(index, mz_array, int_array)`.
- `data::MSIData`: The `MSIData` object.
- `source::ImzMLSource`: The `ImzMLSource`.
# Returns
- `nothing`
# Example
```julia
_iterate_spectra_fast_impl(data, 1) do mz, intensity
# Process the spectrum here
end
```
""" """
function _iterate_spectra_fast_impl(f::Function, data::MSIData, source::ImzMLSource, indices_to_iterate::Union{AbstractVector{Int}, Nothing}) function _iterate_spectra_fast_impl(f::Function, data::MSIData, source::ImzMLSource, indices_to_iterate::Union{AbstractVector{Int}, Nothing})
if isempty(data.spectra_metadata) if isempty(data.spectra_metadata)
@ -1366,6 +1648,23 @@ Internal implementation of the fast iterator for `.mzML` files. It iterates
through each spectrum, decodes the Base64 data on the fly, and calls the through each spectrum, decodes the Base64 data on the fly, and calls the
provided function. It bypasses the `GetSpectrum` cache to avoid storing provided function. It bypasses the `GetSpectrum` cache to avoid storing
all decoded spectra in memory. all decoded spectra in memory.
# Arguments
- `f::Function`: A function to execute for each spectrum, with the signature `f(index, mz_array, int_array)`.
- `data::MSIData`: The `MSIData` object.
- `source::MzMLSource`: The `MzMLSource`.
# Returns
- `nothing`
# Example
```julia
_iterate_spectra_fast_impl(data, 1) do mz, intensity
# Process the spectrum here
end
```
""" """
function _iterate_spectra_fast_impl(f::Function, data::MSIData, source::MzMLSource, indices_to_iterate::Union{AbstractVector{Int}, Nothing}) function _iterate_spectra_fast_impl(f::Function, data::MSIData, source::MzMLSource, indices_to_iterate::Union{AbstractVector{Int}, Nothing})
# This implementation is for mzML. To improve disk I/O, we can reorder the read # This implementation is for mzML. To improve disk I/O, we can reorder the read
@ -1430,6 +1729,17 @@ number of available threads (`Threads.nthreads()`).
- `data`: The `MSIData` object. - `data`: The `MSIData` object.
- `indices_to_iterate`: An optional `AbstractVector{Int}` specifying a subset of - `indices_to_iterate`: An optional `AbstractVector{Int}` specifying a subset of
spectrum indices to iterate over. If `nothing`, it iterates over all spectra. spectrum indices to iterate over. If `nothing`, it iterates over all spectra.
# Returns
- `nothing`
# Example
```julia
_iterate_spectra_fast(data, 1) do mz, intensity
# Process the spectrum here
end
```
""" """
function _iterate_spectra_fast(f::Function, data::MSIData, indices_to_iterate::Union{AbstractVector{Int}, Nothing}=nothing) function _iterate_spectra_fast(f::Function, data::MSIData, indices_to_iterate::Union{AbstractVector{Int}, Nothing}=nothing)
# Default to all indices if not specified # Default to all indices if not specified
@ -1518,6 +1828,22 @@ end
Parallel version that processes all spectra in the dataset. Parallel version that processes all spectra in the dataset.
Convenience wrapper for the indexed version. Convenience wrapper for the indexed version.
# Arguments
- `f::Function`: A function to execute for each spectrum, with the signature `f(index, mz_array, int_array)`.
- `data::MSIData`: The `MSIData` object.
# Returns
- `nothing`
# Example
```julia
_iterate_spectra_fast_parallel(data, 1) do mz, intensity
# Process the spectrum here
end
```
""" """
function _iterate_spectra_fast_parallel(f::Function, data::MSIData) function _iterate_spectra_fast_parallel(f::Function, data::MSIData)
all_indices = collect(1:length(data.spectra_metadata)) all_indices = collect(1:length(data.spectra_metadata))
@ -1529,6 +1855,22 @@ end
Thread-safe version that processes all spectra in the dataset. Thread-safe version that processes all spectra in the dataset.
Convenience wrapper for the indexed version. Convenience wrapper for the indexed version.
# Arguments
- `f::Function`: A function to execute for each spectrum, with the signature `f(index, mz_array, int_array)`.
- `data::MSIData`: The `MSIData` object.
# Returns
- `nothing`
# Example
```julia
_iterate_spectra_fast_threadsafe(data, 1) do mz, intensity
# Process the spectrum here
end
```
""" """
function _iterate_spectra_fast_threadsafe(f::Function, data::MSIData) function _iterate_spectra_fast_threadsafe(f::Function, data::MSIData)
all_indices = collect(1:length(data.spectra_metadata)) all_indices = collect(1:length(data.spectra_metadata))

9
src/MSI_src.jl
View File

@ -84,6 +84,15 @@ using Setfield # For immutable struct updates
Opens a .mzML or .imzML file and prepares it for data access. Opens a .mzML or .imzML file and prepares it for data access.
This is the main entry point for the new data access API. This is the main entry point for the new data access API.
# Arguments
- `filepath::String`: Path to the .mzML or .imzML file to open.
- `cache_size::Int`: Number of spectra to cache in memory.
- `spectrum_type_map::Union{Dict{Int, Symbol}, Nothing}`: Optional mapping of spectrum indices to types.
# Returns
- An `MSIData` object.
""" """
function OpenMSIData(filepath::String; cache_size=300, spectrum_type_map::Union{Dict{Int, Symbol}, Nothing}=nothing) function OpenMSIData(filepath::String; cache_size=300, spectrum_type_map::Union{Dict{Int, Symbol}, Nothing}=nothing)
local msi_data local msi_data

10
src/MzmlConverter.jl
View File

@ -842,6 +842,16 @@ Main workflow function to convert a .mzML file to an .imzML file.
* `target_file`: Path for the output .imzML file (the .ibd will be named accordingly). * `target_file`: Path for the output .imzML file (the .ibd will be named accordingly).
* `img_width`: width dimention for the creation of the x axis * `img_width`: width dimention for the creation of the x axis
* `img_height`: height dimention for the creation of the y axis * `img_height`: height dimention for the creation of the y axis
# Returns
- `true` if the conversion was successful.
- `false` if the conversion failed.
# Example
```julia
ImportMzmlFile("test.mzML", "test.txt", "test.imzML")
```
""" """
function ImportMzmlFile(source_file::String, sync_file::String, target_file::String; img_width::Int=0, img_height::Int=0) function ImportMzmlFile(source_file::String, sync_file::String, target_file::String; img_width::Int=0, img_height::Int=0)
if !isfile(source_file) if !isfile(source_file)

28
src/ParserHelpers.jl
View File

@ -28,8 +28,18 @@ end
Reads a stream line-by-line until a line matches the provided regex. Reads a stream line-by-line until a line matches the provided regex.
# Arguments
- `stream::IO`: The stream to parse the metadata from.
- `regex::Regex`: The regex to match.
# Returns # Returns
- A `RegexMatch` object if a match is found, otherwise throws an error. - A `RegexMatch` object if a match is found, otherwise throws an error.
# Example
```julia
find_tag(open("test.mzML"), r"<spectrum")
```
""" """
function find_tag(stream, regex::Regex) function find_tag(stream, regex::Regex)
while !eof(stream) while !eof(stream)
@ -47,8 +57,18 @@ end
Retrieves an attribute's value from an XML tag string. Retrieves an attribute's value from an XML tag string.
# Arguments
- `source::String`: The string to parse the attribute from.
- `tag::String`: The tag to match.
# Returns # Returns
- A `RegexMatch` object containing the attribute and its value. - A `RegexMatch` object containing the attribute and its value.
# Example
```julia
get_attribute("<spectrum index=\"1\">", "index")
```
""" """
function get_attribute(source::AbstractString, tag::String = "([^=]+)") function get_attribute(source::AbstractString, tag::String = "([^=]+)")
# Construct the regex pattern string # Construct the regex pattern string
@ -65,10 +85,16 @@ struct. It reads accessions to determine the data format (e.g., `Float32`),
compression status (`zlib`), and axis type (m/z vs. intensity). compression status (`zlib`), and axis type (m/z vs. intensity).
# Arguments # Arguments
- `stream`: An IO stream positioned at the start of the `cvParam` block. - `stream::IO`: An IO stream positioned at the start of the `cvParam` block.
# Returns # Returns
- A `SpecDim` struct populated with the parsed configuration. - A `SpecDim` struct populated with the parsed configuration.
# Example
```julia
configure_spec_dim(open("test.mzML"))
```
""" """
function configure_spec_dim(stream) function configure_spec_dim(stream)
axis = SpecDim(Float64, false, 1, 0, UNKNOWN) # Add UNKNOWN as default mode axis = SpecDim(Float64, false, 1, 0, UNKNOWN) # Add UNKNOWN as default mode

283
src/imzML.jl
View File

@ -17,6 +17,21 @@ Core Functions:
# #
# ============================================================================ # ============================================================================
"""
parse_imzml_header(stream::IO)
Parses the header of an imzML file.
# Arguments
- `stream::IO`: The stream to parse.
# Returns
- `instrument_meta::InstrumentMetadata`: The instrument metadata.
- `param_groups::Dict{String, SpecDim}`: The parameter groups.
- `img_dims::Vector{Int32}`: The image dimensions.
"""
function parse_imzml_header(stream::IO) function parse_imzml_header(stream::IO)
# Initialize all return values and temporary variables # Initialize all return values and temporary variables
instrument_meta = InstrumentMetadata() instrument_meta = InstrumentMetadata()
@ -149,6 +164,14 @@ end
axes_config_img(stream) axes_config_img(stream)
Determines the storage order of the m/z and intensity arrays. Determines the storage order of the m/z and intensity arrays.
# Arguments
- `stream::IO`: The stream to parse.
# Returns
- `param_groups::Dict{String, SpecDim}`: The parameter groups.
""" """
function axes_config_img(stream::IO) function axes_config_img(stream::IO)
param_groups = Dict{String, SpecDim}() param_groups = Dict{String, SpecDim}()
@ -176,6 +199,14 @@ end
get_spectrum_tag_offset(stream) get_spectrum_tag_offset(stream)
Calculates the character offset within a `<spectrum>` tag, ignoring attribute values. Calculates the character offset within a `<spectrum>` tag, ignoring attribute values.
# Arguments
- `stream::IO`: The stream to parse.
# Returns
- `offset::Int`: The character offset.
""" """
function get_spectrum_tag_offset(stream::IO) function get_spectrum_tag_offset(stream::IO)
offset = position(stream) offset = position(stream)
@ -197,6 +228,20 @@ end
get_spectrum_attributes(stream, hIbd) get_spectrum_attributes(stream, hIbd)
Reads metadata to determine the byte offsets and data types for reading spectra. Reads metadata to determine the byte offsets and data types for reading spectra.
# Arguments
- `stream::IO`: The stream to parse.
- `hIbd::IO`: The IBD file handle.
# Returns
- `mz_first::Int`: The index of the first array
- `array_lengths::Vector{Int}`: The lengths of the arrays.
- `offset_mz::Int`: The byte offset of the mz array.
- `offset_intensity::Int`: The byte offset of the intensity array.
- `data_type_mz::DataType`: The data type of the mz array.
- `data_type_intensity::DataType`: The data type of the intensity array.
""" """
function get_spectrum_attributes(stream::IO, hIbd::IO) function get_spectrum_attributes(stream::IO, hIbd::IO)
# Look for position x # Look for position x
@ -230,7 +275,19 @@ function get_spectrum_attributes(stream::IO, hIbd::IO)
end end
# Helper function to read an entire spectrum block from the stream """
read_spectrum_block(stream::IO)
Helper function to read an entire spectrum block from the stream.
# Arguments
- `stream::IO`: The stream to parse.
# Returns
- `lines::Vector{String}`: The lines of the spectrum block.
"""
function read_spectrum_block(stream::IO) function read_spectrum_block(stream::IO)
lines = String[] lines = String[]
in_spectrum = false in_spectrum = false
@ -265,7 +322,20 @@ function read_spectrum_block(stream::IO)
end end
end end
# Helper function to parse coordinates from a spectrum block """
parse_coordinates(spectrum_block::String)
Helper function to parse coordinates from a spectrum block.
# Arguments
- `spectrum_block::String`: The spectrum block to parse.
# Returns
- `x::Int32`: The x-coordinate.
- `y::Int32`: The y-coordinate.
"""
function parse_coordinates(spectrum_block::String) function parse_coordinates(spectrum_block::String)
x = Int32(0) x = Int32(0)
y = Int32(0) y = Int32(0)
@ -279,7 +349,20 @@ function parse_coordinates(spectrum_block::String)
return x, y return x, y
end end
# Helper function to parse spectrum mode from a spectrum block """
parse_spectrum_mode(spectrum_block::String, global_mode::SpectrumMode)
Helper function to parse spectrum mode from a spectrum block.
# Arguments
- `spectrum_block::String`: The spectrum block to parse.
- `global_mode::SpectrumMode`: The global mode.
# Returns
- `mode::SpectrumMode`: The spectrum mode.
"""
function parse_spectrum_mode(spectrum_block::String, global_mode::SpectrumMode) function parse_spectrum_mode(spectrum_block::String, global_mode::SpectrumMode)
if occursin("MS:1000127", spectrum_block) if occursin("MS:1000127", spectrum_block)
return CENTROID return CENTROID
@ -290,7 +373,26 @@ function parse_spectrum_mode(spectrum_block::String, global_mode::SpectrumMode)
end end
end end
# Helper function to parse binary data arrays """
parse_binary_data_arrays(spectrum_block::String, spectrum_start_line::String,
default_mz_format::DataType, default_intensity_format::DataType,
mz_is_compressed::Bool, int_is_compressed::Bool)
Helper function to parse binary data arrays from a spectrum block.
# Arguments
- `spectrum_block::String`: The spectrum block to parse.
- `spectrum_start_line::String`: The start line of the spectrum.
- `default_mz_format::DataType`: The default data type for mz arrays.
- `default_intensity_format::DataType`: The default data type for intensity arrays.
- `mz_is_compressed::Bool`: Whether mz arrays are compressed.
- `int_is_compressed::Bool`: Whether intensity arrays are compressed.
# Returns
- `array_data::Vector{Tuple{Bool, Int32, Int64, Int64}}`: The array data.
"""
function parse_binary_data_arrays(spectrum_block::String, spectrum_start_line::String, function parse_binary_data_arrays(spectrum_block::String, spectrum_start_line::String,
default_mz_format::DataType, default_intensity_format::DataType, default_mz_format::DataType, default_intensity_format::DataType,
mz_is_compressed::Bool, int_is_compressed::Bool) mz_is_compressed::Bool, int_is_compressed::Bool)
@ -342,8 +444,32 @@ function parse_binary_data_arrays(spectrum_block::String, spectrum_start_line::S
return current_array_data return current_array_data
end end
"""
parse_imzml_spectrum_block(stream::IO, hIbd::Union{IO, ThreadSafeFileHandle}, param_groups::Dict{String, SpecDim},
width::Int32, height::Int32, num_spectra::Int32,
default_mz_format::DataType, default_intensity_format::DataType,
mz_is_compressed::Bool, int_is_compressed::Bool, global_mode::SpectrumMode)
# Unified function to parse spectrum metadata Unified function to parse spectrum metadata
# Arguments
- `stream::IO`: The stream to parse.
- `hIbd::Union{IO, ThreadSafeFileHandle}`: The IBD file handle.
- `param_groups::Dict{String, SpecDim}`: The parameter groups.
- `width::Int32`: The width of the image.
- `height::Int32`: The height of the image.
- `num_spectra::Int32`: The number of spectra.
- `default_mz_format::DataType`: The default data type for mz arrays.
- `default_intensity_format::DataType`: The default data type for intensity arrays.
- `mz_is_compressed::Bool`: Whether mz arrays are compressed.
- `int_is_compressed::Bool`: Whether intensity arrays are compressed.
- `global_mode::SpectrumMode`: The global mode.
# Returns
- `spectra_metadata::Vector{SpectrumMetadata}`: The spectrum metadata.
"""
function parse_imzml_spectrum_block(stream::IO, hIbd::Union{IO, ThreadSafeFileHandle}, param_groups::Dict{String, SpecDim}, function parse_imzml_spectrum_block(stream::IO, hIbd::Union{IO, ThreadSafeFileHandle}, param_groups::Dict{String, SpecDim},
width::Int32, height::Int32, num_spectra::Int32, width::Int32, height::Int32, num_spectra::Int32,
default_mz_format::DataType, default_intensity_format::DataType, default_mz_format::DataType, default_intensity_format::DataType,
@ -412,6 +538,23 @@ function parse_imzml_spectrum_block(stream::IO, hIbd::Union{IO, ThreadSafeFileHa
return spectra_metadata return spectra_metadata
end end
"""
load_imzml_lazy(file_path::String; cache_size::Int=100)
The main entry point for loading an imzML file.
It works by locating the .ibd file and parsing the .imzML file
Lazily loads the data and then creating a MSIData object with the
parsed information acquired by the helper functions.
# Arguments
- `file_path::String`: The path to the imzML file.
- `cache_size::Int`: The size of the cache.
# Returns
- `msi_data::MSIData`: The MSI data.
"""
function load_imzml_lazy(file_path::String; cache_size::Int=100) function load_imzml_lazy(file_path::String; cache_size::Int=100)
println("DEBUG: Checking for .imzML file at $file_path") println("DEBUG: Checking for .imzML file at $file_path")
if !isfile(file_path) if !isfile(file_path)
@ -534,6 +677,32 @@ function load_imzml_lazy(file_path::String; cache_size::Int=100)
end end
end end
"""
parse_compressed(stream::IO, hIbd::Union{IO, ThreadSafeFileHandle}, param_groups::Dict{String, SpecDim},
width::Int32, height::Int32, num_spectra::Int32,
default_mz_format::DataType, default_intensity_format::DataType,
mz_is_compressed::Bool, int_is_compressed::Bool, global_mode::SpectrumMode)
Helper function to parse spectrum metadata from an imzML file.
# Arguments
- `stream::IO`: The stream to parse.
- `hIbd::Union{IO, ThreadSafeFileHandle}`: The IBD file handle.
- `param_groups::Dict{String, SpecDim}`: The parameter groups.
- `width::Int32`: The width of the image.
- `height::Int32`: The height of the image.
- `num_spectra::Int32`: The number of spectra.
- `default_mz_format::DataType`: The default data type for mz arrays.
- `default_intensity_format::DataType`: The default data type for intensity arrays.
- `mz_is_compressed::Bool`: Whether mz arrays are compressed.
- `int_is_compressed::Bool`: Whether intensity arrays are compressed.
- `global_mode::SpectrumMode`: The global mode.
# Returns
- `spectra_metadata::Vector{SpectrumMetadata}`: The spectrum metadata.
"""
function parse_compressed(stream::IO, hIbd::Union{IO, ThreadSafeFileHandle}, param_groups::Dict{String, SpecDim}, function parse_compressed(stream::IO, hIbd::Union{IO, ThreadSafeFileHandle}, param_groups::Dict{String, SpecDim},
width::Int32, height::Int32, num_spectra::Int32, width::Int32, height::Int32, num_spectra::Int32,
default_mz_format::DataType, default_intensity_format::DataType, default_mz_format::DataType, default_intensity_format::DataType,
@ -690,6 +859,13 @@ end
Finds the intensity of the most intense peak within a mass tolerance window. Finds the intensity of the most intense peak within a mass tolerance window.
This optimized version uses binary search for efficiency. This optimized version uses binary search for efficiency.
# Arguments
- `mz_array::AbstractVector{<:Real}`: The m/z array.
- `intensity_array::AbstractVector{<:Real}`: The intensity array.
- `target_mass::Real`: The target mass.
- `tolerance::Real`: The tolerance.
# Returns # Returns
- The intensity (`Float64`) of the peak if found, otherwise `0.0`. - The intensity (`Float64`) of the peak if found, otherwise `0.0`.
""" """
@ -727,6 +903,13 @@ end
Extracts an image slice for a given m/z value without plotting. Extracts an image slice for a given m/z value without plotting.
This is a performant function that iterates through spectra once. This is a performant function that iterates through spectra once.
# Arguments
- `data::MSIData`: The MSIData object.
- `mass::Real`: The target mass.
- `tolerance::Real`: The tolerance.
- `mask_path::Union{String, Nothing}`: The path to the mask.
# Returns # Returns
- A `Matrix{Float64}` representing the intensity slice. - A `Matrix{Float64}` representing the intensity slice.
""" """
@ -811,6 +994,13 @@ end
Extracts multiple image slices for a given list of m/z values in a single pass. Extracts multiple image slices for a given list of m/z values in a single pass.
This is a highly performant function that iterates through the full dataset only once. This is a highly performant function that iterates through the full dataset only once.
# Arguments
- `data::MSIData`: The MSIData object.
- `masses::AbstractVector{<:Real}`: The target masses.
- `tolerance::Real`: The tolerance.
- `mask_path::Union{String, Nothing}`: The path to the mask.
# Returns # Returns
- A `Dict{Real, Matrix{Float64}}` mapping each mass to its intensity slice matrix. - A `Dict{Real, Matrix{Float64}}` mapping each mass to its intensity slice matrix.
""" """
@ -916,6 +1106,19 @@ end
Generates and saves a plot of a single image slice for a given m/z value. Generates and saves a plot of a single image slice for a given m/z value.
This function closely imitates the logic of the original `GetSlice` but uses This function closely imitates the logic of the original `GetSlice` but uses
the modern `MSIData` access patterns and robust peak finding. the modern `MSIData` access patterns and robust peak finding.
# Arguments
- `msi_data::MSIData`: The MSIData object.
- `mass::Real`: The target mass.
- `tolerance::Real`: The tolerance.
- `output_dir::String`: The output directory.
- `stage_name::String`: The name of the stage.
- `bins::Int`: The number of bins.
- `mask_path::Union{String, Nothing}`: The path to the mask.
# Returns
- `fig::Figure`: A figure with the data plotted with colors applied.
""" """
function plot_slice(msi_data::MSIData, mass::Real, tolerance::Real, output_dir::String; function plot_slice(msi_data::MSIData, mass::Real, tolerance::Real, output_dir::String;
stage_name="slice_mz_$(mass)", bins=256, mask_path::Union{String, Nothing}=nothing) stage_name="slice_mz_$(mass)", bins=256, mask_path::Union{String, Nothing}=nothing)
@ -1115,6 +1318,16 @@ end
norm_slices_hist(slices, bins; prob=0.98) norm_slices_hist(slices, bins; prob=0.98)
Normalizes a set of image slices based on a shared histogram range. Normalizes a set of image slices based on a shared histogram range.
# Arguments
- `slices`: The image slices.
- `bins`: The number of bins.
- `prob`: The target cumulative probability.
# Returns
- `norm_img`: The normalized image slices.
- `mass_bounds`: The bounds for each mass.
""" """
function norm_slices_hist(slices, bins; prob=0.98) function norm_slices_hist(slices, bins; prob=0.98)
n_files, n_masses = size(slices) n_files, n_masses = size(slices)
@ -1145,7 +1358,17 @@ end
Linearly scales the intensity values in a slice to a specified number of levels. Linearly scales the intensity values in a slice to a specified number of levels.
The output is an array of `UInt8` values. The output is an array of `UInt8` values.
This is a modernized version of `IntQuantCl`. This is a modernized version of the original `IntQuantCl`. by Nacho
# Arguments
- `slice`: The image slice.
- `levels`: The number of levels.
- `mask_matrix`: The mask matrix.
# Returns
- `quantized`: The quantized image slice.
- `bounds`: The bounds for each mass.
""" """
function quantize_intensity(slice::AbstractMatrix{<:Real}, levels::Integer=256; mask_matrix::Union{BitMatrix, Nothing}=nothing) function quantize_intensity(slice::AbstractMatrix{<:Real}, levels::Integer=256; mask_matrix::Union{BitMatrix, Nothing}=nothing)
local max_val local max_val
@ -1198,6 +1421,16 @@ end
Reduces the number of points in a spectrum for faster plotting, while preserving peaks. Reduces the number of points in a spectrum for faster plotting, while preserving peaks.
It divides the m/z range into `n_points` bins and keeps only the most intense point from each bin. It divides the m/z range into `n_points` bins and keeps only the most intense point from each bin.
# Arguments
- `mz`: The m/z values.
- `intensity`: The intensity values.
- `n_points`: The number of points.
# Returns
- `mz`: The downsampled m/z values.
- `intensity`: The downsampled intensity values.
""" """
function downsample_spectrum(mz::AbstractVector, intensity::AbstractVector, n_points::Integer=2000) function downsample_spectrum(mz::AbstractVector, intensity::AbstractVector, n_points::Integer=2000)
if isempty(mz) || length(mz) <= n_points if isempty(mz) || length(mz) <= n_points
@ -1375,6 +1608,15 @@ end
Saves an 8-bit indexed image as a BMP file. Saves an 8-bit indexed image as a BMP file.
This is a modernized version of `SaveBitmapCl`. This is a modernized version of `SaveBitmapCl`.
# Arguments
- `name`: The name of the file.
- `pixMap`: The pixel map.
- `colorTable`: The color table.
# Returns
- `nothing`
""" """
function save_bitmap(name::String, pixMap::Matrix{UInt8}, colorTable::Vector{UInt32}) function save_bitmap(name::String, pixMap::Matrix{UInt8}, colorTable::Vector{UInt32})
# Get image dimensions # Get image dimensions
@ -1443,6 +1685,14 @@ end
generate_palette(colorscheme, n_colors=256) generate_palette(colorscheme, n_colors=256)
Generates a BMP-compatible UInt32 color palette from a ColorScheme. Generates a BMP-compatible UInt32 color palette from a ColorScheme.
# Arguments
- `colorscheme`: The colorscheme.
- `n_colors`: The number of colors.
# Returns
- `palette`: The color palette.
""" """
function generate_palette(colorscheme, n_colors=256) function generate_palette(colorscheme, n_colors=256)
palette = Vector{UInt32}(undef, n_colors) palette = Vector{UInt32}(undef, n_colors)
@ -1461,6 +1711,27 @@ end
# Viridis color palette (256 colors) - defined as constant # Viridis color palette (256 colors) - defined as constant
const ViridisPalette = generate_palette(ColorSchemes.viridis) const ViridisPalette = generate_palette(ColorSchemes.viridis)
"""
generate_colorbar_image(slice_data::AbstractMatrix, color_levels::Int, output_path::String,
bounds::Tuple{Float64, Float64};
use_triq::Bool=false, triq_prob::Float64=0.98,
mask_path::Union{String, Nothing}=nothing)
Generates a colorbar image for a given slice of data.
# Arguments
- `slice_data`: The slice data.
- `color_levels`: The number of color levels.
- `output_path`: The output path.
- `bounds`: The bounds for the colorbar.
- `use_triq`: Whether to use triq.
- `triq_prob`: The probability of using triq.
- `mask_path`: The path to the mask.
# Returns
- `fig::Figure`: A figure with the colorbar.
"""
function generate_colorbar_image(slice_data::AbstractMatrix, color_levels::Int, output_path::String, function generate_colorbar_image(slice_data::AbstractMatrix, color_levels::Int, output_path::String,
bounds::Tuple{Float64, Float64}; bounds::Tuple{Float64, Float64};
use_triq::Bool=false, triq_prob::Float64=0.98, use_triq::Bool=false, triq_prob::Float64=0.98,

31
src/mzML.jl
View File

@ -19,6 +19,24 @@ const DATA_FORMAT_ACCESSIONS = Dict{String, DataType}(
"MS:1000523" => Float64 "MS:1000523" => Float64
) )
"""
parse_instrument_metadata_mzml(stream::IO)
Parses instrument metadata from an mzML file.
# Arguments
- `stream::IO`: The stream to parse the metadata from.
# Returns
- `InstrumentMetadata`: The instrument metadata.
# Example
```julia
parse_instrument_metadata_mzml(open("test.mzML"))
```
"""
function parse_instrument_metadata_mzml(stream::IO) function parse_instrument_metadata_mzml(stream::IO)
println("DEBUG: Starting mzML instrument metadata parsing...") println("DEBUG: Starting mzML instrument metadata parsing...")
# Initialize with default values from the InstrumentMetadata constructor # Initialize with default values from the InstrumentMetadata constructor
@ -328,6 +346,19 @@ end
Finds the index offset in an mzML file by reading from the end. Finds the index offset in an mzML file by reading from the end.
Optimized version with better buffer management. Optimized version with better buffer management.
# Arguments
- `stream::IO`: The stream to parse the metadata from.
# Returns
- `Int64`: The index offset.
# Example
```julia
find_index_offset(open("test.mzML"))
```
""" """
function find_index_offset(stream::IO)::Int64 function find_index_offset(stream::IO)::Int64
file_size = filesize(stream) file_size = filesize(stream)

1
start_MSI_GUI.jl
View File

@ -14,6 +14,7 @@ if get(ENV, "GENIE_ENV", "dev") != "prod" && !isfile(manifest_path)
Pkg.gc() Pkg.gc()
elseif get(ENV, "GENIE_ENV", "dev") != "prod" elseif get(ENV, "GENIE_ENV", "dev") != "prod"
@info "Manifest.toml found. Skipping Pkg.instantiate() for faster boot. Delete Manifest.toml if you need to re-instantiate." @info "Manifest.toml found. Skipping Pkg.instantiate() for faster boot. Delete Manifest.toml if you need to re-instantiate."
ENV["GENIE_ENV"] = "prod"
end end
using Genie using Genie