Insights into the transcriptional and translational mechanisms of linear organellar chromosomes in the box jellyfish Alatina alata (Cnidaria: Medusozoa: Cubozoa)

ABSTRACT

Background: In most animals, the mitochondrial genome is characterized by its small size, organization into a single circular molecule, and a relative conservation of the number of encoded genes. In box jellyfish (Cubozoa, Cnidaria), the mitochondrial genome is organized into 8 linear mito-chromosomes harboring between one and 4 genes each, including 2 extra protein-coding genes: mt-polB and orf314. Such an organization challenges the traditional view of mitochondrial DNA (mtDNA) expression in animals. In this study, we investigate the pattern of mitochondrial gene expression in the box jellyfish Alatina alata, as well as several key nuclear-encoded molecular pathways involved in the processing of mitochondrial gene transcription. Results: Read coverage of DNA-seq data is relatively uniform for all 8 mito-chromosomes, suggesting that each mito-chromosome is present in equimolar proportion in the mitochondrion. Comparison of DNA and RNA-seq based assemblies indicates that mito-chromosomes are transcribed into individual transcripts in which the beginning and ending are highly conserved. Expression levels for mt-polB and orf314 are similar to those of other mitochondrial-encoded genes, which provides further evidence for them having functional roles in the mitochondrion. Survey of the transcriptome suggests recognition of the mitochondrial tRNA-Met by the cytoplasmic aminoacyl-tRNA synthetase counterpart and C-to-U editing of the cytoplasmic tRNA-Trp after import into the mitochondrion. Moreover, several mitochondrial ribosomal proteins appear to be lost. Conclusions: This study represents the first survey of mitochondrial gene expression of the linear multi-chromosomal mtDNA in box jellyfish (Cubozoa). Future exploration of small RNAs and the proteome of the mitochondrion will test the hypotheses presented herein.

KEYWORDS:

• Alatina alata
• linear mitochondrial DNA
• mitochondrial gene expression
• mitochondrial polymerase
• tRNA C-to-U editing

Background

The mitochondrion, largely responsible for cellular respiration, is thought to have originated through a single endosymbiotic event of an α–proteobacterium¹ very early in the evolution of eukaryotes (though see Pittis and Gabaldòn (2016)² for an alternative hypothesis). Subsequent gene transfer to the nucleus of the host cell and loss of genes from the proto-mitochondrion³ has rendered the typical mitochondrial genome a small molecule that encodes only a fraction of the genes necessary for the functioning and maintenance of the mitochondrion. In most animals (Metazoa), mitochondrial DNA (mtDNA) comprises a single, circular and compact molecule of ∼17 kbp in size, encoding 13 (+2/−2) protein-coding genes involved in oxidative phosphorylation, 2 rRNA (rRNA) genes, and up to 25 tRNA (tRNA) genes involved in the translation of mitochondrial-encoded proteins. While gene order in the mt genome is highly variable, most protein-coding and ribosomal genes are interspersed with tRNA genes. The replication of the typical circular mtDNA is performed through the formation of a D-loop, where the replication of the two strands is out of sync.⁴ The “tRNA punctuation model” proposes that mtDNA is transcribed into one or two polycistronic precursor transcripts (pre-mRNA molecules) that are then excised at the tRNA sites,⁵ releasing mostly single-gene (monocistronic) mRNAs, rRNA and tRNA transcripts.⁶

In the early branching metazoan phylum Cnidaria (sea anemones, corals, hydras, hydroids, and jellyfish), the mt genome encodes for only one or 2 tRNA genes,^7-17 challenging the tRNA punctuation model of gene expression. Species from the cnidarian class Anthozoa (the clade containing reef-building corals, sea anemones, soft corals, and their relatives) possess relatively large, circular mitochondrial genomes in which protein-coding genes are interspaced with intergenic regions (IGRs) that might play a role in mitochondrial mRNA maturation (e.g.^18,19). The mt genome in the cnidarian subphylum Medusozoa (e.g. hydroids and jellyfish) is organized into mono-, bi-, or multi-linear chromosomes.^7-13,20 Medusozoan mtDNA is very compact with a few small IGRs.^12,15 In hydrozoans, Kayal et al. (2015)¹⁵ suggested that both rRNA genes and stemloop-like secondary structures formed by IGRs could be recognized by the nucleolytic processing mechanism, producing mature often multi-cistronic mRNA transcripts.

Cubozoans (box jellyfish) have mtDNA made up of eight small linear chromosomes.^12,13,20 Mitochondrial chromosomes range from 4 to 6 kbp in length and harbor between one and four genes, flanked by long and conserved telomere-like inverted repeats (IRs) with opposite orientation thought to be involved in the maintenance of the linear mito-chromosomes.¹³ In addition to the 13 protein-coding genes common in Metazoa, most medusozoan mtDNA (excluding hydroidolinan hydrozoans) encodes for two extra protein-coding genes (mt-polB and orf314). Mohammadou et al (2004)²¹ have found that mt-polB has affinities with polymerases from the B family. Kayal et al. (2012)¹² have suggested that orf314 might be acting as a terminal protein by binding the 5′ single-stranded end of the linear chromosomes as described for the liverwort Marchantia polymorpha²² and some yeast species.^23,24 mt-polB and orf314 are present in all medusozoan subclasses and lack internal stop codons; in cubozoans mt-polB and orf314 are found together on an independent linear mito-chromosome. The maintenance of these genes through long evolutionary times (acquired before the divergence of Medusozoa and maintained in several medusozoan clades) has been used as evidence to suggest that they have functional roles in the mitochondrion.^12,13 Additional evidence on the expression of these genes is needed to confirm and refine this hypothesis and to assess specific functions.

Here, we surveyed recently released RNA-seq data from the box jellyfish Alatina alata (SRR1952741) and assembled the near complete multi-chromosomal mitochondrial genome, as well as several mitochondrial-targeted genes involved in mtDNA expression, including but not limited to mitochondrial aminoacyl-tRNA synthases (mt-aaRS), polymerases, and ribosomal proteins. Our analysis offers insights into the expression of mtDNA genes and provides additional evidence about the functional roles of mt-polB and orf314.

Material and methods

We obtained a recently published transcriptome assembly for Alatina alata²⁵ as well as 2 publicly available short read archive (SRA) datasets from the National Center for Biotechnology Information (NCBI; ncbi.nlm.nih.gov/sra/): SRR217435 containing SOLID DNA-seq reads and SRR1952741 containing ILLUMINA HiSeq 2500 RNA-seq reads. The RNA-seq library for the data set used in this study was prepared using the TruSeq mRNA enrichment protocol, where oligo-dT beads capture RNA molecules with poly(A) tail, before random reverse transcription is performed.²⁵ Tissue on which the transcriptome assembly was based did not include any gut tissue, lessening concerns about contamination from other organisms. We mapped the SOLID DNA-seq data to a published mt genome of Alatina alata (JN642330-JN642344) using Bowtie2²⁶ with –end-to-end and –very-sensitive options to estimate individual chromosome and gene coverage.

Fig. S1 details the bioinformatics workflow we used in this study. In short, we captured and assembled mitochondrial contigs from the assembled transcriptome using the BLAST plugin in Geneious v.9.²⁷ We filtered raw reads from the ILLUMINA HiSeq RNA-seq data with Trimmomatic v.0.32²⁸ and mapped filtered reads to the assembled mitochondrial chromosomes using Bowtie2, with default settings. We used Geneious to visualize and manually check the assemblies, and to build consensus sequences. We annotated individual chromosomes following the reference mitogenome and predicted the origin and end of transcription. Finally, we mapped trimmed RNA-seq reads back to the mitochondrial consensus sequences to estimate individual chromosome and gene coverage with Bowtie2 using –end-to-end and –very-sensitive options (Fig. S1).

We also investigated several genes involved in the expression and function of the mitochondrion from the RNA-seq data: the glutamyl-tRNA amidotransferase subunits (GatA-C) that are absent in other cnidarian classes,²⁹ the mitochondrial ribonucleases complex P (MRPP1-3) and Z (RNase Z) responsible for the cleavage of the 5′ and 3′ ends of the tRNA, the mitochondrial poly(A) RNA polymerase (MTPAP) that creates the 3′ poly(A) tails of mitochondrial transcripts, the full set of aminoacyl-tRNA synthetases (aaRSs), the mitochondrial DNA-directed RNA polymerase (POLRMT), and mitochondrial ribosomal proteins involved in transcription. First, we identified homologs of the genes of interest in several species from the Uniprot database (http://www.uniprot.org/; last accessed 3/5/2016; the list of query sequences is found in Table S1). We then used the BLAST plugin in Geneious to search the transcriptome assembly for A. alata for putative homologs. All BLAST hits with a max e-value of 1e-1 were retained. Retained hits were were aligned and assembled when possible into consensus contigs. The resulting contigs were individually checked against the NCBI protein database (https://blast.ncbi.nlm.nih.gov/Blast.cgi) using blastx to remove potential bacterial contaminants (the top 10 or more hits of bacterial origin) and genes not involved in mitochondrial expression (eukaryote hits with the majority reporting other functions). When the identity of the gene was dubious (not a majority reporting a single function), the translated coding region was blasted against the Uniprot database. Incomplete transcripts were extended, when possible, to include parts of the 5′ and 3′ untranslated regions (UTRs) using a combination of the mirabait script from MIRA v.4,³⁰ MITObim v.1.8,³¹ and the overlap-layout-consensus assembler implemented in Geneious. The presence of putative mitochondrial targeting-signal peptides was assessed using the online versions of the programs MitoProt II,³² TargetP v1.1,³³ Predotar³⁴ and Mitofates.³⁵ Given that the bioinformatics tools for signal peptide prediction are trained on only a handful of model organisms, we considered a signal peptide present when probabilities were >0.5 for at least 2 of the 4 prediction programs. All identified mitochondrion-targeted genes were deposited in the FigShare digital repository and can be accessed at https://dx.doi.org/10.6084/m9.figshare.3380848.

Since two out of 3 components of the mitochondrial ribonucleases complex P (MRPP1-2) are part of large gene families, we used a phylogenetic approach to confirm the identifications of all 3 mitochondrial ribonucleases complexes (MRPP1-3). This approach also allowed us to test the accuracy of our BLAST-based method described above. We first downloaded the sequences that matched “trna methyltransferase 10” and “reviewed:yes” for MRPP1 (240 sequences), “3 hydroxyacyl coa dehydrogenase” and “reviewed:yes” for MRPP2 (282 sequences), and “mitochondrial ribonuclease protein p” for MRPP3 (136 sequences) from the Uniprot database (last accessed 4/5/2016). We then aligned each set with the corresponding sequence identified from Alatina alata using the online program MAFFT v.7³⁶ with the L-INSi algorithm. Gblocks v0.91b³⁷ was used on mrrp2 with gap allowed at all positions; the MRRP1 and MRRP3 alignments were not filtered since Gblocks would remove most (MRPP1) or all (MRPP3) the data given the high level of sequence variability within the alignments. Phylogenetic analyses were performed under a Maximum Likelihood framework using RAxML v.8.2.7³⁸ with the PROTGAMMAWAG model of sequence evolution (100 independent searches for the best tree; 1000 bootstraps). All alignments and trees are deposited in the FigShare digital repository and can be accessed at https://dx.doi.org/10.6084/m9.figshare.3380848.

We obtained the full set of sequences belonging to the polymerase α superfamily (DNA_pol_A PF00476) from the Pfam database³⁹ and created a local HMM database using HMMER v3.1b1.⁴⁰ We searched the transcriptome assembly using HMMER and isolated contigs matching the gamma polymerase gene (polG). We used Geneious to assemble the positive matches into a single contig and checked the assembly by mapping the filtered reads against the polG contig with Bowtie2 (see above).

Results and discussion

Expression of mt genes in Alatina alata

We mapped SOLID DNA-seq data from the SRR archive (SRR217435) to published mitochondrial chromosomes obtained from these data (JN642330-JN642344¹³) and found that all eight mito-chromosomes have similar but unequal (variance 18%) read coverage. That is, no chromosome is found to have considerably more coverage than the others (Fig. 1A; Table S2). Interestingly, we found the genes cox1, cox2, mt-polB, nad4 and orf314 having the highest coverage in the DNA-seq SOLID data (Table 1). The variation in the coverage between various mito-chromosomes in the DNA-seq SOLID data does not appear to correspond to the functional requirements of mitochondrial genes. For instance, despite the key role played by mt-rRNAs, rnl and rns display no noticeable deviation from the average read coverage, in stark contrast with RNA-seq coverage (Fig. 1A, Table 1). It is very likely that all eight mito-chromosomes are in equimolar proportion in the mitochondrion of A. alata.

Figure 1. Mitochondrial read coverage values in the mtDNA of Alatina alata. (A) Nucleotide coverage of the unique regions of individual mitochondrial genes from the DNA-seq and RNA-seq data. (B) Nucleotide coverage of the unique regions of mitochondrial chromosomes from the DNA-seq and RNA-seq data. The Y axis represents mean coverage in kbp.

Table 1. Read coverage of the uniq regions of the mitochondrial genes in the DNA-seq and RNA-seq datasets.

	RNA-seq to shortened gene				DNA-seq to shortened gene
mt gene	size (bp)	read^* (x10^3)	nucl^# (x10^3)	read cov^§	nucl cov^$	size (bp)	read (x10^3)	nucl (x10^3)	read cov	nucl cov
nad2	1199	2.0	196.7	2	164	1198	81.5	4076.2	68	3403
nad5	1591	2.4	237.2	1	149	1591	102.3	5115.2	64	3215
rns	896	108.8	10890.6	121	12155	896	56.4	2819.6	63	3147
nad6	537	1.3	132.8	2	247	537	40.6	2030.9	76	3782
nad3	351	0.8	79.1	2	225	351	26.0	1301.0	74	3707
nad4L	291	0.6	55.8	2	192	291	19.9	995.2	68	3420
nad1	967	2.1	215.0	2	222	969	69.1	3458.2	71	3569
cox2	719	9.2	920.8	13	1281	698	68.0	3402.5	97	4875
atp8	210	2.5	253.3	12	1206	210	17.5	874.9	83	4166
atp6	705	9.5	950.3	13	1348	705	59.5	2977.2	84	4223
cox3	679	10.2	1019.2	15	1501	679	55.9	2797.3	82	4120
cox1	1535	30.0	3001.6	20	1955	1535	168.3	8417.8	110	5484
rnl	1471	243.4	24238.7	165	16478	1473	129.8	6489.1	88	4405
nad4	1358	4.4	442.6	3	326	1360	129.0	6453.4	95	4745
orf314	294	0.2	24.6	1	84	294	26.9	1342.7	91	4567
polB	957	1.3	126.0	1	132	1002	96.3	4815.2	96	4806
cob	1003	6.3	628.9	6	627	1003	66.1	3305.0	66	3295

read

* : raw count of read; nucl

^# : raw count of nucleotide; read cov

^§ : read coverage; nucl cov

^$ : nucleotide coverage.

We recovered nearly complete mitochondrial chromosomes from the RNA-seq reads, including part of the terminal inverted repeats (IRs). Compared to the reference mito-genomes, 10 bp are missing in the RNA-based assemblies, which we define as the origin of transcription (OT). The OT is part of a 30 bp region present on all mito-chromosomes with putative stem-loop structures predicted by the DNA folding program from the mfold Web Server with default parameters⁴¹ (Fig. 2). The complete conservation of the OT and its putative secondary structure support its role as the origin of transcription for the mito-chromosomes. Smith et al (2012)¹³ described a gene-conversion model for the presence of partial duplication at the subtelomeric regions of all mito-chromosomes. The level of conservation of subtelomeric regions varies among mito-chromosomes: the coding regions of seven out of the eight mito-chromosomes start with the 5′end of rnl, with the highest conservation found in nad2 (Fig. S2). A small (∼5 bp) sequence conservation was found between the 5′end of rnl and rns.¹³

Figure 2. Secondary structure of the putative origin of translation (OT) found on all 8 chromosomes of the mt genome of Alatina alata using the DNA folding program from the mfold Web Server with default parameters. The greyed area is present at the DNA level but absent from the RNA-seq data.

The pattern of gene expression corresponds to the organization of the genes on each chromosome, with no break within mitochondrial contigs, as described in other medusozoans,¹⁵ where ≥10 bp stem-loop structures might be responsible in mitochondrial pre-mRNA processing. Such a pattern suggests that mito-chromosomes are expressed into single per-chromosome transcripts (Fig. 3A) starting at the OT region a few base pairs before the IR upstream of the coding region, and ending 161 bp within the IR downstream of the coding region. Average coverage ranged from ∼130 to ∼16,800 nucleotides per base of the unique regions of the genes (the large regions of genes not duplicated on other chromosomes), and the highest values were observed for the mito-chromosomes harboring the rRNA genes (Fig. 1B; Table 1). High-level representation of mt-rRNA in the RNA-seq data is expected given that ribosomes are essential in the expression of the mitochondrial genes. Interestingly, the 2 chromosomes containing cox1-3 and atp6 and 8 genes are also highly expressed in the RNA-seq data (Fig. 1B; Table S2).

Figure 3. Gene expression in the mitochondrion of Alatina alata. (A) Pattern of gene expression and their coverage in the DNA-seq and RNA-seq data sets for a model mito-chromosome. Graphs of read coverage do not correspond to the true values, but rather approximate graphical representations for visualization purposes. CR = putative control region; IGR = inter-genic region; EF? = putative elongation factor involved in the expansion of the rRNA; RNase? = putative enzyme responsible for recognizing and excising stemloop structures found at the 3′end of rRNA and 5′ of tRNA genes. Colored arrows correspond to protein coding genes (green), ribosomal gene (red), tRNA gene (purple), and inverted repeats or IR (cyan). (B) Nuclear genes involved in the maturation of the mitochondrial tRNA-Met: RNase P complex made of three subunits MRPP1-3 responsible for excising extra nucleotides at the 5′end of the tRNA; RNase Z responsible for excising extra nucleotides at the 3′end of the tRNA; CCAse responsible for adding CCA to the 3′end of the tRNA. (C) Import and editing of cytoplasmic tRNA-Trp, where the CCA anticodon is edited via an as yet unknown C-to-U editing enzyme (?). (D) Composition of the mitochondrial ribosome large (39S) and small (28S) subunits. Secondary structures of the large (rnl) and small (rns) rRNA subunits are in red; numbers correspond to individual mitochondrial ribosomal proteins (MRPs) with identified (in black) or lacking (in red) mitochondrial-targeting signal peptides. Number in bold is the newly identified mitochondrial ribosome small subunit S39.

Conservation of subtelomeric regions suggests that transcription can start on any mito-chromosome and an as yet unknown elongation factor is required to recognize the correct gene (Fig. 3A). Mitochondrial ribosomal proteins or rRNA genes might act as elongation factors in the transcription. Mito-chromosome 2 harbors rns along with four protein-coding genes (nad1, 3, 4L and 6), and displays the highest coverage of RNA reads (mean coverage ∼4 kbp; Fig. 1B). Read coverage is not evenly distributed along this mito-chromosome: the protein-coding genes have read coverage ∼55 times lower than rns (Table 1). The dramatic drop in read coverage is observed at the 3′end of rns, which can fold into a stem-loop structure (Fig. S2). We propose that the secondary structure at the 3′end of rns is recognized by some factor(s) and helps regulate transcription by detaching the RNA POLRMT from its template.

We found no read spanning beyond one-fourth (161 out of 702 bp) of the IR downstream of the coding region of the chromosomes. Studies have shown polyadenylation of both rRNA and mRNA in mammalian mitochondria.⁴² In the transcriptome of Alatina alata, we identified and annotated several complete sequences (including both the 5′ and 3′ UTR regions) homologous to the poly(A) polymerase (PAP) gene, one of which (MTPAP PIP5K1A) likely contains an N-terminal mitochondrial-targeting signal peptide (Table 2). Assuming similar polyadenylation is taking place in the mitochondria of all animal species, we posit that the recognition site for MTATP in Alatina is about 1/4^th of the way into the IR of the mitochondrial chromosomes.

Table 2. Read coverage of some mitochondrial-targeted genes in the RNA-seq data set.

gene	CDS size^¶	mt-signal	mt-sign size^⌘	MitoProtII	TargetP	Predotar	MitoFates
polG	3429	+	20/35	0.99	0.88	0.93	0.79
CCAse1	1494	−	33	1.00	0.37	0.02	0.30
ELAC2a	>1134	?	?	?	?	?	?
ELAC2b	>1612	−	39	0.9945	0.443	0.11	0.126
MRRP1	1134	+	15/−	0.86	0.39	0.51	0.71
MRRP2	774	−	14/69	0.41	0.65	0.13	0.27
MRRP3	1941	+	−/48	0.97	0.88	0.75	0.04
MTPAP PIP5K1A	2160	+	16/9	0.90	0.68	0.63	0.20
POLRMT	3681	+	47/18	0.98	0.79	0.07	0.28
PTCD1	1527	−	−/−	0.1953	0.183	0.07	0.464

Presence of mt-signal is assumed when at least two prediction software display probabilities >0.50 (in bold). CDS size

^¶ : size of the coding regions found in the transcriptome assembly; mt-sign size

^⌘ : number of amino acids composing the mt-targeting signal peptide as predicted by MitoProt II v1.01 and TargetP v1.1. Probability of mitochondrial targeting were obtained using MitoProt II v1.01, TargetP v1.1, Predotar v.1.03 and MitoFates.

A functional mitochondrial-encoded type B DNA polymerase

We recovered the complete sequence of the polG gene whose product is responsible for mtDNA replication, including the mitochondrial-targeting signal peptide (Table 2). The product of this gene is similar in size to its counterparts in other animals. Medusozoan mito-genomes (excluding hydroidolinan hydrozoans) harbor two additional protein-coding genes orf314 and mt-polB. Our previous studies^12,13 have suggested that orf314 and mt-polB play roles in the maintenance and replication of the mitochondrial telomeres: the product of mt-polB has affinities with members of the DNA-directed DNA polymerase type B found in plant and fungal mitochondria, and viruses; orf314 was hypothesized to encode for a DNA-binding protein that could sit at the 5′ single-stranded end of mitochondrial telomeres to prevent their degradation by exonucleases and may be involved in the maintenance of the integrity of the mito-chromosomes during replication.^7,12 It has been posited that the linearization of medusozoan mtDNA originated in association with the incorporation of a linear plasmid containing both mt-polB and orf314.¹² The presence of these genes on a separate mito-chromosome in A. alata was already suggested as evidence of positive selection for their maintenance.¹³ We recovered the chromosome containing these genes from the RNA-seq data with per-gene coverage values similar but below that of the other protein coding genes (Table 1 and S2), indicating that these genes are expressed and polyadenylated (i.e., present in the mRNA pool). These data provide further evidence that both orf314 and mt-polB are functional in the mitochondrion.

Near-parallel loss of mt-tRNA and mt-aaRS in Cubozoa suggests tRNA editing

The only complete mitochondrial genome for Cubozoa published to date, that of Alatina alata (JN642330-JN642344), encodes one single mitochondrial tRNA for methionine (mt-tRNA^Met or trnM),¹³ suggesting that the other tRNAs are imported from the cytoplasm. We recovered the complete set of cytoplasmic aminoacyl-tRNA synthetases (nc-aaRS) but only one mitochondrial aminoacyl-tRNA synthetase (mt-PheRS) from the RNA-seq data (Table 3). This pattern is in accord with the evolutionary scenario of parallel loss of mt-tRNAs and mt-aaRS, excluding mt-PheRS, which is made of a single protein as opposed to nc-PheRS composed of 2 subunits encoded by 2 separate genes.^29,43 Three full transcripts (Glu/ProRS, IleRS and ValRS) did not harbor a recognizable 5′ mitochondrial-targeting signal peptide (Table 3), while the sequence for HisRS was missing the 5′ end where the mt-targeting signal peptide is expected to be found. The absence of mt-targeting signal peptides in Glu/ProRS, IleRS and ValRS suggests some post-transcriptional modification that might be responsible for the transport of these nc-aaRSs into the mitochondrion, though the possibility of a deviant targeting signal not identified by the prediction programs cannot be ruled out.

Table 3. List of aminoacyl-tRNA synthetases (aaRS) found in the RNA-seq data of Alatina alata.

name	CDS size^¶	mt-signal	mt-sign size^⌘	MitoProtII	Target P	Predotar	MitoFates
A	2985	+	29/29	0.97	0.89	0.81	0.99
C	2343	+	16/100	0.62	0.56	0.47	0.22
D	1803	+	72/30	0.99	0.81	0.43	0.69
E/P	4371	−	NA	0.04	0.11	0.01	0.00
mt-F	1356	+	37/37	0.99	0.76	0.02	0.46
nc-Fa	1512	−	26/−	0.13	0.08	0.00	0.01
nc-Fb	1782	−	NA	0.25	0.06	0.00	0.00
G	2283	+	44/35	0.99	0.59	0.72	0.60
H	>1578	?	NA	NA	NA	NA	NA
I	3543	−	NA	0.01	0.06	0.00	0.02
K	1893	+	35/42	0.82	0.61	0.68	0.61
L	3294	+	31/22	0.99	0.47	0.52	0.72
M	2970	+	28/12	0.98	0.94	0.72	0.92
N	1896	+	39/−	0.97	0.50	0.10	0.83
Q	2427	+	26/26	0.97	0.84	0.51	0.78
R	2175	+	55/56	1.00	0.90	0.90	1.00
S	1617	+	14/27	0.98	0.93	0.66	0.62
T	2172	+	50/50	0.99	0.72	0.29	0.84
V	3288	−	NA/−	0.14	0.50	0.38	0.79
W	1413	+	27/47	0.99	0.76	0.41	0.89
Y	1797	+	57/58	0.99	0.89	0.86	0.58

Aminoacyl-tRNA synthetases are listed according to their one letter name. Presence of mt-signal is indicated when at least two prediction software display probabilities >0.50 (in bold). mt-F: monomeric mitochondrial aaRS for phenylalanine; nc-Fa and nc-Fc: α and β subunits of the dimeric cytoplasmic aaRS for phenylalanine. CDS size

^¶ : size of the coding regions found in the transcriptome assembly; mt-sign size

^⌘ : number of amino acids composing the mt-targeting signal peptide as predicted by MitoProt II v1.01 and TargetP v1.1. Probabilities of mitochondrial targeting were obtained using MitoProt II v1.01, TargetP v1.1, Predotar v.1.03 and MitoFates.

Previous studies have suggested that the loss of mt-tRNAs preceded that of the corresponding mt-aaRSs in Cnidaria⁴³; cases where an mt-tRNA is present in the absence of the corresponding mt-aaRS can be explained by recruitment of the mt-tRNA for other functions, such as a structural component of the mitochondrial ribosome.²⁹ In mammals, where mt-tRNA^Val has been found as a structural part of the mitochondrial ribosome, the tRNA gene (trnV) is encoded between the small and large ribosomal subunits (rns and rnl respectively) in the mt genome, leading to a polycistronic transcript containing all three RNAs.^44,45 In demosponge mtDNA, trnG (mt-tRNA^Gly) and trnV are the most common tRNA genes found between rns and rnl, but never trnR (mt-tRNA^Arg),⁴⁶ yet mt-ArgRS is one of the five mt-aaRSs (the others being mt-CysRS, mt-LysRS, mt-ThrRS and mt-ValRS) that seem to have been lost in this clade.²⁹ This pattern cannot explain the role of these tRNAs as structural parts of the ribosome as shown in mammals.⁴⁴ The presence of mt-tRNA^Met on three out of the eight mitochromosomes of A. alata concurrent with the absence of mt-MetRS in the nucleus is an interesting pattern that cannot be explained by the role of this tRNA as a structural component of the mitochondrial ribosome.

We investigated several key genes responsible for mitochondrial tRNA processing, namely RNase P (for 5′-end maturation), RNase Z (for 3′-end maturation), and CCase (responsible for adding CCA to the 3′-end of the tRNA)(Table 2). In humans, RNase P is a protein complex composed of three subunits: MRPP1, one out of at least three homologs of a tRNA methyltransferase 10 (C homolog); MRPP2, a 3-hydroxyacyl-CoA dehydrogenase; MRPP3, the mitochondrial ribonuclease P protein 3.^47,48 We identified homologs to both MRPP2 and MRPP3 (the latter with a mitochondrial targeting signal peptide) but not to MRPP1 (Table 2). Instead, our phylogenetic analysis identified a tRNA methyltransferase 10 A homolog (Fig. S3) with a mitochondrial targeting signal peptide, while human homologs do not harbor such a signal peptide (data not shown). We therefore decided to label the transcript tRNA methyltransferase 10 homolog A as MRPP1 (Table 2). Phylogenetic analyses of gene trees for transcripts annotated MRPP1-3 confirm they being homologous to MRPP1-3 genes derived from other organisms (Figs. S3-5). Lopez Sanchez et al. (2011) recently showed that MRPP1 and MRPP3 are sufficient for processing the 5′end of mitochondrial tRNAs in humans.⁴⁹

The same authors and others identified the zinc phosphodiesterase ELAC protein 2 as RNase Z responsible for processing the 3′end of mitochondrial tRNAs.^49,50 We recovered two transcripts that encode for proteins similar to ELAC2: the shorter transcript (ELAC2a) displayed more similarity (in the form of higher matching base and alignment length BLAST results in the Uniprot database) with homologs from other taxa but was incomplete, missing the 5′ terminal; the longer transcript (ELAC2b) lacked a mitochondrial targeting signal peptide (Table 2). Finally, we identified one CCAse homolog that lacked an acceptable mitochondrial targeting signal peptide (Table 2).

The presence of mitochondrial targeting signals on the transcripts of several mitochondrial tRNA-processing genes combined with the high expression of mt-tRNA^Met (Table S3) suggests that trnM requires a functional tRNA structure in the mitochondrion of this species (Fig. 3B). Given the very distinct bacterial-like translation initiation system in the mitochondrion compared to the eukaryote-like one for the cytosol, the loss of trnM is thought to be nigh impossible.⁵¹ Yet, the loss of mt-MetRS requires nc-MetRS, whose transcript harbors an mt-targeting signal peptide (Table 3) to recognize and aminoacetylate mt-tRNA^Met in the mitochondrion. Alternatively, trnM could play a role as initiation of translation in the mitochondrion.

Both Octocorallia and Cubozoa have lost both the mitochondrial tRNA and the mitochondrial aminoacyl-tRNA synthetase for trpyptophan (mt-tRNA^Trp and mt-TrpRS, respectively).²⁹ In the absence of both mt-tRNA^Trp and mt-TrpRS, tryptophan residues are processed by a combination of cytoplasmic nc-tRNA^Trp and nc-TrpRS. In the ctenophore Pleurobrachia bachei, the loss of mt-tRNA^Trp has resulted in the reassignment of the TGA codon to serine.^29,52,53 In cubozoan mtDNA, tryptophan residues are encoded by two codons TGG and TGA at a relatively similar proportion (Table S3). In octocorals, TGG is the main codon used for tryptophan, the exceptions seemingly resulting either from sequencing errors or present in the extra protein coding gene homolog to bacterial mutS⁵⁴(Table S3). In the cytoplasm, nc-tRNA^Trp only recognizes TGG, while TGA is a stop codon. In the absence of mRNA editing in Alatina (proportion of TGA was similar in RNA and DNA-seq data; Table S3), we suggest that the anticodon of nc-tRNA^Trp(CCA) is edited in the mitochondrion to recognize and aminoacylate TGA codons (anticodon UCA) through a C → U editing mechanism (Fig. 3C). Anticodon C to U editing has been reported in the mitochondrion of Trypanosoma allowing the recognition of TGA as tryptophan in the absence of an mt-tRNA^Trp(UCA) gene,⁵⁵ and in marsupials where some of the mt-tRNA^Asp(GCC) are edited to mt-tRNA^Asp(GUC).⁵⁶ Alternatively, addition factors (RNA and/or proteins) would be required to guide TGA recognition by a native nc-tRNA^Trp. Future studies focusing on short RNAs should allow testing our hypotheses.

The mitochondrial ribosomal proteins complex

Ribosomes, essential pieces of the translation machinery, are complex molecules made of two subunits (large and small), each composed of RNA (of different origin for the cytoplasm and the mitochondrion) and many proteins. In animals, all mitochondrial proteins are encoded by the nucleus, while two RNA subunits (12S and 16S for small and large RNA subunits, respectively) are encoded by the mitochondrial DNA (rns and rnl). We used the dataset of a recent study on Holozoa⁵⁷ to investigate ribosomal protein coding genes in the transcriptome of A. alata. We recovered the same cnidarian set of large (40 gene transcripts) and small (29 gene transcripts) subunits of the mitochondrial ribosome, with a few exceptions, plus the unassigned MRP63 (Fig. 3D; Table 4).

Table 4. List of mitochondrial ribosomal proteins (MRPs) found in the RNA-seq data of Alatina alata.

MRPs	CDS size^¶	mt-signal	mt-sign size^⌘	MitoProtII	TargetP	Predotar	MitoFates
L2	966	+	22/69	0.8928	0.894	0,41	0.883
L3	1131	+	46/39	0.9912	0.831	0,42	0.540
L4	909	+	47/26	0.9969	0.884	0,66	0.989
L9	681	+	34/9	0.7681	0.560	0,50	0.655
L10	651	+	51/109	0.9916	0.647	0,50	0.800
L11	459	−	−/−	0.3073	0.216	0,44	0.039
L12	597	+	46/46	0.9672	0.729	0,22	0.620
L13	510	−	−/103	0.3097	0.758	0,30	0.068
L14	576	+	111/27	0.9548	0.759	0,47	0.615
L15	795	−	28/−	0.6694	0.400	0,15	0.488
L17	486	+	−/13	0.7122	0.718	0,21	0.450
L18	690	+	45/45	0.6141	0.842	0,55	0.698
L19	984	−	−/−	0.6420	0.190	0,09	0.175
L20	381	+	71/38	0.9999	0.876	0,64	0.878
L21	705	+	45/71	0.9902	0.890	0,55	0.818
L22	612	+	30/30	0.9904	0.924	0,94	0.996
L23	471	+	−/5	0.8558	0.812	0,12	0.500
L24	738	+	21/119	0.9780	0.641	0,68	0.365
L28	519	+	33/15	0.9809	0.616	0,77	0.906
L29/L47	597	+	28/64	0.9374	0.872	0,52	0.694
L30	429	−	−/−	0.6505	0.272	0,21	0.013
L32	318	+	13/50	0.5893	0.466	0.26	0.660
L33	378	−	−/15	0.4681	0.079	0,00	0.026
L34	264	+	23/41	0.9761	0.787	0,88	0.950
L35	480	+	47/109	0.9755	0.928	0,53	0.796
L36	216	+	50/39	0.9912	0.717	0,35	0.794
L37a	1053	+	−/35	0.9135	0.675	0,09	0.148
L37b	1143	+	11/36	0.9610	0.794	0,52	0.380
L38	975	+	−/29	0.5728	0.753	0,01	0.915
L40	540	−	17/−	0.9278	0.151	0,48	0.006
L41	303	+	32/42	0.8645	0.595	0,63	0.871
L43	423	−	−/−	0.9212	0.262	0,22	0.377
L44	1002	+	19/12	0.9418	0.706	0,67	0.751
L45	801	+	35/26	0.9767	0.858	0,47	0.062
L46	909	+	11/13	0.7244	0.624	0,32	0.840
L48	600	+	13/45	0.9404	0.655	0,52	0.775
L49	468	+	25/44	0.9936	0.907	0,46	0.196
L52	345	−	35/−	0.8568	0.391	0,47	0.354
L54	357	+	19/19	0.9452	0.631	0,26	0.829
L55	420	+	−/6	0.8509	0.777	0,44	0.482
S2	810	+	21/36	0.9895	0.944	0,77	0.945
S3/S24	495	+	46/28	0.9719	0.815	0,73	0.750
S5	1281	+	14/−	0.8815	0.448	0,23	0.828
S6	363	−	−/−	0.2581	0.094	0,05	0.132
S7	711	+	15/−	0.8070	0.748	0,83	0.427
S9	1170	+	38/63	0.9000	0.643	0,01	0.914
S10	747	+	37/48	0.9680	0.839	0,03	0.976
S11	648	+	61/61	0.9905	0.490	0.39	0.613
S12	606	+	153/83	0.9986	0.770	0,81	0.799
S14	486	+	14/53	0.6913	0.816	0,41	0.751
S15	657	+	14/101	0.6887	0.664	0,45	0.719
S16	375	+	37/10	0.9752	0.846	0,81	0.553
S17	279	−	−/−	0.9399	0.144	0,10	0.028
S18A	393	+	26/17	0.9935	0.839	0,59	0.985
S18B	546	+	38/31	0.9879	0.959	0,40	0.983
S18C	501	+	38/55	0.8529	0.713	0,22	0.433
S21	270	+	40/28	0.6619	0.622	0,06	0.135
S22	1068	+	45/29	0.9663	0.799	0,57	0.248
S23	567	+	15/6	0.6029	0.458	0,11	0.668
S25a	492	+	−/52	0.5498	0.425	0,07	0.526
S25b	483	+	−/−	0.7713	0.468	0,64	0.111
S26	771	+	27/19	0.9902	0.905	0,58	0.810
S28	498	−	33/−	0.8869	0.353	0,26	0.387
S29	1215	+	37/38	0.9914	0.824	0,22	0.480
S31	2307	+	30/48	0.8050	0.781	0,86	0.500
S33	312	+	31/38	0.8297	0.589	0,02	0.477
S35	1047	−	39/−	0.6656	0.327	0,00	0.159
S36	376	+	21/20	0.9792	0.802	0.57	0.650
S39	2043	+	47/40	0.8781	0.832	0.49	0.944
MRP63	354	+	27/−	0.9871	0.444	0,89	0.305

Presence of mt-signal is indicated when at least two prediction software display probabilities >0.50 (in bold). CDS size

^¶ : size of the coding regions found in the transcriptome assembly; mt-sign size

^⌘ : number of amino acids composing the mt-targeting signal peptide as predicted by MitoProt II v1.01 and TargetP v1.1. Probabilities of mitochondrial targeting were obtained using MitoProt II v1.01, TargetP v1.1, Predotar v.1.03 and MitoFates.

From the large ribosomal subunit, L1, L16, L27, L50, L53 and L56 were not recovered. From the small ribosomal subunit, only S30 was not recovered. We also found that the subunits L11, L13, L15, L19, L30, L33, L40, L43 and L52 from the large subunit, and S6, S17, S28 and S35 from the small subunit had low probability (≤0 .5 for at least three out of the four prediction software used) for possessing mitochondrial-targeting signal peptides (Table 4). Sheel et al. (2014) also found L50 and L56 to be absent from the hydrozoan Hydra magnapapillata.⁵⁷ Assuming the lack of these genes reflects their loss in the cubozoan A. alata, the loss of L56 could be a synapomorphy of Medusozoa (Hydra and Alatina), while the loss of L1, L16, L27 and S30 appear to be autapomorphies to Alatina.

Interestingly, we found two versions of subunits L37 and S25, suggesting gene duplication, though an assembly artifact is conceivable. Both candidates for subunit L37 had recognizable mitochondrial-targeting signal peptides but the longer transcript (L37b) displayed higher probability of export to the mitochondrion, while the shorter transcript (L37a) showed more similarity (in the form of higher matching base and alignment length BLAST results in the Uniprot database) to homologs from other taxa (Table 4). Similarly, while we identified a mitochondrial-targeting signal peptide on the two candidate transcripts for subunit S25, the longer transcript (S25b) showed more similarity to homologs from other taxa (Table 4). We also found the transcript for subunit S31 to be far larger in size than homologs from other metazoans, with an addition of ∼840 nucleotides (∼280 amino acids) on the 5′-end of the coding DNA sequence (CDS). The function of the extra amino acids in S31 was not determined.

Lastly, we identified another mitochondrial ribosomal protein subunit (S39) with similarities to pentatricopeptide repeat-containing protein 3 previously not found by Sheel and Hausdorf (2014).⁵⁷ Using the BLAST search tool on the Uniprot database, we identified homologs from Nematostella vectensis (Cnidaria: Anthozoa; predicted protein v1g212337), Hydra vulgaris (Cnidaria: Hydrozoa; pentatricopeptide repeat-containing protein 3, mitochondrial T2M872) and Amphimedon queenslandica (Porifera: Demospongiae; uncharacterized protein LOC100633076) with mitochondrial targeting signal peptides (data not shown).

Conclusion

Cubozoans have one of the most unusual mtDNA organizations in Metazoa, where genes are distributed on eight linear chromosomes with long terminal inverted repeats. We explored the recently published RNA-seq data from the box jellyfish Alatina alata to better understand the mechanisms involved in the expression of linear mitochondrial chromosomes. We propose a pattern of mitochondrial gene expression where each mito-chromosome is transcribed into a single transcript starting at a well conserved origin of transcription (OT; 30 bp) and ending with the polyadenylation site about 1/4^th of the way into the inverted repeat (IR) downstream of the coding regions of the mito-chromosomes. We captured the two extra protein-coding genes mt-polB and orf314 in the RNA-seq data with coverage levels similar to other mitochondrial protein-coding genes, providing further evidence for the functionality of these genes. We also found that the evolution of mitochondrial translational machinery in A. alata displays similarities to what has been documented in other non-bilaterians.²⁹ Finally, we suggest that the mt-tRNA^Met is likely a functional tRNA that could be recognized and aminoacetylated by the nc-MetRS, and that an imported nc-tRNA^Trp is edited in order to recognize TGA codons as tryptophan. Future studies exploring the repertoire of small RNAs and the proteome of the mitochondrion will be able to test these hypotheses.

Disclosure of potential conflicts of interest

No potential conflicts of interest were disclosed.

Acknowledgments

Some of this work was performed using resources of the Laboratories of Analytical Biology at the Smithsonian National Museum of Natural History. EK and BB were supported by a Smithsonian Institution Peter Buck postdoctoral fellowship. We also want to thank two anonymous reviewers for providing constructive criticisms to earlier versions of this manuscript.