Simulation of pattern of gene flow in Canjerana fragments in the Brazilian Atlantic Rainforest for evaluating genetic..

Gene flow is important for the conservation of genetic resources to allow connectivity of geographically isolated populations and which genetic variability is reduced. Gene movement is a function of flow rate and model. Understanding how gene flow occurs can contribute to the conservation and selection of priority populations that could benefit from an eventual intervention. Simulation softwares allow making inferences about past events based on current datasets or predict future phenomena under real genetic scenarios. Adverse phenomena can be predicted and actions can be taken to avoid them. The aim of this study was to identify a model and the gene flow rates that could explain genetic structure of eight forest fragments of Cabralea canjerana in development in the Brazilian Atlantic Rainforest. To do this, simulations were performed with the EASYPOP software using a microsatellite marker dataset obtained for the species by Melo and collaborators, in 2012, 2014 and 2016. We tested five models and nine migration rates and we selected the model that produced values closer to those previously obtained for them. Criteria used for selection were the observed and expected heterozygosity and the Wright’s F Statistics obtained in the simulations. The gene flow model selected was the isolation by distance model that used a rate of 0.1. We observed high levels of genetic differentiation among the fragments as result of their reproductive isolation. To allow homogenization of the allelic frequencies through gene flow, the solution would be to create ecological corridors with the aim of connecting distant fragments.


INTRODUCTION
In general, most conservation programs are focused on reducing gene flow barriers between isolated populations that carry a small fraction of the total variability in relation to continuous populations. Populations of a species can be reproductively isolated by a physical barrier (e.g. mountains) or by fragmentation, which is often caused by several types of human activities including construction of roads, railways, cities, and expansion of farmland. In small fragments, events like genetic drift and inbreeding may occur, which can reduce genetic variability (PRIMACK & RODRIGUES, 2001;BROQUET et al., 2010). Serrote et al. Over time, populations with low genetic variability have low evolutionary potential and limited ability to respond to environmental changes in long term. While a large area habitat can hold a single large population, it is possible that none of its fragments can maintain a subpopulation for an extended period (PRIMACK & RODRIGUES, 2001).
Population connectivity and redistribution of genetic variability among the individuals allow the homogenization of allelic frequencies, thus providing the viability of populations with sizes below those required to control for the adverse effects of drift and inbreeding. Thus, small isolated populations start to behave as a large panmitic population with greater genetic variability and this reduces the risk of its extinction. Furthermore, the restoration of gene flow homogenizes the allelic frequencies, thereby maintaining genetic variability and polymorphism (SLATKIN, 1985;ZUCCHI et al., 2003;PINHO & HEY, 2010).
Improving our understanding of gene flow (model and rate) between geographically isolated populations is important for planning conservation strategies, including the selection of priority populations that would benefit from an eventual intervention. Wright's island model is a conventional model used in population genetics to evaluate gene flow between populations and allows the estimation of the actual number of migrants (Nm) by means of the F ST statistic for an n set of populations (WRIGHT, 1965). However, this model does not have realistic assumptions for all populations, such as the equilibrium between migration and genetic drift, if there is equal gene exchange between populations, and if all populations have equivalent sources of migrants. Also, the estimates obtained with this model do not reflect the contemporary variation in gene exchange between populations or current changes in the dispersal process.
Due to these limitations, there is need for more accurate methods for unraveling the models and gene flow rates for specific populations in order to ensure eventual success of genetic conservation programs. With advances in computer technology, simulation softwares were developed to make inferences about past events from actual data, as well as those that perform simulations forward in time, allowing predictions of future phenomena using real genetic scenarios (YUAN et al., 2012;DALQUEN et al., 2012). In this simulation research we focalized Cabralea canjerana ssp. canjerana (Meliaceae), a dioecious forest species that is considered model to Brazilian Atlantic Rainforest conservation studies (KLEIN, 1979). The Atlantic Rainforest is the Brazilian biome most affected by environmental fragmentation. Today's the biome presents itself as a mosaic composed of few relatively large areas. Thus, the forest fragments of diverse sizes and shapes, assume fundamental importance for the perenniality of the Atlantic Rainforest biome (ZAÚ, 1998).
The objective of the present study was to identify the model and the gene flow rate that best explain the genetic structure of eight fragments of C. canjerana using data from microsatellite markers.

Simulations
The EASYPOP Version 2.0.1 software (BALLOUX, 2001) was used to simulate different models and flow rates, which were selected based on how well they could explain the genetic structure of eight C. canjerana fragments obtained by six microsatellite markers from MELO (2012), MELO et al. (2014) and MELO & FRANCESCHINELLI (2016). In these population fragments studied were sampled: (#1) 27, (#2) 33, (#3) 31, (#4) 33, (#5) 17, (#6) 11, (#7) 39, and (#8) 30 individuals. All of the populations were located in the Atlantic Rainforest biome in the Environmental Protection Area Fernão Dias, which is in the south of the state of Minas Gerais, in Brazil ( Figure 1 and Table 1). This area has suffered in the last 60 years a great deforestation for extractivism and potato planting. The Environmental Protection Area Fernão Dias has 180,073 hectares and was created in 1997 as an environmental compensation plan for the duplication of BR 381 (MELO, 2012).

Simulation parameters
In this study, we tested the following migration models: one-dimension stepping stone, island, hierarchical stepping stone ('contact zone'), hierarchical island, and spatial (isolation by distance). For the hierarchical stepping stone ('contact zone') model, two groups were considered with four fragments in each group. The first group consisted of 1, 2, 3, and 4; and the second group was composed of 5, 6, 7 and 8 fragments. The contact zone between the groups was composed of fragments 4 and 5. In the hierarchical island model, two archipelagos were composed of four islands (fragments) each. The archipelago for each fragment was the same composition as the previous model groups and the distribution of the fragments used was the same as that in the hierarchical stepping stone model.
In the spatial model simulation, two dimensions based on the real geographic position of the fragments were used. The average distance of pollen dispersion was 127 m (MELO & FRANCESCHINELLI, 2016). For the seeds, the same dispersal value was used because, according to SLATKIN (1985), seed dispersal distance estimates are similar to those of pollen.
The following parameters were considered in our simulations: a dioic species, migration rates between 0.1 and 0.9, six loci, mixed mutation Singlestep Model (SSM) with a proportion of K allele Model (KAM) mutation events (with a ratio of 0.1) and 25 possible allelic states, maximum variability of the initial population, and a mutation rate of 0.00005, which is in accordance with the mean observed in nature (PIRES et al., 2011). The generation number for each simulation was 100. For each model and migration rate combination 100 independent replicates were conducted. We selected the model that showed values of observed and expected heterozygosity and Wright's F-Statistics closest to those previously obtained with microsatellite markers by MELO (2012)   Serrote et al.

Model selection
Based on observed and expected heterozygosity average values of 0.69 and 0.72, respectively, obtained by MELO et al. (2014) we pre-selected the following models for gene flow: the one-dimensional stepping stone model with migration rates of 0.5 and 0.6; island model, which has a migration rate of 0.1; the hierarchical stepping stone model with a migration rate 0.2; and the spatial model, which has a migration rate of 0.1. Using these models, we obtained statistics for genetic structure similar to those obtained with molecular markers ( Table 2).
The exclusion criteria used in the model selection was pollen dispersal distance; the average distance obtained from the microsatellite markers was 127.33 m, and a range of 42 -1,040 m. Given that physiological barriers may prevent gene exchange between fragments, the statistical pairwise difference (R ST ) was assessed. The data obtained by MELO & FRANCESCHINELLI (2016) showed a strong correlation with the Mantel test (r= 0.836, p<0.0001) ( Table 3). It was, therefore, considered that fragments with a distance of more 1,040 m could not exchange genes with each other. The fulfillment of the first law of geography (TOBLER, 1970) was also observed, with less genetic differentiation occurring among closer fragments and vice versa.
Based on the pollen dispersal distance and geographic arrangement of the population fragments, the fragment groups 1, 2, 3, and 4 were exchange genes with each other. This exchange was also observed between fragments 5 and 6, but not in the groups that were over 1,040 m apart (MELO & FRANCESCHINELLI, 2016). Thus, only the spatial model of isolation by distance with migration rate of 0.1 could explain the pattern of gene flow between the eight C. canjerana fragments studied ( Table 2). The difference of statistics in migration rate 0.1 in relation to others is significant, by t test at 5% error probability. Spatial model and migration rate are the most suitable for the genetic characterization of these population fragments.

Characterization of the fragments
The gene flow rate between C. canjerana fragments was 0.1, as indicated in the isolation-bydistance model. Only proximal fragments could  exchange genes between them, if they were limited to the maximum dispersion distance of 1,040 m; however, a migration rate of 0.1 underestimates the actual movement of pollen and seeds, since it considers only the successful migrants. According to ENNOS (1994), the estimation of gene flow in plant populations is complex because it is done in two ways -through both pollen and seeds. Gene flow via pollen involves the dispersal of pollen of one population and the successful fertilization in the flowers of another population; while gene flow via seeds involves the dispersal of seeds of one population and their successful establishment in a new population (RIDLEY & MALLORY-SMITH, 2015). A complete description of gene flow between plant populations should ideally include the evaluation of the relative importance of both pollen and seeds as gene flow agents. This analysis is very complex and highlights the value of simulation models, which can be based on the genetic structure that is actually observed in the population.
Due to species population fragmentation, both the observed and expected heterozygosity decreased over 100 generations (Figure 2A). This reduction in heterozygosity occurred due to the loss of alleles by genetic drift. Simulations performed by STEFENON & COSTA (2012) with EASYPOP software indicated an increase in the inbreeding coefficient over generations in small populations compared to larger populations, resulting in lower genetic variability. However, random mating contributed to an observed heterozygosity that was similar to that expected from the Hardy-Weinberg equilibrium in all generations. This result confirms the importance of intra-population gene flow for the conservation of genetic variability within these fragments. A study of Luehea divaricata fragments in the Pampa biome showed that the largest fraction of genetic variability was distributed within the fragments due to high crossing rates (NAGEL et al., 2015).
Another effect of gene flow is the reduction of genetic divergence between populations.
To study this effect in the C. canjerana fragments, we plotted F ST values with migration rates of 0.1 and 0.9 ( Figure 2B). There was a differentiating coefficient near zero when the migration rate was 0.9, resulting in homogeneous fragments. This indicates that fragmentation of the population was not strongly affected by inbreeding, which is attributed almost exclusively to two-parent crossings within fragments. A migration rate of 0.1 caused high levels of inbreeding due to reproductive isolation that limited the gene exchange between fragments. Genetic differentiation was more pronounced in the first generations, as observed in the stabilization after the thirtieth generation. This delayed stabilization suggested that the loss of genetic variability in all fragments did not allow genetic differentiation to continue growing.
The decrease in heterozygosity suggested the same behavior in the future, requiring the adoption of interventional measures that reduce the risk of population extinction. Similarly, the simulation using a higher migration rate allowed us to predict how much genetic variability was lost due to the low migration rate, further confirming the need for the development of strategies for the conservation of these genetic fragments.
Assuming that all fragments of the first group (1, 2, 3, 4) exchange genes between each other, as a function of their geographical distance and the correlation between the geographical distance and the value of R ST , we proposed the creation of ecological corridors that would connect the fragments 7 and 8, 6 and 7, 5 and 8, and 6 and 8 as a potential intervention plan. This plan could also include a contact area between fragments 4 and 5, which would follow the hierarchical island model for gene flow, and form a second island between the fragments 5, 6, 7 and 8. A more economical but less effective alternative could include ecological corridors only among the 4 and 5, 6 and 7, and 7 and 8 fragments, in which case the gene flow would be governed under the mixed model of islands and stepping stones.

Implications for conservation
The studied area has a deforestation history of more than 60 years, which led to its fragmentation.
In addition, C. canjerana does not have long-distance pollen dispersers. The main pollen dispersers are the moths of the Lepidoptera order, which have short flight behavior, pollinating nearby trees (PIZO & OLIVEIRA, 1998). The dispersion of fruits and seeds is zoocorical, mainly by birds, attracted by the orange color of the aril (PIZO, 1997). As a result, high levels of genetic differentiation were observed among fragments due to their reproductive isolation characterized by low migration rates. This condition propitiates the occurrence of inbreeding and genetic drift that lead to progressive loss of genetic variability and, ultimately, population extinction.
Therefore, the findings in this study allowed a better management of the area toward its conservation by enabling gene flow among fragments, thus increasing the genetic variability.

CONCLUSION
The eight fragments of C. canjerana studied have high levels of genetic differentiation due to reproductive isolation (the gene flow model is isolation by distance), and because a gene flow rate equivalent to 0.1 is not sufficient for the homogenization of allelic frequencies. In this condition these fragments demand special attention to reduce their vulnerability to adverse events. We suggest the creation of ecological corridors in order to increase the flow rate, connect the distant fragments, and provide homogenization of the fragments via gene flow.