Digital Droplet PCR-Based Environmental DNA Tool for Monitoring Cryptocaryon irritans in a Marine Fish Farm from Hong Kong

: The adoption of new investigative strategies based on environmental DNA (eDNA) can be used to monitor parasites, associated bacterial microbiomes, and physical-chemical parameters in ﬁsh farms. In this study, we used the economically important and globally distributed ﬁsh ciliate parasite Cryptocaryon irritans as a model to understand the parasite abundance and potential drivers of its presence in marine ﬁsh farms. Environmental (rainfall) and physical-chemical (temperature, oxygen, salinity, pH) data collected from a marine ﬁsh farm in Hong Kong were analyzed together with the eDNA approach targeting C. irritans abundance based on digital droplet PCR and 16S metagenomics to determine associations and triggers between parasites and speciﬁc bacterial groups. Rainfall and temperature demonstrated positive associations with high abundance of C. irritans (eDNA) at the studied marine ﬁsh cage farm. However, rainfall was the only parameter tested that demonstrated a signiﬁcant association with parasite eDNA, indicating that the raining season is a risky period for ﬁsh farmers in Hong Kong. Coraliomargarita was the bacterial genus with the most signiﬁcant relationship with low abundance of C. irritans in water. Understanding the environmental triggers of ciliate parasites propagation and associated bacterial microbiome could elucidate new insights into environmental control, microbial management, and promote the reduction of chemical use in marine ﬁsh farms.


Introduction
The aquatic environment in marine fish farms is complex and interactive. These farms are exposed to multiple aquatic microbes which propagate and co-exist with fish [1][2][3]. Recent studies have started to unveil close relationships between the bacterial microbiomes pathogens and parasites in fish farms before animals show disease symptoms [28]. As disease outbreaks are not solely caused by the presence of a pathogen/parasite, but by its combination with the susceptibility of the animals (health status, microbiome community), and the quality of their environment (physical-chemical variables), it is possible to use quantification of pathogens/parasites eDNA to predict diseases outbreaks in aquaculture before sick animals present themselves.
Droplet digital PCR (ddPCR) is a novel tool to quantify the absolute number of specific DNA sequences at very low concentrations [25,29,30]. Replicates of conventional qPCR cannot resolve the "concentration plateau" problem due to the sample variability at low DNA concentrations, commonly found in eDNA studies, yet ddPCR circumvents these limitations due to the increased sensitivity and dilution of PCR inhibitors through partitioning a single sample into thousands of droplets which run in parallel [26,30]. The nanodroplet PCR reactions are independent, single amplification events. The ddPCR technology has been shown to have higher accuracy in abundance estimation, producing more precise, reproducible, and statistically significant results [31,32]. Recently, ddPCR has also been applied to quantify water and foodborne noroviruses concentrations from water samples [26]. In the case of C. irritans infections, only the trophont or "feeding" stage is the most observable on fish [13]. However, the characteristic white spots may not always be obvious to farmers before mass mortalities are observed [33]. Hence, eDNA method coupled with ddPCR could potentially allow for accurate and early quantification of parasite loads in farming systems without its visible presence.
In the present study, a ddPCR assay for the detection C. irritans from water samples (eDNA) was developed to determine its absolute abundance in a marine fish farm from Hong Kong. Additionally, in order to understand the potential drivers behind the C. irritans abundance in the environment, we assessed its relationship with physical-chemical and environmental parameters, as well as the associated bacterial microbiome. The application of this novel approach demonstrates its potential as a promising environmental monitoring and management tool for mitigating risks of disease outbreaks in commercial fish farms. This important information can also be adopted by environmental and governmental agencies for policy development and regulation of new aquaculture sites and environmental protection.

Farm Selection, Fish Species, and Parasite Identification
A marine fish farm located at Three Fathoms Cove, Hong Kong reporting yearly mortalities with C. irritans was adopted as the location for this study. Cages contained pompano (Trachinotus blochii), Chinese bahaba (Bahaba taipingensis), and rabbitfish (Siganus fuscescens). T. blochii, which was more susceptible to C. irritans infections, were transferred across the net-pen cages, while B. taipingensis was added during the sampling period. This is a common practice adopted by farmers in Hong Kong and China to help reduce infections (moving fish around may have a "dilution" factor on infections) caused by parasites and other pathogens [34]. In this site, severe mortalities of juvenile T. blochii caused by C. irritans were reported by the farmer and diagnostics confirmed through the ambulatory of veterinary services of the Jockey Club College of Veterinary Medicine and Life Sciences from City University of Hong Kong (Figure 1a).

Water Sample Collection
Water samples were collected weekly between January and August 2020 from five net-pen cages located on the commercial fish farm described in 2.1 (Figure 1b,c). One liter of water sample was collected from each of the five sea cages approximately 50 cm below the water surface using plastic bottles which were previously sterilized with 10% bleach. To ensure there was no cross-contamination among water samples during sampling and DNA extraction, negative controls (field blank, Figure 1d) were deployed. In each sampling event, five field blanks consisted of the same plastic bottle used in the field sampling with one liter of tap water each and briefly submerged into the farm cages. The field blanks were filtered the same way as the water sample and included in each extraction batch. All samples were stored on ice and transported immediately for processing at the laboratory for Infectious Diseases and Public Health Department at City University of Hong.

Environmental and Physical-Chemical Parameters Measurement
A series of physical-chemical and environmental parameters were collected by the In-Situ's Aqua TROLL 500 Multiparameter Sondes and HydroMace data loggers (In-Situ Inc., Waltham, MA, USA), together with rain gauges installed at the farm study site. Temperature, dissolved oxygen, salinity, turbidity, pH, ORP (oxidizing or reducing potential), and rain levels were collected in real-time via the In-Situ's data services platform every five minutes during the entire 7-month study period. Data were transmitted from the farm site via telemetry to HydroVu and was accessed on-line in real time.

Environmental DNA (eDNA) Filtration, Extraction, and Parasite Assay
Each one liter water sample was pumped through a 3 µm, 47 mm diameter cellulose nitrate membrane (GE Whatman) using a portable vacuum pump system (Grover Scientific, Rosslea, Australia). Filter membranes were stored in a 1.5 mL Eppendorf tube at −20 • C for DNA extraction. All equipment (vacuum pump, forceps and extension pipes) were cleaned using 10% bleach for three minutes and then rinsed with tap water for five minutes between each filtration to prevent cross-contamination. The genetic material from half of each filter membrane was extracted using a CTAB (cetyl trimethyl ammonium bromide) DNA extraction protocol described in Bastos Gomes et al. [35] and resuspended in 50 µL of 1X TE buffer. The DNA quality was checked by a 0.8% agarose gel and the DNA was quantified with NanoDrop spectrometer (Thermofisher Scientific, Waltham, MA, USA).
A specific C. irritans TaqMan assay (Thermofisher Scientific, USA; Table 1) was developed based on the SSU-rDNA sequence of C. irritans obtained from infected Asian Seabass (Lates calcarifer) from an Australian marine farm (data not published) to quantify the absolute copy numbers in water samples by ddPCR.

Droplet Digital PCR (ddPCR) Assay
The working area (bench and fume hood) was sterilized with 10% bleach and 70% ethanol before every assay set up. All equipment (pipettes, filter tips, ultra-pure water, tube racks, tubes, well plates) were UV-sterilized for 30 min prior to use. DNA extracts were diluted ten-fold with ultra-pure distilled water (Invitrogen, USA). ddPCR was performed on a QX200 Droplet Digital PCR (Bio-Rad Laboratories, Inc., Hercules, CA, USA) according to the manufacturer's instructions. In summary, ddPCR reactions consisted of 10 µL 2X ddPCR Supermix for Probes (Bio-Rad Laboratories, Inc., USA), 1 µL of TaqMan assay (TaqMan, Thermofisher Scientific, USA), and 9 µL of diluted DNA sample. For ddPCR negative controls, the 9 µL DNA samples were replaced by the same volume of ultra-pure distilled water. The positive control corresponded to 5.85 × 10 13 copies/µL (synthetic sequence of the SSU-rDNA gene from C. irritans based on isolates from infected barramundi L. calcarifer (Supplementary Table S1)) was diluted six times using a 100-fold serial dilution to achieve 58.5 copies/µL. Sixty-two field blanks (around 25% (62 samples) of all extracted samples) were randomly selected for ddPCR to detect any cross-contamination during sample processing. An 8-well cartridge was loaded with 20 µL of the reaction mixture and 70 µL of oil for droplet generation, in which 40 µL of droplet emulsion was pipetted into a 96-well plate. The plate was sealed by a pierceable foil heat seal, and then amplified in C1000 Touch Thermocycler (Bio-Rad Laboratories, Inc., USA) using the following conditions: activation at 95 • C for 10 min, 40 cycles of 94 • C for 30 s and 60 • C for 1 min, followed by enzyme deactivation at 98 • C for 10 min and held at 4 • C for 30 min. The droplets were read by a QX200 droplet reader and the concentration (absolute count of copies/µL) was analyzed by the QuantaSoft Analysis Pro Software (Bio-Rad Laboratories, Inc., USA). Droplets were considered positive when the background threshold limit (the cut-off value of positive droplets) was above 2000 (fluorescence unit). The quality control of reactions was defined as >16,000 droplets in each well and <5 positive droplets in the negative controls. The copy number of DNA was assumed to have a Poisson distribution; hence the lower limit of detection (LLOD) was interpreted by at least five positive droplets and the upper limit of detection (ULOD) by at least five negative droplets [30,36]. The total copy number presented in the one liter water sample (x) was determined by the given data of absolute count (in target copies/µL), then calculated relative to the total volume of ddPCR reaction (20 µL) and the starting volume of samples before DNA extraction from half of a filter membrane (two of 0.5 L). The formula was carried out as:

Library Preparation for Bacterial Microbiome Analysis
The 16S metabarcoding library preparation was carried out at Groken Bioscience Ltd. (Hong Kong). The V3 and V4 regions of the 16S gene were amplified using the 16S Amplicon PCR primers (341F (5 -CCTAYGGGRBGCASCAG-3 ) and 806R (5 -GGACTACNNGGG-TATCTAAT-3 )). Amplicon PCR reactions were carried out with Phusion High-Fidelity PCR Master Mix (Bio-Rad Laboratories, Inc., USA). The reaction volumes and cycling conditions followed the manufacturer's protocol. Amplicons of 450-550 bp were selected by a 2% agarose gel and purified with Qiagen Gel Extraction Kit (Qiagen, Germany). The sequencing libraries were generated using NEBNext Ultra DNA Library Prep Kit for Illumina, following manufacturer's protocols and index codes were attached. The library quality was assessed on the Qubit 2.0 Fluorometer (Thermo Scientific, Waltham, MA, USA) and Agilent Bioanalyzer 2100 system. The library was then normalized, pooled, and denatured to be sequenced on a NovoSeq Illumina platform, and 250 bp paired-end reads were generated.

Statistical Analysis and Bioinformatics
The association between the observed variables, and the abundance of C. irritans in water, was analyzed by a principal component analysis (PCA). Abundance of C. irritans eDNA data was Ln transformed to conform to a normal distribution prior to analysis (Kolmogorov-Smirnov test, p > 0.05). The original set of five variables of interest were reduced into two principal components for which eigenvalues were >1 and represented as a two-dimensional plot. To explore the strength, direction and significance of relationships between variables, a correlation matrix based on Spearman rho (r) bivariate analyses was generated.
A heatmap was also generated to explore the relationship between different genera with high and low C. irritans abundance. Genera were normalized with DESeq2 s variance stabilizing transformation, and differential abundance testing performed using a likelihood ratio test on the C. irritans abundance (high vs. low). The transformed abundances for taxa significantly affected by C. irritans were then used to create a matrix with mean abundances for high and low C. irritans abundance, with the matrix then applied to the heatmap function with clustering by rows (taxonomy).
Both pre-analytical bioinformatics and statistical analysis were conducted using R Studio Version 3.6.1 [37] with a pipeline adapted from Workflow for Microbiome Data Analysis: from raw reads to community analyses [38,39]. DADA2 [40] was used for quality filtering and trimming, demultiplexing, denoising, and taxonomic assignment (using the SILVA Database). The microDecon package [41] was used to remove homogenous contamination from samples using blanks. For statistical analysis, a phyloseq object was created using the package Phyloseq [42] with taxa filtered by prevalence (threshold = 0.01) and agglomerated at the genus level.
A generalized logistic regression model (binomial) was used to explore the effect of both the microbiome and environmental parameters on C. irritans abundance (eDNA) (high vs. low) using the package lme4 [43]. Negative-binomial modelling using DESeq2 was used, with its built-in variance stabilizing transformation and likelihood ratio test, to identify genera associated with parasite abundance (eDNA) and that could subsequently be included in the regression model [44]. Both the resulting DESeq2-normalised counts for physical-chemical and environmental parameters were centered and scaled to avoid convergence issues, and multicollinearity assessed, with collinear variables subsequently removed. The resulting fixed factors: dissolved oxygen (mg/L), rain (mm), salinity (PSU), NS4 marine group (genus), Coraliomargarita (genus), and Arcobacter (genus), as well as the random factor date were included in the model C. irritans (eDNA)~(6 Parameters) + (1|Date).
Backwards selection [43][44][45][46] was then implemented to simplify the model by comparing Akaike's Information Criterion (AIC) scores between regression models and removing predictors that were not contributing to the model. The final model was Cryptocaryon (eDNA)~rain + dissolved oxygen + NS4 marine group + Coraliomargarita + (1|Date). The significance of the fixed effects variables in this final model was then assessed using analysis of deviance (Type II Wald Chi-square test) from the car package [45].

Results
Two main components of a principal component analysis (PCA) explained 71% of the variances observed among physical-chemical and environmental parameters monitored on farm (Figure 2). A positive association was observed between dissolved oxygen and pH (Figure 2; quadrant I). A positive association was also found among abundance of C. irritans in water (eDNA), rainfall, and water temperature (Figure 2; quadrant II), which were inversely related to salinity (Figure 2; quadrant III), demonstrating that levels of C. irritans were higher during more intensive rain and when temperatures were higher in Hong Kong (summer period in the north hemisphere). Cryptocaryon irritans loads increased between March and May 2020. Figure 3 summarizes the seasonality of C. irritans abundance, physical-chemical and the environmental parameters (dissolved oxygen, water temperature, rainfall and salinity) during the seven months period of water collection. C. irritans abundance (eDNA) followed mainly the trends of increasing temperature and rainfall and decreased dissolved oxygen levels. C. irritans loads were particularly high in May 2020, when temperatures surpassed 25 • C, dissolved oxygen readings were consistently dropped below 6 mg/L and rainfall was intensifying in Hong Kong (Figure 3).
The association among the absolute abundance of C. irritans (eDNA), water quality, environmental parameters, and microbial composition in sea cages from Hong Kong were assessed between January and August 2020. Table 2 summarizes the output for the generalized logistic regression model (GLM) exploring the effects of all biotic (bacterial microbiome) and abiotic variables (physical-chemical and environmental parameters) on high vs. low C. irritans load (eDNA) (based on mean values from seven months of sampling) in water samples, whereas Figure 4 represents the mean (±SD) values of each predictor relative to parasite load.  Among the physical-chemical and environmental parameters tested on the farm, rainfall was the only parameter which demonstrated a significant positive association with the abundance of C. irritans (eDNA) ( Figure 4D), which was also the highest when the abundance of the parasite was high in water. In regards to bacteria identified in the study, the association between high and low parasite load was only significant for Coraliomargarita, for which the relative abundance was high when C. irritans eDNA was low ( Figure 4A).
The top twenty most abundant bacterial genus for both high and low C. irritans abundance during the seven months period can be observed on the heatmap of Supplementary Figure S1. The relationship between the most significant bacterial genus and high and low abundance of C. irritans (eDNA) are represented in the heatmap of Figure 5. As corroborated by the GLM, Coraliomargarita genus had a relative low abundance when C. irritans was high in the environment during the sampling period of this study.

Discussion
The present study adopted an innovative approach to monitor the marine parasitic protozoan Cryptocaryon irritans based on eDNA-ddPCR approach in a sea cage fish farm from Hong Kong between January and August 2020. Additionally, we coupled parasite abundance with the bacterial microbiome (eDNA-NGS of 16S gene), physical-chemical and environmental parameters (using sensors data in 5-min intervals) recorded in the same location to understand the drivers of parasite propagation. Cryptocaryonosis has been reported as one of the common parasitic diseases in Hong Kong fish farms [47]. Unfortunately, the participant farm in our study did not monitor fish parasite prevalence or abundance as part of their management practices, which limited our ability to link our physical-chemical and environmental parameters, bacterial microbiome, and C. irritans abundance in water (eDNA) with fish infections or mortalities. The fact that the farm physical-chemical parameters (e.g., dissolved oxygen, temperature) were constantly oscillating along the seven months of sampling, mortalities triggered by these changes together with high loads of C. irritans in water could have happened and not being recorded by farmer workers. The lack of fish health monitoring and data records in fish farms seems to be a common trend observed in small scale fish farms of Hong Kong, and South East Asia in general. This fact reinforces the need for the development of indirect ways to identify disease risk (early warning systems) in Asian aquaculture farms.
The difficulty in managing parasite outbreaks in open sea cage fish farms [10,20] could potentially be mitigated by targeted eDNA monitoring strategies associated with physical-chemical and environmental parameters, such as those employed in the present study. In fact, researchers who looked at the ciliate diversity in relation to environmental factors showed that abiotic parameters including water temperature, salinity, pH, and dissolved oxygen are key factors affecting ciliate richness, evenness, and community composition [48]. A recent study of Atlantic salmon (Salmo salar) has demonstrated the role of abiotic parameters in pathogen eDNA collected from seawater similar to our study [49]. Our proposed model (which we acknowledge to be in its exploratory stages) could provide a framework of supportive methods to better understand the complex intricacy between biotic and abiotic factors influencing ciliate dynamics in marine environments. Such framework has strong potential to support the development of warning systems to assist farmers to potentially predict the risk parasites (such as C. irritans) outbreaks, prevent future fish losses and environmental impacts.
Another important aspect of adopting environmental DNA-based monitoring in fish farms is its non-destructive sampling approach [28,50,51]. Traditionally parasites burden is monitored in fish farms by handling animals, use of anesthetics, and preventive bath treatment, which are stressful for fish and can predispose animals getting sick after any of these procedures [52]. Using eDNA-based monitoring tool could bring benefits for farmers with less animal handling, such as lower mortalities due to stress (related to handling). This less aggressive approach could reduce the cost of sampling [53], reduce the use of chemicals, improve animal welfare [54,55], drive positive changes of the bad public perception towards fish farms (related to animal welfare), and support farms with sustainability certification.
In this study a significant relationship between the C. irritans abundance and rainfall (p < 0.05) was found. Drastic changes in salinities can certainly affect parasite propagation and have a direct association with rainfall. During heavy rainfall periods the salinity measurement in fish farms could decrease to 15 or 10 ppt, promoting cytolysis of trophonts and tomonts parasite stages which can significantly cause reduce outbreaks of cryptocaryonosis [17,20,21]. However, the observed salinity changes in the present study were potentially not sufficient to cause a disruption in C. irritans development and propagation. Although salinity was inversely related with rainfall, salinity of seawater in the farm area decreased only about 3 ppt overall (i.e., from 34.5 to 31.5 ppt) in May/June periods of heavier rainfall and higher C. irritans abundance. One factor which may have contributed to higher ciliate loads during this rainy period can be explained by a disturbance of the aquatic environment and the sea bed (only five meters deep in the studied farm) through increased nutrient input from surrounding hills and resuspension of organic sediments, which could predispose the propagation of ciliates such as C. irritans. Another potential factor for higher ciliate loads identified in periods of heavier rainfall would be directly related to immunity of fish (hosts) farmed within the sampled area. A previous study on the European sea bass (Dicentrarchus labrax) and sea bream (Sparus auratus) reported that fish are weakened during physiological stresses induced by intense precipitation events and become more susceptible to parasites, predisposing fish mortalities [22].
The present study recorded a steady increase in water temperature from March (~20 • C) to June (~28 • C), which coincides with periods of higher levels of C. irritans in water. While in our model temperature was not the most significant parameter to be associated with parasite abundance (perhaps due to only seven months of water sampling), temperature is an important general abiotic factor to consider for the development of any C. irritans warning system. The development of trophonts and the division/release of the infective tomites are promoted by warmer temperatures, with optimal growth at 23-30 • C [17,21]; whereas C. irritans can remain dormant in water below 15 • C [19]. Moreover, there is evidence that prevalence of C. irritans can vary widely with seasons as it is temperature dependent (24~27 • C) and its life cycle can be accelerated at high temperatures (especially over 27 • C) [14,56]. This fact is particularly important considering the predicted changes in global environmental conditions due to climate change. Increase in water temperatures has caused serious cryptocaryonosis in southern China from spring to autumn [57], which is in agreement with what has been observed through the eDNA approach in the present study.
In relation to the associations observed between bacterial communities, Coraliomargarita genus was significantly associated to the low abundance of C. irritans in water. Based only on our data set, it is difficult to determine the exact reasons for this association as this is the first time that this bacterial genus has been reported in association with ciliate parasites from marine environments. Some coral reef fish seems to harbor Coraliomargarita during larval stage in their gastrointestinal tract, indicating a specific niche for this genus [58,59]. This genus could be a common bacterial group found in the waters around the studied farm (shallow area, five meters deep) considering its high organic content and the fact that the farm location has experienced previous algal bloom events (personal communication, farm owner). Uncultured Coraliomargarita genus apparently grows preferentially on high organic compounds which are important indicators of increased phytoplankton concentration and algal blooms, for which this genus was reported to account for 1.15-5.51% of the bacterial communities in marine cage farms [60,61]. The association of Coraliomargarita with algal blooms which are usually triggered by physical-chemical and environmental parameters should be further studied as it may clarify any link of this bacterial group and the abundance of C. irritans in marine environments. Furthermore, future studies should integrate data on nutrient concentration in water, turbidity, chlorophyll-a concentration, or even some other indicators of organic matter to clarify the potential link between this bacterial genus and parasites.
Among the most significant bacterial genera associated with the high and low abundance of C. irritans, Arcobacter and NS4 marine group were the ones associated with high levels of C. irritans in water. Higher C. irritans was also associated with relative high levels of NS4 marine group in water. NS4 marine group is another indicator of increased phytoplankton concentration and algal blooms aforementioned [61]. The intensive aquaculture in coastal marine cages caused the dominance of NS4 marine group in Vietnam due to its strong correlations with inorganic nutrients such as ammonium and phosphate (good predictor of algae proliferation) [62]. The Arcobacter genus on the other hand, similarly to Coraliomargarita, is a common bacterial group found in coral reef fish larvae bodies [59]. This genus is one of the food-borne and water-borne human pathogens associated with aquaculture animals representing a serious concern for public health worldwide [63]. Previous studies have showed pathological infections of rainbow trout (Oncorhynchus mykiss) and zebrafish (Danio rerio) with Arcobacter species [64,65]. Given that information related to the functional roles of these two bacterial genera is limited [66], our study may help to expand the still limited knowledge base required to implement microbiome management with the goal of predicting the risk of outbreaks in fish farms

Conclusions
This is the first study to adopt an eDNA-ddPCR-based method to monitor abundance of C. irritans in marine waters of Hong Kong. Digital droplet PCR has demonstrated to be an important new tool in clinical microbiology to quantify pathogens [67]. Assays based on ddPCR technology can be more sensitive as a screening tool of species-specific parasites in marine environments, since the end-point quantification method is less affected by inhibitors in water samples compared to conventional quantitative PCR [68][69][70]. Adopting eDNA quantification based on a ddPCR assays in combination with models of physical-chemical and environmental parameters together with the marine bacterial composition could provide important insights about aquatic animal health and microbial management for the aquaculture industry. While this study reports the first exploratory insights into the relationship between a parasitic ciliate, bacterial communities, and physicalchemical/environmental data analysis from a marine fish farm, this information could be further explored as the basis for the development of environmental policies by government agencies. Furthermore, the proposed model could support companies or governmental agencies on mapping out new aquaculture sites with a focus on environmental protection and sustainability.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/d13080350/s1. Figure S1: Heatmap clustered by rows showing the DESeq2 normalized abundance of the top twenty most abundant bacterial genus for both high and low (based on mean values) C. irritans abundance (eDNA) between January and August 2020 from a marine sea cage fish farm from Hong Kong. Table S1: Synthetic sequence of the positive control targeting 18S SSU-rDNA gene from C. irritans for quantitative PCR and digital PCR.