Homopolymer tract organization in the human malarial parasite Plasmodium falciparum and related Apicomplexan parasites

Genome-wide surveys of homopolymer tract frequencies in eukaryotic genomes reveal a striking enrichment of long poly dA.dT tracts in intergenic regions (IGR) compared to expected frequencies of random tracts of equivalent base composition [1–3]. Utilizing techniques such as MAINE-FAIRE (micrococcal nuclease-assisted isolation of nucleosomal elements and formaldehyde assisted isolation of regulatory elements) or chromatin immunoprecipitation, both coupled to high throughput sequencing approaches, reveal a relative depletion of nucleosomes over poly dA.dT tracts of increasing length [2, 4–6]. The shorter helical repeat distance within poly dA.dT tracts, along with a narrower and deeper minor grove with a defined spine of hydration, appear to energetically disfavour the necessary remodeling of these tracts for nucleosome binding (for review see [7]). Over recent years, high-resolution maps of nucleosome positioning in a number of eukaryotes have revealed that these poly dA.dT tracts represent a canonical feature of an intrinsic nucleosome positioning code within DNA sequences [4, 5, 8–11]. A second canonical feature, a 10 bp periodicity of AA/TT/TA dinucleotide repeats, provides an opposing role by promoting sharp bending of DNA necessary for wrapping around nucleosomes. Thus, the absence of dinucleotide repeats, together with the presence tracts of poly dA.dT, provides for an ordered, and modifiable, nucleosome landscape where nucleosome depleted regions (NDR) typically act as barriers between well-ordered arrays of nucleosomes organized over exonic sequences [5, 7, 9, 10]. Critically, these flanking 5’ and 3’ NDR, located in the immediate IGR, provide key sites for the regulation of gene expression. 5’ NFR facilitating access of specific transcription factors and the basal transcriptional apparatus to the genomic template, with the more recently described 3’ NDR offering a region for RNA polymerase II (RNAPolII) complex disassembly over transcription termination sites as well as a site for initiation of antisense transcription [4, 11].

Analysis of chromosomes 2 and 3 of Plasmodium falciparum, the aetiological agent of the most virulent form of human malaria, reveals higher than expected frequencies of over-long polydA.dT tracts in IGR [1, 3]. Whilst few studies have explored the function of poly dA.dT tracts in the control of gene expression in this parasite, those that have support a role for these tracts in altering absolute levels of gene expression [12, 13]. Genome-wide profiling of nucleosome positioning during the intraerythrocytic (IE) stage of development indicates variation in positioning subject to temporal developmentally-linked control [14–16]. Maximal levels of nucleosome occupancy are apparent in the final stages of IE development with minimum nucleosome occupancy levels coincident with S-phase and the highest levels of transcriptional activity [16]. P. falciparum is unusual in that the majority of its genome is maintained in a euchromatic state, with no evidence of the highest orders of chromatin packing as a result of the apparent absence of histone H1 (for reviews see [17, 18]). Given recent findings that suggest the spatial organization of chromatin within the nucleus affects gene expression [19], it appears likely that dramatic macro-scale rearrangement of chromatin accompanied by micro-scale nucleosome rearrangements over IGR play a major role in directing the cascade of developmentally-linked mRNA steady state levels during intraerythrocytic schizogony. Spatial variations in nucleosome occupancy in P. falciparum are generally apparent over key expression landmarks [16, 20]. Nucleosomes are preferentially distributed over exonic sequences whilst intergenic regions immediately flanking exonic sequences are relatively nucleosome depleted [14–16]. To date, however, no evidence of a spatial correlation between poly dA.dT tracts and the positioning of these gene-flanking NFR has been demonstrated. Moreover, whilst in other organisms the 5’ flanking NDR typically contain the transcription start site, this is unlikely in P. falciparum. Precise mapping of transcriptional start sites is challenging within the extreme AT nucleotide content of P. falciparum intergenic regions (typically exceeding 80-90%), although a recent study performed by us suggests that these are likely located some 600–1350 bp upstream of the start of the coding sequence [21–23], at least 300-400 bp beyond the mapped 5’ NDR border.

The global impact on human and animal health imposed by the protist parasites of the Apicomplexan phylum, of which P. falciparum is perhaps the best-known member, has driven a sustained effort to sequence their genomes over recent years. Thus, complete annotated genomes are available for: (i) the human malaria parasites P. falciparum, P. vivax and P. knowlesi as well as murine models for this disease (P. berghei and P. yoelii), (ii) the closely related hematozoans of the order Piroplasmarida Theileria spp. and Babesia bovis, aetiological agents of the bovine diseases tropical theileriosis, East Coast fever and babesiosis, (iii) three genomes for Cryptosporidium spp., aetiological agents of a life-threatening diarroheal disease in immunocompromised individuals and (iv) two additional coccidian parasites of the sub-order Eimeriorina, Toxoplasma gondii and Neurospora caninum, responsible for abortive diseases in sheep and cattle, respectively (see Figure 1 for cladogram of these organisms) [20, 24–33]. These resources offer an opportunity to perform a detailed comparative analysis of homopolymeric tract organisation within IGR distinct from previous comparative studies that have emphasised evolutionary diversity in the organisms investigated. In a recent report exploring the size and organisation of IGR within the Apicomplexan phylum, we demonstrated a consensus gene-spacing rule that is shared between the moderately compact genomes in this phylum despite the huge variation in the sizes of their genomes (8-63Mbp) [22]. That is, the size of IGR reflects the nature of the core transcriptional activity over the IGR; group A IGRs, containing two promoters (head-to-head flanking genes), are larger than group B IGRs, which contain one promoter and one terminator (head-to-tail flanking genes), and these in turn are larger than type C IGRs that contain two terminators (tail-to-tail flanking genes). For those parasite genomes with moderately compact gene density (c. 2.0-4.8kbp/ORF), we showed that irrespective of the actual mean sizes of the IGR there is a consensus 3:2:1 ratio in the median size of types A:B:C IGR.

This study was designed to address two aims. First, provide a comparative analysis of homopolymer tract frequency, size and spatial organisation in what may be considered functionally comparable regions of flanking IGR in these evolutionary-related pathogenic organisms. The second aim then correlates the spatial organisation of poly dA.dT tracts in P. falciparum with available nucleosome mapping data to explore a role for poly dA.dT tracts in directing nucleosome positioning over NDR associated with key expression landmarks.

Simple sequence repeats, such as homopolymeric tracts, are important molecular tools in the investigation of genetic disease and evolution. As such, significant efforts have been invested in understanding the forces that shape their frequency, length, spatial distribution and function [1, 3, 35–37]. It is apparent that there is a dynamic balance between stochastic seeding and expansion of homopolymer tracts, their mutation, and non-stochastic selective pressures (e.g. intrinsic DNA factors such as triplet code, nucleosome positioning information and initiation of transcription), that together direct their evolution and persistence. Two leading theories for the seeding and expansion of homopolymer tracts invoke either slipped-strand replication errors or the action of transposable elements such as retrotransposons [38–41]. The evidence for slipped-strand replication errors directing tract expansion is initially conceptually compelling. In general, representation of homopolymeric tracts increases logarithmically when tract lengths are equal to, or exceed, approximately 7 bases – this figure apparently representing a minimal thermodynamic threshold for slipped-strand errors during replication [1]. This theory, however, was challenged in a study of homopolymeric tract frequencies in 27 organisms which showed that there is no single threshold tract-length for homopolymer overrepresentation [3]. Instead, this work demonstrated that the threshold tract-length for overrepresentation is dependent upon the base composition (Fraction_i). Thus, for polydA.dT tracts, as Fraction_A,T increases, the threshold tract length for overrepresentation of these tracts similarly increases. This study of homopolymer tract frequency in Apicomplexan protozoa agrees with this observation, although our study shares with the Zhou et al. study [3] a similar limitation in range of average genomic AT content of organisms (between 50-85%). Whether this same correlation between threshold tract-length for polydA.dT overrepresentation and Fraction_A,T exists in more GC rich organisms needs to be established to explore the universality of this observation and its implication for homopolymer tract-length expansion.

The second theory hypothesises that polydA.dT tracts are seeded from polyadenylated transcripts of retrotransposons. This study of Apicomplexan protozoa presents a challenge, although perhaps not insurmountable, to this hypothesis. Whilst transposable elements are ubiquitously found in the genomes of metazoan eukaryotes, and even within lineages of protozoa closely related to Apicomplexans, there is a lack of overwhelming evidence for such transposable elements in the extant genomes of the Apicomplexan organisms investigated here [42–45]. Apicomplexans share the features of relatively small and moderately compact genomes, and are characterised by gene loss and expansion of species-specific gene families – features attributed to adaptation to a parasitic life cycle and virulence within the host. Critically, whilst there is extensive synteny within families of these parasites (e.g. comparing P. falciparum with P. knowlesi), there is minimal synteny between more distantly related members of the Apicomplexan lineage (e.g. comparing Plasmodium with Cryptosporidium) [42]. This evidence of dramatic genomic rearrangements within ancestral lineages, much more so than in other eukaryotic lineages of similar evolutionary age, coupled with evidence of widespread intron loss during ancestral Apicomplexan genome evolution, have led to the suggestion that retrotransposons were present in Apicomplexan ancestral lineages and may therefore represent the original source of polydA.dT tracts in modern genomes. That said, Plasmodium spp. and Theileria spp. share a recent common ancestor, and show diametrically opposed patterns of homopolymeric tract arrangement in their modern genomes. Thus, any model based on the action of retrotransposons seeding homopolymer tracts within ancestral Apicomplexan lineages would also have to account for the apparent absence of overrepresented homopolymeric tracts in modern Theileria spp. genomes.

Previous comparative studies of homopolymer tract frequency and/or positioning have emphasized phylogenetic diversity in the organisms studied, quite distinct to this study limited to a single phylum [1, 3]. However, this study reports features of homopolymer tract frequency, length and distribution common to these previous comparative analyses. Specifically, the correlation between base composition and length of sequence investigated with homopolymer tract representation, thresholds of overrepresentation and degree of overproportionment. Critically, however, analysis of representation and proportion of polydA.dT and polydG.dC tracts within the IGR reveals intriguing, and novel, lineage-specific patterns within Apicomplexans. Whilst members of the Plasmodium and Cryptosporidium families show a greater level of overrepresentation and overproportionment in polyd.A.dT tracts compared to polydG.dC tracts, a pattern generally shared with other eukaryotes, we find a reversal of this pattern in the coccidians and complete absence in the piroplasmida. The apparent deficit of homopolymeric tracts in the piroplasmida cannot be simply accounted for by the relatively short length of the IGR investigated in these organisms (150-550 bp), as short IGR (200-650 bp) are also analysed here in Cryptosporidium spp. [22]. Further, whilst the coccidians have the most GC-rich genomes of the Apicomplexans investigated here, and may therefore be expected to have longer polydG.dC tracts, both humans and mice share a similar base composition and clearly show higher levels of overrepresented and overproportioned polydA.dT tracts in their IGR [3]. Resolution of the impact of these distinct patterns of tract organisation will require nucleosome-positioning data currently not available for these organisms. Of note, however, is that whilst these organisms show distinct patterns of tract organisation in their IGR, the preferential spatial arrangement of polydA.dT tracts immediately flanking ORF is conserved. This organisation presumably confers a comparable role in nucleosome barrier placement at these sites, and suggests that canonical features of an intrinsic nucleosome positioning code are similar in all Apicomplexans.

In many eukaryotes, the 5’ NDR contains the transcription start site. Whilst accurate transcription start site mapping data are not currently available for P. falciparum, a similar placement of the transcription start site in the 5’ NDR is unlikely. P. falciparum transcripts contain significant untranslated regions (c. 800–1800 bases), with recent modelling performed by us indicating that some 75-80% of this is accounted for by the 5’UTR – suggesting that transcription start sites typically lie some 600 to 1350 bp upstream of the translational start site, much further than the 200-300 bp of the 5’ NDR [22]. In the absence of mapped transcriptional start sites, we instead correlated nucleosome positional mapping data and the spatial arrangement of polydA.dT with bioinformatic predictions of core promoters [20, 46]. Some caution needs to be applied to the interpretation of these data, as structural effects of homopolymer tracts on DNA bending are incorporated into the bioinformatic algorithm used to predict core promoters [46]. However, it is noteworthy that the centre of the NDR located over the most confidently predicted core promoters correlates with a peak of polydT on the sense strand, the 3’ end of which lies immediately upstream of the likely transcription start site. This arrangement is very similar to that recently described in the similarly AT rich Dictyostelium discodieum[47]. Here, the authors speculate that the RNA polymerase II preinitiation complex exploits the relative instability of the polyrU.dT RNA:DNA hybrid in abortive transcription, a mechanism key for directing a programme of coordinated temporal gene expression. It is therefore noteworthy that a not dissimilar programme of temporal global upregulation of RNAPolII pre-initiation complex (PIC) processivity is found in P. falciparum[48]. During intraerythrocytic schizogony, RNAPolII PIC is constitutively associated with promoter regions, and presumably undergoes repeated rounds of abortive transcription prior to an onset of increased RNAPolII processivity in the early trophozoite stage, some 18–22 h after the invasion of the host erythrocyte [49]. The molecular basis of this increase in RNAPolII processivity remains to be determined, although conserved features of the T-regulatory loop and/or the heptad repeats within the carboxy-terminal domain of the P. falciparum RNAPolII large subunit would suggest regulation through kinase phosphorylation activity at these sites [50].

Placement of poly dA.dT tracts in the centre of NDR organised over P. falciparum promoters supports their role as intrinsic factors in determining nucleosome positioning in these regions. In budding yeast, poly dA.dT tracts have also been demonstrated to act as extrinsic factors in this same role by modulating interactions of promoter DNA with transcription factors, the RNAPolII complex and chromatin remodelling machinery to provide additional context, and directionality, for distinct promoter classes. The data presented here represent the average nucleosome positioning data, and whilst sufficient to demonstrate the intrinsic role, lacks context to explore any extrinsic role. That such a role may exist in P. falciparum is suggested by observed temporal variation in nucleosome positioning over promoters and the preferential positioning of non-canonical histone containing nucleosomes over IGR [14, 15, 51]. Resolution would require additional data correlating temporal gene expression patterns with binding of components of the transcriptional apparatus and/or chromatin remodelling machinery.

A key finding of this study was the highly overrepresented and overproportioned polydA.dT tracts in the IGR of the Plasmodium spp. compared to the other Apicomplexans investigated. Significantly, analysis of homopolymer tract organisation over the open reading frames in Plasmodium spp. reveals that polydA.dT are similarly highly overrepresented and overproportioned in these regions. This organisation of polydA.dT over open reading frames is quite distinct from the other Apicomplexans investigated here, and eukaryotic genomes in general [1, 3]. The presence of polyA tracts within P. falciparum open reading frames is well documented; leading to a pronounced codon usage bias and inclusion of low complexity amino acid stretches (particularly polylysine and polyasparagine) within the protein products. This same bias in codon usage is also evident throughout Plasmodium spp., similarly trending towards the inclusion of polyA tracts, although not to the extent of the extremely AT-rich P. falciparum. These low complexity regions (LCR) have been implicated in a variety of regulatory mechanisms including control of translational efficiency, a “smokescreen” of immunodominant regions assisting in the evasion of the host adaptive immune response, or as a protective response to heat-shock [52–56]. Reconciling these diverse adaptive and mechanistic roles for low complexity regions has recently been facilitated by their classification into three distinct groups based on their AT-content and heterozygosity [56]. Poly dA.dT would form part of the polyN group, exhibiting both high AT content and high levels of heterozygosity. No clear role for the polyN class of LCR exists, although speculation for a role in translational efficiency or an equally plausible neutral role has been made. For the neutral role, some form of selection pressure would be required to promote and/or retain these polyN LCR. This pressure may exist if we consider a role for poly dA.dT tracts in the formation of origins of replication in Plasmodium spp., as have been implicated in both budding and fission yeast [4, 6], during the critical process of infection and colonisation of the female mosquito following a blood meal on an infected vertebrate host. Plasmodium spp. are haploid for the majority of their life-cycle, with diploid stages only present in early stages of development in the mosquito midgut. Gametogenesis starts in erythrocytes and is completed when the macro- (female) and micro- (male) gametocytes are taken up during a mosquito blood-meal. The final maturation step for the microgametocyte, termed exflagellation, is a remarkable process taking some 20 minutes during which the microgametocyte undergoes three rounds of replication to generate eight flagellated microgametes [57–59]. This rapid replicative cycle would require each of the three rounds of DNA replication to take place over a 3–5 minute window, which, assuming a processivity of the order of 1-2kbp/min for eukaryotic replication [60], would require an exceptionally high density of origins of replication throughout the genome. As some 50% of Plasmodium spp. genomes encode protein, with a gene density of between 2.6-4.6 kb/gene, these origins of replication likely cannot be restricted solely to IGR, and may therefore provide a novel explanation to account for the significant overrepresentation and overproportionment of polydA.dT observed in both ORF and IGR of Plasmodium spp. genomes.

Unlike previous comparative studies that emphasise evolutionary diversity in their determination of features of homopolymeric tract representation, length and spatial organisation, this study is instead restricted to a single phylum of unicellular parasites that all share moderately compact genomes. We describe features of polydA.dT tract organisation within this phylum that support a canonical role as intrinsic regulators of nucleosome positioning as well as findings that support nucleotide fraction as a key determinant in the thermodynamic threshold for tract expansion. Critically, we also present evidence for a novel lineage-specific organisation of homopolymer tract organisation in both intergenic and genic compartments along with evidence that overproportionment of homopolymer tract length may be dependent on the available intergenic space into which expand. Given the general lack of specific transcription factors in Apicomplexan genomes, a dynamic programme of nucleosome binding and rearrangement likely plays a significant role in the temporal and absolute regulation of transcription. Our observations relating to polydA.dT spatial arrangement and enrichment likely reflects the impact of chromatin structure and function in shaping the genomes of these parasites important for human and animal health.