The Southern Bluefin Tuna Mucosal Microbiome Is Influenced by Husbandry Method, Net Pen Location, and Anti-parasite Treatment

Aquaculture is the fastest growing primary industry worldwide. Marine finfish culture in open ocean net pens, or pontoons, is one of the largest growth areas and is currently the only way to rear high value fish such as bluefin tuna. Ranching involves catching wild juveniles, stocking in floating net pens and fattening for 4 to 8 months. Tuna experience several parasite-induced disease challenges in culture that can be mitigated by application of praziquantel (PZQ) as a therapeutic. In this study, we characterized the microbiome of ranched southern Bluefin Tuna, Thunnus maccoyii, across four anatomic sites (gill, skin, digesta, and anterior kidney) and evaluated environmental and pathological factors that influence microbiome composition, including the impact of PZQ treatment on microbiome stability. Southern bluefin tuna gill, skin, and digesta microbiome communities are unique and potentially influenced by husbandry practices, location of pontoon growout pens, and treatment with the antiparasitic PZQ. There was no significant relationship between the fish mucosal microbiome and incidence or abundance of adult blood fluke in the heart or fluke egg density in the gill. An enhanced understanding of microbiome diversity and function in high-value farmed fish species such as bluefin tuna is needed to optimize fish health and improve aquaculture yield. Comparison of the bluefin tuna microbiome to other fish species, including Seriola lalandi (yellowtail kingfish), a common farmed species from Australia, and Scomber japonicus (Pacific mackerel), a wild caught Scombrid relative of tuna, showed the two Scombrids had more similar microbial communities compared to other families. The finding that mucosal microbial communities are more similar in phylogenetically related fish species exposes an opportunity to develop mackerel as a model for tuna microbiome and parasite research.

Aquaculture is the fastest growing primary industry worldwide. Marine finfish culture in open ocean net pens, or pontoons, is one of the largest growth areas and is currently the only way to rear high value fish such as bluefin tuna. Ranching involves catching wild juveniles, stocking in floating net pens and fattening for 4 to 8 months. Tuna experience several parasite-induced disease challenges in culture that can be mitigated by application of praziquantel (PZQ) as a therapeutic. In this study, we characterized the microbiome of ranched southern Bluefin Tuna, Thunnus maccoyii, across four anatomic sites (gill, skin, digesta, and anterior kidney) and evaluated environmental and pathological factors that influence microbiome composition, including the impact of PZQ treatment on microbiome stability. Southern bluefin tuna gill, skin, and digesta microbiome communities are unique and potentially influenced by husbandry practices, location of pontoon growout pens, and treatment with the antiparasitic PZQ. There was no significant relationship between the fish mucosal microbiome and incidence or abundance of adult blood fluke in the heart or fluke egg density in the gill. An enhanced understanding of microbiome diversity and function in high-value farmed fish species such as bluefin tuna is needed to optimize fish health and improve aquaculture yield. Comparison of the bluefin tuna microbiome to other fish species, including Seriola lalandi (yellowtail kingfish), a common farmed species from Australia, and Scomber japonicus (Pacific mackerel), a wild caught Scombrid relative of tuna, showed the two Scombrids

INTRODUCTION
Bluefin tuna is one of the highest value fish in the world. In 2014, the three species of bluefin including Atlantic Bluefin Tuna (ABT), Pacific Bluefin Tuna (PBT), and Southern Bluefin Tuna (SBT) had a combined total value of USD $610-660 million for fishers with a total value of USD $2-2.5 billion at final sale (Macfadyen, 2016). To meet high demands, bluefin tunas are ranched in open ocean net pens, or pontoons, which currently represents approximately 20% of total global bluefin production. In Australia, the SBT wild-caught fishery began in 1949 (Serventy, 1956), peaking in 1982 at an annual catch of 21,000 tons (Geen and Nayar, 1989); after this peak catch restrictions were imposed on the fishery. Ranching of SBT involves the collection and transfer of wild juvenile tuna into static pontoon enclosures where they are reared and fattened for up to 8 months. Over this period tuna are fed whole baitfish, which results in a near doubling of fish biomass. Ranching began in Port Lincoln Australia in 1990 and in 2017-2018 had an annual production value of AUD$166 million (François et al., 2010;Ellis and Kiessling, 2016;EconSearch, 2019). One of the primary challenges to promoting and ensuring long term success of the industry is to maintain and improve fish health.
When cultured in open ocean pontoons, fish are exposed to free-living microbes, including opportunistic pathogens and infectious parasites (Nowak, 2007). Several parasites have been identified in ranched SBT including Miamiensis avidus (scuticociliate), Caligus sp. (copepod "sea lice"), Cardicola spp. (blood flukes), and Hexostoma thynni (gill fluke) (Nowak et al., 2003). Digenean blood flukes, specifically Cardicola spp., infect all three species of bluefin tuna during ranching which can lead to morbidities and mortalities causing significant losses for producers (Balli et al., 2016;Ellis and Kiessling, 2016). Among the most destructive parasites for SBT are blood flukes that include two species, Cardicola forsteri and C. orientalis (Cribb et al., 2000;Colquitt et al., 2001). In captivity, infections peak approximately 2 months post-relocation to pontoons with SBT having on average 27 flukes per fish (Aiken et al., 2006). Blood flukes produce eggs which are deposited in the gills (up to 9 million eggs per individual) resulting in morbidity and mortality (Colquitt et al., 2001;Bullard and Overstreet, 2002;Shirakashi et al., 2012b). Currently, praziquantel (PZQ) is the only therapeutic used to treat C. forsteri and C. orientalis infections. PZQ is administered orally through feeds by injecting 75-150 mg/kg body weight into sardines which are then fed to SBT (Hardy-Smith et al., 2012). Treatment with PZQ (15 mg/kg body weight) has a positive impact on adults SBT health by eradicating adult flukes and has shown to be effective in PBT juveniles as well (Shirakashi et al., 2012a). Lower doses (3.75-7.5 mg/kg) have also shown positive effects for adult fluke eradication with organ elimination after 24 h (Ishimaru et al., 2013). In addition to therapeutic treatments, pontoon location and alterations to husbandry practices, such as ranching at greater depths, can also lower infection outcomes with a combined impact of decreasing mortalities from 15 to <1% (Kirchhoff et al., 2011).
Marine fish harbor specialized microbiomes within the gill, skin, and gastrointestinal tract that are distinct from the surrounding seawater and influenced by environmental variables (Minich et al., 2019a). These microbial communities play important roles in mucosal barrier defenses to protect against pathogen invasion and maintain host health. Manipulation and optimization of the microbiome of farmed fish through prebiotic, probiotic, or synbiotic interventions are active areas of research to promote fish health and resilience (Gomez et al., 2013;Kelly and Salinas, 2017). Understanding how the fish mucosal microbiome is impacted by the environment and pathological status of the host, e.g., parasite infection load, is largely unknown yet may be important for the development of therapeutics and/or diagnostics (Li et al., 2016). Secondary infections with bacteria or viruses can also be enhanced during parasitic infections but little is known about blood flukes specifically (Kumon et al., 2002;Boxaspen, 2006;Lhorente et al., 2014;Novak et al., 2016). Since blood fluke egg deposition in the gills leads to pathology and respiratory problems, evaluating microbiome changes between infected and noninfected individuals may reveal dysbiosis signatures that further influence host health. Moreover, the presence of a parasite infection itself may lead to changes in the mucosal bacteria communities of the fish which could, in turn, be used to predict, diagnose, and monitor parasite infections. In Atlantic salmon, skin infections by sea lice resulted in decreased richness of skin bacteria which was further driven by an enrichment of potential pathogenic microbes including Vibrio spp., Pseudomonas spp., and Tenacibaculum spp. (Llewellyn et al., 2017). Understanding the mechanisms and interactions between infection, treatment, host response, and the microbiome, the collective "pathobiome" (Bass et al., 2019), are thus important for reducing the impact of infections.
In this study, ranched SBT were sampled to characterize the microbial diversity associated with mucosal body sites, including gill, skin, and gut, providing the first assessment of microbiome diversity in this ecologically and commercially important fish species. PZQ treatment has previously been shown to reduce Cardicola spp. parasite egg abundances in the gills (Power et al., 2019). To better establish the impact of antiparasitic treatment on microbiome composition, we analyzed mucosal microbiome composition as a function of parasite infection and PZQ treatment. We hypothesized that prevalence of Cardicola eggs in the gill would correspond to an altered microbial community. In addition, we hypothesized that among the body sites evaluated, the gastrointestinal tract would have the highest probability of being impacted by the oral treatment of PZQ delivered through feeds. During the course of a harvest event, a total of 65 SBT across six different pontoons from five spatially separate companies were sampled. Assessment of blood fluke prevalence, gill egg densities, and fish biometrics (length and weight) were collected and analyzed with respect to mucosal microbiome composition. Lastly, we compare the mucosal microbiome of SBT to another Scombrid species to evaluate its utility as a possible model organism for marine fish microbiome studies. A meta-analysis of the SBT microbiome compared to two other commercially important fish species suggests the Pacific chub mackerel, Scomber japonicus, may be a suitable model to study microbiome dynamics in marine aquaculture species, including bluefin tuna.

Experimental Animals and Sampling
Southern Bluefin Tuna "SBT, " Thunnus maccoyii, reared in ocean pontoons (also known as ranching) in Port Lincoln, Australia, were opportunistically sampled during a harvest event from July 10 to 20 2018. A total of 65 fish were sampled across six pens from five total companies. Between 9 and 12 fish were sampled from each pen. Fish were measured for total length (mm) and mass (grams) and Fulton's condition factor "Kfactor" was calculated (Fulton, 1904;Froese, 2006). Fish (n = 45) from four of the six pens were treated with the antiparasitic drug, anthelmintic drug praziquantel (PZQ) orally (injected into baitfish) at a dose of 30 mg/kg bodyweight 4-5 weeks after transfer to pontoons in a single treatment over 2 consecutive days by the prescribing veterinarian (APVMA Permit Number-Per 85738) (Power et al., 2019) while fish in the other two pens did not get treated with PZQ (n = 20). Presence of adult flukes in heart and egg densities within the gill were measured per fish as previously described (Aiken et al., 2006;Power et al., 2019). Comparison of the SBT microbiome to other fish samples, including Pacific chub mackerel (MKL) (Minich et al., 2020b) and yellowtail kingfish (YTK), was performed computationally within Qiita . All MKL samples were collected from the wild while the YTK were farmed. Details for fish sample collections are listed in the metaanalysis section below.

Microbiome Assessment
Southern Bluefin Tuna microbiome samples were collected from four body sites (gill, skin, digesta, and anterior kidney) from each of the 65 fish for a total of 260 samples. Gill and skin communities were collected by swabbing approximately 25 cm 2 surface area using wooden shaft cotton swabs. For digesta and anterior kidney samples, approximately 100 mg of fecal matter or tissue was collected with a swab and placed in a 1.5 ml tube. To preserve microbiome integrity, all microbiome samples were stored in 95% ethanol for 4-6 h at room temperature while on the boat and then later transferred to a −20 • C freezer where samples were stored until DNA extraction (approximately 2 weeks later). We chose 95% ethanol over RNAlater because it has a superior performance in preservation of microbiome integrity, meaning there is the least bias when using 95% ethanol (Song et al., 2016). In addition 95% ethanol is readily available globally and has been demonstrated to successfully store samples for both bacteria (Vogtmann et al., 2017;Knight et al., 2018) and viruses while not interfering with the extraction (Minich et al., 2020a).
Molecular methods outlined in the Earth Microbiome Project (earthmicrobiome.org) were used to process samples. Specifically, swabs of either gill, skin, or digesta stored in 95% ethanol were transferred to a 2 ml bead beating tube at the beginning of DNA extraction. Numerous large scale studies have shown that direct extraction of swab heads are an effective way to collect and process microbiome samples (Thompson et al., 2017;McDonald et al., 2018). DNA was extracted in single tubes to avoid well-to-well contamination (Minich et al., 2019b) using the MoBio PowerSoil kit. Multiple titration replicates of 10-fold serial dilutions of positive controls (n = 21) (Escherichia coli) were included to enable detection of background contaminates and determine the limit of detection of the method using the Katharoseq method (Minich et al., 2018b,c). Miniaturized (5 ul) PCR reactions were used to amplify genomic DNA (200 nl) (Minich et al., 2018a) using 16S rRNA V4 515/806 EMP primers Walters et al., 2016). Equal volumes of amplicons (1 ul) were pooled across samples and processed through the Qiagen PCR cleanup kit and then sequenced on a MiSeq 2 × 250 bp run (Caporaso et al., 2012). Sequencing analysis was performed using QIIME 2 (Bolyen et al., 2019, 2) and Qiita .

Statistics and Meta-Analysis
Samples were trimmed to 150 bp and then demultiplexed and processed through the deblur  pipeline to generate unique sOTUs (sub-Operational Taxonomic Units) also referred to as ASVs (Amplified Sequence Variants). We use sOTUs terminology here. Samples were rarified to 1000 reads as determined by Katharoseq cutoff. All downstream analyses were performed using the rarified biom table. Two sOTUs were removed from the dataset. The first was the sample used as a positive control (Enterobacteriaceae: TACGGAGGGTGC AAGCGTTAATCGGAATTACTGGGCGTAAAGCGCACGCAG GCGGTTTGTTAAGTCAGATGTGAAATCCCCGGGCTCAAC CTGGGAACTGCATCTGATACTGGCAAGCTTGAGTCTCGT AGAGGGGGGTAGAATTCCAGG) and the second was a known PCR mastermix contaminant (Eisenhofer et al., 2018) (Pseudomonas veronii: TACAGAGGGTGCAAGCGTTAA TCGGAATTACTGGGCGTAAAGCGCGCGTAGGTGGTTTGT TAAGTTGGATGTGAAATC CCCGGGCTCAACCTGGGAACT GCATTCAAAACTGACTGA CTAGAGTATGGTAGAGGGTGG TGGAATTTCCTG) which we observe disproportionally in low biomass positive controls. The sOTU table was annotated using both Greengenes and the SILVA SSU 1.2.11 database using a minimum identity with query sequence (0.9) and default parameters 10 neighbors per query sequence (Supplementary Table S1; McDonald et al., 2012;Quast et al., 2013;Yilmaz et al., 2014). While both annotations are included in the table  (Supplementary Table S1), the Greengenes annotation is used for subsequent analyses and visualizations. Alpha diversity was calculated using richness (total unique sOTUs), Shannon evenness, and Faith's Phylogenetic diversity (Shannon, 1948;Whittaker, 1972;Faith, 1992). Alpha diversity comparisons were calculated using non-parametric (Kruskal and Wallis, 1952) test with multiple comparisons done using Benjamini-Hochberg correction with a 0.05 FDR cutoff (Benjamini and Hochberg, 1995). Pairwise comparisons were performed using a non-parametric Mann-Whitney test (Mann and Whitney, 1947). Statistical comparisons were performed within PRISM 8.0. Beta diversity was calculated using both weighted and unweighted UniFrac distances (Lozupone and Knight, 2005;Lozupone et al., 2010). To test which metadata categories or variables were associated with beta diversity, multivariate statistical testing was done using ADONIS which is a modified version of PERMANOVA (Anderson, 2001). Specifically, the following categorical variables were assessed for their impact on the microbiome: sample type, company, and PZQ treatment. Following those results, body sites were independently assessed for the impacts of PZQ treatment, company, parasite eggs in gill, parasite flukes in heart, tuna condition factor, tuna length, and tuna mass. Differentially abundance measures were calculated within Calour using the rank mean DS-FDR function to compare two groups with a 0.1 FDR (Jiang et al., 2017;Xu et al., 2019).
A metanalysis was performed to test the hypothesis that phylogenetically related fish species have greater similarity in microbial community composition. To enable this comparison, we used publicly available datasets in Qiita . Both alpha and beta diversity were compared. Note, only datasets which used the Earth Microbiome Project standardized protocols (same extraction and primer set) were used in the comparison as both DNA extraction and primer choice lead to biases in analyses (Brooks et al., 2015;Fouhy et al., 2016). To do this analysis, publicly available microbiome data from gill, skin, and digesta microbiomes of Pacific chub mackerel (MKL), Scomber japonicus, (Qiita ID 11721, ERP2664334) and yellowtail kingfish (YTK), Seriola lalandi, were combined with SBT using the Qiita database. Mackerel samples were wild fish sampled throughout 2017 caught off the Scripps Institution of Oceanography Pier in the Eastern Pacific Ocean, San Diego. The yellowtail kingfish was sampled from pontoons in the Western Pacific Ocean, NSW Australia. Mackerel are phylogenetically closer to SBT in the Scombridae family (mackerels, bonitos, and tunas) (Collette and Nauen, 1983;Collette et al., 2001), whereas yellowtail are within the Carangidae family (jacks and pompanos) (Swart et al., 2015). For specimens collected in this study, YTK are however similar to SBT in terms of trophic level, geographic location of sampling, and farmed rather than wild caught.

Sampling Design and Results
Southern Bluefin Tuna were sampled during annual harvest across 10 days (July 10 2018 to July 20 2018) in Port Lincoln, Australia. A total of six pens were sampled among five different companies (Figure 1a). Two of the pens (n = 20 fish) were not treated with the anthelmintic praziquantel (PZQ) while the other four pens (n = 45 fish) were treated. Four anatomical sites were sampled for microbiome analysis (Figure 1b) including gill, skin, digesta, and anterior kidney. Fish biometrics including fork length (mm), mass (grams), and condition factor "K-factor" were measured for 62 fish whereas data from three fish from July 15th (host_subject_ID: SBT_51, SBT_57, and SBT_58) were not recorded (Figure 1). Tuna ranged in length from 920 to 1180 mm with a median length of 1070 mm (Figure 1c). Tuna ranged in mass from 16200 to 35000 grams with a median mass of 25950 grams (Figure 1d). The K-factor ranged from 1.81 to 2.73 with a median of 2.10 (Figure 1e). Blood fluke counts and parasite egg counts in the gill were recorded for all 65 fish. Blood fluke counts ranged from 0 (30 fish) to 6 with a median and mean intensity of 1 (Figure 1f). Parasite egg counts in the gill ranged from 0 (4 out of 65) to 3.94, with a median of 0.45 and mean intensity of 0.846 eggs per mm filament (Figure 1g).

Microbiome Analysis
A total of 260 host-associated SBT microbiome samples were processed through the pipeline from 65 unique fish across four body sites including gill, skin, digesta, and anterior kidney. In addition, 21 positive controls were included to determine the limit of detection of the pipeline. After applying the Katharoseq formula, the limit of detection of the pipeline, where 50% of the reads of a positive control map to the known control, was 168 reads. When the 90% threshold is applied, as recommended by the Katharoseq method, the limit of detection was 405 reads indicating that any samples with at least 405 reads could be included (Supplementary Figure S1). To be conservative on read depth, we included samples with at least 1000 reads and thus rarified to 1000 reads. At this sequencing depth, a total of 98 samples out of the 260 passed QC ( Table 1). The majority (93.8%) of anterior kidney samples failed and thus were excluded from downstream analyses due to low successful sample size. A final table of 94 samples with 794 unique sOTU features was annotated (Supplementary Table S1). The Good's coverage analysis demonstrated that we were capturing the majority of microbial diversity even with just 1000 reads as the median values for gill, skin, and digesta were 0.997, 0.994, and 0.997, respectively (Supplementary Figure S2). Digesta samples had the highest success rate (38 out of 65) followed by skin (33 out of 65) and gill (23 out of 65).
The microbiomes of the three body sites were compared across the fish. For alpha diversity measures, richness differed across body sites (P = 0.0243, KW = 7.435) with gill and skin communities having a higher richness than digesta (P < 0.05) FIGURE 1 | Sampling design of ranched SBT from Port Lincoln, Australia. (a) Fish were sampled from a total of six pens spanning five companies. (b) Upon harvest, four body sites were sampled and stored in 95% EtOH for the microbiome including the gill, skin, digesta (hindgut), and anterior kidney. Fish biometrics across all samples collected for (c) fork length (mm: median, IQR), (d) mass (grams: median, IQR), and (e) calculated condition factor (median, IQR). Blood fluke parasite counts from (f) heart (flukes per heart: median, IQR), and (g) gill (eggs per mm gill filament: median, IQR). (Figure 2A). Shannon diversity differences across body sites was even stronger (P < 0.0001, KW = 29.97) with both gill and skin having higher evenness than digesta samples (P < 0.0001) ( Figure 2B). However, phylogenetic diversity did not differ across body sites ( Figure 2C). Upon comparing beta diversity, samples were most strongly influenced by body site location for both Unweighted UniFrac ( Figure 2D) (Adonis: P = 0.001, F = 5.731, R 2 = 0.112) and Weighted UniFrac ( Figure 2E) (P = 0.001, F = 13.065, R 2 = 0.223) ( Table 2). Microbes sampled from the SBT mucosal sites represented 20 different phyla, with digesta being enriched in Tenericutes and Spirochaetes relative to other phyla (Figure 2F). At the order level, the Tenericutes were primarily Mycoplasmatales. The Spirochaetes were primarily Brevinematales. Both the Vibrionales and Pseudomonadales were also prevalent across digesta samples within the Gammaproteobacteria (Supplementary Figure S3). Skin and gill microbes were enriched in Proteobacteria, Cyanobacteria, Firmicutes, and Bacteroidetes relative to other phyla ( Figure 2F). Skin and gill samples were primarily enriched in Vibrionales and Pseudomonadales (Supplementary Figure S3). Microbiome variation among body sites was strongest with Weighted UniFrac. Company (geographic location) and PZQ treatment was also significant but less so than body site variation ( Table 2). Subsequent analyses were therefore performed on each body site independently ( Table 2).

Praziquantel Treatment and the Microbiome
The impact of PZQ treatment on the microbiome was evaluated using alpha and beta diversity measures. Treatment with PZQ resulted in digesta samples having lower microbial richness (Mann-Whitney: P = 0.025, U = 98) (Figure 3A), Shannon evenness (Mann-Whitney: P = 0.0328, U = 101) (Figure 3B), and phylogenetic diversity (Mann-Whitney: P = 0.003, U = 76) ( Figure 3C). Although not significant, likely due to low sample size, all measures of microbial diversity in the gill were lower when fish were treated with PZQ. To determine if particular sOTUs were associated with PZQ treatment in the various body sites, we applied the rank mean DS-FDR method for each body site comparing samples which were either treated or not treated with PZQ. From the skin comparison, a sOTU within the Mycoplasmataceae family (Tenericute phyla) was more enriched in skin samples of fish which did not get treated with PZQ ( Figure 3D). For digesta samples, a total of 8 sOTUs were also associated or enriched in fish which were PZQ naïve. This included two Pseudomonas sOTUs, two Acinetobacter sOTUs, Brevundimonas, Delfita, and a sOTU within the Brevinemataceae family ( Figure 3E). All sOTUs classifications were verified with both Greengenes and SILVA databases. The Greengenes Brevinemataceae annotated sOTU was "unclassified" in the SILVA database, with the closest phylogenetic similarity to uncultured Spirochaetes sampled from fish guts, identified from a BLAST search (DQ340184.1, HE586962.1).

Cross Species Microbiome Comparison
To understand the SBT microbiome in relation to other fish species, we performed a metanalysis comparing the SBT mucosal microbiome to a geographically (Australia) and trophically similar (tertiary carnivore) fish species, yellowtail kingfish (YTK), Seriola lalandi (also farmed), and a phylogenetically similar yet geographically and trophically dissimilar fish Pacific chub mackerel (MKL), Scomber japonicus (from the wild). Both MKL and SBT are within the same family, Scombridae (tunas, mackerel, and bonito). The gill, skin, and digesta samples were similarly sampled across all species by the same researcher and further processed using the same molecular methods. Wild S. japonicus were sampled from the Eastern Pacific Ocean in San Diego CA as part of a fish microbiome time series study (Minich et al., 2019a). To verify that 1000 reads was sufficient to interpret and make conclusions from the data, we also compared samples rarified at 5,000 and 10,000 reads. Higher sampling depth resulted in fewer samples being compared, but the trends remain the same indicating that 1000 reads is sufficient for comparisons while enabling the highest number of samples to be compared (Supplementary Figure S4). Body sites were independently compared across all samples for both alpha and beta diversity using richness, Faith's Phylogenetic Diversity, Unweighted UniFrac and Weighted UniFrac (Figure 4).
(P < 0.0001, KW = 68.88, n = 146), even more so than gill microbial diversity (Figures 4A,B). On the skin, for both richness and phylogenetic diversity, a gradient was observed with SBT being lowest followed by MKL and then YTK. Digesta richness and phylogenetic diversity was lower in SBT compared to MKL (P < 0.0001) while MKL was also lower than YTK (P < 0.001) ( Figure 4A). For phylogenetic diversity, SBT was lower than MKL (P < 0.0001) for all body sites and lower than YTK for both gill and skin but not digesta (P < 0.0001) (Figure 4B).
Beta diversity was assessed using Unweighted UniFrac, which gives rare taxa an equal weight, and Weighted UniFrac, which weights taxa based on their relative abundance in the sample. For each unique body site (gill, skin, and digesta), all SBT samples were compared to MKL and YTK samples rarified at 1000 reads. For Unweighted UniFrac, there was no significant difference in microbial diversity in the gill (Figure 4C) whereas for skin, distances to SBT were lower for MKL than YTK (Mann-Whitney P < 0.0001, U = 633804). Digesta distances of YTK and SBT were lower as compared to MKL (Mann-Whitney P < 0.0001, U = 789291) ( Figure 4C). For Weighted UniFrac, distances of all three body sites of MKL to SBT were lower than YTK (gill: Mann-Whitney P < 0.0001, U = 341508, difference = 0.080; skin: Mann-Whitney P < 0.0001, U = 395081, difference = 0.135; digesta: Mann-Whitney P < 0.0001, U = 1704289, difference = 0.063) (Figure 4D).

DISCUSSION
Parasite infections cause significant economic complications in marine aquaculture systems, particularly in high value species such as tuna (Shinn et al., 2015). In this study we set out to describe how the fish mucosal microbiome is associated by parasitic infection and treatment with praziquantel (PZQ) across three body sites including the gill, skin, and digesta of Southern Bluefin Tuna (SBT). Results indicate that the microbiome composition is most explained by body site location, followed by geographic location (company) and lastly PZQ treatment. The microbiomes at each body site were associated with company and PZQ treatment, but were not associated with blood fluke infection (measured as the number of adult blood flukes in the heart, or parasitic egg density in the gills). Lastly, in a meta-analysis, we show that mackerel may have a more similar microbial community to tuna. Since both can be infected by the same parasite and also have a similar microbiome, we suggest the feasibility of using a mackerel as a model for SBT for understanding mucosal microbiome community dynamics. Further experiments would be needed to validate the model. Our findings indicate that husbandry practices, including geographic location, influence the mucosal microbiomes of ranched SBT, and that parasite infection may not have a significant impact on fish microbiomes.

Microbiome Composition of Body Sites
Our overall DNA sequencing success rate was rather low for a typical microbiome study (gill = 35.4%, skin = 50.8%, digesta = 58.5%) and we attribute this to either poor field preservation or that the samples were too low of biomass. In the field, we used 95% ethanol to preserve samples, and although 95% ethanol was shown to best preserve the microbiome community for field collection in human stool samples, it can also inhibit particular steps of molecular assays including DNA extraction and PCR, and we would instead recommend further optimization of this preservation method in the future for fish or instead use cryopreservation, e.g., dry ice (Kemp et al., 2006;Demeke and Jenkins, 2010;Song et al., 2016). A second explanation is that the samples were low biomass meaning there was very few microbial cells. Since we used column cleanups, it's likely that our limit of detection was between 50,000 and 100,000 cells, thus if a sample had less than this it would potentially fail sequencing. To improve this, we recommend using magnetic bead cleanup methods in the future (Minich et al., 2018b). Since 95% ethanol has been shown to work successfully for storage of RNA viral samples (Minich et al., 2020a), its most likely that sample dropout in this particular study was likely due to low biomass of the samples.
Sample type was the strongest predictor of the microbiome in the dataset. Mucosal environments of fish, including the gill, skin, and gut, are inhabited by commensal microbial communities which can have both negative and positive impacts on fish health (Gomez et al., 2013;Beck and Peatman, 2015). Recent advances in genomics methods including transcriptomics, microbiome, and proteomics have enabled the characterization and monitoring of these communities in relation to host response and environmental changes (Salinas and Magadán, 2017). Although the anterior kidney is an important organ for immune function (Abelli et al., 1994;Fänge, 1994;Watts et al., 2003;Geven and Klaren, 2017), only four of the 65 samples had detectable microbial DNA. We would expect anterior kidneys to generally lack a microbial community except in cases of systemic infection. In comparison, gill, skin, and digesta samples had a much higher success rate, indicating a rich microbial community. While various studies have focused on fish body site microbiomes independently (primarily gut followed by skin and gill), few have evaluated the cumulative microbiome across multiple body sites for individual fish (Ghanbari et al., 2015;Larsen et al., 2015;Pratte et al., 2018;Minich et al., 2019aMinich et al., , 2020c. Our study highlights how amongst SBT body sites, gills have a high microbial diversity which was not as pronounced in the other fish species (relative to body sites within that species). Gill microbial communities in fish are understudied but have been shown to be influenced by host diet, (Pratte et al., 2018) and environmental conditions such as suspended sediment (Hess et al., 2015). Gill and skin communities, however, were generally more similar and stable as compared to the gut communities in SBT. This stability has also been shown in Pacific chub mackerel, Scomber japonicus (Minich et al., 2019a). Gill and skin communities were primarily enriched in Proteobacteria, while gut communities had higher proportions of Mycoplasmatales within the Tenericutes phylum and Brevinematales within the Spirochaetes phylum. Gill and skin samples were enriched in various sOTUs within Pseudomonadales and Vibrionales both of which can be a source of potential pathogens and thus important to study further for health monitoring (Williams et al., 1974;Austin, 2005). Mycoplasma are a common marine microbe associated with fish guts (Llewellyn et al., 2016). While many host associated genera within Spirochaetes are intestinal pathogens including Treponema, Borrelia, Leptospira, and Brachyspira, little is known about the order Brevinematales (Bellgard et al., 2009;Gupta et al., 2013). Fish mucosal environments including the gill, skin, and digesta can also be influenced by the rearing condition and surrounding environment including the water (Minich et al., 2020c).

Influence of Location or Company on Mucosal Microbiome
When analyzing body sites independently, skin and digesta samples were significantly differentiated by the location of rearing (company). While it is probable that the five companies sampled utilized varying animal husbandry approaches (including diet) which could influence the host-associated microbial communities, it is equally likely that the pontoon location had an effect. Depth of ocean floor has been shown to positively reduce parasite infection rates in SBT as the sediment and intermediate host, the polychaete Longicarpus modestus, is further away from the fish and current velocities are higher (Cribb et al., 2011;Kirchhoff et al., 2011). Sediment is known to impact microbial communities of the water column, so it is likely this would additionally impact the fish microbiome (Hess et al., 2015). Another possible explanation for a company effect would be that biofouling on farming structures could differ based on location in the ocean and net cleaning frequency. Although it is not known if elevated biofouling leads to increased microbial diversity in fish, it is possible that exposure to different microbes could be greater. During pontoon aquaculture, many benthic organisms at the larval stage will settle on the pen infrastructure, utilizing nutrients from fish feces and excess feed. Many of these organisms are common benthic invertebrates and have been shown to negatively impact aquaculture by consuming oxygen, reducing water flow in the pens, and harboring microbes which can be pathogenic (Cronin et al., 1999;Fitridge et al., 2012;Madin and Ching, 2015). Biofouling communities may change depending on geospatial location of the pontoons which may experience varying oceanographic conditions including currents, wind, and upwelling.
Our study identified Pseudoalteromonas, Psychrobacter, and Vibrio as representing the most abundant genera across the body sites. In culture based studies, Vibrio, Photobacterium, Pseudoaltermonas, Tenacibaculum, and Flavobacteria have all been isolated from the gills of SBT (Valdenegro-Vega et al., 2013). Prominent SBT skin microbial genera identified included Pseudoalteromonas, Vibrio, Acinetobacter, Psychrobacter, and Mycoplasmataceae. Skin microbes including Vibrio, Clostridium, Enterobacter, Klebsiella, and Proteus isolated from other tunas have previously been shown to be important for causing decomposition through histamine production (Yoshinaga and Frank, 1982).

Praziquantel Impacts on the SBT Microbiome
Treatment of PZQ had a moderate effect on microbial communities in SBT. When comparing alpha diversity, digesta samples decreased in richness, Shannon, and Faith's phylogenetic diversity in fish treated with PZQ. When comparing beta diversity, Unweighted UniFrac distances in gill and digesta samples were associated with PZQ treatment which may suggest that rare taxa were more influenced by PZQ. For Weighted UniFrac, only gill and skin samples were influenced, indicating that more abundant microbial taxa were impacted by PZQ treatment. Children successfully cured of schistosomiasis with PZQ had higher Fusobacterium bacteria abundances in the gut microbiome compared to infected children who were not successfully cured (Schneeberger et al., 2018). The associations were not strong compared to other factors, but suggests that PZQ has an impact on the human gut microbiome which could be extrapolated to fish. A caveat of this study is that only one company did not use PZQ thus the effect seen could also be confounded by husbandry differences between these companies. It is also possible that microbes in the gut may metabolize PZQ rendering it less effective or even toxic to the host .

Parasite -Microbiome Interactome
To the best of our knowledge, this is the first study to investigate the impact of blood fluke infection on the mucosal microbiome of a marine fish. However, no significant associations were observed between the microbiome and blood fluke prevalence or intensity in SBT. Hypothetically, parasitic infections can influence the microbial communities of their hosts either directly through grazing or indirectly by promoting secondary infections post trauma. In Atlantic salmon, parasitic copepod sea lice, Lepeophtheirus salmonis, infections were associated with decreased bacterial richness and dysbiosis of the skin microbiome (Llewellyn et al., 2017). In freshwater barramundi farms, parasitic ciliate abundance was associated with changes in the gill microbiome and mortalities (Bastos Gomes et al., 2019). Despite a relatively low incidence of infection in this study, the results suggest that blood fluke infection does not impact the microbiome. Initial studies on blood fluke densities in ranched SBT show that infected fish have on average 27 flukes per fish (Aiken et al., 2006). The median in this study was one blood fluke per fish with a maximum of six, therefore the infection intensity was low compared to previous reports.

Modeling Future SBT Mucosal Microbiome and Parasitome Research
Our study demonstrated that Pacific chub mackerel (MKL), Scomber japonicus, has a similar gill, skin, and digesta mucosal microbiome to SBT and thus is a putative candidate model organism to study the parasitome and mucosal microbiome. Previous studies have shown that platyhelminth parasites are distributed across multiple tuna hosts (Aiken et al., 2007). Further, Caridcola infections occur in other Scombrid species along with other fish families (Nolan et al., 2014). Scomber japonicus has been shown to be successfully treated for skin fluke infections with PZQ (Yamamoto et al., 2011). For microbiome comparisons, digesta samples analyzed by Unweighted UniFrac and alpha measures showed yellowtail kingfish (YTK) was more similar to SBT. Since YTK and SBT are both tertiary carnivores and, in this study, both are farmed, their diet is generally more rich in fish protein than MKL which are small pelagic secondary consumers (Agusa et al., 2007;Hajeb et al., 2010). This association of diet and high trophic level corresponding with low diversity was first observed in mammals, and may also be conserved to some extent in fish (Ley et al., 2008a,b). However, when relative abundance (compositional) measurements are included in the analysis, all MKL body sites were more similar to SBT, which could suggest evidence of phylosymbiosis for highly abundant taxa. For microbiome studies, a link has been demonstrated within plants and animals to suggest that phylogenetically related hosts retain more similar microbiomes (Pollock et al., 2018;Ross et al., 2019;Lim and Bordenstein, 2020). While most of these studies have focused on the gut microbiome, additional body sites have recapitulated this relationship in the skin (Chiarello et al., 2018) but not gills for fish and amphibians. Our study demonstrates how gill, skin, and gut communities of the SBT were more similar to a genetically similar fish within the same family.
Since the MKL microbiome is similar to SBT, it provides an opportunity to explore developing MKL as a model for future SBT research. One of the challenges with marine finfish aquaculture is performing replicated experiments, including disease challenges, that are not feasible to do in a production setting (Salama and Rabe, 2013). Mackerel is an ideal species for aquaculture research as they are globally distributed, not threatened, easy to catch, inexpensive, small, and easy to culture in tanks.

CONCLUSION
Our findings demonstrate that Southern Bluefin Tuna harbor a unique microbiome in each body site. The SBT microbiome is influenced by pontoon location and, to a lesser extent, treatment with the antihelminth therapeutic, praziquantel. No association between blood fluke infection and the microbiome was identified. It is possible, however, that infection intensity was not sufficient or sample size not adequate to identify a relationship. In addition, we showed genetically similar fish (Scombrids) have a more similar microbial community across multiple body sites including gill, skin, and digesta, suggesting possible phylosymbiosis. Based on this finding, the exploration of Pacific chub mackerel as a candidate model organism for studying the microbiome and potential parasitome of SBT is proposed. Tuna are highly valued, and present significant challenges for conducting controlled experiments due to their size and physiology. Future work employing a smaller, highly accessible relative could enable greater gains in research and ultimately, enhanced aquaculture productivity.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories (Qiita ID 11721, ERP2664334). The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

ETHICS STATEMENT
Ethical review and approval was not required for the animal study because all animals were harvested and euthanized as part of the standard harvest by aquaculture company personnel in accordance with standards established by the Australian Southern Bluefin Tuna Association. Microbiome samples were opportunistically obtained from individual specimens post-euthanization.

AUTHOR CONTRIBUTIONS
JM, BN, NB, and EA contributed to the conception and design of the study. JM collected microbiome samples in the field, performed DNA extractions and microbiome processing, analyzed the data, and wrote the first draft of the manuscript. CP performed the analyses on gill egg density and heart fluke counts. MM contributed to analysis and visualization. All authors contributed to the manuscript revision, and read and approved the submitted version.

ACKNOWLEDGMENTS
We thank the University of Tasmania for hosting JM and Royal Melbourne Institute of Technology (RMIT) for enabling sample processing. The Australian Southern Bluefin Tuna Industry Association are acknowledged for hosting the team in Port Lincoln, providing support for the collection of samples, including visits to commercial SBT ranching pontoons.