DOI: 10.5061/dryad.gtht76hvp
title: Data from: The genomics of mimicry: Gene expression throughout development provides insights into convergent and divergent phenotypes in a Müllerian mimicry system

description: A common goal in evolutionary biology is to discern the mechanisms that
 produce the astounding diversity of morphologies seen across the tree of
 life. Aposematic species, those with a conspicuous phenotype coupled with
 some form of defense, are excellent models to understand the link between
 vivid color pattern variations, the natural selection shaping it, and the
 underlying genetic mechanisms underpinning this variation. Mimicry systems
 in which species share a conspicuous phenotype can provide an even better
 model for understanding the mechanisms of color production in aposematic
 species, especially if comimics have divergent evolutionary histories.
 Here we investigate the genetic mechanisms by which mimicry is produced in
 poison frogs. We assembled a 6.02 Gbp genome with a contig N50 of 310 Kbp,
 a scaffold N50 of 390 Kbp, and 85% of expected tetrapod genes. We
 leveraged this genome to conduct gene expression analyses throughout
 development of four color morphs of R. imitator and two color morphs from
 both R. fantastica and R. variabilis which R. imitator mimics. We
 identified a large number of pigmentation and patterning genes
 differentially expressed throughout development, many of them related to
 melanophores/melanin, iridophore development, and guanine synthesis. We
 also  identify the pteridine synthesis pathway (including genes
 such as qdpr and xdh) as a key driver of the variation in color between
 morphs of these species, and identify several plausible candidates for
 coloration in vertebrates (e.g., cd36, ep-cadherin, perlwapin). Finally,
 we hypothesize that keratin genes (e.g., krt8) are important for producing
 different structural colors within these frogs.
description: <em>G</em><em>ene expression sample
 preparation:</em>                 Samples were
 prepared differently for the mimic (<em>R. imitator</em>) and
 the model species (<em>R. fantastica </em>and<em> R.
 variabilis</em>). During the course of our work, we discovered that
 there were multiple groups approaching the same questions using collected
 samples from different species, but at slightly different timepoints. In
 light of this, we chose to combine our efforts into a single manuscript in
 an attempt at making broader inferences. We acknowledge this, and as a
 result, the data in this manuscript are analyzed in a manner concordant
 with these differences. <em>Ranitomeya
 imitator:</em> The initial breeding stock of
 <em>Ranitomeya imitator</em> was purchased from Understory
 Enterprises, LLC (Chatham, Canada). Frogs used in this project represent
 captive-bred individuals sourced from the following wild populations going
 clockwise from top left in Figure 1: Tarapoto (green-spotted), Sauce
 (orange-banded), Varadero (redheaded), and Baja Huallaga (yellow-striped).
 Tadpoles were sacrificed for analyses at 2, 4, 7, and 8 weeks of age. We
 sequenced RNA from a minimum of three individuals at each time point from
 the Sauce, Tarapoto, and Varadero populations (except for Tarapoto at 8
 weeks), and two individuals per time point from the Huallaga population.
 Individuals within the same time points were sampled from different family
 groups (Appendix Table 1). Tadpoles were anesthetized
 with 20% benzocaine (Orajel), then sacrificed via pithing. Whole skin was
 removed and stored in RNA later (Ambion) at -20°C until RNA extraction.
 Whole skin was lysed using a BeadBug (Benchmark Scientific, Sayreville,
 NJ, USA), and RNA was then extracted using a standardized Trizol protocol.
 RNA was extracted from the whole skin using a standardized Trizol
 protocol, cleaned with DNAse and RNAsin, and purified using a Qiagen
 RNEasy mini kit. RNA Libraries were prepared using standard poly-A tail
 purification with Illumina primers, and individually barcoded using a New
 England Biolabs Ultra Directional kit as per the manufacturer’s protocol.
 Individually barcoded samples were pooled and sequenced using 50 bp paired
 end reads on three lanes of the Illumina HiSeq 2500 at the New York Genome
 Center. <em>Ranitomeya fantastica</em>
 and<em> Ranitomeya variabilis</em>: We set
 up a captive colony in Peru (see Appendix) consisting of between 6 and 10
 wild collected individuals per locality. We raised the tadpoles on a diet
 consisting of a 50/50 mix of powdered spirulina and nettle, which they
 received 5 times a week. Tadpoles were raised individually in 21oz plastic
 containers, within outside insectaries covered with 50% shading cloth, and
 water change was performed with rainwater.  Three tadpoles per stage (1,
 2, 5, 7, and 8 weeks after hatching; see Appendix Table 1) were fixed in
 an RNAlater (Ambion) solution. To do so, tadpoles were first euthanized in
 a 250 mg/L benzocaine hydrochloride bath, then rinsed with distilled water
 before the whole tadpole was placed in RNAlater and stored at 4°C for 6h
 before being frozen at -20°C for long-term storage. Before RNA extraction,
 tadpoles were removed from RNA later and the skin was dissected off. Whole
 skin was lysed using a Bead Bug, and RNA was then extracted using a
 standardized Trizol protocol. RNA libraries were prepared using standard
 poly-A tail purification, prepared using Illumina primers, and
 individually dual-barcoded using a New England Biolabs Ultra Directional
 kit. Individually barcoded samples were pooled and sequenced on four lanes
 of an Illumina HiSeq X at NovoGene (California, USA). Reads were paired
 end and 150 base pairs in length.
 <em>Differential gene expression:</em>
                 We indexed our new <em>Ranitomeya imitator
 </em>genome using STAR version2.7.10a_alpha_220601 (Dobin et al.
 2013). We aligned our reads to our genome using STAR
 version2.7.10a_alpha_220601 (Dobin et al.
 2013), allowing 10 base mismatches (--outFilterMismatchNmax 10),
 a maximum of 5 multiple alignments per read (--outFilterMultimapNmax 5),
 and discarding reads that mapped at less than 50% of the read length
 (--outFilterScoreMinOverLread 0.5). We then counted aligned reads using
 htseq-count version 2.0.2 (Anders, Pyl, and
 Huber 2015).
 <strong>References:</strong>
 Anders, Simon, Paul Theodor Pyl, and Wolfgang Huber. 2015.
 “HTSeq—a Python Framework to Work with High-Throughput Sequencing Data.”
 Bioinformatics  31 (2): 166–69. Dobin, Alexander,
 Carrie A. Davis, Felix Schlesinger, Jorg Drenkow, Chris Zaleski, Sonali
 Jha, Philippe Batut, Mark Chaisson, and Thomas R. Gingeras. 2013. “STAR:
 Ultrafast Universal RNA-Seq Aligner.” Bioinformatics  29 (1):
 15–21.
description: # The genomics of mimicry: gene expression throughout development provides
 insights into convergent and divergent phenotypes in a Müllerian mimicry
 system
 [https://doi.org/10.5061/dryad.gtht76hvp](https://doi.org/10.5061/dryad.gtht76hvp) ## Description of the data The data in this dataset are gene count files of gene expression in the skin for individual frogs. Files are named using this scheme: "Genus abbreviation_specific epithet_morph name_age in weeks_unique identifier.gene.counts". So for example, file "R_fantastica_banded_8_M15_088_L2.gene.counts" refers to an *R. fantastica *individual of the banded population that was 8 weeks old. Files are tab-delimited and contain two columns, one for the transcript ID from our annotation, and a second of the number of counts of that particular gene. Genomic data are on the European Nucleotide Archive under accession PRJEB28312. The methods describe how data was collected and processed. ## Code Code is available on [Zenodo.](https://zenodo.org/doi/10.5281/zenodo.4758403)