Loss of SNAI1 induces cellular plasticity in invasive triple-negative breast cancer cells.

The transcription factor SNAI1 mediates epithelial-mesenchymal transition, fibroblast activation and controls inter-tissue migration. High SNAI1 expression characterizes metastatic triple-negative breast carcinomas, and its knockout by CRISPR/Cas9 uncovered an epithelio-mesenchymal phenotype accompanied by reduced signaling by the cytokine TGFβ. The SNAI1 knockout cells exhibited plasticity in differentiation, drifting towards the luminal phenotype, gained stemness potential and could differentiate into acinar mammospheres in 3D culture. Loss of SNAI1 de-repressed the transcription factor FOXA1, a pioneering factor of mammary luminal progenitors. FOXA1 induced a specific gene program, including the androgen receptor (AR). Inhibiting AR via a specific antagonist regenerated the basal phenotype and blocked acinar differentiation. Thus, loss of SNAI1 in the context of triple-negative breast carcinoma cells promotes an intermediary luminal progenitor phenotype that gains differentiation plasticity based on the dual transcriptional action of FOXA1 and AR. This function of SNAI1 provides means to separate cell invasiveness from progenitor cell de-differentiation as independent cellular programs.

Introduction

Breast cancer (BRCA), the most widespread tumor among women, is a heterogeneous disease with diverse pathological features, molecular signatures and clinical outcomes [1]. Inter-tumoral heterogeneity originates from the distinct mammary epithelial cell types that serve as the cell of origin for the accumulation of oncogenic processes [1]. Alternatively, mechanisms of de- or trans-differentiation of mammary epithelial cells eventually cause phenotypic and molecular heterogeneity in breast tumors among patients [2, 3]. Accordingly, breast cancers are classified as tumors expressing estrogen receptor-α (ERα/ESR1), progesterone receptor (PGR) and/or epidermal growth factor (EGF) family receptor HER2/ERBB2, which emanate from luminal epithelial progenitors and are classified as luminal-A or -B. Triple-negative breast cancers (TNBC) do not express the above three receptors and are subdivided into basal-like A and -B, mesenchymal (or claudin-low), mesenchymal stem-like, luminal androgen receptor and immunomodulatory [2, 4–7]. Interestingly, the majority of TNBCs, basal-like BRCA, exhibit gene expression profiles which are closely related to luminal progenitors and mammary stem cells, while oncogenic transformation of luminal cells generates tumors with basal-like characteristics [8-11].

A process of de- or trans-differentiation within the tumorigenic mammary tissue is the epithelial-mesenchymal transition (EMT) [12, 13]. Cells undergoing EMT remodel their secreted extracellular proteins and express lower numbers of E-cadherin, tight junctional components and cytokeratins, accompanied by higher numbers of N-cadherin, vimentin and certain integrins [13, 14]. A given breast tumor is composed of mixtures of epithelial, mesenchymal and intermediate phenotypes, with evidence of mesenchymal-epithelial transition (MET) [15, 16]. This heterogeneity explains the recent TNBC sub-classification explained above [5, 7, 17]. The TNBC features are reproduced well in many BRCA cell lines [18, 19], such as the mesenchymal and metastatic MDA-MB-231, originally derived from a patient pleural effusion [20]. The plasticity changes observed in breast cancer have been explained by the differential activities of signaling and transcriptional programs. A prominent signaling pathway in this respect is the transforming growth factor β (TGFβ) pathway that drives the EMT in every cell type examined so far [14]. TGFβ signaling is directly coupled to the transcriptional induction of a cohort of transcription factors that promote EMT (EMT-TFs), such as SNAI1, SNAI2, ZEB1, ZEB2, TWIST1, and TWIST2, many of which (SNAI1, ZEB1, ZEB2) associate with TGFβ signal transducers of the SMAD family and control their activity [14].

The current study focuses on the EMT-TF SNAI1, an evolutionarily conserved zinc finger transcription factor, that together with SNAI2/SLUG, forms a small family of EMT inducers during embryogenesis [21, 22]. A central function of SNAI1 is the transcriptional repression of the E-cadherin (CDH1) gene [23, 24], of tight junctional genes [25] and of the fructose-1,6-bisphosphatase gene that controls the rate of glycolysis [26]. Genome-wide chromatin association studies have revealed multiple genes to which SNAI1 binds in breast and colorectal cancer cells [27, 28]. In ovarian cancer cells, SNAI1 function was linked to cell-cell and cell-matrix adhesion [29]. Silencing of SNAI1 in MDA-MB-231 cells revealed loss of invasiveness in vitro, accompanied by reduced tumor growth and metastatic potential in xenografted mice [30, 31]. SNAI1 knockout in additional TNBC cells (HCC38 and Hs578T) indicated the generation of intermediary quasi-epithelial phenotypes, explained by the remaining and compensatory expression of SNAI2 in the same cells [28, 32]. Yet, these studies could not identify possible links between the altered TNBC phenotype and key transcription factors that define mammary epithelial cell differentiation. Loss of SNAI1 may lead to lineage plasticity that could explain some of the phenotypes in heterogeneous breast cancers. This hypothesis is examined in the present study that links SNAI1 knockout to basal-luminal plasticity.

Materials and methods

Cells, reagents and treatments

Human TNBC-basal MDA‐MB‐231 cells obtained from ATCC (HTB-26, LGC Standards, UK) and Hs578T cells obtained from Dr. Paraskevi Heldin were cultured in Dulbecco’s modified Eagle’s medium (DMEM; Sigma-Aldrich AB, Stockholm, Sweden). Luminal BRCA MCF-7 and ZR-75-1 cells obtained from Dr. Paraskevi Heldin were cultured in DMEM and in Roswell Park Memorial Institute (RPMI)‐1640 (Sigma-Aldrich AB), respectively. All culture media were supplemented with 10% fetal bovine serum (FBS; Biowest, Almeco A/S, Esbjerg, Denmark), 100 U/ml penicillin and 100 μg/ml streptomycin (Sigma-Aldrich AB). All cells were kept in a humidified incubator at 37 °C and 5% CO2. All cell lines were frequently analyzed and found free of mycoplasma and authenticated using STR analysis. Cells were starved for 18 h in serum-free DMEM or RPMI and treated for different time periods with recombinant human TGFβ1 (PreproTech EC Ltd, London, UK) at a final concentration of 5 ng/ml, as indicated. Cells were treated with the indicated concentrations of enzalutamide (MDV3100, Selleckhem, Houston, USA) in DMEM or RPMI containing 10% charcoal-stripped FBS.

CRISPR/Cas9 knockout

MDA-MB-231 cells were transfected with CRISPR/Cas9 and HDR plasmids targeting SNAI1 (Santa Cruz Biotechnology Inc., CA, USA). A pool of three plasmids, each encoding one guide RNA targeting exon-1, exon-2 or exon-2/intron-2 junction, was co-transfected with the Cas9 nuclease and HDR plasmids. Two days post‐transfection, cells were selected with 4 μg/ml puromycin (Merck/Millipore, Stockholm, Sweden) and single-cell colonies were expanded. Knockout clones were validated using immunoblotting and quantitative RT-PCR to analyze the mRNA levels of each of the three SNAI1 exons.

Transient transfection with siRNA

Cells were transfected with siRNAs (20 nM; ON-TARGETplus Non-targeting Control Pool (D-001810-10), ON-TARGETplus Human FOXA1 set of four individual siRNAs (LQ-010319-00-0002, Dharmacon/VWR, Stockholm, Sweden) using Lipofectamine RNAiMAX according to the manufacturer’s instructions (Thermo Fisher Scientific, Stockholm, Sweden).

Plasmid transfections

SNAI1-KO cells were transfected using X-tremeGENE 360 transfection reagent (Merck/Millipore) with the peGFPNI SNAI1 plasmid kindly provided by Dr. Antonio Garcia de Herreros (IMIM-Hospital del Mar, Barcelona, Spain). After 48 h, cells expressing high GFP were sorted from the pool of the transfected cells by flow cytometry in 96-well plates and single-cell colonies were grown in the presence of 2.0 mg/ml geneticin (Thermo Fisher Scientific) for two weeks, generating stable SNAI1 overexpression clones.

3D culture of human mammospheres

Mammosphere cultures of MDA-MB-231 and ZR-75-1 cells were performed using the hanging drop method by seeding 10 000 cells per drop in complete growth medium with 20% methylcellulose for 48 h. Mammospheres were collected and resuspended in MEM-F12 (Sigma-Aldrich AB) supplemented with 25 ng/ml epidermal growth factor, 25 ng/ml basic fibroblast growth factor, 1×B27 (Thermofisher Scientific) and 2% methylcellulose (Sigma-Aldrich AB). One mammosphere was placed per well in 96-well round bottom ultra-low attachment microplates (CORNING, Wiesbaden, Germany) and cultured for 5 days. Phase-contrast pictures were acquired every 24 h for 5 days. The cross-section area of the mammospheres was calculated using the Image-J software (National Institutes of Health, Bethesda, MD, USA).

Cell proliferation

For proliferation assessment, 50,000 MDA-MB-231 or ZR-75-1 cells per well were seeded in 6-well plates, and cell counting was performed on day five after seeding using the Luna automated cell counter (Logos Biosystems, France). The number of viable cells was used for further analysis.

Cell viability

Cytotoxicity to decreasing serial dilutions of doxorubicin (D1515; 2 μM–7.8 nM), paclitaxel (T7402; 1 μM–3.9 nM), both from Sigma-Aldrich AB and to enzalutamide (MDV3100, Selleckhem; 50–0.1 μM) was monitored at day 2 and 5 by PrestoBlue HS reagent, following the manufacturer’s instructions (Thermo Fisher Scientific). Dimethyl-sulfoxide (DMSO) served as a vehicle. Cells treated with enzalutamide were cultured in charcoal-stripped FBS.

Extreme limiting dilution assay (ELDA)

Cells were seeded in low-attachment 96-well plates (CORNING) in decreasing serial dilutions (100-1 cells/well), in 200 μl of serum-free MEM-F12 (Sigma-Aldrich AB) supplemented with 25 ng/ml epidermal growth factor, 25 ng/ml basic fibroblast growth factor, 1×B27 (Thermofisher Scientific). Six technical replicates for each cell plating density were created. On day 10, the number of wells containing spheres was recorded and analyzed using the online ELDA analysis program (http://bioinf.wehi.edu.au/software/elda) [33].

Invasion of mammosphere cells

Mammospheres of MDA-MB-231 and ZR-75-1 cells generated as described above, were collected after 48 h and resuspended in a collagen I solution (1.7 mg/ml, PureCol, Advanced BioMatrix, Inc., San Diego, CA, USA) in DMEM or RPMI. One mammosphere was embedded per well, and phase-contrast images were acquired at the onset of embedding and after 24 and 48 h. Invasive growth was calculated as the area occupied by cells between the invaded area and the core mammosphere using Image-J (National Institutes of Health).

Zebrafish extravasation assay

Staging and embryo production of Tg(Fli1:EGFP) zebrafish (Danio rerio), whose vasculature is marked in green, were conducted as described [34] and embryos were maintained at 34 °C. Empirical determination of sample size that provided power for discrimination between conditions was used and for this reason 200 embryos were injected per condition in order to reach a final number of viable embryos of more than 100. No randomization method was applied and the microinjector was blinded to the groups of injected cancer cells. MDA-MB-231-WT and SNAI1-KO cells were stained with 4 ng/μl CM-Dil Dye (ThermoFisher Scientific) for 30 min at 37 °C. At 48 h post-fertilization, approximately 400 CM-Dil Dye-labeled cells were loaded into borosilicate glass capillary needles (1 mm O.D. × 0.78 mm I.D, Harvard Apparatus, Holliston, MA, USA), and were injected into the duct of Cuvier of dechorionized embryos, anesthetized with 0.003% tricaine (Sigma-Aldrich AB), and mounted on a 10-cm Petri dish coated with 1% agarose gel, using a Pneumatic Picopump and a manipulator (WPI, Stevenage, UK). Based on pre-established criteria, only microscopy-verified, correctly injected and viable zebrafish were retained at 34 °C and imaged automatically (ImageXpress Nano, Molecular Devices, USA) 24 h post-implantation; cells extravasated from the circulation at the posterior part of 150 zebrafish per cell line were counted.

Immunoblotting

Total cellular proteins were extracted in lysis buffer (20 mM Tris-HCl, pH 7.5, 1% nonidet P‐40, 150 mM NaCl) supplemented with a complete protease inhibitor cocktail (Roche Diagnostics Scandinavia AB, Bromma, Sweden). The lysates were sonicated for 5 min at 4 °C (30 sec ON/30 sec OFF), cleared by centrifugation, and protein concentration was measured by Bradford assay (Sigma-Aldrich AB). Equal amounts of protein (40 μg) were subjected to sodium dodecyl sulfate polyacrylamide gel electrophoresis and transferred to a PVDF membrane using a Bio-Rad wet transfer unit (Bio-Rad Laboratories Inc., Sundbyberg, Sweden). Filters blocked with 5% bovine serum albumin (BSA; Saveen Werner, Limhamn, Sweden) were incubated with primary antibodies and horseradish peroxidase-conjugated anti-mouse or anti-rabbit secondary antibodies (Thermofisher Scientific) as listed in Table S7, followed by enhanced chemiluminescence assays using the Millipore kit (Merck/Millipore). The original, uncropped immunoblots are listed as the last Supplementary figure.

Immunofluorescence microscopy

MDA-MB-231 and ZR-75-1 cells fixed in 3.7% formaldehyde for 15 min, were incubated in 0.1 M glycine for 45 min and permeabilized in 0.5% Triton X-100 for 10 min, blocked in 5% BSA/PBS for one h at room temperature and incubated with primary antibodies (Table S7) in 1% BSA/PBS overnight at 4 °C. Samples were incubated with Alexa-Fluor-488 or Alexa-Fluor-546 secondary antibodies (Thermofisher Scientific) at a dilution of 1:500 in PBS for 1 h at room temperature. Phalloidin (1:500 in PBS) staining lasted for either 1 h at room temperature (2D) or overnight at 4 °C (3D). The nuclei were stained afterward with 4′,6-diamidino-2-phenylindole (DAPI, Sigma-Aldrich AB) at a dilution 1:1000 in PBS for 10 min at room temperature. Coverslips mounted with Fluoromount-G (Southern Biotech, AH Diagnostics, Solna, Sweden) were examined by a Nikon Eclipse 90i fluorescence microscope (Nikon Corp., Tokyo, Japan) or a Zeiss LSM700 confocal microscope (Carl Zeiss AB, Stockholm, Sweden). Five to 6 random pictures were taken with a 10× or 20× objective at the same exposure time.

Chromatin immunoprecipitation

Cells were harvested (20 × 106) by trypsinization followed by PBS wash at room temperature and crosslinked using 10 ml of 1% formaldehyde (dissolved in PBS) for 10 min at room temperature and kept on gentle shaking. The reaction was quenched with 0.125 M final glycine concentration for 5 min, followed by centrifugation at 1000 × g at 4 °C for 10 min. The crosslinked pellets were washed once with cold PBS followed by a centrifugation step at 1000 × g at 4 °C for 10 min. The cell pellet was resuspended in 1 ml of SDS lysis buffer (0.1% SDS, 0.5% Triton X-100, 20 mM Tris-HCl, pH 8, and 150 mM NaCl, 1 mM PMSF, and 1× protease inhibitors) and incubated on ice for 30 min with continuous pipetting. Then, the chromatin of the lysed pellet was sheared using a Bioruptor for a total of 30 cycles (30 s ON, 30 s OFF at high pulse). The insoluble components were removed by centrifugation at 13,000 × g at 4 °C for 10 min, while 5% of the total lysate served as input. Cleared lysates were incubated with the respective antibodies (FOXA1 or IgG) for immune precipitation overnight at 4 °C. An aliquot of 4 μg of antibody was used per 1 mg of lysate protein. We allowed the immune complexes to bind to Dynabeads Protein A (ThermoFisher Scientific) for 3 h at 4 °C with gentle rotation and then the immune complexes bound to beads were separated using magnetic precipitation followed by two washes of 750 µl low salt buffer (0.1% SDS, 1% Triton-X 100, 2 mM EDTA, 20 mM Tris-HCl, pH 8, 150 mM NaCl, 0.5 mM PMSF and 1× protease inhibitors) for 5 min at 4 °C with gentle rotation. The immune complexes were then washed with 750 µl of high salt buffer (0.1% SDS, 1% Triton-X 100, 2 mM EDTA, 20 mM Tris-HCl, pH 8, 500 mM NaCl, 0.5 mM PMSF and 1× protease inhibitors) for 5 min at 4 °C with gentle rotation. To elute immune-precipitated material from the beads, we added 400 μl of elution buffer (50 mM Tris-HCl, pH 8, 20% SDS, 10 mM EDTA and 0.5 mM PMSF) and incubated at 55 °C for 30 min with agitation. A de-crosslinking step was performed by incubating the beads with 300 µl of proteinase K buffer (100 mM NaCl, 10 mM Tris-HCl pH 7, 0.5% SDS, 10 µg/ml Proteinase K) at 65 °C overnight. DNA isolation was performed using the QIAquick PCR purification kit (QIAGEN). Sequences of oligonucleotide primers used in ChIP assays are listed in Table S8. The data were plotted as fold-enrichment of the specific ChIP signal relative to the IgG control and presented as average values with standard deviations of four biological replicates, each with triplicate technical repeats.

RNA analysis

One μg of total RNA extracted with Nucleospin RNA-plus Kit (Macherey-Nagel, AH Diagnostics, Solna, Sweden) and quantified on a NanoDrop 2000 (Thermofisher Scientific) was reverse transcribed using the iScript cDNA synthesis kit (Bio-Rad Laboratories Inc.) according to the manufacturer’s protocol. Sequences of oligonucleotide primers used in RT-qPCR assays are listed in Table S8. Specific target gene expression was normalized to the reference gene GAPDH.

Real-time PCR analysis

Real-time PCR of cDNA or chromatin DNA (ChIP analysis) samples was performed on a Bio-Rad CFX96 cycler (Bio-Rad Laboratories Inc.) using the qPCRBIO SyGreen 2× Master Mix (PCR Biosystems, London, UK) and primers (Table S8). Specific amplicon numbers were normalized to a reference gene (usually GAPDH). DNA amplification was calculated based on the ΔΔCt method and plotted as average values with standard deviations of at least three biological replicates, each with triplicate technical repeats.

RNA-seq analysis

Five hundred nanograms total RNA extracted using ReliaPrep™ RNA Cell Miniprep System (Promega, Madison, WI, USA), were checked using the Agilent-2100 Bioanalyzer System and subjected to library preparation utilizing the TruSeq stranded total RNA Gold library preparation kit with RiboZero Gold treatment and unique dual indexes according to the manufacturer’s instructions (Protocol # 1000000040499, Illumina Inc., San Diego, CA, USA). The paired-end reads of 150 bp were generated in an S4 flowcell with v1.5 sequencing chemistry on a NovaSeq-6 000 platform (Illumina Inc.). Read quality was investigated by FastQC tool v0.11.9, adapters and low-quality reads were trimmed by trimmomatic tool v0.36 [35], and trimmed reads were aligned against the reference human genome (GRCh38) using STAR aligner v2.7.2b with two-pass mode [36]. After checking alignment quality by SAMStat tool v1.5.1 [37], high-quality (Q = 30) aligned reads were annotated and quantified against the gencode comprehensive gene annotation release 38 (GRCh38.p13) using featureCounts of subread package v2.0.0 [38]. Differential gene expression analysis was performed by DESeq2 Bioconductor package in R [39], thus calculating normalized counts per million reads. We considered cut-off criteria of log2 fold-change ±2 and false discovery rate (FDR) < 0.05 for differential expression, visualized in RStudio v1.4.1717 with R v4.0.5. All primary RNA-sequencing data are available at GEO (accession number GSE210870; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE210870).

Pathway-enrichment, network analysis and ChIP-sequencing annotation

We investigated enriched hallmarks and biological processes based on transcriptome profiles utilizing the gene set enrichment analysis (GSEA) tool and the molecular signature database MSigDBv6 [40, 41]. Pathway- and transcription factor-enrichment analysis of differentially expressed genes was based on the g:Profiler package in R [42], and network analysis on EnrichmentMap plugin-v3.3.3 [43] and visualization via Cytoscape-v3.8.2 [44]. We annotated ChIP-sequencing BED files using the annotatePeaks.pl program integrated into HOMER v4.11 suite [45].

TCGA and Cancer Cell Line Encyclopedia (CCLE) data analysis

The log2 transformed transcript per kilobase million (TPM) values of the queried genes were obtained from the Pan-Cancer Atlas analysis of the Cancer Genome Atlas (TCGA) [46] available at cBioportal [47]. We examined the expression distribution of the queried genes in different subtypes and compared the median expression difference using the Wilcoxon rank-sum test (U test). Multinomial logistic regression analysis was based on the neural network-based nnet package in R. TCGA samples (n = 945) formed a training (n = 609) and a validation set (n = 336), before testing the multivariate model in the training set, and investigating the predictive power of the model in the validation set. Statistical significance of the variables was based on a two-tailed Z test and a confusion matrix assessed performance of the model. For binomial logistic regression analysis, we classified TCGA samples into basal and non-basal tumors, and implemented a generalized linear model to fit the data in R. Statistical significance of the variables was assessed by ANOVA testing, and performance of the predictive model by ROCR package in R [48]. For cell line data, we downloaded TPM pseudo-counts of 61 BRCA cell lines with the metadata from the CCLE [49], available at the depmap portal [50]. For the compatibility with CCLE data, we re-aligned and calculated the log2 pseudo-counts TPM of MDA-MB231-WT and SNAI1-KO cells using the RSEM tool [51]. Heatmap visualization and hierarchical clustering were done using the hclust algorithm in pheatmap package in R.

Survival analysis

Relapse-free survival (RFS) and distant metastasis-free survival (DMFS) in basal BRCA patients were predicted using the KM plotter database [52]. Samples were stratified, compared using auto-select cut-off and p-value was calculated by log-Rank test.

Data analysis and statistics

Results express mean values of at least three independent experiments (biological repeats) as explained in the figures. The precise number of replicates is indicated in every figure legend. The number of technical and biological replications was determined by the efficiency of the cell-based assay used. Error bars represent standard error of the mean (SEM). We selected the appropriate statistical method based on sample content and variation within each group of data that was compared. The variance was similar between the groups that have been compared. Accordingly, two-group comparisons were performed using either two-tailed unpaired Student’s t-test or Wilcoxon rank-sum test (U test). Statistical significance is represented by p-values *p ≤ 0.05, **p ≤ 0.01, ***p ≤ 0.001. Additional statistical methods are described in the previous method sections.