High phylogenetic utility of an ultraconserved element probe set designed for Arachnida

12 Pages·2016·0.7 MB·English

MolecularEcologyResources(2016) doi:10.1111/1755-0998.12621
High phylogenetic utility of an ultraconserved element probe
set designed for Arachnida
JAMES STARRETT,*1,2 SHAHAN DERKARABETIAN,*†2 MARSHAL HEDIN,*
ROBERT W. BRYSON JR.,‡§ JOHN E. MCCORMACK§ andBRANT C. FAIRCLOTH¶
*DepartmentofBiology,SanDiegoStateUniversity,5500CampanileDrive,SanDiegoCA92182,USA,†DepartmentofBiology,
UniversityofCalifornia,Riverside,900UniversityAvenue,RiversideCA92521,USA,‡DepartmentofBiologyandBurke
MuseumofNaturalHistoryandCulture,UniversityofWashington,4331MemorialWayNortheast,SeattleWA98195,USA,
§MooreLaboratoryofZoology,OccidentalCollege,LosAngelesCA90041,USA,¶DepartmentofBiologicalSciencesandMuseum
ofNaturalScience,LouisianaStateUniversity,BatonRougeLA70803,USA
Abstract
Arachnidaisanancient,diverseandecologicallyimportantanimalgroupthatcontainsanumberofspeciesofinterest
formedical,agriculturalandengineeringapplications.Despitetheirimportance,manyaspectsofthearachnidtreeof
life remain unresolved, hindering comparative approaches to arachnid biology. Biologists have made considerable
effortstoresolvethearachnidphylogeny;yet,limitedandchallengingmorphologicalcharacters,aswellasadearthof
genetic resources, have hindered progress. Here, we present a genomic toolkit for arachnids featuring hundreds of
conservedDNAregions(ultraconservedelementsorUCEs)thatallowtargetedsequencingofanyspeciesinthearach-
nid tree of life. We used recently developed capture probes designed from conserved regions of available arachnid
genomestoenrichasampleoflocifrom32diversearachnids.Sequencecapturereturnedanaverageof487UCEloci
forallspecies,witharangefrom170to722.PhylogeneticanalysisoftheseUCEsproducedahighlyresolvedarachnid
tree with relationships largely consistent with recent transcriptome-based phylogenies. We also tested the phyloge-
neticinformativenessofUCEprobeswithinthespider,scorpionandharvestmanorders,demonstratingtheutilityof
thesemarkersatshallowertaxonomicscalesandsuggestingthattheselociwillbeusefulforspecies-leveldifferences.
This probe setwillopenthedoortophylogenomic andpopulationgenomicstudiesacrossthearachnidtreeoflife,
enablingsystematics,speciesdelimitation,speciesdiscoveryandconservationofthesediversearthropods.
Keywords: Araneae, arthropods, Opiliones, phylogenomics, population genomics, Scorpiones, target enrichment, ultra-
conservedelements
Received28July2016;revisionreceived24September2016;accepted4October2016
arachnids have been described, with spiders and mites
Introduction
ranking among the most diverse of all animal orders
Arachnidaisanextremelyancientanddiversearthropod (Harvey2002;Zhang2011).Yet,morethanhalfofarach-
lineage, including conspicuous taxa such as spiders, nidspeciesdiversityremainstobediscovered(Chapman
scorpions, mites and ticks. The oldest known arachnid 2009). Arachnida also contains a number of medically
fossilsconsistmostlyofscorpionsandextincttrigonotar- importantvenomousordisease-vectorspeciesandmany
bids from the Silurian period, 443–416 MYA (million important agricultural mite pests, while spiders are of
years ago; Laurie 1899; Jeram et al. 1990; Dunlop 1996; particular interest to biologists and engineers for the
Dunlop et al. 2008). More than 110 000 species of strong and elastic silk fibres they produce. Despite the
attention arachnids have received for their ecological
Correspondence:ShahanDerkarabetian,DepartmentofBiology, importanceandpracticalutilitytohumans,phylogenetic
San Diego State University, 5500 Campanile Drive, San Diego relationships among and within many arachnid orders
CA92182,USA;E-mail:[email protected]
remainuncertain.Attherootofthisproblemisthatmor-
phological characters are limited and difficult to inter-
1Present address: Department of Biological Sciences,
pret(Shultz2007),andgenomicresourcesforthisgroup
AuburnUniversity,101RouseLifeSciences,AuburnAL
are sparse. Adding to the difficulty is the uncertainty in
36849,USA
the rooting of the arachnid tree, with fossil,
2Equalcontribution.
©2016JohnWiley&SonsLtd
2 J. STARRETT ET AL.
morphological and molecular data recovering drastic 2016). We then demonstrated the phylogenetic utility of
discrepanciesinearlyarachnidrelationshipsandplacing this probe set at multiple evolutionary timescales, from
traditionalnonarachnidchelicerates(e.g.Limulus)within thosespanninghundredsofmillionsofyearstospecies-
Arachnida (Wheeler & Hayashi 1998; Masta et al. 2009; leveldivergencesoflessthan10millionyears.Werecon-
Regieret al.2010;Sharmaet al.2014).Genomewidephy- structed well-supported phylogenies between and
logenetic markers are essential for resolving deep rela- within arachnid orders and demonstrated that these
tionships within Arachnida and for helping to uncover UCE baits may also be useful when reconstructing
thenumerousarachnidspeciesthatawaitdiscovery. species-levelrelationships.
Ultraconserved elements (UCEs) provide one poten-
tial target for a universal set of genomic markers for
Materialsandmethods
arachnids and would allow researchers to collect geno-
mic information from diverse taxa across the arachnid
Ultraconservedelementlibraryconstructionand
tree of life. UCEs are segments of DNA that are highly
enrichment
conserved across divergent taxa (Bejerano et al. 2004)
andarethoughttoregulateand/orenhancegeneexpres- Ultraconserved element capture was tested using the
sion (Alexander et al. 2010). UCEs can also be used as principal arachnid bait set (Faircloth 2016) on 32 arach-
anchors to target and retrieve variable DNA sequence nidsamplesrepresentingsixorders(Table 1).Theorders
flankingthecoreUCEregions.ThisflankingDNAshows selected span the root of Arachnida (Regier et al. 2010;
atrendofincreasinggeneticvariabilityasdistancefrom Sharma et al. 2014). Within three orders (Araneae, Opil-
the core UCE increases (Faircloth et al. 2012). UCEs are ionesandScorpiones),taxawereselectedfrommajorlin-
ideal markers for molecular systematics for several rea- eages,includingacross theinternal rootfor each(Hedin
sons. Whereas transcriptomes require high-quality RNA et al.2012;Bondet al.2014;Sharmaet al.2015).Toassess
as input, the UCE protocol only requires DNA, and UCE variability between closely related taxa, we
enrichmentscanbeperformedusingrelativelylowstart- included two turret spiders from the Antrodiaetus riversi
ing DNA concentrations. This allows the method to be complex (Hedin et al. 2013), two congeneric harvestmen
extendedtosmall-bodiedtaxa,largecollectionsofspeci- (Briggsus pacificus and B. bilobatus) and a second bark
mens preserved for Sanger-based DNA research and scorpion Centruroides sculpturatus to complement the
even to ‘standard’ museum specimens with varying published C. exilicauda genome. Voucher specimens are
levelsofDNAdegradation(McCormacket al.2015;Blai- deposited in the San Diego State University Terrestrial
mer et al. 2016). Homology between UCE loci across Arthropod Collection (SDSU_TAC). In the case of small
divergent taxa is also easy to assess because the core arachnids (e.g. mites and ticks), whole specimens were
UCEregionoftendisplays>95%sequencesimilarity,and used in extractions, but vouchers from the same locality
UCE cores are rarely duplicated (Derti et al. 2006). aredepositedintheSDSU_TAC.
AlthoughcoreUCEregionsshowreducedsequencevari- Genomic DNA was extracted from legs or whole
ation, UCE-flanking DNA shows levels of phylogenetic specimens using the DNeasy Blood and Tissue Kit
informativeness equal to or exceeding that of tradition- (Qiagen, Valencia, CA) following the manufacturer’s
allyusedprotein-codingmarkers(Gilbertet al.2015). protocol. DNA concentrations were determined with a
Ultraconserved elements have been successfully used Qubit fluorometer (Life Technologies, Inc) and run out
in phylogenetic studies at multiple taxonomic levels, ona1%agarosegeltoassessquality.Theloweststarting
including recent divergences (~5 MYA, within species; DNA quantity for processing was 108 ng (Siro acaroides,
Smith et al. 2014; Harvey et al. 2016; Manthey et al. 2016) a small-bodied harvestmen typically <2 mm in length),
andancient divergences (e.g. all amniotes; Faircloth et al. while most samples started with approximately 500 ng
2012). Originally developedfor use intetrapods, thevast (Table 1). Samples with high molecular weight DNA
majority of UCE phylogenetic studies have been con- were fragmented with a QSonica Q800R sonicator for
ducted on vertebrate taxa, including Amniota (Faircloth 10–12 cycles of 20 s on/20 s off, resulting in fragments
et al. 2012), mammals (McCormack et al. 2012), birds predominantly in the range of 300–1000 bp. DNA from
(McCormacket al.2013;Smithet al.2014;Meiklejohnet al. sample AR014 (Neomolgus littoralis) was partially
2016), reptiles (Crawford et al. 2012, 2015; Streicher et al. degradedandnotsonicated.
2016) and fish (Faircloth et al. 2013). More recently, UCE Libraries were prepared with the KAPA Hyper Prep
probe sets have been developed for use in insect taxa, Kit (Kapa Biosystems), with a generic SPRI substitute
includingHymenoptera(Fairclothet al.2015),Coleoptera, used for bead clean-up steps (Rohland and Reich 2012;
Diptera,LepidopteraandHemiptera(Faircloth2016). Glenn et al. 2016). Universal adapters were ligated onto
Here, we performed an in vitro test of a recently end-repaired and A-tailed DNA fragments. Each adap-
developed bait set targeting arachnid UCEs (Faircloth ter-ligatedlibrarywasamplifiedina50-lLtotalreaction
©2016JohnWiley&SonsLtd
UTILITY OF AN ARACHNID UCE PROBE SET 3
Table1 Sampleinformationandsequencingstatistics
DNA Rawread Reads
Order Species Voucher (ng) Pool pairs passQC %PassQC Contigs UCEloci
Acariformes Neomolguslittoralis AR014 500 G 633368 595582 94.03 26192 170
Amblypygi Damondiadema AR010 500 F 1623520 1484780 91.45 26761 442
Araneae Antrodiaetusriversi AR17 494 A 1277026 1216925 95.29 39879 576
Araneae Antrodiaetusriversi MY1253 360 A 1240653 1171503 94.43 28151 571
Araneae Calisogasp. MY4085 500 A 1283059 1209706 94.28 37800 664
Araneae Diguetiasignata G0580 500 C 2109935 1988677 94.25 76172 504
Araneae Euagruschisoseus MY4295 500 A 1672688 1592368 95.20 46272 674
Araneae Habronattustarsalis HA0955 500 C 1016446 957507 94.20 30683 484
Araneae Hebestatistheveneti MY4392 500 A 1327914 1257588 94.70 34651 722
Araneae Hypochiluspococki H0595 500 C 1143368 1081293 94.57 32241 552
Araneae Kukulcaniasp. G0551 500 C 1032516 984284 95.33 32448 466
Araneae Liphistiusmalayanus MIS1 500 A 832843 781359 93.82 26169 557
Araneae Megahexurafulva MY4378 500 A 1680701 1602923 95.37 60357 565
Araneae Nesticusgertschi N0806 456 B 1354456 1238198 91.42 19918 543
Araneae Oecobiusnavus GAR4 280 C 1596815 1523543 95.41 39711 539
Araneae Porrhothelesp. MY0857 500 A 1121499 1046279 93.29 33797 628
Opiliones Bishopellalaciniosa OP0569 500 D 1486862 1385611 93.19 42527 329
Opiliones Briggsusbilobatus OP3559 216 E 1148872 1093882 95.21 26325 406
Opiliones Briggsuspacificus OP3625 320 D 1355200 1267279 93.51 50366 399
Opiliones Fumontanadeprehendor OP0623 432 D 1265798 1174121 92.76 28841 356
Opiliones Leiobunumcalcar OP1089 500 D 1315511 1242038 94.41 47403 346
Opiliones Metanonychus OP3704 247 D 1313727 1227411 93.43 43275 377
n.nigricans
Opiliones Sabaconcavicolens OP1518 500 D 990843 907152 91.55 18060 293
Opiliones Siroacaroides OP3383 108 D 1260358 1179969 93.62 33373 369
Scorpiones Bothriuruskeyserlingi AR021 500 I 1238117 1183837 95.62 23438 500
Scorpiones Centruroidessculpturatus AR015 500 H 3913560 3750267 95.83 100622 706
Scorpiones Diplocentruspeloncillensis AR016 500 H 533530 507517 95.12 14225 402
Scorpiones Hadrurusarizonensis AR017 500 H 734122 700561 95.43 17206 565
Scorpiones Paravaejovisspinigerus AR018 500 H 1035158 990302 95.67 22872 478
Thelyphonida Mastigoproctusgiganteus AR011 500 F 926004 871445 94.11 22195 431
Mean 1315482.3 1240463.57 94.22 36064.33 487.13
DNAreferstostartingquantitythatwasprocessed.VoucherspecimensarehousedintheSanDiegoStateUniversityTerrestrialArthro-
podCollection(SDSU_TAC).
volume, which consisted of 15 lL of adapter-ligated Target Enrichment of Illumina Libraries v. 1.5 protocol
DNA, 1X KAPA HiFi HotStart ReadyMix and 5 lM of (http://ultraconserved.org/#protocols). Custom TruSeq
each Illumina TruSeq dual-indexed primer (i5 and i7) adapter blockers (Glenn et al. 2016) and standard
with modified 8-bp indexes (Glenn et al. 2016). Amplifi- MYbaits blockers were annealed to 147 ng/uL library
cationconditionswere98 °Cfor45 s,followedby10–12 pools, followed by hybridization to the master arachnid
cyclesof98 °Cfor15 s,60 °Cfor30 sand72 °Cfor60 s bait set (Faircloth 2016). Hybridizations were performed
and then a final extension of 72 °C for 60 s. After bead at65 °Cfor24 h.Afterhybridization,librarypoolswere
clean-up of amplified libraries, equimolar amounts of bound to Dynabeads MyOne Streptavidin C1 magnetic
libraries were combined into 1000 ng total pools. Pool beads(LifeTechnologies)forenrichment.Weperformed
combinations ranged from 1 to 8 individual libraries with-beadPCRrecoveryoftheposthybridizationenrich-
(Table 1). Due to a low concentration following library ments in a 50-lL reaction volume consisting of 15 lL
preparation and amplification, the pool containing sam- enriched DNA, 1X KAPA HiFi HotStart ReadyMix and
pleAR014(Neomolguslittoralis)includedaloweramount 5 lM each of TruSeq forward and reverse primers.
oflibrary(108.9 ng)forthisspecimencomparedwiththe Amplification conditions were 98 °C for 45 s, followed
othertwolibrariesthatweincluded(333 ngeach). by 18 cycles of 98 °C for 15 s, 60 °C for 30 s and 72 °C
Targetenrichmentoflibrarieswasperformedusinga for 60 s and then a final extension of 72 °C for 5 min.
MYbaits custom kit (MYcroarray, Inc.) following the Following PCRrecovery, libraries werequantified using
©2016JohnWiley&SonsLtd
4 J. STARRETT ET AL.
a Qubit fluorometer and diluted to 5 ng/lL. We per- GTRGAMMA model. ML analyses were conducted on
formed qPCR quantification of enriched library pools, threematricesforeachconcatenateddataset,whichcon-
and prepared a 10 lM mix of each pool at equimolar sisted of taxon coverages of 90%, 70% and 50% for each
ratios.Wesequencedthelibrarypoolusingapartialrun locus. Topologies and support scores based on 70% and
ofpaired-end150-bpsequencingonanIlluminaNextSeq 90% taxon coverage were nearly identical to those from
(GeorgiaGenomicsFacility). 50% taxon coverage, and we show matrix statistics for
only the 50% data sets and trees for the 50% and 90%
data sets. For Araneae, Opiliones and Scorpiones, we
Readprocessing,contigassemblyandmatrixcreation
conductedcoalescent-basedanalyseswithASTRALversion
RawreaddatawereprocessedusingthePHYLUCEpipeline 4.10.8 (Mirarab et al. 2014; Mirarab & Warnow 2015).
(Faircloth 2015). Adapter removal and quality control Here, trees for each locus were generated in RAXML with
trimming were conducted using the ILLUMIPROCESSOR 500 bootstrap replicates from the 90%, 70% and 50%
wrapper (Faircloth 2013) using default values. Reads taxoncoveragedatasets.Allanalyseswereconductedon
wereassembledusingTRINITYversionr2013-02-25(Grab- alate2015iMacwitha4-GHzInteli7processor,withthe
herret al.2011).Contigsfromallsampleswerematched exceptionofcontigassemblywithTRINITY,whichwasrun
toprobesusingminimumcoverageandminimumiden- ona12coreCentOSlinuxmachinewith48GBofRAM.
tity values of 65. We additionally extracted UCE loci in
silico from available arachnid genomes and Limulus
Results
polyphemus.UCElociwerealignedusingMAFFT(Katoh&
Standley 2013) and trimmed with GBLOCKS (Castresana Sequencingresults,assemblystatisticsandthenumberof
2000;Talavera&Castresana2007)asimplementedinthe UCElocirecoveredarepresentedinTable 1forsamples
PHYLUCE pipeline. Multiple data sets were created for sequencedinvitrointhisstudyandinTableS1(Support-
downstream analyses. First, a ‘2perArachnid’ data set ing information) for published genomes used in silico.
wascreatedthatcontaineduptotworepresentativesam- On average, we produced 1 315 482 rawreads per sam-
ples from each arachnid order sequenced, plus Limulus ple, with an average of 1 240 464 reads (94.2%) passing
as an outgroup. This included UCE data extracted from quality control. Assemblies resulted in an average of
both novel samples and previously deposited genomes. 35 884 contigs per sample. The number of UCE loci
Representative samples were chosen to span the root recovered from newly sequenced samples varied
node of the relevant arachnid orders and based on the between 170 and 722 (average = 487), while the average
number of UCE loci recovered. Because the purpose of recovered in silico from published genomes was 675,
our study was not to reconstruct arachnid phylogeny, including 555 UCEs from the outgroup Limulus. The
we did not include all arachnid orders. Second, we average number of recovered UCE loci differed among
assembled a ‘UCEsample’ data set that included only orders (unpaired two-tailed t-test, t = 6.4308, df = 38,
samples newly sequenced for this study. Third, three P-value<0.0001),beinghigheringroupsfromwhichthe
individual data sets were created that included all sam- probes were designed (Araneae, Parasitiformes, Scorpi-
plesfromwithintheordersAraneae,OpilionesandScor- ones: 605 loci) versus those that were not used in probe
piones.Foreachoftheseindividual-order datasets,two design (Acariformes, Amblypygi, Opiliones, Thely-
matrices were created, one including the amblypygid phonida: 376 loci). We tested whether the number of
Damon as the outgroup and a second without an out- samples included in a hybridization pool influenced the
group. Data matrices without outgroups were used for numberofUCElocirecoveredandfound nocorrelation
determining matrix statistics, with the number of parsi- (R2 = 0.03),althoughtheOpiliones-onlypoolsrecovered
mony-informative characters computed using PAUP* 4.0 the fewest loci (Fig. S1, Supporting information). Matrix
(Sinauer Associates, Inc.). Finally, to assess species-level statisticsarepresentedfora50%taxoncoveragedataset
utility of UCEs, three congeneric data sets were created inTable 2.Werecoveredidenticaltopologiesandnearly
for Antrodiaetus (turret spiders), Briggsus (Briggs’ identicalsupportscoresformatriceswith70%taxoncov-
harvestmen)andCentruroides(barkscorpions). erage(notshown)despiteanapproximate37%decrease
inlocusnumber.Acrossthefinalmatrices,anaverageof
589UCElociwereincludedinthe50%matrices,thelow-
Phylogeneticanalysis
est being the Opiliones + OUT data sets. Matrix lengths
Data matrices including an outgroup taxon were also varied from 105 337 bp in the ‘UCEsample’ matrix to
subjected to phylogenetic analyses using RAXML HPC v8.0 452 309 bp in the Scorpiones matrix. The average per-
(Stamatakis 2014), implementing the rapid bootstrap centageofparsimony-informativecharacterswas28.99%.
algorithm (Stamatakis et al. 2008) plus ML tree For shallow timescale comparisons, Briggsus (Briggs’
search option (-f a), 200 bootstrap replicates and the harvestmen) resulted in the lowest number of UCE loci
©2016JohnWiley&SonsLtd
UTILITY OF AN ARACHNID UCE PROBE SET 5
Table2 Matrixstatistics
50%taxoncoverage
Dataset N Loci Length Meanlength Min–maxlength PI %PI
2perArachnid 11 602 199810 331.91 113–853 60468 30.3
UCEsample 32 510 105337 206.54 80–746 49362 46.9
Araneae+OUT 20 686 169134 246.55 92–833 67644 40.0
Araneae 19 724 181915 251.26 92–833 71786 39.5
Opiliones+OUT 9 435 158314 363.94 97–893 33263 21.0
Opiliones 8 381 153187 402.07 97–914 29807 19.5
Scorpiones+OUT 8 627 318178 507.46 145–1136 65227 20.5
Scorpiones 7 749 452309 603.88 99–1141 64892 14.4
Mean 589.3 217273 364.2 55306.1 29.0
PI,parsimony-informativesites.
recovered(Table 3;292loci,totallengthof172 562,aver- relatively few published arachnid genomes. Here, we
agelocus length 591),but contained the highest number demonstrate the utility of a UCE probe set targeting
of variable sites (25 524—14.8%). Conversely, the around 1000 loci that works for all arachnid species
Centruroides (bark scorpions) comparison recovered the tested (and likely also with chelicerate outgroups).
highest number of UCE loci (585 loci, total length of Across all samples, which include the deepest diver-
583 454,averagelocuslength997.4),butthelowestnum- gencewithinArachnida(Regieret al.2010;Sharmaet al.
berofvariablesites(7660—1.3%).Theproportionofvari- 2014), we recovered 510 UCE loci with 49 362 parsi-
able loci for all comparisons was above 0.9, with an mony-informative sites. For data sets within scorpions,
averagenumberofvariablesitesperlocusrangingfrom spiders and harvestmen, which each span root nodes
14to87(Table 3). estimatedathundredsofmillionsofyearsofdivergence
Maximum-likelihood analyses of concatenated data (Hedinet al.2012;Bondet al.2014;SharmaandWheeler
sets demonstrated the utility of the arachnid-specific 2014), we recovered 749, 724 and 381 loci, respectively.
UCE loci in resolving relationships among and within Wedidnotfindasignificantcorrelationbetweenlibrary
arachnid orders (Figs 1–3). For the limited deep-level poolsizeandaveragenumberoflocirecovered(Fig.S1,
sample data set (i.e. ‘2perArachnid’), bootstrap support Supporting information), andthus,the lower numberof
was 100% for all nodes except the node uniting Scorpi- loci recovered for Opiliones likely reflects the fewer
ones to Araneae + Pedipalpi (=Amblypygi + Thely- genomic resources for this group and their distant rela-
phonida) (Fig. 1). For the Araneae data set, most nodes tionship to the taxa from which the probes were devel-
were fully supported (BS = 100%), while two nodes oped.
within Entelegyne spiders had bootstrap values of 97%
(Fig. 2). Within Opiliones, all nodes received bootstrap
Phylogeneticutility
support ≥99% (Fig. 3), and within Scorpiones, all nodes
werefullysupportedexceptforonenodewithbootstrap Recent arachnid phylogenomic data sets based on
support of 69% (Fig. 3). For the ASTRAL analyses of transcriptomes have provided better-resolved and well-
Araneae, Opiliones and Scorpiones, most relationships supportedphylogeniescomparedwithpriormorphologi-
were consistent across analyses of different matrix cal and genetic work using few loci (Hedin et al. 2012;
completeness, and with concatenated results (Fig. S3, Bond et al. 2014;Fern(cid:1)andez et al. 2014; Sharma et al.2014,
Supporting information). We elaborate on notable dis- 2015;Garrsionet al.2016).However,obtainingsequencecover-
crepanciesinthediscussion. age across taxa from transcriptomes requires high-quality RNA
and a consistent expression pattern, which may be difficult to
obtain for many nonmodel taxa. Additionally, UCEs can be
Discussion
obtainedfrommuseumspecimensorotherpotentiallydegraded
samples(McCormacket al.2015;Blaimeret al.2016).
ArachnidUCEs
Our reconstructed phylogenies based on UCEs
The development of genetic markers for arachnids has (Figs 1–3, S2 and S3, Supporting information) demon-
lagged behind that of many groups due, in part, to the strate the utility of this probe set for resolving relation-
ancientdivergenceswithinthegroup(atleast400MYA; ships at a wide range of divergence levels, from the
Rehm et al. 2012; Rota-Stabelli et al. 2013) and the deepest bifurcations in the arachnid tree to shallower
©2016JohnWiley&SonsLtd
6 J. STARRETT ET AL.
Table3 Congenericmatrixstatistics
Total Mean Min–max Variable #Polymorphic Proportion Averagevariable
Dataset Loci length length length sites %variable loci ofpolymorphic sitespervariablelocus
Antrodiaetus 480 329110 685.65 129–1248 8389 2.55 476 0.99 17.62
Briggsus 292 172562 590.97 144–1190 25524 14.79 292 1.00 87.41
Centruroides 585 583454 997.36 243–1164 7660 1.31 527 0.90 14.54
Fig.1 Maximum-likelihood phylogeny based on 50% minimum taxon coverage ‘2perArachnid’ data set, including representative
imagesofincludedarachnidorders.Nodeshave100%bootstrapsupport,unlessotherwiseindicated.LettersA–Fcorrespondtorepre-
sentativeimagesofincludedarachnidorders.Photographcredit:MHedin(A,C,D),MErbland(B),WESavary(E,F).[Colourfigure
canbeviewedatwileyonlinelibrary.com].
©2016JohnWiley&SonsLtd
UTILITY OF AN ARACHNID UCE PROBE SET 7
Araneae
Damon AR010
Liphistius MIS1 Mesothelae
A
Megahexura MY4378
t
y
p
Antrodiaetus MY1253 o M
id y
Antrodiaetus AR17 e g
a a
lo
Euagrus MY4295 A m
v o
Porrhothele MY857 ic r
u p
Acanthoscurria la ha
PRJNA222716 rio e
Calisoga MY4085 id
e
a
Hebestatis MY4392
Hypochilus H595 HH
ya
pp
Kukulcania G551 olo
c
LGoCxAo_sc0e0l1e1s88405.1 hilidgyna Ara
0.06 Diguetia G580 aee ne
+ o
Habronattus HA955 m
o
Stegodyphus E r
n p
97 GCA_000611955.2 te h
Oecobius GAR4 le ae
g
97 y
Nesticus N806 n
a
Parasteotoda e
GCA_000365465.1
Latrodectus
GCA_000697925.1
Fig.2 Maximum-likelihoodphylogenybasedon50%minimumtaxoncoverageAraneae+OUTdataset.Nodeshave100%bootstrap
support,unlessotherwiseindicated.BoxesatnodesindicatenodesthatdisagreewiththeGarrsionetal.(2016)transcriptomestudy.
[Colourfigurecanbeviewedatwileyonlinelibrary.com].
divergences within genera. In our limited, deep-level phylogenies (Hedin et al. 2012; Bond et al. 2014; Fern(cid:1)an-
sampledataset(i.e.‘2perArachnid’),werecoveredwell- dez et al. 2014; Sharma et al. 2015; Garrsion et al. 2016)
supported monophyletic groups forallorders for which and a traditional PCR approach (Dimitrov et al. 2016).
more than one individual was included (Fig. 1). Addi- Within Araneae, relationships that are well supported
tionally, weobtained well-supportedsisterrelationships with transcriptome data (see Fig. 2 of Garrsion et al.
between Thelyphonida and Amblypygi (Pedipalpi), and 2016)arealsorecoveredwithUCEdatainbothconcate-
aPedipalpisistergrouptoAraneae,consistentwithrela- nated and coalescent analyses under different missing
tionships based on morphology, Sanger sequence data dataschemes(Fig. 2,S2andS3,Supportinginformation).
and transcriptome data (Giribet et al. 2002; Shultz 2007; These relationships and recovered clades include a split
Regieret al.2010;Sharmaet al.2014).Thelowsupportin betweentheAtypoideaandAvicularioideawithinmyga-
ourphylogenyfortheplacementofScorpionesmaybea lomorphs,andasplitbetweenHaplogynaeandEnteleg-
result of sparse sampling, or it may reflect the difficulty ynae within araneomorphs. We also recover a
indeterminingthepositionofthisgroupduetoapossi- well-supported Paleocribellate (Hypochilus) + Filistatidae
blerapidradiationdeepinthearachnidtreeoflife. (Kukulcania) sister group, a relationship only recently
Phylogenetic analyses of the UCE data for more hypothesized based on transcriptome data (Bond et al.
densely sampled groups (Araneae, Opiliones, Scor 2014; Garrsion et al. 2016). Opisthothelae (Mygalomor-
piones; Figs 2 and 3) also recover topologies that are phae + Araneomorphae) is recovered in all analyses
mostly congruent with recent transcriptome-based except the ASTRAL analysis of the 90% complete data set,
©2016JohnWiley&SonsLtd
8 J. STARRETT ET AL.
Opiliones Fig.3 Maximum-likelihood phylogenies
based on 50% minimum taxon Opil-
Damon AR010
iones+OUTandScorpiones+OUTdata
Eupnoi sets.Nodeshave100%bootstrapsupport,
Leiobunum OP1089 unless otherwise indicated. Boxes at
Dyspnoi Palpatores nodes indicate nodes that disagree with
Sabacon OP1518 recenttranscriptomestudies(Hedinetal.
2012; Sharma et al. 2015). [Colour figure
Siro OP3383 Cyphophthalmi canbeviewedatwileyonlinelibrary.com].
Bishopella OP569
99
99
L
Fumontana OP623 a
n
ia
Metanonychus OP3704 t
o
r
e
0.08 Briggsus OP3625 s
Briggsus OP3559
Scorpiones
Damon AR010
Mesobuthus

GCA_000484575.1 B
u
Centruroides th
GCA_000671375.1 id
a
Centruroides AR015
Diplocentrus AR016
Paravaejovis AR018 I
u
r
id
a
Hadrurus AR017
0.08
69
Bothriurus AR021
whichresultedinaMesothelae + Mygalomorphaesister before meaningful comparisons can be made between
relationshipwithlowsupport(Fig.S3,Supportinginfor- thesedifferentdataclasses.
mation).Entelegynerelationshipswerestableandhighly Within Opiliones, analyses of concatenated data
supportedacrossourconcatenatedanalysesunderdiffer- under all missing data schemes resulted in strong sup-
ent missing data schemes (Figs 2, S2, Supporting infor- port for a Palpatores clade (Eupnoi + Dyspnoi) and a
mation). These relationships were not consistent with Laniatoresclade(Figs 3,S2andS3,Supportinginforma-
those obtained from coalescent analyses where relation- tion). These clades have also received strong support in
shipshadlowsupportandwerenotconsistentacrossthe analysesoftranscriptomedata(Hedinet al.2012;Sharma
50/70 and 90% complete data sets (Fig. S3, Supporting et al.2014)andafive-genedatasetcombined withmor-
information). Support is generally lower for deep nodes phological characters (Sharma & Giribet 2014). With a
within Entelegynae than that observed for much of the single exception, all missing data schemes for both con-
restofthespiderphylogenybasedontranscriptomedata catenated and coalescent analyses of UCEs resulted in
(Garrsion et al. 2016). The conflict in relationships and strong support for a sister relationship between
low support in both UCEs and transcriptomes indicate Cyphophthalmi (Siro) and Laniatores. This relationship
that more comprehensive taxon sampling is needed conflicts with the well-accepted Phalangida clade
©2016JohnWiley&SonsLtd
UTILITY OF AN ARACHNID UCE PROBE SET 9
(Palpatores + Laniatores), which has received strong
support based on transcriptome data (Hedin et al. 2012; es 0.6
s UCEsample
Sharma et al. 2014). However, the 90% ASTRAL analysis ba 0.5
doesrecoverPhalangidawithlowsupport(Fig.S3,Sup- ant
porting information). Themisplacement ofSirodoes not vari 0.4
appeartobeanartefactofmissingdatagiventheconsis- y of 0.3
c
tency and strong support across our missing data n 0.2
e
u
schemes.Wesuggestthatthesparsetaxonomicsampling q
e 0.1
andthedistantavailableoutgrouparelikelyculpritsfor Fr
0.0
theconflictinrelationships.
−400 −200 0 200 400
Within Scorpiones (Figs 3, S2 and S3, Supporting
Distance from center of alignment (bp)
information),phylogeneticanalysesofUCEsrecoverthe
same split between parvorders Buthida and Iurida as 0.5
s
inferred from transcriptome data (Sharma et al. 2015). ase 0.4 Antrodiaetus riversi (10.76%)
However,thewell-supportedpositionofBothriurusdeep nt b
within Iurida in our phylogeny conflicts with its more aria 0.3
basalpositioninIuridabasedontranscriptomeanalyses v
of
(Sharma et al. 2015). Notably, the 90% and 50% concate- y 0.2
c
nated trees yield identical topologies and similar sup- en
u 0.1
port, indicating that the relationships and support are eq
Fr
robust to missing data when using concatenation (or at
0.0
leastnotanartefactofmissingdata).AllASTRALanalyses
produceanidenticaltopology(Fig.S3,Supportinginfor- 0.5
mation), which differs from the concatenated phyloge- s
e
s
nies in the relationships within Iurida. The short ba 0.4 Briggsus spp. (14.81%)
internodes, likely due to an ancient and rapid radiation, nt
a
account for these incongruences between data sets and ari 0.3
v
previousanalyses. of
y 0.2
Comparison of UCE sequences across close relatives nc
e
within spiders, harvestmen and scorpions indicates qu 0.1
e
that our probe set will have applications beyond deep- Fr
0.0
level phylogenetics. The UCE loci show promise for
−500 0 500
use in species-level phylogenetics and species delimita-
Distance from center of alignment (bp)
tion because the flanking regions show increasing vari-
ability as distance from the core UCE increases
s
e 0.025
(Fig. 3). UCE enrichment produced 480 loci with more s
a Centruroides spp. (9.41%)
b
than 8000 variable sites across two members of the nt 0.020
a
Antrodiaetus riversi complex, 292 loci with more than ari
25 000 variable sites across Briggsus species and 585 of v 0.015
loci with more than 7000 variable sites across two Cen- ncy 0.010
truroides species (Table 3). These patterns of UCE vari- ue
q
ability are correlated with divergence in mitochondrial Fre 0.005
cytochrome oxidase I (Fig. 4). The proportion of poly-
0.000
morphic loci (minimum 0.9) and average variable sites
−600 −400 −200 0 200 400 600
per polymorphic locus (minimum 14) seen across these
Distance from center of alignment (bp)
pairs is greater than those reported in the study of
Smith et al. (2014; maximum 0.77 and 3.2, respec- Fig.4 Variability increases as distance from the core UCE
tively), which demonstrated the utility of UCEs for increases. Data points with no variability were removed. For
species-level phylogenetic and species delimitation Antrodiaetus(turretspiders)andBriggsus(harvestmen)compar-
isons,thescaleofx-axisisidentical.Numbersinparenthesesare
analyses in birds. UCEs may be advantageous for shal-
uncorrected COI divergence values. In Antrodiaetus, an outlier
low-level studies in Arachnida because other tech-
wasremovedforbetterviewing.CytochromeoxidasesubunitI
niques, such as RADseq, are susceptible to locus
(COI)distancesforCentruroides(barkscorpions)werecalculated
dropout in arachnid species with relatively deep levels fromGenBankAccessionnos.AY995831.1(C.sculpturatus)and
of divergence (Bryson et al. 2016; Derkarabetian et al. AY995833.1(C.exilicauda).
©2016JohnWiley&SonsLtd
10 J. STARRETT ET AL.
2016). Our UCE probe set will also provide a comple- ChapmanAD(2009)NumbersofLivingSpeciesinAustraliaandtheWorld,
mentary/alternative resource to the recently developed 2ndedn.AustralianBiodiversityInformationServices,Canberra.
Crawford NG, Faircloth BC, McCormack JE, Brumfield RT, Winker K,
spider anchored-enrichment loci (Hamilton et al.
GlennTC(2012)Morethan1000ultraconservedelementsprovideevi-
2016b), which were used to produce a comparable dencethatturtlesarethesistergroupofarchosaurs.BiologyLetters,8,
number of loci (455) and a well-resolved phylogeny 783–786.
for the North American tarantula genus Aphonopelma CrawfordNG,ParhamJF,SellasAB,etal.(2015)Aphylogenomicanaly-
sisofturtles.MolecularPhylogeneticsandEvolution,83,250–257.
(Hamilton et al. 2016a).
DerkarabetianS,BurnsM,StarrettJ,HedinM(2016)Populationgenomic
Thearachnid-specificUCEshaveutilityandphyloge- evidence for multiple refugia in montane-restricted harvestmen
neticinformativenessatalllevelsofArachnida,spanning (Arachnida, Opiliones, Sclerobunus robustus) from the southwestern
UnitedStates.MolecularEcology,25,4611–4631.
extremely ancient divergence times between orders (at
DertiA,RothFP,ChurchGM,WuC(2006)Mammalianultraconserved
least 400 MYA) to more recent congeneric divergences
elements are strongly depleted among segmental duplications and
(<10 MYA for Antrodiaetus, Hedin et al. 2013). Thus, copynumbervariants.NatureGenetics,38,1216–1220.
thesemarkershavethepotentialtobeanintegralcompo- Dimitrov D, Benavides LR, Arnedo MA etal. (2016) Rounding up the
usual suspects: a standard target-gene approach for resolving the
nent of future comparative studies. Arachnids are an
interfamilialphylogeneticrelationshipsofecribellateorb-weavingspi-
extremelydiversegroup,andyetitisestimatedthatless
ders with a new family-rank classification (Araneae, Araneoidea).
than half of arachnid species have been formally Cladistics,doi:10.1111/cla.12165(inpress).
described (Chapman 2009).Species delimitation increas- DunlopJA(1996)AtrigonotarbidarachnidfromtheUpperSilurianof
Shropshire.Palaeontology,39,605–614.
ingly relies on phylogenomic data that can be sampled
DunlopJA,TetlieOE,PrendiniL(2008)ReinterpretationoftheSilurian
consistentlyacrossdiversetaxa(Leach(cid:1)eet al.2014;Smith scorpionProscorpiusosborni(Whitfield):integratingdatafromPalaeo-
et al. 2014; Rannala 2015; Harvey et al. 2016). The recov- zoicandRecentscorpions.Palaeontology,51,303–320.
ery of homologous UCEs across Arachnida makes these FairclothBC(2013)Illumiprocessor:aTrimmomaticwrapperforparallel
adapterandqualitytrimming.doi:10.6079/J9ILL.
lociavaluableresourcefordiscoveryofnewspeciesand
FairclothBC(2015)PHYLUCEisasoftwarepackagefortheanalysisof
for inferring phylogenetic relationships in understudied conservedgenomicloci.Bioinformatics,32,786–788.
arachnidgroups. FairclothBC(2016)Identifyingconservedgenomicelementsanddesign-
inguniversalprobesetstoenrichthem.bioRxivpreprintdoi:10.1101/
077172,bioRxiv077172.
Acknowledgements FairclothBC,McCormackJE,CrawfordNG,HarveyMG,BrumfieldRT,
GlennTC(2012)Ultraconservedelementsanchorthousandsofgenetic
For help during laboratory work, we thank Carl Oliveros and markersspanningmultipleevolutionarytimescales.SystematicBiology,
Rich Harrington. John Klicka provided other assistance. Tim 61,717–726.
Burkhart, Ken MacNeil, Peter Clausen and Cor Vink provided FairclothBC,SorensonL,SantiniF,AlfaroME(2013)Aphylogenomic
perspectiveontheradiationofray-finnedfishesbasedupontargeted
assistancewithspecimencollection.Livespecimenphotographs
sequencingofultraconservedelements(UCEs).PLoSONE,8,e65923.
were provided by Mardon Erbland and Warren E. Savary.
doi:10.1371/journal.pone.0065923.
Funding for this research was provided by NSF Grant DEB-
FairclothBC,BranstetterMG,WhiteND,BradySG(2015)Targetenrich-
1354558 to MH. Startup funds from LSU and DEB-1242260 to mentofultraconservedelementsfromarthropodsprovidesagenomic
BCFsupportedthiswork. perspective on relationships among Hymenoptera. Molecular Ecology
Resources,15,489–501.
Fern(cid:1)andezR,HormigaG,GiribetG(2014)Phylogenomicanalysisofspi-
dersrevealsnonmonophylyoforbweavers.CurrentBiology,24,1772–
References
1777.
Garrsion NL, Rodriguez J, Agnarsson I etal. (2016) Spider phyloge-
AlexanderRP,FangG,RozowskyJ,SnyderM,GersteinMB(2010)Anno-
nomics:untanglingthespidertreeoflife.PeerJ,4,e1719.doi:10.7717/
tatingnon-codingregionsofthegenome.NatureReviewsGenetics,11,
peerj.1719.
559–571.
GilbertPS,ChangJ,PanCetal.(2015)Genome-wideultraconservedele-
BejeranoG,PheasantM,MakuninIetal.(2004)Ultraconservedelements
ments exhibit higher phylogenetic informativeness than traditional
inthehumangenome.Science,304,1321–1325.
genemarkersinpercomorphfishes.MolecularPhylogeneticsandEvolu-
BlaimerBB,LloydMW,GuilloryWX,BradySG(2016)Sequencecapture tion,92,140–146.
andphylogeneticutilityofgenomicultraconservedelementsobtained
GiribetG,EdgecombeGD,WheelerWC,BabbittC(2002)Phylogenyand
frompinnedinsectspecimens.PLoSONE,11,e0161531.doi:10.1371/
systematic position of opiliones: a combined analysis of chelicerate
journal.pone.0161531.
relationshipsusingmorphologicalandmoleculardata1.Cladistics,18,
BondJE,GarrisonNL,HamiltonCA,GodwinRL,HedinM,AgnarssonI 5–70.
(2014) Phylogenomics resolves a spider backbone phylogeny and
GlennTC,NilsenR,KieranTJetal.(2016)AdapteramaI:universalstubs
rejectsaprevailingparadigmfororbwebevolution.CurrentBiology,
andprimersforthousandsofdual-indexedIlluminalibraries(iTru&
24,1765–1771.
iNext).bioRxiv,doi:10.1101/049114.
Bryson RW, Savary WE, Zellmer AJ, Bury RB, McCormack JE (2016)
GrabherrMG,HaasBJ,YassourMetal.(2011)Full-lengthtranscriptome
Genomic data reveal ancient microendemism in forest scorpions
assembly from RNA-Seq data without a reference genome. Nature
across the California Floristic Province. Molecular Ecology, 25, 3731– Biotechnology,29,644–652.
3751.
Hamilton CA, Hendrixson BE, Bond JE (2016a) Taxonomic revision
Castresana J (2000) Selection of conserved blocks from multiple
of the tarantula genus Aphonopelma Pocock, 1901 (Araneae, Mygalo-
alignmentsfortheiruseinphylogeneticanalysis.MolecularBiologyand
morphae, Theraphosidae) within the United States. ZooKeys, 560,
Evolution,17,540–552. 1–340.
©2016JohnWiley&SonsLtd