PARP1 expression in soft tissue sarcomas is a poor‐prognosis factor and a new potential therapeutic target

Soft tissue sarcomas (STSs) are aggressive tumors with few efficient systemic therapies. Poly(ADP‐ribose) polymerase‐1 (PARP1) inhibitors represent an emerging therapeutic option in tumors with genomic instability. The genomics of STSs is complex in more than half of cases, suggesting a high level of inherent DNA damage and genomic instability. Thus, STSs could be efficiently targeted with PARP inhibitors. Promising preclinical results have been reported, but few data are available regarding PARP1 expression in clinical samples. We examined PARP1 mRNA expression in 1464 clinical samples of STS, including 1432 primary tumors and 32 relapses, and searched for correlations with clinicopathological features, including metastasis‐free survival (MFS). Expression was heterogeneous across the samples, not different between primary and secondary tumors, and was correlated to PARP1 DNA copy number. In the 1432 primary tumors, the ‘PARP1‐high’ samples were associated with younger patients, more frequent locations at the extremities, superficial trunk and head and neck, more leiomyosarcomas and other STSs and less liposarcomas and myxofibrosarcomas, more grade 3, more high‐risk CINSARC tumors, and more ‘chromosomically instable’ tumors. They were associated with shorter MFS, independently of other significant prognostic features, including the CINSARC signature. We found a strong involvement of genes overexpressed in the ‘PARP1‐high’ samples in cell cycle, DNA replication, and DNA repair. PARP1 expression refines the prediction of MFS in STSs, and similar expression exists in secondary and primary tumors, supporting the development of PARP1 inhibitors.


Introduction
Soft tissue sarcomas (STSs) are rare, severe, and heterogeneous tumors including many different pathological subtypes (Casali et al., 2018). Surgery is the main treatment in early stages, but more than 40% of operated patients will ultimately experience metastatic relapse and die. The survival benefit of adjuvant anthracycline-based chemotherapy remains unproven, likely in part because of the absence of accurate prognostic features and predictors of response to chemotherapy. Identification of new prognostic features such as the promising CINSARC gene expression signature (Chibon et al., 2010) is warranted. In patients with metastatic disease not amenable to curative-intent surgery, the first-line systemic treatment involves palliative chemotherapy that has very little change over the three past decades and remains based on doxorubicin. After intolerance or failure, the second-line therapies include chemotherapies (ifosfamide, dacarbazine, trabectedin, eribulin) and targeted therapy (pazopanib), but the results remain disappointing. Clearly, the improvement of systemic therapies and identification of new prognostic and therapeutic targets are crucial.
An emerging therapeutic option in oncology concerns PARP inhibitors (Lim and Tan, 2017). Poly (ADP-ribose) polymerase-1 (PARP1) is a nuclear chromatin-associated protein involved in several biological processes including cell proliferation, apoptosis, malignant transformation, transcriptional regulation, and DNA repair. It is essential to the base excision repair of DNA single-strand breaks (SSBs). In response to DNA damage, PARP1 senses and binds to DNA nicks and breaks, resulting in activation of catalytic activity, causing poly(ADP)ribosylation of PARP1 itself, as well as other acceptor proteins such as histones and topoisomerases. This modification stimulates the recruitment and activity of other components of DNA repair pathways (Ame et al., 2004). In its absence, DNA SSBs accumulate and degenerate to DNA double-strand breaks, which are not appropriately repaired if the BRCA pathway is deficient or dysfunctional. This is thought to explain the exquisite sensitivity to PARP inhibitors of tumors with BRCA inactivation, a concept called 'synthetic lethality' (Bryant et al., 2005;Farmer et al., 2005). Today, different PARP inhibitors are marketed or in development in advanced solid tumors (Lim and Tan, 2017) with BRCAness features.
Although STSs do not have a characterized defect in BRCA1/2, their genomics is complex in more than half of the cases, suggesting genomic instability and eventual possible deficiency in DNA damage repair and a high level of inherent DNA damage, as recently reported for leiomyosarcomas (Chudasama et al., 2018). Thus, STSs could be efficiently targeted with PARP inhibitors to drive cells to synthetic lethality. In sarcomas, promising preclinical data have been reported, notably in Ewing sarcoma and in STSs (Chudasama et al., 2018;Laroche et al., 2017;Stewart et al., 2014;Vormoor and Curtin, 2014). Recently, the safety of the combination of trabectedin chemotherapy and olaparib PARP inhibitor in second-line or furtherline has been shown in patients with advanced sarcomas (Grignani et al., 2018), with a promising 18% partial response in patients with STS. The response rate and progression-free survival were higher in patients with high PARP1 tumor expression, confirming their preclinical findings (Pignochino et al., 2017).
However, very few data are available in the literature regarding PARP1 expression in clinical STS samples. To our knowledge, only two studies are available and concern only 91 malignant peripheral nerve sheath tumors (MPNSTs) (Kivlin et al., 2016) and 112 STSs (Kim et al., 2016). To fill this gap, we examined PARP1 mRNA expression in a series of 1464 clinical samples of STS, including 1432 primary tumors and 32 relapses, and searched for correlations with clinicopathological features, including metastasis-free survival (MFS).

Soft tissue sarcoma samples and data sets
We retrospectively gathered clinicopathological and gene expression data of clinical STS samples from 16 public data sets (Baird et al., 2005;Barretina et al., 2010;Beck et al., 2010;Cancer Genome Atlas Research Network, 2017;Chibon et al., 2010;Detwiller et al., 2005;Gibault et al., 2011;Gobble et al., 2011;Hajdu et al., 2010;Henderson et al., 2005;Nakayama et al., 2007;Nielsen et al., 2002;Renner et al., 2013;Skubitz et al., 2012;West et al., 2005;Ylipaa et al., 2011). Sets and raw data were collected from the National Center for Biotechnology Information (NCBI)/GenBank GEO and ArrayExpress databases and authors' web sites (Table S1). The selection of data sets was based according to the availability of clinical and expression data, including PARP1 expression measurement. Samples had been profiled using DNA microarrays or RNASeq. The pooled data set contained a total of 1464 clinical samples of primary STS, including 1432 primary STSs and 32 STS relapses. These relapse samples were included in order to compare the PARP1 mRNA expression level between primary tumors and relapse samples, since these later will be the first candidates to PARP inhibitors in their clinical development. We also collected PARP1 DNA copy number, DNA methylation, and DNA mutational data of 224 STS primary tumors profiled in the The Cancer Genome Atlas (TCGA) data set (Cancer Genome Atlas Research Network, 2017) using SNP array and whole-exome sequencing, respectively.

Gene expression data analysis
The pre-analytic processing first included normalization of each data set separately, by using robust multichip average (Irizarry et al., 2003) with the nonparametric quantile algorithm for the raw Affymetrix data and quantile normalization for the available processed non-Affymetrix microarray data. Normalization was done in R using Bioconductor and associated packages. Then, we mapped hybridization probes across the different technological platforms as reported (Bertucci et al., 2014). When multiple probes mapped to the same GeneID, we retained the one with the highest variance in each data set. We log 2 -transformed the already normalized TCGA RNAseq data. Next, the batch effects were corrected across the 16 studies using z-score normalization. Briefly, for each expression value in each study separately, all values were transformed by subtracting the mean of the gene in that data set divided by its standard deviation, mean and standard deviation (SD) being measured on leiomyosarcoma samples. We applied to each data set separately two gene expression signatures: the CIN-SARC signature (Chibon et al., 2010) and the Carter's chromosomal instability signature (Carter et al., 2006).
To decipher the biological pathways associated with PARP1 expression in STSs, we applied a supervised analysis to expression profiles of the 224 TCGA samples (learning set) to search for genes differentially expressed between the 'PARP1-high' vs 'PARP1-low' classes (cut-off defined as the median expression level across all samples). We used a moderated t-test with empirical Bayes statistic included in the limma R packages. False discovery rate (Hochberg and Benjamini, 1990) was applied to correct the multiple testing hypothesis: The significant genes were defined by p < 1%, q < 1%, and fold change superior to |1.59|. The robustness of the resulting gene list was tested in the validation set of 1208 remaining samples (592 'PARP1-low' samples and 616 'PARP1-high' samples) by computing for each tumor a metagene-based prediction score defined by the difference between the 'metagene PARP1-up' (mean expression of all genes upregulated in the 'PARP1-high' class) and the 'metagene PARP1-low' (mean expression of all genes upregulated in the 'PARP1-low' class). This score was then compared between the 'PARP1-high' and 'PARP1low' samples. Ontology analysis of the resulting gene list was based on Gene ontology (GO) biological processes of the Database for Annotation, Visualization and Integrated Discovery (DAVID; http://david.abcc. ncifcrf.gov/).

Statistical analysis
Correlations between the PARP1 expression-based classes (low vs high) and the clinicopathological factors were calculated with Student's t-test for the continuous variables and the Fisher's exact test for the binary variables. Our primary endpoint, MFS, was calculated from the date of diagnosis until the date of metastatic relapse. The follow-up was measured from the date of diagnosis to the date of last news for event-free patients. Survival was calculated using the Kaplan-Meier method, and curves were compared with the logrank test. Univariate and multivariate analyses were done using Cox regression analysis (Wald test). The variables tested in univariate analysis included the PARP1-based classification (low vs high), patients' age and gender, pathological type, grade, and tumor size, depth, tumor site, and the CINSARC risk (high vs low). Multivariate analysis incorporated all variables with a P-value inferior to 5% in univariate analysis. The likelihood ratio (LR) tests were used to assess the prognostic information provided beyond that of the CINSARC signature, assuming a chi-square distribution. Changes in the LR values (LR-ΔX 2 ) measured quantitatively the relative amount of information of one model compared with another. All statistical tests were two-sided at the 5% level of significance. Statistical analysis was done using the survival package (version 2.30) in the R software (version 2.15.2) (R Foundation for Statistical Computing, Vienna, Austria). The paper was written in accordance with reporting recommendations for tumor marker prognostic studies (REMARK) criteria (McShane et al., 2005). were available. Their characteristics are summarized in Table 1. The median patients' age was 63 (range, 2-93) years. The sex ratio was balanced, with 49% of females. The most frequent anatomical sites were extremities, followed by internal trunk; 84% of tumors were deeply seated, below or through the superficial fascia. As expected, the most frequent pathological types were liposarcomas, leiomyosarcomas, and undifferentiated sarcomas. The median pathological tumor size on the operative specimen was 9 cm. Most of the samples were F ed eration Nationale des Centres de Lutte Contre le Cancer (FNCLCC) grade 3, and 47% of samples were classified as high risk according to the CINSARC signature and 56% as chromosomically instable according to the Carter's signature.

PARP1 expression in soft tissue sarcomas
PARP1 mRNA expression was heterogeneous through the 1432 primary tumors with a three-decade range of values, and it was similar between primary and secondary tumors (Fig. 1A). We searched for correlation between PARP1 expression and DNA alterations at different levels (copy number, methylation, and mutation of PARP1 gene on chromosome 1) that were simultaneously annotated in the 224 TCGA samples. No sample showed PARP1 methylation or mutation. By contrast, 35 samples (16%) showed copy number alteration (CNA), including two deletions (homozygous losses: log 2 ratio>|0.5|), 22 heterozygous losses, nine gains, and two amplifications (log 2 ratio>|1|). There was a positive correlation between PARP1 mRNA expression and DNA copy number, with increasing expression from samples with deletion, then loss, then no CNA, then gain, and finally amplification (P = 2.82E-08, ANOVA; Fig. 1B).

PARP1 expression and clinicopathological characteristics
Within the 1432 primary tumors, and when compared with the 'PARP1-low' class, the 'PARP1-high' class was associated (Table 1) with younger patients' age (P = 2.18E-04, Student's t-test) and with (Fisher's exact test) more frequent tumor locations at the extremities, superficial trunk and head and neck and less internal trunk locations (P = 1.90E-02), pathological subtypes with more leiomyosarcomas and other STSs and less liposarcomas and myxofibrosarcomas (P = 2.39E-11), higher pathological grade 3 (P = 1.57E-03), high-risk CINSARC class (P = 3.62E-05), and higher Carter's signature-based chromosomal instability (P = 4.56E-24). There was no correlation with patients' gender, depth location, and pathological tumor size.

PARP1 expression and associated biological processes
To further explore the biological alterations associated with the PARP1 expression status, we compared the whole-genome expression profiles of the 'PARP1-high' (N = 89) and 'PARP1-low' (N = 135) primary tumor samples in the TCGA data set (Fig. 3A). We identified 530 genes differentially expressed, including 359 genes overexpressed and 171 genes underexpressed in the 'PARP1-high' class (Table S2). The robustness of this gene signature was confirmed in the pool of all other independent sets (1208 primary tumors) by using a metagene-based prediction score (Fig. 3B): The score was higher in the 'PARP1-high' samples than in the 'PARP1-low' samples (P = 2.56E-43, Student's t-test).
Ontology analysis (Table S3, Fig. 3C) showed strong involvement of genes overexpressed in the 'PARP1high' samples in cell cycle, chromosome segregation, DNA replication, and DNA repair.

Discussion
The need for new therapeutic and/or prognostic targets is crucial in STSs. Because of the promising therapeutic value of PARP inhibitors in oncology and the paucity of data in the literature, we analyzed PARP1 mRNA expression in 1432 previously untreated operated STS samples and 32 relapses. We showed that higher expression was an independent negative prognostic factor for MFS of patients with primary tumors. To our knowledge, this is by far the largest study analyzing PARP1 expression in STSs. PARP1 tumor expression was heterogeneous between samples. The analysis of the 224 TCGA primary tumor samples profiled at both RNA and DNA levels revealed a positive correlation between PARP1 mRNA expression and DNA copy number; however, gain/amplification was not the sole mechanism of high expression, which was also found in tumors without such alterations. Of note, expression was similar in primary and secondary tumors. No PARP1 mutation was reported in those 244 primary tumors, nor in the 1215

5-year
Log-rank test, P = 2.55E-11 metastatic samples of the GENIE AACR database (data not shown). The frequency of PARP1 amplification was also very low in the TCGA primary tumors (0.8%) and the GENIE metastatic samples (0.2%; data not shown). This wide range of expression values within the 1432 operated primary tumors provided opportunity to search for correlations with clinicopathological features. Our analysis was based on discrete values using the median PARP1 expression level across the 1432 samples as cut-off, but similar correlations were found with continuous values. An optimal expression cut-off able to stratify patients on MFS was measured by means of ROC analysis at 11.62 in the learning set (N = 340 samples), very close to our median cut-off measured at 11.64 (data not shown). It was validated in the validation set (N = 338 samples) with significant MFS difference between the 'PARP1-high' vs 'PARP1-low' classes (data not shown). The concordance rate between the two classifications in the whole series was very high, equal to 98.5%. Correlations existed between the two PARP1 classes ('high' and 'low') and patients' age, tumor site, pathological type and grade, the CINSARC classification and a chromosomal instability signature, with younger patients in the 'PARP1-high' class, more frequent tumor locations at the extremities, superficial trunk and head and neck, more leiomyosarcomas and other STSs and less liposarcomas and myxofibrosarcomas, more pathological grade 3, more high-risk CIN-SARC tumors, and more 'chromosomically instable' tumors. Such association with adverse prognostic features was confirmed in univariate analysis with shorter MFS in the 'PARP1-high' class. This negative prognostic value remained significant-suggesting independence -in multivariate analysis. It was notably independent from the pathological type. Interestingly, the prognostic value was observed in liposarcomas and undifferentiated sarcomas, but not in leiomyosarcomas. The fact that one marker has a prognostic value different according to the pathological type of STM is not surprising given the big intertype differences. In liposarcomas, the PARP1-based classification was associated with the pathological subtypes (well differentiated/dedifferentiated, myxoid, pleomorphic), but its prognostic value persisted in multivariate analysis including these later (data not shown). In leiomyosarcomas, CINSARC was associated with the PARP1-based classification and with MFS, whereas strikingly PARP1-based classification had no prognostic value. In undifferentiated sarcomas, none clinicopathological prognostic variable was associated with the PARP1-based classification that was itself associated with MFS. The analysis of larger series of samples per pathological type is required to better understand such differences. For several decades, efforts have been made to improve the prognostic classification of STSs; different tumor cell-intrinsic molecular parameters have been proposed, mainly related to cell cycle such as CINSARC (Chibon et al., 2010), as well as parameters related to immune microenvironment, such as PDL1/CD274 expression (Bertucci et al., 2017). The on the Akaike information criterion retained the CIN-SARC/PDL1/PARP1 combination as the best MFS predictive model (P = 5.99E-11; data not shown).

A B
Whether the prognostic value of PARP1 expression reflects the metastatic risk and/or the response to eventual adjuvant chemotherapy deserves analysis of a larger and informative series of patients. Indeed, information about delivery or not of adjuvant chemotherapy was available for only 374 out of 678 patients (29 with and 345 without chemotherapy). This is a bias of retrospective analyses. The analysis of MFS in the 345 chemotherapy-untreated cases showed a trend (P = 0.067, log-rank test) for longer MFS in the 'PARP1-low' class (71% 5-year MFS, 95% CI: 64-78) than in the 'PARP1-high' class (60% 5-year MFS, 95% CI: 50-73). No MFS analysis could be done in the 29 patients treated with adjuvant chemotherapy because of the small series size, thus impeding interaction analysis.
Considering the multiple functions of PARP1 including cell proliferation and DNA repair, it was not surprising to find high expression associated with shorter MFS, as already reported in other cancers (Goncalves et al., 2011;Li et al., 2016;Mahe et al., 2015;Michels et al., 2015;Murnyak et al., 2017). Regarding STSs, to our knowledge, only two studies have described PARP1 expression in clinical samples (Kim et al., 2016;Kivlin et al., 2016). All were based on IHC and tissue microarrays and used different antibodies, scoring systems and cut-offs. In the MDA Cancer Center study (Kivlin et al., 2016), 91 MPNSTs and 24 neurofibromas were analyzed: Overall, MPNST samples had higher levels of PARP1 expression than neurofibromas; no correlation with clinicopathological features was reported except with survival, which was nonsignificantly higher in MPNSTs with high vs low expression. The authors concluded that ubiquitous expression pattern of PARP1 supports the use of PARP inhibitors in MPNSTs, but that larger series of samples needed to be analyzed. In the Korean study (Kim et al., 2016), 112 STSs, representing 17 different pathological types, were tested. As found in our study, leiomyosarcomas, undifferentiated sarcomas, and other types were more frequently PARP1-positive than liposarcomas and myxofibrosarcomas; high expression was associated with higher pathological grade and higher mitotic count; and PARP1 expression was an independent negative prognostic feature regarding event-free survival and diseasespecific survival. Of course, the role of such an overexpression in STS initiation or progression, if any, remains to be elucidated. As expected, given the role of PARP1, many genes identified in our supervised analysis as overexpressed in the 'PARP1-high' samples were associated with cell proliferation, which could in part explain such poor prognostic value. However, the PARP1 expression prognostic value remained independent from the CINSARC signature, possibly reflecting the impact of a reaction against genetic instability and DNA damages. Indeed, several genes involved in DNA repair were overexpressed in the 'PARP1-high' samples, as previously reported in STSs (Kim et al., 2016). Such association of genes overexpressed in the 'PARP1-high' samples with ontologies representing known functions of PARP1 protein provides indication that increased PARP1 mRNA expression in STS is likely associated with increase in its biological activity and thus its protein expression. The correlation between mRNA and protein expression is also corroborated by the finding of similar clinicopathological correlations of PARP1 expression at the protein level in the Korean series (Kim et al., 2016) and the mRNA level in our present series.

Conclusion
We showed that PARP1 mRNA expression is heterogeneous in STS and associated with metastatic relapse independently from the other prognostic features, including the proliferation-associated CINSARC signature, the most robust prognostic signature reported to date in STSs. The strength of our study lies in the size of the series (the largest series of tumors reported to date regarding analysis of PARP1 expression), its originality (the first one describing PARP1 mRNA expression in STSs), and the biological and clinical relevance of PARP1 expression and its independent prognostic value. Limitations include its retrospective multicentric nature and associated biases such as lack of information about the time interval between imaging used during follow-up, and possible heterogeneity across patients, the heterogeneity with several different STS pathological types, and a limited number of cases in certain types. No overall survival analysis could be done because of the lack of information both quantitative and qualitative. The analysis of larger series, retrospective, then prospective, is warranted to confirm our observation and to assess each pathological type independently. If such prognostic value is confirmed, PARP1 expression might refine the prediction of metastatic relapse and improve our ability to tailor adjuvant chemotherapy. Given this unfavorable prognostic value, STS patients with high level of PARP1 expression would warrant a more aggressive treatment plan, which might include PARP1 inhibitors possibly associated with trabectedin given the predictive value of high PARP1 expression (Grignani et al., 2018;Pignochino et al., 2017) or with other DNA-damaging agents. We found similar expression level in secondary vs primary tumors. Even if analysis of larger series of metastatic samples is warranted, our results support the ongoing development of PARP inhibitors in STSs. In the future, it will be important not only to test whether PARP1 mRNA expression can predict the clinical response to PARP1 inhibitors or DNA-damaging agents, but also to validate our findings at the protein level using IHC that remains more convenient for use in clinical routine.

Supporting information
Additional supporting information may be found online in the Supporting Information section at the end of the article. Fig. S1. MFS in patients with different STS pathological types according to PARP1 expression. (A) Kaplan-Meier MFS curves in 256 patients with liposarcoma, according to the PARP1-based classification ('PARP1low' and 'PARP1-high' classes). (B) Similar to A, but in 202 patients with undifferentiated sarcoma. (C) Similar to A, but in 149 patients with leiomyosarcoma. Fig. S2. MFS in patients with STS according to PARP1 expression, CINSARC signature, and PDL1 expression. Kaplan-Meier MFS curves in 470 patients with STS, informative for the three variables: PARP1 expression (high and low), CINSARC signature (highrisk and low-risk), and PDL1 expression (high and low). The PDL1 legend and the colors in the table to the right of the figure define the eight patients groups. Table S1. List of soft tissue sarcoma data sets included. Table S2. List of 530 genes differentially expressed between the 'PARP1-high' and 'PARP1-low' sample classes. Table S3. Ontology analysis of the 530 genes differentially expressed between the 'PARP1-high' and 'PARP1-low' sample classes.