Abstract
Marjolein J.A. Weerts1, Anieta M. Sieuwerts1, Marcel Smid1, Maxime P. Look1, John A. Foekens1, Stefan Sleijfer1, John W.M. Martens1
1Department of Medical Oncology and Cancer Genomics, Erasmus MC Cancer Institute, Erasmus University Medical Center, Rotterdam, The Netherlands
Correspondence to:
Marjolein J.A. Weerts, email: [email protected]
Keywords: breast cancer, mitochondrial DNA, mtDNA content, epithelial-to-mesenchymal transition, distant metastasis
Received: November 27, 2015 Accepted: March 28, 2016 Published: April 11, 2016
ABSTRACT
Reduced mitochondrial DNA (mtDNA) content in breast cancer cell lines has been associated with transition towards a mesenchymal phenotype, but its clinical consequences concerning breast cancer dissemination remain unidentified. Here, we aimed to clarify the link between mtDNA content and a mesenchymal phenotype and its relation to prognosis of breast cancer patients. We analyzed mtDNA content in 42 breast cancer cell lines and 207 primary breast tumor specimens using a combination of quantitative PCR and array-based copy number analysis. By associating mtDNA content with expression levels of genes involved in epithelial-to-mesenchymal transition (EMT) and with the intrinsic breast cancer subtypes, we could not identify a relation between low mtDNA content and mesenchymal properties in the breast cancer cell lines or in the primary breast tumors. In addition, we explored the relation between mtDNA content and prognosis in our cohort of primary breast tumor specimens that originated from patients with lymph node-negative disease who did not receive any (neo)adjuvant systemic therapy. When patients were divided based on the tumor quartile levels of mtDNA content, those in the lowest quarter (≤ 350 mtDNA molecules per cell) showed a poorer 10-year distant metastasis-free survival than patients with > 350 mtDNA molecules per cell (HR 0.50 [95% CI 0.29–0.87], P = 0.015). The poor prognosis was independent of established clinicopathological markers (HR 0.54 [95% CI 0.30–0.97], P = 0.038). We conclude that, despite a lack of evidence between mtDNA content and EMT, low mtDNA content might provide meaningful prognostic value for distant metastasis in breast cancer.
INTRODUCTION
Mitochondria play a role in many cellular processes including oxidative phosphorylation, redox homeostasis, controlling calcium levels for regulation of signal transduction pathways, and the intrinsic apoptotic pathway [1]. The mitochondria contain their own genome – mtDNA – encoding their own translational machinery and 13 crucial proteins for the oxidative phosphorylation system. Related to energy needs, numerous mtDNA molecules may exist in a single cell. This number is not only dependent on the amount of mitochondria per cell but also on the number of mtDNA molecules per mitochondrion. Broad ranges in mtDNA content have been reported, from a few molecules in embryonic and pluripotent stem cells [2, 3] up to several thousands in subcutaneous adipocytes [4] or cardiac myocytes [5]. The cell-specific mtDNA content is assumed to be fairly stable under physiological conditions but can be altered by stress such as exogenous toxins [6], viral infection [7] and by genetic mutations [8]. The effects of changes in mtDNA content are illustrated in several mtDNA depletion syndromes [9], which are all characterized by impaired energy production.
Several studies examined mtDNA content in the context of cancer but so far no clear picture has emerged. In preclinical models, depletion of mtDNA yielded both increased and decreased in vitro tumorigenic phenotypes [10–17]. The in vivo findings using mouse xenografts are indecisive as well, as both gain and loss of tumorigenic potential upon mtDNA depletion has been reported [16–21]. Additionally, contradictory findings have been described for mtDNA content in human tumor specimens compared to their healthy counterparts in multiple cancer types (as reviewed in [22, 23]).
With regard to breast cancer, the impact of the mtDNA content on phenotype, prognosis and drug response has been investigated in several studies. Lower mtDNA content is observed in approximately 70% of breast cancer specimens when compared to their surrounding normal epithelium [24–31]. There are indications that low mtDNA content in breast cancer may yield a more aggressive phenotype and altered therapy responses. First, depletion of mtDNA in in vitro models affects the mRNA and protein expression levels of several genes involved in epithelial-to-mesenchymal transition (EMT) [12, 14]. The transition towards the mesenchymal phenotype has been implied as an essential mechanism contributing to cancer dissemination [32]. Consequently, low mtDNA content as a marker for the mesenchymal phenotype potentially identifies tumor aggressiveness. Second, a link between reduced mtDNA content and resistance to anti-estrogen regimens has been established in in vitro models [33]. Nevertheless, no association between estrogen receptor status and mtDNA content was observed in breast tumors [24–29]. Also, reduced mtDNA content was linked to a shift in drug response for breast cancer cell lines [17, 24, 34]. An in vitro reduction in mtDNA content revealed increased sensitivity to cisplatin [17] and doxorubicin [24], but also decreased sensitivity to vincristine, paclitaxel and – in contrast to a previous study – doxorubicin [34]. In a small patient cohort, low mtDNA content was associated with longer disease-free survival in patients receiving adjuvant chemotherapy, whereas this was not the case for patients not receiving adjuvant treatment [24]. Few additional studies reported on breast cancer patient disease free- or overall survival in relation to tumorous mtDNA content [25–27]. However, these studies had either relatively small sample sizes or no information about treatments administered, the mtDNA content determination methods varied, and results were inconclusive.
Here, we further explore the putative link between mtDNA content and prognostic features in breast cancer. In a broad panel of human breast cancer cell lines the link between mtDNA content and a mesenchymal phenotype was studied by correlating it with expression levels of EMT-related genes and with the intrinsic subtypes of breast cancer [35, 36]. In a well-defined patient cohort of primary breast tumor specimens [37], tumor mtDNA content was examined in relation to expression levels of EMT-related genes, to the intrinsic subtypes, as well as to established clinicopathological variables. Primarily, in our cohort of primary breast cancer patients with lymph node-negative disease who did not receive any (neo)adjuvant systemic therapy, we examined the prognostic value of mtDNA content using distant metastasis-free survival as the main endpoint.
RESULTS
mtDNA content in breast cancer cell lines and primary tumor specimens
In total, we analyzed DNA extracts from 42 breast cancer cell lines and 207 primary tumor specimens. Multiplex real time quantitative PCR (qPCR) targeting a nuclear-encoded and a mitochondrial-encoded gene combined with array-based copy number changes of the nuclear-encoded gene to correct for sample specific somatic variation at the reference locus was used to obtain the mtDNA content in the DNA extracts of these samples. Inter-assay variability of the multiplex qPCR assay was monitored using the calibration curves taken along in each run (n = 7). Amplification in the calibration curve samples was linear between 0.16 and 16 ng DNA per reaction with mean efficiencies and standard error of 97.6 ± 4.4% for nuclear encoded HMBS and 91.5 ± 5.2% for mitochondrial encoded MT-TL1. Copy number variation of the nuclear encoded HMBS gene was observed in 39% of the breast cancer cell lines including 1 with homozygous loss, 12 with heterozygous loss and 3 with gain, and in 14% of the primary tumor specimens including 20 with heterozygous loss and 10 with gain. Because of a homozygous HMBS loss in SUM1315MO2, this cell line was excluded from further analysis. Because of absence of HMBS qPCR signal amplification in three primary tumor specimens, these samples were excluded from further analysis as well. The median mtDNA content and interquartile ranges (IQR) in the 41 breast cancer cell lines and in the 204 primary breast tumor specimens were respectively 489 (IQR 360) and 462 (IQR 294) mtDNA molecules per cell.
mtDNA content and the mesenchymal characteristics
In vitro reduction of mtDNA content has been linked to changes in expression of the EMT-related genes CDH1 [12, 14], CDH2 [14], ESRP1 [14], FN1 [14], MMP9 [14], SNAI1 [14], SNAI2 [14], TGFB1 [12], TGFBR1 [12], TWIST1 [14] and VIM [12, 14]. To address whether a more mesenchymal phenotype is a physiological characteristic linked to low mtDNA content [12, 14], we analyzed the relation between mtDNA content and the RNA expression levels of genes related to EMT. Expression data for the above mentioned genes were available for 40 of the 41 breast cancer cell lines and all 204 primary breast tumor specimens. Expression data of TGFBR1 was excluded because the probe gave expression levels close to background noise. Correlation between gene expression levels and mtDNA content did not exceed a correlation coefficient ρ of 0.35, and we could not demonstrate statistical significance after correction for multiple testing (all P > 0.027, Supplementary Table S1) in the breast cancer cell lines. In our primary breast tumor specimens, correlation between mtDNA content and the expression of ESRP1 (ρ = 0.25, P < 0.001), SNAI1 (ρ = 0.23, P < 0.001) and TGFB1 (ρ = 0.18, P < 0.01) was statistically significant after correction for multiple testing (Supplementary Table S1). To further explore the link between mtDNA content and EMT, we analyzed the association between mtDNA content and the intrinsic subtypes of breast cancer, which have been assigned with epithelial or mesenchymal characteristics [36, 38, 39]. Comparisons between the intrinsic subtypes for both the breast cancer cell lines as well as the primary tumor specimens did not show differences in mtDNA content among the subtypes (P > 0.05) (Table 1).
Table 1: mtDNA content in the intrinsic breast cancer subtypes
| Subtype | n (%) | mtDNA content (IQR) | P-value |
|---|---|---|---|---|
Breast cancer cell lines | Basal | 5 (12.5%) | 269 (149) | 0.1† |
ERBB2 | 7 (17.5%) | 620 (521) | ||
Luminal | 19 (47.5%) | 518 (359) | ||
Normal | 9 (22.5%) | 489 (142) | ||
Primary breast tumor specimens | Basal | 65 (31.8%) | 454 (287) | 0.8† |
ERBB2 | 35 (17.2%) | 566 (351) | ||
Luminal A | 56 (27.5%) | 423 (224) | ||
Luminal B | 40 (19.6%) | 514 (343) | ||
Normal | 8 (3.9%) | 377 (286) |
Median mtDNA content [number of mtDNA molecules per cell] with interquartile range (IQR) for each group and corresponding probabilities (P value) for equal distribution using Kruskal-Wallis one-way analysis of variance (†).
Association of tumor mtDNA content with established prognostic clinicopathological variables
In our patient cohort, we analyzed tumor mtDNA content in relation to patient age at diagnosis, menopausal status, tumor size, histological grade, estrogen receptor status, progesterone receptor status and ERBB2 amplification (Table 2). Because the currently used conventional histological grade (modified Bloom-Richardson) was not available for nearly 25% of our cohort – samples originated from multiple hospitals and from time periods when histological grading according to Bloom-Richardson was not common – we included a molecular grading system shown to be equivalent, the qRT-PCR genomic grade index (GGI) [40]. There were no statistically significant associations between the tumor mtDNA content and age at diagnosis, menopausal status, estrogen receptor status, progesterone receptor status or ERBB2 amplification status (P > 0.05). However, tumors smaller than 2 cm had statistically significant lower mtDNA content (median 421 mtDNA molecules per cell) compared to tumors larger than 2 cm (median 514 mtDNA molecules per cell) (Mann-Whitney P = 0.019). In addition, tumor mtDNA content varied between the GGI groups (Kruskal-Wallis P = 0.028), with the highest mtDNA content in the GGI group representing poorly differentiated high grade 3 tumors (median 523 mtDNA molecules per cell). However, we did not observe a significant trend in tumor mtDNA content across the GGI groups (Cuzick’s test for trend P = 0.066).
Table 2: Association between clinicopathological variables and mtDNA content
Variable | Group | n (%) | mtDNA content (IQR) | P-value |
|---|---|---|---|---|
Age at diagnosis | ≤ 40 | 21 (10.3%) | 491 (532) | 0.21† |
> 40–55 | 88 (43.1%) | 433 (273) | ||
> 55–70 | 64 (31.4%) | 466 (259) | ||
> 70 | 31 (15.2%) | 546 (490) | ||
Menopausal status | Pre | 99 (48.5%) | 427 (323) | 0.15# |
Post | 105 (51.5%) | 500 (280) | ||
Tumor size | ≤ 2 cm | 99 (48.5%) | 421 (280) | 0.019# |
> 2 cm | 105 (51.5%) | 514 (382) | ||
Genomic Grade Index | 1 | 35 (17.2%) | 440 (225) | 0.028† |
2 | 59 (32.4%) | 410 (260) | ||
3 | 103 (50.5%) | 523 (389) | ||
Estrogen receptor status | Negative | 87 (42.7%) | 483 (328) | 0.12# |
Positive | 115 (56.4%) | 424 (290) | ||
Progesterone receptor status | Negative | 97 (47.5%) | 480 (335) | 0.073# |
Positive | 96 (47.1%) | 413 (281) | ||
ERBB2 amplification | Negative | 169 (82.8%) | 454 (287) | 0.46# |
Positive | 29 (14.2%) | 463 (385) |
Number of patients and corresponding median mtDNA content [number of mtDNA molecules per cell] with interquartile range (IQR) for each group and corresponding probabilities (P value) for either equal distribution using Mann-Whitney U test (#) or Kruskal-Wallis one-way analysis of variance (†). Due to missing values the numbers of samples per variable do not always add up to 204.
Distant metastasis-free survival and primary tumor mtDNA content
Finally, we studied in our patient cohort the prognostic value of tumor mtDNA content with respect to the length of distant metastasis-free survival. All included breast cancer patients presented as lymph node-negative and did not receive any (neo)adjuvant systemic treatment. The distribution of mtDNA content in our cohort was skewed and could not be normalized by transformation (Skewness and Kurtosis test P < 0.05). To assess tumor mtDNA content for the length of metastasis-free survival in our exploratory analysis, we first divided the cohort based on mtDNA content quartiles in four groups (Q1–Q4). Because the patients within the first quarter Q1 presented a different rate of metastasis-free survival compared to the other quarters (Q2–Q4) (Supplementary Figure S1), we divided the cohort in two patient groups of low mtDNA content (Q1 with mtDNA content ≤ 350 mtDNA molecules per cell) versus the rest (Q2–Q4 with mtDNA content > 350 mtDNA molecules per cell). To visualize the length of metastasis-free survival as a function of the levels of tumor mtDNA content (Q1 vs Q2–Q4) we used the Kaplan-Meier survival analysis method (Figure 1). Patients in the low mtDNA content group Q1 showed a higher metastasis probability (log-rank P = 0.047). In univariate and multivariable Cox regression analysis including only the 186 patients with no missing values (Table 3), patients in the Q2–Q4 mtDNA content group showed a longer distant metastasis-free survival compared to patients in the low mtDNA tumor content group (univariate: HR 0.50, 95% CI: 0.29–0.87, P = 0.015; multivariable: HR 0.54, 95% CI: 0.30–0.97, P = 0.038).
Figure 1: Kaplan-Meier curve showing probability of distant metastasis-free survival as a function of tumor mtDNA content of 204 patients (60 events). Numbers of patients at risk at 24 month time intervals are indicated.
Table 3: Univariate and multivariable analyses for distant metastasis-free survival in lymph node-negative patients who did not receive any (neo)adjuvant systemic therapy
Variable | Group | n (%) | Univariate | Multivariate | ||
|---|---|---|---|---|---|---|
Hazard ratio (95% CI) | P value | Hazard ratio (95% CI) | P value | |||
Age at diagnosis | ≤ 40 | 19 (10.2%) | 1 |
| 1 |
|
> 40–55 | 81 (43.5%) | 0.40 (0.19–0.88) | 0.022 | 0.34 (0.15–0.76) | 0.009 | |
> 55–70 | 59 (31.7%) | 0.42 (0.19–0.95) | 0.038 | 0.27 (0.08–0.92) | 0.037 | |
> 70 | 27 (14.5%) | 0.39 (0.14–1.05) | 0.062 | 0.25 (0.06–0.97) | 0.045 | |
Menopausal status | Pre | 91 (48.9%) | 1 |
| 1 |
|
Post | 95 (51.1%) | 0.94 (0.55–1.61) | 0.8 | 1.55 (0.54–4.44) | 0.4 | |
Tumor size | ≤ 2 cm | 87 (46.8%) | 1 |
| 1 |
|
> 2 cm | 99 (53.2%) | 1.06 (0.62–1.83) | 0.8 | 0.92 (0.54–1.64) | 0.8 | |
Genomic Grade Index | 1 | 33 (17.7%) | 1 |
| 1 |
|
2 | 56 (30.1%) | 1.62 (0.63–4.18) | 0.3 | 1.43 (0.54–3.76) | 0.5 | |
3 | 97 (52.2%) | 2.65 (1.03–6.83) | 0.043 | 2.52 (0.95–6.66) | 0.063 | |
Progesterone receptor status | Negative | 96 (51.6%) | 1 |
| 1 |
|
Positive | 90 (48.4%) | 0.74 (0.38–1.45) | 0.4 | 0.88 (0.43–1.81) | 0.7 | |
ERBB2 amplification | Negative | 159 (85.5%) | 1 |
| 1 |
|
Positive | 27 (14.5%) | 1.39 (0.70–2.78) | 0.3 | 1.45 (0.70–2.97) | 0.3 | |
mtDNA content | ≤ 350 | 48 (25.8%) | 1 |
| 1 |
|
> 350 | 138 (74.2%) | 0.50 (0.29–0.87) | 0.015 | 0.54 (0.30–0.97) | 0.038 | |
Number of patients and corresponding hazard ratio for distant metastasis-free survival with its 95% confidence intervals (CI) and corresponding probabilities for equal risk (P value) for each group. Analyses were stratified for estrogen receptor status and limited to the 186 patients (54 events) with no missing values.
DISCUSSION
Many contradictions about the physiological consequences of reduced mtDNA content exist in the literature. A critical reduction in mtDNA content compromises mitochondrial functioning with downstream effects. Subsequent changes in cellular processes such as aerobic respiration, calcium homeostasis or the intrinsic apoptotic pathway could in turn impact tumorigenic properties. Previous findings have pointed towards a link between low mtDNA content and breast cancer aggressiveness but the exact association remains uncertain. Here, to elucidate its potential as a prognostic marker, the putative relation between tumor mtDNA content and mesenchymal features or distant metastasis-free survival in breast cancer was explored. Using a quantitative PCR approach, the mtDNA content of 41 breast cancer cell lines and 204 primary breast tumor specimens was obtained. A correction for copy number variations of the nuclear-encoded reference locus (HMBS) minimized bias due to tumor-related genomic aberrations in the obtained number of mtDNA molecules per cell. Furthermore, the quantitative mtDNA target in our current assay lies outside of the common deletion region [41] and therefore it is likely that we measure a mixture of functional and dysfunctional mtDNA molecules.
Previous in vitro studies reported induction of EMT and stem-cell features upon depletion of mtDNA [12, 14, 19]. In the panel of breast cancer cell lines – homologous cell populations – no relation between mtDNA content and expression levels of genes involved in EMT could be demonstrated. In the cohort of primary breast tumor specimens – more heterogeneous cell populations – we find a positive but weak correlation (ρ ≤ 0.25) between mtDNA content and ESRP1, SNAI1 and TGFB1. In in vitro studies, a reduction in mtDNA content resulted in decreased ESRP1 protein levels but, contradictory to our findings, increased SNAI1 mRNA expression or TGFB protein expression. In addition, we could not demonstrate a difference in mtDNA content between the intrinsic subtypes in our cell line panel nor in the cohort of primary breast tumor specimens. Mesenchymal properties have been attributed to the basal and normal-like subtypes, whereas the luminal subtypes are generally epithelial [36, 38, 39]. Accordingly, the evaluated EMT-related genes were commonly highly statistically significant related to each other and differentially expressed between the intrinsic subtypes within our cell line panel and the cohort of primary breast tumor specimens (Supplementary Table S1). Apart from a true lack of association between mtDNA and EMT-features in breast cancer, there are several other reasons which may explain the absence of this association in our data set. To understand the physiological effects of mtDNA content, previous studies suggesting a relation between mtDNA and EMT often used cell lines artificially depleted of mtDNA, termed rho0 clones [42]. The endogenous mtDNA content of the cell lines and primary tumor specimens in our study is a few hundred molecules per cell, which is still orders of magnitude higher than of the rho0 clones. Since the extent of mtDNA reduction is of importance in gaining tumorigenic properties, as demonstrated in glioblastoma models [43], perhaps the mtDNA content in our data set is not at the critically low level necessary to induce a transition towards a mesenchymal phenotype. Alternatively, low mtDNA levels may be important during the process of EMT, but might be restored to normal after the transition is accomplished. Despite the unknown exact reason, we conclude that in our data set mtDNA content is not related to the molecular features connected to a mesenchymal-like phenotype.
To address a possible relation between mtDNA content and aggressive behavior in vivo, we analyzed primary breast tumor mtDNA content and prognosis in a cohort of 204 breast cancer patients. Notably, the mtDNA content in our primary tumor specimens is only an estimate, representing not only a heterogeneous tumor cell population but also non-neoplastic cells incorporated in the tumor specimen. However, because no evidence for a relation between mtDNA content and tumor infiltrating lymphocytes [44] was observed (Supplementary Table S2) and stromal content was minimized (Materials and Methods), we estimate the contribution of non-neoplastic cells to the final mtDNA content to be minimal. A few associations between mtDNA content and clinicopathological variables have been reported in previous studies, albeit never consistently [24–29]. These studies included either a low number of study participants or heterogeneous groups regarding treatment regimen or disease stage, making interpretation difficult. In this study, we included a population of lymph node-negative primary breast cancer patients who did not receive any (neo)adjuvant systemic treatment. In this patient group, lower mtDNA content was observed in tumors smaller than 2 cm across compared to tumors larger than 2 cm across. Previous studies could not demonstrate such a difference [27, 28] or reported lower mtDNA content in tumors over 5 cm across compared to smaller tumors [26]. Larger tumors presumably underwent more cell divisions potentially resulting in additional replication-induced mtDNA damage [45], which in turn might require additional compensatory mtDNA molecules to maintain proper mitochondrial functioning. It is also plausible that a hypoxic environment in larger tumors reduces mtDNA content as suggested previously [26]. However, in our primary tumor cohort no relation was observed between mtDNA content and hypoxia-related gene expression [46] as surrogate for the hypoxic state of the tumor (Supplementary Table S2). In addition, our results show a relation between mtDNA content and GGI, a gene expression-based identifier of the histological grade of tumors [40], with the highest grade representing poorly differentiated tumors showing higher mtDNA content. However, we could not demonstrate a conventional significant trend between mtDNA content and GGI. The relation between mtDNA content and histological grade has been reported before [26]. An increase in mtDNA content occurs in early S-phase of the cell cycle [47], and we attribute the relation between grade 3 tumors and higher mtDNA content to the high-proliferative nature of these higher grade tumors. Nevertheless, we note that the median difference in mtDNA content for both tumor size and GGI is only 20%, making a substantial biological consequence of these associations less likely.
Importantly, our cohort is highly suitable to study the prognostic value of mtDNA content for distant metastasis-free survival because all included patients presented with lymph node-negative disease and did not receive any (neo)adjuvant systemic treatment. The size of our cohort did not allow for separate analyzes for estrogen receptor-negative and -positive tumors, which show different proportionality over time (test of proportional hazards assumption P = 0.016). Therefore, stratification for estrogen receptor status was applied in all proportional hazard analyses. After adjustment for established prognostic clinicopathological variables, we observed a prognostic effect for mtDNA content. The patients with the 25% lowest mtDNA content (≤ 350 mtDNA molecules per cell) showed a significant unfavorable prognosis with shorter time to metastasis compared to patients with higher mtDNA content. One previous study reported on low tumor mtDNA content corresponding to a higher risk of death [26]. However in that study, no clear information was provided about treatments administered, disease stage at diagnosis and other clinical variables included in their statistical analysis. This renders interpretation and comparison with that previous study difficult. Interestingly, low mtDNA content predicted for a favorable response to anthracycline treatment in a small patient cohort [24]. It is plausible that cells with low mtDNA content are susceptible to such regimen, because damage of mtDNA in cells containing fewer mtDNA molecules can affect mitochondrial functionality more effectively. In our cohort we could study the prognostic value of mtDNA content independent of treatment regimen.
To conclude, we demonstrate a link between particularly low mtDNA content and metastatic potential in breast cancer, which appears to be independent of the mesenchymal phenotype. Low primary tumor mtDNA content potentially identifies patients with unfavorable prognosis but at the same time might predict therapeutic efficacy of DNA-damaging treatment regimen in this group. Larger cohorts of uniformly treated patients are necessary to validate these results and to further unravel the clinical relevance of mtDNA content determination in cancer.
MATERIALS AND METHODS
Study cohort and sampling
We employed a panel of 42 breast cancer cell lines (including BT20, BT474, BT483, BT549, CAMA1, DU4475, EVSAT, HCC1937, Hs578T, MCF7, MDAMB134VI, MDAMB157, MDAMB175VII, MDAMB231, MDAMB330, MDAMB361, MDAMB415, MDAMB435s, MDAMB436, MDAMB453, MDAMB468, MPE600, OCUBF, OCUBM, SKBR3, SKBR5, SKBR7, SUM102PT, SUM1315MO2, SUM149PT, SUM159PT, SUM185PE, SUM190PT, SUM225CWN, SUM229PE, SUM44PE, SUM52PE, T47D, UACC812, UACC893, ZR751 and ZR7530 [36, 44]). In addition, DNA extracts from fresh frozen primary breast tumor specimens from an earlier study [37] were selected from our bio-bank at the Erasmus MC. The study was approved by the medical ethics committee of the Erasmus MC (MEC 02.953) and conducted in accordance to the code of conduct of Federation of Medical Scientific Societies in the Netherlands. Whenever possible, we adhered to the Reporting Recommendations for Tumor Marker Prognostic Studies (REMARK) [48]. Patient selection criteria have been described before [49] and include lymph node-negative primary breast cancer with local treatment but no systemic (neo)adjuvant therapies. Our selection (Supplementary Figure S2) was based on availability of genotypic data and gene expression data from the primary tumors (n = 337) and availability of uniformly extracted DNA (see below) (n = 250). Next, specimens with a tumor cell percentage below 50% were excluded to minimize skewed values due to stromal cell contamination (n = 38). In addition, five patient samples were ineligible in retrospect and excluded. Thus, mtDNA content was examined in a total of 207 patients. Patients’ follow-up involved examinations every 3 months for the first two 2 years, every 6 months for years 3–5, and every 12 months from year 5 onwards. Estrogen receptor and progesterone receptor status were determined as described before [50]. Evaluation of ERBB2 amplification via RNA expression levels and the qRT-PCR Genomic Grade Index were determined as described before respectively [51] and [40].
DNA extraction
DNA was extracted from cultured cell lines using the DNeasy Blood & Tissue kit (Qiagen, Venlo, Netherlands) according the suppliers’ protocol. We selected the DNA which was previously extracted from cryostat sections of the primary tumor tissues based on uniformity in extraction procedure (QIAamp DNA mini kit (Qiagen) as described before [37]). DNA extracts were quantified using the Qubit dsDNA HS assay kit (Life Technologies, Carlsbad, United States of America) and all samples were diluted to a concentration of 0.2 ng/μL DNA prior to mtDNA content analysis.
Copy number analysis
Copy number variation of the nuclear encoded HMBS gene – which served as a reference to obtain mtDNA content – was obtained from our previously described microarray data (Gene Expression Omnibus database accession numbers GSE10099 [37] and GSE41308 [52]). The breast cancer cell lines were genotyped on the Genome Wide Human SNP Array 6.0 (Affymetrix, Santa Clara, United States of America), the primary tumor specimens on the GeneChip Human Mapping 100K SNP Array (Affymetrix).
mtDNA content
Mitochondrial DNA content was determined in duplicate runs using a multiplex quantitative PCR targeting the nuclear HMBS gene (chr11q23.2-qter) and the mitochondrial MT-TL1 (chrMT 3212–3319). Primers targeting the nuclear encoded HMBS gene (forward 5′-TGAGGCGGATGCAGATAC-3′ and reverse 5′- CCCACCCACGGTAGTAATTC-3′ (Life technologies)) yielded a 201 bp amplicon quantitatively detected using a CY5 labeled probe (5′-[CY5]TATCAGCCAAGCCTCCGAAC[BHQ2]-3′ (Sigma Aldrich, St. Louis, United States of America)). Primers targeting the mitochondrial encoded MT-TL1 (forward 5′-CACCCAAGAACAGGGTTTGT-3′ and reverse 5′-TGGCCATGGGTATGTTGTTA-3′ (Life Technologies)) yielded a 108 bp amplicon quantitatively detected using a HEX labeled probe (5′-[HEX]TTACCGGGCTCTGCCATCT[BHQ1]-3′, (Sigma Aldrich)) [53]. Reactions included 1x Absolute QPCR Mix containing SYBR Green and ROX (AB-1163 Life Technologies) in the presence of 100 nM mtDNA primers, 360 nM nDNA primers and 100 nM probes. The 45-cycle PCR was carried out at a 62°C annealing temperature and probe fluorescence was monitored using ROX, HEX, CY5 and FAM filters on Mx3000P or Mx3005P qPCR systems (Agilent Technologies, Waldbronn, Germany). Quantification cycle values (Cq [dRN]) were obtained using the adaptive baseline approach (MxPro v4.10) up to cycle 35 with fixed fluorescence thresholds at 0.004 dRn. Performance of singleplex PCR and multiplex PCR runs was comparable (Supplementary Table S3). Performance of the assay at variable ratios of artificial HMBS (289 bp linear: ACA GAC GGG GTC CTT TCA TTC GAG GCT GGG CTG AGG CGG ATG CAG ATA CGG CCC CTT TGG GAA GAC ACG TTC CAC TTT TGA TTC ATA GGA GAG AGT ATC AGC CAA GCC TCC GAA CTG CAC ACA AAC GTC TTA GAA GTG CGC CTT CTT TTT GTG TTA TAG TGG TCT CCC AGC CAC AGC CAA CGC TCC AAG TCC CCA GCT GTG ACA CAC CTA CTG AAT TAC TAC CGT GGG TGG GAG GCC GCC GTG GGC CTT TCC ATT ACG AGC CTG CTT GCC GAG CCC TGG GCT TGT GCA C) and artificial MT-TL1 (180 bp cloned in circular 2374 bp pMA-T vector: TAT CAT CTC AAC TTA GTA TTA TAC CCA CAC CCA CCC AAG AAC AGG GTT TGT TAA GAT GGC AGA GCC CGG TAA TCG CAT AAA ACT TAA AAC TTT ACA GTC AGA GGT TCA ATT CCT CTT CTT AAC AAC ATA CCC ATG GCC AAC CTC CTA CTC CTC ATT GTA CCC ATT CTA ATC GCA ATG GCA) was linear (Supplementary Table S4). DNA input of the breast cancer cell lines and the primary breast tumor specimens was standardized for 1 ng DNA per reaction. A calibration curve containing a pool of DNA isolates from independent fresh frozen tumors was taken along as internal control to monitor inter assay variation. Obtained Cq values were used to calculate the ratio of mitochondrial DNA opposed to nuclear DNA by the relative quantitation method (2^ΔCq [54]). Multiplying this ratio by the copy number of HMBS (obtained as described above) resulted in the number of mtDNA molecules per cell as mtDNA content.
Gene expression analysis
Gene expression data of the cell lines was obtained from our previously described triplicate microarray data (Gene Expression Omnibus database accession number GSE41313 [52]) on the Human Genome HT_HG-U133_Plus_PM GeneChip 96-well arrays (Affymetrix). Data of all breast cancer cell lines were available with the exception of SUM225CWN. Gene expression data of the primary breast tumor specimens was obtained from our previously described microarray data (Gene Expression Omnibus database accession number GSE2034 [50] and GSE5327 [55]) on the Human Genome HG-U133a GeneChip 96-well arrays (Affymetrix). Subtype classification was based on expression of the intrinsic gene set defined by Perou et al. [56]. Cell line DU4475 could not be classified to a subtype group and was therefore excluded from the intrinsic subtype analysis. For individual genes, levels based on log2 transformed distances to the geometric mean for each probe set were obtained for probe IDs 201131_s_at (CDH1), 203440_at (CDH2), 219121_s_at (ESRP1), 210495_x_at (FN1), 203936_s_at (MMP9), 219480_at (SNAI1), 213139_at (SNAI2), 203085_s_at (TGFB1), 213943_at (TWIST1) and 201426_s_at (VIM). The tumor infiltrating lymphocyte classification as low TIL and high TIL was based on the immune signature probe-set by Massink et al. [44]. Classification as low hypoxia response and high hypoxia response was based on the expression of the hypoxia-related gene signature as described before by Chi et al. [57].
Statistical analyses
All analyses included the average mtDNA content obtained from the duplicate analysis for each individual sample. Data distribution was tested using Skewness-Kurtosis tests for normality. Numerical correlations between RNA expression levels and mtDNA content were investigated using the Spearman rank correlation and corrected for multiple testing using the false discovery rate controlling procedure [58]. Categorical comparisons of the intrinsic subtypes or grouped clinical variables and mtDNA content were employed using either Mann-Whitney U-tests (two groups) or Kruskal-Wallis one-way analysis of variance (multiple groups). When appropriate, we performed Cuzick’s test for trend across ordered categorical variables. Kaplan-Meier survival plots and log-rank tests were used to assess the differences in time to distant metastasis between mtDNA content groups. Proportional hazard analyses for distant metastasis-free survival were performed using Cox proportional-hazards regression methods. We stratified for estrogen receptor status and censored for 10 years clinical follow-up (most patients are redirected to their general practitioner at that point in time) to maintain proportionality (test of proportional hazards assumption using the Schoenfeld residuals P > 0.05). Univariate analysis was done on the individual clinicopathological variables, multivariable analysis included all clinicopathological variables and mtDNA content. All statistical tests were two-sided, and P values smaller than 0.05 were considered as statistically significant. Clinical variables were statistically analyzed in Stata version 13.1 (StataCorp LP, College Station, United States of America). Other analyses were performed using Spotfire 7.0.0 (TIBCO, Palo Alto, United States of America).
ACKNOWLEDGMENTS
The authors thank Antoinette Hollestelle for providing the DNA extracts from cultured cell lines.
CONFLICTS OF INTEREST
Authors declare that they have no conflicts of interest.
FUNDING
This work was funded by a research grant from Philips Research.
REFERENCES
1. Wallace DC. Mitochondria and cancer. Nat Rev Cancer. 2012; 12:685–98.
2. A Dickinson KYY, Donoghue J, Baker MJ, Kelly RDW, McKenzie M, Johns TG, St. John JC. The regulation of mitochondrial DNA copy number in glioblastoma cells. Cell Death and Differentiation. 2013; 20:1644–1653.
3. Facucho-Oliveira JM, Alderson J, Spikings EC, Egginton S, St John JC. Mitochondrial DNA replication during differentiation of murine embryonic stem cells. J Cell Sci. 2007; 120:4025–34.
4. Gahan ME, Miller F, Lewin SR, Cherry CL, Hoy JF, Mijch A, Rosenfeldt F, Wesselingh SL. Quantification of mitochondrial DNA in peripheral blood mononuclear cells and subcutaneous fat using real-time polymerase chain reaction. J Clin Virol. 2001; 22:241–247.
5. Miller FJ, Rosenfeldt FL, Zhang CF, Linnane AW, Nagley P. Precise determination of mitochondrial DNA copy number in human skeletal and cardiac muscle by a PCR-based assay: lack of change of copy number with age. Nucleic Acids Res. 2003; 31.
6. Mansouri A, Gaou I, De Kerguenec C, Amsellem S, Haouzi D, Berson A, Moreau A, Feldmann G, Letteron P, Pessayre D, Fromenty B. An alcoholic binge causes massive degradation of hepatic mitochondrial DNA in mice. Gastroenterology. 1999; 117:181–90.
7. Casula M, Bosboom-Dobbelaer I, Smolders K, Otto S, Bakker M, de Baar MP, Reiss P, de Ronde A. Infection with HIV-1 induces a decrease in mtDNA. J Infect Dis. 2005; 191:1468–71.
8. Moraes CT, Shanske S, Tritschler HJ, Aprille JR, Andreetta F, Bonilla E, Schon EA, DiMauro S. mtDNA depletion with variable tissue expression: a novel genetic abnormality in mitochondrial diseases. Am J Hum Genet. 1991; 48:492–501.
9. El-Hattab AW, Scaglia F. Mitochondrial DNA depletion syndromes: review and updates of genetic basis, manifestations, and therapeutic options. Neurotherapeutics. 2013; 10:186–98.
10. Kulawiec M, Safina A, Desouki MM, Still I, Matsui S, Bakin A, Singh KK. Tumorigenic transformation of human breast epithelial cells induced by mitochondrial DNA depletion. Cancer Biol Ther. 2008; 7:1732–43.
11. Singh KK, Kulawiec M, Still I, Desouki MM, Geradts J, Matsui S. Inter-genomic cross talk between mitochondria and the nucleus plays an important role in tumorigenesis. Gene. 2005; 354:140–146.
12. Naito A, Cook CC, Mizumachi T, Wang M, Xie CH, Evans TT, Kelly T, Higuchi M. Progressive tumor features accompany epithelial-mesenchymal transition induced in mitochondrial DNA-depleted cells. Cancer Sci. 2008; 99:1584–8.
13. Moro L, Arbini AA, Yao JL, di Sant’Agnese PA, Marra E, Greco M. Mitochondrial DNA depletion in prostate epithelial cells promotes anoikis resistance and invasion through activation of PI3K/Akt2. Cell Death Differ. 2009; 16:571–583.
14. Guha M, Srinivasan S, Ruthel G, Kashina AK, Carstens RP, Mendoza A, Khanna C, Van Winkle T, Avadhani NG. Mitochondrial retrograde signaling induces epithelial-mesenchymal transition and generates breast cancer stem cells. Oncogene. 2014; 33:5238–50.
15. Pelicano H, Xu RH, Du M, Feng L, Sasaki R, Carew JS, Hu YM, Ramdas L, Hu LM, Keating MJ, Zhang W, Plunkett W, Huang P. Mitochondrial respiration defects in cancer cells cause activation of Akt survival pathway through a redox-mediated mechanism. J Cell Biol. 2006; 175:913–923.
16. Yu M, Shi Y, Wei X, Yang Y, Zhou Y, Hao X, Zhang N, Niu R. Depletion of mitochondrial DNA by ethidium bromide treatment inhibits the proliferation and tumorigenesis of T47D human breast cancer cells. Toxicol Lett. 2007; 170:83–93.
17. Cavalli LR, VarellaGarcia M, Liang BC. Diminished tumorigenic phenotype after depletion of mitochondrial DNA. Cell Growth Differ. 1997; 8:1189–1198.
18. Morais R, Zinkewich-Peotti K, Parent M, Wang H, Babai F, Zollinger M. Tumor-forming ability in athymic nude mice of human cell lines devoid of mitochondrial DNA. Cancer Res. 1994; 54:3889–96.
19. Magda D, Lecane P, Prescott J, Thiemann P, Ma X, Dranchak PK, Toleno DM, Ramaswamy K, Siegmund KD, Hacia JG. mtDNA depletion confers specific gene expression profiles in human cells grown in culture and in xenograft. BMC Genomics. 2008; 9:521.
20. Tan AS, Baty JW, Dong LF, Bezawork-Geleta A, Endaya B, Goodwin J, Bajzikova M, Kovarova J, Peterka M, Yan B, Pesdar EA, Sobol M, Filimonenko A, et al. Mitochondrial Genome Acquisition Restores Respiratory Function and Tumorigenic Potential of Cancer Cells without Mitochondrial DNA. Cell Metab. 2015; 21:81–94.
21. Guo JH, Zheng L, Liu WY, Wang XS, Wang ZM, Wang ZH, French AJ, Kang DC, Chen L, Thibodeau SN, Liu WG. Frequent Truncating Mutation of TFAM Induces Mitochondrial DNA Depletion and Apoptotic Resistance in Microsatellite-Unstable Colorectal Cancer. Cancer Res. 2011; 71:2978–2987.
22. Yu M. Generation, function and diagnostic value of mitochondrial DNA copy number alterations in human cancers. Life Sci. 2011; 89:65–71.
23. Cook CC, Higuchi M. The awakening of an advanced malignant cancer: an insult to the mitochondrial genome. Biochim Biophys Acta. 2012; 1820:652–62.
24. Hsu CW, Yin PH, Lee HC, Chi CW, Tseng LM. Mitochondrial DNA Content as a Potential Marker to Predict Response to Anthracycline in Breast Cancer Patients. Breast J. 2010; 16:264–270.
25. Yu M, Zhou YL, Shi YR, Ning LS, Yang Y, Wei XY, Zhang N, Hao XS, Niu RF. Reduced mitochondrial DNA copy number is correlated with tumor progression and prognosis in Chinese breast cancer patients. IUBMB Life. 2007; 59:450–457.
26. Bai RK, Chang J, Yeh KT, Lou MA, Lu JF, Tan DJ, Liu H, Wong LJ. Mitochondrial DNA content varies with pathological characteristics of breast cancer. J Oncol. 2011; 2011:496189.
27. Tseng LM, Yin PH, Chi CW, Hsu CY, Wu CW, Lee LM, Wei YH, Lee HC. Mitochondrial DNA mutations and mitochondrial DNA depletion in breast cancer. Genes Chromosomes Cancer. 2006; 45:629–38.
28. Fan AX, Radpour R, Haghighi MM, Kohler C, Xia P, Hahn S, Holzgreve W, Zhong XY. Mitochondrial DNA content in paired normal and cancerous breast tissue samples from patients with breast cancer. J Cancer Res Clin Oncol. 2009; 135:983–9.
29. Barekati Z, Radpour R, Kohler C, Zhang B, Toniolo P, Lenner P, Lv Q, Zheng H, Zhong XY. Methylation profile of TP53 regulatory pathway and mtDNA alterations in breast cancer patients lacking TP53 mutations. Hum Mol Genet. 2010; 19:2936–46.
30. Mambo E, Chatterjee A, Xing M, Tallini G, Haugen BR, Yeung SC, Sukumar S, Sidransky D. Tumor-specific changes in mtDNA content in human cancer. Int J Cancer. 2005; 116:920–4.
31. McMahon S, LaFramboise T. Mutational patterns in the breast cancer mitochondrial genome, with clinical correlates. Carcinogenesis. 2014; 35:1046–54.
32. Thiery JP. Epithelial-mesenchymal transitions in tumour progression. Nat Rev Cancer. 2002; 2:442–54.
33. Naito A, Carcel-Trullols J, Xie CH, Evans TT, Mizumachi T, Higuchi M. Induction of acquired resistance to antiestrogen by reversible mitochondrial DNA depletion in breast cancer cell line. Int J Cancer. 2008; 122:1506–11.
34. Yu M, Shi Y, Wei X, Yang Y, Zang F, Niu R. Mitochondrial DNA depletion promotes impaired oxidative status and adaptive resistance to apoptosis in T47D breast cancer cells. Eur J Cancer Prev. 2009; 18:445–57.
35. Sorlie T, Perou CM, Tibshirani R, Aas T, Geisler S, Johnsen H, Hastie T, Eisen MB, van de Rijn M, Jeffrey SS, Thorsen T, Quist H, Matese JC, et al. Gene expression patterns of breast carcinomas distinguish tumor subclasses with clinical implications. Proc Natl Acad Sci U S A. 2001; 98:10869–74.
36. Hollestelle A, Peeters JK, Smid M, Timmermans M, Verhoog LC, Westenend PJ, Heine AA, Chan A, Sieuwerts AM, Wiemer EA, Klijn JG, van der Spek PJ, Foekens JA, et al. Loss of E-cadherin is not a necessity for epithelial to mesenchymal transition in human breast cancer. Breast Cancer Res Treat. 2013; 138:47–57.
37. Zhang Y, Martens JWM, Yu JX, Jiang J, Sieuwerts AM, Smid M, Klijn JGM, Wang YX, Foekens JA. Copy Number Alterations that Predict Metastatic Capability of Human Breast Cancer. Cancer Res. 2009; 69:3795–3801.
38. Prat A, Parker JS, Karginova O, Fan C, Livasy C, Herschkowitz JI, He X, Perou CM. Phenotypic and molecular characterization of the claudin-low intrinsic subtype of breast cancer. Breast Cancer Res. 2010; 12:R68.
39. Sarrio D, Rodriguez-Pinilla SM, Hardisson D, Cano A, Moreno-Bueno G, Palacios J. Epithelial-mesenchymal transition in breast cancer relates to the basal-like phenotype. Cancer Res. 2008; 68:989–97.
40. Toussaint J, Sieuwerts AM, Haibe-Kains B, Desmedt C, Rouas G, Harris AL, Larsimont D, Piccart M, Foekens JA, Durbecq V, Sotiriou C. Improvement of the clinical applicability of the Genomic Grade Index through a qRT-PCR test performed on frozen and formalin-fixed paraffin-embedded tissues. BMC Genomics. 2009; 10.
41. Cortopassi GA, Arnheim N. Detection of a Specific Mitochondrial-DNA Deletion in Tissues of Older Humans. Nucleic Acids Res. 1990; 18:6927–6933.
42. King MP, Attardi G. Human cells lacking mtDNA: repopulation with exogenous mitochondria by complementation. Science. 1989; 246:500–3.
43. Dickinson A, Yeung KY, Donoghue J, Baker MJ, Kelly RD, McKenzie M, Johns TG, St John JC. The regulation of mitochondrial DNA copy number in glioblastoma cells. Cell Death Differ. 2013; 20:1644–53.
44. Massink MP, Kooi IE, Martens JW, Waisfisz Q, Meijers-Heijboer H. Genomic profiling of CHEK2*1100delC-mutated breast carcinomas. BMC Cancer. 2015; 15:877.
45. Ju YS, Futreal CA, Gerstung M, Martincorena I, Nik-Zainal S, Ramakrishna M, Davies HR, Papaemmanuil E, Gundem G, Shlien A, Bolli N, Behjati S, Tarpey PS, et al. Origins and functional consequences of somatic mitochondrial DNA mutations in human cancer. Elife. 2014; 3.
46. Chi JT, Wang Z, Nuyten DS, Rodriguez EH, Schaner ME, Salim A, Wang Y, Kristensen GB, Helland A, Borresen-Dale AL, Giaccia A, Longaker MT, Hastie T, et al. Gene expression programs in response to hypoxia: cell type specificity and prognostic significance in human cancers. PLoS Med. 2006; 3:e47.
47. Trinei M, Berniakovich I, Pelicci PG, Giorgio M. Mitochondrial DNA copy number is regulated by cellular proliferation: A role for Ras and p66 (Shc). Biochim Biophys Acta. 2006; 1757:624–630.
48. McShane LM, Altman DG, Sauerbrei W, Taube SE, Gion M, Clark GM. Statistics Subcommittee of the NCIEWGoCD. Reporting recommendations for tumor marker prognostic studies (REMARK). J Natl Cancer Inst. 2005; 97:1180–4.
49. Smid M, Wang YX, Zhang Y, Sieuwerts AM, Yu J, Klijn JGM, Foekens JA, Martens JWM. Subtypes of breast cancer show preferential site of relapse. Cancer Res. 2008; 68:3108–3114.
50. Wang Y, Klijn JG, Zhang Y, Sieuwerts AM, Look MP, Yang F, Talantov D, Timmermans M, Meijer-van Gelder ME, Yu J, Jatkoe T, Berns EM, Atkins D, et al. Gene-expression profiles to predict distant metastasis of lymph-node-negative primary breast cancer. Lancet. 2005; 365:671–9.
51. van Agthoven T, Sieuwerts AM, Meijer-van Gelder ME, Look MP, Smid M, Veldscholte J, Sleijfer S, Foekens JA, Dorssers LC. Relevance of breast cancer antiestrogen resistance genes in human breast cancer progression and tamoxifen resistance. J Clin Oncol. 2009; 27:542–9.
52. Riaz M, van Jaarsveld MT, Hollestelle A, Prager-van der Smissen WJ, Heine AA, Boersma AW, Liu J, Helmijr J, Ozturk B, Smid M, Wiemer EA, Foekens JA, Martens JW. miRNA expression profiling of 51 human breast cancer cell lines reveals subtype and driver mutation-specific miRNAs. Breast Cancer Res. 2013; 15:R33.
53. Bai RK, Wong LJC. Simultaneous detection and quantification of mitochondrial DNA deletion(s), depletion, and over-replication in patients with mitochondrial disease. J Mol Diagn. 2005; 7:613–622.
54. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(T)(-Delta Delta C) method. Methods. 2001; 25:402–408.
55. Yu JX, Sieuwerts AM, Zhang Y, Martens JW, Smid M, Klijn JG, Wang Y, Foekens JA. Pathway analysis of gene signatures predicting metastasis of node-negative primary breast cancer. BMC Cancer. 2007; 7:182.
56. Perou CM, Sorlie T, Eisen MB, van de Rijn M, Jeffrey SS, Rees CA, Pollack JR, Ross DT, Johnsen H, Akslen LA, Fluge O, Pergamenschikov A, Williams C, et al. Molecular portraits of human breast tumours. Nature. 2000; 406:747–52.
57. Chi JT, Wang Z, Nuyten DS, Rodriguez EH, Schaner ME, Salim A, Wang Y, Kristensen GB, Helland A, Borresen-Dale AL, Giaccia A, Longaker MT, Hastie T, et al. Gene expression programs in response to hypoxia: cell type specificity and prognostic significance in human cancers. PLoS Med. 2006; 3:e47.
58. Benjamini Y, Hochberg Y. Controlling the False Discovery Rate - a Practical and Powerful Approach to Multiple Testing. J R Stat Soc Series B Stat Methodol. 1995; 57:289–300.
