A Predictive 7-Gene Assay and Prognostic Protein Biomarkers for Non-small Cell Lung Cancer.

PURPOSE: This study aims to develop a multi-gene assay predictive of the clinical benefits of chemotherapy in non-small cell lung cancer (NSCLC) patients, and substantiate their protein expression as potential therapeutic targets. PATIENTS AND METHODS: The mRNA expression of 160 genes identified from microarray was analyzed in qRT-PCR assays of independent 337 snap-frozen NSCLC tumors to develop a predictive signature. A clinical trial JBR.10 was included in the validation. Hazard ratio was used to select genes, and decision-trees were used to construct the predictive model. Protein expression was quantified with AQUA in 500 FFPE NSCLC samples.
RESULTS: A 7-gene signature was identified from training cohort (n = 83) with accurate patient stratification (P = 0.0043) and was validated in independent patient cohorts (n = 248, P < 0.0001) in Kaplan-Meier analyses. In the predicted benefit group, there was a significantly better disease-specific survival in patients receiving adjuvant chemotherapy in both training (P = 0.035) and validation (P = 0.0049) sets. In the predicted non-benefit group, there was no survival benefit in patients receiving chemotherapy in either set. The protein expression of ZNF71 quantified with AQUA scores produced robust patient stratification in separate training (P = 0.021) and validation (P = 0.047) NSCLC cohorts. The protein expression of CD27 quantified with ELISA had a strong correlation with its mRNA expression in NSCLC tumors (Spearman coefficient = 0.494, P < 0.0088). Multiple signature genes had concordant DNA copy number variation, mRNA and protein expression in NSCLC progression.
CONCLUSIONS: This study presents a predictive multi-gene assay and prognostic protein biomarkers clinically applicable for improving NSCLC treatment, with important implications in lung cancer chemotherapy and immunotherapy.

Introduction

Lung cancer is the leading cause of cancer-related deaths in the world, and non-small cell lung cancer (NSCLC) accounts for almost 80% of lung cancer deaths [1]. The heterogeneous nature of lung cancer makes it a very difficult disease to treat. Major histology of NSCLC includes lung adenocarcinoma and squamous cell lung carcinoma. Surgical resection is the major treatment for early stage NSCLC. However, about 22–38% of stage I NSCLC patients will develop tumor recurrence within five years following the surgery [2]. It is therefore important to select early stage NSCLC patients for more aggressive treatment. While adjuvant chemotherapy of stage II and stage III disease has resulted in 10–15% increased overall survival [3], the prognosis for early stage NSCLC remains poor [4]. Currently, there are no clinically available molecular assays to predict the risk for tumor recurrence and the clinical benefits of chemotherapy in NSCLC patients.

Immunotherapy has rapidly gained attention of oncologists as an effective and less toxic treatment than chemotherapy in patients with advanced lung cancers [[5], [6], [7], [8]]. A recent study used paired single cell analysis to compare normal lung tissue and blood with tumor tissue in stage I NSCLC, and found that early-stage tumors had already begun to alter the immune cells in their microenvironment [8]. These results suggest that immunotherapy could potentially be used to treat early stage lung cancer patients. However, predictive biomarkers of immunotherapy are not well established except PD-1 or PD-L1, and it is unlikely that a single marker is sufficient.

High-throughput technologies, such as microarray and RNA-seq, promise the discovery of novel biomarkers from genome-scale studies. The FDA conducted a systematic evaluation and suggested continued usefulness of legacy microarray data and established microarray biomarkers and predictive models in the forthcoming RNA-seq era [9]. However, several disadvantages have limited the application of high-throughput techniques in routine clinical tests, including costs, reproducibility, and data analyses [10]. Compared with microarray and RNA-seq, quantitative real-time RT-PCR (qRT-PCR) is more efficient, consistent, and able to measure gene expression over a greater dynamic range. It requires only a small sample and can be modified to quantify gene expression in formalin-fixed and paraffin-embedded (FFPE) tissues [11]. The combined use of real-time qRT-PCR with high-throughput analysis can overcome the inherent biases of the high-throughput techniques and is emerging as the optimal method of choice to translate genome research into clinical practice [12]. The protein expression validation of the identified mRNA biomarkers could substantiate their ultimate functional involvements in disease, and may lead to the discovery of potential proteomic biomarkers in abundant FFPE samples for broader applications in community hospitals.

DNA microarray-based studies identified gene expression-based NSCLC prognostic [13] and predictive biomarkers [14,15]. A qRT-PCR based 14-gene assay by Kratz et al. [16] is prognostic of non-squamous NSCLC outcome in FFPE tissues and is ready for wide-spread clinical applications. However, this 14-gene assay is limited to non-squamous NSCLC and is not shown to be predictive of the clinical benefits of chemotherapy.

In this study, a combined analysis of genome-wide transcriptional profiles and qRT-PCR was utilized to develop a multi-gene assay both prognostic of NSCLC outcome and predictive of the benefits of chemotherapy. Patient cohorts from multiple hospitals in the US and JBR.10 data [14] were used to validate this multi-gene assay. Protein expression of the identified biomarkers was also evaluated in patient tissue samples and correlated with the mRNA expression and DNA copy number variation to substantiate their functional involvement and potential as therapeutic targets in chemotherapy and immunotherapy, in addition to companion tests.

Materials and Methods

Patient Samples

Clinical characteristics of patient cohorts used in qRT-PCR assays is summarized in Table 1. All NSCLC patients were staged I, II, or IIIA at the time of diagnosis. Tumor tissues were collected in surgical resections and were snap-frozen at −80 °C until used for RNA extraction. Tumor cell content was above 50% for qRT-PCR assays. Those with missing AJCC staging information, missing histology, death within 30 days of resection or from other disease conditions were excluded from further analysis. A total of 122 NSCLC patient samples were obtained from Case Western Reserve University (CWRU) Comprehensive Cancer Center. Total RNA of good quality was extracted from 89 tumor specimens. Good quality RNA from 101 lung adenocarcinoma tumor specimens was obtained from University of Michigan (UM) Comprehensive Cancer Center, with detailed description of patients, tissue specimens and mRNA quality check provided in [17]. A total of 65 NSCLC tumor specimens from NorthShore University HealthSystem Kellogg Cancer Center and 49 specimens from West Virginia University Cancer Institute [Mary Babb Randolph Cancer Center (MBRCC)] generated good quality mRNA. The tissue collection in this study was approved by an Institutional Review Board (IRB) at each institution.

Table 1

Clinical information of non-small cell lung cancer patient cohorts collected for the qRT-PCR analysis.

		CWRU (n = 89)	MBRCC (n = 49)	UM (n = 101)	NorthShore (n = 65)
Age	Mean (Std error)	70.11 (0.94)	66.70 (1.25)	67.04 (0.96)	69.64 (1.02)
	<60	15 (15.15%)	7 (14.29%)	28 (27.72%)	7 (10.77%)
	≥60	84 (84.85%)	39 (79.59%)	73 (72.28%)	48 (73.85%)
	Missing		3 (6.12%)		10 (15.38%)
Sex	F	52 (52.53%)	23 (46.94%)	53 (52.48%)	34 (52.31%)
	M	47 (47.47%)	26 (53.06%)	48 (47.52%)	21 (32.31%)
	Missing				10 (15.38%)
Smoking	Current	43 (43.43%)	1 (2.04%)		Yes	60 (92.31%)
	Former	40 (40.40%)	3 (6.12%)
	Never	8 (8.08%)			No	5 (7.69%)
	Passive	1 (1.01%)
	Other	1 (1.01%)
	Missing	6 (6.06%)	45 (91.48%)
AJCC stage	I	46 (46.46%)	27 (55.10%)	59 (58.42%)	46 (70.77%)
	II	46 (46.46%)	16 (32.65%)	16 (15.84%)	15 (23.08%)
	III	6 (6.06%)	6 (12.25%)	26 (25.74%)	4 (6.15%)
	Missing	1 (1.01%)
Chemotherapy	Yes	29 (29.29%)	27 (55.10%)	24 (23.76%)	28 (40.03%)
	No	52 (52.53%)	20 (40.82%)	77 (76.24%)	36 (55.38%)
	Missing	13 (13.13%)	2 (4.08%)		1 (1.54%)
Histology	Adenocarcinoma	65 (65.66%)	27 (55.10%)	101 (100%)	43 (66.15%)
	Squamous	27 (27.27%)	14 (28.57%)		11 (16.92%)
	Other	7 (7.07%)	8 (16.33%)		6 (9.23%)
	Missing	5 (5.05%)			5 (7.69%)
Differentiation	Well	5 (5.05%)		28 (27.72%)	20 (30.77%)
	Moderate	44 (44.44%)			4 (6.15%)
	Moderate to Poorly	4 (4.04%)		39 (38.61%)	22 (33.85%)
	Poorly	35 (35.35%)		34 (33.66%)	17 (26.15%)
	Missing	11 (11.11%)			2 (3.08%)
Tumor Grade	1	5 (5.05%)	3 (6.12%)		20 (30.77%)
	2	44 (44.44%)	18 (36.73%)		19 (29.23%)
	3	36 (36.36%)	22 (4.90%)		21 (32.31%)
	Other	3 (3.03%)
	Missing	11 (11.11%)	6 (12.25%)		5 (7.69%)

RNA Extraction, and Quality and Concentration Assessments

Total RNA was extracted from snap-frozen tumor tissues using a RNeasy mini kit according the manufacturer's protocol (Qiagen, USA), followed by elution in 30 μL of RNase-free water and storage at −80 °C. The quality and integrity of the RNA, the 28S to 18S ratio, and a visual image of the 28S and 18S bands were evaluated on the 2100 Bioanalyzer (Agilent Technologies, CA). RNA assessed as having good quality from 304 tumor samples was included for further analysis. The RNA concentration of each sample was assessed using a Nanodrop-1000 Spectrophotometer (NanoDrop Tech, Germany).

Generation of Complementary DNA (cDNA)

The reverse transcriptase polymerase chain reaction was used to convert the high-quality single-stranded RNA samples to double-stranded cDNA, using an Applied Biosystems GeneAmp® PCR 9700 machine (Foster City, CA). For standardization across all samples, one microgram of RNA was used to generate cDNA.

Real-time Quantitative RT-PCR Low-density Arrays

Real-time qRT-PCR assays of independent patient cohorts of NSCLC tumor samples were used to further select biomarkers to form a multi-gene assay from prognostic genes identified from microarray data in our previous studies [[18], [19], [20], [21]]. The identified prognostic genes were initially validated with multiple independent NSCLC microarray data publically available [[18], [19], [20], [21]]. Based on the validation results, 160 prognostic genes and three housekeeping genes were included in the qRT-PCR experiments. The three housekeeping genes were 18S, UBC, and POLR2A due to their confirmed constant mRNA expressions across samples [18].

We analyzed 337 tumor samples with good RNA quality using TaqMan microfluidic low-density array (LDA) plates on an ABI 7900HT Fast RT-PCR instrument (Applied Biosystems). Total RNA samples were analyzed on an Agilent 2100 Bioanalyzer RNA 6000 Nano LabChip. The report was generated by the SDS2.3 software (Applied Biosystems). In the report, the number of cycles required to reach threshold fluorescence (Ct) and ∆CT for each sample relative to the control gene defines the expression pattern for a gene. The gene expression data were further analyzed using the 2− method [22].

Statistical and Computational Analysis

Prognostic biomarkers were evaluated with Cox proportional hazard model. Hazard ratio was used in the evaluation of prognostic performance of biomarkers. If a biomarker gives a hazard ratio > 1, it means that patient samples predicted as high risk are more likely to have a poor outcome. In the evaluation of genes in qRT-PCR assays, ∆C was used as a covariate in Cox model. If a gene as a hazard ratio > 1, it means that down-regulation of this gene is associated with a poor outcome and up-regulation of this gene is associated with a good outcome in NSCLC patients; otherwise, if a gene has a hazard ratio < 1, it means that down-regulation of this gene is associated with a good outcome and up-regulation of this gene is associated with a poor outcome in NSCLC patients. During the evaluation, UBC (Hs00824723_m1) was chosen as the house keeping gene to normalize gene expression. The CWRU cohort was used as the training set, and seven genes were selected to form a prognostic classifier based on decision trees. These seven genes are ABCC4 (Hs00988717_m1), CCL19 (Hs00171149_m1), SLC39A8 (Hs00223357_m1), CD27 (Hs00154297_m1), FUT7 (Hs00237083_m1), DAG1 (Hs00189308_m1), and ZNF71 (Hs00221893_m1). The 7-gene prognostic model was validated with independent patient cohorts (UM, MBRCC, and NorthShore). In Kaplan-Meier analysis, log-rank tests or Wilcoxon tests were used to assess the difference in probability of survival of different prognostic groups. All the analyses were performed with packages in R or SAS unless otherwise specified.

Validation on Clinical Trial JBR.10

Data from JBR.10 was obtained from NCBI Gene Expression Omnibus with accession number GSE14814 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=gse14814). A total of 133 non-small cell lung cancer samples were profiled for gene expression using Affymetrix 133A platform [14]. Patients were all in early stage (I or II). Patient samples assayed in the same batch with consecutive accession numbers ranging from GSM370913 to GSM371002 (n = 90) were used in the validation of the 7-gene signature. Among these patient samples, those who died from other disease conditions were excluded from further analysis. ABCC4 (203196_at), CCL19 (210072_at), CD27 (206150_at), DAG1 (205417_s_at and 212128_s_at), FUT7 (210506_at and 217696_at), SLC39A8 (209266_s_at, 209267_s_at, 216504_s_at, and 219869_s_at), and ZNF444 (218707_at and 50376_at) were used in validating the qRT-PCR based multi-gene assay. For a gene with multiple probe sets, the one with the highest expression value (yielding the clearest signal) in each sample was chosen to represent the gene expression. ZNF71 was not available in the GSE14814 dataset. ZNF444 was chosen to replace ZNF71 to validate the qRT-PCR results, because both ZNF444 and ZNF71 are at locus NC_000019.10 in Chromosome 19 and belong to zinc finger protein family. To be compatible with the ∆Ct values in qRT-PCR data, log2 transformed microarray data was used in the analysis, and the expression values of UBC minus those of selected probes were used in the normalization of the microarray data.

Tissue Microarrays (TMA)

Samples from 2 retrospective collections of lung cancer were examined in TMA format from Yale University Pathology Archives; Cohort A (YTMA 250 [n = 298]) and Cohort B (YTMA 79 [n = 202]). TMAs consisted of 0.6 mm cores in 1 (Cohort A) and 2 fold (Cohort B) redundancy. TMAs were prepared according to standard methods. Cohort A comprises 314 serially collected NSCLC who underwent surgical resection of their primary tumor between 2004 and 2011. Cohort B comprises of 202 serially collected NSCLC patients who underwent surgical resection of their primary tumor between 1988 and 2003. All tissue was used after approval from the Yale Human Investigation Committee protocol #9505008219, which approved the patient consent forms or in some cases waiver of consent. The actual number of samples analyzed for each study is lower, due to unavoidable loss of tissue or the absence or limited tumor cells in some spots as is commonly seen in TMA studies. NSCLC patients in stage I, II, and IIIA were included in the analysis. Those who died with no evidence of disease were excluded from further analysis.

Quantitative Immunofluorescence

FFPE whole-tissue sections, tissue microarrays (TMAs) and cell pellets were processed at Yale Cancer Center/Pathology Tissue Microarray Facility with details provided in Supplementary File 1.

Primary antibodies were followed by incubation with Alexa 546–conjugated goat anti-mouse secondary antibody (Life Technologies) diluted 1:100 in rabbit EnVision reagent (Dako) for 1 h. ZNF71 signal was amplified with Cy5-Tyramide (Perkin Elmer) for 10 min, and then nuclei were stained with 0.05 mg DAPI in BSA-tween for 10 min. Slides were mounted with ProlongGold (Life Technologies). Two TBS-T and one TBS wash was performed between each step after the primary antibody.

Immunofluorescence was quantified using automated quantitative analysis (AQUA) Fluorescent images of DAPI, Cy3 (Alexa 546-cytokeratin), and Cy5 (ZNF71) for each TMA spot were collected. Image analysis was carried out using the AQUAnalysis software (Navigate Biopharma Inc.), which generated an AQUA score for each compartment by dividing the sum of target pixel intensities by the area of the compartment in which the target is measured. AQUA scores were normalized to the exposure time and bit depth at which the images were captured, allowing scores collected at different exposure times to be directly comparable. Specimens with <5% tumor area per region of interest were not included in AQUA analysis for not being representative of the corresponding tumor specimen.

Enzyme-Linked Immunosorbent Assay (ELISA)

A total of 38 NSCLC patient tissue samples were selected for ELISA assays, including 29 tumor resections of lung adenocarcinoma and squamous cell lung cancer and 9 matched adjacent normal lung tissue samples. The DuoSet ELISA Development Systems from R&D Systems (Minneapolis, MN; catalog number: DY382-05) were used for quantifying protein expression of T-Cell Activation Antigen CD27 (CD27)/Tumor Necrosis Factor Receptor Superfamily, Member 7 (TNFRSF7) in NSCLC patient tissue samples, according to manufacturer's protocol. The ELISA assay results were quantified using the Synergy H1 Hybrid Multi-Mode Microplate Readers from BioTek Instruments, Inc. (Winooski, VT). Samples that yielded a positive OD values were included for further analysis. Statistical analysis was done using a two-sample t-test assuming unequal variances. The concordance between CD27 mRNA and protein expression was evaluated with Spearman correlation coefficient.

Results

A 7-Gene NSCLC Prognostic and Predictive Assay

The NSCLC prognostic biomarkers identified with hybrid feature selection models [18,19] and molecular network approach [20,21] in our previous studies were validated with multiple independent microarray datasets. Based on the validation results in microarray data, 160 genes were selected for assays using low-density microfluidic qRT-PCR arrays. Among 160 genes analyzed in the qRT-PCR assays, a 7-gene signature was identified from training cohort obtained from Case Western Reserve University (CWRU; n = 83). Details of the decision tree based 7-gene prognostic and predictive model are provided in Supplementary Fig. 1. In the training cohort (CWRU), the 7-gene model stratified patients into two prognostic groups with significantly different disease-specific survival (P = 0.0043; Fig. 1A). Moreover, in the 7-gene assay predicted chemotherapy benefit (high-risk) patient group, there was a significantly prolonged disease-specific survival (P = 0.035; Fig. 1B) in adjuvant chemotherapy treated patients (ACT) compared with the observation group (OBS) who did not receive any chemotherapy. Specifically, the 30 months survival rate was <0.4 in the high-risk patients who did not receive chemotherapy (the OBS group), and the 30 months survival rate was 100% (5/5) in patients receiving adjuvant chemotherapy (the ACT group). In contrast, there was no survival benefit in receiving chemotherapy (P = 0.31; Fig. 1C) in the 7-gene assay predicted non-benefit (low-risk) group. Consistent prognostic and predictive results were confirmed in the validation set (n = 248), including NSCLC patients from another three hospitals (UM, MBRCC, and NorthShore) as well as a clinical trial JBR.10 [14] (Fig. 1D–F). In the validation set, the 7-gene signature generated significant prognostic stratification (P < 0.0001; Fig. 1D). In the predicted benefit (high-risk) patient group, there was a significant prolonged disease-specific survival in the ACT group compared with the OBS group (P = 0.0049; Fig. 1E). Specifically, the 5-year survival rate was 70.9% (39/55) in the high-risk patients who received adjuvant chemotherapy (the ACT group), whereas the 5-year survival rate was 45.8% (22/48) in high-risk patients who did not receive adjuvant chemotherapy (the OBS group). In contrast, in the predicted non-benefit (low-risk) group, there was no survival benefit in the ACT group compared with the OBS group (P = 0.46, Fig. 1F). It is noteworthy that in the predicted non-benefit (low-risk) group, patients who received adjuvant chemotherapy (ACT) had a worse post-surgical survival in the long term compared with those who did not receive any chemotherapy (OBS) in both training and validation sets (Fig. 1C and F). These results further corroborate the 7-gene model prediction of non-benefit that patients would suffer from unnecessary cytotoxicity side-effects of chemotherapy instead of benefiting from it. Overall, these results demonstrate that the 7-gene assay is both prognostic of NSCLC clinical outcome and predictive of the benefits from chemotherapy.

Fig. 1

Kaplan-Meier analyses of the 7-gene model. Patient stratification in training cohort CWRU (A), CWRU high-risk group (B), CWRU low-risk group (C), validation set (D), validation set high-risk group (E), and validation set low-risk group (F). ACT: Adjuvant chemotherapy group; OBS: observation group without chemotherapy. The validation set includes patient cohorts from MBRCC, UM, JBR.10, and Northshore. The 7-gene signature stratified patients into high-risk and low-risk groups in both training (A) and validation (D) sets. In the high-risk groups from training (B) and validation (E) sets, there were significant survival benefits in patients receiving adjuvant chemotherapy (the ACT group) compared with those who did not receive any chemotherapy (the OBS group). In the low-risk groups from training (C) and validation (F) sets, there were no significant survival benefits in patients receiving adjuvant chemotherapy (the ACT group) compared with those who did not receive any chemotherapy (the OBS group). P values were assessed with log-rank tests.

The chemoresponse prediction for specific therapeutic agents was examined in the identified 7 biomarkers. In particular, gene expression of ATP binding cassette subfamily C member 4 (ABCC4) was predictive of chemoresistance in patients receiving carboplatin, cisplatin, and Taxol, with under-expressed mRNA (higher ∆Ct) value associated with significantly decreased hazard ratio of death from disease and tumor recurrence (Table 2). In patients treated with carboplatin plus Taxol, using ∆Ct value of ABCC4 in Cox model, the hazard ratio of death from disease of was 0.43 (95% CI: [0.208, 0.888], P = 0.02) and the hazard ratio of recurrence was 0.343 (95% CI: [0.122, 0.968], P = 0.04), both statistically significant. In patients treated with Taxol, the hazard ratio of death from disease of ABCC4 ∆Ct value was 0.403 (95% CI: [0.194, 0.834], P = 0.01, Cox model) and the hazard ratio of recurrence was 0.48 (95% CI: [0.253, 0.912], P = 0.02, Cox model), both statistically significant. In patients treated with either carboplatin plus Taxol, carboplatin plus Taxotere, cisplatin plus Taxotere, or cisplatin plus Taxol, the hazard ratio of death from disease of ABCC4 ∆Ct values was borderline significant (hazard ratio: 0.528 [0.271, 1.028], P = 0.06, Cox model) and the hazard ratio of recurrence was significant at 0.545 (95% CI: [0.298, 0.998], P = 0.049, Cox model; Table 2). The expression of fucosyltransferase 7 (FUT7) was predictive of chemosensitivity to carboplatin, with under-expressed mRNA (higher ∆Ct value) associated with significantly increased hazard ratio of death from disease (hazard ratio: 1.605 [1.058, 2.435], P = 0.026, Cox model; Table 2). The expression of zinc finger protein 71 (ZNF71) was also predictive of chemosensitivity in patients treated with either carboplatin plus Taxol, carboplatin plus Taxotere, cisplatin plus Taxotere, or cisplatin plus Taxol, with a significant hazard ratio of death from disease 1.986 (95% CI: [1.001, 3.938], P = 0.049, Cox model; Table 2). Solute carrier family 39 member 8 (SLC39A8) was predictive of chemoresistance to Taxol, with a borderline significant hazard ratio of recurrence 0.584 (95% CI: [0.33, 1.03], P = 0.06, Cox model; Table 2). The expression of SLC39A8 was also predictive of chemoresistance to Alimta, with a borderline significant hazard ratio of recurrence 0.49 (95% CI: [0.219, 1.098], P = 0.08, Cox model; Table 2).

Table 2

Predictive biomarkers of chemoresponse in non-small cell lung cancer. Hazard ratios were computed with Cox proportional hazard model using ∆Ct values in qRT-PCR assays.

Genes	Chemotherapeutic agents	Hazard ratio of death from disease with 95% CI	Hazard ratio of recurrence with 95% CI	Chemosensitive/resistant
ABCC4	Carboplatin + Taxol	0.43 [0.208, 0.888]⁎	0.343 [0.122, 0.968]⁎	Chemoresistant
	Taxol	0.403 [0.194, 0.834]⁎	0.48 [0.253, 0.912]⁎	Chemoresistant
	Carboplatin + Taxol	0.528 [0.271, 1.028]#	0.545 [0.298, 0.998]⁎	Chemoresistant
	Carboplatin + Taxotere
	Cisplatin + Taxotere
	Cisplatin + Taxol
FUT7	Carboplatin	1.605 [1.058, 2.435]⁎	–	Chemosensitive
ZNF71	Carboplatin + Taxol	1.986 [1.001, 3.938]⁎		Chemosensitive
	Carboplatin + Taxotere
	Cisplatin + Taxotere
	Cisplatin + Taxol
SLC39A8	Taxol	–	0.584 [0.33, 1.03]#	Chemoresistant
	Alimta	–	0.49 [0.219, 1.098]#	Chemoresistant

Hazard ratio significant at P < 0.05.

Hazard ratio borderline significant at P < 0.08.

The 7-gene NSCLC prognostic and predictive signature is involved in cell to cell signaling and interaction, inflammatory response, and cellular movement in Ingenuity Pathway Analysis (Qiagen, Redwood City, CA). Based on the molecular network of the 7 NSCLC biomarkers (Supplementary Fig. 2A), these identified biomarkers have interactions with major inflammatory and cancer signaling hallmarks such as TNF, PI3K, NF-κB, and TGF-β. The top pathways involving the 7 signature genes and their interaction partners are nNOS signaling in skeletal muscle cells, CD27 signaling in lymphocytes, and agrin interactions at neuromuscular junction (Supplementary Fig. 2B). The 7-gene signature identified in this study does not overlap with the NSCLC gene signatures reported in previous studies [13,[15], [16], [17],[23], [24], [25]].

Protein Expression of ZNF71 is Prognostic of NSCLC Outcome

To substantiate the functional involvement of the identified 7 signature genes, protein expression of these biomarkers was evaluated with immunohistochemistry (IHC). Based on the IHC results, biomarkers with staining of good quality in FFPE NSCLC tumor tissues were further quantified with AQUA. Protein expression of ZNF71 was identified as prognostic of NSCLC outcome in two TMA cohorts (Fig. 2A). Details of the AQUA assays and results are provided in Supplementary File 1. Based on the quantitative AQUA scores representing ZNF71 protein expression levels in tumor tissues, a cutoff point was defined for patient prognostic stratification in training cohort YTMA250 (n = 145). Specifically, when loge-transformed ZNF71 AQUA scores were greater than or equal to 7.9, patients had significantly better disease-specific survival (P = 0.021) than those with a lower ZNF71 protein expression level (Fig. 2B). This cutoff was further validated with significant patient stratification (P = 0.047, Fig. 2C) in an independent cohort YTMA79 (n = 46). Higher protein expression of ZNF71 is significantly associated with better patient survival, which is concordant with its mRNA results in multiple independent patient cohorts and its observed association with chemosensitivity in Taxol (Taxotere) plus platinum-based treatment in NSCLC patients (Table 2). These results indicate that ZNF71 is a prognostic protein biomarker and might be a potential therapeutic target of NSCLC. Furthermore, ZNF71 had a 7% (19/271) of loss of DNA copy number in tumors vs. normal lung tissues in a NSCLC patient cohort from Starczynowski et al. [26] (n = 271; Supplementary Table 1). These results suggest the concordance in the loss of DNA copy number, down-regulated mRNA and protein expression of ZNF71 in lung cancer progression.

Fig. 2

Kaplan-Meier analyses of ZNF71 protein expression quantified by AQUA. ZNF71 immunofluorescence images of different expression levels in TMA. Patients were stratified into two groups based on ZNF71 AQUA scores. Patients with loge(ZNF71 AQUA Score) ≥ 7.9 had a low-risk and those with loge(ZNF71 AQUA Score) < 7.9 had a high-risk for tumor metastasis in training cohort YTMA250 (B) and validation cohort YTMA79 (C). P values were assessed with Wilcoxon tests.

Concordant mRNA and Protein Under-expression of CD27 in NSCLC Progression

The protein expression level of CD27 was quantified with ELISA assays in NSCLC tumor tissues (n = 29) and normal adjacent lung tissues (n = 9). Spearman correlation coefficient between mRNA and protein expression of CD27 is 0.494 (P < 0.0088; Fig. 3A) in tumor tissues. CD27 had an average protein expression of 599.06 pg/mL in low-risk patients with a better disease-specific survival, and an average protein expression of 245.5 pg/mL in high-risk patients with a poorer disease-specific survival in ELISA assays. CD27 had significant under-expression in high-risk patients vs. low-risk patients at mRNA level with a fold-change of 0.17 (P < 0.00001) and a fold-change of 0.41 (P < 0.02) at protein level (Fig. 3B). CD27 had an average protein expression of 191 pg/mL in normal lung tissues. CD27 had significant protein over-expression in NSCLC tumor vs. normal tissues with a fold-change of 2.56 (P < 0.025), while mRNA expression in tumor vs. normal tissues was not significantly different (Fig. 3B). The over-expressed CD27 protein in NSCLC tumors is concordant with an observed 4% (11/271) of gain or amplification of DNA copy number in tumors vs. normal lung tissues in the NSCLC patient cohort from Starczynowski et al. [26] (n = 271; Supplementary Table 1). Overall, these results demonstrate that CD27 had concordant under-expression at both mRNA and protein levels in NSCLC patients with a poor outcome and a greater chance of tumor recurrence and metastasis. The overexpressed CD27 protein level in NSCLC tumor vs. normal lung tissues indicates that CD27 regulation in tumorigenesis and metastatic processes is different. Our results confirm the potential role of CD27 as a target in lung cancer immunotherapy [27,28].

Fig. 3

Comparison of mRNA and protein expression of CD27 in NSCLC patient samples. (A) Scatterplot with regression line for CD27 mRNA (relative quantity) in qRT-PCR and protein expression (pg/mL) in ELISA assays of 29 NSCLC tumor resections. RQ: relative quantity, measured as 2− values in qRT − PCR with UBC as the control gene. R: Spearman correlation coefficient. (B). Comparison of CD27 fold-change in NSCLC vs. normal lung tissues and high-risk vs. low-risk NSCLC tumors in qRT-PCR and ELISA assays. High-risk NSCLC patients had a poor survival outcome and low-risk NSCLC patients had a good survival outcome. Bar plot shows mean + SE. *: P < 0.05.

Discussion

Lung cancer is the second most common cancer in both men and women, and remains the highest cancer-related mortality with a death rate higher than colon, prostate, and breast cancer combined. Currently, there is no clinically available multi-gene assay to prognosticate and predict the benefits of chemotherapy in NSCLC patients for improved personalized treatment. Immunotherapy is more effective and less toxic than chemotherapy in advanced lung cancers [[5], [6], [7], [8],29,30], and recent studies show promise of immunotherapy in early stage lung cancer patients [8]. Nevertheless, predictive biomarkers and therapeutic targets of immunotherapy are not well established.

There were abundant publically available microarray data generated in NSCLC patient tissues. Although microarray platforms are phasing out, the legacy data and biomarkers identified in microarray platforms are still useful in the RNA-seq era [9]. However, high-throughput platforms such as microarrays and RNA-seq are not suitable for routine clinical tests. Validation of biomarkers identified from high-throughput technologies with qRT-PCR emerges as the most promising experimental protocol for developing multi-gene assays for clinical applications.

NSCLC prognostic biomarkers were identified with hybrid feature selection models [18,19,31] and molecular network approach [20,21] in our previous studies. The hybrid feature selection models [18,19,31] contain multiple layers of gene selection algorithms in the process of biomarker identification. This scheme takes advantage of different algorithms in different stages of gene shaving, in order to identify the gene signatures with the optimal performance. The molecular network approach [20,21] constructs genome-scale co-expression networks in good-prognosis and poor-prognosis patient groups separately, and compares the network structures of these two patient groups to identify disease-specific network modules. Next, genes with concurrent co-expression with multiple major lung cancer signaling hallmarks were pinpointed from disease-specific network modules for further gene signature identification. This approach embedded biological relevance into biomarker identification. The signature genes identified with these sophisticated approaches were validated with multiple independent publically available microarray datasets. Genes with consistent expression patterns in multiple validation sets were included in qRT-PCR assays. The 7-gene signature identified in qRT-PCR assays were prognostic and predictive of chemoresponse in patient cohorts from multiple hospitals and JBR.10.

The identified 7 signature genes have interactions with major inflammatory and cancer signaling hallmarks including TNF, PI3K, NF-κB, and TGF-β (Supplementary Fig. 2A). Multiple signature genes are potential targets in cancer immunotherapy. Specifically, reduction of DAG1 may increase susceptibility of muscle fibers to necrosis [32]. A study shows that DAG-1 cells are resistant to TNF-α and IFNγ-induced apoptosis, with implications in bladder cancer progression and resistance to immunotherapy [33]. CD27 is part of TNF receptor family, and overexpression of CD27 induces NF-κB activation involving signaling transduction of TNF receptor-associated factors [34]. CD27 was also reported as a potential target of cancer immunotherapy [27,28]. The synergy between PD-1 blockade and CD27 stimulation for CD8+ T-cell driven anti-tumor immunity was reported recently [35], indicating the therapeutic potential of CD27 in neoadjuvant PD-1 blockade in resectable lung cancer [36]. The zinc finger protein ZNF71 is induced by TNF-α [37] and ZNF71 SNP was found to be associated with asthma in human serum [38]. CCL19 is regulated by multiple NF-κB and INF family transcription factors in human monocyte-derived dendritic cells [39]. ABCC4 is associated with multiple drug resistance in cancer [40] and smooth muscle cell proliferation [41], and interacts with PI3K in cancer prognosis and drug resistance [42]. Our results on ABCC4 in Table 2 are consistent with its functional role and reported drug resistance. FUT7 interacts with TNF-α in human bronchial mucosa [43] and its induction at sites of tumor cell arrest is involved in metastasis [44]. NF-κB was reported to regulate expression of the zinc transporter SLC39A8 [45]. Indirect interactions between TGF-β and SLC39A8 are involved in tumorigenesis [46] and fibrogenic response [47].

The 7-gene signature identified in this study does not overlap with the NSCLC gene signatures reported in recent studies [15,16,[23], [24], [25]]. However, several biomarker genes identified in this study belong to the same families or functional categories as the biomarkers identified in [[14], [15], [16]]. In particular, FUT7 from the current study and FUT3 from Kratz et al. [16] are both fucosyltransferase and involved in metabolism. In the 12-gene prognostic and predictive signature from Tang et al. [15], two genes belong to the same family or share similar functions as the 7-gene signature. Specifically, SLC35A5 from Tang et al. [15] and SLC39A8 from this study both belong to solute carrier superfamily, and ATPase Phospholipid Transporting 8A1 (ATP8A1) from Tang et al. [15] and ATP Binding Cassette Subfamily C Member 4 (ABCC4) from this study are both involved in energy metabolism. The 15-gene prognostic and predictive gene signature of JBR.10 [14] also contains two genes that share similar functions as the 7-gene signature. ATPase Na+/K+ Transporting Subunit Beta 1 (ATP1B1) from Zhu et al. [14] and ABCC4 from this study are again involved in energy metabolism, and ZNF236 from Zhu et al. [14] and ZNF71 identified in this study both belong to zinc finger protein family. Overall, the 7-gene signature presented in this study and two previous gene signatures from Zhu et al. [14] and Tang et al. [15] are all prognostic of NSCLC outcome and predictive of the benefits of chemotherapy. These three gene signatures all contain a biomarker related to ATP activities and energy metabolism. Other shared gene families between the 7-gene signature and these two signatures include zinc finger protein and solute carrier superfamily. The 7-gene signature and the practical prognostic gene assay for non-squamous NSCLC by Kratz et al. [16] both contain biomarkers from fucosyltransferase family. These common gene families shared by the NSCLC gene signatures with promise for clinical utility might be functionally involved in tumor metastasis with implications in lung cancer therapy.

The protein expression of the identified 7 signature genes was also validated in this study. In particular, ZNF71 protein expression quantified with AQUA was a prognostic biomarker in two NSCLC patient cohorts (n = 191). Higher mRNA and protein expression of ZNF71 is both associated with good prognosis, and ZNF71 mRNA is predictive of chemosensitivity in Taxol (Taxotere) plus platinum-based treatment in NSCLC patients. These results demonstrate that ZNF71 mRNA and protein expression can both be used in prognostication of NSCLC in clinical applications and ZNF71 may be a therapeutic target. CD27 had highly correlated mRNA and protein expression, with significant under-expression in poor prognostic (high-risk) NSCLC patients. CD27 mRNA and protein expression could potentially be used as a biomarker and target in lung cancer immunotherapy. Protein expression of CCL19 was also confirmed with ELISA in NSCLC tumor and adjacent normal tissues. CCL19 protein was under-expressed in NSCLC tumor tissues compared with normal lung tissues, with no statistically significant difference (results not shown). CCL19 also had lower protein expression in poor-prognosis (high-risk) NSCLC patients compared with good-prognosis (low-risk) patients, with no statistically significant difference (results not shown). The trend of CCL19 protein expression was qualitatively concordant with its mRNA expression that higher expression of CCL19 is associated with good prognostic outcome of NSCLC. CCL19 had a 12.5% (34/271) of a loss of DNA copy number in tumors vs. normal lung tissues in the NSCLC patient cohort from Starczynowski et al. [26] (n = 271; Supplementary Table 1), which suggests a loss of DNA copy number and down-regulated mRNA and protein expression of CCL19 in NSCLC progression. In our previous integrated DNA copy number and gene expression regulatory network analysis of NSCLC metastasis, CCL19 is a driver gene and CD27 expression is modulated by CCL19 in squamous cell lung cancer patients with good prognosis [48]. Together with the molecular network reported in the literature (Supplementary Fig. 2A), the interaction between CCL19 and CD27 could be through PI3K and NF-κB complexes. In addition, FUT7 and DAG1 had concordant loss or deletion of DNA copy number (Supplementary Table 1) and down-regulated gene expression in NSCLC progression (Table 2 and Supplementary Fig. 1).

Conclusions

This study presents a 7-gene predictive assay based on qRT-PCR to improve NSCLC treatment in clinics. This 7-gene assay provides accurate prognostication and prediction of the clinical benefits of chemotherapy in multiple patient cohorts from the US hospitals and the clinical trial JBR.10. The 7-gene assay is enriched in inflammatory response. The protein expression of ZNF71 is prognostic of NSCLC outcome in two independent patient cohorts, which is concordant with its mRNA expression. These results demonstrate that ZNF71 is a prognostic protein biomarker and a potential therapeutic target of NSCLC. The protein expression of CD27 was strongly correlated with its mRNA expression in NSCLC tumor tissues, and could potentially serve as a biomarker and target of immunotherapy in lung cancer. Multiple signature genes had concordant DNA copy number variation, mRNA and protein expression in NSCLC progression. The results presented in this study are important for precision therapy in NSCLC patients with implications in developing new therapeutic strategies to combat this deadly disease.

Funding

This work is supported by NIH Grants R01/R56LM009500 (Guo), P20RR16440 ARRA Supplement (Guo), and R01ES021764 (Guo). The WV-INBRE (IPA license) is supported by NIH Grant P20GM103434. The content is solely the responsibility of the authors and does not necessarily represent the official views of the NIH.