Effects of the Temporal Scale of Observation on the Analysis of Aquatic Invertebrate Metacommunities

The development of metacommunity theory has boosted the implementation of numerous empirical tests with field data, mostly focused on the role of spatial and environmental gradients on metacommunity organization. These studies showed an important dependence of the results on the observational scale considered, i.e., spatial grain, sampling spacing, and extent. However, few works deal with time per se as a component explaining metacommunity structure, even when data from periodic sampling are available. We suggest adding time explicitly to metacommunity analysis, but taking into account that the temporal scale of observation could affect the estimation of the relative influence of environment, space, and time, as previously recorded for spatial scale variation. Here, we analyze temporal scale dependence using simulated and empirical metacommunities of aquatic invertebrates. The effects of the study duration (i.e., temporal extent) were stronger when most metacommunity variation occurred along the temporal axis, so that local communities were spatially homogenized under high dispersal rates. Contrarily, dispersal limitation and niche differentiation (depending on the spatio-temporal structure of the environment) kept constant the spatial heterogeneity of the metacommunity, reducing the temporal variation and the importance of the temporal scale of observation. Our results highlight the importance of the temporal scale chosen for the analysis of metacommunity dynamics and emphasizes the temporal perspective of metacommunities, suggesting novel and interesting avenues in this research program.


INTRODUCTION
Ecological communities are assembled by a complex interaction of processes, such as niche-related selection, dispersal, and ecological drift (Vellend, 2010;Leibold and Chase, 2017). Our perception of these processes strongly depends on the observational scale considered (Soininen et al., 2011;Viana and Chase, 2019), but most studies have treated this issue from a spatial point of view, suggesting a similar behavior for temporal scales. Nevertheless, the effects of variable temporal scales remain largely unknown (Korhonen et al., 2010;Tomašových and Kidwell, 2010;Dornelas et al., 2014).
The relative importance of selection, dispersal and drift vary with the spatial extent of the study following a generalized conceptual model (Leibold et al., 2004;Leibold and Chase, 2017). At a small spatial extent, environmental conditions can be mostly homogeneous. Then, all the species in the metacommunity may have similar environmental niches and only stochastic dynamics (i.e., ecological drift) foster metacommunity variation (Neutral Theory archetype, NT; Hubbell, 2001). However, despite the small spatial extent, localities could be environmentally heterogeneous, sorting species with different niches. Then, high dispersal (facilitated by short distances among sites) could maintain populations at suboptimal environmental conditions, hindering the deterministic responses of the species to the environment (Mass Effects archetype, ME; Mouquet and Loreau, 2003). At an intermediate spatial scale, environmental heterogeneity can increase, and niche filtering then originates spatial differences in community composition, while dispersal rates still allow the movement of species to reach (without surplus) their potentially suitable localities (Species Sorting archetype, SS; Chase and Leibold, 2003). At a larger spatial extent, dispersal limitation may spatially restrict the distribution of species, creating spatial dissimilarities independently of the environmental conditions and species' niches (Patch Dynamics archetype, PD; Leibold et al., 2004), and biogeographic effects.
Under a temporal focus, local temporal turnover can increase with the study duration (i.e., temporal extent; Wu and Li, 2006) just because of stochastic fluctuations, as NT and PD (even ME) assume (Leibold and Chase, 2017). However, environmental changes can also determine the temporal turnover due to species sorting, potentially reaching rates higher than those predicted by stochastic dynamics (Dornelas et al., 2014). At the metacommunity level, long-term studies have shown either stable (e.g., Azeria and Kolasa, 2008;Huttunen et al., 2018;Lindholm et al., 2020) or changing metacommunities through time, due to environmental fluctuations or disturbances observed (e.g., Datry et al., 2016;Sarremejane et al., 2017;Cañedo-Argüelles et al., 2020). Nevertheless, the spatio-temporal scale of observation in each study could explain these differences, depending on whether or not the study duration covered those relevant environmental fluctuations (Korhonen et al., 2010).
Previous studies estimated the relative relevance of different ecological processes by means of partitioning the metacommunity variation in species composition between the effects of environmental and spatial variables (e.g., Cottenie, 2005). With this approach, we assume that the fraction of the species variation explained by environmental variables is associated with selection (SS) and the rest of the variation (included pure space and unexplained fractions) is due to neutral dispersal and drift (NT and PD, even ME). Spatial variables can capture any spatial variation in the metacommunity, although we usually choose those spatial variables that are associated with broad-scale patterns . These spatial patterns may be generated by dispersal limitation (i.e., poor connected localities follow different stochastic dynamics; pure spatial fraction) and/or be associated with spatial patterns of the environment (the shared fraction between environment and space). Some studies have also considered introducing temporal variables in variation partitioning analyses (e.g., Anderson and Cribble, 1998;Muylaert et al., 2000;Padial et al., 2014). As with spatial variables, the temporal variables can explain either community changes related to (stochastic) fluctuations and biological cycles or those affected by a temporally-structured environment (the shared fraction between environment and time). Considering the relevance of spatio-temporal variables in these analytical methods is crucial, since metacommunities are to be considered as dynamic structures, with both spatial and temporal variation.
Here, we used simulated metacommunities to analyze the effects of study duration on the relative importance of environmental, spatial, and temporal variables. We expected that (1) the effects of the variation in study duration should be relevant when the metacommunity varies along the temporal axis (high temporal variation; Figures 1A,B). This could be a consequence of high dispersal rates, diluting the spatial variation through time. In this context, local communities may experience a process of mixing by the high number of colonizers, temporally synchronizing the metacommunity . On the other hand, we envisage that (2) the study duration should be less relevant when the spatial structure of the metacommunity is maintained through time (low temporal variation; Figures 1C,D). Dispersal limitation can create permanent differences among localities and therefore spatial patterns would acquire more relevance (with high spatial variation; Figure 1D). Similarly, environmental conditions could also create spatial or temporal patterns in those metacommunities depending on the environmental variation, increasing the importance of the temporal scale when the environmental conditions vary through time. These two predictions can be considered two extremes of a continuum between unique temporal or spatial variation. Additionally, we analyze empirical data sets of aquatic invertebrates from a temperate setting to compare our expectations as resulting from the simulations with the field data.

Model for Simulated Metacommunities
The objective of our simulations was to emulate the temporal dynamics of metacommunities constrained by different conditions of selection and dispersal (Figure 2A). To this aim, we obtained species matrices (i.e., locality × species) at several points in time for each simulated metacommunity. Then, we sampled the species matrices changing the study duration but keeping constant the number of temporal sampling points. Consequently, we increased the time lag between sampling events as we extended the study duration (Figures 2B,C).
We used a model based on previous metacommunity simulations (Gravel et al., 2006;Sokol et al., 2017;Thompson et al., 2020). These models allow the observation of simulated metacommunity dynamics along different gradients of selection and dispersal. The model used had two parts differentiated in the simulation routine. First (Equation 1), species disperse among localities and then (Equation 2), the species compete at each locality depending on their relative abundances and the environmental conditions, determining the species composition for the next round. In the Equations 1 and 2, N ij (t) is the abundance of individuals of the species i in the site j at time t, N ij (t + 1) d is this abundance after dispersal (i.e., after Equation 1), and N ij (t + 1) c after environmental filtering and competition (i.e., after Equation 2). M is the total number of localities (being k any other than j) and S the regional species richness (where g are different species than i).
(2) In step 1, dispersal rates were determined by the proportion of potential emigrants (a; Loreau et al., 2003). We randomized the number of emigrants (aN ij (t)) following a Poisson distribution (Thompson et al., 2020). Following Equation 1, all the localities received the same number of immigrants of each species. Additionally, we also considered an unequal distribution of emigrants relying on the distance among localities (Gravel et al., 2006;Sokol et al., 2017;Viana and Chase, 2019), increasing dispersal resistance with the distances among sites (see Supplementary Material). For simplicity, we only show here the results in a scenario where emigrants could reach all the localities with the same probabilities (i.e., very low distance limitation).
In step 2 (Equation 2), local dynamics were influenced by the abundance of each species after dispersal (N ij (t + 1) d ), the intrinsic growth rate (r i ) and competition (α). The intrinsic growth rate (r i ) relied on the performance of each species (λ i ) to the environmental conditions (E), which characterized the strength of the niche differentiation (Tilman, 2004;Gravel et al., 2006). The performance was defined by a Gaussian distribution, where µ i is the optimal environmental value and σ i is the niche breadth (Equation 3). We randomly modified the resulting r i λ i (E)N ij (t + 1) d following a Poisson distribution in each step, adding stochasticity to the model (Hubbell, 2001;Adler et al., 2007). We also considered intra-(α intra ) and interspecific (α inter ) competition in each locality (equation 2; Chesson, 2000;Adler et al., 2007).

Scenarios and Simulation Routine
The model ran in simulated landscapes of 10 localities with 10 species during 3,000 timesteps. Localities were randomly distributed in a square of 100 units of side. We simulated one environmental variable ranging from 0 to 1 for each landscape, which was spatially and temporally autocorrelated based on an exponential covariance model (Thompson et al., 2020). We set a spatial and temporal scale of autocorrelation of 50 and 500 units, respectively (see Supplementary Material). We generated different landscapes for each simulation.

Dispersal (a)
Constant frequency B C A FIGURE 2 | Summary of the simulated scenarios across niche breadth (i.e., selection; σ ) and dispersal rate (a) gradients (A) and a summary of the temporal sampling procedure in simulated (B) and empirical (C) metacommunities. In the conceptual models, polygons, and colors represent localities with different suitable conditions for an hypothetical species depending on its niche. The lines and arrows between the polygons show the connectivity among them.
The initial species matrix composition was randomly configured with a mean of 1 individual of each species per locality in a Poisson distribution. During the first 200 timesteps the environmental variable was kept constant in each locality (maintaining the spatial variability). In this period, we added each 10 steps more individuals following the same procedure as in the starting species matrix. After that, the environmental variable varied spatially and temporally and we did not add more individuals. We later removed the first 1,000 timesteps (200 + 800 steps) for further analyses to avoid the effect of the initial conditions on the results (Thompson et al., 2020).
We considered two levels of niche differentiation in the simulations (Figure 2A). In both scenarios, we set the same intrinsic growth rate for all the species (r = 10). The niche differentiation was determined by the niche breadth (σ i ) and it was the same for all the species in each metacommunity (Gravel et al., 2006;Sokol et al., 2017;Viana and Chase, 2019). In the neutral conditions (σ = 10), niche breadths were so wide to overlap species niches and cover the whole environmental gradient with the highest fitness (consequently all species have the same fitness at all the environmental values; Hubbell, 2001). Otherwise, with niche differentiation (σ = 0.5), niches were not totally overlapped and fitness was different for each species relying on their optimal environmental value (µ i ) and the environmental conditions at each site and time. We also assumed stable local coexistence in both scenarios. For this purpose, in the neutral scenario with equivalent fitness among species, stabilizing processes were weak (α intra = 0.050; α inter = 0.048). But we reinforced stabilization in the niche scenario due to established differences in fitness (α intra = 0.050; α inter = 0.028; Adler et al., 2007).

Empirical Metacommunities
We used three databases of aquatic invertebrates from previous works, which were here reanalyzed with other purposes. The first database had 10 localities in a stream sampled on 11 occasions for 19 months (the first months were sampled seasonally and the last ones monthly; Mezquita et al., 1999;Rueda et al., 2002). The other two databases correspond to two groups of interdunal ponds, with 9 and 13 localities, respectively, sampled monthly during 12 months (Valls et al., 2013;Rueda, 2015). All the data originate from the eastern Iberian Peninsula, with an intraannual fluctuation of temperatures and precipitations typical of Mediterranean regions (i.e., mild winters, hot and dry summers, and rainfalls concentrated in autumn months).
We selected four groups of arthropods (Ostracoda, Odonata, Coleoptera, and Diptera), each one corresponding to a taxonomic group with similar trophic features for the encompassed species, which potentially compete in a region for similar resources (Hubbell, 2001). Ostracods are benthic microcrustaceans which are predominantly omnivorous and they spread passively among ponds (Mesquita-Joanes et al., 2012). The other three groups have active dispersal with winged adults, although they have different dispersal capabilities and body sizes (Schmidt-Kloiber and Hering, 2015). Dragonflies and damselflies (Odonata) are predators, as well as the coleopterans (here, only family Dytiscidae). However, these two groups have different prey preference, attack strategies and life histories. Finally, dipterans (only families Stratiomyidae, Psychodidae, Ephydridae, and Dixidae) present aquatic larvae mostly consuming detritus. The criteria for selecting these groups were based on their abundance and species richness in the datasets. Overall, we analyzed nine datasets, separating the invertebrate groups in each landscape setting (3 in each landscape). More information of the empirical data is available in the Supplementary Material.

Sampling Metacommunities
We sampled the simulated and the empirical metacommunities across a gradient of study durations (i.e., temporal extents), maintaining the number of temporal sampling points or frequency (therefore varying the lag between samples; Figures 2B,C). In the simulated data, we sampled the metacommunities varying the temporal extent from 20 to 2,000 timesteps ( Figure 2B), fixing the frequency at 10 times. In a similar way, empirical metacommunities were sampled at different study durations maintaining the number of temporal points at a frequency of 6 times ( Figure 2C). After being sampled, the simulated and the empirical data have the same arrangement in matrices, and we applied the same analytical routine for both types of datasets, as follows.

Metacommunity Variation Partitioning
The empirical and the simulated data (after sampling) were analyzed by means of variation partitioning between environmental, spatial and temporal variables (Borcard et al., 1992;Anderson and Cribble, 1998;Peres-Neto et al., 2006). Environmental variables in the empirical dataset were previously log-transformed, and later we carried out a PCA to use the principal components as orthogonal environmental variables with all the environmental variability measured. In the simulated data, we used directly the simulated environmental variable (E).
We used Moran Eigenvector Maps (MEMs; Dray et al., 2006) to model spatio-temporal variables (see Supplementary Material). For this purpose, we created a three-dimensional network with the sampling points and the links among them, distributing them across space (x and y axes) and time (z axes). The links among samples were decided in two ways. In the spatial axes, the samples of each temporal point were connected following a Gabriel graph criterion, but never across time (Legendre and Legendre, 2012). In the temporal axis, we linked independently all the temporal points of each site separately, uniquely connecting each sample with the previous and the next one in the same site (Legendre and Gauthier, 2014). Therefore, we did not establish temporal connections between different sites. All the links were row standardized depending on the number of links of each point (e.g., 1 link = 1, 2 links = 0.5, and 0.5, 3 links = 0.3, 0.3, and 0.3, and so on), obtaining the weights of each link and calculating the MEMs .
With this method, the resulting MEMs modeled spatiotemporal patterns. We separated them in pure spatial and pure temporal MEMs by means of a two-way ANOVA, using each MEM as a response variable and spatial sites (space) and temporal points (time) as categorical explanatory variables (based on Legendre et al., 2010;Legendre and Gauthier, 2014).
We considered spatial MEMs those with significant space, and temporal MEMs with significant time. Some of them could be significant for both space and time or for any of them. These MEMs were added into the two groups (spatial and temporal), and eventually, variation partitioning can attribute the explanation associated to these MEMs as the shared fractions among space and time (spatio-temporal fraction).
The variation partitioning procedure was based on RDA to calculate the fraction of variation of the species matrix explained by environmental, spatial, and temporal variables as R 2 Legendre and Legendre, 2012). RDA presents some important limitations which could affect the explained fractions , although it is one of the most used in metacommunity analysis. For example, RDA uses linear regressions to obtain the R 2 , although the species response to an environmental gradient could be nonlinear (niches are usually Gaussian as in our simulation). In order to minimize this issue, we added a quadratic term for the environmental variables. In RDA, species matrices were transformed with the Hellinger method (i.e., the squareroot-transformation of the relative abundance; Legendre and Gallagher, 2001). Additionally, we performed the same analyses by means of RDA without the quadratic term, so as with CCA and dbRDA (see Supplementary Material).
Variation partitioning has further limitations. Spurious correlations between environmental and spatio-temporal variables can occur due to spatio-temporal autocorrelation of the environment (Smith and Lundholm, 2010). This artificially increases the fractions (as R 2 ) shared between the environment and the spatio-temporal variables, whereas this spurious correlation is not associated with the environment. We corrected R 2 using a method based on MSR (Wagner and Dray, 2015), creating replicates of an environmental matrix with the same autocorrelation properties but removing their relationship with the species matrix (as a null model; Clappe et al., 2018). Previously to the variation partitioning, we selected subsets of each set of explanatory variables (environmental and all the MEMs with significant positive Moran's indices, before separating them in spatial and temporal) by means of forward selection, with a double stopping criterion (Blanchet et al., 2008a). We applied a Bonferroni correction for multiple testing to the p-values of the selected variables, i.e., those with adjusted p-values lower than 0.05.

RESULTS
In the simulated metacommunities, the explained fractions (adjusted R 2 ) varied across selection and dispersal levels (Figure 3). Generally, the effects of the temporal scale increased as we incremented the dispersal rates (a), for both niche breadth levels (σ ). These temporal scale effects raised the relevance of the temporal variables, and reduced the fractions explained by other factors. The shared fraction between space and time was negligible in all the scenarios, indicating pure temporal or spatial patterns in the metacommunity structure due to the design of the simulations (e.g., equal distribution of emigrants among the localities).
With neutral conditions (σ = 10; Figures 3A-C, top subpanels), spatial variables dominated at low dispersal rates (a = 0), whereas temporal variables were prevalent at high dispersal rates (a = 0.4). At intermediate dispersal rates (a = 0.02), the pattern was a combination of the two extremes.
Dispersal determined the main dimension (spatial or temporal) of change of the metacommunity. We can better understand this behavior if we observe the temporal dynamics of a single random species in each of the ten sites separately (Figures 2A-C, bottom subpanels). At low dispersal, spatial variation was high and maintained through time due to neutral dynamics. Therefore, space was relevant independently of the study duration. This spatial variation disappeared as we increased the dispersal rates (i.e., synchronizing the localities), until all the variation was concentrated in the temporal axis. However, the relevance of the temporal variables was not constant, increasing logarithmically with the study duration (Figures 3B,C).
When selection was included in the simulations by establishing narrower niches (σ = 0.5; Figures 3D-F) environmental effects appeared, and spatial and temporal variables exhibited a similar behavior to the neutral scenarios. However, the relevance of the environmental variable also depended on dispersal rate, accounting for larger fractions of metacommunity variation at intermediate dispersal rates (a = 0.02). Generally, the shared fraction between environment and space was higher at short study duration, whereas the shared fraction between environment and time slightly increased with the study duration (a pattern similar to pure spatial and temporal fractions). With low dispersal (a = 0), the species dynamics showed non-occupied suitable localities (λ ≥ 0.9; see the orange line for N = 0 in Figure 3D, bottom subpanel). In these conditions, the environmental variable was less relevant than the spatial variables, which captured this spatial variation in species distribution. At intermediate dispersal rates (a = 0.02), the species abundances fitted the environmental suitability through time, according to a constant dominant relevance of the environmental variable (Figure 3E). At high dispersal rates ( Figure 3F) the represented species (lower subpanel) survived even under unsuitable conditions. However, its abundances were higher in the suitable than in the unsuitable conditions, showing an important role for environmental effects at least at short and intermediate study durations. We can also observe a synchronization of the local populations, increasing the relevance of the temporal effects as we extended the study duration.
When exploring the empirical data on aquatic invertebrate metacommunities, explanatory variables accounted only for a relatively low proportion of the metacommunity variation, compared with the simulated scenarios, and their relative role depending on time extent varied widely among groups (Figure 4). Generally, pure environmental fractions were more relevant than the others, explaining about 5-25% of the metacommunity variation, particularly for the ponds (Figures 4B,C). In these data, the spatio-temporal position of the samples were more complex than in the simulations, and we can observe some relevance of combined spatiotemporal effects (i.e., the shared fraction between space and time). Extending the analyses of the empirical data at different temporal scales, while keeping the number of temporal sampling points fixed, we can observe changes in the variation partitioning results, but the effects of the study duration were unclear (Figure 4). Extending the study duration seems to generally increase the relevance of the environmental variables and decrease the spatio-temporal fractions. However, the opposite was true in some cases, and anyway these changes were relatively very low compared with the simulations.

DISCUSSION
Our results show that the temporal scale of observation affects our perception of the main processes contributing to metacommunity organization. Previous studies already  Ponds 1 (B), and Ponds 2 (C) in a gradient of study durations from 6 to 19 months (A) or 6-12 months (B,C), with six temporal sampling points. In the Stream (A), we displayed ostracods (left), odonates (center), and coleopterans (right). In the Ponds 1 (B), we represented dipterans (left), odonates (center), and coleopterans, and in the Ponds 2 (C), ostracods (left), dipterans (center), and coleopterans (right).
highlighted the influence of spatial scale on the inference of these processes when analyzing communities sampled once (e.g., Condit et al., 2002;Heino et al., 2017;Viana and Chase, 2019), or comparing the spatial effects between several temporal points of the same metacommunity (e.g., Langenheder et al., 2012;Fernandes et al., 2014;Castillo-Escrivà et al., 2017). Other studies have focused on the effects of temporal scales on temporal turnover, but at one spatial point (e.g., Korhonen et al., 2010;Tomašových and Kidwell, 2010). Here, we integrated both perspectives, taking into account the spatio-temporal structure of the data as a whole and testing the influence of sampling period extent on the results of variation partitioning.
The simulations pointed out how dispersal and selection play a key role on metacommunity variation in both the spatial and temporal dimensions. Low dispersal rates allow maintaining the spatial variation through time, because local dynamics are spatially independent without a flux of organisms (Koelle and Vandermeer, 2004). In these cases, long study durations present the same results as short study durations in the simulations, suggesting to focus our sampling design only across space when dispersal is very limited (PD archetype), taking into account ranges of distributions and dispersal abilities of the organism studied (Wiens, 1989). On the other side, high dispersal reduces the spatial variation and temporally synchronizes all the local communities Gonzalez et al., 2009;Pandit et al., 2013). This synchronization implies that the whole metacommunity variation goes to the temporal axis (as expected with dispersal surplus in ME, and NT archetypes), and the longer the study duration, the higher the influence of time needed to capture these effects (e.g., ecological drift, progressive dispersal) on metacommunity variation.
The distribution of the environmental variation through time and space can also determine the importance of an adequate spatial or temporal scale of observation in nicheconstrained metacommunities (SS archetype; Korhonen et al., 2010;Viana and Chase, 2019). In our simulations, the environmental variation was predominantly spatially distributed and the localities offered the whole environmental gradient at all the timesteps (even though the local environment was changing through time). Consequently, species always had available localities with suitable conditions when a locality became unsuitable; the key is an adequate dispersal to reach the optimal sites. Notwithstanding this environmental heterogeneity, extending the study duration affected the relevance of the shared fractions between environment and space-time, decreasing the fraction shared with space and slightly increasing the fraction overlapped with time. However, if the environmental variability would be expanded through time (as under climate change scenarios; Thompson et al., 2015), the temporal extent of observation would have a more important effect on detecting such environmental relevance.
In the empirical metacommunities, we observed relatively large differences when modifying the temporal scales of observation, usually increasing pure environmental effects and decreasing spatial pure and shared fractions. However, there was a large variation in this response among groups of organisms and landscape settings. According to the results of the simulations, the studied empirical metacommunities most probably rely on low to intermediate dispersal rates, combined with strong environmental filtering. However, the low explained proportion of metacommunity variation also suggests that finescale stochastic dynamics might have an important role in the empirical settings. Such low explained proportions are indeed very common in metacommunity studies (Cottenie, 2005), and may originate, among other sources, from the sampling design and methods used. First, the spatial and temporal scales could not be adequate to fully capture metacommunity variation, as we find in the simulations, emphasizing the importance of the sampling design (Viana and Chase, 2019). Second, we might have not acquired some relevant environmental information, failing to fully explain niche-related patterns and therefore increasing the unexplained fraction. Third, the statistical methods used might not be the best for the data being analyzed (e.g., linear models for non-linear environmental responses), although recent studies have developed new and promising methods (Clappe et al., 2018;. Moreover, the large unexplained variation in empirical data demonstrates multiple and complex variation sources largely disregarded in the simulations, such as methodological errors when sampling or processing the samples, the influence of rare species (Magurran and Henderson, 2003) or trophic interactions (Guzman et al., 2019;García-Girón et al., 2020).
Previous empirical studies analyzing spatio-temporal patterns of metacommunities (Muylaert et al., 2000;Ysebaert and Herman, 2002;Padial et al., 2014) highlighted the relevance of the environment, even when environmental variables were temporally structured. According to our results, this also suggests a dependence on the temporal scale of observation, particularly whenever the environmental heterogeneity increases with the study duration. However, these studies only discussed pure spatial and temporal patterns, whereas we added explicitly mixed spatio-temporal patterns and the results could be different because of this combination (although some of the studies had indeed detected an overlap between space and time; Padial et al., 2014). Nevertheless, the spatio-temporal MEMs used here have some issues, such as their origin on symmetric connections, which may not be adequate for time (unlike AEMs, which consider an asymmetric relationship between points; Blanchet et al., 2008b).
In our empirical metacommunities, the shared fraction between space and time was important in some cases, indicating that the spatio-temporal variation in natural metacommunities might be more complex than that shown by our simulations. For example, ephemeral ponds do not disappear and appear always at the same time, creating spatio-temporal patterns in the data. The sampling design in these cases should take into account both spatial and temporal scales together. Therefore, the sampling design should be based on a previous knowledge of the temporal fluctuations of the landscape, this being important to find natural references of stable (such as tropical systems or stable interstitial zones; Hubbell, 2001;Dumas, 2002) and unstable landscapes (such as temporary ponds or intermittent rivers; Castillo-Escrivà et al., 2017;Cid et al., 2020) to test the theoretical predictions. This also encourages the use of more complex landscapes in the simulations for specific cases, taking into account different types of environmental variation (e.g., mosaic, gradient; Viana and Chase, 2019), the type of spatial distribution of localities (e.g., regular, random; Henriques-Silva et al., 2015), and regular or unexpected temporal events (e.g., droughts, seasonality; Tonkin et al., 2017).
We applied a basic model for the simulations that assumed that all the organisms of each metacommunity had the same population parameters. The study of metacommunity dynamics requires more complex models, considering the survival of the individuals through time (as a kind of temporal dispersal, such as diapause or lethargy) or the population structure (e.g., eggs, juveniles, adults). Simulating or monitoring metacommunities with different life-history strategies (Verberk et al., 2008) may allow comparisons at several temporal scales of observation among groups that differ in development time, synchronization in the reproduction, type of reproduction (e.g., semelparity or iteroparity) or with a population age-structure with changing dispersal capabilities among stages (e.g., juvenile and adult dimorphism). This type of comparison among organisms with different traits (e.g., body sizes and dispersal capabilities) has been a frequent approach in spatial studies of metacommunities (e.g., Soininen et al., 2011;Astorga et al., 2012;De Bie et al., 2012).
Another limitation of our simulations was to consider the same life-history traits (e.g., fecundity, survival, life cycle, dispersal capabilities) for all species. Natural communities are more complex, and variable life-history strategies within the community may allow the coexistence of species with similar trophic requirements (Amarasekare, 2003). For example, metacommunities could be composed by good competitors (with long life cycles, iteroparous, low dispersal rates, and niche dependence) and good colonizers (with short development times, semelparous, high dispersal rates, and low niche dependence; Chave et al., 2002). The prevalence of a trait could determine the relevance of the assembly processes in a metacommunity at different observational scales, whereas rare species may drive an increase in the amount of unexplained variation (Magurran and Henderson, 2003). Additionally, priority effects can also strongly influence the relevance of temporal factors in metacommunities (Fukami, 2015). As a process that is basically temporal, we can emulate historic and priority effects considering resistant stages to adverse transient environments (Wisnoski et al., 2019) and/or setting different inter-specific competition parameters for the modeled species (Thompson et al., 2020). Despite all these limitations, our model provides a general view of the effects of the temporal scale of observation on understanding metacommunity dynamics and methodological problems related to its study.
We focused on the temporal scale of observation using only one spatial scale, but the distribution of the metacommunity variation in the spatial and the temporal dimensions depends on both scales (increasing the extent of one of them decreases the relative relevance of the other). It is important to predict at which spatial and temporal scales the metacommunity variation is balanced in both dimensions, allowing to establish equivalences between spatial and temporal units (Adler and Lauenroth, 2003). We could find the equivalence by analyzing the data, but we need to develop deductive methods to estimate this equivalence taking into account the community structuring processes. This equivalence between space and time is crucial to plan research or conservation projects.

CONCLUSIONS
The inference of metacommunity structuring processes is influenced by the temporal scale of observation, depending on the distribution of the metacommunity variation between the spatial and the temporal dimensions. In our simulations, the variation was accumulated in the spatial dimension when dispersal rates were low, and in this case, an increase of the study duration did not have any effect on the inference of the underlying processes. On the other hand, high dispersal rates synchronized all the local communities, reducing the metacommunity variation to the temporal dimension. In this case, the study duration influenced our estimation of the metacommunity organization processes, logarithmically increasing the role of time-related effects (dispersal movements, ecological drift) when increasing the temporal extent of the study. In addition, we found the temporal scale effects to depend on the distribution of the environmental variation, which can generate more spatial or temporal variation in the metacommunity. In the empirical data, we observed only a slight influence of increasing temporal scale, mostly producing an increase of the role of environmental effects and decrease of spatial effects. This suggests metacommunity variation was rather distributed along the spatial axes and/or that the temporal scales considered were not adequate (maybe too short) to observe an increase in the role of time. However, the empirical data is expected to be more complex than the simulations, possibly including many unmeasured spatial, environmental, and stochastic effects. The present study highlights the importance of selecting an adequate observational scale to study or assess natural metacommunities, and the necessity to develop better methods to model and analyze spatio-temporal dynamics of natural systems.

DATA AVAILABILITY STATEMENT
All datasets generated for this study are included in the article/Supplementary Material.

AUTHOR CONTRIBUTIONS
FM-J, JR, and AC-E conceived the study. JR carried out the fieldwork and the identification of the organisms of the empirical data. JR and FM-J provided the empirical datasets. AC-E and FM-J designed the simulations. AC-E programmed the code and analyzed the data. AC-E wrote the first draft of the paper, with further revisions and contributions by FM-J. All authors approved the final version.

FUNDING
Partial funding has been provided by the METACOM-SET project (code CGL2016-78260-P), funded by the Spanish Ministry of Economy and Competitiveness (Agencia Española de Investigación, and FEDER).