Worldwide co-occurrence analysis of 17 species of the genus Brachypodium using data mining

Implementation of the algorithm in Python

The methodology presented in Silva et al. (2016) was implemented in Python 2.7 followed by the Algorithm 1 (Silva et al., 2016), (https://github.com/simonorozcoarias/co-ocurrence_analysis). Two parameters were needed to create a transactions file, the geo-located data with the following structure: species_name, latitude, longitude and distance in meters. The first step for the algorithm was to create a folder named “DE_results” with the output files (transactions.csv, all_rules.txt, positive_rules.txt, negative_rules.txt), using all records of the data file. The distance between all records was calculated, one by one, based on its latitude and longitude using the formula of Haversine (Chopde & Nichat, 2013). If the distance was lower than a parameter (1, 5, 10, or 20 km), a transaction was created and stored in the transaction file. Each transaction was composed by the name of the species involved in it, separated by a comma.

The second step was to define the basic values of parameters such as minimum support, minimum confidence, negative minimum lift, and positive minimum lift (Table 2), which lift is the measure of importance of a rule. Next, rules were generated through the a priori algorithm presented in Agrawal & Srikant (1994), using the transaction file previously created. These rules are computed from the data and, unlike the if-then rules of logic, association rules are probabilistic in nature. This algorithm created a set of rules, where each rule contained antecedent (if) and consequent (then) elements and a confidence value, that is, if the algorithm creates the following association rule species X ⇒ species Y with confidence 0.5 indicates that species X and Y are usually found together and the confidence value 0.5 show us that if you find species X there is a 50% chance to find species Y. In addition, the frequency of each element of the rule was calculated using the Chi-squared function from the scipy.stats python package (Woods & Amasino, 2015). This function calculated the Chi-squared and p-value for the rule. We calculated support and lift values from definitions presented in Silva et al. (2016). We considered a rule as positive if it had a p-value equal to ∼0, and a rule as negative if it had a p-value equal to ∼1. Positive and negative rules were stored in the all-rules file.

In this study was used the following set of measures to evaluate the quality of the rules generated by the algorithms of association rules: support, confidence, lift (Han & Kamber, 2006), Chi-square (Hahsler, Grün & Hornik, 2005), and p-value (Liu, Zhang & Wong, 2011). The first two measures of support and confidence were used to define the species’ pairs and groups; the third (lift) assesses the association type (positive or negative), and the Chi-square and p-value measures consider the degree of independence of the species. For example, considering two species, B. retusum and B. phoenicoides, the support is the probability P of transactions with both species and is defined as support (B. retusum → B. phoenicoides) = P (B. retusum ∪ B. phoenicoides). The confidence is defined as the frequency with which items are found in the transaction B. retusum containing B. phoenicoides, and is defined as the conditional probability confidence (B. retusum → B. phoenicoides) = P (B. retusum | B. phoenicoides). The lift is the measure of the importance of a rule and can be defined by P (B. retusum ∪ B. phoenicoides)/(P (B. retusum) ∗ P (B. phoenicoides)) (Silva et al., 2016). The Chi-squared of a rule may be computed directly from the values of confidence, support, and lift (interest) of the rule in question (Alvarez, 2003). In the case of association rules the p-value of a rule, R = (B. retusum → B. phoenicoides), is defined as the probability of observing R or a rule more extreme than R, given that the two sides of R are independent. If a rule R has low p-value, then R has a low chance of occurring if its two sides are independent (Liu, Zhang & Wong, 2011).

Creation of the dataset for the evaluation of co-occurrence

The previously implemented algorithm was executed with seven different datasets, which contained a diverse amount of records and species (Table 3). We chose to analyze different combinations of species to ensure that the co-occurrences established have consistency and were not affected by the number of species involved or by the number of records. Whereby, we used two, three, and four species (datasets 1, 2, and 3) with the aim of identifying the minimum number of species possible for implementing the algorithm and to discover positive rules, considering that the dataset contained a number of diverse records and distributions.

In total we used 17 species of the genus Brachypodium, downloaded from the GBIF data repository (http://www.gbif.org) or obtained from bibliography sources (López-Alvarez et al., 2015; Schippmann, 1991). All datasets can be downloaded at https://github.com/simonorozcoarias/co-ocurrence_analysis/tree/master/raw_data. All of the data used in this study was validated and filtered according to each species and its native distribution (Table 1). In addition, each dataset was tested with four different distances (1, 5, 10, and 20 km), due to the cosmopolitan distribution of the genus Brachypodium, since some species have very specific distributions and it is necessary to establish their co-occurrence at different distances.

Exhaustive analysis of positive rules found in datasets with 15 and 17 species

Given the importance of analyzing the coexistence of a large number of species, datasets 6 and 7 acquire relevance in the study, since they bring us a little closer to a real situation in nature. For dataset 6, taking into account a distance of 20 km, we found only three species present in the rules; B. hybridum in five rules, B. phoenicoides and B. retusum in six rules. On the contrary, no positive rules were generated with the other distances. Nonetheless, by analyzing dataset 7 we found that using distances of 1, 5, and 10 km, the rules generated were composed of only four species, but with 20 km the rules were composed of five species, where the new species corresponds to B. rupestre (Fig. 1). Also, we analyzed the role for each species in their positive rules (as antecedent or consequent) using the 20 km distance for datasets 6 and 7 (Fig. 2).

We observed that as the tested distance increased so did the complexity of rules. Using distances such as 1, 5, and 10 km rules with only two species were generated, but with 20 km there appeared rules that were composed of three species. After creating the rules, the algorithm filtered them using the lift value in order to classify them as positive, negative, independent, or negligible. We found a relation of co-occurrences between some species contained within positive rules generated, and compared those relations with non-filtered generated rules.

Finally, we found only one difference between species that were contained in datasets 6 and 7, B. hybridum appeared just in the first. Using the number of appearances of each species in the created transactions and its relation with the number of total transactions, we found the importance that the amount of registries have in the creation of positive rules and its impact on other transactions.

Species distribution models quantify the relationship between species and their environments without considering the possible biotic interactions (Pollock et al., 2014). Therefore, not all the features that influence species co-occurrence will be captured by environmental variables. Accordingly, other alternatives to quantify co-occurrence can provide insights into the underlying causes of similarities and dissimilarities in distributions among species.

In the natural habitats and distributions of the genus Brachypodium, it is evident that some species are co-occurring. For instance, in the case of annual species such as B. distachyon, B. stacei, and B. hybridum (López-Alvarez et al., 2015) or in the perennial Eurasian species where B. pinnatum, B. rupestre, and B. sylvaticum grow in mesic to humid open grasslands and forests (Catalán et al., 2016; Scholthof et al., 2018). This kind of algorithms is very useful, because it enables us to mathematically check the co-occurrence of biological species. This study allowed us to determine the relationships between the different species of the genus Brachypodium using various data sets. We observed that it is necessary to obtain a broad number of generated transactions composed of various species in order to create rules for relations, since with few transactions the p-value has no statistical significance (Table 4). Therefore, the algorithm ignores this relation in the positive rule creation process. Table 4 clearly shows this case with datasets 1, 2, and 3.

We found that the coexistence of a small group of species can be evaluated and correlated with deeper studies of ecological niche models, in dataset 2, the coexistence of B. distachyon, B. stacei, and B. hybridum was of 20% for a distance of 10 km. These results agree with the data reported in the study of López-Alvarez et al. (2015), where the calculated area of B. hybridum was 1,464,910 km² and the overlap area of B. distachyon and B. stacei was of 294,041 km², representing a coexistence of 20%. Using a distance of 20 km, we discovered that B. distachyon was found where B. stacei was present in 31% of the cases. It is possible to verify this result taking into account data from López-Alvarez et al. (2015), where the overlapping areas of B. distachyon and B. stacei in addition to the overlapping area of the three species (B. distachyon, B. stacei, and B. hybridum) is 29.6%, as the area where the two species are coexisting. Thus, our results are consistent with those reported by previous studies.

However, by adding more species into the analysis, we wanted to test the potential of the algorithm to estimate the possible correlations that occur between species in their natural habitat. Therefore for the datasets 4, 5, 6, and 7, the rule B. retusum vs. B. phoenicoides was generated with a support of 0.5 in the first three and with one in the last, demonstrating that these two species are coexisting in a large extension of their distribution, and that its correlation is so high that it is not affected by the amount of additional species being tested in the datasets, helping us to determine that when there is a real correlation between species, these will be reflected in the support of their association. Different researchers have reported that both species can be found together in North Africa and in Western Europe, according to their geographic distribution (Catalan et al., 2015; Schippmann, 1991); even Khan & Stace, (1999), reported the association of B. phoenicoides with B. retusum rather than with B. pinnatum; which was not found in this study either.

The rule B. sylvaticum vs. B. pinnatum was found in datasets 4, 5, and 7 with a coexisting rate of 50% and confidence values of 0.84 and 0.99 in the opposite direction; recent studies of Díaz-pérez et al., (2018) revealed alleles associated with genomes of these species are present in the allotetraploid B. phoenicoides, therefore find coexisting these two species can corroborate and help explain these results. In the same way, genomes of B. pinnatum and B. sylvaticum together with B. arbuscula were potentially involved in the origins of seven allopolyploid core perennial species: B. phoenicoides, B. kawakamii, B. madagascariense, B. retusum, B. flexum, and B. bolusii (Díaz-pérez et al., 2018). Others correlations found in 20 km for datasets 4, 5 and 7, among these species can be explained perfectly, either because they are their parents, or because it is a child coexisting with their parent, for example, B. phoenicoides with B. pinnatum or B. retusum with B. sylvaticum; we also found association of B. rupestre with B. sylvaticum and B. pinnatum that for Díaz-pérez et al., (2018) this specie have alleles associated to the genome of B. glaucovirens, but also to B. sylvaticum, B. pinnatum, and B. arbuscula.

Analyzing the correlation of all the species (dataset 7), only four of them presented positive co-occurrence relationships up to a distance of 10 km (B. phoenicoides, B. retusum, B. pinnatum, and B. sylvaticum), while B. rupestre appears only in 20 km, then it is not common to find B. Rupestre near other Brachypodium species in a distance lower than 10 km. Species with a large amount of records impacted the co-occurrence analysis, and only species that appeared in almost 15% of all transactions were considered with statistical significance to create positive rules by the algorithm. Figure 3 presents the case of B. hybridum, where in dataset 6 it appeared in 15.9% of the transactions, while in dataset 7 it was 3.9%; here, we must bear in mind that dataset 6 did not present the two species with the widest distribution, which can play an important role in the number of correlations that are established. Hence, B. hybridum lost statistical importance and this is the reason why it was not included in the positive rules created by dataset 7. We found an interesting behavior by B. hybridum as this species always appeared as antecedent in the generated rules (Fig. 2A), thus demonstrating its importance in creating co-occurrence relationships. After 20 km, complex rules were generated, which included more than two species in more than 50% of all rules, such as the case of dataset 7 (see Fig. 4). This is possible because the size of the study area increases and more species can be found within, thus creating more transactions. This was observed in the case of B. sylvaticum that has co-occurrence relations with many species such as B. pinnatum, B. rupestre, B. retusum, and B. phoenicoides (see Fig. 5), due to its broad distribution in Europe, Asia, and North of Africa. When B. sylvaticum and B. pinnatum were not included in the analysis, we found close relations between B. hybridum, B. phoenicoides, and B. retusum in a distance of 20 km, and almost 40% of those co-occurrences included three species, suggesting their importance (Fig. 5A).

Finally, we found that species with a restricted distribution (endemic) like B. arbuscula and B. boissieri, or particularly present in a continent or place as is the case of B. flexum, B. bolusii from African or B. kawakami from Taiwan, B. madagascariense from Madagascar, and B. mexicanum from American, they did not present correlations with other species, which was to be expected, helping us to corroborate that the statistics are correct.

This study can play an important role in the knowledge of the associations of the species, in a level of experimental ecology, improving our capacity to predict correlations in species from big spatial data. For example, Swenson & Jones (2017) clearly demonstrate that the bioinformatics and data mining techniques are vital to analyze large volumes of biological data that involves of plant ecology.