Metabolic Specialization and Codon Preference of Lignocellulolytic Genes in the White Rot Basidiomycete Ceriporiopsis subvermispora

Ceriporiopsis subvermispora is a white-rot fungus with a high specificity towards lignin mineralization when colonizing dead wood or lignocellulosic compounds. Its lignocellulose degrading system is formed by cellulose hydrolytic enzymes, manganese peroxidases, and laccases that catalyze the efficient depolymerization and mineralization of lignocellulose. To determine if this metabolic specialization has modified codon usage of the lignocellulolytic system, improving its adaptation to the fungal translational machine, we analyzed the adaptation to host codon usage (CAI), tRNA pool (tAI, and AAtAI), codon pair bias (CPB), and the number of effective codons (Nc). These indexes were correlated with gene expression of C. subvermispora, in the presence of glucose and Aspen wood. General gene expression was not correlated with the index values. However, in media containing Aspen wood, the induction of expression of lignocellulose-degrading genes, showed significantly (p < 0.001) higher values of CAI, AAtAI, CPB, tAI, and lower values of Nc than non-induced genes. Cellulose-binding proteins and manganese peroxidases presented the highest adaptation values. We also identified an expansion of genes encoding glycine and glutamic acid tRNAs. Our results suggest that the metabolic specialization to use wood as the sole carbon source has introduced a bias in the codon usage of genes involved in lignocellulose degradation. This bias reduces codon diversity and increases codon usage adaptation to the tRNA pool available in C. subvermispora. To our knowledge, this is the first study showing that codon usage is modified to improve the translation efficiency of a group of genes involved in a particular metabolic process.


Introduction
The main carbon source synthesized through photosynthesis which plays a central role in the carbon cycle of the planet is lignocellulose. Lignin, one of the compounds of lignocellulose, is a recalcitrant, aromatic, and amorphous polymer that protects lignocellulose from microbial attack. A small group of filamentous fungi from the basidiomycete phylum is unique in its ability to efficiently degrade lignocellulose [1]. Together they are collectively known as white rot fungi, developing an enzymatic machinery that allows degradation of the three main components of lignocellulose: lignin, cellulose, and hemicellulose. White-rot fungi mineralize lignin as a strategy to access cellulose and

Calculation of Relative Synonymous Codon Usage (RSCU) and Codon Adaptation Index (CAI)
The Relative Synonymous Codon Usage (RSCU) and Codon Adaptation Index (CAI) were calculated with the Emboss program [30] based on the frequency of codon usage of C. subvermispora.

Determination of Adaptation to the tRNA Pool
The adaptation to the tRNA pool of genes present in C. subvermispora was determined by calculating the values of tAI and AAtAI [31]. We calculated tAI as established by dos Reis [24] (Equations (3) and (4)), estimating the relative abundance of tRNAs from the number of copies of each tRNA gene in the C. subvermispora genome. The number of copies of each tRNA gene was determined with the tRNAscan-SE program [25] in Linux, which was used to analyze the unmasked assembled C. subvermispora genome.
W i is the relative adaptiveness of the ith codon to the tRNA pool, n i is the number of tRNA isoacceptors that recognize the ith codon, tGCN ij is the number of copies of the jth tRNA gene that recognizes the ith codon, and S ij is the selective constraint in the efficiency of codon-anticodon pairing.
The adaptation of a gene to the tRNA pool is calculated according to Equation (4), in which w i is defined as the ratio between W i and W max (W i /W max ), ikg is the codon defined by the kth triplet of gene g, and l g is the length of gene g in codifying codons.
AAtAI g = l g k=1 AA w ikg 1/l g (5) AA w i is defined as the ratio between W i and W AAmax (W i /W AAmax ), where W AAmax is the highest value of Wi among codons coding for the same amino acid. ikg is the codon defined by the kth triplet of gene g, and l g is the length of gene g in codifying codons. As AAtAI is calculated from Equation (5), which is similar to the calculation of CAI, we used the EMBOSS software [30], entering Wi data in replacement of the frequency of codon usage. Wi was calculated using the procedure described for tAI in Equation (3).

Determination of effective Number of Codons (Nc)
The number of effective codons for the C. subvermispora genes was calculated using the CodonW program (http://codonw.sourceforge.net/).

Phylogenetic Analysis
Multialignment between sequences of tRNA genes was performed using ClustalW [32]. The parameters were set up to align sequences using IUB as substitution matrix. The evolutionary history was inferred using the Neighbor-Joining method [33]. Interior branch test with a bootstrap of 1000 was used to analyze confidence of the tree [34]. The evolutionary distances were computed using the Maximum Composite Likelihood method [35]. The rate of variation among sites was modeled with a γ distribution (shape parameter = 1). Multialignment and evolutionary analyses were conducted in MEGA5 [36].

Graphs and Statistical Methods
The program SigmaPlot 11 was used for graphs and statistical tests. The significance of the differences or correlations among the data groups obtained was evaluated with the Rank Sum Test non-parametric test for comparing two groups and the non-parametric Spearman Rank Order test for correlations, using a p-value < 0.05 as a cutoff.

C. subvermispora tRNAs
Genome analysis of C. subvermispora by tRNAScan-SE identified a total of 192 tRNAs in 32 scaffolds ( Figure 1). About 72% of the tRNA genes presented introns ( Table 1). The scaffold with the highest number of tRNA genes was scaffold 1, which contains 20 copies of various tRNA genes (Table S2). The tRNA with the highest number of gene copies corresponded to tRNAs charging glycine, with 17 gene copies distributed in eight scaffolds (scaffolds 1, 2, 5, 7, 9, 10, 13, and 20). The tRNAs for cysteine and tryptophan presented the lowest number of gene copies, each with three gene copies in three scaffolds (Table 1).

tRNA Abundance and Codon Usage in C. subvermispora
The expansion of certain families of tRNAs in the genome of C. subvermispora could be the result of an evolutionary pressure to increase their expression. In some organisms such as E. coli and yeast, the number of copies of a tRNA gene is proportional to the abundance in the genome of the decoded codon [24]. This proportionality is explained because during translation process, aminoacyl-tRNAs that decode frequently used codons have a higher rate of consumption. To sustain adequate translation, cells must balance synthesis and consumption rates of aminoacyl-tRNAs. The increase in copy number of a gene that encodes a tRNA that recognizes a frequently used codon, allows increasing expression of this tRNA and its aminoacylated form to balance its rate of consumption.
To test if the expansion of certain tRNA families in C. subvermispora is related to the increment of the specific codons, we assessed whether there is a correlation between the frequency of codon usage and the amount of tRNA genes that decode these most highly used codons. We identify a positive correlation between these parameters (ρ = 0.406, p = 0.0016, n = 61). Dos Reis et al. show that the Relative Adaptiveness to the tRNA pool (w) which takes in account that codons can be recognized by anticodons with perfect or imperfect match (wobble codon-anticodon recognition rules) with different affinities is a better parameter to measure the adaptation of a codon to their decoding tRNAs than the absolute number of tRNA [24]. When the frequency of codon usage was correlated with the Relative Adaptiveness to the tRNA pool (w), an increased correlation was observed (ρ = 0.459, p = 2.2 × 10 −4 , n = 61) ( Figure 2).

Phylogenetic Analysis of C. subvermispora tRNA Genes
To determine whether the high number of tRNA genes charging the same amino acid corresponds to related gene copies, a phylogenetic reconstruction and evolutionary distance calculation were performed using the tRNA sequences identified by the tRNAscan-SE program. Phylogenetic reconstruction indicated that most tRNA genes that code for the same amino acid group together, with the exception of the tRNA genes that charge arginine, valine, and alanine ( Figure 1), which form two groups in each case. For tRNAs loading arginine, Group I comprises genes presenting anticodons with the WCG sequence. In contrast, group II presents a YCK anticodon consensus sequence. In tRNAs that load valine, Group I corresponds to tRNA genes without introns, while group II includes all valine tRNA genes containing introns. Genes coding for group I of alanine tRNAs exhibit anticodons with the consensus sequence YGC, whereas, the anticodon sequence is AGC for group II (Figure 1). The tRNA genes 66 and 155 also show a different pattern. tRNAScan prediction indicates that the amino acid loaded by tRNA155 should be serine, however, the sequence of this gene grouped with threonine charging tRNA genes. Additionally, the tRNA66 gene is expected to load isoleucine, though this gene does not group with isoleucine charging tRNAs ( Figure 1). To identify tRNA genes that are repeated, the evolutionary distance between different tRNA genes was calculated. Genes with values of evolutionary distance equal to 0.000 were selected. This analysis identified 15 tRNA genes that are repeated between two and ten times. The group of tRNA genes that showed the greatest expansion corresponds to those tRNAs carrying glutamic acid and glycine. One glycine tRNA gene is repeated twice and the other is repeated ten times. The tRNA gene for glutamic acid is also repeated ten times (Table S3). These genes are scattered along the genome, with the exception of tRNA genes 91, 92, 94, 95, 96, and 97, which code for glycine and are located at adjacent positions.

tRNA Abundance and Codon Usage in C. subvermispora
The expansion of certain families of tRNAs in the genome of C. subvermispora could be the result of an evolutionary pressure to increase their expression. In some organisms such as E. coli and yeast, the number of copies of a tRNA gene is proportional to the abundance in the genome of the decoded codon [24]. This proportionality is explained because during translation process, aminoacyl-tRNAs that decode frequently used codons have a higher rate of consumption. To sustain adequate translation, cells must balance synthesis and consumption rates of aminoacyl-tRNAs. The increase in copy number of a gene that encodes a tRNA that recognizes a frequently used codon, allows increasing expression of this tRNA and its aminoacylated form to balance its rate of consumption.
To test if the expansion of certain tRNA families in C. subvermispora is related to the increment of the specific codons, we assessed whether there is a correlation between the frequency of codon usage and the amount of tRNA genes that decode these most highly used codons. We identify a positive correlation between these parameters (ρ = 0.406, p = 0.0016, n = 61). Dos Reis et al. show that the Relative Adaptiveness to the tRNA pool (w) which takes in account that codons can be recognized by anticodons with perfect or imperfect match (wobble codon-anticodon recognition rules) with different affinities is a better parameter to measure the adaptation of a codon to their decoding tRNAs than the absolute number of tRNA [24]. When the frequency of codon usage was correlated with the Relative Adaptiveness to the tRNA pool (w), an increased correlation was observed (ρ = 0.459, p = 2.2 × 10 −4 , n = 61) ( Figure 2).
Synonymous codons are not used equally in an organism. As one tRNA can decode several synonymous codons with different affinities, the expansion of some tRNA families in C. subvermispora may be related to the preferential use of certain synonymous codons in coding regions of C. subvermispora.
To assess this hypothesis, we correlated RSCU values with Relative Adaptiveness to the tRNA pool. Non-statistical correlation was observed, in part because w values are normalized with respect to the tRNA with the highest number of genes able to decode it (W i /W max ) and not with respect to the pool of tRNAs that decode the complete set of synonymous codons. When W values were normalized with respect to the total amount of tRNAs that decode a set of synonymous codons, a strong correlation with the RSCU values (ρ = 0.628, p = 0, n = 61) was observed. This suggests that among synonymous codons, those highly represented in the C. subvermispora genome tend to be decoded by those tRNAs with a high number of gene copies. synonymous codons with different affinities, the expansion of some tRNA families in C. subvermispora may be related to the preferential use of certain synonymous codons in coding regions of C. subvermispora. To assess this hypothesis, we correlated RSCU values with Relative Adaptiveness to the tRNA pool. Non-statistical correlation was observed, in part because w values are normalized with respect to the tRNA with the highest number of genes able to decode it (Wi/Wmax) and not with respect to the pool of tRNAs that decode the complete set of synonymous codons. When W values were normalized with respect to the total amount of tRNAs that decode a set of synonymous codons, a strong correlation with the RSCU values (ρ = 0.628, p = 0, n = 61) was observed. This suggests that among synonymous codons, those highly represented in the C. subvermispora genome tend to be decoded by those tRNAs with a high number of gene copies.

Relationship between Gene Expression Level, Codon Bias, and Translational Efficiency in C. subvermispora
The bias in codon usage and adaptation to the tRNA pool modulates translational efficiency. Thus, highly expressed genes tend to use codons that are over-represented in the genome, which in turn present greater availability of tRNAs [24]. To determine whether this relationship exists in C. subvermispora, expression levels of C. subvermispora genes were correlated with their adaptation values to the tRNA pool and to codon bias. Adaptation to the tRNA pool was evaluated using two indexes: (i) tAI, which measures adaptation of the tRNA pool compared to the relative amount of each tRNA gene, and (ii) AAtAI, which evaluates whether a gene preferentially uses the most abundant tRNA charging a particular amino acid. Codon bias was analyzed by calculation of CAI (Codon Adaptation Index), and CPB (Codon Pair-Bias), to assess if gene expression is correlated with bias in usage of codons or of codon pairs. Codon bias also was evaluated using Nc value to determine if the transcription level is associated with a decrease in the diversity of codons used. We analyzed expression levels determined by RNAseq published by Hori et al. 2014 which was obtained from C. subvermispora grown on Ball-Milled Aspen medium [17]. CAI, tAI, AAtAI, and CPB showed high degrees of correlation among them (Table S4), however, these indicators exhibit very low correlation coefficients with expression levels reported by Hori. The high correlation between the different indicators of codon bias and translational efficiency indicates that the frequency with which codons are used in C. subvermispora is related to the abundance in the genomes of the tRNAs that decode these codons.

Relationship between Gene Expression Level, Codon Bias, and Translational Efficiency in C. subvermispora
The bias in codon usage and adaptation to the tRNA pool modulates translational efficiency. Thus, highly expressed genes tend to use codons that are over-represented in the genome, which in turn present greater availability of tRNAs [24]. To determine whether this relationship exists in C. subvermispora, expression levels of C. subvermispora genes were correlated with their adaptation values to the tRNA pool and to codon bias. Adaptation to the tRNA pool was evaluated using two indexes: (i) tAI, which measures adaptation of the tRNA pool compared to the relative amount of each tRNA gene, and (ii) AAtAI, which evaluates whether a gene preferentially uses the most abundant tRNA charging a particular amino acid. Codon bias was analyzed by calculation of CAI (Codon Adaptation Index), and CPB (Codon Pair-Bias), to assess if gene expression is correlated with bias in usage of codons or of codon pairs. Codon bias also was evaluated using Nc value to determine if the transcription level is associated with a decrease in the diversity of codons used. We analyzed expression levels determined by RNAseq published by Hori et al. 2014 which was obtained from C. subvermispora grown on Ball-Milled Aspen medium [17]. CAI, tAI, AAtAI, and CPB showed high degrees of correlation among them (Table S4), however, these indicators exhibit very low correlation coefficients with expression levels reported by Hori. The high correlation between the different indicators of codon bias and translational efficiency indicates that the frequency with which codons are used in C. subvermispora is related to the abundance in the genomes of the tRNAs that decode these codons.

Transcriptional Response to Growth on Ball-Milled Aspen (BMA), Codon Bias, and Translational Efficiency
Growth of C. subvermispora in natural environments is dependent on wood. In 2012, Fernandez-Fueyo et al. [8] reported microarray experiments that compared gene expression of C. subvermispora grown on glucose and on Ball-Milled Aspen (BMA) as carbon sources. Saline media with BMA has been used as a laboratory medium that mimics growth on wood to analyze expression of genes that are transcriptionally regulated by growth on wood. To analyze if genes regulated by conditions that mimic growth on wood, such as BMA have a different adaptation to the tRNA pool or codon bias, we used the microarray data published by Fernandez-Fueyo [8] and defined four groups of genes: group A corresponds to genes where expression was reduced at least 2 times with a p-value lower than 0.05. Group B includes all genes which showed increased expression of at least 2 times with a p-value lower than 0.05. Group C corresponds to all genes with non-significant differences (p > 0.05) and group D contains all genes with low changes in expression (<2 fold) that are statistically significant. (p < 0.05). Our results show that group B has lower values of Nc and higher values of CAI, tAI, AAtAI, and CPB. Groups A, C, and D show non-significant differences among them. This implies that genes induced by wood preferentially use a reduced set of codons that are better adapted to the tRNA pool present in C. subvermispora (Figure 3). grown on glucose and on Ball-Milled Aspen (BMA) as carbon sources. Saline media with BMA has been used as a laboratory medium that mimics growth on wood to analyze expression of genes that are transcriptionally regulated by growth on wood. To analyze if genes regulated by conditions that mimic growth on wood, such as BMA have a different adaptation to the tRNA pool or codon bias, we used the microarray data published by Fernandez-Fueyo [8] and defined four groups of genes: group A corresponds to genes where expression was reduced at least 2 times with a p-value lower than 0.05. Group B includes all genes which showed increased expression of at least 2 times with a pvalue lower than 0.05. Group C corresponds to all genes with non-significant differences (p > 0.05) and group D contains all genes with low changes in expression (<2 fold) that are statistically significant. (p < 0.05). Our results show that group B has lower values of Nc and higher values of CAI, tAI, AAtAI, and CPB. Groups A, C, and D show non-significant differences among them. This implies that genes induced by wood preferentially use a reduced set of codons that are better adapted to the tRNA pool present in C. subvermispora (Figure 3). When we correlated the CAI, tAI, AAtAI, Nc, and CPB indexes with the ratio between expression in BMA and glucose culture medium, a statistically significant correlation was observed. Positive correlations were found with almost all indicators used (CAI, tAI, AAtAI, and CPB), the exception was Nc that showed a negative correlation. The higher correlation was identified in CAI and tAI indexes, and in genes that showed significant differences for expression in BMA saline medium compared to expression in glucose-supplemented saline medium (Table 2). This increase in correlation coefficients can be explained if growth on lignocellulose exerts pressure on codon usage When we correlated the CAI, tAI, AAtAI, Nc, and CPB indexes with the ratio between expression in BMA and glucose culture medium, a statistically significant correlation was observed. Positive correlations were found with almost all indicators used (CAI, tAI, AAtAI, and CPB), the exception was Nc that showed a negative correlation. The higher correlation was identified in CAI and tAI indexes, and in genes that showed significant differences for expression in BMA saline medium compared to expression in glucose-supplemented saline medium (Table 2). This increase in correlation coefficients can be explained if growth on lignocellulose exerts pressure on codon usage of genes involved in the metabolization of this carbon source, thereby selecting those codons that increase the translational efficiency of these genes.  Interestingly, when genes from Group B were sorted according to their codon usage adaptation values or to the tRNA pool, we found that genes coding for ribosomal proteins presented the highest CPB, tAI, and AAtAI values (Table S5). An increase in the expression of ribosomal proteins may lead to improved ribosome biogenesis, which in turn increases the overall protein biosynthetic capacity, as observed in yeast growing in rich media [37,38]. Thus, exposure of C. subvermispora to wood or lignocellulose might lead to an increase of the overall translation rate enhancing the synthesis of proteins related to lignocellulose metabolisms.

Translational Efficiency and Codon Bias in Lignocellulolytic Genes
Genome analysis of C. subvermispora indicates that this organism presents an expansion of the number of genes directly related to the mineralization and hydrolysis of lignocellulose. The genome contains 16 annotated genes of ligninolytic peroxidases (13 Manganese peroxidase, one Versatile peroxidase, one lignin peroxidase and one generic peroxidase [9]), seven genes coding for laccases, and 14 genes coding for proteins containing a cellulose-binding domain. Moreover, it also shows an expansion of the auxiliary enzymes required for lignin degradation with four genes for cellobiose dehydrogenases (CDH), five genes for ∆-12 dehydrogenases, four genes for ∆-9 dehydrogenases, five genes for Aryl-alcohol oxidase, four genes for Methanol oxidases, two genes for Aryl alcohol dehydrogenases, three genes for copper radical oxidases, and 14 for glucose methanol choline oxidoreductase [7,8]. Genes belonging to the same family exhibit differential expression, which might reflect that they serve slightly different functions in the mineralization/hydrolysis process of lignocellulose [8]. To determine whether genes of the same family present a similar bias of codon usage and of codon adaptation or adaptation to the tRNA pool, we compared the values of Nc, tAI, AAtAI, CAI, and CPB of these genes. These values were also normalized with respect to the mean and standard deviation of the values obtained for all genes encoded in C. subvermispora (Z-values) [39]. This normalization was applied to identify if and how values for lignocellulolytic genes differ with respect to these same values in genes not directly related to the lignocellulose degradation processes (Table 3).
Genes encoding manganese peroxidase generally show above-average adaptation values to the tRNA pool and codon usage (Z-value > 0). The manganese peroxidase gene with the highest level of expression (transcript ID: 50297) proved to be the most adapted to the tRNA pool, with tAI values that are more than two standard deviations from the mean tAI values of other C. subvermispora genes (Z-tAI > 2). We also observed that the most highly induced manganese peroxidase gene in BMA medium (transcript ID:129418) also shows a high adaptation value to the tRNA pool (Z-tAI = 1.619). Genes encoding laccases showed a similar trend, as the single gene that significantly changed its expression levels after growth in BMA medium (transcript ID: 130783) showed an above average value of tAI (Z-tAI = 0.461).
Among genes encoding for proteins with cellulose binding domains, which includes cellulases and other enzyme that bind and hydrolyzes polymers related to cellulose, the gene encoding for cellulase GH7-CBM1 (transcript ID: 148588) showed the highest adaptation value to the tRNA pool (Z-tAI = 1.756) and high expression level under growth with glucose as the sole carbon source. This gene also increased its expression in BMA medium. Additionally, several genes coding for cellulases GH7-CBM1 and GH5-CBM1 showed Nc values below 40, indicating a strong bias in the use of synonymous codons (Table 3). Cellulases from the GH12 family also showed high CAI values. However, in this group of genes no association between expression levels or adaptation to the tRNA pool or any other index used in this work was found.  Genes for CDH show little bias in codon usage (Nc~50). Only the CDH gene which is induced by BMA (transcript ID: 84792) shows an Nc value of 48 and slightly higher than mean (Z > 0) adaptation values to the tRNA pool ( Table 3). The ∆12-dehydrogenase genes showed a weak codon bias, with the exception of gene with transcript ID: 124050; this particular gene showed a Nc value of 38 and an above average adaptation to the tRNA pool. This gene also showed strong expression in glucose-containing medium, but its expression was not modified in BMA medium. The opposite behavior was observed in the 9-dehydrogenase genes, where the gene with the greatest induction in BMA (transcript ID: 129048) showed a lower Nc value (41.46), suggesting a strong bias in codon usage (Table 3). Regarding auxiliary genes encoding for Aryl-alcohols oxidase, Methanol oxidase, Aryl-alcohol dehydrogenase and copper radical oxidases, those induced by BMA show CAI and tAI values higher than the average (Z-values > 1). An exception were the genes encoding for Methanol oxidase (transcript ID 151964) and Aryl-alcohol dehydrogenase (transcript ID: 126785) which were not induced by BMA but shown Z-values > 1.5 (Table S6).
When genes related to mineralization and digestion of lignocellulose were arranged according to their adaptation values to codon usage or the tRNA pool, genes with greater adaptation values encoded ligninolytic peroxidases, and proteins with cellulose-binding domains. Interestingly, manganese peroxidases are more adapted to the tRNA pool, while proteins with cellulose-binding domains showed higher adaptation to codon usage.
High CAI values with tAI values lower than the optimal should indicate the use of frequent codons with low availability of tRNAs to decode them. A reduction of the speed of translation has been associated with the use of rare codons or codons decoded by tRNAs with low availability. This reduction of speed allows the proper folding of the nascent protein. Proteins containing cellulose-binding domains show a complex structure, where a proper folding of at least two domains should require some reduction in the speed of translation [40,41]. Ribosome profiling experiments could help to define if these proteins require coordination between ribosome translocation and protein folding.

Discussion
The development of massive sequencing technologies, bioinformatics sequence analysis and synthetic biology have established that synonymous mutations, far from being silent, play an important role in the fine-tuning of protein synthesis efficiency and the role of different functional forms [42,43]. The first clues or hints of the importance of synonymous mutations for translational efficiency arose from the identification of bias in the use of synonymous codons that were detected in highly expressed ribosomal proteins from E. coli, Bacillus subtilis, and S. cerevisiae [22,44,45]. In viruses, genes coding for highly required proteins, such as those of the virus capsid or nucleoprotein, show higher adaptation values to codon usage of the host than other viral genes [39,46]. The generation of synthetic genes has revealed that the amount of synthesized protein can vary over thousand fold solely by changing the composition of synonymous codons [47]. Synthetic viruses, constructed from highly virulent viruses, show reduced replication and very low host mortality when their codon usage is changed to codon combinations present in low frequency in the host genome [48]. Further, the arrangement of codons within a gene is not random and is related to the proper folding of the nascent proteins and also to the proper recycling of ribosomes [31]. For example, the low adaptation level of the FRQ gene of Neurospora crassa to codon usage is essential to maintain the rhythm of the circadian clock [49]. Continuously increasing evidence supports the fact that bias in the use of synonymous codons plays an active role in the fine-tuning of protein production. The relationship between the use of synonymous codons and translation efficiency lies in the abundance of cognate tRNAs. Highly used codons tend to have a greater number of tRNAs that recognize them [50]. In this work, we have identified 192 genes encoding tRNAs in the white-rot basidiomycete C. subvermispora, a number that is similar to the number of tRNA genes present in other fungi such as Aspergillus fumigatus (178 genes) and Schizosaccharomyces pombe (186 genes) (http://gtrnadb.ucsc.edu/). Moreover, C. subvermispora tRNA genes present the interesting feature that each type of tRNA gene (i. e. that loads the same amino acid), groups in a different clade, indicating that the current pool of tRNAs present in the C. subvermispora monokaryotic strain B has arisen from gene duplication processes, and probably not from horizontal transfer or recombination. This result is consistent with the absence of a known C. subvermispora sexual stage [12]. Interestingly, we identified an expansion of the tRNA genes coding for amino acids glycine and glutamine. The high similarity of their sequences indicates that this expansion may have occurred recently. The expansion of this set of tRNA genes could be the result of the presence of some unidentified SINE elements, which are scored as tRNAs by the tRNAscan program [51]. Another alternative is that this expansion is guided by the need to synthesize large quantities of glycine-and glutamic acid-rich proteins, requiring an enrichment of this set of tRNA genes. This strategy is similar to that used by bacteriophages to increase their rate of protein synthesis in a new host, that involves carrying tRNA genes for those codons that are most frequent in bacteriophage genes [52]. Recently it was discovered that HIV also uses a similar strategy, specifically packaging a pool of tRNAs whose codons are poorly represented in the human genome [53]. Interestingly, the addition of glutamic acid to the culture medium increases the production of cellulases by the brown rot fungus Fomitopsis sp. RCK2010 [54], and increases the production of manganese peroxidases and laccases in some white rot fungi [55]. On the other hand, the amino acid glycine is a precursor of heme biosynthesis [56], an important cofactor present in manganese peroxidases and cytochromes, genes that are highly abundant in C. subvermispora.
We also found a strong bias in codon usage, codon pairs, and adaptation to the tRNA pool of C. subvermispora genes involved in lignocellulose degradation. Similar bias has been described in lignin peroxidase from Phaenorochaete chrysosporium [57]. Since a bioinformatic identification of these genes can only be applied to those directly involved in the mineralization/digestion of lignocellulose, such as manganese peroxidase type enzymes, laccases, and cellulases, we devised a functional definition, whereby genes involved in lignocellulose degradation are those induced in the presence of BMA, a model substrate of wood (Group B). Using this functional definition, we were able to detect biases for all evaluated parameters. In turn, the lack of a strong correlation between gene expression levels and adaptation to the tRNA pool or codon bias of C. subvermispora grown in a medium containing glucose or BMA, is consistent with this postulate. A similar correlation between gene function, levels of transcriptional induction, and adaption to the translational machinery is expected when C. subvermispora grown on wood or other media that resemble its natural substrates. Future experiments need to be addressed to test this hypothesis. Studies of phylogenetic reconstruction indicate that the lignocellulose-degrading machinery may have arisen during the Paleozoic [58]. It is reasonable to pose that efficient use of this carbon source requires a metabolic adaptation that involves the whole organism, in addition to the emergence of new types of enzymes. Therefore, it is expected that genes encoding for enzymes that are part of the general metabolism of C. subvermispora also respond to exposure of lignocellulose, as observed in microarray experiments reported by Fernandez Fueyo et al. 2012 [8] were 293 genes modify significantly its expression at least 2 fold (205 from Group B and 88 from Group A). Among them, we can find genes related to the Krebs cycle, cell division, and translation (Table S5).
Expansion of cellulase, laccase, and manganese peroxidase gene families in the C. subvermispora genome permitted that each of these genes developed differential expression; however, the precise contribution of each of these genes in the digestion and mineralization of lignocellulose is not clear. In bacteria, an increase of the copy number of a particular gene is a strategy for increasing its expression. This genomic adaptation explains the expansion of genes directly related to the digestion and mineralization of lignocellulose. However, the complex expression pattern of these genes suggests that this gene expansion process is more intricate than simply a response to the need of producing more enzymes. The high bias in codon usage, adaptation of the tRNA pool, and codon diversity, in some genes directly related to the digestion and mineralization processes of lignocellulose together with the increase in expression of ribosomal proteins in BMA, suggest that C. subvermispora also uses increased translational efficiency as an additional strategy to increase the production these specific set of proteins. Thus, the increase in copy number may be linked to the generation of a diverse array of enzymes that can process a wide range of substrates. In support of this hypothesis, it has been shown that C. subvermipora manganese peroxidases present different kinetic parameters and substrate specificities [8,9].

Conclusions
Our results suggest that lignocellulose degradation by C. subvermispora has modified the genome structure of this fungus, changing the bias in codon usage, the tRNA gene pool, and codon diversity in genes that are induced in the presence of wood substrates, in order to optimize the production of these proteins. This strategy may be particularly useful in slow-growing organisms, such as C. subvermispora that cannot increase the production of enzymes by increasing cell mass. To our knowledge, this study is the first example to show metabolic adaptation to a particular ecological niche by modification of the genetic structure of an organism in favor of a selective increase of the translational efficiency of genes involved in metabolizing specific substrates that determine its adaptation to a particular environment.
Supplementary Materials: The following are available online at http://www.mdpi.com/2073-4425/11/10/1227/s1, Table S1: CPS values of pair-codon present in C. subvermispora, Table S2: tRNA encoded in the Draft genome of C. subvermispora, Table S3: tRNA repeated in C. subvermispora, Table S4: Correlation coefficient among indicators. Table S5: Values of CAI, tAI, AAtAI, Nc, and CPB in C. subvermipora genes expressed in Salt medium with BMA or Glucose. Table S6: Codon bias and translation efficiency of genes encoding for enzymes auxiliary in the lignocellulose degradation. Funding: This research was in part funded by USACH, grant number 021971GM_DAS to M.T., and supported by DI-ULagos to AG.