General
All chemicals were obtained from Sigma-Aldrich unless specified otherwise.
Vesicle isolation
A. californica (150–200 g) were obtained from the National Resource for Aplysia and kept in a 14 °C aquarium filled with Instant Ocean (Aquarium Systems). Animal euthanasia was performed in accordance with the AVMA Guidelines for the Euthanasia of Animals: 2020 Edition (Section S6.3.1.1 Noninhaled agents for immersion). Animals were anesthetized using a 333 mM MgCl2 solution injected into the body cavity (50% volume/body).
For DCV isolation, the AG was isolated by manual dissection, and placed into a microcentrifuge tube containing 1 ml of ASW. Due to the holocrine release mechanism, DCV secretion from the AG was induced by gentle trituration with a polypropylene Pasteur pipette, releasing intact DCVs into the ASW solution. Next, 100 μl of 500 mM ammonium acetate solution was prespotted onto an ITO-unpolished float glass slide, Surface Resistivity = 70–100 Ω (Delta Technologies). Then, 50 μl of vesicle-containing ASW solution was pipetted into the 100 μl of ammonium acetate solution on the ITO-glass slide for a total volume of 150 μl. The resulting vesicle solution was then rinsed with 500 mM ammonium acetate with simultaneous aspiration of the solution, leaving DCVs seeded across the slide. All steps were visually monitored using an inverted microscope. Each slide held DCVs from three animals (biological replicates) and a total of three slides (technical replicates) were prepared. For LV isolation, the procedure used for the DCVs was repeated using red hemiduct tissue, with the other steps being the same as for the DCV isolation; 123 LVs were analyzed from the same three biological replicates used for the DCV isolation.
Vesicle imaging
Brightfield images were acquired on an Axio Imager M2 (Zeiss) equipped with an AxioCam ICc 5, a ×0.63 camera adapter and a transmitted light visible light-light emitting diode lamp. Images were acquired in mosaic mode using a ×10 objective with 30% tile overlap. The resulting tiles were stitched before exporting in TIFF-file format using ZEN 2.0 Pro edition (Zeiss) software.
Scanning electron microscopy images were acquired using an FEI Quanta FEG 450 environmental scanning electron microscope (FEI). Images were acquired using an accelerating voltage of 10 kV, dwell time of 10 µs and working distance of 6.6 mm.
Image-processing for vesicle recognition in microMS
We utilized biological image analysis software, including ImageJ26 and CellProfiler27, for vesicle pixel recognition on the brightfield microscopy image. Importantly, other image-processing software can be used for this process as well. CellProfiler pipeline modules were created to selectively mask the DCVs and LVs on the brightfield microscopy image to remove the background and retain just information on the vesicles of interest as the foreground. Masking the identified vesicles creates a binary image marking the pixel locations of vesicles on the glass slide. The microMS software was used to translate the pixel locations on microscopy images to physical coordinates of the mass spectrometer’s stage. A 200-µm-distance filter was then applied in microMS (removing vesicles located closer than 200 µm to each other from the target list). Alternatively, ImageJ was used as a simple thresholding strategy by measuring the difference between the vesicle (foreground) and background pixel intensity, which allows a threshold to be set, leaving only objects of interest for recognition in microMS.
Matrix application
Matrix deposition for MALDI MS analysis was performed using a glass sublimation apparatus (Wilmad-LabGlass) filled with 2,5-dihydroxybenzoic acid as the MALDI matrix. The slide was attached to the cold-finger and vacuum was created using a rotary vane pump (Edwards Vacuum, model E2M30). The sublimation apparatus was placed on a sand bath preheated to 150 °C for 8 min. Matrix deposition was followed by recrystallization using 5% methanol. Recrystallization was performed using a 100 × 15-mm2 polystyrene Petri dish as a recrystallization chamber. The ITO-glass slide was attached to the top of the recrystallization chamber and a filter paper (Whatman Grade 1 Qualitative Filter Paper, Thermo Fisher Scientific) was wetted with 1 ml of 5% methanol. The chamber was sealed using tape and placed in an oven at 85 °C for 1.5 min. After removal of the recrystallization chamber from the oven, the slide was immediately removed from the recrystallization chamber and allowed to dry in a nitrogen chamber until analysis.
MALDI MS measurements
High-throughput single-DCV and single-LV analyses were performed on a SolariX XR 7T Fourier-transform ion cyclotron resonance mass spectrometer equipped with an APOLLO II dual MALDI/ESI source (Bruker) using an m/z range of 150–4,500. Data were acquired at 1 M giving a mass resolution of 107,000 at m/z 535 and 19,070 at m/z 3,922, yielding a transient length of 0.721 s. The instrument was operating in positive-mode using a Smartbeam-II UV laser (Bruker) set to ‘Ultra mode’, which yields a 100-µm-diameter laser footprint. Each MALDI acquisition consisted of two accumulations comprised of 400 laser shots each, at a frequency of 1,000 Hz. DCV and LV stage coordinates and geometry files were generated using microMS as previously described8.
Data preprocessing
Data preprocessing was performed using Compass Data Analysis 4.4.2 (Bruker) and MATLAB 2018b (MathWorks). Internal calibration using the exact mass of AG peptides was performed using AGPB1 [71–81] (m/z 1,221.6878), AGPA1 [22–34] (m/z 1,396.7225), AGPA1 [156–173] (m/z 1,908.8761), AGPA2 [33–69] (m/z 3,922.9478) and AGPB1 [85–118] (m/z 4,031.0053). Peak picking and peak export for statistical analysis in MATLAB 2018b were set to a signal-to-noise ratio of 5 with a relative intensity threshold of 0.01%. A nonuniform bin width was used for mass spectral alignment. For DCV data analysis, mass features were truncated at m/z 1,100 for downstream data analysis of only AG peptides. For LV and DCV vesicle classification, mass features in m/z 500–1,100 were selected for downstream data analysis. Internal calibration using AG peptides was not performed for the LV and DCV classification tasks. Using the target list provided by microMS, the pixel coordinates were used to find and crop the locations of individual vesicles across the original microscopy image. The corresponding mass spectra were matched with the appropriate vesicle for visual evaluation of corresponding single vesicles.
CX decomposition and statistical analysis
Unsupervised approaches such as principal component analysis (PCA) and, more generally, matrix decomposition or factorization, are valuable data analysis tools to enable field-specific interpretation of high-dimensional datasets. However, the interpretation is limited due to the complicated eigenspace obtained from PCA, deterring our further understanding of the feature space of the data matrix. CX matrix decomposition is designed to obtain a low-rank approximation of the data in terms of actual rows or columns14. Given an m × n data matrix A, the algorithm decomposes it into an m × c matrix C and a c × n matrix X, where C is expressed by c number of column vectors of the original data. The statistical leverage scores are used to rank and select the columns of C from A, which can be obtained by
$$l_j = mathop {sum }limits_{i = 1}^k v_{ji}^2$$
where lj is the leverage score for the jth column/feature, vi is the right singular vectors obtained by the singular value decomposition and k is the rank to be selected. X is then determined by minimizing the error:
$$mathop {{min }}limits_X left| {left| {A – CX} right|} right|_F$$
The reconstruction error evaluation and the rank k selection are provided as Supplementary Fig. 8. Based on the evaluation, the top 200 features are selected to form the columns of the matrix C. k-means clustering, and the Wilcoxon rank-sum test, were performed using the Python-based open source package SCANPY28. The stacked violin plots, shown in Fig. 2d, of the normalized signal intensities represent the identified peak features across the three clusters, with the vesicle type (or cluster) assignments obtained by k-means clustering. The y axis is the root-mean-squared-normalized peak intensity and each violin contains the distribution of the normalized intensities of the corresponding features in the x axis.
Peptide sequencing by LC–MS/MS
Peptide extracts (n = 3) were obtained by manually grinding entire AG tissue in 500 μl of acidified methanol followed by evaporation and reconstitution of each extract in 0.1% formic acid. A nanoElute (Bruker) ultra-high-pressure nano-flow chromatography system was coupled to a trapped ion mobility–quadrupole time-of-flight mass spectrometer (timsTOF Pro, Bruker) with a CaptiveSpray nano-electrospray ion source (Bruker) equipped with an external column oven. Mobile phases A and B were water with 0.1% formic acid (v/v) and acetonitrile with 0.1% formic acid (v/v), respectively. Samples were loaded onto a precolumn peptide trap (Acclaim PepMap 100 C18, 1 × 5 mm2, 5-µm particle size, Thermo Fisher Scientific) using solvent A for off-line desalting. Next, the trap was placed in-line with the analytical column and peptide separation was performed at 40 °C with a uniform flow of 300 nl min−1 on a C18 ReproSil AQ column (Bruker FIFTEEN, P/N no. 1842621: 150 mm × 75 µm, 1.9-µm particle size, pore size 120 Å) equilibrated at 2% B. A linear gradient of solvent B was applied as follows: 2–10% within 5 min, 10–50% in the next 115 min, followed by a washing step at 95% B and re-equilibration, during which data collection was not performed. The mass spectrometer was operated in parallel accumulation–serial fragmentation (PASEF) mode for peptide sequencing. The mass range for the precursor ion was set to m/z 100–1,700, ion mobility 1/K0 range 0.6–1.6 V s cm−2. Fragmentation was performed with 10 PASEF scans, cycle time of 1.1 s, during which collision energy varied linearly between 20 and 59 eV depending on precursor 1/K0 value within the set range. Active dynamic exclusion of precursor ions was set to 0.4 min.
Peptide identification by bioinformatics
MS raw files were processed with PEAKS Online29 (Bioinformatics Solution) using the DeNovo, database (DB) and post-translational modification (PTM) protocols, sequentially. Peptide sequence tags obtained by the DeNovo process (80% average local confidence score cut-off) were searched against the Aplysia RefSeq database available from the NCBI (GCF_000002075.1). The protein database was filtered to include proteins up to 1,000 amino acids long. Search parameters for DB included: parent mass error tolerance 20.0 ppm, fragment mass error tolerance 0.03 Da, no enzyme, digest mode unspecific. Next, PTMs were identified by searching the data for amidation, acetylation (K) and acetylation N terminus, pyro-glutamylation from (Q) and (E), phosphorylation (STY) and half-disulfide bridge. The false discovery rate was determined by decoy fusion method and the threshold set to 1% for peptides. Proteins with −log10 P > 20 and at least one unique peptide are reported.
Vesicle classification through machine learning
We adapted a previously described machine learning strategy6 for vesicle type classification to differentiate between DCVs and LVs based on their lipid contents. Features in the m/z range of 500–1,100 were used for the classification task. Gradient boosting trees were trained with a threefold validation; in each fold of the model performance, the metrics were computed to obtain classification accuracies, confusion matrices and receiver operating characteristic curves (Supplementary Fig. 11). The most contributing features to the classification task were selected via SHAP, a game theory approach for model explanations, and were obtained through a Python implementation of SHAP. A total of 97 features with nonzero mean SHAP values from the output of trained models were selected and 36 were putatively annotated by searching against the LIPID MAPS30 database with a 7-ppm tolerance. Models were then retrained with the SHAP-selected features as well as the annotated features to verify the discriminative ability of those features. t-SNE using the cosine distance was used to visualize the lipid differences between DCVs and LVs in a low-dimensional space.
Reporting Summary
Further information on research design is available in the Nature Research Reporting Summary linked to this article.

