Timber tracking of Jacaranda copaia from Amazon Forest using DNA fingerprint

1. Introduction

Great part of trade timber harvest around the world from natural forest has illegal origin [1,2,3] and even the legal origin comes from unsustainable forest logging [4,5,6,7]. This is a fact, especially for tropical forests, as Amazon Forest, where greatest part world plant and animal biodiversity lives [8]. Both illegal and legal logging contributes to biodiversity loss, decreasing the human resources potential to find medicinal and future timber sources. The opening of roads and highways within the Amazon Forest is linked to the felling of trees and the export of illegally cut timber, as easy road access contributes to entering and leaving the forest, with about 95% of deforestation occurring up to 5,5 km of roads [9]. Illegal logging is also an economic problem for the legal market, as extremely low-priced wood competes with legal logging, where costs are higher [10,11].

Due to that, many international lows have been established in countries around de world to avoid illegal timber importations, as for example, the timber regulation of the US (Lacey Act), EU (No. 995/2010), and Australia (prohibition bill), which made it illegal to import timber and its derived products (furniture, musical instruments, etc) originated from illegal extraction [11]. These rules require that importing companies declare the country and specific place of origin, names of the species of all plants contained in their products and guarantee that these wood products have legal origin and have been extracted in a sustainable way and in accordance with the laws of the country of origin. In addition to legislation targeting illegal timber trade, the Convention on International Trade in Endangered Species- CITES [12,13] allows import only from natural sources and plantations of species listed in Appendix II. In Brazil, in 2022 the Normative Instruction (No. 8, of March 25, 2022) was implemented, which “Establishes the procedures for authorizing the export of wood products and by-products of native species of natural or planted forests, aiming to complement, relatively control the export of loads of native wood within the scope of the Brazilian Institute for the Environment and Renewable Natural Resources (IBAMA), Normative Instruction No. 21, of December 24, 2014, Normative Instruction No. 17, of December 1, 2021, and Ordinance No. 8, of January 3, 2022.” and an Ordinance (No. 8, of January 3, 2022, amended by Ordinance No. 46, of January 6, 2022), which “Establishes, within the scope of IBAMA, the Single Consent Platform of Brazil - Pau Brasil for use in foreign trade activities involving products and by-products of biodiversity”. Three lists are found in the country, according to Decree nº 3.607, of September 21, 2000, indicating the species that are already considered threatened with extinction, and those that appear in the second list, which includes species that are not currently threatened, but that are at serious risk unless trade traffic is strictly observed and laws and regulations are strictly followed. In the third annex are the species that the Brazilian agency itself has listed, so that there is greater care in the exploration, which is restricted or even blocked for better export control. For example, Jacaranda spp. are listed in the first annex [13], so its timber and derived products originated from legal and sustainable forest logging in natural forest can be exported.

The biggest problem about the illegality of the trade is in the documentation of the origin of the wood, having falsified certificates and documentation and with missing or false information [14]. According to the Brazilian Federal Police, in 2021, 90% of the wood extracted from the Amazon Forest had an illegal origin [15,16]. Such timber tracking methods are required, but methods based only on documentation are sensitive to manipulation and forgery. Therefore, exporting companies and institutions responsible for controlling the origin of imported timber, such as customs, federal police and Interpol need reliable tools to prove and confirm the declared origin of wood and its derivative products, traded internationally.

To timber tracking the species, country and specific site origin of imported timber, many methods have been tested, such as chemical analysis or mass spectrometry species timber differences [17,18,19,20], stable isotopes [21], near-infrared spectroscopy [22], wood anatomy [23,24], and DNA fingerprint [1,2,3,11,25,26,27,28,29,30,31,32,33]. Results of different methods have been shows strong potential for species determination, country and site origin, in special using DNA fingerprint [1,2,3,11,25,26,27,28,29,30,31,32,33]. Methods such as wood anatomical, isotopic, and spectrometric methods are limited to all species, country and site origin determination, due to variations in tissue type, individual sample age, individuals and population genetic differences or environmental influences on timber composition [31]. Thus, molecular methods that allow correct identification of tree species and tracking of timber origin are essential for controls on the legality of timber by public authorities, industry, and trading companies.

In the current work we investigated the use of the DNA fingerprint method to track the intensive and high value wood of the Neotropical pioneer tree Jacaranda copaia (Aubl.) D. Don. (Bignonaceae). The genus Jacaranda sp. present about 53 species, all Neotropicals, most of them are mainly found in the Cerrado around the Amazon Forest, where J. copaia is the only species of this genus widely distributed in the Amazon [34]. The species occurs from northern to western South America, from Belize to Bolivia, Brazil, French Guiana, and Peru [34]. In Brazil, the species is found in the states of Acre, Amapá, Amazonas, Maranhão, Mato Grosso, Pará, Rondônia, and Roraima, and populations usually have more than one tree per hectare [34]. This is a fast-growing species, with an average annual increase in diameter at breast height (DBH) of 2.05 cm and height of 1.98 m, and in adulthood, trees have a straight stem, reaching 106 cm in DBH and 45 m in height [34,35]. Its wood is light and used for furniture [35,36]. The species is hermaphroditic, self-incompatible and about 40 species of bees, butterflies and hummingbird wasps were detected as potential pollinators, although Euglossa spp. and Centris spp. bees was detected as the main pollinators in the Tapajós National Forest, Brazil [37]. Reproduction is primarily mediated by outcrossing, ranging from 0.942–0.993 [38,39]. The fruits can have up to 250 seeds and the winged seeds are dispersed by the wind [37]. Mating between related trees has been reported and attributed to the fact that populations have intrapopulational spatial genetic structure, ranging 50–500 m, due to the short-distance seed dispersal (up to 100 m), but also due short pollen dispersal distance, ranging from 34–90 m [39].

Here, samples from 43 populations of J. copaia collected in Bolivia, Brazil, French Guiana and Peru were used to quantify genetic diversity and population structure, and to test the ability of chloroplast, mitochondrial and nuclear single nucleotide polymorphism (SNPs) markers to determine the geographic location of origin of wood of the species harvested from the Amazon Forest.2. Materials and methods

2.1. Sampling

Were collected cambium or leave samples from 832 trees from 43 natural populations (2–31 individuals) in the Amazon rainforest of four countries (Bolivia, Brazil, French Guiana and Peru) and all trees sampled were georeferenced with GPS usage (Table 1, Figure 1). All samples were stored in a labelled plastic bag with silica gel. The samples were collected by the Institut National de la Recherche.Agronomique- INRA together with the forest authorities (Office National des Forêts, ONF) in seven populations in French Guiana (2–30 individuals per site); In Brazil, samples were collected in national forests, extractive reserves, ecological stations, and national parks with the support of Chico Mendes Institute of Biodiversity (Brazilian governmental institution), totaling 12 sites (4–31 individuals per site); In Peru, samples were collected in national forests, extractive reserves, ecological stations, national parks, farms, and forest concessions, totaling 19 sites (2–30 individuals per site); In Bolivia samples were collected in farms and forest concessions from five sites (5–29 individuals per site). However, due to the low sample size in five populations (2 trees), for genetic analysis, these individuals were grouped with the closest population (Table 1), as well as we also divided the Brazilian samples in West (six populations) and East (six populations) origin and Peru in North (nine populations) and South (eight populations) origin. The minimum distance among sampled trees within populations was 50 m and the distance among populations ranged from 23–2648 km. After sampling, all collected plant material was stored and dried in silica gel. All samples were registered in a database at the Thünen Institute (SampleDataBase, Grosshansdorf, Germany).

2.2. DNA Extraction and SNPs Analysis

Cambium and leaf samples collected in Brazil were sent the Laboratory of Population Genetics and Forestry of São Paulo State University in Ilha Solteira, Brazil (UNESP-FEIS) for DNA isolation. Samples collected in French Guyana, Peru and Bolivia were sent to Thünen Institute facilities in Großhansdorf, Germany, for DNA isolation. DNA isolation from leaf and cambium was carried out according to Dumolin et al. [40]. The samples were screened for 128 SNP and INDEL markers using the MassARRAY^® iPLEX^TM genotyping, where 113 were nuclear SNPs (nSNPs), 11 cloroplastidial SNPs (CpSNPs), and four mitochondrial SNPs (MtSNPs), all selected for genetic tracking analysis due to show a minimum amplification rate of 95% [41].

2.3. Genetic Diversity and Population Differentiation

The genetic diversity for nSNP was determined for each population and country by the total number of alleles ($K$), percent of polymorphic loci ($P\%$), observed ($H_o$), and expected ($H_e$) heterozygosity. The mean fixation index ($F$) was estimated to quantify the inbreeding within each population or country and the statistical significance of the $F$ values was determined permuting alleles among individuals. For CpMtSNP, the genetic diversity was determined for each population by $K$, percent of polymorphic loci ($P$%), and $H_e$. We also estimated the percent of polymorphic loci ($P\%$), observed for all nCpMtSNP loci. The genetic differentiation ($F_{st}$) was estimated among all populations, populations within country, among countries, and pairwise populations for all nCpMtSNP markers. These analyses were carried out using the GDANT 1b software [42]. The pairwise $F_{st}$and spatial distance among populations was used to investigate the isolation by distance (IBD) gene dispersal, using the Spearman correlation coefficient (ρ). We also assessed the IBD by the estimate the spatial genetic structure (SGS) based on calculation of the coancestry coefficient (${\theta }_{ij}$) described in Loiselle et al. [43], between mean pairs of individuals within 10 distance-classes (0–198, 199 –441, 442–655, 656–860, 861–1068, 1069–1409, 1410–1838, 1839–2770, 2771–2478, and 2479–2646 km) determined using the same number of pairs per classes, and using the SPAGEDI 1.5 software [44]. The statistical significance of the average ${\theta }_{ij}$ of each distance class was derived by comparing the limits of the confidence interval at 95% probability for the average ${\theta }_{ij}$ for each distance class, estimated permuting (1000 times) genotypes between distances classes, using the SPAGEDI software.

2.4. Bayesian Clustering Analysis

Bayesian clustering method [45] implemented in the software STRUCTURE v.2.3.4 [46] was used for the 128 SNPs makers to check for the number of genetic groups/populations. We set the length of burn-in and Markov chain Monte Carlo to 10.000 and tested K values from 2–10 for 30 times. We used the admixture model with correlated allele frequencies. The optimal number of genetic clusters was estimated with the ΔK method [46]. For each optimal K, data from the 30 STRUCTURE runs was permuted with CLUMPP v.1.1.2 [47] to obtain the final Q values for each individual. For determine the optimal number of genetic clusters was determined with the ΔK method described by Evanno et al. [46] using the software CLUMPAK [48].

2.5. Genetic Assignment

Bayesian method [49] implemented in GeneClass 2.0 [50] was used to assign group (population) and individuals to its population and country. Both group and individual assignment were caried out for all nCpMtSNP loci, and the most likely group determined by the highest score by the Bayesian criteria was used as an indicator of the power of the markers to compute the proportion of correctly assigned groups or individuals in self-assignments tests [51]. Here the individuals of the reference data were self-classified to the sampled groups (populations and countries) using the leave-one-out approach [52]. We also estimated the overall, ≥ 80% and ≥ 95% score rate of group and individuals assigned to correct origin population.

3. Results

3.1. Genetic Diversity

The total number of alleles ($K$) among populations for 113 nSNPs ranged from 122–226 alleles (mean of 183) and for 15 CpMtSNPs from 14–19 alleles, mean of 15.9 alleles (Table 2). Percent of polymorphic loci ($P\%$) for 113 nSNPs was higher (ranged among populations from 23–88.5%) than for 15 CpMtSNPs (0–26.7%) and 128 nCpMtSNPs (20.3–78.1%, mean of 55.9%). The observed heterozygosity ($H_o$) ranged from 0.023–0.343 (mean of 0.178), the expected heterozygosity ($H_e$) for nSNPs ranged among populations from 0.028–0.348 (mean of 0.204), where for CpMtSNPs ranged from 0–0.16 (mean of 0.024). The values of $H_e$ for nSNPs were high (≥ 50% of the maximum value for biallelic loci, 0.5) for 21 populations, where the highest values were observed for Brazil and Bolivia populations. The mean intrapopulation fixation index ($F$) ranged from -0.086–0.295 and was significantly (P< 0.05) higher than zero in six of the 38 populations, suggestion inbreeding. At countries level, for nSNPs, $K$ and $P\%$were highest in Brazil for nSNPs (226, 100%, respectively) and CpMtSNPs (30, 100%, respectively), and lowest in French Guiana (nSNPs= 200, 76.1%, respectively; CpMtSNPs= 18; 20%, respectively). Bolivia presented the highest $H_o$ (0.319) and $H_e$ (0.359) values and lowest $F_{IS}$ (0.113), where French Guiana presented the lowest $H_o$ (0.095) and $H_e$ (0.192) values and the highest $F_{IS}$ (0.506). The $H_e$ for CpMtSNPs was also highest in Brazil (0.251) and lowest in French Guiana (0.013). Spearman correlation (ρ) for nSNPs was significantly higher than zero (ρ= 0.36, P= 0.026) between sample size ($n$) and total number of alleles per population ($K$), between $K$ and observed heterozygosity ($H_o$: ρ= 0.675, P= 0), and between $K$ and expected heterozygosity ($H_e$: ρ= 0.691, P= 0), showing that high sample sizes increases allele diversity and populations with higher number of alleles presented higher heterozygosity’s (Supplementary Materials, Table S1).

3.2. Population Differentiation

To determine differences between nuclear SNPs, chloroplastidial and mitochondrial markers (CpMtSNPs) in estimating genetic differentiation ($F_{st}$) between samples, we compared the results between all 128 nCpMtSNPs, 113 nSNPs, and 15 CpMtSNPs (Table 3). The results showed that the $F_{st}$ among all populations and countries was higher for CpMtSNPs (0.942 and 0.695, respectively) than for all nCpMtSNPs (0.484 and 0.295, respectively) and nSNPs (0.415 and 0.233, respectively). Among populations within Brazil, Peru and Bolivia countries, the $F_{st}$ was also higher for CpMtSNPs (0.974, 0.741, and 0.735, respectively) than for all nCpMtSNPs (0.3, 0.466, and 0.142, respectively) and nSNPs (0.224, 0.456, and 0.107, respectively). In contrast, among populations within French Guiana, the $F_{st}$ was higher for nCpMtSNPs (0.12) and nSNPs (0.117) than CpMtSNPs (0.004). The results also showed that the estimates of $F_{st}$ for nCpMtSNPs, nSNPs, and CpMtSNPs was generally highest within Peru and Brazil, where the number of sampled populations was highest (Peru= 17, Brazil= 12), indicating that the number of sample populations is important to detect high population differentiation within countries. The highest $F_{st}$ for CpMtSNPs and nCpMtSNPs than for nSNPs indicated that the combination of 113 nuclear and 15 cytoplasmatic SNP markers increase the capacity of detect genetic differences among populations.

3.3. Isolation by Distance

The Spearman correlation coefficient ($\rho$) between pairwise $F_{st}$ and spatial distance between populations, based on 128 nCpMtSNP loci was significantly (P< 0.001) higher than zero ($\rho$= 0.506), indicating an isolation by distance (IBD) genetic pattern (Figure 2). The estimate of pairwise coancestry for nSNP and nCpMtSNP loci between individuals within ten distance classes was significantly higher than zero up to the distance of 655 km, where for CpMtSNP loci was significant up to 1068 km, and non-significant or significantly lower than zero in the other distance classes (Figure 3).

3.4. Bayesian Cluster Analysis

The results of the STRUCTURE analysis showed that based on delta K, the best representation of the data is obtained with a K of 4 or 9 different genetic groups (Figure 3). With 5 groups, there is less mixture of gene pools within individuals. With K= 5, 85% of all individuals had a gene pool membership coefficient above 0.9, while for K= 4, only 76% of all individuals had gene pool membership coefficients above 0.9. The STRUCTURE analysis based on CpMtSNPs showed a clear distinction between French Guiana and North Brazil populations (red) from North Peru, South Brazil (6BW-An; 7B-Ja), East Bolivia (34Bo-Ri, 35Bo-De) (green and white), and from South Peru (29P-In; 30P-Md, 31P-Ca, 32P-Am, 33P-Ta), Southwest Brazil (9BW-Xa, 10B-Co) and Bolivia (36Bo-Vi, 37Bo-Ch; 38Bo-Ma) (blue) (Figure 3A). The STRUCTURE analysis based on nSSR, there are a few locations with mixtures of individuals from different groups (e.g. Flona do Jamari (7BW-Ja), Resex Tapajos-Arapins (14BE-Ar), and ESEC Maraca (5BW-Ma) in Figure 3B). The samples from French Guiana are very similar among them-selves, and created one genetic group (K2 in Figure 3B). We found one genetic group (K3) which we only found in two Brazilian populations (8BW-Jr. 11BE-Am) and another genetic group (K1) distributed more in the middle of the Amazon basin (Figure 3B).

3.5. Genetic Assignment

The grouped sample assignment test for all nCpMtSNPs was able to self-assignment of 100% of population (score of 100%) to the correct origin country and population (Table 4). The individual sample assignment test for all nCpMtSNPs was able to self-assignment of 30.5–100% (mean of 85.7%) of individuals to the correct origin population (score ranging from 59.3–100%). The rate of correct individual assignment for scores ≥ 80% and ≥ 95% were 71 and 62.3%, respectively. The assignment test for individuals for was able to self-assignment among 56.4–98.6% (mean rate of 98.7%) of country origin (mean score ranging from 74.9–99.2%). The rate of correct individual assignment to the correct country source with scores ≥ 80% and ≥ 95% ranged from 49.2–97.6% and 31.9–97.6%, respectively. In the cases of wrong self-assignments, the scores ranged from 45.5–100%, where the distance between the correct to assigned population ranged from of 29–591 km.

Spearman correlation (ρ) for all SNPs was significantly higher than zero (P≤ 0.015) between correct country assignment and $H_o$ and $H_e$ (ρ: ranging from 0.809–0.905), as well as between all SNPs and $K$(ρ= 0.74, P= 0.04) (Supplementary material, Table S2). Spearman correlation for all SNPs and nSNPs was also significantly higher than zero (P≤ 0.026) between correct population assignment and $K$, $H_o$, and $H_e$ (ρ: ranging from 0.36–0.691) and significantly lower (P≤ 0.003) than zero between correct population assignment and $F$ for all SNPs (ρ= -0.473) and nSNPs (ρ= -0.504).

4. Discussion

4.1. Genetic Diversity

The ilegal timber trade, both species-specific and illegally sourced, has become a major problem for tropical forests. Control of the chain of custody for timber originating from the Amazonia Forest, as well as timber from other biomes and regions of the planet, is urgent to verify the species, countries, and place of origin declared in the transport and export documentation and restrict or even irradicate illegal logging. Our study shows the potential of DNA fingerprinting to track the country and population origing for timber, as well as to follow and verify the chain of custodia for J. copaia timber products.

Our study displayed only moderate levels of genetic diversity for the nSNP and low levels for CpMtSNP markers for J. copaia. However, the observed ($H_o$) and expected ($H_e$) heterozygosity were generally moderate and high, respectively, and are within the reported pattern identified for other Neotropical trees using nSNPs markers, where $H_o$has been observed ranging from 0.017–0.39 and $H_e$ ranging from 0.02–0.371 (Table S4). The genetic diversity of J. copaia was especially low in some North Peru populations. The $H_e$ values were also higher for nSNPs than for CpMtSNPs due to the higher polymorphism of nSNP loci.

The estimate of intrapopulation fixation index ($F$) indicated inbreeding in some populations. Inbreeding has also been observed for other Neotropical trees using SNPs markers (Table S4). However, because J. copaia is self-incompatible [37], and our within population samples were taken from geographically distant trees, the observed inbreeding is very likely an artefact of the Wahlund effect [53] due to mixtures of samples from different subpopulations.

4.2. Population Genetic Differentiation

The presence of population genetic differentiation and intrapopulation spatial genetic structure (SGS), or the occurrence of isolation by distance patterns (IBD), is key to assigning timber to different origins [11,29,54,55,56], including country, population, and regions within countries. High genetic differentiation among the different genetic groups increases the success of genetic assignment [54]. Our results showed strong genetic differentiation among all populations and between countries, as well as SGS and a pattern of IBD. The $F_{st}$ for all SNPs, nSNPs, and CpMtSNps was generally high, especially for CpMtSNps, indicating that SNP markers have a strong capacity to determine the correct timber population of origin. The genetic differentiation was higher among all populations than among countries. Within countries for all SNPs markers, $F_{st}$ was highest among the Brazilian and Peru population. These results indicate strong potential for assigning timber origins between populations and countries.

4.3. Genetic Assignment and Practical Application

The results indicate a high power to correctly assign groups at both levels: between countries and populations (100%). Our results confirm that this approach offers the possibility of high levels of success in the grouping of timber samples based on their origin populations. This success can be attributed to the general wide genetic differentiation among countries and populations. Even when the differentiation was only moderate ($F_{st}$< 0.2) (French Guiana, Bolivia), differences in allele frequencies between analyzed groups were enough to produce high scores for the correct origin. However, it is important to emphasize that this is true when the reference data set contains genetic data of the timber’s original population. Therefore, the success in determining the origin of J. copaia timber, as well as other species, greatly depends on the reference data set developed by the authorities involved in controlling the origin of timber (reference populations). This data set should have genetic data from the population where harvesting is legally practiced.

The results for individual assignment test for all SNPs and nSNPs were lower than group assignment, but also showed a high power of correct assignment among countries and populations. Peru and its regions (North and South) showed a low rate of assignment. According to the Spearman correlation (ρ) for countries and populations for assignment of all SNPs and nSNPs versus $K$, $H_o$, and $H_e$ for population level for assignment of nSNPs versus $F$, the $H_o$ and $H_e$, following by $K$ are the indices with the most influence on the success of individual assignment tests. The Peru North population generally showed the lowest $H_o$ and $H_e$ values as well as the highest levels of inbreeding, which can explain the low assignment scores. In general, the results for the individual assignment test for all SNPs and nSNPs indicate that this approach has the power to determine the specific site of timber origin. In using CpMtSNP loci markers, it is necessary to develop a large number of loci to improve both group and individual assignment tests. Similar results of better group that individual assignment tests were reported for Hymenaea sp. [1,2] and Handroanthus sp. [57]. For Swietenia macrophylla, based on nSSR loci, the assignment test was higher at the country (82%) than population (71%) level [11]. For the Malaysian Gonystylus bancanus tree, the self-assignment rate using a set of 16 nSSRs was lower (55%) than that observed here at the population level [58]. For SNP data of Entandrophragma cylindricum, the assignment at the country level ranged from 66–74%, depending on the assignment method [53]. Many other studies on tropical, African, and European tree species have been developed using microsatellite and SNPs markers to track timber and the main conclusion is that, due to the presence of SGS and genetic differentiation, DNA fingerprinting is the most effective tool to track the country and population of origin [54,55,56,57,58,59,60,61,62,63,64,65].

Finally, we recommend that the timber sector add such genetic controls as an independent audit beyond the paper-based proof of the chain of custody [11]. It is important to note that the power of the genetic reference data to detect false declarations reaches 100% if more than one sample with the same declaration is tested.

4.4. Conclusion

The genetic differentiation ($F_{st}$) among all populations, countries, and regions within countries was generally high, especially based on CpMtSNP loci. Furthermore, there is a strong isolation by distance pattern among populations, favoring the tracking of group or individual samples to the correct site. For self-assignment tests, we were able to obtain 100% accuracy to the country, population, and region origin for all samples using all SNPs and nSNPs. Our results show that the use of all SNP or nSNP markers are suitable to correctly determine the country and population of origin for J. copaia timber, thus offering a very useful tool for customs and local and international police. The J. copaia reference database of our study represents a robust assignment tool available to timber companies or government agencies to test and validate origin declarations. We recommend the use of the method described herein for other native tropical species, since it is highly effective in identifying the origin of wood, thus helping the police and other relevant agencies in the definition of illegally deforested areas, as well as unsustainable extraction.