Signatures of selection for bonamiosis resistance in European flat oyster (Ostrea edulis): New genomic tools for breeding programs and management of natural resources

Abstract The European flat oyster (Ostrea edulis) is a highly appreciated mollusk with an important aquaculture production throughout the 20th century, in addition to playing an important role on coastal ecosystems. Overexploitation of natural beds, habitat degradation, introduction of non‐native species, and epidemic outbreaks have severely affected this important resource, particularly, the protozoan parasite Bonamia ostreae, which is the main concern affecting its production and conservation. In order to identify genomic regions and markers potentially associated with bonamiosis resistance, six oyster beds distributed throughout the European Atlantic coast were sampled. Three of them have been exposed to this parasite since the early 1980s and showed some degree of innate resistance (long‐term affected group, LTA), while the other three were free of B. ostreae at least until sampling date (naïve group, NV). A total of 14,065 SNPs were analyzed, including 37 markers from candidate genes and 14,028 from a medium‐density SNP array. Gene diversity was similar between LTA and NV groups suggesting no genetic erosion due to long‐term exposure to the parasite, and three population clusters were detected using the whole dataset. Tests for divergent selection between NV and LTA groups detected the presence of a very consistent set of 22 markers, located within a putative single genomic region, which suggests the presence of a major quantitative trait locus associated with B. ostreae resistance. Moreover, 324 outlier loci associated with factors other than bonamiosis were identified allowing fully discrimination of all the oyster beds. A practical tool which included the 84 highest discriminative markers for tracing O. edulis populations was developed and tested with empirical data. Results reported herein could assist the production of stocks with improved resistance to bonamiosis and facilitate the management of oyster beds for recovery production and ecosystem services provided by this species.


| INTRODUC TI ON
The flat oyster Ostrea edulis has been one of the most appreciated seafood products since the Roman times and represented an important aquaculture species in Europe throughout the 20th century.
Additionally, oyster reefs play an important ecosystem role providing multiple benefits such as water filtration, food and habitat for other animals, shoreline stabilization, carbon burial, nutrient regeneration, and coastal fisheries (Grabowski & Peterson, 2007). Ostrea edulis habitats are endangered in Europe and have been identified as a priority for protection and restoration in the European MPAs (OSPAR Commission, 2011).
Overexploitation of natural beds, introduction of non-native oysters, adverse effects of climate change and specific epidemic outbreaks, such as the "shell disease" caused by the fungus Ostracoblabe implexa (Alderman & Jones, 1971), the protist Marteilia refringens (Grizel et al., 1974) and the protozoan parasite Bonamia ostreae (Pichot, Comps, Tigé, Grizel, & Rabouin, 1980), have determined the decline of flat oyster wild production down to marginal figures in many European areas (a peak of 32,995 t in 1961 to 3,120 t in 2016; http://www.fao.org/fishe ry/stati stics/ global-produ ction/ es). Bonamiosis, endemic in many areas of flat oyster distribution, is currently the main concern for flat oyster production and the main contributory factor to its decline (Engelsma, Culloty, Lynch, Arzul, & Carnegie, 2014); therefore, effective control of this parasitism is considered essential for the recovery of the species. Successful breeding programs aiming to obtain resistant strains have been accomplished in several oyster species including different pathogens, such as the ostreid herpesvirus 1, Haplosporidium nelsoni, Perkinsus olseni, Roseobacter crassostreae, Marteilia sydneyi or B. ostreae, by selecting survivors after long exposure in heavily affected areas (Beattie, Davis, Downing, & Chew, 1988;Dégremont, Nourry, & Maurouard, 2015;Dove, Nell, & O'Connor, 2013;Ford & Haskin, 1987;Frank-Lawale, Allen, & Degremont, 2014;Lynch, Flannery, Hugh-Jones, Hugh-Jones, & Culloty, 2014;Ragone Calvo, Calvo, & Burreson, 2003). The procedure usually involves the evaluation of populations in the field to identify resistant individuals to be used as breeders for the next generation and can take up to 4-5 years. This has limited its efficacy due to the inability to maintain constant selection pressure when disease exposure is low or absent. Accordingly, detection of quantitative trait loci (QTL) and markers associated with resistance to bonamiosis could greatly accelerate the selection of breeders making the process much more efficient. Such resistance markers have been successfully found for some diseases in mollusks (He, Yu, Bao, Zhang, & Guo, 2012;Meistertzheim et al., 2014;Nie, Yue, & Liu, 2015;Nikapitiya et al., 2014;Normand et al., 2014;Raftos, Kuchel, Aladaile, & Butt, 2014;Schmitt, Santini, Vergnes, Degremont, & de Lorgeril, 2013), and particularly, several QTL for resistance to bonamiosis have been identified in flat oyster (Harrang et al., 2015;Lallias et al., 2009). These studies typically used low-medium density genetic maps and a small number of families and were able to identify some QTLs explaining up to 17% of the phenotypic variance in the trait (Harrang et al., 2015;Lallias et al., 2009). However, this information, although useful, should be considered preliminary since it is far from elucidating the genetic architecture of bonamiosis resistance in the oyster populations of the Atlantic area, considering the low number of analyzed families. Moreover, the detected QTL should be considered as suggestive, since they were identified at chromosome level (p < 0.05) and in only one or two families.
Recently, a gene expression study using a flat oyster hemocyte oligo-microarray (Pardo et al., 2016) enabled the identification of specific genes and enriched pathways differentially expressed (DE) between two flat oyster stocks with differential susceptibility to bonamiosis: a naïve one (NV) from a free-bonamiosis area and a long-term affected one (LTA), which had shown some resistance to B. ostreae (da Silva, Fuentes, & Villalba, 2005;Villalba, da Silva, & Fuentes, 2007). The LTA stock showed a substantially more intense and rapid response to infection than the NV, with the activation of several signal transduction pathways early after challenge with B. ostreae, likely representing an important factor for the resistance observed (Ronza et al., 2018). Furthermore, cell-extracellular matrix (ECM) interactions and proteases inhibitors were essential in the defense response of LTA oysters, while the broad down-regulation of histone-related genes in NV oysters suggested a fundamental role of these proteins during the infection process. The set of DE genes detected in this study and the proper oyster database, where a substantial number of gene-linked markers are available (Pardo et al., 2016), represent a useful repository to build upon the identification of genetic markers associated with candidate genes linked to bonamiosis resistance.
Marine species typically show low structuring due to high effective population sizes and high gene flow facilitated by the absence of physical barriers, which lead to high genomic homogenization across populations (Vandamme et al., 2014;Vilas et al., 2015). However, habitat shifts (ecotones) and oceanic currents/ fronts have been demonstrated to act as practical barriers (Blanco-González, Knutsen, & Jorde, 2016;Vera et al., 2016), and particularly, natural selection in response to abiotic and biotic factors has been demonstrated to actively shaping specific genomic regions related to adaption even in a context of high gene flow (Vandamme et al., 2014;Vilas, Bouza, Vera, Millán, & Martínez, 2010;Vilas et al., 2015). The ability to distinguish between neutral and adaptive genetic variation has become a central issue in evolutionary biology, as it would allow the understanding of population structure in both historical/demographic K E Y W O R D S Bonamia ostreae, candidate genes, disease resistance, divergent selection, genetic traceability, Ostrea edulis, SNP array and adaptive terms (Bernatchez, 2016;Nielsen, Hemmer-Hansen, Larsen, & Bekkevold, 2009), thus being an essential information for sustainable fisheries management. In turn, genetic diversity in the wild represents the raw material for the foundation of aquaculture broodstocks and consequently a reference to identify selection signatures for targeted traits in farmed populations which can be detected via genome scanning (e.g., Liu et al., 2017). Several European flat oyster populations have been subjected to the selective pressure of bonamiosis for many generations, while other oyster populations have never been in contact with the parasite. This scenario raises the hypothesis that specific genomic regions related to parasite resistance have been under selection and therefore modified by this selective pressure. In fact, increase of resistance to B. ostreae due to natural selection associated with long exposure of oyster populations to this parasite was supported by field experiments (Elston, Kent, & Wilkinson, 1987;da Silva et al., 2005). Consequently, a genetic comparison of several LTA vs. NV oyster populations originating from different sources and geographical regions, thus being essentially considered as independent replicates, could aid to identify those genomic regions subjected to selection.
In the current study, we characterized and compared two sets of three LTA vs. three NV flat oyster populations collected throughout the European Atlantic area and involving the three genetically differentiated regions previously reported (Vera et al., 2016) to identify candidate genomic regions related to bonamiosis resistance by using two different strategies: (a) a panel of validated genetic markers from DE candidate genes between LTA and NV populations; (b) a large panel of SNPs covering the whole genome at a high density (1 SNP per 80 kb on average). A consistent set of 28 closely linked genetic markers from a putative single genomic region was identified which is suggestive of a strong QTL differentiating both population sets.  Figure 1). These three flat oyster genetic regions seem to show high connectivity due to larval advection according to their low differentiation (global F ST = 0.0079), and this low but significant differentiation has been suggested to be driven by oceanic fronts present in the region (Vera et al., 2016). Three of the oyster beds, Ortigueira (ORT), Quiberon (QUI), and Rossmore (ROS), had been exposed to B. ostreae since the early 1980s (McArdle, McKiernan, Foley, & Jones, 1991;Pichot et al., 1980;Polanco, Montes, Outon, & Melendez, 1984)

| DNA extraction and SNP validation and genotyping
Total DNA was extracted from gill tissue using the e.Z.N.A. E-96 mollusk DNA kit (OMEGA Bio-tech), following the manufacturer recommendations. Two different strategies were followed for single nucleotide polymorphisms (SNPs) selection and genotyping to identify genetic markers associated with resistance to bonamiosis. The first one was focused on the search for SNPs associated with candidate (denoted "CAND" markers) differentially expressed genes (DEGs) between LTA vs. NV oysters (Ronza et al., 2018). The flat oyster transcriptome database available for the species (Pardo et al., 2016) and the list of the 837 DEG previously reported by Ronza et al. (2018) were the starting point for SNPs selection. This selection was performed according to technical features (SNP flaking regions-100bp-without additional SNPs) and the functional relevance of the genes related to bonamiosis response (Ronza et al., 2018). Selected SNPs were validated on a MassARRAY platform (Sequenom) at the USC node of the Spanish National Centre of Genotyping (CeGen ISCIII). The distinct mass of the extended primer identifies the SNP allele. MALDI-TOF mass spectrometry analysis in an Autoflex spectrometer was used for allele scoring. All the 199 oysters sampled were genotyped for CAND SNPs (Table 1).

| Genetic diversity and population structure
The mean number of alleles per locus (Na), observed (Ho), and expected (He) heterozygosity, departure from Hardy-Weinberg F I G U R E 1 Geographical situation of the Ostrea edulis locations in the present study. Dotted lines and roman numbers show the OSPAR regions. The three population genetic clusters previously described by Vera et al. (2016) are also represented (orange circle: Spanish cluster; blue ellipse: Irish/British/French cluster; gray ellipse: Dutch/Danish cluster). Locations long-term exposed to Bonamia are indicated in red, while naïve locations are indicated in green. Location codes are shown in Table 1 TBay LRy LIM I II III V IV ORT QUI ROS equilibrium (HWE) and the estimators of Wright (1951) F-statistics at each locus were calculated using GENEPOP v4.0 (Rousset, 2008) and ARLEQUIN v3.5 (Excoffier & Lischer, 2010). Linkage disequilibrium was checked with exact tests to look for genotyping association using GENEPOP v4.0, and they were specifically applied to check linkage among loci related to bonamiosis resistance.
Pairwise F ST values between oyster beds were calculated with ARLEQUIN v3.5 using 10,000 permutations to test for significance.
Genetic subdivision was explored using an approach that does not require defining populations a priori. In particular, a model-based clustering approach, implemented in STRUCTURE v2.3.4 (Pritchard, Stephens, & Donnely, 2000), which constructs a number of genetic clusters (K) from a collection of individual multilocus genotypes and estimates for each individual the fraction of its genome that belongs to each cluster, was applied. A burn-in period of 10,000 steps and 100,000 Monte Carlo replicates for different values of K ranging from 1 to 7 was used. For each K tested, five replicates were run to assess their reliability. After checking that K = 1 was not the most probable scenario (see Results), the STRUCTURE results were reanalyzed with the program STRUCTURE HARVESTER v0.3 (Earl & vonHoldt, 2012) to estimate the best K value between 2 and 7 following the ΔK method described by Evanno, Regnaut, and Goudet (2005). The software CLUMPP v1.1.2 (Jakobsson & Rosenberg, 2007)

| Divergent selection for bonamiosis resistance: outlier detection and gene mining
Different statistical strategies were used to identify outlier loci subjected to selection. BAYESCAN v2.01, a method based on a Bayesian approach to identify outliers using F ST values under an inland model (Foll & Gaggiotti, 2008), was run using 20 pilot runs, prior odds value of 10, with 100,000 iterations, burn-in of 50,000 iterations and a sample size of 5,000. Two different scenarios were tested with BAYESCAN: identifies outliers using the FDIST F ST method using a maximum likelihood approach (Beaumont & Nichols, 1996), but has the advantage of incorporating a priori information regarding population structure. This method was applied to investigate outlier loci under the two scenarios also tested with BAYESCAN: (a) ARL IND and (b) ARL POOL , but additionally including a third one, ARL HIER , where a hierarchical island model was applied grouping oyster beds by their resistance to bonamiosis. All ARLEQUIN runs were performed with 25,000 simulations, 10 groups, and 100 demes per group. Strategies to identify outliers under selection can produce type I (false positive) or type II (false negative) errors.
BAYESCAN follows a more conservative approach than ARLEQUIN (Narum & Hess, 2011), and accordingly, loci with a false discovery rate (FDR, q-value) <0.05 with BAYESCAN were considered consistent outliers. ARLEQUIN has the tendency to produce a higher proportion of false positives, and therefore, loci at p < 0.01 were considered consistent outliers, while those at p < 0.05 suggestive outliers. Markers shared (or excluded) by the different methods and scenarios were assumed to be the most reliable ones (see Results).
All sequences including divergent outlier SNPs were compared against the O. edulis database (Pardo et al., 2016)

| Testing population discriminative and individual allocation capability of SNPs
Outlier SNPs putatively related to other factors than bonamiosis resistance were identified by discarding, from the total outliers detected with all strategies, those detected in the LTA vs. NV comparisons, putatively associated with bonamiosis since this was the factor used for grouping the set of populations analyzed covering the whole Atlantic area. This subset of outlier SNPs was used to estimate their capability for discriminating populations and for allocating individuals to their populations. These outliers were ordered according to global F ST from higher (more discriminant) to lower (less discriminant) values, and then, subsets of markers were set up by progressively adding one by one SNPs following the F ST rank. The USEPOPINFO option of the program STRUCTURE 2.3.4 was used to define the six reference populations studied. Then, for each SNP subset and individual tested, the admixture coefficient q was obtained for each of the six reference populations. With this information, we estimated the proportion of individuals assigned to their population according to the highest q (among the six calculated) and the highest one was considered the correct allocation. Using this information, we also estimated the average q value for each population and subset of SNPs tested.
Finally, an assignment test using the Bayesian method of Rannala and Mountain (1997) implemented in GENECLASS 2.0 (Piry et al., 2004) was performed to check for the power of the SNPs to assign individuals to their populations.

| SNP genotyping and genetic diversity
Among the 837 DEGs identified by Ronza et al. (2018), 411 contained reliable SNPs in the flat oyster database and 390 of them were consistently annotated. After filtering using criteria

| Signals of selection for CAND markers
We looked for outliers under selection by comparing the six oys-

| Signals of selection for ARRAY markers
The number of outliers detected with BAYESCAN was much lower

| Genomic information of candidate outliers: linkage disequilibrium (LD) and gene mining
Considering the density of the flat oyster SNP panel, we hypothesized that some outlier loci may be located at the same genomic region associated with divergence between LTA and NV populations. This would be especially true if these regions (QTLs) showed a strong effect on bonamiosis resistance, determining a significant selective sweeping around the favorable mutations.
Thus, we checked for linkage disequilibrium (LD) among the 22 outlier POOL-S loci considering all oyster samples as a single population, but also within each of the six populations evaluated. There is a trade-off when using these strategies; using the whole dataset, LD is increased due to admixture, and thus, we have lower precision for pinpointing associated genomic regions. However, it has the advantage of increasing sample size and, thus, the statistical power to detect associations; however, testing LD in populations at HWE enables fine mapping of the associated genomic regions, but it has the disadvantage of decreasing sample size and therefore statistical power. Due to the lower structure of flat oyster populations and the small sample size (~16) in our study, the first strategy should show higher resolution. Additionally, we performed the same analysis with 50 randomly selected loci from the whole 14,065 loci to be taken as an unlinked dataset reference. Using the 22 outlier POOL-S , a total of 206 significant LD tests were detected out of 231 performed (p < 0.05; 89.2%), 173 after Bonferroni correction (p < 0.00022; 74.9%), which strongly supports linkage among several of these loci (Table S3). When using the reference of 50 randomly selected unlinked loci in the whole population data, 56 significant LD tests were detected at p < 0.05 out of 903 performed (6.2%) and 3 after Bonferroni correction (0.3%). A LD pairwise matrix was constructed for integrating data, and strong linkage signals were observed among 19 loci (161,992, 162,137, 163,443, 163,478, 165,368, 168,827, 168,879, 172,878, 174,144, 174,273, 179,850, 183,021, 183,584, 183,747, 184,858, 186,592, 186,978, 199,647, and 203,903), putatively constituting a single QTL related to flat oyster resistance to bonamiosis (Table 3).
Even the three remaining loci (168,138, 174,715, and Oed_rep_ c2446_3082e) showed weak signals of linkage with the main outlined cluster, suggesting location in the vicinity. This observation was also confirmed at population level, and although the figures were notably lower than those with the whole dataset (average LD at p < 0.05: ~ 40%; average after Bonferroni correction: ~12%) (Table 4) due to the lower statistical power, the proportion of LD was still much higher than the unlinked reference and involved the most consistent set of loci detected with the whole data as a single population (Table S4). We also checked for LD among the 87 outlier POOL-L loci being the figures lower than with the 22 outlier POOL-S SNP panel, although again much higher than those expected by chance. A total of 1,861 significant tests were detected over 3,741 (p < 0.05; 49.7%) and 633 after Bonferroni correction (16.9%) ( Table   S5). The most consistent SNPs detected using the whole dataset as a single population were also detected at population level as in the previous case with the 22 outlier POOL-S SNP panel (Table S6) However, when this information was integrated with the most consistent 19 outlier POOL-S loci outlined before, they showed strong linkage signals (Table S7)

| Testing outlier loci in different grouping scenarios: bonamiosis resistance vs. other environmental factors
We compared pairwise F ST values between the six oyster populations using the panels of 22 and 87 outlier POOL loci (Table 5) with regard to the whole 14,065 SNP dataset, basically representing the neutral background (Table 6) LTA and NV groups. Bayesian clustering with the 22 loci and 87 outlier POOL-L loci panels rendered in both cases K = 2 (Ln P TA B L E 3 Consistent pairwise linkage disequilibria involving the most strict 22 outlier POOL-S loci related to bonamiosis resistance using the whole population data. Significant Bonferroni p < 0.00013, 0.00013 < p < 0.001, 0.001 < p < 0.01, and 0.01 < p < 0.05 values are shown in red, yellow, green, and blue color, respectively   (Table S9).

| A practical tool for allocating individuals to their populations
The sampling design in our study was thought to identify SNP outliers related to bonamiosis resistance following a population genomics approach. Therefore, despite signals for selection related TA B L E 6 Pairwise F ST values between populations using the whole SNP dataset (below diagonal) and the most consistent 324 outlier IND loci (above diagonal) F I G U R E 3 Population structure of flat oyster beds analyzed with STRUCTURE (left part) and Discriminant Analysis of Principal Components (DAPC) (right part) including the 22 most consistent outlier POOL-S loci regarding bonamiosis resistance (upper part) and the 87 outlier POOL-L loci following a more relaxed criterion including all BY POOL (24), the Oed_rep_c2446_3082e CAND marker, and those shared by ARL POOL AND ARL HIER at p < 0.01 (62) (bottom part). For STRUCTURE analysis, each vertical bar represents one individual and its color proportion the posterior probability for assignment of each individual to the different clusters (K = 2) inferred by the program. Population codes are shown in Table 1

| D ISCUSS I ON
European flat oyster production is severely threatened by B. ostreae in the European Atlantic area, and despite the connectivity among oyster populations mainly through larvae advection (Vera et al., 2016), the presence of this parasite is apparently restricted to specific Atlantic areas (Engelsma et al., 2014). We took advantage of the restricted parasite distribution and the low genetic structure of flat oyster in the Atlantic area and selected three LTA and three NV beds from the three different genetic units previously reported (Vera et al., 2016) trying to identify genomic regions showing signatures of selection related to bonamiosis resistance. Genomic signatures for divergent or balancing selection have been reported despite high gene flow between populations Vandamme et al., 2014;Vilas et al., 2015). We assumed as a working hypothesis that outlier loci for divergent selection between both oyster groups could underlie QTLs for resistance to the parasite and that selective sweeps around these regions would determine significant differentiation over the neutral background for a set of markers considering the marker density of the screening performed.
Similar approaches have been carried out by comparing domestic vs.  Table 1. For DAPC analysis, weight of retained discriminant analysis (DA) eigenvalues representing >90% of the variance is shown wild strains or wild populations to identify candidate genomic regions related to selection in crops (Fawcett et al., 2013;Price et al., 2018) or animals (Bradic, Teotonio, & Borowsky, 2013;Hohenlohe et al., 2010;Prado, Vera, Hermida, Blanco, et al., 2018). that populations long-term exposed to the parasite are more resistant than naïve ones (Elston et al., 1987;da Silva et al., 2005) and that breeding programs respond to selection for bonamiosis resistance, which demonstrates the existence of underlying genetic variation for resistance to B. ostreae Naciri-Graven, Martin, Baud, Renault, & Gerard, 1998). Furthermore, some QTLs were detected in preliminary studies at family level (Harrang et al., 2015;Lallias et al., 2009). We applied two complementary strategies to look for genetic markers associated with resistance in a population scenario. On one hand, we focused on makers linked to candidate immune genes (CAND) showing differential expression between LTA F I G U R E 5 Potential of the 324 outliers IND loci to discriminate flat oyster populations. The outliers IND SNPs were ranked from higher to lower F ST and, using a progressive higher number of SNPs of decreasing resolution, checked the proportion of individuals correctly allocated to their population using the admixture coefficient q-value (a) and the mean q-value for each of the populations studied (b) and NV oysters previously reported by Ronza et al. (2018), for which 37 out of 77 markers were validated using Sequenom technology.
This is a rather low validation success regarding previous reports using a similar approach with anonymous markers (Cruz et al., 2017;Siccha-Ramírez et al., 2018) and could be related to the use of mature mRNA sequences as a reference for primer design or the highly polymorphic nature of bivalve species causing interference in probe sequence binding (Pino-Querido et al., 2015). The lack of the whole assembly of the flat oyster genome, which represents an important limitation as discussed below, precluded the location of introns, and thus the adequate positioning of primers for some loci. Among the 37 validated loci, we detected two CAND SNPs evidencing signals for divergent selection between LTA and NV oyster groups, one as-  (Gutierrez et al., 2017) to perform a high-resolution genomic screening with primarily anonymous markers. We identified a set of 21 highly confident SNP markers differentiating flat oyster populations with all strategies applied and 65 additional ones including other loci at p < 0.01 mostly with ARL pooling strategies.
By combining CAND and ARRAY SNP information, we defined a set of 22 very consistent outlier loci (1 CAND and 21 ARRAY markers) related to divergent selection between LTA and NV groups taking as reference the neutral background for Atlantic oyster outlined before. Remarkably, most of these markers showed strong LD among them, apparently targeting to a specific genomic region in the flat oyster genome, suggestive of a large QTL as the primary determinant of resistance to bonamiosis. When a less strict, although consistent 87 SNP set was used by adding another 65 SNPs, another cluster of nine loci demonstrated to be linked to that outlined before, providing additional support for the existence of a single major QTL.
We looked for annotation of the outlier ARRAY loci using the flat oyster database and the C. gigas genome, and significant matches and annotation of candidate markers associated to bonamiosis resistance were poor, particularly in the most interesting genomic region identified. Clearly, the lack of the whole genome flat oyster assembly and the low quality of the C. gigas genome limited successful annotation, gene mining, and functional enrichment, which strongly recommending the urgency of a draft reference genome for O. edulis.
This genome would greatly facilitate the identification of responsible genes and mutations for achieving a more resistant strain to be used for recovering production and appropriate management of flat oyster resources in the Atlantic area. , a fish species whose dispersal mechanism is mainly related to larval advection .
In this study, we used a population genomics approach to identify putative selective sweeps harnessing QTLs related to bonamiosis resistance in flat oyster, the main threat for its production and viability.
The combination of candidate gene and genomic screening strategies enabled preliminary ascertainment of the genomic architecture of resistance to bonamiosis, where a putative major QTL appears to explain most of the divergence between LTA and NV populations. Other minor QTLs scattered across the genome, including two candidate genes associated with enriched pathways and functions activated against the parasite, would also explain the divergence between both population sets. The lack of a flat oyster reference genome precluded a high-resolution analysis genes and genetic variants associated with resistance, which remains an urgent task for the near future. In any case, the information obtained provides a set of markers to be evaluated as tools to enhance selective breeding programs and for appropriate management and recovery of flat oyster beds.

ACK N OWLED G EM ENTS
This work was funded by the OYSTERECOVER project (FP7-SME-2008-2-243583) from the European Community's Seventh Framework Programme, the European Regional Development's funds (FEDER), and Xunta de Galicia local government (GRC2014/010, R2014/046). The development and provision of the medium-density SNP array for oysters was supported by Biotechnology and Biological