Identification of a Master Regulator of Differentiation in Toxoplasma

Abstract

Toxoplasma gondii chronically infects a quarter of the world's population, and its recrudescence can cause life-threatening disease in immunocompromised individuals and recurrent ocular lesions in the immunocompetent. Acute-stage tachyzoites differentiate into chronic-stage bradyzoites, which form intracellular cysts resistant to immune clearance and existing therapies. The molecular basis of this differentiation is unknown, despite being efficiently triggered by stresses in culture. Through Cas9-mediated screening and single-cell profiling, we identify a Myb-like transcription factor (BFD1) necessary for differentiation in cell culture and in mice. BFD1 accumulates during stress and its synthetic expression is sufficient to drive differentiation. Consistent with its function as a transcription factor, BFD1 binds the promoters of many stage-specific genes and represents a counterpoint to the ApiAP2 factors that dominate our current view of parasite gene regulation. BFD1 provides a genetic switch to study and control Toxoplasma differentiation and will inform prevention and treatment of chronic infections.

Keywords: Toxoplasma gondii; bradyzoite; chronic infection; differentiation; master regulator; single-cell RNA-sequencing.

Graphical abstract

Figure 1

A Genetic Screen Identifies a Putative Regulator of T. gondii Differentiation (A) Construction of a differentiation reporter strain that constitutively expresses RFP and Cas9 and conditionally expresses mNeonGreen (mNG) under the regulation of the bradyzoite-specific BAG1 promoter. SAG2Y is a bradyzoite-specific surface marker. Images were taken after 48 h of growth under unstressed or alkaline-stressed conditions. Scale bar is 10 μm. (B) Transfection and selection for a gRNA targeting the surface antigen SAG1 resulted in gene disruption in 98% of the resulting population. The mean was plotted for n = 3 biological replicates; 92–102 vacuoles were scored for each replicate; ^∗∗∗∗p < 0.0001 by Student’s two-tailed t test. (C) Percentage of alkaline-stressed reporter parasites expressing mNG, quantified by FACS. The mean ± SD was plotted for n = 2–4 biological replicates. (D) RNA sequencing and differential expression (DE) analysis identified 1,311 genes as significantly upregulated (green) and 933 genes as significantly downregulated (red) in bradyzoites (adjusted p < 0.001), with 1,240 and 700 changing 2-fold or more, respectively (dotted lines). The analysis was based on n = 3 independent experiments. (E) Screening and analysis workflow. The log₂-transformed fold changes from the input library to the final unstressed or stressed mNG⁺ populations are defined as fitness or differentiation scores, respectively. (F and G) Fitness and differentiation scores at the gene level following screening L1 (F) or L2 (G). (H) Fitness and differentiation scores for individual gRNAs in L2. See also Figure S1 and Data S1 and S2.

Figure S1

Stage-Specific RNA-Seq and Forward-Genetic Screening in a Differentiation Reporter Strain, Related to Figure 1 (A) Sample FACS plots for the bradyzoite reporter strain at 24 or 48 h of growth under unstressed or stressed conditions. 10,000 events per plot. Gates for mNG⁺ and mNG^– populations are highlighted. (B) Principal component analysis of stage-specific RNA-seq replicates. (C–D) Comparison of differentially expressed genes identified in this study and a previously published dataset (Fritz et al., 2012). Color assigned by adjusted p-value < 0.001 (C) or rank of base mean expression as calculated by DESeq2 (D). (E) Efficiency of construct integration into the bradyzoite reporter strain. Plaquing efficiency was compared between parasites transfected and selected for integration of a plasmid encoding a gRNA against SAG1 and the pre-transfection population to obtain a viability-normalized integration rate. Mean ± SD is plotted for n = 3 independent experiments. (F) Generation of a BFD1-deficient reporter strain by introducing a frameshift using Cas9. The bradyzoite reporter strain was transfected with a plasmid encoding a gRNA targeting the first exon of TGME49_200385 (BFD1) to isolate a strain with a frameshift mutation (BFD1^frameshift), resulting in a premature stop codon after amino acid 251.

Figure S2

An Updated BFD1 Gene Model and Phylogeny, Related to Figure 2 (A) Updated gene model and protein sequence of TGME49_200385. Sequencing of cDNA confirmed changes to the gene model after the 5^th exon annotated on ToxoDB v. 42, which change the reading frame of the latter third of the protein. DNA-binding domains (SM00717) highlighted in blue. (B) Neighbor-joining phylogenetic tree of SANT/Myb-like DNA-binding domains (SM00717) present in representative apicomplexan genomes, along with human c-Myb and CDC5L. Clades containing c-Myb and CDC5L are highlighted in blue and orange, respectively. Alignment performed using ClustalW. Scale bar is substitutions per site.

Figure 2

BFD1 Is a Nuclear Factor Necessary for Differentiation in Cell Culture (A) Neighbor-joining tree showing the phylogenetic relationship of the concatenated Myb-like domains from BFD1 and its closest homologs in other apicomplexans and humans. Tissue-cyst-forming species are indicated (green circles). Bootstrap values for 1,000 trials are displayed. Scale bar is substitutions per site. (B) Diagram of BFD1 and human c-Myb highlighting the Myb-like domains (blue). The DNA-binding repeats of BFD1 are similar to the second and third repeats of c-Myb. (C) Generation of ΔBFD1 and BFD1^WT or BFD1^ΔMYB complemented parasites. To create ΔBFD1 parasites, the endogenous coding sequence was replaced with a fluorescent cassette. The knockout was complemented with a wild-type (WT) (BFD1^WT) or DNA-binding deficient (BFD1^ΔMYB) Ty-tagged allele at the endogenous locus. (D) Plaque assays of indicated strains grown under unstressed conditions for 14 days. Scale bar is 1 cm. (E) Representative vacuoles after 48 h of alkaline stress. FITC-labeled Dolichos biflorus lectin (DBL) specifically stains differentiated vacuoles. Ty was stained with BB2 (magenta), and DNA was stained with Hoechst (blue). Scale bar is 10 μm. (F) Quantification of differentiation in WT, knockout, and complemented parasites following 48 h of alkaline stress, 48 h of compound 1 treatment, or occurring spontaneously under unstressed conditions in the same time frame. The mean ± SD was plotted for n = 3–8 biological replicates, with percentage of DBL positive vacuoles calculated from at least 100 vacuoles per replicate. ^∗∗∗∗p < 0.0001, ^∗p < 0.05, Student’s one-tailed t test. See also Figure S2.

Figure 3

BFD1 Is Necessary for Formation of Brain Cysts in Mice (A) Timeline of mouse infections. Groups of CD-1 female mice were inoculated i.p. with 500 tachyzoites per animal from each strain or mock inoculated with PBS. Cyst formation was assayed in moribund animals starting 2 weeks post-infection and in all surviving animals at 5 weeks post-infection. (B) Mean normalized weights of animals in each group. Graph represents mean ± SEM for all surviving animals at a given time point. Graphs are for n = 5 mock-inoculated mice and n = 15 for each parasite strain. (C) Survival curve of animals in (B). (D) Representative cysts from WT and ΔBFD1::BFD1^WT-infected animals. The cyst wall was stained with DBL (green) and individual parasites with anti-CDPK1 (magenta). Scale bar is 20 μm. (E) Cyst burden per animal, denoting those sacrificed before (open circles) or after (closed circles) 5 weeks of infection. Cysts per brain were estimated from counting four blinded replicates, with a limit of detection of 56–71 cysts per brain, depending on the volume of the sample analyzed. Mean is plotted for each group. ^∗∗p < 0.01, Student’s one-tailed t test. See also Figure S3.

Figure S3

Virulence and Brain Cyst Formation by ΔBFD1 Parasites in CBA/J Mice, Related to Figure 3 (A–B) In groups of 5, female CBA/J mice were inoculated with 100 or 2,000 tachyzoites i.p. of WT, ΔBFD1, or ΔBFD1::BFD1^WT and surviving animals were sacrificed 2 weeks post-infection to assay brain cyst formation (A). Cyst burdens were estimated by counting 4 blinded samples from each animal. Mean ± SD is plotted with each dot representing an animal; ^∗∗∗∗p < 0.0001, Student’s one-tailed t test (B). (C–E) In groups of 5, female CBA/J mice were inoculated with 500 or 10,000 tachyzoites i.p of WT, ΔBFD1, or ΔBFD1::BFD1^WT. Starting at 3 weeks post-infection, brains were isolated from moribund animals, and at 5 weeks post-infection all surviving animals were sacrificed (C). Survival curve of animals infected with 10,000 (dotted lines) or 500 (solid lines) tachyzoites (D). Brain cyst burden of moribund or sacrificed animals, estimated by counting 4 blinded samples from each animal. Mean ± SD is plotted with each dot representing an animal (E).

Figure S4

Profiling Toxoplasma Differentiation at Single-Cell Resolution, Related to Figure 4 (A–B) Distribution of UMIs (A) or unique genes detected (B) across single cells from indicated samples and time points. Pre-processing quality control cutoffs required a minimum of 200 and a maximum of 10,000 UMIs. (C) Percentage of UMIs corresponding to ribosomal genes. Pre-processing quality control cutoffs allowed a maximum of 40% rRNA reads. (D) The majority of variance between cells is driven by cell-cycle and stage-specific genes. Plotting the 18 principal components (PCs) determined to be statistically significant by permutation analysis. The first three PCs explain 66.4% of variance at the 72 h time point. (E) Pearson correlations of cell embeddings in PCs 1–18 to cell scores for G1, S/M, or bradyzoite-specific gene signatures. (F) Violin plots of expression of early bradyzoite marker genes in wild-type parasites after 24, 48 or 72 h of growth under unstressed or stressed conditions. (G) UMAP visualization as in Figure 4E colored by expression of CST1 or other early bradyzoite markers. (H) Endogenous tagging of two early bradyzoite markers shows localization to the cyst wall. Cultures were fixed and stained after 72 h under alkaline stress. Scale bar is 10 μm. (I) BFD1 is not stage-specific but is modestly overrepresented in bradyzoite-containing clusters. UMAP visualization as in Figure 4E colored by expression of BFD1. (J) UMAP visualizations colored by expression of marker genes identified as specifically upregulated in bradyzoites during the S/M phase of the cell cycle.

Figure 4

BFD1 Is Required to Express Bradyzoite Genes and Initiate Differentiation (A) Clustering of unstressed parasites of both genotypes after 72 h of growth (3,149 cells total) visualized by uniform manifold approximation and projection (UMAP). (B) UMAP from (A) shaded by score for expression of known G1 or S/M-specific gene sets. (C) Average expression profiles of all genes differentially expressed within each cluster identified in (A) across a microarray dataset of synchronized tachyzoites (Behnke et al., 2010). Colors correspond to those used in (A). (D) Proportion of cells in G1 (0 and 1) or S/M (2, 3, 5, and 6) clusters. (E) Cells shaded by expression of the canonical stage-specific genes SAG1 and BAG1 following clustering of all parasites from all time points, genotypes, and growth conditions. UMAP visualization is downsampled to 500 cells from each combination of time point, genotype, and growth condition (6,000 cells total). (F) UMAP as in (E), with cells highlighted by sample of origin for WT (green) or ΔBFD1 (blue) parasites. (G) Distribution of cell scores for the expression of genes highly upregulated in bradyzoites. (H) UMAP of all WT and ΔBFD1 parasites from unstressed and stressed cultures at 72 h time point. (I and J) UMAP projection from (H), shaded by scores for expression of highly upregulated bradyzoite genes as in G, or S/M-specific genes as in (B) (I) or colored by sample of origin (J). (K) Representative WT and ΔBFD1 vacuoles at 72 h post-alkaline stress. GAP45 is a marker for the inner membrane complex. Scale bar is 10 μm. (L) Proposed model of cell-cycle progression for WT and ΔBFD1 parasites under the different treatments. See also Figure S4 and Data S3.

Figure 5

BFD1 Is Sufficient to Drive Differentiation in the Absence of Stress (A) Normalized counts for selected genes in unstressed or stressed stage-specific RNA-seq replicates from Figure 1D. (B) Representative WT and ΔBFD1::BFD1^WT-Ty vacuoles after 48 h of growth under unstressed or stressed conditions. Ty was immunostained with BB2 (magenta). Scale bar is 5 μm. (C) Immunoblot of lysates from unstressed or stressed ΔBFD1::BFD1^WT-Ty parasites for the presence of BFD1-Ty or CDPK1. The expected molecular weight of BFD1 is 262.5 kDa. (D) Generation of a conditionally stabilized BFD1 overexpression strain. Addition of the stabilizing-ligand Shield-1 inhibits degradation of BFD1. (E) Representative images of ΔBFD1/DD-BFD1-Ty parasites grown in standard medium with vehicle (left) or 3 μM Shield-1 (right) added. After 4 days of growth, parasites grown in vehicle lyse the host cell monolayer, while parasites treated with Shield-1 continue to grow intracellularly. Scale bar is 20 μm. (F) Timeline to assay differentiation following induction of BFD1 expression. After 4 h, standard medium was replenished adding vehicle or to 3 μM Shield-1. After an additional 48 h, differentiation was quantified by DBL positivity. (G) Representative vacuoles of ΔBFD1/DD-BFD1-Ty parasites grown for 48 h in standard medium with vehicle or 3 μM Shield-1 and stained for Ty or DBL. Scale bar is 10 μm. (H) Quantification of (G). Mean plotted for n = 3–4 biological replicates, with DBL-positive vacuoles counted from at least 100 vacuoles per replicate. ^∗∗∗∗p < 0.0001, Student’s one-tailed t test. (I) RNA sequencing and differential expression (DE) analysis identified 1,998 genes as significantly upregulated (green) and 2,110 genes as significantly downregulated (red) in ΔBFD1/DD-BFD1-Ty parasites grown for 48 h in 3 μM Shield-1 compared to vehicle alone (adjusted p < 0.001). Of these, 585 and 655 genes changed 2-fold or more (dotted lines). (J) Comparison of significantly regulated genes observed under alkaline stress or BFD1 expression (adjusted p < 0.001, 2-fold change or greater; points below this threshold were not plotted). Pearson correlation is indicated. (K) Overlap between differentially expressed genes plotted in (J). See also Figure S5 and Data S4.

Figure S5

Conditional or Transient Overexpression of BFD1 Is Sufficient to Induce Differentiation in the Absence of Stress, Related to Figure 5 (A) Distribution of lengths of 5,699 previously annotated Toxoplasma 5′ UTRs. By manual annotation, the 5′ UTR of BFD1 is 2,709 bp (green arrow), placing it in the 98^th percentile. (B) Diagram of the BFD1 locus indicating RNA-seq read density from tachyzoites (blue) and indicating the position of upstream open reading frames (uORFs) in the 5′ UTR, and the BFD1 coding sequence. (C–D) Constructs and experimental workflow for transient expression of BFD1. Epitope-tagged cDNA versions of wild-type BFD1 (BFD1^WT) or BFD1 lacking its DNA-binding domains (BFD1^ΔMYB) are under the regulation of the TUB1 promoter (C). Parasites were immuno-labeled for Ty (magenta) and for differentiation with FITC-conjugated DBL (green) and differentiation in WT or ΔBFD1 parasites 48 h after transient overexpression of BFD1^WT or BFD1^ΔMYB was quantified. Scale bar is 10 μm. Ty⁺ vacuoles were identified and then scored for DBL positivity as shown in representative images. Mean ± SD is plotted for n = 2 independent replicates, 17–61 vacuoles counted per replicate; ^∗p-value < 0.05, ^∗∗p-value < 0.01; Student’s one-tailed t test (D). (E) Representative images of WT, ΔBFD1, or ΔBFD1/DD-BFD1-Ty parasites grown in standard medium supplemented with vehicle or 3 μM Shield-1. After 4 days of growth, all host monolayers had been lysed except for the Shield-1-treated ΔBFD1/DD-BFD1-Ty cultures in which parasites continued replicating intracellularly. Scale bar is 10 μm. (F) Principal component analysis of RNA-seq of WT, ΔBFD1, or ΔBFD1/DD-BFD1-Ty parasites grown in media containing 3 μM Shield-1 or vehicle alone.

Figure 6

BFD1 Binds to Transcriptional Start Sites of Differentially Regulated Genes (A) Overview of the CUT&RUN protocol. Antibody-directed MNase activity preferentially creates short fragments inclusive of the region bound by the protein of interest. These short fragments diffuse out of the nucleus and can then be enriched and sequenced. (B) Distribution of BFD1 peaks ordered according to their quality score, plotted relative to the nearest transcriptional start site (TSS). Each row is a single gene with one or more associated BFD1 peaks. Genes are categorized according to their differential expression in stage-specific RNA sequencing. (C) MEME or HOMER motifs significantly enriched at BFD1 binding sites. (D) Genes with BFD1 binding motifs within a −443 to +494 bp window of their TSS tend to be upregulated following alkaline stress. Distribution of fold changes for genes in each category. ^∗∗∗p < 0.001, ^∗∗∗∗p < 0.0001, Mann-Whitney test. (E) Representative loci indicating binding of BFD1 upstream of many of the most differentially expressed genes and its own transcriptional start site. See also Data S5.