Nora H Prior
email: nhprior@umd.edu

Marie S A Fernandez

Hédi A Soula

Clémentine Vignal

Water restriction influences intra-pair vocal behavior and the acoustic structure of vocalisations in the opportunistically breeding zebra finch (Taeniopygia guttata)

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

Temperate seasonally-breeding songbirds experience profound fluctuations in singing associated with the changing seasons [START_REF] Slagsvold | Bird song activity in relation to breeding cycle, spring weather, and environmental phenology[END_REF][START_REF] Catchpole | Bird Song: Biological Themes and Variations[END_REF]. More broadly, there is growing evidence from many species which suggests there are widespread changes to auditory-vocal communication systems throughout the year: including variation in the type and composition of songs and calls [START_REF] Marler | The voice of the chaffinch and its function as a language[END_REF][START_REF] Brenowitz | Plasticity of the adult avian song control system[END_REF][START_REF] Marler | Bird calls: their potential for behavioral neurobiology[END_REF][START_REF] Maddison | Soft song during aggressive interactions: seasonal changes and endocrine correlates in song sparrows[END_REF][START_REF] Lohr | The function of song types and song components in Grasshopper Sparrows (Ammodramus savannarum)[END_REF], acoustic structure of song syllables [START_REF] Maddison | Soft song during aggressive interactions: seasonal changes and endocrine correlates in song sparrows[END_REF], and auditory perception [START_REF] Lucas | Seasonal variation in avian auditory evoked responses to tones: a comparative analysis of Carolina chickadees, tufted titmice, and white-breasted nuthatches[END_REF][START_REF] Bass | Functional aspects of song learning in songbirds[END_REF][START_REF] Caras | Modulation of Peripheral and Central Auditory Processing by Estrogens in Birds[END_REF].

Whereas, a higher amount of courtship song during breeding seasons can easily be described as an example of breeding behavior coupled to the appropriate environmental context, the functional significance of these other changes in auditory-vocal communication systems remain largely unknown.

For opportunistically breeding species, it is unclear how environmental condition alone induces any changes in auditory-vocal communication systems, such as is seen in seasonallybreeding species. The Australian zebra finch is an opportunistically breeding species. Zebra finches time breeding bouts based on many environmental factors such as food and water availability, the presence of green grasses, humidity and day length (Perfito 2010). However, zebra finches within the same population, under good environmental conditions, are not necessarily all in the same breeding state (Perfito 2010). Furthermore, in the wild, zebra finch pairs do not always breed synchronously [START_REF] Zann | The zebra finch: a synthesis of field and laboratory studies[END_REF]). In the laboratory, experimental water restriction has reliably been used as a method for manipulating male and female breeding readiness [START_REF] Perfito | Tonic activation of brain GnRH immunoreactivity despite reduction of peripheral reproductive parameters in opportunistically breeding zebra finches[END_REF][START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF][START_REF] Prior | Neuroendocrine regulation of long-term pair maintenance in the monogamous zebra finch[END_REF]. The physiological effects of water restriction on males and females are somewhat varied; however, both sexes do experience similar shifts in neuroendocrine state as what has been described for seasonally breeding species (e.g. brain and circulating steroid levels) (Perfito 2010;[START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF][START_REF] Perfito | Seasonal-like variation in song control system volume of wild zebra finches[END_REF][START_REF] Prior | Neuroendocrine regulation of long-term pair maintenance in the monogamous zebra finch[END_REF].

While zebra finches will engage in courtship behaviors throughout the year, the fact that they coordinate breeding bouts at the level of the pair raises the question of whether shifts in vocal dynamics, that proceed breeding, could be used to coordinate breeding attempts. Indeed there is some evidence that vocal behavior may shift with environmental condition in zebra finches, as well as evidence that there are season-like shifts in neuroendocrine state centrally, which could mediate season-like behavioral variation (Prior andSoma 2015:Perfito et al., 2015).

For wild zebra finches, there is evidence that undirected song is more common during nonbreeding periods (Dunn and Zann 1996a;Dunn and Zann 1996b). Additionally, during periods of active nesting zebra finches produce a synchronized vocal exchange (duet) that is used to coordinate bi-parental care [START_REF] Elie | Vocal communication at the nest between mates in wild zebra finches: a private vocal duet?[END_REF][START_REF] Boucaud | Vocal negotiation over parental care? Partners adjust their time spent incubating based on their acoustic communication at the nest[END_REF][START_REF] Boucaud | Acoustic communication in zebra finches signals when mates will take turns with parental duties[END_REF]). More broadly, there is a shift in call types used during nesting [START_REF] Gill | Patterns of call communication between group-housed zebra finches change during the breeding cycle[END_REF], and some call types are used almost exclusively during nesting [START_REF] Elie | Vocal communication at the nest between mates in wild zebra finches: a private vocal duet?[END_REF][START_REF] Gill | Patterns of call communication between group-housed zebra finches change during the breeding cycle[END_REF]).

Here we tested the effect of water restriction, as an environmental manipulation that is known to impact breeding, on three aspects of zebra finch vocal behavior: (1) overall vocal activity, (2) patterns of vocal exchanges, and (3) the acoustic structure of vocalizations. We manipulated breeding readiness using established experimental water restriction paradigms. In order to compare the effect of water restriction to previous studies, we quantified circulating testosterone in males and circulating corticosterone levels in males and females. Despite that zebra finches are gregarious, they actually appear to spend the majority of their life with their pair bonded mate [START_REF] Zann | The zebra finch: a synthesis of field and laboratory studies[END_REF][START_REF] Mccowan | The size and composition of social groups in the wild zebra finch[END_REF]; thus, we chose to focus on describing the effects of water restriction on the vocal behavior of well-established zebra finch pairs. Additionally, in order to control for social condition during recordings, we recorded vocal behavior while pairs were visually isolated (Perez et al. 2015b). This forced all intra-pair communication to be through the vocal-acoustic domain and not visual or tactile.

Materials and Methods

Subjects

Twenty-eight paired adult domesticated zebra finches (>120 days old) were used for this study. All pairs were bonded for >1 year and had been allowed to breed together for several breeding cycles. Thus, only successfully-bonded pairs were used. Each pair was housed in their own cage (dimensions 50.8 × 40.5 × 30.5 cm). The colony was maintained on a 14:10 h light: dark cycle with temperatures between ~24-26 °C. All individuals had ad libitum access to seed, cuttlefish bone, and grit in addition to supplemental green vegetables once per week. Note that over the course of the study one individual died, fifteen pairs were recorded in the control condition (CON), while only fourteen of these pairs were also recorded in the water restricted condition (WR).

Experimental Design

Water availability was manipulated in order to influence reproductive condition. Pairs were slowly water restricted over the course of six weeks from 5mL down to 1mL per individual per week. This protocol is similar to what has previously been published [START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF][START_REF] Prior | Acute and chronic effects of an aromatase inhibitor on pair-maintenance behavior of water-restricted zebra finch pairs[END_REF]. Briefly, at the start of water restriction, two empty water towers were placed on each cage. Initially, 2mL of water was administered via pipette into each water tower five days a week. The amount of water administered was gradually decreased to .5mL twice a week in each water tower.

Vocal behavior of each pair was recorded three times during Water Restriction (WR) and three times during Water ad libitum (CON) resulting in six recordings/pair (each recording ~4 hours). From these recordings, we quantified the overall vocal activity (song and call rate), the quality of vocal exchanges (male song composition and patterns of intra-pair calling exchanges), and the acoustic structure of vocalisations (Figure 1). To control for potential confounds of habituation to the recording session, the order of water restriction was counterbalanced across pairs: nine pairs were recorded first under the CON and then under the WR condition (Cohort 1), and six pairs were recorded first under the WR and then under the CON condition (Cohort 2).

For both cohorts there was a break (37 days for cohort 1 and 16 days for cohort 2) between conditions. Importantly, this amount of time (>2 weeks) has previously been shown to allow recovery from water restriction [START_REF] Vleck | Reproduction in zebra finches: hormone levels and effect of dehydration[END_REF]. In order to assess the effect of water restriction on systemic hormone levels, we collected blood samples from the brachial vein and quantified circulating corticosterone and testosterone levels from each individual 0-3 days following the completion of recordings for each condition (CON and WR). Thus, there were at least 2 weeks between blood sample collection and the acoustic recordings.

Recording vocal behavior

To identify effects of water restriction on acoustic communication, we recorded intra-pair vocal behavior while the male and female of each pair were visually isolated. Since the only means of contact between partners during visual separation is acoustic, this paradigm allowed us to disentangle effects of our environmental manipulation on acoustic communication from more general pairing behaviors. Previous research has shown that intra-pair vocal dynamics are more stereotyped while visually isolated (Perez et al. 2015b). While visually-isolating pairs introduces a perturbation that could itself impact behavior, this isolation removes confounds that could result from effects of water restriction on other sensory systems involved in communication dynamics.

Each pair was moved from their home cage to a sound attenuation chamber (cage dimensions 63.5 × 33.0 × 30.5 cm) the day before recordings (~14:00 h). The male and female were physically separated in the testing chamber by placing each bird in individual but adjacent cages. Each cage contained two perches, a seed cup and a water tower, which was either filled or empty depending on water condition. Initially pairs remained in visual contact to facilitate habituation to the testing chamber. After ~ 3 h, an opaque partition was placed between the two cages (Figure 1). Behavior was not recorded until the following day which ensured behavior was not recorded immediately after the stress of handling and moving. In total pairs had at least 15 h to habituate to the testing set up.

The following morning, passive recordings of each pair were made using a digital recorder (SongMeters SM2, 16-bit, 44 kHz sampling rate; Wildlife Acoustics Inc., Concord, MA, U.S.A.) with omnidirectional microphones (SMX-II, Wildlife Acoustics Inc.; omnidirectional, flat frequency response 20 Hz-20 kHz, sensitivity 36 ± 4 dB, 0 dB 1∕4 1 V/pa at 1 kHz). One microphone was positioned above each partner's cage. More specifically, the microphone was positioned downwards above the cages in order to minimize amplitude variations due to the position of the bird. Each recording was ~4 h long from ~08:00-12:00, the period starting at lights on. There were six recordings in total (Recording number (Nb), 1-6), with ~12 h of recordings/pair/water condition.

Scoring Vocal Behavior

All vocalisations were extracted from recordings using in-house software, as previously described [START_REF] Elie | Dynamics of communal vocalizations in a social songbird, the zebra finch (Taeniopygia guttata)[END_REF]Perez et al. 2015b) (see supplemental information).

Songs

Songs were identified manually (see supplemental information). In total we recorded 205 songs from 13 males during CON and 28 songs from 8 males during WR. Characteristics of each song were manually scored (N.H.P.) using Praat (V 5.3.56, 2014) [START_REF] Boersma | PRAAT: a system for doing phonetics by computer[END_REF].

Song bouts were defined as a single rendition of song, separated by at least 1 sec. Because there were drastically fewer WR songs and most males did not sing on each recording date we analyzed song rate by calculating a single song rate per male for each water condition (CON and WR).

In order to determine the effect of water restriction on the composition and acoustic structure of song, we further processed song in two ways. First, global song elements were extracted manually, as concise sequences, in Praat. More specifically, the following sequences were extracted: an introductory note sequence as well as every full and partial motif (Figure 2A).

Any additional inserted elements (e.g. distance calls) were not extracted for analysis. Second, the following four individual notes were extracted: the first and last introductory notes, as well as two song syllables: one from the first half of the motif (early syllable) and the second from the second half of the motif (late syllable) (Figure 2B). For a given male, the same two song syllables were extracted from each song bout. We preferentially selected harmonic, 'stacked' syllables to facilitate consistent labeling of the note boundaries. We used the extracted global song elements to generate four characteristics of song: (1) mean duration of full motifs, (2) the ratio of the duration of full motifs to the duration of the total song bouts (the duty cycle), (3) the mean inter-sequence interval (between the introductory note sequence of first motif and subsequent motifs) and (4) the ratio of the duration of the introductory note sequence to the total duration of motifs in the corresponding song bout.

Additionally, analyses of song composition were performed on the full and partial motifs (see supplemental information). The extracted introductory notes and song syllables were used for acoustic structure analysis (see supplemental information).

Calls

In addition to song, we extracted all the calls during a recording. Our initial observations suggested that pairs engaged in call and response periodically throughout the recordings and that during quieter periods there was typically only one individual calling. Thus, in order to focus on intra-pair calling, we defined calling 'bursts', where calling activity was elevated (overall call rate per pair was >10% above the average for the recording) (see supplemental information). By isolating bursts specifically for analysis, we increased the likelihood that we were examining parameters of true vocal exchanges. We quantified burst rate (# bursts/ hour) and call rate per individual (for the full recording and within bursts specifically), as well as two parameters of calling exchanges: (1) temporal coordination (cross-correlation), and (2) the predictability in the pattern of turn-taking (% Markov fit), independent of the timing (see supplemental information).

In highly coordinated pairs, we expected that the delay of the male's (M) response to his female (F) partner would be the same as that F's response to the M (comparison of M-F and F-M calling delay). Furthermore, we expected higher consistency in the delay times between calls (delay times of the M answer to the F, and the F answer to the M). These would be examples of increased temporal coordination.

Previous work has shown visually isolated pairs are more Markovian in their alternation of calls (Perez et al. 2015b). More Markovian patterns mean that the identity of a given caller depends on the identity of the previous caller. In this way, the pattern of turn taking is assessed independently of the calling delays. This phenomenon may be analogous to how speech patterns in humans change between talkers when they are face-to-face in contrast to when they are on the phone. More Markovian exchanges are likely a result of the birds' increased attention to the acoustic channel (Perez et al. 2015b). Together, the cross-correlation and % Markov fit provide a description of the pattern of intra-pair vocal dynamics.

Additionally, we randomly selected a subset of 50 calls from each recording for each individual to use for the analysis of acoustic structure (function: "random", a uniform random selection in python). These calls were then manually examined, and only cleanly recorded calls were kept for further analysis (calls that did not include noises from wing/body movements, the cage, or other vocalizations from the partner). On average there were 121 ± 32 (Mean ± SD) calls per individual per condition. Final sample sizes are: N=1991, Female CON; N=1608, Female WR; N=1683, Male CON; and N=1670, Male WR. The majority of calls were 'stack' calls [START_REF] Maat | Zebra finch mates use their forebrain song system in unlearned call communication[END_REF][START_REF] Gill | Patterns of call communication between group-housed zebra finches change during the breeding cycle[END_REF]) (see supplemental information).

This subset of extracted calls was used to quantify the effect of water condition on the acoustic structure of calls. We used two complementary methods of analysis: (1) full-spectrum analysis (Mouterde et al. 2014a;[START_REF] Elie | The vocal repertoire of the domesticated zebra finch: a datadriven approach to decipher the information-bearing acoustic features of communication signals[END_REF], and (2) quantification of pre-defined parameters of the spectral envelope (Seewave package [START_REF] Sueur | Seewave, a free modular tool for sound analysis and synthesis[END_REF])) implemented in R software (v. 3.2.3, 2015, R Foundation for Statistical computing).

Analysis of acoustic structure

We used a standard two-tier approach (full spectrum analysis and pre-defined acoustic features) to describe the effects of water restriction on the acoustic structure of vocalisations. The two approaches we used answer slightly different questions. First, we asked whether water condition affected the acoustic structure of vocalisations in any way, i.e. "Is there any difference?". In order to comprehensively analyze the spectrogram, we used an in-house program (written in R) to summarize all the information in the entire spectrogram by using all the data points of the spectrograms as raw measures in a PCA (Mouterde et al. 2014b;Elie and Theunissen 2016) (see supplemental information). The full spectrum analysis allows for an unbiased assessment of how WR affected the spectrogram (there is no a priori selection of specific acoustic elements); however, the principal components cannot be easily related back to specific acoustic features. Thus, we also asked "are there acoustic differences in traditionally defined parameters?". For this analysis we used Seewave package in R [START_REF] Sueur | Seewave, a free modular tool for sound analysis and synthesis[END_REF] to quantify the effect of water condition on four pre-defined acoustic features: an estimate of call loudness (root-mean-square RMS of the temporal envelope) and three parameters of the spectral envelope (median freq (kHz), stDev (kHz), and IQR (kHz)) (see supplemental information).

Circulating testosterone and corticosterone levels

For each individual, blood samples were collected during both treatments after acoustic recordings were taken. Note that this means there was a minimum of 2 weeks between blood sample collection and a subsequent recording. Both partners of a pair were bled at the same time by two experimenters (M.A.S.F and N.H.P). All blood samples were collected between 11:00-13:00 h and within 3 min of opening the cage door to catch individuals (Mean ± SEM: 55 ± 6 sec). Blood samples were centrifuged to obtain plasma (10 min at 10,000 x g), which was stored at -20°C until further processing. Circulating concentrations of steroids were assayed using enzyme immunoassay kits (Cayman Chemical Co., Ann Arbor, MI: Testosterone, 582701; Corticosterone, 500655) (see supplemental information).

Statistics

All statistical analyses were carried out in R (v. 3.2.3, 2015, R Foundation for Statistical computing). We used generalized linear-mixed models (LMMs; function lmer from the lme4 Package). The design of our study included three potential factors, Water Condition (WR and CON), Recording number (1-6: indicating the order recordings were made in), and Cohort (Cohort 1, CON first and Cohort 2, WR first). We used Recording Nb as a continuous variable in our models. We did not include Cohort as a factor. For pair-level dependent variables, Pair was included as a random factor. For individual-level dependent variables, Individual was included as a random factor, and Sex was included as a between-subjects factor. Prior to interpretation, we checked the validity of each model by plotting the distribution of the residuals. Data were transformed as necessary. Model summaries are provided in Supplementary Table 1-7.

Results

WR effects on circulating corticosterone and testosterone

Consistent with prior research, water restriction had no effect on circulating corticosterone levels in males or females (Males, mean ± SEM, CON: 1.64 ± 0.43 ng/mL; WR:

1.17 ± 0.14 ng/mL; χ 2 (1) <0.01, P=0.987; Females, mean ± SEM, CON: 1.34 ± 0.23 ng/mL; WR: 1.43 ± 0.35 ng/mL; χ 2 (1) =1.48, P=0.225). Additionally, circulating testosterone levels were low for males regardless of water condition and there was no main effect of water restriction on circulating testosterone (mean ± SEM, CON: 0.41 ± 0.10 ng/mL; WR: 0.47 ± 0.15 ng/mL; χ 2 (1) =0.40, P=0.527). The low circulating level of testosterone for males in both treatment groups suggests that the effect of water restriction was smaller in the current study than what has been previously reported [START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF].

WR decreases overall vocal activity

Overall, in response to water restriction, zebra finches vocalized less. This effect was clearest on male song rate during water restriction (χ 2 (1) =6.70, P=0.009) (Figure 3A). In addition, there was a small, but statistically significant decrease in call rate (calculated from the entire recordings) for both males and females in response to water restriction (Females, mean ± SEM, CON: 3.03 ± 0.43 calls/min; WR: 2.67 ± 0.29 calls/min. Males, mean ± SEM, CON: 3.56 ± 0.57 calls/min; WR: 2.92 ± 0.34 calls/min) (Water Condition χ 2 (1)=6.66, P=0.010. Sex χ 2 (1) =0.87, P=0.350. Recording nb χ 2 (1) =0.09, P=0.763. Water Condition × Sex χ 2 (1) =0.24, P=0.877).

As we were primarily interested in intra-pair calling behavior, we also quantified the effect of water restriction on calling during 'bursts', periods of elevated calling activity. There was no main effect of water restriction on burst rate (number of bursts/ hour of recording) (Water Condition χ 2 (1) =0.65, P=0.420) (Figure 3B); however, as with call rate during the entire recording, water restriction caused a decrease in male and female call rate within bursts (Water Condition χ 2 (1) =4.53, P=0.033; Sex χ 2 (1) =0.01, P=0.940; Water Condition × Sex χ 2 (1) =0.02, P=0.881) (Figure 3C). Throughout the course of the experiment, the burst rate increased from an average of 4.14 ± 0.31 bursts/ h in recording 1 to 6.01 ± 0.47 bursts/ h in recording 6 (mean ± SEM) (Recording nb χ 2 (1) =19.29, P<0.001); however, call rate within bursts did not change over the course of the experiment (Recording nb χ 2 (1) =0.04, P=0.850).

WR influenced song composition

Water restriction did not change the duration of male full motifs (Water Condition, χ 2 (1) =2.88, P=0.090) (Figure 4A), but significantly affected other global song elements. More specifically, when water restricted, males produced songs with a decreased proportion of full motifs (duration of full motifs/total song bout (s)) (Water Condition, χ 2 (1) =4.46, P=0.037) (Figure 4B), an increased duration of inter-element intervals (Water Condition, χ 2 (1) =5.81, P=0.016) (Figure 4C), and an increased proportion of introductory notes (duration of introductory notes/ full motifs (s)) (Water Condition, χ 2 (1) =7.18, P=0.007) (Figure 4D).

WR does not affect F-M calling dynamics

In contrast to male song, there was no significant effect of water restriction on the patterns of intra-pair calling (temporal coordination or turn taking). First, there was no effect of water restriction on our measure of temporal coordination (cross correlation of F-M and M-F calling delay) (Figure 5A). The cross-correlation demonstrates that there was consistent temporal delay in M-F calling which was unaffected by water condition. Calls were answered in ~200ms (regardless of whether the male or female was answering). Additionally, there was no effect of water restriction on our estimate of turn-taking within a pair (% Markov Fit) (Water Condition: χ 2 (1) =3.03, P=0.081; Recording Nb: χ 2 (1) =0.07, P=0.789, Figure 5B).

WR has sex specific and vocalisation specific effects on acoustic structure

WR has effects on acoustic parameters of song elements

In order to describe the effect of Water Condition on the acoustic structure of songs, we used two complementary analyses. From the full spectrum analysis, we found no effect of Water Condition on song elements, but some effects on individual notes and song syllables (Table 1).

More specifically, for the last introductory note PC1 and PC2 were significantly affected by water restriction. Additionally, for the early song syllable PC3 was affected by Water Condition.

Consistent with the full spectrum analysis, we saw effects of Water Condition on predefined acoustic features. Again, there was no effect of Water Condition on full motifs (Table 1;Table 2). There was a significant effect of Water Condition on amplitude and interquartile-range (IQR) for partial motifs (IQR was lower during WR and RMS was higher: Table 1;Table 2). Similar to our full spectrum analysis, we saw significant effects of Water Condition on the last introductory note and early song syllable (Table 1; Table 2). More specifically, the last introductory note was at a lower amplitude when produced during water restriction (Table 1; Table 2). Additionally, the early song syllable had an increased IQR (Table 1; Table 2).

Sex specific effects of water restriction on the acoustic parameters of calls

As with song, we used two complementary methods to examine the effect of Water Condition on the structure of calls (predominately short, stack-like calls). First, we conducted the full spectrum analysis on all of the calls. There was a main effect of Water Condition on PC3 (Table 3). However, there was also a highly significant Water Condition × Sex for PC1, PC2, and PC3 (Table 3). Therefore, we also ran the full spectrum analysis on female and male calls separately. From these results, we see a clear effect of Water Condition on female calls (PC1, Table 3) and male calls (PC1 and PC2, Table 3).

To complement the full spectrum analysis, we also examined the effect of Water Condition on specific pre-defined acoustic features. For each of our acoustic parameters (RMS, median frequency, SD, and IQR), we found highly significant Water Condition × Sex interactions (Table 3 andTable 4). WR had different effects on males and females (Table 3). For female calls, WR decreased amplitude of calls (Table 3,Table 4), but increased Median Frequency (Table 3,Table 4), StDev (Table 3,Table 4), and IQR (Table 3,Table 4). In contrast, WR increased the amplitude of male calls (Table 3,Table 4), while decreasing the Median Frequency (Table 3, Table 4).

Discussion

There are widespread changes in both vocal production and auditory perception in seasonally breeding species, which raises the question of whether auditory-vocal communication in opportunistically breeding species is also affected by reproductive states and/or environmental conditions [START_REF] Bass | Functional aspects of song learning in songbirds[END_REF]Alward et al. 2017b). Our results suggest that zebra finch vocal production is affected by environmental condition. More specifically, we describe several effects of water restriction on male and female vocal behavior: water restriction (1) decreased the vocal activity of both males and females (calling and male song), (2) shifted the composition of male song so that it more closely resembled undirected song, and (3) influenced the acoustic structure of all vocalizations in sex-and vocalisation-specific ways which we discuss below. Based on previous research, we would anticipate that these changes are discriminable and could be salient to partners and other conspecifics. The fact that these effects are vocalization and sex specific raises interesting questions about how WR is impacting the acoustic structure of vocalizations. In addition to neuroendocrine changes associated with decreased breeding readiness, there may be context-specific changes in socio-sexual motivation in response to a lower environmental state.

Assuming that these are acoustic indicators of breeding readiness, this raises the question of whether this information is used by the partner and other conspecifics. In order to determine if these effects are specific to the pair-directed behavior, it would be useful to identify the effects of water restriction on vocal dynamics for other social relationships and social conditions as well.

Effects of water restriction on vocal activity

Seasonally breeding bird species experience profound neural and vocal plasticity related to breeding season [START_REF] Bass | Functional aspects of song learning in songbirds[END_REF]. In our current study, the decrease in vocal output is similar to the decrease in song output seen in temperate-breeding songbirds. This is not altogether surprising considering zebra finches have similar shifts in neuroendocrine state associated with breeding as seasonally-breeding species. For zebra finches, changes in reproductive and or breeding state are associated with changes in circulating steroid levels [START_REF] Perfito | Opportunism at work: habitat predictability affects reproductive readiness in free-living zebra finches[END_REF][START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF][START_REF] Gill | Patterns of call communication between group-housed zebra finches change during the breeding cycle[END_REF][START_REF] Prior | Sex steroid profiles in zebra finches: Effects of reproductive state and domestication[END_REF], brain steroid levels [START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF], and nuclei volume for brain regions associated with song production (HVC and area X) [START_REF] Perfito | Seasonal-like variation in song control system volume of wild zebra finches[END_REF]. While we did not see changes in circulating testosterone levels in males here, we did not measure other circulating or brain steroids. Furthermore, in zebra finches, as in seasonally-breeding birds, sex steroids can have significant effects on male song [START_REF] Arnold | The effects of castration and androgen replacement on song, courtship, and aggression in zebra finches (Poephila guttata)[END_REF][START_REF] Prior | Context-dependent effects of testosterone treatment to males on pair maintenance behaviour in zebra finches[END_REF]. In temperate, seasonally-breeding species, the decrease in song output during non-breeding periods is generally assumed to be the result of changes in behavioral context and neuroendocrine condition, which together result in the absence of breeding-condition song signals. Similarly, we interpret the effect of water restriction on song output in male zebra finches to reflect the changes in neuroendocrine state induced by water restriction.

The relationship between calling behavior and singing is not well understood, and thus it is particularly interesting that we find some evidence that water restriction decreases calling activity in addition to song. This effect could be interpreted as a decrease in (1) general vocal output, (2) attention or focus on the auditory channel, or (3) sociosexual motivation. There is evidence in zebra finches that the same brain nuclei in the song control system that regulate song, also regulate aspects of calling behavior [START_REF] Maat | Zebra finch mates use their forebrain song system in unlearned call communication[END_REF][START_REF] Benichov | The forebrain song system mediates predictive call timing in female and male zebra finches[END_REF]).

Therefore, it is possible that the effect of water restriction on specific brain regions could influence call and song output similarly. Another possibility is that the motivation to call and sing are related. The fact that water restriction did not impact the quality of calling exchanges between zebra finch mates suggests that water restriction did not have an impact on the quality of the pair bond and/or motivation to stay bonded. This is consistent with what is known about pair bonding in zebra finches [START_REF] Zann | The zebra finch: a synthesis of field and laboratory studies[END_REF][START_REF] Prior | Neuroendocrine regulation of long-term pair maintenance in the monogamous zebra finch[END_REF]. In this study, we chose to examine behavior under a relatively controlled social context (visual separation) and to focus on pair bonded birds. Therefore, it is unclear to what extent this effect of water restriction would extend to vocal behavior in other social contexts. In seasonally-breeding song sparrows, for example, there is a decrease in song output in the winter, yet males sing at the same rate during simulated territorial intrusions [START_REF] Maddison | Soft song during aggressive interactions: seasonal changes and endocrine correlates in song sparrows[END_REF]). In zebra finches, which are non-territorial and form life-long pair bonds, it would be particularly interesting to determine whether calling activity is similarly affected by water restriction within other social relationships (e.g. familiar and novel birds) and social contexts (e.g. traveling, foraging, etc).

Effect of water restriction on acoustic structure

Several lines of more recent research have described plasticity in the acoustic structure of calls depending on ontogeny [START_REF] Villain | Parental influence on begging call structure in zebra finches (Taeniopygia guttata): evidence of early vocal plasticity[END_REF], social context [START_REF] Hile | Call convergence within groups of female budgerigars (Melopsittacus undulatus)[END_REF], and internal physiological state [START_REF] Perez | The acoustic expression of stress in a songbird: does corticosterone drive isolation-induced modifications of zebra finch calls?[END_REF]Perez et al. 2015a;[START_REF] Perez | Corticosterone triggers high-pitched nestlings' begging calls and affects parental behavior in the wild zebra finch. Behavioral Ecology:arw069 Perfito N[END_REF]. Based on our current knowledge of zebra finch perception, we have every reason to expect that zebra finches can hear very small differences in the acoustic structure of their harmonic calls. Amongst birds, zebra finches are particularly sensitive to small spectro-temporal changes in harmonic complexes [START_REF] Lohr | Detection of changes in timbre and harmonicity in complex sounds by zebra finches (Taeniopygia guttata) and budgerigars (Melopsittacus undulatus)[END_REF][START_REF] Dooling | Auditory temporal resolution in birds: discrimination of harmonic complexes[END_REF][START_REF] Dooling | Auditory temporal resolution in the Zebra Finch (Taeniopygia guttata): A model of enhanced temporal acuity[END_REF][START_REF] Lohr | The discrimination of temporal fine structure in calllike harmonic sounds by birds[END_REF][START_REF] Prior | Acoustic fine structure may encode biologically relevant information for zebra finches[END_REF], therefore, while the effects of water restriction that we report here are small (e.g. 100 hz in median frequency), it is reasonable to assume these small differences would be perceptible to zebra finches.

One of the simplest explanations for seasonal plasticity in acoustic structure is that this plasticity is due to changes in syrinx size mediated by circulating testosterone levels (Suthers et al. 1999). However, the differences we report are unlikely to be related to shifts in syrinx size. First, we did not see changes in circulating testosterone levels in the periphery. Furthermore, we would expect that changes in syrinx volume would result in a change in fundamental frequency across the entire song [START_REF] Cynx | Testosterone implants alter the frequency range of zebra finch songs[END_REF][START_REF] Alward | Pleiotropic control by testosterone of a learned vocal behavior and its underlying neuroplasticity[END_REF]. Alternatively, changes in brain sex steroid levels alone within key nuclei of the song control system can impact spectral features of song syllables (Alward et al. 2017a;Alward et al. 2017b). It is likely that zebra finches experience significant changes in brain steroid levels independent of peripheral changes [START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF], and it is possible, although it has not been tested, that there are also effects of WR on other neuromodulatory systems as well as sex steroid receptors which could explain these results.

Here we only looked at stack-like calls and male song, thus it is unclear if these changes in acoustic structure extend across call types. Regardless of the mechanism, the changes in acoustic structure could serve as a functional cue of reproductive condition.

Functional significance of flexibility of calls

It may seem surprising that the effect of water restriction on calls differed between males and females. The results of our full spectrum analysis suggest that WR has a significant effect of water restriction on both male and female calls. However, for each of our pre-defined acoustic features, the effect of water restriction was in the opposite direction for males and females.

Female water restricted calls were lower amplitude, higher frequency (median frequency) and more variable (standard deviation and inter-quartile range). In contrast, male water restricted calls were higher amplitude and had a lower frequency (median frequency). From an ecological perspective, it makes sense that females may be more sensitive to environmental perturbations (Perfito 2010;[START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF], and it is possible that female reproductive condition may be a bottleneck on the timing of breeding bouts. Furthermore, water restriction has a more profound effect in females, including significantly reducing the size of the female ovary and oviduct. In contrast, water restriction has only a small effect on male testis size [START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF][START_REF] Perfito | Tonic activation of brain GnRH immunoreactivity despite reduction of peripheral reproductive parameters in opportunistically breeding zebra finches[END_REF]. Although in the current study, the effects of water restriction on male song suggest that there was a shift in the central neuroendocrine state of male zebra finches. Thus, one interpretation of the sex differences in the effect of water restriction on calls is that water restriction had differing effects on the neuroendocrine state of males and females. Further research describing the effect of water restriction and reproductive readiness on vocal dynamics across social conditions and contexts would help clarify what functional relevance this difference may have.

Importantly, the stack-like call that we investigated here is not specific to nesting periods [START_REF] Gill | Patterns of call communication between group-housed zebra finches change during the breeding cycle[END_REF]. Thus, regardless of the mechanisms resulting in changes in acoustic structure, we would expect that environmental and/or reproductive condition dependent flexibility in stacklike calls would precede other early breeding behaviors such as identifying a nesting site and building a nest. This raises the hypothesis that such vocal flexibility could play a functional role in coordinating partners' breeding attempts.

Functional significance of the effect of water restriction on male song

This is not the first report suggesting that breeding condition affects male song in zebra finches. Previous paradigms using water restriction have reported similar decreases. In [START_REF] Prior | Effects of water restriction on reproductive physiology and affiliative behavior in an opportunistically-breeding and monogamous songbird, the zebra finch[END_REF] we did not find a statistically significant effect of water restriction on male song during a partner separation and reunion paradigm; however, qualitatively males did spend less time singing. Additionally, [START_REF] Rashotte | Influence of food and water availability on undirected singing and energetic status in adult male zebra finches (Taeniopygia guttata)[END_REF] also reported a significant effect of brief (2 day) water deprivation on song.

Beyond just a decrease in overall song output, we might also expect a change in the type of song and/or performance. Zebra finches produce two types of song, directed and undirected, which are composed of the same stereotyped motif but differ in more subtle aspects such as the number of introductory notes and overall stereotypy of the performance [START_REF] Sossinka | Song types in the zebra finch Poephila guttata castanotis[END_REF]. Performances of directed song, or courtship song, involve the male facing the female and often involve other aspects of courtship dance [START_REF] Zann | The zebra finch: a synthesis of field and laboratory studies[END_REF]. In contrast, undirected song can be performed anywhere, and the male appears to be singing to himself [START_REF] Zann | The zebra finch: a synthesis of field and laboratory studies[END_REF]. In our current study, these pairs were always visually isolated. Therefore, all song produced under this behavioral paradigm would technically be classified as undirected song. However, taken together, the effects of water restriction on song that we report here are broadly consistent with a shift from directed to undirected zebra finch song. More specifically, water restriction increased inter-sequence interval and the proportion of introductory notes as would be expected if there were a shift from directed to undirected song [START_REF] Woolley | Social context-induced song variation affects female behavior and gene expression[END_REF][START_REF] Rajan | Behavioral and neural signatures of readiness to initiate a learned motor sequence[END_REF]. Additionally, our full spectrum analysis revealed that water restriction had the most profound effect on the last introductory note and the early song syllable, again similar to what would be expected from undirected song [START_REF] Rajan | Behavioral and neural signatures of readiness to initiate a learned motor sequence[END_REF].

A key difference between undirected and directed song is the motivational state of the male. Likewise, we would expect a similar shift in motivational state (away from courtship/breeding behaviors) following water restriction. Together these results raise the hypothesis that the changes in song are related to shifts in sexual motivation following a reduced state of breeding readiness.

Conclusion

Our current results contribute to several lines of evidence suggesting there are widespread seasonal changes in audio-vocal communication [START_REF] Catchpole | Bird Song: Biological Themes and Variations[END_REF][START_REF] Bass | Functional aspects of song learning in songbirds[END_REF].

Our work supports the notion that this vocal plasticity extends to non-seasonally-breeding species as well as to calling behavior. Furthermore, as our data suggests there is vocal plasticity in response to environmental conditions, this raises the question of whether information on breeding-state is present within zebra finch vocalisations and how zebra finches may use this information to coordinate breeding attempts.

Table 1: Effect of water restriction on the first three components of the full spectrum analysis and pre-defined acoustic parameters. Note: To assess the effect of water restriction on the acoustic structure of male song we used mixed linear models (e.g. lmer(acoustic parameter ~ Water Condition + (1|Ind)). Because so few males sang, we only included water condition (not Recording nb) in these models. Bolded values indicate statistical significance (P<0.05). DF=1 for water condition.

Full Spectrum

Pre-Defined acoustic parameters summarizes the parameters of vocal behavior we used as dependent variables and the analyses conducted. Briefly, timestamped calls and male song were extracted, and we used these vocalisations to quantify three characteristics of vocal activity (male song rate, burst rate, and call rate). Additionally, we quantified the effect of water restriction on the composition of male song (e.g. global song elements such as introductory note sequence, full motifs and partial motifs) and on the quality of the partners calling exchanges (cross-correlation of male-female calls, and % Markov Fit of patterns of turn taking). Finally, we used two methods (full spectrum analysis and pre-defined acoustic features) to assess the effects of water restriction on the acoustic structure of vocalisations.

Figure 1 :

 1 Figure 1: Schematic of our general experimental paradigm. The vocal behavior of established

Figure 2 :

 2 Figure 2: Spectrogram of a male song bout. Examples of extracted song elements are labeled:

Figure 3 :

 3 Figure 3: Effect of water restriction on vocal activity: (A) song rate (per hour). Note that only

Figure 4 :

 4 Figure 4: Effect of water restriction on song composition (A) duration of full motifs (s) (B) ratio

Figure 5 :

 5 Figure 5: Effect of water restriction on coordination of calling exchanges. (A) Coordination

but very small effect of WR on partial motifs

			Axis 1		Axis 2		Axis 3		RMS		Med (kHz)	SD (kHz)	IQR (kHZ)
			ChiSq	P	ChiSq	P	ChiSq	P	ChiSq	P	ChiSq	P	ChiSq P	ChiSq	P
	Global	Full Motif 0.02	0.891	0.25	0.615	0.72	0.397	2.90	0.088	3.03	0.082	0.15	0.696	0.33	0.565
	Song
	Elements	Partial	0.07	0.791	0.10	0.757	0.05	0.827	3.85	0.049	0.16	0.689	0.07	0.798	4.31	0.038
		Motif
			Q1: Is there evidence of any difference? NO		Q2: Are there acoustic differences in traditionally defined parameters? Yes,
	statistically significant, Individual First IN 0.07 0.798 1.26 0.261 0.12 0.726 0.48 0.487 0.26	0.609	0.22	0.641	3.47	0.063
	Notes	Last IN	11.50	<0.001 2.46	0.117	6.67	0.010	5.68	0.017	0.51	0.477	1.71	0.191	0.36	0.546
		Early	2.36	0.125	1.71	0.191	4.87	0.027	3.08	0.080	1.11	0.293	2.14	0.143	3.85	0.049
		Song
		Syllable
		Late Song	1.26	0.261	2.17	0.141	0.27	0.605	0.05	0.831	2.50	0.114	2.78	0.096	0.79	0.373
		Syllable
			Q1: Is there evidence of any difference? YES in the last	Q2: Are there acoustic differences in traditionally defined parameters? YES
			Introductory note and early song syllable

in the last Introductory note and early song syllable

Table 2 :

 2 Effect of water restriction (Mean ± SEM) on pre-defined acoustic parameters of song.

	Note: Bolded values indicate statistical significance (P<0.05).

Table 3 .

 3 Effect of water restriction, sex, and recording number on the first three components of the full spectrum analysis and pre-

	defined acoustic parameters.

Q2: Are there acoustic differences in traditionally defined parameters? YES (

for all four parameters there is a Water Condition X Sex interaction) Female

		Water Cond	14.68	<0.001	0.92	0.338	0.03	0.854	52.42	<0.001	121.90	<0.001	12.10	<0.001	44.61	<0.001
	Calls
		Record. Nb	55.97	<0.001	1.67	0.196	0.007	0.935	21.80	<0.001	46.80	<0.001	7.73	<0.001	1.78	0.180
			Q1: Is there evidence of any acoustic difference? YES		Q2: Are there acoustic differences in traditionally defined parameters? YES (all four
									parameters)
	Male	Water Cond	4.91	0.027	4.98	0.026	1.70	0.195	30.36	<0.001	24.74	<0.001	3.15	0.080	1.41	0.236
	Calls	Record. Nb	11.16	<0.001	7.64	0.006	7.23	0.007	0.80	0.371	1.07	0.301	0.37	0.542	2.11	0.150
			Q1: Is there evidence of any acoustic difference? YES		Q2: Are there acoustic differences in traditionally defined parameters? YES (RMS and
									median freq)

Table 4 .

 4 Effect of water restriction (Mean ± SEM) on pre-defined frequency parameters of male and female calls.

	678
	679
			RMS		Med (kHz)		SD (kHz)		IQR (kHZ)	Not
			CON	WR	CON	WR	CON	WR	CON	WR	e:
		Female Calls	10.31 ± 0.05	9.91 ± 0.06	3.03 ± 0.01	3.15 ± 0.01 1.82 ± 0.005 1.83 ± 0.005 2.56 ± 0.02 2.69 ± 0.02	Bol
		Male Calls	10.01 ± 0.05	10.44 ± 0.06 3.12 ± 0.03	3.02 ± 0.01 1.86 ± 0.005 1.85 ± 0.005 2.73 ± 0.02 2.70 ± 0.02	ded
										valu
	685	es indicate statistical difference (P<0.05) between CON and WR conditions.
	686
	687

Acknowledgements

We would like to thank Colette Bouchut, Nicolas Boyer, Amanda Cocquelet, Karen Tronchere, Dr. Ingrid Boucaud, Dr. Emilie Perez, Dr. Avelyne Villain, and the entire ENES lab for help with animal husbandry, data collection, and discussion of statistical analyses. For feedback on the manuscript, we thank Dr. Matthew D. Taves, Dr. Benjamin A. Sandkam, Adam Fishbein, and Dr. Robert J. Dooling. This work was conducted in accordance with the French and European animal care legislation under authorization number 42-218-0901-38 SV 09 (ENES Lab, Direction Départementale des Services Vétérinaires de la Loire).

Funding

Funding including: A French Agence Nationale de la Recherche grant (project 'Acoustic Partnership') to C.V., an IUF grant (Institut Universitaire de France) to C.V., a joint NSF/ANR e CRCNS grant 'AuComSi' for M.S.A.F. and H.A.S., and a post-doctoral fellowship from Université Jean Monnet to N.H.P.

Conflict of Interest.

The authors have no conflict of interest to declare.

Ethical Approval

Informed Consent

Not applicable.

Condition + scale(Recording nb) + (1|Ind)). Bolded values indicate statistical significance (P<0.05). DF=1 for water condition, sex 676 and recording number. 677


Christian Lécot
email: christian.lecot@univ-smb.fr

Pierre L'ecuyer
email: lecuyer@iro.umontreal.ca

El Rami

Ali Haddad
email: rami.haddad@usj.edu.lb

Tarhini

Haddad

Ali Tarhini
email: ahtarhini@yahoo.fr

Quasi-Monte Carlo simulation of coagulation-fragmentation

Keywords: Quasi-Monte Carlo method, Coagulation equation, Fragmentation equation, Stochastic particle method, Low-discrepancy sequence 2000 MSC: 11K45, 65C05, 82C80

Introduction

Coagulation-fragmentation models have applications in many domains of science and technology: aerosol dynamics, polymerization, and combustion processes, for example [START_REF] Aldous | Deterministic and stochastic models for coalescence (aggregation and coagulation): a review of the mean-field theory for probabilists[END_REF]. M. von Smoluchowski modelled the coagulation of particles in a medium where only binary collisions occur: the numbers of particles of different sizes obey an infinite system of ordinary differential equations [START_REF] Smoluchowski | Versuch einer mathematischen Theorie der Koagulationskinetik kolloider Lösungen[END_REF]. H. Müller rewrote the equations in terms of an integro-differential equation for the time evolution of the particle size density function [START_REF] Müller | Zur allgemeinen Theorie der raschen Koagulation[END_REF]. S. K. Friedlander considered the result of binary collisions occuring simultaneously with splitting of single particles in two new particles [START_REF] Friedlander | On the particle size spectrum of a condensing vapor[END_REF]. This gives the coagulation-fragmentation equation for the density c(x, t) of particles of size (or mass) x at times t (see [START_REF] Dubovskiǐ | Existence, uniqueness and mass conservation for the coagulation-fragmentation equation[END_REF]): x 0 F (x -y, y)c(x, t)dy, x > 0, t > 0,

with the initial condition c(x, 0) = c 0 (x) (x > 0), which is a given nonnegative function. The coagulation kernel K(x, y) describes the rate of formation of a particle of size x + y by coagulation of two particles of size x and y; the fragmentation kernel F (x, y) is describing the rate of formation of particles of size x and y by fragmentation of a particle of size x+y. Both kernels are assumed to be non-negative and symmetric. For an integer m ≥ 0, the moment of order m is defined by: M m (t) := +∞ 0

x m c(x, t)dx. The total number of particles, M 0 (t) may decrease by coagulation or increase by fragmentation, while the total mass of particles, M 1 (t) remains constant if there is neither gelation (formation of particles of infinite size by coagulation: see [START_REF] Aldous | Deterministic and stochastic models for coalescence (aggregation and coagulation): a review of the mean-field theory for probabilists[END_REF]) nor shattering (formation of zero-size particles by fragmentation: see [START_REF] Mcgrady | Shattering" transition in fragmentation[END_REF]). We denote M 1 := M 1 (0).

Existence, uniqueness, boundedness, and positiveness of the solutions of coagulation-fragmentation equation have been considered in [START_REF] Dubovskiǐ | Existence, uniqueness and mass conservation for the coagulation-fragmentation equation[END_REF]. Exact solutions are only known for pure coagulation (F = 0) or pure fragmentation (K = 0) equation and for specific kernels and initial data, so efforts are being pursued to obtain accurate numerical solutions. We focus on stochastic (or quasi-stochastic) strategies. Several Monte Carlo (MC) methods have been developed concurrently. In [START_REF] Eibeck | Approximative solution of the coagulationfragmentation equation by stochastic particle systems[END_REF] the density is approximated using a particle system with a variable particle number. In contrast to this direct simulation scheme, where the particles represent the number density, a mass flow scheme was introduced in [START_REF] Babovsky | On a Monte Carlo scheme for Smoluchowski's coagulation equation[END_REF]: the particles represent the mass density so that the number of particles is kept constant throughout the simulation. The article [START_REF] Goodson | Simulation of coalescence and breakage: an assessment of two stochastic methods suitable for simulating liquid-liquid extraction[END_REF] presents a mass flow scheme for the coagulation-fragmentation equation. In [START_REF] Zhao | Analysis of four Monte Carlo methods for the solution of population balances in dispersed systems[END_REF] time-driven MC are opposed to event-driven MC. In timedriven simulations a time step is chosen, and all events are implemented within that step; in event-driven MC, first an event is selected to occur, and time is advanced by an appropriate increment. It is found that event-driven methods generally provide better accuracy; but time-driven algorithms are more suitable in cases where the equation is to be solved within a larger process simulator that performs explicit integration in time; see [START_REF] Kotalczyk | A time-driven constant-number Monte Carlo method for the GPU-simulation of particle breakage based on weighted simulation particles[END_REF] for recent numerical experiments.

The aim of the present work is to extend to the coagulation-fragmentation equation the time-driven constant-number method introduced in [START_REF] Lécot | A quasi-Monte Carlo method for the coagulation equation[END_REF][START_REF] Lécot | A quasi-Monte Carlo scheme for Smoluchowski's coagulation equation[END_REF] for pure coagulation. The method uses quasi-random numbers (low-discrepancy point sets) in place of pseudo-random numbers. The basic notations and concepts of quasi-Monte Carlo (QMC) methods are given in [START_REF] Niederreiter | Random Number Generation and Quasi-Monte Carlo Methods[END_REF]. We denote I := [0, 1). For a dimension s, the Lebesgue measure is denoted λ s . For U = {u 0 , . . . , u N -1 } ⊂ I s and for a Borel set B ⊂ I s the local discrepancy is For many MC schemes, it is possible to develop corresponding QMC algorithms by replacing the pseudo-random numbers with quasi-random numbers. In [START_REF] Lécot | A quasi-Monte Carlo method for the coagulation equation[END_REF][START_REF] Lécot | A quasi-Monte Carlo scheme for Smoluchowski's coagulation equation[END_REF], it was shown that it is convenient to take special measures in order to benefit from the great uniformity of quasi-random points: as time advances, the simulation particles are renumbered according to size at every time step.

D N (B, U ) := 1 N 0≤k<N 1 B (u k ) -λ s (B),
The present paper is restricted to the presentation and analysis of a QMC (deterministic) scheme, which is able to produce an approximation of the density c(x, t) for any (x, t): this may be of interest in some cases (we think of inverse problems). The remainder of the paper is organized as follows. In Section 2, the QMC algorithm is presented. We establish in Section 3 a convergence result for the discrepancy of the numerical particles relative to the exact mass distribution. Results of computational experiments are shown in section 4.

The quasi-Monte Carlo scheme

The case of pure coagulation was analyzed in our previous paper [START_REF] Lécot | A quasi-Monte Carlo method for the coagulation equation[END_REF]. If F is not identically null, we assume that, for every x > 0, the function y ∈ (0, x) → F (x -y, y) is not identically null. We multiply Eq. (1) by x/M 1 . If we write x = x -y + y in the first integral, it can be split up into two integrals, which are equated by a change of variable; we apply the same transformation to the last integral and we use the symmetry of K and F (see [START_REF] Babovsky | On a Monte Carlo scheme for Smoluchowski's coagulation equation[END_REF][START_REF] Eibeck | Stochastic particle approximations for Smoluchowski's coagulation equation[END_REF][START_REF] Goodson | Simulation of coalescence and breakage: an assessment of two stochastic methods suitable for simulating liquid-liquid extraction[END_REF]). By introducing the mass density function f (x, t) := xc(x, t)/M 1 (which is a probability density function), we obtain the mass-flow equation:

for x > 0, t > 0, ∂f ∂t (x, t) = x 0 K(x -y, y)f (x -y, t)f (y, t)dy - +∞ 0 K(x, y)f (x, t)f (y, t)dy + +∞ 0 F (x + y, x)f (x + y, t)dy - x 0 F (x, y)f (x, t)dy, (2)
where K and F are the modified coagulation and fragmentation kernels respectively defined by: K(x, y) := M 1 K(x, y)/y (for x, y > 0), and F (x, y) := yF (x -y, y)/x (for x > y > 0). We denote f 0 (x) := xc 0 (x)/M 1 the initial data. We introduce a weak formulation of Eq. ( 2), so we define a set of test functions. If R * + := (0, +∞), then a function σ : R * + → R + is said to be simple if its image is a finite pointset of R + ; let S(R * + ) be the set of all measurable simple functions on R * + . We multiply Eq. ( 2) by a simple function σ ∈ S(R * + ) and we integrate over R * + . We write the first (double) integral on the right-hand side of the new equation (by using variables (y, z), where z := x-y):

+∞ 0 +∞ 0 K(z, y)f (z, t)f (y, t)σ(y+z)
dydz and we write the third double integral on the same side (by using variables (x, z), where z := x + y):

+∞ 0 z 0 F (z, x)f (z, t)σ(x)dxdz.
By renaming then integration variables in the first and third double integrals, we obtain a first weak formulation of the mass-flow equation:

d dt +∞ 0 f (x, t)σ(x)dx = +∞ 0 +∞ 0 K(x, y)f (x, t)f (y, t)(σ(x + y) -σ(x))dydx + +∞ 0 x 0 F (x, y)f (x, t)(σ(y) -σ(x))dydx
Before doing the numerical integration, we convert the integral over [0, x] to an integral over [0, 1]. For x ≥ y > 0, let

F(x, y) := y 0 F (x, w)dw, Ḟ(x) := F(x, x), F(x, y) := F(x, y) Ḟ(x) .
For x > 0, let F x be the function y ∈ [0, x] → F(x, y) ∈ [0, 1], with F x (0) := 0. We assume that this function is strictly increasing; let F -1 x be the inverse function. In the integral over [0, x], we change variable y to u := F x (y), and we obtain another weak formulation of the mass-flow equation:

d dt +∞ 0 f (x, t)σ(x)dx = +∞ 0 +∞ 0 K(x, y)f (x, t)f (y, t)(σ(x + y) -σ(x))dydx + +∞ 0 Ḟ(x)f (x, t) 1 0 σ • F -1 x (u)du -σ(x) dx. (3)
We suppose that K and Ḟ are bounded and we set K ∞ := sup x,y>0 K(x, y), Ḟ∞ := sup x>0 Ḟ(x). For b ≥ 2 and m ≥ 1, we put N := b m : this is the number of numerical particles used for the simulation. We need a low discrepancy sequence for the time evolution: U = {u 0 , u 1 , . . .} ⊂ I 3 . For n ∈ N, we write: U n := {u p : nN ≤ p < (n + 1)N }. We make the following assumptions: U is a (t, 3)-sequence in base b (for some t ≥ 0) and

π 3 2 (U n ) is a (0, m, 2)-net in base b (where π 3 2 is the projection defined by π 3 2 (x 1 , x 2 , x 3 ) := (x 1 , x 2 )).

Initialization

We choose a set X 0 = {x 0 0 , . . . , x 0 N -1 } ⊂ R * + of N particles such that the initial mass probability f 0 (x)dx is approximated by the probability distribution:

f 0 (x) := 1 N 0≤k<N δ(x -x 0 k ),
where δ(x -ξ) is the Dirac measure at a point ξ ∈ R * + . This reverts to generating N samples from the density function f 0 .

Time discretization

A fixed time step ∆t is chosen such that ∆t( K ∞ + Ḟ∞ ) < 1. We set t n := n∆t and we denote f n (x) := f (x, t n ). We suppose that a set X n = {x n 0 , . . . , x n N -1 } ⊂ R * + of N particles has been computed so that

f n (x) := 1 N 0≤k<N δ(x -x n k )
approximates (in a sense to be made precise below) the exact mass probability f n (x)dx. The approximation of the solution at time t n+1 is calculated as follows.

(i) Renumbering the particles Particles are relabeled at the beginning of the time step by increasing mass:

x n 0 ≤ x n 1 ≤ • • • ≤ x n N -1 .
This sorting guarantees theoretical convergence: see Section 3 below.

(ii) Coagulation/fragmentation

We define an auxiliary probability measure g n+1 by using an explicit Euler scheme in time for Eq. ( 3):

+∞ 0 g n+1 (x)σ(x) = 1 N 0≤k<N 1 - ∆t N 0≤ <N K(x n k , x n ) -∆t Ḟ(x n k ) σ(x n k ) + ∆t N 0≤k<N 1 N 0≤ <N K(x n k , x n )σ(x n k + x n ) + Ḟ(x n k ) 1 0 σ • F -1 x n k (u)du .
Since we replace the derivative with a finite difference quotient, the measure g n+1 approximates f n+1 (x)dx, but it is not a sum of Dirac measures. We recover this kind of approximation if we use a QMC quadrature rule. Let

R k, := [k/N, (k + 1)/N ) × [ /N, ( +
Γ n+1 σ (u) := 0≤k, <N 1 R k, (u 1 , u 2 ) 1 -χ n k, (u 3 ) -ϕ n k (u 3 ) σ(x n k ) +χ n k, (u 3 )σ(x n k + x n ) + ϕ n k (u 3 )σ • F -1 x n k 1 -u 3 ∆t Ḟ(x n k ) (for u = (u 1 , u 2 , u 3 ) ∈ I 3 ), which is such that +∞ 0 g n+1 (x)σ(x) = I 3 Γ n+1 σ (u)du.
We determine f n+1 by performing a QMC quadrature in

I 3 : +∞ 0 f n+1 (x)σ(x) = 1 N nN ≤p<(n+1)N Γ n+1 σ (u p ).
The calculation on a time step may be summarized as follows. If u ∈ [0, 1), let k(u) := N u . Then, for every p with nN ≤ p < (n + 1)N , we have:

x n+1 k(u p,1 ) =            x n k(u p,1 ) + x n k(u p,2 ) if u p,3 < ∆t K x n k(u p,1 ) , x n k(u p,2 ) , F -1 x n k(u p,1 ) 1 -u p,3 ∆t Ḟ(x n k(u p,1 ) ) if 1 -∆t Ḟ x n k(u p,1 ) ≤ u p,3 , x n k(u p,1 )
otherwise.

For every p, the numbers u p,1 and u p,2 select particles; the particle k(u p,1 ) has for coagulation partner the particle k(u p,2 ) and the coagulation probability is p c := ∆t K(x n k(u p,1 ) , x n k(u p,2 ) ); the fragmentation probability of particle k(u p,1 ) is p f := ∆t Ḟ(x n k(u p,1 ) ). Then u p,3 selects an event and in the case of fragmentation, it determines the size of the new particle. The algorithm is a slight modification of the method devised for pure coagulation and uses threedimensional quasi-random numbers as well. This is interesting, because it is well known that QMC methods perform better in low dimensional spaces. Nevertheless, the extension of the convergence analysis is not straightforward, because fragmentation introduces new linear terms (the coagulation terms are quadratic in c).

Convergence analysis

We prove a convergence result for the QMC scheme previously described. First, we need to adapt the basic tools of QMC methods to the algorithm.

Let f be a probability density function on R * + . If z > 0, then σ z is the indicator function of (0, z). We define the local discrepancy of the set

X = {x k : 0 ≤ k < N } ⊂ R * + relative to f by: D N (z, X; f ) := 1 N 0≤k<N σ z (x k ) - +∞ 0 σ z (x)f (x)dx.
The star discrepancy of X relative to f is:

D N (X; f ) := sup z>0 |D N (z, X; f )|.
The error of the scheme at time t n is defined to be D N (X n ; f n ). The concept of variation of function in the sense of Hardy and Krause can be extended to a function φ defined on R * s + and is denoted by V (φ). The Koksma inequality can be generalized as follows (the proof follows the general outline of the proof of the classical Koksma inequality given in [START_REF] Zaremba | Some applications of multidimensional integration by parts[END_REF]). Proposition 1. Let f be a probability density function over R * + . If φ has bounded variation V (φ) on R * + , then for any

X = {x k : 0 ≤ k < N } ⊂ R * + , 1 N 0≤k<N φ(x k ) - +∞ 0 φ(x)f (x)dx ≤ V (φ)D N (X; f ).
We introduce the following intermediate terms.

• The local truncation error :

ε n z := 1 ∆t +∞ 0 (f n+1 (x) -f n (x))σ z (x)dx - +∞ 0 +∞ 0 K(x, y)f n (x)f n (y)(σ z (x + y) -σ z (x))dydx - +∞ 0 Ḟ(x)f n (x) 1 0 σ z • F -1 x (u)du -σ z (x) dx.
• The additional error : e n z = e n z,K + e n z,F , where

e n z,K := +∞ 0 +∞ 0 K(x, y)f n (x)f n (y)(σ z (x + y) -σ z (x)) - +∞ 0 +∞ 0 K(x, y)f n (x)f n (y)(σ z (x + y) -σ z (x))dydx, e n z,F := +∞ 0 Ḟ(x)f n (x) 1 0 σ z • F -1 x (u)du -σ z (x) - +∞ 0 Ḟ(x)f n (x) 1 0 σ z • F -1 x (u)du -σ z (x) dx.
• The QMC integration error :

d n z := 1 N nN ≤p<(n+1)N Γ n+1 σz (u p ) - I 3 Γ n+1 σz (u)du.
We have the recurrence formula:

D N (z, X n+1 ; f n+1 ) = D N (z, X n ; f n ) -∆tε n z + ∆te n z + d n z .
The local truncation error is bounded as follows.

Lemma 1. If for every x > 0, the function t → f (x, t) is twice continuously differentiable over (0, T ) and if f, ∂f ∂t , ∂ 2 f ∂t 2 are integrable over R * + × (0, T ), then, for t n+1 ≤ T ,

|ε n z | ≤ +∞ 0 t n+1 tn ∂ 2 f ∂t 2 (x, t) dtdx.
The part of the additional error associated to coagulation is bounded in our previous paper (see Lemma 13.4 of [9]). Lemma 2. If for every y > 0, the function K(•, y) : x ∈ (0, +∞) → K(x, y) is of bounded variation V ( K(•, y)) and sup y>0 V ( K(•, y)) < +∞, and if for every x > 0, the function

K(x, •) : y ∈ (0, +∞) → K(x, y) is of bounded variation V ( K(x, •)) and sup x>0 V ( K(x, •)) < +∞, then |e n z,K | ≤ sup y>0 V ( K(•, y)) + sup x>0 V ( K(x, •)) + 3 K ∞ D N (X n ; f n ).
We have the following bound for the part of the additional error associated to fragmentation.

Lemma 3. If for every y > 0, the map F(•, y) : x ∈ [y, +∞) → F(x, y) is of bounded variation V (F(•, y)) and sup y>0 V (F(•, y)) < +∞, then |e n z,F | ≤ sup y>0 V (F(•, y)) + Ḟ∞ D N (X n ; f n ).
Proof. We have

e n z,F = +∞ 0 f n (x) x 0 F (x, y)(σ z (y) -σ z (x))dy - +∞ 0 f n (x) x 0 F (x, y)(σ z (y) -σ z (x))dy.
If we define a function: ϕ(x) := x 0 F (x, y)(σ z (y) -σ z (x))dy (for x > 0), then

e n z,F := 1 N 0≤k<N ϕ(x n k ) - +∞ 0 ϕ(x)f n (x)dx.
Using the generalized Koksma inequality leads to: (•,y)) + Ḟ∞ , the result of the lemma follows.

|e n z,F | ≤ V (ϕ)D N (X n ; f n ). Since V (ϕ) ≤ sup y>0 V (F
For bounding the QMC integration error d n z , we use the following result (see Lemma 3.4 of [START_REF] Niederreiter | Point sets and sequences with small discrepancy[END_REF]): Lemma 4. Let U be a (t, m, s)-net in base b. For every elementary interval J ⊂ I s-1 and for every u s ∈ I, we have:

|D b m (J × (u s , 1), U )| ≤ b t-m .
Lemma 5. If K is of bounded variation V HK ( K) in the sense of Hardy and Krause, and if for every y > 0, the function F(•, y) : x ∈ [y, +∞) → F(x, y) is of bounded variation V (F(•, y)) and sup y>0 V (F(•, y)) < +∞, then

|d n z | ≤ (2 + c K ∆t)b -(m-t)/3 + (1 + c F ∆t) (m-t)/2 ,
where

c K := 4V HK ( K) + 3 K ∞ and c F := sup y>0 V (F(•, y)) + sup x≥y>0 F(x, y). ( 4
)
Proof. The function Γ n+1 σz is the indicator function of some subset P n z of I 3 , thus d n z = D N (P n z , U n ). We have the disjoint union: P n z = (P n z,0 \ P n z,1 ) ∪ P n z,2 ∪ P n z,F , where

P n z,0 := 0≤k, <N x n k <z R k, × I, P n z,1 := 0≤k, <N x n k <z R k, × (I n k, ∪ J n k ), P n z,2 := 0≤k, <N x n k +x n <z R k, × I n k, , P n z,F := 0≤k, <N R k, × J n k ∩ u ∈ I : F -1 x n k 1 -u ∆t Ḟ(x n k ) < z .
We also have the disjoint union: P n z = P n z,0 \ P n z,1 ∪ P n z,2 ∪ P n z,F , where

P n z,1 := 0≤k, <N x n k <z R k, × I n k,
and

P n z,F := 0≤k, <N x n k ≥z R k, × (1 -∆tF(x n k , z) , 1).
Consequently, d n z can be split up:

d n z = d n z,K + D N (P n z,F , U n ), where d n z,K := D N (P n z,0 , U n ) -D N (P n z,1 , U n ) + D N (P n z,2 , U n ).
In [START_REF] Lécot | A quasi-Monte Carlo method for the coagulation equation[END_REF] (Lemma 13.6) we have proved:

|d n z,K | ≤ (2 + c K ∆t)b -(m-t)/3
. We now establish a bound for |D N (P n z,F , U n )|. By denoting I k := [k/N, (k + 1)/N ), we have:

P n z,F = 0≤k<N x n k ≥z I k × I × (1 -∆tF(x n k , z), 1).
We set σ c z := 1 -σ z and we introduce the function:

φ n z (u) := 0≤k<N 1 I k (u)F(x n k , z)σ c z (x n k ), u ∈ I,
where 1 I k is the indicator function of

I k . Then P n z,F = {u ∈ I 3 : u 3 > 1 -∆tφ n z (u 1 )}. Let d 1 be an integer such that d 1 ≤ m -t: its value is determined hereafter. For a 1 ∈ N with 0 ≤ a 1 < b d 1 , we define: I a 1 := [a 1 b -d 1 , (a 1 + 1)b -d 1 ) and P n z,F := 0≤a 1 <b d 1 I a 1 × I × (1 -∆t inf I a 1 φ n z , 1), P n z,F := 0≤a 1 <b d 1 I a 1 × I × (1 -∆t sup I a 1 φ n z , 1), ∂P n z,F := 0≤a 1 <b d 1 I a 1 × I × [1 -∆t sup I a 1 φ n z , 1 -∆t inf I a 1 φ n z ].
We have

P n z,F ⊂ P n z,F ⊂ P n z,F and P n z,F \ P n z,F ⊂ ∂P n z,F . Consequently, D N (P n z,F , U n ) -λ 3 (∂P n z,F ) ≤ D N (P n z,F , U n ) ≤ D N (P n z,F , U n ) + λ 3 (∂P n z,F ).
Since the I a 1 × I are disjoint elementary intervals in base b, an application of Lemma 4 yields: max

|D N (P n z,F , U n )|, |D N (P n z,F , U n )| ≤ b d 1 +t-m . Besides, λ 3 (∂P n z,F ) = ∆t b d 1 0≤a 1 <b d 1 sup I a 1 φ n z -inf I a 1 φ n z .
If we introduce the function

ψ z (x) := F(x, z)σ c z (x) (for x > 0), then φ n z (u) = ψ z (x n k(u)
). We define the sets:

E n a 1 := [x n a 1 b m-d 1 , x n (a 1 +1)b m-d 1 -1 ].
Since the particles are numbered by increasing size, we get:

u ∈ I a 1 ⇒ x n k(u) ∈ E n a 1 . Consequently we have: sup I a 1 φ n z -inf I a 1 φ n z ≤ sup E n a 1 ψ z -inf E n a 1 ψ z .
Since

0≤a 1 <b d 1 sup E n a 1 ψ z -inf E n a 1 ψ z ≤ V (ψ z ) ≤ V (F(•, z)) + sup x>0 F(x, z),
we obtain:

|D N (P n z,F , X n )| ≤ b d 1 +t-m + c F ∆tb -d 1 . By choosing d 1 = (m - t)/2 , we get |D N (P n z,F , U n )| ≤ (1 + c F ∆t)b -(m-t)/2
. In conjunction with the bound for d n z,K , this gives the final result. By combining the results of Lemmas 1, 2, 3 and 5, we obtain an upper bound for the error D N (X n ; f n ); it resembles the final result of [START_REF] Lécot | A quasi-Monte Carlo method for the coagulation equation[END_REF] and is similarly established. Since N = b m , we have 1/b (m-t)/3 = O(1/N 1/3 ) and 1/b (m-t)/2 = O(1/N 1/2 ). Proposition 2. If (1) for every x > 0, the function t → f (x, t) is twice continuously differentiable over (0, T ) and f, ∂f ∂t , ∂ 2 f ∂t 2 are integrable over R * + × (0, T ), (2) K is of bounded variation V HK ( K) in the sense of Hardy and Krause, (3) for every y > 0, the function F(•, y) : x ∈ [y, +∞) → F(x, y) is of bounded variation V (F(•, y)) and sup y>0 V (F(•, y)) < +∞, then

D N (X n ; f n ) ≤ e ctn D N (X 0 ; f 0 ) + ∆t +∞ 0 tn 0 e c(tn-t) ∂ 2 f ∂t 2 (x, t) dtdx + 2 ∆t + c K 1 b (m-t)/3 + 1 ∆t + c F 1 b (m-t)/2 e ctn -1 c
where t n = n∆t, c K , c F are given by Eq. ( 4) and

c := sup x>0 V ( K(x, .)) + sup y>0 V ( K(., y)) + 3 K ∞ + sup y>0 V (F(•, y)) + Ḟ∞ .
This upper bound does not converge as fast as we would like, but is nevertheless a worst-case deterministic bound. The O(N -1/3 ) convergence of the bound is due to QMC integration of indicator functions in 3D; the O(1/∆t) is due to the summation of integration errors at every time step. The same result would have been obtained if the same QMC pointset has been used at every time step (which leads in practice to a method where the error increases with the number of time steps). A more complicated error analysis would take into account that an infinite QMC sequence is used, but in this case, the QMC integration of indicator functions would be in 4D and leads to O(N -1/4 ) convergence. The present bound only guarantees that the error can be arbitrarily small, if ∆t is small enough and N is large enough. It does not aim to provide an optimal time step for a given N . Numerical experiments show that the method converges faster than O(N -1/2 ) and the error decreases with ∆t.

Numerical results

The efficiency of the QMC scheme is tested empirically. To do this, approximate solutions are computed in a case where analytical solutions for the moments are available. A quantity of interest when simulating liquidliquid extraction problems is the Sauter mean diameter of drops, which is calculated with the second and third moments. The QMC solutions are also compared with those given by the MC scheme adapted from the algorithm described in [START_REF] Babovsky | On a Monte Carlo scheme for Smoluchowski's coagulation equation[END_REF]. We consider the case where coagulation and fragmentation occur simultaneously with the following kernels: K(x, y) = 1 and F (x, y) = 2/(x + y), which is considered in [START_REF] Goodson | Simulation of coalescence and breakage: an assessment of two stochastic methods suitable for simulating liquid-liquid extraction[END_REF] (since K and F are unbounded, the hypotheses of Proposition 2 are not satisfied). Then Eq. ( 1) becomes:

∂c ∂t (x, t) = 1 2 x 0 c(x -y, t)c(y, t)dy -c(x, t) +∞ 0 c(y, t)dy +2 +∞ x 1 y c(y, t)dy -c(x, t), x > 0, t > 0. (5)
We restrict ourselves to monodisperse initial condition: c(x, 0) = δ(x -1). By multiplying Eq. ( 5) by x m and using integration over (0, +∞), we get differential equations for the moments M m (t), which can be solved to obtain:

M 2 (t) = 3 -2e -t/3 and M 3 (t) = 18 -36e -t/3 + 19e -t/2 .

We compute the solution up to time T = 1.0 with N particles (N varying from 2 7 to 2 20 ) and P time steps (P varying from 2 0 × 100 to 2 4 × 100).

All QMC computations use a (1, 3)-sequence in base 2 constructed by H.

Niederreiter [START_REF] Niederreiter | Random Number Generation and Quasi-Monte Carlo Methods[END_REF]. The moment M m (t n ) at time t n can be approximated, using the MC or QMC scheme, by:

M m,N,P,n := 1 N 0≤k<N (x n k ) m-1 . (6)
First we compute the absolute error in the moment estimation at time αT , for α = 1/100, 1/10 and 1:

E α m,N,P := |M m (αT ) -M m,N,P,αP |.
Figure 1 represents the log-log plots (base 2) of these absolute errors, for different values of N and P . In all cases, the QMC simulations give better re- sults (for the same discretization parameters) than those given by the Monte Carlo scheme. The error seems stable on the time interval [1/100, 1/10], so that it makes sense to compare the methods on the whole time interval [0, 1]: the mean absolute error in the moment estimation is calculated as follows: Figure 2 shows the log-log plots (base 2) of the mean absolute error in the estimation of the moments of order 2 and 3, for different values of N and P . In both cases, the QMC simulations show faster convergence than the Monte Carlo scheme.

Conclusion

In this paper, we have presented and analyzed a new algorithm for the approximation of the continuous coagulation-fragmentation equation. A sample of N numerical particles is used to simulate the behavior of the system as a whole. Time is discretized into P time steps and the mass density is approximated by a sum of N Dirac measures. Three-dimensional quasi-random points are used to change the particles masses according to the dynamics of the equation.

A deterministic error bound is established. The accuracy of the algorithm is assessed through comparison with exact values, in a test case where analytical results are available. The numerical approximations given by the QMC scheme converge to the exact results when P and N are large enough. The errors given by the QMC simulation are always smaller than the errors given by the MC method using the same P and N .

As announced, we have restricted the study to a deterministic scheme. Randomized QMC would be an interesting alternative (see [START_REF] L'ecuyer | A randomized quasi-Monte Carlo simulation method for Markov chains[END_REF][START_REF] L'ecuyer | Sorting methods and convergence rates for Array-RQMC: Some empirical comparisons[END_REF]), but requires a different analysis: this will be the subject of a forthcoming paper.

0 KF

 0 -y, y)c(x -y, t)c(y, t)dy-+∞ (x, y)c(x, t)c(y, t)dy + +∞ 0 (x, y)c(x + y, t)dy -1 2

where 1 B

 1 denotes the indicator function of B. The discrepancy of U is D N (U ) := sup J |D N (J, U )|, where the supremum is taken over all subintervals J ⊂ I s . The star discrepancy of U is D N (U ) := sup J |D N (J , U )|, where J runs through all subintervals of I s with a vertex at the origin. QMC methods are versions of MC methods, where the pseudo-random samples are replaced with low-discrepancy point sets. Powerful methods of constructing low-discrepancy point sets are based on the theory of (t, m, s)-nets and (t, s)sequences. For an integer b ≥ 2, an elementary interval in base b is an interval of the form s i=1 [a i b -d i , (a i + 1)b -d i ), with integers d i ≥ 0 and 0 ≤ a i < b d i . If 0 ≤ t ≤ m are integers, a (t, m, s)-net in base b is a point set U of N = b m points in I s such that D N (J, U ) = 0 for every elementary interval J in base b with measure b t-m . If b ≥ 2 and t ≥ 0 are integers, a sequence u 0 , u 1 , . . . of points in I s is a (t, s)-sequence in base b if, for all integers n ≥ 0 and m > t, the points u p with nb m ≤ p < (n + 1)b m form a (t, m, s)-net in base b.

  1)/N ) ⊂ I 2 be an elementary interval in base b: we denote by 1 R k, its indicator function. Let χ n k, be the indicator function of I n k, := [0, ∆t K(x n k , x n )) and let ϕ n k be the indicator function of J n k := [1 -∆t Ḟ(x n k ), 1). To σ ∈ S(R * + ) corresponds the indicator function:

Figure 1 :

 1 Figure 1: Error in the second (left) and third (right) moment estimation at time 1/100 (top) 1/10 (middle) and 1 (bottom), as a function of N (number of particles varying from 2 7 to 2 20 ), for P (number of time steps) between 100 and 1 600. Log-log plots of the MC (dashed line) versus QMC (solid lines) results.

  hT /100) -M m,N,P,hP/100 |.

Figure 2 :

 2 Figure 2: Mean error in the second (top) and third (bottom) moment estimation as a function of N (number of particles varying from 2 7 to 2 20 ), for P (number of time steps) between 100 and 1 600. Log-log plots of the MC (dashed line) versus QMC (solid lines) results.