Community assembly in Nothobranchius annual fishes: Nested patterns, environmental niche and biogeographic history

Abstract The assembly of local communities from regional species pools is shaped by historical aspects of distribution, environmental conditions, and biotic interactions. We studied local community assembly patterns in African annual killifishes of the genus Nothobranchius (Cyprinodontiformes), investigating data from 168 communities across the entire range of regionally co‐existing species. Nothobranchius are small fishes associated with annually desiccating pools. We detected a nested pattern of local communities in one region (Southern Mozambique, with Nothobranchius furzeri as the core and dominant species), but no nestedness was found in the second region (Central Mozambique, with Nothobranchius orthonotus being the dominant species). A checkerboard pattern of local Nothobranchius community assembly was demonstrated in both regions. Multivariate environmental niche modeling revealed moderate differences in environmental niche occupancy between three monophyletic clades that largely co‐occurred geographically and greater differences between strictly allopatric species within the clades. Most variation among species was observed along an altitudinal gradient; N. furzeri and Nothobranchius kadleci were absent from coastal plains, Nothobranchius pienaari, Nothobranchius rachovii, and Nothobranchius krysanovi were associated with lower altitude and N. orthonotus was intermediate and geographically most widespread species. We discuss implications for ecological and evolutionary research in this taxon.


| INTRODUCTION
The principles of biological community assembly are of central importance to ecologists (Connor, Collins, & Simberloff, 2013;Connor & Simberloff, 1979;Diamond, 1975;Gaston et al., 2000;Götzenberger et al., 2012;Henriques-Silva, Lindo, & Peres-Neto, 2013;Holyoak, Leibold, & Holt, 2005;Hubbell, 2001;Parr, 2008;Ulrich & Gotelli, 2007). The species occurrence is contingent upon its evolutionary and biogeographic history (historical aspects of distribution), environmental characteristics (habitat suitability) and the presence of other species (biotic interactions). The composition of local communities is restricted by diversity of the regional species pool, setting the upper limit to their alpha diversity, and is shaped by a set of "assembly rules" related to the ecological niche and dispersal abilities of potentially co-occurring species (i.e., meta-community) (Diamond, 1975;Gaston et al., 2000;Henriques-Silva et al., 2013;Hubbell, 2001). The assembly patterns identify the non-random component in local community composition by comparing certain parameters of an observed dataset with the same parameters in multiple randomized datasets (Connor & Simberloff, 1979;Götzenberger et al., 2012).
Nestedness and negative species co-occurrence are two main patterns of local community assembly. The nested pattern describes the positive co-occurrence of certain species, with species-poor communities being predictable subsets of species-richer communities.
Negative species co-occurrence tends to produce checkerboard patterns where some species combinations are absent, while the distribution of other species can still be positively associated (Diamond, 1975). Historical, biotic and abiotic factors can all drive nested and checkerboard patterns (Connor et al., 2013;Ulrich & Gotelli, 2007).
To understand which processes underpin the observed distributional patterns, examination of biotic and abiotic species-occurrence associations is necessary.
The concept of ecological niche, describing species' response to the distributions of resources and potentially co-existing taxa, serves as a basis for assessing the ecological similarities and differences among species (Barve et al., 2011). The ecological niche characterizes the conditions for species to persist and the degree of niche overlap between potentially co-existing species determines the pattern of their local co-occurrence, including environmental aspects of the ecological niche. However, checkerboard patterns can arise for reasons other than inter-specific competition. Difference in habitat preferences, recent geographic speciation, and constraints on dispersal can all lead to negative co-occurrence between species (Connor et al., 2013).
Meta-communities (sets of local communities linked by potential dispersal) and their local subsets vary widely in species richness.
Species-rich communities contributed most to our current understanding of the assembly of local communities (Holyoak et al., 2005).
Species-poor communities in well-defined habitat patches are more amenable to accurate census. However, very low species diversity leads to low statistical power of some standard analytical procedures that have been developed for species-richer communities (Barve et al., 2011;Götzenberger et al., 2012;Henriques-Silva et al., 2013). For taxonomically defined, species-poor communities, the assembly rules can be easily dissected to specific combinations of species and directly related to particular environmental factors.
We studied local community assembly patterns in a clade of African annual fishes. Annual killifish of the genus Nothobranchius (Cyprinodontiformes, Nothobranchiidae) are small fishes that are tightly associated with ephemeral pools in the East African savannah (Wildekamp, 2004). The pools only exist for a short period during and after the rainy season, and annual killifishes are the only teleost fishes adapted to these temporary habitats. Every year, fish hatch from desiccation-resistant eggs when the pools fill with rainwater (Polačik, Donner, & Reichard, 2011) and complete their life cycle during 1-11 months Reichard, 2015), with eggs surviving the dry period buried in the sediment (Cellerino, Valenzano, & Reichard, 2016;Reichard, Cellerino, & Valenzano, 2015).
Several Nothobranchius species often inhabit the same pool (Reichard, 2015;Reichard, Polačik, & Sedláček, 2009;Watters, 2009;Wildekamp, 2004). When Nothobranchius species co-exist in a pool, limited space does not permit any substantial microhabitat separation. The diet of co-occurring species is largely overlapping (Polačik, Harrod, Blažek, & Reichard, 2014;; there is no apparent trophic specialization. Other teleost fishes only invade temporary pools in the event of substantial flooding, except those in the flat floodplains of the lower reaches of major rivers, where co-existence with non-annual fishes is more frequent. In such cases, Nothobranchius fishes are constrained to particular microhabitats within a floodplain wetland matrix, not being able to sustain stable populations in direct competition with other teleost fishes (Reichard, 2015). Other vertebrates in temporary pools are Protopterus spp. lungfishes and tadpoles (Reichard, Polačik, Blažek, & Vrtílek, 2014), though their ecological relationship to Nothobranchius fishes is not known.
Phylogeographic analysis of all three complexes, based on nuclear and mitochondrial genetic markers, has revealed strong spatial structure of populations, indicative of very low dispersal ability (Bartáková et al., 2013(Bartáková et al., , 2015. All study species are illustrated on Figure 1.
In the present study, we concentrated on understanding the finerscale processes that have led to particular local community assembly.
Specifically, we asked the following questions. First, is the distribution of the three clades geographically congruent and is the distribution of species within each clade allopatric? Second, are local Nothobranchius communities composed of a predictable subset of regionally present species and is the assembly of local communities driven by positive or negative associations among species? Third, is the presence of particular species in local communities explained by local environmental factors?
We used presence/absence data from nine expeditions to the region (between 2008 and 2015) to reconstruct the distribution of each clade and species, and the patterns of their co-occurrence and nestedness. We then compared the relative abundance of co-occurring species from a subset of pools and investigated the factors responsible for the observed patterns by analyzing species-habitat associations.
We used multivariate environmental niche modeling (Broennimann et al., 2012;Warren, Glor, & Turelli, 2008) to compare niche overlap between the clades and between species within the clades. We then used univariate linear models to identify the key local habitat characteristics important for the presence of particular species. Nothobranchius presence). Repeated sampling of the same sites usually recovered the same Nothobranchius communities within particular pools, though in 20 cases, repeated sampling over years revealed an additional species. All suitable sites were targeted, but we generally omitted pools in close proximity to each other (and likely connected during flooding). Pools directly connected to large rivers were not sampled on the basis of previous knowledge (Wildekamp, 2004) that they do not include Nothobranchius fishes. Stagnant pools within the riverbed of smaller temporary streams were included in the sampling design, as they often conform to annual fish habitat characteristics (Reichard, 2015). Over the years, we aimed to sample all accessible regions in Southern and Central Mozambique and the adjacent part of Zimbabwe that were deemed likely to contain Nothobranchius fishes.

| Fish collection
To collect fish, a 5-mm-mesh dip net with a triangular metal frame (45 × 45 cm) on a 1.5-m wooden pole was used in most sites.
Typically, 15-40 hauls were performed at each site but more hauls were performed when fish density was low. Fewer hauls were taken if the site was too small to accommodate 15 hauls and more hauls were sometimes taken when we aimed to collect a large number of individuals for a different purpose (e.g., population genetic analysis).
In larger sites, dip net sampling was supplemented by sampling using a seine net (length 2.7 m, depth 0.7 m, mesh size 4 mm), especially when we aimed for a larger sample size. The mesh size used retained adult Nothobranchius unselectively and there was no apparent speciesspecific bias in the probability of capture. Sampling effort variation within and among years had a negligible effect on species detection.
Fish were unambiguously identified to species in the field. Most individuals were returned to the pools, but in some cases a subsample was taken to the laboratory for further research.

| Habitat data
We estimated a set of habitat variables for each site. GPS reading (longitude, latitude) was taken during each collection, along with altitude (to the nearest 1 m) using Garmin nüvi 550 (Garmin Ltd.).

F I G U R E 1 Adult males of all study species
Over the first 5 years of data collection (2008-2012), a set of habitat variables was taken according to the sampling protocol. variables that were collected (e.g., maximum water depth, water temperature and water conductivity) varied considerably both diurnally and seasonally and were not useful in predicting the presence of particular species. Further details on habitat data collection are given in Reichard et al. (2009Reichard et al. ( , 2014.
Given that habitat data were collected over separate years and were subject to intra-seasonal variability, we averaged continuous and

| Analysis of co-occurrence
A pairwise co-occurrence between species with overlapping ranges was quantified using the probabilistic model of species co-occurrence on presence-absence data (Griffith, Veech, & Marsh, 2015). The analysis compared the observed and expected frequencies of cooccurrence between the pair of species. The expected frequency is based on a random distribution of each species, independent of the other species. Statistically significant over-and under-representations of particular pairwise associations were interpreted as positive and negative co-occurrence. We quantified relative species abundance as the relative proportion of individuals of a particular species in a Nothobranchius community for a given site, calculating the relative abundance of each species for each combination of species separately. The co-occurrence analysis was conducted in R 2.14.2 (R Foundation for Statistical Computing, Vienna, Austria), using the cooccur library (Griffith et al., 2015). We used Venn diagrams to visualize proportional species co-occurrence.

| Nestedness analysis
Nestedness was analyzed using the number of decreasing fills (NODF) (Almeida-Neto, Guimarães, Guimarães, Loyola, & Ulrich, 2008) as the nestedness metric. This index can range from 0 (completely random) to 100 (perfectly nested). Statistical significance was tested by comparing the NODF values with 999 random assemblages created according a null model. Choosing the appropriate null model is the most controversial aspect of statistical inference in nestedness analysis (Ulrich, Almeida-Neto, & Gotelli, 2009). The equiprobable null model (that preserves total number of species occurrences in the original matrix, but does not preserve margin totals) is truly random, but criticized for being prone to type I errors (Ulrich et al., 2009;Wright, Patterson, Mikkelson, Cutler, & Atmar, 1998). Fixed-fixed null models (that preserve both margin totals) are suggested as an alternative, but often fail when matrices are extremely nested, because there are too few possible matrix rearrangements (Almeida-Neto et al., 2008). This limitation is especially apparent in small-size matrices, as in our dataset. As a compromise, we chose a null model that constrained only column totals (proportional to species relative abundances; Randnest, or c0) (Jonsson, 2001), and allowed more row rearrangements in our species-poor matrix. We did not use matrix temperature to estimate nestedness (Atmar & Patterson, 1993), because this method has been criticized for being prone to type I errors (Ulrich et al., 2009;Wright et al., 1998). The nestedness analysis was conducted using the vegan library (Oksanen et al., 2012).

| Multivariate niche modeling
Niche overlap between clades and between species within clades was measured using the niche similarity test (Warren et al., 2008). This analysis compares realized niche overlap against a series of overlaps calculated from randomized datasets to examine whether similarity between the niches of two taxa is different from that expected by chance (Warren et al., 2008). Using this test, it is also possible to compare environmental niches of species with non-overlapping geographic ranges (Broennimann et al., 2012). The latter approach, which we have undertaken, measures niche overlap along two-dimensional space determined by multivariate analysis gradients, weighting densities of species occurrences by densities of environmental factors along these gradients. Analogically, we used the niche equivalency test that examines "niche conservatism in the strictest sense" (Broennimann et al., 2012;Warren et al., 2008). The null hypothesis for this analysis postulates that niches occupied by two entities are identical.
Environmental gradients were constructed using principle component analysis (PCA

| Univariate analyses of environmental predictors of species presence
The same set of predictors as in the PCA was used to test the univariate response of particular species presence to each environmental variable. Only sites from within the geographic range of particular species were selected for analysis. Within the R-clade, three sites without Nothobranchius were deleted because it cannot be determined whether they belonged to the range of N. pienaari or N. rachovii.
Species presence was modeled using a Bernoulli generalized linear model (GLM) with a polynomial effect for altitude (to model constraints at both extremes of the gradient) and linear effects for substrate and turbidity, and pool size, littoral, submersed and Nymphaea vegetation that were all converted into three categories. Before applying statistical models a data exploration was undertaken, following the protocol described in Ieno and Zuur (2015). The data were examined for outliers in the response and explanatory variables, homogeneity and zero inflation in the response variable, collinearity between explanatory variables, and the nature of relationships between the response and explanatory variables. We present full models for four study species but backward selection of the minimal adequate model produced concordant results. The number of sites sampled in the N. rachovii and N. krysanovi ranges was too low to provide reliable estimates from GLM analysis.

| Overview and geographic patterns
Overall, 374 pools were sampled and 168 (45%) contained at least one Nothobranchius orthonotus was also the numerically dominant species and, at sites where all three species co-existed, N. pienaari and N. kadleci occurred at similar abundances. In two-species communities with either N. kadleci or N. pienaari, quantitative dominance of N. orthonotus was relatively smaller than in three-species communities (Figure 4b).  Figure 5d).

| Ecological niche modeling
At the level of clades, the niche equivalency and niche similarity tests detected a difference in environmental niche between the F-and R-clades. The O-clade niche was similar to both the R-and F-clade niches using either test (Table 2a) (Table 1c). There were too few N. krysanovi populations to permit analysis for this species, despite a superficial similarity between the N. rachovii and N. krysanovi niches (Figure 5g-i).

| Nothobranchius furzeri
Detailed habitat data were available from 154 pools within the potential range of this species, with 68 N. furzeri populations found (56% of F I G U R E 2 Distribution maps of three study clades and species co-existence. The presence of (a) Nothobranchius furzeri (yellow) an

Nothobranchius kadleci (dark red), (b) Nothobranchius orthonotus (violet), and (c) Nothobranchius pienaari (black), Nothobranchius rachovii (blue) and
Nothobranchius krysanovi (red) populations at particular sites. Small empty points represent sampled sites where population of a given clade was absent. (d) A combined dataset illustrating geographic aspect of species co-occurrence in local Nothobranchius communities, where sites with single, two and three co-existing species are illustrated. An inset details the region with a parapatric distribution of N. rachovii and N. pienaari. The maps were created in R environment (R Core Team 2015) [packages sp (Bivand, Pebesma, & Gomez-Rubio, 2013), rgdal (Bivand, Keitt, & Rowlingson, 2016), maps (Becker, Wilks, Brownrigg, Minka, & Deckmyn, 2016) and GIStools (Brunsdon & Chen, 2014)]. The altitudinal gradient for Mozambique was downloaded from http://www.diva-gis.org/gdata F I G U R E 3 Quantification of species co-existence in local Nothobranchius communities. Venn diagrams visualizing proportional species co-occurrence (a) within Nothobranchius furzeri range, (b) within Nothobranchius kadleci range and (c) outside the range of F-clade. Diagrams were constructed using Venny 2.0.2 (Oliveros, 2007(Oliveros, -2015 F I G U R E 4 Relative species abundance. Quantitative estimates of relative species abundance were calculated as relative proportion of individuals of a particular species in a Nothobranchius community for a given site. Co-existence between (a) Nothobranchius furzeri and its congeners, (b) Nothobranchius kadleci and congeners, and (c) Nothobranchius orthonotus and Nothobranchius pienaari are illustrated. Mean values with 95% confidence intervals were calculated from a set of all separate species ratios. The number of communities used for quantitative estimates is given for each combination F I G U R E 5 Niche similarity among and within clades. Two-dimensional visualization of niches along environmental gradient detected by the principle component analysis (PCA). (a-c) clade-specific niches with (d) correlation of environmental predictors with PCA axes, and species-specific niches for (e) Nothobranchius furzeri, (f) Nothobranchius kadleci, (g) Nothobranchius pienaari, (h) Nothobranchius rachovii and (i) Nothobranchius krysanovi. Gray shading indicates density of occurrence along two PCA gradients. The solid and dashed contours indicate 100% and 50% of available environment, respectively the pools). The presence of N. furzeri was most significantly related to altitude (Table 3), with populations absent from coastal areas (no occurrence below 24 m a.s.l.) and rarely recorded at high altitude ( Figure 6). Nothobranchius furzeri was more common at sites with a soft substrate and positively associated with littoral, submerged and Nymphaea vegetation (Table 3).

| Nothobranchius kadleci
Twenty-one N. kadleci populations were found across 78 pools investigated within the potential range of this species (27%). As with N. furzeri, the N. kadleci was absent from coastal areas (a single population present at 16 m a.s.l., others at 26 m a.s.l. and higher) and highaltitude sites ( Figure 6). Nothobranchius kadleci were also more likely to occur at sites with higher water turbidity (Table 3).

| Nothobranchius orthonotus
Nothobranchius orthonotus populations were found in 93 of 258 pools sampled within its range (36%). This species was present in coastal plains (lowest recorded altitude 2 m a.s.l.) and the probability of its presence declined above altitudes of approximately 100 m a.s.l.
( Figure 6). This species was positively associated with pools that contained more littoral vegetation (Table 3).

| Nothobranchius pienaari
The species was found in 45 of 202 pools sampled within its range (22%). There was a negative linear relationship between N. pienaari presence and altitude (Figure 6), though N. pienaari was found as high as 181 m a.s.l. N. pienaari populations were positively associated with T A B L E 2 Observed niche overlaps and their statistical significance. The number of pools used to construct the niche for each taxon (N) and values of Schoener's D between environmental niches occupied by Nothobranchius clades and species within the clades (D) Statistically significant results are highlighted in bold. Altitudinal effect is visualized in Figure 6. All significant relationships between Nothobranchius population presence and vegetation are positive.

N D
sites that contained littoral and Nymphaea vegetation, and tended to have more solid substrate (Table 2).

| DISCUSSION
Interactions between local and regional processes are central in shaping local community assembly (Barve et al., 2011;Connor & Simberloff, 1979;Connor et al., 2013;Diamond, 1975;Gaston et al., 2000;Götzenberger et al., 2012;Henriques-Silva et al., 2013;Hubbell, 2001;Ulrich & Gotelli, 2007). The diversity of the regional species pool sets the upper limit to potential alpha diversity and the ability of species to disperse within the region enables local communities to diverge.

Different local environmental conditions and species interactions then
shape the particular composition of local communities. Previous research has documented a strong role of historical vicariance on largescale distributional patterns of Nothobranchius fishes. The genus is separated into four major allopatric lineages (Dorn, Musilová, Platzer, Reichwald, & Cellerino, 2014), of which the Southern lineage, consisting a total of six species assigned to three monophyletic clades, was investigated in the present study throughout its entire distribution.
Phylogeographic analysis of all three complexes, based on nuclear and mitochondrial genetic markers, has revealed strong spatial structure of populations, indicative of very low dispersal ability (Bartáková et al., 2013(Bartáková et al., , 2015. In the present study, an altitudinal gradient was demonstrated to play a role in geographic distribution of the main clades (and, consequently, of species within the clades). The F-clade, with N. furzeri and N. kadleci, was absent from coastal plains, while the R-clade (N. pienaari, N. rachovii, N. krysanovi) was underrepresented (though not absent) at relatively higher altitude. This pattern may be related to historical dispersal or environmental conditions. We think that a major role for historical dispersal (or a lack of it) is unlikely. While the source, We suggest a major role for local conditions on the lack of F-clade species in coastal plains. Importantly, the altitudinal gradient in the study area coincides with a gradient of aridity where high-altitude sites are significantly drier, with lower total rainfall, higher evapotranspiration rates, and lower predictability of rains (Terzibasi Tozzini et al., 2013;Vrtílek & Reichard, 2016b). The coastal plains are more humid, with temporary pools inundated for a longer period Terzibasi Tozzini et al., 2013). The flat topography of the coastal region likely connects isolated habitats with permanent water bodies more frequently than pools at a higher altitude. Hence, coastal plain sites are more likely to be invaded by non-annual fishes (Reichard, 2015) and N. furzeri and N. kadleci may be more sensitive to the presence of non-annual fishes than the other Nothobranchius species. The exclusion of F-clade species from coastal plains cannot be directly associated with increased salinity; water conductivity in most coastal plain pools is within the range recorded from higher-altitude sites (Reichard et al., 2009)  . Nevertheless, the more humid coastal region receives more seasonal and non-seasonal rainfall Terzibasi Tozzini et al., 2013). It is possible that the eggs of N. furzeri and N. kadleci are more prone to non-seasonal rains, triggering their hatching outside the appropriate season and compromising population persistence. This possibility is amenable to experimental testing.
The upper limit of distribution of Nothobranchius species is constrained by altitude; the presence of all species was rare above 150 m a.s.l. This lack of Nothobranchius populations is apparently related to at an altitude of almost 400 m a.s.l. (Watters, Wildekamp, & Cooper, 2009) Given extensive research effort in the area, it is likely that these represent isolated populations (Shidlovskiy et al., 2010), located at the border with Mozambique. It remains to be investigated whether highaltitude sites represent ancient populations remaining from a period of more suitable environmental conditions (i.e., in regions currently too dry to sustain more viable populations that were frequent in that region in the past), ancient populations that acted as sources to colonize lower-altitude sites, or recent colonization events.
We demonstrated that local Nothobranchius communities showed the most frequent species are also the most numerous. This positive abundance-occupancy relationship is frequently reported from other communities (Gaston, 2003 was visited repeatedly and the qualitative community composition was stable across years. Nothobranchius orthonotus absence at this site is, therefore, unlikely related to an inadequate sampling. It is interesting that the species negatively associated in two-species communities were frequently co-existing in three-species communities. We suggest that the observed checkerboard pattern is related to environmental characteristics rather than biotic interactions between Nothobranchius species. F-clade and R-clade species were associated with most dissimilar positions on environmental niche gradients, especially along the vegetation-turbidity axis, and their environmental niches overlap in the segment of N. orthonotus preference ( Figure 5). The range of environmental characteristics suitable for co-existence of F-clade and R-clade species was, therefore, always suitable also for N. orthonotus.
The absence of N. orthonotus is some of these pools was likely a result of stochastic local extinctions (or a lack of colonization); species from all three clades co-existed frequently where habitat conditions were favorable for all three species (Figure 3a,b). Hence, the rarity of local communities with exclusive co-existence of F-clade and R-clade species was apparently mediated by environmental conditions rather than inter-specific incompatibilities. In conclusion, we suggest that environmental conditions rather than biotic interactions likely shaped the assembly of local Nothobranchius communities. We acknowledge, however, that regional species pool in our study was small, limiting inferences from our analyses on the nestedness and checkerboard pattern to our particular study system.
Within clades, no sympatric co-existence of sister species was found. Nothobranchius furzeri and N. kadleci are known to produce viable and fertile hybrids in the laboratory (Ng'oma, Groth, Ripa, Platzer, . It is currently assumed that allopatric speciation is the single driver of Nothobranchius diversification (Dorn et al., 2014). However, Nothobranchius genomes are evolving at a fast rate (Reichwald et al., 2015;Valenzano et al., 2015), with intra-specific differentiation in the sex determination system in N. furzeri . Divergence in sex determination systems and sex-linked genes is known to drive speciation across several taxa (Kitano et al., 2009;Qvarnström & Bailey, 2009). Interestingly, species of the Rclade differ cytogenetically (Shidlovskiy et al., 2010), and it is possible that a combination of allopatric diversification and major cytogenetic incompatibilities during secondary sympatry may promote speciation at least in some Nothobranchius lineages. The parapatric distribution between N. pienaari and N. rachovii, where their populations were located within a single coastal swamp system as close as 5 km apart, is generally consistent with this scenario. If post-mating barriers are indeed fully developed in this species pair (Shidlovskiy et al., 2010), it remains to be investigated whether a lack of premating barriers or ecological exclusion interferes with the local co-existence of the two species. Nothobranchius species are certainly a promising system in which to study prezygotic, postzygotic and ecological barriers to speciation (Ng'oma et al., 2014;Sedláček, Baciaková, & Kratochvíl, 2014;Valenzano et al., 2015).
Other lineages of annual killifishes are distributed in the Neotropics.
While the distribution of all annual killifishes is restricted to ephemeral pools (Polačik & Podrabsky, 2015), such pools vary in their connectivity with permanent water systems. South American annual fishes co-exist with a diverse community of non-annual fishes (Lanés et al., 2016;Loureiro et al., 2015;Nico & Thomerson, 1989). Unlike in the African annual fishes, particularly good understanding of community assembly exists for communities of annual Austrolebias fishes in Uruguay (García et al., 2009;Loureiro et al., 2015). Neotropical Austrolebias species are of a comparable size to the African Nothobranchius and co-exist in communities comprising up to five species (Laufer et al., 2009;Loureiro et al., 2015). The distribution of individual species is best explained by past allopatric fragmentations and subsequent range expansion involving secondary contacts (García et al., 2009), likely driven by recent (Holocene) repeated marine transgressions (Sprechman, 1978).
Stable species co-existence was suggested to be mediated by interspecific differences in developmental time (Laufer et al., 2009), morphological specializations to diversified feeding niches (Costa, 2009;Keppeler et al., 2015) and spatial segregation of preferred microhabitats (Loureiro et al., 2015). At a regional scale, significant checkerboard pattern with a strict lack of co-occurrence between certain species pairs was detected (Loureiro et al., 2015). Austrolebias species inhabit pools closely related to the active stream alluvia (Loureiro et al., 2015), making it easier to disperse at least in the downstream direction. Yet, there are strict regional limits to species occurrences (Loureiro et al., 2015;Volcan, Gonçalves, Lanés, & Guadagnin, 2015), suggesting that competitive interactions restrict potential range expansions and coexistence of related Austrolebias species. The ranges of individual species vary from extensive to locally endemic (Loureiro et al., 2015;Volcan et al., 2015), being analogous to the range size variation seen in Nothobranchius (Reichard, 2015). Overall, historical and local effects interact to shape Austrolebias community assembly and abrupt limits of species distributions in a similar way to the African Nothobranchius species.
In conclusion, we found large co-occurrence of all three clades across Southern and Central Mozambique, with F-clade being absent in the region north of the Zambezi River. We demonstrated negative and positive co-occurrence of particular Nothobranchius species, including a checkerboard pattern with the lack of exclusive co-existence of the F-and R-clades species in two-species communities. We argue that, at the level of clades of closely related species, co-occurrence patterns arise from specific associations with environmental characteristics rather than from biotic interactions and environmental factors associated with altitudinal and precipitation gradients apparently drive geographic distribution. The species with most widespread populations had also the highest relative abundance in local communities. Sister species never co-existed. This suggests a major role for geographic isolation in speciation and potential for competitive exclusion in co-existence of sister species. Given that other aspects of Nothobranchius ecology are relatively well known (Cellerino et al., 2016), this group of fish presents a useful model for evolutionary and ecological research.