In support of local lobstermen and the New York State Department of Environmental Conservation, American lobsters were collected from WLIS in July 2019, 50 miles south of Montauk in October 2019, and ELIS in August 2020 (Figure 1A; Supplementary Table S1). The lobsters were immediately placed in a cooler with ice packs to keep them alive and delivered to the lab within 2 h. Healthy and diseased lobsters were kept separately to avoid cross-contamination. Weight and length measurements were taken, and the lobsters were then dissected to collect various tissues, such as the carapace, green gland, hepatopancreas, intestine, and testis. Following a gentle wash with phosphate-buffered saline (PBS) at pH 7.4, the samples were stored in RNAlater stabilization solution (Qiagen) at −80°C until microbial genomic DNA extraction. The length from the rear of the eye socket to the end of the carapace was an average of 25.4 cm (23.5–27.0 cm) and the weight an average of 632.4 g (528.5–763.5 g). Lobsters with shell lesions regardless of severity were classified as diseased (ESD, n = 8), and lobsters with no apparent signs of carapace lesions were considered healthy (HTH, n = 10). Metadata for sample information is available in Supplementary Table S1.

Microbial genomic DNA was extracted from both shells and internal organs of lobsters using the DNeasy Blood & Tissue kit (Qiagen) according to the manufacturer’s protocols, with a minor adjustment in the pretreatment of the shells. Briefly, shells in RNAlater stabilization solution were cut into small square pieces, placed in a microcentrifuge tube in 320 μl of PBS, and incubated on ice for 1 h to help detach bacteria from the shells. The incubated shells in PBS were disrupted and homogenized using a handheld homogenizer (Bel-Art). The resulting homogenizing solution, which contains bacteria from the surface of shells as well as bacteria present in the shells, was subsequently used for DNA extraction. For internal organs, 10–25 mg of the hepatopancreas, intestine, green gland, and testis were used for the extraction of DNA according to the manufacturer’s instructions. The DNA concentration and purity from all samples were measured by a NanoDrop One^C spectrophotometer (Thermo Scientific).

16S rRNA V3-V4 amplicon sequencing library preparation and Illumina MiSeq sequencing were conducted at GENEWIZ, Inc (South Plainfield, NJ, United States). Briefly, microbial genomic DNA for each sample was normalized to 20 ng/μl using a Qubit 2.0 Fluorometer (Invitrogen) and was used to generate amplicons using a MetaVx™ Library Preparation kit (GENEWIZ). The V3-V4 region of the 16S rRNA gene was amplified using the forward and reverse primers (341 Forward: CCTACGGRRBGCASCAGKVRVGAAT, 806 Reverse: GGACTACNVGGGTWTCTAATCC). Second limited-cycle amplification was performed to add multiplexing indices (barcode) and Illumina sequencing adapters. DNA libraries were validated by Agilent 2100 Bioanalyzer (Agilent Technologies) and quantified by Qubit 2.0 Fluorometer (Invitrogen). DNA libraries were multiplexed and loaded on an Illumina MiSeq instrument according to manufacturer’s instructions (Illumina). Sequencing was performed using the 2×250 paired-end configuration, and image analysis and base calling were conducted by the MiSeq Control Software embedded in the Illumina MiSeq instrument.

Two carapace samples from ESD and one from HTH were excluded from the analysis due to sequencing failure. Detailed information on the sequencing data statistics and quality is provided in Supplementary Table S2. 16S rRNA data analysis was completed with the QIIME2 package (2022.v8) in Python where forward and reverse reads were joined and assigned to samples based on their barcodes (Bolyen et al., 2019). After removing the barcode and primer sequences using qiime2/q2-cutadapt, quality filtering on joined sequences was performed. Sequences were discarded if they did not meet the following the parameters: Sequence length > 200 bp, min quality = PHRED 4, quality-window = 3, and min length fraction = 0.75. Sequences that passed this quality filtering were used to perform the next bioinformatics analysis. The QIIME package was used to perform operational taxonomic units (OTUs) picking, and the sequences were clustered into OTUs using the VSEARCH plugin with a clustering threshold of 99% sequence identity against the pre-clustered Silva 138 database (min merge length = 250; Rognes et al., 2016). A taxonomic category was assigned to each OTU with a confidence threshold of 0.7 using the Ribosomal Database Program (RDP) classifier. To standardize the OTUs data matrix across all samples produced by 16S rRNA gene sequencing, rarefaction curves were calculated with the satisfied depth (>155) using the QIIME (v1.9.1.) alpha-rarefaction curve analysis package. For the α-diversity analysis, diversity indices (ACE, Chao1, Shannon, and Simpson) were computed for each sample and compared between ESD and HTH using the QIIME plug-in. To assess similarities and differences between the carapace microbiota of ESD and HTH at the genus level, a principal coordinate analysis (PCoA) based on Bray-Curtis, weighted UniFrac, and unweighted UniFrac distances was conducted and visualized using ggplot2 in the R package v4.1.3 (Wickham, 2016). PCoA plots for microbiota obtained from the green gland, hepatopancreas, intestine, and testis were performed using Bray-Curtis distance and plotted using PAST3 (v4.03; Hammer et al., 2001). Venn diagrams showing the number of core genera shared by all lobsters within each group, the number of bacterial genera with an average relative abundance of ≥0.1%, and the number of bacterial genera present in any of the samples from each group were created using Bioinformatics & Evolutionary Genomics.¹ We plotted the 10 most abundant bacterial genera within each body site as bar graphs or stacked bar graphs using GraphPad Prism v8.4.2 (GraphPad Software, San Diego, United States). Functional profiles of carapace microbiota were predicted based on 16S rRNA genes using Phylogenetic Investigation of Communities by Reconstruction of Unobserved States (PICRUSt) v1.1.4 (Langille et al., 2013), and were analyzed in the Statistical Analysis of Metagenomic Profiles (STAMP) v2.1.3 (Parks et al., 2014) to identify and visualize functional genes (level 3 KEGG orthology) that were statistically significant between groups. Hierarchical clustering analysis of five tissues based on the genus abundance data was conducted using the hclust function with the ward.D2 method in the R package v4.2.2 (Gregory et al., 2021).

To compare carapace microbiota between groups, we conducted a non-parametric permutational multivariate analysis of variance (PERMANOVA) test on the Bray-Curtis, weighted UniFrac, and unweighted UniFrac distances using pairwiseAdonis (v0.4) in the R package (v4.1.3). We also analyzed the similarity of microbiota obtained from the green gland, hepatopancreas, intestine, and testis between ESD and HTH within each body site using PERMANOVA on the Bray-Curtis distance in PAST3 (v4.03). The relative abundance of bacterial taxa was reported as means with standard errors of the mean (SEM), and differences between ESD and HTH were assessed using the Mann–Whitney U test in the R package (v4.1.3). The abundance of KEGG orthology gene families between groups was analyzed using the two-sided, non-parametric White’s t-test in STAMP (v2.1.3). We considered differences with a p-value ≤0.05 as statistically significant.

American lobsters were obtained from three regions (Figure 1A): Eastern Long Island Sound (ELIS), Western Long Island Sound (WLIS), and 50 miles south of Montauk. Lobsters with ESD were only found in ELIS and offshore. Healthy lobsters showed no apparent signs of lesions on the carapace (Figure 1B), while lobsters with ESD mostly showed pits or erosion over the entire body (Figure 1C). There was one severe case with ulcers that removed the carapace and exposed underlying connective tissues to the external environment (Figure 1D). This may show the possibility that carapace bacteria can penetrate to the internal organs through shell lesions.

Forty samples were sequenced on the Illumina MiSeq, generating an average of 54,118 high-quality reads per sample with >Q20 (97%) and > Q30 (94%). The average number of non-chimeric reads filtered through the QIIME2 pipeline was 17,384, which were used to cluster bacterial OTUs through the SILVA 138v 16S rRNA gene reference database with a 70% classification confidence threshold. A total of 4,954 OTU classifications were generated, of which 1,022 were classified at the genus level across 40 samples. The OTU classification data also provided a framework for understanding microbial distribution and diversity by lobster tissue or ESD (Supplementary Table S2).

To determine species diversity in a sample, we examined the α-diversity indices including the ACE, Chao1, Shannon, and Simpson (Supplementary Figure S1). The ESD group showed 150.98 ± 66.09 for ACE, 141.15 ± 65.64 for Chao1, 4.91 ± 0.55 for Shannon, and 0.93 ± 0.03 for Simpson. The HTH group presented 227.20 ± 106.73 for ACE, 211.51 ± 101.94 for Chao1, 5.27 ± 0.90 for Shannon, and 0.94 ± 0.05 for Simpson. There were no significant differences (Mann–Whitney U test) in α-diversity indices between ESD and HTH. To determine the similarity and dissimilarity of carapace microbiota between the ESD and HTH groups, we performed the β-diversity analyses at the genus level using PCoA based on the Bray-Curtis (Figure 2A) and weighted UniFrac (Figure 2B) distances. The first two axes of Bray-Curtis PCoA explained 18.4 and 15.9% of variation in bacterial community structure, respectively, with significant separation of carapace microbiota between ESD and HTH (PERMANOVA, p = 0.018). The first two axes of weighted UniFrac PCoA explained 34.5 and 17.9% of variation in bacterial community structure, respectively, with significant separation of carapace microbiota between ESD and HTH (PERMANOVA, p = 0.024). Meanwhile, unweighted UniFrac PCoA showed no significant separation (PERMANOVA, p = 0.18) between ESD and HTH (Supplementary Figure S2). The data suggest that changes in relative abundances of microbial taxa are likely attributed to distinct carapace microbiota in ESD.

We analyzed and compared the carapace microbiota structure between lobsters with ESD and healthy lobsters. Lobsters with ESD had a significantly higher abundance (p = 0.04) of the phylum Bacteroidota (formerly Bacteroidetes) in their carapace microbiota, whereas healthy lobsters exhibited a higher abundance (p = 0.05) of the phylum Pseudomonadota (formerly Proteobacteria). Detailed information on the relative abundance of bacteria at each taxonomic level is provided in Supplementary Table S3. To further reveal carapace microbiota structure associated with ESD, we examined the composition of carapace microbiota at the genus level (Figure 3A). The number of bacterial genera that were found in any of the lobsters of each group was 367 in ESD and 465 in HTH, where 313 genera were shared in both groups. The number of bacterial genera that were found in any of the lobsters of each group with an average relative abundance of ≥0.1% was 64 in ESD and 81 in HTH, where 44 were shared between groups. Core genera refer to bacterial genera shared by all of the lobsters within each group. The number of core genera was 18 in ESD and 14 in HTH, where 4 were shared between groups. Core genera in ESD and HTH are listed in Supplementary Table S4. We also investigated the relative abundance of bacteria genera in the carapace of ESD and HTH. The top 10 most abundant bacterial genera in all samples included Aquimarina, Halocynthiibacter, Tenacibaculum, Loktanella, Maritimimonas, Cocleimonas, Vibrio, Cohaesibacter, Maribacter, and Perspicuibacter (Figure 3B). Aquamarina (13.5%), was the most abundant bacterial genus in ESD, and Halocynthiibacter (4.6%) was the most prevalent in HTH. There was no significant difference in the relative abundance of bacterial genera between ESD and HTH, except for Aquimarina, which was significantly higher (p < 0.01) in ESD than in HTH. To identify the function of carapace microbiota associated with ESD, we used PICRUSt, which predicts abundance of functional genes (i.e., KEGG) based on 16S rRNA gene data. The PICRUSt analysis of all carapace samples revealed the presence of 170 gene families in KEGG. Among them, 25 showed significant differences (p ≤ 0.05) between ESD and HTH (Figure 4), with 12 being more abundant in ESD and 13 more abundant in HTH. Notably, gene families involved in amino acid metabolism, including histidine, tyrosine, phenylalanine, and phosphonate and phosphinate metabolism, were enriched in ESD, whereas gene families involved in lipid metabolism, such as glycerophospholipid metabolism and primary bile acid biosynthesis, were enriched in HTH. Additionally, pathways related to genetic information processing, including base excision repair, sulfur relay system, aminoacyl-tRNA biosynthesis, and ribosome biogenesis in eukaryotes, were found to be more abundant in HTH than in ESD. Taken together, these data demonstrate that the structure and function of the carapace microbiota associated with ESD are distinct from those of healthy carapace microbiota.

To assess whether shell infection affects microbiota in the internal organs, we examined microbiota at the genus level in the green gland, hepatopancreas, intestine, and testis of both ESD and HTH. The lobsters showed different bacterial communities depending on their body site (Supplementary Figure S3). PCoA plots based on the Bray-Curtis distance showed no significant differences between ESD and HTH in the microbiota of the green gland, hepatopancreas, intestine, and testis (Supplementary Figure S4). We then investigated whether the three most abundant core carapace bacteria in ESD (Aquimarina, Halocynthiibacter, and Tenacibaculum) were present in the internal organs of lobsters. In lobsters with ESD (Figure 5A), Aquimarina and Halocynthiibacter were detected in the green gland, hepatopancreas, and testis, but not in the intestine. Tenacibaculum was present in the green gland, hepatopancreas, and intestine, except for the testis. The results of hierarchical clustering analysis, based on genus abundance data from five tissues, indicate that the carapace microbiota was more closely related to the green gland than other organs (Figure 5A). In contrast, in lobsters from HTH (Figure 5B), Aquimarina and Halocynthiibacter were absent in all internal organs tested, and Tenacibaculum was only detected in the green gland. Hierarchical clustering analysis revealed that the carapace microbiota was more similar to testis microbiota (Figure 5B). These findings suggest that carapace bacteria associated with ESD may have the ability to invade the internal tissues, particularly the green gland.