Imbalance of the intestinal virome caused by a conditional deletion of the vitamin D receptor from epithelial cells, Paneth cells, or myeloid cells in murine models

Background Vitamin D receptor (VDR) is highly conserved in mammals, and its deciency is associated with various diseases (e.g., cancer, infection, and chronic inammation). Prior research has described the effect of VDR on bacteria; however, little is known regarding the effects of VDR on viruses. We hypothesize that VDR is a regulator of the virome and virus-bacterial interactions. We conditionally deleted VDR from intestinal epithelial cells (VDR ΔIEC ), Paneth cells (VDR ΔPC ), and myeloid cells (VDR ΔLyz ) in mice. Feces were collected for shotgun metagenomic sequencing and metabolite proling. To test the functional changes, we evaluated pattern recognition receptors (PRRs) and analyzed microbial metabolites.

infection [14], and hepatitis C virus infection [16]. An epidemiological study showed an inverse association between vitamin D levels in the serum and infections of the upper and lower respiratory tracts [17]. Vitamin D de ciency also contributes to the pathogenesis of HIV infection by negatively modulating innate and adaptive immune responses [18,19]. At the gene level, a meta-analysis determined that the VDR gene polymorphism FokI is consistently associated with host susceptibility to infection by respiratory syncytial virus [20]. Gene variation in VDR was also reported to correlate with persistent hepatitis B infection in clinical patient samples [21]. Consistent with the epidemiological data, activation of the vitamin D/VDR pathway in response to respiratory syncytial virus infection occurs in lung cells via TLR3 signaling, which results in the induction of the CAMP gene [22,23].
The microbiota consist of bacteria, viruses, fungi, multicellular parasites and archaea. The virome is a collection of nucleic acids, both RNA and DNA that compose the viral community associated with a particular ecosystem of microbiota. The virome includes eukaryotic viruses, endogenous retroviruses, bacterial viruses (i.e., bacteriophages), and archaeal viruses and is one of the least understood components of the microbiota [24]. Dysbiosis of the microbiome not only leads to intestinal in ammatory and infectious diseases but also to diseases beyond the gastrointestinal tract [5,25]. At present, we are beginning to understand the in uence of VDR on microbial homeostasis, which is critical in various diseases [5,[26][27][28][29]. However, the effects and mechanisms of VDR on the virome have not been fully elucidated.
In this study, we hypothesize that host factors (e.g., VDR status in speci c tissues) are regulators of the intestinal virome, virus-bacterial interactions, and microbial metabolites. Speci cally, we used novel animal models and statistical and bioinformatic tools and models to understand host factors and aspects of metabolites and the microbiome, including bacteria and viruses, in vivo. We conditionally deleted VDR from intestinal epithelial cells (VDR ΔIEC ), Paneth cells (VDR ΔPC ), and myeloid cells (VDR ΔLyz ) of mice. Fecal samples were collected for microbial DNA extraction and shotgun metagenomic sequencing. From a functional perspective, we investigated the virus infection-associated metabolites and virome-related PRRs of colonic epithelial cells. Understanding the role of VDR in altering the multikingdom aspects of the microbiome, including bacteria and viruses, may help to elucidate the mechanisms by which VDR regulates microbial homeostasis in terms of health and disease.

Experimental animals
We conditionally deleted vitamin D receptor from intestinal epithelial cells, Paneth cells, and myeloid cells of the mice by crossbreeding VDR LoxP mouse and cre mouse. The VDR LoxP mice were formerly developed by Dr. Geert Carmeliet [30]. VDR ΔIEC , VDR ΔPC and VDR ΔLyz mice were obtained by crossing the VDR LoxP mice with villin-cre, Defa6-cre and Lyz-cre mice, respectively. Villin-cre and Lyz-cre mice were purchased from Jackson Laboratories (Bar Harbor, ME, USA). DEFA6-cre mice were from Dr. Richard Blumburg [31].
Conditional knockout (KO) mice were derived from heterozygous mating pairs such that wild-type and conditional KO mice came from the same litter. The same breeding method was used for the VDR LoxP mice.
The mice aged 6 to 8 weeks were chosen from each group and cohoused until the experiments were performed. All the animals were housed in the Biologic Resources Laboratory (BRL) at the University of Illinois at Chicago (UIC) and utilized in accordance with the UIC Animal Care Committee (ACC) and O ce of Animal Care and Institutional Biosafety (OACIB) guidelines. The animal work was approved by the UIC O ce of Animal Care (ACC 15231,17-218, and 18-216).

Real-time quantitative PCR
Total RNAs were extracted from colon epithelial cells of four genotype mice of VDR LoxP , VDR ΔIEC , VDR ΔPC and VDR ΔLyz , using TRIzol reagent (Life technologies, Carlsbad, CA, USA). RNAs were rst reversetranscribed into cDNA with the iScript cDNA synthesis kit (Bio-Rad, Hercules, CA, USA) according to the manufacturer's manual. The CTFX 96 Real-time system (Bio-Rad, Hercules, CA, USA) and SYBR Green Supermix (Bio-Rad, Hercules, CA, USA) were used for the quantitative real-time PCR with the RT-cDNA reaction products according to the manufacturer's directions. All expression levels were normalized to βactin levels of the same sample. Percent expression was calculated as the ratio of the normalized value of each sample to that of the corresponding untreated ones. All real-time PCR reactions were performed in triplicate. Optimal primer sequences target for PRRs were obtained from Primer Bank [32] (http://pga.mgh.harvard.edu/primerbank/) as listed in Table 1.

Fecal sample collection and shotgun metagenomic sequencing
We used whole-genome shotgun sequencing to sequence fecal samples. Fresh fecal samples from each group (VDR LoxP : male n=3 and female n=7; VDR ΔIEC : male n=5 and female n=5; VDR ΔPC : male n=5 and female n=5; VDR ΔLyz : male n=5 and female n=5) were collected from the colon and placed into the sterile tubes. The samples were kept at low temperature with dry ice and were sent to the University of Illinois at Chicago Research Resources Center for genomic sequencing. The DNAs of samples were extracted using DNeasy Power Fecal Kit (Qiagen, Hilden, Germany) based on manufacturer's instructions with a slight modi cation. The samples were heated at 65°C for 10 min before homogenizing with FastPre-24 5G bead-beating device (MP Biomedicals, Solon, OH, USA) at 6 m/s for 40 s. These genomic DNA was fragmented into relatively small pieces (generally 250-600 bp fragments) prior sequencing. Sequencing was performed using a Illumina HiSeq system, as in our previous publications [31] [33]. Following quality checking, ltering the reads, removing noisy sequences, metagenomic assembly [34], the resulting assemblies were ltered to exclude contigs shorter than 1,000 nucleotides and all remaining set of DNA reads (or contigs) were classi ed with Centrifuge [35], an e cient metagenomic classi er capable of indexing the entirety of nucleotide (nt). Next, the taxonomic annotation of each contig (genes) was obtained by searching for the comprehensive NCBI GenBank non-redundant nucleotide database (as described in https://merenlab.org/2016/06/18/importing-taxonomy) ([36]). . After removing the identical sequences with ≥99% identity of each other to make it nonredundant, 289,629,756 readable sequences were yielded in 40 studied samples with an average of 7,240,744 reads per sample. Of these sequences, a total of 100,013,480 reads were taxonomic alignments with an average of 2,500,337 reads per sample.
Metabolite sample preparation and metabolomic data analysis Sample preparation and process. As in our pervious report [33], metabolite pro ling was performed on Metabolon platform (Metabolon, Inc., Durham, NC, USA). Brie y, fecal samples from mouse groups (VDR LoxP : male n=6 and female n=10; VDR ΔIEC : male n=8 and female n=9; VDR ΔLyz : male n=5 and female n=5) were maintained at -80 o C and were prepared using the automated MicroLab STAR® system from Hamilton Company. Several recovery standards were added prior to the rst step in the extraction process for QC purposes. The proteins were precipitated with methanol under vigorous shaking for 2 min (Glen Mills GenoGrinder 2000) followed by centrifugation aim to remove protein, dissociate small molecules bound to protein or trapped in the precipitated protein matrix, and to recover chemically diverse metabolites. Next, samples were placed on a TurboVap® (Zymark) brie y to remove the organic solvent. The sample extracts were stored overnight under nitrogen before preparation for analysis. The resulting extract was divided into ve fractions for analysis by four methods using Ultrahigh Performance Liquid Chromatography-Tandem Mass Spectroscopy (UPLC-MS/MS) (see details in the bioinformatic analysis below) and one for backup.
Bioinformatic analysis of metabolomic data. Bioinformatic analysis was performed using the Metabolon informatic system (the LAN backbone, and a database server running Oracle 10.2.0.1 Enterprise Edition), which consists of four major components, including the Laboratory Information Management System (LIMS), the data extraction and peak-identi cation software, data processing tools for QC and compound identi cation, and a collection of information interpretation and visualization tools for use by data analysis. The hardware and software foundations for these informatic components were the LAN backbone and a database server running Oracle 10.2.0.1 Enterprise Edition.
Data extraction and compound identi cation. After raw data extraction, peak-identi cation and QC process using Metabolon's hardware and software, compounds were identi ed by comparison to library entries of more than 3300 currently commercially available puri ed standard compounds or recurrent unknown entities. Metabolon maintains a library based on authenticated standards that contains the retention time/index (RI), mass to charge ratio (m/z), and chromatographic data (including MS/MS spectral data) on all molecules present in the library. Furthermore, biochemical identi cations are based on three criteria: retention index within a narrow RI window of the proposed identi cation, accurate mass match to the library +/-10 ppm, and the MS/MS forward and reverse scores between the experimental data and authentic standards.
Metabolite quanti cation and data normalization. Peaks were quanti ed using area-under-the-curve. A data normalization step was performed to correct variation resulting from instrument inter-day tuning differences as necessary or for purposes of data visualization. Additionally, biochemical data were normalized to other factors (e.g., cell counts, total protein as determined by Bradford assay, and osmolality) to account for differences in metabolite levels due to differences in the amount of material present in each sample.

Statistical analysis
All the tests as speci ed in related analysis performed in study were two-sided. A P-value ≤0.05 was considered statistically signi cant. For a large number of tests, the P-values were adjusted to account for false positives using the False Discovery Rate (FDR) [37] and the q-value (FDR-corrected p-value) ≤0.05 was reported to be signi cant. For all the data, we ran Shapiro-Wilk test to determine whether parametric ANOVA or non-parametric Kruskal-Wallis test is used to group comparisons.
Alpha (within-sample) diversity (e.g., Shannon diversity) measures the number (richness) and/or distribution (evenness) of species within a single sample, while beta (between-sample) diversity (e.g., Bray-Curtis dissimilarity) measures the differences of microbial composition in one sample compared to another [38]. For the metagenomic sequencing data, the raw read counts were normalized to Counts Per Million (CPM). Shannon diversity was used to examine the diversity of the gut microbiome in conditional VDR KO and control mice and the differences between groups were tested with Kruskal-Wallis test.
Bray-Curtis dissimilarity index was used to detect the differences or dissimilarities among studied groups. We rst performed principal coordinate analysis (PCoA) to visualize the Bray-Curtis dissimilarities among groups. Next, we performed the nonparametric PERMANOVA to evaluate whether the group (VDR ΔIEC , VDR ΔPC , and VDR ΔLyz mice vs. VDR LoxP mice) has a signi cant effect on overall gut microbiota composition. Next, we performed pairwise PERMANOVA using RVAideMemoire package. Finally, we conducted analysis of similarity (ANOSIM), another nonparametric procedure, based on a permutation test for rank dissimilarities among between-and within-groups to con rm the hypothesis testing results from those of PERMANOVA and obtained the pairwise-comparison results. The Spearman correlation analysis of viruses and bacteria was performed via the Hmisc package [39].
The statistical analysis of metabolites was performed using Welch's two-sample t-test on log transformed data. Welch's two-sample t-test is used to test whether two unknown means are different from two independent populations, which enables unequal variances and has an approximate t-distribution with degrees of freedom estimated using Satterthwaite's approximation. The log 2 fold-changes were reported as the effect sizes in comparisons between groups. The statistical analyses for microbiome and metabolomic data were performed using R [40], R packages of ampvis2, microbiome, phyloseq and vegan [41], as well as ggplot2 [42] and ggpubr packages [43].

Results
Overall virome community composition in conditional VDR knockout mice In all fecal samples, 4,048 viral species were identi ed. Distinct viral compositions and communities, which differed in both diversity and composition, were present between the control mice and the conditional VDR knockout mice ( Figure 1). We presented the 10 most abundant viral species for individual samples in each group ( Figure 1). The most abundant species in all the groups were Vibrio phage JSF5, Vibrio phage JSF6, bovine viral diarrhea virus 1, Escherichia coli O157 typing phage 7, human alphaherpesvirus 1, Lactobacillus prophage Lj771, in uenza A virus, Lactobacillus phage KC5a, and viruses incertae sedis unidenti ed phage ( Figure 1A). When individually analyzing the species, Vibrio phages JSF5 and JSF6 were signi cantly more abundant in all three VDR KO mouse models compared to control VDR LoxP mice. Escherichia coli O157 typing phage 7 was found to be signi cantly more abundant in VDR ΔPC mice, compared to the control mice. Meanwhile, Lactobacillus phage phiadh was considerably less abundant in VDR ΔLyz mice compared to control VDR LoxP mice ( Figure 1B) (P<0.01).
Altered diversity of the virome in conditional VDR knockout mice Shannon diversity is commonly used to characterize species diversity in a community [38]. We found that the three conditional VDR knockout groups had lower Shannon diversities than the control mice. The Shannon diversity at the viral species level was signi cantly lower in VDR ΔIEC , VDR ΔPC , and VDR ΔLyz mice than in VDR LoxP control mice (P=0.05, P=0.05, and P=0.04, respectively) ( Figure 2A).
The Bray-Curtis dissimilarity index was used in this study to measure the dissimilarities of samples. We rst performed PCoA and found viral dissimilarities between conditional VDR knockout mice and control mice ( Figure 2B). VDR ΔIEC mice partially overlapped with control VDR LoxP mice, whereas VDR ΔPC and VDR ΔLyz mice were completely separated from VDR LoxP mice. The group differences can be explained by a total of 37.2% (29.3% + 7.9%) of the variations among the animals.
Next, we performed nonparametric PERMANOVA to evaluate whether VDR status impacts the overall intestinal viral pro le. The sequential test "Group and Gender" showed that the dissimilarities among groups were signi cantly different (P=0.011). Because the overall dissimilarity among groups was signi cantly different, we performed a pairwise PERMANOVA and found that the Bray-Curtis dissimilarities of viruses of the VDR ΔIEC , VDR ΔPC , and VDR ΔLyz mice were signi cantly different from the dissimilarity in the VDR LoxP mice. Furthermore, the differences in dissimilarity among groups were con rmed by ANOSIM (analysis of similarity), where the rank dissimilarities between and within groups were signi cantly different (P=0.001). VDR ΔIEC , VDR ΔPC , and VDR ΔLyz mice had higher dissimilarities (i.e., lower similarity) than VDR LoxP mice ( Figure 2C).
Altered abundance of the virome in conditional VDR knockout mice We observed that a total of 12 viral species were differential in the conditional VDR knockout mice compared to the control VDR LoxP mice, of which 6 had q-values <0.001, 2 had q-values <0.01 and 4 had qvalues <0.05 ( Figure 3A).  Figure 3A).

Sex-based differences in the gut virome altered by VDR status
To investigate the impact of sex on the alteration of the virome community in the intestines of the studied mice, we illustrated signi cantly different virus species abundances in male and female mice with logarithmic fold-changes and q-values ( Figure 3B). Seven species were detected to be signi cantly altered in the female VDR ΔIEC /VDR LoxP comparison but not when comparing males of the same two groups, the altered species included 2 enriched species, Vibrio phages JSF5 and BVDV1, and 5 depleted species, Lactobacillus prophage phiadh, Lactobacillus prophage KC5a, Lactobacillus prophage Lj771, macacine alphaherpesvirus 1 and Catovirus CTV1. However, only Vibrio phage JSF5 was enriched in female VDR ΔPC mice compared to control VDR LoxP mice ( Figure 3B), which is the same as in the group-factor analysis. Five species in the female VDR ΔLyz /VDR LoxP comparison were found to be signi cantly differential (q<0.05), including 2 enriched species:, Vibrio phage JSF5 and bovine viral diarrhea virus 1, and 3 depleted species, Lactobacillus prophage phiadh, Cherry green ring mottle virus and Lactobacillus prophage KC5a, while in males, only Vibrio phage JSF5 was enriched and only Lactobacillus prophage phiadh was depleted ( Figure 3B). Overall, more altered viral species in all conditional VDR knockout mice were found in the female mice than in male mice, which indicated the impact of VDR status on sex differences.
VDR deletion led to signi cantly enriched Vibrio phages Vibrio phages target Vibrio cholerae bacteria, which can secrete cholera toxin and cause watery diarrhea in patients [44,45]. In this study, we found a markedly enriched abundance of Vibrio phage JSF5 and Vibrio phage JSF6 in VDR ΔIEC , VDR ΔPC , and VDR ΔLyz mice compared with control VDR LoxP mice ( Figure 1). Moreover, differential analysis found that Vibrio phage JSF5 was signi cantly enriched in all three conditional VDR-knockout mice (q<0.01) compared to the control mice, and Vibrio phage JSF6 was also more enriched in VDR ΔLyz mice in comparison to VDR LoxP mice (q<0.01) (Figure 3). When the sex factor was included in the analysis, the fold-changes of Vibrio phage JSF5 were found to be signi cant in both female and male VDR ΔLyz mice but only in female VDR ΔIEC and VDR ΔPC mice ( Figure 3).

Altered bacterial abundance in conditional VDR knockout mice
We reasoned that the abundance of bacteria in the intestines of the mice should be altered due to the dysbiosis of bacteriophages. In this study, we found that Vibrio cholerae, the host of Vibrio phage JSF5 and Vibrio phage JSF6, which were largely enriched in conditional VDR knockout mice ( Figure 1B), was detected at low concentrations in knockout mice compared with control mice ( Figure 4A). Signi cantly lower abundance of this bacterium was found in VDR ΔLyz mice, which matches the nding that both Vibrio phage JSF5 and Vibrio phage JSF6 were depleted in VDR ΔLyz mice. A similar situation was also found in Lactobacillus gasseri, the bacterial host of Lactobacillus prophage phiadh and Lactobacillus prophage KC5a ( Figure 1B; Figure 4A). Furthermore, as the target host of Escherichia coli O157 typing phage 7, E. coli was found to have the opposite phage species abundance in the mice studied here ( Figure 1B; Figure 4A). The altered bacterial abundance in the microbial community supports our ndings of virome changes.
To further evaluate the altered bacteria, we showed the altered bacterial species with q-values <0.1 in differential analysis in our three mouse models ( Figure 4B). We found that Haemophilus ducreyi and Mesorhizobium huakuii were signi cantly depleted in VDR ΔIEC mice compared with control VDR LoxP mice ( Figure 4B). Meanwhile, ve bacterial species were signi cantly depleted, and three bacterial species were signi cantly enriched in VDR ΔPC mice, compared with the control. These depleted species were Haemophilus ducreyi Kushneria konosiri, Microbacterium sp. LKL04, Isosphaera pallida, and Actinomyces radingae. Three enriched bacterial species were Bacteroides uniformis, Faecalibaculum rodentium, and Cutibacterium acnes (q<0.01. Figure 4B). Furthermore, we found that nine bacterial species were signi cantly depleted, and three bacterial species were enriched in VDR ΔLyz mice. These depleted bacterial species were Bi dobacterium pseudolongum, Bi dobacterium choerinum, Bi dobacterium animalis, Bordetella pseudohinzii, Haemophilus ducreyi, Clostridium perfringens, Streptomyces peucetius, and Bacteroides acidifaciens. The three enriched bacterial species were Ralstonia solanacearum, Chroococcidiopsis thermalis, and Cutibacterium acnes ( Figure 4B). The less abundant Haemophilus ducreyi is a gram-negative bacterium and causative agent of genital ulcer disease chancroid [46] and was detected in the three conditional VDR knockout mice. Similar to viruses, more bacterial species were changed in VDR ΔPC and VDR ΔLyz mice compared to the control mice ( Figure  4B).

Correlation of viral and bacterial alterations in conditional VDR KO mice
Bacteria and viruses are essential for protective, metabolic, and physiological functions. We found that the viral (consisting mostly of bacteriophages) and bacterial species abundances were altered in conditional VDR KO mice. To investigate the interactions between bacteria and viruses, we performed a correlation analysis of viruses and bacteria. All the bacterial and viral species with signi cant fold changes in the differential analysis (q<0.05) were included in the correlation analysis ( Figure 4C and Supplement Table 1).
In VDR ΔIEC mice, bovine viral diarrhea virus 1 was signi cantly negatively correlated with the bacteria Mesorhizobium huakuii (correlation coe cient (r)=-0.63; P=0.05) ( Figure 4C). M. huakuii induces the formation of nitrogen-xation nodules on its host plant Astragalus sinicus and has been assigned to a new biovariety based on its host range and taxonomic characteristics [47]. M. huakuii isolates were also found to have endotoxic activity against lipopolysaccharides [48]. In VDR ΔPC mice, the Vibrio phage JSF5 was dramatically positively correlated with the bacteria Cutibacterium acnes (r=0.98; P<0.0001) ( Figure  4C). Cutibacterium acnes, formerly known as Propionibacterium acnes, is an anaerobic, aerotolerant, bacillus-shaped bacterium. It is ubiquitously found as a commensal on the surface of skin in areas rich in oleic and palmitic acids [49]. Moreover, variable results for the association between C. acnes and ulcerative colitis were found [50].  Figure 4C). In addition, we found a positive correlation between Vibrio phage JSF5 and Cutibacterium acnes in VDR ΔPC mice and VDR ΔLyz mice. Taken together, these data indicate the critical role of viral and bacterial interactions in intestinal microbial homeostasis with the support of VDR.

VDR status altered the expression of PRRs in colonic epithelial cells
To examine the impact of virome dysbiosis, virus-related receptors in the colon were evaluated ( Figure 5). TLR3 and TLR7 are transmembrane PRRs located in endosomes that recognize nucleic acids and mediate cell extrinsic virus recognition [51]. We observed that expression of TLR3 and TLR7 was upregulated in the conditional VDR knockout mice compared with the control mice, with a signi cant difference in the VDR ΔLyz group vs. the control group (P<0.05) ( Figure 5).
NLRs characterized by the presence of a conserved NOD motif, comprise a large receptor family, including NOD1, NOD2, and NLRPs [52,53]. NLRs are activated not only in response to viruses but are also important modulators of other virus sensing pathways [51]. Here, we observed upregulated expression of NOD1 and NLRP6 and signi cantly increased NOD2 RNA in VDR ΔLyz mice compared with control VDR LoxP mice (P<0.05) ( Figure 5). C-type lectin receptors (CLRs) also mediate cell-extrinsic sensing of speci c viruses by binding viral glycans [54]. Compared with the control VDR group, the expression of CLEC4L, one of the CLRs, was signi cantly upregulated in all three conditional VDR knockout mice, especially in VDR ΔIEC and VDR ΔPC mice ( Figure 5). Overall, the signi cant alterations of PRRs in conditional VDR knockout mice suggest the in uence of VDR on intestinal homeostasis and expression of PRRs.

VDR status altered metabolites related to bacteriophage infection in feces
Microbiota-derived metabolites are chemical messengers that elicit a profound impact on host physiology. Thus, we investigated bacterial metabolites that are related to bacteriophage infection ( Figure   6). We found that glycolysis glucose was signi cantly decreased in VDR ΔIEC mice compared with the control (P<0.05), while ribulose/xylulose and xylose were signi cantly increased in VDR ΔLyz mice compared with the control (P<0.05) ( Figure 6A). Moreover, most of the long-chain fatty acids were signi cantly increased in both VDR ΔIEC and VDR ΔLyz female mice compared to control mice, while some fatty acids were only increased in VDR ΔIEC mice. 10-Hydroxystearate, which is related to phage infection, was decreased in VDR ΔLyz mice and increased in female VDR ΔIEC mice compared to control mice (P<0.05) ( Figure 6B). Metabolic alterations of nucleotides were also detected in both conditional knockout mice, such as uridine in VDR ΔIEC mice and 3-ureidoisobtyrate in VDR ΔLyz mice (P<0.05) ( Figure 6C). Similarly, we observed amino acid alterations in the feces of conditional VDR knockout mice when compared to those of the control mice. For instance, we observed decreased phage infection-related serine in both VDR ΔIEC and VDR ΔLyz female mice and decreased cellular virus infection-related glutamine in VDR ΔIEC and VDR ΔLyz female mice compared to the control ( Figure 6D).

Discussion
In the current study, we report that conditionally deleted VDR from intestinal epithelial cells, Paneth cells, and myeloid cells signi cantly altered the virome pro le and virome-bacterial interactions in mice. The changes in Vibrio phages, Lactobacillus phages, and E. coli typing phage were signi cantly related to VDR status. Seven more virus species were signi cantly changed in VDR ΔLyz mice compared to control mice, while BVDV1 was signi cantly enriched in VDR ΔIEC mice. The altered viral abundance was greater in female mice than in male mice. Bacteriophages modulate their hosts directly by affecting their mortality and horizontal gene transfer or indirectly by impacting target/host bacteria and hence may alter the microbiome. After conditionally deleting VDR, the virome and other aspects of the microbiota were changed, which further led to microbial dysbiosis. From a functional perspective, the signi cant changes in PRRs and metabolites related to infection suggest the in uence of VDR on intestinal virus homeostasis.
Vibrio phages JSF5 and JSF6 are organisms hosted by Vibrio cholera, which can secrete cholera toxin and cause watery diarrhea in patients [44,45]. The Lactobacillus phages phiadh and KC5a are members of the viral family Siphoviridae and can preferentially infect the bacteria Lactobacillus gasseri, which is a common component of intestinal mucosae and plays an important role in modulating the gut immune system [55]. These observations were consistent with our results that VDR deletion in myeloid cells had a more severe in uence on the intestinal virome balance. Similarly, the interactions in the microbial community are important to maintain homeostasis and host health. For instance, viruses can bind to bacterial products (e.g., lipopolysaccharide or HBGA (human blood group antigen)-like substances), increase virome stability and protect virome from physical stresses [56].
We found that the Lactobacillus phage KC5a and E. coli O157 typing phage 7, which belong to the family Myoviridae, and Lactobacillus phage phiadh, which belongs to the family Siphoviridae, were altered in our conditional VDR knockout mouse models compared with the control mice. A study indicated that the CD susceptibility gene ATG16L1 was phenocopied in mice infected with murine norovirus [57][58][59]. Meanwhile, we have demonstrated that VDR transcriptionally regulates ATG16L1, and VDR status may be a determinant of IBD risk through its actions on ATG16L1 [26,60]. However, these results may support a model in which bacteriophages contribute to the development of bacterial dysbiosis associated with IBD. Bacteriophages modulate their hosts directly by affecting their mortality and horizontal gene transfer and by further altering the components of the intestinal community and contributing to dysbiosis. The virome could be a candidate biomarker for human IBD [61] because CD and UC were associated with a signi cant expansion of Caudovirales bacteriophages. The phage families of the Caudovirales order, including the Siphoviridae, Myoviridae and Podoviridae, were found to be altered in murine colitis and human IBD patients [62]. Therefore, further studies are needed to elucidate the important and real role of the intestinal virome in the development and progression of IBD.
Paneth cells play an integral role in shaping the microbiome and host defenses. The absence of VDR in Paneth cells impairs antimicrobial function, affects microbial assemblage, and increases susceptibility to colitis and infection [26,31]. The altered virome in VDR ΔPC mice further indicated the important role of VDR and Paneth cells in shaping the intestinal microbiome and secreting antimicrobial peptides or metabolites. Paneth cells may indirectly impact the intestinal microbiota or the outcome of a viral infection by adjusting the population of bacteria. Paneth cells also play a key role in host defense by sensing microorganisms through TLRs [63]. It is clear from our data that Paneth cell dysfunction leads to dysbiosis and a compromised epithelial barrier. Our previous reports demonstrated that VDR could promote healthy microbial metabolites and a healthy microbiome to prevent obesity [33]. Furthermore, it has been reported that increased serum 25(OH)D was associated with increased bene cial bacteria, such as Parabacteroides, which were suppressed in IBD patients and were altered in an obesity mouse model [33,64]. We have reported that human VDR gene variation determines the abundance of Parabacteroides [28]. In this study, we also found a decreased abundance of Parabacteroides in three VDR knockout mice (VDR ΔPC , VDR ΔIEC and VDR ΔLyz ) compared with VDR LoxP controls. This further con rms the critical role of VDR in shaping the microbiome at the genetic level.
Eukaryotic viruses, such as adenoviruses, hepatitis B virus, hepatitis C virus and human immunode ciency virus [65], are also present in the intestinal viromes of some individuals, which indicates the potential infectious capability of these viruses in the host. BVDV1 can infect mice and cause signi cant histopathological damage [66]. The signi cantly increased abundance of intestinal BVDV1 in VDR conditional knockout mice suggested the importance of VDR in viral infection and survival inside the host. Furthermore, the association between vitamin D de ciency and the pathogenesis and course of HIV disease has also been recognized recently. Infants born from HIV-infected women with vitamin D de ciency are at an increased risk of infection and a decreased survival rate [67]. In addition, the VDR polymorphism FokI was evaluated using a meta-analysis to be consistently associated with susceptibility to infection to respiratory syncytial virus [20].
The innate immune response is nonspeci c, is the rst line of defense against infectious agents and initiates antigen presentation, including responses to viral infection. PRRs and related pathways involved in intestinal virus sensing, such as TLRs, are mainly mediated by cell-extrinsic virus recognition [68]. We found signi cant alterations in PRRs, including upregulated expression of TLR3, TLR7 and NOD2 in VDR ΔLyz mice and increased expression of CLEC4L in VDR ΔIEC and VDR ΔPC mice, suggesting the in uence of VDR on intestinal virus homeostasis. Our data further suggest that the impacts of tissue-speci c VDR deletion on intestinal receptors are different. For instance, VDR deletion in myeloid cells would signi cantly increase TLR3, TLR7 and NOD2, while VDR deletion in intestinal epithelial cells would more likely increase CLEC4L levels. TLRs are known to be affected by VDR in monocytes and epidermal keratinocytes [12]. In mice administered an antiviral drug cocktail, the depletion of gut viral entities could increase susceptibility to DSS via TLR3 and TLR7 signaling, with a nal increased level of Caudovirales, which is observed in IBD patients [61,69]. Viral dysbiosis in the intestine is consistent with our reports that VDR deletion leads to a higher risk of infection and chronic in ammation [27,31]. Vitamin D metabolites have long been known to support innate antiviral effector mechanisms, including induction of antimicrobial peptides and autophagy. Laboratory data relating to effects of vitamin D on host responses to SARS-CoV-2 speci cally are scarce [1,2].
The diversity and abundance of bacteriophages have been proposed to affect bacterial communities in the gut [70]. We found a signi cant correlation between the virome and bacteria altered by VDR status. The enrichment of bacteriophages and related bacteria were consistent, such as Vibrio phages and their host Vibrio cholerae. Bi dobacterium animalis, one of the bacteria widely used for probiotics in clinical trials [71,72], was negatively correlated with BVDV1 and Vibrio phage JSF6. The impaired protection of B. animalis in the intestinal barrier may further enhance infection with pathogenic viruses, such as BVDV1 and Vibrio phage JSF6. Similarly, Bacteroides acidifaciens, which was found to have strong protective effects against colitis [73], was found to be less abundant and thus negatively correlated with BVDV1.
Meanwhile, Bordetella pseudohinzii showed a positive correlation with Vibrio phage JSF5. As a new member of the genus Bordetella, B. pseudohinzii was rst identi ed and isolated from laboratory-raised mice and has now been detected in mouse facilities worldwide [74]. The infected mice presented elevated numbers of neutrophils in bronchoalveolar lavage uid and in ammatory signs in histopathology, although no obvious clinical symptoms were shown [75]. These causative agents may induce more severe in ammation through cooperation with intestinal microbes. Similarly, the abundances of Vibrio phage JSF5 and Bordetella pseudohinzii [50], were changed in VDR ΔPC and VDR ΔLyz mice compared to control mice.
Microbial metabolites are important players in diverse cellular processes and functions. We found marked changes in virus-related metabolites and pathways, such as fatty acid metabolism. It has been suggested that most viruses require lipids or intermediates of lipid synthesis to replicate, and many of them actively induce lipid metabolic pathways to sustention a favorable replication environment [76].
Pathogenic viruses use lipid droplets as a platform for viral replication, and many non-oncogenic viruses are related to various metabolic alterations during infection, such as glycolysis, nucleotide synthesis, fatty acid biosynthesis and glutaminolysis [65,77]. This nding may explain our results that many longchain fatty acid and fatty acid metabolites were increased in VDR ΔIEC and VDR ΔLyz mice. Although among the most abundant microbes in the gut, phages are also among the least understood [78]. Using gnotobiotic mice, Hsu et al. found that phage predation not only directly impacts susceptible bacteria but also leads to cascading effects on other bacterial species via interbacterial interactions. Moreover, the shifts in the microbiome caused by phage predation have a direct consequence on the intestinal metabolome as revealed by metabolomic pro ling [79]. We also found that the levels some phage infection-related metabolites were altered in the conditional VDR knockout mice compared to the control. For example, VDR deletion decreased fecal serine amino acids, which is consistent with a previous study [79]. Even though various aspects of host central carbon metabolism have been shown to be related to virus infection, several viruses were also found to increase the consumption of key nutrients such as glucose and glutamine and converge on similar metabolic pathways for anabolism [65], which was different from our results. It should be mentioned that all these reports are based on host cell metabolic data, while our analysis is based on intestinal microbial metabolites. Furthermore, it is noteworthy that the precise metabolic changes induced by speci c viruses are often context-dependent and can vary even within the same viral family or largely depend on the individual host [65]. However, the regulatory complexity of viral metabolism in chronic diseases is an area of investigation in the future.

Conclusions
The virome includes viruses that infect host cells, virus-derived elements in the genome, and viruses that infect the broad array of other types of microorganisms that inhabit the host [80]. Dysbiosis is associated with the treatment e ciency of vitamin D supplements [81,82]. Our data demonstrate that conditional VDR knockout causes changes in the abundance and diversity of the virome and functional changes in viral intestinal receptors, which may further induce intestinal dysbiosis and the risk of infection. We also found that the related fecal metabolites were altered in the mice with tissue-speci c deletion, which further con rmed the consequences of dysbiosis. The marked alteration of gut viruses (especially bacteriophages) by VDR status may aid the development of intestinal phage-bacteria or intestinal virusbacteria therapy against pathobionts [56,83]. Our study lls the knowledge gaps of how the virome is affected by VDR in a tissue-speci c manner. The physiological relevance of these changes will be assessed in the future in digestive disease and infectious models. Notably, there is a growing body of evidence suggesting that VDR activation has a regulatory role in mutualistic intestinal virome-host interactions, and more information on the interactions between the sensing of vitamin D and VDR is needed [84]. Availability of data and material:

List Of Abbreviations
data and material will be available by request.

Competing interests:
There are no competing interests. The authors declare no con ict of interest. Funding: We would like to acknowledge support from the NIDDK grants R01 DK105118 and R01DK114126, VA  Relative viral abundance at the species level (family/f_; s_; the unidenti ed species were named with super level and other) is shown with the top 10 species, and less abundant species were grouped as "others". Species were colored using the key as listed on the right side of the gure. Each bar represents an individual mouse (n=10 each group). (B) The virus differential abundance at count per million of the top 10 species in four studied mouse genotypes are illustrated with different colors. Statistical analysis in each group was compared to control VDRLoxP mice. Mean ± SD, n=10 per group; ** P-value < 0.01, * Pvalue < 0.05, Kruskal-Wallis test by ranks.  comparisons are shown with fold-changes (logFCs) and q-values in three comparisons between the conditional VDR knockout groups (VDRΔIEC, VDRΔPC, and VDRΔLyz) and control mice (VDRLoxP). Both fold-change and q-value were colored using the key as listed on the right side of the gure, n=10 each group.

Figure 4
Altered bacterial abundance and correlations with viral alteration. (A) Bacterial species related to the altered bacteriophages were shown with bacterial differential abundance at count per million, including  Altered pattern recognition receptors in the colon of VDR conditional knockout mice. Expression of virus translocation-related PRRs in the intestine, including TLR3, TLR7, NOD1, NOD2, NLRP6, and CLEC4L, was performed with colon epithelial cells using real-time quantitative PCR. Mean ± SD, n=3; ** P-value < 0.01, * P-value < 0.05, one-way ANOVA. acids (D). Mean ± SD, n=5-10; * P-value indicated as the color bar on the right and together with the number inside the blocks, two-way ANOVA.

Figure 7
The interrelations of the host, bacteria, and virus. In the host, microbiota, including both bacteria and viruses, could be affected by immune activities, diet, health status and genetic background (e.g., vdr gene). In turn, the homeostasis of bacteria and viruses is essential to keep the host healthy. Moreover,