Symbiont Genomic Features and Localization in the Bean Beetle Callosobruchus maculatus

The bacterial community of the bean beetle C. maculatus was characterized throughout development using Illumina amplicon sequencing of the bacterial 16S rRNA gene. A total of 480,536 high-quality sequences were generated across 20 samples representing eggs, larvae, and adults. These sequences were subsequently binned to 425 amplicon sequence variants (ASVs) (100% sequence similarity). After taxonomic assignment, 42 ASVs were removed because they were classified as chloroplasts (23) or Archaea (5). The 5 ASVs assigned to Archaea belonged to the YLA114 phylum and were detected in just one egg sample, comprising 5% of all reads within the sample. The rarity of these microbes within bean beetles, in addition to their presence in nonscarabaeoid beetles being unusual, suggests that Archaea may be a contaminant in this sample.

The remaining 383 ASVs contained an average of 24,129 reads per sample, representing 95.4% of the original sequences (see Table S1 in the supplemental material). These ASVs were assigned to 91 bacterial species spanning 81 genera (97% and 95% ASV similarities, respectively). Eggs more commonly associate with Proteobacteria, followed by Actinobacteria, Firmicutes, and Bacteroidetes, whereas larvae and adults largely harbor Firmicutes and, to a lesser extent, Proteobacteria (Table S2). The bacterial communities of larvae and adults, compositionally simple, are predominantly represented by Staphylococcus, Enterococcus, and Enterobacter (Fig. 1). The high relative abundance of Staphylococcus in larvae and adults is especially notable, ranging between 99.5% and 97.6% of the total sequences per sample, respectively (Fig. 1). This is consistent with culture-dependent approaches pointing to Staphylococcus as the most common bacterial genus consistently isolated from C. maculatus (21, 24).

Several distance metrics were calculated to assess the variation in bacterial community structure between sexes, body sites (whole body versus digestive tract), and developmental stages. Beta diversity analyses based on weighted UniFrac distance matrices (Table S3) revealed that the microbiota composition of the beetle host does not differ between sexes (pseudo-F = 1.01 and P = 0.44 by permutational multivariate analysis of variance [PERMANOVA]). Likewise, no variation in community composition was observed between body sites (pseudo-F = 0.21 and P = 0.45 by PERMANOVA), confirming that the gut bacterial community is indeed representative of the adult beetle as a whole, as has been described for other insect groups (25). However, there is significant variation between developmental stages. The most significant variation that we observe compositionally is that for eggs relative to larvae and adults (pseudo-F = 13.72 and q = 0.03, and pseudo-F = 45.86 and q = 0.006, respectively, by pairwise PERMANOVA). The latter two do not differ from each other (pseudo-F = 1.35 and q = 0.29 by pairwise PERMANOVA). Overall, these results are consistent with those obtained from the other distance matrices tested, Bray-Curtis and Jaccard (Table S3).

An analysis of the “core microbiome” (i.e., those ASVs shared by 100% of all samples) revealed that just one ASV is maintained throughout the developmental cycle of C. maculatus (Fig. 2). Phylogenetic placement (see below) classified this ASV as Staphylococcus gallinarum. While the transmission route was not directly elucidated, the presence of this ASV in every sample across all developmental stages suggests that S. gallinarum is vertically transmitted from mother to offspring in the bean beetle.

The differences in bacterial diversity observed throughout development are most severe between the egg and larval stages of C. maculatus. While this may stem from uneven sampling depth, it could also reflect the development of the gut as a habitat to house a specialized community adapted to the unique physicochemical conditions forming along the digestive tract (26). These conditions can vary substantially in both pH and oxygen availability (25, 27), potentially sieving environmental microbes conditioned to the aerobic surfaces on the egg. The compositional stability observed posteclosion in C. maculatus may also reflect a specific and consistent immune profile in both larvae and adults (26), which can confer tolerance to commensal or beneficial members (28, 29) while eliminating environmental or pathogenic strains (30).

As 16S rRNA gene profiling revealed that the bacterial community inhabiting the digestive tract of C. maculatus is representative of the beetle as a whole, we aimed to localize the symbiont along the insect’s gut using fluorescence in situ hybridization (FISH). Performed on representative semithin sections of the digestive system, FISH revealed that bacterial colonization in the gut is largely restricted to the foregut, evidenced by the clear, dual fluorescence of eubacterial and S. gallinarum-specific probes in this region but not in the midgut or hindgut (Fig. 3). Within the foregut, symbiont localization is exclusively extracellular along the epithelial lining, with limited bacterial presence in the lumen. The overlapping fluorescence of the S. gallinarum-specific probe (Cy5) and the general eubacterial probe (Cy3) is consistent with 16S rRNA gene profiling of the digestive tract revealing that the gut bacterial community is largely dominated by the firmicute.

Whole-genome sequencing of S. gallinarum isolated from the beetles’ digestive tract was carried out to understand the biological role of the symbiont within C. maculatus’ gut. A hybrid assembly of long MinION and short Illumina reads yielded a contiguous and circularized genome with general features that are consistent with those of other members of the Staphylococcus genus (31). The bacterial symbiont harbors one chromosome (2.89 Mb; 33% GC content) (Fig. 4A and B) and one plasmid (46 kb; 28% GC content). The direction of gene transcription is asymmetrical with respect to the vertical axis of the genome map (Fig. 4A). The change in direction is observed around half past 6, which is supported by GC skew values. This observation suggests that the replication termination site of C. maculatus-associated S. gallinarum is located around the 195° position on the genome map and not directly opposite the replication origin, resulting in two replichores of different lengths. This lack of symmetry is common in bacteria and often correlates with the enrichment of the leading replicating strand in essential genes, which are also more conserved than genes in the lagging strand (32). In Staphylococcus sp., the leading strand can additionally contain horizontally acquired genes coding for traits specific to each strain within a species (33).

Staphylococcus strains, historically characterized as human pathogens, are classified by their ability to coagulate blood into coagulase-positive and -negative strains. Whereas Staphylococcus aureus normally presents a symmetrical GC skew across the vertical axis of the genome map, other strains within the genus, such as S. saprophyticus, S. lugdunensis, and other coagulase-negative strains, do not (34). S. gallinarum is known to be a coagulase-negative species (35), and while a coagulase test was not performed on C. maculatus-associated S. gallinarum, this bacterium presents the expected GC skew of a coagulase-negative strain. Unlike their coagulase-positive counterparts, negative strains usually lack aggressive virulence properties and maintain symbiotic interactions with their hosts, colonizing surfaces and mucous membranes of animals (36).

To resolve the phylogenetic placement of C. maculatus-associated S. gallinarum relative to representative Firmicutes, a phylogenomic tree was constructed on the basis of 49 marker genes (Table S4). This confirmed the placement of the beetle’s symbiont within the Staphylococcus genus and demonstrated that it is most closely related to Staphylococcus gallinarum, forming a well-supported monophyletic clade (Fig. 4C). Staphylococcus gallinarum was initially isolated and described from a poultry sample (35) and has since been isolated from the skin and respiratory tracts of a number of mammal (37) and insect (23, 38, 39) herbivores. Interestingly, it is yet to be described in the absence of a host.

The genome of C. maculatus-associated S. gallinarum harbors 2,753 putative protein-coding genes, with an average gene size of 892 bp, resulting in a coding percentage of 84.9%. Functional annotation assigned specific functions to 70% of the genes according to the Clusters of Orthologous Groups (COG) database, whereas the rest were assigned to hypothetical proteins or lacked a known function. The genome is equipped with 7 copies of the ribosomal operon containing the 3 structural rRNA genes (5S, 16S, and 23S) and 60 tRNA genes assigning all 20 amino acids. Most genes in the chromosome are involved in metabolic functions (38%), followed by genes involved in information storage and processing (18%) and cellular processes and signaling (11%). Within the metabolic functions, “carbohydrate transport and metabolism” is the most abundant gene category.

Alongside other leaf beetles (Chrysomelidae), C. maculatus is a member of the speciose Phytophaga clade of herbivorous beetles. Recent genomic insights into beetle adaptation and diversification revealed that the evolution of specialized herbivory across the coleopteran order, and the origin of the Phytophaga, coincided with the genomic integration of microbial genes encoding plant cell wall-degrading enzymes (PCWDEs) (40, 41). These enzymes, typically glucoside hydrolases (GHs), upgrade the digestive physiology of herbivorous beetles to process a diet rich in recalcitrant polymers such as cellulose, hemicellulose, and pectin (40, 41).

In encoding a subset of pectinases (GH 28) and mannanases (GH 5_10) endogenously, C. maculatus is unique, relative to other members of the Phytophaga, in lacking all of the horizontally acquired cellulases belonging to the GH families 9, 45, and 48 (40–42). While C. maculatus-encoded mannanases can partially degrade cellulose alongside galactomannan (glucose is an isomer of mannose) (43), most herbivorous beetles deploy multiple cellulolytic enzymes to synergistically monomerize the polysaccharide. As symbiont acquisition can offset the loss of endogenous PCWDEs (44–46), we explored whether S. gallinarum might play a similar role for the bean beetle by encoding cellulases to complement the insect’s repertoire of endogenous enzymes. In cataloging symbiont-encoded PCWDEs, our annotation revealed 35 enzymatic families (Table S5) putatively involved in the breakdown of plant cell wall polymers. Notably, this included a pair of endoglucanases (GH 6) predicted to internally cleave the β-1,4-glycosidic bonds in cellulose.

Genome annotation revealed the presence of metabolic pathways for the biosynthesis of five B vitamins (Fig. 5A): thiamine (B₁), riboflavin (B₂), pantothenic acid (B₅), biotin (B₇), and folic acid (B₉). The pathway leading to the biosynthesis of nicotinic acid (B₃) from aspartate as a precursor was incomplete, lacking the first three enzymes. Finally, the biosynthetic pathway leading to pyridoxine (B₆) was completely absent. Many seed-based diets are deficient in B vitamin content, necessitating supplementation from specialized endosymbionts (47) or gut bacteria (16). Symbiont-mediated B vitamin supplementation has been demonstrated in grain-feeding beetles, including Lasioderma serricorne, Stegobium paniceum (48), and S. oryzae (49). Similarly, firebugs rely on dietary supplements from their gut-associated Coriobacteriaceae symbionts to offset nutritional limitations arising from specializing on cottonseeds deficient in thiamine, riboflavin, pantothenic acid, and pyridoxine (16).

Likewise, metabolic pathways involved in the production of essential and aromatic amino acids are present in the genome of C. maculatus-associated S. gallinarum. In addition to complete pathways for all essential amino acids (except arginine and lysine), S. gallinarum encodes a complete shikimate pathway for the biosynthesis of chorismate and prephenate, both precursors of tyrosine (Fig. 5B). Tyrosine is required for the biosynthesis of melanin and catecholamines, which are essential in the tanning and sclerotization of the insect cuticle (50). Given that the shikimate pathway is frequently absent in animals, beetles feeding on deficient diets often rely on bacteria providing tyrosine precursors to build the hardened elytra that are characteristic of these insects (18, 51, 52). Preventing desiccation in dry environments through a strongly sclerotized cuticle, these types of associations may be particularly critical for grain beetles such as S. oryzae, O. surinamensis (19, 20), and, potentially, C. maculatus.

In summary, we show that the simple bacterial community previously described for C. maculatus is stable across development and that the microbial assembly along the insect’s digestive tract is representative of the insect’s microbiome as a whole. In characterizing the metabolic features of its sole core member, S. gallinarum, we highlight several putatively symbiotic functions toward upgrading the physiology and nutritional ecology of its herbivorous host. Granivorous insects rely on a multitude of metabolic functions from their obligate symbionts, ranging from nutritional supplementation (16, 18–20) to the degradation of recalcitrant or toxic plant compounds (17). The annotation of a complete biosynthetic pathway for the production of tyrosine by S. gallinarum is notable given recent findings implicating the symbiont-produced amino acid as an important feature defining numerous beetle-bacterial symbioses, including grain-feeding species belonging to the Silvanidae and Curculionidae families (18, 20, 51–53). As a precursor for the biosynthesis of both melanin and catecholamines, tyrosine is central to the tanning and hardening of beetle cuticles (18, 20, 51–53).

Pertaining to the digestive ecology of C. maculatus, symbiont acquisition may compensate for the loss of several endogenous cellulases from the insect’s genome (42, 43, 54). Relative to other members of the Phytophaga, the bean beetle no longer encodes cellulolytic enzymes belonging to GH families 9, 45, and 48 (41, 42). Alongside pectinases and xylanases, the acquisition of these digestive enzymes is cited as a key innovation that spurred the evolution of herbivory in beetles (41). Our annotation of complementary cellulases belonging to the GH 6 family in the genome of S. gallinarum is consistent with the cellulolytic metabolism of the firmicute. By potentially underlying a digestive role for its host, symbiont-produced endoglucanases would complement a range of mannanases and pectinases endogenously retained by C. maculatus. Such a division of labor is demonstrated in other Chrysomelidae subfamilies where the loss of insect pectinases is offset following the independent acquisition of foregut-associated pectinolytic symbionts (44–46).

Toward determining the exact mutualistic role of S. gallinarum, future efforts will establish a targeted and efficient symbiont-clearing method, coupled with the development of a chemically defined artificial diet. Elucidating the metabolic importance of the symbiont will further our understanding of an economically important insect species and possibly shed light on novel avenues for the control of the bean beetle.

A laboratory culture of C. maculatus was established using beetles collected in Columbus, OH. The insects were reared in plastic containers (20 by 13 by 12 cm) at room temperature and supplied with cowpea seeds (V. unguiculata) ad libitum. The grains and rearing cages were autoclaved to minimize exposure to environmental microbes.

To characterize the bacterial community associated with C. maculatus during development and to investigate the heritability of the insect’s bacterial associates, a cohort of beetles from four different females was established. One replicate per line was collected for every sample type (i.e., larvae, adult females and adult males, and eggs), resulting in a total of four replicates per developmental stage. Eggs were obtained by allowing F₁ females to oviposit on autoclaved glass beads for 24 h. All eggs laid within 48 h by each female were pooled and constituted a single replicate to ensure sufficient DNA extraction yields. For larvae and whole adults, entire individuals were used to extract genomic DNA. To assess whether the gut bacterial symbionts of C. maculatus are representative of the community associated with the beetle as a whole, we additionally sampled isolated guts. Prior to dissection, animals were submerged in 0.1% sodium dodecyl sulfate (SDS) and rinsed with sterile water to remove surface contaminants. A continuous gut (foregut, midgut, and hindgut) was dissected under sterile water and placed into individual Eppendorf tubes for DNA extraction.

Insects were euthanized by freezing at −20°C for 15 min. Prior to extraction, egg, larva, and adult samples were surface sterilized with ethanol, subsequently flash-frozen with liquid nitrogen, homogenized with a pestle, and resuspended in 0.5 ml of sterile water. DNA was extracted using a phenol-chloroform protocol according to the instructions at www.beanbeetle.org (accessed January 2021). Subsequently, the extracted DNA was used for sequencing with the Illumina MiSeq platform and paired-end 250-bp technology. Sequencing was performed at Mr. DNA (Shallowater, TX, USA). The V4 variable region of the bacterial 16S rRNA gene was sequenced using the Earth Microbiome Project primer pair 515F and 806R (55). Barcodes were attached to the forward primer. Amplification was performed using the HotStarTaq Plus master mix kit (Qiagen, USA) under the following PCR conditions: 3 min at 94°C, followed by 28 cycles of 30 s at 94°C, 40 s at 53°C, and 1 min at 72°C, with a final 5-min elongation step at 72°C. Amplicons were run in a 2% agarose gel to determine amplification success. Samples were purified using calibrated AMPure XP beads, and purified amplicons were used to construct Illumina DNA libraries.

The quality of raw sequences was controlled by performing q25 trimming of the ends of reads with Trimmomatic (56). Forward and reverse sequences were joined in Qiime v1 (57) (script join_paired_ends.py) and reoriented where necessary. Sequences shorter than 200 bp were removed. Processing of high-quality reads was performed using Qiime v2 2017.4 (58). Joined reads were demultiplexed using the q2-demux plug-in, and forward and reverse primers were removed with q2-cutadapt. Reads were further subjected to a quality-filtering step using DEBLUR (59) to denoise data by setting –p-trim-length at 273 bp. Amplicon sequence variants (ASVs) were aligned using mafft (60). In order to build a phylogenetic tree, the alignment was filtered to eliminate highly variable positions that added noise to the tree, which was constructed using fastTree (61) via q2-phylogeny, and rooted at the midpoint. In order to assess the quality of our sequencing depth, our samples were rarefied to 3,000 sequences per sample. After rarefaction, one egg sample was observed to have too few reads and was eliminated from the data set. Taxonomy was assigned to the ASV table using q2-feature-classifier (62) classify-sklearn against the Greengenes database released in May 2013 (99% sequence similarity) (63). ASVs corresponding to Archaea, mitochondria, and chloroplasts were eliminated. In order to facilitate visual representation of the bacterial community, ASVs were grouped according to their taxonomy assignment at 95% sequence similarity (genus level).

Three beta diversity metrics (weighted UniFrac [64], Bray-Curtis dissimilarity, and the Jaccard similarity coefficient), as well as principal coordinates, were calculated using q2-diversity after rarifying samples (subsampled without replacement) to 3,000 sequences per sample. We statistically tested whether ASV association with C. maculatus was dependent on life stage, sex, or body site (whole body versus gut) by independently analyzing the matrices with PERMANOVA. Where necessary, the resulting P values were adjusted to q values for multiple comparisons using Benjamini-Hochberg false discovery rate correction (65). P values of <0.05 and q values of <0.05 were considered significant.

The ASV table was employed to calculate the “core” microbiome, defined as all the ASVs present in 100% of all samples in each life stage. Taxa shared by all developmental stages were visualized through a Venn diagram generated with the free online tool Venny (66).

To localize the bacterial symbiotic community within C. maculatus, fluorescence in situ hybridization (FISH) was performed on semithin section preparations as previously described by Weiss and Kaltenpoth (67). Whole beetles, fixed in 4% paraformaldehyde (PFA) in phosphate-buffered saline (PBS), were embedded in Technovit 8100 (Heraeus Kulzer, Wehrheim, Germany), and semithin sections (8 μm) were generated on a rotary microtome (Leica RM 2245) using a glass blade. Hybridization was achieved by incubation for 16 h at 50°C with hybridization buffer (0.9 M NaCl, 0.02 M Tris-HCl [pH 8.0], 0.01% SDS) containing 50 μl/ml of the probe Staph_Callosobruch-Cy5 (GGACTACGGGGGTATCTAATCC) targeting the Callosobruchus-associated Staphylococcus symbiont (500 nM), the general eubacterial probe EUB-3888-Cy3 (GCTGCCTCCCGTAGGAGT) (500 nM), as well as 5 μg/ml 4′,6-diamidino-2-phenylindole (DAPI) for counterstaining of host cell nuclei. Slides were incubated for 20 min in 50°C preheated washing buffer (0.1 M NaCl, 0.02 M Tris-HCl [pH 8.0], 0.01% SDS, 5 mM EDTA), followed by washing with distilled water (dH₂O) for another 20 min. Upon drying at room temperature, Vectashield was applied prior to observation using a Zeiss (Jena, Germany) AxioImager.Z2 fluorescence microscope with Apotome.2.

To better understand the possible symbiotic role of C. maculatus-associated S. gallinarum, the microbe was isolated on culture plates. Insects were euthanized by freezing at −80°C for 5 min. Samples were surface sterilized by submerging single individuals in 10% bleach for 3 s and sterile water for 10 s, followed by submersion in 70% ethanol for 5 s and a final rinse with sterile water for 10 s. Individual beetles were placed in Eppendorf tubes containing 450 μl of a sterile saline solution (Sigma). Samples were homogenized using a sterile pestle, and debris was pelleted by centrifugation (20 s at 8,000 rpm). One hundred microliters of this mixture was plated on each of the following media and incubated at 30°C: nutrient agar, eosin methylene blue, blood agar, and phenylethyl alcohol agar. Single colonies were further streaked onto LB agar plates to obtain pure cultures that were subsequently stored in glycerol stocks until further use.

Pure isolates were identified via colony PCR using the eubacterial primer pair 27F-1429R (68). The temperature profile was as follows: an initial step of 95°C for 10 min was followed by 36 cycles consisting of 30 s at 95°C, 30 s at 55°C, and an elongation step of 1.5 min at 72°C, followed by 4 min at 72°C. PCR products were visualized in a 1.5% agarose gel. Successful amplicons were sent to Eurofins for sequencing. Taxonomy was assigned to each isolate by comparing each sequence to the NCBI database using BLAST (69).

The genome of a C. maculatus-associated S. gallinarum isolate from C. maculatus’ digestive tract was sequenced in order to explore its general features and metabolic capabilities. Genomic DNA was extracted from the isolate using the DNeasy blood and tissue DNA extraction kit (Qiagen), according to the manufacturer’s protocol for Gram-positive bacteria. An unsheared DNA library was prepared with an SQK-LSK109 ligation sequencing kit (Oxford Nanopore Technologies). DNA was size selected (enrichment for reads of >3 kb) during the ligation step by washing the sample with long-fragment buffer (Oxford Nanopore Technologies). DNA was sequenced in-house with the Oxford Nanopore MinION sequencer, and base calling was performed with MinKNOW 19.06.8 (Oxford Nanopore Technologies) with the “fast” algorithm. The flow cell used was FLO-FLG001 (Oxford Nanopore Technologies). Genome sequencing resulted in 88,949 reads, an N₅₀ of 13,447, and an average read length of 6,857.4, as determined by the “assembly-stats 1.0.1” script, which was run with default parameters (https://github.com/sanger-pathogens/assembly-stats [accessed March 2021]). No quality control (QC), error correction, or adapter trimming was performed. A replicate sample was also sequenced with the Illumina HiSeq platform utilizing 150-bp paired-end technology at a depth of ∼80 million reads (Omega Bioservices, GA, USA). An Illumina DNA library was prepared using the Kapa HyperPlus kit (Roche) according to the manufacturer’s protocol. A hybrid assembly was carried out by combining long Nanopore with short Illumina reads using Unicycler v0.4.8 (70), which automatically circularized the genome and rotated it to DnaA. The genome was confirmed to be complete and visualized by employing Bandage (v0.8.1) with default parameters. Annotation of contigs was performed using PROKKA 1.11 (71). Gene annotations were clustered into Clusters of Orthologous Groups (COG) categories and KO (KEGG Orthology) pathways using eggNOG (72, 73). Carbohydrate-active enzymes were predicted using DIAMOND within the dbCAN2 meta server (74). GC skew values were calculated using the script GC_content.pl by Damien Richard (https://github.com/DamienFr/ [accessed January 2021]), using default parameters. Genomic data were visualized with Circos (75).

Phylogenetic inference was performed using a set of 49 universal and highly conserved genes defined by COG and shared by members of the Firmicutes phylum (see Table S3 in the supplemental material). Extracted sequences were aligned using MUSCLE (76), ahead of concatenation. Maximum likelihood (ML) tree inference, under the G+T+I model, was performed using GARLI 2.0 as implemented in the Cipres Science Gateway (77), with 1,000 bootstrap replicates.

Genomic sequencing data generated for this article are available in the NCBI and the National Center for Biotechnology Information databases under BioProject accession number PRJNA694570. Amplicon sequencing data and sample metadata are available at https://figshare.com/articles/dataset/16s_rRNA_gene_bacterial_community_profiling_files/14330867. Genome annotation data (gff) are available at https://figshare.com/articles/dataset/gff_file_for_S_gallinarum_s_genome/14340344.