1 INTRODUCTION

House dust mites (HDM) are considered the most important source of indoor allergens worldwide, leading to perennial rhinitis and asthma, and worsening atopic dermatitis.1 One of the major producers of house dust allergens and the most widespread species is Dermatophagoides pteronyssinus. The main route of exposure to HDM allergens is via their fecal particles, which become easily airborne.2 Allergy to HDM is currently treated by allergen-specific immunotherapy (AIT) using natural allergen extracts produced from mite-derived fractions obtained in large-scale commercial cultures. These extracts are very complex since they include biological components from both the mites and their associated microbiomes.3, 4 Noteworthy, mite cultures commonly contain immunoactive endotoxins of bacterial origin such as liposaccharide (LPS) that are systematically controlled as per pharmacological certifications.5, 6

During the last decade, viral metagenomics hand in hand with breakthrough mining techniques on high throughput sequencing data are revealing the extraordinary diversity and ubiquity of viruses in the biosphere.7 A paradigmatic example is the unprecedented number of novel virus species being discovered in arthropods, showing genomic structures and phylogenetic relationships exceeding by far our previous knowledge on invertebrate RNA viruses.8 In the case of mites, arachnids of the Acari sub-class also including ticks, viral metagenomic research is still in an incipient stage. Some notable exceptions are studies carried out on species of economically important impact, such as the honeybee ectoparasite Varroa destructor,9, 10 the phytophagous pest Tetranychus urticae,11 and, very recently, the medically important mites D. pteronyssinus, Dermatophagoides farinae, and Tyrophagus putrescentiae.12

Expanding our understanding on the viruses associated to HDM is of medical interest. Commercial HDM cultures may be infected with non-human infecting viruses that make their way to the finished pharmaceutical products used to diagnose and treat allergy. This could have unknown consequences for both the mite culture and the exposed human subject, which should be considered in terms of safety and/or efficacy. In the case of mites, viruses could affect the culture performance or the production of allergens; in turn, in humans, they might trigger immunological or inflammatory responses. Although not expected to show pathogenic effects on humans, some components of the mite-derived virions could show antigenic properties or elicit innate immune responses that further interact with the allergic response, as has been demonstrated for pathogenic viruses.13, 14 Interestingly, the inclusion of inactive viral-like particles as an adjuvant in AIT products has been suggested due to their immunomodulatory effects.15, 16 Additionally, the study of mite-pathogenic viruses could eventually be promising to develop new biocontrol methods for HDM.

Herein, we report the discovery and characterization of seven RNA viruses in D. pteronyssinus. We detect the presence of viral RNA in different in-house laboratory colonies by direct molecular methods. We localize picornavirus-like particles in the digestive system and fecal pellets of the mite by electron microscopy. Additionally, we demonstrate that fecal-borne viruses can be a way of horizontal transmission between mites, and we gain insight on the effects of viral infection over mite population growth. Finally, by in silico screening and direct molecular methods, we confirm the presence of homolog sequences to these viruses in other sources beyond our laboratory colonies, including different laboratory stocks of D. pteronyssinus and other Astigmata mites, samples from house dust containing wild populations of D. pteronyssinus, and in pharmaceutical products used to diagnose HDM allergy.

2 MATERIALS AND METHODS 2.1 Mite cultures

Our D. pteronyssinus stock colony derives from a commercial culture17, 18 and was established in our laboratory in 2012. Inbred D. pteronyssinus colonies were obtained in a former inbreeding project by initial mating of adult virgin couples from the stock colony followed by consecutive brother-sister sib-mating.19 Standard culture conditions and sampling of body and fecal fractions are described as supporting information (Text S1, Figure S1).

2.2 Environmental samples, reference standards, and commercial products

House dust environmental samples were vacuum-collected in mattresses of HDM-sensitized subjects in Northwestern Spain, a region where the prevalence of D. pteronyssinus is well documented.20 Samples from Galicia were collected at the Santiago de Compostela area, and samples from Asturias were collected at Infiesto and Villamayor municipalities. The occurrence of D. pteronyssinus ​in house dust samples was determined by specific molecular methods based on the amplification of ribosomal DNA (rDNA) (Text S1). United States Food and Drug Administration (FDA) reference D. pteronyssinus extracts E10-Dp, E11-Dp, and E12-Dp were provided by the Center for Biologics Evaluation and Research (CBER, Silver Spring, USA). The World Health Organization and International Union of Immunological Societies (WHO/IUIS) D. pteronyssinus international standard 82/518 was provided by the National Institute for Biological Standards and Control (NIBSC, Potters Bar, UK). The following commercial D. pteronyssinus skin-prick testing products from four international suppliers available in the Spanish market were analyzed (in alphabetic order): ALK-Abelló (Madrid, Spain; batch A7438_CT01), LETIPharma (Tres Cantos, Spain; batch 20985617), Roxall Medizin (Zamudio, Spain; batch TCDPT 191115), and Stallergens-Greer (Barcelona, Spain; Alyostal batch 9717109). Note that the results have been shown in the manuscript by using a code in a random order instead of the name of the company (ie, company 1–4).

2.3 In silico analysis of viruses

Viral RNA sequences were identified in unmapped reads resulting from the alignment of in-house D. pteronyssinus WMC RNAseq libraries (four replicates) against the D. pteronyssinus genome assembly ASM190122v2 (RefSeq accession GCF_001901225.1).21 Specific methods to assemble and characterize viral sequences, together with the screening of public Astigmata transcriptomes and phylogenetic analyses, are detailed as supporting info (Text S1).

2.4 Detection of viral RNA by reverse transcription and PCR

The presence of viral RNA in WMC, growth medium, mite body fraction, fecal fraction, reference standards, commercial products, and environmental samples was assessed by RT-PCR. Detailed methods for RNA extraction, reverse transcription and PCR, including to the amplification of a picornavirus (DerpV3) negative (antigenomic) strand, are described as supporting information (Text S1, Table S1).

2.5 Virus visualization by electron microscopy

The localization of virus particles in D. pteronyssinus tissues was studied by Transmission Electron Microscopy (TEM). Virion morphology was inspected by negative staining, and fecal pellets were visualized by Scanning Electron Microscopy (SEM). See Text S1 for detailed methods.

2.6 Transmission of viruses to D. pteronyssinus

An inbred virus-free D. pteronyssinus colony (free of the seven identified viral genotypes as determined by RT-PCR) was inoculated with virus particles derived from the fecal material of the stock colony. The inoculated colony was established and characterized by RT-PCR to detect viral acquisition. Additional culture bioassays were conducted to determine the effects of viruses over mite growth. See Text S1 for methodological details.

3 RESULTS 3.1 Virus discovery

Our analysis identified seven RNA viral contigs within a transcriptome built from unmapped reads obtained from the alignment of D. pteronyssinus WMC RNAseq libraries against the NCBI´s RefSeq genomic assembly of the species.21 BLASTx against the NCBI GenBank´s “RNA viruses” database showed protein identities ranging from 22.8 to 99.1% (BLASTx best hits; Table S2). The identified viral contigs were named Dermatophagoidespretonyssinus virus 1–7, DerpV1 to DerpV7 (Genbank accessions MW355885 to MW355891). Six out of the seven contigs showed similarity to (+)ssRNA viruses: DerpV1/2/3, that were related to the order Picornavirales; DerpV4, to the family Togaviridae; DerpV5, to the families Tombusviridae and Nodaviridae; and, DerpV6, that could not be related to any viral family. In addition, one contig, DerpV7, could be associated to (–)ssRNA viruses of the family Chuviridae. See Figure 1 for a representation of their genomic structure and conserved domains. Further analyses for each viral contig regarding sequence similarity and genomic features are detailed in Text S2, Figure S2, Table S2 and S3. Interestingly, the alignment of our D. pteronyssinus RNAseq libraries from total RNA of WMC against the D. pteronyssinus RefSeq genome assembly resulted in low average alignment rate of 48% (Table 1). The alignment of the resulting unmapped reads to a fasta file containing our seven viral sequences revealed that approximately 75% of these corresponded to viruses. Overall, only 12% of the reads did not align either to the D. pteronyssinus genome or to viruses, hence at least 40% of the total reads corresponded to viruses. As the RNAseq libraries were obtained from WMC, including diet and fecal pellets, we additionally analyzed the presence of viruses on RNAseq libraries obtained from purified bodies of females, males, nymphs, and larvae. Virus-matching reads were found in all developmental instars and adults, but the percentage of reads corresponding to viruses on those libraries was lower (13–30%) than on WMC (Table 1). Finally, the coverage (ie, abundance) of each virus was estimated by total-read-base-counts (ie, the sum of per base read depths) using BedCov, being DerpV3 the most abundant RNA virus in all samples followed by DerpV2, DerpV1, and DerpV4; the less abundant were DerpV5, DerpV6, and DerpV7 (Table S4).

Genomic schematic representation of the RNA viruses identified in this study. (A) DerpV1 (GenBank accession MW355891); (B) DerpV2 (MW355890); (C) DerpV3 (MW355889); (D) DerpV4 (MW355888); (E) DerpV5 (MW355887); (F) DerpV6 (MW355886); (G) DerpV7 (MW355885). Positive strands are represented (5’ to 3’). Arrows indicate open-reading frames (ORF); directions correspond to the orientation of the respective gene. Putative conserved domains are shown as rectangles, an asterisk denotes protein prediction could not established based on structure (Phyre2), but only on InterPro database search. Black squares show poly(A) tails. Empty circles at terminus indicate linearization of a circular molecule. Black and empty triangles show positions of forward and reverse PCR primers, respectively. Acronyms: CAP (mRNA-capping region of the RdRP); HEL (viral RNA helicase); PEP (peptidase C3); RdRP (RNA-dependent RNA polymerase); VMT (viral methyltransferase); VP# (capsid viral proteins)

TABLE 1. Presence of viral-like RNAseq reads in D. pteronyssinus in-house transcriptomes RNAseq librariesa Average total reads/library Alignment rates (%)b To Dp RefSeq assembly Unmapped reads to viruses Reads unmapped to Dp or viruses (%)c Reads corresponding to viruses (%)c WMC (n = 4) 30,929,331 47.5 ± 1.3 75.1 ± 0.8 12.2 ± 0.1 40.4 ± 1.3 Females (n = 3) 24,313,936 72.6 ± 0.6 48.6 ± 0.9 14.4 ± 0.2 13.1 ± 0.5 Males (n = 3) 25,108,996 51.1 ± 4.6 61.4 ± 3.6 18.7 ± 0.3 30.2 ± 4.3 Nymphs (n = 3) 24,677,562 70.6 ± 3.1 54.7 ± 4.8 13.1 ± 0.1 16.3 ± 3.1 Larvae (n = 3) 22,281,857 64.6 ± 2.7 62.1 ± 3.3 13.4 ± 0.3 22.0 ± 3.0 Note Alignment rates were determined, first, by aligning raw RNAseq reads to the RefSeq genomic assembly ASM190122v2 (accessionGCF_001901225.1), and, second, after aligning the resulting unmapped reads to the viral RNA sequences identified in this study (DerpV1-7). Figures are averages ±SEM. “Dp” denotes D. pteronyssinus; “viruses” refers to the fasta file containing all seven viral contigs from this study. a RNAseq libraries were obtained using total RNA extracted from stock whole mite cultures (WMC; n = 4), and from pooled samples of purified bodies at different developmental stages: females (n = 3); males (n = 3); nymphs (n = 3), larvae (n = 3). b "Alignment rate" denotes the percentage of RNAseq reads being mapped to a reference after HISAT2 alignment. c Calculated percentages of RNAseq reads with respect to the total reads. 3.2 Molecular detection of viral RNA in D. pteronyssinus laboratory colonies

A multiplex RT-PCR method was designed for the simultaneous detection of six of the viruses, DerpV1 to 6. For virus DerpV7, a singleplex PCR strategy was adopted due to its undetectable amplification by multiplex PCR. As expected, all seven viruses were detected using these methods on RNA from WMC, containing bodies and fecal pellets (Figure S3). In a next step, the method was applied to detect viruses in RNA samples from purified bodies from our stock colony, and from twelve other D. pteronyssinus colonies established previously after a sib-mating inbreeding program started out from different mite couples of the stock.19 Interestingly, none of the inbred colonies contained all seven viruses identified on the stock colony, indicating that the isolation of mite couples during the inbreeding program had influenced viral transmission. DerpV1 was detected in colonies deriving from the original couple “I”; DerpV1/4 in colonies from “O”; DerpV2/3/4 in colonies from “C”; and no viral RNA was detected in colonies from couples “J” and “M” (Table S5). DerpV7 could not be detected in samples of purified bodies from any of the colonies, including the stock colony. The detection of viral RNA in fecal samples by RT-PCR was only achieved for DerpV3 after applying a modified singleplex protocol (10-fold less RT template) potentially due to inhibition from the fecal RNA extract (data not shown). In addition, for DerpV3, the most abundant virus, the intermediate replicative form (antigenomic RNA, negative strand) could also be detected in samples from WMC and purified bodies from the stock colony using a RT-PCR strategy of high specificity (Figure S4). Finally, attempts failed to detect any of the viruses in RNA extracts from the culture growth medium.

3.3 Detection of virus particles by electron microscopy

Virus particles of about 30 nm in diameter and near-hexagonal/octagonal shape were detected in the digestive tract of D. pteronyssinus, while they could not be found in other tissues, such as the nervous system, muscle, fat body, or gonads. They were specifically located in the cytoplasm of epithelial cells of the midgut. The distribution of viral particles could vary from scattered in the cytosol, including microvilli, enclosed in vesicle-like structures, or forming dense paracrystalline arrays often associated to host cellular membranes (Figure 2A,B,C). Virus particles could also be found at high density in the gut lumen, within non-excreted fecal pellets, specifically into what appeared to be cell remnants shed from the gut epithelium (Figure 2D,E,F). In addition, viral particles were detected by negative staining electron microscopy in fecal pellets directly purified from WMC (Figure 2G,H). No other morphologically distinct viral particles, potentially corresponding to different virus species, were identified.

Detection of virus particles in the digestive system of D. pteronyssinus by electron microscopy. (A) Histological section showing digestive epithelial tissue with the occurrence of virus particles; high dense pockets of particles are indicated with arrowheads and rectangles (close-ups available in panels B, C). (B) Virus particles scattered in the cytoplasm or in association to vesicle-like enclosures. (C) Paracrystalline arrays of virus particles. (D) Histological section of an adult female showing the posterior midgut bearing a fecal pellet (top-right quadrant). The fecal pellet encloses cell remnants containing virus particles. The area corresponding to the rectangle is magnified in panel E; note it was obtained from an adjacent histological cut. (E) Digestive cell remnant heavily infected with virus particles, rectangle close-up in panel F. (F) Virus particles arranged in paracrystalline arrays, note membrane-like structures delimiting the area of infection. (G) Scanning electron micrograph of purified fecal pellets. (H) Negative staining electron micrograph of virus particles released from fecal pellets; white arrowheads indicate individual virions, note their polygonal contour. BL, basal lamina; CR, digestive cell remnant carrying virus; Cu, exocuticle; EC, gut epithelial cell; FP, fecal pellet; IT, intermediate tissue; Lu, gut lumen; Ov, ovipositor; PM, peritrophic matrix; Nu, nucleus; Mt, mitochondrion; Mu, visceral muscle; Mv, microvilli. Scales: A = 5 µm; B = 200 nm; C = 200 nm; D = 40 µm; E = 3 µm; F = 500 nm; G = 20 µm; H = 50 nm

3.4 Transmission of RNA viruses and effects on mite population growth

In order to assess the potential of feces as a route of transmission of viruses to other mites, fecal pellets purified from the virus-containing stock colony were used to inoculate mites separated from the inbred colony J-1–1–1–1, which was free of the seven viruses identified in this study (ie, virus-free; Table S5). After two and six culture cycles since the initial fecal spike (81 and 268 days, respectively), the presence of viral RNA in the inoculated and the original inbred colonies was assessed by multiplex RT-PCR, results are shown in Figure 3. As expected, none of the studied viruses was detected in the original non-inoculated J-1–1–1–1 colony. However, positive detection of viral RNA was observed for DerpV1/3/5/6 by multiplex PCR on inoculated J-1–1–1–1 WMC, as well as on separated bodies from this colony, showing the same results after two and six culture cycles. The presence of DerpV2, but not DerpV4 or DerpV7, was further confirmed by singleplex RT-PCR in all samples from the inoculated colony. Additionally, the overall effect of the identified viruses on the fitness of D. pteronyssinus was evaluated by comparing the growth of small-scale populations from adults of the original virus-free inbred colony and from the inoculated colony. Mite population numbers after 3 weeks increased over 10-fold for both colonies, but no significant difference was recorded between cultures on neither total mites (p = .169; Figure 4A), nor the percentage of laid eggs (p = .066; Figure 4B). Likewise, no significant difference was found on the duration of the standard culture cycle in flasks (p = .452; Figure 4C).

Fecal-mediated transmission of RNA viruses from the stock D. pteronyssinus colony into the inbred colony J-1–1–1–1. Six virus species (DerpV1-6) were screened by multiplex RT-PCR and agarose gel electrophoresis. RNA samples were extracted from both whole mite cultures (WMC) or purified mite bodies (PMB). J-1–1–1–1 colonies were sampled at the 2nd and 6th culture cycle after the initial inoculation with fecal pellets from the stock colony; both non-inoculated (N-FI) and inoculated (FI) cultures were analyzed. A no template control was included (NTC). The molecular-size marker was GeneRuler™ 100 bp Plus. Amplicon sizes are (in decreasing order; use lane “WMC” as reference): 716 bp (DerpV5), 622 bp (DerpV3), 549 bp (DerpV2), 468 bp (DerpV1), 372 bp (DerpV4), 271 bp (DerpV6)

Impact of the inoculation of viruses in the population growth of a D. pteronyssinus colony. The inbred colony J-1–1–1–1, reported as virus-free by molecular methods, was inoculated with fecal-derived viruses from the stock colony, virus-inoculated. The growth parameters of both colonies were compared. (A) Total mite population after 3 weeks of small-scale culture (n = 14 microtubes). (B) Percentage of non-hatched eggs in small-scale populations after 3 weeks (n = 14). (C) Duration of the standard-scale culture cycle (in 2 g flasks; n = 9 cycles; first cycle not included). Columns indicate mean values; errors bars indicate standard errors of the mean; means were statistically compared by unpaired t test (n.s. denote “not significant,” p > .5)

3.5 Detection of viral RNA in other sources

Nine environmental samples obtained from house dust were screened. First, the occurrence of D. pteronyssinus was detected by rDNA amplification, showing detectable levels in seven out of the nine samples. Then, viral presence was assessed by RT-PCR, showing two positive samples for at least one virus, including DerpV1 and DerpV2 (Table 2; Figure S5). Additionally, four commercial allergy diagnostic products and four international reference standard extracts from D. pteronyssinus were analyzed by RT-PCR. Remarkably, all the pharmaceutical products and the three CBER/FDA standards were positive for at least one virus (Table 2; Figure S5). Likewise, RNA from all the virus species identified in this study, except for DerpV7, could be detected in at least one of the commercial extracts or international standards. The two more prevalent viruses were DerpV1 and DerpV3, with six and five positives out of the eight analyzed samples, respectively. Direct sequencing of amplification products for DerpV1/2/3 from all different origins revealed that each individual sample exhibited a distinct variant on the amplified region (468, 549, and 622 bp, respectively; Table S1, Figure 1), with sequence identities to our assembled viral sequence ranging from 95.3–97.0%, 95.1–95.3%, and 92.1–94.2%, respectively (Figure S6,S7,S8). In order to further explore the presence of viruses similar to those in this study, public transcriptome databases were additionally screened revealing the occurrence of homologous sequences to DerpV1/3/4/6/7 in two other D. pteronyssinus laboratory stock cultures (accessions SRX2026624, SRX4482484)22, 23 (Table S6). Besides, similar sequences to DerpV1/3/6 were also found in RNAseq data from the HDM D. farinae and the storage mite T. putrescentiae, albeit their coverage was insufficient to assess viral similarity between host species, except for one case. Alignment of RNAseq reads suggested the occurrence of a DerpV1-like virus in T. putrescentiae (approximately 88% of DerpV1´s sequence with 91% identity). Finally, no sequence similarity was found in other mite species of the order Astigmata.

TABLE 2. Detection of D. pteronyssinus-derived viral RNA in allergy international standards, commercial pharmaceutical products, and house dust samples by RT-PCR International reference standards Commercial pharmaceutical products House dust samples

CBER/FDA

E10-Dp

CBER/FDA

E11-Dp

CBER/FDA

E12-Dp

WHO/IUIS

82/518

Comp.

1

Comp.

2

Comp.

3

Comp.

4

Asturias

1

Asturias

2

Asturias

3

Galicia

1

Galicia

2

Galicia

3

Galicia

4

Galicia

5

Galicia

6

DerpV1 +b +b +b +b +b +b +b +b DerpV2 +b +b DerpV3 +b +b +b +b +b DerpV4 + DerpV5 + + DerpV6 + + + DerpV7 rDNA Dpa n.a. n.a. n.a. n.a. n.a. n.a. n.a. n.a. + + + + + + + a The presence of D. pteronyssinus (Dp) in house dust samples was assessed by the specific amplification of ribosomal DNA (rDNA). N.a. denotes "not analyzed". b RT-PCR products were analyzed by Sanger sequencing with forward and reverse specific primers. 4 DISCUSSION

The virome of D. pteronyssinus resolved in this study comprised seven different viral contigs, all of them corresponding to RNA viruses. RNA sequences associated to the expression of putative DNA viruses have not been identified in our analysis. Most of the identified viruses show relationship with arthropod-infecting taxa and/or their closest virus genotypes are associated to arthropods or invertebrates, but not to mammals. Five out of the seven viral contigs showed >95% protein sequence identity to viral accessions from a similar D. pteronyssinus meta-transcriptomic survey published during the preparation of this manuscript.12 The similarity of DerpV6 and DerpV7 to known viruses in databases was very low, suggesting they are putatively novel (extended details in Text S2). The three most abundant viruses, DerpV1/2/3, were related to the order Picornavirales. Among them, DerpV1 and DerpV3 showed capsid domains with structural homology to those from formal members of the picornavirus family Dicistroviridae, as well as phylogenetic relationship to informally classified dicistro-like viruses. However, their monocistronic genomic organization differs from all currently assigned dicistroviruses, which are bicistronic with separate non-structural and structural ORFs.