Mathilde Jay

Anne Régnier

Anaïs Dasnon

Killian Brunet

Marie Pelé
email: marie.pele@iphc.cnrs.fr

The light is red: uncertainty behaviours displayed by pedestrians during illegal road crossing

Keywords: illegal road-crossing, pedestrian, uncertainty, hesitation, decision-making, culture

Road accidents involving pedestrians are a reality of urban life. Pedestrian risk is now well known and documented from the perspective of drivers. However, pedestrian behaviour plays a central role in road accidents, notably in terms of illegal road crossing at signalized intersections. This study focuses on pedestrians crossing illegally at a signal light, and specifically investigates uncertainty behaviour, also referred to as hesitation, which occurs when a pedestrian slows down or stops his/her crossing movement then (1) abandons the crossing by returning to the kerb or (2) accelerates to cross the road more quickly. We sought to understand the causes of this behaviour in France and Japan, two countries where interesting differences have already been demonstrated in the way pedestrians behave. The results show a longer period of uncertainty for pedestrians in Japan compared to France.

Japanese pedestrians also hesitated longer when they were alone. This study demonstrates a tendency to speed up if there are a number of pedestrians already crossing the road, but abandoning behaviours were more frequently observed than acceleration. This study confirms that pedestrians may misevaluate the moment to cross and hesitate when they realise that they have made a mistake, thus increasing the risk of an accident. These results could help to find solutions that prevent illegal and dangerous road-crossing behaviours.

INTRODUCTION

In 1896, Bridget Driscoll was killed by a motor vehicle during a demonstration in the United Kingdom; her death is the first recorded case of pedestrian road fatality [START_REF] Porter | Reducing road traffic. First fatal car crash in Britain occurred in 1898[END_REF].

The number of road traffic accidents has since continued to rise, with pedestrians representing 22% of world road traffic death (World Health Organization, 2015, WHO/NMH/NVI/15.6).

Despite the efforts made by governments to enhance road safety practices, there are still too many accidents involving pedestrians, most of which occur in urban agglomerations. In France, collisions with pedestrians represent 14% of lethal road accidents, and 67% of these deaths occur in the urban environment (ONISR, French Observatory of Road Safety, 2017).

These deaths are all the more inacceptable given their largely predictable and preventable nature. Indeed, risks to pedestrians are now well known and documented in terms of driver behaviour (for speeding, see the review by Rosén et al., 2011; for drink driving see [START_REF] Hingson | Epidemiology and consequences of drinking and driving[END_REF] to road infrastructures (see [START_REF] Bunn | Traffic calming for the prevention of road traffic injuries: systematic review and meta-analysis[END_REF][START_REF] Hamed | Analysis of pedestrians' behavior at pedestrian crossings[END_REF] and vehicle design (see [START_REF] Desapriya | Do light truck vehicles (LTV) impose greater risk of pedestrian injury than passenger cars? A meta-analysis and systematic review[END_REF].

Another key factor of pedestrian safety is the behaviour of the pedestrians. [START_REF] Zhuang | Pedestrians' crossing behaviors and safety at unmarked roadway in China[END_REF] showed that 66% of pedestrians crossing outside signalized intersections in China did not check for oncoming vehicles before crossing, and 16.1% did not even do so whilst crossing. Among pedestrians who saw an approaching vehicle after stepping off the kerb, 40.6% stopped crossing, 11.4% stepped backwards and 31.9% accelerated to complete the crossing. In this example, more than 40% of pedestrians demonstrated uncertainty behaviours.

The authors observed many instances where pedestrians were using cell phones or listening to music while walking or even while crossing roads. The risky pedestrian behaviour described in this case is also problematic in the context of illegal road-crossings at signalized intersections (Diaz, 2002;Preusser et al., 2002;Evans and Norman, 1998;Evans and Norman, 2003;[START_REF] Rosenbloom | Crossing at a red light: Behaviour of individuals and groups[END_REF]King et al., 2009;Xu et al., 2013).

Pedestrian behaviours vary according to the socio-demographic traits of individuals.

Men show more risky behaviours than women, who display more compliance and conformity [START_REF] Osman | Conformity or compliance? A study of sex differences in pedestrian behaviour[END_REF][START_REF] Tom | Gender differences in pedestrian rule compliance and visual search at signalized and unsignalized crossroads[END_REF]. These behaviours can also vary according to countries and social conventions [START_REF] Sueur | Different risk thresholds in pedestrian road crossing behaviour: a comparison of French and Japanese approaches[END_REF]. For instance, French pedestrians crossed against the red light much more often than their Japanese counterparts (41.9% vs 2.1%) [START_REF] Pelé | Cultural influence of social information use in pedestrian road-crossing behaviours[END_REF]. As unsafe behaviours, pedestrian infractions can be classified as either errors or violations (Reason, 1990). While an error is unintentional and associated with individual information, violation is deliberate and associated with social factors (Reason, 1990). Social information (gained from surrounding people) is an important factor in collective human phenomena and followership [START_REF] Faria | Collective behavior in road crossing pedestrians: the role of social information[END_REF][START_REF] Boos | Leadership in Moving Human Groups[END_REF]. For instance, Dyer and colleagues (2009) showed that a small percentage of informed individuals (5%) was sufficient to guide a larger group of uninformed people to a target place. At a signalised intersection, some pedestrians do not check the light colour but simply follow when other pedestrians have already started to cross [START_REF] Pelé | Decision-making processes underlying pedestrian behaviours at signalised crossings: Part 2. Do pedestrians show cultural herding behaviour[END_REF]. Trusting wrong or unreliable information and following someone who is crossing without first checking the light colour (personal information) leads to an increased risk of accidents and injuries [START_REF] Rosenbloom | Crossing at a red light: Behaviour of individuals and groups[END_REF]. Indeed, a conflict between these two types of information (social vs personal) may lead the pedestrian to change his/her behaviour while already crossing. Our study focuses on this specific behavioural change in pedestrians crossing illegally at a signal light. We have called this uncertainty behaviour (also referred to as hesitation). We define this specific behaviour as the moment the pedestrian slows down or stops his/her crossing movement then either (1) abandons the crossing by turning back to the kerb or (2) accelerates to cross more quickly.

Both of these reactions are potentially dangerous for the pedestrian, who may be hit by vehicles driving on inner or outer road lanes.

To our knowledge, this is the first study to investigate the factors (sociodemographic and environmental) affecting uncertainty behaviour. This specific behaviour is a clear signal that the pedestrian has misevaluated the situation, and it deserves in-depth analysis. We also sought to identify any factors that could affect the abandon and acceleration rates observed in pedestrians. This was achieved by studying the behaviours of pedestrians crossing streets in two different countries: Japan and France. We have already demonstrated interesting differences in the way pedestrians behave in these two countries [START_REF] Sueur | Different risk thresholds in pedestrian road crossing behaviour: a comparison of French and Japanese approaches[END_REF][START_REF] Pelé | Cultural influence of social information use in pedestrian road-crossing behaviours[END_REF]. First, there was a higher frequency of infractions in France than in Japan. In view of this finding, we predict a greater relative number of uncertainty behaviours in French pedestrians. Women and the elderly are expected to show more uncertainty behaviours and the longest uncertainty times in comparison to other categories of pedestrians when crossing at the red light. Indeed, women appear to focus more on the behaviour of other pedestrians at the time of crossing (i.e. prioritizing social information) [START_REF] Tom | Gender differences in pedestrian rule compliance and visual search at signalized and unsignalized crossroads[END_REF] while older pedestrians (> 65 years old) appear to prioritize their course of travel over other environmental information (Tapiro et al., 2016a). Men and the youngest pedestrians are expected to accelerate after hesitation. At a pedestrian crossing red light, combining the number of pedestrians already crossing with the pedestrian departure order (rank) should increase the rate of illegal crossings and subsequently the rate of uncertainty behaviours, as pedestrians show mimetic process and trust in social information [START_REF] Faria | Collective behavior in road crossing pedestrians: the role of social information[END_REF][START_REF] Bode | Information use by humans during dynamic route choice in virtual crowd evacuations[END_REF]. We also expect any pedestrians who are using their phones or accompanied by another person to show more uncertainty behaviours than other pedestrians, as these factors are known as distractors [START_REF] Neider | Pedestrians, vehicles, and cell phones[END_REF]Tapiro et al., 2016b;[START_REF] Tapiro | The Effect of Environmental Distractions on Child Pedestrian's Crossing Behavior[END_REF].

Moreover, the probability of abandoning instead of accelerating should increase with the number of lanes on the road, as the risk and the perception of risk increase with the crossing distance.

MATERIALS & METHODS

Study sites

Pedestrians were observed at three sites in Strasbourg, France and at four sites in Nagoya, Japan (see details in Table 1). On each site, the speed of vehicles was limited to 50km/h. The speed of pedestrians when crossing the road was estimated by scoring the speed of 20 random pedestrians at each site. It was not significantly different from one site to another (permutation test for independent samples: maxT=2.22, p=0.168). At some sites, vehicles approaching the intersections were allowed to enter the road despite the green light for pedestrians; these drivers were aware that crossing pedestrians had priority and drove much more slowly than vehicles that were heading straight on. However, the driver of an approaching vehicle may be less careful if pedestrians cross at a red light, as the former has the right of way. Pedestrians therefore run a much higher risk of being injured when crossing at a red light. There was no button for pedestrians to trigger the green light at any of the sites.

Data scores

Data were scored over a six-day period for each site, for one hour per day during working days and between 12:00 and 17:00 to ensure that data included different groups of pedestrians (workers, families with children, retirees, etc.) but also excluded movements generated by social events (tourism, festivals, concerts, etc.). Observation time then varied within locations to get a wider pool of pedestrians. This scoring duration is sufficient to provide a large dataset [START_REF] Faria | Collective behavior in road crossing pedestrians: the role of social information[END_REF][START_REF] Sueur | Different risk thresholds in pedestrian road crossing behaviour: a comparison of French and Japanese approaches[END_REF][START_REF] Yannis | Pedestrian gap acceptance for midblock street crossing[END_REF]. A video camera was set up in order to score the light colour in locations that were chosen to ensure that crossing pedestrians were visible at all times. We used the behavioural sampling method [START_REF] Altmann | Observational Study of Behavior: Sampling Methods[END_REF] to score the behaviours of pedestrians in one direction only, i.e. those recorded by the camera. Pedestrians were not informed about the purpose of the study to avoid any modification of their behaviour. As both cities are touristic, pedestrians are accustomed to seeing tourists taking pictures or videos. Our camera was placed on the other side of the road to the starting point of pedestrians. It was not placed directly in front of the crosswalk but was located slightly to one side to make it less of a distraction for pedestrians before and during crossing. None of the pedestrians who noticed the camera changed their direction. To avoid influencing pedestrian behaviour, we purposely did not take any other equipment such as counters or pocket PCs. Any observations of road-crossing that were hampered by a visual obstacle (i.e. a car or a truck in front of the video camera) and the observations occurring immediately before and after it were removed from the analyses. We also removed data in which cyclists and tourists were among the pedestrians. We differentiated tourists from citizens without difficulty, as they were mostly in large groups accompanied by a guide, or carried specific equipment (camera or travel guide) and could be dressed differently to local citizens.

Research ethics

This study only involves anonymous observations and anonymous data scores. We followed the ethical guidelines for research of both countries and received ethical approval from both the Institut Pluridisciplinaire Hubert Curien, Strasbourg, France and the Primate Research Institute, Kyoto University, Japan to conduct our observations. As all data were anonymous, focal individuals were given sequential numerical identities according to the day of observation, the crossing concerned and the time of the road crossing. To avoid influencing pedestrians, the experimenters did not approach them. However, pedestrians had the opportunity to learn about the study by an informative medium in their language (French or Japanese) if requested after seeing the camera. If they wished to do so, they also had the possibility to obtain an email address and phone number to contact the corresponding author at a later date. Neither of the two Japanese pedestrians nor the three French ones who asked for information refused to participate in the study.

Data analyses

Uncertainty behaviours were analysed from video recordings by two teams of two coders (AR with AD and MJ with KB, each team coding for French and Japanese videos), using BORIS software [START_REF] Friard | BORIS: a free, versatile open-source event-logging software for video/audio coding and live observations[END_REF] and demonstrating all coder reliability of 80%. We only scored and analysed uncertainty behaviours observed at the red light on pedestrian crossings. However, in order to have a rate of uncertainty behaviour per illegal road-crossings, we also scored the number of illegal road-crossings behaviours in France and in Japan, regardless of the occurrence of uncertainty behaviour.

We defined the start of a crossing as the moment the pedestrian lifted one of his/her feet to cross. The uncertainty behaviour was defined as the moment the pedestrian stopped his/her crossing movement, or slowed down. This hesitation was then followed by one of two possible behaviours, i.e. scenario 1, the abandon of the road crossing, which is characterized by returning to the kerb, or scenario 2, acceleration, characterized by an acceleration of the crossing movement (Figure 1). These two possible behaviours following hesitation will hereafter be referred to as "subsequent behaviours". We never observed a pedestrian who continued to cross at the same speed after hesitating, with speed estimated according to the number of steps per second.

Using these behavioural sequences, we defined uncertainty time as the time comprised between the beginning of the uncertainty behaviour and the beginning of the following behaviour (abandon or acceleration).

To determine the environmental and sociodemographic factors affecting uncertainty behaviour and its possible consecutive behaviours, we scored the following variables for each focal pedestrian:

-gender (male or female), -age, estimated at 10-year intervals from 0-9, 10-19 [ . . . ] to 70-89.

-country (France or Japan), -the number of lanes to be crossed, -whether or not the focal pedestrian was using their phone, -whether or not the pedestrian was accompanied (i.e. the pedestrian was talking with another person or walking very closely and synchronously with her/him, arriving at the road and crossing it at the same time), -the total number of pedestrians involved in the road-crossing event, -the departure order of pedestrians, where the first pedestrian to leave is ranked as 1, the second as 2, and so on. No uncertainty behaviour was noted after the sixth rank.

Statistical analyses

We first used chi-square tests to assess differences between countries and between genders. To evaluate the effect of age on uncertainty, we calculated the percentage of uncertainty behaviours (i.e. pedestrians as one behaviour per pedestrian) in each age category.

Two samples were used for this purpose: (1) a random sample of 1,278 pedestrians (N=551 for France and N=727 for Japan, with 2x10 minutes scorings for each site on two different videos) and (2) our sample of 118 pedestrians crossing at the red light and showing uncertainty in each country (N=80 for France and N=38 for Japan, see details for each site in Table 1). For each sample, we divided the number of individuals per age category by the total number of uncertainty behaviours (*100). We then calculated the differences between the two percentages in each age category: negative values indicated that fewer individuals showed uncertainty behaviours than the average for a given age category, while positive values indicated that more individuals showed uncertainty behaviours than the average. A Kolmogorov-Smirnov test was realised to compare the distributions of the hesitation frequency according to departure order. For this test, the mean departure order rank ± standard deviation per uncertainty behaviours is given per country.

A generalized linear model (GLM) was used to test the influence of the independent variables on our response variables. We first analysed the uncertainty time using a Gamma law (link = log) with the following independent variables: the age and the gender of the hesitating pedestrian, whether s/he was using a phone, whether s/he was accompanied by at least one other person, the country, the number of lanes, the number of waiting pedestrians, the order of pedestrian departure, the [country*gender] interaction and the [country*accompanied] interaction. Another possible biological and viable interaction would have been [country*age], but this was not tested due to insufficient data per age category. We then analysed the probability of abandon or acceleration using a Binomial law (abandon noted as 1, acceleration noted as 0) and with the following dependent variables: the age and the gender of the hesitating pedestrian, whether s/he was using a phone, whether s/he was accompanied by at least one other person, the country, the number of lanes, the number of waiting pedestrians and the order of pedestrian departure. For each GLM, we ran multi-model inferences to rank candidate models according to (i) their respective Akaike information criterion after correction for small sample sizes (AICc) and (ii) normalized Akaike weights (AICw) [START_REF] Burnham | Multimodel inference understanding AIC and BIC in model selection[END_REF]. ΔAICc is the difference in AICc between a given model and the model with the lowest AICc. AICw indicates the probability that a given model will be the best among candidate models [START_REF] Burnham | Multimodel inference understanding AIC and BIC in model selection[END_REF][START_REF] Burnham | AIC model selection and multimodel inference in behavioral ecology: some background, observations, and comparisons[END_REF][START_REF] Johnson | Model selection in ecology and evolution[END_REF].

All models were considered, from the null model (0 variable) to the full models (all tested variables). Models were ranked and compared according to their AICc, with the lowest AICc indicating the best model. Models with a difference of AICc (ΔAICc) <4 were considered equally possible candidates and their statistics (estimate, z-score and p-value) were averaged according to Bartoń's method (2013). The null model was also included as a possible candidate but was never among the models with lowest AICc. Averaged model coefficients were obtained for models with a ΔAICc < 4. Model inference and averaging were carried out with the R package 'MuMIn' [START_REF] Bartoń | MuMIn: multi-model inference. R package version. 1[END_REF].

A Pearson correlation test was used to assess whether the relative number of abandons and accelerations increased or decreased significantly according to the pedestrian's departure order (this analysis was carried out until the sixth rank, as no pedestrians were observed hesitating after this rank order).

Statistical analyses were performed using Rstudio 1.0.143 (Team RStudio, 2015) and the significance level was set at 0.05.

RESULTS

a) Factors affecting uncertainty behaviour

We observed 80 uncertainty behaviours for 1,464 infractions at Strasbourg (i.e.

5.46%) whilst 38 uncertainty behaviours were observed for 369 infractions at Nagoya (10.29%), which is different (chi²=11.4, p<0.001). The gender of the pedestrian had no effect on the likelihood of hesitating while crossing against a red light in either country (chi²<0.48, p>0.490). For the age variable, a zig-zag product on the likelihood of hesitating was observed on the graph for Japanese pedestrians. This is probably due to the small dataset with uncertainty behaviours. The age graph for French pedestrians seemed more coherent, with the highest number of uncertainty behaviours recorded for 20-year-old pedestrians and the lowest for their 60-year-old counterparts (Figure 2). The tendency to hesitate (K-s test, W = 6074, pvalue = 0.06178) was observed in the first pedestrians to cross the road in France (average±stdv=1.8±1.2), whilst the same tendency was shown by pedestrians following individuals who had already started crossing in Japan (average±stdv=2.2±1.6). We observed ten hesitating individuals using their phones (8.4% of hesitating individuals), which is no higher than the rate of individuals using their phones whilst crossing at a red pedestrian light (8.6%).

The uncertainty time GLM showed that this dependent variable is influenced by the country (Table 2, z=2.134, p=0.033), by being accompanied (z=2.526, p=0.012) and by the [country*accompanied] interaction (z=2.087, p=0.037). Uncertainty time appears to be longer in Japan, and longer when the pedestrian is alone in both countries, but more importantly this last effect is enhanced in Japan (i.e. a stronger effect was observed in Japan, with longer uncertainty times when alone). Other factors (see Table 2) did not affect the time of uncertainty.

b)

Factors affecting subsequent behaviours

After the 118 scored uncertainty behaviours, pedestrians showed 90 abandons (76.3%) and 28 accelerations (23.7%). Abandoning a crossing therefore appears then more probable than accelerating. Surprisingly, the percentage of accelerations is higher in Japan compared to France (34.2% and 18.7%, respectively). The GLM for the types of subsequent behaviours (abandoning or accelerating, Table 3) showed that only the number of lanes [2x1] is significant (z=3.445, p=0.001): the abandon rate increased with the number of lanes, whilst the acceleration rate decreased (Figure 3). Due to our small sample size, we did not obtain a significant effect for [2x2] and [2x3] lanes (where only abandon behaviours were observed).

Concerning the departure order (z=1.682, p=0.093), a Spearman correlation test on the percentage of abandons/accelerations according to the pedestrian departure order revealed that the relative number of abandons decreases with the departure order, whilst the relative number of accelerations increases (rhoabandon=-0.86, rhoacceleration=0.86, p=0.026, Figure 4).

DISCUSSION

This study highlights the factors influencing the uncertainty behaviour of pedestrians crossing against the red light and the subsequent behaviours, namely abandoning or accelerating the crossing. Japanese citizens are known to respect social rules more than

French or even western citizens [START_REF] Benedict | The Chrysanthemum and the Sword: Patterns of Japanese Culture[END_REF][START_REF] Komiya | A cultural study of the low crime rate in Japan[END_REF][START_REF] Markus | Culture and the Self -Implications for Cognition, Emotion, and Motivation[END_REF][START_REF] Ito | Socio-cultural backgrounds of Japanese interpersonal communication style[END_REF]. These cultural differences are not explained by differences in environmental factors but rather by the emergence of social conventions (Centolla and Baronchelli, 2015). A previous study showed that only 2.1% of crossings at signalized intersections are illegal in Japan, compared to 41.9% in France [START_REF] Pelé | Cultural influence of social information use in pedestrian road-crossing behaviours[END_REF]. However, only 5% of French pedestrians showed uncertainty behaviours when crossing against the red light, compared to 10% of the Japanese pedestrians. It is important to note that it was impossible to carry out observations of a road with 3*2 lanes in Strasbourg to compare with our observations in Nagoya. This limitation could decrease the ratio of illegal crossings and the ratio of uncertainty behaviours in Japan compared to France.

The uncertainty time is also longer in Nagoya (Japan, 1.55±1.45 s) than in Strasbourg (France, 1.22±0.97 s). These two elements -higher rate and longer time in Japan -lead us to suggest that Japanese pedestrians who crossed illegally mostly did so by following other pedestrians. Indeed, hesitating seems to occur in Japan when a pedestrian follows others already crossing against the red light (rather than due to inattention from using their phone, see explanations below) then realises the misevaluation. [START_REF] Pelé | Cultural influence of social information use in pedestrian road-crossing behaviours[END_REF] have already described this phenomenon in Japanese pedestrians. Indeed, pedestrians show a strong mimetism by following social information when they cross the road, sometimes inadvertently crossing against the red light [START_REF] Faria | Collective behavior in road crossing pedestrians: the role of social information[END_REF][START_REF] Pelé | Cultural influence of social information use in pedestrian road-crossing behaviours[END_REF]. This could explain why surprised pedestrians hesitate when they realise that they are crossing illegally. The information they obtain after starting to cross might be personal (by checking the light colour or seeing a car arriving) or social (noticing that other pedestrians do not cross). These conclusions may be confirmed by the fact that accompanied pedestrians in Japan showed a lower uncertainty time than pedestrians who crossed alone: the accompanied pedestrian could be therefore helped by social information. These results and previous findings [START_REF] Pelé | Cultural influence of social information use in pedestrian road-crossing behaviours[END_REF]) confirm a strong inter-country difference in the causes of illegal road-crossings and the causes of uncertainty behaviour. [START_REF] Pelé | Cultural influence of social information use in pedestrian road-crossing behaviours[END_REF] showed that whilst most pedestrians deliberately cross against the red light in France, they appear to do so accidentally in Japan.

Many studies have shown an age and gender effect on crossing behaviours [START_REF] Tom | Gender differences in pedestrian rule compliance and visual search at signalized and unsignalized crossroads[END_REF][START_REF] Holland | The effect of age, gender and driver status on pedestrians' intentions to cross the road in risky situations[END_REF][START_REF] Oxley | Differences in traffic judgements between young and old adult pedestrians[END_REF]. However, no such effects were observed in our study for the uncertainty behaviour rate, or for the uncertainty time and its possible subsequent behaviours. We suggest that this absence of effect may be due to two consecutive indirect filters resulting from our hesitating pedestrian population and linked to the data set: we studied uncertainty behaviours in a population of individuals who showed such behaviours (compared to those who did not show them). In turn, this population was part of a larger population of individuals: those crossing at a red light (compared to those who crossed at the green light). We assume that these two consecutive filters resulting from our dataset could lead us to solely observe pedestrians presenting a certain level of risk taking.

For instance, women and/or elderly pedestrians, known to take fewer risks than men or adolescents [START_REF] Sueur | Different risk thresholds in pedestrian road crossing behaviour: a comparison of French and Japanese approaches[END_REF][START_REF] Faralla | Risk Taking and Social Exposure[END_REF]Wilson and Dally, 1985) did not cross at the red light or rarely did so, thus could bias the dataset and providing non-significant results

for age and/or sex effect on uncertainty. We also expected men to display more accelerations than women, but this gender effect does not appear to be statistically significant, again probably due to the aforementioned filtering. We checked this filtering hypothesis in our dataset: this was checked with Japan pedestrians but not with French ones. Indeed, we observed a change in gender-ratio between total crossings and illegal crossings and between illegal crossings and uncertainty behaviours in Japan but this was not the case with France (see Table 4). Considering a possible age effect, we expected more uncertainty behaviours from the youngest (10-year-old) and oldest (more than 60-year-old) pedestrians than in other age categories. However, the results were inconsistent in Japan (Figure 2a) and were unexpected in France (Figure 2b), with 20-year-old French pedestrians showing the highest rate of uncertainty behaviours and 60-year-old pedestrians displaying the lowest rate of uncertainty behaviours. This result may be explained by our low sample size, but we cannot exclude the possible effect of decreasing cognitive abilities with age [START_REF] Brown | Driving and dementia: a review of the literature[END_REF][START_REF] Farias | Degree of discrepancy between self and other-reported everyday functioning by cognitive status: Dementia, mild cognitive impairment, and healthy elders[END_REF][START_REF] Iosa | Development and decline of upright gait stability[END_REF][START_REF] Salthouse | The processing-speed theory of adult age differences in cognition[END_REF], responsible for a lack of reaction or absence of uncertainty behaviour in the elderly when in risky situations [START_REF] Holland | The effect of age, gender and driver status on pedestrians' intentions to cross the road in risky situations[END_REF][START_REF] Oxley | Differences in traffic judgements between young and old adult pedestrians[END_REF].

Whilst we did not observe an effect of environmental variables on the uncertainty time, we did observe an effect on the probability of an acceleration or abandon behaviour after uncertainty. A higher number of lanes resulted in a decrease in the probability that a pedestrian would accelerate and an increase in the probability that they would abandon the crossing. Only abandon scenarios were displayed by pedestrians illegally crossing roads with six lanes (three lanes in each direction). As the number of lanes increases, the number of vehicles and subsequently the risk of accidents increases. Thereby, abandoning the crossing appears to be a safer behaviour for pedestrians than accelerating. Indeed, we decrease our perception abilities when rushing, and we thus increase the probability of an error of judgement. This processing-speed theory [START_REF] Salthouse | The processing-speed theory of adult age differences in cognition[END_REF] -where accuracy of decisions decreases with speed [START_REF] Pelé | Decision-making theories: linking the disparate research areas of individual and collective cognition[END_REF][START_REF] Nagengast | Risk sensitivity in a motor task with speed-accuracy trade-off[END_REF][START_REF] Sueur | Risk should be objectively defined: reply to Zentall and Smith[END_REF] has been well demonstrated for driving speed but can also be confirmed in pedestrians [START_REF] Aarts | Driving speed and the risk of road crashes: A review[END_REF][START_REF] Kloeden | Reanalysis of Travelling Speed and the Risk of Crash Involvement in Adelaide South Australia[END_REF][START_REF] Kloeden | Travelling Speed and Risk of Crash Involvement on Rural Roads[END_REF].

Our results also show that the probability of abandoning decreases with the departure order of the pedestrian (i.e. with the number of pedestrians who have already crossed or started crossing) while the probability to accelerate increases. This mimetic effect is due to conformity bias [START_REF] Cialdini | Social influence: Compliance and conformity[END_REF], and increases the risk of being injured or killed. Referred to as "sheep" or "herd" behaviour, this copying behaviour is well known in human beings in domains such as financial markets, fashion, purchasing or crowd behaviour [START_REF] Banerjee | A Simple Model of Herd Behavior[END_REF][START_REF] Bikhchandani | A Theory of Fads, Fashion, Custom, and Cultural Change as Informational Cascades[END_REF][START_REF] Chen | Herd behavior in purchasing books online[END_REF][START_REF] Hirshleifer | Herd Behaviour and Cascading in Capital Markets: a Review and Synthesis[END_REF][START_REF] Mawson | Understanding Mass Panic and Other Collective Responses to Threat and Disaster[END_REF][START_REF] Raafat | Herding in humans[END_REF]. It seems to be deeply rooted in human behaviour due to the gregariousness and sociality of the human species, and can be observed in many other social animal species [START_REF] Couzin | Self-organization and collective behavior in vertebrates[END_REF][START_REF] King | Where Next? Group Coordination and Collective Decision Making by Primates[END_REF][START_REF] Krause | Living in Groups[END_REF].

However, it may also reflect a sense of safety in numbers; car drivers are more likely to see a group and are hence more likely to slow down or stop, even if they have right of way [START_REF] Harrell | Factors influencing pedestrian cautiousness in crossing streets[END_REF][START_REF] Rosenbloom | Crossing at a red light: Behaviour of individuals and groups[END_REF].

Finally, contrary to what we expected and what has been described in the literature (McKnight and McKnight, 1993;[START_REF] Neider | Pedestrians, vehicles, and cell phones[END_REF]Strayer et al., 2003), we did not observe a decrease in attention when pedestrians were using a phone. Only ten hesitating individuals were observed using their phone (8.4% of the sample size). This does not seem no higher than the rate of individuals using a phone whilst crossing against the red light (8.6%), and these pedestrians did not show a higher uncertainty time than other pedestrians. The rate recorded in our study is lower than the 20% of distracted pedestrians found by Bungum and colleagues (2005). Again, this result may be due to the low sample size caused by the double filtering (see above).

Overall, combined with the findings of previous studies [START_REF] Pelé | Cultural influence of social information use in pedestrian road-crossing behaviours[END_REF], these results indicate the existence of two categories of rule breakers among pedestrians, namely those who deliberately do so (as described in [START_REF] Pelé | Cultural influence of social information use in pedestrian road-crossing behaviours[END_REF]) and those who do so inadvertently (shown here through uncertainty behaviour). This second category appears to represent the most vulnerable pedestrians, because they seem to misevaluate a safety crossing.

Additional signals could be set up to provide more adequate protection of pedestrians at accident-prone intersections. Sound signals have already been set up to help impaired or blind people [START_REF] Tauchi | Audible traffic signal for visually impaired persons using multiple sound outputs[END_REF][START_REF] Tauchi | A new audible traffic signal that helps blind pedestrians to cross intersections by providing navigation clues[END_REF]. In our case, we could imagine less pleasant signals such as the new water spray systems that surprised rule-breaking pedestrians in China's Hubei Province (BBC News). Further research will then be needed to measure the success rate of such innovative systems in saving lives.

TABLES

Figure 1 :

 1 Figure 1: Temporal scheme of the different behaviours scored in this study.

Figure 2 :

 2 Figure 2: Differences in the percentages of 118 pedestrians showing uncertainty behaviours (N= 80 for France and N= 38 for Japan) and 1,278 pedestrians of a sample size (N=551 for France and N= 727 for Japan) in a. Japan and b. France. Positive values showed that more pedestrians in a given age category showed uncertainty behaviour than expected; negative values showed that fewer pedestrians in a given age category showed uncertainty behaviour than expected.

Figure 3 :

 3 Figure 3: Percentage of abandons (black bars) and accelerations (grey bars) according to the number of lanes.

Figure 4 :

 4 Figure 4: Rate of abandonment (black bars) and acceleration (grey bars) according to the departure order of the pedestrian.

Table 1 :

 1 Details concerning the studied sites in Strasbourg (France) and Nagoya (Japan). A total of 80 and 38 uncertainty behaviours were observed in France and in Japan, respectively, on the lines "Number of uncertainty behaviours counted during 6 hours of observation".

	Pictures of the different intersections are available in the Supplementary Material section of
	Pelé et al. 2017, via this link: https://doi.org/10.6084/m9.figshare.c.3677176.v1
		Strasbourg-France
	Sites	Train Station	Pont des Corbeaux	Place Broglie
	Coordinates	48.584474, 7.736135	48.579509, 7.750745	48.584559, 7.748628
	Number of lanes in each direction	1	2	1
	Mean pedestrian flow per hour	667	612	850
	Mean road-crossing speed (m.s-1 )	0.96±0.05	1.11±0.29	1.01±0.16
	Data collection dates	02/07-07/07/2014 01/10-25/10/2014 15/02-09/03/2015
	Number of uncertainty behaviours
	counted during 6 hours of	36	23	21
	observation
	Percentage of illegal crossings	58.29%	32.92%	44.02%
		Nagoya-Japan
	Sites	Train Station	Maruei	Excelco	Osu-Kannon
	Coordinates	35.170824, 136.884328	35.168638, 136.905740	35.166891, 136.907284	35.159316, 136.901697
	Number of lanes in each direction	3	one-way street	1	1
	Mean pedestrian flow per hour	480	645	869	814
	Mean road-crossing speed (m.s-1 )	1.10±0.22	1.15±0.21	0.98±0.21	1.07±0.18
	Data collection dates	13/06-05/07/2011		27/01-05/02/2015
	Number of uncertainty behaviours
	counted during 6 hours of	3	15	4	16
	observation
	Percentage of illegal crossing	8.34%	14.53%	2.11%	8.86%

Table 2 :

 2 Results of the model selection (see statistical analyses section) with the uncertainty time (in seconds) as a dependent variable. Considering the number of lanes, the reference chosen is [1x1].

	Independent variables	Estimate Std. Error	Adjusted SE	z value Pr(>|z|)
	(Intercept)	0.869	0.281	0.284	3.064 0.002
	Accompanied [Yes]	-0.733	0.288	0.290	2.526 0.012
	Country [France]	-0.629	0.292	0.295	2.134 0.033
	Accompanied [Yes]* Country [France]	0.660	0.313	0.316	2.087 0.037
	Number of waiting pedestrians	-0.014	0.014	0.014	1.033 0.302
	Gender [Man]	-0.105	0.191	0.193	0.543 0.587
	Number of lanes [2x1]	0.368	0.286	0.289	1.275 0.202
	Number of lanes [2x2]	0.478	0.323	0.326	1.465 0.143
	Number of lanes [2x3]	-0.224	0.487	0.493	0.455 0.649
	Using a phone [Yes]	-0.157	0.260	0.263	0.595 0.552
	Age	-0.003	0.005	0.005	0.547 0.584
	Departure Order	-0.005	0.054	0.054	0.099 0.921
	Gender [Man]*Country [France]	-0.200	0.303	0.306	0.652 0.514

Table 3 :

 3 Results of the model selection (see statistical analyses section) with the type of behaviour following uncertainty as dependent variable (abandon = 1, acceleration = 0).

	Independent variables	Estimate Std. Error Adjusted SE z value Pr(>|z|)
	(Intercept)	-0.574	1.067	1.075	0.534 0.593
	Number of lanes [2x1]	2.882	0.829	0.837	3.445 0.001
	Number of lanes [2x2]	0.625	0.874	0.883	0.708 0.479
	Number of lanes [2x3]	17.015	1377.357	1392.471	0.012 0.990
	Departure Order	-0.303	0.178	0.180	1.682 0.093
	Age	0.023	0.018	0.018	1.273 0.203
	Gender [Man]	0.648	0.550	0.556	1.166 0.244
	Number of waiting pedestrians	-0.054	0.050	0.051	1.072 0.284
	Country [France]	0.670	0.832	0.841	0.797 0.425
	Using a phone [Yes]	0.490	0.984	0.995	0.492 0.623
	Accompanied [Yes]	-0.179	0.558	0.563	0.319 0.750

Table 4 :

 4 Relative percentage of men and women in France and Japan according to the total crossing (at red and green light), the illegal crossing (red light) and the uncertainty behaviours.

			Total	Illegal	Uncertainty
			crossing	crossing	behaviours
	France	Men	44.4	44.4	46.25
		Women	55.6	55.5	53.75
	Japan	Men	42.2	59	42.1
		Women	57.7	41	57.9

ACKNOWLEDGEMENTS

We thank Caroline Bellut, Elise Debergue, Charlotte Gauvin, Anne Jeanneret, Thibault Leclere, Lucie Nicolas, Florence Pontier and Diorne Zausa for their help collecting data, and Cédric Sueur for his help with the analyses. We are grateful to Kunio Watanabe and Hanya Goro (Primate Research Institute, Kyoto University) for their help in authorizing data collection and to Joanna Munro for English language editing. We also thank the three anonymous referees for their constructive comments, which helped to improve this paper.


Mathieu Fauvel
email: mathieu.fauvel@inra.fr

Mailys Lopes
email: ys.lopes@inra.fr

Titouan Dubo
email: dubo@hotmail.fr

Justine Rivers-Moore
email: justine.rivers-moore@inra.fr

Pierre-Louis Frison
email: pierre-louis.frison@univ-mlv.fr

Nicolas Gross
email: nicolas.gross@inra.fr

Annie Ouin
email: annie.ouin@ensat.fr

Prediction of plant diversity in grasslands using Sentinel-1 and -2 satellite image time series

Keywords: Satellite image time series, Sentinel-1&-2, Statistical learning, Grasslands, taxonomic diversity

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

Introduction

Natural and semi-natural grasslands cover 22% of the European agricultural land surface [START_REF] Bengtsson | Grasslands-more important for ecosystem services than you might think[END_REF]. Grasslands are today one of the most endangered ecosystems due to land use change, agricultural intensification, and abandonment [START_REF] Pärtel | Biodiversity in temperate european grasslands: origin and conservation[END_REF]. Grasslands host a unique biodiversity, which support the provision of key ecosystem services such as carbon storage, food production crop polination, pest regulation and contribute to landscape scale amenities such as landscape scenic view. Permanent grasslands are crucial to maintain biodiversity in many agricultral landscape [START_REF] Watkinson | Grasslands, grazing and biodiversity: editor's introduction[END_REF][START_REF] Habel | European grassland ecosystems: threatened hotspots of biodiversity[END_REF], and particularly insect diversity: they provide host plants, nectar, pollen [START_REF] Ockinger | Semi-natural grasslands as population sources for pollinating insects in agricultural landscapes[END_REF] as well as nest sites [START_REF] Carrié | Bee diversity in crop fields is influenced by remotely-sensed nesting resources in surrounding permanent grasslands[END_REF][START_REF] Potts | Abiotic and biotic factors influencing nest-site selection by halictus rubicundus, a ground-nesting halictine bee[END_REF]).

Yet, despite recent regulation in favor of grasslands, agriculture intensification, abandonment and urbanization generate a decrease of both grasslands surface area and plant diversity [START_REF] O'mara | The role of grasslands in food security and climate change[END_REF][START_REF] Newbold | Has land use pushed terrestrial biodiversity beyond the planetary boundary? a global assessment[END_REF] leading to a loss of biodiversity (plants, insects, birds) and its related ecosystem services. Therefore, there is an urgent need to monitor grasslands over large extents in order to assess their current state in terms of area covered but also to estimate their plant diversity and their carrying capacity for insect populations. While the area information can be extracted from land cover databases (e.g. CORINE Land Cover [START_REF] Heymann | CORINE land cover: Technical guide[END_REF]) diversity and carrying capacity data are scarcely available. Though, spatial location, extent and quality of grasslands are crucial for their conservation per se but also to produce maps of related ecosystem services in agricultural areas. For instance, mobile agents based ecosystem services [START_REF] Kremen | Pollination and other ecosystem services produced by mobile organisms: a conceptual framework for the effects of land-use change[END_REF] such as pollination and pest regulation rely on the presence of a mosaic of crops and semi-natural elements among which grasslands are crucial and represent the largest surface area.

Plant diversity is usually assessed by botanical surveys of the grasslands and indices related to the diversity in species of a community are commonly computed (richness, Shannon and Simpson diversity [START_REF] Magurran | Measuring Biological Diversity[END_REF]). Functional diversity can also be derived from the botanical survey through functional types, i.e., set of species exhibiting similar attributes having similar effect on ecosystems [START_REF] Diaz | Vive la difference: Plant functional diversity matters to ecosystem processes[END_REF]. Field surveys, sometimes called ground-truth in the remote sensing community, require significant human and material resources, the knowledge of the assessor and a sampling strategy, which make them expensive and time consuming [START_REF] Magurran | Measuring Biological Diversity[END_REF]. Ecological surveys are thus limited in spatial extent and in temporal frequency, limiting grassland monitoring to a local scale and usually over a short period of time. Although field surveys provide valuable and high quality data at a point scale, they cannot easily be upscaled while taking into account the landsape heterogeneity. Therefore, field surveys alone cannot address the needs to monitor grassland biodiversity over large spatial extents and other techniques relating to field surveys should be considered.

Among such complementary techniques, remote sensing has received an increasing attention in the last decades [START_REF] Turner | Sensing biodiversity[END_REF][START_REF] Buck | Using information layers for mapping grassland habitat distribution at local to regional scales[END_REF]. However, grasslands have been rarely studied in the remote sensing literature compared to other land covers like crops or forest [START_REF] Newton | Remote sensing and the future of landscape ecology[END_REF]. Most of the studies focusing on grasslands have agronomic applications, such as estimating biomass productivity and growth rate [START_REF] Gu | Developing a 30-m grassland productivity estimation map for central nebraska using 250-m modis and 30-m landsat-8 observations[END_REF][START_REF] Li | Monitoring and modeling spatial and temporal patterns of grassland dynamics using time-series modis ndvi with climate and stocking data[END_REF] or deriving biophysical parameters like the Leaf Ara Index (LAI), the fraction of Photosynthetically Active Radiation (fPAR) or the fraction of Vegetation Cover (fCOVER) [START_REF] Darvishzadeh | Lai and chlorophyll estimation for a heterogeneous grassland using hyperspectral measurements[END_REF][START_REF] He | Reflectance measures of grassland biophysical structure[END_REF].

Regarding ecological applications, Corbane et al. ( 2013) mapped semi-natural grasslands types using two RapidEye images, while [START_REF] Stenzel | Identification of high nature value grassland with remote sensing and minimal field data[END_REF] used five RapidEye images to classify high nature value grasslands. Species richness of herbaceous ecosystems such as grasslands, prairies and savannah was assessed using airborne hyperspectral imagery in [START_REF] Möckel | Airborne hyperspectral data predict fine-scale plant species diversity in grazed dry grasslands[END_REF], [START_REF] Gholizadeh | Detecting prairie biodiversity with airborne remote sensing[END_REF] and [START_REF] Oldeland | Does using species abundance data improve estimates of species diversity from remotely sensed spectral heterogeneity?[END_REF], respectively. Airborne hyperspectral imagery has also been used to map groups of species [START_REF] Schuster | Grassland habitat mapping by intraannual time series analysis -comparison of rapideye and terrasar-x satellite data[END_REF], flower types [START_REF] Landmann | Application of hyperspectral remote sensing for flower mapping in african savannas[END_REF] or pollination types [START_REF] Feilhauer | Mapping pollination types with remote sensing[END_REF].

Schematically, the different approaches to analyse grassland with remote sensing data can be divided into three categories: 1. Techniques that use vegetation indices or spectral reflectance as explanatory variables [START_REF] Rapinel | Evaluation of sentinel-2 time-series for mapping floodplain grassland plant communities[END_REF][START_REF] Féret | Detecting the phenology and discriminating mediterranean natural habitats with multispectral sensors-an analysis based on multiseasonal field spectra[END_REF] 2. Techniques that use the spatial heterogeneity of vegetation indices as explanatory variables [START_REF] Goodin | Seasonality of finely-resolved spatial structure of ndvi and its component reflectances in tallgrass prairie[END_REF][START_REF] Tuanmu | A global, remote sensing-based characterization of terrestrial habitat heterogeneity for biodiversity and ecosystem modelling[END_REF], 3. Techniques that use the "spectral variation hypothesis" (SVH) to extract explanatory variables [START_REF] Lopes | Spectro-temporal heterogeneity measures from dense high spatial resolution satellite image time series: Application to grassland species diversity estimation[END_REF][START_REF] Hall | Spectral heterogeneity of quickbird satellite data is related to fine-scale plant species spatial turnover in semi-natural grasslands[END_REF].

In terms of remote sensing data, the majority of the studies having ecological schemes were conducted using very high spatial resolution data (≤ 5 meters/pixel), e.g., hyperspectral data issued from an airborne sensor or multispectral satellite data such as RapidEye. Although a very high spatial resolution might seem necessary to estimate grasslands plant diversity as plant are typically small in grasslands, these types of acquisitions are limited in time and in space because of their current cost. In addition, medium spatial resolution satellite image time series (e.g, MODIS or SPOT-VEGETATION) have been found not appropriate for small and heterogeneous elements, such as grasslands. In particular, [START_REF] Ali | Satellite remote sensing of grasslands: from observation to management[END_REF] have stated that the moderate spatial resolution "precludes their use for field-scale application in many countries". As a consequence, the grassland biodiversity monitoring from remote sensing data is nowadays limited to "low temporal resolution & small spatial cover & high spatial resolution" [START_REF] Giménez | Determination of grassland use intensity based on multi-temporal remote sensing data and ecological indicators[END_REF][START_REF] Wang | Mapping foliar functional traits and their uncertainties across three years in a grassland experiment[END_REF][START_REF] Gholizadeh | Detecting prairie biodiversity with airborne remote sensing[END_REF][START_REF] Lopatin | Mapping plant species in mixed grassland communities using close range imaging spectroscopy[END_REF] or "high temporal resolution & large spatial cover & coarse spatial resolution" [START_REF] Ge | Modeling alpine grassland cover based on modis data and support vector ma-chine regression in the headwater region of the huanghe river, china[END_REF][START_REF] Schmidt | Mapping spatio-temporal dynamics of the cover and management factor (c-factor) for grasslands in switzerland[END_REF][START_REF] Yang | Modeling grassland above-ground biomass based on artificial neural network and remote sensing in the three-river headwaters region[END_REF][START_REF] Barrett | Assessment of multi-temporal, multi-sensor radar and ancillary spatial data for grasslands monitoring in ireland using machine learning approaches[END_REF][START_REF] Si | Mapping spatio-temporal variation of grassland quantity and quality using meris data and the prosail model[END_REF].

The launch of Sentinel-1&-2 satellites, issued from the European Union Copernicus and European Space Agency program, has enabled continuous and regular monitoring of small vegetated areas over large spatial extents thanks to their high spatial resolution (10m per pixel) and their frequent revisit (every 5 days for Sentinel-2) [START_REF] Drusch | Sentinel-2: ESA's optical high-resolution mission for GMES operational services[END_REF][START_REF] Defourny | Near real-time agriculture monitoring at national scale at parcel resolution: Performance assessment of the sen2-agri automated system in various cropping systems around the world[END_REF]. Sentinel-1&-2 satellite image time series (SITS) open new possibilities to bridge the gap between high spatial resolution and high temporal resolution. First, their high temporal resolution enables the monitoring of plant phenology. Because plant communities differ in their temporal behavior, the phenological diversity measured by these sensors could be used as a proxy to estimate the plant diversity. Second, they provide complementary information about the grassland cover: optical data are influenced by chemical composition of the elements while radar data are sensitive to the geometric configuration of the elements [START_REF] Joshi | A review of the application of optical and radar remote sensing data fusion to land use mapping and monitoring[END_REF]. Then, being available for free, they can be used to process large areas and complement field surveys at a reduced cost. Yet, recent works using Sentinel-2 images are still performed on small areas (a few square kilometers) [START_REF] Rapinel | Evaluation of sentinel-2 time-series for mapping floodplain grassland plant communities[END_REF]. Although interesting, studies over such small areas do not have to cope with a high diversity of grasslands in terms of composition neither with the spatial variability of meteorological, topographical conditions and farmers pratices. Hence, it is unlikely that classification accuracies reported in such cases can be reached on larger spatial areas. Furthermore, they do not cope with the spectro-temporal variability of Sentinel satellite image series.

In this paper, we use Sentinel's optical and radar SITS to estimate plant taxonomic and functional diversity in grasslands over a a large spatial extent (40 000 km 2 ) for one growing season. Fifteen months of acquisitions of Sentinel-1&-2 images were used to predict at the pixel level plant diversity indices (richness, Shannon and Simpson diversity indices) and functional diversity indices related to pollination (flower diversity, insect dependence) of grasslands located in farmed landscapes in the southwest of France. Various linear and non-linear algorithms were investigated to cover different kinds of possible statistical relationships between the biodiversity indices and the SITS.

Materials

Study site

The study area is part of the Long-Term Socio-Ecological Research site "Pyrénées-Garonne" located in southwest France near the city of Toulouse (43 • 17 0 N, 0 • 54 0 E), see Figure 1. This hilly area of around 900 km 2 is characterized by a mosaic of crops, small woods and grasslands. It is dominated by mixed crop-livestock farming. Grasslands provide food for cattle by grazing and/or producing hay or silage. They range from mono-specific, annual grasslands sown with rye-grass (improved with mineral fertilizing and mown up to three times a year) to permanent, semi-natural grasslands composed of spontaneous plant species (not fertilized and mown once a year). Grasslands are mainly located on steep slopes, whereas annual crops are in the valleys on the most productive, often irrigated lands. The climate is sub-Atlantic with sub-Mediterranean and mountain influences (mean annual temperature, 12.5 • C; mean annual precipitation, 750 mm). Botanical surveys were conducted in 83 plots from distinct grasslands between April and May 2018 after the flowering and before the mowing. The average distance between two plots was 13.4 km with a standard deviation of 9.1 km and the minimum and the maximum distances were 0.08 km and 44.37 km, respectively. Each plot was composed of 5 quadrats of 1 square meter area (Figure 2). The center of each quadrat corresponded to the theoretical center of a Sentinel-2 pixel. Note that since the absolute geolocation is below 11 meter at 95.5% confidence (ESA, 2019), a possible mis-registration between Sentinel-2 data and GNSS can be observed for few pixels at some dates. After a ground checking and getting owners agreement, the position of the central quadrat was controlled with a high-accuracy GNSS handhelds (Trimble Geo 7X Ter-raSync). The grasslands were digitized using the agricultural Land Parcel Information System "Registre Parcellaire Graphique" (RPG) [START_REF] Cantelaube | Le registre parcellaire graphique: des données géographiques pour décrire la couverture du sol agricole[END_REF].

Field botanical surveys

Plants were identified at the species level and the abundance of each species was estimated by its cover percentage value in the one square meter area, according to 12 bins, see Table 1. Vegetation height was measured in two locations in each quadrat and averaged. A total of 415 quadrats (leading to 415 Sentinel pixels), belonging to 83 grasslands, was recorded and used in this study.

Three indices were computed using vegan package of R 3.2.3 [START_REF] Oksanen | Vegan: Community ecology package.ordination methods, diversity analysis and other functions for community and vegetation ecologists[END_REF][START_REF] Team | R: A language and environment for statistical computing[END_REF]: plant richness S q , defined as the number of plant species in quadrat q, the Shannon diversity index

H (q) = - Sq s=1 p qs ln(p qs ),
where p qs is the relative abundance of the species s in quadrat q and Simpson index of diversity (refers to Simpson index in the paper for simplicity) ,

L(q) = 1 - Sq s=1 p 2 qs .
The Shannon index measures heterogeneity taking into account the number of species and the relative abundance of each species, while the Simpson index is the inverse of a dominance index reflecting how much the community is dominated (or not) by few species. We also consider vegetation attribute related to pollination. Following [START_REF] Clough | Density of insect-pollinated grassland plants decreases with increasing surrounding land-use intensity[END_REF], Insect Dependence was assessed for each plant species, based on pollination syndrome (insects, wind, water..) and its frequency. The values were: 1 if insects are the major pollination vector, 0.5 if insects and an other vector are the major pollination vector, 0 if insects are not an important pollination vector. The Reward Index combines the amount of nectar and pollen which values run from: 0 for "none", 1 for "little", 2 for "present", 3 for "plenty". The sum of nectar and the amount of pollen score lead to values from 0 to 6 for the Reward Index of the species. Averaged at the community level using the FD package Laliberté and Legendre (2010) of R 3.2.3, the diversity and the abundance of these two vegetation attributes may reflect how species impact community and ecosystems. Finally, the Color Diversity was derived from the Shannon index,

CD(q) = - Jq j=1 r qj ln(r qj ),
where r qj is the sum of the relative abundance of the species whom major color of flowers is j in quadrat q and J q is the number of major colors in quadrat q.

The correlation matrix between the biodiversity indices is reported in Figure 3.

Satellite data

Both the optical and radar data were processed using the iota 2 software [START_REF] Inglada | iota2a25386[END_REF], built on the Or-feoToolbox (OTB Development Team, 2018).

Optical data

Four Sentinel-2 tiles of Level 2A were downloaded from the Theia Land Data Centre 1 , corresponding to all available acquisitions between August 26, 2017 and December 03, 2018 2 . Only the bands natively at 10 m/pixel were used in the experiments (blue, green, red and nearinfrared). Surface reflectance time series were produced using the MAJA (Multi-sensor atmospheric correction and cloud screening-ATCOR Joint Algorithm) processing chain developed by the CNES-CESBIO and DLR. It includes orthorectification, atmospheric correction, clouds and shadows detection [START_REF] Baetens | Validation of copernicus sentinel-2 cloud masks obtained from maja, sen2cor, and fmask processors using reference cloud masks generated with a supervised active learning procedure[END_REF]. All the acquisitions were re-sampled and gap-filled (to handle clouds and shadows) onto the same set of dates (every 10 days, starting from 2017-08-29 and ending in 2018-12-02), as proceeded in [START_REF] Inglada | Assessment of an operational system for crop type map production using high temporal and spatial resolution satellite optical imagery[END_REF]. The final SITS contains 47 dates, in the visible and near infrared bands, corresponding to a total of 188 spectro-temporal features. The normalized difference vegetation index (NDVI) and the first principal components, accounting for 99% of the total variance, were extracted as alternative features [START_REF] Tupin | Remote Sensing Imagery[END_REF].

Figure 4 presents the reconstructed NDVI (in black) of the five quadrats for four plots as well as their corresponding detected shadows and clouds.

Radar data

Sentinel-1 Level-1 GRD (Ground Range Detected) time series recorded in interferometric wide swath mode were used. Four orbits corresponding to the Sentinel-2 extent were downloaded from the French Sentinel collaborative ground segment 3 , corresponding to the same temporal period as Sentinel-2 time series. GRD time series were calibrated to γ 0 , orthorectified and filtered using a multitemporal spatial filter [START_REF] Quegan | Filtering of multichannel sar images[END_REF]. The parameters of the spatio-temporal filter were a spatial window of size 5×5 pixels and all the dates were used for the temporal filters. VV and VH polarizations were processed as well as the ratio VV/VH. The features were re-sampled onto the same set of dates (every 10 dates starting from 2017-08-27 to 2018-11-30).

Figure 4 presents the reconstructed γ 0 (VV) polarization (in blue) for descending orbits for the same four plots.

Methods

For this work, a set of radiometric features commonly used with optical and radar data is investigated. Table 2 shows the various combinations of features used in the regression algorithms 4 .

1 http://www.theia-land.fr/en/presentation/products 2 Images tagged as composed of more than 90% of cloudy pixels by the MAJA processing chain are not processed by the data center.

3 https://peps.cnes.fr/rocket/#/home 4 Dynamic habitat indices, such as in [START_REF] Hobi | A comparison of dynamic habitat indices derived from different modis products as predictors of avian species richness[END_REF], as well as VV/VH ascending/descending features alone were tested, but the results were too bad and are not reported in the paper for the sake of clarity. 117

2 0 1 7 -0 9 2 0 1 7 -1 0 2 0 1 7 -1 1 2 0 1 7 -1 2 2 0 1 8 -0 1 2 0 1 8 -0 2 2 0 1 8 -0 3 2 0 1 8 -0 4 2 0 1 8 -0 5 2 0 1 8 -0 6 2 0 1 8 -0 7 2 0 1 8 -0 8 2 0 1 8 -0 9 2 0 1 8 -1 0 2 0 1 8 -1 1 2 0 1 8 -1 2 0.0 0.2 0.4 0.6 0.8 1.0 -15 -10 -5 0 NDVI Plot 943_2647, Shannon index=0 γ O (V V )

Regression techniques

Five conventional machine learning algorithms were used in the experiments with different set-up [START_REF] Lary | Machine learning in geosciences and remote sensing[END_REF]: Linear regression (with ridge and lasso regularization) (Hastie et al., 2009, Chap. 3), K-Nearest Neighbors (K-NN) (Hastie et al., 2009, Chap. 13), Kernel Ridge Regression (KRR) (Murphy, 2013, Chap. 14), Random Forest (RF) (Hastie et al., 2009, Chap. 15) and Gaussian Process (GP) [START_REF] Rasmussen | Gaussian Processes for Machine Learning[END_REF]. K-NN, KRR, RF and GP are non-parametric and non-linear algorithms that have been used successfully in many remote sensing applications, see for instance [START_REF] Maxwell | Implementation of machine-learning classification in remote sensing: an applied review[END_REF]. In all the experiments, the Scikit-learn python library was used [START_REF] Pedregosa | Scikit-learn: Machine learning in Python[END_REF]. Hyperparameters of ridge, lasso, K-NN, RF and KR methods were fitted through crossvalidation (Hastie et al., 2009, Chap 7) and those for GP were found by maximizing of the likelihood [START_REF] Rasmussen | Gaussian Processes for Machine Learning[END_REF].

Statistical accuracy assessment

The accuracy of the prediction of a biodiversity variable value y was assessed using the coefficient of determination (r 2 ) [START_REF] Draper | Applied regression analysis[END_REF]:

r 2 = 1 - nv i=1 (y i -ŷi ) 2 nv i=1 (y i -ȳ) 2
where n v is the number of validation sample, y i is the true biodiversity variable value associated with the i th sample, ŷi is its corresponding predicted value and ȳ is the mean value of the validation sample.

In order to correctly estimate r 2 , the k-fold cross validation was used. Cross-validation directly estimates the generalization error of a given method (Hastie et al., 2009, Chap. 7): The data is split into K almost equal-size folds, and the model is fitted on K -1 folds. The prediction error, here r 2 , is computed on the remaining unseen fold. This process is done for all folds k = 1, 2, . . . , K. The average error over the K folds is the estimated generalization error. However, in our experimental setting, the 5 quadrats corresponding to one plot are spatially correlated and this must be taken into account when constructing the K folds of the cross-validation [START_REF] Hammond | Optimistic bias in classification accuracy assessment[END_REF]. It is proposed in this work to perform a grouped cross validation by splitting the samples according to their plot membership [START_REF] Pohjankukka | Estimating the prediction performance of spatial models via spatial k-fold cross validation[END_REF][START_REF] Le Rest | Spatial leave-one-out cross-validation for variable selection in the presence of spatial autocorrelation[END_REF]: quadrats from a same plot are all included in one fold and not shared among several folds. Hence, in the experiments, a five fold-grouped cross validation was performed, where the first four folds contain 16 plots (80 quadrats) and the last fold contains 19 plots (95 quadrats). The same splits were used for all experiments, i.e., the error assessment was done using the same folds whatever the methods, data sources and indices. Figure 5 shows the whole process.

RS

Results

The Appendix C shows the results for the prediction of each variable. For clarity, only the best results for each variable corresponding to an average r2 > 0.3 are presented in the core of the paper, Table 3. 3). The black line is the identity line, i. e., when ŷ = y; r2 = 0.45. The colormap corresponds to the confidence interval σ of the prediction.

Biodiversity indices prediction accuracy

Of all indices, four indices were predicted with an estimated r2 above 0.3 and three indices above 0.4. The best predictions are obtained for the taxonomic diversity indices, i.e., Richness, Shannon and Simpson. These indices are strongly correlated (Figure 3), but they are also correlated with Color Diversity, which is well predicted too. The two others indices, Insect Dependence and Reward Index are not correlated with the others ones and are less accurately predicted.

The scatter plot of the prediction for the best result is displayed in Figure 6, and the second and third best results are displayed in Appendix B, figures B.12 and B.13. They were built during the spatial cross validation. They also show the confidence interval of the Random Forest predictions for each plot (see the colormap). Their main tendency is that lower values of indices is overestimated (above the black line) while high values are underestimated (below the black line).

Impact of the regression algorithm

Linear method with ridge regularization performs the worst in terms of prediction accuracy. Lasso and K-NN provide in-between results, sometimes Lasso works better, e.g., for Shannon, sometimes K-NN works better, e.g., for Simpson and sometimes they perform equally bad than the other methods. Among linear models, Lasso provides the best results, but it does not reach the best accuracy and the corresponding estimated r2 are usually low. Non linear regression methods, such as RF, GP and kernel ridge, provide overwhelmingly the best results.

From Table 3, RF is the method that provides the best overall results for three indices with an estimated r2 above 0.4. Kernel methods follow closely. Gaussian process performs slightly better on average than Kernel ridge. Regarding the influence of the kernel used in the GP, no significant differences are observed. For instance, Shannon indices were predicted with a r2 of 0.39±0.20 for the Rational quadratic kernel and of 0.38±0.19 for the RBF kernel, see Table C.6.

From the different results, no best method clearly appears. Depending on the indices, they may perform almost equally well for some indices while for others indices one method clearly outperforms the other ones. For instance, for Richness in table C.4, the estimated highest r2 are 0.32 for the GP with a rational quadratic kernel and 0.34 for the Kernel Ridge with a RBF kernel. However, for Simpson (table C.7), RF provides the best prediction by far with a r2 of 0.45 while GP RBF and GP RQ reach only 0.31.

Influence of the data source

From the tables in Appendix C, the optical data alone leads in general to better results than the radar data alone, see as example the Color Diversity results. Radar features did not lead to results close to the optical ones, for indices predicted with an r2 significantly higher than 0.

For the optical features, the NDVI and the PCA features are less informative than the original set of spectrotemporal features and the prediction accuracy is in general lower. However, combining the NDVI and spectrotemporal features can improve the prediction performance marginally for some indices, as it can be seen for Shannon index, table C.6. Surprisingly, using only the red and infrared bands -those used to compute the NDVI -results in highest, or close to highest, accuracy for most of the biodiversity indices, see for instance in Table 3.

Combining optical and radar data has positive influence on the prediction accuracy for some indices. For instance in Table 3, the four last indices were best predicted using a full stack of S1 and S2 features. Again, applying PCA does not lead to an increase of the prediction accuracy.

Influence of the spatial auto-correlation

For the three best estimated indices, the prediction accuracy was also estimated using a conventional cross-validation, i.e., allowing plots to be split into several folds. Figure 7 shows the prediction accuracy with and without the spatial auto-correlation. It is clear from the results that the quality of the prediction is overestimated when spatial auto-correlation is not taken into account. For instance for the Simpson index, the estimated r2 is approximately 35% lower when the spatial cross-validation is used.

Simpson

Large scale map

Prediction of the Simpson index at the pixel level was done for all the grassland parcels identified from the Land Parcel Information System (RPG) over the four Sentinel-2 tiles with the data and method that resulted in the highest prediction accuracy (i.e., R-IR satellite data and RF algorithm, see Table 3). The confidence interval was also computed. Figure 8 shows an extract of the predicted Simpson value in a square of 1 km 2 centered on the plot 10_470. It also shows a very high spatial resolution aerial image from 2016 over the same extent as well as the confidence interval associated with the prediction of the Simpson index obtained with the Random Forest. Extracts for all plots are provided in the supplementary materials.

From the prediction, spatial heterogeneity in terms of Simpson value (and its associated confidence value) can be seen in the parcel, corresponding to intra-parcel heterogeneity of species distribution that can be captured with a spatial resolution of 10m/pixel. Furthermore, inter-parcel variability can be observed, even for adjacent fields.

The distribution of the predicted Simpson values over the entire area is reported in Figure 9 for the most representative types of grasslands identified on the RPG (i.e., fallows, permanent grassland, long rotation meadow, temporary grassland, pastoral surface). These types accounted for approximately 80% of the total grasslands area of the four tiles. The distribution of the Simpson index values exhibit strong differences w.r.t the grassland classes. For instance, permanent grassland, long rotation meadow and temporary grassland have in average the highest values while young fallow and pastoral surface have the lowest values. Interestingly, young and old fallow present significant differences in terms of distribution, young fallow having lower Simpson indices (i.e., lower plant diversity) than old fallow.

Discussion

By using dense S2 and S1 SITS, this work contributed to the remote sensing of grasslands biodiversity using new generation high resolution SITS. It follows the works of [START_REF] Féret | Detecting the phenology and discriminating mediterranean natural habitats with multispectral sensors-an analysis based on multiseasonal field spectra[END_REF] and [START_REF] Rapinel | Evaluation of sentinel-2 time-series for mapping floodplain grassland plant communities[END_REF] but rather than performing classification of grasslands types, the prediction of biodiversity indices and plant traits was considered as in [START_REF] Wang | Mapping foliar functional traits and their uncertainties across three years in a grassland experiment[END_REF]. Our study differed from the previously cited works by two important points. First, the area considered in our study site is much larger and thus accounts for more variability in the indices and in the remote sensing data, making the statistical learning problem more difficult but more realistic. Second, the spatial auto-correlation was taken into account when estimating the prediction error but leads to a lowest value (and rather more correct, i.e., less biased) of the prediction accuracy [START_REF] Le Rest | Spatial leave-one-out cross-validation for variable selection in the presence of spatial autocorrelation[END_REF].

Prediction performance for biodiversity indices

Prediction of diversity measures were more accurate for diversity indices (Simpson and Shannon) than for species richness indices (e.g. Richness) as found in [START_REF] Oldeland | Does using species abundance data improve estimates of species diversity from remotely sensed spectral heterogeneity?[END_REF]. The Simpson index gives more weight to species with higher proportions while the richness is only based on species presence/absence. The presence of rare species (about 1% of cover) in the under storey will influence species richness but would be hardly detected by Sentinel sensor while dominant species are more detectable. Dominance indices such as Simpson are thus more indicated to estimate plant diversity with remote sensors. Functional diversity measures (Color Diversity, Insect Dependence and Reward Indices) showed a r2 above 0.3, which is promising as their link to optical or radar signal are less clear than species diversity. Even if flower cover never exceeds 20% of the grassland area, it may have been detected by the satellite sensors. Temporal variability of the signal may be linked to the different phenophase (period of flowering) in the grasslands. The link between diversity indices functionally related to pollination and satellite signals is truly interesting and opens new research frontiers. In particular, the prediction of land carrying capacity for pollinators at large extent based on actual relation between signal and functionality rather than expertise on suitability for pollinators based on satellite derived land cover is of most interest.

The spatial resolution of Sentinel images (10m×10m per pixel) might be too coarse to capture fine variation of plant composition. [START_REF] Gholizadeh | Detecting prairie biodiversity with airborne remote sensing[END_REF] showed a loss For the fallow classes or long rotation meadow, number in brackets indicates their ages. For the pastoral surfaces, the string in brackets indicates the most predominant resources in the parcels. The densities have been estimated by a non-parametric kernel density estimator with a Gaussian kernel and using the Scott rules for the selection of the bandwidth [START_REF] Scott | Multivariate Density Estimation: Theory, Practice, and Visualization[END_REF].

of correlation between the spectral feature and α-diversity index (Shannon) for simulated hyperspectral image with a spatial resolution above 2m per pixel. Nonetheless, they used only one date and captured little information about phenology. [START_REF] Feilhauer | Assessing floristic composition with multispectral sensors-a comparison based on monotemporal and multiseasonal field spectra[END_REF] found that spectral coverage of the entire solar-reflective domain is the most important characteristic to successfully assess floristic variation in grasslands, while multi-seasonal data did not improve it. These results obtained from simulated data contradict ours obtained with satellite imagery where the best prediction accuracies were reached using only the red and infrared bands and SITS. Our prediction accuracies are equivalent to the ones of [START_REF] Möckel | Airborne hyperspectral data predict fine-scale plant species diversity in grazed dry grasslands[END_REF] in predicting the inverse Simpson index in grasslands using airborne hyperspectral data. Their approach based on spectral reflectance reached a r2 of 0.45 using 245 wavebands and 0.40 using 35 wavebands. Our results suggest that dif-ferences in phenology can be used as a proxy for plant diversity when spatial resolution and spectral resolution are not high enough.

Similar ideas were considered by [START_REF] Rapinel | Evaluation of sentinel-2 time-series for mapping floodplain grassland plant communities[END_REF] who assumed that differences in phenology can be accounted for using SITS. However, in large scale setting, clouds and shadows are unavoidably present in the optical images, particularly during the key period for plant flowering (April-June), resulting in missing values that are reconstructed during the pre-processing. Figure 10 shows the proportion of masked pixels corresponding to the 415 quadrats: for the period between 2018-05-21 and 2018-06-20, the information was missing because of a too high level of clouds proportion (see section 2.3.1) and thus the corresponding dates were interpolated [START_REF] Inglada | Assessment of an operational system for crop type map production using high temporal and spatial resolution satellite optical imagery[END_REF]. Hence, differences in terms of phenology during that period of time were not captured by the optical sensors and could explain why the indices based on flower composition only were not accurately predicted. This also suggests that other drivers influencing the botanical composition are captured by the optical SITS, such as agricultural practices [START_REF] Moog | Comparison of species composition between different grassland management treatments after 25 years[END_REF].

By comparison with our previous works (Lopes et al., 2017a,b), in which biodiversity indices were surveyed at the field scale, here the sampling was done at the pixel scale, with a spatial match between the sensor acquisition and the field work. This sampling strategy has a clear positive influence on the predictive accuracy, since the estimated r2 for the Shannon index was no more than 0.17 in (Lopes et al., 2017a). The heterogeneity of grasslands, in particular semi-natural grasslands, can be high and thus the spatial measurement unit of botanical composition should match as much as possible the sensor's spatial resolution. These findings show that the remote sensing-based estimation of plant diversity in grasslands should be done at the pixel level and not at the field level.

SITS pre-processing

In very cloudy situations, SAR data, not affected by clouds, provide additional information that can help for prediction [START_REF] Clerici | Fusion of sentinel-1a and sentinel-2a data for land cover mapping: a case study in the lower magdalena region, colombia[END_REF]. Indeed, the combination of Sentinel-1 and Sentinel-2 SITS improved the prediction accuracy for some indices, e.g. Richness and Insect Dependence, but not for all. Sentinel-1 data can complement Sentinel-2 information for the learning algorithm, but SAR SITS alone are not sufficient. SAR data has been shown to provide relevant information for grassland management practices monitoring and cutting practices in particular [START_REF] Voormansik | Observations of cutting practices in agricultural grasslands using polarimetric sar[END_REF][START_REF] Hill | Relating radar backscatter to biophysical properties of temperate perennial grassland[END_REF][START_REF] Dusseux | Contribution of radar images for grassland management identification[END_REF]. Hence, although SAR data was obviously not expected to be sensitive to botanical composition, it could have been sensitive to biomass, providing information about variations in plant height that could be associated with variation in species. Yet, S1 data was not found correlated to plant height in this study. However, the aforementioned works were conducted on a reduced geographical area and the variability in SAR acquisitions, e.g. due to different orbit and incident angles, were negligible in comparison with our context. Another possible explanation could be the smoothing effect of the spatiotemporal filter that alters the spatial resolution of S1-SITS and thus reduces the sensibility of such filtered data. This suggests in future studies to investigate a more sophisticated pre-processing of Sentinel-1 for a joint use with Sentinel-2 data.

In this work, the full stack of temporal acquisition was used, and the irregular temporal acquisitions between various Sentinel tiles was handled using linear interpolation [START_REF] Inglada | Assessment of an operational system for crop type map production using high temporal and spatial resolution satellite optical imagery[END_REF]. This leads to a large amount of features, as displayed in Table 2, that can complexify the statistical learning step. In such situation, it is common to reduce the dimensionality, i.e., the number of features used. In this work, a conventional principal component analysis was used but it resulted in a loss of prediction accuracy. More advanced feature extraction have been investigated, such as partial least square or forward feature reduction, but they did not improve the prediction accuracy. As [START_REF] Rapinel | Evaluation of sentinel-2 time-series for mapping floodplain grassland plant communities[END_REF], we found that using the full set of Sentinel-2 bands (restricted here to spectral bands natively at 10m per pixel) works well in most of cases. Yet, using only red and near-infra red channels, with all the temporal features, provided very good results too. It may indicate that ignoring the spectro/backscatter or temporal structure of the Sentinel SITS when performing the dimension reduction is not appropriate and that further investigation is needed to extract relevant information from the data without loosing significant details.

Importance of the spatial auto-correlation

Spatial auto-correlation was scarcely taken into account in previous studies on grasslands while it can have a strong influence on the estimated prediction accuracy [START_REF] Hammond | Optimistic bias in classification accuracy assessment[END_REF]. Results obtained is this work are similar to those reported in [START_REF] Pohjankukka | Estimating the prediction performance of spatial models via spatial k-fold cross validation[END_REF]. This suggests that working on a large spatial area that allows the construction of a validation set that is spatially uncorrelated to the training set is an essential step to assess the quality of the results. Fortunately, with the increasing availability of high resolution SITS it will not be a critical issue in the future. However, for small spatial coverage data, such as hyperspectral or UAV ones, a careful statistical analysis should be done to prevent optimistic results.

Large scale prediction

The large scale map we produced highlights the grassland intra-parcel heterogeneity in terms of plant diversity, underlining the limits of predicting biodiversity indices at the field scale [START_REF] Lopes | Spectro-temporal heterogeneity measures from dense high spatial resolution satellite image time series: Application to grassland species diversity estimation[END_REF]. Such a map can be used to assess the relationships between environmental factors and plant diversity at fine scale.

Part of the prediction associated with the large confidence interval observed in the images can be explained by geographical factors. First, the edges of the parcels could correspond to mixed pixels and therefore exhibit atypical values. Second, individual trees and hedgerows are sometimes included in the field. Their spectral-temporal features do not correspond to inputs from which the system is learned and unconfident predictions are done. Looking at the confidence interval thus could inform on the reliability of the predicted values, at the pixel scale.

The decreasing average predicted Simpson values from the permanent, long rotation and temporary meadows observed on figure 9 are in accordance with the importance of long term no tillage and low intensity practices to keep grassland biodiversity at high level [START_REF] Klimek | Plant species richness and composition in managed grasslands: The relative importance of field management and environmental factors[END_REF]. However, the temporary grasslands present a distribution of predicted Simpson values similar to the permanent ones. This pattern may be explained by the fact that the number of reported permanent grasslands by farmers is lower than the actual number (they are declared as "temporary" instead). This is because permanent grasslands have tillage restriction, while temporary grasslands do not. As a consequence, a large number of declared temporary grasslands are, actually, permanent grasslands.

Perspectives

Perspectives for this work are two-fold. First, field work for different grassland areas (mountains, wetland ...) should be conducted to enrich the learning database and to assess the robustness of the prediction at the country scale. Second, different stratification strategies should be investigated, for instance according to their types or their eco-climatic conditions. Third, additional data set, such as DEM or surface temperature, should be used jointly for the prediction in order to capture more information about grasslands. This would require specific methodological developments to properly use such a heterogeneous data set.

Conclusions

The prediction of grasslands biodiversity indices using Sentinel-1&-2 times series was investigated in this work. We have shown that abundance-based indices based on the dominance of some species (Simpson) can be predicted over a large area with a coefficient of determination above 0.4 using Sentinel-2 SITS.

This study indicates that the temporal information contained in the SITS can compensate for the limited spectral information and spatial resolution compared to hyperspectral imagery. Another finding is that the spatial autocorrelation tends to bias the estimation of the prediction accuracy and thus, should be taken into account during validation.

The contribution of this work does not lie in selecting the best set of dates nor the best machine learning regression methods but in assessing the potential of Sentinel-1 and -2 dense times series as explanatory variables for the prediction of grasslands biodiversity indices. We believe that our results, using large scale data with various agricultural practices for different meteorological and topographic conditions, demonstrate the capacity of such data to monitor grasslands from an ecological viewpoint. In particular, intra-parcel variability was highlighted in this work and can be monitored over large areas.

Figure 1 :

 1 Figure 1: Study area. The region in the black frame is a zoom of the red area at the bottom of the image. The blue circles indicate the sampled plots in the grasslands. See section 2.2 for a full description. Data: 1:10m Natural Earth I (https://www.naturalearthdata.com); BD ORTHO R 50 cm IGN (background aerial photography)

10mFigure 2 :

 2 Figure 2: Field sampling: a plot (light red area) is made of 5 quadrats (red squares), 4 of them being located at 10 m from the center one. The light blue squares represent the corresponding Sentinel-2 pixel grid.

Figure 3 :

 3 Figure 3: Correlation matrix of the biodiversity variables using the Spearman's rank correlation coefficient.

Figure 4 :

 4 Figure 4: Reconstructed S2-NDVI (in black) and S1-VV polarization (in blue) with a 10 days step for different plots. For each plot, the 5 pixels corresponding to each quadrat are plotted. The red vertical lines indicate the dates of the detected clouds and shadows for the raw S2 pixels, before the gap-filling.

Figure 5 :

 5 Figure 5: Spatial cross validation processing flow. For clarity, only the data assignment between the folds and the training step for the estimation of r25 is shown. This process is performed for every biodiversity variable.

Figure 6 :

 6 Figure 6: Scatter plot of the predicted Simpson values obtained with RF and R-IR data (see table3). The black line is the identity line, i. e., when ŷ = y; r2 = 0.45. The colormap corresponds to the confidence interval σ of the prediction.

Figure 7 :

 7 Figure7: Influence of the spatial auto-correlation on the estimation of the prediction accuracy. "Spatial-CV " indicates that the spatial auto-correlation was taken into account during the fold construction of the cross validation, see section 3.2 for details.

Figure 8 :Figure 9 :

 89 Figure 8: Simpson index around plot 132_1699. (a) Very high spatial resolution aerial image from 2016, (b) Prediction of the Simpson value for pixels corresponding to grasslands in the national agricultural Land Parcel Information System and (c) Confidence interval value associated to the prediction.

Figure 10 :

 10 Figure 10: Proportion of free (in yellow) and masked (in gray) pixels for each temporal acquisition.

Table 1 :

 1 Plant cover percentage bins.

	Coverage estimation Assigned value
	(%)	(%)
	[1, 2.5]	1
	[2.5, 7.5]	5
	[7.5, 11.5]	10
	[11.5, 20]	15
	[21, 30]	25
	[31, 40]	35
	[41, 50]	45
	[51, 60]	55
	[61, 70]	65
	[71, 80]	75
	[81, 90]	85
	[91, 99]	95
	100	100

Table 2 :

 2 List of inputs used for the prediction. The last column provides the number of variables for each input.

	Acronym	Description

Table 3 :

 3 Estimated accuracy for the best predicted variables (r 2 > 0.3) and the corresponding standard deviation (σ r 2 ) of the 5-fold spatial cross validation.

	Data	Variable	Method	r2	σr 2
	R IR	Simpson	RF	0.45 0.13
	R IR	Shannon	RF	0.43 0.13
	S2	Color Diversity	RF	0.40 0.08
	S1 S2	Richness	KRidge	0.34 0.15
	S1 S2	Insect Dependence KRidge	0.32 0.15
	S1 S2	Flowers Richness	KRidge	0.32 0.16
	S1 S2	Reward Index	GP RBF	0.32 0.20

2017-09 2017-10 2017-11 2017-12 2018-01 2018-02 2018-03 2018-04 2018-05 2018-06 2018-07 2018-08 2018-09 2018-10 2018-11 2018-12 0.0

Acknowledgement

This work was supported by the SEBIOREF (National Institute for Agronomic Research and Occitanie Région) project which funded part of the field work. Catherine Bataillon, Sylvie Ladet, François Calatayud, Wilfried Heintz, Xavier Henry and Jérôme Willm are gratefully acknowledged for the field work in Spring 2018. We also acknowledge farmers for allowing full access to their grasslands. Mathieu Fauvel is thankful to Jordi Inglada, Vincent Thierion and Arthur Vincent for their support and help with iota 2 software. The authors would also like to thank the reviewers and the handling editor whose comments and suggestions improved this paper.

 3). The black line is the identity line ŷ = y; r2 = 0.40. The colormap corresponds to the confidence interval σ of the prediction.