Gonéri Le Cozannet
email: g.lecozannet@brgm.fr

Thomas Bulteau

Bruno Castelle

Roshanka Ranasinghe

Guy Wöppelmann

Jeremy Rohmer

Nicolas Bernon

Déborah Idier

Jessie Louisor

David Salas-Y-Mélia

Quantifying uncertainties of sandy shoreline change projections as sea level rises

Sandy shorelines are constantly evolving, threatening frequently human assets such as buildings or transport infrastructure. In these environments, sea-level rise will exacerbate coastal erosion to an amount which remains uncertain. Sandy shoreline change projections inherit the uncertainties of future mean sea-level changes, of vertical ground motions, and of other natural and anthropogenic processes affecting shoreline change variability and trends. Furthermore, the erosive impact of sea-level rise itself can be quantified using two fundamentally different models. Here, we show that this latter source of uncertainty, which has been little quantified so far, can account for 20 to 40% of the variance of shoreline projections by 2100 and beyond. This is demonstrated for four contrasting sandy beaches that are relatively unaffected by human interventions in southwestern France, where a variance-based global sensitivity analysis of shoreline projection uncertainties can be performed owing to previous observations of beach profile and shoreline changes. This means that sustained coastal observations and efforts to develop sea-level rise impact models are needed to understand and eventually reduce uncertainties of shoreline change projections, in order to ultimately support coastal land-use planning and adaptation.

Since 1870, sea level has been rising, mainly due to the melting of land-ice and ocean expansion caused by anthropogenic climate warming [START_REF] Church | A 20th century acceleration in global sea-level rise[END_REF][START_REF] Church | Sea Level Change[END_REF][START_REF] Dangendorf | Reassessment of 20th century global mean sea level rise[END_REF][START_REF] Dieng | New estimate of the current rate of sea level rise from a sea level budget approach[END_REF][START_REF] Chen | The increasing rate of global mean sea level rise during 1993-2014[END_REF] . While the most immediate impact of sea-level rise is increased coastal flooding hazards [START_REF] Nicholls | Sea level rise and its impact on coastal zones[END_REF][START_REF] Woodworth | Evidence for century-timescale acceleration in mean sea levels and for recent changes in extreme sea levels[END_REF][START_REF] Hallegatte | Future flood losses in major coastal cities[END_REF][START_REF] Hinkel | Coastal flood damage and adaptation costs under 21st century sea level rise[END_REF][START_REF] Sweet | From the extreme to the mean: Acceleration and tipping points of coastal inundation from sea level rise[END_REF][START_REF] Wahl | Understanding extreme sea levels for broad-scale coastal impact and adaptation analysis[END_REF][START_REF] Vitousek | Doubling of coastal flooding frequency within decades due to sea level rise[END_REF][START_REF] Moftakhari | Compounding effects of sea level rise and fluvial flooding[END_REF] , there are significant concerns regarding shoreline retreat as well [START_REF] Nicholls | Sea level rise and its impact on coastal zones[END_REF][START_REF] Leatherman | Sea level rise shown to drive coastal erosion[END_REF][START_REF] Stive | How important is global warming for coastal erosion?[END_REF][START_REF] Moore | Complexities in barrier island response to sea level rise: Insights from numerical model experiments, north carolina outer banks[END_REF][START_REF] Hinkel | A global analysis of erosion of sandy beaches and sea level rise: An application of DIVA[END_REF][START_REF] Cazenave | Sea level rise and its coastal impacts[END_REF][START_REF] Ranasinghe | Assessing climate change impacts on open sandy coasts: A review[END_REF][START_REF] Vitousek | Can beaches survive climate change?[END_REF][START_REF] Jiménez | Impacts of sea level rise-induced erosion on the Catalan coast[END_REF][START_REF] De Winter | Sensitivity analysis of climate change impacts on dune erosion: case study for the Dutch Holland coast[END_REF][START_REF] Enríquez | Changes in beach shoreline due to sea level rise and waves under climate change scenarios: application to the Balearic Islands (western Mediterranean)[END_REF][START_REF] Monioudi | Assessment of island beach erosion due to sea level rise: the case of the Aegean archipelago (Eastern Mediterranean)[END_REF][START_REF] Luijendijk | The State of theWorld's Beaches[END_REF][START_REF] Mentaschi | Global long-term observations of coastal erosion and accretion[END_REF] . In particular, beaches backed by sandy deposits are receiving particular attention for several reasons: first, they represent 31% of the world's ice-free coasts [START_REF] Luijendijk | The State of theWorld's Beaches[END_REF] ; second, they are potentially highly sensitive to sea-level changes [START_REF] Leatherman | Sea level rise shown to drive coastal erosion[END_REF][START_REF] Stive | How important is global warming for coastal erosion?[END_REF][START_REF] Moore | Complexities in barrier island response to sea level rise: Insights from numerical model experiments, north carolina outer banks[END_REF][START_REF] Monioudi | Assessment of island beach erosion due to sea level rise: the case of the Aegean archipelago (Eastern Mediterranean)[END_REF] ; third, it has been estimated that at least 24% and up to 70% of the world beaches are already under chronic erosion, although with large regional and local differences [START_REF] Luijendijk | The State of theWorld's Beaches[END_REF][START_REF] Bird | Coastline Changes. A Global Review[END_REF] ; finally, beaches are both valuable for tourism and as buffer zones during extreme events such as storms.

Sandy shoreline change projections along a given coast need to consider sediment losses and gains caused by a number of hydro-sedimentary processes acting at different timescales [START_REF] Cowell | The coastal-tract (part 2): applications of aggregated modeling of lower-order coastal change[END_REF] . Coastal change is driven by a myriad of processes interacting with one another at a wide range of spatial and temporal scales, making the use of process-based coastal evolution modelling difficult at operational levels. Hence, instead of attempting to quantify the impact of each process causing sediment transport, coastal adaptation practitioners have pragmatically relied on extrapolations of past observations in order to anticipate future shoreline changes. In the absence of human interventions, estuaries or other major sediment sources or sinks, this results in the following equation:

∆ = + . + + L L n Tx Lvar S (1)
r r slr 0 where L r0 and L r are the shoreline positions in a cross-shore direction at the initialization and after n years, respectively, Tx is the linear trend over multi-decadal to centennial timescales, and Lvar represents the seasonal, inter-annual and decadal modes of variability of shoreline changes. Neither Tx nor Lvar is related to a single physical process: for example, the effects of longshore gradients in sediment transport are included in Lvar when they result in evolution and cycles occurring at interannual to decadal timescales, while there are included in Tx if their manifestation involves longer timescales [START_REF] Antonilez | A multiscale climate emulator for long-term morphodynamics (MUSCLE-morpho)[END_REF] . Finally, the term ∆S slr quantifies the impacts of relative sea-level changes, that is, those due to ocean thermal expansion, land water and ice contributions, and vertical ground motions (uplift or subsidence) [START_REF] Stammer | Causes for contemporary regional sea level changes[END_REF] . This is denoted "coastal impact model" hereafter.

In the area of coastal prospective planning, a comprehensive description of the uncertainties of future shoreline positions is required in order to avoid maladaptation [START_REF] Magnan | Addressing the risk of maladaptation to climate change[END_REF] . Shoreline change projections based on equation (1) need to consider the uncertainties of future sea-level rise under different climate forcings [START_REF] Church | Sea Level Change[END_REF][START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] , of local vertical ground motions [START_REF] Wöppelmann | Vertical land motion as a key to understanding sea level change and variability[END_REF] , of other modes of variability of shoreline positions (Tx and Lvar), as well as uncertainties of coastal impact models (∆S slr ). The latter source of uncertainty has not been quantified so far.

We estimate the uncertainties of coastal impact models (∆S slr ) by considering the difference between two existing approaches: the Bruun rule [START_REF] Bruun | Sea level rise as a cause of shore erosion[END_REF] and the Probabilistic Coastline Recession model [START_REF] Ranasinghe | Estimating coastal recession due to sea level rise: beyond the Bruun rule[END_REF][START_REF] Li | Probabilistic modelling of extreme storms along the Dutch coast[END_REF][START_REF] Wainwright | Moving from deterministic towards probabilistic coastal hazard and risk assessment: Development of a modelling framework and application to Narrabeen Beach, New South Wales, Australia[END_REF][START_REF] Dastgheib | Regional scale risk-informed land-use planning using probabilistic coastline recession modelling and economical optimisation: East coast of Sri Lanka[END_REF][START_REF] Verschuur | Implications of uncertain Antarctic ice sheet dynamics for managing future coastal erosion: a probabilistic assessment[END_REF] . The Bruun rule is the most commonly used and historical approach to assess sea-level rise impacts on shorelines [START_REF] Bruun | Sea level rise as a cause of shore erosion[END_REF][START_REF] Davidson-Arnott | Conceptual model of the effects of sea level rise on sandy coasts[END_REF][START_REF] Karunarathna | Multi-timescale morphological modelling of a dunefronted sandy beach[END_REF] . It assumes a landward translation of the beach profile as sea level rises. The Probabilistic Coastline Recession model (PCR) is a recently introduced approach that quantifies sediment losses at the dune toe during storms, as well as the nourishment of the dune by aeolian sediment transport processes between storms [START_REF] Ranasinghe | Estimating coastal recession due to sea level rise: beyond the Bruun rule[END_REF] . Over multi-decadal timescales, the superimposition of unchanged storms with rising mean sea levels results in more frequent and larger sediment losses in the PCR model. The Brunn rule and the PCR models are not only based on different assumptions regarding the physical processes guiding the response of sandy shorelines to sea-level rise, but they also provide different results [START_REF] Ranasinghe | Estimating coastal recession due to sea level rise: beyond the Bruun rule[END_REF][START_REF] Toimil | Managing coastal erosion underclimate change at the regional scale[END_REF] . Today, both coastal impact models are equally difficult to validate due to the scarcity of coastal data and the complexity of the hydrosedimentary processes involved [START_REF] Cooper | Sea level rise and shoreline retreat: time to abandon the Bruun Rule[END_REF][START_REF] Cowell | Management of uncertainty in predicting climate-change impacts on beaches[END_REF] . Faced with this structural uncertainty [START_REF] Rougier | Model and data limitations: the sources and implications of epistemic uncertainty[END_REF] , stakeholders concerned with coastal adaptation generally lack the relevant observations to validate one particular model, and may therefore assign an equal confidence to both, following the maximum entropy principle [START_REF] Mishra | Assigning Probability Distributions to Input Parameters of Performance Assessment Models[END_REF] . Hence, we consider the difference between the Bruun and PCR models as a first order measure of the coastal impact model error.

Equation (1) requires coastal data, which are not available at the global scale. Hence, we implement the approach in the sandy coast of Aquitaine in southwestern France (Supplementary Materials 1, 2 and 3), where observations of shoreline positions [START_REF] Castelle | Spatial and temporal patterns of shoreline change of a 280-km high-energy disrupted sandy coast from 1950 to 2014: SW France[END_REF] and of other metrics allow to minimize uncertainties caused by the lack of data, thus avoiding mixing them with uncertainties due to the variability of the different processes (Supplementary Material 1, section 1). We select four coastal sites because they are representative of the broad range of variability and trends in shoreline positions in the sandy coast of Aquitaine (Supplementary Materials 1, 3 and 4): sites #1 and #2 are stable today, whereas sites #3 and #4 are eroding. The variability around this trend reaches +/-21 m at site #3, approximately +/-9 m at site #2 and #4 and only +/-3.2 m at site #1. Furthermore, site #2 is located close to a permanent Global Navigation Satellite System (GNSS) station, providing a precise evaluation of vertical ground motions (subsidence or uplift) and of its contribution to relative sea level changes. The uncertainties of vertical ground motions at site #2 are estimated from the SONEL database [START_REF] Santamaria-Gomez | Uncertainty of the 20th century sea level rise due to vertical land motion errors[END_REF] , while those at sites #1, #3 and #4 are based on a global analysis of coastal vertical ground motions based on GNSS measurements and Global Isostatic Adjustment (GIA) modeling. Hence, site #2 provides a testbed to appraise the usefulness of GNSS instrumentation in reducing uncertainties of shoreline change reconstructions and projections.

We compute shoreline change projections for the four coastal sites in Aquitaine by means of a Monte-Carlo procedure, within which probabilistic input parameters (Supplementary Material 2) are propagated through equation (1) [START_REF] Vitousek | Can beaches survive climate change?[END_REF][START_REF] Anderson | Doubling of coastal erosion under rising sea level by mid-century in Hawaii[END_REF] . We estimate past regional sea-level changes and their uncertainties using tide-gauge records and permanent GNSS stations providing estimates of coastal vertical ground motions [START_REF] Santamaria-Gomez | Uncertainty of the 20th century sea level rise due to vertical land motion errors[END_REF] (Fig. 1). We use the future sea-level rise projections provided by Kopp et al. [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] that we correct from local vertical ground motions (Fig. 1; see Methods, subsection 1). These probabilistic sea-level projections [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] are essentially consistent with the 5 th Assessment Report of the International Panel on Climate Change (IPCC) [START_REF] Church | Sea Level Change[END_REF] , which remains today the reference source of information from the perspective of coastal adaptation practitioners (see the discussion section). For example, for RCP 8.5 by 2100, the 17 th and 83 rd quantiles of the probabilistic projections of Kopp et al. [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] indicate a global mean sea level rise of 0.62 to 1 m. This compares well with the likely range of the IPCC projections for the same time horizon (0.52 to 0.98 m), corresponding to a 66-100% probability according to the IPCC uncertainty language, or roughly 66% probability according to a later communication of the IPCC authors [START_REF] Church | Sea-level rise by 2100[END_REF] . The other terms in equation (1) (Tx and Lvar) are evaluated empirically, based on past observations of L r obtained through various data (e.g. historical charts, satellite and aerial photographs, shoreline surveys) available from the Aquitaine Coastal Observatory (Supplementary Materials 1 and 2) for each coastal impact model.

We evaluate the uncertainties of coastal impact models against all other sources of uncertainties using a global sensitivity analysis [START_REF] Saltelli | Global sensitivity analysis: the primer[END_REF] (see Methods, subsection 3). This latter procedure provides quantitative estimates of the contribution of each uncertain parameter to the variability of shoreline change projections, while addressing explicitly interactions between uncertain parameters involved in equation (1). The global sensitivity analysis involves a large number of simulations: the PCR model would be run approximately 100,000 times each year to reduce errors of sensitivity indices to less than 1%. As the computation time of the PCR model is prohibitive for such a number of simulations, we use a surrogate model of the PCR model to reduce the computation time (see Methods, subsection 4 and Supplementary Material 5).

Our probabilistic approach delivers a range of contrasting shoreline change projections. We show that by 2100, depending on the location of interest, the statistical uncertainties of sea-level projections account for up to approximately half of that due to the choice of a particular coastal impact model (section 2). However, there remain residual uncertainties, the probability of which is difficult to quantify, and which could bring future shoreline positions outside the range of probabilistic projections shown in section 2 (see section 3). Despite these residual uncertainties, this study shows that the uncertainties of coastal impact models should be considered in future global studies aiming at quantifying sea-level rise impacts and their uncertainties.

Results

Shoreline change projections. This subsection presents the shoreline change reconstructions (1807 to 2014) and projections (2000 to 2200) obtained from the propagation of probabilistic uncertainties through equation (1) at the coastal site #1 (for other sites, see Supplementary Materials 6 to 11). We provide shoreline change projections and reconstructions based on equation (1) implementing the Bruun rule (Fig. 2) and the surrogate PCR model (Fig. 3). For consistency, the effects and uncertainties of past sea-level changes have been included in past shoreline change positions as assumed in equation 1, although they are expected to remain small [START_REF] Castelle | Spatial and temporal patterns of shoreline change of a 280-km high-energy disrupted sandy coast from 1950 to 2014: SW France[END_REF][START_REF] Cozannet | Approaches to evaluate the recent impacts of sea level rise on shoreline changes[END_REF] .

At each time-step from 1807 to 2200, Figs 2 and 3 provide the median shoreline positions along a cross-shore transect. In addition, the uncertainties of shoreline change projections are conveyed through percentile levels. Thus, the probability for a shoreline position to remain between the 17 th -83 rd percentile levels (most intense gray, green, blue or red colors) is 66%. Note that according to these projections, there is a small probability that shorelines advance seaward in the future. This is because sea-level rise projections used here [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] neither exclude a decrease in the rate of sea-level rise after some decades, nor a drop in sea-level rise itself in the case of RCP 4.5 and RCP 2.6 (Fig. 1). The uncertainties of the shoreline change projections are very large (50 to 100 m by 2200), especially for the simulations based on the Bruun rule (Fig. 2). These uncertainties are analyzed in the following section by means of the global sensitivity analysis.

All median projections shown in Figs 2 and 3 involve shoreline erosion, which corresponds to the superimposition of a multi-decadal eroding trend at site #1 with the effects of sea-level rise. Retreat values are larger for scenarios with highest radiative forcing: median values reach 40 m and 100 m for RCP 4.5 and RCP 8.5 by 2100 with the projections based on the Bruun rule. Projections based on PCR are less sensitive to sea-level rise, so that shoreline retreat values and the related statistical uncertainties are approximately five times smaller in Fig. 3 compared to Fig. 2 (note similar order of magnitudes in Supplementary Materials 6 to 11 for sites #2, #3 and #4). This agrees well with previous studies showing that by 2100, the Bruun estimate lies in the range of 4-40% exceedance probability with respect to the corresponding approach based PCR estimates [START_REF] Ranasinghe | Estimating coastal recession due to sea level rise: beyond the Bruun rule[END_REF][START_REF] Toimil | Managing coastal erosion underclimate change at the regional scale[END_REF] .

For coastal adaptation practitioners, such results mean that at similar coastal sites where buildings and infrastructures are located close to the coast, the modeling approach based on the Bruun rule implies large shoreline retreats for high emission scenarios in a few decades from now, superimposing a trend of 0.5 m/yr or more to current shoreline retreat rates. Such a change would be associated with expensive coastal protection measures or the relocation of numerous assets at risk. Conversely, RCP 2.6 and the PCR model generate less erosive trends, so that no specific adaptation to sea-level rise-induced shoreline erosion would be required. In all cases, a strong increase in shoreline erosion is not expected before 2050, which gives time to plan adaptation. [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] ; Reconstruction: see Methods, subsection 1). 1). This index is commonly interpreted as the expected proportion of the total variance of the shoreline change that would be removed if we were able to learn the true value of the uncertain parameter. It is commonly used to rank the importance of model parameters according to their impact for the variability of the model outcome [START_REF] Saltelli | Global sensitivity analysis: the primer[END_REF] .

Variance-based global sensitivity analysis.

The variance of shoreline change projections is partitioned differently depending on the period of time considered (Fig. 4): uncertainties in the long-term trends (not including sea-level rise, that is, Tx in equation ( 1)) have a larger impact in the past (before 1950), whereas for present days, uncertainties are overwhelmingly due to the current modes of shoreline change variability, associated with the random nature of waves and currents at timescales ranging from days to decades through seasons (Lvar in equation ( 1)). Whatever the coastal site and the period of time considered, the uncertainties of the regional sea-level reconstructions can be neglected with little impact to shoreline change projections, and so can the uncertainties of shoreface slopes estimations. Where vertical ground motions are known precisely from a permanent GNSS station, their impact on shoreline change projections is negligible (site #2), whereas they can account for up to almost 20% of the variance by 2050 at site #1. Hence, while all sites considered are exposed to quite similar waves and tide conditions, the relative importance of each source of uncertainty varies along the coast depending on the local processes and the knowledge available.

If we knew which coastal impact model was correct, the variance of shoreline change projections could be reduced by 20 to 40% by 2100 depending on the coastal site of interest. Indeed, Fig. 4 shows that for sites #1, #2 and #4, the impact of this source of uncertainty grows rapidly over the coming decades, then peaks at 40% during the second half of the 21 st century, and finally decreases after 2100, as sea-level projections are becoming more uncertain and account for a larger fraction of the variance of shoreline change projections. For site #3, the uncertainties of shoreline change trends and variability are much larger, thus reducing the fraction of variance that can be attributed to the choice of a coastal impact model. However, whatever the coastal site considered, the uncertainty due to coastal impact models accounts for twice the impact of the uncertainties of probabilistic sea-level projections (yellow area in Fig. 4).

It could be argued that although the Bruun rule and the PCR models are two common alternatives, they could both be wrong. However, if the true response of shorelines to future sea level rise is outside the range predicted These projections include uncertain shoreface slopes, vertical ground motions, sea-level changes, shoreline and change variability from event-scale to inter-decadal timescales as well as an uncertain multi-decadal trend (see Methods). Note that sea-level projections used here [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] consider that a sea-level drop is possible (although very unlikely) beyond 2100 (Fig. 1). Hence, shoreline changes below the median can bend downwards beyond 2100. by the two coastal impact models, the relative importance of this structural uncertainty will even be greater than estimated in Fig. 4.

To summarize, the main sources of uncertainties of shoreline reconstructions and projections depend on the local coastal settings and vary with time, with future sea-level rise and the choice of a particular coastal impact model becoming prominent beyond 2050. This statement relies on state-of-the-art probabilistic uncertainty quantification exercises, which, nevertheless, still contains some inherent limitations (see section 3). It is more than likely that our understanding of uncertainties will evolve as new knowledge becomes available, leading to new, better assessments and rankings of uncertainty sources. Meanwhile, the different sources of uncertainties shown in Fig. 4 can be addressed differently, either by trying to reduce them using new observations wherever possible (e.g., of shoreline changes and vertical ground motions) [START_REF] Cazenave | Monitoring coastal zone changes from space[END_REF] , or by selecting appropriate decision-making frameworks to minimize the risks of maladaptation despite uncertainties [START_REF] Magnan | Addressing the risk of maladaptation to climate change[END_REF][START_REF] Hinkel | Sea level rise scenarios and coastal risk management[END_REF] .

Discussion: Residual Uncertainties

Shoreline change projections presented here assume that the uncertainties of each input parameter can be described using probability distributions (based on observations in Supplementary Material 2). The uncertainties of shoreline change projections compare well with those of idealized energetic sandy shorelines exposed to sea-level changes close to the global average [START_REF] Cozannet | Uncertainties in sandy shorelines evolution under the Bruun rule assumption[END_REF] . Note that the majority of coasts are expected to be affected by sea-level changes close to or slightly higher than the global average [START_REF] Carson | Coastal sea level changes, observed and projected during the 20th and 21st century[END_REF] . However, in addition to the uncertainties presented in Fig. 4, there are additional unknowns due to a lack of confidence in the probability distributions themselves [START_REF] De Vries | How to interpret expert judgment assessments of 21st century sea level rise[END_REF][START_REF] Bakker | Sources and implications of deep uncertainties surrounding sea level projections[END_REF][START_REF] Stephens | Applying principles of uncertainty within coastal hazard assessments to better support coastal adaptation[END_REF] . These unknowns are called residual uncertainties hereafter, following the terminology of Robinson et al. [START_REF] Robinson | Accounting for Residual Uncertainty: Updating the Freeboard Guide[END_REF] . They are listed in Table 1.

Sea-level rise projections presented here implicitly assign a very low probability to Antarctic marine ice-cliff instabilities [START_REF] Deconto | Contribution of Antarctica to past and future sea level rise[END_REF] , although the probability that such instabilities occur during the 21 st century remains unquantified so far [START_REF] Kopp | Evolving Understanding of Antarctic Ice-Sheet Physics and Ambiguity in Probabilistic Sea-Level Projections[END_REF] . Sea-level projections considering marine ice-cliff instabilities are significantly higher than those used in this study, as their median attains 1.8 m as early as 2100 [START_REF] Kopp | Evolving Understanding of Antarctic Ice-Sheet Physics and Ambiguity in Probabilistic Sea-Level Projections[END_REF][START_REF] Bars | A high-end sea level rise probabilistic projection including rapid Antarctic ice sheet mass loss[END_REF] . Our results reflect the fact that according to the sea-level projections of Kopp et al. [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] and those provided by the IPCC [START_REF] Church | Sea Level Change[END_REF] , future sea-level rise remains a slow process characterized by a large inertia and involving multi-centennial timescales [START_REF] Clark | Consequences of twenty-first-century policy for multi-millennial climate and sea-level change[END_REF] . Replacing the current reference projections 2,32 by those derived from DeConto and Pollard 61 would change coastal impact assessments [START_REF] Wong | Impacts of Antarctic fast dynamics on sea level projections and coastal flood defense[END_REF] , so that this residual source of uncertainty would deserve specific attention from an adaptation perspective. Furthermore, the selection of a particular climate change scenario has little impact to future shoreline positions within the probabilistic framework presented in this study (brown area in Fig. 4). This reflects the fact that sea level projections used here [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] are not differentiated enough according to climate change scenarios to create large differences in shoreline change projections given the uncertainties of other processes (Fig. 1). In these projections, the dynamic contribution of ice in Antarctica is largely independent from greenhouse gas emissions due to the lack of understanding of processes taking place today [START_REF] Ritz | Potential sea level rise from Antarctic ice-sheet instability constrained by observations[END_REF] . However, sea level rise will continue for millennia even if we keep climate warming below the target of the Paris agreement, and the rates of sea level changes over the next centuries will strongly depend on current mitigation policies [START_REF] Clark | Consequences of twenty-first-century policy for multi-millennial climate and sea-level change[END_REF][START_REF] Rockström | A roadmap for rapid decarbonization[END_REF] .

Vertical land motions are subject to large residual uncertainties as well: indeed, the analysis of eleven years of continuously recorded GNSS data at Cap Ferret suggests that site #2 is affected by a subsidence of -1.2 ± 0.6 mm/yr [START_REF] Santamaria-Gomez | Uncertainty of the 20th century sea level rise due to vertical land motion errors[END_REF] , which is qualitatively in agreement with the independent estimate from supplementary levelling data (Source: Aquitaine Coastal Observatory). However, it is unsure that the pointwise information of the permanent GNSS station located at Cap Ferret is representative of the nearby area. White areas indicate interactions between parameters, corresponding to shoreline positions, which can be only reached if at least two uncertain parameters deviate from their mean (see Methods). The graph reads as follows: for site #1, by 2200, uncertainties in regional sea-level rise projections (yellow) account for approximately 30% of the variance of shoreline change projections.

Source of uncertainty Uncertainties quantification in this study

Residual uncertainties (not quantified in this study)

Future sea-level rise Probabilistic regional sea-level rise projections [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] Possibility of rapid melting of ice-sheets [START_REF] Deconto | Contribution of Antarctica to past and future sea level rise[END_REF][START_REF] Wong | Impacts of Antarctic fast dynamics on sea level projections and coastal flood defense[END_REF] Vertical ground motions Geodetic uncertainties at the permanent GNSS at Cap Ferret [START_REF] Santamaria-Gomez | Uncertainty of the 20th century sea level rise due to vertical land motion errors[END_REF] Representativeness of the GNSS records (linear extrapolations in time and space)

Shoreline change variability and trends

Seasonal, interannual and multi-decadal (~50 years) shoreline evolutions Structure of equation ( 1); Possible alteration of natural variability modes and or self-organized patterns [START_REF] Coco | Patterns in the sand: From forcing templates to self-organization[END_REF] Uncertainties of the shoreface and upper shoreface slopes

Observed variability of shoreface and upper shoreface slopes

The actual variability may be larger and the slopes may exceed their current variability in the future

Shoreline evolution modeling framework Structural uncertainties (2 different modeling approaches) Limited confidence in both modeling frameworks

Uncertainties of the PCR model Statistical uncertainties, reflecting the variability of the response of the PCR model to virtual time series of events (Supplementary Material 1 Section 2; Supplementary Fig. 5).

Structural uncertainties: choices made in the different modules of the PCR model, especially the dune erosion module.

Uncertainties of the Bruun rule

Observed variability of shoreface slopes Limitations of the Bruun rule [START_REF] Cooper | Sea level rise and shoreline retreat: time to abandon the Bruun Rule[END_REF][START_REF] Ranasinghe | Rising seas and retreating coastlines[END_REF] Table 1. Residual uncertainties of probabilistic shoreline change projections.

The terms in equation ( 1) contain some uncertainty as they are estimated based on limited observations of shorelines, and foreshore/upper foreshore slopes (Table 1 and Supplementary Material 1). However, the structure of equation ( 1) itself can be discussed because it assumes that only sea-level rise will modify the sediment budgets, that the system will remain mostly unaffected by human interventions, and that potential changes of other factors such as waves remain within the error bars of Tx and Lvar. As noted in the introduction, an alternative would be to aggregate the contribution of each process known to cause net multi-decadal shoreline changes. These would include the effects of persistent longshore sediment transport gradients, decreased inner shelf sediment supply to the shoreface, sediment losses to the submarine canyon south of Aquitaine and where the longshore drift terminates, and aeolian transport (the impact of reduced sediments from rivers is believed to be negligible in the 4 sites investigated here). In Aquitaine, as in many coastal sites worldwide, these processes remain poorly quantified, suggesting that the actual uncertainties associated with Tx might be larger than those presented in Supplementary Material 3. However, it should be noted that the retrospective simulations of past shoreline positions in the 19 th century (in gray), which are known from historical maps, are very sensitive to errors in Tx. If we had used erroneous multi-decadal trends, this would results in past shoreline positions hardly reconcilable with this independent piece of evidence. This makes us confident that the multi-decadal trends used in this work are close to the ground truth. Furthermore, current modes of shoreline change variability could be altered by climate change [START_REF] Charles | Climate change impact on waves in the Bay of Biscay, France[END_REF][START_REF] Vousdoukas | Projections of extreme storm surge levels along Europe[END_REF] . However, climate change models do not indicate that a modification of storm patterns (intensity, trajectories) should be expected in the Bay of Biscay, so that the impacts on wave regimes and surges are expected to remain small compared to those of sea-level rise. Hence, we believe that using equation ( 1) with the values of Tx and Lvar presented in this study is the best possible approach given the present knowledge on the study site. Rare events such as the unusual 2013/2014 sequence of winter storms suggest that more research is needed in this area [START_REF] Castelle | Impact of the winter 2013-2014 series of severe western Europe storms on a double-barred sandy coast: Beach and dune erosion and megacusp embayments[END_REF][START_REF] Masselink | Extreme wave activity during 2013/2014 winter and morphological impacts along the Atlantic coast of Europe[END_REF] .

Coastal impact model uncertainties are quantified here using the differences between two well-established models. However, there is no guarantee that this metrics is appropriate to appraise the real variability of possible impacts of sea-level rise on sandy shorelines. Among the two models, the PCR model results, owing the inherent probabilistic approach adopted in the model, incorporate the statistical uncertainty reflecting the variability of profile response to virtual time series of storm events. However, structural uncertainties, due to choices made during the implementation of the PCR model itself, are not quantified. These include uncertainties due to choices made while generating virtual time series of events (Supplementary material 4, section 3), while computing total water level estimates at the coast [START_REF] Melet | Under-estimated wave contribution to coastal sea-level rise[END_REF] , and when computing dune recovery and erosion [START_REF] Larson | An analytical model to predict dune erosion due to wave impact[END_REF] . Based on our knowledge in Aquitaine and on our experience in developing and applying the PCR model elsewhere in the world, we feel that the dune erosion module is potentially the largest source of structural uncertainties [START_REF] Rougier | Model and data limitations: the sources and implications of epistemic uncertainty[END_REF] of the PCR model outcome in our particular application. In particular, the dune erosion module [START_REF] Larson | An analytical model to predict dune erosion due to wave impact[END_REF] assumes that the dune toe evolves parallel to the slope of the upper shoreface, which remains within the same range over time.

In this article, we concentrated the analysis on quantifiable uncertainties (using the term introduced by Walker et al. [START_REF] Walker | Defining uncertainty: a conceptual basis for uncertainty management in model-based decision support[END_REF] ) without attempting to quantify residual uncertainties. However, we not only treated parametric uncertainties (i.e. uncertainties related to the choice of the values in the model parameters), but also model uncertainties, i.e. structural uncertainties (uncertainties related to the selection in the most appropriate model structure given our problem). The latter uncertainty source is, to our best knowledge, rarely addressed in the literature and falls usually under the category of residual uncertainties.

Conclusion

Despite limitations, our results show that coastal managers in charge of identifying and implementing appropriate adaptation decisions on sandy beaches need to consider uncertainties of sea-level rise together with those of coastal impact models. In fact, the structure of the two coastal impact models available today implies that slopes used to quantify the impacts of sea level rise may vary from the upper foreshore slopes (as in the surrogate PCR model) to the entire foreshore slopes, i.e., measured over a transect starting at the "depth of closure" and ending at the dune toe. This means that for most of the beaches worldwide, structural uncertainties due to the selection of the coastal impact model can be expected to have large impacts in shoreline change projections.

This statement is especially valid at multi-decadal to centennial timescales, which are the most relevant for land-use planning. Hence, our results call for more research in the area of coastal impact models in order to support coastal adaptation. Over the past decades, sea-level rise has been accelerating [START_REF] Dieng | New estimate of the current rate of sea level rise from a sea level budget approach[END_REF][START_REF] Chen | The increasing rate of global mean sea level rise during 1993-2014[END_REF] , and a further acceleration is expected to take place without mitigation of climate change. At the same time, coastal evolution modeling is progressing toward appropriate complexity approaches assimilating coastal observations [START_REF] Vitousek | Can beaches survive climate change?[END_REF][START_REF] French | Appropriate complexity for the prediction of coastal and estuarine geomorphic behaviour at decadal to centennial scales[END_REF] . We speculate that the impacts of sea-level rise on sandy shorelines should be increasingly observable in the coming decades, so that global, long term, repeated, accurate and precise observations of shoreline positions will be especially relevant to support the development of more trustworthy coastal impact models.

Methods

Sea-level reconstructions and projections. Sea-level projections are those of Kopp et al. [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] (Fig. 1). These projections provide the probability of future sea-level rise at La Rochelle. By analyzing tide gauge records, Kopp et al. [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] estimate a subsidence at La Rochelle to be -0.55 +/-0.52 mm/yr. However, this subsidence is within the error bars. Consequently, the null hypothesis of a stable location cannot be ruled out at the 95% confidence level. Furthermore, the permanent GNSS station co-located with the La Rochelle tide gauge indicates that the tide gauge is stable too [START_REF] Santamaria-Gomez | Uncertainty of the 20th century sea level rise due to vertical land motion errors[END_REF] . Hence, we removed the subsidence estimation in the projections of Kopp et al. [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] . We assume that the probabilistic sea-level rise projections at La Rochelle are applicable in Aquitaine once corrected from the local vertical ground motions.

To reconstruct past sea levels in the Bay of Biscay, we averaged yearly records of 15 stations available in the Permanent Service for Mean Sea Level (PSMSL) corrected for local vertical land motions. We rejected 5 additional stations because they displayed anomalous trends (Pointe Saint Gildas, Le Verdon, Pasajes, Santander 2 and Gijon 2). For each station, the local rates of sea-level rise and their uncertainties were computed using a forward-backward Kalman filter [START_REF] Corriou | Commande des procédés[END_REF][START_REF] Visser | A review of trend models applied to sea level data with reference to the "accelerationdeceleration debate[END_REF] as in Rohmer and Le Cozannet [START_REF] Rohmer | Dominance of the mean sea level in the high-percentile sea levels time evolution with respect to largescale climate variability: a Bayesian statistical approach[END_REF] . Local vertical land motions were obtained either using permanent GNSS stations from the SONEL database [START_REF] Santamaria-Gomez | Uncertainty of the 20th century sea level rise due to vertical land motion errors[END_REF] , or, in the absence of GNSS station, using global isostatic adjustment models [START_REF] Peltier | Global glacial isostasy and the surface of the ice-age Earth: the ICE-5G (VM2) model and GRACE[END_REF] available at tide gauge in Jevrejeva et al. [START_REF] Jevrejeva | Nonlinear trends and multiyear cycles in sea level records[END_REF] . Where no GNSS information was available, where the GNSS station was located too far from the tide gauge, or where the GNSS records displayed large step discontinuities or a clear non-linear behavior, we assigned a Gaussian uncertainty of 2 mm/yr (1-sigma) to the vertical land motion value, as obtained from the histogram frequency distribution of all trends in the SONEL database corrected from the global isostatic adjustment. Note that larger subsidence or uplift values would be detectable in the tide gauge sea-level time series. Finally, we computed past mean sea levels and their uncertainties using a weighted least square regression. The method used here assumes that vertical land motions are linear over the timescales considered, and that all stations measure the same signal of mean sea level in the Bay of Biscay.

We assume that vertical land motions at the site #2 are those measured by the GNSS station of the Cap-Ferret. In other sites, as no permanent GNSS with sufficiently long records are available, we modelled the uncertainties due to vertical ground motions through a centered normal distribution with a standard deviation of 2 mm/yr, as obtained from the histogram analysis of all trends computed from the coastal permanent GNSS stations in the SONEL database [START_REF] Wöppelmann | Vertical land motion as a key to understanding sea level change and variability[END_REF] .

Coastal impacts models. The Bruun rule quantifies the shoreline retreat in response to sea-level rise as follows [START_REF] Bruun | Sea level rise as a cause of shore erosion[END_REF] :

∆ ξ α = ∆ S tan( ) (2)
slr where ξ ∆ is the cumulated rise of sea level and α tan( ) is the foreshore slope from the depth of closure to the top of the upper shoreface, usually close to 1% [START_REF] Nicholls | Assessing erosion of sandy beaches due to sea level rise[END_REF] .

We implement a variant of the PCR model adapted to the Aquitaine coast and present the results for site #2. Following Ranasinghe et al. [START_REF] Ranasinghe | Estimating coastal recession due to sea level rise: beyond the Bruun rule[END_REF] , we use the Larson et al. [START_REF] Larson | An analytical model to predict dune erosion due to wave impact[END_REF] formula to compute sediment losses during storms. Furthermore, we parametrize the dune recovery rate in the PCR model using observations in order to minimize this source of uncertainty [START_REF] Ranasinghe | Estimating coastal recession due to sea level rise: beyond the Bruun rule[END_REF][START_REF] Li | Probabilistic modelling of extreme storms along the Dutch coast[END_REF][START_REF] Wainwright | Moving from deterministic towards probabilistic coastal hazard and risk assessment: Development of a modelling framework and application to Narrabeen Beach, New South Wales, Australia[END_REF] : we estimate gains at the dune toe under calm weather conditions at 25 ± 15 m 3 /m/ year, based on an analysis of 9 dune profile campaigns of the Aquitaine Coastal Observatory from 2007 to 2015. Finally, we incorporate tides as suggested in Larson et al. [START_REF] Larson | An analytical model to predict dune erosion due to wave impact[END_REF] . The model inputs are virtual time series of events (surge level, significant wave height, peak period and peak direction, tidal levels, event duration, spacing between events) capturing the statistical dependence between variables, their seasonality and event grouping, which are generated using waves and surges hindcasts [START_REF] Charles | Climate change impact on waves in the Bay of Biscay, France[END_REF][START_REF] Paris | des surcotes-décotes dans le Golfe de Gascogne. 1 ère édition des Journées REFMAR les 18-19juin 2013[END_REF] (see Supplementary Material 1, section 2).

Global sensitivity analysis.

The Variance-Based Sensitivity Analysis (VBSA) quantifies the contribution of input variables and parameters to the variance of the outcomes of a model 51 (see Chu-Agor et al. [START_REF] Chu-Agor | Exploring vulnerability of coastal habitats to sea level rise through global sensitivity and uncertainty analyses[END_REF] ; Wong and Keller [START_REF] Wong | Deep Uncertainty Surrounding Coastal Flood Risk Projections: A Case Study for New Orleans[END_REF] , Le Cozannet et al. [START_REF] Cozannet | Uncertainties in sandy shorelines evolution under the Bruun rule assumption[END_REF][START_REF] Cozannet | Evaluating uncertainties of future marine flooding occurrence as sea level rises[END_REF] for applications of global sensitivity analysis in the area of coastal impacts of sea-level rise). Let us define f as the model computing future shoreline change positions. Considering the n-dimensional vector X as a random vector of independent random variables X i (i = 1,2, …, n) (Table 1), then the output Y = f(X) is also a random variable (as a function of a random vector). VBSA determines the part of the total unconditional variance Var(Y) of the output Y resulting from each input random variable X i . The partial and total variances of Y are assessed based on the functional analysis of variance decomposition of f 86 , into summands of increasing dimension (provided that f can be integrated). Each of these terms can be evaluated through multidimensional integrals, which can be approximated through Monte-Carlo-based algorithms.

On this basis, the Sobol' indices (ranging between null and unity) can be defined as:

= | ( ) S Var E Y X Var Y [ ( )] (3)
i i j

The first-order indices S i are referred to as "the main effects of X i " and can be interpreted as the expected proportion of the total variance of the output Y (i.e. representing the uncertainty in Y) that would be removed if we were able to learn the true value of X i . This index provides a measure of importance (i.e. a sensitivity measure) useful when ranking, in terms of importance, the different input parameters. As a general principle, we attempted to reduce the subjectivity in the representation of our knowledge by relying on the maximum entropy principle for probability law selection [START_REF] Mishra | Assigning Probability Distributions to Input Parameters of Performance Assessment Models[END_REF] . We use the sequential algorithm of Saltelli et al. [START_REF] Saltelli | Variance based sensitivity analysis of model output. Design and estimator for the total sensitivity index[END_REF] , using the R implementation of the Jansen formula [START_REF] Saltelli | Variance based sensitivity analysis of model output. Design and estimator for the total sensitivity index[END_REF][START_REF] Jansen | Analysis of variance designs for model output[END_REF] , and a total of 200,000 model evaluations per time step, which allows to reduce the errors of the first-order indices below 1%.

Surrogate of the PCR model. Under stable sea levels, the PCR models behaves as follows: depending on the frequency and intensity of storms and on the initial conditions, the shoreline, identified as the dune toe, moves around an equilibrium shoreline position, toward which the PCR model converges after a transitional phase of 20 to 30 years. As sea level rises, high water levels during storms become more frequent and induce a retreat of the dune toe parallel to the slopes of the upper shoreface [START_REF] Larson | An analytical model to predict dune erosion due to wave impact[END_REF] , so that sediment blown to the dune by aeolian processes cannot compensate losses and a new equilibrium position is found for the dune toe [START_REF] Ranasinghe | Estimating coastal recession due to sea level rise: beyond the Bruun rule[END_REF] . The Figure provided in Supplementary Material 5 displays this equilibrium response as a function of sea-level change rates. It shows that its value is close to the amount of sea-level rise divided by the slope of the upper shoreface (See section 4 in Supplementary Material 1). The same equilibrium response is obtained for sea-level scenarios following step functions or parabolic curves, as in the sea-level projections used here [START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] . This response, which is here illustrated in the case of site #2, can be explained by the basic principles of the dune erosion model, which assumes that the dune toe evolves parallel to the upper shoreface [START_REF] Larson | An analytical model to predict dune erosion due to wave impact[END_REF] , while the slope of the upper shoreface itself remains within the same range over time (see Supplementary Material 2). While other dune evolution models may deliver different responses, the more complex model SBEACH has provided similar responses so far [START_REF] Ranasinghe | Estimating coastal recession due to sea level rise: beyond the Bruun rule[END_REF] .

Figure 1 .

 1 Figure 1. Sea level reconstructions and projections used in this study (source of projections: Kopp et al.[START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] ; Reconstruction: see Methods, subsection 1).

  Figure 4 partitions the variance of shoreline change projections by displaying the main effects (i.e. 1 st order Sobol' indices) of the uncertain parameters in equation (

Figure 2 .

 2 Figure 2. Reconstructions and projections of shoreline positions using the Bruun rule exemplified at site #1 (See Supplementary materials 1, 3 and 4), provided in the form of median and percentile levels at different timeframes from 1807 to 2200. The reference median shoreline position is arbitrarily set to 0 by 2000, with negative values corresponding to shoreline accretion (seaward) and positive to shoreline retreat (landward).These projections include uncertain shoreface slopes, vertical ground motions, sea-level changes, shoreline and change variability from event-scale to inter-decadal timescales as well as an uncertain multi-decadal trend (see Methods). Note that sea-level projections used here[START_REF] Kopp | Probabilistic 21st and 22nd century sea-level projections at a global network of tide-gauge sites[END_REF] consider that a sea-level drop is possible (although very unlikely) beyond 2100 (Fig.1). Hence, shoreline changes below the median can bend downwards beyond 2100.

Figure 3 .

 3 Figure 3. Same as Fig. 2, with shoreline change projections and reconstructions using the PCR coastal impact model.

Figure 4 .

 4 Figure 4. Variance-based global sensitivity analysis of the shoreline change model response as a function of time, for the four selected sites in Aquitaine (see Supplementary Material 1). For each date considered, the curves indicate the fraction of the variance of shoreline change projections that could be removed if input parameters were known (see main text). The effect of interactions between parameters is indicated as well.White areas indicate interactions between parameters, corresponding to shoreline positions, which can be only reached if at least two uncertain parameters deviate from their mean (see Methods). The graph reads as follows: for site #1, by 2200, uncertainties in regional sea-level rise projections (yellow) account for approximately 30% of the variance of shoreline change projections.

Scientific REPORTS | (2019) 9:42 | DOI:10.1038/s41598-018-37017-4

Scientific REPORTS |

Acknowledgements

This study was supported by the Climate Change projects of the "Aquitaine Coastal Observatory" and by the Risk Prevention Department of BRGM in Orléans, and additional support from the ECLISEA project (European advances on CLImate services for coasts and SEAs) funded through the ERA4CS framework (European Research Area for Climate Services). BC funded by project SONO (ANR-17-CE01-0014) from the Agence Nationale de la Recherche (ANR). RR is supported by the AXA Research Fund and the Deltares Strategic Research Programme "Coastal and Offshore Engineering". The SONEL data assembly center is acknowledged for providing comprehensive access to GNSS data and metadata. We thank Carlos Oliveros, Cyril Mallet, Ali Dastgheib, Aurélie Maspataud, Franck Desmazes and Edward Anthony for useful advices and discussions. We thank Robert Kopp for making his data available. We thank three reviewers for their comments, which allowed to improve this paper.

Data Availability Statement

Results and codes will be provided in an open archive.

Author Contributions

Additional Information

Supplementary information accompanies this paper at https://doi.org/10.1038/s41598-018-37017-4.

Competing Interests:

The authors declare no competing interests.

Publisher's note: Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


B Tuffin

How to obtain a "valid" estimation?

S-valued regenerative process X = (X (t) : t ≥ 0)

Goal: Compute α = E[T ] ,
where

T = inf{t ≥ 0 : X (t) ∈ A}
is the hitting time of subset A Regeneration times 0 = Γ(0) < Γ(1) < • • • , with iid cycles ((τ (k), (X (Γ(k -1) + s) : 0 ≤ s < τ (k)) : k ≥ 1) τ (k) = Γ(k) -Γ(k -1), length of the kth regenerative cycle Ratio expression:

α = E[T ∧ τ ] P(T < τ ) . α = E[T ; T < τ ] + E[τ + T -τ ; T > τ ] = E[T ; T < τ ] + E[τ ; T > τ ] + E[T -τ ; T > τ ] = E[T ∧ τ ; T < τ ] + E[T ∧ τ ; T > τ ] + E[T -τ | T > τ ] P(T > τ ) = E[T ∧ τ ] + α(1 -P(T < τ ))
B. Tuffin (Inria) Hitting times MCM 2019 6 / 44

Regenerative simulation W (k) = inf{t ≥ 0 : X (Γ(k -1) + t) ∈ A} first hitting to A after regeneration Γ(k -1)

I (k) = I(W (k) < τ (k)) with I(•) is the indicator function Definition (Ratio estimator) α(n) = (1/n) n k=1 [W (k) ∧ τ (k)] (1/n) n k=1 I (k)
.

Proposition (Central Limit Theorem)

n 1/2 [α(n) -α] ⇒ σ 2 N (0, 1)
as n → ∞, where

σ 2 2 = E[(T ∧τ ) 2 ] p 2 -2α E[T I(T <τ )] p 2 + α 2 p .
Rare events: hitting A rarely occurs before τ

Denominator p in α = E[T ∧τ ] P(T <τ ) a small probability =⇒ requires an acceleration technique Fraction β of cycles used to estimate the numerator with crude MC Fraction 1 -β to estimate the denominator with a variance reduction technique

Inefficiency of crude Monte Carlo for the denominator

Compute the denominator/probability

p = E[1 [T <τ ] ] << 1 n iiid Y i Bernoulli r.v.:
1 if the event is hit and 0 otherwise. To get a single occurence, we need in average 1/p replications (10 9 for p = 10 -9 ), and more to get a confidence interval. In most cases, you will get (0, 0) as a confidence interval. n Ȳn Binomial with parameters (n, p) and the confidence interval is Ȳn -

c β p(1 -p) √ n , Ȳn + c β p(1 -p) √ n . Relative half width c β σ/( √ np) = c β (1 -p)/p/n → ∞ as p → 0. For a given relative error RE , the required value of n = (c β ) 2 1 -p RE 2 p , inversely proportional to p.
Two main families of techniques:

Splitting (also called subset simulation) and Importance Sampling.

B. Tuffin (Inria) Hitting times MCM 2019

Robustness properties

In rare-event simulation models, we often parameterize with a rarity parameter > 0 such that µ = E[Y ( )] → 0 as → 0.

An estimator Y ( ) is said to have bounded relative variance (or bounded relative error) if σ 2 (Y ( ))/µ 2 ( ) is bounded uniformly in .

Interpretation: estimating µ( ) with a given relative accuracy can be achieved with a bounded number of replications even if → 0.

Weaker property: asymptotic optimality (or logarithmic efficiency) if lim →0 ln(E[Y 2 ( )])/ ln(µ( )) = 2.

Stronger property: vanishing relative variance: σ 2 (Y ( ))/µ 2 ( ) → 0 as → 0. Asymptotically, we get the zero-variance estimator.

Other robustness measures exist (based on higher degree moments, on the Normal approximation, on simulation time...). Importance Sampling (IS)

Let Y = h(X ) for some function h where Y obeys some probability law P. IS replaces P by another probability measure P, using

E [Y ] = h(x)dP(x) = h(x) dP(x) d P(x) d P(x) = Ẽ [h(x)L(x)]
L = dP/d P likelihood ratio, Ẽ is the expectation associated to probability law P.

Required condition: d P(x) = 0 when h(x)dP(x) = 0.

Unbiased estimator:

1 n n i=1 h(X i )L(X i ) with (X i , 1 ≤ i ≤ n) i.i.d;
copies of X , according to P. Goal: select probability law P such that

σ2 [h(X )L(X )] = Ẽ[(h(X )L(X )) 2 ] -µ 2 < σ 2 [h(X )].
B. Tuffin (Inria) Hitting times MCM 2019

Example

We want to estimate the probability that a random variable exceeds T (area in grey under the density f (t)).

0 1 2 3 4 5 6 T f (t) t density f (t)
Reminder: the probability to be in an interval [a, b] is the measure of the area

Importance sampling

Sample according to another density f increasing the probability to be > T .

0 1 2 3 4 5 6 T f (t) f (t) t density f (t) B. Tuffin (Inria)
Hitting times MCM 2019

Importance sampling

Sample according to another density f increasing the probability to be > T .

Rare set reached!

0 1 2 3 4 5 6 T f (t) f (t) t density f (t) B.
(x) = P[X 0 = x], initial probabilities up to a stopping time T µ(x) = E x [Y ].
IS replaces the probabilities of paths (x 0 , . . . , x n ),

P[(X 0 , . . . , X T ) = (x 0 , . . . , x n )] = π 0 (x 0 ) n-1 j=1 P(x j-1 , x j ), by P[(X 0 , . . . , X T ) = (x 0 , . . . , x n )] st Ẽ[T ] < ∞.
For convenience, the IS measure remains a DTMC, replacing P(x, z) by P(x, z) and π 0 (x) by π0 (x).

Then L(X 0 , . . . ,

X T ) = π 0 (X 0 ) π0 (X 0 ) T -1 j=1 P(X j-1 , X j ) P(X j-1 , X j ) .

Zero-variance IS estimator for Markov chains simulation

Restrict to an additive (positive) cost

Y = T j=1 c(X j-1 , X j ) For hitting proba: c(x, z) = 1 if z ∈ A, 0 otherwise, µ(x) ≡ p(x)
For hitting time: c(x, z) avg time in x.

Is there a Markov chain change of measure yielding zero-variance?

We have zero variance with P(x, z) = P(x, z)(c(x, z) + µ(z))

w P(x, w )(c(x, w ) + µ(w )) = P(x, z)(c(x, z) + µ(z)) µ(x) .
Implementing it requires knowing µ(x) ∀x ∈ S, the quantities we wish to compute.

Outline

A short tutorial on rare-event simulation for reg. systems System with c types of components. X = (X 1 , . . . , X c ) with X i number of up components.

IS
1: state with all components up.

Failure rates are O(ε), but not repair rates. Failure propagations possible.

System down (in A) when some combinations of components are down.

Goal: compute µ(1) ≡ p(1) with p(y ) probability to hit A before 1 starting from y (denominator of the ratio est. of MTTF)

Simulation using the embedded DTMC. Failure probabilities are O(ε) (except from 1). How to improve (accelerate) this?

Existing method: ∀y = 1, increase the probability of the set of failures to constant 0.5 < q < 0.9 and use individual probabilities proportional to the original ones (SFB), or uniformly (BFB).

Failures not rare anymore. BRE property verified for BFB. Recall the zero-variance approximation: P(x, z) = P(x, z)(c(x, z) + p(z))

w P(y , w )(c(x, w ) + p(w )) The idea is to approach p(y ) by the probability p(y ) of the path from y to A with the largest probability Intuition: as → 0, we get a good idea of the probability.

Proposition

Bounded Relative Error proved (as → 0) in general. Even Vanishing Relative Error if p(y ) contains all the paths with the smallest degree in .

Other simple version: approach p(y ) by the (sum of) probability of paths from y with only failure components of a given type. Gain of several orders of magnitudes + stability of the results with respect to the literature.

HRMS: numerical illustrations

Comparison of BFB and Zero-Variance Approximation (ZVA). c = 3 types of components, n i of type i failure rates ε, 1.5ε, and 2ε 2 , repair rate 1 System is down whenever fewer than two components of any one type are operational.

n i ε µ 0 BFB est ZVA est BFB σ 2 ZVA σ 2 0.001 2.6 × 10 -3
2.7 × 10 -3 2.6 × 10 -3 6.2 × 10 -5 2.2 × 10 -8 0.01 1.8 × 10 -7 1.9 × 10 -7 1.8 × 10 -7 6.3 × 10 -11 2.0 × 10 -14 0.001 1.7 × 10 -11

1.8 × 10 -11 1.7 × 10 -11 8.8 × 10 -19 1.2 × 10 -23 0.1 6.0 × 10 -8 4.8 × 10 -8 6.0 × 10 -8 8.1 × 10 -10 1.6 × 10 -10 0.001 3.9 × 10 -28

(1.8 × 10 -40 ) 3.9 × 10 -28 (3.2 × 10 -74 ) 1.4 × 10 -55

Crude estimators of MTTF

Notations for an S-valued regenerative process X = (X (t) : t ≥ 0)

Compute α = E[T ],
where T = inf{t ≥ 0 : X (t) ∈ A} is the hitting time of subset A

Crude estimators of MTTF

Notations for an S-valued regenerative process X = (X (t) : t ≥ 0)

Compute α = E[T ], where T = inf{t ≥ 0 : X (t) ∈ A} is the hitting time of subset A Regeneration times 0 = Γ(0) < Γ(1) < • • • , with iid cycles ((τ (k), (X (Γ(k -1) + s) : 0 ≤ s < τ (k)) : k ≥ 1) τ (k) = Γ(k) -Γ(k -1), length of the kth regenerative cycle W (k) = inf{t ≥ 0 : X (Γ(k -1) + t) ∈ A} first hitting to A after regeneration Γ(k -1) I (k) = I(W (k) < τ (k)) with I(•) is the indicator function

Crude estimators of MTTF

Notations for an S-valued regenerative process X = (X (t) : t ≥ 0)

Compute α = E[T ], where T = inf{t ≥ 0 : X (t) ∈ A} is the hitting time of subset A Regeneration times 0 = Γ(0) < Γ(1) < • • • , with iid cycles ((τ (k), (X (Γ(k -1) + s) : 0 ≤ s < τ (k)) : k ≥ 1) τ (k) = Γ(k) -Γ(k -1), length of the kth regenerative cycle W (k) = inf{t ≥ 0 : X (Γ(k -1) + t) ∈ A} first hitting to A after regeneration Γ(k -1) I (k) = I(W (k) < τ (k)) with I(•) is the indicator function Ratio expression: α = E[T ∧ τ ] P(T < τ ) . Definition Direct estimator α 1 (m) = 1 m m j=1 T (j). Ratio estimator α 2 (n) = (1/n) n k=1 [W (k) ∧ τ (k)] (1/n) n k=1 I (k)
.

(Known) Central limit theorems If p = P(T < τ ) > 0:

Proposition (Direct estimator)

m 1/2 [α 1 (m) -α] ⇒ σ 1 N (0, 1)
as m → ∞, where

σ 2 1 = α 2 + E[(T ∧ τ ) 2 ] p -2α E[T I(T < τ )] p .
Proposition (Ratio-based estimator)

n 1/2 [α 2 (n) -α] ⇒ σ 2 N (0, 1)
as n → ∞, where

σ 2 2 = E[(T ∧ τ ) 2 ] p 2 -2α E[T I(T < τ )] p 2 + α 2 p .
Asymptotic result with a decreasing sequence of reachable sets Theorem (Both estimators asymptotically identical)

Assume E[τ 3 ] < ∞. If t b p b → ∞ as b → ∞, then we have that as b → ∞, √ t b p b αi,b (t b ) E[T b ] -1 ⇒ E[τ ]N (0, 1), i = 1, 2, and √ t b p b α1,b (t b ) E[T b ] - α2,b (t b ) E[T b ] ⇒ 0.
B. Tuffin (Inria) Hitting times MCM 2019

Efficient Regenerative IS estimators extensively studied.

Question: What about the direct estimator? Can its combination with IS yield an efficient estimator?

We will play with the toy example:

0 1 2 2 1 with embedded DTMC 0 1 2 1 /(1 + ) 1/(1 + ) E (T ) = ∞ n=0 (n + 1) 1 2 + 1 1 + 1 1 + n 1 + = 1 + 3 2 2 E [(T ) 2 ] = ∞ n=0 (n + 1) 2 1 2 + 1 1 + 2 1 1 + n 1 + = (2 + )(1 + 3 ) 2 4(1 + ) 4 E (N) = ∞ n=0 (2 + 2n) 1 1 + n 1 + = 2(1 + )
with N:# transitions in a run.

Failure biasing

Change the probability of making a failure transition to be ρ, independent of

0 1 2 1 ρ 1 -ρ Ẽ [(T L) 2 ] = E [(T ) 2 L] = ∞ n=0 (n + 1) 2 1 2 + 1 1 + 2 1 1+ n 1+ 2 (1 -ρ) n ρ Converging sum iff 1/((1 + ) 2 (1 -ρ)) < 1, i.e., ρ small enough ρ < 1 - 1 (1 + ) 2 = 2 -3 2 + o( 2 ).

Failure biasing

Change the probability of making a failure transition to be ρ, independent of

0 1 2 1 ρ 1 -ρ Ẽ [(T L) 2 ] = E [(T ) 2 L] = ∞ n=0 (n + 1) 2 1 2 + 1 1 + 2 1 1+ n 1+ 2 (1 -ρ) n ρ Converging sum iff 1/((1 + ) 2 (1 -ρ)) < 1, i.e., ρ small enough ρ < 1 - 1 (1 + ) 2 = 2 -3 2 + o( 2 ). But Ẽ (N) = ∞ n=0 (2 + 2n)(1 -ρ) n ρ = 2 ρ .
The average simulation time for a single run will increase to infinity as → 0!

Zero-variance approximation

For a CTMC with transition matrix (P x,y ) x,y ∈S , if E ,x expectation starting from x,

Px,y = P x,y 1/λ(x) + E ,y (T )

E ,x (T )
yields an estimator with variance zero. On our toy example, the only probability we can change is from 1

0 1 2 1 ρ 1 -ρ ρ = 1 + 1 1+ + 0 1+2 2 2 = 2 3 (1 + ) 2 (1 + 2 )
yields variance 0.

B. Tuffin (Inria) Hitting times MCM 2019 34 / 44

Zero-variance approximation

For a CTMC with transition matrix (P x,y ) x,y ∈S , if E ,x expectation starting from x,

Px,y = P x,y 1/λ(x) + E ,y (T )

E ,x (T )
yields an estimator with variance zero. On our toy example, the only probability we can change is from 1

0 1 2 1 ρ 1 -ρ ρ = 1 + 1 1+ + 0 1+2 2 2 = 2 3 (1 + ) 2 (1 + 2 )
yields variance 0.

But the estimation takes on average longer time, 2 ρ = Θ( -3 ), as gets closer to zero. An approximation of the zero-variance IS can be inefficient, producing an unbounded work-normalized relative variance.

Figure

  Figure: Original probabilities Figure: Probabilities under IS/BFB

  Let α1,b (c) and α2,b (c) be the estimators obtained after c units of computational time To hope for consistency and CLTs, we need a computational budget t b for which t b p b → ∞ as b → ∞

  Instead of counting 1 each time we are > T and look at the average value for each sample value t i , we count 1(t i > T ) f (t i ) f (t i ) (ratio of heights under densities at t i ) and look again at the average value ⇒ unbiased estimation: the true probability is estimated.

			f (t i )
	density f (t)
			f (t i )
	3	3.2	3.4	3.6	3.8	4
			t
	Tuffin (Inria) B. Tuffin (Inria)		Hitting times Hitting times			MCM 2019 MCM 2019	14 / 44 16 / 44

IS for a discrete-time Markov chain (DTMC) {X j , j ≥ 0} Y = h(X 0 , . . . , X T ) function of the sample path with P = (P(x, z)) x,z∈S transition matrix, π 0

  application: simulation of highly reliable Markovian systems

	Highly Reliable Markovian Systems (HRMS)
	Mean Time To Failure (MTTF) estimation by simulation: direct or
	regenerative estimator?
	Crude estimations
	Comparison of crude estimators
	Importance Sampling estimators
	Quantiles and tail-distribution measures
	Definitions
	Exponential approximation and associated estimators
	Numerical examples
	B. Tuffin (Inria)	Hitting times	MCM 2019	20 / 44

i.e., the proportion of points is the probability under the new density does not correspond to the grey area, but to the blue one.

Crude direct and ratio-based estimators are asymptotically equivalent (in two asymptotic contexts)

When IS is used, the regenerative expression is rather advised.B. Tuffin (Inria)Hitting times MCM 2019

Zero-variance approximation

Use a heuristic approximation μ(•) and plug it into the zero-variance change of measure instead of µ(•) P(y , z) = P(y , z)(c(y , z) + μ(z))

w P(y , w )(c(y , w ) + μ(w ))

Main ideas

Glynn, Nakayama & T., WSC, 2017 Two potential estimators:

Direct estimator: repeat experiments up to failure of the system, and compute the average value Literature, regenerative estimator: expresses the MTTF as a ratio of quantities over regenerative cycles Two potential estimators:

Direct estimator: repeat experiments up to failure of the system, and compute the average value Literature, regenerative estimator: expresses the MTTF as a ratio of quantities over regenerative cycles Question: Is there a reason why the regenerative estimator is used? Which one is "better"?

Two potential estimators:

Direct estimator: repeat experiments up to failure of the system, and compute the average value Literature, regenerative estimator: expresses the MTTF as a ratio of quantities over regenerative cycles Question: Is there a reason why the regenerative estimator is used? Which one is "better"?

Contributions

Crude (direct and regenerative) estimators are asymptotically similar in performance, in rare event settings

For Importance Sampling estimators, the regenerative one yield a efficient estimator when the crude can not.

Estimators α 1 (m) and α 2 (n) are actually very similar

If N(j) = inf{k > N(j -1) :

Is an estimator more efficient than the other?

B. Tuffin (Inria) Hitting times MCM 2019 29 / 44

Question: which estimator is "more efficient"?

Estimators α 1 (m) and α 2 (n) are actually very similar

If N(j) = inf{k > N(j -1) :

Is an estimator more efficient than the other? Two asymptotic settings

Failure biasing

Change the probability of making a failure transition to be ρ, independent of

Hitting times MCM 2019

Zero-variance approximation

For a CTMC with transition matrix (P x,y ) x,y ∈S , if E ,x expectation starting from x,

Px,y = P x,y 1/λ(x) + E ,y (T )

yields an estimator with variance zero.

Discussion on the impact of the approximation

, we retrieve a variance zero.

Discussion on the impact of the approximation

, we retrieve a variance zero. For ρ = 3 (approximation of good asymptotic order), the variance is Θ( -2 ), but the work-normalized relative variance is unbounded due to the computational time.

Discussion on the impact of the approximation For ρ = 2 3 (1+ ) 2 (1+2 ) , we retrieve a variance zero. For ρ = 3 (approximation of good asymptotic order), the variance is Θ( -2 ), but the work-normalized relative variance is unbounded due to the computational time.

For ρ = 2 3 (exact first-order term), the variance is Θ(1), which is better but still not sufficient to yield a bounded work-normalized variance.

Discussion on the impact of the approximation For ρ = 2 3 (1+ ) 2 (1+2 ) , we retrieve a variance zero. For ρ = 3 (approximation of good asymptotic order), the variance is Θ( -2 ), but the work-normalized relative variance is unbounded due to the computational time.

For ρ = 2 3 (exact first-order term), the variance is Θ(1), which is better but still not sufficient to yield a bounded work-normalized variance.

Much better than an exact first-order approximation is required. Hard to obtain in practice.

B. Tuffin (Inria)

Hitting times MCM 2019

Conclusions on MTTF estimation

We have compared two standard estimators of the MTTF for regenerative processes a direct one expressed as the average of simulated times to failure one making use of the regenerative structure

Outline

A short tutorial on rare-event simulation for reg. systems

Basic idea

Let F be the cumulative distribution function of T Goal: For fixed 0 < q < 1, estimate the q-quantile (0 < q < 1)

and the conditional tail expectation (CTE)

.

Hitting times MCM 2019

Decomposition

Using τ i = Γ i -Γ i-1 and M the number of first cycles not reaching A

M geometric r.v. with P(M = k) = p(1 -p) k where p = P(T < τ ).

Recall that the regenerative structure of X allows to express

Asymptotic regimes/exponential approximation

Introduction of a rarity parameter Assumption: p ≡ p → 0 as → 0.

Ex HRMS: Probability of reaching a failed state before coming back to the initial (perfectly working) state goes to 0 with failure rates Ex GI/G/1 queue: considering a receding set of states (number of customers)

Theorem (Known result)

The scaled hitting time T /α converges weakly to an exponential: for each x ≥ 0, P (T /α ≤ x) → 1 -e -x as → 0.

B. Tuffin (Inria) Hitting times MCM 2019 40 / 44

Quantile and CTE estimators based on the exponential approximation

From

we get

Using the ZVA efficient estimator α of α, we get


Mathieu Bourdeau

Serge Iyock

David Werner

Elyes Nefzaoui

Sensor networks for existing residential buildings performance evaluation

Enhancing residential building energy efficiency has become a critical goal to take up current challenges of human comfort, urbanization growth and the consequent energy consumption increase. In a context of integrated smart infrastructures, sensor networks offer a relevant solution to support building energy monitoring, operation and prediction. The amount of accessible data also opens new prospects to better consider key parameters such as human behavior and to lead to more efficient energy retrofit of existing buildings. However, sensor networks planning and deployment is a particularly complex task facing many challenges and affecting the performances of such a promising solution. The present paper highlights watchpoints from the implementation of an instrumentation solution for the study and evaluation of a deep energy retrofit of existing collective residential buildings. Observations will be grouped in four categories. They will be illustrated with actual field situations and discussed to provide potential solutions for efficient future sensor networks deployment.

Introduction

Residential buildings are one of the major energy consumers and greenhouse gas (GHG) emitters, with 38.1% of the final energy consumption and 36% of the GHG emissions in Europe [START_REF] Ademe | Climat, air et énergie |Chiffres-clés[END_REF]. Efforts have been made to reduce the impact of residential buildings on energy consumption with various codes, standards and thermal regulation mandatory compliance for designs of new buildings [START_REF] Harmonie | 1974-2020 : l'évolution de la réglementation thermique[END_REF], [START_REF]Standard 90.1-2013 -Energy Standard for Buildings Except Low-Rise Residential Buildings (SI Edition)[END_REF]. However, given the slow turnover of the building stock in European countries, the effect of these regulations is limited [START_REF] Insee | Le parc de logements en France au 1 er janvier 2017[END_REF]. Hence, existing building retrofit turns into a priority [START_REF]Le Plan de rénovation énergétique de l'habitat (PREH)[END_REF]. On the other hand, ambitious massive smart sensor networks plans [START_REF]Décret n° 2010-1022 du 31 août 2010 relatif aux dispositifs de comptage sur les réseaux publics d'électricité en application du IV de l'article 4 de la loi n° 2000-108 du 10 février 2000 relative à la modernisation et au développement du service public de[END_REF]- [START_REF]Smart Metering deployment in the European Union | JRC Smart Electricity Systems and Interoperability[END_REF] have been launched in recent years to enhance buildings energy efficiency. Smart sensor networks have been largely deployed in buildings for energy monitoring and operation [START_REF] Fan | Unsupervised data analytics in mining big building operational data for energy efficiency enhancement: A review[END_REF], [START_REF] Bourdeau | Buildings energy consumption generation gap: A post-occupancy assessment in a case study of three higher education buildings[END_REF] or building energy consumption forecasting [START_REF] Bourdeau | Modeling and forecasting building energy consumption: A review of data-driven techniques[END_REF], [START_REF] Jain | Forecasting energy consumption of multi-family residential buildings using support vector regression: Investigating the impact of temporal and spatial monitoring granularity on performance accuracy[END_REF]. Studies have used sensors to highlight and characterize the link between energy efficiency and inhabitants' behavior [START_REF] Pisello | Human-based energy retrofits in residential buildings: A cost-effective alternative to traditional physical strategies[END_REF]- [START_REF] Nefzaoui | Coévolutions des systèmes techniques et des comportements : le cas de la rénovation énergétique[END_REF], identified as one of the main source of energy performance gaps [START_REF] De Wilde | The gap between predicted and measured energy performance of buildings: A framework for investigation[END_REF]. Data collection solutions have also been proposed to supervise building energy retrofits [START_REF] Cstb | ISABELE mesure la performance énergétique intrinsèque d'un[END_REF]- [START_REF] Ben | Energy retrofit and occupant behaviour in protected housing: A case study of the Brunswick Centre in London[END_REF]. However, instrumentation solutions deployment is a complex task, facing many challenges and difficulties. Thus, it directly impacts on the quality of the analyses and the efficiency of related energy savings measures [START_REF] Pisello | Human-based energy retrofits in residential buildings: A cost-effective alternative to traditional physical strategies[END_REF], [START_REF] Jankovic | Lessons learnt from design, off-site construction and performance analysis of deep energy retrofit of residential buildings[END_REF], [START_REF] Calì | Energy performance gap in refurbished German dwellings: Lesson learned from a field test[END_REF]. In this context, we present a case study of an instrumentation solution deployed in existing residential buildings. The instrumentation takes part in a larger project that aims to study the impact of deep energy retrofit measures on building energy consumption. The present paper is organized as follows: Section 2 presents the project case study and Section 3 describes details of the instrumentation solution. Then a feedback on encountered watch points for sensor networks implementation is proposed in Section 4. Watch-points will be illustrated in the presentation with examples of field situations. They will be discussed to help provide feedback on sensor network deployment.

The case study

In the present project, a group of three existing social residential buildings -63 apartments over a 4,600 m² living area -is considered for field experimentation. Buildings were built in 1974 and are located in Paris (France) eastern suburb. They only had minor retrofit actions and are now planned to undergo deep energy retrofit measures over an eighteen-month period in 2019-2020. Retrofit measures will be implemented in occupied site as apartments will be occupied by tenants during the whole retrofit period.

Instrumentation methodology 3.1.

Instrumentation plan Instrumentation is performed in two steps and includes more than 280 sensors. The first step covers communal areas and entire building level and involves 10% of the total number of sensors. The second step covers 10 households. Step 1 sensors have been installed and data collection is operational. Step 1 serves as a test to gain knowledge on sensors operation, and to identify the key issues and watch points. Step 2 is currently being implemented. The sensor network aims to provide a four-level instrumentation scale (residence, building, apartment and room level) with a particularly fine time-step granularity (from one-minute to onehour time-step). Types of collected data can be divided into five categories: 1) Outdoor environment quality (OEQ) is provided by a weather station measuring outdoor temperature, humidity, rain level, solar irradiation, wind speed and direction 2) Indoor Environment Quality (IEQ) provided by sensors in communal areas as well as apartments which measure indoor temperature, humidity, luminosity, CO2 level and indoor surface temperature of exterior walls at one-hour time-step 3) Inhabitants' behavior characterized with presence and windows opening detection, 4) energy consumption which is monitored through electricity, natural gas, heating and domestic hot water consumption metering at one-minute time-step and 5) the building thermal insulation through walls surface temperatures and heat fluxes. All sensors communicate using LoRa technology, on operated or private communication networks, exception made of two types of electricity sensors directly communicating data with GPRS network. Therefore, gateways have been installed to enhance the network communication capabilities. All data are stored on an FTP server with one data file per sensor per day.

Description of the instrumentation management process

The overall instrumentation of the present project can be divided into six different stages, spread out over an estimated timespan of one year (started in October 2018 and expected to end in September 2018). All the six stages have a crucial impact on the success of the instrumentation and require different knowledge and skills. First is the accurate definition of the instrumentation specifications (1) to meet the needs of the project. In parallel, a review of service providers offering adapted solutions is done (2) and apartments are recruited to participate to the instrumentation study [START_REF]Standard 90.1-2013 -Energy Standard for Buildings Except Low-Rise Residential Buildings (SI Edition)[END_REF]. Once specifications are organized in a clear framework, they are submitted to identified service providers [START_REF] Insee | Le parc de logements en France au 1 er janvier 2017[END_REF]. Several rounds of discussions allow modifications on specifications to adapt to the reality of the market and for service providers to refine their offer to meet the needs of the project to the best of their capabilities. Every round dismisses several service providers to select the most relevant one. Once the quote is signed, sensors can be ordered, delivered and installed [START_REF]Le Plan de rénovation énergétique de l'habitat (PREH)[END_REF]. Finally, once the sensor network is functional comes the hand-over of the instrumentation solution (6).

Watch-points from sensor network implementation

During the planning and deployment of the first phase of the building instrumentation, several watchpoints and issues have arisen. These are reported in the present paper and are grouped in five categories: (1) onsite installation conditions and environment, (2) characteristics and purpose of instrumentation, (3) service provider and (4) project tracking and management. The experimentation site directly impacts on the challenges and success of building instrumentation. Indeed, a main feature of the present project is to collect data on existing and occupied buildings. Therefore, identification and verification of the facilities to be instrumented is critical. Moreover, characteristics of the instrumentation environment can impact on the sensor communication quality and installation feasibility, inducing additional costs and delays. A particular attention should also be given in the technical characteristics of the instrumentation solution. Specifications on the data to be collected is impacting on many related factors such as the choice of an adapted data collection method (wired or wireless), the communication networks, the requirements on the sensor operation (time-step, adaptability to the environment, lifespan). Also, because of the purpose of the study -a research project on inhabited building energy performance -the solution has to comply to other constraints such as being non-intrusive and ensuring personal data protection. The service provider for the instrumentation solution has also a significant impact on the execution of the sensor network deployment plan. Indeed, it may be difficult to find a service provider ensuring the provision of relevant sensors, their integration in a data collection platform and installation onsite. It is often necessary to contact more than one provider to have the whole range of services required. It becomes even more complex when custom services are needed, such as specific sensor characteristics or sensor network maintenance. Finally, a rigorous project tracking and management process is critical for all six stages described. Many details may induce important delays, additional costs and other issues. The four categories of watchpoints will be illustrated with various examples and situation encountered during the implementation of the first instrumentation phase of the present project. Discussions will also aim to provide a guideline for sensor network deployment for future projects.