Circadian dynamics of the teleost skin immune-microbiome interface
Microbiome volume 9, Article number: 222 (2021)
Circadian rhythms of host immune activity and their microbiomes are likely pivotal to health and disease resistance. The integration of chronotherapeutic approaches to disease mitigation in managed animals, however, is yet to be realised. In aquaculture, light manipulation is commonly used to enhance growth and control reproduction but may have unknown negative consequences for animal health. Infectious diseases are a major barrier to sustainable aquaculture and understanding the circadian dynamics of fish immunity and crosstalk with the microbiome is urgently needed.
Here, using rainbow trout (Oncorhynchus mykiss) as a model, we combine 16S rRNA metabarcoding, metagenomic sequencing and direct mRNA quantification methods to simultaneously characterise the circadian dynamics of skin clock and immune gene expression, and daily changes of skin microbiota. We demonstrate daily rhythms in fish skin immune expression and microbiomes, which are modulated by photoperiod and parasitic lice infection. We identify putative associations of host clock and immune gene profiles with microbial composition. Our results suggest circadian perturbation, that shifts the magnitude and timing of immune and microbiota activity, is detrimental to fish health.
The substantial circadian dynamics and fish host expression-microbiome relationships we find represent a valuable foundation for investigating the utility of chronotherapies in aquaculture, and more broadly contributes to our understanding of the role of microbiomes in circadian health of vertebrates.
Circadian rhythms—endogenous daily cycles in physiological and behavioural processes—are a ubiquitous phenomenon to life. Living organisms are adapted to anticipate the daily variations in light, temperature, or food availability driven by the relentless 24-h rotation of Earth. Circadian rhythms are orchestrated by ‘clock genes’ driving transcriptional-translational autoregulatory feedback loops , which are transduced to temporally coordinate biological activities. Immune functions are energetically costly  and often highly rhythmic, enabling organisms to mount their most efficient response at times when risk of infection or injury is highest [3,4,5]. Conversely, immune factors and infections can affect expression of molecular clocks [6,7,8] and subsequent rhythmic phenotypes [9, 10]. Disruption of normal circadian cycles can impact immune functioning [11, 12] and may increase disease risks .
A primary function of immune systems is to protect the host from invading pathogenic microbes. However, animals are invariably colonised by a suite of microorganisms—their ‘microbiome’—which span the spectrum of symbiosis from mutualists to opportunistic pathogens. In vertebrates, it is increasingly apparent that immune systems and microbiomes are intricately linked, together mediating homeostasis and influencing disease outcomes [14, 15]. Intriguingly, microbiomes may also be rhythmic, exhibiting diurnal fluctuations in community composition and activity . In studies of the mammalian gut, it has been demonstrated that not only does expression of host clock genes shape microbiome rhythms  but also disruption of microbial rhythms in turn impacts host circadian functioning .
Aquaculture is the world’s fastest growing food sector, but infectious disease is the principle barrier to sustainability  and a multi-billion-dollar problem for the global industry . Whilst understanding of fish microbiomes is still in its infancy compared to mammalian systems, there is rapidly growing interest in their role for fish nutrition, health and disease resistance [21,22,23,24,25]. Photoperiod manipulation is commonly used in fish farms, with extended day lengths, and, in the extreme, constant light, to promote increased growth rates, or control maturation and reproduction [26,27,28]. Fish are thought to have a decentralised clock, with cells from multiple tissues expressing circadian genes [29, 30], self-sustained rhythmicity and light responsiveness (see  for review). In common with higher vertebrates, fish appear to exhibit circadian rhythmicity in certain immune factors [29, 30, 32,33,34,35]. Therefore, extreme lighting regimes may have profound implications for fish health and response to infection. Moreover, there are indications that infection and/or stress may impact expression of fish circadian clocks [36, 37]. Currently, the extent to which light manipulation practices contribute to disease in aquaculture is unknown. More fundamentally, the daily dynamics of the fish immune-microbiome interface is yet to be explored. Uncovering the effects of infection and photoperiod on fish immune and microbiome rhythms will be pivotal for both aquaculture disease mitigation strategies, and a broader understanding of the role of holobiont chronobiological interactions for animal health.
Here, using rainbow trout (Oncorhynchus mykiss) as a model, we combine 16S rRNA gene metabarcoding and direct mRNA quantification methods to simultaneously characterise the circadian dynamics of skin clock and immune gene expression, and daily changes of skin microbiota. We compare circadian rhythms of host clock and immune gene expression and microbial community composition in healthy fish under regular light-dark cycles (12 h light, 12 h dark, hereafter 12:12 LD) with those in fish experimentally infected with the ectoparasite crustacean Argulus foliaceus and/or raised under constant light (24:0 LD). In addition, we assess rhythmicity in the functional potential of rainbow trout skin microbiomes and establish host expression-microbiome association networks.
Photoperiod impacts host responses to infection
Photoperiod (12:12 LD vs 24:0 LD) had no significant impact on growth of juvenile rainbow trout over the 16-week trial period (Supplementary Figures 1a and 1b). However, a significantly higher number of Argulus lice survived 7 days post-inoculation on fish maintained in constant light conditions (t115 = −8.418, P = 1.23 × 10−23, Supplementary Figure 1c). To examine overall immune responses to Argulus infection, we grouped fish from all timepoints, and contrasted expression of 27 genes from innate and adaptive immune pathways between treatment groups (12:12 LD control, 12:12 LD infected, 24:0 LD control, 24:0 LD infected). Infected rainbow trout had significantly higher expression of 24 immune genes (89%) under 12:12 LD, whereas only 14 (52%) were significantly higher in infected fish compared to healthy controls under constant light (Fig. 1). Two genes (c3 and tgfb) were significantly reduced by infection in both light conditions (Fig. 1). Expression levels were broadly similar amongst infected groups, although upregulation of the pro-inflammatory interleukins il4 and il6 was lower under constant light (Fig. 1). Conversely, comparisons of healthy (unchallenged) fish under 12:12 LD and 24:0 LD revealed a substantial difference in immune expression profiles, with unchallenged fish under constant light exhibiting elevated expression levels in 21 genes (78%), more similar to both infected groups in most immune genes (Fig. 1).
Circadian rhythmicity of host expression is altered by infection and photoperiod
Under 12:12 LD, core and accessory vertebrate clock genes exhibited significant circadian rhythmicity in healthy rainbow trout skin (Fig. 2, Supplementary Table 1, Supplementary Figure 2). Many of these genes are also found to be expressed rhythmically in fish raised in constant light (Fig. 2, Supplementary Table 1, Supplementary Figure 2) and when fish are placed into ‘free-running’ (constant dark, DD) conditions (Supplementary Figure 3, Supplementary Table 1). However, overall expression levels of clock genes are elevated in the absence of light cues (Fig. 2, Supplementary Figure 2), except for timeless (suppressed expression in 24:0 LD). In addition, bmal2, clock1b, per1, and rora exhibited a significantly different phase of expression in constant light (Supplementary Table 1, Fig. 2, Supplementary Figure 2).
Argulus lice infections had variable impacts on the expression levels and rhythmicity of the clock genes. When contrasted with healthy control groups, some gene rhythms were dampened in infected fish (i.e., significantly reduced amplitude; 12:12 LD clock3, 24:0 LD per1), rendered arrhythmic (cry2 in 24:0 LD) and/or phase-shifted (bmal1 in both light treatments, cry1 and per1 in 12:12 LD, clock3 in 24:0 LD). Rhythms of clock gene expression in infected fish under the two photoperiod treatments did not differ in amplitude. But, bmal1, clock1b, clock3, cry1, per1, per2, rory, and timeless had significantly different phases of expression between infected fish under 12:12 LD and those raised in constant light. In addition, bmal2, clock3, csnk1d, per2 and reverbb had increased rhythm mesors (average expression) in 24:0 LD, whilst timeless was suppressed (Supplementary Table 1, Fig. 2, Supplementary Figure 2).
Significant rhythmicity in expression was found in both innate and adaptive immune markers (Supplementary Table 1, Supplementary Figures 4 and 5), with a substantial proportion remaining rhythmic under free-running (DD) conditions (Supplementary Figure 6). The cathelicidins (cath1, cath2), igd, il17a, and tbx21, whilst rhythmic in healthy fish under 12:12 LD, were arrhythmic in fish maintained in constant light (Supplementary Table 1, Supplementary Figures 4 and 5). Of the immune genes rhythmic in healthy fish under both light conditions, the innate markers chi, hamp and nos2, and the adaptive markers cd4, cd8a, foxp3b, igm, igt, tcrb and tgfb had significantly different mesors; with the exception of nos2, all were more highly expressed in 24:0 LD. However, some of these more highly expressed genes (cd4, foxp3b, hamp, igt, tgfb), and others with similar expression levels between photoperiods (il4, tlr9), were phase-shifted in constant light (Supplementary Table 1, Supplementary Figures 4 and 5).
Fewer immune genes were rhythmically expressed in infected fish: 76% and 67% of rhythmic genes found in healthy fish were also rhythmic in the 12:12 LD and 24:0 LD infected groups respectively. Under 12:12 LD, the vast majority (94%) of the immune genes assayed with rhythmicity in both healthy and infected fish exhibited higher mesors in the infected group. In contrast, only 57% of immune genes with rhythms in healthy and infected fish in 24:0 LD had different expression levels (Supplementary Table 1). Only tbx21 had a significantly altered amplitude in rhythm, with a higher amplitude in infected fish at 12:12 LD compared to both healthy 12:12 LD fish and infected fish in constant light. Argulus infection also shifted the phase of expression of mhcii under 12:12 LD and c3, nos2 and igt in 24:0 LD (Supplementary Table 1).
Argulus infection impacts skin mucus microbiome communities
After read pre-processing, error correction, chimaera removal and filtering, a total of 1037 amplified sequence variants (ASVs) were found across all samples. Rarefaction curves confirmed a minimum read depth of 2000 was sufficient to reach saturation of diversity in rainbow trout skin (Supplementary Figure 7a). Background water samples were distinct from fish groups (Supplementary Figure 7b) and had a significantly higher alpha diversity (average Faith’s phylogenetic distance: Control 12:12 LD = 7.10, Infected 12:12 LD = 7.63, Control 24:0 LD = 5.59, Infected 24:0 LD = 7.63, Water = 20.67, Supplementary Table 2, Supplementary Figure 7c). Contrasts of alpha diversity amongst fish samples revealed that the microbiomes of healthy fish under constant light were significantly less diverse than all other groups (all pairwise Kruskal-Wallis tests P < 0.001, Supplementary Table 2). Multivariate permutational analysis of beta diversity indicated significant compositional differences amongst all groups (Supplementary Figure 7b, Supplementary Table 3). Metagenomic sequencing of a subset of samples resulted in six high quality (completion > 90%, contamination < 5%) and five moderate quality (> 60% completion, < 5% contamination) metagenomic assembled genomes (MAGs, Table 1).
The skin microbiome communities—determined by 16S rRNA profiling—in all groups were dominated by Proteobacteria, with Pseudomonadaceae and Burkholderiaceae accounting for over 50% of the communities in all groups and timepoints (Fig. 3). Taxonomic assignment of metagenomic reads revealed similar profiles, including dominance of the genus Pseudomonas (Supplementary Figure 8). The most abundant MAG was determined to be a Pseudomonas species from the P. fluorescens lineage , most closely related to P. chloraphis (Supplementary Figures 9 and 10). Wilcoxon rank-sum testing and DESeq2 on 16S data both revealed substantial differences in the relative abundances of microbial taxa between healthy and lice-infected fish (Fig. 4). At the higher taxonomic levels, healthy fish under both light treatments had a greater proportion of Actinobacteria and Firmicutes lineages, whilst both infected fish groups had increased Bacterodia lineages (Fig. 4A). At the genus level, many Gammaproteobacteria were more abundant in both infected groups (e.g. Aeromonas, Perlucidibaca, Undibacterium, Fig. 4B). Bacteroidia genera, including several Chryseobacterium, Flectobacillus and Flavobacterium ASVs, were also increased in infected fish, with Flavobacterium accounting for some of the highest fold-changes in abundance (Fig. 4B). One MAG was phylotyped as Flavobacterium sp. (most closely related to the fish pathogen F. columnare, Supplementary Figure 11), which was significantly more abundant in infected fish (Supplementary Figure 9, Supplementary Table 4). Full lists of 16S rRNA differentially abundant taxa are provided in Supplementary Table 5.
16S rRNA functional prediction (using inferences based on reference genomes similarity ) of microbiomes revealed putative differences in the activity of microbial communities amongst healthy and infected fish. LefSe analyses indicated pathways enriched in healthy fish groups were predominantly degradative classes including amino acid, aromatic compound and carbohydrate degradation (Supplementary Table 6). In contrast, functional enrichment of lice-infected fish microbiomes was dominated by biosynthetic pathways in both light conditions, particularly those involved in cofactor, carrier and vitamin biosynthesis (Supplementary Table 6). Similarly, functional analyses of metagenomic data revealed substantial differences between healthy and infected fish (Supplementary Table 7). Like the 16S rRNA results, we found increased abundance of functions related to lipid, nucleotide and thiamine processing in infected fish, and increased degradation-related (e.g. urea cycle) and energy generation (e.g. fermentation) functions in healthy fish. Overall, in 16S rRNA functional predictions, a greater number of pathways were identified as differentially abundant between healthy and infected fish in 24:0 LD, suggestive of a greater disruption in microbiota functional potential due to parasitic infection in fish maintained under constant light.
Circadian rhythmicity of skin microbiota and association with host gene expression
Circadian rhythmicity in 16S rRNA relative abundance was apparent in 49 skin bacteria genera in one or more of the treatment groups (Supplementary Table 8, Fig. 5). We also found significant differences in abundance of several MAGs between the two circadian timepoints used for full metagenomic analyses (Supplementary Figure 9, Supplementary Table 4), supporting profound daily shifts in community composition.
Of the 41 genera rhythmic in both healthy and infected fish at 12:12 LD, 17 (41.5%) had significantly different mesors. In contrast, 23 (60.5%) had significantly different mesors when comparing healthy and infected fish under constant light. Perlucidibaca, Undibacterium, and Rhodoferax had significantly greater rhythm amplitudes in infected fish under both light treatments. In addition, Flectobacillus, Alkanibacter and an unassigned Burkholderiaceae genus had higher rhythm amplitudes in infected 12:12 LD fish, whilst Duganella had higher amplitude in 24:0 LD infected fish only. Under 12:12 LD, lice infection significantly altered rhythm phases of seven bacteria genera (unknown Rhizobiaceae, unknown Rickettsiales, Deefgea, Massilia, unknown Neisseriaceae, unknown Chitinophagales and Legionella). Pseudoclavibacter was the only genus found to have altered rhythm phase in 24:0 LD healthy vs infected comparisons.
Visualisation of the timings of peak abundances of rhythmic taxa indicated no clear phylogenetic patterns (e.g. rhythmic Proteobacteria genera peak abundances were spread across the circadian cycle, Fig. 5A). However, when considering the rhythms of the functional potential of the microbiome communities, we found evidence of temporal patterns (Fig. 6). In healthy fish under 12:12 LD, the majority of rhythmic biosynthetic (e.g. heme b, L-lysine and isoprene biosynthesis) and energy generation (e.g. glycolysis, TCA cycle) functions peaked in the first hours of light, whilst degradation function peaks were found primarily in dark hours. In contrast, in infected fish under 12:12 LD, rhythmic biosynthetic and energy generation functions predominantly peak in abundance towards the end of the dark period, whilst degradation pathways peaked just before dark. Supporting these findings, analysis of metagenomic sequences revealed substantial differences in abundance of bacterial functions between the onset of light and dark periods in both healthy and infected fish (Supplementary Table 7).
Constant light conditions also appeared to shift the broad temporal patterns of function abundances. In healthy fish under 24:0 LD, many biosynthetic pathways (e.g. l-valine, heme b and enterobactin biosynthesis) peaked at a circadian time similar to those in the 12:12 LD group. However, we also found a large cluster of biosynthetic pathways peaking later in the day (e.g. fatty acid, spirillozanthin and coenzyme M biosynthesis). In infected fish under 24:0 LD, biosynthetic pathway rhythms were more dispersed, with peaks spread around the majority of the 24-h cycle. For degradation and generation of energy pathways in both healthy and infected fish under 24:0 LD, we found multiple clusters of peak abundances around the 24 h cycle, rather than a single predominant cluster as in 12:12 LD conditions (Fig. 6).
We used co-occurrence network analyses to assess associations of host gene expression and their microbiomes, using betweeness centrality scores and number of connections (degrees) to identify influential genes and bacteria genera [40, 41]. In healthy 12:12 LD fish, there was a high level of connectivity within host immune and clock genes, and within microbial taxa (Fig. 7). Links across the gene expression and bacteria subnetworks were primarily via the rhythmically expressed clock genes clock1b, clock3, bmal1, rora and csnk1d. However, expression of the toll-like receptors tlr2 and tlr9 were significantly associated with abundance of Bacillus and Enhydrobacter, respectively. In contrast, networks of infected fish under 12:12 LD revealed a higher level of connectivity between host expression and bacteria (Fig. 7). The immune markers cd4 and tcrb, and the clock gene reverbb were most influential in terms of their betweeness centrality scores and number of significantly associated microbial taxa (Fig. 7).
In contrast to 12:12 LD, clock genes were less influential (in terms of centrality) in gene-microbe networks for uninfected fish under constant light (Supplementary Figure 12). However, several immune genes (igd, ifng, nos2, hamp, tcrb, foxp3b) were significantly associated with one or more bacteria genera. Tcrb was most influential by betweeness centrality (expression positively correlated with Janthinobacterium and negatively with Flavobacterium), whilst ifng was linked to the highest number of taxa (Escherichia-Shigella, Pseudomonas, Varioivorax, Stenotrophomonas and Pseudoclavibacter). Similar to 12:12 LD contrasts, the network of infected fish under 24:0 LD showed a higher level of connectivity between host gene expression and microbiota compared to the healthy network (Supplementary Figure 12), with the immune markers cd8a and tcrb found to be the most influential genes (in terms of number of associations with taxa and centrality score).
We demonstrate the daily dynamics of immune expression and microbiome composition in fish skin and show ectoparasite infection and constant light—a commonly used environmental condition in aquaculture—can significantly alter circadian rhythms of immunity and microbiota, which may be detrimental for host disease resistance. In addition, network analyses of host gene expression and their microbiomes reveal clock expression, and T cell populations are likely key in shaping the skin host-microbiome interface of teleosts. Our examination of the skin circadian immune response to infection under extreme photic regimes are directly relevant to fish culture practices; fish peripheral tissues are thought to have entrainable, light-responsive clocks , which may make them particularly susceptible to negative health consequences from constant lighting as used in aquaculture.
Over our trial period, we found no significant difference in the growth of rainbow trout fry maintained under 12:12 LD and constant light (24:0 LD) when fish were provided equivalent food rations. However, when challenged with Argulus lice, their ability to clear infection was significantly altered by photoperiod. Under constant light, rainbow trout had a significantly higher lice burden 1 week after inoculation, indicating a reduced ability to mount an effective immune response. These findings are consistent with previous studies showing extended day length increases ectoparasite susceptibility and altered expression in specific immune genes in sticklebacks . Immune profiles in uninfected fish showed elevated levels of expression in both innate and adaptive pathways under constant light. When infected with lice, rainbow trout under both photoperiods showed similar patterns of immune gene responses, except for the interleukins il4 (mediator of Th2 differentiation) and il6 (key to initiate inflammation) which were expressed at lower levels in constant light. Early inflammatory responses and subsequent initiation of Th2 processes are thought to be critical to resistance of crustacean ectoparasites in salmonids . Taken together, chronic elevation of the immune gene expression—which may result in immune exhaustion  or other immunopathologies —and reduced ability to mount effective responses key to lice resistance suggest rearing of fish in the absence of light cues are likely to be detrimental for health.
The impact of photoperiod on overall magnitude of immune gene activation is not the only factor important to parasite resistance; the rhythmicity and the appropriate timing of immune activity (i.e. when fish are maximally vulnerable to pathogen attack) may also be key to pathogen defences. Under regular light-dark cycles, we show rainbow trout skin is highly rhythmic in expression of the core vertebrate clock genes and many immune genes in both innate and adaptive pathways. In essence, we find the highest expression of pro-inflammatory markers (e.g. il6, il17a) at the onset of the light period and peaks in anti-microbial peptide genes (e.g. cathelecidins) mid-light phase, whilst immunoglobulin and T cell markers were highest during dark hours. The timing of different facets of immune systems are considered to have evolved to offer hosts greatest protection from invading pathogens when most likely to encounter them, whilst avoiding energetically inefficient and potentially immunopathological risk of continual immune activation . We found that constant light resulted in arrhythmic expression of genes involved in mucosa anti-microbial (e.g. cathelecidins, igd, il17a) and Th1 (tbx21) responses. Furthermore, genes with phase-shifted expression rhythms in constant light were dominated by those involved in T cell differentiation and regulation (e.g. cd4, foxp3b, il4, tgfb). Loss of synchrony between host immunity and parasite activity and/or immune evasion rhythms is very likely to be detrimental for host fitness and survival . Our results indicate that this is a factor in the reduced clearance of lice in fish reared in constant light. Clearly, the impacts of light cycle perturbation, be it intentional such as in aquaculture or unintentionally due to light pollution , must be more carefully considered for animal health.
The primary function of fish skin mucus is as a protective barrier and it hosts diverse communities of microbes  which are thought to contribute to protection from microbial pathogens via competitive and/or antagonistic activities [50, 51]. Whilst pathogenic taxa occur mostly at low levels in healthy teleost microbiomes, their proliferation is a common signal of microbiome perturbation and dysbiosis . Argulus lice infestations are commonly observed alongside bacterial, fungal or viral infections . Here, using evidence from both 16S rRNA profiling and metagenomic sequencing, we demonstrate significant reorganisation of bacterial communities and their potential functional activities in rainbow trout skin when infected with A. foliaceus, including notable increases in abundance of genera associated with infectious disease [54, 55]. Flavobacterium spp. were found to have marked increases in infected fish, and the Flavobacterium MAG (closely related to F. columnare, the agent of columnaris disease) was more abundant in Argulus-infected samples. Fish lice may elicit host immune profiles and/or destabilise skin microbiota communities resulting in reduced ‘colonization resistance’ , or be direct vectors [56, 57]. Further research into the microbiota of Argulus and other fish ectoparasites, and their pathogen vectoring capabilities, will be valuable for understanding their role in coinfection dynamics. Intriguingly, rainbow trout raised in constant light had a significantly lower microbiome diversity and, when challenged with Argulus, exhibited greater shifts in taxonomic composition, and possibly functional potential, compared to fish under regular light-dark regimes. Given the growing body of evidence for the importance of ‘healthy’ microbial communities  for effective host homeostasis and disease resistance [59, 60], characterising circadian disruption to microbiomes is important for understanding animal disease risks.
We demonstrate significant daily dynamicity in the skin microbiome of rainbow trout; a substantial proportion of bacteria genera exhibit rhythmic changes in relative abundance, suggesting a temporal structure to microbiome functional activity. Parasitic infection appears to perturb microbiome composition, and our data predict shifts in the timings of peak biosynthetic, degradative and energy generation pathway activity in the microbial community. Understanding of the functional importance to the host of commensal microbiota in teleost skin is still in its infancy , and 16S rRNA predictive functional analyses are only indicative of actual microbial activity . However, we also find strongly contrasting taxonomic and functional profiles between the two timepoints used for full-scale metagenomic analysis, which further supports daily cycles in microbiome composition and functions. Temporal metatranscriptomic and/or metaproteomic profiling will be an important means to build upon our results and to confirm functional cycles and their potential significance to host health. Moreover, the impact of circadian disruption on other aspects of fish biology, such as hormone levels , diet/feeding patterns and behaviour  and their subsequent influence on microbiomes, requires further research. Nevertheless, as interest builds towards the utility of microbiome engineering strategies to promote health and productivity in aquaculture [23, 52, 64], we propose that a chronobiological understanding of fish microbiomes may be crucial for their effectiveness. The daily rhythms of both fish host immunity and their microbiome communities, for example, could be critical to uptake and establishment of probiotics treatments. Chronotherapeutics—the timed application of treatments and vaccines —in human medicine holds great promise for improving efficacies but is yet to be given full consideration for managed animal health.
In the mammalian gut—by far the most studied host-microbiome interface—there is a complex interplay between immune factors that shape microbial communities, and, conversely, microbiota profoundly affecting immune system development and maintenance [14, 15]. Mammal gut microbiome daily rhythms may themselves play a role in host circadian health [66, 67]. However, in other tissues, and particularly for non-mammalian vertebrates, host immune-microbiome connectivity and circadian dynamics remain poorly understood. For teleosts, there is evidence that macrophages  and adaptive immune components (e.g. T cells  and immunoglobulins ) may be key to mucosal microbiome composition. Our study is the first to present an integrated analysis of skin microbiomes with a broad set of immune and circadian clock gene expression profiles in fish. We found genes of the core secondary feedback loops (e.g. bmal, clock, rora, csnk1d) that define the vertebrate molecular clock to be strongly associated with microbial taxa relative abundances in uninfected rainbow trout under 12:12 LD, yet these direct clock-microbe associations were largely absent in constant light. Similarly, mice faecal microbiota composition appears closely linked to bmal1, with knock-outs resulting in arrhythmicity and altered abundance of microbial taxa . Our results suggest this arm of the biological clock may be pivotal to orchestrating changes in mucosal microbiomes across vertebrates. However, we also find perturbation of microbial communities via ectoparasite infection reconfigures the connectivity of host expression and microbiota. In both photoperiods, lice-infected fish immune-microbe networks show a greater level of connectivity between host immune gene expression and microbial taxa compared to uninfected individuals. In particular, our results indicate T cell markers to be central to this host-microbiome interface during ectoparasite infection. Under 12:12 LD, we find the T helper cell gene cd4 to be strongly linked to microbiome composition, whilst in constant light the cytotoxic T cell marker cd8a appears to be more influential to microbiome-immune associations. For teleost fish, the ratios and distributions of T cell populations are not well defined [71, 72], although CD4+ and CD8+ subsets appear to have different roles in pathogen defence . Our results suggest their relative importance to shaping fish mucosal microbiomes, or vice versa, warrant further investigation. Disentangling the directionality of the associations, we find via controlled manipulations of host immune cell populations [74, 75], clock gene expression  and microbiota  will undoubtedly be key to advancing the concept of circadian holobiont health.
Our study demonstrates the complex daily interaction of fish immune expression and microbiomes, which are impacted by photoperiod and infection status. There is rapidly growing recognition for the detrimental impacts of circadian rhythm perturbation in human medicine , though little attention has been paid to the implications for animal health. In an industry that heavily utilises light manipulation, contemporary aquaculture practices may be significantly exacerbating current disease issues. We provide here an important resource for furthering efforts to integrate chronobiology into animal disease mitigation strategies. In addition, as artificial light at night (i.e. light pollution) encroaches on ever greater proportions of the world’s ecosystems [78, 79], we propose that circadian disruption may have as yet undiscovered implications for health and disease dynamics in wild animal populations.
Experimental design and sample collection
Juvenile female triploid rainbow trout fry (O. mykiss, 10-day post-yolk sac absorption, n = 500) were obtained from a commercial hatchery (Bibury Trout Farm, UK). Fry were visually and microscopically determined free of parasitic infections upon arrival and maintained in a re-circulating aquaculture system (RAS) in Cardiff University (water temperature 12 ± 0.5 °C, pH 7.5 ± 0.2). The rainbow trout were randomly assigned to duplicate tanks (45 × 60 × 60 cm, 150 L) under one of two photoperiod conditions; 12:12 LD (lights on at zeitgeber time 0; ZT0, off at ZT12) or 24:0 LD (constant light, 24:0 LD). Each tank was individually illuminated with a full-spectrum white LED bar (80 lux at surface) and surrounded with blackout material to ensure no disturbance from ambient light. Fish were fed with a commercial trout feed (Nutraparr, Skretting, UK) ad libitum at ZT2-3 and ZT9-10 daily. Water oxygen saturation (> 90%), ammonia (< 0.02 mg/L), nitrite (< 0.01 mg L−1) and nitrate (< 15 mg L−1) were maintained within an appropriate range.
After 1 month acclimation to light conditions, 130 fish from each light treatment were individually isolated in 1 L plastic containers. Half of the fish from each light treatment (n = 65 per treatment) were individually inoculated with ten Argulus foliaceus metanauplii (24-h post-hatching). Argulus metanauplii were obtained from eggs of wild-caught adult pairs (sourced from Risca Canal, Newport), maintained at Cardiff University. Egg strings were collected and hatched under laboratory conditions according to Stewart et al. (2018). Inoculations were performed at ZT4-5. Fish were individually held in a glass container with 50 ml of tank water and 10 metanauplii added. Fish were observed until all lice had attached (within 2 min) and then returned to their 1 L container. Control fish (those not inoculated with Argulus lice) were also held for 2 min in 50 ml of water to control for handling stress. Water in all individual containers were changed daily, feeding continued on schedule outlined above, and light conditions maintained at same intensity, spectrum and duration as during acclimation period. The remaining fish were maintained in the RAS system. Once a week, 30 random fish per light treatment were weighed (g) and measured (standard length, SL in cm) for 16 weeks to monitor growth rates. General linear models of standard length and weight, including photoperiod and sampling day, were used to assess differences in growth between light treatments. All procedures were performed under Home Office project licence PPL 303424 with full approval of Cardiff University Animal Ethics committee.
One week after inoculation, sampling of fish was performed over a 48-h period to encompass two full circadian cycles. Starting at ZT0 (lights on in 12:12 LD treatment), every 4 h, five fish from each condition (12:12 LD control, 12:12 LD Argulus-infected, 24:0 LD control, 24:0 LD Argulus-infected) were euthanised using an overdose of tricaine methanesulfonate (MS222, 500 mg L− 1) according to Home Office Schedule 1. At timepoints during dark periods in 12:12 LD treatment, fish were handled and euthanised in dim red light. Immediately after euthanasia, infected fish were visually inspected to quantify number of lice surviving and the lice removed to ensure they were not included in tissue samples. Welch’s two sample T test was used to determine difference in infection load (number of Argulus) between light treatments. All sampled fish were weighed (g) and measured (standard length, SL in cm). Skin swabs (MWE MW-100) were rubbed along the entire lateral body surface five times each side and immediately frozen at −80 C to preserve skin mucus microbiota for DNA extraction. All skin from immediately posterior to opercula to the caudal peduncle was dissected using sterile forceps, preserved in RNAlater (Invitrogen), and stored at −80 °C until RNA extraction. All dissections for each timepoint were performed within an hour window. At each timepoint-treatment combination, 10 ml of water from all containers was pooled and frozen at −80 °C to provide background controls for skin microbiome analyses. To test for endogenous expression rhythms, an additional 65 uninfected fish maintained at 12:12 LD were individually isolated and held in constant darkness (DD). After 24 h, starting at ZT0, five fish every 4 h were sampled as above.
RNA extraction, gene expression quantification and analyses
Total RNA was individually extracted from each skin sample using RNeasy Mini kits (Qiagen). RNA was quantified using Qubit Broad Range RNA assays (ThermoFisher Scientific). mRNA expression patterns in the skin were measured by Nanostring analysis, following manufacturer’s guidelines, at Liverpool Centre for Genomic Research. The nCounter PlexSet oligonucleotide and probe design was performed at NanoString Technologies (NanoString Technologies) for 48 genes, including four housekeeping genes (Supplementary Table 9). The oligonucleotide probes were synthesised at Integrated DNA Technologies. Titration reactions were performed according to supplier’s instructions with RNA inputs between 250 and 700 ng to determine the required RNA amount for hybridization reaction. Six hundred nanograms of total RNA per sample was used for PlexSet hybridization reaction for 20 h according to manufacturer’s instructions.
Samples were processed on a nCounter MAX prep station (NanoString Technologies) and cartridges were scanned in a generation II nCounter Digital Analyzer (NanoString Technologies). RCC files (nCounter data files) were used for data analysis. RCC files were imported into the NanoString nSolver 4.0 analysis software and raw data pre-processing and normalisation was performed according to manufacturer’s instructions for standard procedures (positive normalisation to geomean of top 3 positive controls, codeset content normalisation using housekeeping genes hprt1, polr1b, polr2i, and codeset calibration with the reference sample). The housekeeping gene rplp0 and aanat2 expression were not detected and excluded from analyses.
To assess overall differences in immune responses to infection under the different light treatments, pairwise t tests comparing normalised expression of immune genes were performed in R (version 4.0.3). To detect rhythmicity in expression of clock and immune genes, empirical JTK Cycle (eJTK_cycle ) analyses were applied with a set period of 24 h, a phase search every 4 h from ZT0 to ZT20, and an asymmetry search every 4 h from ZT4 to ZT20. FDR-corrected empirical p values less than 0.1 were considered moderately rhythmic [81,82,83], and less than 0.05 strongly rhythmic. CircaCompare  was used to estimate rhythmic genes’ peak expression time, mesor and amplitude, and to statistically contrast rhythms.
DNA extraction, 16S rRNA gene metabarcoding and metagenomic analyses
DNA was extracted from skin swabs using Qiagen DNeasy Blood and Tissue kits according to Gill et al.  to maximise lysis of microbiome community and DNA recovery. PCR amplification of the 16S rRNA V4 region, using 515F and 806R primers, was performed in triplicate for each DNA extract, pooled and prepared for Illumina MiSeq sequencing according to Caporaso et al.  (1-step PCR 16S amplification and incorporation of Illumina adapters/indexes). Gel electrophoresis was used to estimate concentrations for pooling individual amplicon libraries. Negative controls for extractions and PCR, and mock community positive (HM-783D, BEI Resources) controls were included for sequencing. Libraries were sequenced using a 2 × 250 bp Illumina MiSeq v2 run at the Cardiff University School of Biosciences Genomics Hub.
Paired-end demultiplexed Illumina sequencing reads were imported into the Quantitative Insights Into Microbial Ecology 2 (QIIME2) . Sequences were then quality filtered, dereplicated, chimaeras identified and paired-end reads merged in QIIME2 using DADA2  with default settings (--p-trunc-len-f 225, --p-trunc-len-r 196, --p-max-ee-f/r 2, --p-trunc-q 2, minimum overlap = 12 bp, no mismatch). Classification of amplicon sequence variants (ASVs) was performed using a Naïve Bayes algorithm trained using sequences representing the bacterial V4 rRNA region available from the SILVA database (https://www.arb-silva.de/download/archive/qiime; Silva_132), and the corresponding taxonomic classifications were obtained using the q2-feature-classifier plugin in QIIME2. The classifier was then used to assign taxonomic information to representative sequences of each ASV. Following rarefaction analysis, samples with less than 2000 sequences were excluded from further analyses. QIIME2 was used to analyse alpha (Kruskal-Wallis pairwise tests of Faith’s phylogenetic distance) and beta (pairwise PERMANOVA) diversity measures. ASVs were filtered to exclude those assigned to eukaryotes or eukaryotic organelles and include ones with at least 100 copies in at least two samples. The QIIME2 output data were imported in RStudio (version 1.3.959) with the Bioconductor package phyloseq , for subsetting, normalising and plotting of the data. No reads passed data pre-processing and filtering in extraction blanks and PCR negative controls. The mock community positive control profile matched the manufacturer’s expected composition and relative abundances (data not shown).
Differential abundance of ASVs between healthy and infected fish in both light treatments was determined using DESeq2 , with FDR-corrected p values less than 0.05 considered significant. Differential abundances of all taxonomic levels were also determined and visualised using MicrobiomeAnalyst  heat trees  using default settings. We inferred the microbial gene content from the taxa abundance using PICRUSt2 . We used LefSe analyses to identify group differences in the inferred gene abundance of MetaCyc pathways, using the online galaxy server (https://huttenhower.sph.harvard.edu/galaxy/). LDA scores > 2.0 were considered significant. Rhythmicity of microbial genera and MetaCyc pathway abundances were determined following the same methods as gene expression (see above). To determine potential associations of host gene expression and the microbiome, Spearman correlation tests (R package Hmisc rcorr function) were performed including only genera found in at least 50% of samples in each treatment group. Corrected p values (using qvalue R package FDR correction) of less than 0.05 were considered significantly correlated. Correlation networks were visualised using gephi  (with Force Atlas2 algorithm) and influential nodes determined using degree centrality scores and number of connections (degrees).
A subset of 12 samples (3 healthy and 3 infected individuals under 12:12 LD from two timepoints; onset of light and onset of dark, randomly selected) were prepared for full metagenomic sequencing. To reduce host DNA, aliquots of swab extracts (the same used for 16S rRNA profiling) were prepared using the NEBNext Microbiome Enrichment Kit (New England Biolabs) according to manufacturer instructions. The microbial-enriched DNA (6-24 ng input per sample) was then prepared for sequencing using the Illumina DNA Prep Kit (Illumina) according to manufacturer instructions. Libraries were quantified and pooled equimolarly, based on Qubit HS DNA assays. Library quality was checked using Agilent Tapestation D1000 assays. Indexed libraries were sequenced on a high-throughput 2 × 150 bp NextSeq 550 run at Cardiff University School of Biosciences Genomics Hub.
Adapter removal and quality trimming of raw metagenomic reads was performed using Trimmomatic v0.39  with default parameters except increased quality score thresholds of Q30. Duplicate reads were removed using seqkit  and paired reads repaired using BBtools . Deduplicated reads were filtered against the rainbow trout genome  using minimap2 . Paired and single filtered reads from all samples (78.85 Gb, average 6.57 Gb per sample) were used for metagenome co-assembly in Megahit v1.2.9  using meta-sensitive pre-set. Sample reads were individually mapped to assembled contigs using bowtie2  with default parameters, and contigs were binned using Metabat2 v2.15 . Metagenomic assembled genomes (MAGs) were checked for quality using CheckM  and refined using Anvi’o 7 . MAG phylotyping was performed using Amphora2 , with phylotype assigned to lowest taxonomic level with at least 75% agreement in markers. Species trees, using all single-copy orthologs, were generated for selected MAGs (bin 3 Flavobacterium sp., bin 11 Pseudomonas sp.), against other named species within their genus with complete genomes (Flavobacterium = 33, Pseudomonas = 92), using Orthofinder . Differential MAG abundance between treatment groups was determined using DESeq2  with FDR-corrected p value threshold < 0.05. Trimmed reads were also aligned against the NCBI nr protein database using Diamond  and imported into Megan 6  for taxonomic profiling and functional assignments. Raw counts of SEED subsystems were imported into R for normalisation and differential abundance testing using DESeq2 (FDR-corrected p value threshold < 0.05).
Availability of data and materials
16S rRNA gene sequence data are submitted to the NCBI Sequence Read Archive (SRA) under accession PRJNA694669. Metagenomic sequence data are submitted to NCBI SRA under accession PRJNA757611. All other datasets analysed during the current study and their associated scripts are available from the corresponding author on request.
Dunlap JC. Molecular bases for circadian clocks. Cell. 1999;96(2):271–90. https://doi.org/10.1016/S0092-8674(00)80566-8.
Lochmiller RL, Deerenberg C. Trade-offs in evolutionary immunology: just what is the cost of immunity? Oikos. 2000;88(1):87–98. https://doi.org/10.1034/j.1600-0706.2000.880110.x.
Wang W, Barnaby JY, Tada Y, Li H, Tör M, Caldelari D, et al. Timing of plant immune responses by a central circadian regulator. Nature. 2011;470(7332):110–4. https://doi.org/10.1038/nature09766.
Scheiermann C, Kunisaki Y, Frenette PS. Circadian control of the immune system. Nat. Rev. Immunol. 2013;13(3):190–8. https://doi.org/10.1038/nri3386.
Curtis AM, Bellet MM, Sassone-Corsi P, O’Neill LAJ. Circadian clock proteins and immunity. Immunity. 2014;40(2):178–86. https://doi.org/10.1016/j.immuni.2014.02.002.
Okada K, Yano M, Doki Y, Azama T, Iwanaga H, Miki H, et al. Injection of LPS causes transient suppression of biological clock genes in rats. J. Surg. Res. 2008;145(1):5–12. https://doi.org/10.1016/j.jss.2007.01.010.
Cavadini G, Petrzilka S, Kohler P, Jud C, Tobler I, Birchler T, et al. TNF-α suppresses the expression of clock genes by interfering with E-box-mediated transcription. Proc. Natl. Acad. Sci. 2007;104(31):12843–8. https://doi.org/10.1073/pnas.0701466104.
De Leone MJ, Hernando CE, Romanowski A, Careno DA, Soverna AF, Sun H, ... & Yanovsky MJ. Bacterial infection disrupts clock gene expression to attenuate immune responses. Current Biology. 2020;30(9):1740–47.
Shirasu-Hiza MM, Dionne MS, Pham LN, Ayres JS, Schneider DS. Interactions between circadian rhythm and immunity in Drosophila melanogaster. Curr. Biol. 2007;17(10):R353–5. https://doi.org/10.1016/j.cub.2007.03.049.
Marpegán L, Bekinschtein TA, Costas MA, Golombek DA. Circadian responses to endotoxin treatment in mice. J. Neuroimmunol. 2005;160(1-2):102–9. https://doi.org/10.1016/j.jneuroim.2004.11.003.
Castanon-Cervantes O, Wu M, Ehlen JC, Paul K, Gamble KL, Johnson RL, et al. Dysregulation of inflammatory responses by chronic circadian disruption. J. Immunol. 2010;185(10):5796–805. https://doi.org/10.4049/jimmunol.1001026.
Adams KL, Castanon-Cervantes O, Evans JA, Davidson AJ. Environmental circadian disruption elevates the IL-6 response to lipopolysaccharide in blood. J. Biol. Rhythms. 2013;28(4):272–7. https://doi.org/10.1177/0748730413494561.
Touitou Y, Reinberg A, Touitou D. Association between light at night, melatonin secretion, sleep deprivation, and the internal clock: health impacts and mechanisms of circadian disruption. Life Sci. 2017;173:94–106. https://doi.org/10.1016/j.lfs.2017.02.008.
Thaiss CA, Zmora N, Levy M, Elinav E. The microbiome and innate immunity. Nature. 2016;535(7610):65–74. https://doi.org/10.1038/nature18847.
McDermott AJ, Huffnagle GB. The microbiome and regulation of mucosal immunity. Immunology. 2014;142(1):24–31. https://doi.org/10.1111/imm.12231.
Zarrinpar A, Chaix A, Yooseph S, Panda S. Diet and feeding pattern affect the diurnal dynamics of the gut microbiome. Cell Metab. 2014;20(6):1006–17. https://doi.org/10.1016/j.cmet.2014.11.008.
Liang X, Bushman FD, FitzGerald GA. Rhythmicity of the intestinal microbiota is regulated by gender and the host circadian clock. Proc. Natl. Acad. Sci. 2015;112(33):10479–84. https://doi.org/10.1073/pnas.1501305112.
Leone V, Gibbons SM, Martinez K, Hutchison AL, Huang EY, Cham CM, et al. Effects of diurnal variation of gut microbes and high-fat feeding on host circadian clock function and metabolism. Cell Host Microbe. 2015;17(5):681–9. https://doi.org/10.1016/j.chom.2015.03.006.
Stentiford GD, Sritunyalucksana K, Flegel TW, Williams BAP, Withyachumnarnkul B, Itsathitphaisarn O, et al. New paradigms to help solve the global aquaculture disease crisis. PLoS Pathog. 2017;13(2):e1006160. https://doi.org/10.1371/journal.ppat.1006160.
Lafferty, K. D. et al. Infectious diseases affect marine fisheries and aquaculture economics. (2015).
Tarnecki AM, Burgos FA, Ray CL, Arias CR. Fish intestinal microbiome: diversity and symbiosis unravelled by metagenomics. J. Appl. Microbiol. 2017;123(1):2–17. https://doi.org/10.1111/jam.13415.
Ghanbari M, Kneifel W, Domig KJ. A new view of the fish gut microbiome: advances from next-generation sequencing. Aquaculture. 2015;448:464–75. https://doi.org/10.1016/j.aquaculture.2015.06.033.
Perry WB, Lindsay E, Payne CJ, Brodie C, Kazlauskaite R. The role of the gut microbiome in sustainable teleost aquaculture. Proc. R. Soc. B. 2020;287(1926):20200184. https://doi.org/10.1098/rspb.2020.0184.
Wong S, Waldrop T, Summerfelt S, Davidson J, Barrows F, Kenney PB, et al. Aquacultured rainbow trout (Oncorhynchus mykiss) possess a large core intestinal microbiota that is resistant to variation in diet and rearing density. Appl. Environ. Microbiol. 2013;79(16):4974–84. https://doi.org/10.1128/AEM.00924-13.
Ingerslev H-C, Strube ML, Jørgensen LG, Dalsgaard I, Boye M, Madsen L. Diet type dictates the gut microbiota and the immune response against Yersinia ruckeri in rainbow trout (Oncorhynchus mykiss). Fish Shellfish Immunol. 2014;40(2):624–33. https://doi.org/10.1016/j.fsi.2014.08.021.
Biswas AK, Seoka M, Tanaka Y, Takii K, Kumai H. Effect of photoperiod manipulation on the growth performance and stress response of juvenile red sea bream (Pagrus major). Aquaculture. 2006;258(1-4):350–6. https://doi.org/10.1016/j.aquaculture.2006.03.048.
Rad F, Bozaoğlu S, Gözükara SE, Karahan A, Kurt G. Effects of different long-day photoperiods on somatic growth and gonadal development in Nile tilapia (Oreochromis niloticus L.). Aquaculture. 2006;255(1-4):292–300. https://doi.org/10.1016/j.aquaculture.2005.11.028.
Berrill IK, Porter MJR, Smart A, Mitchell D, Bromage NR. Photoperiodic effects on precocious maturation, growth and smoltification in Atlantic salmon, Salmo salar. Aquaculture. 2003;222(1-4):239–52. https://doi.org/10.1016/S0044-8486(03)00125-X.
Onoue T, Nishi G, Hikima J, Sakai M, Kono T. Circadian oscillation of TNF-α gene expression regulated by clock gene, BMAL1 and CLOCK1, in the Japanese medaka (Oryzias latipes). Int. Immunopharmacol. 2019;70:362–71. https://doi.org/10.1016/j.intimp.2019.02.004.
Zhang P, Yu C, Sun L. Japanese flounder (Paralichthys olivaceus) Bmal1 is involved in the regulation of inflammatory response and bacterial infection. Aquaculture. 2020;525:735330. https://doi.org/10.1016/j.aquaculture.2020.735330.
Frøland Steindal IA, Whitmore D. Circadian clocks in fish—what have we learned so far? Biology (Basel). 2019;8(17).
Binuramesh C, Michael RD. Diel variations in the selected serum immune parameters in Oreochromis mossambicus. Fish Shellfish Immunol. 2011;30(3):824–9. https://doi.org/10.1016/j.fsi.2011.01.004.
Ángeles Esteban M, Cuesta A, Rodríguez A, Meseguer J. Effect of photoperiod on the fish innate immune system: a link between fish pineal gland and the immune system. J. Pineal Res. 2006;41(3):261–6. https://doi.org/10.1111/j.1600-079X.2006.00362.x.
Lazado CC, Skov PV, Pedersen PB. Innate immune defenses exhibit circadian rhythmicity and differential temporal sensitivity to a bacterial endotoxin in Nile tilapia (Oreochromis niloticus). Fish Shellfish Immunol. 2016;55:613–22. https://doi.org/10.1016/j.fsi.2016.06.040.
Taira G, Onoue T, Hikima J, Sakai M, Kono T. Circadian clock components Bmal1 and Clock1 regulate tlr9 gene expression in the Japanese medaka (Oryzias latipes). Fish Shellfish Immunol. 2020;105:438–45. https://doi.org/10.1016/j.fsi.2020.07.009.
Ellison AR, Uren Webster TM, Rey O, Garcia de Leaniz C, Consuegra S, Orozco-terWengel P, et al. Transcriptomic response to parasite infection in Nile tilapia (Oreochromis niloticus) depends on rearing density. BMC Genomics. 2018;19(1):723. https://doi.org/10.1186/s12864-018-5098-7.
Ellison AR, Uren Webster TM, Rodriguez-Barreto D, de Leaniz CG, Consuegra S, Orozco-terWengel P, et al. Comparative transcriptomics reveal conserved impacts of rearing density on immune response of two important aquaculture species. Fish Shellfish Immunol. 2020;104:192–201. https://doi.org/10.1016/j.fsi.2020.05.043.
Nikolaidis M, Mossialos D, Oliver SG, Amoutzias GD. Comparative analysis of the core proteomes among the Pseudomonas major evolutionary groups reveals species-specific adaptations for Pseudomonas aeruginosa and Pseudomonas chlororaphis. Diversity. 2020;12(8):289. https://doi.org/10.3390/d12080289.
Douglas GM, Maffei VJ, Zaneveld JR, Yurgel SN, Brown JR, Taylor CM, ... & Langille MG. PICRUSt2 for prediction of metagenome functions. Nature Biotechnology. 2020;38(6):685–88.
Greer R, Dong X, Morgun A, Shulzhenko N. Investigating a holobiont: microbiota perturbations and transkingdom networks. Gut Microbes. 2016;7(2):126–35. https://doi.org/10.1080/19490976.2015.1128625.
Greenblum S, Turnbaugh PJ, Borenstein E. Metagenomic systems biology of the human gut microbiome reveals topological shifts associated with obesity and inflammatory bowel disease. Proc. Natl. Acad. Sci. 2012;109(2):594–9. https://doi.org/10.1073/pnas.1116053109.
Whiting JR, Mahmud MA, Bradley JE, MacColl ADC. Prior exposure to long-day photoperiods alters immune responses and increases susceptibility to parasitic infection in stickleback. Proc. R. Soc. B. 2020;287:20201017.
Braden LM, Koop BF, Jones SRM. Signatures of resistance to Lepeophtheirus salmonis include a TH2-type response at the louse-salmon interface. Dev. Comp. Immunol. 2015;48(1):178–91. https://doi.org/10.1016/j.dci.2014.09.015.
Saeidi A, Zandi K, Cheok YY, Saeidi H, Wong WF, Lee CYQ, et al. T-cell exhaustion in chronic infections: reversing the state of exhaustion and reinvigorating optimal protective immune responses. Front. Immunol. 2018;9:2569. https://doi.org/10.3389/fimmu.2018.02569.
Graham AL, Allen JE, Read AF. Evolutionary causes and consequences of immunopathology. Annu. Rev. Ecol. Evol. Syst. 2005:373–97.
Westwood ML, O’Donnell AJ, de Bekker C, Lively CM, Zuk M, Reece SE. The evolutionary ecology of circadian rhythms in infection. Nat. Ecol. Evol. 2019;3(4):552–60. https://doi.org/10.1038/s41559-019-0831-4.
Reece SE, Prior KF, Mideo N. The life and times of parasites: rhythms in strategies for within-host survival and between-host transmission. J. Biol. Rhythms. 2017;32(6):516–33. https://doi.org/10.1177/0748730417718904.
Longcore T, Rich C. Ecological light pollution. Front. Ecol. Environ. 2004;2(4):191–8. https://doi.org/10.1890/1540-9295(2004)002[0191:ELP]2.0.CO;2.
Ross AA, Hoffmann AR, Neufeld JD. The skin microbiome of vertebrates. Microbiome. 2019;7(1):1–14. https://doi.org/10.1186/s40168-019-0694-6.
Pérez-Sánchez T, Balcázar JL, García Y, Halaihel N, Vendrell D, de Blas I, et al. Identification and characterization of lactic acid bacteria isolated from rainbow trout, Oncorhynchus mykiss (Walbaum), with inhibitory activity against Lactococcus garvieae. J. Fish Dis. 2011;34(7):499–507. https://doi.org/10.1111/j.1365-2761.2011.01260.x.
Balcázar JL, Vendrell D, de Blas I, Ruiz-Zarzuela I, Gironés O, Múzquiz JL. In vitro competitive adhesion and production of antagonistic compounds by lactic acid bacteria against fish pathogens. Vet. Microbiol. 2007;122(3-4):373–80. https://doi.org/10.1016/j.vetmic.2007.01.023.
Llewellyn MS, Boutin S, Hoseinifar SH, Derome N. Teleost microbiomes: the state of the art in their characterization, manipulation and importance in aquaculture and fisheries. Front Microbiol. 2014;5:207.
Walker PD, Flik G, Bonga SEW. The biology of parasites from the genus Argulus and a review of the interactions with its host. Host-Parasite Interact. 2004;55:107–29. https://doi.org/10.4324/9780203487709-6.
Adikesavalu H, Patra A, Banerjee S. Sarkar, A. & Abraham, T. J. Phenotypic and molecular characterization and pathology of Flectobacillus roseus causing flectobacillosis in captive held carp Labeo rohita (Ham.) fingerlings. Aquaculture. 2015;439:60–5. https://doi.org/10.1016/j.aquaculture.2014.12.036.
Loch TP, Faisal M. Emerging flavobacterial infections in fish: a review. J. Adv. Res. 2015;6(3):283–300. https://doi.org/10.1016/j.jare.2014.10.009.
Jakob E, Barker DE, Garver KA. Vector potential of the salmon louse Lepeophtheirus salmonis in the transmission of infectious haematopoietic necrosis virus (IHNV). Dis. Aquat. Organ. 2011;97(2):155–65. https://doi.org/10.3354/dao02414.
Ahne W. Argulus foliaceus L. and Piscicola geometra L. as mechanical vectors of spring viraemia of carp virus (SVCV). J. Fish Dis. 1985;8(2):241–2. https://doi.org/10.1111/j.1365-2761.1985.tb01220.x.
Lloyd-Price J, Abu-Ali G, Huttenhower C. The healthy human microbiome. Genome Med. 2016;8(1):1–11. https://doi.org/10.1186/s13073-016-0307-y.
Harris EV, de Roode JC, Gerardo NM. Diet–microbiome–disease: investigating diet’s influence on infectious disease resistance through alteration of the gut microbiome. PLoS Pathog. 2019;15:e1007891.
Cani PD, van Hul M, Lefort C, Depommier C, Rastelli M, Everard A. Microbial regulation of organismal energy homeostasis. Nat. Metab. 2019;1(1):34–46. https://doi.org/10.1038/s42255-018-0017-4.
Fiore CL, Jarett JK, Steinert G, Lesser MP. Trait-based comparison of coral and sponge microbiomes. Sci. Rep. 2020;10(1):1–16. https://doi.org/10.1038/s41598-020-59320-9.
Brüning A, Kloas W, Preuer T, Hölker F. Influence of artificially induced light pollution on the hormone system of two common fish species, perch and roach, in a rural habitat. Conserv. Physiol. 2018;6:coy016.
Sánchez-Vázquez FJ, Tabata M. Circadian rhythms of demand-feeding and locomotor activity in rainbow trout. J. Fish Biol. 1998;52(2):255–67. https://doi.org/10.1111/j.1095-8649.1998.tb00797.x.
Pérez-Sánchez T, Ruiz-Zarzuela I, de Blas I, Balcázar JL. Probiotics in aquaculture: a current assessment. Rev. Aquac. 2014;6(3):133–46. https://doi.org/10.1111/raq.12033.
Ballesta A, Innominato PF, Dallmann R, Rand DA, Lévi FA. Systems chronotherapeutics. Pharmacol. Rev. 2017;69(2):161–99. https://doi.org/10.1124/pr.116.013441.
Pearson JA, Wong FS, Wen L. Cross talk between circadian rhythms and the microbiota. Immunology. 2020;161(4):278–90. https://doi.org/10.1111/imm.13278.
Gibbs JE, Butler TD. Circadian host-microbiome interactions in immunity. Front. Immunol. 2020;11:1783.
Earley AM, Graves CL, Shiau CE. Critical role for a subset of intestinal macrophages in shaping gut microbiota in adult zebrafish. Cell Rep. 2018;25(2):424–36. https://doi.org/10.1016/j.celrep.2018.09.025.
Brugman S, Schneeberger K, Witte M, Klein MR, van den Bogert B, Boekhorst J, et al. T lymphocytes control microbial composition by regulating the abundance of Vibrio in the zebrafish gut. Gut Microbes. 2014;5(6):737–47. https://doi.org/10.4161/19490976.2014.972228.
Xu Z, Takizawa F, Casadei E, Shibasaki Y, Ding Y, Sauters TJC, et al. Specialization of mucosal immunoglobulins in pathogen control and microbiota homeostasis occurred early in vertebrate evolution. Sci. Immunol. 2020;5(44):eaay3254. https://doi.org/10.1126/sciimmunol.aay3254.
Takizawa F, Dijkstra JM, Kotterba P, Korytář T, Kock H, Köllner B, et al. The expression of CD8α discriminates distinct T cell subsets in teleost fish. Dev. Comp. Immunol. 2011;35(7):752–63. https://doi.org/10.1016/j.dci.2011.02.008.
Kelly C, Salinas I. Under pressure: interactions between commensal microbiota and the teleost immune system. Front. Immunol. 2017;8:559. https://doi.org/10.3389/fimmu.2017.00559.
Sukeda M, et al. Innate cell-mediated cytotoxicity of CD8+ T cells against the protozoan parasite Ichthyophthirius multifiliis in the ginbuna crucian carp, Carassius auratus langsdorfii. Dev. Comp. Immunol. 2020;115:103886.
Moore JC, Mulligan TS, Yordán NT, Castranova D, Pham VN, Tang Q, et al. T cell immune deficiency in zap70 mutant zebrafish. Mol. Cell. Biol. 2016;36(23):2868–76. https://doi.org/10.1128/MCB.00281-16.
Sugimoto K, Hui SP, Sheng DZ, Nakayama M, Kikuchi K. Zebrafish FOXP3 is required for the maintenance of immune tolerance. Dev. Comp. Immunol. 2017;73:156–62. https://doi.org/10.1016/j.dci.2017.03.023.
Alifu Y, Kofuji S, Sunaga S, Kusaba M, Hirayama J, Nishina H. The light-inducible genes Per2, Cry1a, and Cry2a regulate oxidative status in zebrafish. Biol. Pharm. Bull. 2021;44(8):1160–5. https://doi.org/10.1248/bpb.b21-00432.
Luna GM, Quero GM, Kokou F, Kormas K. Time to integrate biotechnological approaches into fish gut microbiome research. Curr. Opin. Biotechnol. 2022;73:121–7. https://doi.org/10.1016/j.copbio.2021.07.018.
Hölker F, Moss T, Griefahn B, Kloas W, Voigt CC, Henckel D, et al. The dark side of light: a transdisciplinary research agenda for light pollution policy. Ecol. Soc. 2010;15(4). https://doi.org/10.5751/ES-03685-150413.
Tidau S, Smyth T, McKee D, Wiedenmann J, D’Angelo C, Wilcockson D, Ellison A, Grimmer AJ, Jenkins SR, Widdicombe S, Queirós AM, Talbot E, Wright A, Davies TW. Marine artificial light at night: An empirical and technical guide. Methods Ecol Evol. 2021;12:1588–1601. https://doi.org/10.1111/2041-210X.13653
Hutchison AL, Maienschein-Cline M, Chiang AH, Tabei SMA, Gudjonson H, Bahroos N, et al. Improved statistical methods enable greater sensitivity in rhythm detection for genome-wide data. PLoS Comput Biol. 2015;11(3):e1004094. https://doi.org/10.1371/journal.pcbi.1004094.
Wang Y, et al. A proteomics landscape of circadian clock in mouse liver. Nat. Commun. 2018;9:1–16.
Lafaye G, Desterke C, Marulaz L, Benyamina A. Cannabidiol affects circadian clock core complex and its regulation in microglia cells. Addict. Biol. 2019;24(5):921–34. https://doi.org/10.1111/adb.12660.
Cui P, et al. Identification of human circadian genes based on time course gene expression profiles by using a deep learning method. Biochim. Biophys. Acta (BBA)-Mol Basis Dis. 2018;1864:2274–83.
Gill C, van de Wijgert JHHM, Blow F, Darby AC. Evaluation of lysis methods for the extraction of bacterial DNA for analysis of the vaginal microbiota. PLoS One. 2016;11(9):e0163148. https://doi.org/10.1371/journal.pone.0163148.
Caporaso JG, Lauber CL, Walters WA, Berg-Lyons D, Huntley J, Fierer N, et al. Ultra-high-throughput microbial community analysis on the Illumina HiSeq and MiSeq platforms. ISME J. 2012;6(8):1621–4. https://doi.org/10.1038/ismej.2012.8.
Bolyen E, Rideout JR, Dillon MR, Bokulich NA, Abnet CC, al-Ghalith GA, et al. Reproducible, interactive, scalable and extensible microbiome data science using QIIME 2. Nat. Biotechnol. 2019;37(8):852–7. https://doi.org/10.1038/s41587-019-0209-9.
Callahan BJ, McMurdie PJ, Rosen MJ, Han AW, Johnson AJA, Holmes SP. DADA2: high-resolution sample inference from Illumina amplicon data. Nat. Methods. 2016;13(7):581–3. https://doi.org/10.1038/nmeth.3869.
McMurdie PJ, Holmes S. phyloseq: an R package for reproducible interactive analysis and graphics of microbiome census data. PLoS One. 2013;8:e61217.
Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-Seq data with DESeq2. Genome Biol. 2014;15(12):550. https://doi.org/10.1186/s13059-014-0550-8.
Dhariwal A, Chong J, Habib S, King IL, Agellon LB, Xia J. MicrobiomeAnalyst: a web-based tool for comprehensive statistical, visual and meta-analysis of microbiome data. Nucleic Acids Res. 2017;45(W1):W180–8. https://doi.org/10.1093/nar/gkx295.
Foster ZSL, Sharpton TJ, Grünwald NJ. Metacoder: An R package for visualization and manipulation of community taxonomic diversity data. PLoS Comput. Biol. 2017;13(2):e1005404. https://doi.org/10.1371/journal.pcbi.1005404.
Bastian, M., Heymann, S. & Jacomy, M. Gephi: an open source software for exploring and manipulating networks. in Proceedings of the International AAAI Conference on Web and Social Media vol. 3 (2009).
Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30(15):2114–20. https://doi.org/10.1093/bioinformatics/btu170.
Shen W, Le S, Li Y, Hu F. SeqKit: a cross-platform and ultrafast toolkit for FASTA/Q file manipulation. PLoS One. 2016;11(10):e0163962. https://doi.org/10.1371/journal.pone.0163962.
Bushnell, B. BBTools. https://sourceforge.net/projects/bbmap/.
Berthelot C, et al. The rainbow trout genome provides novel insights into evolution after whole-genome duplication in vertebrates. Nat. Commun. 2014;5:1–10.
Li H. Minimap2: pairwise alignment for nucleotide sequences. Bioinformatics. 2018;34(18):3094–100. https://doi.org/10.1093/bioinformatics/bty191.
Li D, Liu C-M, Luo R, Sadakane K, Lam T-W. MEGAHIT: an ultra-fast single-node solution for large and complex metagenomics assembly via succinct de Bruijn graph. Bioinformatics. 2015;31(10):1674–6. https://doi.org/10.1093/bioinformatics/btv033.
Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat. Methods. 2012;9(4):357–9. https://doi.org/10.1038/nmeth.1923.
Kang DD, Li F, Kirton E, Thomas A, Egan R, An H, et al. MetaBAT 2: an adaptive binning algorithm for robust and efficient genome reconstruction from metagenome assemblies. PeerJ. 2019;7:e7359. https://doi.org/10.7717/peerj.7359.
Parks DH, Imelfort M, Skennerton CT, Hugenholtz P, Tyson GW. CheckM: assessing the quality of microbial genomes recovered from isolates, single cells, and metagenomes. Genome Res. 2015;25(7):1043–55. https://doi.org/10.1101/gr.186072.114.
Eren AM, Esen ÖC, Quince C, Vineis JH, Morrison HG, Sogin ML, et al. Anvi’o: an advanced analysis and visualization platform for ‘omics data. PeerJ. 2015;3:e1319. https://doi.org/10.7717/peerj.1319.
Wu M, Scott AJ. Phylogenomic analysis of bacterial and archaeal sequences with AMPHORA2. Bioinformatics. 2012;28(7):1033–4. https://doi.org/10.1093/bioinformatics/bts079.
Emms DM, Kelly S. OrthoFinder: phylogenetic orthology inference for comparative genomics. Genome Biol. 2019;20(1):1–14. https://doi.org/10.1186/s13059-019-1832-y.
Buchfink B, Xie C, Huson DH. Fast and sensitive protein alignment using DIAMOND. Nat. Methods. 2015;12(1):59–60. https://doi.org/10.1038/nmeth.3176.
Huson DH, Auch AF, Qi J, Schuster SC. MEGAN analysis of metagenomic data. Genome Res. 2007;17(3):377–86. https://doi.org/10.1101/gr.5969107.
We thank the University of Liverpool Centre for Genomic Research and Cardiff University School of Biosciences Genome Hub for their assistance in data generation. We also thank members of the Cable research group for their assistance in animal husbandry.
This study was funded by a BBSRC Discovery Fellowship awarded to AE (BB/R010609/1).
Ethics approval and consent to participate
All procedures were performed under Home Office project licence PPL 303424 with full approval of Cardiff University Animal Ethics committee.
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
. Gene rhythm summary
. mic rhythm summary
Average A) standard length and B) weight of trout (±1 S.E.) over 16-week growth trial under 12:12 LD (orange) and 24:0 LD (yellow). C) Boxplots of number of Argulus foliaceus lice infecting fish 7 days post-inoculation. Supplementary Figure 2. Mean expression (± 1 S.E.) of accessory clock genes of uninfected (cyan) and Argulus-infected (orange) rainbow trout maintained at 12:12 LD (left) and 24:0 LD (LL, right). Expression is normalised counts of mRNA copies detected via Nanostring nCounter. Curves denote cosinor waveform fitted using CircaCompare. Grey shading indicates time periods in darkness (grey dashing indicates equivalent 12:12 LD light transitions on LL plots). Supplementary Figure 3. Mean expression (± 1 S.E.) of clock genes of rainbow trout under 12:12 LD and DD (free-running, constant darkness). Expression is normalised counts of mRNA copies detected via Nanostring nCounter. Curves denote cosinor waveform fitted using CircaCompare. Grey shading indicates time periods in darkness (grey dashing indicates subjective day-night transition in DD). Supplementary Figure 4. Mean expression (± 1 S.E.) of innate immune genes of uninfected (cyan) and Argulus-infected (orange) rainbow trout maintained at 12:12 LD (left) and 24:0 LD (LL, right). Expression is normalised counts of mRNA copies detected via Nanostring nCounter. Curves denote cosinor waveform fitted using CircaCompare. Grey shading indicates time periods in darkness (grey dashing indicates equivalent 12:12 LD light transitions on LL plots). Only genes with significant rhythm in one or more groups shown. Supplementary Figure 5. Mean expression (± 1 S.E.) of adaptive immune genes of uninfected (cyan) and Argulus-infected (orange) rainbow trout maintained at 12:12 LD (left) and 24:0 LD (LL, right). Expression is normalised counts of mRNA copies detected via Nanostring nCounter. Curves denote cosinor waveform fitted using CircaCompare. Grey shading indicates time periods in darkness (grey dashing indicates equivalent 12:12 LD light transitions on LL plots). Only genes with significant rhythm in one or more groups shown. Supplementary Figure 6. Mean expression (± 1 S.E.) of immune genes of rainbow trout under 12:12 LD and DD (free-running, constant darkness). Expression is normalised counts of mRNA copies detected via Nanostring nCounter. Curves denote cosinor waveform fitted using CircaCompare. Grey shading indicates time periods in darkness (grey dashing indicates subjective day-night transition in DD). Supplementary Figure 7. A) Rarefaction plots of detected amplified sequence variants (ASVs) by sampling depth. B) NMDS ordination of microbiome profiles. C) Alpha diversity plots by treatment group. Argulus_12 = Argulus-infected rainbow trout under 12:12 LD, Argulus_24 = Argulus-infected rainbow trout under 24:0 LD, Control_12 = healthy rainbow trout under 12:12 LD, Control_24 = healthy rainbow trout under 24:0 LD. Supplementary Figure 8. Comparison of rainbow trout skin microbiome taxonomic profiles between 16S rRNA metabarcoding (taxonomic assignment of ASVs in Qiime2) and metagenomic sequencing (taxonomic assignment of reads in Megan). Supplementary Figure 9. Left: Heatmap of MAG relative abundance (scaled to median), yellow indicates increased abundance, blue indicates reduced abundance. Supplementary table 4 provides significant differences among treatment group. Right: Examples of MAGs found to be differentially abundant between infection status and/or timepoint. Supplementary Figure 10. Species tree of Pseudomonas genomes generated by OrthoFinder. MAG from current study highlighted in bold. Supplementary Figure 11. Species tree of Flavobacterium genomes generated by OrthoFinder. MAG from current study highlighted in bold. Supplementary Figure 12. Co-occurrence networks of microbial genera (pink) and host gene expression (orange = clock, green = immune, blue = corticotropin) in healthy (top) and Argulus-infected (bottom) trout under 24:0 LD. Node and label size scaled to degree centrality score. Label colour denotes rhythmicity (black = rhythm FDR p-value <0.05, grey = rhythm FDR p-value >0.05). Connection colour indicates association (grey = positive, red = negative, determined by Spearman correlation tests) and connection width scaled to correlation strength (thicker lines denote a higher correlation coefficient).
Supplementary Datafile 1
About this article
Cite this article
Ellison, A.R., Wilcockson, D. & Cable, J. Circadian dynamics of the teleost skin immune-microbiome interface. Microbiome 9, 222 (2021). https://doi.org/10.1186/s40168-021-01160-4
- Circadian rhythm
- Clock gene expression
- Parasite infection