Deep mutational scanning studies suggest that synonymous mutations are typically silent and that most exposed, nonactive-site residues are tolerant to mutations. Here, we show that the ccdA antitoxin component of the Escherichia coli ccdAB toxin-antitoxin system is unusually sensitive to mutations when studied in the operonic context. A large fraction (∼80%) of single-codon mutations, including many synonymous mutations in the ccdA gene shows inactive phenotype, but they retain native-like binding affinity towards cognate toxin, CcdB. Therefore, the observed phenotypic effects are largely not due to alterations in protein structure/stability, consistent with a large region of CcdA being intrinsically disordered. E. coli codon preference and strength of ribosome-binding associated with translation of downstream ccdB gene are found to be major contributors of the observed ccdA mutant phenotypes. In select cases, proteomics studies reveal altered ratios of CcdA:CcdB protein levels in vivo, suggesting that the ccdA mutations likely alter relative translation efficiencies of the two genes in the operon. We extend these results by studying single-site synonymous mutations that lead to loss of function phenotypes in the relBE operon upon introduction of rarer codons. Thus, in their operonic context, genes are likely to be more sensitive to both synonymous and nonsynonymous point mutations than inferred previously.
Deep mutational scanning studies suggest that synonymous mutations are typically silent and that most exposed, nonactive-site residues are tolerant to mutations. Here, we show that the ccdA antitoxin component of the Escherichia coli ccdAB toxin-antitoxin system is unusually sensitive to mutations when studied in the operonic context. A large fraction (∼80%) of single-codon mutations, including many synonymous mutations in the ccdA gene shows inactive phenotype, but they retain native-like binding affinity towards cognate toxin, CcdB. Therefore, the observed phenotypic effects are largely not due to alterations in protein structure/stability, consistent with a large region of CcdA being intrinsically disordered. E. coli codon preference and strength of ribosome-binding associated with translation of downstream ccdB gene are found to be major contributors of the observed ccdA mutant phenotypes. In select cases, proteomics studies reveal altered ratios of CcdA:CcdB protein levels in vivo, suggesting that the ccdA mutations likely alter relative translation efficiencies of the two genes in the operon. We extend these results by studying single-site synonymous mutations that lead to loss of function phenotypes in the relBE operon upon introduction of rarer codons. Thus, in their operonic context, genes are likely to be more sensitive to both synonymous and nonsynonymous point mutations than inferred previously.
Mutations that lead to changes in the protein’s amino acid sequence are often expected to alter the protein structure, stability, and/or activity and are referred to as missense mutations. Unlike missense mutations, synonymous mutations alter only the DNA (and mRNA) sequence without affecting the amino acid sequence because of the degeneracy of the genetic code (Crick ). Synonymous mutations were classically believed to be silent in their effects on protein function. While this is often true, there have been several pioneering studies that show prominent effects of synonymous mutations on protein function (Komar ; Kimchi-Sarfaty ; Sander ). Synonymous mutations are also reported to be associated with over 50 human diseases (Sauna and Kimchi-Sarfaty 2011). Some of the mechanisms known to be involved in the altered phenotypes, observed upon introduction of multiple synonymous mutations, are due to the effects on protein expression and folding (Buhr ; Rodnina 2016), transcription (Zhao ), as well as mRNA stability due to changes in secondary structure or altered mRNA decay rates (Presnyak ; Hanson and Coller 2018). Such effects are often found to be brought about by the bias due to unequal codon usage (Sharp ; Karlin ; Quax ) and tRNA abundance (Gorochowski ) as well as effects on translation initiation (Li ), elongation or termination (Rodnina 2016).Most deep mutational scanning (DMS) studies involving parallel high-throughput investigation of phenotypic fitness in large number of single-site mutants of proteins have revealed synonymous mutations to be largely neutral (Hietpas ; Jiang ; Melamed ; Roscoe ; Wu ), with the exception of studies such as in β-lactamase where a very small number of synonymous mutants (∼2%) display <50% of wild-type (WT) activity (Firnberg ) and in Saccharomyces cerevisiae Hsp90 where one out of fifteen synonymous mutations tested shows significant loss of fitness (Fragata ). Several low throughput studies have investigated phenotypic fitness effects of small sets of synonymous and nonsynonymous (missense) single-nucleotide or single-codon mutants and attempted to elucidate molecular mechanisms responsible for the observed patterns of distribution of fitness effects (DFEs) (Sanjuán ; Carrasco ; Domingo-Calap ; Lind ; Peris ; Cuevas ; Schenk ; Bailey ; Agashe ; Kristofich ; Lebeuf-Taylor ). These studies either use adaptive evolution experiment to screen for emerging beneficial mutants or use site-directed mutagenesis (SDM) to construct a small subset of single-site mutants to study their fitness. Several single-site synonymous mutants have been identified to exhibit lower fitness than WT with reduced growth rates, while others have higher fitness with improved growth rates in various bacterial genes under selection pressure (Lind ; Schenk ; Bailey ; Agashe ; Kristofich ). A number of studies in viruses also identify small numbers of synonymous mutants to be deleterious or even lethal (Sanjuán ; Carrasco ; Domingo-Calap ; Wu ), as analyzed and detailed in an useful review on synonymous mutants in laboratory evolution experiments (Bailey ). However, most of these DFE studies show deleterious synonymous mutations to have smaller effects than nonsynonymous mutations. A study of DFE in two ribosomal proteins reveals a large fraction of mutants to be weakly deleterious and fitness effects of synonymous and nonsynonymous mutants to be comparable (Lind ). While some of the studies cite altered expression levels (higher transcript levels) or mRNA stability alterations causing the observed synonymous mutational effects, the molecular basis of altered fitness in synonymous mutants has not been elucidated in most of the above described cases. For nonsynonymous mutations, the large effects of amino acid changes on protein structure, stability and function may often overshadow and obscure the mutational effects on mRNA transcript levels, structure, stability, and translation.Mutational effects on protein function have been studied extensively in mono-cistronic proteins or in isolated, heterologous expression systems. Studies of bacterial genomes however reveal a high occurrence of operonic systems comprising of functionally related protein coding genes proximally spaced on a DNA stretch, under a shared promoter (Jacob and Monod 1961; Yanofsky 1981; Balázsi ; Price ; Zhang ). Operonic genes are typically separated by less than 20 base pairs (Moreno-Hagelsieb and Collado-Vides 2002). Operons exhibit characteristic tightly regulated mRNA expression from a common promoter, followed by translation of the gene products from the same mRNA molecule. The translation re-initiation efficiency of genes clustered in operons is often dependent on the space between genes (Chemla ). Consequently, the DNA and mRNA sequence features are expected to play an important role in gene expression and regulation in operons. Moreover, transcription and translation are known to be spatiotemporally coupled in prokaryotes (Byrne ; Stent 1964; Miller ). In tightly regulated operonic systems, downstream functions are often acutely dependent on the relative levels of gene products of the operonic genes. Thus, such operonic systems provide a sensitive and relatively unexplored readout to study phenotypic effects of synonymous mutations.Bacterial type II toxin–antitoxin (TA) systems, where the antitoxin and the toxin proteins are expressed as part of a single operon, serve as an attractive model to study co-regulation of gene expression as well as transcriptional and translational coupling (Masachis ). These TA systems usually code for an upstream labile antitoxin protein responsible for negating the toxicity of the downstream toxin protein, the latter can cause growth arrest or cell death (Goeders and Van Melderen 2014; Deter ). The Escherichia coli F-plasmid borne ccd operon encodes a TA system that comprises of a labile 72 residue CcdA antitoxin, that prevents killing of the bacterial cells by binding to the 101-residue toxin component, CcdB under the control of a single auto-regulated promoter (Miki ; Tam and Kline 1989). Both genes are co-expressed in low amounts in F-plasmid bearing E. coli cells, and their expression is autoregulated at the level of transcription (Afif ). If the cell loses F-plasmid, the labile CcdA is degraded by the ATP dependent Lon protease, releasing CcdB from the complex to act on its DNA gyrase target. Free CcdB poisons gyrase by forming a CcdB-gyrase-DNA ternary complex and induces double-stranded breaks in DNA (Bahassi ; Van Melderen 2002). This eventually leads to cell death (Van Melderen ). Therefore, mutations which disrupt CcdA antitoxin function or lower the CcdA:CcdB ratio in vivo, can lead to cell death. We therefore carried out saturation mutagenesis studies of the antitoxin ccdA, from the F plasmid encoded CcdAB TA system in E. coli to understand how point mutations, especially synonymous mutations, might affect gene function and regulation in such an operonic context.We find that the labile antitoxin CcdA in the operonic context shows considerably higher mutational sensitivity than most proteins studied to date (Gupta and Varadarajan 2018). Several synonymous point mutations in CcdA lead to loss of antitoxic function and significant reduction in cell survivals, in a manner dependent on the E. coli codon preference. Our data suggests that reduction in optimality of the codons upon mutation results in decreased levels of CcdA protein, likely by altering the translation efficiency. The reduced CcdA:CcdB protein ratio in vivo in such mutants, results in heightened toxic effects on bacterial cells. Moreover, the low CcdA:CcdB ratio also results in upregulation of ccd transcription, further amplifying the effect of mutation, resulting in cell death. Improved translation initiation of the downstream toxin CcdB, upon mutations in the ccdA gene, especially in the ribosome binding site (RBS) of CcdB, also appears to contribute to increased toxicity in cells. The ccd and likely most other TA operons are highly sensitive genetic circuits that can be used to probe effects of mutations on gene function, in vivo. We observe that a large fraction (∼80%) of single codon variants of the ccdA gene (both synonymous as well as nonsynonymous mutations) results in partial or complete loss of antitoxic function, these are referred to as inactive phenotypes.
Results
Unexpectedly High Mutational Sensitivity in the ccdA Antitoxin Gene, in its Native Operonic Context
In the ccd operonic context, mutations in the 72 residue antitoxin CcdA that impair its antitoxic function are expected to promote cell death (fig. 1). To understand the consequences of mutations in ccdA on bacterial growth, the ccd operon consisting of the promoter, ccdA and ccdB genes (fig. 1) was cloned into the pUC57 plasmid. Restriction digestion sites introduced in the construct for ease in cloning resulted in a few changes with respect to the WT sequence of the ccdAB operon, at certain positions in the noncoding region upstream of the CcdA gene (supplementary fig. S1, Supplementary Material online). A site-saturation mutagenesis (SSM) library of CcdA was created in this operonic construct using NNK primers (supplementary fig. S2, Supplementary Material online). The library comprised of both synonymous and nonsynonymous mutants of CcdA, along with ∼10% WT population arising owing to the cloning strategy (see Materials and Methods). In order to retrieve the full ccdA library in its operonic context, it was cloned and constructed in an E. coli strain (Top10 gyr) resistant to CcdB toxicity. This Top10 gyr strain has a point mutation in DNA-gyrase that prevents CcdB binding and toxicity (see Supplementary Methods online, Supplementary Material online). Following selection in an E. coli (Top 10) strain that is sensitive to the CcdB toxicity, the operonic ccdA mutant library was subjected to deep sequencing. The enrichment score (ES) of each CcdA mutant was computed from the available sequencing reads after selection (obtained from plasmid isolated from the sensitive strain) versus before selection (obtained from plasmid isolated from the resistant strain), relative to corresponding reads for WT CcdA (supplementary fig. S2 and Supplementary Methods online, Supplementary Material online). The ESseq score is thus defined as the relative ratio of fraction of mutant reads obtained after selection to that obtained before selection, normalized to WT reads. Normalization to WT reads allows facile estimation of the enrichment of the individual mutants with respect to WT. Low enrichment is indicative of reduced survival of a mutant with respect to WT. The ESseq score can therefore be used as a measure of the relative activity of mutants, which in turn is related to but not quantitatively identical to their fitness, as discussed below.
High mutational sensitivity of ccdA antitoxin gene in its operonic context in Escherichia coli (A) A schematic overview of the CcdAB TA module and the molecular mechanisms by which CcdB toxin triggers growth arrest in cells (top), how antitoxin CcdA rescues growth arrest by binding and sequestering CcdB (middle), and how mutations in CcdA can affect the antitoxic functions in the CcdAB TA system and promote CcdB mediated cell toxicity (bottom). (B) Schematic representation of ccdAB operon cloned in pUC57 vector. (C) Distribution of the mutational effects for synonymous and nonsynonymous single-site ccdA mutants, in the operonic context. ESseq values for the ccdA mutants in the library are inferred through deep sequencing. A low ESseq value is indicative of low antitoxic activity of the corresponding ccdA mutant, relative to WT (ESseq = 1, log10ESseq = 0). The median values for synonymous and nonsynonymous mutant subsets are shown in red and blue dashed lines, respectively. WT ESseq value is shown in black dashed line. The fractional frequency of mutants = number of mutants with a particular range of ESseq/total number of mutants. (D) Distribution of the inferred fitness values for synonymous and nonsynonymous single-site ccdA mutants. The WT fitness and median values for synonymous and nonsynonymous mutant subsets are shown in dashed lines. (E) Phenotypic landscape for different residues in CcdA inferred from deep sequencing data. For each residue position in CcdA, the fraction of mutants showing inactive (ESseq ≤ 0.1), marginally active (0.1 < ESseq < 0.7), active (0.7 ≤ ESseq ≤ 1.5) and hyperactive phenotypes (ESseq > 1.5) is plotted. Only positions which have data for more than 10 mutants available have been plotted here. A large number of mutants across the length of ccdA display a loss of a function phenotype.
High mutational sensitivity of ccdA antitoxin gene in its operonic context in Escherichia coli (A) A schematic overview of the CcdAB TA module and the molecular mechanisms by which CcdB toxin triggers growth arrest in cells (top), how antitoxin CcdA rescues growth arrest by binding and sequestering CcdB (middle), and how mutations in CcdA can affect the antitoxic functions in the CcdAB TA system and promote CcdB mediated cell toxicity (bottom). (B) Schematic representation of ccdAB operon cloned in pUC57 vector. (C) Distribution of the mutational effects for synonymous and nonsynonymous single-site ccdA mutants, in the operonic context. ESseq values for the ccdA mutants in the library are inferred through deep sequencing. A low ESseq value is indicative of low antitoxic activity of the corresponding ccdA mutant, relative to WT (ESseq = 1, log10ESseq = 0). The median values for synonymous and nonsynonymous mutant subsets are shown in red and blue dashed lines, respectively. WT ESseq value is shown in black dashed line. The fractional frequency of mutants = number of mutants with a particular range of ESseq/total number of mutants. (D) Distribution of the inferred fitness values for synonymous and nonsynonymous single-site ccdA mutants. The WT fitness and median values for synonymous and nonsynonymous mutant subsets are shown in dashed lines. (E) Phenotypic landscape for different residues in CcdA inferred from deep sequencing data. For each residue position in CcdA, the fraction of mutants showing inactive (ESseq ≤ 0.1), marginally active (0.1 < ESseq < 0.7), active (0.7 ≤ ESseq ≤ 1.5) and hyperactive phenotypes (ESseq > 1.5) is plotted. Only positions which have data for more than 10 mutants available have been plotted here. A large number of mutants across the length of ccdA display a loss of a function phenotype.Relative fitness of mutants of microbial genes is routinely estimated from their competitive growth rates (Lind ; Zwart ; Lebeuf-Taylor ) or from the change in ratio of sequencing reads of mutant versus WT over time (Hietpas ; Roscoe ; Fragata ; Flynn ). Laboratory evolution experiments of variant populations over several generations are unsuitable for accurate estimation of fitness effects in deleterious and lethal mutants. In an all-against-all competition over long timescales, mutants with higher growth rates tend to enrich exponentially over time and wipe out the slow growing or growth–defective variants completely from the population. This is beneficial in adaptive evolution experiments where the aim is to screen for improved variants or new function. However, in studies that aim to identify and investigate the fitness effects of both deleterious and beneficial variants, growth rate based screening can be disadvantageous. Unlike many previous studies that have investigated DFE of nonessential protein coding genes (Lind ; Bailey ; Lebeuf-Taylor ), CcdA antitoxic function is essential for survival of cells containing the CcdAB operon. We observed that a large number of ccdA mutants that showed an inactive phenotype in our plate-based screen, failed to grow in liquid culture when propagated from colonies, indicating a near complete loss-of-function. This prevents estimation of relative fitness effects of ccdA mutants from a conventional growth based assay in liquid culture. Therefore we designed a phenotypic screen of the ccdA mutant library that relies on quantifying mutants in terms of the colonies that appear after transformation. Plasmid pools used for deep sequencing were directly purified from colonies scraped from the plates.The deep sequencing was performed in three biological replicates, which showed a high correlation (Pearson correlation coefficient = 0.98 approximately; supplementary fig. S3, Supplementary Material online) among themselves. The ESseq was calculated and assigned only to mutants having a minimum of twenty reads in the resistant strain, for a given replicate. The ESs averaged over three replicates were used to classify mutants based on their activity. The WT CcdA has an ESseq value of 1. Of 2272 possible NNK single codon CcdA mutations (71 positions × 32 codons), scores could be assigned to 1,528 mutants. A drastic reduction in cellular growth in the sensitive strain (after selection samples) indicated a substantial prevalence of mutants exhibiting an inactive phenotype in the library. We observed large fractions of both synonymous and nonsynonymous single site mutants in ccdA to reduce antitoxin activity in the operonic context (fig. 1). CcdA mutants with lower enrichment are expected to have decreased antitoxic activity, thus failing to rescue CcdB mediated cell death. The low incidence of active mutations demonstrates the extreme sensitivity of the ccdAB operonic construct to mutations in ccdA. We classified mutants having ESseq values ≤ 0.1 as inactive. Surprisingly, ∼80% of all CcdA mutants have ESseq values ≤ 0.1 thus displaying an inactive phenotype (fig. 1and). A phenotypic assay was conducted for a subset of CcdA mutants by individual transformation in resistant (control) and sensitive (selection experiment) strains, and dilution plating of the single mutants. While the estimated colony forming units (CFUs) for all constructs were similar to that of WT in the resistant strain, the CFU/mL varied over a range of 106 across mutants in the sensitive strain and showed good agreement with the deep sequencing derived mutational scores (supplementary fig. S4, Supplementary material online; Table 1). Phenotypic assays of individual single mutants with ESseq values less than 0.1 confirmed that they show significant growth defects (significantly lower numbers of colonies after transformation) in the sensitive strain. Even when a few colonies are obtained for these inactive mutants, these fail to grow upon inoculation into liquid media or re-streaking on agar plates. We found equal fractions (∼80%) of both synonymous (n = 62) and nonsynonymous (n = 1,466) ccdA mutants to be highly inactive, with ESseq values ≤ 0.1. The distribution of ESseq scores for ccdA synonymous mutants (median value of 0.03) is typically higher than that for nonsynonymous mutants (median value of 0.004) (fig. 1). To compare our results with previous investigations of fitness distributions of synonymous versus nonsynonymous mutants (Bailey ), we have converted ESseq scores to fitness (w) (see Methods), using a previously described formula (Bailey ; Lebeuf-Taylor ). The median fitness values for synonymous and nonsynonymous ccdA mutants are 0.9 and 0.86, respectively (fig. 1) and are lower than the WT fitness of 1. In contrast to the previous investigations of mutational effects on microbial growth rate in liquid culture, the current study measures the growth of mutants in terms of colonies appearing on plates after transformation. As discussed above, several ccdA mutants with ESseq < 0.1 and corresponding w < 0.93 fail to grow in liquid culture. On the other hand, all previously studied microbial gene variants with reported fitness values could be successfully grown in liquid culture for several generations (Lind ; Lebeuf-Taylor ; Bailey ) implying that the identified inactive ccdA mutants are more deleterious than these previously documented cases of fitness defects. Therefore, the fitness (w) measure in the current study provides the upper limit of fitness estimations in ccdA mutants and should not be directly compared with those in the previous studies.
Table 1.
CcdA Mutant Phenotypes Inferred From Deep Sequencing and Individual Transformation.
Mutant Name[a]
Mutant Codon
Description
Phenotype[b]
ESseq
Highest Culture Dilution With Visible Growth on Plate[c]
T8
ACT
Synonymous
Hyperactive
2.482
105
V9
GTG
Synonymous
Hyperactive
16.505
106
S11
TCT
Synonymous
Inactive*
0.011
105
S11G
GGT
Nonsynonymous
Inactive
0.002
101
S13
AGT
Synonymous
Inactive
0.066
101
L16
CTG
Synonymous
Active
1.286
106
L17
CTG
Synonymous
Inactive
0.066
102
V22
GTG
Synonymous
Inactive
0.009
103
S25
TCG
Synonymous
Inactive
0.002
102
V28
GTG
Synonymous
Inactive
0.011
102
S29
AGT
Synonymous
Inactive
0.044
101
T30S
TCG
Nonsynonymous
Moderately inactive*
0.598
101
T31
ACG
Synonymous
Inactive
0.007
102
R38
AGG
Synonymous
Inactive
0.007
103
M52I
ATT
Nonsynonymous
Hyperactive
2.062
104
V55
GTT
Synonymous
Hyperactive*
1.624
103
V55
GTG
Synonymous
Moderately inactive
0.486
103
V55S
TCG
Nonsynonymous
Inactive
0.017
101
R57
AGG
Synonymous
Hyperactive
1.765
105
R57A
GCG
Nonsynonymous
Hyperactive
8.718
105
F58L
CTG
Nonsynonymous
Active
0.73
104
G63L
CTG
Nonsynonymous
Inactive
0.002
101
S64
AGT
Synonymous
Inactive
0.002
102
A66P
CCT
Nonsynonymous
Hyperactive
5.416
106
N69L
CTG
Nonsynonymous
Inactive
0
101
R70
CGG
Synonymous
Hyperactive
9.673
106
D71
GAT
Synonymous
Inactive
0.007
103
D71R
AGG
Nonsynonymous
Inactive
0
101
Mutants are labeled as WT amino acid identity_ residue position_mutated amino acid identity in case of nonsynonymous variants and as WT amino acid identity_ residue position in case of synonymous variants.
Mutant phenotype is assigned based on the deep sequencing derived ESseq values. These results correlate well with growth of individual mutants on plate, except in case of some mutants, namely S11_TCT, T30S_TCG and V55_GTT where growth on plates and deep sequencing results do not agree. These cases have been marked with an *.
Highest dilution showing growth is obtained from experiments involving individual transformation of the mutant plasmids and plating of serial dilutions to validate deep sequencing results. Spotting of serial dilutions of culture correlated well with plating results. WT ccdA construct grows till highest culture dilution of 104.
CcdA Mutant Phenotypes Inferred From Deep Sequencing and Individual Transformation.Mutants are labeled as WT amino acid identity_ residue position_mutated amino acid identity in case of nonsynonymous variants and as WT amino acid identity_ residue position in case of synonymous variants.Mutant phenotype is assigned based on the deep sequencing derived ESseq values. These results correlate well with growth of individual mutants on plate, except in case of some mutants, namely S11_TCT, T30S_TCG and V55_GTT where growth on plates and deep sequencing results do not agree. These cases have been marked with an *.Highest dilution showing growth is obtained from experiments involving individual transformation of the mutant plasmids and plating of serial dilutions to validate deep sequencing results. Spotting of serial dilutions of culture correlated well with plating results. WT ccdA construct grows till highest culture dilution of 104.
Inactive Synonymous Mutants and Codon Specific Effects in ccdAB Operonic Construct
The majority of the ccdA single-site synonymous mutants show severe growth defects (fig. 2 and Table 2). The growth defects are inferred in terms of the number of colonies formed (and thus deep sequencing reads) by each mutant on plates following transformation in sensitive Top10 E. coli cells and overnight growth at 37°C (see Materials and Methods). These synonymous inactive mutants are dispersed throughout the length of the CcdA protein, with mutations at only a few specific positions at the N- and C-terminal regions showing ESseq scores greater than 1. We could confidently assign ESseq scores to only 62 synonymous mutants having adequate representation in the library in terms of the numbers of reads. Since we used NNK (K = G/T) primers during mutagenesis and library preparation steps, we only have 50% of all 64 codons at each position represented in the library. These 62 variants bear single codon mutations at 36 positions out of the 72 residues in CcdA. Forty-nine of these 62 available synonymous mutants have ESseq scores less than 0.1 and can therefore be classified as highly inactive. To our knowledge, such a high fraction of growth-defective single-site synonymous variants has not been previously documented in any deep mutational scan or SDM study. The central stretch of CcdA (residues 20–55) was especially sensitive to synonymous mutations, with all documented synonymous mutations in this region being highly inactive. Synonymous mutations at the N terminal and C-terminal regions of CcdA resulted in both active and inactive variants. We also observed that most synonymous active mutants were hyperactive in comparison to WT with ESseq > 1.5.
Inactive phenotype of synonymous mutants and codon specificity effects for nonsynonymous mutants are nonuniformly distributed across the length of ccdA. (A) Distribution of log10ESseq scores for all 62 available single synonymous mutants as a function of residue position. The WT ESseq score of 1 (log10ESseq = 0) is depicted as a black dashed line. Loss of function mutations with low ESseq values are distributed across the length of the ccdA gene. A few mutants show WT like or higher ESseq values and are largely clustered at the terminal regions of the gene. The phenotype of each synonymous mutant of CcdA is listed in Table 2. (B) Averaged pairwise absolute differences in the ESseq scores amongst codons as a function of residue position in CcdA. The differences between ESseq values of each pair of synonymous codons encoding the same mutant amino acids are calculated. The average of the absolute differences of all mutant codon pairs for a residue position is estimated as follows: Average , where a and b are different codons encoding a missense mutation, and n is the total number of such codons at a position that have data for fitness scores. Regions close to the termini show the greatest variation in mutational effects, whereas the central region is largely populated by loss of function mutations.
Table 2.
Mutational Sensitivities of all Available Synonymous Mutant Codons at Each Residue Position of CcdA Library in Operonic Context.
Residue Position
Synonymous Mutant Codons
8
ACG
ACT
9
GTG
10
GAT
11
AGT
TCG
TCT
12
GAT
13
AGT
TCG
TCT
16
CTG
CTT
17
CTG
CTT
TTG
19
GCT
22
GTG
GTT
24
ATT
25
AGT
TCG
TCT
26
GGG
27
CTT
TTG
28
GTG
GTT
29
AGT
TCG
TCT
30
ACG
ACT
31
ACG
ACT
36
GCG
37
AGG
CGG
38
AGG
CGG
39
CTT
TTG
40
AGG
CGG
45
AAG
47
GAG
55
GTG
GTT
57
AGG
CGT
60
GAG
62
AAT
63
GGG
GGT
64
AGT
TCG
66
GCG
67
GAT
69
AAT
70
CGG
CGT
71
GAT
The codons (synonymous with respect to WT sequence of CcdA) have been colored based on the deep sequencing phenotypic activity (ESseq scores) in operonic context. Dark green depicts hyperactivity (ESseq > 1.5), while light green denotes WT like activity (0.7 ≤ ESseq ≤ 1.5). Yellow and red indicate slightly inactive mutants (0.1 < ESseq < 0.7) and highly inactive mutants (ESseq ≤ 0.1), respectively.
Inactive phenotype of synonymous mutants and codon specificity effects for nonsynonymous mutants are nonuniformly distributed across the length of ccdA. (A) Distribution of log10ESseq scores for all 62 available single synonymous mutants as a function of residue position. The WT ESseq score of 1 (log10ESseq = 0) is depicted as a black dashed line. Loss of function mutations with low ESseq values are distributed across the length of the ccdA gene. A few mutants show WT like or higher ESseq values and are largely clustered at the terminal regions of the gene. The phenotype of each synonymous mutant of CcdA is listed in Table 2. (B) Averaged pairwise absolute differences in the ESseq scores amongst codons as a function of residue position in CcdA. The differences between ESseq values of each pair of synonymous codons encoding the same mutant amino acids are calculated. The average of the absolute differences of all mutant codon pairs for a residue position is estimated as follows: Average , where a and b are different codons encoding a missense mutation, and n is the total number of such codons at a position that have data for fitness scores. Regions close to the termini show the greatest variation in mutational effects, whereas the central region is largely populated by loss of function mutations.Mutational Sensitivities of all Available Synonymous Mutant Codons at Each Residue Position of CcdA Library in Operonic Context.The codons (synonymous with respect to WT sequence of CcdA) have been colored based on the deep sequencing phenotypic activity (ESseq scores) in operonic context. Dark green depicts hyperactivity (ESseq > 1.5), while light green denotes WT like activity (0.7 ≤ ESseq ≤ 1.5). Yellow and red indicate slightly inactive mutants (0.1 < ESseq < 0.7) and highly inactive mutants (ESseq ≤ 0.1), respectively.The mutational effects observed in case of synonymous mutations led us to investigate codon specific effects amongst nonsynonymous CcdA mutants. We found several nonsynonymous mutations at multiple positions having distinct phenotypes for different codons that code for the same mutant amino acid. Since NNK codons were used for mutagenesis and preparation of the SSM library, we could only investigate phenotypes of mutated codons having G/T at the 3rd nucleotide position. Therefore, the phenotypic differences between codons could be investigated for a subset of mutated amino acid residues that are encoded by more than one codon, namely Alanine, Glycine, Leucine, Arginine, Proline, Threonine, Serine, and Valine. We also found that phenotypic differences amongst codons in nonsynonymous mutants were heightened in the N- and C- terminal positions of CcdA, while all mutations in the CcdA central stretch were inactive irrespective of the identity of the mutated codons as well as amino acids (fig. 2). These results indicate the observed effects are not due to introduction of a mutated amino acid residue, but due to changes in the DNA sequences. The CcdA library under study comprises of a very small number of active mutants, which are both synonymous and nonsynonymous and most are substantially more active in comparison to the WT construct (figs. 1 and 2). These active mutants were populated at the N- and C-terminal residues, where mutations led to the most diverse range of activities (fig. 2).
Growth Defects Caused by ccdA Mutations Are Not Due to CcdB Binding Defects
The growth defects caused by ccdA mutations are only observed in the Top 10 strain that is sensitive to CcdB toxicity and not in the Top10 gyr strain (supplementary fig. S4, Supplementary Material online), where the mutated Gyrase renders the strain resistant to CcdB mediated growth arrest. This clearly indicates that the growth defect displayed by ccdA mutations is caused by decreased inhibition of CcdB toxin in the sensitive strain. CcdA is known to rescue growth arrest by binding to CcdB with very high affinity, thus preventing CcdB from binding and poisoning DNA gyrase (Maki ; De Jonge ). Decreased inhibition could arise because of mutations affecting the folding of CcdA, mutations at residues in direct contact with CcdB, or because the total amount of expressed CcdA has been altered by mutation. In order to examine whether mutations significantly affect CcdA conformation, we characterized the ability of CcdA mutants to bind CcdB, using yeast surface display (YSD) methodology, where CcdA-CcdB binding can be decoupled from the toxic effect of CcdB on bacterial growth. CcdA was fused to S. cerevisiae Aga2p surface protein, displayed on the cell surface and its binding to purified CcdB toxin protein was measured using flow cytometry (fig. 3). Although CcdA shows high mutational sensitivity throughout its length based on the phenotypic study in the operonic construct in E. coli (fig. 3), when YSD CcdA libraries were constructed from plasmid DNA isolated from the sensitive and resistant E. coli strains by cloning in pETcon vector (fig. 3), the YSD study revealed only minor differences in CcdB binding between these two libraries (fig. 3). In addition, studies of individual mutants using YSD reveal that all CcdA synonymous mutants and WT protein bound equivalently to CcdB (fig. 3), irrespective of the corresponding distinct mutant phenotypes in the native operonic context in E. coli. Synonymous mutations have previously been suggested to affect translation rate and co-translational folding, thus altering protein activity (Krasheninnikov ; Komar ; Sander ). Such folding defects are unlikely in the case of CcdA, since its C-terminal domain, involved in CcdB binding, is intrinsically disordered and remains natively unfolded in unbound conditions. A more extensive DMS study of CcdA mutants in the YSD system indicates that mutations that affect the CcdB-binding affinity of CcdA primarily occur at residues directly involved in contact with CcdB (Chandra ). While YSD fails to replicate the in vivo conditions of E. coli and coupling of transcription-translation to protein function, the binding assay on yeast surface clearly shows that most mutations do not alter the inherent CcdB-binding function of the intrinsically disordered CcdA molecule. Further, we find no discernable difference in the distribution of phenotypes or the ESseq values between mutations at CcdB interacting and CcdB noninteracting residues in the C-terminal CcdB binding domain of CcdA, classified based on the available structure of CcdA C-terminal domain bound to dimeric CcdB (De Jonge ) (supplementary fig. S5, Supplementary Material online). These observations indicate that mutational effects on binding affinity to CcdB contribute insignificantly to the mutant phenotypes observed for ccdA in its native operonic context.
Observed growth defects in CcdA mutants are typically not caused by loss of CcdB binding activity. (A) Average log10ESseq values as a function of residue position. The error bars represent the associated standard deviation and the dashed line indicates the WT ESseq value of 1. (B) Schematic representation showing (top) a monocistronic fusion construct of ccdA gene for (bottom) yeast surface display and subsequent probing of surface expression and CcdB binding using flow cytometry (see Materials and Methods). (C) The CcdA mutant libraries derived before and after selection in E.coli strains were cloned in the YSD vector, surface displayed and probed for surface expression and CcdB binding using FACS. Double plots and histograms are shown in the upper and lower panel respectively. The resistant strain derived library shows only marginally lower CcdB binding than the sensitive strain derived library, indicating that the number of CcdA mutants having significantly lower CcdB binding activity in comparison to WT CcdA is very low (population marked with circle). The lower panel shows the CcdB-binding histograms of the two libraries studied using flow cytometry, which are not significantly different. (D) Binding of individual CcdA mutants to WT CcdB probed by yeast surface display. CcdB binding affinity was assessed by titrating the cells displaying an individual CcdA mutant with various concentrations of biotinylated CcdB, followed by binding to streptavidin conjugated AlexaFluor-633. Labeled cells were analyzed on a BD FACSARIA III. The binding affinity was calculated by fitting median fluorescence intensity to a one site binding model. All single-site synonymous mutants of CcdA (showing significant variation in phenotypes in the operonic study, see Table 1) showed comparable binding affinity to WT (see Table 3). The double mutant F65G/W44L was taken as a negative control that shows impaired CcdB binding.
Observed growth defects in CcdA mutants are typically not caused by loss of CcdB binding activity. (A) Average log10ESseq values as a function of residue position. The error bars represent the associated standard deviation and the dashed line indicates the WT ESseq value of 1. (B) Schematic representation showing (top) a monocistronic fusion construct of ccdA gene for (bottom) yeast surface display and subsequent probing of surface expression and CcdB binding using flow cytometry (see Materials and Methods). (C) The CcdA mutant libraries derived before and after selection in E.coli strains were cloned in the YSD vector, surface displayed and probed for surface expression and CcdB binding using FACS. Double plots and histograms are shown in the upper and lower panel respectively. The resistant strain derived library shows only marginally lower CcdB binding than the sensitive strain derived library, indicating that the number of CcdA mutants having significantly lower CcdB binding activity in comparison to WT CcdA is very low (population marked with circle). The lower panel shows the CcdB-binding histograms of the two libraries studied using flow cytometry, which are not significantly different. (D) Binding of individual CcdA mutants to WT CcdB probed by yeast surface display. CcdB binding affinity was assessed by titrating the cells displaying an individual CcdA mutant with various concentrations of biotinylated CcdB, followed by binding to streptavidin conjugated AlexaFluor-633. Labeled cells were analyzed on a BD FACSARIA III. The binding affinity was calculated by fitting median fluorescence intensity to a one site binding model. All single-site synonymous mutants of CcdA (showing significant variation in phenotypes in the operonic study, see Table 1) showed comparable binding affinity to WT (see Table 3). The double mutant F65G/W44L was taken as a negative control that shows impaired CcdB binding.
Table 3.
Binding Affinities of CcdA Mutants and WT Calculated Using Titration of Yeast Surface Displayed CcdA Molecules With Purified CcdB.
CcdA
Kd (nM)
Phenotype
WT
0.24 ± 0.09
Active
L17_CTG
0.4 ± 0.1
Inactive
S25_TCG
0.24 ± 0.07
Inactive
T8_ACT
0.31 ± 0.12
Active
L16_CTT
0.32 ± 0.06
Active
The dissociation constants (K) were obtained by fitting the fluorescence intensity titration traces to a single-site ligand binding model. In the operonic context in E. coli, the single codon synonymous mutants L17_CTG and S25_TCG mutants display an inactive phenotype, while T8_ACT and L16_CTT have an active phenotype based on the deep sequencing studies.
Binding Affinities of CcdA Mutants and WT Calculated Using Titration of Yeast Surface Displayed CcdA Molecules With Purified CcdB.The dissociation constants (K) were obtained by fitting the fluorescence intensity titration traces to a single-site ligand binding model. In the operonic context in E. coli, the single codon synonymous mutants L17_CTG and S25_TCG mutants display an inactive phenotype, while T8_ACT and L16_CTT have an active phenotype based on the deep sequencing studies.We also studied the distribution of mutational phenotypes (for nonsynonymous mutations) at different classes of functionally important residues in the N-terminal domain of CcdA, identified based on the available NMR structure of the DNA bound dimeric CcdA N terminal domain (Madl ). We found that mutations in the DNA binding residues in the CcdA N-terminal domain have no distinguishable difference in terms of the growth phenotype, when compared to other noninteracting exposed residues (supplementary fig. S5, Supplementary material online). Surprisingly, mutations at buried positions in the N terminal domain show significantly reduced activity and more severe growth defects, especially in the case of mutations to nonaliphatic residues, when compared to any other classes of CcdA residues (supplementary fig. S5, Supplementary Material online). This suggests that mutations that can disrupt the N-terminal buried core of CcdA can induce increased protein degradation or proteolysis, ultimately producing the observed inactive phenotypes.
E. coli Genomic Codon Usage Contributes to Growth Phenotypes of Single Codon Mutants in CcdA
The frequency of different degenerate codons (for any amino acid) in the genome of an organism varies significantly. We investigated the relationship between the phenotypes of mutants and relative codon usage (RCU) (in the E.coli genome) of the introduced codon. The distribution of RCU values amongst active and inactive CcdA mutants was significantly different (fig. 4). It has been suggested that having a higher percentage of frequent codons in a gene leads to higher amounts of expressed protein, across various kingdoms of life (Sharp ; Shields and Sharp 1987; Sharp and Devine 1989; Karlin ). Decreased levels of CcdA protein in vivo will result in the presence of unbound CcdB molecules that will poison bacterial DNA Gyrase, resulting in cell death. While there is a large fraction of inactive mutants in the ccdA single mutant library, we observe that ∼70% of ccdA mutants displaying an active phenotype harbor frequent codons (based on E. coli genome) with similar or higher codon usage frequency (CUF) relative to the WT codon. The observed high mutational sensitivity in the central region of the ccdA gene might result from the high prevalence of frequent codons in this region (supplementary fig. S6, Supplementary Material online). Thus, most mutations at these positions lead to the introduction of rarer codons, possibly leading to lower CcdA protein expression and subsequent cell death. On the other hand, the distribution of relative tRNA abundance values among active and inactive mutants was similar (fig. 4). There was also no statistically significant difference between inactive and active mutant distributions for other sequence features such as relative GC (Guanine-Cytosine) content and number(s) of base changes upon mutation (fig. 4and).
Relation between DNA sequence features and growth phenotypes of CcdA mutants. All available mutants in the CcdA library (N = 1,528) were divided into active (N = 314) and inactive (N = 1,214) classes based on an ESseq score cut-off of 0.1. Mutant distribution was plotted with respect to various sequence based features, namely RCU (A), relative tRNA abundance (B), relative GC content (C), and number of base changes upon mutation (D) calculated for mutant codon with respect to WT codon (See Methods). The distribution of sequence feature scores for inactive and active mutants have been plotted as box plots, where the internal horizontal line depicts the median, the box the interquartile range (IQR), edges of the box the first and third quartile (Q1 and Q3), the whiskers minimum (Q1-1.5IQR) and maximum (Q3 + 1.5IQR) of the distribution and the black circles the outliers. Statistical significance was assessed by a Mann-Whitney Rank Sum Test. A P-value > 0.05 indicates that there is no statistically significant difference between the distributions of the indicated parameter for active and inactive mutants, at a 95% confidence interval. Only the distribution of RCU is significantly different across the inactive and active class of mutants. The mutant distribution with respect to RCU was also plotted for ccdA synonymous mutants (N = 62) (E), nonsynonymous mutants (N = 1,466) (F) as well as the mutants showing codon specific effects (N = 550) (G), all of which showed significant differences in distributions of active and inactive classes. CcdA nonsynonymous mutants encoded by different codons (for mutated amino acids encoded by more than one codon) showing different phenotypes for the different codons are referred to as mutants showing codon specific effects. (H) Effects of mutations in the CcdB RBS site in ccdA gene on phenotype in the native operonic context. The free energy of binding to the aSD sequence on the ribosome, denoted as ΔGbind (CcdB RBS) was calculated for each of the possible single codon mutations at the CcdA residue positions 70 and 71 that constitute the RBS of the downstream CcdB gene. Calculations were done using the RNAsubopt program in the Vienna RNA package. Distribution of ΔGbind (CcdB RBS) values for inactive and active mutants at residue positions 70 and 71 of ccdA are significantly different. Statistical significance was assessed by a Mann-Whitney Rank Sum Test.
Relation between DNA sequence features and growth phenotypes of CcdA mutants. All available mutants in the CcdA library (N = 1,528) were divided into active (N = 314) and inactive (N = 1,214) classes based on an ESseq score cut-off of 0.1. Mutant distribution was plotted with respect to various sequence based features, namely RCU (A), relative tRNA abundance (B), relative GC content (C), and number of base changes upon mutation (D) calculated for mutant codon with respect to WT codon (See Methods). The distribution of sequence feature scores for inactive and active mutants have been plotted as box plots, where the internal horizontal line depicts the median, the box the interquartile range (IQR), edges of the box the first and third quartile (Q1 and Q3), the whiskers minimum (Q1-1.5IQR) and maximum (Q3 + 1.5IQR) of the distribution and the black circles the outliers. Statistical significance was assessed by a Mann-Whitney Rank Sum Test. A P-value > 0.05 indicates that there is no statistically significant difference between the distributions of the indicated parameter for active and inactive mutants, at a 95% confidence interval. Only the distribution of RCU is significantly different across the inactive and active class of mutants. The mutant distribution with respect to RCU was also plotted for ccdA synonymous mutants (N = 62) (E), nonsynonymous mutants (N = 1,466) (F) as well as the mutants showing codon specific effects (N = 550) (G), all of which showed significant differences in distributions of active and inactive classes. CcdA nonsynonymous mutants encoded by different codons (for mutated amino acids encoded by more than one codon) showing different phenotypes for the different codons are referred to as mutants showing codon specific effects. (H) Effects of mutations in the CcdB RBS site in ccdA gene on phenotype in the native operonic context. The free energy of binding to the aSD sequence on the ribosome, denoted as ΔGbind (CcdB RBS) was calculated for each of the possible single codon mutations at the CcdA residue positions 70 and 71 that constitute the RBS of the downstream CcdB gene. Calculations were done using the RNAsubopt program in the Vienna RNA package. Distribution of ΔGbind (CcdB RBS) values for inactive and active mutants at residue positions 70 and 71 of ccdA are significantly different. Statistical significance was assessed by a Mann-Whitney Rank Sum Test.We also studied the distribution of these DNA sequence features among the active and inactive variants for the subsets of synonymous mutations, nonsynonymous, as well as for the nonsynonymous substitutions which showed codon specific effects. Inactive mutations show significantly lower RCU values than active mutants in all subsets consistently (fig. 4–).Codon adaptation index (CAI) is a commonly used parameter that describes the degree of codon bias for a whole gene or sequence and helps to infer the relative adaptiveness of the sequence (Sharp and Li 1987). While it is routinely used in studies that investigate functional changes in sequences where multiple rarer or more frequent codons have been introduced, a single-site mutation (change of one codon) does not produce a large change in the CAI. CcdA WT gene has a CAI of 0.76 and CAI calculated for single codon changes are in the range of 0.7–0.78. Such small changes are not expected to show any distinguishable effects on relative fitness of mutants. We find that distribution of relative CAI values (mutant CAI/WT CAI) across inactive and active mutants fail to show any significant difference in the Mann–Whitney test (supplementary fig. S7, Supplementary Material online). On the other hand, RCU uses the relative fraction of usage of the synonymous codon encoding a particular amino acid and also takes into account the most common codon (for the WT amino acid). Estimating phenotypic effects due to altered synonymous codon usage fraction at the particular site using RCU helps to reveal that codon usage bias does indeed play a role in affecting the fitness in single-site mutants as described earlier (fig. 4and–).
Strength of the CcdB RBS Determines CcdA Mutant Phenotype
A primary sequence element that is known to be a strong determinant of the efficiency of translational coupling is the Shine Dalgarno (SD) sequence. Early studies by Das and Yanofsky in the Trp operon showed that translation initiation at any start site, located near a functional stop codon of the upstream gene may be influenced by the strength sequence and location of the SD region and its spacing from the start codon (Das and Yanofsky 1989). The RBS for the ccdB gene lies within the CcdA coding region (ccdA gene), primarily spanning CcdA residues 70–71. Therefore, mutations in these residues may affect the relative CcdA and CcdB levels by changing the strength of the RBS for CcdB.An analysis of all available codon substitutions at residues 70 and 71 of CcdA from the deep sequencing data indicated that the predicted strength of the RBS for CcdB expression was a significant contributor to the phenotype of these mutants (fig. 4). Mutants at residues 70 and 71 displaying inactive phenotype were predicted to improve the CcdB RBS strength. The molecular basis of one such mutation, D71R was examined. Introduction of the AGG arginine codon at the 71st position led to highly inactive phenotype, but no inactivity was observed for the other arginine codons. Therefore, it is unlikely that the inactive phenotype observed for the AGG codon is due to effects on protein stability or protein-protein interaction. It is known that the consensus sequence for the RBS, having the maximal strength is 5′ AGGAGG 3′ (based on the corresponding complementary sequence of 16S rRNA, 5′CCUCCU 3′). Therefore, it is likely that a change from the suboptimal RBS for CcdB present in the ccdAB operon (5′AGGGAC 3′), to a consensus RBS (as is predicted in case of D71R_AGG mutation), decreases the CcdA: CcdB ratio in the cell, there by leading to cell death. This prediction was consistent with data from quantitative proteomics (fig. 5).
Decreased CcdA:CcdB protein ratio leads to loss of function phenotype for inactive synonymous CcdA mutants and also increases ccdAB mRNA levels by feedback autoregulation of the operon. (A) Relative levels of CcdA and CcdB peptides in E. coli (Top 10 gyr strain) lysates were determined for WT as well as selected single synonymous mutants of CcdA using a quantitative proteomics approach. All synonymous inactive mutants studied have a decreased CcdA:CcdB protein ratio. The inactive nonsynonymous mutant D71R which is predicted to improve CcdB translation due to a stronger CcdB RBS, also displays a decreased CcdA: CcdB protein ratio relative to WT. (See Table 1 for phenotypic scores of mutants.) F indicates the value for F-plasmid. Mutants showing active and inactive phenotypes in the operonic study are marked in green and red, respectively. (B) Inactive synonymous mutants show higher levels of ccdAB mRNA, while active mutants have lower or WT comparable levels of ccdAB mRNA in Top10 gyr E. coli cells. The mRNA levels of ccdAB specific transcripts were determined using q-RTPCR. ccdA and ccdB regions were both independently quantitated. Mean Ct represents the threshold cycle for amplification obtained from duplicate samples from two different experiments. Fold change in ccdA and ccdB transcript levels are with respect to the WT Ct values. Error bars indicate estimated standard deviation of the measurements. The mutants are labeled as WT amino acid identity_residue position_mutated amino acid identity in case of nonsynonymous mutants and as WT amino acid identity_residue position in case of synonymous mutants. (C) Schematic model for CcdA:CcdB ratio dependent autoregulation of ccdAB operon transcription as described previously (Vandervelde ).
Decreased CcdA:CcdB protein ratio leads to loss of function phenotype for inactive synonymous CcdA mutants and also increases ccdAB mRNA levels by feedback autoregulation of the operon. (A) Relative levels of CcdA and CcdB peptides in E. coli (Top 10 gyr strain) lysates were determined for WT as well as selected single synonymous mutants of CcdA using a quantitative proteomics approach. All synonymous inactive mutants studied have a decreased CcdA:CcdB protein ratio. The inactive nonsynonymous mutant D71R which is predicted to improve CcdB translation due to a stronger CcdB RBS, also displays a decreased CcdA: CcdB protein ratio relative to WT. (See Table 1 for phenotypic scores of mutants.) F indicates the value for F-plasmid. Mutants showing active and inactive phenotypes in the operonic study are marked in green and red, respectively. (B) Inactive synonymous mutants show higher levels of ccdAB mRNA, while active mutants have lower or WT comparable levels of ccdAB mRNA in Top10 gyr E. coli cells. The mRNA levels of ccdAB specific transcripts were determined using q-RTPCR. ccdA and ccdB regions were both independently quantitated. Mean Ct represents the threshold cycle for amplification obtained from duplicate samples from two different experiments. Fold change in ccdA and ccdB transcript levels are with respect to the WT Ct values. Error bars indicate estimated standard deviation of the measurements. The mutants are labeled as WT amino acid identity_residue position_mutated amino acid identity in case of nonsynonymous mutants and as WT amino acid identity_residue position in case of synonymous mutants. (C) Schematic model for CcdA:CcdB ratio dependent autoregulation of ccdAB operon transcription as described previously (Vandervelde ).
Mutations in CcdA Lead to Altered Levels of CcdA Proteins and Therefore Modified CcdA:CcdB Ratio in Cell
Codons for inactive CcdA mutants have significantly lower RCU in comparison to active mutants, suggesting that the observed mutational phenotypes might result from alterations in translational rate of the CcdA or CcdB gene products. To further probe protein levels in mutants, we carried out a proteomics study. The extremely low levels of CcdA and CcdB proteins in vivo cannot be accurately characterized using classical methods like SDS-PAGE or Western Blotting. The ratio of CcdA: CcdB in WT was found to be slightly larger than 1 in the proteomics study (fig. 5). Interestingly, all the five inactive CcdA mutants tested show a lower CcdA:CcdB ratio (fig. 5). Introduction of rarer codons thus appears to reduce the translation of CcdA, decreasing the CcdA:CcdB ratio of proteins in the cell, causing more severe growth arrest. The inactive D71R_AGG mutant (CcdB RBS site mutation in ccdA) also exhibits a lower CcdA:CcdB ratio, likely due to increased translation of CcdB protein. Active mutant R57A and CcdAB operon in F-plasmid show comparable protein ratio relative to WT.
Increased Levels of Operonic mRNA in Inactive Mutants Confirm Decrease in CcdA:CcdB Ratio
The five inactive and one hyperactive mutants of CcdA tested in the proteomic study have relatively reduced and elevated levels of CcdA proteins relative to WT respectively. The ccdAB expression at the transcription level is autoregulated by the CcdA-CcdB protein complex in a CcdA:CcdB ratio dependent manner (Tam and Kline 1989; Afif ). When the CcdA:CcdB protein level ≥ l, elongated multimeric CcdA-CcdB complexes with alternating dimeric CcdA2 and CcdB2 are formed which bind to the ccdAB operator region with high affinity, repressing gene expression. On the other hand, when CcdA:CcdB < 1 (protein level), hetero-hexameric complexes of one CcdA dimer bound to two CcdB dimers (CcdB2-CcdA2-CcdB2) are formed that have low binding affinity to the operator region and leads to a derepressed state of operon, allowing transcription to proceed (Vandervelde ).To study the possible effects of CcdA sequence changes resulting from single codon mutations on gene expression, we quantified ccdAB operon specific RNA using RT-qPCR using both ccdA and ccdB gene specific primers. Relative amounts of ccdA and ccdB encoding transcripts were near identical in WT and were not significantly changed in case of most mutations (ratio of ccdA and ccdB specific mRNA levels ≈ 1) (fig. 5), as expected for transcription from an operon. However, significant increases in the mRNA levels for both the ccdA and ccdB amplicons were observed in case of all inactive mutants (fig. 5). This is consistent with the lowered CcdA:CcdB protein levels in the inactive mutants observed in the previous section. A lowered CcdA:CcdB protein ratio will result in de-repression of the operon, leading to upregulation in mRNA production (fig. 5). Analysis of whole genome RNA-seq and ribosome profiling data available for the E. coli genome indeed indicate that while the mRNA levels for most type II TA systems were similar for the toxin and the antitoxin consistent with our observations in ccdAB, the antitoxin protein levels can be upto two fold higher than the toxin under normal growth conditions (Deter ). This indicates that differential translation of the ccdA and ccdB genes and/or differential proteolysis of the two proteins play a vital role in maintaining the TA protein levels in vivo. Altered DNA sequences affecting the relative levels of the TA proteins can thus be detrimental to functions of TA systems.
Investigation of Possible Effects of Ribosomal Pausing on CcdA Mutant Phenotype
An analysis of bacterial genomes has suggested that SD like sequences within the coding region of genes could be potential ribosomal pause sites and is therefore under negative selection (Li ; Mohammad ). To examine if mutations in CcdA generate such potential pause sites, the difference in interaction energies with the E. coli ribosomal anti-SD sequence between single synonymous mutants and the WT sequence was calculated using the RNAsubopt program in the Vienna RNA package (Gruber ). The data indicate that the inactive and active synonymous mutations show no distinguishable difference in anti-SD interaction (supplementary fig. S8, Supplementary Material online). The extent of ribosome pausing at each codon and the associated effects on in vivo translation speed have been estimated previously (Chevance ). In the present study, we did not observe any significant differences in distribution of these relative ribosome stalling values for different codons between synonymous inactive and active mutants of CcdA (supplementary fig. S8, Supplementary Material online).
Effect of Single Codon Synonymous Mutations in Another TA Module, RelBE
To test if other TA modules are also sensitive to synonymous mutations, the RelBE system was used as a test case. Unlike the plasmidic ccdAB system, relBE is natively a chromosomal TA system, expressed as a part of the relBEF operon (Cherny ). The RelE toxin is a site-specific, ribosome dependent mRNA endonuclease. Expression of free RelE leads to bacterial growth inhibition. Mild overexpression of RelE is also reported to increase persister cell frequency (Fasani and Savageau 2013). While unrelated in sequence, the RelB protein from E. coli is similar to the CcdA antitoxin in both fold and thermodynamic properties, and inhibits the action of RelE by forming a tight complex. The RelBE complex acts as a transcriptional repressor and autoregulates its own expression (Cherny ). Homologs of this system are found in many bacterial species including pathogens like Haemophilus influenzae. An analysis of ribosome profiling data for the E. coli K12 MG1655 strain indicates that the translation efficiency for the relB gene is 3.6-fold higher than the relE gene (Deter ), although both the codon and tRNA adaptation indices are similar for both genes, similar to what we observed in the CcdAB system.To study mutational effects, relBE operon was cloned in pUC57 plasmid. We studied five single codon synonymous mutations, distributed across the length of the relB gene in a native relBE operonic context, including both three and two base substitutions. These mutations were to both more frequent and rarer codons. A control mutation at residue 77 of relB that strengthens the SD sequence of RelE was also studied. This is expected to lead to an increase in RelE levels and inhibition of cell growth (fig. 6). Using a similar screen of cell death versus growth as in the case of CcdAB system, we monitored the growth of single synonymous mutants in the WT (E. coli BW25113) as well as in the toxin deleted (E. coli BW25113 ΔrelE) strain. The two relB variants with rarer mutated codons tested here had a decreased growth in the WT strain (fig. 6and). We also found that two mutations of relB with higher codon usage with respect to WT codon display improved growth phenotype. The growth defective D77 variant involves introduction of a rarer codon besides strengthening the RBS. The relB mutant phenotypes indicate the existence of a relatively general phenomenon wherein a single-codon change in DNA sequence upon synonymous mutations can have observable phenotypic effects for TA operon systems, where small perturbations in the system can compromise cellular fitness.
Phenotypic effects of select single synonymous substitutions in relB antitoxin gene of the RelBE operon. (A) List of various synonymous mutations made in the relB gene in the native operon context along with their RCU upon mutation with respect to the WT codons. A single base substitution at the D77 residue of RelB which increases the strength of the SD sequence for downstream toxin RelE was used as a control mutant. D77 mutant also involves introduction of a rarer codon. (B and C) Growth of single synonymous mutants in E. coli BW25113, sensitive strain. Following transformation, ten-fold dilutions were plated (B) as well as spotted (C) on LB-Agar plates. The number of colonies formed in the dilution plate with ∼50–100 colonies were counted for each mutant and CFU/mL was calculated and reported. (B) The D77 synonymous mutant, expected to increase the levels of the downstream RelE toxin, displays an inactive phenotype. RelB synonymous mutants with rarer codons display decreased growth, while frequent codon bearing mutants exhibit higher growth relative to WT. The error bars indicate standard error. ‘*’ indicate significant differences from the WT (*** P value < 0.001, paired t-test).
Phenotypic effects of select single synonymous substitutions in relB antitoxin gene of the RelBE operon. (A) List of various synonymous mutations made in the relB gene in the native operon context along with their RCU upon mutation with respect to the WT codons. A single base substitution at the D77 residue of RelB which increases the strength of the SD sequence for downstream toxin RelE was used as a control mutant. D77 mutant also involves introduction of a rarer codon. (B and C) Growth of single synonymous mutants in E. coli BW25113, sensitive strain. Following transformation, ten-fold dilutions were plated (B) as well as spotted (C) on LB-Agar plates. The number of colonies formed in the dilution plate with ∼50–100 colonies were counted for each mutant and CFU/mL was calculated and reported. (B) The D77 synonymous mutant, expected to increase the levels of the downstream RelE toxin, displays an inactive phenotype. RelB synonymous mutants with rarer codons display decreased growth, while frequent codon bearing mutants exhibit higher growth relative to WT. The error bars indicate standard error. ‘*’ indicate significant differences from the WT (*** P value < 0.001, paired t-test).
Discussion
In type II TA systems, the toxin and antitoxin can form higher order TA complexes which bind to the TA promoter repressing its transcription. Conditional co-operativity (T:A protein ratio dependent co-operative interaction of TA complex with the operator DNA) is known to be a dominant mechanism in regulating the transcription of most TA systems (Cataudella ). Given the higher rate of degradation for the antitoxin due to its disordered nature, the antitoxin needs to be produced at a higher rate than the toxin, for the cells to be viable. However, mechanisms that can lead to differential levels of the toxin and the antitoxin protein from a common mRNA are not well understood. While protein levels can in principle be regulated at both the transcriptional as well as the translational level, RNA-seq data from E. coli indicate that mRNA levels for the toxin and the antitoxin are not significantly different, as is expected for genes encoded by a single polycistronic mRNA (Deter ). Ribo-seq data available for TA systems indicate that the protein synthesis rates (calculated from ribosome densities) for the antitoxin are at least two fold higher than the toxin, for six TA pairs for which sufficient data is available (Deter ). These observations indicate that there are indeed nucleotide sequence dependent features, which allow regulation at the level of translation, leading to differential synthesis rates for the two genes in the operon. Maintaining a perfect balance between the antitoxin and toxin proteins is essential for proper functioning of TA modules. Small changes in gene expression caused by altered DNA sequences, both in coding or noncoding regions can potentially be detrimental to this balance and result in observable phenotypes. In the present study, we examined mutational effects in the bacterial ccdA antitoxin gene on phenotype, in its TA operonic context, using saturation mutagenesis coupled to deep sequencing. We observe significantly high fractions (∼80%) of both synonymous and nonsynonymous mutants to show reduced antitoxic functions in the operonic context, leading to reduced survival of the bacterial host cells. An extensive review by Bailey et al (Bailey ) summarizes the DFEs of synonymous versus nonsynonymous mutations for twelve different proteins and viral genomes. Studies investigating DFEs in bacterial and yeast proteins have revealed that most mutants (>95%) retain at least 50% fitness (growth-rate based) with respect to WT (Lind ; Schenk ; Lebeuf-Taylor ). Similar results are seen in a comprehensive study of Hsp90 (Flynn ). Other studies revealed DFE of missense mutants to be bimodal. For bacterial TEM-1 (Firnberg ), ∼25% missense mutations showed loss-of-function and another ∼30% showed reduced fitness. For yeast Hsp90, ∼50% missense mutations showed 10–50% reduced fitness (Hietpas ). In case of various viral genes/genomes 20–40% mutants show lethal (loss-of-function) and 10–35% show deleterious (reduced fitness) phenotypes (Carrasco ; Domingo-Calap ; Peris ; Wu ). Notwithstanding these variations in DFE for missense mutants, all previous studies indicate most synonymous mutants to be near neutral with the exception of cases where only a small number of synonymous mutants displays a range of reduced (20–70%) fitness (Carrasco ; Lind ; Cuevas ; Wu ). For ccdA, we observe that 80% of both synonymous and nonsynonymous ccdA mutants are observed to have ESseq value ≤ 0.1 (ten-fold lower colonies obtained following transformation into Top10 E.coli strain relative to WT) in the present study. The major peak in the distribution of ESseq scores in ccdA mutants is shifted towards lower values relative to the WT. While ESseq is not identical to conventional fitness measures, it is still apparent that the fraction of synonymous (as well as nonsynonymous) mutations showing significantly reduced activity is higher than those seen in the other systems, some of which are summarized by Bailey et al. Most importantly, we observed that most inactive ccdA mutants displaying ESseq < 0.1 produced 10–1,000-fold lower number of colonies than WT ccdA when transformed in the CcdB-sensitive E.coli strain Top10 (but similar numbers of colonies when transformed in toxin resistant strain Top10 gyr to WT). In the previous DFE studies involving bacterial/yeast growth-rate based phenotypic investigation, only mutants with comparable number of colonies to ancestor/WT sequence retrieved after transformation were used for further DFE investigations (Lind ; Bailey ; Lebeuf-Taylor ). This indicates that the ccdAB operonic system under study is apparently more sensitive to mutations than those previously studied. We also found the ccdAB operonic system under study was more sensitive to nonsynonymous mutations than synonymous mutations in the ccdA gene. A reduced severity in fitness defects in case of synonymous relative to nonsynonymous mutations is also observed in few other cases (Fragata ; Bailey ), though the overall mutational sensitivity in these studies was found be much lower than observed for CcdA.The ccd operon used in the present experiments differs slightly in sequence from the one on F-plasmid (supplementary fig. S1, Supplementary Material online), resulting in a slightly reduced CcdA: CcdB protein ratio (fig. 5). This can be attributed to the presence of a suboptimal ccdA RBS in the noncoding region (supplementary fig. S1, Supplementary Material online) in our WT CcdAB construct. Small decreases in the ratio caused by mutations have dramatic phenotypic effects. This results from the very high toxicity of CcdB, the feedback regulation present in the system and the additional amplification resulting from use of a high copy number pUC plasmid to house the ccd operon. The readout here is a function of the CcdA:CcdB ratio. While it is very challenging to measure absolute levels of either protein, we are able to measure the ratio using proteomics. The data in figure 5 show, consistent with the observed phenotypes, that the CcdA:CcdB ratio is <1 for the four inactive synonymous mutants tested, whereas it is ∼1.4 in the WT operon. For reasons that are currently unclear, in the ccd operonic construct we have used, mutations that enhance the CcdA:CcdB ratio and result in higher values of ESseq (indicating beneficial variants) are rare, in contrast to other systems described in Bailey et al (Bailey ) and Flynn et al (Flynn ). The present ccd construct amplifies effects of point mutations and allows detection of changes that would be missed in other systems where small changes in expression level do not result in observable phenotypes. Strong mutant phenotypes were also found for single-site synonymous variants of the relB gene in RelBE TA operon, also cloned in a high copy number plasmid. Single codon changes in relB produced observable changes in phenotype comparable to those in the ccdAB system. The current study therefore showcases the utility of TA operons for providing general insights into the plethora of mechanisms by which point mutations can affect phenotype.The present studies reveal that most single codon substitutions, including synonymous mutations along the length of ccdA lead to a loss of function phenotype. Many nonsynonymous mutations also show codon specific phenotypes, meaning that significantly different phenotypes are observed upon introduction of different codons that encode the same mutated amino acid residue. These studies demonstrate that there can be strong selection for specific synonymous codons. This complicates interpretation and use of dN/dS ratios for studying selection pressure especially in operons. The observed dependence of the mutant phenotype on RCU of mutant with respect to WT suggests that single codon mutations in CcdA likely influence the translation efficiency of the CcdA, thus altering the relative levels of CcdA and CcdB proteins in the cell. A high frequency of frequent codons in the central region of ccdA gene suggests that the observed high mutational sensitivity at this region is because any codon introduced upon mutation is expected to be rarer relative to the WT codon. At the ccdA gene extremities there is larger variation in codon usages values (supplementary fig. S6, Supplementary Material online). Hence mutated codons at these regions may be more or less frequent relative to WT. Altered DNA sequence can directly affect ribosome assembly, ribosome stalling, rate of translation, or promote altered mRNA structure, which in turn can have multifaceted consequences on protein production. The quantitative proteomics study reveals that for loss of function mutants, the in vivo CcdA:CcdB protein level ratio is lower than that of WT, resulting in cell death. Codon usage has previously been found to be a determinant of protein translation efficiency, accuracy as well as kinetics (Ikemura 1985; Sørensen ; Quax ; Lyu ). Most prior studies do not characterize the change in protein levels associated with single synonymous codon substitutions. However, total functional protein levels have been found to increase in case of beneficial single-codon synonymous mutants in TEM-1 by upto 6–7-fold (Zwart ) and in ADAMTS13 enzyme by 1.4-fold (Hunt ), without detected effects on transcription levels and protein conformation respectively. An earlier study of an operonic system has revealed that synonymous mutations alter relative protein levels of operonic gene products upto 3-fold, causing altered growth rates (Kristofich ). In vivo translation rates of glutamic acid codons (GAA and GAG) were measured to be significantly different (and correlated to respective codon usage), despite both codons being read by same tRNA (Sørensen and Pedersen 1991). In FAE and CAT enzymes, protein levels were found to be reduced in a number of variants containing multiple synonymous substitutions across the genes (Amorós-Moya ; Agashe ). A study of introduction of a small number of rarer codons in the most highly expressed genes in E. coli also suggested that codon usage of a gene can affect the translation rate of gene with altered codons as well as of other genes due to perturbation of cellular tRNA supply (Frumkin ).The prevalence of large phenotypic variations amongst CcdA mutants exclusively at the N-terminal (first 20 residues) and C-terminal residues (56–72) might also be a result of mutational effects on translation initiation of ccdA and ccdB genes respectively. Mutations near the ccdA gene terminus, at residue positions 70 and 71 (that also constitute the RBS of CcdB) appear to affect translational efficiency of the downstream CcdB gene. Mutants at these positions which are predicted to show improved ribosome binding exhibit a lower in vivo CcdA:CcdB protein ratio and higher toxicity in E. coli. Mutations that decrease the CcdA:CcdB ratio result in de-repression of the operon, this is confirmed by the observed increased ccd mRNA expression levels in these mutants. The elevated operonic mRNA level is expected to amplify the mutational effect on differential translation efficiency of CcdA and CcdB.Earlier studies which observe a strong influence of codon rarity and codon pair bias on protein level (Quax ; Boël ) typically used over-expressed proteins. In natural contexts, many bacterial genes are organized in operons to ensure co-regulation of genes involved in similar pathways and functions. In the special case of TA operons, fine tuning of the relative amounts of antitoxin and toxin is essential and directly linked to bacterial fitness. Given the abundance and diversity of these TA systems in many bacterial species, including the pathogen Mycobacterium tuberculosis, it is reasonable to assume that the cell must have evolved multiple mechanisms to ensure their stringent spatial and temporal regulation. We find that these relative protein levels can be perturbed through changes in DNA sequences. Although effects of synonymous mutations are often attributed to co-translational folding defects and altered activity of protein, we find that most of the inactive ccdA mutants are not CcdB binding defective, based on YSD studies. We also observed that single mutations in DNA interacting residues of the CcdA N-terminal domain did not affect growth phenotype more drastically than mutations in other residues in CcdA, indicating no significant loss of DNA binding in such mutants. In contrast however, the heightened inactivity was observed in case of mutations at buried positions in the structured N-terminal domain. This suggests probable disruption of the hydrophobic core in the N-terminal domain, leading to increased proteolysis and decreased CcdA protein level. Besides the amino acid substitution related effects observed in the N-terminal domain, we also observe large phenotypic variation in the ccdA activity amongst a subset of missense mutants coded by different codons in both the N- and C- terminal regions of CcdA that cannot be explained by the codon usage bias or tRNA abundance. Interestingly, these effects occur in spans of about 40 bases which is on the order of the footprint of a single ribosome on mRNA in prokaryotes (Mohammad ). Thus, it is possible that the high diversity in phenotype for synonymous codons exclusively at the two termini of ccdA may be caused by altered translation initiation rates for CcdA and CcdB in the operon with point mutations at the N- and C-terminus positions of the ccdA gene respectively.Previous studies on prokaryotic translation also suggest mRNA structure and accessibility near the translation start site play important roles in translation initiation, a rate limiting step of translation (Studer and Joseph 2006; Gualerzi and Pon 2015; Mustoe ). Though we find no significant correlation of the phenotype with change in GC content or number of base changes upon mutations in these regions, we suspect that local mRNA structural changes upon mutations in the 5′ end of the ccdA gene affect the translation initiation of CcdA while those in the 3′ end of the ccdA gene affect the translation initiation of CcdB, leading to the large codon-specific phenotypic variation observed at both termini of the ccdA gene. While the RNA structures of the ccdAB transcript predicted by the MFold (Zuker 2003) and RNAfold (Gruber ) web servers show that several single codon mutations can affect the predicted mRNA structural architecture of ccdA mRNA and its stability, no consistent difference is observed amongst the structures for mutants displaying neutral, inactive and hyperactive phenotypes. Since mRNA stability, transcription and translation are coupled in prokaryotes, to understand mutational effects of mRNA structure on protein levels and activity, one must determine the in vivo mRNA structure of the mutants which is beyond the scope of the present study.
Materials and Methods
Cloning of WT and Saturation Mutagenesis Library of CcdA
A 986 bp region comprising of the complete ccdAB operon along with its upstream and downstream regulatory regions (based on F-plasmid sequence deposited in GenBank) was synthesized at GenScript (USA) and cloned in the EcoRV site of the pUC57 plasmid vector. A few mutations were introduced in this cloned ccdAB construct relative to F-plasmid sequence, in order to introduce restriction enzyme sites to facilitate cloning of ccdA or ccdB mutant libraries (supplementary fig. S1, Supplementary Material online). A single-site saturation mutagenesis (SSM) library of ccdA in this operonic construct was then constructed using inverse PCR methodology (Jain and Varadarajan 2014) (see Supplementary Methods online, Supplementary Material online).
Sample Preparation for Deep Sequencing and Data Processing
The plasmid library purified from the resistant Top10 gyr strain was considered as the initial unselected library. This was then transformed into the sensitive Top10 strain for selection. We used high efficiency electro-competent cells of Top10 and Top10 gyr strains, having equal efficiencies of 109 CFU/µg of pUC57 plasmid DNA. The detailed method for construction of the Top10 gyr strain is available in Supplementary Methods, Supplementary Material online. Approximately equal numbers of transformants (∼108 in number) from both strains obtained after 14 h of growth at 37°C on LB-Agar plates containing 100 μg/mL Ampicillin, were scraped into 10 mL LB media and directly used for plasmid purification using a Thermo Scientific Plasmid purification kit for both the unselected and selected libraries. These transformations and subsequent plasmid isolation procedures were performed in three biological replicates. The CcdA gene was PCR amplified from each purified plasmid sample, with primers (Sigma-Aldrich) containing condition specific, six-base long Multiplex IDentifier (MID) tags (barcodes). The resulting 270 bp long PCR products containing the full ccdA gene, were pooled, gel-band purified, and sequenced using Illumina Sequencing, on the NovaSeq6000 Platform, at Macrogen, South Korea. The overall quality of the sequencing data was assessed using the FASTQC software. Further analysis was performed using our in-house Perl scripts. Reads were separated into ‘bins’ based on their MID tags. The downstream primer sequence was used to identify forward and reverse reads in each bin. A Phred score cutoff of 20 and a minimum read length cutoff of 75 were used to filter low quality reads. Reads were converted into FASTA format and aligned with the ccdA gene sequence using the WATER program of the EMBOSS package (Rice ; Carver and Bleasby 2003). Default values were used for all the parameters except the Gap Opening Penalty, which was increased to 20. For the present CcdA DMS library, we have used a paired end sequencing platform. Forward and reverse reads of the same read pair were merged together if the overlapping regions were identical, and all mismatches between forward and reverse reads were discarded. Only those reads which cover the entire ccdA gene length were considered for further analysis. Reads with insertions, deletions and multiple mutations were omitted, and only mapped single mutants were used for further data analyses. The total number of usable reads including the WT reads in the three replicates were 201938, 478,296, and 100,205, respectively. Since NNK primers were used for construction of mutant library, 2.7% of all constructs are expected to be WT sequence. However, deep sequencing of unselected library (as well as Sanger sequencing of 15 colonies) revealed ∼10% of the constructed ccdA library to be WT. We suspect the additional 7.3% arises from carry-over of the WT plasmid used as template in mutagenesis PCR reactions that has not been completely removed by the DpnI digestion and subsequent clean up steps.We used ×200 theoretical coverage for each mutant for deep sequencing and have achieved high average read counts of each mutant in the unselected library (mean = 263 and median = 101). The average quality score for this final dataset is 36.76 and majority reads have Q score > 35. Assuming a Q score of 35, the probability of a sequencing error is 1/3,162. Assuming that all four bases are equally likely to be assigned an error at a given position, the probability of erroneously observing a specific base substitution at a specific position in a read is 1/12,649. For a given error probability and total number of reads, the probability of getting at least the observed number of reads for a specific mutant codon can be calculated assuming a binomial distribution. This binomial cumulative distribution function can be approximated by a Poisson cumulative distribution function. For 100,000 total reads and a probability P < 0.05, the minimum number of reads is 13 for a single base substitution when Q35 score is taken. This means that any single-site mutant (with 1 nucleotide change with respect to WT) with 13 or more reads has a probability of <0.05 of being observed by chance. In present analysis, we have discarded all mutants with less than 20 reads in the unselected library. This read cutoff, combined with high average read number, ensures the accuracy of the ESseq values and the correctness of phenotypic assignments from the DMS data.
Assignment of Mutant Phenotypes Based on the Deep Sequencing Data
Read numbers for mutants at all 71 positions (2–72) in CcdA were analyzed. Mutants with <20 reads in the resistant strain were not considered for analysis. The total number of reads in different conditions was calculated. Read numbers for each mutant at a given condition were normalized to the total number of reads in that condition. Mutational enrichment was defined as the ratio of the normalized reads after selection to the normalized reads before selection, for a given mutant. We further normalized the mutant scores with respect to the WT fitness score to obtain the ES.Survival score of a mutant iWhere total # of reads refers to all single mutants and WT reads identified in a sample.ES of a mutant i, ESseq is a measure of the enrichment (growth advantage) of a mutant relative to the WT construct in the sensitive strain. The ESseq scores were calculated for all mutants in the library that had >20 reads in the resistant strain, separately for each replicate. ESseq values used in the work are the average of the three replicates.For each residue in the CcdA protein, reads for all 32 possible mutant codons (NNK codons) were analyzed in the sensitive versus resistant strain, and phenotypes were assigned in terms of the ESseq scores. Positions showing codon specific mutation effects were analyzed separately.
Fitness (w) Estimation
Fitness was calculated as previously described (Bailey ; Lebeuf-Taylor ) using the formula,The number of generations was taken as 37 and calculated using a lag time of 1.4 h and doubling time of 20.2 min as reported for E. coli growth on plates (Fujikawa and Morozumi 2005).
In Vivo Activity of Individual Single-site CcdA Mutants
Selected single-site synonymous mutants of CcdA (in pUC57-ccdAB construct) were designed and synthesized from GenScript and sequence confirmed (see Table 1). 40 ng of each plasmid was individually transformed into the sensitive and the resistant strains and plated on LB-amp agar plates in ten-fold serial dilutions and grown at 37°C for 16 h. The number of CFUs for the transformants in Top10 strain obtained in different dilution plates was counted to confirm their in vivo activity. For easy readout, serially diluted transformation mix samples were also spotted using 3 μL volumes per spot to visualize differences among mutant phenotypes.
Monitoring the Relative Levels of CcdA and CcdB Using a Quantitative Proteomics Approach
The relative levels of CcdA and CcdB proteins in vivo were quantified from cell lysates using mass spectrometry based quantitative proteomics. All the experiments were carried out in two or more technical replicates except for F plasmid sample. Detailed Methodology is available in the Supplementary Methods, Supplementary Material online. The ratio of CcdA: CcdB protein levels were determined for each replicate and the mean values with standard deviations have been reported.
Estimation of ccd Specific mRNA Levels Using qRT-PCR
For RNA isolation and quantification, 3 mL of culture was grown for each of the mutants of ccdA transformed in E. coli Top10 gyr strain. Cells were grown to saturation under shaking conditions at 37°C, 180 rpm, pelleted and total RNA was extracted by the RNAsnap method (Stead ). Chromosomal and plasmidic DNA was removed by treatment with 2 units of DNase1 (New England Biolabs) for 2 h at 37°C, followed by Sodium Acetate-ethanol precipitation of RNA. RNA was quantified by nanodrop spectrophotometric estimation, and quality was assessed by agarose gel electrophoresis (samples with A260/A280 = 2, were used for further studies), prior to downstream processing. 2 µg of total RNA was taken in a sterile, RNase-free microcentrifuge tube and 0.5 µg of random hexamers (Sigma-Aldrich) was added to serve as primer for cDNA synthesis along with DEPC (Sigma-Aldrich) treated double autoclaved MilliQ water. The mixture was heated to 60°C for 10 min to melt secondary structure within the template. The tube was cooled immediately on ice to prevent secondary structure from reforming. To this mix, 200 units of SuperScript III Reverse Transcriptase (Invitrogen) were added, along with dNTPs, RT Reaction Buffer and 25 units of Ribonuclease Inhibitor following the manufacture’s protocol. The reaction mix was incubated for 60 min at 37°C. cDNA was directly used for quantitative PCR (Q-PCR) analysis. Q-PCR was set up with ccdA as well as ccdB gene specific primers (250 nM) using the cDNA template and Bio-Rad iQ SYBR Green Supermix (×1) in 20 µL total reaction volume, in a Bio-Rad iQ5 machine. The thermocycle parameters used for the Q-PCR were initial 95°C for 10 min followed by 95°C for 30 s, 56°C for 30 s and 72°C for 30 s in 40 repeat cycles. 16S rRNA was used as the internal positive control for total RNA quantification. Reactions with no reverse transcriptase as well as no template were used as a negative control.Threshold cycle of amplification in PCR reaction was analyzed and automatically provided by the Bio-Rad iQ5 Optical System Software Version 2.1. All reactions were carried out in triplicates. The negative control experiments had Ct values in the range of 29–31.
Computational Prediction of mRNA Secondary Structure
The initial 150 bp of the ccdAB transcript were submitted to mFOLD (Zuker 2003) and RNAfold (Gruber ) for prediction of secondary structure. All energy parameters were set at default values. Detailed output was obtained in the form of structure plots with reliability information, single strand frequency plots, and energy dot plots. Local mRNA secondary structural elements encompassing the region of synonymous inactive mutations were further analyzed.
Calculation of Parameters Describing Various DNA Sequence Features
The CUF of E.coli K12 strain (Nakamura ) (based on 14 coding gene sequences) was used in the study since the Top10 strain is a derivative of the K12 strain. These values of CUF in K12 was also found to correlate very well (r = 0.98) with the E. coli genome CUF table at the GenScript website (https://www.genscript.com/tools/codon-frequency-table). The CUF fraction was used to calculate the RCU for mutants, using the following formula:The CAI, CAI for WT CcdA sequence and all single-site mutants were calculated using a Python implementation of CAI and 274 highly expressed genes from E. coli K12 as reference (Lee 2018).The tRNA abundance raw data (Dong ) depicts the abundance of tRNA molecules in an E. coli cell that decode one or more codons. We assigned the tRNA abundance values to the respective codons that they decode, and calculated the relative tRNA abundance values associated with each codon using a similar normalization as described above for calculation of RCU.The ribosome stalling data (Chevance ) describes the observed averaged pausing of ribosome on each codon at the A-site for all genes in E. coli. The ribosomal pausing was estimated by introducing each codon individually into a 6X-His leader sequence of the histidine operonic system and quantified by downstream HisD-LacZ fusion protein (β-lactamase) activity. For the current work, the relative ribosome stalling for each mutation was calculated by dividing the mutant codon pause value by the WT codon pause value.The relative GC content for mutations was calculated by dividing the fractional GC content of mutant codon by that of WT codon. The number of base changes is the number of positions in a 3 bp codon where the mutant codon varies from the WT codon.To investigate if there is any significant difference in the distribution of these sequence parameters across the inactive and active classes of mutants, we chose to use the nonparametric Mann–Whitney test since our data is ordinal and does not follow a normal distribution.
Calculation of Anti-SD Interaction Energies for Synonymous Mutants in ccd mRNA
The difference in the interaction energy between single synonymous mutants in ccdA and the WT sequence, with the consensus anti-SD (aSD) sequence (5′ CCUCCUAU 3′), was calculated for a window of eight nucleotides using the RNAsubopt program from the RNA Vienna package 2.4.3 (Gruber ). The calculations were performed for ten such windows encompassing the mutant codon. The average energy difference for each mutant across these ten windows was compared for all the available synonymous active and inactive mutants. Such aSD like sequences are known ribosomal pause sites (Li ).The binding energies for different mutants at the RBS of CcdB in the ccdA gene (residue positions 70 and 71) to the aSD sequence (5′ CCUCCUAU 3′) was also calculated similarly, to look for possible ccdA mutants (at residue positions 70 and 71) that can alter the translation initiation of CcdB.
YSD of CcdA Library and Single Mutants to Probe CcdB Binding
The CcdA libraries, individually recovered both from the resistant and the sensitive strains, were amplified from the pUC57-ccdAB plasmid vector and recombined into the YSD vector, pETcon (Addgene plasmid # 41,522) between NdeI and XhoI sites, using yeast homologous recombination in S. cerevisiae EBY100 strain. Positive clones, sufficient in number to cover the library diversity, were obtained on selective SDCA plates. The transformant pools obtained were grown in liquid SDCAA media for 48 h and stored in aliquots of 109 cells per ml in SDCA media containing 20% glycerol at −70°C. The same methodology was used for cloning WT and individual synonymous mutants of CcdA into pETcon. The CcdA coding region was cloned as a C-terminal fusion to the yeast Aga2p protein with a C-terminal c-myc tag under control of the GAL promoter. The library was displayed on the yeast surface as described (Chao ). Cells were induced and grown in SGCA media containing 2% galactose at 30°C for 16 h. The surface expression of CcdA was monitored by incubation of 106 cells with anti-c-myc antibody raised in chicken (1:400 dilution) that binds to the Myc tag at the C-terminus of CcdA. Anti-chicken IgG conjugated AlexaFluor-488 (1:300) was used as the secondary antibody. CcdB binding to the surface expressed CcdA was assessed by incubating cells with a fixed concentration of 2 nM biotinylated CcdB (or varying concentrations between 0.1 pM and 200 nM for titration experiments) and binding of streptavidin conjugated AlexaFluor-633 secondary antibody (1:2,000) to biotinylated CcdB was monitored. All the primary and secondary antibodies were obtained from Invitrogen. Labeled cells were analyzed on a BD FACSARIA III. The titration curves were fit to a single ligand binding model to obtain the dissociation constants as described previously (Rathore ; Chandra ). Double plots for mean fluorescence intensities for both expression and binding for the CcdA library as well as for single synonymous mutants were analyzed to monitor differences with respect to the WT.
Monitoring the Effect of Single Synonymous Mutations in the RelBE Operon
Design of Mutants
The WT sequence for the RelBE operon from E.coli K12MG1655 strain was retrieved from the NCBI database. The entire transcription unit of 910 bp (Bech ) including another gene relF, which is expressed as a part of the same operon, was cloned in the EcoRV site of the pUC57 vector at GenScript. Five different constructs having synonymous mutations in the relB gene were synthesized with a stop codon in the relE toxin gene to ensure normal cell growth even after inactivation of the antitoxin, for ease of cloning. Three base synonymous substitutions are only possible with serine codons. Three base synonymous substitutions with varying codon preferences were made in three of the serine codons, S3, S16, and S28, all at the N-terminal half of the relB coding sequence. One representative two base substitution at position R7 was also tested. A single base substitution at the D77 residue was made to increase the strength of the SD sequence for relE and was used as a positive control mutant. To obtain an easy phenotypic readout, the introduced stop codon in the relE gene was reverted, using the inverse PCR strategy described earlier (Jain and Varadarajan 2014).
Phenotype Testing
Due to absence of a strain that is resistant to toxin action, we employed a toxin deletion strain (E. coli BW25113 ΔrelE) for the propagation of mutants. By virtue of this deletion, the WT chromosomal operon is de-repressed, leading to an excess of antitoxin in the cell, relative to the WT strain (E. coli BW25113 relBE). The inverse PCR products for individual mutants were phosphorylated and ligated, as described previously. Ligation mixes were individually transformed into the E. coli BW25113 ΔrelE strain, to recover the plasmid. Single clones were isolated, and sequence confirmed. Equal amounts of the mutant plasmids (∼30 ng) were transformed in E. coli BW25113 WT strain (having chromosomal relBE genes intact) to test their phenotype, transformation mixes were plated on LB-amp agar plates in ten-fold serial dilutions and grown at 37°C for 16 h. The number of CFUs for the transformants in E. coli BW25113 WT strain obtained in different dilution plates was counted. Serially diluted transformation mix samples were also spotted using 2 μL volumes per spot to visualize differences among mutant phenotypes.Click here for additional data file.
Authors: Natalie De Jonge; Abel Garcia-Pino; Lieven Buts; Sarah Haesaerts; Daniel Charlier; Klaus Zangger; Lode Wyns; Henri De Greve; Remy Loris Journal: Mol Cell Date: 2009-07-31 Impact factor: 17.970
Authors: Ryan Hunt; Gaya Hettiarachchi; Upendra Katneni; Nancy Hernandez; David Holcomb; Jacob Kames; Redab Alnifaidy; Brian Lin; Nobuko Hamasaki-Katagiri; Aaron Wesley; Tal Kafri; Christina Morris; Laura Bouché; Maria Panico; Tal Schiller; Juan Ibla; Haim Bar; Amra Ismail; Howard Morris; Anton Komar; Chava Kimchi-Sarfaty Journal: Int J Mol Sci Date: 2019-11-15 Impact factor: 5.923
Authors: Julia M Flynn; Ammeret Rossouw; Pamela Cote-Hammarlof; Inês Fragata; David Mavor; Carl Hollins; Claudia Bank; Daniel Na Bolon Journal: Elife Date: 2020-03-04 Impact factor: 8.140