BACKGROUND
the methylation of dna is involved in the control of chromatin folding, protein complex assembly, and gene expression.
it is a normal process of cellular development, aging, paternal/maternal genetic imprinting, and x-chromosome inactivation  <cit> .
dna methylation is the only known heritable modification of dna, and is aberrant in cancer, involving both excessive and insufficient methylation of dna regions of cytosine–guanine bonds , which then alters the expression of certain genes  <cit> .
the aberrant hypermethylation of dense cpg sites in many cancers occurs at tumor suppressor genes in a non-random, tumor-specific pattern  <cit> .
given the possible role of aberrant dna methylation in cancer initiation and progression, a significant effort has been directed to identify dna methylation biomarkers in cancer and to use such markers to predict therapeutic responses and disease prognosis  <cit> .
a notable example is the association between the silencing of the o6-methylguanine-dna methyltransferase  gene and the improved response to temozolomide and radiation in patients with glioblastoma  <cit> .
others include the association of the methylation of pten, igfbp- <dig>  and various dna repair enzymes  with either resistance  or sensitivity  to chemotherapy or radiation .
for lung cancer, dna methylation of genes may be useful for assessing cancer risk based on the analysis of the sputum of smokers  <cit> , and may be prognostic in early-stage lung cancer  <cit> .
to better understand the dna methylation changes associated with gene expression in lung cancer, we determined the methylation status of more than  <dig>  cpg sites across  <dig>  genes at selective promoter regions in  <dig> non-small cell lung cancer  cell lines using the illumina humanmethylation <dig> beadchip.
we integrated this dataset with the gene expression profiles of the same panel of cell lines to identify the genes with expression levels that were most affected by dna promoter methylation.
we sought to determine the biologic significance of genes that were significantly repressed in association with methylation in silico, as well as experimentally, both in vitro and in vivo.
RESULTS
identifying genes aberrantly associated with dna methylation in nsclc cell lines
for each of the  <dig> nsclc cell lines examined, we analyzed the relationship between the degree of methylation  of a particular gene and the expression level of that gene.
we sought to identify genes that were significantly repressed in association with dna methylation, which we subsequently refer to as srams, a term first used in the study of breast cancer  <cit> .
we computed the spearman correlation of the relationship between methylation and gene expression for all the cell lines, assuming that genes negatively regulated by dna promoter methylation have a negative spearman rho value.
this approach would only identify differentially regulated genes and not genes that are fully methylated or unmethylated in these cell lines.
we first tested our approach by looking at the spearman correlation for cdkn2a, a gene well known to be frequently downregulated by promoter methylation in nsclc  <cit> .
in the illumina assay, four probes assessed four cpg sites within the cdkn2a promoter.
we found that cdkn2a was differentially methylated in the various cell lines, and the degree of methylation was significantly and inversely correlated with gene expression .
this relationship was preserved in each of the four probes examined .figure  <dig>
differential genes are significantly repressed in association with methylation  in nsclc cell lines.
a) example of the relationship between gene expression and methylation status for a single probe interrogating a single cpg site within the promoter region of cdkn2a.
b) frequency histogram of the spearman rho values for all probes that correspond to either cpg or non-cpg island sites.
red curves are theoretical gaussian representation of the distribution of rho given the number of probes for each gene set and the actual mean and standard deviation values; blue curves represent the empirical distribution from the data; dashed burgundy lines represent normal distribution curves for comparison; and the arrow identifies the point where the two curves cross and an excess of probes are seen with a rho ≤ - <dig> .
c) relationship of rho to the distance from the transcriptional start site  for cpg and non-cpgi probes.
most of the probes  that are candidate srams cluster within 500 bp of the tss.
we expanded our analysis to the rest of the probes.
we looked at the distribution of rho values for probes located within or outside of cpg islands .
rho values of cpgi probes had a mean distribution that was significantly shifted to the left of  <dig> compared to that of non-cpgi probes .
for both sets of probes, an excess was seen in the experimental dataset at the rho cutoff of - <dig>  that was more than expected when compared to a random distribution .
this was not observed for probes that had a positive correlation with gene expression.
there were also no differences regarding whether the probes corresponded to cpgi or non-cpgi genes .
greater than 80% of these probes lie within 500 bp of the transcription starting site , corresponding to sites under transcriptional control by dna methylation .
we found  <dig> probes that correspond to  <dig> unique srams, genes with expression levels that were inversely correlated with the degree of dna methylation .
we found the same pattern as that of cdkn2a for other genes known to be regulated by methylation , as well as for many novel genes not previously known to be regulated by methylation .
we validated a set of srams using pyrosequencing and quantitative real-time polymerase chain reaction  in  <dig> cell lines .figure  <dig>
validation of the integrative analysis for select genes.
three genes with cpg islands were evaluated.
top row is the relationship between array methylation data and illumina gene expression; middle row provides pearson correlation of the array methylation data with the pyrosequencing methylation status; and bottom row shows the relationship of methylation status by pyrosequencing with the quantitative real-time pcr results for the expression of these genes.
srams differentiate tumor from normal lung and distinguish tumor subtypes
as aberrant dna methylation is one of the defining characteristics of cancer, srams should reflect alterations that discriminate cancer from normal tissue and reflect the underlying tumor biology.
we evaluated the methylation status of srams in a set of  <dig> lung tumors from a japanese cohort, hereafter termed the kyoto tumor set, which included  <dig> unmatched normal lung tissues arrayed using the illumina humanmethylation <dig> beadchip.
since the cell line sram dataset was generated from the humanmethylation <dig> platform, we took  <dig> probes  that overlapped between the two illumina platforms for the analysis.
we found that srams discriminated nsclc from normal lung tissue, with clusters of probes that were either predominantly methylated in normal tissue or methylated in tumors .
srams also distinguished squamous carcinoma from adenocarcinoma, the two major histological subtypes of nsclc, as well as lung adenocarcinomas with or without an egfr mutation .
these observations were corroborated using sram expression profiles of lung adenocarcinomas included in the director’s challenge consortium for the molecular classification of lung adenocarcinoma  <cit> .
srams preferentially separated egfr mutants from wild-type tumors .
in the kyoto tumor set, patients that separated into these clusters had significant survival differences, similar to that seen for the director’s challenge tumor set .figure  <dig>
nsclc srams are differentially methylated between normal and tumor tissues segregate into molecular pathologic subtypes.
a) heatmap of the hierarchical clustering of  <dig> matched probes for srams in the kyoto set of  <dig> tumors and compared to  <dig> unmatched normal lung tissue samples.
b) using matched  <dig> probes from the humanmethylation <dig> to the humanmethylation <dig> platforms, unsupervised clustering of  <dig> tumors from the kyoto tumor set was conducted, with tumor clustering along lines of histology , egfr mutation status, and tumor grades.
c-d) unsupervised clustering based on gene expression of the matched gene regulatory module  genes on the director’s challenge dataset of  <dig> lung adenocarcinomas.
the grms were able to segregate the tumors in the director’s challenge into two main groups, with egfr mutations and well-to-moderately differentiated tumors in group  <dig>  and wild-type egfr and moderately-to-poorly differentiated tumors in group  <dig>  e) sram segregates nsclc tumors into prognostic groups in the kyoto and director’s challenge tumor sets.
branch point c <dig> vs c <dig> mostly segregates squamous from adenocarcinoma, without a difference in prognosis; whereas branch point c <dig> vs. c <dig> is highly prognostic because of the separation of egfr mutant from wild-type tumors.
for the director’s challenge tumor set, the two clusters that preferentially separate egfr mutant and wild-type tumors are also prognostic for survival.
srams and the epithelial-to-mesenchymal transition  phenotype
to further define the functional significance of genes negatively regulated by dna promoter methylation, we performed a gene set enrichment analysis  using the spearman correlation rho value as the ranking variable  <cit> .
the top gene sets that included genes negatively regulated by dna promoter methylation had been obtained from diverse cancer datasets and were involved in various biologic functions, such as cell migration in bladder cancer and resistance to gefitinib in nsclc, and included genes downregulated with e-cadherin knockdown in human breast mammary epithelial  cell lines, genes differentially expressed in metastatic melanoma, genes methylated in glioblastoma and pancreatic cancer, and genes differentially regulated in luminal vs. basal/mesenchymal breast cancer cells .
we postulated that biologic functions identified by the gsea were related to the process of emt, a developmental and adaptive cellular process that has been associated with resistance to cancer therapies and regulation of metastasis  <cit> .
we performed a hierarchical cluster analysis of the cell lines using the srams, and overlaid the relative protein expression of e-cadherin, a cell adhesion molecule that is downregulated during emt and which plays a key role in the signaling and regulation of emt  <cit> , along with the expression of emt-related genes zeb <dig>  vim, twist <dig>  fn <dig>  cdh <dig>  and cdh <dig>  we found that the identification of srams clustered cells into two groups: “epithelial-like ” cells, which expressed high e-cadherin levels; and “mesenchymal-like ” cells, which expressed low e-cadherin levels .
to identify genes within the srams that distinguish e cells from m cells, we performed a wilcoxon rank sum test.
out of the  <dig> unique genes,  <dig>  were differentially methylated between e and m cells, with a p-value <  <dig> .
we called these genes emt-srams.
a majority of emt-srams  were methylated in m cells compared to e cells .
a few of the emt-srams, including cdh <dig>  were already known to be epigenetically regulated during the process of emt  <cit> .
however, the epigenetic regulation of the remaining emt-srams had not been established, including that of spint1
<cit> , syk
<cit> , mst1r
<cit>  and esrp1
<cit> .
interestingly, one of these genes is axl, which encodes a receptor tyrosine kinase that is associated with emt and with resistance to egfr tyrosine kinase inhibition  <cit> .
we found that the expression of axl was inversely correlated with promoter methylation , and its promoter methylation was positively correlated with e-cadherin expression , indicating that axl may be epigenetically regulated during emt.figure  <dig>
specific sram methylation patterns define epithelial  and mesenchymal  cell types.
a) hierarchical clustering of srams for a subset of cell lines  with rppa data for e-cadherin  segregates cell lines into high and low ecad groups.
representative genes that are preferentially unmethylated  or methylated  in the mesenchymal cells are shown.
illumina gene expression levels of emt markers , twist <dig>  fibronectin , cdh <dig>  and cdh1) are overlaid to confirm the e vs m identity of the cell lines.
b) wilcoxon rank sum test was used to identify the srams that are differentially methylated between e and m cells, which we call the emt-srams.
the identified  <dig> probes  were used for hierarchical clustering of the cell lines.
emt-srams that are predominantly methylated in m cells are called m-srams, and those methylated in e cells are called e-srams.
c) emt-sram association with emt network genes using curated network analysis.
closest second neighbor analysis of genes associated with core emt factors  finds numerous associations with genes within the emt-sram set.
red spheres represent srams that are methylated in m cells; green spheres are those methylated in e cells; white spheres represent linking genes.
to further determine whether emt-srams define an epigenetically regulated group of genes relevant to emt, we compared our gene list to a variety of publically available gene expression signatures from studies that experimentally manipulated cells to undergo emt.
specifically, we examined datasets from prostate cancer cells transduced with snai1
<cit> , a nsclc epithelial cell line  transduced with snai <dig>  zeb <dig> or stimulated with tgf-beta  <cit> , and hbme cells transduced with twist
<cit> .
we applied hypergeometric statistics to test the significance of overlap between the entire dataset of srams  and the emt-srams .
we found that the srams significantly overlapped with the genes in the published datasets.
on average, there was a 16% overlap between the srams and published emt gene signatures.
emt-srams had greater overlap with the published emt gene signatures, ranging from 18% in the breast cancer cell line signature to 48% in the h <dig> cells treated with tgf-beta .
the majority of overlapping genes were downregulated, which is consistent with our observation of a greater proportion of genes being methylated in m cell lines.
there appeared to be a greater overlap between genes from h <dig> cells and the srams or the emt-srams than either the prostate or breast cancer datasets.
this suggests that while there could be a common core group of genes for both the srams and the emt-associated srams, the tissue-specific srams may reflect the underlying biology of the tissue of origin.
using a curated network extracted from pathway commons  <cit> , which consists of  <dig>  genes and over a million biological interactions, we performed a network analysis using markers of emt, including  <dig> transcription factors and major regulators of emt , and three recognized markers of emt , and the emt-sram gene list.
we were able to map 73/ <dig> emt-srams.
the pair-wise shortest distance between markers of emt with emt-srams was calculated using dijkastra’s algorithm  <cit> : 51/ <dig> emt-srams  in the network had the shortest path lengths of  <dig> and  <dig>  connecting at least one of the  <dig> markers of emt .
this implies that a large proportion of the emt-srams are functionally related to regulators of emt.
using hypergeometric statistics, we created a core emt subnetwork  that accounted for all biological interactions between markers of emt-linking genes , and emt-srams .
while  <dig> % of the genes in the entire curated network were cited in the literature as being related to emt,  <dig> % of the genes in the subnetwork were cited thusly .
this indicates that there was a significant enrichment for known markers of emt in our subnetwork, as well as novel interactions.
taken together, these findings demonstrate that a subset of nsclc srams is related to emt, and that epigenetic regulation through dna methylation may play an important role in the process of emt.
emt-related genes are epigenetically regulated
to determine whether genes associated with emt were epigenetically regulated, we created an emt gene set based on the e-cadherin differential expression  in our panel of nsclc cell lines.
a total of  <dig> genes were ≥ 4-fold differentially expressed between these two groups and thus were included in the emt gene set.
we then identified the publicly available gene expression profiles of nsclc cell lines before and after treatment with the hypomethylating agent 5-azacytidine  from a single study  <cit> .
we performed a gsea using the a <dig> mesenchymal cell line.
we found the emt gene set to be negatively enriched when genes were ranked according to the rho value between gene expression and the average beta value in our integrative analysis .
we found the emt gene set to be positively enriched when genes were ranked according to the log <dig> fold-change before and after treatment with 5aza  .
we found statistical significance  in the overlap between genes that were overexpressed in response to 5aza and genes that were regulated by dna promoter methylation and were represented in the emt gene set enrichment.figure  <dig>
srams induced with 5aza are strongly associated with m cell lines.
a) gsea of srams in a <dig>  finds a substantial proportion of genes to have a negative association with gene sets that are repressed during cell migration and emt states; this relationship is readily reversed by 5aza treatment.
b) data from shames et al.
<cit>  demonstrates more m-srams being induced  than suppressed in m cells, which is not apparent in the e cells.
a similar trend is noticed in the e cells for the e-srams, although not significant.
c) bayesian statistical method to determine the posterior probability that more genes increase  expression after treatment with 5aza by the type of cell line, using either an uninformative prior  or a conservative prior .
since the m cells have proportionally more srams compared to e cells, we next determined whether 5aza was upregulating a larger number of genes in m cells  compared to e cells .
we performed a gsea that included the emt gene set, with genes ranked according to the log <dig> fold-change before and after treatment with 5aza, in all  <dig> nsclc cell lines  <cit> .
while the enrichment of the emt gene set was significant for all but one e cell , the number of genes induced by 5aza with a fold-change ≥  <dig> and the enrichment scores tended to be higher in m cells  compared to e cells  .
we confirmed these results with the dataset from heller et al.
<cit>  where two of the three cell lines  overlapped with that in the shames et al.
dataset.
we performed gsea analysis of these two lines and compared the 5aza and 5aza + trichostatin  treated cells compared to untreated cells.
we found similar enrichment of the emt gene set for the a <dig> cell line but not for the e cell type h <dig> .
to further test the hypothesis that m-related srams are more likely than e-related srams to exhibit increased expression in m cells after 5aza treatment, we used a bayesian statistical method .
we assigned each gene in each cell line to one of three categories , as defined by a change in expression that was at least a 2-fold change in response to 5aza.
we calculated the posterior probability that more genes increase  expression after treatment with 5aza, by the type of cell line, using either an uninformative prior or a conservative prior that assumes most genes do not change expression.
we found that on average, the probability of induced gene expression in the m cells is high for the m-related sram category of genes, and is high for the e cells in the e-related sram category of genes .
these results corroborate the in silico data that indicated that more genes were methylated in the m cells and were inducible upon 5aza treatment.
association between sram methylation and resistance to erlotinib in wild-type egfr/krasnsclc cell lines
since the srams segregate the cell lines into e and m cell types, and emt has been implicated in erlotinib resistance  <cit> , we determined whether the two clusters of nsclc cell lines defined by srams  could be associated with erlotinib resistance.
since egfr mutation is known to be associated with sensitivity to erlotinib, and kras mutations are associated with resistance, we focused on the group of wild-type egfr/kras cell lines.
we found that the half maximal inhibitory concentration  for erlotinib was significantly higher in cell lines that segregated to clusters with methylated srams compared to those that segregated to clusters with unmethylated srams .figure  <dig>
extent of methylated srams associated with erlotinib response in cell lines and patients in battle- <dig> trial.
a) hierarchical clustering of emt-srams in  <dig> wild-type  egfr nsclc cell lines, with erlotinib ic <dig> levels overlying individual cell lines.
g <dig> and g <dig> represent the two cell lines clustered by sram methylation.
kras mutant cells are represented by asterisks.
b) erlotinib ic <dig> values of g <dig> vs g <dig> clusters for the combined cell lines, and the wt egfr/kras cells.
g <dig> cluster had significantly higher ic <dig> values than g <dig> cluster.
error bars indicate average and 95% confidence interval.
c) hierarchical clustering of expression of srams in the panel of cell lines with rppa data for e-cadherin  clusters of e and m cells.
specific gene markers that relate to emt are on the right.
d) hierarchical clustering based on first principle component analysis  of sram expression signature in the tumor biopsies from patients in battle- <dig> trial clusters tumors into groups that resemble the separation seen in the cell lines.
e) eight-week disease control  of patients treated with erlotinib or sorafenib based on the first pca separation.
egfr mutant tumors were excluded.
sram expression profiles in tumors segregate along e vs m expression patterns and erlotinib response
we next determined whether the sram profile was associated with clinical benefit in patients.
since we found no available clinical datasets with dna methylation profiling and erlotinib drug response information, we analyzed the gene expression profile for srams in baseline tumor biopsies collected in patients with stage iv nsclc that were included in the clinical trial biomarker-integrated approaches of targeted therapy for lung cancer elimination   <cit> .
gene expression profiles were available for  <dig> patients treated with erlotinib , and for  <dig> patients treated with sorafenib , all of whom were evaulable for the primary endpoint of the trial, which was 8-week disease control.
interestingly, when srams were used to cluster tumor samples collected at baseline in battle- <dig>  a similar pattern was observed in vitro .
in order to summarize the effect of the srams, we computed the first principal component in each sample.
among the patients treated with erlotinib who had wild-type egfr, the first principal component was lower in the patients who had a clinical status of disease control  .
no difference was observed in patients treated with sorafenib .
when we further analyzed the  <dig> patients treated with erlotinib who had both egfr and kras wild-type tumors, we found that the first principle component trended lower in patients with a clinical status of 8-week disease control .
CONCLUSIONS
in this study, we carried out an integrative analysis of gene expression and promoter cpg methylation profiles in nsclc cell lines to identify a group of genes that were strongly associated with the degree of dna methylation.
we found that a subset of these epigenetically regulated genes reside within pathways related to epithelial and mesenchymal states of cancer.
many of the genes that distinguish e versus m cells in nsclc significantly overlap with the genes that are suppressed during the emt processes found in other cancer systems, such as prostate and breast cancers.
we found that many of the genes were methylated in m cells, and that 5aza upregulated genes that were silenced by dna methylation in these cells compared to the epithelial counterpart.
furthermore, we found that srams could predict erlotinib resistance in the nsclc cell lines, and the gene expression signature of the srams predicts erlotinib response in a cohort of patients with nsclc who were treated with erlotinib.
our data suggest that the srams of various tumors may hold biomarkers that are predictive of response to drug therapy.
transcriptional regulation involves the interplay of a host of epigenetic factors, such as chromatin remodeling, histone modification, non-coding rnas, and dna methylation.
methylation of cpg sites promotes the recruitment of methyl-binding factors, which serves as a scaffold for protein complexes that contribute to nucleosome compaction, chromosomal condensation and transcriptional silencing  <cit> .
while a promoter cpgi is classically linked to gene expression, the definition of what constitutes a cpgi is fairly arbitrary and based on computational prediction.
studies have shown that a large number of cg-rich areas within regions of the promoter and first exon that are not cpgi’s play important roles in transcriptional regulation  <cit> .
we found that for probes corresponding to both cpg and non-cpg sites, a significant excess was seen with a rho cutoff of less than - <dig>  than expected from normal distribution, suggesting the expression of genes with this type of inverse relationship with dna methylation corresponds to the srams.
we also demonstrated from this data that there is an equal likelihood for srams to correspond to promoters with cpgi’s as to those without cpgi’s.
more than 80% of the cpg sites that are related to transcriptional regulation fall within 500 bp upstream or downstream of the tss, which is consistent with what has been reported in the literature.
indeed, there were approximately  <dig> % of the promoter cpg sites to correspond to srams, and  <dig> genes with differential expression based on promoter cpg methylation in the nsclc cell line to correspond to ~ <dig> % of the genome.
this is somewhat similar to what has been reported in the literature regarding the percentage of genes induced after 5aza treatment in multiple cancer types.
however, for one study in nsclc cell lines, the number of srams we found appears to be somewhat lower  <cit> .
this can be explained by the fact that 5aza not only exerts gene induction through hypomethylation of normally hypermethylated promoters, but also has secondary effects, such as the dna damage response  <cit> .
furthermore, our approach identified only genes that are differentially methylated and expressed between cell lines, and not genes that are uniformly suppressed in cancer cells by dna methylation, such as hic1
<cit>  and sox17
<cit> , which may be induced by 5aza.
we attempted to identify genes whose differential methylation patterns could help cluster cells according to functional subgroups; therefore, we believe our approach identifies these differentially regulated genes.
indeed, as proof of the validity of our method, we identified many genes that are known to be differentially methylated across various cancer types , in addition to many novel genes not previously known to be regulated by dna methylation.
we found that nearly 20% of the lung cancer srams are related to genes that distinguish epithelial from mesenchymal cell types, which is similar to what was recently reported  <cit> .
we compared our emt-sram gene list with that from published gene expression datasets derived from various cell systems treated to undergo emt, and found significant overlap in the gene sets.
this was particularly apparent in the downregulated genes, suggesting that a significant proportion of genes that are downregulated during the emt process may be regulated by dna methylation.
using gene interaction network analysis, we found that almost 70% of the genes in the emt-sram panel are one to two gene neighbors away from six of the key factors that play critical roles in the emt process .
a great majority of the methylated genes in the emt-sram panel are seen in mesenchymal cells.
this is consistent with data showing that the emt state induces epigenetic alterations in numerous genes in a non-random fashion.
in breast cancer, basal-like/mesenchymal breast cancers exhibit patterns of methylation in genes such as twist and the estrogen receptor gene that are not found in the luminal type of cells  <cit> .
e-cadherin, or cdh <dig>  a gene involved in cell adhesion and signaling that plays a central role in emt, is suppressed by dna methylation during permanent and irreversible emt.
this is seen for ras-transformed breast cancer cells stimulated with serum.
emt induced in breast cancer cells by treatment with isolated tgf-beta was found to be reversible when the transforming growth factor was removed  <cit> .
however, when non-transformed aml <dig> hepatocytes were induced to undergo emt with tgf-beta, dna methylation was not observed, even though abundant epigenetic reprogramming is known to occur in specific chromatin regions throughout the genome  <cit> .
these data support the proposition that dna methylation of emt-related factors is a non-random event that depends on the transient versus permanent state of emt.
while similar sets of genes may be suppressed during transient emt, this is not due to dna methylation but rather is regulated by chromatin factors.
however, when subject to certain long-term environmental conditions as well as oncogenic factors, these emt-related genes become progressively and irreversibly methylated  <cit> .
cells that acquire these characteristics are locked into an epithelial-like or mesenchymal-like state, which is seen in established cell lines of various origins and also in primary tumors.
tumors that are locked into mesenchymal-like states may be responsible for intrinsic resistance to drug therapies, such as erlotinib.
the mesenchymal-associated srams could unveil important pathways involved in the adaptive regulation of genes directly involved in drug sensitivity, such as the expression of axl in erlotinib resistance.
we found that 5aza treatment induced many of the genes that were silenced by dna methylation, and that there was a preferential induction of genes in the mesenchymal cells, although a great majority of the genes were uninducible, as seen in breast cancer  <cit> .
erlotinib therapy in combination with epigenetic agents may be a promising avenue to help reverse resistance to egfr tyrosine kinase inhibitors in wild-type egfr cells, as shown in a study of cell lines, where egfr silencing by dna methylation contributed to gefitinib resistance that was reversed with decitabine treatment  <cit> .
while current uses of decitabine are mainly reserved for hematologic malignancies, the future use of this drug in solid tumors may require a combination of epigenetic agents, such as histone deacetylase inhibitors with decitabine, or selective small molecular dna methyltransferase inhibitors.

 the methylation of dna is involved in the control of chromatin folding, protein complex assembly, and gene expression. it is a normal process of cellular development, aging, paternal/maternal genetic imprinting, and x-chromosome inactivation. dna methylation is the only known heritable modification of dna, and is aberrant in cancer, involving both excessive and insufficient methylation of dna regions of cytosine–guanine bonds, which then alters the expression of certain genes. the aberrant hypermethylation of dense cpg sites in many cancers occurs at tumor suppressor genes in a non-random, tumor-specific pattern. given the possible role of aberrant dna methylation in cancer initiation and progression, a significant effort has been directed to identify dna methylation biomarkers in cancer and to use such markers to predict therapeutic responses and disease prognosis. a notable example is the association between the silencing of the o6-methylguanine-dna methyltransferase gene and the improved response to temozolomide and radiation in patients with glioblastoma. others include the association of the methylation of pten, igfbp- and various dna repair enzymes with either resistance or sensitivity to chemotherapy or radiation. for lung cancer, dna methylation of genes may be useful for assessing cancer risk based on the analysis of the sputum of smokers, and may be prognostic in early-stage lung cancer. to better understand the dna methylation changes associated with gene expression in lung cancer, we determined the methylation status of more than cpg sites across genes at selective promoter regions in non-small cell lung cancer cell lines using the illumina humanmethylation beadchip. we integrated this dataset with the gene expression profiles of the same panel of cell lines to identify the genes with expression levels that were most affected by dna promoter methylation. we sought to determine the biologic significance of genes that were significantly repressed in association with methylation in silico, as well as experimentally, both in vitro and in vivo.  identifying genes aberrantly associated with dna methylation in nsclc cell lines for each of the nsclc cell lines examined, we analyzed the relationship between the degree of methylation of a particular gene and the expression level of that gene. we sought to identify genes that were significantly repressed in association with dna methylation, which we subsequently refer to as srams, a term first used in the study of breast cancer. we computed the spearman correlation of the relationship between methylation and gene expression for all the cell lines, assuming that genes negatively regulated by dna promoter methylation have a negative spearman rho value. this approach would only identify differentially regulated genes and not genes that are fully methylated or unmethylated in these cell lines. we first tested our approach by looking at the spearman correlation for cdkn2a, a gene well known to be frequently downregulated by promoter methylation in nsclc. in the illumina assay, four probes assessed four cpg sites within the cdkn2a promoter. we found that cdkn2a was differentially methylated in the various cell lines, and the degree of methylation was significantly and inversely correlated with gene expression. this relationship was preserved in each of the four probes examined.figure differential genes are significantly repressed in association with methylation in nsclc cell lines. a) example of the relationship between gene expression and methylation status for a single probe interrogating a single cpg site within the promoter region of cdkn2a. b) frequency histogram of the spearman rho values for all probes that correspond to either cpg or non-cpg island sites. red curves are theoretical gaussian representation of the distribution of rho given the number of probes for each gene set and the actual mean and standard deviation values; blue curves represent the empirical distribution from the data; dashed burgundy lines represent normal distribution curves for comparison; and the arrow identifies the point where the two curves cross and an excess of probes are seen with a rho ≤ -. c) relationship of rho to the distance from the transcriptional start site for cpg and non-cpgi probes. most of the probes that are candidate srams cluster within 500 bp of the tss. we expanded our analysis to the rest of the probes. we looked at the distribution of rho values for probes located within or outside of cpg islands. rho values of cpgi probes had a mean distribution that was significantly shifted to the left of compared to that of non-cpgi probes. for both sets of probes, an excess was seen in the experimental dataset at the rho cutoff of - that was more than expected when compared to a random distribution. this was not observed for probes that had a positive correlation with gene expression. there were also no differences regarding whether the probes corresponded to cpgi or non-cpgi genes. greater than 80% of these probes lie within 500 bp of the transcription starting site, corresponding to sites under transcriptional control by dna methylation. we found probes that correspond to unique srams, genes with expression levels that were inversely correlated with the degree of dna methylation. we found the same pattern as that of cdkn2a for other genes known to be regulated by methylation, as well as for many novel genes not previously known to be regulated by methylation. we validated a set of srams using pyrosequencing and quantitative real-time polymerase chain reaction in cell lines.figure validation of the integrative analysis for select genes. three genes with cpg islands were evaluated. top row is the relationship between array methylation data and illumina gene expression; middle row provides pearson correlation of the array methylation data with the pyrosequencing methylation status; and bottom row shows the relationship of methylation status by pyrosequencing with the quantitative real-time pcr results for the expression of these genes. srams differentiate tumor from normal lung and distinguish tumor subtypes as aberrant dna methylation is one of the defining characteristics of cancer, srams should reflect alterations that discriminate cancer from normal tissue and reflect the underlying tumor biology. we evaluated the methylation status of srams in a set of lung tumors from a japanese cohort, hereafter termed the kyoto tumor set, which included unmatched normal lung tissues arrayed using the illumina humanmethylation beadchip. since the cell line sram dataset was generated from the humanmethylation platform, we took probes that overlapped between the two illumina platforms for the analysis. we found that srams discriminated nsclc from normal lung tissue, with clusters of probes that were either predominantly methylated in normal tissue or methylated in tumors. srams also distinguished squamous carcinoma from adenocarcinoma, the two major histological subtypes of nsclc, as well as lung adenocarcinomas with or without an egfr mutation. these observations were corroborated using sram expression profiles of lung adenocarcinomas included in the director’s challenge consortium for the molecular classification of lung adenocarcinoma. srams preferentially separated egfr mutants from wild-type tumors. in the kyoto tumor set, patients that separated into these clusters had significant survival differences, similar to that seen for the director’s challenge tumor set.figure nsclc srams are differentially methylated between normal and tumor tissues segregate into molecular pathologic subtypes. a) heatmap of the hierarchical clustering of matched probes for srams in the kyoto set of tumors and compared to unmatched normal lung tissue samples. b) using matched probes from the humanmethylation to the humanmethylation platforms, unsupervised clustering of tumors from the kyoto tumor set was conducted, with tumor clustering along lines of histology, egfr mutation status, and tumor grades. c-d) unsupervised clustering based on gene expression of the matched gene regulatory module genes on the director’s challenge dataset of lung adenocarcinomas. the grms were able to segregate the tumors in the director’s challenge into two main groups, with egfr mutations and well-to-moderately differentiated tumors in group and wild-type egfr and moderately-to-poorly differentiated tumors in group e) sram segregates nsclc tumors into prognostic groups in the kyoto and director’s challenge tumor sets. branch point c vs c mostly segregates squamous from adenocarcinoma, without a difference in prognosis; whereas branch point c vs. c is highly prognostic because of the separation of egfr mutant from wild-type tumors. for the director’s challenge tumor set, the two clusters that preferentially separate egfr mutant and wild-type tumors are also prognostic for survival. srams and the epithelial-to-mesenchymal transition phenotype to further define the functional significance of genes negatively regulated by dna promoter methylation, we performed a gene set enrichment analysis using the spearman correlation rho value as the ranking variable. the top gene sets that included genes negatively regulated by dna promoter methylation had been obtained from diverse cancer datasets and were involved in various biologic functions, such as cell migration in bladder cancer and resistance to gefitinib in nsclc, and included genes downregulated with e-cadherin knockdown in human breast mammary epithelial cell lines, genes differentially expressed in metastatic melanoma, genes methylated in glioblastoma and pancreatic cancer, and genes differentially regulated in luminal vs. basal/mesenchymal breast cancer cells. we postulated that biologic functions identified by the gsea were related to the process of emt, a developmental and adaptive cellular process that has been associated with resistance to cancer therapies and regulation of metastasis. we performed a hierarchical cluster analysis of the cell lines using the srams, and overlaid the relative protein expression of e-cadherin, a cell adhesion molecule that is downregulated during emt and which plays a key role in the signaling and regulation of emt, along with the expression of emt-related genes zeb vim, twist fn cdh and cdh we found that the identification of srams clustered cells into two groups: “epithelial-like ” cells, which expressed high e-cadherin levels; and “mesenchymal-like ” cells, which expressed low e-cadherin levels. to identify genes within the srams that distinguish e cells from m cells, we performed a wilcoxon rank sum test. out of the unique genes, were differentially methylated between e and m cells, with a p-value < . we called these genes emt-srams. a majority of emt-srams were methylated in m cells compared to e cells. a few of the emt-srams, including cdh were already known to be epigenetically regulated during the process of emt. however, the epigenetic regulation of the remaining emt-srams had not been established, including that of spint1, syk, mst1r and esrp1. interestingly, one of these genes is axl, which encodes a receptor tyrosine kinase that is associated with emt and with resistance to egfr tyrosine kinase inhibition. we found that the expression of axl was inversely correlated with promoter methylation, and its promoter methylation was positively correlated with e-cadherin expression, indicating that axl may be epigenetically regulated during emt.figure specific sram methylation patterns define epithelial and mesenchymal cell types. a) hierarchical clustering of srams for a subset of cell lines with rppa data for e-cadherin segregates cell lines into high and low ecad groups. representative genes that are preferentially unmethylated or methylated in the mesenchymal cells are shown. illumina gene expression levels of emt markers, twist fibronectin, cdh and cdh1) are overlaid to confirm the e vs m identity of the cell lines. b) wilcoxon rank sum test was used to identify the srams that are differentially methylated between e and m cells, which we call the emt-srams. the identified probes were used for hierarchical clustering of the cell lines. emt-srams that are predominantly methylated in m cells are called m-srams, and those methylated in e cells are called e-srams. c) emt-sram association with emt network genes using curated network analysis. closest second neighbor analysis of genes associated with core emt factors finds numerous associations with genes within the emt-sram set. red spheres represent srams that are methylated in m cells; green spheres are those methylated in e cells; white spheres represent linking genes. to further determine whether emt-srams define an epigenetically regulated group of genes relevant to emt, we compared our gene list to a variety of publically available gene expression signatures from studies that experimentally manipulated cells to undergo emt. specifically, we examined datasets from prostate cancer cells transduced with snai1, a nsclc epithelial cell line transduced with snai zeb or stimulated with tgf-beta, and hbme cells transduced with twist. we applied hypergeometric statistics to test the significance of overlap between the entire dataset of srams and the emt-srams. we found that the srams significantly overlapped with the genes in the published datasets. on average, there was a 16% overlap between the srams and published emt gene signatures. emt-srams had greater overlap with the published emt gene signatures, ranging from 18% in the breast cancer cell line signature to 48% in the h cells treated with tgf-beta. the majority of overlapping genes were downregulated, which is consistent with our observation of a greater proportion of genes being methylated in m cell lines. there appeared to be a greater overlap between genes from h cells and the srams or the emt-srams than either the prostate or breast cancer datasets. this suggests that while there could be a common core group of genes for both the srams and the emt-associated srams, the tissue-specific srams may reflect the underlying biology of the tissue of origin. using a curated network extracted from pathway commons, which consists of genes and over a million biological interactions, we performed a network analysis using markers of emt, including transcription factors and major regulators of emt, and three recognized markers of emt, and the emt-sram gene list. we were able to map 73/ emt-srams. the pair-wise shortest distance between markers of emt with emt-srams was calculated using dijkastra’s algorithm: 51/ emt-srams in the network had the shortest path lengths of and connecting at least one of the markers of emt. this implies that a large proportion of the emt-srams are functionally related to regulators of emt. using hypergeometric statistics, we created a core emt subnetwork that accounted for all biological interactions between markers of emt-linking genes, and emt-srams. while % of the genes in the entire curated network were cited in the literature as being related to emt, % of the genes in the subnetwork were cited thusly. this indicates that there was a significant enrichment for known markers of emt in our subnetwork, as well as novel interactions. taken together, these findings demonstrate that a subset of nsclc srams is related to emt, and that epigenetic regulation through dna methylation may play an important role in the process of emt. emt-related genes are epigenetically regulated to determine whether genes associated with emt were epigenetically regulated, we created an emt gene set based on the e-cadherin differential expression in our panel of nsclc cell lines. a total of genes were ≥ 4-fold differentially expressed between these two groups and thus were included in the emt gene set. we then identified the publicly available gene expression profiles of nsclc cell lines before and after treatment with the hypomethylating agent 5-azacytidine from a single study. we performed a gsea using the a mesenchymal cell line. we found the emt gene set to be negatively enriched when genes were ranked according to the rho value between gene expression and the average beta value in our integrative analysis. we found the emt gene set to be positively enriched when genes were ranked according to the log fold-change before and after treatment with 5aza. we found statistical significance in the overlap between genes that were overexpressed in response to 5aza and genes that were regulated by dna promoter methylation and were represented in the emt gene set enrichment.figure srams induced with 5aza are strongly associated with m cell lines. a) gsea of srams in a finds a substantial proportion of genes to have a negative association with gene sets that are repressed during cell migration and emt states; this relationship is readily reversed by 5aza treatment. b) data from shames et al. demonstrates more m-srams being induced than suppressed in m cells, which is not apparent in the e cells. a similar trend is noticed in the e cells for the e-srams, although not significant. c) bayesian statistical method to determine the posterior probability that more genes increase expression after treatment with 5aza by the type of cell line, using either an uninformative prior or a conservative prior. since the m cells have proportionally more srams compared to e cells, we next determined whether 5aza was upregulating a larger number of genes in m cells compared to e cells. we performed a gsea that included the emt gene set, with genes ranked according to the log fold-change before and after treatment with 5aza, in all nsclc cell lines. while the enrichment of the emt gene set was significant for all but one e cell, the number of genes induced by 5aza with a fold-change ≥  and the enrichment scores tended to be higher in m cells compared to e cells. we confirmed these results with the dataset from heller et al. where two of the three cell lines overlapped with that in the shames et al. dataset. we performed gsea analysis of these two lines and compared the 5aza and 5aza + trichostatin treated cells compared to untreated cells. we found similar enrichment of the emt gene set for the a cell line but not for the e cell type h. to further test the hypothesis that m-related srams are more likely than e-related srams to exhibit increased expression in m cells after 5aza treatment, we used a bayesian statistical method. we assigned each gene in each cell line to one of three categories, as defined by a change in expression that was at least a 2-fold change in response to 5aza. we calculated the posterior probability that more genes increase expression after treatment with 5aza, by the type of cell line, using either an uninformative prior or a conservative prior that assumes most genes do not change expression. we found that on average, the probability of induced gene expression in the m cells is high for the m-related sram category of genes, and is high for the e cells in the e-related sram category of genes. these results corroborate the in silico data that indicated that more genes were methylated in the m cells and were inducible upon 5aza treatment. association between sram methylation and resistance to erlotinib in wild-type egfr/krasnsclc cell lines since the srams segregate the cell lines into e and m cell types, and emt has been implicated in erlotinib resistance, we determined whether the two clusters of nsclc cell lines defined by srams could be associated with erlotinib resistance. since egfr mutation is known to be associated with sensitivity to erlotinib, and kras mutations are associated with resistance, we focused on the group of wild-type egfr/kras cell lines. we found that the half maximal inhibitory concentration for erlotinib was significantly higher in cell lines that segregated to clusters with methylated srams compared to those that segregated to clusters with unmethylated srams.figure extent of methylated srams associated with erlotinib response in cell lines and patients in battle- trial. a) hierarchical clustering of emt-srams in wild-type egfr nsclc cell lines, with erlotinib ic levels overlying individual cell lines. g and g represent the two cell lines clustered by sram methylation. kras mutant cells are represented by asterisks. b) erlotinib ic values of g vs g clusters for the combined cell lines, and the wt egfr/kras cells. g cluster had significantly higher ic values than g cluster. error bars indicate average and 95% confidence interval. c) hierarchical clustering of expression of srams in the panel of cell lines with rppa data for e-cadherin clusters of e and m cells. specific gene markers that relate to emt are on the right. d) hierarchical clustering based on first principle component analysis of sram expression signature in the tumor biopsies from patients in battle- trial clusters tumors into groups that resemble the separation seen in the cell lines. e) eight-week disease control of patients treated with erlotinib or sorafenib based on the first pca separation. egfr mutant tumors were excluded. sram expression profiles in tumors segregate along e vs m expression patterns and erlotinib response we next determined whether the sram profile was associated with clinical benefit in patients. since we found no available clinical datasets with dna methylation profiling and erlotinib drug response information, we analyzed the gene expression profile for srams in baseline tumor biopsies collected in patients with stage iv nsclc that were included in the clinical trial biomarker-integrated approaches of targeted therapy for lung cancer elimination. gene expression profiles were available for patients treated with erlotinib, and for patients treated with sorafenib, all of whom were evaulable for the primary endpoint of the trial, which was 8-week disease control. interestingly, when srams were used to cluster tumor samples collected at baseline in battle- a similar pattern was observed in vitro. in order to summarize the effect of the srams, we computed the first principal component in each sample. among the patients treated with erlotinib who had wild-type egfr, the first principal component was lower in the patients who had a clinical status of disease control. no difference was observed in patients treated with sorafenib. when we further analyzed the patients treated with erlotinib who had both egfr and kras wild-type tumors, we found that the first principle component trended lower in patients with a clinical status of 8-week disease control.  in this study, we carried out an integrative analysis of gene expression and promoter cpg methylation profiles in nsclc cell lines to identify a group of genes that were strongly associated with the degree of dna methylation. we found that a subset of these epigenetically regulated genes reside within pathways related to epithelial and mesenchymal states of cancer. many of the genes that distinguish e versus m cells in nsclc significantly overlap with the genes that are suppressed during the emt processes found in other cancer systems, such as prostate and breast cancers. we found that many of the genes were methylated in m cells, and that 5aza upregulated genes that were silenced by dna methylation in these cells compared to the epithelial counterpart. furthermore, we found that srams could predict erlotinib resistance in the nsclc cell lines, and the gene expression signature of the srams predicts erlotinib response in a cohort of patients with nsclc who were treated with erlotinib. our data suggest that the srams of various tumors may hold biomarkers that are predictive of response to drug therapy. transcriptional regulation involves the interplay of a host of epigenetic factors, such as chromatin remodeling, histone modification, non-coding rnas, and dna methylation. methylation of cpg sites promotes the recruitment of methyl-binding factors, which serves as a scaffold for protein complexes that contribute to nucleosome compaction, chromosomal condensation and transcriptional silencing. while a promoter cpgi is classically linked to gene expression, the definition of what constitutes a cpgi is fairly arbitrary and based on computational prediction. studies have shown that a large number of cg-rich areas within regions of the promoter and first exon that are not cpgi’s play important roles in transcriptional regulation. we found that for probes corresponding to both cpg and non-cpg sites, a significant excess was seen with a rho cutoff of less than - than expected from normal distribution, suggesting the expression of genes with this type of inverse relationship with dna methylation corresponds to the srams. we also demonstrated from this data that there is an equal likelihood for srams to correspond to promoters with cpgi’s as to those without cpgi’s. more than 80% of the cpg sites that are related to transcriptional regulation fall within 500 bp upstream or downstream of the tss, which is consistent with what has been reported in the literature. indeed, there were approximately % of the promoter cpg sites to correspond to srams, and genes with differential expression based on promoter cpg methylation in the nsclc cell line to correspond to ~ % of the genome. this is somewhat similar to what has been reported in the literature regarding the percentage of genes induced after 5aza treatment in multiple cancer types. however, for one study in nsclc cell lines, the number of srams we found appears to be somewhat lower. this can be explained by the fact that 5aza not only exerts gene induction through hypomethylation of normally hypermethylated promoters, but also has secondary effects, such as the dna damage response. furthermore, our approach identified only genes that are differentially methylated and expressed between cell lines, and not genes that are uniformly suppressed in cancer cells by dna methylation, such as hic1 and sox17, which may be induced by 5aza. we attempted to identify genes whose differential methylation patterns could help cluster cells according to functional subgroups; therefore, we believe our approach identifies these differentially regulated genes. indeed, as proof of the validity of our method, we identified many genes that are known to be differentially methylated across various cancer types, in addition to many novel genes not previously known to be regulated by dna methylation. we found that nearly 20% of the lung cancer srams are related to genes that distinguish epithelial from mesenchymal cell types, which is similar to what was recently reported. we compared our emt-sram gene list with that from published gene expression datasets derived from various cell systems treated to undergo emt, and found significant overlap in the gene sets. this was particularly apparent in the downregulated genes, suggesting that a significant proportion of genes that are downregulated during the emt process may be regulated by dna methylation. using gene interaction network analysis, we found that almost 70% of the genes in the emt-sram panel are one to two gene neighbors away from six of the key factors that play critical roles in the emt process. a great majority of the methylated genes in the emt-sram panel are seen in mesenchymal cells. this is consistent with data showing that the emt state induces epigenetic alterations in numerous genes in a non-random fashion. in breast cancer, basal-like/mesenchymal breast cancers exhibit patterns of methylation in genes such as twist and the estrogen receptor gene that are not found in the luminal type of cells. e-cadherin, or cdh a gene involved in cell adhesion and signaling that plays a central role in emt, is suppressed by dna methylation during permanent and irreversible emt. this is seen for ras-transformed breast cancer cells stimulated with serum. emt induced in breast cancer cells by treatment with isolated tgf-beta was found to be reversible when the transforming growth factor was removed. however, when non-transformed aml hepatocytes were induced to undergo emt with tgf-beta, dna methylation was not observed, even though abundant epigenetic reprogramming is known to occur in specific chromatin regions throughout the genome. these data support the proposition that dna methylation of emt-related factors is a non-random event that depends on the transient versus permanent state of emt. while similar sets of genes may be suppressed during transient emt, this is not due to dna methylation but rather is regulated by chromatin factors. however, when subject to certain long-term environmental conditions as well as oncogenic factors, these emt-related genes become progressively and irreversibly methylated. cells that acquire these characteristics are locked into an epithelial-like or mesenchymal-like state, which is seen in established cell lines of various origins and also in primary tumors. tumors that are locked into mesenchymal-like states may be responsible for intrinsic resistance to drug therapies, such as erlotinib. the mesenchymal-associated srams could unveil important pathways involved in the adaptive regulation of genes directly involved in drug sensitivity, such as the expression of axl in erlotinib resistance. we found that 5aza treatment induced many of the genes that were silenced by dna methylation, and that there was a preferential induction of genes in the mesenchymal cells, although a great majority of the genes were uninducible, as seen in breast cancer. erlotinib therapy in combination with epigenetic agents may be a promising avenue to help reverse resistance to egfr tyrosine kinase inhibitors in wild-type egfr cells, as shown in a study of cell lines, where egfr silencing by dna methylation contributed to gefitinib resistance that was reversed with decitabine treatment. while current uses of decitabine are mainly reserved for hematologic malignancies, the future use of this drug in solid tumors may require a combination of epigenetic agents, such as histone deacetylase inhibitors with decitabine, or selective small molecular dna methyltransferase inhibitors. 
methylation of the candidate genes was associated with erlotinib resistance in wild-type egfr cell lines. the expression profiles of the candidate genes were associated with 8-week disease control in patients with wild-type egfr who had unresectable non-small cell lung cancer treated with erlotinib, but not in patients treated with sorafenib. we integrated the gene expression profiles and data of gene promoter methylation for a large panel of non-small cell lung cancer cell lines, and identified candidate genes with expression levels that were inversely correlated to the degree of dna methylation. we used gene set enrichment analysis of the genes ranked by the degree of correlation between gene expression and dna methylation to identify gene sets involved in cellular migration and metastasis. our unsupervised hierarchical clustering of the candidate genes segregated cell lines according to the epithelial-to-mesenchymal transition phenotype. our results demonstrate that the underlying biology of genes regulated by dna methylation may have predictive value in lung cancer that can be exploited therapeutically. dna methylation is associated with aberrant gene expression in cancer, and has been shown to correlate with therapeutic response and disease prognosis in some types of cancer. we sought to investigate the biological significance of dna methylation in lung cancer. dna methylationepithelial-mesenchymal transitionerlotiniblung cancerissue-copyright-statement© the author 2014 genes related to the epithelial-to-mesenchymal transition, such as axl, esrp hoxb and spint1/ were among the nearly 20% of the candidate genes that were differentially methylated between epithelial and mesenchymal cells. greater numbers of genes were methylated in the mesenchymal cells and their expressions were upregulated by 5-azacytidine treatment. we found these candidate genes to be differentially methylated in normal lung tissue versus non-small cell lung cancer tumors, and segregated by histologic and tumor subtypes.
 dna methylation is associated with aberrant gene expression in cancer, and has been shown to correlate with therapeutic response and disease prognosis in some types of cancer. we sought to investigate the biological significance of dna methylation in lung cancer.  we integrated the gene expression profiles and data of gene promoter methylation for a large panel of non-small cell lung cancer cell lines, and identified candidate genes with expression levels that were inversely correlated to the degree of dna methylation. we found these candidate genes to be differentially methylated in normal lung tissue versus non-small cell lung cancer tumors, and segregated by histologic and tumor subtypes. we used gene set enrichment analysis of the genes ranked by the degree of correlation between gene expression and dna methylation to identify gene sets involved in cellular migration and metastasis. our unsupervised hierarchical clustering of the candidate genes segregated cell lines according to the epithelial-to-mesenchymal transition phenotype. genes related to the epithelial-to-mesenchymal transition, such as axl, esrp hoxb and spint1/ were among the nearly 20% of the candidate genes that were differentially methylated between epithelial and mesenchymal cells. greater numbers of genes were methylated in the mesenchymal cells and their expressions were upregulated by 5-azacytidine treatment. methylation of the candidate genes was associated with erlotinib resistance in wild-type egfr cell lines. the expression profiles of the candidate genes were associated with 8-week disease control in patients with wild-type egfr who had unresectable non-small cell lung cancer treated with erlotinib, but not in patients treated with sorafenib.  our results demonstrate that the underlying biology of genes regulated by dna methylation may have predictive value in lung cancer that can be exploited therapeutically. electronic supplementary material the online version of this article contains supplementary material, which is available to authorized users. keywords dna methylationepithelial-mesenchymal transitionerlotiniblung cancerissue-copyright-statement© the author 2014 