VapC21 Toxin Contributes to Drug-Tolerance and Interacts With Non-cognate VapB32 Antitoxin in Mycobacterium tuberculosis

The prokaryotic ubiquitous Toxin-antitoxin (TA) modules encodes for a stable toxin and an unstable antitoxin. VapBC subfamily is the most abundant Type II TA system in M. tuberculosis genome. However, the exact physiological role for most of these Type II TA systems are still unknown. Here, we have comprehensively characterized the VapBC21 TA locus from M. tuberculosis. The overexpression of VapC21 inhibited mycobacterial growth in a bacteriostatic manner and as expected, growth inhibition was abrogated upon co-expression of the cognate antitoxin, VapB21. We observed that the deletion of vapC21 had no noticeable influence on the in vitro and in vivo growth of M. tuberculosis. Using co-expression and biophysical studies, we observed that in addition to VapB21, VapC21 is also able to interact with non-cognate antitoxin, VapB32. The strength of interaction varied between the cognate and non-cognate TA pairs. The overexpression of VapC21 resulted in differential expression of approximately 435 transcripts in M. tuberculosis. The transcriptional profiles obtained upon ectopic expression of VapC21 was similar to those reported in M. tuberculosis upon exposure to stress conditions such as nutrient starvation and enduring hypoxic response. Further, VapC21 overexpression also led to increased expression of WhiB7 regulon and bacterial tolerance to aminoglycosides and ethambutol. Taken together, these results indicate that a complex network of interactions exists between non-cognate TA pairs and VapC21 contributes to drug tolerance in vitro.


INTRODUCTION
Tuberculosis (TB) with an incident rate and deaths of 10 million and 1.4 million individuals, respectively, was the leading cause of death worldwide among infectious diseases in 2018 (Glaziou et al., 2018). M. tuberculosis is a highly successful pathogen because of its ability to subvert host antimicrobial pathways and persist in host tissues (Ehrt and Schnappinger, 2009). However, the exact molecular switches that enable M. tuberculosis to slow down metabolism and enter into dormant or latent state still remains unknown. Toxin-antitoxin (TA) systems, initially referred as addiction systems, are auto-regulatory operon encoding for a labile antitoxin and stable toxin (Lobato-Marquez et al., 2016;Page and Peti, 2016;Harms et al., 2018;Holden and Errington, 2018;Kang et al., 2018). The toxin-mediated growth arrest is mostly bacteriostatic, reversible and regulated by the expression levels of cognate antitoxins (Muthuramalingam et al., 2016). The characterized toxins most likely inhibit bacterial growth by targeting an essential cellular processes such as protein synthesis or cell wall synthesis or cell division or DNA replication (Lobato-Marquez et al., 2016;Harms et al., 2018). TA systems have been broadly characterized into 6 different groups based on the nature and neutralization mechanisms of antitoxins (Schuster and Bertram, 2013;Page and Peti, 2016). Among these, Type II is the most well characterized TA family where the antitoxin negates the activity of cognate toxins by forming proteinprotein complexes. The genome of M. tuberculosis encodes for a large repertoire of Type II TA systems such as mazEF, relBE, higBA, parDE, and vapBC (Pandey and Gerdes, 2005;Ramage et al., 2009;Akarsu et al., 2019;Tandon et al., 2019b). Using inducible expression systems, it has been reported that ectopic expression of the majority of these toxins inhibits bacterial growth in a bacteriostatic manner (Gupta, 2009;Ramage et al., 2009;Winther et al., 2016;Agarwal et al., 2018). Among these, MazF toxins cleaves mRNA in a sequence specific but ribosome independent manner and contribute cumulatively to the ability of M. tuberculosis to establish infection in host tissues (Zhang et al., 2003;Cook et al., 2013;Tiwari et al., 2015). In another study, it has been reported that deletion and overexpression of M. tuberculosis RelE toxins resulted in decreased and increased number of drugtolerant populations, respectively, in vitro (Singh et al., 2010).
Here, we have performed a detailed functional and biochemical characterization of VapBC21 (Rv2757c-Rv2758c) TA system from M. tuberculosis. VapBC21 TA complex is absent in M. smegmatis but present in the genome sequences of members belonging to M. tuberculosis complex (Ramage et al., 2009;Tandon et al., 2019b). Although, the crystal structure of VapC21 is available, the cellular target for VapC21 is still unknown (Jardim et al., 2016). In the present study, we demonstrate that inducible expression of VapC21 inhibited M. smegmatis growth in a bacteriostatic manner. We show that VapC21 is dispensable for M. tuberculosis survival in different stress conditions and in mice model of infection. Using growth assays, we demonstrate that in addition to VapB21, VapB32 coexpression was also able to abrogate the growth inhibition associated with VapC21 expression in M. smegmatis. We also performed surface plasmon resonance (SPR) experiments and size exclusion chromatography-multi angle light scattering (SEC-MALS) to determine the relative stabilities and binding affinities of cognate and non-cognate TA pairs. The transcription profiles obtained upon VapC21 overexpression overlapped with the profiles obtained in M. tuberculosis upon exposure to different stress conditions. VapC21 overexpression also increased the survival of M. smegmatis upon exposure to aminoglycoside and ethambutol. Taken together, we have performed a detailed biochemical and functional characterization of the VapBC21 TA system from M. tuberculosis.

Culture Condition and Generation of Mutant and Complemented Strains
The list of strains and plasmids used in the study are shown in Supplementary Table S1. E. coli strains were cultured in either Luria Bertani Broth (LB) or Terrific Broth (TB) with shaking at 200 rpm at 37 • C. The mycobacterial strains were cultured in Middlebrook 7H9 medium containing 0.2% glycerol, 0.05% Tween-80 and supplemented with 1× ADS at 200 rpm at 37 • C as previously described (Singh et al., 2013). The E. coli and mycobacterial cultures were plated on LB agar and Middlebrook 7H11 Agar supplemented with 1× OADS at 37 • C, respectively. The following antibiotics were added to the medium as and when required; 25 µg/ml kanamycin for both E. coli and mycobacteria, 50 µg/ml ampicillin and 10 µg/ml tetracycline for E. coli, and 150 µg/ml and 50 µg/ml hygromycin for E. coli and mycobacteria, respectively. Unless mentioned otherwise, all chemicals used in the study were procured from Sigma Aldrich, Merck. MIC 99 determination assays of M. tuberculosis strains against various drugs were performed as previously described (Kidwai et al., 2017).
The M. tuberculosis Erdman strain lacking the ribonuclease activity associated with VapC21 was constructed using temperature sensitive mycobacteriophages as previously described (Bardarov et al., 2002). Briefly, approximately 800 bp upstream and downstream regions flanking the vapC21 gene were cloned into cosmid, pYUB854. The recombinant cosmid, pYUB854-vapC21 was PacI digested and packaged into phagemid, phAE159 using MaxPlax TM Lambda Packaging Extract. The recombinant cosmid, phAE159-vapC21 was introduced into electrocompetent cells of M. smegmatis to generate temperature sensitive mycobacteriophages. The mutant strain of M. tuberculosis was constructed using these hightitre mycobacteriophages. The replacement of the VapC21 open reading frame with the hygromycin resistance gene in M. tuberculosis genome was confirmed by Southern blot. For construction of the complemented strain, vapC21 was PCR amplified, cloned into pVV16 under the transcriptional control of the hsp65 promoter. The resulting plasmid, pVV16-vapC21 was introduced into the electrocompetent cells of vapC21 strain.

In vitro Stress and in vivo Experiments
The growth patterns of both wild type and mutant strain were compared by measuring OD 600nm at regular intervals and CFU analysis. The colony morphology and ability of these strains to form biofilms was determined as previously described (Arora et al., 2018;Tiwari et al., 2019). The effect of deletion of vapC21 on the survival of M. tuberculosis was evaluated in the following conditions: oxidative stress (5 mM H 2 O 2 for 3 days), nitrosative stress (5 mM NaNO 2 , pH-5.2 for 3 days), 0.25% SDS for 3 days, 2.5 mg/ml lysozyme for 3 days, nutrient starvation (1× Tris buffered saline, 1× TBS with 0.05% Tween-80 for 7 and 14 days), drugs (isoniazid, rifampicin, and levofloxacin for 14 days) as previously described (Arora et al., 2018;Tiwari et al., 2019). For bacterial enumeration, at designated time points, 10-fold serial dilutions were prepared and plated on Middlebrook 7H11 medium 37 • C for 3-4 weeks.
The animal experiments were performed in accordance with guidelines recommended by the CPCSEA, Govt. of India and appropriate approvals were obtained from the animal ethics committee of Translational Health Science and Technology Institute. Prior to infection, mid-log phase cultures of M. tuberculosis were harvested, washed and single-cell suspensions were prepared. Six to eight week old mice were infected with either wild type or vapC21 or vapC21-CT strains via aerosol route, resulting into implantation of 50-100 bacilli in lung tissues. The lung bacillary loads of infected animals at 4-and 8-weeks post-infection was determined by homogenizing tissues in 2.0 ml of normal saline and 10-fold serial dilutions were plated on Middlebrook 7H11 plates at 37 • C for 3-4 weeks.

Coexpression and Growth Inhibition Studies
For growth inhibition studies in M. smegmatis, vapC21 was PCR amplified and cloned into the anyhydrotetracycline (Atc) inducible integrative vector, pTetR-Int (Ehrt et al., 2005;Agarwal et al., 2018). For coexpression studies, various antitoxins were PCR amplified and cloned into acetamide inducible vector, pLAM12 (van Kessel et al., 2008). The final recombinant constructs were confirmed by DNA sequencing. For growth assays, recombinant M. smegmatis strains were grown until earlylog phase (OD 600nm ∼ 0.2) and expression of toxin and antitoxin was induced by the addition of 50 ng/ml Atc and 0.2% acetamide, respectively. The growth of different strains was determined by either measuring absorbance (OD 600nm ) at regular intervals or by spotting assays. For spotting assays, at 12 h post-induction, an aliquot was removed and 10-fold serial dilutions were spotted on Middlebrook 7H11 plates at 37 • C for 2-3 days.

Cloning, Protein Expression, and Purification
For expression and purification of the VapBC21 complex, vapB21 (Rv2758c) and vapC21 (Rv2757c) were cloned individually into the pET-Duet vector. In order to purify toxin and antitoxins individually, the genes encoding either vapC21 (Rv2757c) or vapB21 (Rv2758c), vapB3 (Rv0550c), vapB4 (Rv0596c), vapB26 (Rv0581), or vapB32 (Rv1113) were PCR amplified and cloned into pET15b. The recombinant BL-21 (λDE3, plysE) strains were grown in TB medium until an OD 600nm ∼ 0.6 and protein expression was induced upon the addition of 1.0 mM isopropyl β-D-1-thiogalactopyranoside. The cultures were induced for 5 h at 37 • C in the case of VapBC21 complex, 16 h at 20 • C for various antitoxins and VapC21 protein. The induced cultures were harvested, resuspended in lysis buffer (10 mM HEPES, pH-8.0, 100 mM NaCl, 100 mM arginine, 10% glycerol and protease inhibitor) and lysed by sonication. The clarified lysates were incubated with 2 ml Ni-Sepharose resin with end-to-end mixing for 4 h and recombinant protein was purified as per manufacturer recommendations. The purified proteins were subjected to 15% Tricine SDS-PAGE, concentrated and stored in storage buffer (10 mM HEPES, 100 mM NaCl, 100 mM arginine, 10% glycerol, 500 mM imidazole) at −80 • C until further use. Further, a Cterminal VapB21 peptide (residues 69-108), synthesized from GeneScript was also used in the study.

NanoDSF, SEC-MALS, and SPR Studies
Thermal unfolding experiments were performed using nanoDSF of different proteins as previously described (Teale and Weber, 1957;Bruce et al., 2019;Chattopadhyay and Varadarajan, 2019;Magnusson et al., 2019). nanoDSF, measures the changes in the fluorescence of the intrinsic fluorophores tryptophan and tyrosine upon protein unfolding in a label free manner. The excitation was performed at 280 nm and emission at 330 and 350 nm and F350/330 was calculated for analysis. Briefly, 10 µM of each protein was filled in capillaries and subjected to heating from 20 • C to 90 • C, with a ramp rate of 1 • C/min. The first derivative of F350/330 was plotted as a function of temperature, to determine the apparent T m for different proteins. For SEC-MALS experiments, 100 µg of VapB antitoxins or VapC toxin or VapBC complexes were injected and separated using pre-equilibrated Superdex-200 analytical gel filtration column connected with in-line UV (SHIMADZU), MALS (mini DAWN TREOS, Wyatt Technology orp) and refractive index detectors (WATERS24614) to determine molecular weight, oligomeric state and aggregation. The equilibration of the instrument was performed in buffer containing 10 mM HEPES, pH-8.0, 100 mM NaCl, 100 mM arginine and 500 mM imidazole. The data was collected and analyzed using ASTRA TM software as previously described (Kesavardhana et al., 2017). The SPR experiments were performed using a Biacore 3000 optical biosensor at 25 • C. The dialyzed VapC21 protein was immobilized on the surface of a CM5 chip using a standard amine coupling kit as per the manufacturer's recommendations. All experiments included an activated and deactivated sensor surface (without VapC21) as negative control. The dialyzed VapB proteins in the concentration range of 25 nM to 5 µM were run across each sensor surface at a flow rate of 30 µL/min in 1× PBS (pH-7.4) running buffer containing 0.005% Tween -20. The regeneration of the sensor surface between runs was performed by washing the sensor chip twice with a solution of 4 M MgCl 2 for 10-30 s at a flow rate 30 µL/min. The binding curve data for each interaction was normalized, fitted to a simple 1:1 Langmuir interaction model using BIA EVALUATION 3.1 software and various kinetic parameters were determined.

RNA Sequencing Experiments
In the case of the overexpression strain, total RNA was isolated from M. tuberculosis harboring either pTetR-int or pTetR-Int-vapC21 grown until an OD 600nm ∼ 0.2 and the expression was induced with 50 ng/ml anhydrotetracycline for 24 h. In order to determine the effect of deletion of vapC21 on M. tuberculosis physiology, parental and mutant strains were grown until midlog phase (OD 600nm ∼ 0.8). For RNA-seq experiments, the bacterial cultures were harvested, washed, and total RNA was isolated using the TRIzol method (Singh et al., 2013). The isolated RNA was subjected to DNase I treatment, cDNA libraries were prepared from rRNA depleted samples and sequenced using a Illumina HiSeq platform at Aggrigenome labs Pvt. Ltd. as previously described (Tiwari et al., 2019). The data obtained from the mutant strain was analyzed as previously mentioned (Tiwari et al., 2019). The quality of the raw data obtained from the overexpression strain was further assessed using FastQCv 0.11.8. Subsequently, the adaptor sequences were removed from QC-passed reads using Trim Galore v0.6.4 and the cleaned reads were subjected to further analysis. The indexing of the reference genomes was performed using Bowtie2 and the trimmed pairedend reads were aligned to either M. tuberculosis Erdman or H 37 Rv genome with TopHat v2.1.1 using default parameter settings for read obtained from mutant and overexpression strain, respectively (Langmead and Salzberg, 2012;Kim et al., 2013). The reference genome was downloaded from ENSEMBL Database using the FTP download. The resulting mapped files in SAM format were converted to BAM format using SAMtools v1.9 (Li et al., 2009). The assembly of transcripts was performed using Cufflinks v2.2.1, the separately assembled transcripts were merged into a cohesive set using Cuffmerge and the differentially expressed transcripts were identified using Cuffdiff (Trapnell et al., 2010). This provided an average expression value for an individual gene in the form of fragment per kilobase of transcript per million mapped reads (FPKM) as the output. The fold change for each transcript was calculated relative to the vector control by taking a ratio of their corresponding FKPM values. The Cuffdiff output of gene counts was subsequently analyzed for differential gene expression using CummeRbund (R-package) software. The volcano plots for the identified transcripts were plotted using the basic plot function in R software.

qPCR Analysis
For qPCR analysis, 1 µg of DNase I treated mRNA was subjected to cDNA preparation using Superscript III reverse transcriptase. The synthesized cDNA was used as a template for qPCR using gene specific primers and SYBR Green mix. The data obtained was normalized to the transcript levels obtained for sigA, housekeeping gene as described previously (Singh et al., 2013).

Network Analysis
The genome-wide knowledge-based protein-protein interaction (PPI) for M. tuberculosis was used as a base network to generate a condition specific network (Mishra et al., 2017). The network is directed in most part, where direction information was obtained from databases such as KEGG and STRING (Kanehisa and Goto, 2000;von Mering et al., 2003). In the network, proteins are represented as nodes while the interactions between proteins are represented as edges. The FKPM values for all the genes were used for assigning weights to corresponding nodes and edges in order to generate a condition-specific PPI network. The node weight (NW) values were calculated for each node in the network corresponding to the fold-change (FC) values using the following equation.
where, i denotes the node in the network. The node weight values were used to calculate the edge weight (EW) values using the following equation.
where, i and j denotes nodes present in an edge.
The shortest paths between all the nodes in the weighted and directed network were computed using Dijkstra's algorithm. The algorithm computes minimum weight shortest paths, in which each path begins from a source node and ends with a sink node, through interacting proteins, choosing the least-cost edge in every step. The edge cost values were used as an input for calculating all vs. all shortest paths in each condition using Zen 1 . This was followed by computing a top-response network using a previously described network analysis pipeline (Sambarey et al., 2013;Sambaturu et al., 2016). Instead of analyzing all the paths, a subnetwork comprising the top-ranking paths were considered to constitute the top-response network. The path score was computed as a summation of the edge weights constituting the path and normalized with the path length. The resultant shortest paths were ranked based on their normalized path score using a percentile approach, with paths having a lower path score given a higher rank. Thereafter, top 0.05 percentile of highest-ranked paths were pooled and viewed as a network. The generated response networks were visualized using Cytoscape 3.7.1 (Shannon et al., 2003).

Drug Persistence Assays in VapC21
Overexpressing Strains of M. smegmatis The concentration of drugs used for the drug tolerance experiments were as follows: 0.78 µM amikacin, 0.2 µM streptomycin and 1.56 µM ethambutol. For drug tolerance experiments, the expression of VapC21 in the overexpression strain was induced by the addition of 50 ng/ml Atc. The induced cultures were diluted to an OD 600nm of 0.2 and subsequently exposed to various drugs. After 12 h post-exposure, an aliquot was removed and CFU enumeration was performed. For bacterial enumeration, 10-fold serial dilutions were prepared and plated on Middlebrook 7H11 plates. Percent survival was calculated from the obtained CFU/ml after incubation with the drug divided by the CFU/ml obtained at time zero.

Statistical Analysis
Differences between groups were determined by paired (twotailed) t test and were considered significant at a P value of < 0.05. GraphPad Prism version 8 (GraphPad Software Inc., CA, United States) was used for statistical analysis and the generation of graphs.

Overexpression of VapC21 Induces Bacteriostasis in M. smegmatis
VapBC family is the most abundant family of Type II TA systems in M. tuberculosis and its genome encodes for approximately 50 VapBC homologs (Ramage et al., 2009;Tandon et al., 2019b). VapC toxins inhibit mycobacterial growth by cleaving either mRNA, rRNA, or tRNA and their activity is neutralized by the levels of their cognate antitoxin (Winther et al., 2013;Cruz et al., 2015;Cintron et al., 2019). In our earlier studies, using the Atc based expression vector, pTetR, we have functionally characterized VapC toxins in M. bovis BCG as either inactive, moderately active or highly active proteins . The toxins were classified as highly active, where no significant increase in absorbance was observed upon induction in liquid cultures. This activity associated with the highly active toxins might be attributed to either increased expression levels or endoribonuclease activity or the essentiality of their cellular targets . In the present study, we have functionally and biochemically characterized VapC21 (Rv2757c) that belongs to the highly active VapC toxins group. Using episomal and integrative inducible vectors, we report that ectopic expression of VapC21 severely inhibited M. smegmatis growth in a bacteriostatic manner and this was abrogated in the presence of its cognate antitoxin, VapB21 (Figures 1A-D). The bacterial counts were reduced by approximately 45-, 615-, and 1500-fold in strains harboring an integrative copy of VapC21 at 9, 15, and 36 h, respectively post-Atc induction in comparison to parental strain ( Figure 1D). In agreement with the absorbance based data, the reduction in viable bacterial counts was reversed by the coexpression of their cognate antitoxin, VapB21 (Figures 1C,D). As expected, LIVE-DEAD cell viability kit revealed that the bacilli overexpressing VapC21 were as viable as the parental strain ( Figure 1E). Further, to determine the effect of VapC21 overexpression on cellular morphology of M. smegmatis, DAPIstaining of parental and the overexpression strain was performed. We noticed that the ectopic expression of VapC21 resulted in the formation of the polar head bulge with condensed nucleoid at 9 h post-Atc induction (data not shown). In contrast, DAPI stained nucleoid was uniformly distributed in vector-only control strain (data not shown). As expected, formation of the polar-head bulge was not observed when both VapC21 and VapB21 were co-expressed in M. smegmatis (data not shown).

Thermal Stability Analysis of the VapBC21 Complex in Comparison to VapB21 and VapC21
For thermal stability assays, VapBC21 complex was purified from the pET-Duet vector, with VapC21 toxin having an aminoterminal 6×-histidine tag. The individually purified VapB21 antitoxin and VapC21 toxin expressed in the pET15b vector both possessed an amino-terminus 6×-histidine tag and carboxy terminus 3×-FLAG tag. The VapBC21 complex and individual proteins were purified using a Ni-NTA purification system. The bound proteins were eluted using an imidazole gradient in the range of 100 mM -900 mM and analyzed on 15% Tricine-SDS PAGE. The purified fractions were pooled, concentrated and subjected to thermal denaturation. Apparent T m for various proteins at 10 µM concentration was determined by measuring changes in intrinsic fluorescence of tryptophan and tyrosine residues as a function of temperature. As expected, purified VapBC21 complex and VapC21 incubated with 20 µM C-terminal VapB21 peptide showed a much higher T m of 73 • C and 62 • C, respectively, as compared to the individual VapC21 toxin, which had a T m of 50 • C ( Figure 1F). In line with published reports, VapB21 was intrinsically disordered and failed to show any proper thermal transition ( Figure 1F). These observations indicate that the complexes are more thermally stable in comparison to the individual toxin and antitoxin.

VapC21 Toxin Is Not Required for Survival of M. tuberculosis Under in vitro Stress Condition
To determine the contribution of VapC21 in survival of M. tuberculosis under different stress conditions, we constructed a vapC21 mutant strain of M. tuberculosis Erdman using temperature sensitive mycobacteriophages (Bardarov et al., 2002). The replacement of VapC21 open reading frame with the hygromycin resistance gene in the genome of the mutant strain was verified by Southern blot (Supplementary Figures S1A,B and data not shown). As shown in Supplementary Figure S1B, the probe hybridized with 2.0 and 3.7 kb fragments, respectively, in lanes corresponding to PvuII digested genomic DNA from parental and mutant strain, respectively. Both parental and vapC21 mutant strain displayed identical growth patterns and no defect was observed until late exponential phase in Middlebrook 7H9 medium (Figures 2A,B). In concordance, the deletion of vapC21 did not alter the colony morphology or biofilm formation of M. tuberculosis (Supplementary Figure S1C). Previously, we have shown that the transcript levels of vapC21 remains unaltered upon exposure to various stress conditions . In agreement, we also observed that the deletion of vapC21 did not impair the ability of M. tuberculosis to survive upon exposure to either oxidative, nitrosative, nutrient starvation, acidic, lysozyme, or detergent stress (Supplementary Figures S1D-G, data not  shown). Further, we also compared the survival of parental and mutant strains upon exposure to drugs with different mechanism of action such as isoniazid, levofloxacin, rifampicin or ethambutol. We observed that both strains displayed comparable MIC 99 values for various drugs evaluated in the study (Supplementary Figure S2A). In agreement, deletion of VapC21 in M. tuberculosis genome does not affect persisters formation in vitro after exposure to isoniazid, rifampicin or levofloxacin (Supplementary Figure S2B).

VapC21 Toxin Is Dispensable for M. tuberculosis Growth in Mice Model of Infection
In order to gain further mechanistic insights into the role of VapC21 in M. tuberculosis physiology, total RNA was isolated from mid-log phase cultures of parental and mutant strains and subjected to RNA-seq analysis. Using a cut-off value of log 2 fold change of 1.0 and P-value < 0.05, we observed that approximately 131 genes were differentially expressed between the two strains. Among DEGs, the expression levels of 50 and 81 transcripts were increased or reduced, respectively ( Figure 2C and Supplementary Table S2). The observed DEGs were further characterized based on their functional category and the majority of DEGs have been annotated as either conserved hypothetical or cell wall associated proteins ( Figure 2D). More detailed analysis of the RNAseq data revealed that the expression of genes adjacent to VapC21, DNA methyl transferases, hsds.1 and hsdM were increased in the mutant strain ( Figure 2C). The transcript levels of ESAT-6 subfamily of small secreted proteins such as esxK and esxL were increased in the mutant strain in comparison to the parental strain ( Figure 2C). The increased expression of hsds.1, hsdM and esxK and esxL in the The volcano plot showing differential expressed genes between parental and mutant strain. The Y -axis and X-axis represents log 10 P-value and log 2 fold change, respectively. The red and blue dots represent DEGs with increased and decreased expression, in the mutant strain, respectively. The transcript levels of genes which remain unchanged between the parental and mutant strain are shown as gray dots. (D) The DEGs identified in panel A are further categorized according to the functional category mentioned in Mycobrowser. The Y -axis represents the number of DEGs for a functional category as mentioned on the X-axis. CW, cell wall and cell processes; LM, lipid metabolism; IMR, intermediary metabolism and respiration; CHP, conserved hypothetical protein; PE/PPE, Pro-Glu/Pro-Pro-Glu proteins; Vir, virulence, detoxification, and adaptation; IP, information pathway; RP, regulatory protein; IS, insertion sequences and phages. (E) The transcript levels of hsdS.1, hsdM, esxL, and esxK was quantified in total RNA isolated from mid-log phase cultures of wild type and vapC21 M. tuberculosis strains by qPCR using gene specific primers. The relative expression of these transcripts were obtained after normalization to levels of sigA, housekeeping gene. The data shown in this panel is mean ± SE of fold change obtained for each transcript from three independent experiments. (F) 6-8 weeks old Female Balb/c mice were infected with either wild type or vapC21 mutant or vapC21 complemented strain via the aerosol route. The lung bacillary loads were determined at 4 and 8 weeks post-infection as described in section "Materials and Methods." The data shown in this panel is mean ± SE of log 10 CFU obtained from 5 animals per group at a given time point. Despite multiple DEGs in the mutant strain, deletion of vapC21 had no effect on the virulence of M. tuberculosis. mutant strain was also validated by qPCR using gene specific primers ( Figure 2E). The transcripts of downregulated genes have been annotated as either conserved hypothetical or proteins involved in intermediary metabolism and respiration and cell wall processes. A few of the repressed transcripts such as nrdZ, Erdman_2888, Erdman_3434, Erdman_1907, Erdman_1908, Erdman_2003 and Erdman_3425 belong to the DosR regulon (Park et al., 2003;Chauhan et al., 2011;Boon and Dick, 2012; Supplementary Table S2). Also, the expression of latency associated antigens such as Erdman_0623, Erdman_1907, Erdman_1908, Erdman_2206 and otsB1 were also reduced in the mutant strain (Schuck et al., 2009;Supplementary Table S2).
We have earlier reported that strains with deletions in either vapBC3 or vapBC4 or vapBC11 or vapC22 are attenuated for growth in comparison to the wild type strain in guinea pigs and mice Deep et al., 2018). The reduced expression of latency associated genes and genes belonging to the DosR regulon suggests that VapC21 might also be important for M. tuberculosis to establish infection in host tissues. In order to investigate the role of VapC21 in M. tuberculosis pathogenesis, the growth patterns of wild type, vapC21 mutant, and vapC21 complemented strain were compared in a murine model of infection (Figure 2F). The aerosol infection of mice resulted in implantation of approximately 50 bacilli in lungs at day 1 post-infection. As shown in Figure 2F, all three strains displayed comparable growth during both acute (4 weeks) and chronic (8 weeks) stage of infection. The bacterial burdens in lung tissues of mice infected with parental, mutant and complemented strain was approximately log 10 4.5 at 8 weeks post-infection ( Figure 2F). These findings suggest that VapC21 individually does not contribute to the survival of M. tuberculosis in lung tissues.

Co-expression Studies to Identify Interactions Between VapC21 With Cognate and Non-cognate Antitoxins
Under normal physiological conditions the antitoxin levels are in excess of toxins and TA complexes function as auto-repressors (Slayden et al., 2018). Further, in addition to cognate pair, interactions also exist between non-cognate antitoxins and toxins (Zhu et al., 2010;Chen et al., 2019). Next, we performed experiments to determine whether VapC21 is able to interact with other non-cognate antitoxins. The interactions between VapC21 and cognate/non-cognate antitoxins was investigated using growth inhibition and rescue experiments in vitro. M. smegmatis mc 2 155 harboring an integrative copy of Atc inducible VapC21 was transformed individually with different pLAM12 constructs carrying an acetamide inducible copy of various antitoxins. The expression of toxin and antitoxin in early-log phase cultures of recombinant strains was induced by the addition of Atc and acetamide, respectively. The growth rescue experiments were performed by measuring OD 600nm at regular intervals. As expected, minimal growth was observed upon induction of toxin expression in M. smegmatis and this was restored upon coexpression of VapC21 along with VapB21 (Supplementary Figure S3A). In our preliminary cross-talk experiments, we observed that co-expression of non-cognate antitoxins VapB17 or VapB26 or VapB29 or VapB30 or VapB32 or VapB36 were able to rescue the growth defect associated with the overexpression of VapC21 (Supplementary Figure S3). However, in our repeat experiments, growth restoration was only observed when VapC21 was co-expressed along with VapB32 ( Figure 3A). In agreement, we also observed that VapB32 was able to abrogate the growth inhibition associated with the expression of VapC21 in our spotting assays (Figure 3B). These experiments suggest that in addition to VapB21, VapC21 might be able to interact with non-cognate antitoxin, VapB32.
inhibition studies, co-expression of VapB3, VapB4, and VapB26 had no effect where as VapB32 overexpression was able to abrogate the growth inhibition activity of VapC21. The noncognate antitoxins were purified as (His) 6 -tagged proteins using Ni-NTA based affinity chromatography. The purified fractions were pooled, concentrated and subjected to SPR analysis. The dissociation constant, K D was calculated for each set of interaction as shown in Table 1. We observed that VapB21 binds to VapC21 with a K D of ∼ 3 nM. In comparison, non-cognate antitoxins VapB4, VapB26, and VapB32 displayed reduced binding to VapC21. The K D values for binding of VapB4, VapB26 and VapB32 with VapC21 was 70, 53, and 18 nM, respectively ( Figure 3C and Table 1). In the case of VapB4, we observed ∼ 10.0-fold decrease in both association and dissociation rates in comparison to VapB21. Also, we observed ∼ 10.0 fold decrease in the association and ∼ 10.0 fold decrease in the dissociation rate in the case of VapB26 and VapB32, respectively, with VapC21 ( Figure 3C and Table 1). In agreement with co-expression data, we did not observe any binding between VapB3 and VapC21. These observations suggest that VapC21 is able to interact with non-cognate VapB antitoxins, in particular VapB32 in vitro.

SEC-MALS Analysis to Determine the Oligomeric States of Various Proteins
The oligomeric states of different purified proteins were analyzed by SEC-MALS under non-denaturing conditions. We observed two to three different peaks in the UV trace for VapBC21, VapB21-Cterminal peptide-VapC21, VapB3-VapC21, VapB4-VapC21, VapB26-VapC21 and VapB32-VapC21 (Figure 4). The molecular weights of various peaks were calculated and are shown in Table 2. The peak 1 of VapBC21 constitutes 87% of the total fraction and the molecular weight of 103 kDa corresponds to the hetero-octameric complex of T 4 A 4 ( Figure 4C). Further, the peak 2 obtained upon incubation of VapC21 with excess of C-terminal VapB21 peptide constitutes 51% of the total fraction ( Figure 4D). The peak corresponded to the hetero-tetrameric form of complex, T 2 A 2 with a molecular weight of 51 kDa ( Figure 4D). However, in the in vitro formed complexes of VapB4 and VapB26 with VapC21, the major peak constitutes 89% and 49%, with a molecular weight of 67 kDa and 51 kDa, respectively (Figures 4H,J). These peaks corresponds to hetero-tetrameric form of the complex, T 2 A 2 . As shown in Figure 4L and Table 2, the peak 2 and peak 3 obtained upon incubation of excess VapB32 with VapC21, constitutes 33% and 27% of the total fraction and represents the molecular weight of 60 kDa and 94 kDa, respectively. The peak 2 and peak 3 corresponds to heterotetrameric (T 2 A 2 ) and hetero-octameric form of the complex (T 4 A 4 ), respectively ( Figure 4L and Table 2). The obtained peak 1 in the case of VapB3, VapB4, VapB21, VapB26, VapB32 and VapC21 represents their dimeric form with an approximate molecular weight of 31, 30, 32, 19, 23, and 39.2 kDa, respectively (Figures 4A,B,E,G,I,K). However, the exact oligomeric status of complexes obtained with SEC-MALS, can differ from what is observed in vivo, owing to the small size of the individual toxins, antitoxins and inability of the column to resolve such small differences in molecular weight.

Transcriptional Response to VapC21 Overexpression in M. tuberculosis
Several studies have shown that overexpression of toxins belonging to TA systems result in transcriptional reprogramming that might enable the bacteria to adapt to different stress conditions (Singh et al., 2010;Deep et al., 2018). We next performed RNA-seq analysis to compare the transcription profiles of parental and VapC21 overexpression strain. For transcription profiling, total RNA was isolated from early-log phase cultures of various strains and subjected to RNA-seq analysis. Using a cut-off of log 2 fold change of ≥ 1.0 or ≤ -1.0 and P-value < 0.05, we observed that the overexpression of VapC21 in M. tuberculosis altered the expression of 445 genes ( Figure 5A and Supplementary Table S3). Among these differentially expressed genes, 215 and 230 transcripts were either upregulated or downregulated, respectively ( Figure 5A and Supplementary Table S3). These DEGs were further annotated according to their functional category as shown in Figure 5B. Among the DEGs with reduced expression, approximately 29% and 23% of the proteins are involved in intermediary metabolism and respiration and cell wall processes, respectively ( Figure 5B). The expression of enzymes involved in lipid metabolism of M. tuberculosis such as pks1, pks2, pks3, pks4, papA3, papA1, ppsA, ppsC, ppsD, ppsE, mmpL8, fadD9, eccCb1, fadE21, eccCa1, eccB1, fadD22, echA21, fadE18, echA1, eccD3, and eccA1 were significantly reduced in the overexpression strain (Supplementary Table S3). The expression of enzymes belonging to either ATP or NADH biosynthesis such as atpE, atpF, atpH, atpC, nuoN, nuoM, nuoL, nuoH, and nuoG were also decreased in the overexpression strain. The transcript levels of esxK, esxL and hsdS.1 which were increased in the mutant strain were observed to be decreased in the overexpression strain. In agreement with earlier reports, transcript levels of non-cognate toxins and antitoxins such as vapC1, mazE3, vapB43, vapB22, vapC15, vapB17, Rv0366c, and vapB15 were also increased in the VapC21 overexpression strain (Supplementary Table S3; Agarwal et al., 2018;Deep et al., 2018). DEGs annotated as regulatory proteins such as sigB, mce2R, whiB1, sigD, furA, whiB7, sigE, clgR were also upregulated in the overexpression strain (Supplementary Table S3). We also observed that the transcriptional response obtained upon VapC21 overexpression shared considerable overlap with the bacterial responses upon exposure to different stress conditions such as in nutrient starvation and enduring hypoxic response (Figures 5C-E). Among the identified DEGs, 90 and 41 transcripts were also differentially expressed in M. tuberculosis during nutrient starvation and enduring hypoxic response, respectively (Figures 5C-E; Betts et al., 2002;Rustad et al., 2008). The expression profile of a subset of DEGs in the VapC21 overexpression strain was confirmed by qPCR using gene specific primers (Figure 6B).

M. tuberculosis VapC21 Toxin Contribute to Persister Formation in vitro
Further, we performed network based analysis on DEGs observed between parental and VapC21 overexpression strain. The master protein-protein interaction network of M. tuberculosis, consists of 3,686 proteins (nodes) and 34,223 molecular interactions (edges) (Mishra et al., 2017). The subnetwork formed by the top-ranked paths consisted of 567 nodes (and 901 edges), of which 113 DEGs were identified in our overexpression RNAseq data ( Figure 6A). Among these, the expression of 105 and 8 transcripts were increased and decreased, respectively.
The subnetwork was analyzed to identify hub nodes (nodes having higher number of interactions with other nodes), of which whib7 (regulatory protein) and rpmG1 (ribosomal protein) were identified as the highest degree hubs. Further, other proteins belonging to functional category of regulatory protein or information pathways such as sigB, sigE, rplS, rplJ, rplM, and rpsN2 were also identified as hub-nodes in our networks ( Figure 6A). In agreement, the transcript levels of genes that are known to be regulated by WhiB7, such as eis (Rv2416c, aminoglycoside acetyltransferase) and Rv1258c (efflux pump) and Rv1988 were also increased in the overexpression strain ( Figure 6C). However, the expression levels of other transcripts for Rv0263c, Rv1257c, Rv1473, Rv2301, and Rv2302 belonging to WhiB7 regulon was comparable in the parental and overexpression strain ( Figure 6C). Previously it has been reported that the expression of whib7 is induced upon response to aminoglycosides, stress conditions and is also associated with intrinsic drug-resistance mechanisms (Morris et al., 2005;Burian et al., 2012;Reeves et al., 2013). This prompted us Degraded Antitoxin Free Antitoxin A 2 (Dimer) T 2 A 2 (Hetero-Tetramer) T 4 A 4 (Hetero-Octamer) Aggregate to investigate whether VapC21 contributes to the formation of drug-tolerant persisters in mycobacteria upon exposure to aminoglycosides. We observed that VapC21 overexpression in M. smegmatis increased the number of amikacin tolerant and streptomycin tolerant persisters by 5.0-fold and 12.0-fold, respectively (Figures 6D,E, *P < 0.05). In addition to WhiB7, the transcript levels of isoniazid inducible operon, iniBAC (Rv0340c-Rv0342c) were also increased in VapC21 overexpression strain. Studies have shown that increased levels of iniBAC operon is associated with the emergence of isoniazid and ethambutol tolerance (Alland et al., 2000;Colangeli et al., 2005). Therefore, we also determined the effect of VapC21 overexpression on emergence of ethambutol tolerant persisters in M. smegmatis. As shown in Figure 6F, overexpression of VapC21 also conferred an increase in the frequency of ethambutol tolerant persisters by 16.0-fold (*P < 0.05). Taken together, these observations suggest that overexpression of VapC21 increased the percentage of bacterial population that survived killing upon exposure to both aminoglycosides and ethambutol.

DISCUSSION
The complex molecular mechanisms of regulatory networks that coordinate mycobacterial adaptation during infection culminating into an active or latent disease and/or subsequent reactivation are poorly understood. M. tuberculosis deploys a multitude of regulatory networks to survive under different stress conditions in host tissues. The repertoire of TA systems is highly conserved in members belonging to the M. tuberculosis complex and a few of these have been biochemically and functionally characterized (Ramage et al., 2009;Tandon et al., 2019b).
In the present study, we have functionally and biochemically characterized the VapBC21 TA system from M. tuberculosis. In concordance with previous reports, we report that overexpression of VapC21 inhibited growth of M. smegmatis in a bacteriostatic manner (Tiwari et al., 2015;Agarwal et al., 2018). The observed growth inhibition upon overexpression of VapC21 could be restored upon coexpression of cognate antitoxin. Live-dead imaging revealed that M. bovis BCG overexpressing VapC21 were viable as reported in the case of overexpression of MazF, RelE, VapCs and ParE toxins. Previous studies have shown that expression of toxins belonging to TA systems results in morphological changes, such as formation of lemon shaped cells (Masuda et al., 2012;Wang et al., 2012). We also observed that overexpression of VapC21 resulted in bulge formation in M. smegmatis while no changes were observed in strains harboring the vector control. Also, in concordance with previous studies, we observed that nucleoid was more localized toward the poles in VapC21 overexpression strain. We hypothesized that morphological changes observed upon VapC21 overexpression could be attributed to the reduced expression of enzymes involved in either peptidoglycan biosynthesis or cell wall synthesis or cell division proteins.
In order to establish a successful infection, M. tuberculosis adapts to a variety of stress conditions that it encounters in the host. Several studies have shown that subsets of TA systems are differentially expressed upon exposure of M. tuberculosis to stress conditions such as oxidative or low oxygen or nutrient limiting conditions (Ramage et al., 2009;Agarwal et al., 2018). These observations indicates that in addition to their postulated role in genome stability, TA systems might work synergistically to enable bacteria to adapt to different stress conditions and persist in host tissues. In order to understand the role of VapC21 in M. tuberculosis physiology, a vapC21 mutant strain was generated using temperature sensitive mycobacteriophages. In agreement, with previously published gene essentiality data, we demonstrate that VapC21 is dispensable for survival of M. tuberculosis in vitro in different conditions. TA systems have been implicated to regulate pathogenesis in Uropathogenic E. coli, S. typhimurium, H. influenzae and M. tuberculosis (Norton and Mulvey, 2012;Ren et al., 2012;De la Cruz et al., 2013;Tiwari et al., 2015;Agarwal et al., 2018;Deep et al., 2018). Previously, we have shown that both MazF and VapC ribonucleases contribute to the ability of M. tuberculosis to establish infection in host tissues. These studies suggested that TA systems such as MazF3, MazF6, MazF9, VapBC3, VapBC4, VapBC11, and VapC22 enable M. tuberculosis to adapt to changes associated with the onset of host adaptive immunity (Tiwari et al., 2015;Agarwal et al., 2018Agarwal et al., , 2020Deep et al., 2018). Despite the reduced expression of genes belonging to the dormancy regulon or latency associated antigens in the mutant strain, we did not observe any significant The overexpression of VapC21 resulted in differential expression profiles of 445 transcripts. The transcripts whose expression were either increased (red dots) or decreased (blue dots) or remained unchanged (gray dots) are shown in the Volcano plot. The Y -axis and X-axis represents log 10 P-value and log 2 fold change, respectively. (B) The DEGs identified in panel A have been categorized as per the functional category mentioned in Mycobrowser. The different functional categories have been described in legend to Figure 2D. The Y -axis represents the number of DEGs for a given functional category and these are mentioned on X-axis.  The transcript levels of DEGs (B) and WhiB7 regulon genes (C) were quantified in early-log phase cultures of parental and VapC21 overexpression M. tuberculosis strain by qPCR using gene specific primers. The relative expression of these transcripts was obtained after normalization to levels of sigA, housekeeping gene. The data shown in this panel is mean ± SE of fold change obtained for each transcript from two independent experiments. (D-F) The expression of toxin was induced in M. smegmatis by the addition of 50 ng/ml Atc for 12 h. The cultures were diluted and exposed to medium containing either amikacin (D) or streptomycin (E) or ethambutol (F) for 12 h. The bacterial enumeration and percentage survival was calculated as described in section "Materials and Methods." The results shown are mean ± SE of percentage survival obtained from three or four independent experiments. Statistically significant differences were obtained for the indicated groups (Paired two-tailed t-test, *P < 0.05). differences in the lung bacillary loads in mice infected with various strains until 8 weeks post-infection. Taken together, these findings indicate that similar to VapC28, RelE1, RelE2, RelE3 and DarTG, VapC21 is also dispensable for M. tuberculosis growth in liquid cultures, stress adaptation and in host tissues (Singh et al., 2010;Agarwal et al., 2018;Zaveri et al., 2020).
Several studies have shown that TA systems constitute an interactome instead of the one-to-one interaction model (Zhu et al., 2010;Chen et al., 2019). For example, Zhu et al. (2010) showed that Rv1991c in addition to Rv1991a also interacts with Rv2801a, Rv0599c, and Rv2595. Also, recently, it has been reported that Rv2063a interacts with both Rv2063 and Rv2801a (Chen et al., 2019). However, the strength of interactions between members of a cognate TA pair was relatively stronger in comparison to those observed in non-cognate TA pairs (Chen et al., 2019). In the present study, we observed that in addition to VapB21, co-expression of VapB32 was able to abrogate growth inhibition associated with overexpression of VapC21 in M. smegmatis. SPR studies revealed the strongest binding of VapC21 was observed with VapB21 followed by non-cognate antitoxins, VapB32, VapB26, VapB4 and VapB3. In agreement, with the reduced binding affinity of VapB26, VapB4 and VapB3 with VapC21 we did not observe any growth restoration upon co-expression of these antitoxins. SEC-MALS analysis revealed that complexes formed between VapC21 and antitoxins VapB21 (cognate) and VapB32 (non-cognate) were hetero-octamers, whereas hetero-tetrameric complexes were formed between VapC21 and non-cognate antitoxins VapB4 and VapB26. Based on these observations, we speculate that TA systems indeed constitute a regulatory network that enables M. tuberculosis to restore growth from stressful conditions. Since these noncognate TA pair interactions are relatively weaker in strength in comparison to the cognate TA pair, we hypothesize that under cellular conditions with reduced levels of cognate antitoxin, these heterologous interactions could facilitate fine tuning the activity of non-cognate toxins in vivo.
The analysis of the transcriptional response to the ectopic expression of VapC21 in M. tuberculosis revealed that the presence of free toxin results in transcriptional reprogramming that overlapped substantially with responses observed upon exposure of bacteria to different stress conditions. The observed pleiotropic response resulting in differential gene expression of 445 transcripts could be attributed to increased transcript levels of regulatory proteins and other non-cognate toxins. Network analysis identified the WhiB7 regulon among the top-response network upon overexpression of VapC21. WhiB7 has been previously reported to be upregulated in response to various antibiotics and is associated with intrinsic drug resistance (Morris et al., 2005;Burian et al., 2012;Reeves et al., 2013). In addition to WhiB7 the transcript levels of another antibiotic response operon, iniBAC, which activates and provides the defense against cell wall inhibiting antibiotics isoniazid and ethambutol were also increased in the overexpression strain (Alland et al., 2000;Colangeli et al., 2005). In agreement, we also observed that overexpression of VapC21 led to an increase in amikacin-tolerant and ethambutol-tolerant persisters in M. smegmatis. Although, the contribution of Type II TA systems in bacterial persistence and drug tolerance are controversial, our results are in agreement with other studies, where overexpression of Type II toxins is associated with increased tolerance of E. coli, M. smegmatis or M. tuberculosis upon exposure to various drugs (Keren et al., 2004;Dorr et al., 2010;Singh et al., 2010;Holden and Errington, 2018;Ronneau and Helaine, 2019;Tandon et al., 2019a).
In summary, we have functionally characterized the VapBC21 TA system from M. tuberculosis. We report that in addition to cognate antitoxin, VapC21 is also able to interact with the non-cognate antitoxin, VapB32. We speculate that cross-talk between non-cognate TA pairs results in a complex regulatory network that fine tunes the activity of toxin and might enable the bacteria to reactivate from a dormant state. VapC21 is nonessential for in vivo growth but is involved in the generation of amikacin-tolerant and ethambutol-tolerant persisters in vitro. These observations suggests that mycobacterial ribonucleases belong to a regulatory network that contributes to disease reactivation and relapse.

DATA AVAILABILITY STATEMENT
The RNA-seq data discussed in the present study has been deposited in NCBI's Gene Expression Omnibus and are accessible through GEO Series accession number GSE154320.

ETHICS STATEMENT
The animal study was reviewed and approved by Translational Health Science and Technology Institute Animal Ethics Committee.  The open reading frame for vapC21 was replaced with the hygromycin resistance gene in the vapC21 mutant strain (MT1 and MT2), using temperature sensitive mycobacteriophages. (B) The replacement of vapC21 with the hygromycin resistance gene in the mutant strain was confirmed by Southern blot. Solid line (P) depicts the region for hybridization of probes. WT represent wild type strain of M. tuberculosis. (C) The colony morphology and biofilm formation of the wild type and vapC21 mutant strain was determined as described in section "Materials and Methods." (D-G) For stress experiments, early-log phase cultures (OD 600nm ∼ 0.2) were exposed to different stress conditions such as oxidative stress (D), nitrosative stress (E), nutrient starvation (F), and 2.5 mg/ml lysozyme (G) as described in section "Materials and Methods." The data shown in panels (D-G) are mean ± SE. obtained from triplicate samples and represents two independent experiments. For drug-tolerance experiment, mid-log phase cultures of various strains were exposed to different drugs for 14 days. For bacterial enumeration, 10.0-fold serial dilutions were prepared and plated on Middlebrook 7H11 at 37 • C for 3-4 weeks. The data shown in this panel is mean ± SE obtained from three independent experiments. FIGURE S3 | Co-expression experiments to determine interaction of VapC21 with cognate and non-cognate antitoxins. For co-expression studies, various pLam12-vapBx constructs were electroporated into M. smegmatis harboring an integrative Atc inducible copy of VapC21. The expression of toxin and antitoxin was induced by the addition of 50 ng/ml Atc and 0.2% acetamide, respectively. The growth of various strains was monitored by measuring OD 600nm and data shown in this figure is representative of two independent experiments. TABLE S1 | List of strains and plasmids used in the present study.