Phosphorylation of yeast transcription factors correlates with the evolution of novel sequence and function

Abstract

Gene duplication is a significant source of novel genes and the dynamics of gene duplicate retention vs loss are poorly understood, particularly in terms of the functional and regulatory specialization of their gene products. We compiled a comprehensive data set of S. cerevisiae phosphosites to study the role of phosphorylation in yeast paralog divergence. We found that proteins coded by duplicated genes created in the Whole Genome Duplication (WGD) event and in a period prior to the WGD are significantly more phosphorylated than other duplicates or singletons. Though the amino acid sequence of each paralog of a given pair tends to diverge fairly similarly from their common ortholog in a related species, the phosphorylated amino acids tend to diverge in sequence from the ortholog at different rates. We observed that transcription factors (TFs) are disproportionately present among the set of duplicate genes and among the set of proteins that are phosphorylated. Interestingly, TFs that occur on higher levels of the transcription network hierarchy (i.e., tend to regulate other TFs) tend to be more phosphorylated than lower-level TFs. We found that TF paralog divergence in expression, binding, and sequence correlates with the abundance of phosphosites. Overall, these studies have important implications for understanding divergence of gene function and regulation in eukaryotes.

Figure 1

(A) Phylogenetic tree showing the predicted evolutionary relationship among major yeast species. Alphabetical letters (A–I) near diverging branches indicate small-scale duplication (SSD) events that are predicted to have occurred during the species divergence. Both SSD and WGD events and the resulting retained genes are as predicted by Wapinski et al. (B) Phosphosites are enriched in WGD and I category duplicates as compared to singleton genes. The number of phosphosites per gene for each duplication event (A–) and WGD was compared to the distribution of phosphosites on singleton genes. The negative log of the resulting p-values of a Wilcoxon signed-rank test is graphed for each category. We indicate the p = 0.05 level with a vertical line. WGD and I category duplicates are phosphorylated significantly above the singleton rate.

Figure 2

(A) Schematic of yeast protein sequence comparison measuring asymmetric sequence evolution of S. cerevisiae paralogs. We compared the amino acid sequences of S. cerevisiae paralogs with their orthologs in K. waltii, a descendent of yeast that did not undergo WGD. We then performed the same comparison for amino acids that are phosphorylated. p refers to the phosphosite conservation rate, a to other amino acid conservation (see Methods). We investigated the difference between |a₁ – a₂| and |p₁ – p₂| as illustrated in the figure. As depicted by the length of the arrows, in some cases one paralog diverges more rapidly than the other from their common ortholog in either amino acid sequence as a whole or only phosphosite sequence. (B) Example alignment comparison. The ASK10 and GCA1 paralogs of S. cerevisiae were aligned to their common K. waltii orthologs Kwal55.20547; a sample stretch of the alignment is depicted. The average amino acid sequence diverged fairly symmetrically in both paralogs so in this case |a₁ – a₂| is small. There are 40 phosphosites between the two paralogs. Those sites are more asymmetric in their divergence: the ASK10 phosphosite conservation rate is 0.22 and the GCA1 phosphosite conservation rate is 0.35. (C) Phosphosites sequence evolution is different among paralogs, whereas nonphosphorylated amino acids evolve similarly. The distributions of the differences between the sequence conservation of the paralogs are plotted for nonphosphorylated amino acids (black) and phosphosites (green). The means of the paralog differences in amino acids conservation and phosphosite amino acid conservation are presented in the table along with the p-value of a Wilcoxon signed-rank test comparing the two distributions.

Figure 3

(A) TFs are enriched among genes retained in WGD and several duplication events prior to WGD. We compared the fraction of TFs among singleton genes to the fraction of TFs for each duplication category. The negative log p-values of the Fisher exact test comparing the fractions are plotted. (B) TFs are more phosphorylated than the average S. cerevisiae gene. The number of phosphosites per gene is graphed for all genes (red) and TFs (blue) in the indicated duplication category.

Figure 4

TF phosphorylation is related to TF regulatory hierarchy. Depicted are the S. cerevisiae TF regulatory network hierarchy levels as defined by Yu and Gerstein. Level 1 TFs are those that directly regulate non-TF genes whereas higher-level TFs are those that tend to regulate the expression of other TFs. The corresponding average number of phosphosites per gene is plotted next to each TF level category. The top level only contains four genes so comparisons of phosphosite distributions of that TF level were not statistically significant and thus not included. The p-values of a t test comparison of level 1 phosphosites per gene distribution to levels 2 and 3 are indicated to the right of the graph. The differences in the means suggest that level 2 and 3 TFs are more phosphorylated on average than level 1 TFs. The t test p-values do not reject this hypothesis.

Figure 5

(A) Distance between gene coexpression vectors correlates with the number of phosphosites for TF paralogs but not for non-TF paralogs. We plotted the distance between the gene coexpression vectors described against the total number of phosphosites for a given paralog pair for pairs where both are TFs (blue) and pairs where neither are TFs (red) (left plot). The correlations between the coexpression distances and the number of phosphosites are presented in the leftmost table along with p-values (t test). The middle panel shows the comparison between phosphorylation number and sequence divergence between the paralogs, and the right panel is the comparison of sequence divergence vs coexpression distance. (B) TF binding divergence correlates with phosphosite number. We calculated the correlation among TF paralog pairs of their binding profiles. We investigated those paralog pairs where both paralogs belonged to the 89 TFs whose intergenic region binding profiles were calculated by Zhu et al. These correlations were plotted against the total number of phosphosites on the TF pair. The more phosphorylated a TF paralog pair, the lower the correlation of binding profiles among the pair (left panel). Similarly, the number of phosphosites correlates with sequence divergence (middle panel) and sequence divergence anticorrelates with the TF binding correlation metric (right panel).