Maps and altitude charts

Along with abundance and distribution, the type of data available to us differs from species to species. For some species, all nest sites are known; for others, we can calculate the density per kilometre square based on the territory mapping surveys. In the case of species that are difficult to detect, we can only estimate the probability of occurrence, but have no information on density. To account for these differences, we designed several map types to represent a species’ occurrence and the changes since 1993–1996. We decided for each species which map type would be most informative. The table in the section «Survey methods» (p. 53-59) indicates which map types were used to depict the situation in 2013–2016 and the changes since 1993–1996.

Density map 2013–2016

For 75 relatively widespread and at least locally abundant species we were able to generate density maps using the territory mapping data.

Based on the number of mapped territories per survey visit and 16 different environmental variables, we estimated the number of territories per species for all kilometre squares in Switzerland and Liechtenstein. We used a binomial mixture model, which corresponds to a Poisson regression model that also incorporates the species' detection probability. This allowed us to account for the fact that not all birds present were detected during the visit.

We calculated the values of the 16 environmental variables for all kilometre squares located fully or partly in Switzerland as well as for the surveyed kilometre squares in neighbouring countries. For seven land-cover variables (1–7 in the table), the percentage area in each kilometre square in Switzerland and Liechtenstein was elicited from the Swiss Land Use Statistics based on the NOAS04 nomenclature. For the 2013–2016 period, we used the 2004–2009 land-use statistics. For kilometre squares beyond the Swiss border (not including Liechtenstein), the values were derived from the CORINE Land Cover data 2012 (data collected in 2011–2012), maintaining consistency with the Swiss Land Use Statistics as far as possible. Five further variables (8–12 in the table) were calculated for 2013–2016 based on the Swiss Topographic Landscape Model or the VECTOR25 data set. For a few square kilometres beyond the Swiss border, the necessary data were digitised based on the current national maps at a scale of 1:25 000. Data concerning nitrogen deposition (13 in the table) were compiled by the company Meteotest for the Federal Office for the Environment; the values are those for 2010. Three variables (14–16 in the table) are based on the digital height model of Switzerland, which covers the surface of Switzerland with a grid at a resolution of 25 m. A value is indicated for each node (e.g. altitude above sea level). All nodes within a kilometre square are averaged to yield the value of the covariable for that square.

It was often not possible to estimate the influence of each variable on a species' density, especially since certain variables were quite strongly correlated. Therefore, variables were omitted in a step-by-step process for each species if their influence was difficult to estimate and, given the species' biology, no effect on density was expected. We used penalised 2D splines to model spatial autocorrelation. This was important because regional differences in a species' population density can occur despite suitable habitat. Including kilometre squares on the other side of the Swiss border allowed us to avoid boundary effects that can occur when using spline functions. Finally, density was set to zero in all kilometre squares above a certain species-specific altitude to avoid artefacts in high-mountain areas. Using the parameter values estimated by the model, we were able to estimate the number of territories for each of the approximately 41 000 kilometre squares in Switzerland and Liechtenstein. For most species, the estimates were unrealistically high in a small number of squares. We adjusted these values downwards for the visual representation. All values above a defined upper limit were reduced to the value of that upper limit. For five species, we set the upper limit manually (Grey Wagtail, Willow Warbler, Wood Warbler, House Sparrow, Italian Sparrow). For the other species, the 99.5 %-quantile of the estimates for all kilometre squares was defined as the upper limit. Squares with an estimated density of less than 0.05 were excluded prior to calculating the 99.5 %-quantile. To improve the visual presentation, the estimated values were slightly smoothed by interpolation.

Density change map

Of the 75 species with a density map for 2013-2016, 70 also feature a density change map showing the change in density from 1993–1996 to 2013–2016. While we were able to account for detection probability when estimating density for 2013–2016, the same was not possible with the 1993–1996 data. This is because the 1993–1996 data set only contained the final number of territories per species and kilometre square, but not the number of territories confirmed by an observation at each visit (hereafter «territory per visit»). Applying the same approach used to produce the density maps, we estimated the density for each kilometre square during 1993–1996 and during 2013–2016, only this time without accounting for detection probability. The data came from the territory mapping surveys in each period. In contrast to the model accounting for detection probability, where the number of detected «territories per visit» was the target variable, we now used the number of detected territories per kilometre square as the target variable.

The same 16 environmental variables used for the density maps served as explanatory variables. For the variables rivers and streams, lake shores, altitude above sea level, exposure and gradient we used the same values for both survey periods. Values for the remaining variables were determined separately for 1993–1996 based on the Swiss Land Use Statistics 1992–1997, the CORINE Land Cover data 1990 and the VECTOR25 data set. Data on nitrogen deposition was provided by Meteotest; they refer to 1990. We used the same set of environmental variables per species for 1993–1996 and 2013–2016. The analysis was performed using a generalised linear model instead of a binomial mixture model.

We also had to take into consideration that an upper limit had been set for territory counts in 1993–1996. For example, once more than ten territory-marking Common Chaffinches had been found, the observer did not have to continue recording this species. To be able to estimate the actual number of territories, we had to incorporate this so-called censoring into the model. To make the two sets of data comparable, the 2013–2016 survey data was censored to match the 1993–1996 approach. Some of the results obtained this way were unsatisfactory and did not correspond well with expert assessments or regional trends generated in the common breeding bird monitoring scheme (MHB). In these cases, we supplemented the censored data from the territory mapping surveys with uncensored data from the MHB scheme. The territory mapping data from 1993–1996 was supplemented with MHB data from the year 2000, the 2013–2016 data with MHB data collected in 2016. This allowed us to improve the estimates of density and density change for several species.

The density estimates for 1993–1996 and 2013–2016 were smoothed separately for presentation. For each kilometre square, we first averaged the density estimates of nine squares in a 3 × 3 km matrix with the corresponding square at its centre and used this value for the central kilometre square. We then calculated the difference between these smoothed values in the 1993–1996 and the 2013–2016 periods.

Distribution change map

We were able to generate distribution change maps for 21 species, most of them scarce. In order to produce an estimate of the species-specific occurrence probability per kilometre square for the 1993–1996 survey period as well, we used data collected in 1993–1996 in territory mapping surveys, complete species lists and individual observations, analysing them with a site-occupancy model.

For 608 of the 2934 territory maps produced in 1993–1996, the results of all survey visits were transferred from paper to digital form. For these areas, repeated records within the same season were now available, allowing us to estimate detection probability using a site-occupancy model. For a large portion of the 1993–1996 surveys, however, we only knew whether a species had been detected at some point during the three visits or not. In order to include this data type in the analysis, we expanded the site-occupancy model as follows: if a species had not been recorded, it was either absent or had remained undetected three times. This occurs at a probability of (1 – p)³, where p is a species’ detection probability during a visit, estimated based on the 608 digitised maps. Using the various types of data available to us, we were thus able to determine the occurrence probability per kilometre square for the 1993–1996 period.

As we had done for the density change maps, we smoothed the occurrence probability estimates using a 3 × 3 km matrix, then calculated the difference between the smoothed values for 1993–1996 and 2013–2016.

Distribution map based on density estimate

For widespread species with large territories, we often had enough records from territory mapping in the kilometre squares to produce a density map. However, the resulting density tended to be too high, because observers had often counted birds whose territories lay mostly beyond the limits of the kilometre square. For Black Kite, Eurasian Buzzard, Common Kestrel, Black Woodpecker and Common Raven as well as for the non-territorial Black Grouse, we therefore decided to convert the density estimates to estimates of occurrence probability using the formula v = 1 – e, where N corresponds to the estimated number of territories per kilometre square and v to the derived occurrence probability. For these six species, we depicted the distribution in 2013–2016 and the change since 1993–1996 based on the data thus transformed.

Distribution map 2013–2016

For 47 species, we created distribution maps that depict the probability of occurrence per kilometre square rather than the density. These species are comparatively rare, and the records collected during the territory mapping surveys in the kilometre squares were insufficient to produce density maps. Data was sourced from the territory mapping surveys as well as from the complete species lists and the individual observations recorded by the Ornithological Institute's volunteer collaborators (registered members of the Ornithological Information Service). As we did not model density, but only occurrence (presence/absence), we reduced the territory mapping data to presence/absence data (values >1 were set to 1). Thus, we only considered whether a species was detected during a visit or not. The complete species lists also gave us information about each species' presence/absence. For the individual observations, we constructed absences as described in Kéry et al.: for example, if a volunteer collaborator reported a Eurasian Sparrowhawk but no Middle Spotted Woodpecker, this was considered an absence record for the Middle Spotted Woodpecker. Only the species in the reporting categories A and B (scarce species) could be treated in this way; regarding species in category C (widespread species), we were confined to data from the territory mapping surveys and the complete species lists (for definitions of reporting categories, see Schmid et al.). Using these three types of data (territory mapping surveys, complete species lists and individual observations) allowed us to combine the advantages of each data set: the territory mapping surveys ensured that geographic distribution was adequate and covered the various habitats, while the species lists and individual observations greatly increased spatial coverage.

For seven wetland species (Little Grebe, Common Little Bittern, Western Water Rail, Common Moorhen, Common Grasshopper-warbler, Savi's Warbler and Great Reed-warbler) we also included data collected during the survey visits for the wetland monitoring scheme as presence/absence data. Because access to some of the most important wetlands is restricted in spring due to their status as protected sites, they are underrepresented in the other data sources.

The same 16 environmental variables used for the density maps again served as explanatory variables. As for the density maps, we pared down the set of environmental variables for certain species, modelled the spatial autocorrelation with penalised 2D splines, and set occurrence to zero in all kilometre squares above a species-specific altitude limit. We chose a site-occupancy model for the analysis. Like in the binomial mixture model, this approach allowed us to include detection probability. Detection probability was estimated separately for each data type (territory mapping surveys, complete species lists and individual observations). Like in the density maps, the spatial unit for modelling was one kilometre square. For each kilometre square in Switzerland, we thus estimated the probability of a species occurring there as a breeding bird.

In most cases, we depicted the occurrence probability estimated by the model. In some cases, we showed the «realised occurrence», which results from a combination of model estimate and raw data. In practice, it means that an area where a species was recorded was always shown as occupied, regardless of whether the model estimated a high occurrence probability or not. If an area had yielded no records, then realised occurrence corresponded to the occurrence probability estimated by the model. On the one hand, we used this approach for species like the Woodlark, where a large part of the population was actually found and recorded by observers. On the other hand, we applied it to various wetland species, because many wetlands are frequently visited by ornithologists, or because the wetland monitoring scheme provides a good level of coverage and the likelihood of detecting the species present is therefore high.

Similar to the density maps, the values were slightly smoothed for presentation.

Point map 2013–2016

The distribution of 77 rare species and colony breeders was represented using a point map. The underlying data set was made up of the results of territory mapping surveys in the kilometre squares, complete species lists and individual observations as well as special surveys (e.g. territory mapping surveys for the wetland monitoring scheme, data from the Rook monitoring programme). Based on this data, we calculated the number of territories or breeding pairs per site and year. A point on the map corresponds to the number of territories or breeding pairs per spatial unit, averaged across the 2013–2016 period. For colonial species, the spatial unit is a colony. In settlements, we generally summarised the data for a whole town or district. In the case of waterbirds, the spatial unit was a lake or pond or a stretch along the shore of a river or large lake. For very rare species, we plotted a point if the species was recorded at a site with a sufficiently high Atlas Code in at least one year during 2013–2016.

Point change map

We produced point change maps for 21 species showing the change in the number of territories or breeding pairs per spatial unit. This was only possible for those species for which we knew the location and average number of breeding pairs in each colony or breeding site for the 1993–1996 survey period. The location of colonies or breeding sites in 1993–1996 was often uncertain. To account for this uncertainty and to improve legibility, sites that were very close together were combined into one point on the point change map. Certain sites depicted on the point map are therefore not displayed separately on the point change map.

Altitude chart

Besides geographical distribution, the 2013–2016 atlas also illustrates each species' distribution along the altitudinal gradient. To this effect, we divided the data shown on the maps into 100-m altitude bands. We summed up the estimated occurrence information per level and divided it by the total for the whole of Switzerland. This gave us the proportion of the total population per altitude level, which we then plotted on the altitude charts.

For many species, we were also able to display the change in altitudinal distribution. For each altitude level, we calculated the difference between the figures in 1993–1996 and in 2013–2016. We then divided this figure by the Swiss total in 2013–2016. This gave us the relative change in numbers per altitude level as a proportion of the total for the 2013–2016 period. The result was displayed on the right hand side of the altitude charts. In this approach, the length of the bars on the two graphs is directly comparable. The altitudinal distribution for 1993–1996 can be elicited simply by subtracting the gains from the current status or by adding the losses.

The change in altitudinal distribution could only be shown for those species for which there was a large amount of data from 1993–1996. Where this was not the case, we provide only the current altitudinal distribution.

Validation

The maps and altitude charts were reviewed by various species experts. If the result was not satisfactory, we tried to improve the estimates by adjusting the selection of environmental variables, giving special consideration to the species' ecological requirements. If no improvement was achieved, we took a more conservative approach. Where it was not possible to produce a satisfactory density map, we generated a distribution map. If the distribution map was inadequate, we used the raw data to represent the distribution at atlas-square resolution. We used this approach to represent the current status in 2013–2016 and the change since 1993–1996.

The models were often unable to correctly estimate a species' upper distribution limit, as the surface area – and with it the amount of available data – decreases rapidly with increasing altitude. Incorrect estimates were particularly apparent in the altitude chart. In such cases we manually set the density or the occurrence probability to zero beyond a certain altitude limit.

Software

We processed the data in R 3.3.2. For statistical analysis, we used JAGS 4.2.x. Data for the point maps were processed in QGIS.