Originally published as JCO Early Release 10.1200/JCO.2005.02.9405 on December 19 2005

This Article Abstract Full Text (PDF) Supplementary Tables Purchase Article View Shopping Cart Alert me when this article is cited Alert me if a correction is posted Services Email this article to a colleague Similar articles in this journal Similar articles in PubMed Alert me to new issues of the journal Save to my personal folders Download to citation manager Rights & Permissions Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Bredel, M. Articles by Sikic, B. I. Search for Related Content PubMed PubMed Citation Articles by Bredel, M. Articles by Sikic, B. I.

Tumor Necrosis Factor-–Induced Protein 3 As a Putative Regulator of Nuclear Factor-B–Mediated Resistance to O⁶-Alkylating Agents in Human Glioblastomas

From the Division of Oncology, Center for Clinical Sciences Research, Institute for Computational and Mathematical Engineering, Departments of Neurosurgery, Pathology, and Neurology, Stanford University School of Medicine, Stanford, CA; Ina Levine Brain Tumor Center, Neuro-Oncology Research and Neurosurgery Research, Barrow Neurological Institute of St. Joseph's Hospital and Medical Center, Phoenix, AZ; and the Department of General Neurosurgery, Neurocenter, University of Freiburg, Freiburg, Germany

Address reprint requests to Markus Bredel, MD, PhD, Division of Oncology, Stanford University School of Medicine, 269 Campus Dr, CCSR-1105, Stanford, CA 94305-5151; e-mail: mbredel{at}stanford.edu

	ABSTRACT

TOP ABSTRACT INTRODUCTION PATIENTS AND METHODS RESULTS DISCUSSION Authors' Disclosures of... Author Contributions GLOSSARY REFERENCES

 
PURPOSE: Pre-existing and acquired drug resistance are major obstacles to the successful treatment of glioblastomas.

METHODS: We used an integrated resistance model and genomics tools to globally explore molecular factors and cellular pathways mediating resistance to O⁶-alkylating agents in glioblastoma cells.

RESULTS: We identified a transcriptomic signature that predicts a common in vitro and in vivo resistance phenotype to these agents, a proportion of which is imprinted recurrently by gene dosage changes in the resistant glioblastoma genome. This signature was highly enriched for genes with functions in cell death, compromise, and survival. Modularity was a predominant organizational principle of the signature, with functions being carried out by groups of interacting molecules in overlapping networks. A highly significant network was built around nuclear factor-B (NF-B), which included the persistent alterations of various NF-B pathway elements. Tumor necrosis factor-–induced protein 3 (TNFAIP3) was identified as a new regulatory component of a putative cytoplasmic signaling cascade that mediates NF-B activation in response to DNA damage caused by O⁶-alkylating agents. Expression of the corresponding zinc finger protein A20 closely mirrored the expression of the TNFAIP3 transcript, and was inversely related to NF-B activation status in the resistant cells. A prediction model based on the resistance signature enabled the subclassification of an independent, validation cohort of 31 glioblastomas into two outcome groups (P = .037) and revealed TNFAIP3 as part of an optimized four-gene predictor associated significantly with patient survival (P = .022).

CONCLUSION: Our results offer strong evidence for TNFAIP3 as a key regulator of the cytoplasmic signaling to activate NF-B en route to O⁶-alkylating agent resistance in glioblastoma cells. This pathway may be an attractive target for therapeutic modulation of glioblastomas.

	INTRODUCTION

TOP ABSTRACT INTRODUCTION PATIENTS AND METHODS RESULTS DISCUSSION Authors' Disclosures of... Author Contributions GLOSSARY REFERENCES

 
The prognosis of patients with glioblastoma multiforme has not improved substantially during the last decades, and almost all patients die as a result of their disease. Current treatment approaches are based on radiation therapy and alkylating agent chemotherapy. O⁶-guanine alkylating agents, such as 1,3-bis(2-chloroethyl)-1-nitrosourea (BCNU) and temozolomide (TMZ), are among the most widely used chemotherapeutics in treating glioblastomas because they can efficiently cross the blood-brain barrier. These agents have modest efficacy against glioblastomas.^1,2 A subset of glioblastoma patients demonstrates an initial response, lasting a few or several months and eventually leading to tumor recurrence.

One of the most prominent resistance mechanisms to alkylating agents includes O⁶-methylguanine DNA methyltransferase (MGMT),³ which acts as a suicide enzyme by removing the methyl or chloroethyl damage at the O⁶-position of guanine. Epigenetic MGMT gene silencing via promoter hypermethylation, present in approximately 40% of cases, has been shown to predict outcome in glioblastoma patients treated with BCNU or TMZ.^4,5

The growing awareness that resistance in human cancer is likely regulated by the coordinated alteration of molecular pathways⁶ suggests that many more genes might be involved in the development of resistance phenotypes in glioblastomas than the changes described thus far for a limited number of known resistance genes.³ Resistance of glioblastomas to alkylating agents such as BCNU and TMZ seems to follow a more complex pattern than simple dependence on MGMT levels.^3,7-9

Excessive and prolonged activation of nuclear factor-B (NF-B) has been established as a principal mechanism of tumor chemoresistance, which is primarily mediated by its antiapoptotic activity.^10,11 Some evidence also indicates a link between the NF-B pathway and resistance of glioblastoma cells to O⁶-alkylating agents, and suggests that inhibition of NF-B is a promising means to potentiate the cytotoxic effects of these agents.¹² The NF-B complex consists of a family of heterodimers, of which the p50/p65 heterodimer is the most abundant form. NF-B is active in the nucleus and is inhibited through its sequestration in the cytoplasm by the inhibitors of B (IBs), primarily through the interaction of IB proteins with p65. IB is a target of several well-characterized kinase cascades that activate IB kinases (IKKs), which phosphorylate IB and mark it for degradation via the ubiquination pathway, thereby allowing activation of NF-B. Activated NF-B translocates to the nucleus and binds DNA at B-binding motifs, which initiates gene transcription. Anticancer drugs are known to induce the expression of NF-B target genes through the direct activation of NF-B and the secondary production of NF-B activators.¹¹

There is increasing recognition of the value of comprehensive approaches to the molecular characterization of biologic phenotypes such as drug resistance. We have here utilized an integrated model of glioblastoma resistance to O⁶-alkylating agents and genomics tools to globally explore molecular factors, cellular pathways, and functional interaction networks perturbed during the selection and evolution of drug resistance in glioblastoma cells. Our results highlight a key role of a cytoplasmic signaling cascade that activates NF-B in response to these agents. This pathway reveals the tumor necrosis factor-–induced protein 3 (TNFAIP3) gene as a factor associated with glioblastoma cell resistance.

	PATIENTS AND METHODS

TOP ABSTRACT INTRODUCTION PATIENTS AND METHODS RESULTS DISCUSSION Authors' Disclosures of... Author Contributions GLOSSARY REFERENCES

 
Cell Culture and Selection for Drug Resistant Cell Populations
Primary tumors were assigned a random, two-letter code designation, and the recurrent tumor from the same patient received the same code with the addition of "R" (ME/MER, LX/LXR, DI/DIR). Cell lines were derived from these tumors as described^13,14 and grown in Waymouth MAB 87/3 medium (MAB; Mediatech Inc, Herndon, VA) with 20% fetal calf serum (FCS; Mediatech Inc). Cells from primary and recurrent tumors were selected for resistance to BCNU or TMZ as described.^14,15 Cells were washed with MAB without serum three times; they were then mock-treated using MAB alone, or treated with increasing concentrations of BCNU (2.5, 5.0. 7.5, and 10 µg/mL) or TMZ (2.5, 5.0, 7.5, and 10 µmol/L) in MAB for 1 hour at 37°C with 5% CO₂. Cells were washed and fed with MAB containing 20% serum. The cells were treated for 3 (BCNU) or 5 (TMZ) consecutive days, after which the cells were allowed to grow. Corresponding sets of cells were mock-treated in parallel using MAB media without drug. This step was repeated several times until the resulting cell population was resistant, as evidenced by the absence of cell death after treatment when compared wtih the mock-treated controls. The time required to select for a resistant cell population varied for the different cell lines. Cells from the recurrent tumors had a higher level of intrinsic resistance than cells from the primary tumor. Cells were re-treated with 10 µg/mL BCNU or 10 µmol/L TMZ every eight to 10 passages to maintain the resistant phenotype.

Tumor Specimens and Patients
Thirty-one fresh-frozen glioblastoma specimens were collected and subjected to standard WHO classification.¹⁶ Patients underwent tumor debulking surgery (gross total resection, 84%; subtotal resection, 16%) and were generally treated with an adjuvant regimen that included irradiation (total of 60 Gy) and TMZ (150 to 200 mg per square meter for 5 days during each 28-day cycle). All but four patients were treated with chemotherapy, and one patient also did not receive radiotherapy. For four patients, data on adjuvant therapy were incomplete. Written informed consent was obtained from all patients, and the study was approved by the institutional review board of Stanford University Medical Center (Stanford, CA).

RNA and DNA Preparation
For RNA extraction from cell lines and tumor specimens, cell lysates and samples were homogenized using QIAshredder columns (Qiagen, Valencia, CA) and a rotor-stator homogenizer (Kinematica, Cincinnati, OH), respectively. Total RNA was isolated from cell and tumor homogenates using the RNeasy Mini and RNeasy Lipid Tissue Kits (Qiagen), respectively, and quantified via spectrophotometry. RNA integrity was confirmed using the Agilent 2100 Bioanalyzer (Agilent, Palo Alto, CA). Universal human reference total RNA was purchased from Stratagene (Strategene, La Jolla, CA). Genomic DNA from cell lines was isolated using the DNeasy Tissue Kit (Qiagen), DpnII (New England Biolabs, Beverly, MA) digested, and purified using the QIAquick PCR Purification Kit (Qiagen). Genomic DNA and genomic digest were quantified spectrophotometrically. Human male and female genomic reference DNA was purchased from Promega (Madison, WI).

Microarray-Based Gene Expression Profiling
An indirect, dendrimer-based labeling method¹⁷ was used for microarray hybridization that utilized the Genisphere 3DNA Array 900 labeling system (Genisphere, Hatfield, PA), following the procedural protocol provided by the manufacturer without any modifications. For cDNA synthesis, 3 µg of glioblastoma and universal human reference total RNA were separately reverse transcribed using the Cy5- and Cy3-specific Genisphere primers, respectively, and hybridized together overnight at 65°C to a Stanford human cDNA microarray containing 41,421 cDNA elements, corresponding to 27,290 different UniGene cluster IDs. Microarrays were coated with DyeSaver2 (Genisphere) immediately after the last wash.

Array-Based Comparative Genomic Hybridization
Labeling of digested DNA and microarray-based comparative genomic hybridizations (array-CGH) were performed essentially as described,¹⁸ with slight modifications. Two micrograms of DNA were labeled using random primers (Bioprime Labeling Kit; Invitrogen, Carlsbad, CA). Tumor DNA and reference DNA were fluorescently labeled with Cy5 and Cy3 dye (Amersham Biosciences, Piscataway, NJ), respectively. Tumor DNA was hybridized together with sex-matching reference DNA to the same Stanford human cDNA microarray as mentioned previously.

Data Normalization and Filtering
Microarrays were scanned on a GenePix 4000B scanner (Axon Instruments, Union City, CA). Primary data collection was performed using GenePix Pro 5.1 software. Raw data were deposited into the Stanford Microarray Database. Data were background corrected, filtered using a flag and background filter (1.5 minimal signal-over-background ratio for expression arrays in either channel; 2.5 minimal signal-over-background ratio in the reference channel and regression correlation > 0.6 in both channels for array-CGH), and normalized by the LOWESS normalization function using standardization and normalization of microarray data (SNOMAD) data analysis tools (http://pevsnerlab.kennedykrieger.org/snomad.htm) or the Institute for Genomic Research Microarray Data Analysis System (TIGR-MIDAS) function of the TM4 microarray software suite (www.tigr.org/software/tm4/midas.html). The GoldenPath Human Genome Assembly (http://genome.ucsc.edu, National Center for Biotechnology Information build 34) was used to map fluorescence ratios of the arrayed human cDNAs to chromosomal positions. Chromosomal copy number maps were generated by mean filtering of signal intensity ratios according to 5-megabase (mb) windows moved across the chromosomes in 2.5-mb steps. Gene copy number values were deemed changed as compared with normal human reference DNA if they fell beyond the plus or minus three standard deviations (SDs) range (± 0.2135) of distribution of all signal-intensity ratios of control self-to-self hybridizations. For gene-by-gene integration of gene copy number and gene expression, copy numbers were reported as symmetric three-nearest-genomic-neighbors moving averages.¹⁹ The TreeView software²⁰ was used to display gene expression and gene copy number ratios.

In the cell line model, 9,734 of 37,860 clones with expression in 80% of samples and whose expression levels differed by at least three-fold, in at least one sample, from their mean expression levels across all cell lines were included in downstream statistical analyses. One-class response significance analysis of microarrays (SAM),²¹ which corrects for multiple testing by assigning a false discovery rate–based measure of significance, called q value,²² was utilized after parental transformation of gene expression ratios of resistant sublines to identify genes overexpressed and underexpressed in all in vitro and in vivo BCNU- and TMZ-resistant sublines. Genes identified with q < 0.005 were deemed significant. Nonparametric t-test analysis was performed in R,²³ and was used to allocate additional statistical confidence to clones identified by SAM. Unsupervised hierarchical clustering was performed in Cluster,²⁰ and two-way average linkage clustering was applied, on the basis of Pearson correlation as a distance metric. Principal component analysis based on Pearson correlation was executed in MATLAB (The MathWorks, Natick, MA).

Gene Ontology and Functional Network Analysis
Analyses of gene ontology, canonical pathways, and functional networks were executed using tools from Ingenuity Pathways Analysis (Ingenuity Systems, Mountain View, CA), a Web-delivered application that enables the discovery, visualization and exploration of molecular interaction networks in gene expression data. The gene list identified by SAM, containing GenBank (National Institutes of Health, Bethesda, MD) accession numbers as clone identifiers as well as d-scores, was uploaded into the Ingenuity pathway analysis. Each clone identifier was mapped to its corresponding gene object in the Ingenuity Pathways knowledge base, which represents a proprietary ontology of 300,000 classes of biologic objects spanning genes, proteins, cells and cell components, anatomy, molecular and cellular processes, and small molecules. Semantically consistent pathway relationships are modeled based on a continual, formal extraction from the public domain literature and cover more than 10,300 human genes (www.ingenuity.com/products/pathways_knowledge.html). These mapped focus genes were then used as a starting point for generating biologic networks. A score was computed for each network according to the fit of the original set of significant genes. This score reflects the negative logarithm of the P value that indicates the likelihood of the focus genes in a network being found together as a result of random chance. Using a 99% confidence level, scores 2 were considered significant. Significances for the enrichment of the genes in a network with particular biologic functions or canonical pathways were determined via right-tailed Fisher's exact test with = .05 and the whole database as a reference set. The same computation was used for gene ontology analyses of the initial gene list.

Real-Time Reverse Transcription Polymerase Chain Reaction
Quantitative real-time reverse transcription polymerase chain reaction (qRT-PCR) reactions were performed with the ABI Prism 7900HT Sequence Detection System using SYBR GREEN PCR Master Mix (Applied Biosystems, Foster City, CA). Primers targeting the transcripts of TNFAIP3, NF-B inhibitor (NFKBIA), chromosome 8 open reading frame 4 (C8orf4) and leukemia inhibitory factor (LIF) genes and the glyceraldehyde-3-phosphate dehydrogenase (GAPDH) housekeeping gene were designed with the Primer3 program (http://frodo.wi.mit.edu/cgi-bin/primer3/primer3_www.cgi) and synthesized at the Stanford PAN Facility (for sequences, see Supplementary Table 1, online only). Total RNA was reverse transcribed using the SuperScript first-strand synthesis system and SuperScript II (Invitrogen). Thermocycling for each PCR reaction was carried out in a final volume of 20 µL containing 1 ng of cDNA, forward and reverse primers at 3 µmol/L final concentration, and 1x SYBR GREEN PCR Master Mix. After 10 minutes of initial denaturation at 95°C, the cycling conditions of 40 cycles consisted of denaturation at 95°C for 15 seconds, annealing at 55°C for 30 seconds, and elongation at 72°C for 30 seconds. All reactions were performed in triplicate. Dissociation curve analysis was performed after every run to confirm the primer specificity. Gene quantities were determined using standard curves, constructed by five serial dilutions of RT product of universal human reference RNA (Stratagene), and gene expression levels were reported as ratios of quantities of the target transcript and GAPDH as the reference transcript.

Immunoblotting
The light-enhanced chemiluminescence protocol was used for the detection of specific proteins from total cell lysates prepared using 1x RIPA buffer (150 mmol/L NaCl, 10 mmol/L Tris [pH, 7.2], 0.1% [weight/volume] sodium dodecyl sulfate, 1.0% [volume/volume] Igepal CA-630 [Sigma-Aldrich, St. Louis, MO], 0.5% [weight/volume] sodium deoxycholate, 5 mmol/L EDTA). Blots were exposed to 2 µg/mL of an anti-TNFAIP3 monoclonal antibody (Abcam, Cambridge, MA), recognized by a horseradish peroxidase (HRP) –conjugated goat antimouse secondary antibody (Santa Cruz Biotechnology, Santa Cruz, CA). The ECL-Plus detection system (Amersham Biosciences) was used according to the manufacturer's protocol. Blots were stripped using the Re-Blot Plus mild antibody stripping solution (Chemicon International, Temecula, CA), and reprobed with two independent loading controls including the anti-ß-actin goat polyclonal (Santa Cruz Biotechnology) antibody and the anti-GAPDH mouse monoclonal (Abcam, Cambridge, MA) antibody. A separate gel was stained with the SimplyBlue gel stain (Invitrogen) to ensure equal protein loading. Bands were quantified on an AlphaImager 2200 (Alpha Innotech, San Leandro, CA), and TNFAIP3 expression normalized to loading controls.

NF-B DNA-Binding Activity Assay
Nuclear protein extracts were prepared using the NucBuster Protein Extraction Kit (Novagen, Madison, WI) according to the manufacturer's instructions. DNA-binding activity of NF-B was assayed colorimetrically, utilizing the NoShift Transcription Factor Assay Kit and NoShift NF-B (p65) reagents (both Novagen) according to the manufacturer's instructions. To assess sequence-specific binding activity, 15 µg of sample nuclear extract or 25 µg of HeLa positive control nuclear extract were incubated with various combinations of biotinylated NF-B wild-type dsDNA, specific NF-B competitor dsDNA lacking biotin end labels, and nonspecific, nonbiontinylated dsDNA with a mutant NF-B consensus binding motif. Negative controls consisted of reactions performed in the absence of a binding sequence. HRP-conjugated goat antimouse immunoglobulin G, targeting an anti–NF-B (p65) mouse monoclonal antibody, was used as secondary antibody. All assays were performed in triplicate. Binding activity was measured via colorimetric absorbance at 450 nm on a ThermoMax multiwell spectrophotometer (Molecular Devices, Sunnyvale, CA) using 3,3',5,5'-tetramethylbenzidine as substrate.

Survival Analysis
Overall survival was calculated from the date of tumor diagnosis until death or the last follow-up contact. Data were current as of January 1, 2005. At last follow-up 28% of patients were alive and 72% were dead. Patient subgroups were defined by unsupervised clustering of patients based on gene expression data of the resistance-associated transcripts revealed by SAM in the resistance model and with expression in > 75% of the tumor specimens. Actuarial survival curves between groups were estimated by the Kaplan-Meier product-limit method, and survival distributions between groups were compared using the log-rank test. Univariate and multivariate Cox proportional hazards regression analyses were performed with overall survival as the dependent variable.

	RESULTS

TOP ABSTRACT INTRODUCTION PATIENTS AND METHODS RESULTS DISCUSSION Authors' Disclosures of... Author Contributions GLOSSARY REFERENCES

 
Resistance Signature of Glioblastoma Cells to O⁶-Alkylating Agents
One-class response SAM, following parental transformation, was used to identify gene expression patterns associated with resistance formation to O⁶-alkylating agents. A d-score was assigned to each gene on the basis of change in gene expression relative to the SD of repeated measurements. Permutations of the repeated measurements estimated the q value, a false discovery rate–based measure of significance.²¹ The three sensitive parental sublines DI, LX, and ME were analyzed against the pool of 15 sublines with in vitro resistance to BCNU (DI-B and ME-B), in vivo resistance to BCNU (DIR, LXR, MER), combined in vivo/in vitro resistance to BCNU (DIR-B and MER-B), in vitro resistance to TMZ (LX-T, ME-T) and combined in vivo resistance to BCNU and in vitro resistance to TMZ (DIR-T, LXR-T, MER-T). This analysis revealed a set of 329 transcripts consistently overexpressed (78 clones, 23.7%) or underexpressed (251 clones, 76.3%) in the resistant versus the parental sublines (q < 0.005; Fig 1A; Supplementary Table 2, online only). The well-established resistance marker MGMT³ was among the top-scoring overexpressed transcripts (q < 0.003). High-scoring underexpressed transcripts included the NF-B pathway–modulator TNFAIP3, which encodes the zinc finger protein A20²⁴; the NF-B–inhibiting IB family member NFKBIA; the metastatic colon cancer–downregulated C8orf4²⁵; and the astrocyte differentiation-associated LIF²⁶ (q for all < .003). Figure 1A displays these clones ordered according to d-score significance. Nonparametric t testing in the context of multiple testing was performed to allocate additional confidence to this resistance signature. One hundred forty-one clones passed a P-value threshold < .01 (Supplementary Table 3, online only). Although this analysis substantiated high significances for TNFAIP3 (P < .000001), C8orf4 (P < .000001), NFKBIA (P = .000005), and LIF (P = .000005) in the resistance phenotype, MGMT passed a P value filter < .05 only.

View larger version (47K): [in this window] [in a new window]   Fig 1. Transcriptomic resistance signature of glioblastoma cells to O6-alkylating agents. (A) Heat map displaying the expression of 329 transcripts identified by one-class response significance analysis of microarrays (SAM) to be significantly linked to the resistant phenotype (q < 0.005). Corresponding d scores. Expression levels in the resistant variants are normalized to their corresponding parental cell line (left; green, underexpression; red, overexpression). (B) Unsupervised, two-way average linkage clustering of 141 transcripts, identified by filtering the 329 SAM transcripts using a nonparametric t test in the context of multiple testing with a P value threshold of .01. (C) Dendrogram with cell line labels corresponding to part B. (D) Principal component analysis in the same subset. -R, recurrent; -M, mock; -B, BCNU; -T, TMZ; mock, no drug; TMZ, temozolomide; BCNU, 1,3-bis(2-chloroethyl)-1-nitrosourea; Res, resistant; Sen, sensitive.

Genome-Wide Gene Copy Number Maps
Because gene copy number aberrations (CNAs) have a significant impact on gene expression patterns and also represent common mechanisms of gene activation and inactivation in drug resistance formation, we have mapped gene copy numbers in a high-resolution manner using the same cDNA microarray platform. Figure 2A reports the copy number profile for 33,587 clones mapped along the genome and mean filtered according to 5-mb windows moved across the chromosomes in 2.5-mb steps. Multiple resistance-associated CNAs, including both losses and gains, were found, with the extent of the chromosome involved in the CNAs differing for each chromosome. Although some areas of the profile were similar across most of the samples (eg, chromosome 9q33.2-q33.3), others showed differences based on the individual patient's tumor (eg, chromosome 13), or on whether the cells were from the primary or recurrent tumor (eg, chromosome 7q11.21-q31.31). The left panel of Figure 2A shows the recurrences frequencies of these changes in the resistant variants. Peak recurrences of chromosomal gains were noted for chromosomes 9p23-9p22.3, 10q21.1, 11q14.1, 14q23.2-q31.1, and 17p11.2, as were those of losses for chromosomes 1p35.2-p34.3, 8p23.2-p23.1, 9q33.2-q33.3, 22q11.1-q21.1, and Xp21.2-p21.1.

View larger version (49K): [in this window] [in a new window]   Fig 2. Genomic changes associated with resistance formation of glioblastoma cells to O6-alkylating agents. (A) Highly compressed genome-wide gene copy number map, with 33,587 cDNA clones ordered as a function of genome position (left, chromosome identifiers) and mean filtered according to 5-mb windows moved across the chromosomes in 2.5-mb steps. Gene copy numbers in the resistant variants are normalized to their corresponding parental cell lines. Resistance-associated gene copy number aberrations (CNA) included both gains (red) and losses (green). Corresponding CNA recurrence frequencies (left; dark red, gains; dark green, losses) based on a ± three-standard-deviations threshold of data distribution, plotted for all resistant variants and aligned according to genome order. (B) Coincident alterations in gene copy number and gene expression in a subset of the resistance signature. Copy numbers in resistant cells normalized to their corresponding parental cells are reported as symmetric three-nearest-genomic-neighbors moving averages. Fluorescence ratios indicating genes with less than ± two-fold change in gene dosage in the resistant versus the sensitive sublines have been masked. -R, recurrent; -M, mock; -B, BCNU; -T, TMZ; mock, no drug; TMZ, temozolomide; BCNU, 1,3-bis(2-chloroethyl)-1-nitrosourea; Res, resistant; Sen, sensitive.

View larger version (26K): [in this window] [in a new window]   Fig 3. Candidate resistance genes in glioblastoma cells. (A) Gene expression by microarray and real-time reverse-transcription polymerase chain reaction (qRT-PCR) of four resistance genes. Bar graphs indicate the microarray-assessed gene expression in resistant cells relative to the corresponding parental cells, the heat map correlates the parentally normalized expression between microarray (array-EXP) and qRT-PCR and reports corresponding parentally transformed gene copy numbers (microarray-based comparative genomic hybridization [array-CGH]). Heat maps have been masked (black squares) to only show fluorescent ratios indicating at least ± two-fold changes in resistant cells versus sensitive cells. (B) Mean expression levels of the genes in the resistant variants versus parental cells, error bars indicate the spread of the expression across all resistant variants. (C) Gene expression and gene copy numbers for the tumor necrosis factor-–induced protein 3 (TNFAIP3) regulated replication protein A–interacting protein (RIP) and the established resistance factor O6-methylguanine DNA methyltransferase (MGMT). -B, 1,3-bis(2-chloroethyl)-1-nitrosourea; -M, no drug; T, temozolomide.

View larger version (79K): [in this window] [in a new window]   Fig 4. Functional network analysis of resistance signature. (A) Network map of molecular interactions and subcellular distribution of resistance genes, composed of five significant subnetworks (networks 1, 2, 3, 5, and 6 of Table 1) that demonstrated at least two overlapping nodes. Nodes represent genes, with their shapes representing the functional classes of the gene products, and edges indicate the biologic relationship between the nodes, which include physical and functional interactions. Nodes are color-coded according to their d-score (red, overexpression; green, underexpression). (B) Excerpt of the top-scoring functional network built around nuclear factor-B (NF-B), which included several altered NF-B pathway members and modulators, such as tumor necrosis factor-–induced protein 3 (TNFAIP3) and NF-B inhibitor (NFKBIA), as well as NF-B–interacting genes.

Exploration of the NF-B Candidate Pathway

Candidate Genes
Gene expression levels for the NF-B pathway genes TNFAIP3 and NFKBIA were confirmed by real-time RT-PCR. The relative expression of the two transcripts C8orf4 and LIF, deemed biologically interesting because of their highly significant association with the resistance phenotype, was also confirmed by real-time RT-PCR. C8orf4, which is downregulated in metastatic colon cancer, has been implicated in colon cell differentiation and transforming growth factor ß–induced apoptosis.²⁵ The neuropoietic cytokine family member LIF regulates gliogenesis²⁹ and promotes differentiation of astrocytes.^26,30 Normally, this gene is expressed in glioma cells³¹ and mediates a growth inhibitory effect in these cells.³² Figure 3A correlates the expression levels of these four genes revealed by the microarray and the real-time RT-PCR analyses, and also indicates the gene copy number profiles for these genes. There was high concordance between both analyses in individual sublines (Fig 3A), as well as when the mean transcription levels of all resistant sublines normalized to their corresponding parental cells and their spreads were compared (Fig 3B). Figure 3C interfaces gene expression and copy numbers for RIP and MGMT, indicating the relationship of RIP gene dosage to gene expression in the DI/DIR cell lines, but not LX/LXR and ME/MER. MGMT gene dosage correlates with expression in some, but not all of the cell lines examined.

Correlation of TNFAIP3 Protein and Gene Expression
Because our integrated large-scale screening and inferential statistical approach for resistance factors had highlighted a potentially important role of the NF-B pathway modulator TNFAIP3 in the resistance of glioblastoma cells to O⁶-alkylating agents, we examined whether the reduced expression of the TNFAIP3 transcript in the resistant cells may be reflected by reduced expression of the cognate protein product A20. Figure 5B shows a representative immunoblot in ME cell lines. Figure 5A indicates that the level of A20 protein expression in these cell lines normalized to loading controls. A substantial reduction in protein expression in all resistant variants compared with the sensitive parental subline was noted. Figure 5A also correlates TNFAIP3 protein and gene expression in these cells. In each cell line, the levels of protein product closely mirrored the expression of the corresponding transcript as determined using real-time RT-PCR.

View larger version (22K): [in this window] [in a new window]   Fig 5. Tumor necrosis factor-–induced protein 3 (TNFAIP3) expression and nuclear factor-B (NF-B) activation in resistant glioblastoma cells. (A) TNFAIP3 gene and protein expression in ME cell lines. Expression levels in the resistant variants are indicated as percentages of expression in the parental cells. Gene expression, assessed by real-time reverse transcriptase polymerase chain reaction, is normalized to the glyceraldehyge-3-phosphate dehydrogenase (GAPDH) housekeeping transcript. Protein expression, assessed by immunoblotting, is normalized to loading controls (GAPDH, ß-actin) from two independent blots. (B) Representative blot with TNFAIP3 (approximately 70 KDa), ß-actin (approximately 43 KDa) and GAPDH (approximately 36 KDa) expression. (C) Competitive analysis of nuclear factor-B (NF-B) DNA-binding activity in the same cell lines and HeLa positive control nuclear extract. Black bars, binding activity assessed by biotinylated NF-B wild-type (WT) double-stranded DNA (dsDNA); dark gray bars, WT dsDNA plus nonspecific, nonbiontinylated dsDNA with a mutant NF-B consensus-binding motif (NSC); light gray bars, WT dsDNA plus specific NF-B competitor dsDNA lacking biotin end labels (SC). -R, recurrent; -B, 1,3-bis(2-chloroethyl)-1-nitrosourea; Neg, negative control.

DNA-Binding Activity of NF-B in Sensitive Versus Resistant Cells

Outcome Prediction Model Based On Resistance Signature
We evaluated the clinical impact of our resistance signature in an independent cohort of glioblastomas commonly treated with O⁶-alkylating agents. In order to create a model for predicting survival and response in glioblastomas, we queried the expression status of the 329-transcript resistance signature derived from the cell line model in the gene expression profiles of 31 glioblastomas from a different institution (Stanford University Medical Center, Stanford, CA). One hundred seventy-two transcripts were expressed in at least 75% of the tumors (Supplementary Table 4, online only). Unsupervised hierarchical clustering of the tumors based on these transcripts, applying two-way average linkage clustering based on Pearson correlation as a distance metric, resulted in two major tumor subgroups (group 1, 12 patients; group 2, 19 patients) with distinct gene expression signatures (Fig 6A). Highly correlated expression behavior was observed for two major gene clusters (Fig 6A). One of these clusters demonstrated persistent overexpression in group 1 and persistent underexpression in group 2, and included two coclustering transcripts of TNFAIP3 and one transcript each of NFKBIA, C8orf4, and LIF. The other gene cluster showed overexpression in a subset of tumors in group 2 and included the MGMT transcript. On the basis of these gene expression signatures, we labeled group 1 as a potentially favorable tumor group and group 2 as an unfavorable tumor group. Actuarial survival analysis revealed a significant difference in survival between the two groups (log-rank P = .037; Fig 6C). The overall survival rates at 2 years in groups 1 and 2 were 0.53 and 0.09, and the median survival times were 814 and 412 days, respectively (Fig 6B).

View larger version (48K): [in this window] [in a new window]   Fig 6. Outcome prediction model in 31 glioblastomas. (A) Results of unsupervised hierarchical cluster analysis based on the expression of 172 of 329 resistance-associated transcripts with > 75% expression in the tumors. Two major tumor-classifying gene clusters are denoted by vertical yellow and purple bars. Tumor necrosis factor-–induced protein 3 (TNFAIP3), NF-B inhibitor (NFKBIA), chromosome 8 open reading frame 4 (C8orf4) and leukemia inhibitory factor (LIF), and O6-methylguanine DNA methyltransferase (MGMT) transcripts are indicated by black bars. (B) Kaplan-Meier estimates of overall survival (OS) in all 31 glioblastomas. Dotted lines indicate the 95% CI. (C) Kaplan-Meier estimates of overall survival in the 31 glioblastomas after subdivision into the two groups based on the clustering results of part A (group 1, favorable; group 2, unfavorable; log-rank P = .037). (D) Hierarchical clustering result based on the top-10 genes sufficient in driving the tumor subgrouping into the two main classes in part A (class error rate, 0.03) and associated significantly with OS according to log-likelihood estimate (P = .016). (E) Optimized and simplified predictive model based on the weighted expression of four resistance-associated genes. According to Cox proportional hazards regression analysis, P = .022 for the model as a continuous variable, and log-rank P = .007 for the model as a class based on the two major subgroups defined by unsupervised hierarchical clustering. (F) Partitioning of the tumors into two discrete subgroups (group 1, favorable; group 2, unfavorable) by the four-gene predictor model (see Outcome Prediction Model Based On Resistance Signature section for details; log-rank P = .026).

Because a focused number of genes would be particularly amenable to future target modulation we optimized the model by further minimizing the number of predictive genes on the basis of SAM ranking. Hierarchical clustering according to the top four transcripts (CD44 antigen [CD44], F-box protein 32 [FBXO32], syndecan 1 [SDC1], and TNFAIP3) revealed two major tumor subgroups (14 v 17 tumors), which demonstrated a significantly different outcome in actuarial survival analysis (log-rank P = .007; Fig 6E). Because of the apparent difference in the expression of these genes between the two groups, we performed a Cox proportional hazards regression analysis of the individual genes, with overall survival as dependent variable. TNFAIP3 and CD44 were associated significantly with outcome when considered as individual continuous variables (univariate P = .028 and .032, respectively). We then examined the combined predictive potential of all four genes, using the average expression of the genes weighted for their relative predictive contribution as indicated by the individual Cox scores. As a continuous variable, the combined predictive model performed better in outcome prognostication compared with the individual genes (univariate P = .022). This model remained significantly associated with patient outcome when taking into account patient age (multivariate P = .025), the most important clinical prognosticator for glioblastomas.

Finally, we examined how the model would perform as a class. We used a schema based on median gene expression to subtype the tumors. Here, the expression of each of the four genes was labeled either high or low based on the median level of all tumors and the number of calls was counted. This approach generated a simplified gene expression profile for each tumor, represented by an integer between 0 and 4 for the number of low- and high-expression calls. If a tumor had at least three low calls, it was classified as potentially unfavorable (group 2); if there were one or fewer low calls, the tumor was classified as favorable (group 1); if a tumor had an equal number of low and high calls, the patient's prognosis was not inferred on the basis of the expression data of the four genes ("noninformative" group 3). Such stratification revealed 14 tumors each to fall into groups 1 and 2, and three tumors into group 3. Actuarial survival analysis disclosed a strikingly different outcome between groups 1and 2 (log-rank P = .026), with favorable tumors demonstrating a comparably good prognosis (Fig 6F).

	DISCUSSION

TOP ABSTRACT INTRODUCTION PATIENTS AND METHODS RESULTS DISCUSSION Authors' Disclosures of... Author Contributions GLOSSARY REFERENCES

 
We have identified a distinct genomic signature shared by glioblastoma cells selected for resistance to O⁶-alkylating agents under in vitro and in vivo conditions. Our results highlight the involvement of a cellular pathway of NF-B–mediated resistance to these agents in glioblastoma cells. The contribution of NF-B to anticancer drug resistance has been described in various in vitro and in vivo resistance models.¹⁰ The antiapoptotic activity of NF-B appears to be the most important mode of action mediating the resistance and pro-survival effects of this gene in cancer cells.^10,34 Genotoxic stress resulting from the exposure of tumor cells to O⁶-alkylating agents causes DNA damage and leads to the initiation of apoptosis. NF-B activation abrogates the apoptosis signal in response to these agents.¹⁰ Antagonism of NF-B in malignant gliomas has been shown to render glioma cells more susceptible to BCNU via increased apoptosis,¹² but the DNA damage–induced signaling pathway upstream of IB has not been identified in these cells.

We have found alterations in several NF-B pathway members in glioblastoma cells, which may act synergistically in activating NF-B during resistance formation to O⁶-alkylating agents. The most significant link to resistance to these agents within the NF-B pathway was revealed for TNFAIP3. The consistent downregulation of this gene in our resistance models suggests a potentially important role of this gene in the development of resistance to O⁶-alkylating agents in glioblastoma cells. Protein levels of TNFAIP3 related directly to the expression of the corresponding transcript, as well as to the levels of NF-B–DNA-binding activity of the cells.

The zinc finger protein TNFAIP3³⁵ is a potent inhibitor of NF-B signaling.³⁶ Its mechanism of action involves the cooperative activity of its two ubiquitin-editing domains.²⁴ The amino-terminal domain of TNFAIP3,³⁷ removes lysine-63–linked ubiquitin chains from RIP.²⁴ The carboxy-terminal domain, composed of seven C₂/C₂ zinc fingers,³⁶ then functions as an ubiquitin ligase by polyubiquitinating RIP with K48-linked ubiquitin chains, thereby targeting RIP for proteasomal degradation.²⁴ In addition to consistent downregulation of TNFAIP3 in resistant glioblastoma cells, which frees RIP, we also found a significant upregulation (q < 0.05) of RIP in these cells. This observation supports a cooperative or reciprocal molecular mechanism for these two genes in the resistant phenotype. RIP, which complexes with IB, has been reported to have an essential role in DNA damage–induced NF-B activation, but not in UV-induced NF-B activation.²⁷ RIP-mediated NF-B activation by drug- and irradiation-induced DNA damage is not mediated by autocrine or tumor necrosis factor (TNF) receptor 1 signaling. In line with recent evidence indicating that TNFAIP3 is critical for the regulation of TNF-independent signals that lead to termination of NF-B activity,³⁸ our data provide evidence for an upstream extension of a RIP-mediated signaling cascade that augments NF-B–induced resistance in glioblastoma cells to DNA-alkylating agents. Although the initiation point of this resistance signal remains ambiguous, our data support a model that extends its cytoplasmic pathway to TNFAIP3. We therefore propose a dual mechanism that may synergistically foster the activation of NF-B by these agents. On receiving the nuclear signal in response to DNA damage, upregulation of RIP has been hypothesized to initiate the cytoplasmic signaling that activates NF-B.²⁷ RIP action and its complexing with IB may be enabled and ameliorated by downregulation of TNFAIP3 in the resistant cells, which under normal cell conditions obscures RIP by targeting it to proteasomal degradation. Functional validation will be needed to confirm this candidate resistance pathway in glioblastomas.

Our resistance signature–based outcome predictor model, derived from the cell lines, enabled us to subcategorize an independent cohort of glioblastomas commonly treated with O⁶-alkylating agents into two major groups with apparently different outcomes. On the basis of this link, we reasoned that the observed difference in survival may have been the result of distinct response characteristics of these tumors to therapy rather than different biologic tumor behavior. The inclusion of more genes may make our predictive model perform better in independent validation analyses, but a smaller number of genes would make the model more practical and also amenable to future target modulation approaches. We have therefore reduced our resistance signature to a minimal number of genes for use in constructing a predictive model. We have shown that measurement of the weighted expression of four resistance-related genes was sufficient in predicting patient outcome. This optimized predictor, which included the four transcripts SDC1, CD44, FBXO32, and TNFAIP3, was able to partition glioblastomas into two subgroups according to their survival. The identification of TNFAIP3 as part of this predictor corresponds to its significant association with both resistance formation and NF-B activation in our in vitro and in vivo cell models. It is supported as well by strong evidence for a link between NF-B and glioblastoma cell survival³⁹ and glioblastoma cell resistance to cytotoxic therapy.^12,40

We have found reduced expression for all four genes associated with poor patient outcome. In line with our results, loss of SDC1, a transmembrane type I heparan sulfate proteoglycan and member of the syndecan proteoglycan family, has been linked to unfavorable prognosis of various human malignancies, including squamous cell carcinoma of the head and neck,⁴¹ laryngeal cancer,⁴² poorly differentiated non–small-cell lung carcinoma,^43,44 hepatocellular carcinoma with high metastatic potential,⁴⁵ and gastric cancer.⁴⁶ In addition, SDC1 has been reported as a predictor of chemotherapy efficacy in oral squamous cell carcinoma, with decreased expression in response to cytostatic treatment indicating a poor prognosis.⁴⁷ CD44 has various functions in cell-cell and cell-matrix interactions. Expression of this cell-surface glycoprotein has been linked to increased survival in node-negative, invasive breast cancer,⁴⁸ and to indicate favorable prognosis in epithelial ovarian cancer.⁴⁹ Lack of CD44 expression is also a highly significant factor of poor outcome in neuroblastoma.⁵⁰ CD44v6 has been shown to predict responses to treatment in advanced colorectal cancer.⁵¹ FBXO32,⁵² which constitutes a potential substrate-recognition component of the cell cycle–regulating SKP1-cullin-F-box ubiquitin protein ligase complex and functions in phosphorylation-dependent ubiquitination,^53,54 has not been associated with tumor prognosis and drug resistance so far. ß4GalNAc-T4, which is involved in protein glycosylation,⁵⁵ was the only high-ranking transcript for which increased expression was associated with unfavorable outcome. Although the role of ß4GalNAc-T4 in drug resistance remains enigmatic, altered protein glycosylation has been implicated in tumorigenesis.⁵⁶ A limitation of our study, however, was the relatively small sample size and the not completely uniform treatment of our patients. Assessment of our predictor in a larger, standardized patient population will be necessary to ultimately assign outcome significance to these genes in glioblastomas, and to refine parameters for risk-based stratification.

In summary, our results suggest a role of a cellular pathway that leads to NF-B activation during the emergence of acquired resistance to O⁶-alkylating agents in glioblastoma cells. Although our data indicate the alteration of various members of the NF-B canonical pathway, the endogenous NF-B inhibitor TNFAIP3 was linked most significantly to the resistance phenotype. TNFAIP3 gene as well as protein expression mirrored the level of NF-B activation in these cells. Though it remains unclear how the DNA damage response is linked to the cytoplasm, downregulation of TNFAIP3 may promote the initiation of an RIP-dependent signaling cascade that mediates NF-B–induced cell survival. TNFAIP3 related significantly to patient outcome in a cohort of glioblastomas, and was a member of an optimized four-gene outcome predictor that enabled the subcategorization of these tumors. These observations raise the hope for an amenable target to modulate NF-B–mediated resistance in glioblastomas cells, with the ultimate goal of increasing the efficacy of chemotherapy in patients harboring these challenging malignancies.

	Authors' Disclosures of Potential Conflicts of Interest

TOP ABSTRACT INTRODUCTION PATIENTS AND METHODS RESULTS DISCUSSION Authors' Disclosures of... Author Contributions GLOSSARY REFERENCES

 
The authors indicated no potential conflicts of interest.

	Author Contributions

TOP ABSTRACT INTRODUCTION PATIENTS AND METHODS RESULTS DISCUSSION Authors' Disclosures of... Author Contributions GLOSSARY REFERENCES

 

Conception and design: Markus Bredel, Claudia Bredel, Dejan Juric, Adrienne Scheck, Branimir Sikic Financial support: Markus Bredel, Griffith Harsh, Hannes Vogel, Lawrence Recht, Adrienne Scheck, Branimir Sikic Administrative support: Branimir Sikic Provision of study materials or patients: Griffith Harsh, Hannes Vogel, Lawrence Recht, Adrienne Scheck Collection and assembly of data: Markus Bredel, Claudia Bredel, Dejan Juric, George Duran, Adrienne Scheck Data analysis and interpretation: Markus Bredel, Claudia Bredel, Dejan Juric, Ron Yu Manuscript writing: Markus Bredel Final approval of manuscript: Branimir Sikic

 

	GLOSSARY

TOP ABSTRACT INTRODUCTION PATIENTS AND METHODS RESULTS DISCUSSION Authors' Disclosures of... Author Contributions GLOSSARY REFERENCES

 

Canonical pathway:

A core pathway established for a given molecule in the cell in which molecular interactions occur in a linear and stepwise manner. Although clustering expression data groups functionally related genes, they do not order pathway components according to physical or regulatory relationships. Software is now available for linking significant genes in one's experiments with a world collection of biologic networks created from millions of individually modeled relationships between genes, proteins, complexes, cells, and tissues. The software allows a view of one's data, integrated in biologic networks according to different biologic context and identifies canonical and noncanonical pathways that connect molecules within a biologic network.

CNA (copy number aberration):

Defines the net increase or net decrease in gene copy number (ie, "gene dosage"). It includes extensive low-amplitude and focal high-amplitude changes in gene copy number. Large regional copy number changes—involving chromosomal fragments, chromosomal arms, or whole chromosomes—are typically of low amplitude and commonly referred to as gains and losses. Focal reduction or increment in copy number of a specific gene in the genome without a proportional increase in other—or only few neighboring—genes are commonly denoted as amplification and deletion, respectively. Amplifications can include hundreds of repeat copies of a gene, which may sometimes occur via the excision of a copy of the repeating sequence from the chromosome and its extra chromosomal replication (so-called "double minute chromosomes"). In a deletion, one or both copies—so-called "alleles"—of a gene (which are present at an autosomal locus of the normal human genome) can be lost, commonly designated as hemizygous and homozygous deletion, respectively.

Gene ontology:

Allows for annotating genes and their products with a limited set of attributes, with the three organizing principles being molecular function, biological process, and cellular component. The development of structured, controlled vocabularies (ontologies) that describe gene products in terms of these organizing principles in a species-independent manner is a constantly evolving process.

Microarray-based comparative genomic hybridization (array-CGH):

Array-based comparative genomic hybridization is a method that uses microarrays to probe changes in chromosomal DNA, thereby identifying precise areas in which genetic changes occur in cancer cells.

Molecular interaction network:

Molecular interactions (eg, protein-protein and protein-DNA) made possible as a result of advances in high-throughput technologies are visualized as a network of interactions using software tools that have been devel-oped to visualize and analyze large-scale data. The network is visualized with molecules as nodes and molecular interactions as edges, with "force-directed layout algorithms" used to visualize the molecular interaction network. Some programs used are InterViewer, Pajek, Tulip, and Cytoscape. Increasingly, Web-based application programs (eg, WebInterViewer, Ingenuity Pathway Analysis) are used to produce molecular interaction networks of "good quality," because multiple molecular interaction networks for shared molecules and for interactions shared by all or some of the networks can be searched.

NF-B (nuclear factor kappa B):

NF-B is a transcriptional factor involved in the transcriptional activation of genes that regulate different cellular processes. Its nuclear location is restricted by its interaction with an inhibitor, I-B, which sequesters it in the cytoplasm. When I-B is phosphorylated and degraded in response to different stimuli, NF-B becomes free to enter the nucleus.

O⁶-alkylating agent:

Chemical agents that alkylate (eg, methylate) the O⁶ of guanine in DNA and are used in chemotherapy. Examples are 1,3-bis(2-chloroethyl)-1-nitrosourea (BCNU) and temozolomide (TMZ).

Permutation:

In statistics, the number of ways to count a set of elements. In performing molecular network analysis from microarray experiments, significance analysis of microarrays (SAM) computes a statistic value for each gene, which is a measure of the strength of the relationship between gene expression and the response variable. To do this, SAM uses repeated permutations of the data to determine if the expression of any genes is related significantly to the response.

SAM (significance analysis of microarrays):

A statistical technique using established software that determines the significance in changes of gene expression seen in microarray analysis (eg, cDNA and oligonucleotide microarrays), which measures the expression of thousands of genes and identifies changes in expression between different biologic states. On the basis of changes in gene expression relative to the standard deviation of repeated measurements, SAM assigns a score to each gene. When scores are greater than an adjustable threshold, permutations of repeated measurements are used by SAM to estimate the percentage of such genes identified by chance, the false discovery rate (FDR). In