Comparative transcriptomics reveal tissue level specialization towards diet in prickleback fishes

Beyond a few obvious examples (e.g., gut length, amylase activity), digestive and metabolic specializations towards diet remain elusive in fishes. Thus, we compared gut length, δ13C and δ15N signatures of the liver, and expressed genes in the intestine and liver of wild-caught individuals of four closely-related, sympatric prickleback species (family Stichaeidae) with different diets: Xiphister mucosus (herbivore), its sister taxon X. atropurpureus (omnivore), Phytichthys chirus (omnivore) and the carnivorous Anoplarchus purpurescens. We also measured the same parameters after feeding them carnivore or omnivore diets in the laboratory for 4 weeks. Growth and isotopic signatures showed assimilation of the laboratory diets, and gut length was significantly longer in X. mucosus in comparison to the other fishes, whether in the wild, or in the lab consuming the different diets. Dozens of genes relating to digestion and metabolism were observed to be under selection in the various species, but P. chirus stood out with some genes in the liver showing strong positive selection, and these genes correlating with differing isotopic incorporation of the laboratory carnivore diet in this species. Although the intestine showed variation in the expression of hundreds of genes in response to the laboratory diets, the liver exhibited species-specific gene expression patterns that changed very little (generally <40 genes changing expression, with P. chirus providing an exception). Overall, our results suggest that the intestine is plastic in function, but the liver may be where specialization manifests since this tissue shows species-specific gene expression patterns that match with natural diet. Supplementary Information The online version contains supplementary material available at 10.1007/s00360-021-01426-1.

. A. Respirometer setup showing the chamber submerged in ambient seawater in the closed, recirculating configuration with thermistor and oxygen probe in series. The oxygen consumption of each fish was measured individually in 15minute intervals after a 30 min acclimation period to the respirometry chamber. Following each 15-min interval, the valves of the system were opened manually to the open configuration to exchange with the flow-through, ambient seawater before being closed again for the next measurement period. B. Representative plot showing the oxygen concentrations in the respirometer over time during measurements of fish metabolic rate. The portions with the negative slopes are the measurement periods when the system was closed (15 min intervals). The system was then opened again for several minutes to be flushed by new seawater from the flow-through system, and then closed again to take another measurement on the same fish. This process was repeated three times. The above traces are from an individual of Anoplarchus purpurescens C that weighed 6.65g. Supplemental Figure S2: X. mucosus H wild fish replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced for transcriptomic data. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.

X. mucosus H wild individuals
Supplemental Figure S3: X. mucosus H laboratory omnivore diet fish replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced for transcriptomic data. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.

X. mucosus H lab omnivore individuals
Supplemental Figure S4: X. mucosus H laboratory carnivore diet replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced for transcriptomic data. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.

X. mucosus H lab carnivore individuals
Supplemental Figure S5: X. atropurpureus O wild fish replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced for transcriptomic data. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.

X. atropurpureus O wild individuals
Supplemental Figure S6: X. atropurpureus O laboratory carnivore diet replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced for transcriptomic data. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.

X. atropurpureus O lab carnivore individuals
Supplemental Figure S7: P. chirus O wild fish replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced for transcriptomic data. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.

P. chirus O wild individuals
Supplemental Figure S8: P. chirus O laboratory carnivore diet replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced for transcriptomic data. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.

P. chirus O lab carnivore individuals
Supplemental Figure S9: A. purpurescens C wild fish replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced for transcriptomic data. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.
A. purpurescens C wild individuals Supplemental Figure S10: A. purpurescens C laboratory omnivore diet replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced for transcriptomic data. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.
A. purpurescens C lab omnivore individuals Supplemental Figure S11: A. purpurescens C laboratory carnivore diet replicates. PCA plot to depict the quality check analysis of individual replicates within the same species and diet group and across the three tissues we sequenced. Spheres depict Liver replicates, Squares depict pyloric ceca replicates, and triangles depict mid-intestine replicates.
A. purpurescens C lab carnivore individuals Supplemental Figure S12: PCA plot generated with Batch Quality Check Results. Color indicated batch (blue is first run, orange is second run, green is third run. Samples do not cluster by batch. The standardized Pearson correlation coefficient is 0.87 and the Cramer's V is 0.7, indicating batch does not fully interfere with the signal, with batch 3 showing the most uniqueness because this batch contained the liver samples. All species are represented in each cluster of liver samples (batch 3) on the plot, showing that they are dispersed throughout and not grouping by sample.    Figure S14: An adaptive branch-site random effects likelihood (aBSREL) test for episodic diversification phylogenetic tree constructed for various genes in the pyloric ceca from four prickleback fish species: a) serine protease 27 (PR27), b) tubulin alpha chain (TBA), c) elastase (ELA1). ω is the ratio of nonsynonymous to synonymous substitutions. The color gradient represents the magnitude of the corresponding ω. Branches thicker than the other branches have a p < 0.05 (corrected for multiple comparisons) to reject the null hypothesis of all ω on that branch (neutral or negative selection only). A thick branch is considered to have experienced diversifying positive selection.

A ) B) C)
Supplemental

Supplemental Results: Relative Gene Expression
We used RNA-seq data of the liver, pyloric ceca and mid-intestine to observe the suites of genes that changed with different diets and how species respond to dietary variation. Note that we are only reporting on pathways relevant to digestion and metabolism of specific nutrient classes (Fig. 4, Table 4). If a cluster is not mentioned, yet depicted in the heatmap, then the genes within that cluster were not directly relevant to digestion and nutrient metabolism.

Liver
There were 11 DEGs when comparing wild fish (WF) and laboratory carnivore diet fish (LC) of X. atropurpureus O , out of which 36% were annotated (Supplemental Figure S14).
Cluster 1 (elevated in wild fish) contained genes important for glucose and fatty acid metabolism (Supplemental Table S7). These proteins are important for energy storage, insulin signaling pathway and glucagon signaling pathway. Cluster 2 (elevated in LC fish) contained one unannotated gene.

P. chirus O stands out with 302 DEGs when comparing liver gene expression among WF
and LC P. chirus O , out of which 14% of genes were annotated (Supplemental Figure S15).

Cluster 1 (elevated in wild fish) contains genes involved in cholesterol metabolism
(Supplemental Table S4). Cluster 2 (elevated in LC fish) contains genes for lipid metabolism, fatty acid synthesis, and bile acid biosynthesis.
There are 19 DEGs when comparing WF, laboratory omnivore diet fish (LO) and LC A. purpurescens C , out of which 32% of genes were annotated (Supplemental Figure S16). Cluster 1 (elevated in wild fish) consists of genes for cholesterol homeostasis and genes that play a role in controlling the metabolism of fatty acids, specifically glycerophospholipid metabolism and glycerolipid metabolism (Supplemental Table S4).

Pyloric ceca
There were 1226 DEGs when comparing WF to LC X. atropurpureus O , out of which 68.1% were annotated (Supplemental Figure S17). Wild X. atropurpureus upregulated genes in Cluster 1 (elevated in wild fish), that were involved in digestive processes for chitin degradation, glycolysis, glycogen catabolic process, glycosaminoglycan biosynthesis, bile acid metabolism, proteolysis, carbohydrate metabolic process, collagen metabolic/catabolic process, pentose phosphate pathway, lipid metabolism, and glutamate biosynthetic process (Supplemental Table   S8). LC X. atropurpureus O upregulated genes in Cluster 2 (elevated in LC fish), that were involved in lipid metabolism.
Like in the liver, P. chirus O showed differing DEGs in comparison to the other species, with only 19 DEGs in the pyloric ceca among WF and LC P. chirus O ; 94.7% of the genes were annotated (Supplemental Figure S18). Cluster 1 (elevated in wild fish), featured genes involved in carboxypeptidase activity, proteolysis, and carbohydrate binding (Supplemental Table S8).
There were no Cluster 2 (elevated in LC fish) genes in P. chirus O pyloric ceca.
There were 259 DEGs when comparing WF, LC, and LO A. purpurescens C , out of which 62.5% were annotated (Supplemental Figure S19). Cluster 1 (elevated in wild fish) contained genes involved in endopeptidase/trypsin activity and insulin receptor signaling pathway (Supplemental Table S8). Cluster 2 (wild-omnivore genes) contained genes involved in carbohydrate metabolism, bile acid metabolism, cholesterol catabolism, and fatty acid synthesis.
Cluster 3 (elevated in the lab genes) contains a large amount of genes, although not directly involved in digestion or metabolism.

Mid-intestine
There were 343 DEGs for W and LC X. atropurpureus O , out of which 83.96% were annotated (Supplemental Figure S20). Cluster 1 (elevated in wild fish) contained genes are involved in glycolysis, cholesterol biosynthesis, lipid metabolism, carbohydrate metabolism, protein metabolism, fatty acid biosynthesis, and pentose phosphate pathway (Supplemental Table   S9). Cluster 2 (elevated in LC fish) genes are involved in cholesterol metabolism.
There were 298 DEGs for WF and LC P. chirus O , out of which 90.6% were annotated (Supplemental Figure S21). Cluster 1 (elevated in wild fish) contained 293 of the genes, including those involved in gluconeogenesis, protein deubiquination, creatine metabolic process, and calcium ion transport (Supplemental Table S9). Five genes not directly involved in digestion or nutrient metabolism composed Cluster 2 (elevated in LC fish).
There were 872 DEGs for WF, LC, and LO A. purpurescens C , out of which 82.2% were annotated (Supplemental Figure S22). Cluster 1 (elevated in wild fish) contained genes involved in mannose metabolism, glycogen catabolism, insulin signaling, gluconeogenesis, glucose homeostasis and lipid catabolism (Supplemental Table S9). Cluster 2 (wild-omnivore genes) contained genes involved in bile acid metabolism. Cluster 4 (wild-carnivore genes) contained genes involved in proteolysis. Cluster 5 (carnivore genes) contained genes involved in collagen catabolism. Cluster 6 (omnivore genes) contained genes involved in cholesterol biosynthesis and glucose metabolism. Figure S15: Differential gene expression depicted as a heatmap in the liver of X. atropurpureus O . Yellow indicates elevated relative expression, whereas blue indicates low expression. Each row is a single gene, and genes are clustered in a dendrogram (on left of each heatmap) by similarity of expression patterns. The various clusters of genes are described in Table 4. Each column represents the gene expression in a single tissue from an individual fish, with WF = wild-caught fish, LO = fish fed an omnivore diet in the laboratory (in the case of X. mucosus H and A. purpurescens C ), and LC = fish fed a carnivore diet in the laboratory.

WF WF LC LC
Wild Fish Lab-Carnivore Diet   Figure S24. Correlation matrix of expressed genes in the liver, pyloric ceca, and mid intestine of four prickleback species caught from the wild (WF), or fed omnivore (LO) or carnivore (LC) diets in the laboratory. Clustering of WF are depicted by black bars in between the dendrogram and the correlation matrix, whereas LO fishes are depicted by purple bars, and LC fishes by red bars. All sample names are depicted on the right side (60 individual tissues in total), with symbols for each tissue type used to emphasize clustering. Correlation matrix created with the Trinity toolkit "PtR" and Pearson correlation as sample distances. We emphasize five

Metabolic Rate
Contrary to our expectation that fishes consuming the carnivore diet in the laboratory would have higher metabolic rates than fishes consuming the omnivore diet, routine metabolic rate did not vary among the species or intra-specifically on the different diets, suggesting that body mass is still one of the main determinants of metabolic rate in fishes, and these fishes are all similar in size in comparison to the range of sizes fishes can attain (Gillooly et al. 2001;Clarke and Johnston 1999;Ikeda 2016). More detailed measures of metabolic rate across longer time scales undoubtedly would show differences in Specific Dynamic Action for fishes consuming the different diets in the laboratory (Secor 2009), but we only measured routine metabolic rate.
Given the short period of time over which we measured metabolic rate, our results were possibly influenced by stress and thus, more detailed analyses of metabolic rate in prickleback fishes are