IL-13 receptors as possible therapeutic targets in diffuse intrinsic pontine glioma

Diffuse intrinsic pontine glioma (DIPG) is a universally fatal childhood cancer of the brain. Despite the introduction of conventional chemotherapy and radiotherapy, improvements in survival have been marginal and long-term survivorship is uncommon. Thus, new targets for therapeutics are critically needed. Early phase clinical trials exploring molecularly-targeted therapies against the epidermal growth factor receptor (EGFR) and novel immunotherapies targeting interleukin receptor-13α2 (IL-13Rα2) have demonstrated activity in this disease. To identify additional therapeutic markers for cell surface receptors, we performed exome sequencing (16 new samples, 22 previously published samples, total 38 with 26 matched normal DNA samples), RNA deep sequencing (17 new samples, 11 previously published samples, total 28 with 18 matched normal RNA samples), and immunohistochemistry (17 DIPG tissue samples) to examine the expression of the interleukin-4 (IL-4) signaling axis components (IL-4, interleukin 13 (IL-13), and their respective receptors IL-4Rα, IL-13Rα1, and IL-13Rα2). In addition, we correlated cytokine and receptor expression with expression of the oncogenes EGFR and c-MET. In DIPG tissues, transcript-level analysis found significant expression of IL-4, IL-13, and IL-13Rα1/2, with strong differential expression of IL-13Rα1/2 in tumor versus normal brain. At the protein level, immunohistochemical studies revealed high content of IL-4 and IL-13Rα1/2 but notably low expression of IL-13. Additionally, a strong positive correlation was observed between c-Met and IL-4Rα. The genomic and transcriptional landscape across all samples was also summarized. These data create a foundation for the design of potential new immunotherapies targeting IL-13 cell surface receptors in DIPG.


Introduction
Brainstem gliomas account for up to 20% of all central nervous system tumors in children less than 15 years of age, with a median age at presentation of 6-7 years [1]. Diffuse intrinsic brainstem gliomas (DIPG) comprise 80% of all brainstem gliomas and are typically anaplastic a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 astrocytoma (AA), glioblastoma multiforme (GBM) or Grade 2 lesions [2]. No trials have yet shown benefit of chemotherapy for management of patients with diffuse intrinsic pontine glioma [3]. In the United States, approximately 350 new patients per year are diagnosed with DIPG. Prognosis for patients with DIPG is poor, with a median survival of less than 1 year; fewer than 20% of children are alive at 2 years [4]. Standard therapy consists of conventional local field radiotherapy to a dose of 54-60 Gy for 6 weeks. Without radiotherapy, median survival is approximately 20 weeks [1,5]. Radiotherapy leads to transient improvements in neurological function and improves overall survival by approximately 2-3 months [6]. Multiple studies using a variety of conventional, high-dose and targeted therapies have demonstrated no survival advantage in patients with DIPG [7][8][9][10][11][12][13]. Thus, novel therapies are clearly needed for patients with DIPG.
Early phase clinical trials evaluating the safety and efficacy of the EGFR inhibitors gefitinib and erlotinib have demonstrated modest activity in brainstem gliomas similar to DIPG [14]. Although many previous small molecule inhibitors have demonstrated no significant benefit, an expanding understanding of the mutational landscape and defined therapeutic targets have led us to discoveries of a set of potentially effective small molecule inhibitors. Recently published work from our group prioritized panobinostat [15], an epigenome-modifying small molecule histone deacetylase, via in vitro and in vivo validation on DIPG primary cell cultures.
An alternative approach focuses on antibody drug conjugates (ADCs); the emergence of bivalent antibodies makes targeting both receptors and receptor tyrosine kinases in DIPG a viable therapeutic avenue. Based on our newly-generated data, IL-4, IL-13 and their receptors IL-4Rα, IL-13Rα1, IL-13Rα2 may be reasonable therapeutic targets in DIPG due to differential overexpression on the cell surface of cancerous cells. This is consistent with previous research showing the potential of IL-4 and IL-13Rα2 as potential therapeutic targets in pediatric brain tumors such as glioblastoma [16][17][18][19][20][21][22].
IL-4 and IL-13 tethered exotoxins have already led to clinical trials, e.g., the IL-13PE38QQR trial NCT00880061 at NIH [16,23] and the IL-4::Pseudomonas Exotoxin product PRX321 in development by Protox Therapeutics as NCT00797940 [24]. Similarly, IL-13Rα2 as an immunogen is the basis of a DIPG/high grade glioma vaccine trial [NCT01130077]. This pilot clinical trial evaluated the safety and efficacy of subcutaneous vaccination with glioma-associated antigens (GAAs) and their epitope peptides including IL-13Rα2, EphA2, and Survivin. Five of twenty-six patients experienced inflammatory-associated pseudoprogression, defined as sixmonth progression-free survival following transient clinical progression, with one patient experiencing an event free survival of greater than three years [25,26]. Thus, targeting overexpressed glioma-specific receptors or receptor tyrosine kinases in combination with targeted antibodies may further improve survival. A case report of IL13RA2-targeted chimeric antigen receptor T-cell therapy showed regression of glioblastoma in a human patient [27]. Dual receptor targeting of EGFR and c-Met, possible using certain bivalent antibody approaches [28][29][30], may be reasonable combination therapeutic targets in high grade glial tumors based on expression studies [31][32][33][34][35] from which one DIPG clinical trial employing the EGFR inhibitor, erlotinib, has opened [NCT01182350].
Based on the newly-generated sequencing data, IL-4 and IL-13 axis components were identified as being targets of interest in DIPG. Thus, in this work we explore IL-4 or IL-13 axis components as well as EGFR or c-MET as combination therapy targets in DIPG. We have surveyed archived human DIPG cases (28 by RNA-seq, 17 by immunohistochemistry (IHC)) for IL-4, IL-13, IL-4Rα, IL-13Rα1, IL-13Rα2, EGFR and c-Met expression. The studies herein demonstrate that most DIPG tissues displayed a high content of IL-4 and weak or absent expression of IL-13, yet modest but significant expression of IL-4 and IL-13 was observed at the transcript level. By contrast, IL-13Rα1 and IL-13Rα2 were found to be significantly enriched at the transcript level in DIPG tissues and expression by IHC was found to be high with 41% of these tissues also staining positive for the oncogenes EGFR and c-Met.
RNA-seq identified 3,577 (11.6% of 30,796 genes) differentially upregulated and 3,365 (10.9% of 30,796 genes) differentially downregulated genes in tumor samples relative to normal samples. ) were noted to have gain events in several samples. Frequent 1q gain events have been associated with poor prognosis in a variety of pediatric diseases including neuroblastoma and medulloblastoma [36]. This analysis was generated solely from the Illumina sequencing data, and may not completely capture the mutational landscape of these DIPG samples due to coverage limitations of high throughput sequencing technology. A complete table of genomic variants is provided in S1 Table, and a complete table of gene expression data is provided in S2 Table.

Key mutations identified by Sanger sequencing
Sanger sequencing was performed on VUMC-DIPG-7 and VUMC-DIPG-8 samples independently of the main body of research. Sanger sequencing identified H3F3A (K27M) mutations ( Fig 1A), consistent with other DIPG samples [37]. In previously reported samples, 6/16 (37.5%) samples carry H3F3A mutations; in new samples for this report, 5/22 (22.7%) carry H3F3A mutations. H3F3A mutations were not identified in the Illumina sequencing due to lack of coverage at the H3F3A loci; nonetheless, coverage for other regions was sufficient to identify mutations in other genes.
For reference, we have provided several normal tissue images stained with the different IHC stains, as well as several scored experimental IHC samples (S3 and S4 Figs, respectively).

Discussion
Previous studies have reported that IL-4 blockade exhibits antitumor activity in rodent models of glioma. Follow-up reports elucidated the role of IL-4 and IL-13Rα2 as mediators of aberrant Stat3 signaling driving increased expression of anti-apoptotic genes [38,39]. IL-13Rα1 and IL-13Rα2 are known tumor-associated antigens and recent clinical trials have explored inoculation with IL-13Rα2 and other tumor associated antigens as an immunotherapy [40]. However, the roles of IL-4 and IL-4Rα have not been investigated in DIPG. To this end, we sought to characterize expression of the extended targets of the IL-4/IL-13 signaling axis in DIPG tissues. Of note, exome sequencing did not identify any somatic point mutations or copy number variations in IL-4/IL-13 signaling axis genes. However, consistent with previous reports, RNA-seq gene expression studies revealed significant upregulation of IL-13Rα1 and IL-13Rα2 but no statistically significant differential of IL-4Rα. Immunohistochemistry revealed a high content of IL-4, but not IL-13. However, the receptors for these cytokines (IL-4Rα, IL-13Rα1 and IL-13Rα2) were all expressed in 85% of the cases analyzed. In 41%, expression of these receptors were all labeled as strong and/or the number cells expressing these receptors was more than 25% of the tumor cell population. These results show that IL-4 and IL-13 receptor upregulation is prevalent in DIPG tissues which may be driven by global dysregulation of gene expression mediated by epigenetic mechanisms. Therefore, the potential exists for therapeutic drug targeting of these receptors. The oncogenes EGFR and c-MET were co-expressed in approximately 41% of these gliomas. Additionally, significant correlations with IL-13Rα2 (EGFR) and IL-4Rα and IL-13Rα1 (c-MET), respectively. These results raise the possibility that EGFR may be an epitope worth exploring in combination with IL-13Rα2 immunotherapies. Previous reports have identified the expression of epidermal growth factor variant III (EGFRvIII) in a subset of DIPG tissues with a vaccine pilot study ongoing (NCT01058850) [41,42]. While we did not examine these tissues for expression of EGFRvIII, we found no significant difference in EGFR gene expression in tumor tissues versus matched normal and 41% of cases surveyed expressed EGFR at the protein level. However, recent sequencing efforts in larger sample sizes have shown that recurrent EGFR mutations are infrequent [43]. While the EGFR inhibitors erlotinib and gefitinib have demonstrated activity in this disease, this activity may be due in part to off-target kinase inhibition [44,45]. In summary, we have identified IL-4Rα, IL-13Rα1 and especially IL-13Rα2 as potential therapeutic targets in DIPG. The co-expression of these receptors with EGFR and c-MET may further expand immunotherapy options for this disease.

Patient selection for sequencing experiments
We obtained tumor tissues from sixteen (16)  We also obtained tumor tissues from seventeen (17) new DIPG samples for transcriptome sequencing. Nine (8) samples were provided by Dr. Xiao-Nan Li, four (4) were obtained from OHSU patients, one (1) was obtained from a cell line provided by Dr. Michelle Monje, and

Patient selection for immunohistochemistry experiments
We obtained tumor tissues from six (6) autopsy cases at OHSU, and from two (2) patients that underwent a surgical biopsy of their DIPG tumor at OHSU. Eight (8) samples with matched normal brain tissue samples were provided by collaborator Dr. Xiao-Nan Li of Texas Children's Cancer and Hematology Center. Nine (9) additional cases were provided by Dr. Cynthia Hawkins of Toronto Sick Children's Hospital, in the form of a paraffin block tissue microarray, with an accompanying legend for diagnostic classification of the tissues on the microarray grid. All internal and external cases were morphologically re-evaluated to ensure sample adequacy and diagnostic accuracy. All studies with human tissue were performed with approval of the Oregon Health and Science University (OHSU) institutional review board. Additional information is available in S4 Table. RNA-seq RNA-seq libraries were constructed using the standard TruSeq protocol (Illumina). Briefly, poly(A)+ RNA was isolated from total RNA using oligo-dT attached to magnetic beads. The poly(A)+ was chemically fragmented. Double stranded cDNA was generated using random hexamers as primers for both strands. The cDNAs were blunted-end, then a single A was attached to the 3' end of each strand to promote ligation. Standard Illumina adapters with indices were ligated to the cDNAs. The resulting ligations were amplified using a limited number of rounds of polymerase chain reaction (PCR). Libraries were separated from unincorporated material using AmpPure beads (Beckman Coulter). Libraries were evaluated using the Bioanalyzer (Agilent) and the concentration of each library was determined by real time PCR (Kapa Biosystems) on a StepOnePlus real time PCR machine (Life Technologies). Sequencing was performed on either a HiSeq 2000 (Illumina) or HiSeq 2500 (Illumina). Sequence assembly was performed using either CASAVA (Illumina) or Bcl2Fastq (Illumina).

Somatic point mutation identification by exome capture sequencing
The methodology described below is a minor modification of the approach described in [46]. All captured DNA libraries were sequenced with the Illumina HiSeq in paired-end mode, yielding 80 base pairs from the final library fragments. The reads that passed the chastity filter of Illumina BaseCall software were used for subsequent analysis. Matepairs were pooled and mapped as single reads to the reference human genome (NCBI build 36.1, hg18), excluding unordered sequence and alternate haplotypes, using Bowtie [47], keeping unique best hits, and allowing up to two mismatched bases per read. Likely PCR duplicates, defined as reads with equal match intervals on the reference genome, were removed; and individual basecalls with Phred quality less than Q20 were excluded. A mismatched base (SNV) was identified as a variant when 1) it had at least six reads of support, 2) it constituted at least 10% of the coverage at that position, 3) it was observed on both strands, and 4) it fell within 50 bases of a region targeted for capture. In matched samples (samples with matched normals), a variant was called somatic if 1) there was 8× coverage in the matched normal, 2) it did not occur in the matched normal sample in more than two reads and 2% of the coverage (or 4% of the coverage if the tumor variant fraction was at least 20%), 3) if it had any support in the matched normal, and 4) it was not present in either dbSNP (v137) or the Exome Sequencing Project data set (ESP6500). In unmatched samples, a variant was called probable somatic if it affected the same codon as either a somatic point mutation observed in the matched samples or a somatic mutation from COSMIC (v66). Finally, all somatic/probable somatic variants were screened against the full set of benign samples. Variants were removed from further consideration if they appeared in any benign sample with at least three reads and 10% of the coverage.

Sanger sequencing
Genomic DNA was isolated from tissue, using the QIAamp DNA mini isolation kit (Qiagen, Valencia, CA, USA), and the sequence of interest was amplified by PCR. Primers were described by Wu et al. [37]. PCR products were subsequently sequenced by the dideoxy chaintermination method, using the ABI PrismTM BigDye Terminator kit (Perkin Elmer, Foster City, CA, USA), run on the ABI Prism Genetic Analyser 3100 automatic DNA autosequencer (Perkin Elmer), and analysed with ABI sequence Alignment Editor software.

Annotation
We annotated the resulting somatic mutations using CCDS transcripts wherever possible. If no CCDS transcript was available, we use the coding regions of RefSeq transcripts. HUGO gene names were used.

Exome copy number analysis
Copy number aberrations were quantified and reported for each gene as the segmented normalized log2-transformed exon coverage ratios between each tumor sample and its matched normal as described in [48]. For unmatched samples, we used the average coverage derived from the full set of benign samples as a 'normal pool;' for chromosomes X and Y only, we formed separate normal pools from the XX and XY benign samples. We identified segments as focal if they contained 10 or fewer CCDS gene annotations. To identify a gene as gained or lost we first considered the segmented data, requiring a copy number gain or loss of 30% (ratio > = 1.3 or < = 0.7) to call the gene as gained or lost. In addition, we considered the distribution of copy number estimates for each gene's individual exons. If the mean exonic copy number showed a gain or loss of at least 30% and deviated from the null hypothesis by at least 2.5 s.d., we called the gene gained or lost.

Gene expression in RPKM inferred from RNA-seq data
All transcriptome libraries were sequenced with the Illumina HiSeq in paired end mode. The reads that passed the chastity filter of Illumina BaseCall software were used for subsequent analysis. We trimmed all reads to 85-mers and aligned them to the reference human genome (NCBI build 36.1, hg18), plus a splice junction set including 84 bases on either side of the 2008 Illumina splice junction set, using Bowtie [47] in single read mode, keeping unique best hits and allowing up to two mismatched bases. Matepairs from paired-end runs were pooled and treated as single reads. Next, all of the exons for a single gene were concatenated to form a single 'transcript' for that gene. Reads that mapped to the exons in the gene, as well as reads that mapped to the splice junctions, were remapped to the transcript. We then walked the transcript and summed the coverage at each position, then divided the result by the transcript length × the number of reads in the sample, and then multiplied the result by one million. This method is a modified version of the one described in [49].

DNA fingerprinting using common high-frequency variants
We considered 147 genomic positions that tended to have good coverage in both wholeexome and RNA-seq data and for which at least two alleles were widespread. For each sample, we constructed a 'fingerprint' for those positions with at least six fragments of coverage; an allele was included in the fingerprint if it was seen in at least two fragments of coverage and at least 15% of the total fragments. Two samples were deemed a match if at least 85% of the positions with sufficient coverage in both samples had identical fingerprints. In practice, samples from the same patient matched at more than 90% of positions, whereas samples from different patients matched at fewer than 75% of positions. For this study, all the included samples from the same patient met the 85% match criteria, including all RNA-seq and exome samples. In addition, no samples from different patients exceeded 74% similarity.

Statistical analyses
Bioinformatics and computational methods are described above. For comparison of gene expression levels of the 429 genes associated with the targets of the drug screen, we generated RNA-seq data and then compared tumor tissues to normal tissues to generate significance scores using a t-test with the Benjamini-Hochberg multiple comparison correction.
GraphPad Prism (GraphPad Software Inc., San Diego, CA, USA) was used for all other statistical analyses. Correlation studies were performed using the Pearson product-moment correlation coefficients for select IL-4 signaling axis genes and the oncogenes c-Met and EGFR. Comparisons between groups were determined by Student's t-test. Differences were considered statistically significant when p < 0.05.

Immunohistochemistry
We performed immunohistochemistry using the following primary antibodies, using the man- . Each antibody was first tested with available known positive control normal tissues and tumors, and titrated to obtain an optimal dilution for the best signal/noise ratio. Antigen retrieval procedures were also tested for each antibody, starting with the manufacturer's recommended reagents. The immunohistochemical staining results were interpreted by two pathologists (SHG and VZ). We used a semi-quantitative scoring system for analysis: 0: No labeling, 1: Less than 10% of tumor cells labeled, 2: 10-50% of tumor cells labeled, 3: More than 50% of tumor cells labeled. The labeling of tumor cells, when present, was also visually quantified as weak, medium, or strong. In interleukin receptor studies, cytoplasmic or membrane labeling was taken into account as positive labeling, but nuclear labeling was disregarded, after our antigen adsorption assays revealed that the nuclear staining was not specific. Each assay was conducted with several identically processed positive control tissues and negative controls were included in all assays. Each antibody was performed in all cases. In four cases, due to depletion of the tissue block, we were able to test only some of the seven antibodies. Where a specific antibody testing was not performed due to limited tissue availability, the designation not present (NP) was used to indicate that the assay was not performed.