All zebrafish (Danio rerio) were maintained in accordance with protocols approved by the Institutional Animal Care and Use Committee at Mayo Clinic. Adult zebrafish lines were housed within the Mayo Clinic Zebrafish Core Facility, and mated in false bottom containers to generate offspring for line propagation and experimental purposes. Embryos obtained from individual pair crosses were mixed at 0 dpf for each experiment. The mixes were then divided into groups of 60 fish, which were transferred to 100 × 15 mm petri dishes (Becton, Dickinson and Company) containing 25 ml of 0.5X E2 media [29]. At 1 dpf, all nonviable embryos were removed from each group and the viable embryos were transferred to dishes containing 25 ml of fresh 0.5X E2 media. All embryos were raised in an incubator at 28.5°C with a 14/10-hour light/dark cycle until they were used in the experiments detailed in following Materials and methods subsections.

Zebrafish were obtained from crosses between wild type fish, and raised according to the Materials and methods, Zebrafish husbandry subsection. The samples were collected at 0.25 dpf, 1.0 dpf, 2.0 dpf, 3.0 dpf, 4.0 dpf, 5.0 dpf, 6.0 dpf, and 7.0 dpf, and were treated and stored as previously described [28]. The temporal patterns of eCB gene expression were investigated with quantitative reverse transcriptase-polymerase chain reactions (qRT-PCRs). The samples were processed, ribonucleic acids (RNA) were isolated, complementary deoxyribonucleic acids were synthesized, and qRT-PCRs were performed as previously described [28]. The primers used in the qRT-PCRs were designed to amplify a region of the following eCB genes: cnr1, cnr2, loc793909, abhd4, gde1, napepld, faah, faah2a, and faah2b (S1A Table). Additional primers were designed and used to amplify regions of four selected reference genes previously used in zebrafish gene expression studies: rps6kb1b (Forward 5'-AAATCTCTATGGCGCTCGGACACC-3', Reverse 5'-TGGACTCCTTACACAGCCCGAAATC-3'), eef1a1l1 (Forward 5'-TACAAATGCGGTGGAATCGACAAG-3', Reverse 5'-TCGGCCTTCAGTTTGTCCAACAC-3'), rpl13a (Forward 5'-TCTGGAGGACTGTAAGAGGTATGC-3', Reverse 5'-AGACGCACAATCTTGAGAGCAG-3'), and b2m (Forward 5'-GCCTTCACCCCAGAGAAAGG-3', Reverse 5'-GCGGTTG GGATTTACATGTTG-3') [30–32]. All primers were obtained from Integrated DNA Technologies (Integrated DNA Technologies Inc, Coralville, IA, USA). The obtained data was used to calculate mean expression ± 95% confidence intervals (95% CI) relative to the 5 dpf time point for each target gene. Comparisons between time points were made using a one-way analysis of variance (ANOVA) followed by Sidak's multiple comparisons test. All statistical analyses were performed using Prism 6 software (GraphPad Software, San Diego, CA, USA). To visualize a cross-comparison of all genes’ fold changes in this paper, heatmaps (S4–S11 Figs) were generated using the ‘gplots’ package in R software v3.4.1 (https://www.r-project.org). Fold changes were calculated using the 2^ΔCt formula for each gene pair combination when ΔCt was 0 or a positive number. However, when the ΔCt resulted in a negative value the following formula was used: -2^ΔCt.

Zebrafish in the experimental groups were obtained from crosses between wild type fish, while the zebrafish in the control groups were obtained from outcrosses between heterozygous casz1^mn0001Gt/+ and wild type fish (casz1^+/+) [30]. All fish were raised as described in the Materials and methods, Zebrafish husbandry subsection. The samples were collected at 2 dpf and 4 dpf as previously described [33]. Digoxigenin-labeled probes for mRFP were developed from cDNA clones as previously described [30]. Similarly, primers were obtained from Integrated DNA Technologies and used to develop digoxigenin-labeled probes for faah, faah2a, and faah2b (S1B Table). The in situ hybridization (ISH) experiment was performed using a previously published protocol with the following modifications [33]. In the modified protocol, the 2 dpf samples were treated with proteinase K (Roche) for 20 min and the 4 dpf samples were treated for 40 min. The proteinase K reactions were stopped with 3X 5 min PBT (Bio-Rad Laboratories) washes. The samples were refixed with a 4% PFA solution for 20 min, and then gently shaken during 5X 5 min PBT washes on a rotator (Fischer Scientific). The prehybridization step was completed using a 5 min wash with a 50% hybridization mix / 50% PBT solution, followed by a 4h incubation at 65°C with the hybridization mix. The hybridization step was completed using an overnight incubation at 65°C with hybridization mixes containing probes (1 ng/μl) for each gene of interest. Stringency washes were performed at 65°C to gradually replace the hybridization mixes with SSC buffer containing 0.1% Tween 20 (Bio-Rad Laboratories), which in turn was progressively replaced at room temperature with maleic acid buffer (0.1 M maleic acid, 0.15 M sodium chloride, pH = 7.5) containing 0.1% Tween 20 (MABT) (Bio-Rad Laboratories). The nonspecific antibody binding sites were blocked with 2% Blocking Reagent (Roche) in MABT. The samples were incubated in anti-digoxigenin-AP fab fragments antibody (Roche) diluted at 1:5,000 with blocking buffer. The antibody solution was removed and the samples were washed using 8X 15 min washes with MABT at room temperature. The samples were then equilibrated and stained. A series of single focal plane brightfield images was acquired for each sample using Specimen in a Corrected Optical Rotational Enclosure imaging techniques [34]. The image series were taken on an Axioplan 2 microscope (Carl Zeiss Microscopy) equipped with a Powershot G6 camera (Canon). Each series of images was compiled into a single composite image with Helicon Focus software (Helicon Soft).

The National Center for Biotechnology Information (NCBI) Danio rerio Annotation Release 105 was used to identify the predicted splice variants associated with each eCB gene of interest [35]. The Basic Local Alignment Search Tool (BLAST) was used to identify sequences encoding amino acid residues in faah2a that were putative homologues of the fatty acid amide hydrolase 2 serine-serine-lysine catalytic triad [16,35]. The TALENs were designed to recognize sites conserved in all predicted splice variants by using the Mojo Hand Version 2 software available online at http://www.talendesign.org [36]. The cnr1 TALEN binding sites were positioned around a 15 base pair spacer sequence containing a BstUI restriction enzyme site (S1A Fig). The faah2a TALEN binding sites were positioned around a 13 base pair spacer sequence containing a BsrI restriction enzyme site just upstream of the sequence predicted to encode a catalytic lysine (S1B Fig). Primers were designed to flank these sites so that the fish could be genotyped via restriction fragment polymorphism (RFLP) analyses (S1C Table) [37]. The TALEN vectors were created using the Golden Gate method with a pT3TS-GoldyTALEN destination vector [37,38]. The TALEN mRNAs were synthesized and microinjected as previously described [37]. The fish harboring mutant eCB alleles were outcrossed with a dominant leopard line (Cx41.8^+/tq270) to eventually establish populations of heterozygous F2 mutants [39]. The indels were characterized by sequencing, and the sperm from the F2 populations was cryopreserved [37]. The cnr1 mutant line was assigned the Zebrafish Information Network (ZFIN) designator cnr1^mn49, and the faah2a mutant line was assigned the ZFIN designator faah2a^mn50 [40].

The zebrafish were obtained from in crosses between the F2 heterozygous eCB mutant lines described in the Materials and methods, TALEN-mediated mutagenesis of eCB genes subsection. All viable fish were raised according to the Materials and methods, Zebrafish husbandry subsection. At 3 dpf the fish were individually transferred in 400 μl of fresh 0.5X E2 media to wells on 48-well tissue culture plates (Corning). All of the behavior assays were completed using previously published protocols with the following modifications [41,42]. The 48-well plates were transferred from the incubator to light-box apparatuses 30 min after the onset of the light cycle on 5 dpf. All experiments were performed at 28.5°C, and the light-box apparatuses were configured with custom acrylic templates designed to align two 48-well plates. The apparatuses housed one pair of plates during each assay, with one plate serving as a control plate and the other as an experimental plate. During the hyperosmotic stress assays fish were acclimated to the apparatus for 45 min before the experiment was initiated. The pre-treatment locomotor baseline activity of the fish was filmed for 15 min, and then the treatments were applied to each plate. The control plates were treated with 100 μl of E2 media and the experimental plates were treated with 100 μl of a 500 mM sodium chloride (Sigma-Aldrich) solution prepared in E2 media (+100 mM final sodium chloride concentration). The post-treatment locomotor activity of the fish was filmed for 31 min. During the nicotine assays fish were acclimated to the apparatus for 20 min before the experiment was initiated. The pre-treatment baseline locomotor activity of the fish was filmed for 5 min, and then the treatments were applied to each plate. The control plates were treated with 100 μl of E2 media and the experimental plates were treated with 100 μl of a 250 μM nicotine (Acros Organics) solution prepared in E2 media (+50 μM final nicotine concentration). The post-treatment locomotor responses of the fish were then filmed for 5 min. The locomotor activity of the fish was analyzed at 1 second intervals with MATLAB software (The MathWorks, Natick, MA, USA) to calculate the distances travelled by each fish. After the behavior assays were completed, the individual larval zebrafish were genotyped using the RFLP analysis described in the Materials and methods, TALEN-mediated mutagenesis of eCB genes subsection. The locomotor data was used to calculate means ± 95% CI. Comparisons between groups were made using a two-way analysis of variance (ANOVA) followed by Tukey's honest significant difference test. All statistical analyses were performed using R software. All graphs were generated using R software, and Illustrator CC software (Adobe, San Jose, CA, USA).

The transcript levels of nine eCB genes were analyzed by qRT-PCR using rps6kb1b as the primary housekeeping gene. These levels were assessed relative to 5 dpf, because it corresponded with the stage of development used in our hyperosmotic stress assays. The stage of development significantly affected the expression levels of each eCB gene that was profiled [cnr1: F(7, 16) = 117.70, P < 0.0001; cnr2: F(7, 16) = 58.35, P < 0.0001; loc793909: F(7, 16) = 13.61, P < 0.0001; abhd4: F(7, 16) = 93.30, P < 0.0001; gde1: F(7, 16) = 58.94, P < 0.0001; napepld: F(7, 16) = 123.80, P < 0.0001; faah: F(7, 16) = 91.65, P < 0.0001; faah2a: F(7, 16) = 278.10, P < 0.0001; faah2b: F(7, 16) = 291.20, P < 0.0001]. The transcript levels of eCB receptors encoded by cnr1 and cnr2 increased to a peak at 5 dpf, and then exhibited declines at 6 dpf and 7 dpf (Fig 1A and 1B). The transcript levels of the eCB receptor encoded by loc793909 increased at 2 dpf, and then exhibited an insignificant decrease that was followed by stability through 7 dpf (Fig 1C). The transcript levels of eCB anabolic enzymes encoded by abhd4 and napepld increased through the first 5 dpf, and then exhibited declines at 6 dpf and 7 dpf (Fig 2A and 2C). In contrast, the transcript levels of the eCB anabolic enzyme encoded by gde1 peaked at 1dpf, then declined through 7 dpf (Fig 2B). The transcript levels of eCB catabolic enzymes encoded by faah and faah2a increased through the first 4 dpf, then declined from 5–7 dpf (Fig 3A and 3B). The transcript levels of the eCB catabolic enzyme encoded by faah2b dropped between 0.25 dpf and 1dpf, subsequently increased through 4 dpf, and then declined again through 7 dpf (Fig 3C). The qRT-PCR threshold cycle (Ct) values recorded at 5 dpf for each eCB gene were included in the Supporting Information (S2A Table). The trends in eCB gene expression from 1–7 dpf were largely conserved regardless of the reference gene that was used, however, the stability of reference gene expression was not always conserved at the 0.25 dpf time point (S4–S11 Figs). In addition to the qRT-PCR experiments, ISH was used to investigate the spatiotemporal expression patterns of faah, faah2a, and faah2b. The expression of these genes were not detected by ISH in wild type fish at 2 dpf, however, they were detected at 4 dpf (S2A Fig). The expression of all three serine hydrolases was detected in the intestinal bulb. Additionally, faah and faah2a expression was detected in the liver. Background staining was accounted for by examining mRFP expression in a mix of siblings obtained by outcrossing heterozygous casz1^mn0001Gt/+ fish with wild type (casz1^+/+) fish [30]. The mRFP expression pattern of the GBT0001 line mimics that of casz1, and is only present in fish that harbor the gene-break transposon [30]. In these control fish—which were prepared alongside the experimental groups—no background staining was observed in the negative control wild type siblings while staining was observed in the positive control heterozygous casz1^mn0001Gt/+ siblings (S2B Fig).

TALENs were used to generate zebrafish lines with frameshift mutations in the eCB receptor gene cnr1 and the eCB metabolic enzyme gene faah2a. One pair of TALENs was designed for a target sequence in exon 1 of cnr1 and used to generate a F2 heterozygous mutant line with a 20 base pair deletion (Fig 4, S1A Fig). This mutant line was then used to study how cnr1 regulated stress-associated behavior (Fig 5A). Another pair of TALENs was designed for a target sequence in exon 3 of faah2a and used to generate a F2 heterozygous mutant line with a 10 base pair deletion (Fig 6, S1B Fig). Similarly, this mutant line was then used to study how faah2a regulated stress-associated behavior (Fig 7A). We have previously demonstrated that 5 dpf zebrafish elevate whole-body cortisol levels when challenged with hyperosmotic conditions, and that this neuroendocrine stress response correlated with an increase in stress-associated locomotion [28,42,43]. In the present study, a main affect of treatment was observed in the hyperosmotic stress assays performed with 5 dpf zebrafish obtained from F2 heterozygous eCB mutant in crosses [cnr1: F(1, 663) = 53.75, P < 0.0001; faah2a: F(1, 474) = 54.60, P < 0.0001]. A main affect of genotype was also observed in both of these assays [cnr1: F(2, 662) = 5.307, P < 0.01; faah2a: F(2, 473) = 3.502, P < 0.05]. No significant interaction was observed between treatment and genotype [cnr1: F(2, 659) = 2.114, P > 0.0.5; faah2a: F(2, 470) = 2.516, P > 0.05]. The wild type post-treatment Stress group exhibited a statistically significant increase in locomotor activity relative the wild type post-treatment Control group in each assay (Figs 5B and 7B). The homozygous cnr1 mutant Stress group displayed a significantly higher locomotor response than the homozygous cnr1 mutant Control group and the wild type Stress group (Fig 5B). The homozygous faah2a mutant stress group displayed a locomotor response that was not significantly different from the homozygous faah2a mutant Control group, and that was significantly lower than the wild type Stress group (Fig 7B). A nicotine assay was performed to demonstrate that the attenuated response of the homozygous faah2a mutants was specific to the hyperosmotic stress assay, and not simply attributed to a reduced capacity for locomotion (S3A Fig). Nicotine treatment has been shown to elicit robust locomotor responses from larval zebrafish, and in the present study a main affect of treatment was observed in the nicotine assay [faah2a: F(1, 90) = 81.51, P < 0.0001] [41,44]. No main affect of genotype was observed [faah2a: F(2, 89) = 1.890, P > 0.05], and no significant interaction between treatment and genotype was observed [faah2a: F(2, 86) = 1.224, P > 0.05]. There were no significant differences in the responses of each genotype to nicotine treatment (S3B Fig).