Antoine Adenis

J Bennouna

P L Etienne

E Bogart

E Francois

Etienne Pl E , Bogart

Galais Mp

Ben Abdelghani

Metges Jp

Le Deley Mc

Continuation versus discontinuation of first-line chemotherapy in patients with metastatic squamous cell oesophageal cancer: A randomised phase II trial (E-DIS)

Keywords: oesophageal cancer, chemotherapy, squamous-cell carcinoma, metastatic disease

niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Telephone: +33 467614701; e-mail: antoine.adenis@icm.unicancer.fr # deceased

Introduction

Oesophageal cancer ranks eighth on the list of the most common cancers worldwide, with an estimated 456,000 new cases in 2012 [START_REF] Ferlay | Cancer incidence and mortality worldwide: sources, methods and major patterns in GLOBOCAN 2012[END_REF]. Oesophageal squamous-cell carcinoma (ESCC) represents a majority of all cases of oesophageal cancer globally, and it is particularly prevalent in Eastern Asia, Eastern Africa and South America. ESCC is mainly related to the use of tobacco and alcohol [START_REF] Ferlay | Cancer incidence and mortality worldwide: sources, methods and major patterns in GLOBOCAN 2012[END_REF]. Overall, this cancer is a deadly disease with a mortality-toincidence ratio of 0.88 [START_REF] Ferlay | Cancer incidence and mortality worldwide: sources, methods and major patterns in GLOBOCAN 2012[END_REF]. Fifty per cent of these patients present with synchronous metastases at the time of diagnosis, and most patients who present initially with localised disease eventually develop metachronous metastases [START_REF] Enzinger | Chemotherapy in esophageal cancer[END_REF]. The 3-year survival rate for patients with metastatic oesophageal cancer is less than 1% [START_REF] Bouvier | Management and prognosis of esophageal cancers: has progress been made?[END_REF].

The role of chemotherapy has not yet been fully established in the treatment of metastatic ESCC (mESCC) [START_REF] Adenis | Palliative chemotherapy does not improve survival in metastatic esophageal cancer[END_REF]. Many cytotoxic drugs promote tumor shrinkage [START_REF] Ajani | Activity of taxol in patients with squamous cell carcinoma and adenocarcinoma of the esophagus[END_REF][START_REF] Conroy | Phase II trial of vinorelbine in metastatic squamous cell esophageal carcinoma[END_REF][START_REF] Bleiberg | Randomised phase II study of cisplatin and 5-fluorouracil (5-FU) versus cisplatin alone in advanced squamous cell oesophageal cancer[END_REF][START_REF] Ilson | Phase II trial of paclitaxel, fluorouracil, and cisplatin in patients with advanced carcinoma of the esophagus[END_REF][START_REF] Ilson | Phase II trial of weekly irinotecan plus cisplatin in advanced esophageal cancer[END_REF][START_REF] Conroy | Vinorelbine and cisplatin in metastatic squamous cell carcinoma of the oesophagus: response, toxicity, quality of life and survival[END_REF][START_REF] Lordick | Phase II trial of irinotecan plus docetaxel in cisplatin-pretreated relapsed or refractory oesophageal cancer[END_REF][START_REF] Muro | A phase II study of single-agent docetaxel in patients with metastatic esophageal cancer[END_REF][START_REF] Burkart | A phase II trial of weekly irinotecan in cisplatin-refractory esophageal cancer[END_REF][START_REF] Lee | A phase II study of capecitabine and cisplatin (XP) as first-line chemotherapy in patients with advanced esophageal squamous cell carcinoma[END_REF][START_REF] Wang | A phase II study of oxaliplatin in combination with leucovorin and fluorouracil as first-line chemotherapy in patients with metastatic squamous cell carcinoma of esophagus[END_REF][START_REF] Lee | Capecitabine in combination with either cisplatin or weekly paclitaxel as a first-line treatment for metastatic esophageal squamous cell carcinoma: a randomized phase II study[END_REF][START_REF] Wang | Efficacy and safety of recombinant human lymphotoxin-α derivative with cisplatin and fluorouracil in patients with metastatic esophageal squamous cell carcinoma: A randomized, multicenter, open-label, controlled, phase 2b trial[END_REF] but, until now, there has been no randomised trial that provides unequivocal evidence for a clinical benefit with chemotherapy. For example, the European Organisation for Research and Treatment of Cancer (EORTC) reported the results of a randomised phase 2 trial that compared the activity of cisplatin (CDDP) with or without 5-fluorouracil (FU) in 88 mESCC [START_REF] Bleiberg | Randomised phase II study of cisplatin and 5-fluorouracil (5-FU) versus cisplatin alone in advanced squamous cell oesophageal cancer[END_REF]. Although a greater response rate was observed with the combination, this was achieved without any overall survival (OS) benefit, and at the expense of substantially increased toxicity, including deaths due to toxicity. Recently, the combination of chemotherapy and epidermal growth factor receptor monoclonal antibodies failed to improve the efficacy over chemotherapy alone [START_REF] Lorenzen | Cetuximab plus cisplatin-5-fluorouracil versus cisplatin-5-fluorouracil alone in first-line metastatic squamous cell carcinoma of the esophagus: a randomized phase II study of the Arbeitsgemeinschaft Internistische Onkologie[END_REF][START_REF] Moehler | Cisplatin/5-FU (CF) +/-panitumumab (P) for patients (pts) with non-resectable, advanced, or metastatic esophageal squamous cell cancer (ESCC): An open-label randomized AIO/TTD/BDGO/EORTC phase III trial (POWER)[END_REF]. Guidelines from the European Society of Medical Oncology (ESMO) recommend supportive care or chemotherapy as options in such settings [START_REF] Lordick | Oesophageal cancer: ESMO Clinical Practice Guidelines for diagnosis, treatment and follow-up[END_REF]. Surprisingly, some physicians are convinced that mESCC patients may benefit from chemo. This attitude leads to the delivery of treatments with unproven benefit with side effects and unjustified costs.

The E-DIS trial (ClinicalTrials.gov identifier: NCT01248299) was a randomised discontinuation trial that was offered to mESCC patients who were free from progression after a 6-week selection-phase of FU/platinum-based chemotherapy. This study was designed to estimate the OS of mESCC patients who continued chemotherapy (CT-CONT arm).

Patients and Methods

Patient Population and Study Design

The

Efficacy and Safety Assessments

The primary end point was OS, defined as the time interval from the date of random assignment to the date of death from any cause. Secondary end points included progressionfree survival (PFS), QoL, safety and medical costs. Tumor assessment was performed 6

weeks and 12 weeks after the randomisation date and every 12 weeks thereafter. PFS was assessed according to RECIST 1.1. Quality of life (QoL) was assessed using the EORTC QoL core questionnaire (QLQ-C30) [START_REF] Aaronson | The European Organization for Research and Treatment of Cancer QLQ-C30: a qualityof-life instrument for use in international clinical trials in oncology[END_REF] and with the oesophagus-specific questionnaire (QLQ-OES-18) [START_REF] Blazeby | Clinical and psychometric validation of an EORTC questionnaire module, the EORTC QLQ-OES18, to assess quality of life in patients with oesophageal cancer[END_REF] at baseline and every 6 weeks thereafter until 42 weeks after randomisation. Adverse events (AE) were graded according to the National Cancer Institute

Common Terminology Criteria version 4.0. Patients were observed until death or until 48 months after study entry.

Statistical Considerations

Sample size was calculated to ensure a precision of 12.5% to estimate the 9-month OS rate with a two-sided 85% confidence interval (85% CI) in the CT-CONT arm. Assuming a 9month OS rate of 56% [START_REF] Lorenzen | Cetuximab plus cisplatin-5-fluorouracil versus cisplatin-5-fluorouracil alone in first-line metastatic squamous cell carcinoma of the esophagus: a randomized phase II study of the Arbeitsgemeinschaft Internistische Onkologie[END_REF], 31 patients were required in each arm. The CT-DISC arm served as an internal control. With an anticipated 70% eligibility rate for randomisation after 6 weeks of CT, 88 patients were required in the selection part of the trial in order to randomise 62 patients. In 2013, we found that the actual eligibility rate for the randomised part was only 58%. The protocol was thus amended in order to select 106 patients. Randomisation was performed using a minimization method controlling for the following factors: previous chemotherapy (no vs yes), dysphagia (Atkinson grade 1-2 vs 3-4) and EQ-5D visual analog scale (<40 vs ≥40).

Statistical analyses of efficacy endpoints were performed per randomised arm by the intention-to-treat approach. Survival estimates were calculated per treatment arm using the Kaplan-Meier method from the date of random assignment. An unplanned and exploratory post hoc analysis was performed to estimate OS curves in CT-DISC patients according to whether or not they received post-progression chemotherapy. All randomised patients were included in the QoL analysis. For functional scales, QLQ-C30 and QLQ-OES-18 scores were considered as a definitive deterioration if the score decreased by more than 10 points compared with the score at randomisation, and without later improvement superior to 10 points compared with baseline. For symptom scales, a definitive deterioration was defined as an increase of 10 points or more without subsequent decrease. For each dimension of the QoL questionnaires, the time until definitive deterioration (TUDD) was defined as the time from randomisation to the first observation of a definitive deterioration of the corresponding score or death. Patients alive without reported definitive deterioration were censored at the date of the last follow-up visit. Patients without any QoL questionnaires were censored at randomization [START_REF] Anota | Time to health-related quality of life score deterioration as a modality of longitudinal analysis for health-related quality of life studies in oncology: do we need RECIST for quality of life to achieve standardization?[END_REF]. TUDD was estimated using the Kaplan-Meier method. The impact of treatments on the different dimensions of the QoL was estimated by hazard ratios (HR CT- DISC/CT-CONT ) of QoL deterioration using Cox models. As QoL questionnaires could be missing during follow-up, leading to a possible overestimation of the TUDD, we also performed a sensitivity analysis considering the following imputations in the absence of a definitive deterioration of the QoL score: the date of definitive deterioration was imputed three weeks after the date of the last completion of the QoL questionnaire if the patient died more than two months later, or if the patient was alive without subsequent QoL; patients having completed all QoL questionnaires planned in the study were censored at the date of the last QoL questionnaire plus three weeks.

In the CT-CONT arm, treatment-related AEs were analysed considering the maximum grade per patient and per type of AE.

Results

Between January 2011 and February 2015, 105 patients were recruited in the selection phase of the study. The patient distribution in the trial is illustrated in Figure 1. Ninety patients were eligible for tumor assessment and 69 were found to be disease-controlled with an ECOG performance status of 0 to 2. Among them, 67 patients were randomly assigned, 34 to continue chemotherapy and 33 to discontinue chemotherapy. At the time of the final data extraction, the median follow-up for the randomised patients was 36.9 months.

Patient and treatment characteristics: Baseline characteristics and details of the treatments that were administered before and after randomisation are given in Appendix Table S1 and in Table 1 respectively. The two treatment groups were well balanced with regard to baseline characteristics of randomised patients, and of eligible and treated patients (data not shown).

Among the 31 eligible and treated patients in the CT-CONT arm, 7 (23%) received LVFU2-CDDP and 24 (77%) received FOLFOX. In the CT-DISC group, 24 patients resumed chemotherapy after having progressed, while 8 patients did not (missing data for 1 patient).

Efficacy: At the time of the analysis, 53 deaths had been reported, 51 after disease progression, one due to hepatic encephalopathy 3.5 months, and one from an unknown cause 0.7 months, both in the CT-DISC arm. The estimated 9-month survival rate was 50% (85% CI: 37-62%) and 48% (85% CI: 35-60%) in the CT-CONT and the CT-DISC groups, respectively. The median OS rate was 8.5 months (95% CI: 6.6 to 12 months) and 8.8 months (95% CI: 5.9 to 13.4 months) for CT-CONT and CT-DISC, respectively. The median PFS was 4 months (95% CI: 2.8 to 5.8 months) and 1.4 months (95% CI: 1.4 to 2.7 months)

for CT-CONT and CT-DISC, respectively. Kaplan-Meier survival curves are shown in Figure 2. As illustrated in Figure 3, the median OS was 9.9 months (95% CI: 6.3 to 16.9 months) and 3.5 months (95% CI: 2 to 15.4 months) for the 24 patients (24/33, 72.7%) in the CT-DISC group that received some post-progression chemotherapy and for the 8 patients that did not, respectively.

Quality of Life: QLQ-C30 and QLQ-OES18 domain scores showed no systematic differences between treatment arms at randomisation (data not shown). The median TUDD by treatment arm and Hazard Ratio (HR) for each dimension of both QoL questionnaires are displayed in Table 2 and Figure 4. We observed a trend for a benefit in favor of the CT-CONT arm (HR>1) for most dimensions, with a significant difference for three dimensions of the QLQ-OES18: dysphagia with a median TUDD of 7.3 months (95% CI: 4.2-12.0) and 2.9 months (95% CI: 1.4-4.4) for CT-CONT and CT-DISC, respectively, leading to an estimated HR of 1.98 (95%CI, 1.15-3.4), p-value=0.014; eating with a median TUDD of 7.7 months (95% CI:

5.6-9.5 ) and 2.9 months (95% CI: 2.0-5.9), respectively, with HR=1.75 (1.02-3.02), p-value=0.044; and pain with a median TUDD of 8.1 months (95% CI: 5.6-12.0) and 2.4 months (95% CI: 1.4-3.2), respectively, HR=2.52 (1.43-4.43), p-value=0.001. The results appeared stable in the sensitivity analysis (Appendix Table S3).

Adverse events and other assessments:

There were no treatment-related deaths during the study. As detailed in Appendix Table S2, toxicity was mild and without unexpected AEs. The results of the medical cost analyses will be reported elsewhere.

Discussion

Analysis of the primary end point for this multicenter, randomised, discontinuation phase-2 trial in mESCC patients who continued chemotherapy after having been disease-controlled with a FU/platinum regimen indicated a 50% 9-month survival rate. Patients who had been assigned to CT-DISC had a 48% 9-month survival rate. The median OS and median PFS were 8.5 and 4.0 months in the CT-CONT arm, respectively, and 8.8 and 1.4 months in the CT-DISC arm, respectively. Although one could have anticipated a better OS for patients in the CT-CONT arm that, actually, had been selected on the basis of early control of their disease, our survival results are of the same order of magnitude to those previously reported with FU/platinum regimens in contemporary randomised-studies [START_REF] Wang | Efficacy and safety of recombinant human lymphotoxin-α derivative with cisplatin and fluorouracil in patients with metastatic esophageal squamous cell carcinoma: A randomized, multicenter, open-label, controlled, phase 2b trial[END_REF][START_REF] Lorenzen | Cetuximab plus cisplatin-5-fluorouracil versus cisplatin-5-fluorouracil alone in first-line metastatic squamous cell carcinoma of the esophagus: a randomized phase II study of the Arbeitsgemeinschaft Internistische Onkologie[END_REF][START_REF] Moehler | Cisplatin/5-FU (CF) +/-panitumumab (P) for patients (pts) with non-resectable, advanced, or metastatic esophageal squamous cell cancer (ESCC): An open-label randomized AIO/TTD/BDGO/EORTC phase III trial (POWER)[END_REF]. E-DIS was not designed to detect a superiority of one arm over the other. Nevertheless, the two OS curves are so close and intertwined that one cannot claim a difference in OS. To understand why patients who discontinued chemotherapy had the same OS as patients who continued chemotherapy until progression, we examined the outcome of CT-DISC patients. First, most of these patients (72.7%) received subsequent chemotherapy at progression. As a consequence, the majority of patients in the CT-DISC arm experienced a transient chemotherapy break rather than permanent discontinuation. Second, the shape of the OS curve for patients that resumed chemotherapy appeared better than the OS curve for patients that did not receive post-progression chemotherapy in the CT-DISC arm. Third, the OS curve for patients that resumed chemotherapy resembled those for patients that had been randomly allocated to the CT-CONT arm (Figure 3). Taken together, these features suggest that chemotherapy might have some favorable impact on the OS. However, we acknowledge that our observation is tentative, as the number of patients was low, and as we cannot rule out a selection bias, i.e., a poor performance status may have resulted in some patients not resuming chemotherapy.

In this study, most patients (77%) received FOLFOX and 23% received LVFU2-CDDP. The latter, a biweekly regimen combining CDDP, bolus FU and infusional FU over 2 days, is regarded as a convenient way to deliver FU and CDDP and a safer regimen than the monthly FU-CDDP regimen with FU given as a continuous infusion over 4 or 5 days [START_REF] Bouché | Randomized multicenter phase II trial of a biweekly regimen of fluorouracil and leucovorin (LV5FU2), LV5FU2 plus cisplatin, or LV5FU2 plus irinotecan in patients with previously untreated metastatic gastric cancer: a Federation Francophone de Cancerologie Digestive Group Study--FFCD 9803[END_REF]. Oxaliplatin is a platinum derivative that has a more favorable toxicity profile than CDDP, and FOLFOX has been previously investigated in the treatment of mESCC with comparable efficacy, but with better safety than FU-CDDP [START_REF] Wang | A phase II study of oxaliplatin in combination with leucovorin and fluorouracil as first-line chemotherapy in patients with metastatic squamous cell carcinoma of esophagus[END_REF]. Actually, FOLFOX is becoming popular in daily practice, since it has been shown to be equivalent to FU-CDDP in terms of efficacy in locally advanced settings [START_REF] Conroy | Definitive chemoradiotherapy with FOLFOX versus fluorouracil and cisplatin in patients with oesophageal cancer (PRODIGE5/ACCORD17): final results of a randomised, phase 2/3 trial[END_REF][START_REF] Bascoul-Mollevi | Health-related quality of life results from the PRODIGE 5/ACCORD 17 randomised trial of FOLFOX versus fluorouracil-cisplatin regimen in oesophageal cancer[END_REF].

CT-CONT is associated with a significant delay in the worsening of some major symptoms such as dysphagia, eating and oesophageal pain (Figure 4). The results of the sensitivity analysis performed to control a possible overestimation of the time until definitive deterioration in the primary analysis of QoL were stable. Indeed, the study of the Qol has special relevance when life expectancy is short and when the benefit, if any, of some treatment can be perceived as modest. Consequently, it is satisfactory to observe an apparent benefit of the QoL in patients who continued chemotherapy.

Because of its limitations, E-DIS should be interpreted carefully. First, we acknowledge that part of our intervention was not standardised enough. At the onset of our trial, there was no clear evidence that chemotherapy provided a clinical benefit to mESCC patients, and ESMO guidelines suggested that either the best supportive care or palliative chemotherapy should be considered in this setting [START_REF] Lordick | Oesophageal cancer: ESMO Clinical Practice Guidelines for diagnosis, treatment and follow-up[END_REF]. Consequently, we decided to accept all types of FU/platinum-based regimens as selection treatments. In 2009, at the time when the protocol was being written, we asked trial participants to provide on-demand supportive care to patients regardless of the treatment-arm to which they were assigned. Actually, we cannot certify that this supportive care was standardised adequately to avoid hypothetical differences between the two arms of this open-label trial, leading to a possible risk of systematically over-estimating the net clinical benefit of the control arm, as has been recently

suggested [START_REF] Nipp | Best supportive care in clinical trials: review of the inconsistency in control arm design[END_REF], although we were referring to national recommendations that were supposed to be followed similarly in the two arms of this trial. Second, the main apparent benefit of CT-CONT relies on QoL analysis, which was a secondary endpoint, with multiple comparisons and is subject to possible bias. Finally, E-DIS was designed as a non-comparative trial, therefore, all HRs are exploratory.

In conclusion, chemotherapy until progression provides an OS rate that is numerically similar to chemotherapy interruption in mESCC patients who had been disease-controlled with a 6week selection-course of a FU/platinum regimen. E-DIS provides valuable data to support shared decision-making between physicians and patients regarding CT-CONT/DISC. The HR of definitive deterioration in the CT-DISC group compared to the CT-CONT group. The 95% confidence intervals are estimated in a Cox model for each component separately.

(1) Two patients with an adenocarcinoma were withdrawn from the study after enrolment because they were classified as non-eligible (no measurable disease at study entry). (2) Two patients were not treated due to poor general health status.

(3) Four patients had no tumor evaluation 6 weeks after start of treatment, due to patient's decision, no compliance, concurrent disease, and poor performance status.

Figure 1 -Figure 2 -

 12 Figure 1 -CONSORT diagram. ITT, intention to treat Figure 2 -Kaplan-Meier survival curves for (A) overall survival, and (B) progression-free survival according to treatment arm allocated by randomisation (intention-to-treat).

Figure 3 -

 3 Figure3-Kaplan-Meier estimates of overall survival for patients who were randomly assigned to discontinue chemotherapy (CT-DISC, n=33), according to whether (n=24, blue curve) or not (n=8, yellow curve) they received post-progression chemotherapy (PPC). The dotted curve is the survival curve of patients (n=34) who were randomly assigned to continue chemotherapy (CT-CONT).

Figure 4 -

 4 Figure 4 -Relative treatment effect on the quality of life endpoints (Forest plot of the Hazard Ratio on definitive deterioration of quality of life).

in Appendix Table S1 and inTable 1 .

 1 Patients free from progression after 6 weeks of treatment were assigned equally to

	multicenter E-DIS study selected patients before starting a first-line FU/platinum-based
	chemotherapy. Key selection criteria included histologically confirmed mESCC, measurable
	disease, age greater than 18 years, and an ECOG performance status of 0, 1 or 2. Prior
	chemotherapy was permitted only if it was delivered as a neoadjuvant treatment. The choice
	of the FU/platinum-based regimen was left physician's decision among the following ones:
	FU-CDDP, FU-CDDP-docetaxel, FU-folinic acid-CDDP (LVFU2-CDDP) and FU-folinic acid-
	oxaliplatin combination (FOLFOX). Chemo dosing is detailed either continue the same chemotherapy or to discontinue chemotherapy. In the CT-CONT
	arm, the study treatment continued until there was disease progression, unacceptable
	toxicity, or a patient or physician decision to terminate the treatment. In the CT-DISC arm,
	chemotherapy could be resumed after disease progression. On-demand supportive care was
	offered to any patients. The protocol complied with the recommendations of the 18 th World
	Health Congress (Helsinki, Finland, 1964) and its subsequent amendments, good clinical
	practice guidelines, and other legal requirements. The protocol, including all amendments,
	was reviewed and approved by the CPP Nord-Ouest III ethics committee on August 22th,
	2013. Patients provided written informed consent before enrollment in the study.

The authors have declared no competing interest.

Support and role of the Funding Source

This work was supported by the French Ministry of Health (PHRC-K 2010 19-01). The funding source had no role in the design and conduct of the study; collection, management, analysis and interpretation of the data; preparation or approval of the manuscript and decision to submit the manuscript for publication

Acknowledgments

We are thankful for the contributions from patients, their families and other participating coinvestigators: C. Peugniez, N. Kotecki (Lille, France) and L. Bedenne (Dijon, France). We would also like to thank A. Kramar (Lille, France) for the fruitful discussion we had regarding the statistical plan, Aurélie Thenot-Maricourt (Lille, France) as Project Manager, and the Clinical Research Associates from our institution (Lille, France) for on-site data monitoring.

The Northwest Data Center (CTD-CNO) is acknowledged for managing the data. This was supported by grants from the French National League Against Cancer (LNC) and the French National Cancer Institute (INCa). We also thank the UNICANCER Gastro-Intestinal group and for PRODIGE intergroup for their friendly support. Finally, we are deeply grateful to our colleague and friend, Professor Christophe Mariette (Univ. Lille, Department of Digestive and Oncological Surgery, Claude Huriez University Hospital, Lille, France) who passed away in July, 2017. Without his contribution and support this research project could not have been conducted.

Conflict of interest statement: None declared

TABLE and FIGURE LEGENDS

Table [START_REF] Ferlay | Cancer incidence and mortality worldwide: sources, methods and major patterns in GLOBOCAN 2012[END_REF] -Baseline characteristics of the 67 patients who were included in the randomised part of the study, by treatment group. Appendix Table S1 -Baseline characteristics of the 101 patients registered in the selection part of the study and who started chemotherapy, according to their subsequent participation in the randomised phase of the trial. S2 -Maximal grade per patient and per type of adverse events related to the study treatment in CT-CONT arm. S3 -Quality of life (sensitivity analysis): median time until definitive deterioration (TUDD) by treatment group and Hazard Ratio for each dimension in the quality of life questionnaires (QLQ-C30 and QLQ-OES18).

Appendix Table

Appendix Table

6%

(1) One patient can have several metastatic sites.

(2) FU-CDDP: FU 1,000 mg/m² as a continuous infusion over 96 h, plus CDDP 100 mg/m² day 1 or 2, every 4 weeks.

(3) FU-CDDP-docetaxel: docetaxel 75 mg/m² on day 1, CDDP 75 mg/m² on day 1, and FU 750 mg/m² as a continuous infusion over 120 h, every 4 weeks. (4) LVFU2-CDDP: FU-folinic acid-CDDP, with CDDP 50 mg/m², folinic acid 200 mg/m², and bolus FU 400 mg/m² on day 1 followed by FU 2,400 mg/m² as a continuous infusion over 46 h, every 2 weeks.

(5) FOLFOX: FU-folinic acid-oxaliplatin combination with oxaliplatin 85 mg/m 2 , folinic acid 200 mg/m², and bolus FU 400 mg/m² on day 1 followed by FU 2,400 mg/m² as a continuous infusion over 46 h, every 2 weeks. (6) For 2 patients, known as non-progressive at 6 weeks, detailed information about the tumour response at 6 weeks is not available.


Túlio Honório

Thibault Lemaire

Di Tommaso

Salah Naili

Molecular modelling of the heat capacity and anisotropic thermal expansion of nanoporous hydroxyapatite

Keywords: Hydroxyapatite, Confined water, Thermal expansion, Specific heat, Molecular dynamics, Homogenization

Hydroxyapatite, which is the main mineral phase of mammalian bone, occurs in the form of small bricks of colloidal size organized in a way that leaves room to micro-and mesopores. These pores are filled with an electrolyte and confined fluids are recognized for manifesting different dynamical and structural behaviors when compared to bulk fluids. Research on other nanoporous materials reported that confinement may have repercussions on the effective the thermal properties of these materials. Understanding the physical origin of thermal expansion and heat capacity as a function of the hydroxyapatite porosity is, therefore, crucial to predict the thermo-mechanical behavior of bone. Molecular dynamics simulations of hydroxyapatite nanopores (2 nm ≤ H ≤ 16 nm, where H is the size of the nanopore) in contact with liquid water have been used to determine the effect nanoporosity and water confinement on the heat capacity and thermal expansion of this important biomaterial.

At temperatures corresponding to in vivo conditions, the thermal expansion of water confined in nanopores smaller than 6 nm was is solid-like but becomes liquid-like in larger nanopores. The heat capacity of confined water exhibits a maximum at pore sizes of approximately 7 nm. An up-scaling strategy taking

Introduction

Properties such as thermal expansion and heat capacity characterize the behavior of mammalian bone under varying temperature conditions. In fact, despite it is often assumed that bone is subject to small temperature variations, even subtle differences between the bone surface and deep-body zone temperatures might impact the local metabolism [START_REF] Taylor | Considerations for the measurement of core, skin and mean body temperatures[END_REF]. Moreover, orthopedic [START_REF] Karmani | The thermal properties of bone and the effects of surgical intervention[END_REF][START_REF] Krause | Temperature elevations in orthopaedic cutting operations[END_REF][START_REF] Augustin | Cortical bone drilling and thermal osteonecrosis[END_REF][START_REF] Dolan | Heatshock-induced cellular responses to temperature elevations occurring during orthopaedic cutting[END_REF] and odontological [START_REF] Jacobs | The effect of heat on bone healing. A disadvantage in the use of power tools[END_REF] interventions can to produce significant thermal effects.

For instance, Jacobs et al. [START_REF] Jacobs | The effect of heat on bone healing. A disadvantage in the use of power tools[END_REF] observed temperature rises of at least 5 • C in both graft bed and bone with dental burrs, and Krause et al. [START_REF] Krause | Temperature elevations in orthopaedic cutting operations[END_REF] reported clinical tests during total joint replacement operations in which temperatures exceed 200 • C when irrigation was not used. Temperatures above 45 • C were reported to cause thermal osteonecrosis [START_REF] Augustin | Cortical bone drilling and thermal osteonecrosis[END_REF][START_REF] Dolan | Heatshock-induced cellular responses to temperature elevations occurring during orthopaedic cutting[END_REF]. The recognition of the importance of these thermal effects goes back to Hippocrates' recommendation (500 BC) of frequently plunging the tools into cold water when trephining discs of bone from the skull to prevent the heating of the material [START_REF] Karmani | The thermal properties of bone and the effects of surgical intervention[END_REF][START_REF] Phillips | Greek medicine[END_REF]. Healing following bone surgery may be delayed if the bone cells are injured by thermal effects [START_REF] Eriksson | Temperature threshold levels for heat-induced bone tissue injury: A vital-microscopic study in the rabbit[END_REF].

From another perspective, thermal stimulation of tissues aiming at inducing localized necrosis has been proposed as a therapeutic tool [START_REF] Calttenburg | Thermal properties of cancellous bone[END_REF]. However, relatively few experimental and theoretical studies have been concerned with the thermal properties of bone.

As the main constituent inorganic phase of mineralized tissue found in mammalian bone, tooth enamel and dentin, Hydroxyapatite (HAP), Ca 10 (PO 4 ) 6 (OH) 2 , can be used as a proxy to study the thermal properties of bone. In these tissues, HAP occurs in the form of small 'bricks' of colloidal size [START_REF] Weiner | Organization of hydroxyapatite crystals within collagen fibrils[END_REF][START_REF] Kuhn | Size and Shape of Mineralites in Young Bovine Bone Measured by Atomic Force Microscopy[END_REF] organized in a way that leaves room to pores sizes ranging from 2 to 24 nm [START_REF] Holmes | Gas Adsorption and Surface Structure of Bone Mineral[END_REF] in contact with body fluid. This was for example showed using high resolution transmission electron microscopy (TEM) and electron diffraction as well as 3D stereoscopic TEM [START_REF] Weiner | Organization of hydroxyapatite crystals within collagen fibrils[END_REF], and further confirmed by atomic force microscopy [START_REF] Kuhn | Size and Shape of Mineralites in Young Bovine Bone Measured by Atomic Force Microscopy[END_REF]. Confined water manifests different dynamical and structural behaviors when compared to bulk fluids (e.g. [START_REF] Brovchenko | Interfacial and Confined Water[END_REF]), and research on other porous materials showed that nanoconfinement may have repercussions on their effective thermal properties [START_REF] Etzler | The heat capacity of water in silica pores[END_REF][START_REF] Etzler | Temperature dependence of the heat capacity of water in small pores[END_REF][START_REF] Etzler | The structure and properties of vicinal water: Lessons from statistical geometry[END_REF][START_REF] Xu | Thermal Expansion of Confined Water[END_REF][START_REF] Honorio | Thermal properties of cement-based materials: Multiscale estimations at early-age[END_REF][START_REF] Qomi | Physical Origins of Thermal Properties of Cement Paste[END_REF].

The molecular dynamics (MD) simulation technique is particularly suited to assess the behavior of water confined in nanopores as it allows the quantification of interatomic forces and atomistic phenomena responsible for the anomalous behaviour of water, which otherwise would be difficult to probe experimentally, especially at the temperature and pressure corresponding to in vivo conditions [START_REF] Lemaire | Modelling of the transport in electrically charged porous media including ionic exchanges[END_REF]. Recent works of the authors using MD revealed important aspects of the structure and dynamics of water in HAP-water nanopores [START_REF] Honorio | Anomalous water and ion dynamics in hydroxyapatite mesopores[END_REF][START_REF] Pham | Properties of water confined in hydroxyapatite nanopores as derived from molecular dynamics simulations[END_REF][START_REF] Prakash | Transport properties of water molecules confined between hydroxyapaptite surfaces: A Molecular dynamics simulation approach[END_REF][START_REF] Di Tommaso | Molecular Dynamics Simulations of Hydroxyapatite Nanopores in Contact with Electrolyte Solutions: The Effect of Nanoconfinement and Solvated Ions on the Surface Reactivity and the Structural, Dynamical, and Vibrational Properties of Water[END_REF], including the mechanisms controlling the anisotropic diffusion of water molecules within nanopores [START_REF] Prakash | Anisotropic diffusion of water molecules in hydroxyapatite nanopores[END_REF].

We anticipate that nanoporosity could also have a significant effect on the thermal properties HAP, but despite its importance in biomechanic applications, to the best of the authors' knowledge, no experimental or theoretical data is available regarding the role of nanoporosity on the thermal expansion and heat capacity. In this work, the MD technique has been applied to compute the thermal expansion and heat capacity of bulk HAP and HAP-water nanopores with pore sizes ranging from 2 to 16 nm, and to quantify the effect of water nanoconfinement on the thermal properties of hydroxyapatite. An upscaling strategy based on homogeneization techniques is then used to link the anomalous thermal behaviour in HAP-water nanopores to the effective thermal expansion and effective heat capacity of HAP, as a function of its porosity.

Methods

Heat capacity and thermal expansion from molecular dynamics simulations

At the nanoscopic level, the heat capacity of a material is related to the vibrational modes of atoms (phonons) and the thermal expansion to the anharmonicity of the lattice interaction potential [START_REF] Leibfried | Theory of Anharmonic Effects in Crystals[END_REF]. The heat capacity at constant strain C quantifies how the internal energy U changes with temperature T :

C = ∂U ∂T ij , (1)
where ij are the components of the second-order Lagrangian deformation tensor which is symmetric and associated with a reference volume V 0 . In this paper, we adopt the convention that the subscripts as i, j and k, which vary from 1 to 3, correspond to the Cartesian coordinates x, y and z, respectively (see Fig. 1).

So, the indices 1, 2 and 3 map to x, y and z respectively. For a unit cell, the deformation tensor is defined by:

= 1 2 h -T 0 h T hh -1 0 -I , (2)
where h 0 and h are two tensors associated with the metric of the reference and deformed unit cell, respectively, I is the unit tensor, and the superscript [ ] T denotes the transpose operator. Note that the unit cell is periodically repeated to fill all space whose the form and the volume are completely arbitrary but described by the 3 vectors a, b and c that span the edges of the unit cell. This description can be obtained by arranging the vectors as {a, b, c} to form a [START_REF] Krause | Temperature elevations in orthopaedic cutting operations[END_REF][START_REF] Krause | Temperature elevations in orthopaedic cutting operations[END_REF] matrix h whose the columns are the components of a, b and c (see Fig. 1) in this order.

The heat capacity at constant stress C σ quantifies how the enthalpy changes with the temperature:

C σ = ∂H ∂T σij = ∂(U + V 0 Tr(σ )) ∂T σij , (3)
where σ is the second-order Cauchy stress tensor which is symmetric and whose the components are denoted by σ ij [START_REF] Cagin | Thermodynamic and elastic properties of polyethylene at elevated temperatures[END_REF][START_REF] Lee | Monte Carlo Simulations of Structural Phase Transitions in Metals and Alloys[END_REF][START_REF] Ray | Statistical ensembles and molecular dynamics studies of anisotropic solids. II[END_REF][START_REF] Ray | Molecular dynamics equations of motion for systems varying in shape and size[END_REF] and Tr denotes the trace operator of a tensor. The second-order thermal expansion tensor α, which is symmetric and defined by these components α ij , quantifies how the strain changes with the temperature under a constant stress:

α ij = ∂ ij ∂T σij . (4)
The heat capacities C σ and C , and the elements of the thermal expansion tensor are related by the following expression [START_REF] Wallace | Thermodynamics of Crystals[END_REF]:

C σ -C = V 0 T C ijkl α ij α kl , (5)
where C ijkl is the fourth-order elastic stiffness tensor which is symmetric. In this paper, the Einstein summation convention is used which implies a summation over a set of indexed terms in a formula.

The heat capacities and the thermal expansion tensor have been computed from MD simulations using the fluctuation-dissipation formula. The fluctuation of two dynamical variables A and B is given by δA × δB ens = AB ens -A ens B ens , where the operator δ ens denotes an ensemble average. In the isotension-isothermal (NσT) ensemble, the constant stress heat capacity can be computed from the fluctuations of the enthalpy H [START_REF] Allen | Computer Simulation of Liquids[END_REF]:

C σ = δH 2 N σT kT 2 , ( 6
)
where k is the Boltzmann constant. In the isotension-isothermal (NσT) ensemble, the thermal expansion tensor can then be computed from the cross fluctuations of the strain and the enthalpy [START_REF] Cagin | Thermodynamic and elastic properties of polyethylene at elevated temperatures[END_REF][START_REF] Lee | Monte Carlo Simulations of Structural Phase Transitions in Metals and Alloys[END_REF][START_REF] Ray | Molecular dynamics equations of motion for systems varying in shape and size[END_REF]:

α ij = δ ij × δH N σT kT 2 . ( 7
)
The heat capacity at constant strain C can then be computed using Eq. ( 5).

Simulation details

The atomic structure of HAP was obtained from the hexagonal with P6 3 /m space group cell parameters of Sudarsanan and Young [START_REF] Sudarsanan | Significant precision in crystal structural details. Holly Springs hydroxyapatite[END_REF]. The simulation box is generated by stacking this unit cell 3 times in each in-plane direction and 4 times in out-of-plane direction (see Fig. 1), which results in a total number of atoms in the HAP solid layer of N HAP = 1584 atoms.

For simulation of nanoporous HAP-water systems, HAP platelets were stacked with a surface-to-surface spacing along c-direction varying from 2 to 16 nm. The {0001} surface of HAP, which is the dominant surface in the thermodynamicallystable morphology [START_REF] Mkhonto | A computer modelling study of the effect of water on the surface structure and morphology of fluorapatite: introducing a Ca10(PO4)6f2 potential model[END_REF], was considered as the basal plane. The slit pore was filled with water molecules at an experimental density of 1 g.cm -3 , which is close to the density of SPC/E water at ambient conditions.

The interatomic potential proposed by de Leeuw and co-workers [START_REF] De Leeuw | A computer modelling study of the uptake and segregation of fluoride ions at the hydrated hydroxyapatite (0001) surface: introducing a Ca 10 (PO 4 ) 6 (OH) 2 potential model[END_REF][START_REF] De Leeuw | Molecular-dynamics simulation of MgO surfaces in liquid water using a shell-model potential for water[END_REF] to model crystalline apatite was used for HAP. In this force field, non-bonded interactions are described by Buckingham potentials, phosphate and hydroxyl group bonds as the sum of a Morse and a Coulombic potential, and phosphate bond angles by a harmonic potential. The polarizability of oxygen anion in the phosphate and hydroxyl groups is described by a core-shell model. In HAP-water systems, water were modeled using the extended single-point charge (SPC/E) water model [START_REF] Berendsen | The missing term in effective pair potentials[END_REF], and HAP-water interactions were described using non-bonded potentials derived in previous MD studies [START_REF] Almora-Barrios | Modelling the interaction of a Hyp-Pro-Gly peptide with hydroxyapatite surfaces in aqueous environment[END_REF][START_REF] Hernandez | The effect of water on the binding of glycosaminoglycan saccharides to hydroxyapatite surfaces: a molecular dynamics study[END_REF].

MD simulations were performed using DL POLY 4.05.1 code [START_REF] Todorov | DL poly 3 new dimensions in molecular dynamics simulations via massive parallelism[END_REF]. The systems were first equilibrated in the microcanonical (NVE) ensemble for 1 ns.

Then, a first isotension-isothermal (NσT) run during 1ns was performed to equilibrate the system under a pressure of 1 atm and at a fixed temperature (from 290 to 330 K). The Nosé-Hoover thermostat and barostat, at relaxation time 500 fs for both, were used. Smoothed Particle Mesh Ewald (SPME) method with the acceptable relative error of 10 -6 were employed to compute the long-range electrostatic interactions.

In order to verify the convergence of the estimations of the heat capacity and coefficients of thermal expansion, in the production phase the system was sampled for 2 ns in NσT ensemble (timestep of 1 fs). Figure 2 displays the fluctuations of the enthalpy, volume and deformations according to the in-and out-of-plane directions. To compute the deformations, the reference cell parameters defined via h 0 (see Eq. ( 1)) were determined from the average at zero tension in NσT ensemble (i.e. in a state in which the deformation vanishes [START_REF] Ray | Elastic constants and statistical ensembles in molecular dynamics[END_REF]).

After approximately 1 ns the fluctuations yield smooth Gaussian distributions as shown in Fig. 2. From the estimations of the heat capacity, the isothermal compressibility and the coefficients of thermal expansion, we assume that a simulation run of 1 ns in isotension-isothermal ensemble captures the normal distribution of the thermodynamic quantities (a)-(c) and is enough to converge the results of the computed properties (d)-(f ).

In this study, the pore size H is defined as the projection of c vector on the z-direction minus the effective thickness h HAP of the solid HAP [START_REF] Honorio | Anomalous water and ion dynamics in hydroxyapatite mesopores[END_REF]. Following the definition of effective thickness and pore size used for other layered materials [START_REF] Honorio | Hydration Phase Diagram of Clay Particles from Molecular Simulations[END_REF], h HAP is defined as the distance between the center of the outmost oxygen (the out-most species in HAP surface) in each exposed HAP sur- values reported in the literature. To obtain further validation of our atomistic model, we also computed the components of elastic stiffness tensor using the Parinello-Rahman fluctuation equation [START_REF] Parrinello | Strain fluctuations and elastic constants[END_REF]:

C ijkl = kT V N σT [ δ ij × δ kl N σT ] -1 .
The representation adopted in this paper is called the Voigt's representation in which the new indexes I et J vary in the set {1, ..., 6} such as I = (i, j) and ), and (f ) coefficients of thermal expansion. For the latter, all components of the tensor α, computed from the fluctuations of the deformations with Eq. ( 7)) as well as its trace are reported shown.

This work Other works

C σ (J.kg -1 .K -1 ) 875±72 768 [START_REF] Kijima | Preparation and Thermal Properties of Dense Polycrystalline Oxyhydroxyapatite[END_REF]; 765 (298 K) [START_REF] Egan | High-Temperature Heat Content of Hydroxyapatite[END_REF] 690 (MD 298 K) [START_REF] Cruz | A Molecular Dynamics Study of the Thermodynamic Properties of Calcium Apatites. 1. Hexagonal Phases[END_REF] 632 (MD monoclinic 298 K) [START_REF] Cruz | A Molecular Dynamics Study of the Thermodynamic Properties of Calcium Apatites. 1. Hexagonal Phases[END_REF] α 11 (×10 -5 K -1 ) 0.68±0.50 0.39 [START_REF] Nakamura | Thermal expansion of hydroxyapatite--tricalcium phosphate ceramics[END_REF] α 22 (×10 -5 K -1 ) 1.38±0.81 -

α 33 (×10 -5 K -1 )
1.31±0.32 1.86 [START_REF] Nakamura | Thermal expansion of hydroxyapatite--tricalcium phosphate ceramics[END_REF] (1/3) × Tr(α) (×10 -5 K -1 ) 1.13±1.01 1.82 (monoclinic 298 K) [START_REF] Cruz | Molecular Dynamics Study of the Thermodynamic Properties of Calcium Apatites. 2. Monoclinic Phases[END_REF] C 11 (GPa) 97 140.0 [START_REF] Ching | Ab initio elastic properties and tensile strength of crystalline hydroxyapatite[END_REF] (ab initio); 117.1 [START_REF] Snyders | Experimental and ab initio study of the mechanical properties of hydroxyapatite[END_REF]; 166.7 [START_REF] Hearmon | An Introduction to Applied Anisotropic Elasticity[END_REF]; 115 [START_REF] Lees | Anisotropy in hard dental tissues[END_REF]; 135 [START_REF] Gilmore | Elastic properties of apatites[END_REF] C 22 (GPa) 93 134.8 [START_REF] Ching | Ab initio elastic properties and tensile strength of crystalline hydroxyapatite[END_REF] (ab initio)

C 33 (GPa) 120 174.8 [START_REF] Ching | Ab initio elastic properties and tensile strength of crystalline hydroxyapatite[END_REF] (ab initio); 231.8 [START_REF] Snyders | Experimental and ab initio study of the mechanical properties of hydroxyapatite[END_REF]; 139.6 [START_REF] Hearmon | An Introduction to Applied Anisotropic Elasticity[END_REF]; 125 [START_REF] Lees | Anisotropy in hard dental tissues[END_REF]; 172 [START_REF] Gilmore | Elastic properties of apatites[END_REF] κ (GPa -1 ) 0.013±0.004 0.0118 [START_REF] Ching | Ab initio elastic properties and tensile strength of crystalline hydroxyapatite[END_REF] (ab initio); 0.012 [START_REF] Snyders | Experimental and ab initio study of the mechanical properties of hydroxyapatite[END_REF]; 0.016 [START_REF] Lees | Anisotropy in hard dental tissues[END_REF]; 0.011 [ There is a general good agreement between the properties obtained from our simulations and the experimental values, which corroborates the transferability of the empirical force fields for HAP proposed by de Leeuw and co-workers [START_REF] De Leeuw | A computer modelling study of the uptake and segregation of fluoride ions at the hydrated hydroxyapatite (0001) surface: introducing a Ca 10 (PO 4 ) 6 (OH) 2 potential model[END_REF][START_REF] De Leeuw | Molecular-dynamics simulation of MgO surfaces in liquid water using a shell-model potential for water[END_REF] for the computation of thermal properties.

γ G,ij = V C σ C ijkl α kl . (8)
For hexagonal crystals, there are only two independent Grüneisen parameters [START_REF] Grimvall | Thermophysical Properties of Materials[END_REF]: the in-plane γ G,11 = γ G,22 and the out-of-plane γ G,33 parameters (the basal plane is defined with respect to the {0001} surface). In fact, the results in Fig. 4 (b) shows that the in-plane Grüneisen parameter are approximately equal. On the other hand, the out-of-plane γ G,33 parameter is larger between 290 and 300 K, where there is greater anharmonicity in the z-direction than in-plane directions. At 310 K, the anharmonicity in all direction are similar and between 320 and 330 K the in-plane direction displays greater anharmonicity than the out-of-plane direction. In Fig. 5, the pair distribution functions of HAP computed at 290 and 300 K suggest that the thermal contraction in bulk HAP is associated with a slightly increase of the inter-atomic distances between the hydroxyl oxygens (c), and between the phosphate oxygens (d).

Negative Grüneisen parameters have been observed in layered crystals, due to membrane effects [START_REF] Lifshitz | Thermal properties of chain and layered structures at low temperatures[END_REF], and in amorphous minerals [START_REF] Fabian | Thermal Expansion and Gr\"uneisen Parameters of Amorphous Silicon: A Realistic Model Calculation[END_REF], and they can be attributed to a lateral compression following a strong expansion in the normal direction; this corresponds to Poisson's compression [START_REF] Abdullaev | Grneisen parameters for layered crystals[END_REF]. Non-monotonous thermal expansion coefficients of HAP-β-tricalcium phosphate ceramics have also been reported in experiments [START_REF] Nakamura | Thermal expansion of hydroxyapatite--tricalcium phosphate ceramics[END_REF].

HAP-water nanopores

Figure 6 reports the heat capacities at constant stress, C σ , of nanoporous HAP-water systems as a function of the pore size H, computed at conditions of 310 K and 1 atm. The continuous red line corresponds to results obtained using the mixture rule:

C σ = f HAP C σ,HAP + (1 -f HAP )C σ,W ,
where C σ,HAP = 906 J.kg -1 .K -1 (see Tab. 2) and C σ,SP C/E = 4513 J.kg -1 .K -1

[59] are the bulk values of bulk HAP and SPC/E water (dashed line) at 310 K, respectively, and f HAP is the volume fraction of solid HAP computed from the average pore size (see Section 2.2):

f HAP = 1 -φ = 1 - H || c || ,
where is φ the porosity.

C w σ,conf = C M D σ -f HAP C σ,HAP (1 -f HAP ) , (9)
(dotted blue line in Fig. 6). The heat capacity of confined water reaches a maximum for the pore size of 6 nm. This non-monotonous variation is closed to the to value of 7 nm measured experimentally in nanoporous silica [START_REF] Etzler | The heat capacity of water in silica pores[END_REF][START_REF] Etzler | Temperature dependence of the heat capacity of water in small pores[END_REF] 200 and, according to Etzler and co-workers [START_REF] Etzler | The structure and properties of vicinal water: Lessons from statistical geometry[END_REF], is due to the bimodal distribution of single molecules enthalpies in liquid water [START_REF] Etzler | The heat capacity of water in silica pores[END_REF][START_REF] Etzler | Temperature dependence of the heat capacity of water in small pores[END_REF].

When compared with the bulk, the excess of the heat capacity of confined water is on the order of 1000 J.kg heat capacity of water in nanoporous HAP are not negligible and, therefore, are 205 critical in relevant upscaling modeling approaches.

The coefficients of thermal expansion of HAP-water systems computed at the conditions of 310 K and 1 atm are presented in Fig. 7 as a function of the pore size H. As observed in nanoporous silica [START_REF] Xu | Thermal Expansion of Confined Water[END_REF], the thermal expansion of confined water also depends on the pore size. Here, we also observe that the thermal of-plan thermal expansion typical of a solid close to 10 -5 K -1 . On the other hand, larger pores exhibit the thermal expansion in the out-of-plane direction governed by the thermal expansion of water. In this case, nanoporous HAP 215 exhibits a out-of-plan thermal expansion typical of a liquid like water close to 10 -4 K -1 . This dependence of the thermal expansion on the pore size is likely to impact the macroscopic (effective) thermal expansion of HAP-rich materials and may be used to explain non-linearities in the thermal expansion observed experimentally in porous HAP [START_REF] Fischer | The lattice thermal expansion of hydroxyapatite[END_REF]. In this section, we present an upscaling strategy to determine the effective heat capacity and effective thermal expansion of porous HAP using as input parameters the heat capacities and the thermal expansion coefficients obtained using MD simulations of HAP-water nanopores.

The heat capacity can be extrapolated from other water contents by a rule of mixture because volume (or mass for specific heat capacities) and energy are extensive thermodynamic quantities. However, the significant dependence of the heat capacity on the pore size raises the issues of which value should be adopted in the upscaling procedure. Figure 8 (a) shows the domain of variability of the effective heat capacity of porous hydroxyapatite as a function of the porosity at 310 K and under 1 atm. The the heat capacities of bulk water and of confined water within the 6 nm HAP nanopore were used to define the lower-and upper limits of this domain.

Up-scaling the coefficient of thermal expansion is more complicated than a simple rule of mixture due to the tensorial nature of this property and the dependence on the microstructure. Analytical homogenization techniques present a consistent framework to upscale various intensive quantities and transformation fields resulting of temperature changes [START_REF] Dvorak | On Transformation Strains and Uniform Fields in Multiphase Elastic Media[END_REF]. Here, we adopt the simple Voigt (α ef f V ) and Reuss (α ef f R ) formulas to estimate the effective coefficient of thermal expansion [START_REF] Dewit | Elastic constants and thermal expansion averages of a nontextured polycrystal[END_REF]:

α ef f V = 1 9K V C p,iikl α p,kl Σ , where K V = 1 9 C iijj Σ , (10)
α ef f R = 1 3 Tr(α p ) Σ , (11)
where the subscript p stands for a phase in the composite and the operator Σ stands for a volume average over the representative elementary volume Σ. Voigt and Reuss estimations assume constant strains and constant stresses, respectively. These estimates correspond to the bounds of the effective properties This corroborates the transferability of HAP empirical force fields proposed by de Leeuw and co-workers [START_REF] De Leeuw | A computer modelling study of the uptake and segregation of fluoride ions at the hydrated hydroxyapatite (0001) surface: introducing a Ca 10 (PO 4 ) 6 (OH) 2 potential model[END_REF][START_REF] De Leeuw | Molecular-dynamics simulation of MgO surfaces in liquid water using a shell-model potential for water[END_REF] to compute the thermal properties of this contraction of HAP at temperatures relevant to in vivo conditions.

Nanoconfinement affects the heat capacity and thermal expansion of HAP-water.

HAP-water nanoporous systems have an effective heat capacity that exceeds the estimates obtained using the heat capacity of bulk water. Our predictions show a non-monotonous behavior for the heat capacity of confined water with a maximum value for a pore of approximately 7 nm, which is in agreement with the heat capacities obtained for other nanoporous materials [START_REF] Etzler | The structure and properties of vicinal water: Lessons from statistical geometry[END_REF]. The differences between the heat capacities of bulk water and of water confined in the 7 nm nanopore is very large (1000 J.kg -1 .K -1 ) and could be the origin of the variability of macroscopic heat capacity observed experimentally in HAP-rich materials.

Similarly, the thermal expansion of confined water exhibited a dependence on the confinement. The coefficients associated with the thermal expansion of HAP-water nanoporous system in the out-of-plane direction are close to the thermal expansion of bulk HAP (close to 10 -5 K -1 ) up to nanopores of 6 nm;

for larger pores the out-of-plane thermal expansion coefficients corresponds to the the thermal expansion of bulk water (close to 10 -4 K -1 ). This dependence on the confinement is expected to impact the macroscopic (effective) thermal expansion of HAP-rich materials.

Nanoscale origin of the variability of thermal properties in porous HAP. Homogenization techniques were used to upscale the heat capacity and the coefficients of thermal expansion of HAP-water nanoporous. The anomalous behavior of confined water leads to significant variability at the macroscopic scale in the thermal properties of HAP materials. This finding indicates that nanoscale processes such as the confinement dependency of water properties are crucial to understand and quantify the variability of properties observed in HAP-rich materials as well as other nanoporous materials. Further, since the confined water properties are pore-size related, experimental measurements of heat capacity and coefficient of thermal expansion provide information on the microand mesopore size distributions in HAP.

Conclusion

The heat capacity and thermal expansion of bulk HAP and HAP-water nanoporous systems have been computed by means of MD simulations. These

Fluctuation-dissipation formulaFigure 1 :

 1 Figure 1: Snapshot of configuration of bulk HAP and HAP-water systems with H = 2 nm and H = 6 nm. For reference, the Cartesian coordinates (x, y, z) and crystallographic (a, b, c) frames are also shown.

Figure 2 :

 2 Figure 2: Convergence of thermodynamic quantities in bulk HAP at 310 K and under 1 atm. At left, the fluctuations of the (a) enthalpy, (b) volume, and (c) deformation along the x (11) and z (33) directions. At right, the variation of the (d) heat capacity, (e) isothermal compressibility (obtained from the fluations of the volume by [32]: κ = V 2 N σT kT V

53 ]J

 53 = (k, l) where the indices i, j, k and l vary in the set {1, ..., 3}. The relation between these indices is the following 1 = (1, 1), 2 = (2, 2), 3 = (3, 3), 4 = (2, 3), 5 = (1, 3) and 6 = (1, 2).

Figure 3 22 and C 33 (Figure 3 :

 322333 Figure3shows the temperature dependence of the heat capacity at constant stress C σ , isothermal compressibility κ, and elastic stiffness constants C 11 , C 22 and C 33 (Voigt notation) of bulk HAP under 1 atm. In the range of temperature considered (290 up to 330 K), which correspond to in vivo conditions, these properties can be treated as constant.On the other hand, in the same temperature range the diagonal elements of thermal expansion tensor, α ii , display a significant temperature dependence,

Figure 4 :

 4 Figure 4: Temperature dependence of (a) the coefficients of thermal expansion and (b) the Grüneisen parameters of bulk HAP under 1 atm. The inset shows the principal values α p,11 , α p,22 and α p,33 of the coefficient of thermal expansion.

- 1 .K - 1 .Figure 5 :Figure 6 :

 1156 Figure 5: Radial distribution functions g of (a) P-P, (b) Ca-Ca, (c) Oh-Oh (hydroxyl oxygen) and (d) Op-Op (phosphate oxygen) pairs in bulk HAP at 290 and 330 K both under 1 atm.The insets show the distance of the peaks in the g profiles of oxygen decreasing with the temperature.

210Figure 7 :

 7 Figure 7: Coefficients of thermal expansion of nanoporous HAP-water systems as a function of the pore size H within in vivo (310 K and 1 atm) conditions.

220 3 . 3 .

 33 Upscaling the Heat Capacity and Thermal Expansion of HAP

Figure 8 :

 8 Figure 8: Estimations of the effective (a) heat capacity at constant stress and (b) the coefficients of thermal expansion of saturated HAP-water systems as a function of the porosity at T = 310 K and P = 1 atm conditions.

  thermal properties are crucial to predict the thermo-mechanical behavior of bones and other HAP-rich materials. The simulations were used to assess the effects of temperature in bulk HAP and of confinement in HAP-water thermal properties. This work contributes towards a better understanding of the factors affecting the thermo-mechanical behavior of bones. The identification of nanoscale phenomenon related to the confinement of water as a source of the property variability can be transferred to the study of other nanoporous adsorbing materials. As body fluid and natural water is rich in ions, further work will consider the role of dissolved ions in confined water and their effect on the macroscopic thermal properties of HAP-electrolyte systems.

Table 1 :

 1 Details of HAP-water systems. The number of atoms in HAP solid layer is N HAP

	= 1584.
	H (nm)	H (nm) N w
	2	1.86	455
	3	2.98	682
	4	3.83	910
	6	6.11	1363
	8	7.97	1823
	10	10.09	2276
	12	12.16	2732
	16	16.73	3650

face projected on z-direction plus twice the van der Waals radius of oxygen 135 (r O =0.152 nm). Correspondingly, the average pore size H N σT is computed by: H = c.z N σT -h HAP N σT . Table

1

displays the number of water molecules N w for each target pore size H.

3. Results

140 3.1. Bulk HAP Table

2

reports the heat capacity C σ , the coefficients of the thermal expansion tensor α ii , and the isothermal compressibility κ of bulk HAP computed from MD simulations at 300 K and compared to theoretical and experimental

Table 2 :

 2 Heat capacity Cσ, coefficients of the tensor of thermal expansion α ii and isothermal compressibility κ, and components of stiffness tensor C 11 , C 22 and C 33 (Voigt notation) of bulk HAP at 300 K and under 1 atm. For comparison, results from the literature concerning, except otherwise indicated, experimental evaluation of hexagonal HAP at 300 K are provided.

material. The anisotropic nature of the thermal expansion of bulk HAP was quantified at temperatures corresponding to in vivo conditions, and negative (principal) coefficients of thermal expansion of bulk HAP were observed for temperature exceeding 310 K. This represents the first observation of thermal