2.4 RNA-seq
We terminated the diet treatment at 17 days and each individual was
sacrificed. The interopercle, interopercle-mandibular ligament, and
retroarticular (i.e., IOP-RA complex, Figure 1) was dissected from each
fish for RNA-seq (left-side) or qPCR (right-side). It is to be noted
that this complex is a mix of tissue types (e.g bone, ligament,
epithelial tissue). For RNA-seq, 6 samples from each treatment and
species were stored in Trizol (Invitrogen) at -80°C, homogenized using a
Next Advance Bullet Blender and 5 UFO beads each, and processed using
the phenol-chloroform method of RNA extraction, but did not undergo cDNA
conversion. We standardized each sample to 500ng total RNA in 50uL and
produced libraries using the TruSeq Stranded mRNA Library Prep Kit
(Illumina). Any remaining RNA not used for RNA-seq was stored at -80°C.
Libraries were sequenced at the University of Massachusetts Medical
School Deep Sequencing Core with a HiSeq 4000 with 50 x 50 paired end
reads.
Raw reads from RNA-Seq were assessed using FastQC (Andrews, 2010) and
ends were trimmed accordingly using cutadapt (Martin, 2011). Cleaned
reads were -mapped against the Maylandia zebra genome version UMD2a
(Yates et al., 2020) with Bowtie2 (Langmead and Salzberg, 2012) and a
matrix of read counts was generated from the alignments with
HT-Seq-count (Anders et al., 2015).
We used edgeR (Chen et al., 2014; Robinson et al., 2010) to identify
differentially expressed genes between treatments, as well as those that
showed an additive effect between species and environment. Results were
groundtruthed by visually comparing normalized counts
(counts-per-million; cpm) among treatments. Gene ontology (GO) terms
were assigned using biomaRt (Durinck et al., 2005, 2009) via the
biomartr package (Drost and Paszkowski, 2017) and enrichment analysis
was conducted via topGO (Alexa and Rahnenfuhrer) in the R environment (R
Core Team, 2021) using the weight01 algorithm and Fisher’s exact test.
Tissues from individuals not used for RNA-seq, were stored in Trizol at
-80°C, and homogenized as previously described. We followed the
phenol-chloroform RNA extraction method and converted RNA to cDNA for
gene expression validation purposes. These samples were standardized to
an RNA concentration of 70ng/uL.
qPCR was used to measure gene expression of genes found to be
differentially expressed and/or differentially accessible from the
RNA-seq and ATAC-seq datasets, and those related to bone
development/homeostasis using the comparative CT method. qPCR primer
sequences are listed in the supplementary table (Table S1).