Metastasis-associated MCL1 and P16 copy number alterations dictate resistance to vemurafenib in a BRAF^V600E patient-derived papillary thyroid carcinoma preclinical model

Introduction

The BRAFV600E mutation is the most common genetic alteration in papillary thyroid carcinoma (PTC) and may be associated with progression of PTC to anaplastic thyroid carcinoma (ATC) [1, 2], one of the most lethal human cancers [3, 4]. BRAFV600E is present in about 61% of PTCs as recently confirmed by the PTC TCGA (The Cancer Genome Atlas) and can be considered to be the primary genetic hallmark of PTC [5]. PTC patients harboring BRAFV600E mutation show resistance to radioiodine treatment [4, 6] [7], and have higher rates of recurrence and metastases, and lower survival rates [8-12]. Clearly, new therapeutic options are needed for metastatic and radioiodine-resistant thyroid cancers like BRAFV600E-positive PTC or ATC [13, 14]. Vemurafenib is the first orally available selective inhibitor of BRAFV600E approved by the FDA (Food and Drug Administration) for the treatment of BRAFV600E-melanoma [15, 16]. Although vemurafenib has recently shown promising clinical activity in three patients with metastatic PTC [17], its lack of response with resistance has been frequently observed [18] [19]. To address the unmet clinical need in BRAFV600E-positive metastatic PTC, we have studied the potential mechanism of loss of responsiveness to the BRAFV600E inhibitor vemurafenib, using a strategy referred to as a pre-clinical trial/model. We established short-term primary cell cultures of human thyroid samples, including BRAFV600E-PTC samples derived from metastatic/recurrent PTC, intrathyroidal primary PTC, and human normal thyroid (NT). Here, we describe the development of the first patient-derived BRAFWT/V600E-PTC in vitro and in vivo early intervention pre-clinical model with some similar disease molecular characteristics that are recapitulated. More importantly, this model offers interpretative insight into the concurrent vemurafenib human clinical trials in an independent cohort of patients with metastatic BRAFV600E-PTC, may serve to provide rapid clinical translation of our findings. We aim to investigate somatic copy number alterations (SCNAs) which could be associated with metastatic BRAFWT/V600E-PTC and mechanistically render these carcinomas resistant to the effects of BRAFWT/V600E inhibitors (e.g. vemurafenib) on cell death. We identify high copy number gain of MCL1 (myeloid cell leukemia sequence 1, chromosome 1q) and loss of CDKN2A (P16, chromosome 9p) in metastatic BRAFV600E-PTC cells which are associated with intrinsic resistance to vemurafenib treatment. Collectively, our PTC pre-clinical model suggests that combination of anti-BRAFWT/V600E therapy (e.g. vemurafenib) with inhibitors of pro-survival molecules (i.e. pan-BCL2/MCL1 inhibitors) ameliorates intrinsic resistance to metastatic BRAFV600E-PTC cell death.

Results

Vemurafenib therapy impairs viability of human non metastatic PTC cells in a pre-clinical model of patient-derived PTC harboring BRAFWT/V600E

We have developed the first pre-clinical model of patient-derived PTC with BRAFWT/V600E (Figure 1A) using BRAFWT/V600E inhibitors (i.e. vemurafenib). We established 7 short-term primary cell cultures in vitro of human PTC (which reduce the potential for changes in vitro), including 2 from non angioinvasive and non-metastatic PTC, 2 from angioinvasive PTC with no clinical evidence of neck lymph node (LN) metastasis, 2 from angioinvasive PTC with positive clinical evidence of neck LN metastasis, and, 1 short-term cell culture from a PTC mediastinal LN metastasis (local metastasis) (Suppl. Table 1, Suppl. Figures 1A-B), all from individual patients with near-total thyroidectomy with primary PTC greater than 1.1 cm (pT1b-pT3). For comparison, we also established 4 independent cultures of primary normal thyroid (NT) cells from non-goiterous, non-thyroiditc tissue. All PTC samples were sequenced and showed wild type TP53. Overall, 71.4% (5/7) of the isolated PTC cell batches harbored the heterozygous BRAFWT/V600E mutation (Figure 1B). 14.2 % (1/7) harbored the RET/PTC3 translocation without mutations in BRAF (Suppl. Figure 1E). No mutations included in our genomic sequencing panel were identified in 1 of the 7 PTC samples. Additionally, we have used KTC1 cells, a spontaneously immortalized BRAFWT/V600E cell line established from pleural effusion material (distant metastasis) from a male patient with metastatic and recurrent PTC with wild type TP53 (weak nuclear expression, Suppl. Figure 2) which showed nuclear expression of PAX8 and phospho(p)-ERK1/2 proteins (Suppl. Figure 2). We also used BCPAP cells, with homozygous BRAFV600E, established from a poorly differentiated PTC. The NT cells showed no mutations. Furthermore, in order to assess the morphological growth pattern of these cells, we established a three-dimensional (3D) cell-culture system with Matrigel (extracellular matrix (ECM) in vitro) (Figure 1C). Representative BRAFV600E-PTC cells primarily grew significantly (p=0.02) as numerous cell aggregates (Figure 1C). Furthermore, NT cells engineered to express mutant BRAFV600E grew significantly (p=0.02) as refractile cells compared to NT cells engineered with empty vector (Figure 1C). Representative NT and PTC cells expressed some thyroid differentiation biomarkers and epithelial tumor marker as indicator of tumor purity, and were negative for desmin (mesenchymal marker) (Figure 1D).

We tested the effects of vemurafenib using a dose-response in representative PTC cells with or without BRAFWT/V600E and in NT cells. Ten µM vemurafenib was an effective dose to substantially block the BRAFWT/V600E pathway, specifically reducing pERK1/2 protein expression levels by 98% (IC90) in non-metastatic BRAFWT/V600E-PTC cells as compared to the same cells treated with vehicle (control) within 24 hours and with no substantial changes in pERK1/2 protein expression levels in either PTC or NT cells with BRAFWT (Figure 1E). BRAFWT/V600E-PTC cells morphology was also substantially affected with this latter dose within 24 hours as compared to the vehicle-treated cells (Figure 1F). By contrast, metastatic BRAFWT/V600E-PTC cells were substantially less responsive and more resistant (IC90 for pERK1/2 protein levels was about 30%) to vemurafenib treatment within 24 hours (Figure 1E). The IC50 value (dose-response) for cell viability was 5 µM for vemurafenib in representative non-metastatic BRAFWT/V600E-PTC cells at 48 hours in vitro (Figure 1G) and was significant (p=0.001) as compared to vehicle-treated cells, with no effect on the viability of BRAFWT-PTC or NT cells (Figure 1G). In contrast, 10 µM vemurafenib was an effective (p=0.007) dose able to reduce viability of BRAFWT/V600E-positive metastatic PTC or ATC cells to about 50% vs. vehicle-treated cells (Figure 1G), with no significant effect on the viability of BRAFWT-PTC or NT cells (Figure 1G). Doses greater than 1 µM (e.g. 2 or 10 µM) have also been reported to be effective in BRAFV600E-melanoma or ATC cells in vitro [18] [20]. Moreover, because our primary PTC cells grew as cell aggregates (e.g. spheroids) in culture on the Matrigel, we also investigated the expression of stem-cell markers in PTC and NT cells in vitro (Figure 1H). Interestingly, we found that a sub-population of primary human BRAFWT/V600E-PTC cells expressed substantially about 2.5-fold more PDGFRB (tumor microenvironment-associated pro-angiogenic factor) protein vs. PTC cells with BRAFWT (Figure 1H).

Effects of anti-BRAFV600E therapy in vivo

We have developed the first early intervention pre-clinical mouse trial of human PTC-derived cells with BRAFWT/V600E using vemurafenib (Figure 2A-2B). Immunocompromised mice were orthotopically implanted with the human KTC1 cells derived from a metastatic/recurrent BRAFV600E-positive PTC. These cells were engineered to express green fluorescent protein (GFP). KTC1 orthotopic tumors developed in all mice and were analyzed at 18 days post-tumor cell injection (Figure 2C-2D). Six days after tumor cell injection, mice were randomized for treatment with vemurafenib (early intervention trial). Vemurafenib treatment (tumor size mean=0.22 mm3, SD=0.2, p=0.01) by oral gavage for 12 days resulted in a significant reduction in tumor growth (>90%), compared with control mice treated with vehicle (tumor size mean=4 mm3, SD=1, p=0.01). Also, substantially decreased GFP signal suggested that this reduction was due to suppression of BRAFV600E-positive tumor cell viability and tumor growth (Figure 2C). Importantly, we found that vemurafenib therapy substantially down-regulated nuclear protein levels of cyclin B2, from about 30% of tumor cells in the vehicle-treated tumors to about 2% in vemurafenib-treated tumors (Figure 2F, Supplementary information). Although vemurafenib treatment did not significantly affect protein levels of cleaved caspase-3, marker of apoptosis (Figure 2G), we found a significant (p=0.007) increase of about 2.5-fold human lysosomal endogenous beta-galactosidase (β-gal, senescence marker) expression in the orthotopic tumor cells in the vemurafenib-treated mice (Suppl. Figure 3B) vs. vehicle-treated mice (Suppl. Figure 3A). Because BRAFWT/V600E-positive tumors tend to be more adhesive through ECM in the tumor microenvironment, we performed trichrome staining which showed a substantially abundant amount of the ECM collagen deposition surrounding the BRAFWT/V600E-KTC1 orthotopic tumor, and also intratumoral deposition of collagen; more importantly, vemurafenib anti-BRAFV600E therapy substantially decreased collagen deposition in the orthtotopic tumors, reducing ECM density (Figure 2H). Furthermore, vemurafenib therapy significantly (p=0.007) reduced 2-fold vascular density (Figure 2E) and tumor microenvironment-associated PDGFRB protein levels in the BRAFWT/V600E-KTC1 tumor vs. vehicle-treated mice (Figure 2I, and Supplementary information).

BRAFWT/V600E-PTC cells recruit microvascular endothelial cells and pericytes by regulating pro-angiogenic/metastatic paracrine signaling

We sought to test the hypothesis that BRAFWT/V600E by hyper-phosphorylation of the ERK1/2 triggers PTC lympho-angiogenesis by means of recruitment of human blood and lymphatic microvascular endothelial cells (BEC and LEC, respectively) (Suppl. Figures 4A-4B) and pericyte (Suppl. Figures 4C-4D), which are fundamental cell populations in the tumor microenvironment. We developed a trans-well endothelial cell migration in vitro assay based on PTC- or NT-derived secretome (Suppl. Figure 4A-4B) which revealed that BRAFWT/V600E-PTC cells promoted pericytes (Suppl. Figure 4C, Suppl. information) and endothelial cell migration (Suppl. Figure 4E, Suppl. Figure 4H, Suppl. information), and angiogenesis (i.e. in vitro tubule formation) (Figure 3A, Suppl. Figure 4G, Suppl. information), suggesting activation of potential pro-metastatic paracrine signaling. Tubule formation decreased (p=0.02) 1.5-3.3 fold in the presence of secretome derived from 10 µM vemurafenib-treated metastatic/recurrent BRAFWT/V600E-PTC cells compared to vehicle-treated cells; whereas, it was strongly inhibited from about 4-fold up to 8.6-fold in the presence of vemurafenib-treated secretome derived from non-metastatic BRAFWT/V600E-PTC cells (Figure 3A, Suppl. Figure 4G). LEC tubule formation strongly decreased (9.7-fold) in the presence of 10 µM vemurafenib-treated secretome derived from BRAFWT/V600E-KTC1 cells compared to vehicle-treated secretome derived from KTC1 cells (Suppl. Figure 4G). No statistically significant changes in tubules formation were observed in LEC or BEC grown in the presence of secretome derived from either BRAFWT-PTC cells or NT cells (Figure 3A, Suppl. Figure 4F). Accordingly, we found that some pro-angiogenic and ECM remodeling genes showed moderate or high mRNA expression levels in BRAFV600E-PTC cells compared with BRAFWT-PTC or NT cells (Suppl. Table 4). Vemurafenib treatment reduced their mRNA expression levels compared with BRAFV600E-PTC vehicle-treated cells (p<0.05) (Suppl. information, Suppl. Table 4).

We also used a multiplex ELISA assay that included the most known pro-angiogenic and anti-angiogenic factors. We found that LN metastatic/recurrent BRAFWT/V600E-PTC cells secreted significantly higher protein levels of pro-angiogenic factors (e.g. Leptin, ~2-fold) in vitro compared with BRAFWT-PTC cells (p=0.01) or NT cells (p=0.01) (Figure 3B). Mechanistically, we found that the rate of tubule suppression by the human Leptin blocking antibody was even more potent (4.3-fold decrease, p=0.02) when BEC were grown in the presence of 10 µM vemurafenib-treated secretome derived from primary metastatic/recurrent BRAFWT/V600E-PTC cells (Suppl. Figures 4J), suggesting that inhibition of both BRAFV600E and Leptin pathways was more effective to reduce angiogenesis in vitro. Tumor cell migration/invasion is functionally linked to angiogenic paracrine communication between tumor cells and stromal cells; indeed, we also assessed microenvironment-associated stromal cell types (Figure 3C) which participate in the aberrant behavior of metastatic BRAFV600E-PTC cells using a co-localization approach for PDGFRB and PAX8 (thyrocyte cell marker) proteins; we found that PDGFRB protein levels were significantly higher (2.2-fold increase, p=0.02) in metastatic BRAFV600E-PTC-associated stromal cells (i.e. mesenchymal stem cell–like pericytes) compared to BRAFWT-PTC cells in clinical samples (Figure 3C and Supplementary information). Remarkably, vemurafenib therapy significantly (p=0.01) reduced cell migration and invasiveness by approximately 3.5-fold in non-metastatic BRAFWT/V600E-PTC cells in vitro as compared with vehicle-treated (control) cells (Figure 3D). In contrast, metastatic/recurrent BRAFWT/V600E-PTC cells (Figure 3D and Suppl. Figure 6A-6B, Supplementary information) were less responsive to vemurafenib treatment, with 1.5-fold or 2-fold decrease in migration or invasion. Treatment with vemurafenib did not affect migration or invasion of BRAFWT-PTC and NT cells.

Alterations of MCL1 and P16 somatic copy number in human PTC samples and PTC cell cultures

We have used a new algorithm (see methods) for detecting somatic mutations, insertions, deletions, copy number gain (amplifications), copy number loss, and translocations using a targeted exome sequencing strategy. We found significantly (p=0.0001) that chromosome 1q (Figure 4A-4B, Suppl. Figure 7) showed SCNAs (i.e. copy number gain) of 26 genes (Suppl. Table 3) and chromosome 9p21.3 showed loss of the P16 gene in PTC (Suppl. Table 2, Figure 4G, Suppl. Figure 7). Here, we have focused our attention and validations on MCL1 (a BCL2 protein family member) (one of the 26 genes with copy number gains on chromosome 1), because it is a well-known key-player in cancer progression and metastasis [21], but has not been linked to metastatic BRAFV600E-PTC. To validate this finding, we used quantitative PCR (Figure 4C-4D) using 58 independent samples (including 22 patient-derived independent metastatic PTC clinical samples, Suppl. Table 2, Figure 4C), in particular: 31 PTC patients specimens (23 BRAFV600E-PTC and 8 BRAFWT-PTC), 17 LN metastatic PTC (9 BRAFV600E-LN and 8 BRAFWT-LN), 5 distant metastatic samples from patients with PTC (1 from lungs, 1 from bone and 3 from adrenal glands), and 5 NT (Figure 4C). One DNA control sample from healthy man (human male genomic DNA) was also used (Figure 4C). Additionally, we analyzed genomic DNA from 4 primary PTC cell cultures (3 PTC with BRAFV600E and 1 PTC with BRAFWT), 1 primary LN metastatic BRAFV600E-PTC cell culture, 1 primary NT cell culture and 5 spontaneously immortalized human thyroid cancer cells (i.e. TPC1, BCPAP, KTC1, 8505c and SW1736) (Figure 4C). We found significantly that 34 of 41 (82.9%) BRAFV600E-positive samples vs. 10 of 17 (58.8%) BRAFWT samples (including primary PTC, LN metastatic PTC-derived specimens and primary short-term PTC cell cultures) harbored MCL1 copy number gain (p<0.001), and also as compared to NT tissue samples or healthy man DNA control sample which did not show any MCL1 copy number gain (Figure 4C-4D). Eleven of 12 (91.6%) of the primary BRAFV600E-PTC specimens with a clinical story of LN metastasis and 9 of 9 (100%) BRAFV600E-PTC LN metastatic samples in the neck or mediastinum substantially harbored MCL1 copy number gain. Lesser significant prevalence of MCL1 copy number gain/amplification was found in BRAFWT-PTC LN metastatic samples (5 of 8, 62.5%) as compared to the BRAFV600E-PTC LN metastatic samples (p<0.01). Also, 4 of 5 (80%) PTC samples with distant metastasis showed MCL1 copy number gain (Figure 4C-4D); whereas, NT samples did not show MCL1 copy number gain.

As further validation, high copy number of MCL1 gene also substantially associated with higher MCL1 mRNA expression levels (about 1.5-2.3-fold) in metastatic BRAFV600E-PTC cell cultures (98 mRNA copies/106 18S copies) vs. primary BRAFWT-PTC (64 mRNA copies/106 18S copies) or NT cells (42.7 mRNA copies/106 18S copies) in vitro. MCL1 expression levels showed 104 mRNA copies/106 18S copies by vemurafenib treatment in BRAFV600E-PTC cells and 98 mRNA copies/106 18S copies in the vehicle-treated cells. There was no substantial effect of vemurafenib on MCL1 mRNA expression in BRAFWT-PTC or NT cells. Additionally, immunohistochemistry (IHC) analysis showed substantially strong and diffuse MCL1 protein immunostaining in all neck LN metastatic BRAFV600E-PTC samples vs. primary BRAFV600E-PTC samples, BRAFWT-PTC with or NT samples (Figure 4E), further suggesting that MCL1 SCNAs could be important for metastasis in BRAFV600E-PTC (Figure 4F).

Furthermore, we validated P16 SCNAs (Figure 4G, 4H and 4I) using: (i) 1 control DNA sample from healthy man; (ii) 27 PTC (19 BRAFV600E-PTC and 8 BRAFWT-PTC), 15 LN metastatic PTC (7 BRAFV600E-LN and 8 BRAFWT-LN), 6 distant metastatic samples from patients with PTC (2 from lungs, 1 from bone and 3 from adrenal glands), and 4 NT (Figure 4I) tissue samples; and (iii) 2 primary PTC cell cultures (2 PTC with BRAFV600E and 1 PTC with BRAFWT), 1 primary LN metastatic BRAFV600E-PTC cell culture, and 4 spontaneously immortalized human thyroid carcinoma cells (i.e. TPC1, BCPAP, KTC1, and 8505c). These samples included 21 patient-derived independent metastatic PTC (Suppl. Table 2, Figure 4H). Our validations significantly showed that 33 of 34 (97%) BRAFV600E-positive samples vs. 14 of 18 (77.7%) BRAFWT samples (including primary PTC, LN metastatic PTC-derived specimens and primary short-term PTC cell cultures) loss of the P16 gene (p<0.001), and also as compared to NT tissue samples or healthy man DNA control sample that did not show any loss of P16 (Figure 4G-4I, Suppl. Figure 7). Eight of 9 (88.8%) of the primary BRAFV600E-PTC specimens with a clinical story of LN metastasis, 7 of 7 (100%) BRAFV600E-PTC LN metastatic samples in the neck or mediastinum, and 7 of 7 (100%) of the BRAFWT-PTC LN metastatic samples substantially showed loss of P16. Six of 6 (100%) PTC samples with distant metastasis showed P16 SCNAs compared to NT samples (Figure 4H-4I). In particular, we found P16 heterozygous loss occurred in 13 of 60 (21.6%), and homozygous loss in 39 of 60 (65%) samples analyzed, specifically in: 13 of 19 (68.4%) BRAFV600E-PTC, 6 of 8 (75%) BRAFWT-PTC, 6 of 7 (85.7%) BRAFV600E-PTC with neck LN metastasis, 5 of 8 (62.5%) BRAFWT-PTC samples with neck LN metastasis and 4 of 6 (66.6%) samples with distant metastasis (Suppl. Table 2). Furthermore, IHC analysis showed a substantial loss of P16 protein in the metastatic BRAFV600E-PTC derived KTC1 cells (Figure 4J) vs. primary human NT cells, suggesting that this SCNA may be important for metastatic potential of BRAFV600E-PTC (Figure 4K).

Combination of vemurafenib therapy with the pan-BCL2 inhibitor (obatoclax) increases death of angio-invasive BRAFV600E-PTC cells harboring MCL1 copy number alterations and decreases intrinsic resistance to vemurafenib

Although the lack of response to vemurafenib has been frequently observed in different human cancer cells lines (e.g. melanoma cells) harboring the BRAFV600E mutation, it is still unknown in human metastatic PTC cells. To address this unmet clinical need, we used patient-derived BRAFV600E-positive angioinvasive PTC cells. Importantly, we found that angioinvasive BRAFWT/V600E-PTC cells (i.e. PTC7 cells) with MCL1 copy number gain (Figure 4C) were more resistant to vemurafenib-induced inhibition of cell viability as compared to non-metastatic BRAFWT/V600E-PTC cells with MCL1 neutral copy (Figure 1E-1G). Mechanistically, based on the vemurafenib IC50 calculations (Figure 1G) we found that reduced doses of 2.5 µM vemurafenib or even more effectively 5 µM vemurafenib combined with 100 nM obatoclax significantly ameliorated therapeutic response of angio-invasive BRAFV600E-PTC cells compared to single agent therapy within 48 hours (Figure 5A). This combined treatment resulted in a significant increased cell death (sub-G1 cell accumulation, indicator of apoptosis) detected by flow cytometry analysis by using single agent vemurafenib (13.5-fold) or obatoclax (4.5-fold) treatment compared to vehicle (Figure 5A). More importantly, vemurafenib plus obatoclax combined treatments significantly increased number of BRAFV600E-PTC cells in sub-G1 phase (p<0.05) as compared to single agents treatments (2.2-fold to vemurafenib and 6.8-fold to obatoclax) or to vehicle treatment (30.7-fold) (Figure 5A). Importantly, this combinatorial targeted treatment was also effective to alter BRAFV600E-PTC cell shape (cells became rounded up and detached) as compared to vehicle-treated cells (Figure 5B).

Intriguingly, metastatic BRAFV600E-PTC derived KTC1 cells (Figure 5C) with a lower copy number of MCL1 (+/-) (Figure 4C) were more sensitive to single agents treatments or more effectively to vemurafenib plus obatoclax combined treatments compared to BRAFV600E-PTC cells with a higher copy number of MCL1 (Figure 5C); as result, these cells significantly increased in the sub-G1 phase (apoptosis) (Figure 5C) (p<0.05). In particular, they accumulated in sub-G1 phase 24-fold by vemurafenib (10 µM), 17.6-fold by obatoclax (100 nM), and 97.1-fold by vemurafenib plus obatoclax treatments compared to vehicle (Figure 5C). Furthermore, BRAFV600E-PTC KTC1 cells with MCL1 (+/-) (Figure 5D) showed more substantial alterations in cell shape as compared to metastatic BRAFV600E-PTC cells with a higher copy number of MCL1 (Figure 5B).

Discussion

The pre-clinical trial model for cancer treatment may very well facilitate our ability to combat human cancer [1, 7, 22-27]. Here, we report the first patient-derived pre-clinical model of PTC with BRAFV600E or BRAFWT, or NT using anti-BRAFV600E therapy. Although primary NT or PTC cell survival in culture is limited, the benefit is that primary cells show low likelihood of adopting cellular and molecular changes in vitro.

Our results raise the concern that primary/intrinsic resistance to BRAFV600E inhibitors (i.e. vemurafenib) could occur at or near the acquisition of metastatic capability by BRAFV600E-PTC. Critically, our experimental models identify that genomic co-occurrence of BRAFV600E mutation with either MCL1 copy number gain (chromosome 1q) or P16 (chromosome 9p) loss is functionally associated with metastasis in patients with PTC and patients should be closely monitored for disease recurrence. This paired genomic abnormality could confer primary/intrinsic resistance to vemurafenib therapy, possibly through activation of BRAFV600E parallel signaling pathways regulated by either MCL1 or P16 in metastatic PTC cells. Furthermore, abundance of MCL1 (pro-survival factor) or loss of P16 (tumor suppressor) protein and mRNA levels functionally correlated with their genomic abnormalities. MCL1 is an anti-apoptotic member of the BCL2 family and is known to be focally amplified in about 10.9% of cancers across multiple tissue types (e.g. breast cancer, lung cancer, etc.) and to increase tumor cell survival [21]. MCL1 is also amplified in melanoma cells [28]. Our orthotopic mouse model of tumor harboring the heterozygous BRAFV600E mutation along with MCL1 and P16 copy number alterations showed focal PTC-like nuclear features (i.e. clearing), suggesting that this mouse model may be a helpful translational tool for in vivo pre-clinical testing of anti-BRAFV600E therapy, although these tumor cells (i.e. KTC1 cells) showed in vivo mixed spindled and epithelioid features which are not characteristic of PTC. Vemurafenib therapy did not elicit expression of pro-apoptotic markers in the orthotopic tumor cells; interestingly, tumor cell senescence increased with concomitant tumor size reduction via suppressed pro-survival molecules. Moreover, our analysis of tumor microenvironment both in our orthotopic tumors and in PTC clinical samples provides evidence that metastatic BRAFV600E-PTC with MCL1 and P16 copy number alterations are enriched with PDGFRB-positive stromal cells (i.e. pericyte). Additionally, BRAFV600E-positive PTC cell might synergize with pericytes and sustain the viability via expression and secretion of factors crucial for angioinvasion and angiogenesis. Furthermore, our in vitro data indicate that BRAFV600E might regulate PTC cell migration/invasion and viability, and its inhibition by vemurafenib may induce cell cycle arrest and apoptosis in non-metastatic PTC cells. In contrast, metastatic BRAFV600E-PTC cells with MCL1 copy number gain are substantially less sensitive and responsive to vemurafenib treatment as a single agent. Consistent with our genetic data, we critically show that combined therapy with vemurafenib and pan-BCL2/MCL1 inhibitor decreased viability of angioinvasive BRAFV600E-PTC cells in vitro. This effect could be determined by the inhibition of BCL2 and MCL1 (1q) which are known to inhibit apoptosis and maintain cancer cell survival. Co-treatment with BRAF inhibitors and with inhibitors of BCL2 might overcome intrinsic resistance to BRAF inhibitors (also seen in melanoma cells in vitro and in vivo) [28]. BCL2 inhibitors synergize with MEK inhibitors and induce in vitro death of murine poorly differentiated thyroid carcinoma cells [29]. Moreover, in lung cancer, reduced MCL1 expression sensitized epidermal growth factor receptor mutant non-small cell lung cancers to MEK inhibitors [30]. BCL2 family inhibitors, in combination with a TORC1/2 inhibitor, lead to apoptosis in KRAS- and BRAF-mutants but not wild-type (WT) colorectal cancer cells. This is accompanied by suppression of MCL1 expression in mutant, but not WT, colorectal cancers [31]. However, BCL2 inhibitors may be toxic, and therefore, ascertaining the correct dosage is essential.

Importantly, the recent PTC TCGA data show chromosome 1q amplification in aggressive BRAFV600E-positive PTC [5], suggesting key role for this genomic region in PTC behavior. PTC heterogeneity might affect the sensitivity for detection of copy number alterations (MCL1, P16, etc.) in some cases, with possible epigenetic alterations of MCL1 and P16 expression levels; therefore, we cannot exclude the possibility that the percentage of copy number alterations for these genes might have some variability. However, our results were substantially validated by our integrated in vitro functional assays and in vivo approaches.

Molecular screening for BRAF mutations, as well as MCL1 and P16 copy number alterations might help identify subsets of patients with PTC which have high risk for recurrence and metastasis. They may benefit from combined targeted therapies directed against these alterations. This approach might provide a novel advancement in the therapeutic strategy of aggressive and refractory BRAFV600E-PTC for cases in which surgery is limited by tumor location and extent of disease, or in medically poor surgical candidates. Short-term PTC and NT cell cultures appear to provide a novel predictive model for evaluating therapeutic responses in patients with PTC. Collectively, results from our PTC pre-clinical model provide data that have value beyond predicting the design and outcome of concurrent metastatic BRAFV600E-PTC clinical trials which are using BRAFV600E inhibitors (e.g. vemurafenib therapy) and also identify potentially biomarkers (e.g. MCL1 copy number gain, loss of P16) for aggressive PTC.

Materials and Methods

Histology and immunohistochemistry

Histopathology evaluation of 7 classical type PTC (4 PTC with BRAFV600E and 3 PTC with BRAFWT), 3 neck LN metastatic PTC samples harboring BRAFV600E (n=3) and 3 NT tissues was performed by pathologists (P.S. and V.N.) on hematoxylin and eosin (HE)-stained formalin-fixed paraffin-embedded (FFPE) tissues. For all patients we used discarded and unidentified tissues (Beth Israel Deaconess Medical Center (BIDMC) IRB-approved exemption 4 protocol, Boston). All tissue specimens were fixed with 10% buffered formalin phosphate and embedded in paraffin blocks. These were visualized with an Olympus BX41 microscope and the Olympus Q COLOR 5 photo camera (Olympus Corp., Lake Success, NY, USA). Four µm sections of PTC or BRAFV600E-positive human KTC1 tumor-derived orthotopic mouse tissues (serial sections) were used for IHC. After baking overnight at +37°C, deparaffinization with xylene/ethanol and rehydration was performed. IHC analysis utilized primary antibodies (Supplementary Table 5). The sections, treated by pressure cooker for antigen retrieval (Biocare Medical, Concord, CA), were incubated at 123°C in citrate buffer (Dako Target Retrieval Solution, S1699; DAKO Corp.), cooled and washed with PBS. Antigen retrieval was performed for 60 min at room temperature. The primary antibody was detected using a biotin-free secondary antibody (K4011) (Dako Envision system) and incubated for 30 min. All incubations were carried out in a humid chamber at room temperature. Slides were rinsed with PBS between incubations. Sections were developed using 3,3-diaminobenzidine (Sigma Chemical Co.) as a substrate and were counterstained with Mayer’s hematoxylin [1] [32] [33]. The IHC marker expression was assessed semiquantitatively using the following scoring method: 0 (negative), 1+ (1–10% positive cells), 2+ (11–50% positive cells), and 3+ (more than 50% positive cells). Microvascular density is defined by number of vessels per microscope field showing CD31 (marker of blood endothelial cells) or D2-40 (Podoplanin) (marker of lymphatic endothelial cells) staining. Trichrome staining was performed according to manufacturer instructions (Ventana, USA) in order to assess collagen deposition in the ECM of PTC.

Genotyping oncomap analysis and mass extend sequenom

Genomic DNA was extracted using our previous protocol [21] from overall 13 samples from independent patients, which include: (i) 2 BRAFV600E-PTC, 1 BRAFWT-PTC, 1 mediastinal LN metastatic PTC and 1 NT available FFPE discarded/unidentified samples from 4 independent patients (Beth Israel Deaconess Medical Center (BIDMC) IRB-approved exemption 4 protocol, Boston); (ii) 2 primary PTC cell cultures with BRAFV600E, 1 mediastinal LN metastatic PTC cell culture with BRAFV600E, 1 primary PTC cell culture with BRAFWT and 1 primary NT cell culture; and (iii) 3 spontaneously immortalized human thyroid cancer cells (8505c, BCPAP and KTC1, see cell cultures paragraph for more details). Briefly, 30 µm paraffin sections were lysed in denaturing buffer containing proteinase K (1 mg/ml) (Invitrogen, USA) during overnight incubation at +55°C. DNA was purified using equal volumes of a phenol:chloroform mixture (Invitrogen, USA) and eluted in distilled water. Genomic DNA was quantified using Quant-iT PicoGreen dsDNA Assay Kit (Invitrogen, USA) per manufacturer protocol. 250 ng of genomic DNA was used for mutation detection through Oncomap version 3, which interrogates about 1000 known mutations (including the BRAFV600E mutation) in 112 validated oncogenes and tumor suppressors [34]. All genomic analyses were run in the Center for Cancer Genome Discovery (CCGD, Dana Farber Cancer Institute (DFCI), Harvard Medical School, Boston, MA, USA) [34]. Primers and probes were designed using Sequenom MassARRAY Assay Design 3.0 software, applying default multi-base extension parameters. Whole genome amplification (WGA) was performed using the GenomePlex Complete WGA kit (Sigma, USA) based on chemical fragmentation followed by adapter mediated PCR amplification. Samples were run on the mass spectrometry-based genotyping platform (Sequenom) and analyzed according to current standardized protocols [34]. Sample identity and the possible introduction of artifacts by WGA were evaluated using a 48 Single Nucelotide Polymorphisms (SNPs) panel comparing the pre-WGA to the post-WGA DNA. If ≥3 SNP discrepancies were identified between SNPs found in pre- and post-WGA samples, this sample was discarded. Validations were performed using homogeneous and sensitive Mass Extend (hME) sequenom, a multi-base extension chemistry performed on native unamplified genomic DNA [34]. For all patients, we used discarded/unidentified tissue specimens and consent for genotyping test.

Sanger sequencing and SNaPshot mutational analysis

PCR products were generated and sequenced for the BRAFV600E mutation screening according to previous protocols [35, 36]. SNaPshot analysis covered a broad panel of mutations including the BRAFV600E and was performed according to Dias-Santagata et al. [37].

Targeted exome sequencing

OncoPanel version 2 (OPv2) represents a recent targeted exome sequencing strategy (developed from the CCGD, DFCI, Harvard Medical School) to simultaneously detect mutations, insertions, deletions, translocations and SCNAs. We analyzed 13 samples from independent patients, in particular from: 2 BRAFV600E-PTC, 1 BRAFWT-PTC, 1 mediastinal lymph node (LN) PTC metastasis and 1 NT available FFPE discarded/unidentified samples from 4 independent patients (Beth Israel Deaconess Medical Center (BIDMC) IRB-approved exemption 4 protocol, Boston); as well as we analyzed genomic DNA from 2 primary PTC cell cultures with BRAFV600E, from 1 mediastinal LN metastatic PTC cell culture with BRAFV600E from 1 primary PTC cell culture with BRAFWT and from 1 NT cell culture. Additionally, we analyzed DNA from 3 spontaneously immortalized human thyroid cancer cells (8505c, BCPAP and KTC1 see cell cultures paragraph for more details). Targeted sequencing was achieved by designing RNA baits to capture the exons of 504 genes with relevance to cancer. The bait set was augmented with specific intronic sequences to also detect translocations often involved in cancer. Sequencing libraries were prepared, as previously described in Hettmer et al. [38], starting from 100 ng of genomic DNA. Libraries were quantified by qPCR (Kapa Biosystems, Inc, Woburn, MA) and pooled in equimolar concentrations to 500 ng total and enriched for the OPv2 bait set using the Agilent SureSelect hybrid capture kit. The enriched targeted exon libraries were again quantified by qPCR (Kapa Biosystems, Inc, Woburn, MA) subsequently sequenced in one lane of a Hiseq2000 sequencer (Illumina Inc, San Diego, CA) in a 2× 100 bp pair-end mode. Sequence alignment, demultiplexing and variant calling, including SCNAs and Indels (insertions-deletions) were performed using PICARD (HTSJDK Java library HTSJDK, BROAD Institute, Cambridge, MA, USA), GATK tools, Mutect and IndeLocator as previously described in Hettmer et al. [38]. Common variants were filtered against the 6,500 exome release of the Exome Sequencing Project (ESP) database. Variants represented in either the African-American or European-American and not in COSMIC > 2× were considered to be germline. The ESP filter was only applied to tumor samples that did not have a matched normal DNA sample.

Somatic copy number alteration analysis

Somatic copy number alterations (SCNAs) analysis was performed using Nexus7.1 (BioDiscovery Inc., USA) after calculating the sequencing coverage using PICARD. Coverages were normalized over GC-content using a loess regression and DNA derived from normal and immortalized cultured lympho-blastoid cells by the CEPH (Centre d’Etude du Polymorphisme Humain; http://www.cephb.fr/) as reference. SCNAs were called with the following next generation sequencing (NGS) settings: significant threshold = 1E-6; max contiguous probe spacing or 1000 Kbp; minimum number of probes per segment = 3. SCNA gains has a log ration >0.4 and were called high gain if >0.8. Single copy loss threshold was -0.6 and big loss was -1.2. SCNA’s on chromosomes X and Y were not called.

Translocations analysis

Translocations (rearrangements) were detected by Breakmer, an algorithm that was specifically developed for targeted sequencing by the CCGD (DFCI, Harvard Medical School). Briefly, the method identified the reads that partially aligned to the targeted regions captured and sequenced. The unaligned portion of these ‘soft-clipped’ reads was used to assemble contigs. The contigs were aligned to the reference genome, and structural variants were called based on the alignment. For translocation events, discordantly mapped paired-end reads were identified which flank putative breakpoints.

Validations of MCL1 and P16 SCNAs

Genomic DNA was extracted using our previous protocol [21] from either 58 (for MCL1 copy number assays, see below) or 51 (for P16 copy number assays, see below) FFPE archived clinical samples; in details, for MCL1 copy number assays we used: 31 PTC (23 BRAFV600E-PTC and 8 BRAFWT-PTC), 17 LN metastatic PTC (9 BRAFV600E-LN and 8 BRAFWT-LN), 5 distant metastatic samples from patients with PTC (1 from lungs, 1 from bone and 3 from adrenal glands), and 5 NT; and 1 control sample from healthy man (human male genomic DNA from Promega, USA). For P16 copy number assays we used: 27 PTC (19 BRAFV600E-PTC and 8 BRAFWT-PTC), 15 LN metastatic PTC (7 BRAFV600E-LN and 8 BRAFWT-LN), 6 distant metastatic samples from patients with PTC (2 from lungs, 1 from bone and 3 from adrenal glands), and 4 NT; and 1 control sample from healthy man. Moreover, our SCNAs analyses also included genomic DNA from 4 primary PTC cell cultures (3 PTC with BRAFV600E and 1 PTC with BRAFWT), 1 primary LN metastatic BRAFV600E-PTC cell culture, 1 primary NT cell culture and 5 spontaneously immortalized human thyroid cancer cells (i.e. TPC1, BCPAP, KTC1, 8505c and SW1736). MCL1 and P16 SCNAs analysis was respectively performed in 22 or 21 patient-derived independent metastatic PTC clinical samples (Suppl. Table 2). DNA quality (260/280 ratio) and concentration (ng/µL) were assessed by NanoDrop (ND-1000 Spectrophotometer, Thermo Fischer Scientific Inc, MA, USA); or assessed by Qubit 2.0 fluorometer (Invitrogen-Life Tech., Inc., USA) using the double-stranded DNA broad range assay kit (Invitrogen-Life Tech., Inc., USA). Samples were then diluted in order to generate uniform DNA concentrations (5 ng/µL or 20 ng/μL per replicate/reaction). For the duplex real-time PCR reaction, we used a TaqMan based copy number assay specific for: (i) MCL1 (myeloid cell leukemia sequence 1, gene aliases: BCL2L3, EAT, MCL1-ES, MCL1L, MCL1S, Mcl-1, TM, bcl2-L-3, mcl1/EAT) (NCBI location: Chr.1:150547027-150552214; assay gene location: exon 1; cytoband: 1q21.3; assay reference genome Location: Chr.1:150552073 on NCBI build 37) (Life Technologies, USA, assay ID: Hs01326481_cn, cat#4400291) in 58 samples and (ii) P16 (cyclin-dependent kinase inhibitor 2A or CDKN2A, gene aliases: ARF, CDK4I, CDKN2, CMM2, INK4, INK4A, MLM, MTS-1, MTS1, P14, P14ARF, P16, P16-INK4A, P16INK4, P16INK4A, P19, P19ARF, TP16) (NCBI location: Chr.9:21967751-21994490; assay gene location: exon 2; cytoband: 9p21.3; assay reference genome Location: Chr.9:21974968 on NCBI build 37) (Life Technologies, USA, assay ID: Hs00237642_cn, cat#4400291) in 51 samples for which there was sufficient DNA available. We also used the copy number reference assay human RNAase-P (Life Technologies, USA, cat# 4403326). Briefly, each 5 µL reaction mixture was run in a 384-well plate with 4 replicates per sample each containing 2.5 µL Universal Master Mix II, no UNG, (catalog #4440040, Life Technologies, USA), 0.25 µL copy number reference assay RNAase-P, 0.25 µL copy number assay for MCL-1 or CDKN2A, and 2 µL DNA sample at 2.5 ng/µL. GAPDH (glyceraldehyde-3-phosphate dehydrogenase) was utilized as an additional reference gene. Each 5 µL reaction mixture contained 2.5 µL Universal Master Mix II, no UNG, 0.1 µL of 10 µM forward primer (5’-CCTGACCTGCCGTCTAGAAAA-3’), 0.1 µL of 10 µM reverse primer (5’-CTCCGACGCCTGCTTCAC-3’) according to Li et al. [39], 0.1 µL of 10 µM GAPDH probe (5’-FAM-CTGCCAAATATGATGACATCAAGA-BHQ-3’) (Integrated DNA Technologies, USA), 0.2 µL DNAse/RNAse free water and 2 µL DNA sample at 2.5 ng/µL. The target and reference copy number assays were run in a duplex on Bio-Rad Real-Time PCR System (CFX Connect™ Real-Time PCR Detection System, BioRad, USA) or using the 7900HT Real-Time PCR system (Life Technologies, USA) with the program as 10 minutes polymerase activation at 95°C followed by 40 cycles of 95°C for 15 seconds and 60°C for 1 minute. The Real-Time PCR data were imported to CopyCallerTM Software v2.0 (Life Technologies, USA) for comparative Ct (ΔΔCt) relative quantification analysis using a control sample (Human Male Genomic DNA from Promega, USA) with assumed 2 copies (diploid DNA) of target gene in the following equation: copy number = cnc2–ΔΔCt. Cnc is the copy number of the target sequence (i.e. MCL1) in the calibrator sample (human male genomic  DNA); ΔCt is the mean difference between Ct of target assay (MCL1) and reference assay (GAPDH or RNase-P); and ΔΔCt is the difference between ΔCt for the tested sample and the calibrator sample. We used the following threshold similarly to previous studies [40] [41] for the copy number gain/amplification assay and classified as low-level, moderate-level and high-level amplification, specifically as follows: normal copy number (1.85-2.15 fold), heterozygous deletion (0.90-1.84 fold), homozygous deletions (below 0.9 fold), gain (2.16-2.39 fold), low-level amplification (2.4-3.9 fold), moderate-level amplification (4-9.8 fold) and high-level amplification (greater than 9.8 fold) according to the distributions in reference to normal samples with diploid DNA. Assays were performed at least in duplicate.

Cell cultures

Primary human PTC or NT cells were derived according to procedures used in previous studies [42] [43] from 0.1 to 1 gram excess fresh PTC tissue, or normal fresh thyroid tissue from patients undergoing thyroidectomy (IRB approved for AF). Briefly, we have adopted some changes in the protocol [42] [43] using 1:1 DMEM:Ham’s F12 (Corning, USA) with 10% fetal bovine serum (FBS). We cut tissue into 1 mm3 in a sterile petri dish with cold Hanks Balanced Salt Solution (calcium and magnesium free) (HBSS) (Fischer Scientific, MA, USA). We transferred minced tissue to a 50 ml tube with HBSS. We decanted supernatant red cells. Tissue was digested in a solution with clostridium collagenase (100 U/ml) and Dispase (1 mg/ml) (Stem Cell Technologies, USA) in HBSS or in 1:1 DMEM:Ham’s F12 growth medium with no serum (Corning, USA) for 1.5 hours at 37C in a humidified tissue culture incubator with 5% C02. We filtered the remaining completely digested tissue and pooled supernatants through a 70 micron cell strainer (BD Falcon, USA). We rinsed the strained material with 5 ml of HBSS to maximize yield. This straining step is important to help remove fibroblasts and blood vessels. We spun down the filtered cells at low speed and resuspended them in 1:1 DMEM:Ham’s F12 with 10% FBS plus ampicillin/streptomycin antibiotic and antimycotic solution 100× (Corning, USA) for seeding at about 1×105 cells per mL on plastic petri dishes. PTC cells digest faster than NT cells, and adhere relatively poorly to glass or plastic. Cell viability was generally high and cytokeratin staining confirmed >99% of adherent cells are epithelial using this method. PTC genomic status was determined by genomic sequencing (IRB approved for CN). Thyroid and epithelial cell IHC markers were validated in previous studies [23]. Human microvascular endothelial cells (human blood endothelial cells (BEC) and lymphatic endothelial cells (LEC)) were kindly provided from Dr. Harold F. Dvorak (BIDMC, Harvard Medical School, Boston, USA). Pericytes were obtained from Promo Cell (Heidelberg, Germany). Primary PTC or NT cells were grown in 1:1 DMEM:HAM’s F12 (Corning, USA) supplemented with 10% FBS and ampicillin/streptomycin antibiotic and antimycotic solution 100× (Corning, USA). BEC and LEC were grown in EGM2 bullet kit (Lonza, USA) growth medium with 10% FBS. We used the following spontaneously immortalized human thyroid cancer cells: TPC-1 (PTC, harboring RET/PTC-1 and BRAFWT) was provided by Dr. F. Frasca (University of Catania, Italy); the BCPAP cells which harbor the homozygous BRAFV600E mutation were established from the primary tumor of a 76-year-old woman with poorly differentiated PTC [44] and provided by Dr. G. Damante (University of Udine, Italy); the KTC1 cells which harbor the heterozygous BRAFWT/V600E mutation were established from the metastatic pleural effusion from a recurrent and radioiodine (RAI) refractory PTC in a 60-year-old male patient [45] by Dr. J. Kurebayashi (Department of Breast and Thyroid Surgery Kawasaki Medical School Kurashiki, Japan) and provided by Dr. Rebecca E. Schweppe (University of Colorado, USA); and the 8505c (ATC) cells harboring the hemizygous/homozygous BRAFV600E mutation were purchased from DSMZ (German collection of microorganisms and cell culture, Braunschweig, Germany) [24, 46]. These human thyroid cancer cells were validated by genomic/genotyping analyses (e.g. BRAFV600E mutation analysis, see above) and were grown in high glucose DMEM (Corning, USA) medium supplemented with 10% FBS and plus ampicillin/streptomycin antibiotic and antimycotic solution 100× (Corning, USA). All in vitro assays were performed by growing PTC or NT cells in growth medium supplemented with 0.2% FBS. The human embryonic kidney (HEK) 293T cell line was kindly provided by Dr. Laura E. Benjamin (ImClone Systems, NYC, USA) and grown in DMEM supplemented with 10% FBS and plus ampicillin/streptomycin antibiotic and antimycotic solution 100× (Corning, USA).

Three-dimensional (3D) cell cultures assays

In 3D cell cultures, growth of PTC cells resulted in formation of cell aggregates (e.g. spheroids) that recapitulate some aspects of the in vivo ECM. Using this assay in human breast carcinoma models, it was shown that oncogenes during steps of tumor progression could elicit aggressive cell phenotype [47]. 3D cell cultures were prepared as follows: 4×103 or 20×103 cells were suspended in high glucose DMEM growth media supplemented with 10% FBS and seeded on growth factor-reduced Matrigel (BD Biosciences, USA) in 24-well petri dishes. Phase-contrast images (10× and 20× magnification) were captured using a microscope (Nikon TE300, USA) equipped with a camera (Leica DFC 350 FX, USA).

Vemurafenib preparation

Vemurafenib (PLX4032, RG7204) (Roche, NYC, USA) was dissolved in absolute dimethyl sulfoxide (DMSO, vehicle) (Sigma, USA) to achieve a stock concentration of 10 mM for in vitro assays. Ten mM was diluted to 2 mM in absolute DMSO; then, vemurafenib intermediate doses were diluted in 0.2% FBS high glucose DMEM in order to achieve final concentrations of 0.01 µM, 0.1 µM, 1 µM, 2.5 µM, 5 µM, and 10 µM vemurafenib at 2% DMSO for optimal vemurafenib solubility. Vehicle was 2% DMSO diluted in 0.2% FBS high glucose DMEM. For in vivo studies, a drug suspension was prepared from micro-precipitated bulk powder (MBP) by suspending the drug to a concentration of 25 mg active pharmaceutical ingredient (API)/mL in a 2% solution of hydroxypropylcellulose (vehicle) (Sigma, USA), according to manufacturer instructions. Freshly prepared drug suspensions were stored at 4°C and used within 24-48 hours. Mice were dosed twice daily (8 hours apart) with vehicle alone (control) or with vemurafenib suspensions at 100 mg/kg using an 18G oral gavage needle.

Obatoclax mesylate

Obatoclax [29] also named GX15-070 (Selleckchem, USA) is an inhibitor of the BCL2 proteins family (which it does also include MCL1), was dissolved in absolute DMSO (Sigma, USA) according to the manufactures instructions and then used at 2 nM, 10 nM, 20 nM, 100 nM, 500 nM and 1 µM diluted in 0.2% FBS high glucose DMEM at final 2% DMSO. Combinatorial treatments of vemurafenib plus obatoclax were calculated based on the vemurafenib IC50 values using a decrease of the highest dose (i.e. 10 µM) of vemurafenib, specifically we combined 5 µM or 2.5 µM vemurafenib along with 2 nM, 10 nM, 20 nM, 100 nM, 500 nM or 1 µM obatoclax. All compounds (single agents or combined treatments) were diluted in 0.2% FBS high glucose DMEM at final 2% DMSO. Vehicle was 2% DMSO diluted in 0.2% FBS high glucose DMEM.

Cell viability

Human PTC or NT cells (1×103 cells/well) were cultured in growth medium containing 0.2% FBS (CellGro, USA) in a 96-well sterile culture plate (Thermo Fisher Scientific, USA). Cells were treated with or without various concentrations of vemurafenib for 48 hours. Cell viability was measured using the CellTiter-Glo luminescent cell viability assay kit (Promega, USA). The IC50 (50% maximal inhibitory concentration) for cell viability was determined using the following doses of vemurafenib: 0.01 µM, 0.1 µM, 1 µM, 5 µM and 10 µM.

Cell transfections for lentivurs or retrovirus production

HEK 293T cells (5 ×105) were grown in 60-mm plates and transfected using Fugene-6 (Roche) in OptiMEM (Invitrogen) for 48 hours according to manufacturer instructions.

Over-expression technique

HIV-U6-GFP GL3B lentivirus which expresses GFP (green fluorescent protein) was kindly provided from Dr. Yutaka Kawakami (Division of Cellular Signaling, Institute for Advanced Medical Research, Keio University School of Medicine, Japan). BRAFV600E-pBABE-puro and pBABE-puro (empty vector) retroviral constructs were kindly provided by Dr. W.C. Hahn (Dana-Farber Cancer Institute, Harvard Medical School, Boston) and used for BRAFV600E over-expression studies in NT cells. Either lentiviral or retroviral infections (viral transductions) were performed as follows: HEK 293T cells were seeded in 60-mm dishes and cotransfected the next day with each lentivirus or retrovirus and helper plasmids (i.e., gag-pol and VSV-G were used only for lentiviral infections and provided by Dr. W.C. Hahn). Both pBABE retroviral constructs were similarly packaged using helper plasmid pCL-ampho provided by Dr. W.C. Hahn. Media with progeny virus from transfected HEK 293T cells were collected 48 hours later and filtered with 0.45 μm filters (Millipore) and freshly used to transduce KTC1 cells for 3–6 hours in the presence of 8 μg/mL polybrene (Sigma). All cells used were transduced at a multiplicity of infection of either 50 or 100. Finally, stable transduced cells with HIV-U6-GFP GL3B were selected and sorted by flow cytometry (MoFlo/FACSAria Sorting; Beckman Coulter); and NT cells transduced with pBABE-puro retroviral constructs were treated with puromycin (1 μg/mL) (Sigma-Aldrich) for selection. BRAFV600E over-expression in the cells were confirmed by real-time RT-PCR according to Nucera et al. [1]. All assays were performed in duplicate.

ELISA

PTC or NT cells were cultured at about 80% confluence in 6-well dishes in 1:1 DMEM/Ham’s F12 or high glucose DMEM supplemented with 10% FBS and then treated with 10 μM vemurafenib or vehicle for 24 hrs in 1:1 DMEM/Ham’s F12 or high glucose DMEM supplemented with 0.2% FBS. After 24 hours, the growth medium enriched with PTC or NT cell-derived growth/angiogenic factors (secretome) was collected, removed from dead cell debris by short spin, diluted 1:3, and then secreted VEGFA, FGFβ, EGF, Leptin, TNFα, IL6, IGF-1, TGFβ (protein levels, pg/mL) were determined by ELISA (Enzyme-linked immunosorbent) assay (Signosis, CA, USA) according to manufacturer instructions. TSP-1 protein levels were also determined by ELISA assay (R&D Systems, MN, USA) according to manufacturer instructions. Growth medium was measured by ELISA to subtract background and normalize all samples. Protein levels were also normalized to total protein content (µg/µL).

Cell cycle analysis

PTC-derived cells were cultured at about 60-80% confluence in 6-well dishes in 1:1 DMEM/Ham’s F12 or high glucose DMEM supplemented with 10% FBS then treated with 10 µM vemurafenib, 100 nM obatoclax, combination of 5 µM vemurafenib plus 100 nM obatoclax, or vehicle for 48 hrs in 1:1 DMEM/Ham’s F12 or high glucose DMEM supplemented with 0.2% FBS.. After 48 hours, adherent cells were trypsinized and pelleted with the supernatant. Cell pellets were fixed in pre-chilled (−20°C) ethanol 75%, pelleted at 400 × g at room temperature, suspended in 0.5% PBS/BSA, incubated in 2M HCl 0.5% BSA for 20 min at room temperature, washed with 0.5% PBS/BSA, and centrifuged 5 minutes at 400 × g at room temperature. The cells were suspended in 0.1M sodium borate (pH 8.5) for 2 min at room temperature, washed with 0.5% PBS/BSA, and centrifuged 5 min at 400 × g at room temperature. Finally, the cells were pelleted, washed twice with 0.5% PBS/BSA, and suspended in 500 μL of this solution. Propidium iodide (Sigma) was added to a final concentration of 10 μg/mL with RNase 10 mg/mL (Sigma). Cells were incubated at room temperature for 30 minutes and then analyzed by flow cytometry on a FACS Gallios (Beckman Coulter, Miami, FL, USA).

Apoptosis assay

PTC-derived cells were cultured at about 60-80% confluence in 6-well dishes in 1:1 DMEM/Ham’s F12 or high glucose DMEM supplemented with 10% FBS then treated with 10 µM vemurafenib, 100 nM obatoclax, combination of 5 µM vemurafenib plus 100 nM obatoclax, or vehicle for 48 hrs in 1:1 DMEM/Ham’s F12 or high glucose DMEM supplemented with 0.2% FBS.. After 48 hours, adherent cells were trypsinized and pelleted with the supernatant. Cell pellets were fixed in pre-chilled (−20°C) ethanol 75% for propidium iodide staining and flow cytometry analysis a FACS Gallios (Beckman Coulter, Miami, FL, USA) to evaluate sub-G1 cell population (apoptosis).

Western blotting

Cells were cultured at about 90% confluence in 6-well dishes in 1:1 DMEM/Ham’s F12 or high glucose DMEM supplemented with 10% FBS and then treated with 10 μM vemurafenib or vehicle for 24 hrs in 1:1 DMEM/Ham’s F12 or high glucose DMEM supplemented with 0.2% FBS. After 24 hours, adherent cells were lysed along with the pelleted supernatant. Western blotting assays were performed according to a standard procedure, and the protein lysis buffer was composed of 10 mM Hepes (pH 7.40), 150 mM NaCl, 5 mM EDTA, 1 mM EGTA, 1 mM sodium vanadate, 5 mM sodium flouride, and 1% Triton X-100; protease and phosphatase inhibitors (Pierce) were used for protein extractions [1]. The intensity of the bands was quantified with a densitometer (Quantity One 1-D analysis software, Bio-Rad, USA). The quantity of signal in the vemurafenib lane was divided by the signal of the vehicle (control) lane in the corresponding tubulin or β-actin blot lane. The assessment of the 90% maximal inhibitory concentration (IC90) of ERK1/2 phosphorylation was performed by densitometry comparing phospho-ERK1/2 protein expression levels upon 10 µM vemurafenib treatment vs. vehicle (control) treatment by densitometry. We used the following antibodies: pERK1/2 (cat#9101, Cell Signaling, USA), total-ERK1/2, (cat#9102, Cell Signaling, USA), and beta-actin or tubulin (Sigma, USA).

Cell invasion and migration assays

Cell invasion assays were performed using 24-well BioCoat Matrigel invasion chambers according to manufacturer’s instructions (BD Biosciences, USA). These chambers were used to study the effect of vemurafenib or vehicle on the invasion of PTC or NT cells (8×103 or 25×103 cells/assay) for about 20 hours in culture. Invasion assay was performed in a serum free growth medium and we used 5% FBS as chemo-attractant agent. Migration assay was performed using: 3×103 or 25×103 PTC or NT cells/assay, 6×103 or about 25×103 pericytes/assay, 30-100×103 BEC or LEC/assay loaded into the migration chamber (Corning Incorporated, Corning, NY, USA) and grown for a range of 3-12 hours in culture. PTC or NT cell migration assay was performed in a serum free growth medium and we used 5% FBS as chemo-attractant agent. The secreted pro-angiogenic or anti-angiogenic factors (secretome) derived from the human PTC with or without BRAFWT/V600E or from NT cells treated with vehicle (2% DMSO) or vemurafenib (10 µM) for 24 hours was utilized as chemo-attractant conditioning medium (CM) to induce endothelial cell migration (Supplementary Figure 4A). Migrated or invaded cells were quantified in all assays as described: cells were counted (number of cells/field or number of cells/well) using a 10× or 20× objective, and four fields were chosen per well with two well per each condition.

In vitro angiogenesis assay

In vitro angiogenesis assays were performed as previously described [48] [46]. In brief, BEC or LEC (40-80×103 cells) (with starvation overnight at 1% serum) were suspended in NT or PTC-derived secretome treated with vehicle or vemurafenib (10 µM) for 24 hours and seeded on growth factors-reduced Matrigel (cat#354230, BD Biosciences, USA) in 24-well petri dishes. The secretome derived from the PTC or NT cells was utilized to induce tube (tubule)-like structures formation (Supplementary Figure 4B). After about 6-8 hours of incubation, cells were photographed. The number of tubes was counted using a 10× or 20× objective and four fields were chosen per well with two wells per each condition. Human Leptin blockage was performed using an antibody against Leptin (cat#AF398, R&D Systems, MN, USA).

Quantitative multigene profiling mRNA expression analysis

The expression of EMC and angiogenesis genes was validated by multi-gene transcriptional profiles provide a quantitative view of the expression of many genes [49]. Cells were seeded at about 60-80% confluence in 6-well dishes and grown in 1:1 DMEM/Ham’s F12 or high glucose DMEM supplemented with 0.2% FBS. After 24 hours, RNA isolation was performed by Quiagen columns (Quiagen, USA) following manufacturer protocols. Quantitative multi-gene profiling was performed by absolute quantification using real-time reverse transcriptase PCR (RT-PCR) according to Shih et al. [49]. Primer sequences used for the validation of human genes are reported in Supplementary Table 6. We classified gene expression as ‘low copy number’ if it was below 1 mRNA copy/106 18S copies, ‘moderate copy number’ if it was between 1 and 15 mRNA copies/106 18S copies, and ‘high copy number’ if it was greater than 15 mRNA copies/106 18S copies. Genes showing differences in relative values of comparisons (p values <0.05) were considered differentially expressed.

Orthotopic mouse model

All animal work was done in accordance with federal, local, and institutional guidelines at the Beth Israel Deaconess Medical Center (Boston, MA, USA). Human metastatic KTC1 tumor-derived cells harboring the heterozygous BRAFV600E mutation engineered to express GFP were cultured in 10-cm dishes and grown in high glucose DMEM medium supplemented with 10% FBS, penicillin, streptomycin, and amphotericin at 37°C with 5% CO2 atmosphere. On the day of tumor implantation, the cells were trypsinized, gently centrifuged, and suspended in serum-free high glucose DMEM growth medium to achieve a cell suspension concentration ranged between 3×106 cells and 5×106 cells/10 µL. The cells were kept on ice until implantation. KTC1 cells were orthotopically injected using a Hamilton syringe (Fischer, MA, USA) in the right thyroid of 9 weeks old male NOD SCID gamma mice (strain name: NOD.Cg-Prkdcscid Il2rgtm1Wjl/SzJ; stock number: 005557) (The Jackson Laboratory, Bar Harbor, Maine, USA) (BIDMC IRB, Boston) (n=5 per group) according to our previous experimental procedures [1, 50, 51]. Mice were randomly divided into two equal groups of 5 mice for the purpose of establishing a time course of tumor cell growth and response to the early therapeutic intervention with vemurafenib. Vemurafenib treatment was started 6 days (early intervention trial) after KTC1 tumor cells implantation and performed for about 12 days.

Apoptosis assay in vivo

Apoptosis assay was performed using terminal deoxynucleotidyl transferase dUTP nick end labeling (TUNEL) assay according to manufacturer recommendation (ApopTag peroxidase in situ detection kit) (Millipore Corp., Bedford, MA); apoptotic cells were assessed semi-quantitatively using the following scoring method: 0 negative, 1, 1–10% positive cells, 2, 11–50% positive cells, and 3, more than 50% positive cells.

Statistical analysis

Statistical analysis was carried out using GraphPad Prism 6 software (San Diego, CA, USA) and statistical tools by Microsoft Excel 97-2003 (Boston, MA, USA). Mann-Whitney, one-way ANOVA, T-student and Fisher’s exact tests were used. Data are reported as the averaged value, and error bars represent the standard deviation or standard error mean of the average for each group. Results with p values below 0.05 were considered statistically significant.

Acknowledgments and Funds

Carmelo Nucera (Principal Investigator, Human Thyroid Cancers Preclinical and Translational Research at the Beth Israel Deaconess Medical Center (BIDMC)/Harvard Medical School) was awarded grants by the National Cancer Institute/National Institutes of Health (1R21CA165039-01A1 and 1R01CA181183-01A1), the American Thyroid Association (ATA) and ThyCa:Thyroid Cancer Survivors Association Inc. for Thyroid Cancer Research. Carmelo Nucera was also a recipient of the Guido Berlucchi “Young Investigator” research award 2013 (Brescia, Italy) and BIDMC/CAO Grant (Boston, MA). BIDMC was recipient of a grant by Roche (2011-2012). Roche had no role in the design of this study/manuscript, in the collection, analysis, or interpretation of the data, in the writing and publication of this manuscript.

We thank Drs. Yutaka Kawakami (Keio University, Tokyo, Japan) for HIV-U6-GFP GL3B vector. We thank Drs. Sareh Parangi and Richard Hodin (MGH) for providing some thyroid tumor specimens. We also thank Mr. Bhavik Padmani (Murine Pharmacogenetics Facility), Mr. Neal Smith and Dr. Sushruta Nagarkatti for technical assistance. We thank Drs. Jack Lawler (BIDMC, Harvard Medical School, Boston, USA), John Clohessy (BIDMC, Harvard Medical School), and Florian Karreth (Weill Cornell Medical College, Cancer Center, NYC, USA) for their critical reading of our manuscript.

Conflicts of interests

Authors have not conflicts to disclose for this study/manuscript.

Authors’ contributions

Conception and design: CN

Development of methodology: MD, AH, JNS, ZAA, AF, CS, ECR, GMM, PVH, DDS, CN

Analysis and interpretation of data: MD, PMS, AH, JNS, ZAA, AF, CS, ECR, GMM, PVH, VN, RTB, DDS, PPP, CN

Writing of the manuscript: CN

Review, and/or revision of the manuscript: MD, PMS, AH, JNS, ZAA, AF, CS, ECR, GMM, JK, VN, PVH, RTB, TJG, DDS, PPP, CN

Study supervision: CN

Provision of materials for the study: PMS, AF, AH, JK, VN, MV, TJG, DDS, CN

References

1. Nucera, C., et al., B-Raf(V600E) and thrombospondin-1 promote thyroid cancer progression. Proc Natl Acad Sci U S A, 2010. 107: p. 10649-54.

2. Franco, A.T., et al., Thyrotrophin receptor signaling dependence of Braf-induced thyroid tumor initiation in mice. Proc Natl Acad Sci U S A, 2011. 108: p. 1615-20.

3. Smallridge, R.C., et al., American Thyroid Association guidelines for management of patients with anaplastic thyroid cancer. Thyroid, 2012. 22: p. 1104-39.

4. Ricarte-Filho, J.C., et al., Mutational profile of advanced primary and metastatic radioactive iodine-refractory thyroid cancers reveals distinct pathogenetic roles for BRAF, PIK3CA, and AKT1. Cancer Res, 2009. 69: p. 4885-93.

5. Network, C.G.A.R., Integrated genomic characterization of papillary thyroid carcinoma. Cell, 2014. 159: p. 676-90.

6. Xing, M., BRAF mutation in papillary thyroid cancer: pathogenic role, molecular bases, and clinical implications. Endocr Rev, 2007. 28: p. 742-62.

7. Chakravarty, D., et al., Small-molecule MAPK inhibitors restore radioiodine incorporation in mouse thyroid cancers with conditional BRAF activation. J Clin Invest, 2011. 121: p. 4700-11.

8. Xing, M., et al., Association between BRAF V600E mutation and mortality in patients with papillary thyroid cancer. JAMA, 2013. 309: p. 1493-501.

9. Elisei, R., et al., BRAF(V600E) mutation and outcome of patients with papillary thyroid carcinoma: a 15-year median follow-up study. J Clin Endocrinol Metab, 2008. 93: p. 3943-9.

10. Xing, M., et al., BRAF mutation testing of thyroid fine-needle aspiration biopsy specimens for preoperative risk stratification in papillary thyroid cancer. J Clin Oncol, 2009. 27: p. 2977-82.

11. Alzahrani, A.S. and M. Xing, Impact of lymph node metastases identified on central neck dissection (CND) on the recurrence of papillary thyroid cancer: potential role of BRAFV600E mutation in defining CND. Endocr Relat Cancer, 2012.

12. Xing, M., et al., BRAF V600E and TERT Promoter Mutations Cooperatively Identify the Most Aggressive Papillary Thyroid Cancer With Highest Recurrence. J Clin Oncol, 2014. 32: p. 2718-26.

13. Haugen, B.R. and M.A. Kane, Approach to the thyroid cancer patient with extracervical metastases. J Clin Endocrinol Metab, 2010. 95: p. 987-93.

14. Xing, M., Molecular pathogenesis and mechanisms of thyroid cancer. Nat Rev Cancer, 2013. 13: p. 184-99.

15. Bollag, G., et al., Clinical efficacy of a RAF inhibitor needs broad target blockade in BRAF-mutant melanoma. Nature, 2010. 467: p. 596-9.

16. Bollag, G., et al., Vemurafenib: the first drug approved for BRAF-mutant cancer. Nat Rev Drug Discov, 2012. 11: p. 873-86.

17. Kim, K.B., et al., Clinical Responses to Vemurafenib in Patients with Metastatic Papillary Thyroid Cancer Harboring BRAF(V600E) Mutation. Thyroid, 2013. 23: p. 1277-83.

18. Montero-Conde, C., et al., Relief of feedback inhibition of HER3 transcription by RAF and MEK inhibitors attenuates their antitumor effects in BRAF-mutant thyroid carcinomas. Cancer Discov, 2013. 3: p. 520-33.

19. Sullivan, R.J. and K.T. Flaherty, Resistance to BRAF-targeted therapy in melanoma. Eur J Cancer, 2013. 49: p. 1297-304.

20. Sondergaard, J.N., et al., Differential sensitivity of melanoma cell lines with BRAFV600E mutation to the specific Raf inhibitor PLX4032. J Transl Med, 2010. 8: p. 39.

21. Beroukhim, R., et al., The landscape of somatic copy-number alteration across human cancers. Nature, 2010. 463: p. 899-905.

22. Nehs, M.A., et al., Thyroidectomy with neoadjuvant PLX4720 extends survival and decreases tumor burden in an orthotopic mouse model of anaplastic thyroid cancer. Surgery, 2010. 148: p. 1154-62; discussion 1162.

23. Nucera, C., et al., Targeting BRAFV600E with PLX4720 displays potent antimigratory and anti-invasive activity in preclinical models of human thyroid cancer. Oncologist, 2011. 16: p. 296-309.

24. Nehs, M.A., et al., Late intervention with anti-BRAF(V600E) therapy induces tumor regression in an orthotopic mouse model of human anaplastic thyroid cancer. Endocrinology, 2012. 153: p. 985-94.

25. Zhang, Z., et al., Histone deacetylation of NIS promoter underlies BRAF V600E-promoted NIS silencing in thyroid cancer. Endocr Relat Cancer, 2014. 21: p. 161-73.

26. Chan, C.M., et al., Targeted inhibition of Src kinase with dasatinib blocks thyroid cancer growth and metastasis. Clin Cancer Res, 2012. 18: p. 3580-91.

27. Antonello, Z.A. and C. Nucera, Orthotopic mouse models for the preclinical and translational study of targeted therapies against metastatic human thyroid carcinoma with BRAF or wild-type BRAF. Oncogene, 2013.

28. Haq, R., et al., BCL2A1 is a lineage-specific antiapoptotic melanoma oncogene that confers resistance to BRAF inhibition. Proc Natl Acad Sci U S A, 2013. 110: p. 4321-6.

29. Champa, D., et al., Obatoclax overcomes resistance to cell death in aggressive thyroid carcinomas by countering Bcl2a1 and Mcl1 overexpression. Endocr Relat Cancer, 2014. 21: p. 755-67.

30. Faber, A.C., et al., Differential induction of apoptosis in HER2 and EGFR addicted cancers following PI3K inhibition. Proc Natl Acad Sci U S A, 2009. 106: p. 19503-8.

31. Faber, A.C., et al., mTOR inhibition specifically sensitizes colorectal cancers with KRAS or BRAF mutations to BCL-2/BCL-XL inhibition by suppressing MCL-1. Cancer Discov, 2014. 4: p. 42-52.

32. Nucera, C., et al., Human galectin-3 immunoexpression in thyroid follicular adenomas with cell atypia. J Endocrinol Invest, 2005. 28: p. 106-12.

33. Nucera, C., et al., Maternal thyroid hormones are transcriptionally active during embryo-foetal development: results from a novel transgenic mouse model. J Cell Mol Med, 2010. 14: p. 2417-35.

34. Matulonis, U.A., et al., High throughput interrogation of somatic mutations in high grade serous cancer of the ovary. PLoS One, 2011. 6: p. e24433.

35. Namba, H., et al., Clinical implication of hot spot BRAF mutation, V599E, in papillary thyroid cancers. J Clin Endocrinol Metab, 2003. 88: p. 4393-7.

36. Frasca, F., et al., BRAF(V600E) mutation and the biology of papillary thyroid cancer. Endocr Relat Cancer, 2008. 15: p. 191-205.

37. Dias-Santagata, D., et al., Rapid targeted mutational analysis of human tumours: a clinical platform to guide personalized cancer medicine. EMBO Mol Med, 2010. 2: p. 146-58.

38. Hettmer, S., et al., Mutations in Hedgehog pathway genes in fetal rhabdomyomas. J Pathol, 2013. 231: p. 44-52.

39. Li, J., et al., Antiprimer quenching-based real-time PCR and its application to the analysis of clinical cancer samples. Clin Chem, 2006. 52: p. 624-33.

40. Takei, H., et al., Low-level copy gain versus amplification of myc oncogenes in medulloblastoma: utility in predicting prognosis and survival. Laboratory investigation. J Neurosurg Pediatr, 2009. 3: p. 61-5.

41. Zack, T.I., et al., Pan-cancer patterns of somatic copy number alteration. Nat Genet, 2013. 45: p. 1134-1140.

42. Fischer, A.H., P. Taysavang, and S.M. Jhiang, Nuclear envelope irregularity is induced by RET/PTC during interphase. Am J Pathol, 2003. 163: p. 1091-100.

43. Bond, J.A., et al., In vitro reconstruction of tumour initiation in a human epithelium. Oncogene, 1994. 9: p. 281-90.

44. Fabien, N., et al., Description of a human papillary thyroid carcinoma cell line. Morphologic study and expression of tumoral markers. Cancer, 1994. 73: p. 2206-12.

45. Kurebayashi, J., et al., All-trans-retinoic acid modulates expression levels of thyroglobulin and cytokines in a new human poorly differentiated papillary thyroid carcinoma cell line, KTC-1. J Clin Endocrinol Metab, 2000. 85: p. 2889-96.

46. Shaik, S., et al., SCF(beta-TRCP) suppresses angiogenesis and thyroid cancer cell migration by promoting ubiquitination and destruction of VEGF receptor 2. J Exp Med, 2012. 209: p. 1289-307.

47. Debnath, J., S.K. Muthuswamy, and J.S. Brugge, Morphogenesis and oncogenesis of MCF-10A mammary epithelial acini grown in three-dimensional basement membrane cultures. Methods, 2003. 30: p. 256-68.

48. Ponce, M.L., In vitro matrigel angiogenesis assays. Methods Mol Med, 2001. 46: p. 205-9.

49. Shih, S.C. and L.E. Smith, Quantitative multi-gene transcriptional profiling using real-time PCR with a master template. Exp Mol Pathol, 2005. 79: p. 14-22.

50. Nucera, C., et al., A novel orthotopic mouse model of human anaplastic thyroid carcinoma. Thyroid, 2009. 19: p. 1077-84.

51. Bellelli, R., et al., FOXM1 is a molecular determinant of the mitogenic and invasive phenotype of anaplastic thyroid carcinoma. Endocr Relat Cancer, 2012. 19: p. 695-710.