Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors

ARTICLE IN PRESS Cancer Letters ■■ (2016) ■■–■■

Contents lists available at ScienceDirect

Cancer Letters j o u r n a l h o m e p a g e : w w w. e l s e v i e r. c o m / l o c a t e / c a n l e t

1 2 3 4 5 6

Q2 Original Articles

Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors Q1 Chiara E. Cotroneo a,1,2, Alice Dassano a,1, Francesca Colombo a, Angela Pettinicchio a,

Daniele Lecis b, Matteo Dugo b, Loris De Cecco b, Tommaso A. Dragani a,*, Giacomo Manenti a

7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43

a b

Department of Predictive and Preventive Medicine, Fondazione IRCCS, Istituto Nazionale dei Tumori, Milan, Italy Department of Experimental Oncology and Molecular Medicine, Fondazione IRCCS, Istituto Nazionale dei Tumori, Milan, Italy

A R T I C L E

I N F O

Article history: Received 2 December 2015 Received in revised form 25 February 2016 Accepted 26 February 2016

Q3

Keywords: Animal models Expression quantitative trait loci Gene expression Lung cancer Single nucleotide polymorphisms Transcriptome Pas1

A B S T R A C T

Gene expression modulates cellular functions in both physiologic and pathologic conditions. Herein, we carried out a genetic linkage study on the transcriptome of lung tumors induced by urethane in an (A/J x C57BL/6)F4 intercross population, whose individual lung tumor multiplicity (Nlung) is linked to the genotype at the Pulmonary adenoma susceptibility 1 (Pas1) locus. We found that expression levels of 1179 and 1579 genes are modulated by an expression quantitative trait locus (eQTL) in cis and in trans, respectively (LOD score > 5). Of note, the genomic area surrounding and including the Pas1 locus regulated 14 genes in cis and 857 genes in trans. In lung tumors of the same (A/J x C57BL/6)F4 mice, we found 1124 genes whose transcript levels associated with Nlung (FDR < 0.001). The expression levels of about a third of these genes (n = 401) were regulated by the genotype at the Pas1 locus. Pathway analysis of the sets of genes associated with Nlung and regulated by Pas1 revealed a set of 14 recurrently represented genes that are components or targets of the Ras–Erk and Pi3k–Akt signaling pathways. Altogether our results illustrate the architecture of germline control of gene expression in mouse lung cancer: they highlight the importance of Pas1 as a tumor-modifier locus, attribute to it a novel role as a major regulator of transcription in lung tumor nodules and strengthen the candidacy of the Kras gene as the effector of this locus. © 2016 Published by Elsevier Ireland Ltd.

44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59

Introduction Expression quantitative trait loci (eQTLs) are genomic regions containing genetic variations that influence the transcription of a gene [1,2]. eQTLs can act either in cis, when the genetic variant specifically affects the expression of the allele of a gene located on the same chromosome, or in trans, when the variation alters the expression of both alleles [3]. Cis-acting eQTLs usually contain genetic variants that affect the molecular properties of DNA regulatory elements (e.g. promoters, UTRs, splicing sites, enhancers) located close to, or within, their target genes. Trans-acting loci, usually located in separate areas of the genome from their target genes, contain DNA variants which alter the functional properties or the abundance of diffusible regulatory elements such as transcription factors and microRNAs. However, in genome-wide eQTL studies, the defini-

60 61 62 63 64 65 66

* Corresponding author. Tel.: +39 0223902642; fax: +39 0223902764. E-mail address: [email protected] (T.A. Dragani). 1 These authors contributed equally to the work. 2 Present address: UCD School of Biomolecular and Biomedical Science, University College Dublin, Belfield, Dublin 4, Ireland.

tion of cis- (or local) and trans- (or distal) eQTLs is usually simplified by considering only the physical distance between the variant and the regulated gene(s). eQTL studies can help uncover the functional mechanisms through which disease-associated variants, previously identified in genome-wide association or linkage studies, exert their role in pathogenesis. Indeed, most of these variants map in non-coding regions and therefore are likely to act as regulatory elements [4]. So far, just one genome-wide eQTL study of human lung cancer tissue has been reported, but it focused on genetic regulation acting in cis [5]. Indeed, in humans, the wide genetic heterogeneity among individuals [6] limits the statistical power for the genome-wide detection of eQTLs, hindering the detection of trans-eQTLs which typically have smaller phenotypic effects [3]. The difficulties in detecting eQTLs can be overcome through the use of simpler genetic models, such as inbred mouse strains that are homozygous at almost every genetic locus and that can be crossed to obtain progeny with lower genetic complexity than pedigrees of wild-type animals. In mice, numerous genetic loci are known to modify the development and progression of several cancer types, including lung cancer [7,8]. Except for a few cases where the identification of functional variations in candidate genes suggested a

http://dx.doi.org/10.1016/j.canlet.2016.02.054 0304-3835/© 2016 Published by Elsevier Ireland Ltd.

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88

ARTICLE IN PRESS C.E. Cotroneo et al./Cancer Letters ■■ (2016) ■■–■■

2

89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111

molecular mechanism of action [8–12], the effectors of most of these loci are still to be understood. Regarding the Pulmonary adenoma susceptibility 1 (Pas1) locus, the major locus affecting spontaneous and chemically induced lung tumorigenesis in mouse inbred strains [7], we recently suggested that it exerts its function through the modulation in cis of the expression of the 4A isoform of the cancer-related gene Kras (Kras 4A) in normal lung and in lung tumors [13]. This work was done in (A/J x C57BL/6)F4 intercross mice (hereafter called ABF4 mice) derived from crossing the A/J and the C57BL/6 inbred strains, which carry the susceptibility and resistance alleles of the Pas1 locus, respectively. Here, in the same ABF4 cross, we explored the genome-wide architecture of genetic regulation of gene expression in mouse lung tumor nodules, and identified the cis- and trans-eQTLs acting in this tissue. We observed that the Pas1 locus exerts a pleiotropic effect on transcription in lung tumors, with a large number of transregulated genes not comparable with that of any other locus. By intersecting the set of genes whose levels in lung cancer tissue associated with lung tumor multiplicity with those whose levels were controlled by the genotype at the Pas1 locus, we identified a set of 14 key genes of Pas1 biological function in lung cancer susceptibility.

112

Ethics statement

113 114 115 116 117

All animals received humane care according to the criteria outlined in a protocol approved on December 21, 2006, by the institutional ethical committee for animal use (CESA) at the Fondazione IRCCS Istituto Nazionale dei Tumori.

118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133

In a previous paper [13], we described the generation, treatment and sampling of intercross mice. Briefly, ABF4 mice were treated with a single intraperitoneal injection of urethane (1 g/kg body weight) according to a standard procedure for inducing the development of lung tumors [14]. When the animals were 40 weeks of age, they were killed, the lungs were removed for tumor counting (values of lung tumor multiplicity (Nlung) varied substantially among individuals, ranging from one to 33 in the animals that developed tumors), and a single lung tumor, 1–1.5 mm in diameter, was excised from each for RNA extraction. By histologic analysis of some excised tumors, they were classified as lung adenomas, i.e., the typical mouse lung tumor histotype [15]. Here, we took advantage of existing genotype data, for 142 male mice, regarding 548 informative (polymorphic) non-redundant SNPs dispersed over the whole genome excluding the Y chromosome (average density of coverage, ~ 6.0 Mb/SNP). We also used Nlung data and samples of total RNA from lung tumors of these 142 mice.

Materials and methods

Mouse crosses, biological samples and genotype data

Exon array analysis

Total RNA was amplified using the Low Input Quick Amp WT Labeling kit, ac134 cording to the manufacturer’s instructions (Agilent Technologies, Santa Clara, USA). 135 Fluorescent dye-labeled cRNA was hybridized to SurePrint G3 Mouse Exon 4 × 180K 136 arrays (Agilent); hybridization and washing were performed on Agilent’s microarray 137 platform, according to standard protocols. Microarray images were acquired using 138 an Agilent DNA microarray scanner; raw data were generated using Agilent Feature 139 Extraction software. All microarray data were MIAME compliant and were depos140 141 Q4 ited into the NCBI’s GEO database (http://www.ncbi.nmlm.nih.gov/projects/geo/) with accession number GSE71232. 142 Microarray data were analyzed with a set of customized scripts in R 143 (http://www.r-project.org/). First, raw median probe intensities were pre-processed 144 145 using the limma package [16] of the R/Bioconductor project [17]; then, expression 146 data were summarized at the gene level calculating the geometric mean of probes 147 associated with the same gene. Only probes whose intensity was at least 10% brighter 148 than the 95th percentile of negative control probes in at least 15% of samples were 149 used for gene-level summarization. Finally, we removed genes not listed in the current 150 mouse genome assembly (Mouse GRCm38, according to the Ensembl database [18]) 151 as well as genes now listed as duplicates. After this filtering, performed using the 152 biomaRt R package [19], we obtained a set of expression data for 18,020 genes for 153 each individual of our population. 154 155 eQTL analysis

156 157 158

To identify QTLs controlling gene expression, we did simple interval mapping to assess genetic linkages for all the examined 18,020 genes, using a set of custom R scripts based on the R package qtl [20]. For each examined gene, we obtained a

LOD score curve covering the whole genome; peaks in these curves corresponded to the position of putative eQTLs for that gene on the genetic map (expressed in cM). The genomic coordinates of each peak (in basepairs) were taken as those of the genotyped marker having the highest LOD score among all the markers mapping in proximity of the LOD peak. However, since in this mouse cross each marker identified a genomic region (locus) with an average length of 6 Mb, each of the 548 characterized markers was in linkage disequilibrium with several dozens or hundreds of genetic variants. This means that the observed modulation of gene expression can be due not only to the marker itself, but also to one or more variants in linkage disequilibrium with it. We defined a cis-acting eQTL as a locus that modulates the expression of a gene whose annotated start or end point is ≤ 3 Mb (i.e., ≤ 50% of the average distance between the genotyped markers) from the peak of the LOD score curve. A trans-acting eQTL was initially defined as any statistically significant LOD score peak located further than 3 Mb from the target gene on the same chromosome, or located on another chromosome. For some genes, however, we found multiple peaks, with decreasing LOD scores, on the same chromosome. In most cases, these peaks were partially overlapping and centered about very close positions on the genetic map, suggesting that their presence was caused by only one regulatory locus. However, the low number of meiotic recombinations in the F4 intercross did not allow us to determine if this behavior of the LOD score curves was caused by a single locus or by multiple, very close, eQTLs. Therefore, for every chromosome displaying multiple, statistically significant LOD score peaks, we only considered the one with the highest LOD score, i.e. the peak with the strongest linkage with the examined phenotype. To determine which peaks of a LOD score curve indicated significant linkage, we calculated threshold LOD values by permutation analysis on 62 genes randomly selected among all the genes that had LOD score peaks higher than the genomewide LOD threshold of 3.3 proposed in the literature for mouse crosses [21]. In detail, 31 out of the 62 genes had a LOD peak corresponding to a putative cis-acting eQTL and 31 had a LOD peak corresponding to a putative trans-acting eQTL. After permutation analysis, these 62 genes showed threshold LOD scores peaks (α = 0.05) ranging from 3.8 to 5.7, with a median of 4.0. Therefore, we decided to use a threshold LOD score of 5.0 for all eQTLs; such a threshold is conservative over the proposed genome-wide LOD threshold of 3.3. A circular genome map to visualize inter-chromosomal trans-eQTLs was generated with R scripts for data formatting and Circos software [22].

159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195

Statistical analysis of the transcriptome

196

Statistical association of transcript levels with lung tumor multiplicity (squareroot transformed Nlung values) was tested with custom R scripts by simple linear regression modeling, examining the association of one gene at a time, with Nlung and the expression level of the examined gene treated as continuous variables. The obtained P-values were adjusted for multiple testing following the Benjamini– Hochberg procedure [23]; the threshold for genome-wide significance was placed at an FDR < 0.001.

197 198 199 200 201 202 203

Pathway analysis

204

For genes whose expression level was linked to the genotype in the genomic area surrounding the Pas1 locus, we searched for their participation in annotated biological pathways, looking for cases of functional enrichment, using the Ingenuity Pathway Analysis (IPA, Qiagen; http://www.ingenuity.com/products/ipa) online tool and the Ingenuity Canonical Pathways database. The significance threshold was placed at an FDR < 0.05. The Pas1-regulated genes involved in all the pathways identified by IPA were visually represented as a network using R scripts for data formatting and Cytoscape software [24]. The network was generated by connecting genes (nodes) with a line (edge) when they were involved in at least one common pathway. The size of nodes was scaled proportionally to the total number of IPA-identified pathways in which each gene was involved (Npathways). Nodes were colored in accordance with the direction of association of their expression levels with the number of lung tumor susceptibility (A/J-derived) alleles at the Pas1 locus (red, genes positively associated with the A/J-derived allele; green, genes inversely associated with the A/Jderived allele).

205 206 207 208 209 210 211 212 213 214 215 216 217 218 219

Western blot analysis

220

Tissues from normal lungs and nodules, collected from mice and stored at −80 °C, were finely minced and lysed by boiling 10 minutes in lysis buffer (125 mM Tris HCl pH 6.8, 5% sodium dodecyl sulfate/SDS). After the addition of protease and phosphatase inhibitors, samples were sonicated and centrifuged at 11,000 rpm for 15 minutes at RT. Cleared supernatants were separated by SDS-PAGE using precast 4–12% Bis-Tris NuPAGE gels (Thermo Fisher Scientific) and blotted onto PVDF membranes (Merck Millipore) using the XCell II blot module (Thermo Fisher Scientific). Membranes were saturated in blocking buffer containing Tris-buffered saline (TBS) with 4% BSA for 30 minutes and then incubated overnight with the primary antibodies purchased from Cell Signaling Technology (Phospho-Akt Ser473 #9271; Akt #2920; Phospho-MEK1/2 #9121; MEK1/2 #4694) and Sigma-Aldrich (Phospho-ERK1/2 #M8159; ERK1/2 #M5670; Vinculin #V9131). After 1 h incubation with the appro-

221 222 223 224 225 226 227 228 229 230 231 232

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

ARTICLE IN PRESS C.E. Cotroneo et al./Cancer Letters ■■ (2016) ■■–■■

233 234 235 236 237 238 239 240 241

3

Fig. 1. Four cis-eQTLs with the highest LOD score values. For each of the four cis-regulated genes (Hddc3, Glo1, Btbd9 and BC022687), top panels show the corresponding genome-wide LOD score plots. LOD score peaks indicate the position of each eQTL. Box-plots in the bottom panels show log2-transformed gene expression levels (Y axis) among mice grouped according to genotype at the eQTL marker (X axis). For each marker, “AA” indicates mice homozygous for the A/J-derived allele; “AB” indicates heterozygous mice; “BB” indicates mice homozygous for the C57BL/6-derived allele. For the Hddc3 gene, the number of mice for each genotype at the rs13479363 marker is: AA = 22, AB = 55, BB = 64 (genotype data for one mouse was not available); for the Glo1 and Btbd9 genes, the number of mice for each genotype at the rs3693494 marker is: AA = 24, AB = 73, BB = 44 (genotype data for one mouse was not available); for the BC022687 gene, the number of mice for each genotype at the rs13481655 marker is: AA = 22, AB = 72, BB = 48. The line within each box represents the median; upper and lower edges of each box are 75th and 25th percentiles, respectively; upper and lower bars indicate the highest and lowest values with less than one interquartile range of distance from the edges of each box. Outliers are not displayed to improve the readability of the plot.

242 243 244

priate horseradish peroxidase-conjugated secondary antibody (Sigma-Aldrich), proteins were detected by electrochemiluminescence (ECL) reaction (EuroClone).

245

Data availability

246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278

Supporting Table S1 contains SNPs ID numbers and locations (according to the Mouse Genome Build GRCm38). Supporting Table S2 contains number of lung tumors (Nlung) and genotypes for each of the 142 ABF4 mice. Gene expression data for the 18,020 genes in the ABF4 population were deposited in the Gene Expression Omnibus (GEO) database (www.ncbi.nlm.nih.gov/geo/) with accession number GSE71232.

Results For the current study, we used existing genotype data on 548 informative SNPs from 142 male ABF4 mice [13]. These mice are closely related by descent, as a limited number of genetic recombinations occurred during the six meioses involved in generating the F4 offspring from the two progenitor inbred strains (A/J and C57BL/6). Therefore, large chromosomal regions are in linkage disequilibrium. As a result, we considered the dataset of 548 SNPs to be adequate for genome-wide coverage (excluding the Y chromosome) for this study. Indeed, the average density of coverage was ~6.0 Mb/SNP. The current study also took advantage of the availability of lung tumor RNA from these same mice. Briefly, ABF4 mice had been treated, at 4 weeks of age, with a single urethane injection to induce lung tumors; at 40 weeks of age, one lung tumor per animal was excised for RNA extraction. In the current study, this RNA was hybridized to Agilent SurePrint G3 Mouse Exon 4 × 180K arrays and, after raw data processing and filtering, expression values for 18,020 genes from all the 19 mouse autosomes and the X chromosome were available for study. Numerous transcripts are controlled by local or distant eQTLs in mouse lung tumors To identify loci involved in the modulation of gene transcription in lung tumor tissue, we did simple interval mapping to find those loci where individuals with different genotypes at the cor-

responding marker displayed different expression levels of any of the 18,020 genes, analyzed in a gene-by-gene manner. Using a conservative cutoff (peak LOD > 5.0), we observed that, at some loci, the LOD score curves had complex shapes, with statistically significant leading or trailing sub-peaks along the chromosome (Supporting Fig. S1). Because the small number of recombinations in inbred mice pedigrees limits the resolution power of linkage analysis, we could not establish if these sub-peaks were due to a single regulatory element or to multiple elements located in close proximity. Therefore, to reduce the number of false positives, we applied stringent criteria in selecting cis- and trans-acting eQTLs among the loci having a peak LOD score above the threshold. In particular, we defined a cis-acting eQTL as a locus that modulates the expression of a gene mapping within 3 Mb on the same chromosome. For trans-acting eQTLs, we took the locus with highest peak among all statistically significant LOD peaks on each chromosome, excluding the chromosome where the target gene mapped. Even though such stringent criteria likely generated some false negatives, they were adopted in order to obtain a solid, albeit underestimated, representation of the regulation network of gene expression in mouse lung tumors. Applying these criteria, we found that the transcript levels of 2602 of the 18,020 examined genes (14.44%) were under genetic control in cis or in trans, or both. Indeed, a subset of 156 genes was regulated by both cis- and trans-acting eQTLs. 123 genes were controlled by multiple trans-regulatory elements. Cis-eQTLs in mouse lung tumor tissue In the analysis for cis-acting loci, we identified 1179 cis-eQTLs (Supporting Table S3), involving 363 markers in linkage with the expression of 1179 genes (6.54% of all examined genes). The peak LOD scores ranged from the cutoff of 5.0 to a high of 65.8 (median, 9.6). The four cis-eQTL with the highest LOD score values (regulating the Hddc3, Glo1, Btbd9 and BC022687 genes) are shown in Fig. 1. Upper panels display the genome-wide LOD score plots for each gene, with a distinguishable LOD score peak in correspondence of the po-

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314

ARTICLE IN PRESS C.E. Cotroneo et al./Cancer Letters ■■ (2016) ■■–■■

4

315 316

Table 1 Overview of cis- and trans-acting loci identified in lung tumor tissue of ABF4 mice.

317 318 319 320 321 322 323 324 325

N. of regulated genesa N. of markersb N. of markers in linkage with a single genec N. of markers in linkage with 2–4 genesc N. of markers in linkage with ≥ 5 genesc

cis-eQTLs

trans-eQTLs

1179 (6.54) 363 (66.2) 126 (34.7) 164 (45.2) 73 (20.1)

1579 (8.76) 238 (43.4) 111 (46.6) 87 (36.6) 40 (16.8)

Percentage values (in parenthesis) refer to: athe total number of analyzed transcripts (n = 18,020), bthe total number of informative SNPs (n = 548), cthe total number of markers in linkage with the expression of genes located in cis or in trans.

326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350

sition of each eQTL. The box-plots (Fig. 1, bottom panels) show the differences in gene expression levels among the 142 ABF4 mice, grouped according to the genotype at the eQTL marker. Overall, 126 markers (34.7% of the total markers in linkage with target genes located in cis) were in linkage with a single gene, 164 markers (45.2%) were in linkage with 2–4 different genes, while 73 markers (20.1%) were in linkage with 5 or more genes (Table 1). Altogether, the 363 markers were in linkage with a median of 2 genes (range, 1–34). The markers in linkage with the most genes were rs3696835 (chromosome 17, 34 genes), followed by rs4226520 (chromosome 7, 24 genes), rs8279354 (chromosome 2, 22 genes) and rs3699056 (chromosome 11, 20 genes). The proportion of genes regulated in cis varied among the chromosomes, and ranged from 2.56% of all the genes mapping on chromosome X to 10.70% on chromosome 17 (median, 6.59%). Trans-eQTLs in mouse lung tumor tissue: an emerging role for the Pas1 locus Regarding trans-acting eQTLs, we detected 238 markers in linkage with the expression of 1579 genes (8.76% of all examined genes), and a total of 1715 trans-eQTLs (Supporting Table S4). The peak LOD scores for these trans-eQTLs ranged from 5.0 to 33.2 (median, 5.7). For the four trans-eQTL with the highest LOD score we reported the

genome-wide LOD plots showing the position of the eQTLs (top panels) and the box-plot (bottom panels) showing differences in gene expression levels among mice grouped according to genotype at the marker in correspondence to the eQTL peak (Fig. 2). The circular genome plot in Fig. 3 represents the interchromosomal trans-eQTLs having a peak LOD > 8.0. These involved 51 different markers in linkage with the expression of 139 genes. Considering the genes examined on each chromosome, the proportion found to be regulated in trans ranged from 7.21% (chromosome 6) to 11.86% (chromosome 18), with a median value of 9.50%. The number of markers in linkage disequilibrium with transacting loci per chromosome ranged from 19 (chromosome 19) to 38 (chromosome 1), with a median of 26.5 markers per chromosome. There were 111 markers (46.6% of the total markers in linkage with genes located in trans) in linkage with the expression levels of just one gene, while the remaining 127 were in linkage with multiple genes (Table 1). The median number of genes in linkage per marker was eight (range, 1–378). In particular, a set of markers ranked high for the number of trans-located genes in linkage with their genotype. On chromosome 17, at a distance of around 7 Mb, markers rs3696835 and rs3693494 were in linkage with 80 and 67 genes, respectively, while on chromosome 9, at a distance of around 30 Mb, markers rs13480071 and rs4135590 were in linkage with 46 and 39 genes, respectively. Moreover, on chromosome 12, the gnf12.046.238 marker was in linkage with 22 genes; on chromosome 13, rs3663223 was in linkage with 38 genes and, on chromosome 17, rs3662820 was in linkage with 21 genes. As clearly emerges from Fig. 3, the highest density of transregulatory loci was located on distal chromosome 6 in correspondence of the Pas1 locus, which is a 468 kb region roughly centered around the Kras gene [25]. The associated region, extending from around 140–148 Mb, controlled altogether 857 genes in trans (54.27% of the total trans-controlled genes) and 14 genes in cis. The linked markers were: rs3672808 (located at 139857209 bp), rs3658783 (144315226 bp), Kras-37bp-repeat (145246451 bp), rs6265387 (147254856 bp), and rs3711088 (148311946 bp). Among

351

352 353 354 355 356 357 358 359 360

Fig. 2. Four trans-eQTLs with the highest LOD score values. For each of the four cis-regulated genes (Sfi1, Gm10116, Pisd and Vmn2r74), top panels show the corresponding genome-wide LOD score plots. A LOD score peak is visible in correspondence of the position of each eQTL. Box plots in the bottom panels show log2-transformed gene expression levels (Y axis) among mice grouped according to genotype at the eQTL marker (X axis). For each marker, “AA” indicates mice homozygous for the A/J-derived allele; “AB” indicates heterozygous mice; “BB” indicates mice homozygous for the C57BL/6-derived allele. For the Sfi1 gene, the number of mice for each genotype at the rs4229817 marker is: AA = 41, AB = 70, BB = 31; for the Gm10116 gene, the number of mice for each genotype at the rs8255275 marker is: AA = 14, AB = 64, BB = 64; for the Pisd gene, the number of mice for each genotype at the rs6236348 marker is: AA = 46, AB = 70, BB = 26; for the Vmn2r74 gene, the number of mice for each genotype at the rs6190775 marker is: AA = 28, AB = 72, BB = 42. The lines within each box represents the median; upper and lower edges of each box are 75th and 25th percentiles, respectively; upper and lower bars indicate the highest and lowest values with less than one interquartile range of distance from the edges of each box. Outliers are not displayed to improve the readability of the plot.

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395 396

ARTICLE IN PRESS C.E. Cotroneo et al./Cancer Letters ■■ (2016) ■■–■■

397 398 399

5

Fig. 3. Circular genome plot of 139 inter-chromosomal trans-eQTLs with peak LOD scores > 8 found in lung tumor tissue of ABF4 mice. Arcs connect each eQTL with its target gene(s); colors of the arcs correspond to those of the chromosomes where the eQTLs are located. Labels are the symbols of the 139 genes modulated by the 139 trans-eQTLs. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

400

401 402 403 404 405 406 407 408 409 410 411

these, only the Kras-37bp-repeat marker maps exactly in the Pas1 locus. Nevertheless, all five markers are in tight linkage disequilibrium, due to the small number of recombination events that occurred in our F4 pedigree. Therefore, all of them were considered as a unique haplotype, corresponding to the Pas1 locus. For all of the genes regulated by this locus, individual mice carrying the A/J-derived allele, i.e. the allele linked to lung tumor susceptibility, displayed significantly different expression levels from those carrying the C57BL/ 6-derived allele. Altogether these results show that the Pas1 locus plays a pleiotropic role in regulating gene expression of lung tumor tissue of ABF4 mice.

Pas1 locus modulates several biochemical pathways and regulates a subset of genes associated with lung tumor multiplicity To determine whether the 871 genes whose expression is regulated in cis or in trans by the Pas1 locus are functionally correlated, we used the Ingenuity Pathway Analysis (IPA) tool to search for significant overrepresentations of genes belonging to known biological pathways. With a significance threshold set at FDR < 0.05, this analysis revealed that 151 of the genes were enriched in 41 different pathways (Supporting Table S5). The identified pathways were involved in a wide range of different biological processes but,

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

412 413 414 415 416 417 418 419 420 421 422

ARTICLE IN PRESS 6

423 424 425 426 427 428

C.E. Cotroneo et al./Cancer Letters ■■ (2016) ■■–■■

Fig. 4. Network representation of the 151 genes involved in the 41 Ingenuity canonical pathways overrepresented among the Pas1-regulated genes. Nodes of the network (circles) represent genes involved in at least one of the identified pathways. Green nodes represent genes whose expression levels are inversely associated with the number of lung tumor susceptibility (A/J-derived) alleles at the Pas1 locus; red nodes represent positively associated genes. The size of each node is proportional to the total number of pathways to which the corresponding gene takes part (range, 1–33). Edges (connecting lines) between two genes were drawn every time they took part in the same pathway. The 14 genes highlighted in yellow are those regulated by Pas1 and also associated with Nlung. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

429 430 431 432 433 434 435 436 437 438 439 440

nonetheless, most of them shared a subset of genes. This interconnectivity of the pathways is shown graphically in Fig. 4 where lines (edges) connect genes (circles) involved in common pathways. Each connection represented a functional interaction of a gene (node) with other genes within a single pathway or among multiple pathways. The node degree (number of edges connecting each node of the network to other nodes) indicates that the genes functionally interacted with a median of 27 other genes (range, 1–138; 75th percentile, 37; Supporting Table S6). The size of each node is proportional to Npathways (total number of pathways in which a gene is involved), which ranged from 1 to 33 (median, 2; 75th per-

centile = 4; Supporting Table S6). By selecting those genes whose node degree and Npathways were higher than the 75th percentile, we identified a subset of 26 genes (Mras, Kras, Pik3c2g, Plcb3, Akt3, Prkcq, Plce1, Plch2, Gnao1, Gng2, Gng7, Camk4, Itga4, Itpr2, Nfkbie, Pdgfd, Adcy8, Arhgef15, Rhof, Rnd3, Lef1, Jak2, Jak1, Pde4c, Egfr, Pde3a) that were strongly interconnected with the other nodes of our network and also recurrently found in the 41 Pas1-associated pathways. These 26 genes may reflect not only the downstream effects of the Pas1 locus, but also of other unrelated loci mapping in tight linkage disequilibrium with it. Furthermore, it is also conceivable

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

441 442 443 444 445 446 447 448 449 450 451

ARTICLE IN PRESS C.E. Cotroneo et al./Cancer Letters ■■ (2016) ■■–■■

452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467

Fig. 5. Venn diagram of genes regulated by the Pas1 locus and transcripts whose expression in lung tumor tissue is associated with lung tumor multiplicity of ABF4 mice. The red circle depicts the 871 genes regulated in cis or in trans by the Pas1 locus in lung tumor tissue with peak LOD scores > 5; of these, 381 were upregulated in individuals carrying the A/J-derived allele of the locus (linked to increased lung tumor susceptibility), while 490 were down-regulated. The blue circle depicts the 1124 transcripts whose expression was associated with individual lung tumor multiplicity (Nlung) in a single lung tumor nodule (significance threshold was placed at FDR < 0.001); of these, 607 were positively associated and 517 were negatively associated with Nlung. Shared transcripts are represented in the area of intersection between the two circles. The number of up-regulated genes is shown in bold; the number of “contra-regulated” genes, i.e. genes displaying opposite directions of regulation in the two sets, is shown in regular font (no genes were contraregulated between the two sets); the number of down-regulated genes is underlined. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

468 469 470 471 472 473 474

that, even for those genes controlled by Pas1, only a subset of them is causally involved in lung tumor susceptibility. Therefore, we used an independent approach to investigate how many of the genes targeted by this locus were also associated with this phenotype. First, we analyzed the correlations between the transcript levels of each gene expressed in lung tumor tissue of the ABF4 population (one

7

nodule per mouse) and individual lung tumor multiplicity (Nlung), used as an indicator of lung tumor susceptibility. By linear regression analysis, expression levels of 1124 out of 18,020 expressed genes were significantly associated with Nlung (FDR < 0.001) (Supporting Table S7). Of these transcripts, 607 were positively associated with Nlung, meaning that they were more expressed in animals with higher tumor multiplicity, while 517 were inversely associated, showing lower transcript levels in the tumors of mice with higher Nlung. Then, we compared these 1124 Nlung-associated transcripts with the 871 genes regulated in cis and in trans by Pas1, and we observed an overlapping set of 401 genes (Fig. 5). This set included 189 genes found to be up-regulated and 212 genes found to be down-regulated by both analyses. In no case was the expression of a gene contra-regulated between the two sets: the direction of association with increasing number of lung tumor susceptibility (A/ J-derived) alleles at the Pas1 locus was always concordant with that observed in association with higher lung tumor multiplicity. Of the 26 genes most represented in the Pas1-affected pathways, 14 were also associated with lung tumor multiplicity: Kras, Plcb3, Akt3, Plch2, Itga4, Itpr2, Pdgfd, Rhof, Rnd3, Jak2, Jak1, Pde4c, Egfr, Pde3a. Fig. 6 compares the expression levels of these 14 genes in the 142 ABF4 mice, grouped according to their genotype at the Pas1 locus. Some genes display an allele-dosage effect of the eQTL, meaning that their level of expression in heterozygous animals is intermediate between the two homozygous groups. In other cases, the expression levels in heterozygous animals are similar to either of the two homozygous. Altogether, the 14 genes may constitute a signature for the biological activity of this locus in lung tumorigenesis. Pi3k/Akt and Ras/Erk signaling pathways are activated in lung tumor tissue of ABF4 mice The above identified signature included genes belonging to the Ras–Erk and Pi3k–Akt signaling pathways [26,27]. To test the ac-

475

476 477 478 479 480 481

Fig. 6. Expression levels of the 14 genes, modulated in cis or in trans by the Pas1 locus (peak LOD score ≥5), most represented among the Pas1-affected pathways and associated with lung tumor multiplicity in 142 ABF4 mice. The genotype at the Pas1 locus was assessed by taking the genotype of the Kras-37bp-repeat marker. This marker maps within the Kras gene, which in turn maps approximately in the center of the locus. “AA” stands for mice homozygous for the A/J-derived allele (n = 57); “AB” stands for heterozygous mice (n = 72); “BB” stands for mice homozygous for the C57BL/6-derived allele (n = 13). The line within each box represents the median; upper and lower edges of each box are 75th and 25th percentiles, respectively; upper and lower bars indicate the highest and lowest values with less than one interquartile range of distance from the edges of each box. Outliers are not displayed to improve the readability of the plots.

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516

ARTICLE IN PRESS C.E. Cotroneo et al./Cancer Letters ■■ (2016) ■■–■■

8

517 518 519 520 521 522

Fig. 7. Pi3k–Akt signaling pathway is more activated in tumors of mice carrying the Pas1 susceptible allele than those carrying the resistant one. Normal lung (N) and lung tumor tissue (T) from ABF4 mice homozygous for the A/J- or C57BL/6-derived allele (AA and BB, respectively) and heterozygous mice (AB) were analyzed by western blot for total and phosphorylated Akt (60 kDa), Mek1/2 (45 kDa) and Erk1/2 (42 and 44 kDa, respectively) proteins. Vinculin (124 kDa) was shown as loading control.

523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554

tivation of these pathways, we analyzed by western-blot the phosphorylation status of Mek1/2, Erk1/2 (for Ras–Erk pathway) and Akt (for Pi3k–Akt pathway) proteins both in normal lung and in lung tumors from ABF4 mice carrying the three different genotypes at the Kras-37bp-repeat marker. As expected, all the phospho-proteins examined showed higher levels in lung tumor tissue than in normal tissue, indicating the activation of these pathways in lung tumors (Fig. 7). Of note, pAKT differed also among tumors from mice carrying the A/J derived allele versus mice homozygous for the C57BL/6 derived allele, whereas no differences among genotypes were observed for pMek1/2 and pErk1/2 (Fig. 7). These results indicated the Pi3k–Akt signaling pathway as a main pathway involved in the mechanism of action of Pas1 lung tumor susceptibility. Discussion In this study, carried out in an ABF4 mouse intercross population in which a major locus (Pas1) modulates individual susceptibility to urethane-induced lung tumorigenesis as a monogenic trait [13], we investigated how whole-genome germline variations modulate the transcriptional landscape of lung cancer tissue. The key advantage of this F4 population is that it has three-fold more recombination events than a conventional F2 population of the same size. The higher number of recombination events allows to map genetic loci influencing phenotypes of interest with a better resolution. After assessing the transcript levels of 18,020 genes in a single lung tumor nodule from each mouse, we found that 2,602 genes (14.44%) were under germline control in this tissue. Indeed, by genetic linkage analysis, we identified 1179 and 1579 genes regulated by cis- and trans-eQTLs, respectively (156 genes are under both cis and trans control). The Pas1 locus presented the highest number

of target genes in trans (857 genes) and also regulated 14 genes in cis. Ingenuity Pathway Analysis highlighted a list of 41 biological pathways regulated by Pas1. All these pathways were functionally related, as they shared a subset of 26 genes. Moreover, in the same ABF4 population, we studied the association between the expression levels of the 18,020 transcripts in lung tumor tissue and the total number of lung tumors developed by each mouse (Nlung). This analysis showed that an increase in lung tumor multiplicity was associated with the expression of 1,124 Nlung-associated transcripts. The expression of 401 of these transcripts was also under the genetic control exerted by the Pas1 locus and, among these, 14 of the 26 genes common to the Pas1-regulated pathways were confirmed. The number of cis-eQTLs identified in our study (n = 1179) is comparable to those found so far in different types of human tumors (ranging from around 1000–5300) [5,28,29], including lung adenocarcinoma (1306 cis-eQTLs) [5]. Of note, none of these studies investigated the trans-effects of germline variations on gene expression, due to the low statistical power caused by their smaller effect size [3]. In our study, the use of a mouse pedigree partially overcame this limitation and allowed us to detect the presence of 1715 eQTLs acting in trans. In accordance with their expected smaller effect sizes, LOD scores for these eQTLs were, however, lower than those of cis-acting loci. The observation that the Pas1 locus is linked to 54.27% of the total trans-controlled genes in lung tumor tissue confirms the central role of this locus in mouse lung carcinogenesis, and highlights the strong involvement of Pas1 not only in triggering the susceptibility to urethane-induced lung tumorigenesis, but also in affecting the transcriptional regulation of already developed tumors. Interestingly, this observation indicates that a single genetic locus is enough to induce an unexpectedly large difference in the transcriptome profile in tumors of the same size that develop in closely related individuals. Even though it is conceivable that not all these differences are causally related to individual genetic predisposition to lung carcinogenesis, they nevertheless may alter the biological properties of each tumor, contributing to the well-known phenomenon of interindividual tumor heterogeneity. We speculate that this phenomenon may underlie individual differences in genetic susceptibility to lung cancer and may explain, at least in part, the discrepancies in prognosis and response to therapies commonly observed among patients with the same cancer histotype and clinical stage [30,31]. Further analyses of the germline transcriptional regulation in samples of matched normal lung tissue from this same cross will allow an even better understanding of the genetic architecture of lung cancer susceptibility in mice carrying the Pas1 susceptibility allele. Indeed, these results will allow to discriminate identifying variations affecting gene expression in both normal and tumor lung tissue from those exclusively acting in normal tissue or in neoplastic cells. Similar studies have been already carried out in mouse skin [32]. This work observed a decrease in the number of eQTLs (especially trans-acting loci) in malignant cancer tissue, compared to normal skin, coupled with the existence of a set of tumor-specific eQTLs. In our work, eQTL analysis alone did not allow us to understand the biological role of Pas1-regulated genes and their putative involvement in lung tumor susceptibility. By analyzing the biological pathways overrepresented among the genes controlled by this locus, we identified 41 Pas1-controlled pathways that were involved in a wide set of different biological functions. Among these pathways, we found a set of 26 recurrently represented genes, including membrane receptors (the integrin alpha 4 Itga4 and the epidermal growth factor receptor Egfr), small GTPase proteins (Mras, Kras, Gnao1, Gng2, Gng7, Rnd3) and G-protein associated factors (Arhgef15, Rhof), two Janus kinases (Jak1 and Jak2) along with other signal transduction proteins such as three phospholipase C isoenzymes (Plcb3, Plce1 and Plch2), a catalytic subunit of the

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598 599 600 601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620

ARTICLE IN PRESS C.E. Cotroneo et al./Cancer Letters ■■ (2016) ■■–■■

621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 656 657 658 659 660 661 662 663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684 685 686

phosphatidylinositol-4-phosphate 3-kinase (Pik3c2), the inositol 1,4,5trisphosphate receptor (Itpr2), the homolog 3 of the AKT protein (Akt3) and transcription factors (Nfkbie and Lef1). Of note, all of these genes are components or targets of the Ras–Erk and Pi3k–Akt signaling pathways (reviewed in [26,27]), two interlinked pathways that regulate cell proliferation and survival with well-documented involvement in several types of cancer. Fourteen of these genes were also associated with individual lung tumor multiplicity. By western blot analysis, we observed that the Ras–Erk pathway is equally activated among the three genotypes. Interestingly, Pi3k–Akt pathway activation is differently activated among mice carrying one or two A/J-derived alleles and those homozygous for the C57BL/6-derived allele. This finding suggests that the functionality of the Pas1 locus is mediated by the modulation of the expression of genes belonging to Pi3K–Akt signaling pathway. Of particular note is the Kras gene, which is an upstream regulator of both Ras–Erk and Pi3k–Akt pathways and which maps within the Pas1 locus itself. This gene is often mutated in both mouse and human lung cancer, resulting in the overactivation of its protein product [33–35]. We already reported that the Kras-4A isoform is a strong candidate for the effect of Pas1 on lung tumor susceptibility in this cross, since its expression levels are significantly higher in susceptible individuals [13]. Here, the observation that lung tumors carrying the susceptible Pas1 allele have a differential activation of Pi3k–Akt pathway, in comparison to tumors from resistant individuals (specifically, tumors carrying the A/J-derived allele had higher levels of pAKT than tumors from mice carrying the C57BL/6-derived allele), once again highlights the pivotal role of Kras in ABF4 lung tumors. Considering the differential expression of Kras-4A in lung tumors of susceptible versus resistant mice of this same intercross [13], we suggest that the activation of Pi3k–Akt signaling pathway reflects the specific biological effect of this isoform. In conclusion, our results provide a genome-wide picture of the germline regulation of gene expression in mouse lung cancer. They highlight the importance of Pas1 as a tumor-modifier locus in this cross, and attribute to it a novel role as a major regulator of transcription in lung tumor nodules. Moreover, these results strengthen the candidacy of the Kras-4A isoform as the major effector of Pas1 function. Acknowledgements Valerie Matarese provided scientific editing. This work was funded in part by a grant from Associazione and Fondazione Italiana Ricerca Cancro (AIRC) (http://www.airc.it/) to GM (IG-2014, no. 15797). The funders had no role in the design and conduct of the study, in the collection, analysis, and interpretation of the data, and in the preparation, review, or approval of the manuscript. Conflict of interest All authors declare no conflict of interest.

Q7 Uncited references [36] Appendix: Supplementary material Supplementary data to this article can be found online at doi:10.1016/j.canlet.2016.02.054. References [1] W. Cookson, L. Liang, G. Abecasis, M. Moffatt, M. Lathrop, Mapping complex disease traits with global gene expression, Nat. Rev. Genet. 10 (2009) 184–194.

9

[2] Y. Gilad, S.A. Rifkin, J.K. Pritchard, Revealing the architecture of gene regulation: the promise of eQTL studies, Trends Genet. 24 (2008) 408–415. [3] F.W. Albert, L. Kruglyak, The role of regulatory variation in complex traits and disease, Nat. Rev. Genet. 16 (2015) 197–212. [4] M.T. Maurano, R. Humbert, E. Rynes, R.E. Thurman, E. Haugen, H. Wang, et al., Systematic localization of common disease-associated variation in regulatory DNA, Science 337 (2012) 1190–1195. [5] Q. Li, A. Stram, C. Chen, S. Kar, S. Gayther, P. Pharoah, et al., Expression QTL-based analyses reveal candidate causal genes and loci across five tumor types, Hum. Mol. Genet. 23 (2014) 5294–5302. [6] A. Galvan, J.P. Ioannidis, T.A. Dragani, Beyond genome-wide association studies: genetic heterogeneity and individual predisposition to cancer, Trends Genet. 26 (2010) 132–141. [7] G. Manenti, T.A. Dragani, Pas1 haplotype-dependent genetic predisposition to lung tumorigenesis in rodents: a meta-analysis, Carcinogenesis 26 (2005) 875–882. [8] T.A. Dragani, 10 years of mouse cancer modifier loci: human relevance, Cancer Res. 63 (2003) 3011–3018. [9] K.K. Youngren, D. Coveney, X. Peng, C. Bhattacharya, L.S. Schmidt, M.L. Nickerson, et al., The ter mutation in the dead end gene causes germ cell loss and testicular germ cell tumours, Nature 435 (2005) 360–364. [10] C.A. Ruivenkamp, T. van Wezel, C. Zanon, A.P. Stassen, C. Vlcek, T. Csikos, et al., Ptprj is a candidate for the mouse colon-cancer susceptibility locus Scc1 and is frequently deleted in human cancers, Nat. Genet. 31 (2002) 295–300. [11] Y.G. Park, X. Zhao, F. Lesueur, D.R. Lowy, M. Lancaster, P. Pharoah, et al., Sipa1 is a candidate for underlying the metastasis efficiency modifier locus Mtes1, Nat. Genet. 37 (2005) 1055–1062. [12] M. Iguchi, M. Osanai, Y. Hayashi, F. Koentgen, G.H. Lee, The error-prone DNA polymerase iota provides quantitative resistance to lung tumorigenesis and mutagenesis in mice, Oncogene 33 (2014) 3612–3617. [13] A. Dassano, F. Colombo, G. Trincucci, E. Frullanti, A. Galvan, A. Pettinicchio, et al., Mouse pulmonary adenoma susceptibility 1 locus is an expression QTL modulating kras-4A, PLoS Genet. 10 (2014) e1004307. [14] M.B. Shimkin, G.D. Stoner, Lung tumors in mice: application to carcinogenesis bioassay, Adv. Cancer Res. 21 (1975) 1–58. [15] A.M. Malkinson, Primary lung tumors in mice as an aid for understanding, preventing, and treating human adenocarcinoma of the lung, Lung Cancer 32 (2001) 265–279. [16] M.E. Ritchie, B. Phipson, D. Wu, Y. Hu, C.W. Law, W. Shi, et al., Limma powers differential expression analyses for RNA-sequencing and microarray studies, Nucleic Acids Res. (2015) doi:10.1093/nar/gkv007. [17] R.C. Gentleman, V.J. Carey, D.M. Bates, B. Bolstad, M. Dettling, S. Dudoit, et al., Bioconductor: open software development for computational biology and bioinformatics, Genome Biol. 5 (2004) R80. [18] P. Flicek, M.R. Amode, D. Barrell, K. Beal, K. Billis, S. Brent, et al., Ensembl 2014, Nucleic Acids Res. 42 (2014) D749–D755. [19] S. Durinck, P.T. Spellman, E. Birney, W. Huber, Mapping identifiers for the integration of genomic datasets with the R/bioconductor package biomaRt, Nat. Protoc. 4 (2009) 1184–1191. [20] K.W. Broman, H. Wu, S. Sen, G.A. Churchill, R/qtl: QTL mapping in experimental crosses, Bioinformatics 19 (2003) 889–890. [21] E. Lander, L. Kruglyak, Genetic dissection of complex traits: guidelines for interpreting and reporting linkage results, Nat. Genet. 11 (1995) 241–247. [22] M. Krzywinski, J. Schein, I. Birol, J. Connors, R. Gascoyne, D. Horsman, et al., Circos: an information aesthetic for comparative genomics, Genome Res. 19 (2009) 1639–1645. [23] Y. Benjamini, Y. Hochberg, Controlling the false discovery rate: a practical and powerful approach to multiple testing, J. Roy. Statist. Soc. Ser. B 57 (1995) 289–300. [24] P. Shannon, A. Markiel, O. Ozier, N.S. Baliga, J.T. Wang, D. Ramage, et al., Cytoscape: a software environment for integrated models of biomolecular interaction networks, Genome Res. 13 (2003) 2498–2504. [25] G. Manenti, F. Galbiati, R. Gianni-Barrera, A. Pettinicchio, A. Acevedo, T.A. Dragani, Haplotype sharing suggests that a genomic segment containing six genes accounts for the pulmonary adenoma susceptibility 1 (Pas1) locus activity in mice, Oncogene 23 (2004) 4495–4504. [26] Y. Pylayeva-Gupta, E. Grabocka, D. Bar-Sagi, RAS oncogenes: weaving a tumorigenic web, Nat. Rev. Cancer 11 (2011) 761–774. [27] R. Sever, J.S. Brugge, Signal transduction in cancer, Cold Spring Harb. Perspect. Biol. 5 (2015) doi:10.1101/cshperspect.a006098. [28] Q. Li, J.H. Seo, B. Stranger, A. McKenna, I. Pe’er, T. Laframboise, et al., Integrative eQTL-based analyses reveal the biology of breast cancer risk loci, Cell 152 (2013) 633–641. [29] H. Ongen, C.L. Andersen, J.B. Bramsen, B. Oster, M.H. Rasmussen, P.G. Ferreira, et al., Putative cis-regulatory drivers in colorectal cancer, Nature 512 (2014) 87–90. [30] A. Yoshizawa, N. Motoi, G.J. Riely, C.S. Sima, W.L. Gerald, M.G. Kris, et al., Impact of proposed IASLC/ATS/ERS classification of lung adenocarcinoma: prognostic subgroups and implications for further revision of staging based on analysis of 514 stage I cases, Mod. Pathol. 24 (2011) 653–664. [31] J. Zhang, J. Wu, Q. Tan, L. Zhu, W. Gao, Why do pathological stage IA lung adenocarcinomas vary from prognosis?: a clinicopathologic study of 176 patients with pathological stage IA lung adenocarcinoma based on the IASLC/ ATS/ERS classification, J. Thorac. Oncol. 8 (2013) 1196–1202. [32] D.A. Quigley, M.D. To, I.J. Kim, K.K. Lin, D.G. Albertson, J. Sjolund, et al., Network analysis of skin tumor progression identifies a rewired genetic architecture

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 Q5 704 705 706 707 708 709 Q6 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765 766 767 768 769 770 771 772

ARTICLE IN PRESS 10

773 774 775 776 777 778 779 780 781

C.E. Cotroneo et al./Cancer Letters ■■ (2016) ■■–■■

affecting inflammation and tumor susceptibility, Genome Biol. 12 (2011) [Epub 2011 Jan 18]. R5-2011-12-1-r5. [33] M. You, U. Candrian, R.R. Maronpot, G.D. Stoner, M.W. Anderson, Activation of the ki-ras protooncogene in spontaneously occurring and chemically induced lung tumors of the strain A mouse, Proc. Natl. Acad. Sci. U.S.A. 86 (1989) 3070–3074. [34] F.C. Re, G. Manenti, M.G. Borrello, M.P. Colombo, J.H. Fisher, M.A. Pierotti, et al., Multiple molecular alterations in mouse lung tumors, Mol. Carcinog. 5 (1992) 155–160.

[35] M.S. Lawrence, P. Stojanov, C.H. Mermel, J.T. Robinson, L.A. Garraway, T.R. Golub, et al., Discovery and saturation analysis of cancer genes across 21 tumour types, Nature 505 (2014) 495–501. [36] M. Iguchi, M. Osanai, Y. Hayashi, F. Koentgen, G. Lee, The error-prone DNA polymerase iota provides quantitative resistance to lung tumorigenesis and mutagenesis in mice, Oncogene 33 (2014) 3612–3617.

Please cite this article in press as: Chiara E. Cotroneo, et al., Expression quantitative trait analysis reveals fine germline transcript regulation in mouse lung tumors, Cancer Letters (2016), doi: 10.1016/j.canlet.2016.02.054

782 783 784 785 786 787