Kinome Array Profiling of Patient-Derived Pancreatic Ductal Adenocarcinoma Identifies Differentially Active Protein Tyrosine Kinases

Pancreatic cancer remains one of the most difficult malignancies to treat. Minimal improvements in patient outcomes and persistently abysmal patient survival rates underscore the great need for new treatment strategies. Currently, there is intense interest in therapeutic strategies that target tyrosine protein kinases. Here, we employed kinome arrays and bioinformatic pipelines capable of identifying differentially active protein tyrosine kinases in different patient-derived pancreatic ductal adenocarcinoma (PDAC) cell lines and wild-type pancreatic tissue to investigate the unique kinomic networks of PDAC samples and posit novel target kinases for pancreatic cancer therapy. Consistent with previously described reports, the resultant peptide-based kinome array profiles identified increased protein tyrosine kinase activity in pancreatic cancer for the following kinases: epidermal growth factor receptor (EGFR), fms related receptor tyrosine kinase 4/vascular endothelial growth factor receptor 3 (FLT4/VEGFR-3), insulin receptor (INSR), ephrin receptor A2 (EPHA2), platelet derived growth factor receptor alpha (PDGFRA), SRC proto-oncogene kinase (SRC), and tyrosine kinase non receptor 2 (TNK2). Furthermore, this study identified increased activity for protein tyrosine kinases with limited prior evidence of differential activity in pancreatic cancer. These protein tyrosine kinases include B lymphoid kinase (BLK), Fyn-related kinase (FRK), Lck/Yes-related novel kinase (LYN), FYN proto-oncogene kinase (FYN), lymphocyte cell-specific kinase (LCK), tec protein kinase (TEC), hemopoietic cell kinase (HCK), ABL proto-oncogene 2 kinase (ABL2), discoidin domain receptor 1 kinase (DDR1), and ephrin receptor A8 kinase (EPHA8). Together, these results support the utility of peptide array kinomic analyses in the generation of potential candidate kinases for future pancreatic cancer therapeutic development.


Introduction
In vivo, kinases are heavily trafficked and compartmentalized to microdomains and cellular substructures where they act in concert with other kinases, receptors, and effector proteins to realize intricate signaling mechanisms that give rise to complex cellular behaviors. Despite this complexity, kinases are traditionally studied in isolation: A kinase of interest is typically purified from a sample and the enzymatic activity of the isolated kinase is subsequently interrogated. Kinase isolation is often dependent upon antibody binding and subject to the challenges and expenses of antibody-based purification techniques. Traditional assays remove each kinase of interest from its interacting partners and the physiologic conditions in which it normally operates. By breaking these networks into discrete components and examining individual kinases in isolation, researchers discover initiators of kinase phosphorylation activity and catalogue the peptide targets that are activated or deactivated by direct phosphorylation by a specific kinase. Unfortunately, these purification techniques incentivize researchers to minimize the integrity of a kinase's normal physiologic environment in favor of a clean, easily assayable sample. However, the cleaner a purified kinase sample becomes, the less it can represent that kinase's activity under physiological endogenous conditions. While this reductionism continues to allow researchers to meaningfully investigate discrete components of kinomic signaling networks, the knowledge gained remains incomplete. In this study, we combined PamGene multiplexed kinome activity array data with four bioinformatic pipelines to identify protein kinases responsible for the differential phosphorylation activity observed in patient-derived and commercial pancreatic ductal adenocarcinoma (PDAC) cell lines compared to patient-derived wild-type pancreatic tissue specimens ( Figure 1). Protein kinases are enzymes capable of phosphorylating other proteins to regulate biochemical signaling pathways and modulate cellular behavior [1]. In cancer, many protein kinases are associated with cancer cell initiation, progression, and metastasis, as well as relapse and survival. Small-molecule protein kinase inhibitors represent an increasingly successful therapeutic strategy for a range of cancer types with over 50 FDA-approved protein kinase inhibitors in clinical use [2][3][4] and hundreds more undergoing clinical or preclinical trials. Protein kinase inhibitors improve patient outcomes significantly and protein kinases are, therefore, a popular target for difficult-to-treat malignancies such as pancreatic cancer. Pancreatic cancer patients demonstrate the lowest five-year relative survival rates of any cancer type, with decades of medical research unable to increase these rates beyond 9% [5]. Protein kinase inhibitors erlotinib (OSI-774; Tarceva), everolimus (RAD001; Afinitor), and sumitinib (SU11248; Sutent) are approved for the treatment of pancreatic cancers [4], though the majority of these approvals are for a rare (1% to 2%) subtype of pancreatic cancer known as pancreatic neuroendocrine tumors [6] rather than the far more common PDAC.
Previous, excellent genomic and proteomic studies have characterized PDAC subtypes [7][8][9][10][11]. Different and new PDAC subtypes continue to be defined, with myriad genetic or molecular factors increasing the resolution of the disease spectrum. While the present study could not comprehensively address the complexity and heterogeneity of PDAC tumor subtyping, our results inform growing understanding of the unique protein tyrosine kinase networks that play a role in one or more of the commercial and patient-derived cell lines investigated (Figures 2-4).
To define kinases that are differentially active in one or more PDAC cell lines and identify potentially actionable drug targets, we sought to break away from linear, one-by-one investigations of individual PDAC kinases and consider instead the totality of PDAC kinase networks within each group. To accomplish this, we combined emerging laboratory technologies and contemporary bioinformatic pipelines to identify lead candidate kinases based on peptide phosphorylation signatures. These technologies include the PamStation and the Protein Tyrosine Kinase PamChip. The bioinformatic pipelines include the Kinome Random Sampling Analyzer (KRSA) pipeline, developed by our own laboratory [12][13][14][15][16][17][18][19], and the Upstream Kinase Analysis (UKA) pipeline, which is part of the BioNavigator software tool developed by collaborators at PamGene [20,21]. To passively validate and further contextualize the results of this approach, two additional bioinformatic pipelines were utilized. These additional bioinformatic pipelines include the Post-Translational Modification Signature Enrichment Analysis (PTM-SEA) pipeline developed at the Broad Institute of MIT and Harvard [22], as well as the Kinase Enrichment Analysis Version 3 (KEA3) developed by the Ma'ayan laboratory [23]. Using laboratory equipment designed by PamGene, we measured the relative phosphorylation levels of 198 representative peptide substrates for PDAC cells and wild-type pancreatic tissue before running the resultant data through each of these four self-contained analytical platforms. The analytical results from each pipeline were interpreted individually as well as in combination to maximize the strengths of each pipeline's respective algorithm and reference database set. The results were synthesized to identify kinases likely responsible for the unique phosphorylation signatures of PDAC (Figures 2-4).

PamGene Kinome Activity Profiling Using Protein Tyrosine Kinase PamChip ®
The PamChip is a kind of array. Groups of identical peptide fragments populate a single spot on the chip (Figure 1). The peptide spots contain amino acid sequences of sufficient length and complexity to represent the phosphorylation sites of specific proteins. When these peptide spots are exposed to experimental samples, they become phosphorylated by the activated kinases within that sample. Phosphorylation events correspond precisely with the enzymatic activity of the sample's kinase networks. Fluorescently tagged phosphoantibodies produce a signal proportional to the phosphorylation activity of the kinases contained within the sample. Using specific collections of peptide fragments, the intensities reported by a PamChip provide a signature that can be used to identify the kinases responsible for the observed phosphorylation activity. UKA or KRSA pipelines process differentially phosphorylated peptide fragment identities and signal intensities. For PTM-SEA and KEA3 pipelines to evaluate PamChip signatures, the peptide fragments must first be converted to protein identities. Because some peptide spots on the PamChip represent multiple proteins and because multiple peptide spots may represent different sites on a single protein, this conversion results in a loss of information. The signal intensities on the PamChip-which quantify the degree to which a peptide fragment has been phosphorylated-are also lost when PamChip signatures are converted to the binary "differentially expressed" input categories required for PTM-SEA or KEA3 pipelines. In other words, instead of interpreting "peptide-fragment-1 with signal intensity 10.3," PTM-SEA or KEA3 only process "protein-a." For this reason, we preferentially used the UKA and KRSA pipelines for candidate kinase identification, using the less well-suited (but excellent in their own right) PTM-SEA and KEA3 pipelines for passive validation.
Unabridged comparison of the upstream kinases that each pipeline identified as being responsible for the observed phosphorylation patterns of each pancreatic cancer cell line compared to control wild-type patient-derived pancreatic cells are presented in Figures 2-4 as well as  Supplementary Table S1. The databases associated with each pipeline offer different levels of kinase coverage. This is advantageous for identification of true positives, although the absence of a kinase from a pipeline's identification output cannot be considered a negative indication of that kinase serving as a causative factor in the phosphorylation patterns observed. This approach, therefore, minimizes type I error (i.e., erroneous rejection of a true null hypothesis) and accommodates type II error (i.e., erroneous acceptance of a false null hypothesis). UKA and KRSA bioinformatic pipelines are specifically designed to analyze PamChip kinome activity data. As such, each respective algorithm can integrate multiple PamChip data metrics into its final analytical output. The PTM-SEA and KEA3 pipelines, although capable of excellent analytic activity, are designed for slightly different applications. In their native state, significant programmatic alterations are required to generate results that meaningfully interpret PamChip experimental data. The most significant differences between UKA or KRSA bioinformatic pipelines and PTM-SEA or KEA3 bioinformatic pipelines relate to how these two groups handle multiple dimensions of data and how these two groups natively interpret phosphorylated peptides.

UKA and KRSA Combinatory Analysis
In Table 1, we average the percentile rankings without any concern for false negatives. The goal was to increase positive identification of kinases whose activity is truly a causative factor in the differential phosphorylation patterns observed between tumor and wild-type cells. In Table 2, we attempt to decrease the presentation of false positives by dividing these averages by the number of pipelines employed-in this case, two. In Tables 1 and 2, we also combine the results from the two patient-derived pancreatic cancer cell lines (PDCL15, PDCL5; Patient-Derived) and we combine the results of all pancreatic cancer cell lines (PDCL15, PDCL5, PANC1; All), in order to generate average percentile rankings that inform the general kinase activity of pancreatic cancer. These tables contain only the highest scoring kinases according to UKA and KRSA bioinformatic pipelines. If two or more kinases received identical average percentile ranks or identical weighted average percentile ranks, then rank order was determined arbitrarily. While our complete output data may be found in Table S1, the present tables include only the top 10 highest ranked kinases. Kinases with identical scores were arbitrarily assigned sequential ranks.  PANC1  #1  LCK  SRC  96%  91%  100%  PANC1  #2  DDR2  DDR  96%  96%  PANC1  #3  LYN  SRC  95%  91%  99%  PANC1  #4  SRC  SRC  92%  91%  92%  PANC1  #5  ABL1  ABL  91%  87%  95%  PANC1  #6  TEC  TEC  90%  100%  80%  PANC1  #7  FYN  SRC  90%  91%  88%  PANC1 #8

Expanded PTM-SEA and KEA3 Combinatory Analysis
Tables 3 and 4 expand our analyses to include results obtained through PTM-SEA and KEA3 bioinformatic pipelines. Table 3 provides unweighted average percentile rankings. Table 4 attempts to decrease false positives by weighting these averages according to the number of pipelines which identify a given kinase as being responsible for the observed phosphorylation differences between a pancreatic cancer cell line and patient-derived wild-type pancreas. As above, Tables 3 and 4 also combine results from our two patient-derived pancreatic cancer cell lines (PDCL15, PDCL5; Patient-Derived) as well as from all pancreatic cancer cell lines (PDCL15, PDCL5, PANC1; All).

Identification of Lead Candidate Kinases
Our results confirm the activity of known protein tyrosine kinase-related pathways previously reported as perturbed in pancreatic cancer [24][25][26][27][28]. These results also identify protein tyrosine kinases as yet understudied or unreported in pancreatic cancer. Because our experimental model allows us to maintain the integrity of kinase networks, our data suggest involvement of signaling pathways and regulatory cascades in pancreatic cancer pathophysiology. This study presents evidence in support of the continued development of previously established inhibitory therapeutics that target select protein kinases, such as epidermal growth factor receptor (EGFR) [25], ephrin receptor A2 (EPHA2) [29], and SRC proto-oncogene kinase (SRC) [27], and our experimental results identify new PDAC targets, such as B lymphoid kinase (BLK), lymphocyte cell-specific kinase (LCK), and ABL proto-oncogene 2 kinase (ABL2), which may play a critical role in cancer cell biochemistry or desmoplastic inflammatory cell behavior.
While data gained from reductionist kinase investigations have previously been used to support other complex biochemical studies such as mass spectrometry-based proteomic studies [30], peptide-based kinome array profiling offers unique advantages. Traditionally, peptides identified as kinase targets are probed in vitro and in vivo to examine the behavior of a given kinase. Following that, genetic techniques that knockdown or constitutively express a given kinase further probe the activity of that kinase within normal or experimental biological milieus. But many of these strategies are unable to measure multiple kinases simultaneously while also maintaining kinomic network integrity in complex biological samples. Peptide-based kinome array profiling can accomplish this by overcoming kinase isolation requirements in order to evaluate hundreds of kinases simultaneously. This study was designed as a series of hypotheses generating experiments. Although our experimental designs did not allow us to draw definitive conclusions, our experimental results fall into one of two major categories, which we have termed "reference kinases" and "neoteric kinases." The first category, "reference kinases," represents kinases with well-established roles in human cancer pathophysiology. The second category, "neoteric kinases," represents candidate kinases potentially contributing to PDAC pathology in new or previously understudied ways.

Reference Kinases
Reference kinases include protein tyrosine kinases identified by our study and subsequent bioinformatic analyses that recapitulate previously reported findings in the field of human cancer biology. Kinases in this category provide reference data that passively validate our experimental observations and contextualize our results within the scope of verified kinase discoveries.
One notable reference kinase identified as differentially active in our study of PDAC cells is the EGFR tyrosine kinase. Previously, EGFR has been linked to pancreatic tumor size, advanced clinical staging, and poor survival [28]. The expression frequency of EGFR in human pancreatic carcinomas is reported as 43% [28] and 68.4% in primary invasive ductal carcinoma of the pancreas [31] with elevated expression of EGFR activating ligands also reported. Consistent with these reports, our results show differential EGFR activity in weighted analyses of PDCL5 ( Figure 4, Tables 2 and 4) and aggregated patient-derived cell lines (Table 4). Directionality (increased kinase phosphorylation activity or decreased kinase phosphorylation activity) within each cell line can be gleaned from KRSA's report of the log2-fold change of phosphorylated peptide substrates attributed to each kinase family ( Figure A1, Figure A2, or Figure A3) or by UKA's report of an individual kinase's mean kinase statistic (Table A1). In our study, EGFR demonstrated increased phosphorylation activity in pancreatic cancer compared to control. Inhibition of the EGFR tyrosine kinase improves survival in PDAC animal models [25]. As such, targeted inhibition of tyrosine kinases, including EGFR, is a popular goal of many emerging therapeutic strategies [2]. Erlotinib (OSI-774) is an FDA-approved small-molecule EGFR tyrosine kinase inhibitor for use in pancreatic cancer [3]. Many additional small-molecule tyrosine kinase inhibitors are currently under study or in various stages of clinical trial for their putative role in pancreatic cancer pathophysiology.
The vascular endothelial growth factor receptor (VEGFR) tyrosine kinase family is also heavily implicated in the development of pancreatic cancer [32]. VEGFR-3, also known as fms related receptor tyrosine kinase 4 (FLT4), has been explored as a target for pancreatic cancer therapy [24,33] with significantly upregulated FLT4 expression documented in pancreatic cancer specimens [34,35]. Beyond PDAC, single nucleotide polymorphisms (SNPs) of FLT4 correlate with decreased progression-free survival of patients with gastroenteropancreatic neuroendocrine neoplasms [36]. Our results demonstrate increased FLT4 kinase activity in PDAC cells. In the present study, we identified FLT4 as one of the most differentially active kinases in PDCL5 patient-derived PDAC samples according to average percentile rankings across all pipelines (Table 3). KRSA and UKA directionality metrics ( Figure A3, Table A1) demonstrate increased FLT4 activity in pancreatic cancer compared to wild-type controls.
These results lend strength not only to our experimental and bioinformatic identification of differentially active kinases in pancreatic cancer, but also suggest that kinases identified as among the most strongly differential (e.g., in the top 10) in unweighted average percentile rankings (Table 1 or Table 3) and in the corresponding weighted average percentile rankings ( Table 2 or Table 4) may represent kinases highly likely to contribute to the pathophysiologic processes of PDAC.

Neoteric Kinases
Neoteric kinases represent a second, smaller category of experimental findings that include kinases whose identification in our study and bioinformatic analyses suggest new, hitherto unidentified, or otherwise understudied kinase functionalities in PDAC. Because these kinases are not strictly "novel," we instead call this group "neoteric" in reference to the emerging roles these kinases may play in pancreatic tumor desmoplasia, immune response, and oncometabolism. Beyond the passive validation that our "reference kinase" group provides, this group of "neoteric kinases" provides potentially novel insights. It became evident, after identifying lead candidate kinases, that our data highlight several potential players in unique aspects of PDAC tumor development. Kinases that appear both in our final unweighted average percentile rankings and in our final weighted average percentile rankings are defined as lead candidate kinases and include BLK (Figures 2 and 3, Tables 1 and 2), Fyn-related kinase (FRK) (Figure 3, Tables 1 and 2), Lck/Yes-related novel kinase (LYN) (Figures 2 and 3, Table 1, Table 2, and Table 4 Table 2), discoidin domain receptor 1 kinase (DDR1) (Table S1), and ephrin receptor A8 kinase (EPHA8) ( Table S1). While some of these kinases have been previously identified in kinome or phosphorylome studies of pancreatic cancer (e.g., DDR1 [9], FYN [68]), we classified them as neoteric for the purposes of this discussion because sufficient questions remain as to how these kinases relate to PDAC pathology, treatment, or molecular signaling.
Pronounced deposition of extracellular matrix constituents and the aberrant propagation of fibroblasts characteristic of desmoplasia are common in PDAC tumor microenvironments. Desmoplasia acts as a biophysical barrier contributing to pancreatic cancer therapeutic resistance. Desmoplastic stroma and pancreatic tumor cells interact with one another to elicit complex cellular behaviors with seemingly contradictory roles in PDAC progression [69]. At times pro-tumorigenic [70] and at times anti-tumoral [69], the role of desmoplasia in PDAC remains an active area of study. Our identification of increased HCK, ABL2, DDR1, FYN, and LYN suggests a role for these kinases in the desmoplastic reactions that contribute to the poor survival rates of PDAC patients. HCK overexpression activates fibrotic pathways [71]. ABL2 signaling regulates fibroblast proliferation [72]. DDR1 inhibitors reduce fibrosis in other fibrotic diseases [73][74][75][76]. FYN regulates downstream serine-threonine kinase activities involved in the modulation of fibroblast-epithelial cell interactions and the promotion of organ fibrosis [77,78]. Serotonin promotes fibroblast activation and collagen deposition [79]. LYN mediates pro-tumor serotonin signaling [80].
Desmoplasia serves as the primary source of the cytokines and chemokines that facilitate tumor progression in PDAC [81]. While immunotherapeutic strategies have significantly impacted clinical success in many other human malignancies, pancreatic cancer remains resistant. LCK is an important regulator of immune cell functionality [82]. TEC is a key player in the inflammatory response of pancreatitis [83]. The LCK and TEC kinases identified in the present study may also play a role in the anti-cancer immune response elicited and frequently evaded by PDAC.
The final effector proteins and terminal nodes for many kinase cascades are transcription factors. Our identification of BLK as a differentially active tyrosine kinase in PDAC cells compared to wild-type pancreatic cells presents new insight into the role that the pancreatic and duodenal homeobox 1 (PDX1) transcription factor plays in tumor progression. Overexpression of BLK induces an increase in the PDX1 transcription factor in the cytoplasm and nucleus [84]. The biological functionalities of the PDX1 transcription factor are multitudinous and context dependent. Our group has published extensively on the pro-tumorigenic role of the PDX1 transcription factor in pancreatic cancer progression [85][86][87][88][89][90][91], while other groups have demonstrated antimetastatic [92] and tumor-suppressive effects [93]. Recent evidence shows PDX1 functionality and its multiple-often antagonistic-effects on pancreatic cancer are stage-specific [93,94]. Our data suggest potential mechanistic relationships by which BLK kinase signaling cascades may contribute to PDX1 s multifaceted role in PDAC.
PDX1, originally known as insulin promoter factor 1 (IPF1), also serves as a transcriptional activator for metabolic genes such as insulin and glucose transporter type 2. Intersection between the PDX1 transcription factor, the BLK, and INSR kinases, as well as the role these kinases play in oncometabolic processes, will be examined in future studies. PDAC cells experience extreme deprivation of nutrient and oxygen delivery. Our data also implicate LYN, EPHA8, and FYN kinases as potential actors in oncometabolic PDAC signaling pathways and suggest mechanisms by which these kinases may facilitate oncogenic behavior.
Identification of these kinases in the present study contributes to growing understanding of abnormal fibrotic processes prominent in PDAC, dysregulated transcription factor activity, anti-cancer immune response, and the complex kinomic signaling networks responsible for pathometabolic tumor regulation and nutrient delivery. As novel PDAC subtypes continue to be defined [7,8,10,11,[95][96][97], it is clear that therapeutic strategies for PDAC must take different genetic backgrounds into account. The mutational profiles for the patient-derived cell lines used in this study (Table A2) provide useful insights into which protein tyrosine kinases may serve as effective targets of personalized/precision therapeutic intervention ( Figure A4).

Experimental Design
The experimental design is illustrated in Figure 1. In brief, PDAC epithelial cells and normal pancreatic ductal epithelial cells derived from patients were subjected to kinome array analysis using the PamStation 12 platform. All samples were prepared and assayed sequentially using Tyrosine Kinase PamChips consisting of 196 peptides with known phosphorylation sequences representing over 100 different proteins associated with the activity of upstream kinases. Each sample was assayed in triplicate with the results averaged across three identical kinome array runs.

Cell Lines and Patient-Derived Tissue
We used three different pancreatic cancer cell lines: one commercial cell line and two patient-derived cell lines (PDCL). We used commercial (ATCC CRL-1469) PANC1 cells originating from human pancreatic ductal cells carrying a TP53_R273H mutation and a KRAS_G12D mutation [98]. Two patient-derived cell lines (PDCL5, original name TKCC-05; PDCL-15, original name TKCC-15-Lo) were kindly provided by Andrew Biankin from Wolfson Wohl Cancer Research Centre, UK, with authentication by STR [90]. PDCL5 carries a TP53_G245S mutation and a KRAS_G12V mutation, while PDCL15 carries only a KRAS_G12D mutation (Table A2). These mutational profiles were kindly provided by Andrew Biankin and confirmed by F. Charles Brunicardi and Shi-He Liu. Normal patient-derived pancreatic ductal tissue was harvested from a healthy donor and kindly provided by Camillo Ricordi at the Diabetes Research Institute, University of Miami Miller School of Medicine, under the material transfer agreement. Because traditional two-dimensional cell culture models fail to accurately represent cancer microenvironments, we endeavored to obtain control tissue that more accurately represents physiological conditions. The pancreatic tissue contains ductal cells, acinar cells, and other elements included in the pancreatic microenvironment. While our decision to compare cell lines with wild-type pancreatic tissue may introduce a degree of bias into the study (cell lines and tissue samples contain different cellular contexts and environments), we used identical control wild-type tissue for each cell line comparison. All experiments and procedures were performed in strict compliance with all relevant laws and institutional guidelines.
Cell lines were cultured and lysed 72 h after plating, and tissue samples were processed as previously described [90]. All procedures were performed on ice. Tissue homogenization was performed using a D2400 Homogenizer and 1.5-mm Triple-Pure Zirconium Beads, with five rounds of homogenization and liquid nitrogen cooling to maintain low temperatures and minimize protein degradation. Each round of homogenization consisted of three cycles, with each cycle consisting of 30 s of active homogenization at 7 m/s and 30-second intervals between cycles. Tissue and cell lysate protein extractions were performed using M-PER (mammalian protein extraction reagent) (ThermoFisher, Waltham, MA, USA) and Halt Protease and Phosphatase Inhibitor Cocktails (ThermoFisher). Samples were centrifuged (14,000 RPM, 10 min, 4 • C) before supernatant collection. Total protein concentrations were assayed (Pierce BCA Protein Assay Kit, ThermoFisher) and samples were diluted to 1 µg/µL. All samples were prepared and measured simultaneously. Because freeze-thaw cycles decrease kinase activity [99], multiple aliquots were stored at −80 • C to minimize freeze-thaw cycles, with frozen aliquots used only once for kinome array assays.

Tyrosine Kinase Array
Tyrosine kinase activity was measured with the PamStation 12 instrument (PamGene International, 's-Hertogenbosch, The Netherlands) and PTK (4-well) array PamChips using fluorescently labeled antibodies to detect differential phosphorylation of 196 reporter peptides (including three internal controls) per well. These 196 consensus phosphopeptide sequences were immobilized on porous ceramic membranes. The PamChip wells were blocked with 2% bovine serum albumin (BSA) prior to addition of 1 µg of protein suspended in manufacturer's kinase buffer (PamGene). Next, we added 157 µM adenosine triphosphate (ATP) and FITC-labeled anti-phospho tyrosine antibodies (PamGene) to each well. Homogenized lysates containing active kinases and assay solution were pumped back and forth through PamChip wells in order to facilitate interactions between the active kinases and the 196 immobilized consensus phosphopeptide sequences. Evolve (PamGene) software captured FITC-labeled anti-phospho-antibodies bound to the phosphorylated consensus sequences. Image capture occurred every six seconds for 60 min. After washing, peptide signal intensity was recorded across several exposure times (10, 20, 50, 100, 200 milliseconds). The linear regression slope was calculated in order to provide the peptide phosphorylation intensity signal used in downstream comparative analyses. Signal ratios between pairs of samples were used to calculate fold change (FC) for each peptide. Differential peptide signals greater than or equal to 30% (FC ≥ 1.30 or FC ≤ 0.70) were considered demonstrative of minimum threshold changes in the degree of phosphorylation. This threshold value derived from conservative interpretation of previous literature suggesting even smaller orders of magnitude are sufficiently correlated with biologically relevant signaling changes [16,100,101]. Nonlinear (R 2 values less than 0.90) or undetectable peptides in the post-wash phase were not selected for further analysis. Kinome assays were performed in triplicate with the calculated FC per peptide averaged across three replicates.

Upstream Kinase Identification
Kinase activity corresponded to the degree of consensus peptide phosphorylation as measured using real-time Evolve kinetic image capture software. The raw data generated by the PamStation platform were minimally processed (using threshold changes described above) to generate a list of differentially phosphorylated peptide sequences, which served as input for subsequent bioinformatic analyses. To expand coverage of peptide sequences and maximize identification of candidate upstream kinases, we used four distinct bioinformatic pipelines, each with a semi-overlapping set of databases by which their respective algorithms queried our list of differentially phosphorylated peptide sequences. Because each pipeline also relied upon a unique pipeline-specific algorithm, evaluating the results of one pipeline within the context of the results obtained through the other three pipelines allowed for enhanced perspective. Shared identification of upstream kinases responsible for the observed peptide phosphorylation patterns between pipelines could, therefore, be weighted and integrated into the final analysis. These pipelines include (1) Upstream Kinase Analysis (UKA) from PamGene, (2) Post-Translational Modification Signature Enrichment Analysis (PTM-SEA) from the Broad Institute of MIT and Harvard, (3) Kinase Enrichment Analysis Version 3 (KEA3) from the Ma'ayan laboratory, and (4) Kinome Random Sampling Analyzer (KRSA) developed by our own laboratory [12][13][14][15][16][17]. Meaningful differences between pipelines were compared in the Discussion section, while the methods by which we deployed each pipeline are presented below.

Upstream Kinase Analysis (UKA) Pipeline
Using PamGene's BioNavigator software, we ran our data sets through the Protein Tyrosine Kinase (PTK) Upstream Kinase Analysis (UKA) Knowledge Integration PamApp. These data sets included (1) mean phosphorylated peptide sequences of PDCL15 vs. wild-type, (2) mean phosphorylated peptide sequences of PDCL5 vs. wild-type, and (3) mean phosphorylated peptide sequences of PANC1 vs. wild-type. Sample names served as the factor uniquely defining each observation. The "treatment off chip" factor served as the factor defining the experimental groupings. The default scan rank (4 to 12) and permutation (500) parameters were applied to the analysis, in addition to the default in vitro/in vivo (1), in silico (PhosphoNET) (1), minimal sequence homology (0.9), minimal PhosphoNET prediction score (300), and minimal peptide set (3) parameters. Inclusive percentile ranks were calculated according to the absolute value of UKA's Median Final Score output.

Post-Translational Modification Signature Enrichment Analysis (PTM-SEA) Pipeline
We used the Broad Institute's Single sample Gene Set Enrichment Analysis (ssGSEA) and Post-Translational Modification Signature Enrichment Analysis (PTM-SEA) publicly available repository (https://github.com/broadinstitute/ssGSEA2.0), RStudio Desktop 1.2.5042 (https://rstudio. com), and underlaying R 3.3.1 software environment (https://cran.rstudio.com). We ran our data sets through the PTM-SEA pipeline after modifying the peptide database to include only the peptide sequences present on the PamChip plus all peptide sequences with minimal sequence homologies of 0.9. These data sets include (1) mean phosphorylated peptide sequences of PDCL15 vs. wild-type, (2) mean phosphorylated peptide sequences of PDCL5 vs. wild-type, and (3) mean phosphorylated peptide sequences of PANC1 vs. wild-type. Results were filtered to include only protein tyrosine kinases. PTM-SEA inclusive percentile ranks were determined according to each kinase's respective 1/fdr.pvalue.totalGeoMeanLFC output value.

Kinase Enrichment Analysis Version 3 (KEA3) Pipeline
We used the Ma'ayan laboratory's Kinase Enrichment Analysis Version 3 (KEA3) (https://amp. pharm.mssm.edu/kea3/#) to process our data sets. These data sets include (1) mean phosphorylated peptide sequences of PDCL15 vs. wild-type, (2) mean phosphorylated peptide sequences of PDCL5 vs. wild-type, and (3) mean phosphorylated peptide sequences of PANC1 vs. wild-type. To accommodate the input parameters of KEA3, these peptide sequences were converted to HGNC-approved gene symbols before the data sets were entered into the KEA3 pipeline. Results were filtered to include only protein tyrosine kinases. Average FDR p-values from 0.2, 0.3, and 0.4 LFC cutoff input lists were averaged according to ChengKSIN, PTMsigDB, or PhosDAll database outputs. The resultant ChengKSIN, PTMsigDB, and PhosDAll average values were themselves averaged and −log10 transformations of these averages were used for inclusive percentile ranking calculations.
To accomplish this, we mapped phosphorylation sites within the reporter peptides to individual protein kinases that phosphorylate these sites. To this end, multiple databases were queried including GPS 5.0 (http://gps.biocuckoo.cn), Kinexus Phosphonet (http://www.kinexus.ca), PhosphoELM (http: //phospho.elm.eu.org), and PhosphoSite Plus (https://www.phosphosite.org). In this way, peptide sequences and kinases were matched such that ranked predictions were generated to identify tyrosine kinases most likely to have produced the observed phosphorylation results. Kinases with scores greater than twice the prediction threshold for each phosphorylation site in the GPS 5.0 database were carried forward for downstream analysis. The top five kinase predictions in the Kinexus database were also carried forward for downstream analysis. Additional kinases reported to act on specific phosphorylation sites were identified using PhosphoELM and Phosphosite Plus public databases and carried forward for downstream analysis. Downstream analysis consists of 3000 random sampling iterations in which an equal number of differentially phosphorylated peptide sequences are randomly selected from the total list of 196 phosphopeptide sites on the tyrosine kinase PamChip. Predicted kinases were then mapped to each iteration in order to generate comparative controls to which the experimentally generated kinase lists could be compared. This allowed meaningful approximations of the direction of activity (increased or decreased) and significance for each experimentally identified kinase. Inclusive percentile ranks were calculated according to mean LFC output. KRSA is publicly available at https://github.com/kalganem/KRSA. Because KRSA identifies kinase families rather than individual kinases, all kinases within a family were given identical scores.

Combinatory Analyses
To resolve different output metrics, the results of each respective pipeline were converted to inclusive percentile rankings. These inclusive percentile rankings were aggregated according to cell line and protein tyrosine kinase and averaged per cell line (PANC1, PDCL15, PDCL5) or per cell line group (Patient-derived cell line group: PDCL15 and PDCL5; All group: PANC1, PDCL15, PDCL5). Weighted averages were calculated by dividing average percentile rankings by the number of pipelines under consideration.

Conclusions
This study provides evidence in support of previously reported kinases in human cancer with an emphasis on PDAC. This passive validation supports the strength of ongoing drug development strategies that target protein tyrosine kinases and propounds the utility and accuracy of peptide-based kinomic analytical platforms. Furthermore, our identification and contextualization of candidate or lead candidate kinases responsible for the differential phosphorylation signatures observed between PDAC commercial or patient-derived cell lines compared to wild-type pancreatic patient samples provides evidence of unique kinomic relationships between pancreatic tumor cells and the desmoplastic stromal environments that support tumor progression and cause significant obstacles in pancreatic cancer therapy. Identification of the BLK, HCK, FRK, ABL2, DDR1, LYN, EPHA8, FYN, LCK, and TEC kinases as potentially significant mediators of pancreatic cancer progression and fibrotic/desmoplastic development fits well into established knowledge while also advancing new avenues of drug development and discovery. Additionally, our data provide increased understanding of the relationship between BLK protein tyrosine kinase, PDX1 transcription factor, and pancreatic disease. This study also outlines additional mechanisms by which HCK, ABL2, and DDR1 may play a role in pancreatic cancer and fibrosis. These results also support the role of LYN in oncometabolic processes and posit pathways by which LYN, EPHA8, and FYN may facilitate oncogenic cellular behavior. Lastly, we provide a rationale for continued investigation of the complex interplay between anti-cancer immune response and the activity of LCK and TEC kinases. These findings are summarized in Figure 5, and our companion review piece provides additional information on the role of BLK, HCK, FRK, ABL2, DDR1, LYN, EPHA8, FYN, LCK, and TEC kinases in PDAC and pancreatic cancer desmoplasia [107].   Appendix A Figure A1. Differential phosphorylation levels of peptide sequences attributed to kinase family activity in PANC1 vs. wild-type control. Each red dot in a column represents one peptide sequence whose phosphorylation is performed by that column's kinase family. The y axis reports log2-fold change and the x axis identifies the kinase family. Black dots represent peptides that did not demonstrate differential phosphorylation in PANC1 samples compared to wild-type samples, with the dashed horizontal line representing a positive or negative 0.2 log2-fold change cutoff. Red dots above these lines represent peptides that are more phosphorylated in PANC1 compared to control. Red dots that are below these lines represent peptides that are less phosphorylated in PANC1 compared to control. Figure A2. Differential phosphorylation levels of peptide sequences attributed to kinase family activity in PDCL15 vs. wild-type control. Each red dot in a column represents one peptide sequence whose phosphorylation is performed by that column's kinase family. The y axis reports log2-fold change and the x axis identifies the kinase family. Black dots represent peptides that did not demonstrate differential phosphorylation in PDCL15 samples compared to wild-type samples, with the dashed horizontal line representing a positive or negative 0.2 log2-fold change cutoff. Red dots above these lines represent peptides that are more phosphorylated in PDCL15 compared to control. Red dots that are below these lines represent peptides that are less phosphorylated in PDCL15 compared to control. Figure A3. Differential phosphorylation levels of peptide sequences attributed to kinase family activity in PDCL5 vs. wild-type control. Each red dot in a column represents one peptide sequence whose phosphorylation is performed by that column's kinase family. The y axis reports log2-fold change and the x axis identifies the kinase family. Black dots represent peptides that did not demonstrate differential phosphorylation in PDCL5 samples compared to wild-type samples, with the dashed horizontal line representing a positive or negative 0.2 log2-fold change cutoff. Red dots above these lines represent peptides that are more phosphorylated in PDCL5 compared to control. Red dots that are below these lines represent peptides that are less phosphorylated in PDCL5 compared to control.   [7,9,11,108] Data Repo Appendix B  ABL1  KEA3  NTRK1  NTRK1  KEA3  FLT3  FLT3  KEA3  DDR2  DDR2  KEA3  KIT  KIT  KEA3  PDGFRA  PDGFRA  KEA3  MATK  MATK  KEA3  EPHB3  EPHB3  KEA3  MST1R  MST1R  KEA3  FES  FES  KEA3  FLT4  FLT4  KEA3  SRC  SRC  KEA3  TXK  TXK  KEA3  NTRK3  NTRK3  KEA3  KDR  KDR  KEA3  RET  RET  KEA3  LCK  LCK  KEA3  ABL1  ABL1  KEA3  EPHA2  EPHA2  KEA3  SRMS  SRMS  KEA3  EPHB2  EPHB2  KEA3  FYN  FYN  KEA3  EGFR  EGFR  KEA3  FLT1  FLT1  KEA3 FER FER Table A4. Cont. KEA3  INSR  INSR  KEA3  FGFR4  FGFR4  KEA3  ITK  ITK  KEA3  EPHB1  EPHB1  KEA3  CSF1R  CSF1R  KEA3  PTK6  PTK6  KEA3  CSK  CSK  KEA3  ERBB2  ERBB2  KEA3  NTRK2  NTRK2  KEA3  TYRO3  TYRO3  KEA3  BTK  BTK  KEA3  JAK2  JAK2  KEA3  SYK  SYK  KEA3  LYN  LYN  KEA3  FGFR3  FGFR3  KEA3  PTK2  PTK2  KEA3  FGR  FGR  KEA3  ERBB4  ERBB4  KEA3  YES1  YES1  KEA3  ZAP70  ZAP70  KEA3  JAK3  JAK3  KEA3  MET  MET  KEA3  IGF1R  IGF1R  KEA3  TEC  TEC  KEA3  AXL  AXL  KEA3  ALK  ALK  KEA3  PTK2B  PTK2B  KEA3  PDGFRB  PDGFRB  KEA3  STYK1  STYK1  KEA3  MERTK  MERTK  KEA3  BMX  BMX  KEA3  EPHA3  EPHA3  KEA3  ABL2  ABL2  KEA3  FGFR1  FGFR1  KEA3  EPHA4  EPHA4  KEA3  TYK2  TYK2  KEA3  FRK  FRK  KEA3  FGFR2  FGFR2  KEA3  TNK2  TNK2  KEA3  JAK1  JAK1  KEA3  DDR1  DDR1  KEA3  BLK  BLK  KEA3  HCK  HCK  KEA3  EPHA8  EPHA8  KEA3 TEK TEK Figures 2-4, expanded figure legend. Default outputs from upstream kinase identification pipelines for PDAC cell lines compared to patient-derived wild-type pancreatic tissue. (A) The Kinome Random Sampling Analyzer (KRSA) pipeline identifies differentially active upstream kinases according to mean standard score (Z) values on the x axis, with kinase family on the y axis, and color gradation representing absolute mean standard score. Full red represents an absolute mean standard score above 3 and full white represents an absolute mean standard score below 1. Absolute mean standard scores for each kinase family are calculated by averaging the absolute standard score values of nine standard scores, each derived from input lists consisting of peptide sequences with log-fold change (LFC) in phosphorylation above 0.2, 0.3, or 0.4 in one of our three biological replicates. Each plotted box represents values within and including the first and third percentiles (the 25th and 75th percentiles). Upper or lower whiskers extend from the box to the largest or smallest kinase family standard score that is less than or equal to 150% of the distance between the first and third percentiles, respectively. Outliers defined as greater than 150% of the distance between the first and third percentiles are plotted individually. Dashed lines delineate positive values above 1.5, 1.7, and 2.0 or negative values below −1.5, −1.7, and −2.0. Kinase families appear in descending order of standard score with the most significant differential kinase family activity appearing at the top and bottom of the list. Any zero-value standard score represents a kinase family having identical mean log-fold change values between the phosphorylation levels of representative peptide sequences measured in the PDAC cell line group and the phosphorylation levels of those same representative peptide sequences when measured in the control wild-type patient sample group. (B) The Post-Translational Modification Signature Enrichment Analysis (PTM-SEA) pipeline identifies differentially active upstream kinases according to negative decadic logarithms of the adjusted probability values. These values allow the most significant differentially active kinases to be listed from top to bottom in order of descending significance. Adjusted probability values less than 0.05 are, therefore, represented by x axis values above 1.30. Adjusted probability values correcting for multiple comparisons use the Bonferroni correction method on comparisons between the phosphorylated peptide sequences of the experimental PDAC cell line group and the phosphorylated peptide sequences of the control wild-type patient group in order to determine the extent to which these differences may be attributed to each individual kinase listed on the left. (C) The Kinase Enrichment Analysis Version 3 (KEA3) pipeline also identifies differentially activated upstream kinases according to negative decadic logarithms of the average adjusted probability values and allows the most significant differentially active kinases to be listed from top to bottom in order of descending significance. Adjusted probability values less than 0.05 are, therefore, represented by x axis values above 1.30. Adjusted probability values less than 0.01 are represented by x axis values above 2.

Pipeline Default Standardized
Adjusted probability values correcting for multiple comparisons use the false discovery rate (FDR) correction method to compare the genes containing phosphorylated peptide sequences of the experimental PDAC cell line group with the genes containing phosphorylated peptide sequences of the control wild-type patient group in order to determine the extent to which these differences may be attributed to each individual kinase listed on the left. (D) The Upstream Kinase Analysis (UKA) pipeline identifies differentially activated upstream kinases according to the normalized kinase statistic for change in phosphorylation between the peptide sequences of the experimental PDAC cell line group and the control wild-type patient group. This normalized kinase statistic is calculated by subtracting the mean phosphorylation levels of the peptide sequences associated with a given kinase in the control group from the mean phosphorylation levels of the peptide sequences associated with that same kinase in the experimental group and dividing the resulting value by the standard deviation. A quotient above 0 indicates increased activity of a kinase from the experimental group. A quotient below 0 indicates decreased activity of a kinase from the experimental group. Each graphed point represents the normalized kinase statistic when the UKA algorithm uses only one of the pipeline's available databases. The diameter of a circle increases with the size of the peptide sequences associated with the corresponding kinases in the corresponding database. Color gradation represents specificity scores with full red representing values of two and full black representing specificity scores of zero. Specificity scores are informed by the quotient produced with a dividend corresponding to the number of times that a permuted normalized kinase statistic greatly exceeded the nonpermuted normalized kinase statistic and a divisor corresponding to the number of random permutations, or the quotient produced by the number 1 as a dividend and the number of random permutations as the divisor. Specificity scores represent the negative decadic logarithm of whichever quotient is larger. Individual kinases are listed along the y axis from top to bottom in descending order of median "final score." UKA calculates the "final score" as the summation of the significance score and specificity score of each individual kinase. (E) Comparison of upstream kinase activity as determined by KRSA, PTM-SEA, KEA3, and UKA, each listed on the y axis. Kinase family names appear along the top of the x axis. Individual kinase family members appear along the bottom of the x axis. The comparison includes only kinases identified as one of the top 10 most differentially active upstream kinases by any one of the pipelines. Kinases with identical scores were arbitrarily assigned sequential ranks. However, if two or more kinases tied for the last (10th) position in this list, then the list was extended so as not to arbitrarily exclude an equivalently ranked kinase. The diameter of each black circle corresponds with the percentile in which the labeled kinase appears in the labeled pipeline. The kinase family scores of the KRSA pipeline are repeated for each individual kinase. A white circle demonstrates not only absence from a given pipeline's top 10 list, but also demonstrates absence from a given pipeline's unfiltered exhaustive output.