Transcriptome of Musca domestica (Diptera: Muscidae) Larvae Induced by Bacteria¹

Musca domestica L. (Diptera: Muscidae) larvae are highly adapted to niches in rotting matter, feces, and garbage and often thrive in environments with high levels of bacteria (Liew et al. 2022). Although there are many pathogenic microorganisms on the surface of fly larvae, they are rarely infected (Gold et al. 2022), which can be attributed to their unique immune mechanism (Yoon et al. 2022). The larvae can rapidly produce a variety of immune active substances, including antimicrobial peptides, lysozyme, and lectins (Yee et al. 2022). These active substances kill and inhibit viruses, bacteria, fungi, parasites, and cancer cells (Bavani et al. 2022, Olaniyi et al. 2022). Therefore, the larvae are considered a new type of natural antibiotic that can enhance immunity (Fu et al. 2018). They also are considered a natural, low-cost and sustainable protein source (Fu et al. 2022, Sing et al. 2014, Yang et al. 2019).

In China, house fly larvae are added to feed for livestock, poultry, and aquaculture (Elahi et al. 2020; Li et al. 2018, 2020) and may be used to treat diarrhea caused by pathogenic bacteria, such as Escherichia coli, Salmonella sp., and Staphylococcus aureus (An 2019, Liu et al. 2015, Zhang et al. 2019). Previous studies have shown that the larvae can improve the growth performance of pigs, broilers, and aquatic animals (Khan et al. 2016, Li et al. 2020, Wang et al. 2017a). Adding larvae to feed significantly improves breast muscle growth rate, leg muscle growth rate, and slaughter weight of broilers (Khan et al. 2016). At the same time, larvae added to feeds can improve the blood biochemical indices and immune organ indices of broilers (Song et al. 2021), enhance the immune response ability (Li et al. 2022), and reduce the intestinal population of E. coli (Chen et al. 2013). Furthermore, larvae added to feed can reduce incidences of diarrhea among piglets and treat mastitis in cows (Zou 2007) and pullorum disease in chickens (Liu et al. 2015). In these studies, antimicrobial peptides derived from the larvae are likely important components responsible for the observed responses (Wang et al. 2017b). Antimicrobial peptides are induced and activated in larvae and expressed in large quantities (Wang et al. 2016). The antibacterial peptides induced by E. coli, Salmonella, and S. aureus have broad-spectrum antibacterial effects on Gram-positive bacteria, Gram-negative bacteria, fungi, and viruses (Liu et al. 2015, Wu et al. 2013, Zhou et al. 2013).

Generally, after insects are stimulated by pathogens, some signaling molecules are activated after foreign microbes are detected by recognizing receptors. The main pathways are immunodeficiency (IMD), Janus kinase signal transduction and transcriptional activator (JAK-STAT), and Toll pathway (Zhao et al. 2021). Studies have shown that the humoral immune response of Drosophila melanogaster Meigen is mainly controlled by the Toll and IMD pathways, leading to the production of AMPs (Ertürk-Hasdemir et al. 2009, Shia et al. 2009). In the case of infection by Gram-positive bacteria, the Toll pathway is activated by the peptidoglycan recognition protein SA, SD, and the Gram-positive bacteria binding protein l (GNBPI). The interleukin receptor domain of Toll in the cytoplasm then begins to dimerize and eventually binds to the protein myeloid differentiation primary response gene 88 (Valanne et al. 2011). This complex binds to the adaptor protein Tube, thereby promoting autophosphorylation of protein kinase and phosphorylation and degradation of the nuclear transcription factor (NF-κB) inhibitory protein (cactus). The NF-κB transcription factors (DorsalL or Dif) are then translocated into the nucleus of the cell, which activates AMP transcription in the nucleus (Imler and Hoffmann 2002).

When infected by Gram-negative bacteria, the peptidoglycan recognition protein LC receptor binds to diaminopimelate peptidoglycan and the IMD signaling pathway is activated. IMD is composed of a signal complex with a death domain Fas-associated protein (dFADD) and a homologue of cysteine protease (DREDD), and DREDD is used to cleave IMD. K63 ubiquitinates the cleaved IMD, and the K63 polyubiquitin chain helps the IMD-dFADD-DREDD complex to connect with growth factor β-activating enzymes (TAK1 protein, TAB2 protein), activates the NF-κB inhibitory protein (IκB) kinase, and can degrade IκB to phosphorylate the NF-κB transcription factor (Leulier et al. 2002). After lysis and phosphorylation, the C-terminal portion of Relish remains in the cytoplasm, and the active N-terminal portion is transferred to the cell nucleus, where it activates the transcription of specific AMPs (Ertürk-Hasdemir et al. 2009). Therefore, the induced expression pathway of insect antimicrobial peptides is related to bacterial species (Gram positive or Gram negative). However, the expression pathway of antimicrobial peptides under mixed strain infection conditions needs further study.

Here, we exposed larvae to a mixture of G-positive bacteria (S. aureus and Streptococcus sp.) and Gram-negative bacteria (Escherichia coli and Salmonella sp.) and analyzed the differences of related metabolic pathways induced by mixed bacteria to larval antimicrobial peptides at the transcriptome level under different bacterial densities.

Musca domestica larvae used in these tests were obtained from a colony that was maintained at the Laboratory of Insect Ecological Resources, College of Agriculture, Shihezi University. To rear house fly larvae, we added 6 g of milk powder in the basic medium of 30 g fermented wheat bran + 3 g brown sugar (Jiang and He 2021). Each medium feeds 60 house fly larvae. The colony was maintained on an artificial diet at 25°C and a relative humidity of 60–70%.

Samples of E. coli, Salmonella sp., S. aureus, and Streptococcus sp. were removed from an ultralow temperature freezer maintained at −80°C. After reconstitution, each colony was inoculated into brain-heart infusion broth and cultured in a constant temperature shaker at 37°C for 12 h. Bacteria in each suspension were enumerated and diluted to 1 × 10⁸ colony-forming units (CFUs)/ml for a high dose and 1 × 10⁷ CFUs/ml for a low dose.

Larvae (63 days old) obtained from the laboratory colony were divided into three groups. Two groups of larvae were injected by micro syringe in the abdomens with 0.5 μl of the mixed solution of the four bacteria. One of these groups was injected with the low concentration dose (L) of 1 × 10⁷ CFUs/ml, while the second group was injected with the high concentration dose (H) of 1 × 10⁸ CFUs/ml. The third group served as a control (C) and was injected with the same volume of normal saline. The three groups were held at 25°C for 6 h. Larvae were then collected, immediately placed in liquid nitrogen, and transferred to the −80°C ultralow temperature freezer.

Sample RNA extraction, quality detection, transcriptome sequencing, and statistical analyses were commissioned by Beijing Nuohe Zhiyuan Technology Co., Ltd. Briefly, total RNA was extracted from the cotton samples using TRIzol (Tiangen Biotech [Beijing] Co., Ltd) according to manufacturer’s instructions. RNA purity was checked using the NanoPhotometer® spectrophotometer (Implen USA Inc., Los Angeles, CA). RNA concentration was measured using the Qubit® RNA Assay Kit in the Qubit® 2.0 Flurometer (Life Technologies, Rockville, MD), and RNA integrity was assessed with the RNA Nano 6000 Assay Kit of the Bioanalyzer 2100 system (Agilent Technologies, Santa Clara, CA). All samples had an RNA Integrity Number (RIN) above 6.7. RNA sequencing libraries were prepared with NEBNext® Ultra RNA Library Prep Kit for Illumina® (New England BioLabs, Ipswich, MA) and sequenced on an Illumina Hiseq^TM 2500 platform at an average depth of ∼66 million reads per sample. Raw sequencing reads were quality assessed with FastQ. To pass the initial quality control check, the average Phred score of each base position across all reads had to be at least 30. Reads were further processed by cutting individual low-quality bases and removing adapter and other Illumina-specific sequences with ng-qc using default parameters. HISAT2-2.0.4 was then used to map the trimmed reads to the cotton AD1_ZJU_v2.1 reference genome (Kim et al. 2015, Mortazavi et al. 2008). To quantify gene expression levels, mapped reads were summarized at the gene level using HTSeq v0.6.0 (Anders and Huber 2012). Differential expression analyses were performed with the DESeq2 R package (1.10.1), and Gene Ontology (GO) enrichment analyses were conducted using the clusterProfiler R package (Anders and Huber 2012, Wang et al. 2010). The significance threshold used was an adjusted P value of 0.05 and an absolute fold change of 2 for the differential expression analysis and adjusted P value of <0.05 for GO enrichment analysis (Robinson et al. 2010, Young et al. 2010). We used the clusterProfiler R package to test the statistical enrichment of differential expression genes in the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways (Kanehisa et al. 2008).

Information about the data, for example, the quality value of the original data, was determined. Trimmomatic was used for data processing to obtain clean data of 19.18 Gb. The proportion of Q20 in each sample was >98%, and the proportion of Q30 was >93%. The data could, thus, be used for subsequent analysis. Using Trinity to assemble clean data de novo into transcripts, we obtained 58,701 genes.

A total of 34,666 (59.06%) transcripts was annotated successfully by comparing 58,701 genes obtained by sequencing with each database (Table 1). The number of transcripts annotated in the nucleotide database was the largest (53.4%) followed by that annotated in the nonredundant database (44.51%). The number of transcripts annotated in the KEGG database was the least (5.32%).

We determined that 20,433 genes were similar to those in the GO database (Fig. 1). Genes in the larval transcriptome were divided into 58 branches of molecular function, biological process, and cellular component according to GO function. Among the molecular functions, the binding function involved 11,922 genes. Among the biological process categories, the cell process and metabolic process involved 12,781 and 10,159 genes, respectively. Among the cell component processes, the cell and cell composition involved 13,016 and 12,985 genes, respectively. The expression abundance of other genes varied.

The KEGG metabolic pathways of larval transcripts included a total of 3,120 (5.32%) transcripts (Fig. 2). When annotated, these transcripts were distributed in 33 pathways of cell process, environmental information processing, genetic information processing, metabolic function, and biological system. We found that the metabolic pathways of signal transduction, translation, transport and catabolism, folding classification and degradation, and endocrine system exhibited the largest number of genes.

In the KOG database, 11,420 genes were annotated (Fig. 3), with signal transduction mechanisms (13.68%) ranking first in the 25 categories; followed by general function prediction (12.84%); posttranslational modification, protein turnover, and molecular chaperones (9.52%); transcriptional class (7.67%); translation, ribosomal structure, and biogenesis (7.63%); and amino acid transport and metabolism (6.44%).

The expression of transcripts was calculated, and the differential expression of maggot genes using different induction methods was compared (Fig. 4). Among them, 2,745 genes were upregulated and 779 genes were downregulated in L versus C, whereas H versus C had 1,485 genes upregulated and 3,644 genes downregulated and H versus L had 1,229 genes upregulated and 6,233 genes downregulated.

The biological process of the low concentration dose of the bacterial mixture induced larval responses to the bacteria, defense response to the bacteria, and defense response to external biological stimulation (Fig. 5). Extracellular region, extracellular space, and extracellular space plus other cellular components also were affected, as was the expression of genes related to molecular function, such as serine endopeptidase activity and serine hydrolase activity.

The high concentration dose of the bacterial mixture induced chitin metabolism, amino polysaccharide metabolism, compound metabolism containing glucosamine, and other biological processes (Fig. 5). Extracellular region, extracellular space, extracellular region plus other cellular components, and the molecular functional expressions such as serine endopeptidase activity, chitin binding, and serine hydrolase activity were affected.

The biological processes of chitin metabolism process, glucosamine-containing compound metabolism process, and amino sugar metabolism process induced by the high concentration dose of the bacterial solution were higher than those induced by the low concentration dose of the bacterial solution. They included cellular components such as cytoplasmic, extracellular region, lipid droplets, chitin binding, serine endopeptidase activity, serine hydrolase activity, serine peptidase activity, and protein binding (Fig. 5).

From the enrichment process of different genes in the larvae, it can be found that the most significant pathways of differential gene enrichment following injection with the low concentration dose of the bacterial mixture were pancreatic secretion, neurotrophin signaling pathway, cAMP signaling pathway, and the Toll and IMD signaling pathway. The most obvious pathways of enrichment following injection of the high concentration were ribosome, fatty acid metabolism, glutathione metabolism, and AMP-activated protein kinase (AMPK) signaling pathway. The enrichment of the insulin signaling pathway, HIF-1 signaling pathway, and chemokine signaling pathway in larvae injected with the high concentration dose was more apparent than that induced by injection with the low concentration dose (Fig. 6).

The immune defense mechanism of insects is unique in that they rely mainly on antibacterial substances in the hemolymph to respond to pathogenic microorganisms (Röhrich et al. 2012, Vilcinskas 2013). Previous studies have focused on the differences in the activity of different induced strains against insect antimicrobial peptides (Wang et al. 2016, Zhai et al. 2006) and revealed the differences in the antimicrobial-induced expression of strains from molecular signaling pathways (Ertürk-Hasdemir et al. 2009, Shia et al. 2009). We focused on the difference in the number of strains of mixed bacteria on the induced expression of antimicrobial peptides in larvae of M. domestica. The GO enrichment results showed that the low concentration dose of the bacterial mixture induced the defense response process of the larvae to the bacteria, whereas the high concentration dose of the mixture induced the basic metabolic processes of the larvae. As expected, both concentrations increased the molecular functional expression of serine hydrolase activity that plays an important role in reducing inflammatory response (Nomura et al. 2011) and inhibiting bacterial infection (Surpeta et al. 2021). This finding indicates that serine hydrolase is highly expressed to inhibit bacterial infection and reduce the adverse effects of bacteria on M. domestica larvae.

The expression level of insect antimicrobial peptides is related not only to the type of bacteria inducing the response (Gram positive or Gram negative) (Tu and Wu 2012) but also to the concentration of induced bacteria to which the larvae are exposed (Chen 2015). For example, an increased concentration of E. coli yields a greater concentration of antimicrobial peptide and a greater inhibition zone about the peptides in drop-plate assays (Wu et al. 2013). This finding may be related to the stimulus intensity of pathogens on the immune system.

In this study, the enrichment pathways following injection of the low concentration dose of the bacterial mixture induced M. domestica larval expression of mainly the pancreatic secretion, neurotrophic factor signaling pathway, cAMP signaling pathway, Toll, and IMD signaling pathway. Injection with the high concentration dose induced larval enrichment pathways of mainly the ribosome, fatty acid metabolism, glutathione metabolism, and AMPK signaling pathway. These results indicate that the number of microorganisms in the bacterial solution has a different effect on the pathway expression of infected larvae. In addition, the enrichment of insulin signaling pathway, HIF-1 signaling pathway, and the chemokine signaling pathway in the larvae induced by injection with the high concentration dose of the bacterial mixture was more apparent than that induced by injection with the low concentration dose. This result indicates that the degree of response of the same metabolic pathway is also related to the number of microorganisms to which the larvae are exposed.

When infected with Gram-positive or Gram-negative bacteria, Drosophila show an increased production of antimicrobial peptides through the Toll/IMD signaling pathway (Ertürk-Hasdemir et al. 2009, Shia et al. 2009). This finding indicates that the expression pathways of insect antimicrobial peptides differ with bacterial species or strain. We found that in the case of infection with mixtures of Gram-positive and Gram-negative bacteria, M. domestica larva not only increased the expression of Toll and IMD signaling pathways but also increased the expression of cAMP, AMPK, HIF-1, and other signaling pathways.

AMPK is an enzyme complex that is an important cell metabolic sensor that regulates energy metabolism and also plays a regulatory role in the inflammatory response (Wang et al. 2022). IMD increases intracellular cAMP production by binding to the calcitonin receptor-like site and the modified protein complex activation coreceptor and then activates the AMPK signaling pathway to reduce inflammation (Pang et al. 2016). The HIF-1 signaling pathway is a stress signaling pathway in organisms under hypoxia that increases ATP synthesis and ensures normal metabolism by regulating the function of HIF-1α (Cadenas 2018). This finding indicated that injection with the low concentration dose of the bacterial solution caused larval expression of the Toll and IMD signaling pathways, resulting in the production of antimicrobial peptides in M. domestica larvae, whereas the high dose of bacterial solution produced hypoxia stress on larvae and promoted the basal metabolic function of cells.

Overall, our study revealed the stress changes in M. domestica larvae under pathogen stimulation. From the transcriptome perspective, we found some key pathways, including HIF-1 signaling pathway, AMPK signaling pathway, and Chemokine signaling pathway, which are related to the synthesis and function of antimicrobial peptides. We speculate that these pathways are involved in regulating the antibacterial response of M. domestica larvae to pathogens, and the related mechanisms need further verification and screening to further clarify the molecular mechanism of stress response of larvae caused by pathogens.

This work was supported by the National Natural Science Foundation of China (Grant No. 31660703).