Long-Term Changes in Four Populations of the Spiny Toad, Bufo spinosus , in Western France; Data from Road Mortalities

: Habitat fragmentation is widely recognized as a contributor to the decline of biodiversity, with amphibians one of the key groups impacted. To understand the effects of habitat fragmentation on amphibian populations requires long-term data sets showing population trends. In this paper, road mortalities were employed as proxies to describe long-term numbers of four populations of the spiny toad Bufo spinosus in western France during a 17-year period. Road mortalities were found during all months in all populations but were most frequent during October, November and December, the main migratory period. Large females were found signiﬁcantly more frequently during these migration months, forming 45% of the total sample, compared with their presence from January to September (34.4%). The long-term trends were evaluated using regression analysis of the logarithmic (log e ) transforms of annual counts as dependent variables against year as the independent variables. All coefﬁcients showed no signiﬁcant departure from the 0 hypothetical coefﬁcients, indicative of population stability. This was supported by jackknife analysis, which showed good agreement of the pseudo-regression coefﬁcients with the true equations. Stepwise regression of potential climate impacts on toad numbers suggested rainfall levels in October adjusted to 2- and 3-year lags were involved in driving population change. Road mortality counts were also made during 2020 and 2021 when human movement restrictions were in place due to the COVID-19 pandemic. To estimate the potential impact on this disturbance in the methodology, the Poisson distribution was used to estimate potential differences between what would have been expected counts and the observed counts. The results indicate that the observed mortalities were signiﬁcantly lower than expected in all four populations.


Introduction
Amphibians are now acknowledged as key bio-indicators of ecosystems [1]. However, numerous studies also indicated that they are now one of the most endangered vertebrate groups [2], with habitat loss, climate change, extensive use of agrochemicals, increased UVirradiation and emerging diseases among the factors cited as the causes for the decline [3,4]. Their vulnerability is in part due to many species being dependent on both aquatic and terrestrial environments, which along with a relatively permeable integument renders them particularly susceptible to environmental changes. Identifying the causes of the declines is problematical since amphibian populations often fluctuate over several orders of magnitude [5,6]. Nevertheless, separating long-term declines from short-term fluctuations is a well-known problem [7][8][9][10] and to understand and verify what factors are involved requires a data base of long-term trends, which may require investment in funding as well as time over many years. Additionally, amphibian populations are often in disequilibrium [11], which may mask long-term trends, but for conservation management programs, they are very much needed [12].
A series of methodologies have been used to estimate long or short-term trends [13]. Included are road mortalities, that as well as being employed as proxies for monitoring population changes [13][14][15][16] were also used to estimate changes in, for example, numbers of mammals [17][18][19][20], lizards [21] and snakes [22,23]. The value of the method when applied to amphibians is illustrated in the study of the common European toad Bufo bufo that showed a strong correlation between road mortalities and numbers arriving at breeding ponds [13]. Road mortalities also facilitate detection of amphibians on the open spaces of road surfaces [2,23] and represent independent replicates if carcasses are removed after counting. This avoids double counting and autocorrelation and when applying regression analysis to long time series the equations have 'no memory'-counts for one year do not impact on the following year's counts in a significant way [24].
The spiny toad Bufo spinosus (formerly Bufo bufo in western Europe) is a widespread species of anuran in Europe occupying a range of different habitats, including fragmented landscapes ( Figure 1). Mostly active nocturnally [25], they are explosive breeders and are known for breeding-site fidelity [23,[26][27][28][29][30]. In the study locality, they migrate to large bodies of water in the later months of the year but reproduce in the early months of the following year. In fragmented habitats they are often constrained to crossroads, usually during the hours of darkness that often results in mortalities. Despite reports of population declines in various regions of Europe [13,31,32], it is listed as "least concern" by the IUCN, but as a common species, it has potential to provide insight into long-term population trends and the health of ecosystems [10]. Numerous studies indicated regional differences in population ecology in both B. spinosus and the related B. bufo (e.g., [5,29,[33][34][35][36][37]. A review of long-term changes in amphibian numbers indicated B. spinosus populations were declining at a slightly greater rate than they were increasing (52 % (n = 394) versus 47.5% (n = 360) with n = 4 (0.04%), with populations in the survey showing no change [10]. This suggests a need for additional population data, particularly at the meta-population level and over wider areas of Bufo distribution. [13]. Included are road mortalities, that as well as being employed as proxies for toring population changes [13][14][15][16] were also used to estimate changes in, for ex numbers of mammals [17][18][19][20], lizards [21] and snakes [22,23]. The value of the m when applied to amphibians is illustrated in the study of the common European toa bufo that showed a strong correlation between road mortalities and numbers arri breeding ponds [13]. Road mortalities also facilitate detection of amphibians on th spaces of road surfaces [2,23] and represent independent replicates if carcasses moved after counting. This avoids double counting and autocorrelation and wh plying regression analysis to long time series the equations have 'no memory'for one year do not impact on the following year's counts in a significant way [24].
The spiny toad Bufo spinosus (formerly Bufo bufo in western Europe) is a wide species of anuran in Europe occupying a range of different habitats, includin mented landscapes (Figure 1). Mostly active nocturnally [25], they are explosive ers and are known for breeding-site fidelity [23,[26][27][28][29][30]. In the study locality, th grate to large bodies of water in the later months of the year but reproduce in th months of the following year. In fragmented habitats they are often constrai crossroads, usually during the hours of darkness that often results in mortalities. D reports of population declines in various regions of Europe [13,31,32], it is listed as concern" by the IUCN, but as a common species, it has potential to provide insig long-term population trends and the health of ecosystems [10]. Numerous studie cated regional differences in population ecology in both B. spinosus and the rela bufo (e.g., [5,29,[33][34][35][36][37]. A review of long-term changes in amphibian numbers ind B. spinosus populations were declining at a slightly greater rate than they were i ing (52 % (n = 394) versus 47.5% (n = 360) with n = 4 (0.04%), with populations survey showing no change [10]. This suggests a need for additional population particularly at the meta-population level and over wider areas of Bufo distribution In this paper, a continuous 17-year data set derived from four populations i gion of western France is described. However, unforeseen perturbations in the m ology during 2020 and 2021 due to the COVID-19 pandemic and subsequent chan road traffic levels necessitated that these data be omitted from the long-term an giving a continuous 15-year data set. The following five questions are addressed: In this paper, a continuous 17-year data set derived from four populations in a region of western France is described. However, unforeseen perturbations in the methodology during 2020 and 2021 due to the COVID-19 pandemic and subsequent changes in road traffic levels necessitated that these data be omitted from the long-term analysis, giving a continuous 15-year data set. The following five questions are addressed: (1) What are the general long-term population trends in the study area?

250
(2) Are there differences in long-term numbers between populations?-here defined as subpopulations migrating to different breeding ponds. (3) To what extent do the populations fluctuate on an annual basis? (4) What environmental factors may act as drivers of the populations? This is important since, among other things, rainfall patterns were shown to be influential on amphibians by impacting critical factors such as sex determination [38,39]. (5) What impact did the restrictions on human movement due to the COVID-19 pandemic have on road mortalities during 2020 and 2021?

Study Area
The study area (46 • 27 N; 1 • 53 W) is a fragmented landscape in Vendee, western France ( Figure 2). The locality is dominated by agriculture that had experienced little or no major changes in land use during the study period from 2005 to 2021. The general climate is mild oceanic with mean annual air temperature 13-14 • C. July-August are the hottest and driest months (mean ≈ 22 • C) with high precipitation from October until February. Mortalities for B. spinosis were collected on roads (distance ≈ 26 km) between wetland areas south and east of the town of Lucon (D949/D746) and near the villages of Chasnais (D44/D127), Lairoux (D60) and St Denis du Payre and Grue (D25). The straight-line distance between the main breeding ponds are: D949/D746 to D44/D127, 3.7 km; D44/D127 to D60, 3.2 km; D60 to D25, 3.32 km. Figure 3 shows examples of breeding ponds in the study area. The four populations were treated as distinct populations, in part due to the distance between the breeding ponds and known migratory distances or interactions between subpopulations of B. spinosis, which is a function of the distance between them [26,27,30]. For example, interactions between B. spinosus subpopulations were significantly reduced when the distance was greater than 300 m, with none recorded at 830 m [30]. In Germany, it was observed that certain individuals moved more than 1 km [29], but these distances are lower than the distances between ponds in the present study.
(3) To what extent do the populations fluctuate on an annual basis? (4) What environmental factors may act as drivers of the populations? This is important since, among other things, rainfall patterns were shown to be influential on amphibians by impacting critical factors such as sex determination [38,39]. (5) What impact did the restrictions on human movement due to the COVID- 19 pandemic have on road mortalities during 2020 and 2021?

Study Area
The study area (46°27′ N; 1°53′ W) is a fragmented landscape in Vendee, western France ( Figure 2). The locality is dominated by agriculture that had experienced little or no major changes in land use during the study period from 2005 to 2021. The general climate is mild oceanic with mean annual air temperature 13-14 °C. July-August are the hottest and driest months (mean ≈ 22°C) with high precipitation from October until February. Mortalities for B. spinosis were collected on roads (distance ≈ 26 km) between wetland areas south and east of the town of Lucon (D949/D746) and near the villages of Chasnais (D44/D127), Lairoux (D60) and St Denis du Payre and Grue (D25). The straight-line distance between the main breeding ponds are: D949/D746 to D44/D127, 3.7 km; D44/D127 to D60, 3.2 km; D60 to D25, 3.32 km. Figure 3 shows examples of breeding ponds in the study area. The four populations were treated as distinct populations, in part due to the distance between the breeding ponds and known migratory distances or interactions between subpopulations of B. spinosis, which is a function of the distance between them [26,27,30]. For example, interactions between B. spinosus subpopulations were significantly reduced when the distance was greater than 300 m, with none recorded at 830 m [30]. In Germany, it was observed that certain individuals moved more than 1 km [29], but these distances are lower than the distances between ponds in the present study.

Protocol
Surveying for road B. spinosus mortalities commenced in 2005 on the road networks close to these ponds and were made between four and six times per month from 1 January through to 31 December each year by a single observer on a bicycle at a speed of around 5-10 km/h. Tests of comparison of carcass detection showed that detection increased when using a bicycle compared with motorcar or motorcycle. When found, each individual was measured for snout to vent length (SVL) and removed from the road to avoid double counting. Road traffic volume increased slightly during the study period of 2005 to 2019, with counts showing means of 11.1-60.2 vehicles per hour depending on the road and were greatest on D949 day or night [15].

Statistical Analysis
Long-term population changes were evaluated using regression analysis with annual counts treated as the dependent variables after transformation to logarithms (Loge) and years as the independent variable. This gives, where logeN represents road mortalities, x the year (2005 to 2019), m the regression coefficient and b the y-intercept with ε white noise error, which incorporates measurement error and random variation [24]. Due to a zero count for D25 in 2008, the dependent variable was treated as loge(N + 1). The null hypothesis is that logeN (or loge(N + 1)) is stable when m = 0; significant departures from m indicate population change with positive or negative regression coefficients representing population increase or decrease, respectively. Departures from the 0 regression coefficients were evaluated using t-tests at

Protocol
Surveying for road B. spinosus mortalities commenced in 2005 on the road networks close to these ponds and were made between four and six times per month from 1 January through to 31 December each year by a single observer on a bicycle at a speed of around 5-10 km/h. Tests of comparison of carcass detection showed that detection increased when using a bicycle compared with motorcar or motorcycle. When found, each individual was measured for snout to vent length (SVL) and removed from the road to avoid double counting. Road traffic volume increased slightly during the study period of 2005 to 2019, with counts showing means of 11.1-60.2 vehicles per hour depending on the road and were greatest on D949 day or night [15].

Statistical Analysis
Long-term population changes were evaluated using regression analysis with annual counts treated as the dependent variables after transformation to logarithms (Log e ) and years as the independent variable. This gives, where log e N represents road mortalities, x the year (2005 to 2019), m the regression coefficient and b the y-intercept with ε white noise error, which incorporates measurement error and random variation [24]. Due to a zero count for D25 in 2008, the dependent variable was treated as log e (N + 1). The null hypothesis is that log e N (or log e (N + 1)) is stable when m = 0; significant departures from m indicate population change with positive or negative regression coefficients representing population increase or decrease, respectively.
Departures from the 0 regression coefficients were evaluated using t-tests at n-2 d.f. Tests for normality of the log e transformed annual counts were made using Anderson-Darling a 2 tests. This gave results from a 2 = 0.18 to 0.39 and P values from 0.33 to 0.89; hence, all data were normalised before analysis. The Pearson product moment correlation coefficient was used to compare for similarities in monthly counts.
Jackknifing was used to estimate the potential errors in the true regression coefficients [24]. This nonrandom method produces repeatable results by systematically removing one-year data sets from the sample and re-applying regression analysis. This gives a series of pseudo-m values and enables comparison of the pseudocoefficients with the true coefficients.
Tests for homogeneity of SVL variance between metapopulations were conducted by Leven's test for multiple variances. This test is less sensitive to departures from normality and considers the distances of the observations from their sample medians. The null hypothesis is that SVL variance is in agreement between populations, Ho: δ 2 1 . . . = . . . δ 2 k . ANOVA was then used to test for differences in mean SVLs between populations, with the null hypothesis the means are homogenous. Post hoc tests were conducted with Tukey HSD. Leven's test was also used to test for monthly presence with a null hypothesis of equality of months. Long-term interpopulation correlation of trends between species were assessed with the Pearson correlation coefficient at α = 0.05. To compare the presence of large (≥80 mm SVL) mortalities during the main migration months of October, November and December with their frequencies during January to September, a two-tail z-proportion test was applied at α = 0.05.
To evaluate potential climatic effects on changes in annual counts, stepwise regression using climate data as independent or predictor variables were then compared with annual counts of B. spinosus as the dependent variable. Climatic variables were long-term monthly and total annual rainfall along with annual temperatures during the study period. These data were further analyzed by lagging the climate data to 1, 2 and 3 years previous to annual counts. Initial analysis was by taking the response variable and applying regression analysis to each independent variable for inclusion in the final model [24], with the significance level for independent variable inclusion in the final models α = 0.15. Variables that increased significance of the model were included, those that reduced the significance level were omitted. The regression model has the form where m 1 , m 2 and m 3 . . . are regression coefficients and x 1 , x 2 and x 3 . . . the variables that had a significant influence on changes in N bufo . Independent variable selection is arbitrary, but the method is practical when the number of candidate variables is not large [24].
Data collected during 2020 and 2021 were not used in the long-term data analysis due to the COVID-19 restrictions on people's movement and road traffic volumes on B. spinosus road mortalities. The Poisson distribution was used to evaluate any impact this had on road mortalities. The Poisson equation has the form where the random variable x was defined as the observed number of mortalities during 2020 and 2021, λ the rate parameters as the means of the annual counts from the 15-year time series and e the base natural logarithm. If P (x ≥ λ), this would indicate the observed values of x are reasonable counts and were little influenced by the restrictions. The alternative, P (x < λ), is that the 2020 and 2021 counts are likely due to changes in traffic density and, hence, if included in the analysis would have distorted the long-term trends.
Data analysis was carried out mainly on Minitab 17 along with VassarStats and various other Internet sites.

Monthly Counts
The results here were taken from the 17-year series-that is including 2020 and 2021. Mortalities were found during all months but were negatively skewed towards the latter months in all areas (October, November, December; see Figure 4). A χ 2 goodness-of-fit test with a null hypothesis of equality of month presence (expected = 102.1 toads) versus observed monthly numbers indicated a significant departure from the null model, χ 2 = 922.6 and P < 0.0001. The Pearson correlation coefficient showed strong monthly correlations between populations, r-values from 0.84 to 0.99 and p-values from 0.001 to 0.0001.
Conservation 2022, 2, FOR PEER REVIEW 6 Data analysis was carried out mainly on Minitab 17 along with VassarStats and various other Internet sites.

Monthly Counts
The results here were taken from the 17-year series-that is including 2020 and 2021. Mortalities were found during all months but were negatively skewed towards the latter months in all areas (October, November, December; see Figure 4). A χ 2 goodness-of-fit test with a null hypothesis of equality of month presence (expected = 102.1 toads) versus observed monthly numbers indicated a significant departure from the null model, χ 2 = 922.6 and P < 0.0001. The Pearson correlation coefficient showed strong monthly correlations between populations, r-values from 0.84 to 0.99 and p-values from 0.001 to 0.0001.

Frequency of Large Females
In general (all regions pooled), large (≥80 mm SVL) individuals were found more frequently during the main migration months, which were October, November and December, representing 45% of the total sample, compared with the remaining months from January to September, when they reduced to 34.4%. The difference between months was significant z = 3.72 and P = 0.0001.

Annual Differences in SVL between Populations
Tests for multiple variances and post hoc Tukey HSD indicated SVL variances in the four populations were in agreement (F = 2.13, P = 0.09). However, ANOVA detected differences between means (F = 27.19, P < 0.0001) with Tukey HSD showing all means differed from each other, with the D949/764 showing the largest individuals ( Figure 5; Table 1).

Frequency of Large Females
In general (all regions pooled), large (≥80 mm SVL) individuals were found more frequently during the main migration months, which were October, November and December, representing 45% of the total sample, compared with the remaining months from January to September, when they reduced to 34.4%. The difference between months was significant z = 3.72 and P = 0.0001.

Annual Differences in SVL between Populations
Tests for multiple variances and post hoc Tukey HSD indicated SVL variances in the four populations were in agreement (F = 2.13, P = 0.09). However, ANOVA detected differences between means (F = 27.19, P < 0.0001) with Tukey HSD showing all means differed from each other, with the D949/764 showing the largest individuals ( Figure 5; Table 1).

Long-Term Annual Counts
Regression analysis gave positive values of m in three of the four areas, indic long-term population stability. The population living around D949/D746 showed a ative value of m, but the t-test against the 0 coefficient indicated no significant depa from long-term stability. Regression analysis of pooled data indicated in general p lations were stable during the study period (m = 0.1 ± 0.02, t = 4.79). The Jackknife ana showed good agreement with the true equations but also identified unusually high or numbers that were unexpected given trends at the time. The results are shown an plained in full in Table 2 along with regression lines derived from the equations in F 6A,B.
The means of long-term annual counts ranged from 12 per year (D949/746), (D44/D127), 21.7 (D25) and 24.1 (D60) with F = 1.76, P = 0.16. The post hoc Tukey also showed no differences between inter-population variances in counts with P-va from 0.13 to 0.90 and, hence, long-term population fluctuations were similar.

Long-Term Annual Counts
Regression analysis gave positive values of m in three of the four areas, indicating longterm population stability. The population living around D949/D746 showed a negative value of m, but the t-test against the 0 coefficient indicated no significant departure from long-term stability. Regression analysis of pooled data indicated in general populations were stable during the study period (m = 0.1 ± 0.02, t = 4.79). The Jackknife analysis showed good agreement with the true equations but also identified unusually high or low numbers that were unexpected given trends at the time. The results are shown and explained in full in Table 2 along with regression lines derived from the equations in Figure 6A,B. Table 2. Regression coefficients (m) relating annual mortalities between 2005 and 2019 along with results of the Jackknife analysis (mJK). For m, the ±values are standard errors of the coefficients and the t-tests against the 0 regression coefficients. When m was positive, numbers were generally increasing; when m was negative, numbers were decreasing. The latter was only found for D949/D764, but the result was not significant. For the Jackknife results, means and standard errors are the mean values of the 15-year pseudoregression coefficients with their standard deviation. Right column identifies years when influence functions were detected, indicating unexpected low (L) or high numbers (H) given current trends at the time. All t-tests were set at n − 2 years (15 − 2 = 13 years). See text for further details.   Table 2. See text for further details.   Table 2. See text for further details.

Rainfall, Temperature and Annual Counts
Since the distances between metapopulations were relatively small, climate events were unlikely to be different and so data were pooled. The stepwise regression was applied to a total of 20 data sets and indicated rainfall as a potential influencer on population change when data were lagged to two or three years previous to the year counts. Three data sets were within the P ≤ 0.15 threshold for inclusion in the final model. These were rainfall levels in October adjusted to a 3-year lag (r = 0.79 and P = 0.01), October adjusted to a 2-year lag (r = 0.64 and P = 0.13) and adjusted to a three-year lag in December (r = 0.63 and P = 0.13). The model gave significance at F = 3.47 and P = 0.05. If, however, 'three-year lag in December' was removed, the probability increased to F = 5.61 and P = 0.02, giving a final best-fit model of where N bufo is B. spinosus annual counts. These results suggest that rainfall during October lagged to two and three years previous that from the year of annual counts, which is the beginning of the main migratory movement to breeding sites and may be involved as a key population driver.

Mortalities during 2020 and 2021
The results of the Poisson tests indicated significant departures of the random variables x from the lambda's (x = λ) in all populations, P < 0.0001 to P = 0.04. This confirms the notion that movement restrictions during 2020 and 2021 had a significant effect in lowering mortalities. Estimates of B. spinosus mortalities in the absence of lower traffic volumes are based on the method in Gotelli and Ellison (2004) but here using where the estimate E is derived from λ, the rate parameters in the Poisson equations, P the probability that x > λ, with n as the sample sizes. The values of λ and x and subsequent calculations for each population are fully shown in Table 3. Table 3. Poisson probabilities that road mortality counts during 2020 and 2021 are a true measure of expected mortalities P (λ = x). The rate parameters (λ) are derived from the mean values of the long-term annual counts (15 years) and the random variable x from the observed numbers in each area during 2020 and 2021. The likelihood that the true counts were greater than the rate parameters (x > λ) ranged from 84-100% during 2020 and greater than 99% for 2021.

General Considerations
The results presented here give good insight into questions 1 and 2, indicating that generally all four populations showed asynchronous long-term trends. However, when data from the four populations were pooled, the results showed general stability over the wider area ( Figure 6B). Counts from the D949/D746 area were the most frequent during the early period from 2005 but declined to the 4th rank by 2013 and remained so to the end of the study period. If the populations are operating as metapopulations, this would suggest that asynchronous dynamics may operate as a buffer against the likelihood of simultaneous population crashes and extinctions [40] due to varying degrees of exchange between breeding populations. Hence, although B. spinosus is well known to be philopatric to breeding ponds, some individuals may travel extensively, often over 1 km [29], and potentially colonise areas with depleted numbers. However, interactions between breeding populations apparently declined with increasing distances between breeding ponds, with no interactions above 800 meters, but many individuals (up to 21% of the population) may visit different breeding ponds [30]. Data from the D25 area, for instance, showed relatively sudden and major increases in numbers beginning in 2013 and by the following year the highest counts of the four areas ( Figure 6A). This was mirrored by similar increases in live counts in this area. In contrast, numbers in the D44/D127 and D60 areas were generally stable throughout the study period and, hence, it is conceivable that individual toads from these areas could potentially recolonize the D949/D764 area due to the distances between pond breeding areas (Figure 1), albeit in low numbers given our understanding of toad movements.
What drove the relatively sudden increases around the D25 area? The construction of several new ponds in the locality area during the study period may be important. For example, created ponds have been shown to provide partial mitigation for the loss of the natural breeding areas for amphibians [40][41][42][43]. Game fish species are present in the larger ponds ( Figure 2) and ornamental fish are often introduced into many larger garden ponds (including one next to D127), which is a known associated factor in toad breeding pond preference [2]. For B. spinosus, predatory fish may positively impact on the numbers of alien crayfish Procambarus clarkia that are well-known predators on all stages of amphibian development including adults [44], thus enhancing larvae survivorship. Only one potential breeding pond was established in the other areas during the study period and that was adjacent to the D127, but this was very close to an already major B. spinosus breeding pond with game fish present (Figure 2), while the new pond did not have any game fish. This new pond was quickly colonized by P. clarkia, which in subsequent years impacted the numbers of early amphibian colonizers, with only Rana dalmatina continuing to use it as a breeding site, albeit in low numbers thereafter [16]. Created ponds have also been shown to provide partial mitigation for the loss of the natural amphibian breeding habitat due to differences in hydrologic regimes, vegetation and surrounding terrestrial habitats between created ponds and natural wetlands [42].
Theory predicts, and is supported by studies, that the survival of adults in a population is crucial for the capacity of a population for growth and long-term stability [43,44]. A positive finding of the present study was the long-term presence, as proportions of total annual numbers, of large females in all four samples. Some of these individuals were particularly large (110 mm + SVL) suggesting good age spans and, hence, high reproductive potential. Female common toads (Bufo bufo) in the UK were found to reach life spans of 9 years with males attaining 6 years [43], although in harsh environments females may take 3 to 7 years to reach maturity and many reproduce only once in a lifetime [36)] The finding that in all four populations annual counts fluctuated extensively (question 3) is a natural dynamic in amphibians [7,45], as are local extinctions, the latter event observed in the green frog Pelophylax lessonae in the study locality [6]. Such dynamics are often associated with pond living/breeding biphasic species, where population fluctuation is often determined by changes in the environment [10] with populations likely to be declining most of the time, followed by short periods of increases in numbers [10]. A classic example is a study of B. americanus that initially showed an 18-year trend of long-term population decline but with a rapid recovery during the following 7-year period [46,47]. The data from the present study, albeit with shorter time periods, along with other studies, support these trends [6,14]. Population stability has been shown not to be a key requirement for long-term amphibian persistence, since many species that show high fluctuations also have a high capacity for recovery [10], even from local extinctions [2,6,14].

Comparisons with Sympatric Amphibians
Long-term trends were, in general, in good agreement with population trends in sympatric amphibians recorded in previous studies [6,16]. For example, high numbers of B. spinosus around 2009-2012 were similarly observed in sympatric Triturus marmoratus, Hyla arborea and Pelophylax lessonae [6], with Rana dalmatina also showing comparable long-term trends, albeit based on combined long-term spawn mass and road mortality counts [16]. The high counts recorded during 2009-2012 were followed by a decline, then numbers increased again from approximately 2016. Interspecific similarities were also found in the stepwise regression for climatic data, where lagged rainfall patterns during periods of high movement in November and December were highlighted as potential population influencers [6].

Rainfall and Temperature
The finding of rainfall (albeit lagged rainfall data-question 4) as a potential driver has been shown to affect breeding populations of other species of amphibians by influencing breeding activity and past recruitment [10,11]. Rainfall can alter amphibian habitats, especially water levels, for example small ponds and ditches can dry out early before metamorphosis is complete, while in some years flooding can move spawn to unsuitable habitats (unpublished data for R. dalmatina). Rainfall can also alter water temperatures that in turn may influence sex ratios [38]. The effects of such events on B. spinosus in the study locality are probably much fewer compared with sympatric amphibians given the large water bodies they use for breeding. However, when B. spinosus data were pooled, long-term numbers showed remarkable stability and indeed numbers increased compared with five sympatric species [6,16].

Impact of Lockdown Restrictions
In the present study, lockdowns also impacted sampling frequency, although toad carcasses often have very long duration time on roads [15] and, hence, it is proposed that the observed numbers of carcasses were indeed reduced and not simply an artefact of sampling during lockdown. Therefore, the results from the Poisson distribution from data gathered during the pandemic give an unexpected insight on the significant impact of traffic on toad mortalities during normal times. Previously, researchers have had to rely predominantly on mostly observational approaches to understand road traffic effects on mortalities, but now data from previous years and future years after the anthropause should enable greater understanding of traffic volumes on road mortalities. The data from 2020 and 2021 and the Poisson analysis thus provide an answer to question 5, indicating even short periods of traffic reduction or road closures could be used as part of a management strategy for the conservation of amphibian populations. These may only involve minor changes but could potentially have major benefits for ecosystems and in the longer term, humans [48,49].

Concluding Remarks
Long-term amphibian time series give insight into the range of natural fluctuations that populations experience and present testable hypotheses of long-term trends, including declines. Monthly time series provide additional information on the expected numbers that enter roads outside the main migratory period of these populations, which is important because in other areas of Europe the main migratory period is during early spring, for e.g., [13,34]. This information, together with long-term data series, is key for implementing mitigation measures to at least reduce road mortalities. However, simple times series are unlikely to provide a detailed understanding of the underlying factors involved in amphibian population trends, especially those with complex life cycles. For example, density-dependent regulation is a possible candidate for long-term regulation, but this may occur at the larval, juvenile or adult stage or even all three stages with the larvae interacting in different ways to the terrestrial stages. This can result in highly complex dynamics and wide fluctuations in adult numbers over long time periods. The crucial role amphibians play in many ecosystems can influence the populations of other species, including snakes [48], emphasising the importance of increasing understanding of changes in their population status, especially the nature of their declines. Future data of long-term time series are a critical first step for understanding amphibian population ecology, but measures to reduce road mortalities have already been in place in several European states, for example, in France where underground tunnels for safe road crossings were installed under the A83 autoroute [50] and in the UK where the British Herpetological Society designed and constructed gully pot ladders to enable amphibians to easily escape from being trapped in roadside drains. These are available at very low price and tests have shown that the ladders have the potential to reduce mortalities significantly and have been used by local authorities, road construction companies and developers [51]. However, simple signposts can also be of value, reminding drivers that they are approaching toad and other amphibians' road-crossing areas, which is especially useful when the main migratory periods are known.
Funding: This research received no external funding.

Institutional Review Board Statement:
The study did not require ethical approval since the data was generated from road mortalities of B. spinosus.