Photoperiod Manipulation Affects Transcriptional Profile of Genes Related to Lipid Metabolism and Apoptosis in Zebrafish (Danio rerio) Larvae: Potential Roles of Gut Microbiota

Gut microbiota plays a fundamental role in maintaining host’s health by controlling a wide range of physiological processes. Administration of probiotics and manipulation of photoperiod have been suggested as modulators of microbial composition and are currently undergoing an extensive research in aquaculture as a way to improve health and quality of harvested fish. However, our understanding regarding their effects on physiological processes is still limited. In the present study we investigated whether manipulation of photoperiod and/or probiotic administration was able to alter microbial composition in zebrafish larvae at hatching stage. Our findings show that probiotic does not elicit effects while photoperiod manipulation has a significant impact on microbiota composition. Moreover, we successfully predicted lipid biosynthesis and apoptosis to be modulated by microbial communities undergoing continuous darkness. Interestingly, expression levels of caspase 3 gene (casp3) and lipid-related genes (hnf4a, npc1l1, pparγ, srebf1, agpat4 and fitm2) were found to be significantly overexpressed in dark-exposed larvae, suggesting an increase in the occurrence of apoptotic processes and a lipid metabolism impairment, respectively (p < 0.05). Our results provide the evidence that microbial communities in zebrafish at early life stages are not modulated by a short administration of probiotics and highlight the significant effect that dark photoperiod elicits on zebrafish microbiota and potentially on health.


Introduction
The collection of microorganisms including bacteria, archaea, virus, fungi and protozoa living in different districts of the body as the gastrointestinal tract, skin and mouth give rise to a complex and interconnected ecosystem called microbiota [1]. The vast majority of these microorganisms live in the gastrointestinal tract in a mutually beneficial relationship with the host and their concentration increases moving from the gastric lumen up to the colon/rectum. Gut microbiota plays a key role in maintaining host's health and preventing the insurgence of diseases [2][3][4] and recent studies suggest this role is played via the gut-brain axis [5,6], brain-gut-kidney axis [7], gut-lung axis [8] and others by the production of microbial metabolites that the body would not be able to produce otherwise [9]. Composition of the gut microbiota is usually similar at the phylum level between Danilo Basili and Esmail Lutfi contributed equally to this work.
Electronic supplementary material The online version of this article (https://doi.org/10.1007/s00248-019-01468-7) contains supplementary material, which is available to authorized users. individuals of the same species but diversity and richness of the microbial species may differ as a result of environmental factors, diet, stress, genetics, and other factors [10]. In the last few years, there has been growing evidence that circadian rhythm disorganization, due to environmental cues, has the ability to alter microbiota composition [11,12]. Circadian rhythms are endogenous 24 h rhythmic patterns exhibited by a wide number of organisms, whose main role is to regulate and optimize the functions of cells, organs, systems as well as the animal behaviour [13][14][15]. Thaiss et al. demonstrated that humans and mice intestinal microbiota follows diurnal oscillation in relation to feeding rhythms leading to time-specific compositional and functional profiles across the day [16]. Moreover, humans and mice dysbiosis driven by impaired feeding rhythmicity, as it happens in shift workers and frequent flyers, leads to the disruption of metabolic homeostasis by promoting glucose intolerance and obesity. Recently, Deaver et al. observed Ruminococcus torques, a bacterial species known to negatively affect gut barrier integrity, and Lactobacillus johnsonii, known to help maintaining the intestinal epithelial cell layer, to respectively increase and decrease their abundance in mice undergoing a 4-week period of constant 24 h light [17]. Thus, it is of paramount importance for the microbiota to follow regular diurnal oscillation in order to protect against homeostasis impairment and, consequently, diseases.
The use of probiotic food has been widely studied since the beginning of the twentieth century for its ability to modulate microbiota composition and preserve host's health [18]. According to Merrifield, "a probiotic organism can be regarded as a live, dead or component of a microbial cell, which is administered via the feed or to the rearing water, benefiting the host by improving disease resistance, health status, growth performance, feed utilization, stress response or general vigour, which is achieved at least in part via improving the hosts microbial balance or the microbial balance of the ambient environment" [19]. Probiotic species are able to change the microbial composition of the gut microbiota influencing a wide number of biological processes hence serving as a mean of disease control in aquaculture [20,21]. Among the most widely studied and used microorganisms in probiotic research, we found those belonging to the Bifidobacteria and Lactobacilli genera [22].
In aquaculture, a great effort has focused on identifying ways to ameliorate fish health by both administration of probiotics and modulation of photoperiod. For instance, Ziu et al. fed tilapias with the probiotic Lactobacillus plantarum for 14 days followed by a 3-day suspension and then challenged with Aeromonas hydrophila and found out that suspension of probiotic administration increased susceptibility of the host to A. hydrophila by inducing gut dysbiosis [23]. On the other hand, despite studies investigating photoperiod manipulation have demonstrated its effect on reproductive behaviour and physiology in different fish species [24][25][26], knowledge about microbiota composition following changes in lighting regimes and the underlying biological activity is still poorly known. Moreover, while most of the studies have focused on longterm manipulation of photoperiod, and mainly in adult fish, there is currently a lack of knowledge on whether short-term changes in circadian rhythmicity may affect the physiology of fish at early developmental stages.
The aim of the present study was to establish whether fish microbiota could be modulated by photoperiod manipulation, during the first 24 h since their mouth opening, and whether contemporary administration of beneficial bacteria (probiotics) might be able to affect photoperiod-induced alteration. To this extent, the zebrafish model was used to conduct an in-depth characterization of microbiota composition following disruption of circadian rhythmicity alone or in combination with early administration of the probiotic. By leveraging the power of marker gene-based analysis, we successfully characterized the microbial communities of zebrafish larvae under the different experimental conditions. Moreover, we predicted the functional profile of bacterial communities undergoing different light regimes and validated the findings by qPCR. The results obtained provide new insights about the ability of photoperiod to modulate microbiota composition highlighting the potential effects elicited by exposure to continuous darkness.

Experimental Design
Adult female and male zebrafish (Danio rerio) were purchased from Bologna aquarium (Italy) and acclimatized to the laboratory conditions (27.0 ± 0.5°C under a 12:12-h light/dark cycle). Pairs (seven per tank) were spawned and embryos placed in 1 L tanks under the same laboratory conditions. Temperature was controlled by placing the 1-L tanks within bigger tanks containing the 50 W Magictherm heater (©PRO.D.AC. INTERNATIONAL S.r.l.). At hatching (about 72 h post fertilization, hpf), larvae were divided into 3 groups within the same type of tanks used for the embryos: one group was exposed to a 12:12-h light/dark cycle (LD), one group to 24 h of continuous light (LL) and one group to 24 h of continuous darkness (DD). Illumination was provided by means of fluorescent lamp (36 W, 3350 lm, 4000 K; ©Osram, Germany). In addition, each group received two different treatments: control (C) or probiotic (P), which consisted on the administration of Lactobacillus rhamnosus IMC 501® (C025396A; Synbiotec, Camerino, Italy) via the rearing water at a concentration of 10 6 colonyforming units (CFU) according to previous studies [27,28]. The experiment was set up in triplicates for each condition. After 24 h of exposure, at the same time in the morning, different pools of larvae were euthanized using MS222 (100 mg L −1 ) (Sigma-Aldrich) and stored at − 80°C for high-throughput sequence analysis and gene expression. All procedures involving animals were conducted in accordance with the EU and Italian law on animal experimentation (Directive 2010/63/EU) with no need of requesting ethical approval when larvae within 5 days post fertilization (dpf) are used.

DNA Extraction, PCR and Marker Gene (16S) Amplicon Sequencing
Total DNA was extracted from pools of 100-125 larvae per sample (100 ± 20 mg) using the DNeasy Blood & Tissue Kit (Qiagen) according to manufacturer's instructions. A PCR amplification step was performed to amplify the V4 and V3 variable regions of the 16S rRNA gene using Illumina adapted primer 341F (CCTACGGGNGGCWGCAG) and Illumina adapted barcoded 805R primer (GACTACHVGGGTATCT AATCC) following 16S Metagenomic Sequencing Library Preparation protocol (Illumina, San Diego, CA). Samples and final libraries were quantified and quality tested using the Qubit 2.0 Fluorometer (Invitrogen, Carlsbad, CA). Additionally, libraries were quality checked on Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA). Finally, amplicons were sequenced on the Illumina MiSeq platform run in paired-end mode with 300-bp read length by IGA Technology Services (www.igatechnology.com).

Reads Pre-processing and OTU Assignments
Demultiplexing was performed with CASAVA v. 1.8 and reads not matching indexes or representing the PhiX were removed. Raw sequence reads were processed with the Python package Cutadapt [29] v1.4.2 to remove any residual adapter contamination and quality trimming of paired-end reads was performed using the erne-filter command (Erne v1.4.6, default parameters except min-size = 200) [30]. Reads with a minimum length of 200 bp were retained and analysed with QIIME v1. Briefly, the USEARCH (v8.1.1756, 32-bit) quality filter pipeline was employed to filter chimeric reads, to group replicate sequences, to sort sequences per decreasing abundance, and to finally identify OTUs. OTU picking was achieved applying a minimum pairwise identity threshold of 97%. The most abundant sequence in each OTU was selected to assign a taxonomic classification based on the Greengenes database (v 2013_5) using the RDP classifier (v2.2), clustering the sequences at 97% similarity with a 0.50 confidence threshold. Outliers and singletons were then removed before running downstream analysis.

Microbiota Statistical Analysis
Data analysis was performed within the R statistical environment. Samples were rarefied to the sample with the least reads only for diversity analyses. This choice was driven by the fact that normalizing the samples to account for uneven sampling depth by using rarefaction still represents one of the most promising approaches [31]. The choice of using the sample with the least reads to rarefy the samples was supported by rarefaction curves (Supplementary Material 1). Diversity estimates, rarefaction and principal coordinate analysis (PCoA) were performed using the R package Phyloseq [32]. Statistical differences in alpha diversity were assessed using the ANOVA followed by TukeyHSD post-hoc test with the Benjamini-Hochberg FDR correction. Community composition was analysed using the ADONIS function based on Bray-Curtis distances (R vegan package) [33]. Bar-plots of microbial abundances were drawn firstly taking the average of replicates and then considering taxa whose total abundance across all samples was at least 1%. Differential analysis was performed using raw counts as input into DESeq2 [34,35] and considering a 1% FDR threshold. Functional profiles for 16S rRNA gene sequence data were predicted using the Phylogenetic Investigation of Communities by Reconstruction of Unobserved States (PICRUSt) [36] and analysed using a multi-factorial ANOVA followed by a TukeyHSD post-hoc test within STAMP [37] with a 1% FDR threshold.

qPCR Validation
Functional predictions were validated by means of qPCR. A pool of approximately 40 larvae per sample were homogenized with Precellys Evolution 24 homogenizer coupled to a Cryolis cooler and total RNAs were extracted using TriReagent (Ambion, Alcobendas, Spain), according to the manufacturer's recommendations. RNA was reverse transcribed with the Transcriptor First Strand cDNA synthesis kit (Roche, Sant Cugat del Valles, Spain) and the cDNA obtained was stored at − 20°C. The mRNA transcript levels of key genes drivers of apoptosis (caspases: casp3, casp8 and casp9), circadian rhythm (clocka, clockb and per1a) and lipid biosynthesis and accumulation (hepatocyte nuclear factor 4 alpha, hnf4a; Lanosterol 14a-demethylase, cyp51; Niemann-Pick C1-like 1, npc1l1; fatty acid synthase, fasn; sterol regulatory element-binding transcription factor 1, srebf1; peroxisome proliferator-activated receptor gamma, ppary; 1acylglycerol-3-phosphate O-acyltransferase 4, agpat4; and fat storage-inducing transmembrane protein 2, fitm2) plus two reference genes (beta actin, bactin and acidic ribosomal protein, arp) were examined with a CFX384™ Real-Time System (Bio-Rad, El Prat de Llobregat, Spain). All analyses were performed in triplicate wells using 384-well plates with 2.5 μL itaq SYBR Green Supermix (Bio-Rad, El Prat de Llobregat, Spain), 250 nM forward and reverse specific primers (Supplementary Material 2), and 1 μL diluted cDNA for each sample, in a final volume of 5 μL. The mRNA levels of each target gene analysed were calculated using the Pfaffl method [38], relative to the geometric mean of the two reference genes once demonstrated they were stably expressed by the geNorm algorithm, both implemented in the Bio-Rad CFX Manager 3.1. software. Statistical analysis was initially performed by a two-way ANOVA to assess significance of both photoperiod and probiotic and, for those genes for which probiotic did not elicit any effect, we applied a oneway ANOVA within the control and probiotic conditions alone followed by a Tukey post-hoc test (p value< 0.05). Normality of data was assessed by the Shapiro-Wilk test.

Results
The overarching goal of the present study was to investigate both the ability of probiotic and photoperiod manipulation to modulate the microbiota of zebrafish larvae within the first 24 h after hatching, considered the most sensible window as it represents the time of first opening of the mouth.
High-throughput sequencing following the 16S metagenomic protocol produced 5.75 million paired-end reads 300 bp long, obtaining on average 287,520 reads (min 151,818, max 602,796) for each sample. Average Phred quality score per read was 35. One of the C-LL samples was lost because of the low DNA quality, while data exploratory analysis revealed the presence of two outliers (C-LD and P-LD) that were also removed before running downstream analysis (Supplementary Material 3).

Effects of Both Photoperiod Manipulation and Probiotic Administration on Zebrafish Larvae Microbiota Composition
The ability of either probiotic or photoperiod to affect microbial species richness and evenness was assessed computing alpha diversity using Shannon and Inverse Simpson indexes. Interestingly, a 24-h administration of probiotic occurring straight after hatching was not able to affect any of the indexes while photoperiod affected both metrics (p < 0.01) (Fig. 1).
In order to evaluate relationships among samples and the ability of both probiotic and photoperiod to modulate bacterial composition, a two-dimensional PCoA based on Bray-Curtis [39] distances was performed (Fig. 2). Both LL and DD samples clustered away from LD samples suggesting alteration of lighting regimen may represent a key factor in shaping microbial communities. Moreover, separation between C and P samples within the different photoperiods was not detectable. This finding was also statistically supported by a PERMANOVA analysis, as the photoperiod was the only factor able to affect the beta diversity (p < 0.001).
Overall, the bacterial communities of all treatments were dominated by three main phyla: Proteobacteria, Bacteroidetes and Firmicutes. Moreover, microbial composition of LD samples, regardless of probiotic administration, was made almost entirely of Proteobacteria (98%). The clear separation between LD cycle and DD or LL conditions identified by the PCoA was supported by differences in microbial composition. DD and LL samples were both characterized by a significant reduction of Proteobacteria (Fig. 3). Moreover, the phylum Armatimonadetes, although accounting less than 1% of total abundance, was detected in the LL samples with a 19-fold increase with respect to LD samples, in which was not detected at all (1% FDR).
Most of the OTUs were resolved to either class, order or family level, but for some of them, genus and species taxonomy could be also assigned. At class, order and family level, DD and LL samples showed a more diversified microbial community compared with LD samples, while differences between P and C samples were less pronounced .
A generalized linear model (GLM) was employed to identify genera and species abundance between the different light regimes and as a result of probiotic administration. A considerable number of genera were identified to be differentially abundant between DD-, LL-and LD-treated larvae (Fig. 4a, b). Interestingly, a core microbiota of genera whose abundance was significantly different in both DD and LL samples compared with LD ones was identified.
Genera belonging to this core microbiota whose abundance was increased were Lactobacillus, Segetibacter, Methylobacterium, Mycoplana, Stenotrophomonas, Shinella, Hydrogenophaga, Hypnomicrobium and Ancylobacter, while genera Plesiomonas, Azohydromonas, Limnohabitans and Phascolarctobacterium were found to be less abundant. Moreover, in DD samples, OTUs belonging to the genera Rheinheimera, Pimelobacter, Acidovorax and Brevundimonas were found to be more abundant than in LD samples, while the genera Variovorax, Rikenella, Parabacteroides and Clostridium were significantly less abundant. On the other hand, in LL samples, OTUs belonging to the genera Pirellula, Collinsella, Dorea, Desulforhopalus, Desulfotalea, Desulfosarci na, Paraprevotell a, D esulf obacter, Desulfobulbus, Rikenella, Coprococcus, Fimbriimonas, Tepidimonas, Tepidibacter, Clostridium, Rhodobacter, Haliscomenobacter and Methylotenera were found to increase their abundance compared with LD samples, while abundance of the genera Prostechobacter, Neisseria and Emticicia was significantly reduced. The microbial composition between DD and LL samples was also compared and a total of 22 genera differentially abundant were identified ( Supplementary  Material 7). Surprisingly, when investigating differences in microbial composition between C and P samples, we only identified the genera Candidatus Protochlamydia amoebophila, which could not be accurately identified at the species level, to significantly decrease its abundance in P-treated larvae.
We then set to investigate whether these differences could be observed at the species level. Although we could assign taxonomy at the species level only for a small number of OTUs, we successfully identified species whose abundance was influenced by the different photoperiod exposures regardless of probiotic administration (Table 1).
In DD samples, we identified Methylobacterium organophilum, Pimelobacter simplex, Novosphingobium subterraneum and Acidovorax delafieldii to significantly increase their abundance compared with LD samples while Variovorax paradoxus, Hydrocarboniphaga daqingensis, Rikenella microfusus and Shewanella oneidensis were significantly less abundant. On the contrary, we identified Collinsella aerofaciens, Methylobacterium organophilum, Desulforhopalus singaporensis, Rikenella microfusus, 4 species belonging to the Bacteroides genus (B. massiliensis, B. ovatus, B. uniformis and B. acidifaciens) and 2 species of Fig. 2 Beta diversity. Principal coordinate analysis (PCoA) of 72 hpf zebrafish larvae treated for 24 h with probiotic (P) or not (Control, C) while exposed to different photoperiods: 12:12 h light/dark cycle (LD), continuous light (LL) or continuous darkness (DD). The PCoA analysis shows a clear separation between LD photoperiod and DD and LL photoperiods suggesting the presence of significant changes in the microbial composition. In agreement with the alpha diversity, differences between C and P were not detectable. The PCoA analysis is based on Bray-Curtis distances Fig. 1 Alpha diversity. The plot shows sample species richness according to both the Shannon and the Inverse Simpson indexes of 72 hpf zebrafish larvae treated for 24 h with probiotic (P) or not (Control, C) while exposed to different photoperiods: 12:12 h light/dark cycle (LD), continuous light (LL) or continuous darkness (DD). It is possible to visually detect an increase of alpha diversity measure in LL and DD samples compared to LD samples while discrimination between C and P samples cannot be inferred the Ruminococcus genus (R. bromii and R. gnavus), whose abundance significantly increased in LL samples when compared to larvae exposed to the LD regime. Moreover, L. rhamnosus and Shinella granuli were found to increase their abundance in both DD and LL, while Rikenella microfusus abundance increased in LL samples but decreased in the DD ones when compared to LD samples. Differences in microbial species between DD and LL exposed larvae were also investigated and species belonging to the Bacteroides (B. uniformis, B. massiliensi, B. ovatus and B. acidifaciens) and Ruminococcus (R. gnavus and R. bromii) genera were found to have a significant greater abundance in LL samples (Supplementary Material 8).

Prediction of Bacterial Functional Activity
PICRUSt was employed to predict bacterial functions from the phylogenetic profiles observed. Accuracy of the prediction for each sample was assessed by computing the Weighted Nearest Sequenced Taxon Index (Weighted NSTI) (Supplementary Material 9).
All samples except two had a NSTI value < 0.06 which indicates the good quality of the prediction. The resulting metagenome predictions were categorized in KEGG pathways to identify biological functions potentially modulated by microbial communities of zebrafish larvae exposed to the different photoperiods or administered with probiotic. A total of 23 KEGG pathways were predicted to be modulated by photoperiod-shaped bacterial communities, while none of the KEGG pathways was found to be affected by probiotic administration (Supplementary Material 10).
Interestingly, although significant differences in microbial composition at different taxa levels were identified for both DD and LL samples compared to the LD samples, regardless of probiotic administration, biological activity predictions were only successful for DD samples. DD samples were predicted to be characterized by a greater metabolic potential for lipid biosynthesis, apoptosis and circadian rhythm (Fig. 5).

Experimental Validation of Predicted Biological Activity
Ability of microbial communities arising from different lighting regimes to modulate the predicted biological functions was assessed by investigating expression levels of key target genes (Fig. 6). First, both circadian clock genes which analysed clocka and clockb were found to have a significantly greater expression level in DD samples compared to the other two regimes regardless of probiotic administration, while another clock gene per1a did not show significant differences among any of the conditions. Apoptosis was investigated by quantification of caspase genes (casp3, casp8 and casp9). Casp3 was found to significantly increase its expression in DD samples of larvae not administered with probiotic, while casp8 and casp9 gene expression followed the same trend but was not statistically significant.
Next, ability of DD microbial communities to modulate lipid biosynthetic processes was assessed by quantifying expression levels of key genes involved in lipid metabolism, transport and storage (hnf4a, cyp51, npc1l1, fasn, pparγ, srebf1, agpat4 and fitm2). This goal was achieved by initially performing a two-way ANOVA analysis to investigate effects of both probiotic and photoperiod. Since probiotic was found unable to elicit effects except for the fasn gene, we decided to investigate photoperiod effects within each of the treatments (P and C) by a one-way ANOVA. Expression levels of hnf4a, npc1l1 and srebf1 were found to be significantly increased in DD fish compared to those under the other light regimes both in the presence and absence of probiotic. In addition, cyp51, fasn, pparg, agpat4, ppary and fitm2 were found to significantly increase in DD samples when compared with either LD or LL ones but only in samples not administered with probiotic. Interestingly, fasn was the only gene for which a significant effect elicited by the probiotic could be detected. Surprisingly, in the LL samples none of the genes here examined was found to significantly change its expression level when compared with LD samples regardless of probiotic administration. Fig. 3 Phyla abundances. The stacked barplot shows the phyla abundances of 72 hpf zebrafish larvae treated for 24 h with probiotic (P) or not (Control, C) while exposed to different photoperiods: 12:12 h light/dark cycle (LD), continuous light (LL) or continuous darkness (DD). Replicates were averaged and taxa whose total abundance across all samples was at least 1% were considered

Discussion
In the present study, we investigated the microbial changes occurring in zebrafish larvae undergoing different lighting regimes and/or administered with beneficial bacteria (probiotics) within the first 24 h since their mouth opening (72 hpf), considered the most sensible window for microbial colonization [40,41]. We identified that the administration of the probiotic L. rhamnosus does not significantly affect microbiota composition at the condition tested, while manipulation of photoperiod strongly shapes microbial communities in zebrafish larvae. Moreover, microbial communities undergoing different lighting regimes were successfully predicted to affect the metabolic potential of a wide range of biological pathways. Finally, qPCR analysis demonstrated expression levels of key genes involved in the circadian rhythm, apoptosis and lipid biosynthetic processes to be significantly modulated by the DD regime.
According to Stephens et al., zebrafish larval stage starts at hatching (2-3 dpf) and gut colonization by microbial communities takes place straight after, when the larvae first encounters microbes in the surrounding environment [42]. More specifically, the zebrafish larval mouth opens at approximately 72 hpf and the digestive tract is a continuous tube in connection with the external environment, containing most of the microbiota of the individual [40,41]. We investigated the sensitivity of the zebrafish larvae to probiotics and/or to environmental factors (photoperiod) within the first 24 h since their opening of the mouth. Probiotic administration in 72 hpf zebrafish larvae for 24 h was not able to affect alpha Fig. 4 Genus level analysis. Differential abundance analysis at genus level of 72 hpf zebrafish larvae exposed to different photoperiods: 12:12 h light/dark cycle (LD), continuous light (LL) or continuous darkness (DD). The plot shows genera differentially abundant between a DD and DL samples, and b LL and DL samples, regardless of probiotic administration, along with their fold change. Color coding shows the phyla taxonomic level or beta diversity. Moreover, the only genera found to be significantly less abundant in probiotic-administered samples was the Candidatus Protochlamydia, whose species-level taxonomy could not be accurately identified. Despite the role of Chlamidyae in fish has not been investigated yet, in mammals they are associated with the insurgence of a wide range of pathologies, especially genital infections [43] and pneumonia [44].
The reason why in the current study we have not been able to identify microbial changes in response to probiotic administration could be due to different reasons. Previously, Falcinelli et al. successfully characterized microbial communities changes due to probiotic treatments in 6 dpf larvae after 3 days of treatment [27], whereas in the present study 24 h administration of probiotic to 72 hpf (3 dpf) larvae might not have covered a sufficient window of time for the L. rhamnosus to shape the gut microbial community. However, our inability to capture changes in microbial composition following probiotic administration could be also associated with differences in the gut mucosal composition as previously demonstrated in sea bream (Sparus aurata) by Carnevali and collaborators [45]. In that study, two different probiotic species were administered to sea bream at different developmental stages and found out that the gastrointestinal tract offers variable conditions depending on the stage of development resulting in a diverse microbiota composition. Thus, one possible explanation would be that in the window between 72 and 96 h after fertilization, the zebrafish gut mucosa does not offer favourable conditions for L. rhamnosus to modulate microbial composition. Altogether, these findings suggest that early life stages do not seem to be good targets for acute modulations of gut microbiota by administration of probiotics, to potentially pose benefits to fish health.
Studies investigating the ability of different photoperiod regimes to shape microbiota have already been reported [16,17]. However, to our knowledge, studies focusing on microbial communities' modulation by 24 h dark or light regimes have not been performed. Moreover, we focused on the 24-h post hatching, when microbes first start to colonize the gut. We identified microbiota composition changes at different levels of taxonomy that have been linked to biological activity.
According to Stephens et al. [42], proteobacteria are the most abundant phylum in the microbiota of zebrafish larvae and juveniles, but differences in class composition arise during the different developmental stages. More specifically, γproteobacteria are the most abundant class particularly in zebrafish larvae. Our findings are in agreement with Stephens et al., as Proteobacteria was the most abundant phylum across all the experimental conditions and its abundance was found to significantly decrease following either 24 h of LL or DD conditions. However, Proteobacteria were mainly dominated by β-proteobacteria and, in smaller quantity, by γproteobacteria in larvae undergoing regular LD cycle. This difference is not surprising as the microbiome's structure may be affected by a wide range of factors and there is evidence that microbial composition within individuals belonging to the same species may significantly differ [46,47]. In addition, β-proteobacteria are usually more abundant in freshwater species [48] despite they have been shown to be also important in marine species [49,50].
Significant differences at both genus and species level were observed but the lack of knowledge about the role that each of these taxa play in the maintenance of physiological homeostasis, in human and to a greater extent in fish, limits our understanding to just a small number. Interestingly, a 22-fold increase in the abundance of Ruminococcus species in LL samples compared with LD ones was detected. Despite its role in fish is currently unknown, in humans R. gnavus has been found to cause septic arthritis [51] as well as to correlate with the levels of triglycerides in blood serum [52]. In DD samples, we found species belonging to the Brevundimonas genus to have a 6-fold increase compared to LD samples. This genus has been shown to be one of the most prevalent cause of nosocomial infections [53]. Interestingly, increased abundance of Mycoplana genus was detected in both DD and LL samples in comparison to samples undergoing LD regime. Despite the functional role played by this genus still needs to be investigated, its increased abundance in humans has been associated with multiple sclerosis [54]. These findings suggest that modulation of photoperiod significantly alter microbial communities' composition by creating potentially favourable conditions for species known to be pathogenic in humans.
Alterations of circadian rhythmicity have the ability to affect physiological homeostasis leading to the insurgence of a wide range of diseases [55][56][57][58][59]. In addition, microbial communities play a pivotal role in physiological homeostasis and disruption of microbiome has been associated with the increased occurrence of several diseases as cystic fibrosis, obesity, diabetes, inflammatory bowel disease and chronic obstructive pulmonary disease [60][61][62][63]. In this context, we predicted the metabolic potential of microbial communities arising from the manipulation of normal circadian rhythm and identified lipid biosynthesis, apoptosis and circadian rhythm pathways as the best candidates to be further assessed by qPCR given their association with human diseases. PICRUSt represent a powerful tool to estimate bacterial and archaeal genes present in a microbial community metagenome. However, it does present some disadvantages that include the following: (1) its ability to only predict the portion of the full metagenome targeted by the primers used since the input data is 16S rRNA and eukaryotic or viral contributions to the metagenome cannot be taken into account, (2) its ability to only predict gene families already known and included in the orthology reference used (KEGG in this case) and (3) the fact it is based on evolutionary modelling the gene content of known reference genome; hence, the accuracy of the prediction will depend on the availability of the appropriate references. Despite these limitations, PICRUSt has been widely and successfully used to predict the functional capabilities of a microbial community and here we took advantage of this tool for the same purpose.
Circadian rhythms are endogenous oscillations of about 24 h that regulate organism's physiology and behaviour [64]. They can be entrained by environmental cues called zeitgebers (i.e. light and temperature cycles and periodic food availability) [65]. At a molecular level, the circadian clock mechanism is conserved in vertebrates [15]. Briefly, clock and bmal1 genes form the CLOCK:BMAL1 heterodimer in the cytoplasm that, after translocation into the nucleus, trigger the transcription of target genes as per and cry [66]. Then, a negative feedback is achieved by PER and CRY by forming a Fig. 5 PICRUSt predictions. The plots show the post-hoc results of the ANOVA analysis along with the 95% confidence intervals of 72 hpf zebrafish larvae exposed to different photoperiods: 12:12 h light/dark cycle (LD), continuous light (LL) or continuous darkness (DD). Blue, green and orange bar refer to DD, LD and LL photoperiod, respectively heterodimer complex that, once translocated back into the nucleus, inhibits their own transcription by blocking the CLOCK:BMAL1 complex. Our results suggest an alteration of the circadian rhythmicity in animals undergoing DD photoperiod according to the significant increase induced in the expression level of both clocka and clockb, while the expression of per1a, which plays a key role in the generation of the circadian rhythm, was not found to be affected. LD cycles are needed for the correct onset of behavioural rhythmicity in zebrafish larvae [67]. Moreover, Villamizar and collaborators [68] demonstrated that alteration of those normal LD cycles in this species affects survival and growth, while zebrafish larvae kept under DD conditions died before 18 days post hatching [68], highlighting the detrimental effect of constant darkness. Dekens and Whitmore investigated the expression level of zebrafish clock genes under different light regimes and found out that the expression of per1 loses its rhythmicity in DD conditions, while clock1 (clocka) expression in both LD and Fig. 6 qPCR validation. The boxplots show the expression values of the circadian rhythm, apoptosis and lipid metabolism-related genes analysed in samples of 72 hpf zebrafish larvae treated for 24 h with probiotic (P) while exposed to different photoperiods: 12:12 h light/dark cycle (LD), continuous light (LL) or continuous darkness (DD). Significant differences among fish under different light regimes are indicated by uppercase and lowercase letters for probiotic and control, respectively DD was similar from 1 to 4 dpf [69]. Contrarily to that study, our findings provide evidence that a 24-h exposure to constant darkness in zebrafish larvae at 3 dpf following normal LD cycles is able to significantly alter the expression of clock genes, and to our knowledge this is the first study investigating the effects of such a short-term alteration of normal LD cycles in this species.
Recently, the progress in circadian rhythms research shed light on the circadian regulation of lipid metabolism in mammals [70] and in fish [71]. Clock-controlled genes involved in lipid metabolism, as well as other metabolic pathways, are rhythmically enhanced or repressed by the molecular circadian clock and the loss of clock function is associated with the insurgence of abnormal metabolic phenotypes [72].
Our results revealed hnf4a, cyp51, npc1l, srebf1, fasn, pparγ, agpat4 and fit2, genes involved in lipid metabolism [73][74][75][76][77][78], transport [79] and storage [80,81], to significantly increase their expression in DD samples in comparison to samples undergoing DL and LL photoperiods, suggesting a disruption on lipid turnover and metabolism that could potentially lead to a total body increase of cholesterol and triglycerides. As the role of the gut microbiota in lipid metabolism is now well documented [82,83], these findings suggest that microbial communities associated with a dark regimen could be able to disrupt normal biosynthetic processes affecting organism's health since lipid homeostasis impairment has been associated with the insurgence of a wide range of pathologies [84,85]. Our findings are in agreement with the study of Xie and collaborators [86] and the one of Casado and collaborators [87] where lipid metabolism impairment was observed in rats undergoing different lighting regimes.
Apoptosis, a programmed cell death that involves the genetically determined elimination of cells, is a natural homeostatic mechanism normally occurring during development, aging and different pathologies [88]. Our results show the expression levels of casp3 to significantly increase in DD samples compared to the other regimes suggesting that microbial communities may increase the occurrence of apoptosis. This finding is in agreement with Carballada et al., which identified the presence of apoptotic cells in the epithelium of the epididymis, seminal vesicles, prostate and coagulating gland of the golden hamster (Mesocricetus auratus) following a short-day light regimen (8:16 LD cycle) [89]. A similar study was also performed by Moffatt-Blue and collaborators [90], who observed that a short-day light regimen increased the occurrence of apoptotic follicles in the ovary of Siberian hamsters (Phodopus sungorus). In both studies, the presence of active casp3 was identified by immunodetection. Despite the fact we only analysed gene expression, the changes observed support similar effects in fish. In fact, Dezawa et al. reported increased apoptosis (determined by TUNEL staining) in retinal ganglion cells of carp following exposure to darkness [91].
Overall, the circadian rhythm, lipid metabolism and apoptotic pathways were successfully predicted to be modulated by the microbiota and changes in the expression of key genes driver of the aforementioned biological processes were identified. However, these findings taken together are not conclusive to demonstrate that these observed biological outcomes are truly mediated by the microbiota. Indeed, these outcomes could be the result of a direct effect of circadian disruption secondary to the observed changes in microbiota composition. In this context, the results obtained suggest that the alteration of lipid metabolism and the signal of apoptosis might work through microbial manipulation induced by different light regimes.

Future Directions
Marker gene-based analysis represents a powerful tool to characterize bacterial communities in a given experimental condition. As 16S is shared by all the microorganisms, it is possible to target a wide variety of bacteria. Moreover, given the presence of conserved regions in its sequence, it is possible to easily design primers targeting these regions. Also, since the 16S is one of the most studied and characterized genes, phylogenetic trees are well developed and taxonomic information are easily accessible through public databases. However, this approach just gives indication based on the presence/absence of a given taxa and lacks the ability to put the findings into a functional context. Recently, a few tools with the ability to predict biological activity triggered by microbial communities have been developed (PICRUSt, PAPRICA, tax4fun). Although marker gene-based analysis in combination with these tools gives the opportunity to partially explore the functional capability of a given microbiota, it still lacks the ability to provide a clear understanding of the association between a specific taxon (mainly genera and species) and the observed biological activity. More specifically, mechanistic information on how genera and species modulate biological functions cannot be investigated. In this context, metatranscriptomic approaches represent a more comprehensive methodology to investigate microbial communities and how they modulate an organism's physiology and homeostasis. However, the elevated cost associated with this approach along with the complexity of the data analysis still represents a barrier to the advancement of this field. Our findings represent a very important first step in the characterization of microbial communities arising from different lighting regimes as they provide an overview of specific genera and species in the different conditions. However, we envisage that metatranscriptomic, in addition to studies employing germ-free animals, has the potential to further improve our understanding on the effects of photoperiod manipulation on the organism's health providing the link of the association between microbial communities and organism's physiology.
Data Accessibility Raw sequencing data was deposited as FASTQ files in NCBI Sequence Read Archive (SRA) database under the Bioproject number PRJNA528701.
Authors' Contribution OC conceived and designed the experiment. EL and SF carried out the experiment. DB analysed the data. SB performed the qPCR validation. DB, OC, EL, EC, IN, SF and CB wrote the paper.
Funding Information This research was supported by funds from the MINECO grants (AGL2014-57974-R and AGL2017-89436-R) to both EC and IN, from grant (FAR2018 and FIR2018) to CB and from grant FA2015 to OC.

Compliance with Ethical Standards
Competing Interests The authors declare no competing interests.