Identifying Candidate Druggable Targets in Canine Cancer Cell Lines Using Whole-Exome Sequencing

This article is featured in Highlights of This Issue, p. 1335

Over decades, cancer cell lines have been used as a preclinical model system to delineate histopathology, molecular mechanisms, and therapeutic treatments (1, 2). An extensive genetic catalog of somatic variants across a panel of canine cancer cell lines can provide information regarding oncogene addictions and allow for preclinical assessment of targeted agents in appropriate cell lines. Canine cancer cell lines have been procured for use in understanding veterinary cancers paralleling the use of human cancer cell lines and are an invaluable resource for the field of comparative oncology. Understanding the similarities and differences in putative driving mutations between human and canine cancers allows for a more targeted application of canine cancer clinical trials to human therapeutic development and for more effective utilization of this valuable resource (3).

Comparative oncology has long held an important role in our understanding of the mechanisms of disease pathogenesis with a recent increase in focus on the use of companion animals in translational studies to assess the efficacy of emerging therapies. Two primary areas are responsible for this increased interest in comparative oncology. First, the high failure rate of expensive phase III clinical trials has highlighted the limits of preclinical studies in genetically engineered or immune-compromised animal models to recapitulate the complex biological interactions that drive cancer development and response to therapy. Second, the development of immune-based therapies further illustrates the limitations of preclinical animal models that lack a functional immune system. Consequently, to address this problem, the use of dogs as genetically diverse, immune-competent surrogates for therapeutic development has been explored including DNA sequencing to define the genomic landscape of canine cancers. Recent studies describe the genetic variants present in common canine cancers illustrating the homologies that exist with human cancer (4–7). This study in canine cell lines seeks to bridge the gap between variant identification and functional utility by exploring the impact of putative driver genes on cellular signaling and response to therapy. We have used whole-exome sequencing data from 33 cell lines (10 different cancer types) to gather information on the variants that have been implicated as cancer drivers. The goal of this study was to generate a database of somatic variants and profile the mutational landscape in commonly used canine cell lines. A total of 232 known cancer genes that are implicated in human cancers carried at least one somatic variant among these cell lines. The top mutated cancer gene was TP53 with variations in 36% of cell lines, which is comparable with 42% of human cancers bearing TP53 mutations (8). Cross-comparison to human oncogenic driver/mutational hotspot databases identified 61 somatic variants that are known or putative drivers in human cancer. Functional annotation of all protein-coding somatic variants revealed an enrichment of genes in the MAPK and PI3K–AKT pathways. Functional validation revealed that driving mutations in receptor tyrosine kinases, and members of the MAPK and PI3K–AKT pathways were consistent with activated ERK1/2 and AKT phosphorylation. Sensitivity to the MEK1/2 inhibitor, trametinib, was observed in cell lines with driving mutations that activate the MAPK pathway. Thus, identification of oncogenic variants in canine cancers can be utilized to inform the development of therapeutic strategies.

The 33 canine cancer cell lines utilized in this study were either obtained from researchers at multiple institutions or purchased from the ATCC or established from tumor samples in-house (for details, see Supplementary Table S1A). These cell lines were cultured in DMEM supplemented with 10% FBS, 100 U/mL penicillin, 100 μg/mL streptomycin, 1 mmol/L sodium pyruvate (Thermo Fisher Scientific), 1× MEM nonessential amino acid (NEAA), and 2× MEM vitamin solution (Corning) and were used within 10 passages after thawing. All cell lines were periodically tested for Mycoplasma using a PCR-based method (9). Validation of cell lines was done by short tandem repeat analyses using canine StockMarks Genotyping Kit (Thermo Fisher Scientific) as described previously (refs. 10, 11; Supplementary Materials and Methods). The validation results are cataloged in Supplementary Table S2.

Genomic DNA was isolated from 13 normal canine blood samples and 33 cell lines using Qiagen DNeasy Blood and Tissue Kits. Normal blood samples were taken from dogs representing some of the breeds that are frequently diagnosed with the types of cancers represented by these 33 cell lines. WES of 32 canine cancer cell lines and 13 normal stroma and blood samples was conducted by capturing exonic regions from fragmented genomic DNA using the 43.45 Mb custom Agilent SureSelect XT All Exon V2 (part number: 931198) capture kit. For the Bliley sample, the 53.59 Mb custom Agilent SureSelect XT All Exon V1 (part number: 5190-5452) capture kit was used and sequencing was carried out on HiSeq 2500 (100 bp paired-end reads). The remaining 32 samples were sequenced on an Illumina HiSeq 4000 sequencer to generate 151-bp paired-end reads at the Genomics and Microarray Shared Resource at University of Colorado Anschutz Medical Campus (Aurora, CO).

Read quality was assessed by FastQC and Trimmomatic was used to select high-quality reads (phred score >20) and eliminate adapter sequences. The high-quality reads were mapped against the canine genome (v3.1) using BWA tool. Preprocessing of alignment files for variant calling was done in accordance with GATK best practices (12). Freebayes was used to call variants with a min-alternate-count of 2 and min-alternate-fraction of 0.05. Identified variants were filtered for depth of >10 and QUAL >20 and annotated with SnpEff. For additional methodology, resource version, and references, see Supplementary Materials and Methods and Supplementary Table S3.

Because of a lack of extensive canine population genetic variant datasets, we used the available dbSNP database (build 146), with three external studies (491 dogs) as a resource to remove germline variants from the cell line WES dataset (7, 13, 14). The dbSNP 146 build contains lymphoma somatic variants identified in Elvers and colleagues (2015), which were retained in this study. In addition, the variants obtained from normal blood samples of 13 dogs (collected locally) were also used to eliminate germline variants (Supplementary Table S1B). The normal blood sample variants were called at a minimum alternate count of two to allow for a more stringent selection of somatic variants. Following this germline screen, the resulting list of all somatic SNPs and INDELs from 33 cell lines was designated as level 1 variants (Supplementary Fig. S1). Level 1 includes variants within intergenic, intron, and intragenic regions along with cDNA variants. Level 2 used only protein-coding variants from level 1. These level 2 genes were functionally annotated to identify signaling pathways and gene ontology categories that might play a role in cancer biology (Supplementary Materials and Methods; Fig. 4). In addition, amino acid alterations affecting protein function were analyzed using Sorting Intolerant from Tolerant (SIFT; ref. 15).

Level 3 selects level 2 protein-coding variants for genes listed in the COSMIC database [version 83 (699 genes); Fig. 2], and has been expanded to include 5′ UTR and splice site variants of these genes due to their potential to impact protein expression and function (16). To identify putative driver gene mutations, level 3 cancer gene variants were manually cross-checked against known human oncogenic coding variants using cBioPortal http://www.cbioportal.org/ (Supplementary Materials and Methods).

Validation of selected variants was conducted using Sanger dideoxy sequencing (Genewiz) of amplified products from cell line genomic DNA, microarray, and Western blot analyses. Following amplification of genomic regions, the amplicons were evaluated by gel electrophoresis, isolated, and sequenced using either the forward or reverse amplification primers (Supplementary Table S4). Transcript expression was evaluated using data previously obtained from Affymetrix Canine 1.0 ST microarrays (10). Western blot analysis was used to detect protein levels of NF1, p-ERK1/2, total ERK1/2, p-AKT (S473), and total AKT (Supplementary Materials and Methods).

Canine cell lines were assessed for sensitivity to trametinib, an ATP-noncompetitive inhibitor of MEK1/2. Trametinib is FDA approved as single-agent drug or in combination with dabrafenib for the treatment of BRAF-mutant anaplastic thyroid cancer, melanoma, and non–small cell lung cancer. Following plating in 96-well plates at 1,000 to 6,000 cells/well in complete media, 3–5 wells were treated with serial dilutions of trametinib (Selleck Chemicals, 10 mmol/L stock solutions in DMSO) or DMSO vehicle for controls. Cell proliferation was monitored using percent confluence on the IncuCyte ZOOM Live-Cell Analysis System (Essen BioScience) or a resazurin-based fluorescence assay over a 72-hour period. Vogel drug sensitivity assays were extended to 96 hours because their doubling time is longer at 35.9 hours. Readings at the end of treatment were normalized to time zero and expressed as a percent of the DMSO control. IC₅₀ values were calculated in GraphPad Prism (v7) using a nonlinear regression of trametinib log concentration versus percent of control.

The RMACC Summit supercomputer was used to process canine WES data (17). The raw fastq sequences for canine cancer cell line samples are available from NCBI SRA database under Bioproject: PRJNA503864. The sequences corresponding to normal blood/stroma samples can be obtained from BioProject: PRJNA503860.

The short-read sequences obtained from WES of 33 canine cancer cell lines were analyzed to generate somatic mutation (SNPs and INDELs) profiles. The number of quality reads (phred score >20) ranged from 44 million (CLL1390) to 113 million (Vogel). The mean depth of reads mapping to the CanFam3.1 genome was 181X, with a range from 96X (Bliley) to 268X (Gracie). A total of 66,344 somatic variants were identified across the 33 cell lines. Frequency of mutations varied from 15.79 to129.37 per Mb and on average 62 (±7.5)% of all somatic variants were SNPs (Fig. 1A and B). Across all cell lines, the most frequent single-nucleotide change was C/G > T/A (Fig. 1C), corresponding to an average of 51 (±8.5)% and a range of 32% to 77% of SNPs within each cell line. On average 6 (±1.7)% of all level 1 somatic variants were annotated as synonymous mutations, which were eliminated from downstream analysis (Supplementary Fig. S2A). A total of 8,795 (13.2% of level 1) variants were identified in protein-coding regions of 5,085 genes (Supplementary Table S5). In addition, 1,033 (1.6%) and 51 (0.08%) variants were located in the 5′ UTR and splice regions of genes, respectively. The number of genes with somatic variants in protein-coding regions (level 2) ranged from 80 in STSA-1 to 747 in OSW (Supplementary Fig. S2B). Missense variants were the most prevalent type of mutation (78.2% of total variants), followed by frameshifts (12.7 %), and in-frame insertions and deletions (indel; 6.8%; Fig. 1D). Using the SIFT scoring method, 34.8% of all missense mutations were categorized as deleterious (SIFT < 0.05).

The curated cancer gene census from COSMIC (v83, 699 genes) was used to extract the canine somatic variants in genes that are implicated in cancer, resulting in the identification of 232 genes with 450 variants across 33 cell lines. In addition to protein-coding variants, we incorporated 5′ UTR and splice-region variants. The number of cancer gene variants ranged from 3 in CML-10C2 and STSA-1 to 44 in OSW (Supplementary Fig. S2C). Of these genes, 72 are designated as oncogenes, 80 as tumor suppressors, and 31 function as both. Forty-nine cancer genes were mutated in at least three cell lines (Fig. 2). Thirty-two percent of all mutations in cancer genes were deleterious missense mutations (SIFT < 0.05) and 11.5% were frameshifts. Consistent with human cancers, the most frequently mutated gene was TP53 (12 cell lines) with thirteen missense and two frameshifts variants. Mutations in TP53 were primarily localized in the DNA-binding domain (Fig. 3A). The TP53 mutations were found in 4 of the 9 osteosarcomas, 3 of 4 leukemia/lymphoma lines, the mammary and thyroid carcinomas, both of the mast cell tumor lines, and one of the three hemangiosarcoma lines. For MDM2, the cellular antagonist of TP53, transcript levels were elevated in three additional cell lines and TP53 transcript expression was lost in the one of the leukemia/lymphoma lines (1771) and both histiocytic sarcoma lines, thereby identifying alterations in P53 signaling in over half of the canine cancer cell lines (Supplementary Fig. S3; ref. 18). The mean expression levels of MDM2 in cell lines with mutated TP53 was significantly lower than in cell lines with wild-type TP53 (P = 0.0001, F-statistic). This suggests that the auto-regulatory negative feedback loop of TP53 and MDM2 was altered in cell lines with mutated TP53. Compared with the top 20 most frequently mutated genes in human cancer (https://portal.gdc.cancer.gov/), we have identified an overlap of four genes in addition to TP53. In human cancers, these four genes are ranked at 4th (KMT2D), 7th (ARID1A), 11th (BRAF), and 13th (NF1).

The annotation information for all cancer genes listed in cBioPortal was used to identify 61 variants that are considered candidate driver mutations in 29 known cancer-causing genes (Supplementary Table S6). We have identified at least one driver mutation in 30 of 33 cell lines, including 2 drivers in 4 cell lines, 3 in 7 cell lines, 4 in 2 cell lines, and 8 drivers in the OSW cell line. This provides us with putative canine somatic variants that have been identified as human cancer drivers, either computationally or through direct experimental evidence, for further investigation and potential therapeutic targeting.

The 29 genes with 61 driver mutations were grouped into four functional categories (Supplementary Table S6):

• (i) Chromatin organization (GO:0006325) and chromatin binding (GO:0003682) – Thirteen driver genes were categorized as chromatin organizers, a class which encompasses remodeling, silencing, assembly, and maintenance of chromatin as well as binding to chromatin. The majority of these gene products function as tumor suppressors; however, CIC, KDM6A, KMT2D, and TP53 can also function as oncogenes and PPM1D is considered an oncogene. Frameshifts represent 48.5% of the driver variants identified in these 13 genes, which is significantly higher than the rate in the entire dataset (11.3%), but consistent with loss of tumor suppressor function.

• (ii) Regulation of kinase activity (GO: 0043549) – Kinases and phosphatases are cellular signaling proteins that coordinately regulate a variety of cellular processes. Activating mutations in kinases and their regulators are known to dysregulate tissue growth and cellular proliferation (19). There are 13 genes that have been identified in this category and a majority of them encode proteins with kinase activity (BRAF, ERBB2, KIT, PIK3CA), and GTPase activity (KRAS, NRAS, RAC1). Furthermore, genes encoding phosphatases (PTEN and PTPN11), regulators of the cyclin-dependent kinases (CCND3 and CDKN1A), FBXW7, and the tumor suppressor and GTPase activator of RAS, NF1, were also identified as putative drivers of cancer. Overall 23 (14 missense, 8 frameshifts, and 1 in-frame deletion) variants were identified as driver mutations in genes that regulate kinase activity. Half of these frameshifts were identified in NF1 (CML-6M, Nike, OSW, and STSA-1; Fig. 3B).

• (iii) Small GTPase binding (GO:0031267) and GTPase activity (GO:0003924) – Five genes were identified in this category and four of these genes (BRAF, KRAS, NRAS, and RAC1) overlapped with regulation of kinase activity group (Supplementary Table S6). Included in this category, was the V850fs variant in SRGAP3 identified in the Gracie cell line. Although the cBioPortal database did not provide statistically significant evidence within OncoKB and/or cancer hotspot driver mutations, the most frequently mutated region in SRGAP3 in human cancers was also V850fs. In addition, SRGAP3 acts as a tumor suppressor by negatively regulating RAC1, and low expression of this gene leads to anchorage-independent cell proliferation in breast cancer (20).

• (iv) RNA binding (GO:0044822) – This category included two genes: SPEN, transcriptional repressor, with three frameshifts, and EIF1AX, translation initiation factor, with one missense mutation.

Using Sanger sequencing, we have validated the variants for TP53, KRAS, NRAS, KIT, KMT2D, and PIK3CA. We also used Canine 1.0ST microarrays and Western blot analysis to confirm that cell lines bearing NF1 frameshift variants (OSW, STSA-1, CML-6M, Nike) exhibited decreased transcript and protein expression (Fig. 3C and D). Loss of NF1 protein was greatest in Nike and OSW lines with mutant allelic frequency (AF) near one. Although no NF1 variants were identified in the Gracie cell line, NF1 expression at both mRNA and protein levels was reduced, suggesting a similar decrease in negative regulation of RAS activity by NF1.

To identify mutated genes that play a role in biological processes and signaling pathways associated with cancer beyond those included in the COSMIC database, the bioinformatics resource, DAVID, was used to cluster genes by shared function or domain (see Supplementary Materials and Methods). Using three databases (GO, KEGG, and PFAM), 780 terms were selected with a minimum of one cancer gene and the presence of that term in at least six cell lines. There were eight functional terms that were associated with all 33 cell lines including GO:0044822∼poly(A) RNA binding, GO:0004930∼G-protein coupled-receptor, GO:0006555∼regulation of transcription. Pathways with the most mutated genes were PI3K–AKT (cfa:04151) and MAPK signaling (cfa04010) with 91 and 61 mutated genes across 32 and 29 cell lines, respectively.

We selected 51 functional terms that could be manually clustered into 10 categories associated with cancer: DNA repair, Apoptosis, Cell cycle & proliferation, Radiation response, Chromatin binding/remodeling, Angiogenesis, MAPK kinase signaling and regulation, PI3K–AKT signaling and regulation, Phosphatase activity, and Kinase activity (Supplementary Table S7). A total of 595 genes, including 94 cancer genes, were annotated in these 10 categories. We selected 338 genes (57%) with high impact mutations (tolerated missense were excluded) to plot within the 10 selected functional categories to create a dataset of variants that can be explored as putative drivers in canine cancers (Fig. 4).

Approximately, 81% of cell lines had at least one gene mutation in DNA repair with 67 mutations in 46 genes (66% deleterious and 21% frameshift mutations). The OSW (leukemia/lymphoma) cell line had both the highest number of mutated genes and the highest number of mutated DNA repair genes (7), followed by four each in 1771 (leukemia/lymphoma) and BRMCT (mast cell tumor). Some key mutations identified were in mismatch repair genes: MSH3 (K306fs), MSH6 (R709fs); homologous recombination genes: MCM8 (Q245*), RAD21L1 (D95N); nucleotide excision repair: ERCC4 (S695T, R269S); and a DNA cross-link repair gene: DCLRE1A (Q331*).

We have identified 86 variants in 51 genes that regulate the cell cycle and cellular proliferation. Similar to DNA repair, OSW had the most variants (11), followed by C2 (6), and OSA8 (5). In addition to TP53, NOTCH2, BUB1, CCNB3, LRP2, and XRN2 were frequently mutated genes in this category, each with 3–5 variants (Fig. 4).

There were 28 genes with 48 mutations impacting 85% of cell lines in this category. After TP53 and NF1, the genes with the highest number of variants were BIRC6 (3), and GSDMA (2).

The genes in this functional group regulate cellular responses to radiation. There were 20 genes with 40 variants observed in 24 cell lines. TRIMC13 was mutated in osteosarcoma and transitional cell carcinoma cell lines – HMPOS, Tyler1, and Tyler2 (P19fs). Mutation of TRIM genes has been associated with human blood cancers and may act as either a tumor suppressor or an oncogene (21).

The majority of driver mutations identified play a role in chromatin reorganization. In this study, along with the 21 known cancer genes, we have identified 40 additional genes with 90 variants in this category. The cell lines with the highest number of variants were OSW (14), SB (6), Angus (5), and BRMCT (5). Two additional lysine demethylase genes not included in the COSMIC database were mutated in CML-6M, CMT12 (KDM1A, mutation: A9_A10del), and BRMCT (KDM3A, mutation: S805F). Genes in SWI/SNF complex (e.g., ARID, SMARC genes) were mutated in 30% of the cell line panel.

Multiple components of the MAPK pathway have been implicated in tumor development (22). This category combined four functional terms resulting in 57 genes with 90 variants in 94% of the cell lines. Sixty-seven percent of these genes mapped specifically to the KEGG MAPK signaling pathway (cfa04010; Supplementary Fig. S4; Supplementary Table S8). In addition to the known cancer genes (K/NRAS, BRAF) that regulate MAPK pathway, mutations in RET (R815S), MOS (V292A), PTPRJ (N1148fs, D210fs), and PTPN11 (G503V) were identified in this category.

PI3K–AKT signaling is recurrently altered across a multitude of cancers and over 40 anticancer drugs targeting PI3K–AKT pathway are being tested in different phases of clinical trials (23, 24). We have identified 91 genes in the KEGG PI3K–AKT signaling pathway (cfa04151) with at least one gene mutated in each of the 32 cell lines. The most frequently mutated gene (27% of the cell lines) was PTPRQ (Fig. 4). Fifty-eight percent of genes in this group were mapped to KEGG pathway (cfa04151; Supplementary Fig. S5; Supplementary Table S9). Among the selected groups, PI3K–AKT signaling is the most frequently mutated pathway. However, it is important to note that 17 of the mutated genes encode for extracellular matrix proteins and integrins. Variants in proteins that inhibit AKT, PPP2R3A (R981H), PPP2R2B (D215V), and THEM4 (S171L) were identified in Kinsey, Parks, and Nike, respectively.

Aberrant kinases have been targeted to treat a variety of cancers (25). Eighty-five percent of the cell lines have at least one gene in this group, with BRAF, ALPK2, and OBSCN as the most frequently mutated. OSW and D-17 had the most kinase variants with 19 and 8 genes, respectively. In addition to known cancer genes with altered kinase activity, we have identified variants in cyclin-dependent kinases (CDK9, CDK12, CDK13, CDK14), tyrosine kinases (EPHA2, TNK2, NPR2), and serine/threonine kinases (PRKCD, SIK2, SIK3, BMP2K, PASK, WNK4).

Altered expression of phosphatases can result in dysregulation of multiple signaling pathways, including MAPK, PI3K, Platelet-derived growth factor, Integrin, and EphA signaling pathways (26). Using the four GO and two PFAM functional categories, 25 genes from 23 cell lines were clustered in this group. The cell lines OSW (five), Kinsey (four), and Jones (four) had the highest frequency of variant phosphatases. Twenty protein tyrosine phosphatase variants in genes like PTPRQ, PTPRD, PTPRJ, PTPRO, and PTPRF were identified within 14 cell lines. Five genes that function as dual specificity phosphatases (DUSP16, DUSP9, DUSP5, SSH3, PTEN) carry an alternate allele in at least one of the five cell lines.

Tumor survival under hypoxic conditions is poor; hence, formation of a blood vessel network is essential for solid tumor survival and growth (27). For this functional category, 15 genes with 18 variants were identified across 12 cell lines. Examples of genes identified were PECAM1 (N619I), ENPEP (I346F), and PML (R9fs).

Both the putative driver gene analysis (Supplementary Table S6) and the functional annotation analysis of somatic variants (Fig. 4) identified activation of the MAPK and PI3K–AKT pathways as drivers and potential therapeutic targets in this panel of canine cell lines. As detailed in Supplementary Materials and Methods, we have categorized phosphorylation of ERK1/2 and AKT(S473) as constitutive (equal in the presence/absence of serum), activated (greater phosphorylation in the presence of serum, and low (expression/phosphorylation). On the basis of these criteria, 20 of the cell lines were described as having constitutive and 4 had activated ERK1/2 phosphorylation (Table 1; Supplementary Fig. S6). Constitutive activation was frequently observed in conjunction with activating mutations in BRAF, KRAS, NRAS, and PTPN11, or with inactivating mutations in NF1. Assessment of AKT phosphorylation revealed 13 cell lines exhibiting constitutive AKT activity and 12 exhibiting AKT activation.

Because 20 of the cell lines carried a putative driving mutation in the RTK–RAS–RAF–MAPK pathway and 29 cell lines had either constitutive or activated ERK1/2, we tested for sensitivity to the MEK1/2 inhibitor, trametinib (Fig. 5; Table 1). The IC₅₀ values across this cell line panel ranged from 0.4 nmol/L to 6.8 μmol/L with several cell lines failing to achieve an IC₅₀ at 10 μmol/L, the highest concentration tested. We, however, did not identify any MEK1/2 mutations among the 33 cell lines. Twelve cell lines had an IC₅₀ less than 32 nmol/L and were considered sensitive using the average C_max concentration at steady state (32 nmol/L) as a cutoff (28). Putative driving mutations in NF1, BRAF, NRAS/KRAS, and SRGAP3 were identified in these sensitive cell lines. Two sensitive cell lines had no putative drivers (Kinsey and Cindy). The Kinsey cell line was the most sensitive (IC₅₀ 0.4 nmol/L) and exhibited constitutive ERK1/2 phosphorylation and no activation of AKT signaling. For the 21 cell lines with an IC₅₀ greater than 32 nmol/L, 17 had constitutive or activated AKT. AKT activation could be due to identified mutations in PIK3CA (SB) or PTEN (Jones), or other upstream activating mutations (Table 1). Previous studies have revealed copy number loss of PTEN (29) and microarray analysis revealed significant reduction of PTEN gene expression in 6 of the 19 resistant cell lines (Abrams, CLL1390, CMT12, DH82, HMPOS, and OSA8; ref. 10). Trametinib-sensitive cell lines were not more likely to carry mutated MAPK pathway genes than nonsensitive lines (Fisher exact test, P = 0.27). However, nonsensitive cell lines (including cell lines with MAPK gene variants) had a significantly higher number of TP53 pathway alterations (Fisher exact test, P < 0.05) in comparison to sensitive cell lines (Table 1). Thus, TP53 mutations may contribute to resistance in cell lines predicted to be sensitive to trametinib based on driving mutations in the MAPK pathway.

This study is the first report of next-generation sequence analyses across this panel of 33 canine cancer cell lines. These cell lines have been previously used for preclinical screening of therapeutics to assess their potential for use in the treatment of spontaneous canine cancer (10, 30, 31). We found that the mutational burden of these cell lines varied from 15.79 to 129.37 mutations per MB. These values fell within the described range of variants for human cancers (8). Cell lines with high mutational burden, OSW and Vogel, had deleterious mutations in a number of genes associated with DNA damage and repair including: TP53, RAD51AP1, UIMC1, ZFYVE26, PRKDC, MSH3, and TRRAP (Fig. 4; Supplementary Table S7). Similar to many human cancers, the missense mutational spectrum was dominated by C>T:G>A mutations (32).

To identify putative drivers, we screened for oncogenic mutations in Cancer Census Genes from COSMIC database. From this analysis, the most frequently mutated gene identified was TP53, consistent with the predominance of TP53 mutations in human cancers. We have identified alteration of TP53 signaling in 18 cell lines including TP53 variants, increased expression of MDM2, and decreased expression of TP53 [osteosarcoma 56% (5/9), hemangiosarcoma 67% (2/3), lymphoma/leukemia 100% (4/4)]. In comparison, available WES of canine tumors has identified mutations or copy number loss of TP53 in 83% of osteosarcomas (4), 35% of hemangiosarcomas (6), 15.6% of B-cell lymphomas, and 4.9% of T-cell lymphomas (7). Analysis of canine melanomas identified TP53 mutations or MDM2 amplifications in 43% of canine melanomas (33). We did not identify any TP53 mutations in the five canine melanoma cell lines analyzed in this study.

When directly compared with human variants, 61 putative variants were considered oncogenic (Supplementary Table S6). Included in this list are frameshift or stop-gained mutations in tumor suppressor genes predominantly involved in chromatin binding/organization, loss of PTEN phosphatase (OSW, Jones), NF1, negative regulator of the RAS signaling system, (OSW, STSA-1, CML-6M, Nike; Supplementary Table S6, Fig. 3), and Rho GTPase activating protein SRGAP3 (Gracie). Other oncogenic mutations are conserved variants with homology to activating mutations in human genes such as the V588E variant in BRAF, which is homologous to the V600E mutation, KRAS (E63K, G12D/V), NRAS (Q61R/H), ERBB2 (V659E), KIT (V558_E560del, L575P), PIK3CA (E453K), and PTPN11(G503V).

In canine tumors, the BRAF V588E variant has been identified in canine bladder carcinoma (67%), prostatic carcinoma (80%), pulmonary carcinoma (6%), oral squamous cell carcinoma (11%), melanoma (6%) and melanocytoma (17%), glioma (15%), and peripheral nerve sheath tumors (22%; refs. 34–38). Previous studies exploring canine melanomas have identified the NRAS Q61R mutation in the Jones cell line (33, 39) and indicated that canine melanomas have activating NRAS/KRAS/HRAS mutations in 24% of samples (33). Activating KRAS and NRAS mutations have also been identified in 17% of canine pulmonary carcinomas (38), 58% of acute myelogenous leukemia samples, and NRAS mutations were identified in 14% of acute lymphocytic leukemias (40) while KRAS mutations have been identified in 1 of 14 gastric carcinomas (41), and 4 of 5 pancreatic carcinomas (42). We found activating KRAS or NRAS mutations in 2 of the 6 melanoma lines and 5 of the 33 cell lines including two osteosarcoma cell lines and a bladder carcinoma. The majority of oncogenic variants identified in this study activate the RAS/RAF/MAPK signaling pathway; consequently, we assessed activation of the downstream targets ERK1/2 and AKT in correlation with sensitivity to the MEK1/2 inhibitor trametinib as a sample screen to study targeted therapies in canine cancers. We found that cell lines with activating mutations in the MAPK pathway including BRAF, NRAS/KRAS, and loss of NF1 exhibited constitutive activation of ERK1/2 and were often sensitive to trametinib. In contrast, activating mutations in receptor tyrosine kinases including KIT (mast cell line and melanoma) and ERBB2 (thyroid carcinoma line) were associated with moderate sensitivity to trametinib which may be associated with increased activation of the AKT signaling pathway and TP53 mutations.

Similar to our findings in the mast cell tumor cell lines, activating KIT mutations have been commonly identified in canine mast cell tumors including the L575P variant identified in the BRMCT cell line and a tandem duplication in exon 11 of the C2 cell line (43, 44). A recent report on canine mast cell tumors identified KIT mutations in 29% of tumors. Copy number aberrations were also assessed with deletion of CFA 5, including TP53; gain of CFA 31, including RUNX1; and a small number of subchromosomal regions on CFA 36, occurring significantly more frequently in mutant-KIT tumors (45). KIT variants were also identified in three samples from a panel that included 27 canine melanoma tumors and three canine melanoma cell lines (33). In the same study, a putative tumor suppressor, PTPRJ was mutated in 19% of the samples. We also identified these PTPRJ variants in two of the melanoma cell lines.

Other putative drivers that impinge on the MAPK pathway and modulate trametinib sensitivity include: PTPN11 (SHP2), NF1, RAC1, CIC, and potentially SRGAP3. The phosphatase PTPN11 G503V variant was identified as likely oncogenic in cBioPortal acting through activation of RAS signaling and has been associated with Noonan syndrome (46). Furthermore, next-generation sequencing of human histiocytic sarcomas has identified mutations in PTPN11 (47). Specifically, an E76K gain-of-function mutation has been identified in histiocytic sarcomas in 36.6% of Bernese Mountain Dogs and 8.7% of dogs from other breeds (48). The NF1 protein terminates the RAS active state, thus, frameshift mutations causing loss of functional NF1 would result in constitutive RAS activation (49). RAC1 P29S mutations similar to the P29T mutation identified in the Moresco osteosarcoma cell line have been observed in up to 9% of human melanomas and are associated with reduced sensitivity to RAF and MEK inhibitors and elevated ERK1/2 phosphorylation (50). Similarly, loss of CIC, a transcriptional repressor of the MAPK targets ETV1, 4, and 5, results in activation of MAPK targets independent of ERK1/2 phosphorylation as well as resistance to MAPK pathway inhibitors (51, 52). Loss of the SRGAP2/3 tumor suppressors has been observed in human osteosarcomas and breast cancer where they have also been implicated in metastasis suppression (53). Furthermore, these SRGAP family members may also bind and repress MEKK1 (54) such that their loss may increase ERK1/2 activation. Overall, our data support similar interactions of these mutant proteins with the MAPK pathway that was affirmed by examining the sensitivity of cell lines with mutations in these genes to the MEK1/2 inhibitor, trametinib.

When mutations implicated in the activation of the PI3K–AKT signaling pathway were examined, we identified variants in PTEN (OSW lymphoma and Jones melanoma lines) and PIK3CA (SB hemangiosarcoma). Loss of PTEN expression has previously been described in the HMPOS, Abrams, OSA8, CLL1390, and DH82 cell lines (10, 29). Recent WES of 20 canine hemangiosarcomas identified PIK3CA mutations, (9 cases), TP53 (7 cases), PTEN (2 cases), and a PLCG1 mutation in one case. PIK3CA, PTEN, and PLCG1 were mutually exclusive mutations, with six of the cases bearing the known human PIK3CA H1047R–activating mutation (6). The canine hemangiosarcoma cell lines in this study have shown increased sensitivity to a combined PI3K/mTOR inhibitor. Mutations in PTEN were identified in 9.8% of canine T-cell lymphomas (7). Recurrent genomic losses in PTEN have also been described in canine histiocytic sarcomas and osteosarcoma (55, 56).

In addition to the known cancer genes identified in these cell lines from ten different cancer types, we have cataloged variants from an additional 501 functionally clustered genes. This dataset of variants (Supplementary Tables S5 and S7) can be utilized to screen for novel driver mutations and therapeutic targets. For example, the prevalence of variants in chromatin remodeling genes, members of the SWI/SNF complex, phosphatases, and cell-cycle kinases has been suggested to indicate sensitivity to therapeutic agents that target these regulatory processes (26, 57–59).

In summary, variants identified through WES of this panel of 33 canine cell lines provide a database of somatic mutations that can be explored for their role in the development and progression of canine cancers. While the use of WES largely limits this study to protein-coding sequences, it allows for greater depth of coverage to identify somatic variants and simplifies bioinformatic analysis compared with whole genome sequencing. Most current therapeutic choices in “basket” or “umbrella” trials are based on mutational analysis of protein-coding variants. Thus, in vitro analysis of drug sensitivity in a background of known protein-coding somatic mutations can be used to correlate drug sensitivity to the observed genomic profile. The significance of this study lies not only on generating a searchable dataset of variants (Supplementary Table S5) but also in evaluating drug sensitivity based on driver gene identification (Table 1). Future studies involving correlations of the variants and gene expression profiles in canine cancer cell lines may contribute to our understanding of the mechanisms that contribute to therapeutic sensitivity. Such comparisons between canine and human cancers are only possible when analogous data are available in both model systems. Hence, the significance of this dataset lies in studying cancer by cross-species comparison. These cell lines provide opportunities for preclinical testing and optimization of novel therapeutics prior to conducting clinical trials in spontaneously occurring canine cancers.

No potential conflicts of interest were disclosed.

Conception and design: S. Das, D.L. Gustafson, D.L. Duval

Development of methodology: S. Das, R. Idate, K.E. Cronise, D.L. Duval

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): R. Idate, K.E. Cronise, D.L. Duval

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): S. Das, R. Idate, K.E. Cronise, D.L. Duval

Writing, review, and/or revision of the manuscript: S. Das, K.E. Cronise, D.L. Gustafson, D.L. Duval

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): S. Das, R. Idate, K.E. Cronise, D.L. Gustafson, D.L. Duval

Study supervision: D.L. Duval

We thank Tyler Eike and Kenneth L. Jones for IT support and sequence analysis pipeline development, respectively. We acknowledge the help of Nathaniel Gustafson and Rujuta Idate in screening for driver mutations and pathway analyses. This work was funded in part by The Anschutz Foundation (awarded to D.L. Duval), the Colorado State University Cancer Supercluster (awarded to D.L. Duval), and NIH/NCI (grant no. P30 CA046934, awarded to D.L. Gustafson).

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.