An intronic open reading frame was released from one of group II introns in the mitochondrial genome of the haptophyte Chrysochromulina sp. NIES-1333.

Mitochondrial (mt) genome sequences, which often bear introns, have been sampled from phylogenetically diverse eukaryotes. Thus, we can anticipate novel insights into intron evolution from previously unstudied mt genomes. We here investigated the origins and evolution of three introns in the mt genome of the haptophyte Chrysochromulina sp. NIES-1333, which was sequenced completely in this study. All the three introns were characterized as group II, on the basis of predicted secondary structure, and the conserved sequence motifs at the 5' and 3' termini. Our comparative studies on diverse mt genomes prompt us to propose that the Chrysochromulina mt genome laterally acquired the introns from mt genomes in distantly related eukaryotes. Many group II introns harbor intronic open reading frames for the proteins (intron-encoded proteins or IEPs), which likely facilitate the splicing of their host introns. However, we propose that a "free-standing," IEP-like protein, which is not encoded within any introns in the Chrysochromulina mt genome, is involved in the splicing of the first cox1 intron that lacks any open reading frames.

Introduction

Mitochondrial (mt) genomes can be regarded as a model system for studying intron evolution, as a massive amount of mt genome data including group I and/or group II (gII) introns have been accumulated (for instance, see NCBI Organelle Genome Resource: http://www.ncbi.nlm.nih.gov/genomes/OrganelleResource.cgi?opt=organelle&taxid=2759). Group II introns are found in the genomes of prokaryotes (bacteria and archaebacteria), as well as mitochondria and plastids,, which are derived from an α-proteobacterium and a cyanobacterium, respectively. So far, gII introns have been identified in mt genomes from members of phylogenetically diverse eukaryotic assemblages such as Metazoa,, Jakobida,, Archaeplastida,- Fungi,, Cryptophyta, Haptophyta,, and Stramenopiles.- These gII introns possess features at both the primary and secondary structure levels. At the primary structure level, gII introns possess highly conserved sequence motifs at the 5′ and 3′ ends (i.e., 5′-GTGYG…AY-3′; Y for T or C). At the secondary structure level, we anticipate the transcripts of typical gII introns (intron RNAs) to form a characteristic bulge structure with six stems, so-called domains I to VI. Both primary and secondary structures of intron RNAs are most likely critical for the splicing reaction.

Group II introns can be regarded as mobile genetic elements, which are transmittable between an intron-containing and an intron-lacking loci (intron homing), regardless of their evolutionary distance. The mobility of gII introns are most likely conferred by the proteins encoded within gII introns (intron-encoded proteins or IEPs). Typical IEPs comprise three functionally distinct domains, namely i) reverse transcriptases (RT), ii) domain X, which is also referred to as maturases, and iii) endonucleases (En), although some IEPs were reported to lack En domain., Among the three domains in IEP, RT and En domains are predicted to catalyze reverse transcription of intron RNA and digest the target (intron-lacking) locus, respectively., Domain X may not be responsible for intron mobility, but assists splicing by stabilizing the conformation of intron RNA. Nevertheless, we have known of many “IEP-free” gII introns, and it is difficult to predict the protein factors, which cooperate with a particular IEP-free intron in trans. To our knowledge, there is only a single report that successfully identified the organellar genome-encoded trans factor involved in the splicing of a IEP-free intron.

In this study, we completely sequenced the mt genome of a haptophyte Chrysochromulina sp. NIES-1333, and identified three introns in total, two of those are found in the cox1 gene encoding cytochrome c oxidase subunit 1, and the last one is found in the rnl gene encoding large subunit rRNA. Analyses of the intron sequences suggest that the three introns in the Chrysochromulina mt genome belonged to group II. We identified two open reading frames encoding putative IEPs. Both showed significant sequence similarity to gII intron-hosted IEPs in the mt genome; one is orf627 encoded in the second cox1 intron, and the other is orf584, which is free-standing. Phylogenetic analyses of IEPs and comparisons of intron positions across phylogenetically diverse mt genomes revealed that the Chrysochromulina mt genome shares homologous introns with distantly related mt genomes.

Results and Discussion

Overview of the Chrysochromulina mt genome

The mt genome of Chrysochromulina sp. NIES-1333 was assembled into a circular molecule of 34,291 bp in length with an A + T content of 70.0% (Fig. 1). We identified 16 functionally assignable open reading frames (including those for two IEPs; see below). UGA codons are most likely assigned for tryptophan instead of termination signal (Table S1), as reported previously. This type of mt genetic code was reported in other members of Prymnesiophyceae, one of the two classes of Haptophyta.- We detected 26 tRNA genes and a set of small and large subunits of rRNAs; No 5S rRNA gene was identified. A set of tRNA genes identified in the mt genome is sufficient to translate all amino acid codons except for GGN (N = A, C, G or U) codons for glycine (Table S1). All genes mentioned above were encoded on a single strand. A region with an approximate length of 1.6 Kbp, which contained a single tRNA gene for isoleucine, was found to be duplicated (shown by open arrows in Fig. 1).

Figure 1. Mitochondrial Genome of Chrysochromulina sp. NIES-1333. Protein-coding genes and rRNA genes are represented by boxes. Gray boxes represent two open reading frames, of which amino acid sequences showed significant sequence similarity to intron-encoded proteins. Transfer RNA genes are represented by lines. Introns are shown in dotted lines. Arrows represent duplicated regions.

In terms of gene repertoire, the Chrysochromulina mt genome is fundamentally similar to those of other haptophytes, namely Emiliania huxleyi, Diacronema lutheri (www.bch.umontreal.ca/ogmp/projects/pluth/gen.html‎), Phaeocystis spp., as shown in Table 1.

Table 1. Gene repertoires in haptophyte mitochondrial genomes

	Chrysochromulina sp.	Diacronema lutheri	Emiliania huxleyi	Phaeocystis globosa	Phaeocystis antarica
rnl	Y [1]	Y [1]	Y	Y [1]	Y
rns	Y	Y [1]	Y	Y	Y
rrn5	N	Y	Y	N	N
trn	23 species	22 species	23 species	23 species	23 species
nad1	Y	Y	Y	Y	Y
nad2	N	N	N	Y	Y
nad3	N	N	N	Y	Y
nad4	Y	Y	Y	Y	Y
nad4L	Y	Y	Y	Y	Y
nad5	Y	Y	Y	Y	Y
nad6	Y	Y	Y	Y	Y
cob	Y	Y	Y	Y	Y
cox1	Y [2]	Y [1]	Y	Y	Y
cox2	Y	Y	Y	Y	Y
cox3	Y	Y [1]	Y	Y	Y
atp4	N	Y	Y	Y	Y
atp6	Y	Y [1]	Y	Y	Y
atp8	N	Y	N	Y	Y
atp9	Y	Y [1]	Y	Y	Y
rps3	N	N	Y	Y	Y
rps8	N	N	Y	N	N
rps12	Y	Y	Y	Y	Y
rps14	N	Y	Y	N	Y
rps19	N	Y	N	N	N
rpl14	N	Y	N	N	N
rpl16	Y	Y	Y	Y	Y
dam	N	N	Y	N	N
Others	orf627a  orf538b	orf636c  orf105d	orf104d	N	N

Y, yes; N, no. Numbers of introns are shown in brackets. aEncoded in the second cox1 intron. bFree-standing open reading frame encoding a protein with amino acid sequence similarity to group II intron-encoded proteins. cEncoded in the cox1 intron. dEncodes an uncharacterized protein.

General features of the intron in the Chrysochromulina rnl gene

The rnl gene hosts a single intron encoding no apparent open reading frame (designated as Ch_rnli; Fig. 1). The intron was found to be inserted at the position between the 837th and 838th bases in the Homo sapiens homolog (GeneID: 4550 in NC_012920). Ch_rnli starts with 5′-GTGCG… and ends with …CT-3′, which is similar to the consensus motifs shared among typical gII introns (5′-GTGYG…AY-3′; Y for T or C). Although the intron sequence was too divergent to predict the entire secondary structure, we successfully identified the typical domains, V and VI, which are characteristic secondary structures of gII introns (Fig. S1A), with the aid of MFannot (http://megasun.bch.umontreal.ca/cgi-bin/mfannot/mfannotInterface.pl). Thus, we characterized Ch_rnli as group II.

The homing position of Ch_rnli was found to be identical to those of rnl introns found in a member of Archaeplastida (the red alga, Pyropia haitanensis; NC_017751), two members of Stramenopiles (the brown alga, Pylaiella littoralis, and the diatom, Phaeodactylum tricornutum; NC_003055 and HQ840789, respectively), and a member of Opisthokonta (the fungus, Gigaspora rosea; NC_016985) (Fig. S1A: Note that none of the rnl introns found in the mt genomes of other haptophytes, D. lutheri and Phaeocystis globosa, shared the homing position with Ch_rnli). We predicted the secondary structures of domains V and VI in Ch_rnli and the four introns described above, but detected no apparent homology at the nucleotide sequence level among them (Fig. S1B). Furthermore, Ch_rnli hosts no IEP, which is considered a key aspect to inferring intron evolution. Thus, we hesitate to discuss the evolutionary relationship among Ch_rnli and the introns listed above, solely based on their homing positions.

General features of two introns in the Chrysochromulina cox1 gene

Two introns were found in the Chrysochromulina cox1 gene (Fig. 1). We designate the first and second introns in the cox1 gene as Ch_cox1i1 and Ch_cox1i2, respectively. The two introns commonly start with 5′-GTGCG… and end with …AC-3′, being consistent with the consensus motifs shared among typical gII introns (5′-GTGYG…AY-3′). Ch_cox1i1 was inserted between the second and third letters (phase-2) of the codon corresponding to Phe68 in the H. sapiens cox1 gene, sharing the homing position with the gII introns in cox1 genes of the cryptophyte, Rhodomonas salina, and the diatom, Phaeodactylum tricornatum (Fig. S2C). Ch_cox1i2 was found at phase-2 of the codon corresponding to Phe237 in the H. sapiens cox1 gene, being homologous to those of the gII introns in cox1 genes of the haptophyte D. lutheri and the diatom Ulnaria acus (Fig. S2C). Ch_cox1i2 hosts an intronic open reading frame for an IEP, while Ch_cox1i1 encodes no apparent open reading frame. Both Ch_cox1i1 and Ch_cox1i2 can be folded into the characteristic secondary structures shared among gII introns, albeit with some ambiguity remaining in domain I (indicated as “DI” in Fig. S1B and S1C). All together, we concluded that the two introns belong to group II.

Evolution of Ch_cox1i2 and its IEP

The IEP encoded in the intronic open reading frame of Ch_cox1i2, ORF627, most likely facilitates splicing of the host intron. The ORF627 amino acid sequence showed similarity to other gII intron-hosted IEP sequences deposited in the GenBank database; the top blastp hit was an IEP encoded in the first gII intron of the cox1 gene in the haptophyte D. lutheri (Dl_cox1i), with a 49% sequence similarity and an E-value of 0.0. In both maximum-likelihood (ML) and Bayesian analyses of an IEP alignment (Fig. 2), Chrysochromulina ORF627 formed a clade with two IEPs encoded in cox1 gII introns, namely Dl_cox1i and that of the diatom U. acus (Ua_cox1i) with a ML bootstrap support value (MLBP) of 96% and a Bayesian posterior probability (BPP) of 1.00 (Fig. 2). As we generally believe that gII introns and their IEPs have coevolved, the intimate relationship among the IEPs encoded in Ch_cox1i2, Dl_cox1i, and Ua_cox1i suggests that their host introns are derived from a single ancestral intron bearing an IEP. The single origin of Ch_cox1i2, Dl_cox1i, and Ua_cox1i discussed above is consistent with the fact that the three introns share a homing position (Fig. S2C).

Figure 2. Phylogeny Inferred from 52 Intron-encoded Protein (IEP) Amino Acid Sequences. The IEP alignment was subjected to both maximum-likelihood (ML) and Bayesian methods. As the two methods reconstructed very similar trees, only ML tree is shown here. The tree is rooted by the bacterial sequences. Values at nodes represent ML bootstrap support values greater than 50%. The nodes supported by Bayesian posterior probabilities equal to or greater than 0.95 are highlighted by thick lines. The IEPs encoded in cox1 introns are shaded in orange. The detailed homing positions of cox1 introns are given on the right side of the tree. Codon numbers are based on the Homo sapiens cox1 gene (GenBank accession number, YP_003024028). Free-standing IEPs are highlighted with stars.

The ancestral haptophyte species likely possessed a cox1 gene with a particular gII intron, as Chrysochromulia sp. and D. lutheri are representatives of two major classes, Prymnesiophyceae and Pavlovaphyceae, in Haptophyta, respectively. This scenario suggests that multiple intron losses occurred in the cox1 genes of Emiliania huxleyi,, members of the genus Phaeocystis,, and Isochrysis galbana.

The IEP phylogeny and comparison of homing positions imply that the homologous introns are present in two distantly related branches (haptophytes and diatoms) in the tree of eukaryotes. This sporadic intron distribution can be explained by a scenario incorporating lateral intron transfer. There is an alternative, but less plausible scenario assuming that the cox1 gene in an ancestral organism, which has existed prior to the divergence of major eukaryotic assemblages including diatoms and haptophytes, may have already possessed a gII intron at phase-2 of the codon corresponding to Phe237 in the H. sapiens cox1 gene, and would have been (secondarily) lost in multiple descendants (i.e., ancestral co-occurrence followed by multiple secondary losses). We prefer the scenario incorporating lateral intron transfer to the alternative one, but these scenarios should be reexamined by future studies based on a broader diversity of gII introns (and their IEPs) compared with those considered in the current study.

Link between a free-standing orf584 and an IEP-free Ch_cox1i1

Most IEPs are encoded in intronic open reading frames (as observed in Ch_cox1i2; see above), but a few of those are free-standing in genomes (e.g., orf732 in the Marchantia polymorpha mt genome; highlighted by a star in Fig. 2). Our blast search showed that Chrysochromulina ORF584, which is encoded in a free-standing open reading frame, bore a significant sequence similarity to gII intron-hosted IEPs; The top blastp hit of ORF584 amino acid sequence was an IEP (ORF724) encoded in the first intron of the cox1 gene in the diatom, P. tricornutum (Pt_cox1i1) with a 53% sequence similarity and an E-value of 0.0. ORF584 equips En, RT, and domain X, implying that this protein assists intron splicing. The phylogenetic analyses of the IEP alignment (Fig. 2) recovered a robust affinity between Chrysochromulina ORF584 and ORF724 encoded in Pt_cox1i1 with a MLBP of 100% and a BPP of 1.00 (Fig. 2). This indicates the two proteins were derived from the single ancestral IEP encoded in a gII intron, which is homologous to Pt_cox1i1, the host intron of ORF724. Curiously, Pt_cox1i1 and Ch_cox1i1 appeared to share a homing position (see Fig. S2C). We also noticed that the nucleotide sequence of domain VI in Ch_cox1i1 and that in Pt_cox1i1 are similar to one another (Fig. S1D), although this domain sequences are generally variable among gII introns. The homing position and sequence similarity in domain VI, between Ch_cox1i1 and Pt_cox1i1, suggests that the two introns are homologous to each other. All together, we here propose that ORF584 used to be encoded in Ch_cox1i1, and still assists the splicing of the host intron even after being free-standing secondarily in the current Chrysochromulina mt genome. To the best of our knowledge, co-relation between a particular pair of free-standing IEP and IEP-free introns has been reported only one time prior to this work.

The first intron in the R. salina cox1 gene (Rs_cox1i1) is unlikely to be homologous to Pt_cox1i1 or Ch_cox1i1, although the three introns share the homing position (Fig. S1C). The IEP phylogeny placed the IEP encoded in Rs_cox1i1 in a remote position from the clade of ORF584 and ORF724 (Fig. 2), strongly arguing against the homology between Rs_cox1i1 and Pt_cox1i1/Ch_cox1i1. The homology between Pt_cox1i1 and Ch_cox1i1, which were found in two phylogenetically distantly related species (i.e., a haptophyte and a diatom), can be explained by lateral intron transfer. Nonetheless, we cannot exclude the alternative possibility assuming ancestral co-occurrence followed by multiple secondary losses. We prefer the simplicity of the first scenario incorporating lateral intron transfer, but the alternative scenario should not be ignored before mt genome diversity is sufficiently covered.

Materials and Methods

Cell culture and mt genome sequencing

The haptophyte Chrysochromulina sp. NIES-1333 was purchased from the National Institute for Environmental Study (NIES; 16–2 Onogawa, Tsukuba, Ibaraki 305–8506, Japan). Haptophyte cells were grown in f/2 medium (http://mcc.nies.go.jp/02medium.html#f2) at 20 °C under 14 h light/10 h dark cycles. The cultured cells were harvested by centrifugation. Total DNA and total RNA were extracted from the harvested cells by CTAB buffer and TRIzol (Invitrogen), respectively. Total RNA was used to synthesize cDNA with random hexamers and Superscript II reverse transcriptase (Invitrogen). RNA extraction and cDNA synthesis were conducted following manufacturers’ protocols.

We amplified the entire mt genome by combination of LA PCR with TaKaRa LA Taq DNA polymerase (TaKaRa), genome walking with the GenomeWalker Universal kit (Clonthech), and rolling circle amplification (RCA) with the illustra TempliPhi 100 Amplification kit (GE Healthcare Life Sciences). Amplified DNA fragments of < 3 Kbp-long and those of < 10 Kbp-long were cloned to pGEM T-easy vector (Promega) and pCR-TOPO-XL (Invitrogen), respectively. The short amplicons (< 3 Kbp-long) were sequenced by the Sanger method. 454 pyro-sequencing by the GS-Jr system (454 Sequencing, Roche) was performed on the long amplicons (≥ 10 Kbp). Newbler (454 Sequencing, Roche) was used for de novo assembly of the pyro-sequencing reads. The DNA amplifications and sequencing described above were conducted by following manufactures’ instructions, except the RCA with custom primers instead of random hexamers supplied in the kit, as described in Kamikawa et al. (2014). The custom primers used for the RCA were designed based on the cob and cox3 sequences determined previously. All the DNA sequences obtained were finally assembled into a circular molecule, with an approximate length of 34 Kbp. The complete mt genome sequence was deposited to DDBJ/EMBL/GenBank accession number AB930144.

Genome analyses

Genes encoding proteins and rRNAs were identified by blastx and blastn searches, respectively, against the non-redundant database in National Center for Biotechnology Information (http://blast.ncbi.nlm.nih.gov/Blast.cgi). Transfer RNA-encoding genes were found by using tRNAscan-SE. Independent from the analyses described above, we re-annotated the genome by MFannot (http://megasun.bch.umontreal.ca/cgi-bin/mfannot/mfannotInterface.pl).

During the annotation of the Chrysochromulina mt genome, we noticed that both cox1 and rnl genes are intervened by introns. The precise boundaries of the introns in the mt genome were confirmed by sequencing the corresponding transcripts (cDNAs). The secondary structures of the introns identified in the mt genome were predicted by MFOLD, followed by manual refinement by referring to Toor et al. (2001) and GOBASE database.

Phylogenetic analyses

We found that orf627 and orf584 in the Chrysochromulina mt genome encode IEPs (see above). The conceptual amino acid sequences of the two IEPs were aligned with 46 IEPs encoded in other mt genomes and four bacterial homologs by Muscle. The IEP sequences were retrieved from the GenBank database by referring to pioneering phylogenetic studies (e.g., refs. 17 and 18). After manual refinements and exclusion of ambiguously aligned positions, the final alignment, including 52 IEP sequences with 453 amino acid positions, was used for phylogenetic analyses.

The alignment was subjected to both ML and Bayesian methods using RAxML7.2.6 and PhyloBayes3.3, respectively. We applied LG amino acid substitution model with among-site rate variation approximated by a discrete gamma distribution with four rate categories to ML analyses (LG + Γ + F) of the original alignment and 100 bootstrap replicates. The ML tree was selected by heuristic searches from 10 randomized maximum-parsimony (MP) starting trees. In ML bootstrap analysis, heuristic tree search was performed from a single MP starting tree per replicate. Bayesian analysis was conducted with the LG + Γ+ F model. Two independent Monte Carlo chains were run for 5,800–5,850 cycles, reaching maxdiff value of 0.08353. The first 100 cycles were discarded as “burn-in”; the consensus tree, branch lengths, and BPPs was calculated from remaining trees.