Ncomms7145 PDF
Ncomms7145 PDF
Ncomms7145 PDF
Received 16 Aug 2014 | Accepted 11 Dec 2014 | Published 4 Feb 2015 DOI: 10.1038/ncomms7145 OPEN
1 State Key Laboratory of Veterinary Etiological Biology, Key Laboratory of Veterinary Parasitology of Gansu Province, Lanzhou Veterinary Research Institute,
Chinese Academy of Agricultural Sciences, Lanzhou 730046, Gansu Province, China. 2 Faculty of Veterinary and Agricultural Sciences, The University of
Melbourne, Victoria 3010, Australia. 3 BGI, Shenzhen 518083, China. 4 Department of Veterinary Disease Biology, University of Copenhagen, Copenhagen
2200, Denmark. 5 Institute for Parasitology und Tropical Veterinary Medicine, Freie Universität Berlin, Berlin 14163, Germany. 6 Department of Biochemistry
and Molecular Biology, Monash University, Victoria 3800, Australia. 7 Cancer and Stem Cell Biology, Duke-NUS Graduate Medical School, Singapore 138672,
Republic of Singapore. 8 Genome Institute of Singapore, 60 Biopolis Street, Singapore 138672, Republic of Singapore. 9 Structural Chemistry Program, Eskitis
Institute, Griffith University, Brisbane 4111, Queensland, Australia. 10 HHMI, Division of Biology, California Institute of Technology, Pasadena 91125, California,
USA. * These authors contributed equally to this work. Correspondence and requests for materials should be addressed to X.-Q.Z. (email:
[email protected]) or to P.K.K. (email: [email protected]) or to R.B.G. (email: [email protected]).
P
arasitic worms have a devastating, long-term impact on Table 1 | Summary of the features of the Toxocara canis draft
human health worldwide. For instance, approximately two genome.
billion people are infected with soil-transmitted helminths,
such as hookworms (Necator and Ancylostoma), whipworm
Description
(Trichuris) and the large roundworm (Ascaris), principally in
Genome size (bp) 317,115,901
under-privileged areas of Latin America, Africa and Asia1. The
Number of scaffolds; contigs 22,857; 51,969
global disease burden caused by these parasites is comparable N50 (bp); count42 kb in length 375,067; 2,899
with that of tuberculosis and malaria in disability-adjusted life N90 (bp); count4N90 length 66,363; 938
years1. Ascaris, for example, infects more than one billion people, Genome GC content (%) 40.0
causing nutritional deficiency, impaired cognitive and physical Repetitive sequences (%) 13.5
development, principally in children and, in severe cases, death2. Exonic proportion; including introns (%) 6.8; 49.4
Published information3,4 also shows that a related parasite, Number of putative coding genes 18,596
Toxocara canis (Werner, 1782), infects millions of people in Mean gene size (bp) 8,416
poverty-stricken parts of the USA alone. Toxocariasis, the disease Mean CDS length (bp) 1,156
caused by Toxocara spp., is highly prevalent in many developing Mean exon number per gene 7.4
Mean exon length (bp) 156
countries, where its importance is likely to be seriously
Mean intron length (bp) 1,133
underestimated. Toxocariasis results from zoonotic transmission Coding GC content (%) 47.4
of Toxocara spp. from carnivores, including canids and felids5,6. CEGMA completeness: complete; partial (%) 67.3; 98.4
T. canis of canids is recognized as the main causative agent of
zoonotic disease; this species has a complex life cycle, which can
also involve paratenic hosts such as rodents. In humans,
particularly children, T. canis larvae invade various tissues T. canis
to cause visceral larva migrans, ocular larva migrans, 6037 C. elegans
neurotoxocariasis (including eosinophilic meningoencephalitis) 155
B. malayi
and/or covert toxocariasis7. In addition, clinical and experimental
studies8,9 have indicated an association between T. canis infection P. pacificus 174 59
and allergic disorders, such as asthma, chronic pruritus and 225 37 108 143
A. suum
urticaria. Both the canine and mouse models of T. canis10,11
represent important tools for studies of the biology of the 1925 168 257 5918 1112
parasite, parasite–host interactions and toxocariasis at the 771 1071 436
immuno-molecular level.
To facilitate such investigations, we sequenced and annotated
the 317 Mb draft genome of T. canis and compare it with other Figure 1 | Gene orthology. Venn diagram showing the number of orthologs
nematode genomes. This genome contains at least 18,596 protein- between Toxocara canis and four other nematode species, Ascaris suum,
coding genes, whose predicted products include at least 373 Brugia malayi, Caenorhabditis elegans and Pristionchus pacificus upon pairwise
peptidases, 458 kinases, 408 phosphatases, 273 receptors and comparison.
530 transporters and channels. Notably, the T. canis secretome
(870 molecules) is rich in peptidases proposed to be associated (n ¼ 16,073, 8,550 and 7,696, respectively). We also identified
with the penetration and degradation of host tissues, and an 61 families of DNA transposons (a total of 45,249 distinct
assemblage of molecules suggested to modulate or suppress host- sequences), of which CMC-EnSpm, Novosib and MULE-MuDR
immune responses. This genome should provide a useful resource (n ¼ 19,762, 2,787 and 2,501, respectively) predominated. This
to the scientific community for a wide range of post-genomic richness of families of transposable elements is comparable with
studies of T. canis and could support the development of that of genomes of other parasitic nematodes12,14,15.
new interventions (drugs, vaccines and diagnostics) against
toxocariasis and related nematodiases.
The gene set of T. canis and comparison with other nematodes.
Using transcriptomic data from adult and larval stages of T. canis,
Results de novo predictions and homology-based searching, we identified
Genome assembly and repeat content. We sequenced the 18,596 coding genes, with mean gene, exon and intron lengths of
T. canis genome at 110-fold coverage and produced a final draft 8,416, 156 and 1,133 bp, respectively, and an average of 7.4 exons
assembly of 317 Mb (N50 ¼ 375 kb; Table 1) with a mean GC per gene (Table 1), similar to the findings for our Ascaris gen-
content of 40.0% (5.8% Ns). We detected 98.4% of the 248 core ome12. Compared with the draft genome sequences of A.
eukaryotic genes by Core Eukaryotic Genes Mapping Approach suum12,13, Brugia malayi14, Caenorhabditis elegans16 and
(CEGMA), indicating that the assembly represents a substantial Pristionchus pacificus15, on average, the T. canis genes have
proportion of the entire genome. We estimated a repeat content most sequence similarity to those of A. suum and are significantly
in this draft genome of 13.5% (equating to 42.9 Mb of DNA), longer than in the genomes of the other four nematodes, likely
comprising 1.8% DNA transposons, 4.1% retrotransposons, relating to an increased number and length of exons
1.4% unclassified dispersed elements and 6.4% simple repeats (Supplementary Data 2). Most (67.5%) of the predicted T. canis
(Supplementary Data 1); the overall repeat content is similar genes (Fig. 1) have an ortholog (BLASTp cut-off: 10 5) in A.
to that of the germline genome of the related worm Ascaris suum (n ¼ 11,658; 62.7%), B. malayi (8,395; 45.1%), C. elegans
suum, but considerably higher than its somatic (diminuted) (9,002; 48.4%) or P. pacificus (7,968; 42.8%). A total of 5,918
genome12,13. We identified 72,862 distinct retrotransposon genes are orthologous among all 5 species, 3,557 are shared with
sequences (see Supplementary Data 1) representing at least 68 at least 1 other species of nematode but absent from C. elegans
families (16 LTR, 32 LINE and 20 SINE), with Gypsy, Pao and and 1,925 are shared exclusively with Ascaris to the exclusion of
Copia predominating for LTRs (n ¼ 12,436, 4,859 and 1,302, the other 3 nematodes (Fig. 1). Conversely, 6,037 genes (32.5%)
respectively) and CR1, RTE-RTEand L2 for non-LTRs are unique to T. canis relative to the other 4 nematodes (Fig. 1).
Of the entire T. canis gene set, 5,406 genes (29.1%) have an Also of interest in this context is the panel of channel, pore and
ortholog (r10 8) linked to known biological pathways (Kyoto transporter proteins that we identified here, particularly con-
Encyclopaedia of Genes and Genomes (KEGG); Supplementary sidering that many common anthelmintics bind representatives
Data 3). Mapping to pathways in C. elegans suggested a near of some of these proteins as targets19. We predicted 156 GPCRs
complete complement of KEGG orthology groups (90.6%), also to be encoded in T. canis; these include 111 class A rhodopsin
supporting the CEGMA results. By inference, most T. canis genes family (for example, neuropeptide, serotonin, acetylcholine,
are represented in the present genomic assembly; at least 17,208 dopamine and adrenaline receptors), 21 class B secretin
of all the 18,596 genes predicted here are supported by extensive receptor family (for example, parathyroid hormone, nematode
transcriptomic and/or inferred proteomic data. chemoreception and latrophilin receptors), 10 class C
metabotropic glutamate/pheromone family (for example,
glutamate and gamma-aminobutyric acid (GABA) receptors)
Functional annotation and protein classifications. In total,
and 14 other (for example, Wnt) receptors (Supplementary
14,583 (78.4%) of the protein sequences predicted from the
Data 5). In addition, of the 530 transporters predicted here, we
18,596 coding genes of T. canis were annotated based on the
identified an abundance of those of the solute carrier family
presence of characteristic protein domains (Supplementary
(67.2%), major facilitator superfamily, ABC transporters (10.8%;
Data 4); 12,346 (66.4%) matches were found in C. elegans,
including P-glycoproteins) and major facilitator superfamily
and 10,977 (59.0%), 10,556 (56.8%) and 13,905 (74.8%) had
(3.4%)20,21. In addition, we identified 268 ion-channel proteins
significant matches in the InterProScan, Swiss-Prot and KEGG
(Supplementary Data 5), the majority of which represented the
databases, respectively. Of the 10,977 InterPro matches, 9,121
Cys-loop superfamily (including nicotinic, glycine, aniononic and
(83.1%), 10,494 (95.6%), 1,827 (16.6%) and 412 (3.8%) were in
GABA-A; 39.6%), as well as voltage-gated cation (mainly K þ
Pfam, PANTHER, PRINTS and PIRSF databases, respectively.
channels, including SLO-1; 31.0%) and others related to voltage-
According to the KEGG BRITE hierarchy, predicted proteins
gated cation channels (transient receptor potential; 10.8%), some
represented peptidases (n ¼ 373), kinases (458), phosphatases
of which are known targets for anthelmintic drugs22.
(408), receptors (273), transporters (530), GTPases (127), ion
channels (268) and transcription factors (355); (Supplementary
Data 5), with some proteins inferred to have multiple functions.
Stage-, gut- and gender-enriched transcription profiles. Here
we defined clusters of genes that are significantly differentially
Key enzymes, channels, pore and transporters. We identified transcribed between the T. canis adult and third-stage larvae (L3s
five main classes of peptidases (metallo-, cysteine, serine, threo- from fully embryonated eggs), and between the two sexes of this
nine and aspartic), with the metallo- (n ¼ 165; 44.2%), cysteine nematode, by integrating data from KO enrichment, as well as co-
(107; 28.7%) and serine (60; 16.1%) proteases predominating transcription and functional network analyses (Supplementary
(Supplementary Data 5). The most abundant families in these Data 6–9). In the T. canis life cycle, dioecious adults establish and
classes are the M12 astacins and adamalysins (n ¼ 57), M01 live in the small intestine of the canid host following hepato-
aminopeptidases (19) and M13 neprilysins (13) among the pulmonary migration of parasitic L3s11. By comparison to L3, we
metallo proteases; the C19 ubiquitin-specific proteases (n ¼ 30), showed that pathways enriched in the adult stage involve carbon
C01 papains (that is, cathepsins; 16) and C02 calpain-like metabolism (for example, aldo-2, c04c3.3, enol-1, fbp-1, gpd-3 and
enzymes (16) among the cysteine proteases; and the S01 ‘chy- gpi-1), lysosome activity, transport (ABC transporters),
motrypsins’ (n ¼ 12), S09 prolyl oligopeptidases (12) and S08 xenobiotic metabolism and DNA replication/repair, likely
subtilisins (9) among the serine proteases. These secreted pepti- linked to digestive and reproductive processes in the worm
dases (for example, the M12 metallo, the C01 and C02 cysteine, as (Fig. 2 and Supplementary Data 4, 6 and 7). With reference to the
well as the S01, S08 and S09 serine proteases) are of major adult stage, L3-enriched pathways were mostly linked to neuronal
interest, given their presence in the excretory/secretory (ES) signalling (for example, ckr-2, dop-1, gar-2 and unc-49), and
products of many parasitic helminths, and their central roles in cuticle formation and/or shedding (including col-121, dpy-5 and
tissue degradation and invasion (for example, during migration sqt-1; Fig. 2 and Supplementary Data 4, 6 and 7) which might
and/or feeding) and/or in immune evasion or modulation17,18. ES relate to ‘‘fuzzy coat’’ shedding from L3s during immune attack
peptidases (including aminopeptidases and/or cysteine proteases) in vitro23,24. Interestingly, homologues of various genes known to
are likely to play key roles in these processes in T. canis and might regulate dauer in C. elegans were upregulated (for example, daf-9,
represent important drug or vaccine targets. daf-12, osm-3 and osm-6) or downregulated (two paralogs of mek-
We also identified 458 protein kinases and 408 phosphatases in 2) in the L3 stage25,26 (Fig. 2), which supports the proposal that
T. canis (Supplementary Data 5). The kinome includes a large this stage is in an arrested state of development following the
portion of serine/threonine protein (67.2%) and tyrosine (13.3%), in vitro-hatching of embryonated eggs and maintenance in
as well as a small number of atypical or unclassified kinases culture in the absence of host factors24.
(19.4%). The phosphatome includes predominantly serine- In addition, in the worm, we showed that the intestinal tract of
threonine (81.1%), protein tyrosine (12.7%) and a minority of T. canis was specifically enriched for gene transcripts linked to
other phosphatases. On the basis of the homology to molecules in molecular uptake and degradation via lysosomes (for example,
the KEGG orthology (KO) database, we found 127 GTPases to be vha-2, -4, -6, -8 and -12), the transportation of amino acids,
encoded in the T. canis genome, including 29 large (hetero- sugars, lipids and drugs, the metabolism of xenobiotic metabo-
trimeric) and 98 small (monomeric) G-proteins representing the lism (for example, ugt-22, -34, -45, -49, -63 and cyp14a5), as well
Rab (n ¼ 42), Ras (19), Arf/Sar (16), Rho (16) and Ran (1) as protein digestion (aspartic, cysteine and metallo-peptidases)
families, as well as some unclassified molecules. Examples of these and fatty acid elongation (for example, elo-2, -3, -4, -5, -6 and -9;
include C. elegans homologues eft-1 (gene Tcan_11808), fzo-1 Fig. 2 and Supplementary Data 4, 6 and 7). These processes are
(Tcan_14313), glo-1 (Tcan_03008) and rho-1 (for example, consistent with extensive digestive, absorptive and detoxifying
Tcan_13740; cf. Supplementary Data 5), which play essential functions within the gut of T. canis and similar, in some respects,
roles in embryonic, larval and/or reproductive development and to those suggested for Ascaris27.
are also found in Ascaris12. Therefore, some of these enzymes Given the size of adult T. canis (usually 5–15 cm in length),
might be targets for anti-parasite interventions. we were able to investigate differential transcription between
Lysosome
homologues in any other organisms for which genomic or
Germ cell proliferation
Molecule transporters Embryogenesis
transcriptomic data are publicly available.
Xenobiotics Neuronal signalling
Proteases (1781)
Fatty acid biosynthesis
(7631)
The secretome and molecules for parasite–host interactions.
Proteasome system
Spermatogenesis
We predicted the secretome of T. canis to comprise at least 870
Energy metabolism ES proteins with a diverse array of functions (Supplementary
(1979) Data 10). Most notable are at least 14 peptidases, including 7 SC
Carbon metabolism serine proteases (S1, S9 and S28 families), 3 aspartic proteases (A1
Adults
os 2
os 3
m 6
-2
Lysosomes family), 3 CA/CD cysteine proteases (C1 and C13 families) and 1
da 9
f -1
-
-
m
m
ek
f-
da
ABC transporters
Xenobiotics 1 MA metallo-protease (M14 family), as well as 23 cell-adhesion
2
DNA replication/repair 3 molecules (immunoglobulins, integrins and cadherins), 17 lectins
4
(2656) 1 (C-type) and 6 SCP/TAPS proteins (venom allergen homo-
2
3 logues)28. Many secreted peptidases (representing the
4
1
‘degradome’) and their respective inhibitors have known roles
Neuronal signalling
2
3
in the penetration of tissue barriers and feeding by parasitic
Cuticle formation/shedding worms, including T. canis5,23.
(2833) T. canis larvae invade various tissues, including muscles, brain
C-type and eyes, and cause clinical disease. Such larvae have an
Lectins
Mucins Dauer exceptional ability to evade or block host attack and can survive
Aquaporin L3 for many years in tissues. This ability is associated with the
ANTs deployment of molecules excreted or secreted by the parasite or
released from its surface coat23,24. Here, we predicted 33 proteins
Figure 2 | Stage-, gut- and gender-enriched transcription. Key biological involved in host interactions and/or modulating host immune
pathways or processes in Toxocara canis, for which the gene transcription is responses (Supplementary Data 11). Abundant in O-linked
enriched in adult females (red) and males (blue); in adults and third-stage glycosylations are mucins (n ¼ 7), many of which are likely
larvae (L3; grey); and in the alimentary tract of the adult worms (black box, heavily targeted by IgM antibodies and bound by various pattern-
top left). Molecules enriched in L3s and involved in parasite–host recognition receptors associated with host dendritic cells,
interactions, such as immunomodulation, are also indicated (black box, responsible for inducing a Th2 immune response18.
bottom left; ANT, abundant new transcripts); the number of genes with Immunomodulatory molecules predicted in T. canis also
significantly increased transcription are given in parentheses. The proposed include homologues (E-value cut-off ¼ 10 8) of the B. malayi
arrested or dauer status of the L3 stage of T. canis (from embryontated cystatin CPI-2 (B-cell inhibitor), several transforming growth
eggs) is supported by the transcription profiles (heat-map, right) for gene factor-b and macrophage initiation-factor mimics, numerous
homologues daf-9, daf-12 and osm-3, osm-6 (upregulation; yellow) and two neutrophil inhibitors, oxidoreductases, known to counteract the
mek-2 paralogs (downregulation; orange to red), which are all known to neutrophil oxidative burst, and five close homologues of platelet
promote dauer in C. elegans; the biological replicates (1–4) studied are anti-inflammatory factor-a18. Key examples of T. canis ES
indicated to the right of the heat-map. proteins predicted to be involved in immune evasion include
some ‘hidden’ antigens29 (for example, numerous C-type lectins)
with close homology to vertebrate macrophage mannose or CD23
individual female and male worms. Genes with female-enriched (low affinity IgE) receptors, which mimic host molecules18.
transcription were linked to signal transduction pathways Other representatives include mucins (originating from the
involving mainly neuroactive GABA glycine and acetycholine parasite’s surface coat), phosphatidylethanolamine-binding
receptors, as well as G protein-coupled receptors associated proteins, cathepsins, asparagyl endopeptidases (legumains),
with germ cell proliferation and embryogenesis (Fig. 2 and superoxide dismutases, SCP/TAPS28, olfactomedin, aquaporin,
Supplementary Data 4, 6 and 7). Female-enriched transcription prohibitin and various orphan proteins (including abundant new
related to genes encoding glycosyltransferases of N-acetyl transcripts (ANTs)), identified previously in small-scale
glucose/galactosamine moieties linked to egg shell synthesis in molecular studies of T. canis23,24. Interestingly, although not
oogenesis (for example, gene codes C36A4.4, F07A11.2 and encoded in the genome, ants-3, -5, -30 and -34 were transcribed
F21D5.1), oocyte maturation and embryonic development (for at very high levels in all tissues of one of four female T. canis
example, ceh-13, elt-3, gsk-3, kgb-1, nrh-1, pha-4, rab-11.1, unc-60, studied here and at high levels in all L3 samples, but were
-130 and vab-2), egg laying (for example, lin-29, nsy-1, rhgf-1, transcribed at very low levels or not at all in tissues of the adult
sem-2 and sox-2) and vulva development (for example, arf-1.2, males (Supplementary Data 11). Three-dimensional modelling of
ceh-20, let-23, lin-1 and -61; Fig. 2 and Supplementary Data 4, all predicted ANTs showed that ANT-5 is an RNA-dependent
6 and 7). In contrast, genes with male-enriched transcription were RNA polymerase and that ANT-34 is a RNA helicase, consistent
inferred to be associated with protein degradation, spermatogen- with previous findings30. The absence of the ant genes from the
esis, epinephrine-like hormone regulation (tyramine and octopa- draft genome, the similarity in the structures of ANT-5 and
mine) and energy metabolism (Fig. 2 and Supplementary Data ANT-34 with viral proteins and the inconsistency in their
4, 6 and 7). Key genes with male-enriched profiles encoded transcription in stages, sexes and/or tissues of T. canis suggest
proteasome-related enzyme groups as well as protein-serine/ that ant transcripts are derived from one or more double-
threonine and -tyrosine kinases linked to sperm and spermato- stranded-RNA viruses. This proposal warrants investigation to
genesis. Conspicuous were those linked to the 26S proteasome establish whether ANTs are crucial to the biology of T. canis and
and associated recycling of ubiquitin moieties (for example, pas-1, indeed have a role in regulating transcription30. Overall, the
rpn-1, rps-26, rtp-3 and rpt-5) involved in germline development present genomic and transcriptomic data sets reveal that, on a
and spermatogenesis (for example, cpb-1, fog-3 and spe-6; Fig. 2 global scale, T. canis possesses a major arsenal of ES proteins that
and Supplementary Data 4, 6 and 7). Overall, of all 3,760 genes are involved in manipulating, blocking and/or evading immune
exhibiting gender-enriched transcription, B26% had no responses in host animals. A detailed understanding the roles of
these molecules could pave the way to new intervention strategies, tested and assessed in vivo in experimentally infected animals (for
such as vaccination. example, dogs or mice).
commercial kit (Chemagic STAR DNA Tissue Kit, Perkin Elmer) according to the RNA-seq-, EST- and homology-based gene sets were predicted, which were then,
manufacturer’s protocol. The DNA yield was measured spectrophotometrically together with de novo-predictions, subjected to analysis using the program GLEAN
(Qubit fluorometer dsDNA HS Kit, Invitrogen); DNA integrity was verified by to infer a reference set of genes for T. canis. UTR-regions were removed from
agarose-gel electrophoresis and using a Bioanalyzer (2100, Agilent). Paired-end RNA-seq-based, genome-guided and de novo-assembled transcriptomes using an
(PE) genomic libraries (with 170-, 500- and 800-bp inserts) as well as jumping (J) in-house pipeline. The prediction of the EST data set of T. canis was aided by the
genomic libraries (with 2-, 5-, and 10-kb inserts) were constructed according to NCBI EST database (8 September 2014) using the program BLASTn v.2.2.26.
manufacturer’s instructions (Illumina); to produce sufficient amounts of DNA for Homology-based comparisons were made with the proteomes of A. suum, B.
these latter libraries, 250–500 ng of genomic DNA were subjected to whole-genome malayi, C. elegans and P. pacificus. The reference gene set of T. canis was refined by
amplification using the REPLI-g midi kit (Qiagen), as recommended. All sequen- using gene-prediction models supported by (i) both homology and RNA-seq
cing was carried out on GA II or HiSeq (Illumina; 2 75 or 2 100 reads for evidence (420 mapped reads); (ii) homology-based models only, with no frame-
paired-end libraries, and 2 49 reads for J-libraries), and reads were exported to shift mutations and an alignment length of Z50% of the query sequence; and
the FASTQ format. (iii) RNA-seq models only, but encoding proteins of Z120 amino acids in length.
UTRs from RNA-seq data were added to the gene predictions using TopHat2
v.2.0.7. Finally, genes inferred to encode peptides of Z50 amino acids were
RNA isolation and RNA-seq. RNA was isolated separately from different preserved; genes predicted were represented by their coding and inferred amino
developmental stages and body parts of T. canis (three biological replicates of L3s; acid sequences.
and four replicates of each female and male reproductive tracts; female and male
anterior bodies; female and male guts) employing TriPure (Roche Molecular
Biochemicals). For L3s, packed volumes of 20–50 ml were used, equating to Functional annotation of all predicted protein sequences. First, following the
thousands of larvae. For other stages and components, packed volumes of prediction of the protein-coding gene set for T. canis, each inferred amino acid
100–200 ml from single worms were used. RNA yields were estimated spectro- sequence was assessed for conserved protein domains using InterProScan 5 Web
photometrically (NanoDrop 1000), and integrity was verified using a Bioanalyzer. Services (http://www.ebi.ac.uk/Tools/pfa/iprscan5/) and InterPro 44.0 (http://
RNA-seq was carried out according to manufacturer’s instructions (Illumina). www.ebi.ac.uk/interpro/). Second, amino acid sequences were subjected to BLASTp
v.2.2.26 (E-valuer10 8) using the following protein databases: C. elegans in
WormBase WS240 (www.wormbase.org); Swiss-Prot and TrEMBL (released in
Pre-processing of sequence reads. Genomic and RNA-seq reads derived from March 2013) within UniProtKB; (KEGG release 58); NCBI protein nr (released in
individual libraries were assessed for quality, adaptors removed and sequencing September 2013). On the basis of the KEGG BLAST matches, the KEGG BRITE
errors corrected using an in-house pipeline. Then, any dog (Canis lupus familiaris; hierarchy was used to classify predicted proteins into families; each coding gene
host of T. canis) and bacterial sequences were identified by comparison with the was assessed against the known KO-term BLAST matches using a custom script;
sequences in the ENSEMBL (release 72) and NCBI (8 April 2013) databases, these matches were collected to known protein families in KEGG BRITE defined in
respectively. In brief, all libraries were aligned to the Canis lupus familiaris genome the KEGG Orthology Based Annotation System database. ES proteins were pre-
using SOAPaligner v.2.21 (http://soap.genomics.org.cn/soapaligner.html) and reads dicted using the programs Phobius (http://www.ebi.ac.uk/Tools/pfa/phobius/),
that mapped to the genome removed. To remove bacterial sequences, one million SignalP (http://www.cbs.dtu.dk/services/SignalP/) and TMHMM (http://www.
randomly selected reads per library were aligned (Z80% of read length and Z90% cbs.dtu.dk/services/TMHMM/). The predicted proteins were listed, as were their
sequence identity in overlapping regions) to known bacterial genomes using annotations based on the conserved InterProScan domain matches with respective
BLAST v.2.2.26; low complexity reads were not aligned. gene ontology terms and then (successively) based on their homology matches
(E-value cut-off: r10 8) to proteins in: the Swiss-Prot, C. elegans, the KEGG,
NCBI nr and TrEMBL databases. Any predicted protein without a match
Transcriptome assembly. Pre-processed RNA-seq reads were assembled de novo (E-value cut-off: r10 8) in at least one of these databases was designated an
using the program Trinity (v.r2012-06-08)51 and also employing a genome-guided hypothetical protein. The final annotation for the protein-coding gene set of
approach employing the programs TopHat2 v.2.0.7 (ref. 52) and Cufflinks2 v.2.1.1 T. canis is available for download at http://bioinfosecond.vet.unimelb.edu.au/
(ref. 53). First, to achieve a consensus (non-redundant) set of transcripts, the open Tcanis. Structural modelling was conducted using the program Phyre (http://
reading frames (ORFs) of all transcripts were predicted. Second, a single transcript www.sbg.bio.ic.ac.uk/phyre2/).
with the longest ORF was retained, if another transcript mapped to the same exon
in the reference genome and shared the same reading frame, and if the length of the
exonic (coding) regions covered by both transcripts divided by length of the RNAi pathway prediction. The RNAi effectors were predicted by the program
shortest ORF was 430%. Both genome-guided and de novo-assembled BLAST v2.2.28 (E-valuer10 15) using internal protein database and a database of
transcriptomes for individual stages and tissues of T. canis were used for the effector sequences42. Predictions were verified against the final functional
prediction of genes. annotations for T. canis and with the inferred proteome of A. suum (WormBase;
WS230). For each gene, levels of transcription in different developmental stages
and tissues (gut and reproductive) were normalized using the trimmed mean of
Genome assembly. Pre-processed genomic PE-libraries (insert sizes 170-, 500- M-values scale-normalization method65.
and 800-bp) and J-libraries (insert sizes: 2-, 5 -and 10-kb) were assembled and
scaffolded using the program SOAPdenovo2 v.2.04 (http://soap.genomics.org.cn/
soapdenovo.html) using a k-mer of 43. The GapCloser v.1.10 program within the Differential transcription analysis and integrated enrichment. The analysis of
SOAPdenovo2 suite was used to close gaps in the scaffolded assembly. The com- RNA-seq data representing different developmental stages, body parts and genders
pleteness of the T. canis draft genome assembly was estimated using the program of T. canis was conducted using the program edgeR (http://www.bioconductor.org/
CEGMA and the number of core eukaryotic genes54. packages/release/bioc/html/edgeR.html). Individual pre-processed PE RNA-seq
data sets (three replicates) were mapped to the predicted gene set using TopHat2
v.2.0.7. The numbers of reads that mapped to each gene were enumerated using the
Prediction of repetitive elements. First, repeat sequences were masked using the program SAMtools66. For each replicate of each data set (sample), the resultant
programs LTR_FINDER v.1.0.5 (ref. 55), PILER v.1.0 (ref. 56), RepeatScout v.1.0.5 read counts were used as input data for edgeR. The levels of differential
(ref. 57), RepeatMasker v.open-4.0.3 (ref. 58) and Tandem Repeat Finder v.4.07 transcription were calculated based on pairwise comparisons among samples of
(ref. 59) using a collection of known repeats in Repbase v.17.02 (ref. 60). Then, all T. canis (that is, L3 versus female anterior; L3 versus male anterior; female anterior
transposons were combined and a non-redundant set of representative transposons versus female reproductive tract; male anterior versus male reproductive tract;
produced by selecting the longest transposon (representing each subclass) that female reproductive tract versus male reproductive tract; female anterior versus
overlapped in sequence by 430% with the shortest transposon. female gut; male anterior versus male gut; female gut versus male gut). Using
dispersion factors calculated in edgeR, the genes were considered differentially
transcribed if the fold change compared with the normalized read count data was
Identification and annotation of non-coding regions and protein-coding genes. Z2, and the false discovery rate was r0.05. In addition, enrichment clusters of
The T. canis protein-coding genes were predicted de novo from the repeat-masked significantly differentially transcribed genes were defined (Fisher’s exact test) and
draft genome and also using homology- and RNA-seq-based approaches against an integrated analysis was conducted. Using a custom script, data linked to these
the unmasked genome. The predictions were then merged into a non-redundant gene clusters were enriched from KEGG pathway and KEGG BRITE hierarchy
gene set using the program GLEAN61. The training gene set for the de novo database. Hubs of these gene clusters were identified using functional gene
predictions was inferred from transcriptomes assembled using the program networking (C. elegans; http://geneorienteer.org/) and weighted gene co-expression
Cuffmerge v.2.1.1 (ref. 53). First, for each locus with multiple transcripts, only the network analysis67,68.
longest ORF was retained. Second, protein sequences inferred from the longest
ORFs in the resultant, non-redundant data set were aligned to the NCBI nr
database (released on 4 September 2013) using the program BLASTp. Matches that Orthology. First, orthologs were predicted by pairwise comparison of individual
aligned over 480% of their length both in reference and query and had 475% proteins (predicted from each gene set) among T. canis, A. suum, B. malayi,
identity represented the training set for the de novo gene-prediction programs C. elegans and P. pacificus using the OrthoMCL v.1.4 (ref. 69) and BLASTp (cut-
AUGUSTUS62, GlimmerHMM63 and SNAP64. Using each of these programs, off: 10 5). Results obtained using these two approaches were displayed in a Venn
a gene set was predicted de novo from the assembled draft genome. In addition, diagram.
Prediction of essential molecules. Essential genes/gene products were regulate multiple neuronal circuits in response to sensory cues. Genetics 156,
predicted70, but using reciprocal best BLAST hits to infer orthologous groups. 123–141 (2000).
Lethal phenotypes (WBPhenotype:0000062 and sub-phenotypes (listed in http:// 26. Crook, M. The dauer hypothesis and the evolution of parasitism: 20 years on
www.berkeleybop.org/ontologies/wbphenotype.obo) in C. elegans were identified in and still going strong. Int. J. Parasitol. 44, 1–8 (2014).
WormBase WS240. Network connectivity data for C. elegans were obtained using 27. Wang, Z. et al. Gene expression analysis distinguishes tissue-specific and
WormNet v.2 (http://www.functionalnet.org/wormnet/cgi-perl/WormNet.v2_ gender-related functions among adult Ascaris suum tissues. Mol. Genet.
nga_form.cgi). Genomics 288, 243–260 (2013).
28. Cantacessi, C. et al. A portrait of the "SCP/TAPS" proteins of eukaryotes–
Additional analyses, and use of software for document preparation. We developing a framework for fundamental research and biotechnological
conducted analyses in a Unix environment or Microsoft Excel 2007 using standard outcomes. Biotechnol. Adv. 27, 376–388 (2009).
commands. Bioinformatic scripts required to facilitate data analysis were designed 29. Newton, S. E. & Munn, E. A. The development of vaccines against
using mainly the Python2.6 scripting language and are available via http:// gastrointestinal nematode parasites, particularly Haemonchus contortus.
research.vet.unimelb.edu.au/gasserlab/. Parasitol. Today 15, 116–122 (1999).
30. Callister, D. M., Winter, A. D., Page, A. P. & Maizels, R. M. Four abundant
novel transcript genes from Toxocara canis with unrelated coding sequences
References share untranslated region tracts implicated in the control of gene expression.
1. Hotez, P. J., Fenwick, A., Savioli, L. & Molyneux, D. H. Rescuing the bottom
Mol. Biochem. Parasitol. 162, 60–70 (2008).
billion through control of neglected tropical diseases. Lancet 373, 1570–1575
31. Othman, A. A. Therapeutic battle against larval toxocariasis: are we still far
(2009).
behind? Acta Trop. 124, 171–178 (2012).
2. Crompton, D. W. Ascaris and ascariasis. Adv. Parasitol. 48, 285–375 (2001).
32. Shanmugam, D. et al. Integrating and Mining Helminth Genomes to Discover
3. Hotez, P. J. & Wilkins, P. P. Toxocariasis: America’s most common neglected
and Prioritize Novel Therapeutic Targets. in Parasitic Helminths: Targets,
infection of poverty and a helminthiasis of global importance? PLoS Negl. Trop.
Scceens, Drugs and Vaccines (ed Caffrey, C. R.) 43–45 (Wiley-VCH Verlag Co.
Dis. 3, e400 (2009).
KGaA, 2012).
4. Hotez, P. J., Dumonteil, E., Hefferman, M. J. & Bottazzi, M. E. Innovation for
33. Zhong, W. & Sternberg, P. W. Genome-wide prediction of C. elegans genetic
the ‘bottom 100 million’: eliminating neglected tropical diseases in the
interactions. Science 311, 1481–1484 (2006).
Americas. Adv. Exp. Med. Biol. 764, 1–12 (2013).
34. Lee, I. et al. A single gene network accurately predicts phenotypic effects of gene
5. Gasser, R. B. A perfect time to harness advanced molecular technologies to
perturbation in Caenorhabditis elegans. Nat. Genet. 40, 181–188 (2008).
explore the fundamental biology of Toxocara species. Vet. Parasitol. 193,
35. Campbell, B. E. et al. Atypical (RIO) protein kinases from Haemonchus
353–364 (2013).
contortus - Promise as new targets for nematocidal drugs. Biotechnol. Adv. 29,
6. Macpherson, C. N. The epidemiology and public health importance of
338–350 (2011).
toxocariasis: a zoonosis of global importance. Int. J. Parasitol. 43, 999–1008
36. Lipinski, C. A. Lead- and drug-like compounds: the rule-of-five revolution.
(2013).
7. Rubinsky-Elefant, G., Hirata, C. E., Yamamoto, J. H. & Ferreira, M. U. Human Drug Discov. Today: Technol. 1, 337–341 (2004).
toxocariasis: diagnosis, worldwide seroprevalences and clinical expression of 37. Gaulton, A. et al. ChEMBL: a large-scale bioactivity database for drug
the systemic and ocular forms. Ann. Trop. Med. Parasitol. 104, 3–23 (2010). discovery. Nucleic Acids Res. 40, D1100–D1107 (2012).
8. Pinelli, E., Brandes, S., Dormans, J., Gremmer, E. & van Loveren, H. Infection 38. Campbell, W. C., Fisher, M. H., Stapley, E. O., Albers-Schonberg, G. &
with the roundworm Toxocara canis leads to exacerbation of experimental Jacob, T. A. Ivermectin: a potent new antiparasitic agent. Science 221, 823–828
allergic airway inflammation. Clin. Exp. Allergy 38, 649–658 (2008). (1983).
9. Overgaauw, P. A. & van Knapen, F. Veterinary and public health aspects of 39. Keiser, J. & Utzinger, J. The drugs we have and the drugs we need against major
Toxocara spp. Vet. Parasitol. 193, 398–403 (2013). helminth infections. Adv. Parasitol. 73, 197–230 (2010).
10. Akao, N. Critical Assessment of Existing and Novel Model Systems of 40. Campbell, B. E. et al. Norcantharidin analogues with nematocidal activity in
Toxocariasis. in Toxocara The Enigmatic Parasite (eds Holland, C. V. & Smith, Haemonchus contortus. Bioorg. Med. Chem. Lett. 21, 3277–3281 (2011).
H. V.) 74–85 (CABI Publishing, 2006). 41. de Savigny, D. H. In vitro maintenance of Tococara canis larvae and a simple
11. Schnieder, T., Laabs, E. M. & Welz, C. Larval development of Toxocara canis in method for the production of Toxocara ES antigen for use in serodiagnosis test
dogs. Vet. Parasitol. 175, 193–206 (2011). for visceral larva migrans. J. Parasitol. 61, 781–782 (1975).
12. Jex, A. R. et al. Ascaris suum draft genome. Nature 479, 529–533 (2011). 42. Dalzell, J. J. et al. RNAi effector diversity in nematodes. PLoS Negl. Trop. Dis. 5,
13. Wang, J. et al. Silencing of germline-expressed genes by DNA elimination in e1176 (2011).
somatic cells. Dev. Cell 23, 1072–1080 (2012). 43. McEwan, D. L., Weisman, A. S. & Hunter, C. P. Uptake of extracellular
14. Ghedin, E. et al. Draft genome of the filarial nematode parasite Brugia malayi. double-stranded RNA by SID-2. Mol. Cell 47, 746–754 (2012).
Science 317, 1756–1760 (2007). 44. Xu, M. J. et al. RNAi-mediated silencing of a novel Ascaris suum gene
15. Dieterich, C. et al. The Pristionchus pacificus genome provides a unique expression in infective larvae. Parasitol. Res. 107, 1499–1503 (2010).
perspective on nematode lifestyle and parasitism. Nat. Genet. 40, 1193–1198 45. Chen, N. et al. Ascaris suum: RNAi mediated silencing of enolase gene
(2008). expression in infective larvae. Exp. Parasitol. 127, 142–146 (2011).
16. C. elegans Sequencing Consortium. Genome sequence of the nematode 46. Robinson, L. N. et al. Harnessing glycomics technologies: integrating
C. elegans: a platform for investigating biology. Science 282, 2012–2018 (1998). structure with function for glycan characterisation. Electrophoresis 33, 797–814
17. McKerrow, J. H., Caffrey, C., Kelly, B., Loke, P. & Sajid, M. Proteases in (2012).
parasitic diseases. Annu. Rev. Pathol 1, 497–536 (2006). 47. Sprent, J. F. Observations on the development of Toxocara canis (Werner,
18. Hewitson, J. P., Grainger, J. R. & Maizels, R. M. Helminth immunoregulation: 1782) in the dog. Parasitology 48, 184–209 (1958).
the role of parasite secreted proteins in modulating host immunity. Mol. 48. Carlsgart, J., Roepstorff, A. & Nejsum, P. Multiplex PCR on single
Biochem. Parasitol. 167, 1–11 (2009). unembryonated Ascaris (roundworm) eggs. Parasitol. Res. 104, 939–943 (2009).
19. Kaminsky, R., Mosimann, D., Sager, H., Stein, P. & Hosking, B. Determination 49. Baermann, G. Eine einfache Methode zur Auffindung von Ankylostomum
of the effective dose rate for monepantel (AAD 1566) against adult gastro- (Nematoden) Larven in Erdproben. Geneesk. Tijdschr. Ned.-Indië 57, 131–137
intestinal nematodes in sheep. Int. J. Parasitol. 39, 443–446 (2009). (1917).
20. Wolstenholme, A. J., Fairweather, I., Prichard, R., von Samson-Himmelstjerna, 50. Jacobs, D. E., Zhu, X.-Q, Gasser, R. B. & Chilton, N. B. PCR-based methods for
G. & Sangster, N. C. Drug resistance in veterinary helminths. Trends Parasitol. identification of potentially zoonotic ascaridoid parasites of the dog, fox and
20, 469–476 (2004). cat. Acta Trop. 68, 191–200 (1997).
21. Lespine, A., Menez, C., Bourguinat, C. & Prichard, R. K. P-glycoproteins and 51. Grabherr, M. G. et al. Full-length transcriptome assembly from RNA-Seq data
other multidrug resistance transporters in the pharmacology of anthelmintics: without a reference genome. Nat. Biotechnol. 29, 644–652 (2011).
prospects for reversing transport-dependent anthelmintic resistance. Int. J. 52. Kim, D. et al. TopHat2: accurate alignment of transcriptomes in the presence of
Parasitol. Drugs Drug Resist. 2, 58–75 (2012). insertions, deletions and gene fusions. Genome Biol. 14, R36 (2013).
22. Krücken, J. et al. Anthelmintic cyclcooctadepsipeptides: complex in structure 53. Trapnell, C. et al. Transcript assembly and quantification by RNA-Seq reveals
and mode of action. Trends Parasitol. 28, 385–394 (2012). unannotated transcripts and isoform switching during cell differentiation. Nat.
23. Maizels, R. M. Toxocara canis: molecular basis of immune recognition and Biotechnol. 28, 511–515 (2010).
evasion. Vet. Parasitol. 193, 365–374 (2013). 54. Parra, G., Bradnam, K., Ning, Z., Keane, T. & Korf, I. Assessing the gene space
24. Maizels, R. M., Schabussova, I., Callister, D. M. & Nicoll, G. Molecular Biology in draft genomes. Nucleic Acids Res. 37, 289–297 (2009).
and Immunology of Toxocara canis. in Toxocara The Enigmatic Parasite. (eds 55. Xu, Z. & Wang, H. LTR_FINDER: an efficient tool for the prediction of
Holland, C. V. & Smith, H. V.) 3–17 (CABI Publishing, 2006). full-length LTR retrotransposons. Nucleic Acids Res. 35, W265–W268 (2007).
25. Daniels, S. A., Ailion, M., Thomas, J. H. & Sengupta, P. egl-4 acts through a 56. Edgar, R. C. & Myers, E. W. PILER: identification and classification of genomic
transforming growth factor-beta/SMAD pathway in Caenorhabditis elegans to repeats. Bioinformatics 21(Suppl 1): i152–i158 (2005).
57. Price, A. L., Jones, N. C. & Pevzner, P. A. De novo identification of repeat Melbourne Water Corporation is gratefully acknowledged. N.D.Y. is an NHMRC Early
families in large genomes. Bioinformatics. 21(Suppl 1): i351–i358 (2005). Career Research (ECRF) Fellow. We acknowledge the contributions of all the staff
58. Smit, A. F. A., Hubley, R. & Green, P. RepeatMasker Open-3.0. (1996-2010). members at WormBase (www.wormbase.org). We thank the Howard Hughes Medical
59. Benson, G. Tandem repeats finder: a program to analyze DNA sequences. Institute (HHMI) and the National Institutes of Health (NIH; P.W.S.). We thank
Nucleic Acids Res. 27, 573–580 (1999). BGI-Shenzhen for providing a commercial service for the construction of RNA-seq and
60. Jurka, J. et al. Repbase Update, a database of eukaryotic repetitive elements. genomic libraries, and for sequencing.
Cytogenet. Genome. Res. 110, 462–467 (2005).
61. Elsik, C. G. et al. Creating a honey bee consensus gene set. Genome. Biol. 8, R13 Author contributions
(2007). P.N. collected and prepared T. canis materials. N.D.Y. and R.B.G. isolated genomic DNA
62. Stanke, M., Tzvetkova, A. & Morgenstern, B. AUGUSTUS at EGASP: using and RNA. P.K.K., N.D.Y. and R.B.G. designed and planned the study, which represents a
EST, protein and genomic alignments for improved gene prediction in the major part of P.K.K.’s PhD project. X.W. was responsible for sequencing. H.C. and
human genome. Genome. Biol. 7(Suppl 1): S1.1–8 (2006). P.K.K. undertook the assembly, annotation and analyses, with support from N.D.Y., Q.L.,
63. Majoros, W. H., Pertea, M. & Salzberg, S. L. TigrScan and GlimmerHMM: two J.M., Y.Y., R.S.H. and A.R.J.; P.K.K. and R.B.G. wrote the manuscript, with inputs from
open source ab initio eukaryotic gene-finders. Bioinformatics 20, 2878–2879 other co-authors (N.D.Y., P.N., A.R.J., G.v.S.-H., A.H., P.T., P.W.S., X.F. and X.-Q.Z.).
(2004). R.B.G., P.K.K. and N.D.Y. supervised the project. X.-Q.Z. and R.B.G. raised the funds.
64. Korf, I. Gene finding in novel genomes. BMC Bioinformatics 5, 59 (2004).
65. Robinson, M. D. & Oshlack, A. A scaling normalization method for differential
expression analysis of RNA-seq data. Genome. Biol. 11, R25 (2010). Additional information
66. Li, H. et al. The Sequence Alignment/Map format and SAMtools. Accession codes: Sequence data are accessible via National Center for Biotechnology
Bioinformatics 25, 2078–2079 (2009). Information (NCBI; www.ncbi.nlm.nih.gov) under accession code JPKZ00000000 (ver-
67. Zhang, B. & Horvath, S. A general framework for weighted gene co-expression sion JPKZ00000000.1; GI:734563511; BioProject: PRJNA248777).
network analysis. Stat. Appl. Genet. Mol. Biol. 4, Article 17 (2005).
Supplementary Information accompanies this paper at http://www.nature.com/
68. Langfelder, P. & Horvath, S. WGCNA: an R package for weighted correlation
naturecommunications
network analysis. BMC Bioinformatics 9, 559 (2008).
69. Li, L., Stoeckert, Jr. C. J. & Roos, D. S. OrthoMCL: identification of ortholog Competing financial interests: The authors declare no competing financial interests.
groups for eukaryotic genomes. Genome. Res. 13, 2178–2189 (2003).
70. Doyle, M. A., Gasser, R. B., Woodcroft, B. J., Hall, R. S. & Ralph, S. A. Drug Reprints and permission information is available online at http://npg.nature.com/
target prediction and prioritization: using orthology to predict essentiality in reprintsandpermissions/
parasite genomes. BMC Genomics 11, 222 (2010).
How to cite this article: Zhu, X.-Q. et al. Genetic blueprint of the zoonotic pathogen
Toxocara canis. Nat. Commun. 6:6145 doi: 10.1038/ncomms7145 (2015).
Acknowledgements
This project was funded by the International Science & Technology Cooperation This work is licensed under a Creative Commons Attribution-
Program of China (grant no. 2013DFA31840; X.-Q.Z. and R.B.G.), the Australian NonCommercial-ShareAlike 4.0 International License. The images or
Research Council (ARC) and the National Health and Medical Research Council other third party material in this article are included in the article’s Creative Commons
(NHMRC) of Australia (R.B.G.); it was also supported by a Victorian Life Sciences license, unless indicated otherwise in the credit line; if the material is not included under
Computation Initiative (VLSCI) grant (VR0007; R.B.G.) on its Peak Computing Facility the Creative Commons license, users will need to obtain permission from the license
at the University of Melbourne, an initiative of the Victorian Government. Other support holder to reproduce the material. To view a copy of this license, visit http://
from the Australian Academy of Science, Alexander von Humboldt Foundation and creativecommons.org/licenses/by-nc-sa/4.0/