ORIGINAL RESEARCH published: 27 February 2020 doi: 10.3389/fgene.2019.01406 Single-Nucleotide Polymorphisms (SNP) Mining and Their Effect on the Tridimensional Protein Structure Prediction in a Set of Immunity-Related Expressed Sequence Tags (EST) in Atlantic Edited by: Gen Hua Yue, Salmon (Salmo salar) Temasek Life Sciences Laboratory, Singapore Eva Vallejos-Vidal 1, Sebastián Reyes-Cerpa 2,3, Jaime Andrés Rivas-Pardo 2,3, Kevin Maisey 4, José M. Yáñez 5, Hector Valenzuela 4, Pablo A. Cea 6, Reviewed by: Victor Castro-Fernandez 6, Lluis Tort 1, Ana M. Sandino 4, Mónica Imarai 4* Vladimir M. Milenkovic, and Felipe E. Reyes-López 1* University Medical Center Regensburg, Germany 1 Department of Cell Biology, Physiology and Immunology, Faculty of Biosciences, Universitat Autònoma de Barcelona, Zituo Yang, Barcelona, Spain, 2 Centro de Genómica y Bioinformática, Facultad de Ciencias, Universidad Mayor, Santiago, Chile, National University of 3 Escuela de Biotecnología, Facultad de Ciencias, Universidad Mayor, Santiago, Chile, 4 Centro de Biotecnología Acuícola, Singapore, Singapore Departamento de Biología, Facultad de Química y Biología, Universidad de Santiago de Chile, Santiago, Chile, 5 Facultad de *Correspondence: Ciencias Veterinarias y Pecuarias, Universidad de Chile, Santiago, Chile, 6 Facultad de Ciencias, Universidad de Chile, Mónica Imarai Santiago, Chile monica.imarai@usach.cl Felipe E. Reyes-López Felipe.Reyes@uab.cat Single-nucleotide polymorphisms (SNPs) are single genetic code variations considered one of the most common forms of nucleotide modifications. Such SNPs can be located in Specialty section: genes associated to immune response and, therefore, they may have direct implications This article was submitted to Livestock Genomics, over the phenotype of susceptibility to infections affecting the productive sector. In this a section of the journal study, a set of immune-related genes (cc motif chemokine 19 precursor [ccl19], integrin b2 Frontiers in Genetics (itb2, also named cd18), glutathione transferase omega-1 [gsto-1], heat shock 70 KDa Received: 04 December 2018 Accepted: 24 December 2019 protein [hsp70], major histocompatibility complex class I [mhc-I]) were analyzed to identify Published: 27 February 2020 SNPs by data mining. These genes were chosen based on their previously reported Citation: expression on infectious pancreatic necrosis virus (IPNV)-infected Atlantic salmon Vallejos-Vidal E, Reyes-Cerpa S, phenotype. The available EST sequences for these genes were obtained from the Rivas-Pardo JA, Maisey K, Yáñez JM, Valenzuela H, Cea PA, Unigene database. Twenty-eight SNPs were found in the genes evaluated and Castro-Fernandez V, Tort L, identified most of them as transition base changes. The effect of the SNPs located on Sandino AM, Imarai M and Reyes-López FE (2020) the 5’-untranslated region (UTR) or 3’-UTR upon transcription factor binding sites and Single-Nucleotide Polymorphisms alternative splicing regulatory motifs was assessed and ranked with a low-medium (SNP) Mining and Their Effect on the predicted FASTSNP score risk. Synonymous SNPs were found on itb2 (c.2275G > A), Tridimensional Protein Structure Prediction in a Set of Immunity-Related gsto-1 (c.558G > A), and hsp70 (c.1950C > T) with low FASTSNP predicted score risk. Expressed Sequence Tags (EST) in The difference in the relative synonymous codon usage (RSCU) value between the variant Atlantic Salmon (Salmo salar). Front. Genet. 10:1406. codons and the wild-type codon (DRSCU) showed one negative (hsp70 c.1950C > T) and doi: 10.3389/fgene.2019.01406 two positive DRSCU values (itb2 c.2275G > A; gsto-1 c.558G > A), suggesting that these Frontiers in Genetics | www.frontiersin.org 1 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure synonymous SNPs (sSNPs) may be associated to differences in the local rate of elongation. Nonsynonymous SNPs (nsSNPs) in the gsto-1 translatable gene region were ranked, using SIFT and POLYPHEN web-tools, with the second highest (c.205A > G; c484T > C) and the highest (c.499T > C; c.769A > C) predicted score risk possible. Using homology modeling to predict the effect of these nonsynonymous SNPs, the most relevant nucleotide changes for gsto-1 were observed for the nsSNPs c.205A > G, c484T > C, and c.769A > C. Molecular dynamics was assessed to analyze if these GSTO-1 variants have significant differences in their conformational dynamics, suggesting these SNPs could have allosteric effects modulating its catalysis. Altogether, these results suggest that candidate SNPs identified may play a crucial potential role in the immune response of Atlantic salmon. Keywords: single-nucleotide polymorphism, immune response, synonymous SNP, nonsynonymous SNP, homology modeling, 3D protein structure, molecular dynamics simulation, Salmo salar INTRODUCTION (tirap) (Gowin et al., 2018), c-c chemokine receptor type 5 (ccr5) (Martin et al., 1998; Salkowitz et al., 2003; Fellay et al., Genetic variation occurs within and among populations, leading 2009; Chapman and Hill, 2012), interleukin 6 (il-6) (Zhang et al., to polymorphisms that could be associated with genetic trait or 2012), and il-22 (Zhang et al., 2011), among others. also a phenotype in the presence of an environmental stimulus In species related to aquaculture, SNPs are especially (Brookes, 1999; Rebbeck et al., 2004; Hirschhorn and Daly, important because they may be associated with different 2005). A single-nucleotide polymorphism (SNP) is a single phenotypic traits with economical implications. Therefore, this genetic code variation (i.e., polymorphic). Although multiallelic increase of information has a direct impact on the accuracy of SNPs do exist, the SNPs are usually biallelic (two alternative selection for these traits, improving the rate of genetic gain and bases occur) and require a minimum frequency (>1%) in the production efficiency (Hayes et al., 2007b). The ubiquity of SNPs population (Wang et al., 1998). The SNPs are the most common across the genomes examined to date, has allowed their use as form of variation in the genome and they are extensively used to markers for a wide range of applications including quantitative study genetic differences between individuals and populations. trait locus (QTL) mapping, pedigree analysis, association studies These SNPs may contribute to changes in the genomic sequence, and population genetics, among others. Conversely, whereas the either in the coding (exons), intergenic, or noncoding (introns) effect on gene mutations in mammals has been well documented, region (Dijk et al., 2014; Ahmad et al., 2018). such information in teleost species is still limited. However, SNPs are considered the most useful biomarkers for disease several efforts have been made to provide information regarding diagnosis or prognosis due to their common frequency, ease of the consequences of genetic alterations in immune response- analysis, low genotyping costs, and the possibility to carry out related genes and may influence susceptibility to diseases in fish association studies based on statistical and bioinformatics tools (Kongchum et al., 2011). In this context, SNP variations were (Srinivasan et al., 2016). Thus, SNPs have gained importance as found on il-1b of Cyprinus pellegrini and C. carpio that can be major drivers in disease-association studies in the recent era. In helpful in understanding differential resistance to koi herpesvirus mammals, on the past decade it has been seen an enormous (KHV) and Aeromonas hydrophila, respectively (Jia et al., 2015; progress in identifying hundreds of thousands SNPs to identify Wenne, 2018). On the other hand, three SNPs were identified in associations with complex clinical conditions and phenotypic the leukocyte cell-derived chemotaxin-2 (lect2) gene to be traits associated with hundreds of common diseases (Welter associated with resistance to the big belly disease on Latis et al., 2014; Wijmenga and Zhernakova, 2018). calcarifer (Fu et al., 2014a; Wenne, 2018). Three SNPs in the Furthermore, SNPs may also have a great influence on the mast cell protease 8 (mcp-8) gene were also significantly immune response towards pathogenic challenges and diseases associated with resistance of tilapia to Streptococcus agalactiae outcome, contributing in a range of susceptibility to infections (Fu et al., 2014b; Wenne, 2018). among the individuals. Thus, the SNP may have a protective In Atlantic salmon, some studies have reported SNPs role, may influence the rate of diseases progression or even the associated with a resistance genotype against infectious type of cellular immune response evoked by pathogens (Hill, pancreatic necrosis (IPN) virus (IPNV). IPN is a highly 2001; Skevaki et al., 2015). In this regard, polymorphisms on transmissible disease with worldwide distribution that occurs several immune-related genes have been associated with both at the initial stage of rearing in freshwater and in post- susceptibility to infections including pattern recognition smolts in seawater (Bruno, 2004). Importantly, asymptomatic receptors (prr) and downstream signaling molecules (Skevaki carriers (McAllister et al., 1993), establishment of viral et al., 2015), mannose- binding lectin 2 (mbl2) and toll- persistence (Reyes-Cerpa et al., 2014), IPN-resistance interleukin 1 receptor domain containing adaptor protein phenotype (Reyes-López et al., 2015), and vertical transmission Frontiers in Genetics | www.frontiersin.org 2 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure through eggs (Bootland et al., 1991) have been reported. In this nonsynonymous (nsSNPs) based on the change provoked in matter, several segregating SNP markers linked to a major QTL the predicted amino acid sequence. While for sSNPs a codon associated with resistance against infectious pancreatic necrosis usage analysis was conducted, for those nsSNPs a homology virus (IPNV) in Atlantic salmon from an Scotland commercial modeling analysis was carried out in order to evaluate whether breeding program have been reported (Houston et al., 2012; they could have an effect on the predicted tridimensional protein Wenne, 2018). Similarly, a QTL in Norwegian salmon has been structure. In addition, on those nsSNP in which a significant employed in marker-assisted selection in breeding companies change in the three-dimensional protein structure was observed from Norway and Scotland, which resulted in 75% reduction in by homology modeling, an analysis of such variants over the time the number of IPN-outbreaks in the salmon farming industry. was performed by molecular dynamics (MD) simulation. This This QTL has been located on the SNP-based linkage map and study provide a set of identified candidate SNPs that may help to identified as the epithelial cadherin (cdh1-1) gene with a determine potential correlations between the immunity gene functional involvement in viral attachment and entry of IPNV expression pattern in Atlantic salmon and their response against (Moen et al., 2009; Moen et al., 2015; Wenne, 2018). Based on the pathogens they are exposed under aquaculture conditions. these reports, it seems that the strategy to detect SNPs in immune-related genes could provide a set of candidate polymorphisms that could explain the correlation between the pathogen and the disease phenotype. Particularly, the relevance MATERIAL AND METHODS to identify SNPs in the coding sequence of immunity-related Selection of Genes Modulated genes could explain directly (causal) the variability upon a in Response to IPNV and Sequence specific phenotype evaluated (Carlson et al., 2003; Cheng et al., 2004; Hirschhorn and Daly, 2005). In this context, most of the Cluster Collection The immune-related genes analyzed (cc motif chemokine 19 knowledge about fish immune response is based on large-scale precursor [ccl19], integrin b2 [itb2], glutathione transferase expressed sequence tag (EST) sequencing that has helped to omega-1 [gsto-1], heat shock 70 KDa protein [hsp70], major identify immune-related genes in teleosts. Undoubtedly, the EST histocompatibility complex class I [mhc-I]) were selected based sequencing based on tissues that play a central role on immune on their previously reported expression on splenic leukocytes response contribute to detect gene sequences that are directly isolated from IPNV-infected Atlantic salmon and whose related to host defense functions. The identification of a set of expression was also differentially modulated between the IPN- splenic leukocytes immune-related genes from Atlantic salmon susceptible and IPN-resistant phenotypes (Reyes-López et al., IPNV-infected using EST analysis has been reported (Cepeda 2015). The EST sequences for the previously above-mentioned et al., 2011; Cepeda et al., 2012). Importantly, some of these genes selected candidate genes were downloaded from the Unigene were detected differentially expressed when the IPN-susceptible database (NCBI). The detail regarding all the EST sequences and IPN-resistant phenotypes were compared (Reyes-López analyzed in this study to identify SNPs are indicated on et al., 2015). Thus, the search and identification of SNPs on Supplementary Tables 1–5. these immune-related genes using a data mining strategy may contribute to provide a set of candidate polymorphisms that could help in the progress to establish a link between the possible Data Mining for SNP Identification causes of this differential expression pattern and the IPN- The identification of SNPs by data mining was carried out in phenotype variability. Hayes et al. (Hayes et al., 2007b), order to identify possible functional effects in the sequences described for first time an in silico detection of 2,507 putative related to defense and immune response in Atlantic salmon SNPs in Atlantic salmon from the alignment of 100,866 EST. challenged with IPNV. All the sequences including in the Despite the large number of SNPs identified, there is no gene analysis were first preanalyzed in order to remove any vector directly associated with immune function. sequences or repetitive elements using Cross-Match (Ewing and Therefore, the aim of this study was the identification of SNPs Green, 1998) and RepeatMasker (Tarailo-Graovac and Chen, (in silico) by data mining upon a set of immune-related genes 2009), respectively. The search for SNPs was carried out based on whose expression has been previously reported in response to the multiple sequence alignment including the total number of infection with IPNV in Atlantic salmon. For this purpose, a set of sequences collected for each evaluated gene using the immune-related genes (cc motif chemokine 19 precursor [ccl19], HaploSNPer web-based tool (Tang et al., 2008). From the integrin b2 (itb2, also named cd18), glutathione transferase alignment analysis, each nucleotide variant detected was omega-1 [gsto-1], heat shock 70 KDa protein [hsp70], major considered as putative SNP and it was defined according to the histocompatibility complex class I [mhc-I]) were selected as nucleotide position in the gene sequence. The sequences were target for the in silico SNP search and identification using as then filtered to exclude possible sequencing errors and template the EST sequences for these genes obtained from the noninformative polymorphisms (variants with a frequency Unigene database. Based on their nucleotide sequence, the SNPs lower than 1%). The nonfiltered SNPs were chosen as the most were located in the 5’/3’-UTR or in the translated region. In the probable or reliable SNPs. The nucleotide variations were case of those nucleotide variations located in the translated described according to the nomenclature suggested by Dunnen region, the SNP were classified as synonymous (sSNPs) or and Antonarakis (Den Dunnen and Antonarakis, 2001). Frontiers in Genetics | www.frontiersin.org 3 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure Effect of SNPs on Gene Function scores obtained for the SNP obtained in this study (lower score = In order to determine the SNP location onto the gene region (5’- minimum effect; higher score = maximum effect). Thus, the UTR, coding region, 3’-UTR), the nucleotide sequence was scores obtained for SIFT were ranked from I (tolerated) to IV obtained based on their Unigene annotation from Nucleotide (intolerant); meanwhile the scores obtained for POLYPHEN database (NCBI) and compared them by alignment with BioEdit were ranked from A (benign) to E (probably damaging). The sequence alignment editor (version 7.0.5.3). In addition, the combination of the SIFT and POLYPHEN analysis give a score nucleotide sequence was also used as template to get the whose value was ranked from 1 (minimum effect) to 4 predicted amino acid sequence using ORF Finder tool (NCBI) (maximum effect). The implications of the nonsynonymous to confirm by protein BLAST (NCBI) the annotation for the SNPs reported at protein tertiary structure level was also unigene sequence analyzed. evaluated in silico through homology modeling analysis. The functional impact of the SNP was assessed depending on Homology Modeling the gene region (5’-UTR, coding region, 3’-UTR) on which the In order to evaluate whether the nonsynonymous SNPs reported nucleotide variant was located. For those SNP located in the in this study have an effect on the predicted tridimensional noncoding region (5’-UTR, 3’-UTR) the predictive nucleotide protein structure, a homology modeling analysis was carried out. variant effect upon motifs associated to transcription factor For this, the nucleotide sequence was used as target to obtain the binding sites was evaluated with TFSearch webtool (http:// predicting tridimensional protein structure with the diyhpl.us/~bryan/irc/protocol-online/protocol-cache/ CPHModels-3.0 webserver (Nielsen et al., 2010). The higher TFSEARCH.html) (Heinemeyer et al., 1998). Furthermore, the bit score alignment obtained was chosen as template for GSTO-1 predictive SNP effect onto possible exon splicing enhancer (PDB ID: 1EEM; score = 278; E-value = 4e-75). Alternatively, the [ESEfinder (http://krainer01.cshl.edu/cgi-bin/tools/ESE3/ predicted protein structure for CCL19 (PDB ID: 2HCI; score = esefinder.cgi?process=home) (Cartegni et al., 2003); RESCUE- 55; E-value = 2e-08), ITB2 (PDB ID: 2KCN; score = 42; E-value = ESE (http://hollywood.mit.edu/burgelab/rescue-ese/) 2e-04), HSP70 (PDB ID: 1YUW; score = 1036; E-value = 0) and (Fairbrother et al., 2004)] and exon splicing silencer [FASS-ESS MHC class I (PDB ID: 1KTL; score = 182; E-value = 4e-46) were (http://genes.mit.edu/fas-ess/) (Wang et al., 2004)] was also also obtained. The tridimensional protein structure for each PDB evaluated. Based on these results, the SNP functional effect was ID match was obtained to then be used as template for the predicted according to FASTSNP in order to assign a FASTSNP homology modeling of the above mentioned gene sequences score (Yuan et al., 2006). A FASTSNP score between 0 and 5 was using the MODELLER webtool (Swiss-Model) (Sali and Blundell, assigned to each individually SNP evaluated (representing from 0 1993). The stereochemical quality of the modeled tridimensional to 5 the minimum to maximum functional SNP effect, protein structure was evaluated with Procheck webtool (Swiss- respectively). On the other hand, for those SNP located in the Model) (Laskowski et al., 1996). coding region the nucleotide variant was first individually evaluated with ORF Finder tool (NCBI) in order to compare All-Atom Explicit Solvent MD Analysis the predicted unigene amino acid sequence with the SNP- The protonation state of residues at physiological conditions (pH = containing unigene sequence by BioEdit sequence alignment 7) was assigned using the Propka software included in the Maestro editor (version 7.0.5.3). In the case of the nucleotide variations Suite (Olsson et al., 2011). The models were solvated in TIP3P located in untranslated region (5’-UTR; 3’-UTR) and those did truncated octahedron with an extension of 12 Å over the protein not provoke any change in the predicted amino acid sequence surface, including 4 Cl- atoms to maintain the net charge neutrality (sSNP; synonymous amino acid change), the SNP effect was within the system. Parameters for calculations were derived from evaluated with the FASTSNP decision tree in order to assign a the AMBER ff14SB forcefield (Maier et al., 2015). Energy FASTSNP score (Yuan et al., 2006). In addition, on each sSNP minimization was carried out in four different stages, in each one identified a codon usage analysis was performed according to of them, 5,000 steps of steepest descent followed by 5,000 steps of a Sharp et al. (M.Sharp and Li, 1987). The difference in the relative conjugated gradient were performed. In the first stage, the solute synonymous codon usage (RSCU) value between the variant was fixed with a positional restraint of 500 kcal/mol Å2, and only the codons and the wild-type codon (DRSCU = RSCU mutant - accommodation of the solvent was allowed. Then, the hydrogen RSCU wild-type) was calculated for all the SNPs identified based atoms of the protein were freed from the restraint to allow their on the codon usage database (http://www.kazusa.or.jp/codon) relaxation. In the third stage, the minimization was done imposing a for Salmo salar. On the other hand, in those nucleotide variants lighter restriction of 10 kcal/mol Å2 on the heavy atoms of the in which a change in the amino acid predicted sequence was protein. Lastly, a minimization without restraint was conducted. detected (nsSNP; nonsynonymous amino acid change), the After that, the systems were equilibrated under NVT conditions, nucleotide variation impact was analyzed combining the score heating up from 0 K to 298.15 K in a window of 200 ps and obtained from Sorting Intolerant From Tolerant (SIFT) (applied maintaining the final temperature for 100 ps, using the Langevin to human variant databases) (http://sift.bii.a-star.edu.sg/) (Ng thermostat with a collision frequency of 2 ps-1. Then, the systems and Henikoff, 2003) and POLYPHEN web-based resources were equilibrated for 10 ns under NPT conditions at 298.15 K and 1 (Flanagan et al., 2010). The results were ranked according to atmosphere using the weak-coupling Berendsen barostat. Three the protocol established by Bhatti et al. (Bhatti et al., 2006) with production runs of 100 ns under the NPT conditions with random some modifications in order to homogenize the significance of all velocities seeds were performed for each protein system. All the MD Frontiers in Genetics | www.frontiersin.org 4 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure calculations were performed with periodic boundary conditions nucleotide variations evaluated could be involved in the with a time step of 2 fs, a 10 Å direct space cutoff for PME and regulation of the immune-related genes evaluated. constraining hydrogen atoms with SHAKE algorithm. All the From the total number of SNPs found, the 25% of the simulations were performed using Amber18 with GPU nucleotide substitutions were located in the coding region. At acceleration (Salomon-Ferrer et al., 2013) and the trajectories functional level, three of them (10.71%) were identified as were analyzed using cpptraj (Roe and Cheatham, 2013). synonymous SNP since no variation in the predicted amino acid sequence for itb2 (c.2275G > A), gsto-1 (c.558G > A), and hsp70 (c.1950C > T) was determined. Thus, a low-risk ranking was assigned (FASTSNP score = 1) (Table 2). In terms of codon RESULTS usage, a positive RSCU value was observed for itb2 (c.2275G > A) SNP Identification and gsto-1 (c.558G > A) (Table 3). By contrast, hsp70 (c.1950C > A set of immune-related genes (gsto-1, ccl19, itb2, hsp70, mhc-I) T) showed a negative RSCU value. This antecedents indicate that were chosen based on their previous reported role on the these sSNPs may have an impact in the local rate of translation immune response in Atlantic salmon (Cepeda et al., 2011; elongation. On the other hand, the remaining 14.29% from total Cepeda et al., 2012; Reyes-López et al., 2015). The nucleotide variations were identified in the coding region as identification of SNPs was carried out upon these genes based nonsynonymous SNPs. Importantly, all these nucleotide on data mining analysis. A total of 310 EST sequences obtained variations were only found into the gsto-1 sequence (Table 4). from the Unigene database (NCBI) were analyzed to identify Two of these variations (c.205A > G, c.484T > C) obtained the SNPs on ccl19, itb2, gsto-1, hsp70, and mhc-I. Twenty-eight SNPs highest risk score (SIFT+POLYPHEN score = 5) for were found, broken down into 18 transitions (7 A > G; 11 C > T) nonsynonymous SNPs and whose nucleotide variations and 10 transversions (7 A > C; 1 A > T; 1 C > G; 1 G > T) as the represented the modifications in the predicted amino acid most probable or reliable nucleotide variants (Table 1). From sequence of serine (polar, uncharged R group) by glycine them, 3 SNPs were found for ccl19, 2 SNPs for itb2, 15 SNPs for (nonpolar, flexible and the smallest R group; S26G) and by gsto-1, 4 SNPs for hsp70, and 4 SNPs for mhc-I (Table 1). proline (nonpolar, cyclic, and rigid R group; S119P), respectively. The other two nucleotide substitutions (c.499T > C, c.769A > C) were ranked in the next risk score level SNP Predicted Functionality (SIFT+POLYPHEN score = 4), representing in the predicted From the nucleotide variations identified, the 21.43% of the SNPs amino acid sequence a change of tyrosine (polar and aromatic R were found in the 5’-UTR region and 53.57% in the 3’-UTR group) by histidine (polar and sometimes positively charged region. At untranslated region level, mhc-I showed only imidazole R group; Y124H), and threonine (polar, uncharged R nucleotide modifications at 5’-UTR. In the case of ccl19, itb2, group) by proline (T214P), respectively (Table 4). Taken and hsp70 genes, the modifications were registered on the 3’- together, these results indicate that the SNPs found in the UTR (Table 2). Only in the case of gsto-1 the SNPs were found in gsto-1 could have a relevant impact at functional protein both 5’- and 3’-UTR regions. To evaluate the effect of these SNPs structural level. located in the UTR regions, the FASTSNP decision tree was used in order to assign a FASTSNP score risk for each of the Homology Modeling and All-Atom Explicit nucleotide variations found. The results showed that some Solvent MD Analysis SNPs on 5’-UTR (gsto-1: c.48C > T, c.89A > G; mhc-I: c.17G > In order to evaluate whether the predicted SIFT+POLYPHEN T) and 3’-UTR (ccl19: c.933T > C, c.996C > T; gsto-1: c.1068A > score on the nonsynonymous nucleotide variations for gsto-1 C, c.1177A > C; hsp70: c.2102C > T, c.2139A > G) obtained the have an impact at the protein structure level, a homology modeling minimum FASTSNP score (FASTSNP score = 0), indicating that was performed. To do this, the predicted tridimensional protein these nucleotide variations did not have any predicted structure was obtained, on which the effect of such nucleotide consequence. However, the majority of the SNPs detected on variations determined for gsto-1 were analyzed individually. these regions had a FASTSNP score of 1–3, which indicates that Using the CPHModels-3.0 webserver, the higher score the SNP generates a low to medium impact (Table 2). Based on alignment (score = 278; E-value = 4e-75) obtained for the S. this SNP functional score effect, the results suggest that the salar gsto-1 predicted amino acid sequence was the Homo sapiens TABLE 1 | Summary of single-nucleotide polymorphisms (SNPs) identified in a set of immune-related genes expressed in Salmo salar. UniGene ID Gen Gene acronym No. EST included A>G C>T A>C A>T C>G G>T 1668515 c-c motif chemokine 19 precursor ccl19 44 1 2 – – – – 3025151 cd18 cd18 23 1 – 1 – – – 1183801 glutathione transferase omega-1 gsto-1 138 3 7 5 0 0 0 1247008 heat shock 70 kDa hsp70 53 1 2 1 – – – 3025116 major histocompatibility complex class I mhc-I 52 1 – – 1 1 1 The UniGene ID, number of EST included in the SNP identification per each gene, and the transition (A > G, C > T) and transversion (A > C, A > T, C > G, G > T) nucleotide substitutions are indicated. Frontiers in Genetics | www.frontiersin.org 5 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure TABLE 2 | Evaluation of nucleotide variations located in untranslated region (5’-UTR; 3’-UTR) and synonymous single-nucleotide polymorphisms (SNPs) in the translated (t) region on immune-related genes expressed in Salmo salar. Bold and underlined letter represents the nucleotide variation (transcription factor binding site, and alternative splicing regulatory motif). The arrow behind the sequence on the predicted transcription factor binding site column shows the 5’– > 3’ sense and the motif composition. On the alternative splicing regulatory motif column, normal and italic styles indicate the presence of predicted exon splicing enhancer (ESE) or silencer (ESS), respectively. The FASTSNP score was assigned to each SNP (0: minimum score risk; 5: maximum score risk). NF, not found. ND, not determined. Glutathione-S-transferase omega 1 (GSTO-1; PDB ID: 1EEM) GSTO-1 (Figure 1B) showed high similarity, evidencing that and thus chosen as the best template structure. the salmon predicted modeled protein does not present serious The overall comparison between the human GSTO-1 problems with the structural restrictions dictated by the structure (Figure 1A) and the modeled Atlantic salmon template. Only some local differences were found by Frontiers in Genetics | www.frontiersin.org 6 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure TABLE 3 | Evaluation of relative synonymous codon usage (RSCU) in immune-relevant genes in Salmo salar. Nucleotide accession number Gene Gene acronym sSNP Amino acid Codon change RSCU DRSCU From To From To NM_001165324.1 cd18 cd18 c.2275G>A T ACG ACA 0.48 1.00 0.52 NM_001141472.1 glutathione transferase omega-1 gsto-1 c.558G>A A GCG GCA 0.44 0.84 0.40 BT046112.1 heat shock 70 kDa protein hsp70 c.1950C>T P CCG CCT 1.52 1.04 -0.48 Bold letter represents the nucleotide variations in codon for each synonymous single-nucleotide polymorphism (sSNP). RSCU and DRSCU values were calculated (DRSCU = RSCU mutant - RSCU wild type). TABLE 4 | Effect of nonsynonymous single-nucleotide polymorphisms (SNPs) in glutathione transferase omega-1 coding region of Salmo salar. Nucleotide accession Gene Gene SNP Amino acid SIFT Polyphen SIFT+Polyphen number acronym change score score score NM_001141472.1 Glutathione transferase gsto-1 c.205A>G S26G IV E 5 omega-1 c.484T>C S119P IV E 5 c.499T>C Y124H I D 4 c.769A>C T214P II D 4 SIFT and POLYPHEN scores are indicated from which the risk is obtained to assess whether the amino acid change produces a variation at protein structural level giving a risk of 1–4, being 1 the higher impact risk. S, serine; G, glycine; P, proline; Y, tyrosine; H, histidine; T, threonine. comparing in detail the human and salmon protein model of On the other hand, no appreciable variations in the predicted GSTO-1 at the helixes a4a, a4b, a7, and a8 (Figures 1B, C). To tridimensional S. salar GSTO-1 structure were observed for the further validate the model, a stereochemical evaluation was SNP c499T > C (Y124H substitution; sGSTO-1 Y214H) performed by generating a Ramachandran plot to assess the (Figure 4). By contrast, in the case of the SNP c.769C > C Phy and Psi dihedral angles distribution. The stereochemical (responsible of the modification T214P in the predicted amino quality of the modeled tridimensional protein structure showed acid sequence; sGSTO-1 T214P), a change in the tridimensional that the amino acids of the predicted sGSTO-1 structure were conformation affecting the C-terminus region is observed on the found mainly within the most favored (89.4%) and the a7 helix of the model (Figure 5). Punctually, the amino acid additional allowed (10.1%) energy regions, meanwhile only substitution makes the a7 helix shorter, thus lengthening the a7 the 0.5% were at the disallowed regions (Figure 1D). This helix kink (Figures 5B, C). This region is the special relevance in antecedent indicates the good quality of the predicted sGSTO-1 the protein structure because the residues of the C-terminal end structure obtained. are involved in the formation of the H site, a binding site that The effect of the gsto-1 nonsynonymous nucleotide variations accommodates a hydrophobic motif adjacent to the glutathione detected were individually evaluated by homology modeling of S. binding site (Board et al., 2000). The stereochemical quality of salar GSTO-1. The model of the variant c.205A > G (S26G on the the sGSTO-1 T214P was similar to sGSTO-1 (89.9% of the predicted amino acid sequence; sGSTO-1 S26G) suggests that amino acid in the most favored regions) (Figure 5D). Thus, glycine would increase the conformational freedom on the b2 the predicted salmon GSTO-1 sequence analysis suggests the sheet making it shorter (Figures 2A–C). The stereochemical high-impact of T214P on the protein functionality. quality of the sGSTO-1 S26G protein structure kept the good The predicted protein structure for S. salar CCL19 (PDB ID: quality of the model, with only a slight increase in the number of 2HCI; score = 55; E-value = 2e-08), ITB2 (PDB ID: 2KCN; score = 42; the amino acids in the most favored energy regions (90.3%) E-value = 2e-04), HSP70 (PDB ID: 1YUW; score = 1036; E- compared to the predicted sGSTO-1 structure (Figure 2D). value = 0), and MHC class I (PDB ID: 1KTL; score = 182; E- The SNP c.484T> C (responsible of the modification S119P in value = 4e-46) are also shown (Supplementary Figures 1–4). the predicted amino acid sequence; sGSTO-1 S119P) is located in Since homology modeling provides only restricted the a4 helix kink region (Figure 3). The model with the S119P information about side-chain orientations, three of the GSTO- substitution suggests that proline would makes the helix kink 1 variants that evidenced to be more affected in their structure longer (Figures 3B, C). This modification may generate as (S26G; S119P; T214P) were further studied by all-atom MD consequence that the a4 helix are no longer composed by the simulations. Based on that the homology models showed a good a4a and a4b helix but rather independent two different new a- stereochemical quality, they were further used as initial helices. Compared to the sGSTO-1 predicted structure, the coordinates for MD simulations. In Figure 6 is showed the sGSTO-1 S119P stereochemical quality showed a slight global structure of the GSTO-1 protein and highlights the decrease in the number of the amino acids in the most favored location of each mutation evaluated. Analysis of three replicas energy regions (87.9%) and an increase on those located in the of 100 ns of simulation for each protein model (GSTO-1 wild- additional allowed regions (11.6%) compared to the predicted type, S26G; S119P and T214P) showed that the progression of sGSTO-1 structure (Figure 3D). the root mean square deviation (RMSD) of the Ca atoms as a Frontiers in Genetics | www.frontiersin.org 7 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure FIGURE 1 | Homology modeling for predicted S. salar GSTO-1 based on H. sapiens GSTO-1 structure. (A) Tridimensional H. sapiens GSTO-1 (hGSTO-1) structure. (B) Predicted tridimensional S. salar GSTO-1 (sGSTO-1) structure. The differences on the secondary structure in the sGSTO-1 compared to hGSTO-1 are shown (arrow lines). (C) hGSTO-1 and sGSTO-1 overlay. (D) Ramachandran plot for the predicted sGSTO-1 structure. The amino acid distribution into most favored regions (A,B,L), additional allowed regions (a,b,l,p), generously allowed regions (~a,~b,~l,~p), and disallowed regions (GLU83) is indicated. In the protein structures, a-helices (a1- a8), b-sheets (b1-4) are indicated. The reduced glutathione (GSH) molecule is represented. function of time kept values below 3 Å, suggesting they were protein structure and, in consequence, a potential impairment in stable under the time windows explored (Supplementary the protein functionality. Nevertheless, further studies are Figure 5). needed to unravel the precise biophysical consequences of To assess if the GSTO-1 variants have significant differences these GSTO-1 variants over the function of the enzyme and in their conformational dynamics, we calculated the root mean their consequences at physiological level. square fluctuation (RMSF) over the whole trajectories (300 ns) comparing them with the wild-type trajectories. Importantly, even though the GSTO-1 variants are located far apart within the DISCUSSION protein structure, all of them seem to produce similar effects on conformational flexibility (Figure 7). Compared to the wild-type In this study, we identified a set of SNP located for a group of GSTO-1 (Figure 7A), in the three sGSTO-1 variants explored immune-related genes (ccl19, itb2, gsto-1, hsp70, mhc-I) with [S26G (Figure 7B), S119P (Figure 7C), and T214P (Figure 7D)] differential gene expression in Salmo salar challenged with IPNV the a10 helix located at the C-terminus of the protein is (Cepeda et al., 2011; Cepeda et al., 2012; Reyes-López et al., notoriously more rigid. Moreover, the loop that connects the 2015). These nucleotide variations were localized in both helixes a6 and a7 is significantly rigidized in the sGSTO-1 S26G untranslated regions (5’-UTR; 3’-UTR). In the translated and S119P compared to the wild-type, whereas for the variant region, synonymous (itb2, gsto-1, hsp70) and nonsynonymous sGSTO-1 T214P this loop is only slightly more rigid than the (gsto-1) SNPs were found. The potential impact of the wild-type GSTO-1. The detailed fluctuational profile for each nonsynonymous SNPs was also evaluated based on the GSTO-1 variant compared to the WT is shown on predicted tridimensional S. salar GSTO-1 structure obtained Supplementary Figure 6. Altogether, the differences observed through homology modeling. To our knowledge, currently in the local flexibility of secondary structure elements suggests there are no previous studies in Atlantic salmon in which the that the SNP variations detected have an impact on the GSTO-1 SNP functional effect had been evaluated at predicted Frontiers in Genetics | www.frontiersin.org 8 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure FIGURE 2 | Evaluation of single-nucleotide polymorphism (SNP) c.205A > G effect on the predicted S. salar GSTO-1 (sGSTO-1) structure by homology modeling. (A) Predicted tridimensional sGSTO-1 including the S26G substitution (sGSTO-1 S26G). (B) sGSTO-1 and sGSTO-1 S26P overlay. The punctual amino acid substitution (dotted arrow line) and the region of the secondary structure affected (b2 sheet, circle) are indicated. (C) Enlarged view of the sGSTO-1 S26P substitution shown in (B). sGSTO-1 (transparent) and sGSTO-1 (colored) are shown. The reduced glutathione (GSH) molecule is represented. (D) Ramachandran plot for the predicted sGSTO-1 S26G structure. The amino acid distribution into most favored regions (A,B,L), additional allowed regions (a,b,l,p), generously allowed regions (~a,~b,~l,~p), and disallowed regions (GLU83, ASP97) is indicated. tridimensional protein structure level. These results open the predicted amino acid sequence) appears as an interesting possibility that these candidate SNPs found on genes with alternative to select more accurately those candidate SNPs immune function could be associated to the inter-individual from a set of SNPs identified with a predicted impact. gene expression variability in Atlantic salmon and their Different strategies have been used to identify SNP including response against pathogens they are exposed at aquaculture small (single-strand conformation polymorphism (SSCP), rearing conditions. heteroduplex analyses, random shotgun, polymerase chain The interest in identifying SNPs is given by the large amount reaction (PCR) product sequencing, EST analysis) (Liu and of information that exists in model organisms that relate these Cordes, 2004), and high scale strategies (next generation nucleotide variations with susceptibility to diseases (Mooney se quenc ing (NGS) a nd whole genome s eque ncing) et al., 2010), which may allow to identify a series of candidate (Abdelrahman et al., 2017; Kumar and Kocour, 2017). genes as targets for therapies (Voisey and Morris, 2008). In Particularly, the EST sequence analysis-based SNP salmon, the association between polymorphism and identification is not new although it has not been extensively susceptibility has been suggested (Correa et al., 2017; Holborn explored. In this ambit, the EST sequences available on data et al., 2018). For this reason, the identification of polymorphisms repositories allow to exploit the redundancy existent in these including SNPs in genes of immunological relevance may help to collections making feasible the analysis on transcripts of the establish associations between polymorphism and susceptibility same gene from multiple individuals to identify potential true to diseases affecting the productive sector. In this way, the SNPs minimizing the sequencing errors (Hayes et al., 2007a). identification of SNPs and its analysis through the use of Such volume of information becomes an excellent database for different bioinformatics resources (including those intended for the search of candidate polymorphisms through the use of the evaluation of risk scores and estimating its consequences in bioinformatics approaches (Irizarry et al., 2000; Guryev et al., the predicted three-dimensional protein structure from the 2004). In Atlantic salmon previous reports have been carried out Frontiers in Genetics | www.frontiersin.org 9 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure FIGURE 3 | Evaluation of single-nucleotide polymorphism (SNP) c.484T > C effect on the predicted S. salar GSTO-1 (sGSTO-1) structure by homology modeling. (A) Predicted tridimensional sGSTO-1 including the S119P substitution (sGSTO-1 S119P). (B) sGSTO-1 and sGSTO-1 S119P overlay. The punctual amino acid substitution (dotted arrow line) and the structural region (a4 helix kink, circle) are indicated. (C) Enlarged view of the sGSTO-1 S119P substitution shown in (B). sGSTO-1 (transparent) and sGSTO-1 (colored) are shown. (D) Ramachandran plot for the predicted sGSTO-1 S119P structure. The amino acid distribution into most favored regions (A,B,L), additional allowed regions (a,b,l,p), generously allowed regions (~a,~b,~l,~p), and disallowed regions (GLU83) is indicated. to identify SNPs based on ESTs data (Hayes et al., 2007b) with a synonymous and nonsynonymous modifications. Importantly, high percentage (72.4%) of successfully validated SNPs from the a previous study reported a total of 23.8% of the SNPs found with putative identified SNPs obtained (Hayes et al., 2007a). these functional variations (Hayes et al., 2007b). These Altogether, the search of SNP from the EST database in antecedents suggest that the methodology used in this study Atlantic salmon could contribute to the identification of a set may result in the identification of a set of candidate SNPs using of candidate SNP with a high percentage of success, and also in silico approaches based on information that is currently allows to extract the maximum possible benefit from the large contained from public databases. amount (and often little exploited) of available information in Previous reports have assessed the use of the bioinformatics the data repositories. The use of integrated bioinformatics tools used in our current study (Mohamed et al., 2008; Capasso resources to evaluate the potential consequences of those SNPs et al., 2009; Liu et al., 2009; Wei et al., 2010; Rajasekaran and (in terms of risk factor and predicted protein three-dimensional Sethumadhavan, 2010). Several reported have used FASTSNP as structure) could favor the selection of certain candidate SNPs for a tool to evaluate the SNP significance in genes associated to their validation (both genetic and functional) as well as could disease susceptibility in mammals (Ali Mohamoud et al., 2014; shed lights of their effects from a physiological context point Panda and Suresh, 2014; Phani et al., 2014; Liu et al., 2016; Yu of view. et al., 2017). For SNP identification, we used the FASTSNP According to Hayes et al. (Hayes et al., 2007b), our results decision tree to evaluate whether the nucleotide variation on showed a higher number of transition (A > G, C > T) than noncoding region for ccl19, itb2, gsto-1, hsp70, and mhc-I could transversion (A > C, A > T, C > G, G > T) nucleotide have an predictive effect over gene processes including substitutions. At the same time, our results showed that 25% of transcription binding sites or splicing regulation (Yuan et al., the functional nucleotide variants were associated to 2006). Most of these nucleotide substitutions were ranked with a Frontiers in Genetics | www.frontiersin.org 10 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure FIGURE 4 | Evaluation of single-nucleotide polymorphism (SNP) c.499T > C effect on the predicted S. salar GSTO-1 (sGSTO-1) structure by homology modeling. (A) Predicted tridimensional sGSTO-1 including the Y124H substitution (sGSTO-1 Y124H). (B) sGSTO-1 and sGSTO-1 Y124H overlay. The punctual amino acid substitution (dotted arrow line) is indicated. (C) Enlarged view of the sGSTO-1 Y124H substitution shown in (B). sGSTO-1 (transparent), sGSTO-1 (colored), and the reduced glutathione (GSH) molecule are shown. (D) Ramachandran plot for the predicted sGSTO-1 Y124H structure. The amino acid distribution into most favored regions (A,B,L), additional allowed regions (a,b,l,p), generously allowed regions (~a,~b,~l,~p), and disallowed regions (GLU83, ASP97) is represented. score 1-3, assigned by the predictive change upon transcription relevant role on the antigen presentation and the subsequent factor binding sites. Furthermore, changes affecting predicted immune response activation. In regard to the modulation of exonic splicing enhancer or exonic splicing silencer motif may these genes on Atlantic salmon against pathogens, in IPNV- also compromise the gene function (Fairbrother et al., 2002; infected salmon the expression of all the genes evaluated in this Cartegni et al., 2003; Fairbrother et al., 2004; Wang et al., 2004). current study have been identified and annotated from splenic Importantly, the link between alternative splicing on immune- leukocytes cDNA library obtained from IPNV-infected Atlantic related genes in Atlantic salmon has been previously proposed salmon (Cepeda et al., 2011; Cepeda et al., 2012). In addition, the using IPNV (Maisey et al., 2011), suggesting that SNPs affecting expression for ccl19, itb2, hsp70, and mhc-I was modulated in these regulatory processes could play a relevant role on the host head kidney Atlantic salmon full-sibling families (Reyes-López immune response. Alterations in the immune response et al., 2015). In the case of ccl19, its expression is upregulated in associated to SNP located on untranslated regions have been trout head kidney leucocytes infected with IPNV (Montero et al., reported in mammals including ccl19 (Muinos-Gimeno et al., 2009), suggesting a role in inflammation and immune response 2009; Elton et al., 2010; Cai et al., 2014), itb2 (Bachmann et al., mediated by the recruitment of lymphocytes (Reyes-López et al., 1997; Mueller et al., 2004; Zhao et al., 2013), and mhc-I (Cree 2015). In the same study, we detected a differential expression et al., 2010). By contrast, in teleost there is a clear lack of pattern between IPN-susceptible and IPN-resistant families in information about polymorphisms located on the coding the expression of ccl19 and itb2 (cd18), noting in the susceptible region of these genes, being only reported previously on hsp70 families an abrupt high expression at 1 day post-infection (dpi) (Miichthys miiuy) and particularly in mhc-I of Atlantic salmon to then drop drastically at 5 dpi, meanwhile the expression in the (Grimholt et al., 2002; Grimholt et al., 2003; Miller et al., 2004; IPN-resistant families remained constant (Reyes-López et al., Glover et al., 2007; Wynne et al., 2007), probably because its 2015). In the case of hsp70 and mhc-I a similar upregulation was Frontiers in Genetics | www.frontiersin.org 11 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure FIGURE 5 | Evaluation of single-nucleotide polymorphism (SNP) c.769A > C effect on the predicted S. salar GSTO-1 (sGSTO-1) structure by homology modeling. (A) Predicted tridimensional sGSTO-1 including the T214P substitution (sGSTO-1 T214P). (B) sGSTO-1 and sGSTO-1 T214P overlay. The punctual amino acid substitution (dotted arrow line) and the region of the secondary structure affected (a7 helix, circle) are indicated. (C) Enlarged view of the sGSTO-1 T214P substitution shown in (B). (D) Ramachandran plot for the predicted sGSTO-1 structure. The amino acid distribution into most favored regions (A,B,L), additional allowed regions (a,b,l,p), generously allowed regions (~a,~b,~l,~p), and disallowed regions (GLU83) is indicated. In the protein structures, a-helices (a1- a8), b- sheets (b1-4) are indicated. The reduced glutathione (GSH) molecule is represented. observed in the susceptible and resistant families to IPNV, but may compromise the translation elongation by introducing a differences were found at 5 dpi when expression dropped to basal rarer codon (associated to a negative DRSCU value). This leads to levels in susceptible families but not in those resistant, suggesting the use of a less abundant tRNA that may decrease the rate of that an impaired antigen presentation can contribute to the IPN- local elongation and, in consequence, a lower protein synthesis susceptible phenotype of salmon (Reyes-López et al., 2015). levels and/or protein missfolding. By contrast, the opposite Further studies are needed to elucidate if there is a correlation above-mentioned effect may take place in those nucleotide between the candidate SNPs found from the SNP mining variants with a positive DRSCU value (Sauna and Kimchi- analysis on this current study and IPN-resistance. Sarfaty, 2011). In our study, one negative DRSCU value was Synonymous SNPs have been often called “silent” SNPs assign (hsp70 c.1950C > T), and two positive DRSCU values (itb2 because they do not induce a change in the protein amino acid c.2275G > A; gsto-1 c.558G > A). These values suggest that these composition due to the degenerating genetic code (more than sSNPs may be associated to differences in the local rate of one codon translate the same amino acid). In the last years, elongation. Nonetheless, it is not clear if these nucleotide sSNPs are attracting more attention in mammals since these variants may cause a change at physiological level (Sauna and silent mutations appear to be linked to a large list of diseases Kimchi-Sarfaty, 2011). Further studies are needed to elucidate (Sauna and Kimchi-Sarfaty, 2011). However, in fish no attention the physiological consequences of these polymorphisms. has been paid to these mutations so far. In human, these sSNPs We obtained the score for the nonsynonymous SNPs based can lead to disease by different mechanisms: disrupting splicing on the cutoffs for SIFT and POLYPHEN because they have signals (Parmley and Hurst, 2007) altering regulatory binding previously deemed appropriate for those nucleotide variations sites (Stergachis et al., 2013); and by changing the secondary that in mammals may play a role on infective processes (Xi et al., structure of the mRNA affecting protein expression (Brest et al., 2004; Bhatti et al., 2006). The ranking scheme strategy used in 2011). It has been also described that a change in the codon usage our study has the purpose to evaluate in silico the nucleotide Frontiers in Genetics | www.frontiersin.org 12 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure sequence. GSTO-1 is an enzyme involved in biotransformation of compounds including toxic substances and oxidative stress products, transport of ligands, and regulation of signaling pathways (Burmeister et al., 2008). The GSTO-1 enzyme is responsible for catalyzing the reaction that results from glutathione conjugates generated by GSH in response to damage signals, thus influencing the relationship between reduced glutathione and oxidized glutathione (GSH/GSSG) modulating the expression of cytokines related to the Th1 immune response (Dobashi et al., 2001), thus playing a key role on the host immune response modulation. To date, very few studies have identified the effect of SNP into the gsto-1 sequence (He et al., 2013; Zmorzynski et al., 2015), including teleost species (Liao et al., 2017). The analysis of the four gsto-1 nonsynonymous SNPs found by SIFT and POLYPHEN showed two of them with the second (c.205A > G; c484T > C) and the highest (c.499T > C; c.769A > C) score possible. These results indicate that the protein would be affected on its tridimensional structure. FIGURE 6 | Structure of the GSTO-1 homology model with labeled secondary structure elements. Colored spheres show the position of each Despite the relevance and progress made on the SNP punctual variant: red (S26G), yellow (S119P), blue (T214P). identification, there is a clear lack of studies in which would be proposed predictive models that could help to visualize the consequence of those SNPs located on the translated region and ranked with a high-risk score. In fact, one of the limitations of these tools is there is no reliable proposal about the changes that could take place in the protein structure. Hence, it is not possible to estimate at the structural level if the high impact score for the nucleotide modifications evaluated may have an effect in essential sectors for the protein functionality such as the active site or other fundamental interactions involved in the stabilization of its structure. To carry out this kind of analysis is necessary to characterize the mutant in vitro and determine the protein structure by, for instance, X-ray crystallography. Regrettably, the knowledge about salmon tridimensional protein structures is still very limited, with no protein structural information for GSTO-1. Therefore, in order to evaluate the nonsynonymous SNP effect at protein structure level it is strictly necessary the generation of the predicted tridimensional protein structure using bioinformatics tools. Thus, the homology modeling strategy could be used to predict the tridimensional structure of the protein and, with it, to be able to establish a relationship between the nonsynonymous SNPs- provoked structural protein changes and the relevance of them (SIFT+POLYPHEN score). Although there is no previous antecedents in fish, the SNP effect on the protein structure has been evaluated on genes of immunologic interest (Jang et al., 2017; Majumdar et al., 2017; FIGURE 7 | Differences in conformational dynamics between GSTO-1 wild- Melzer and Palanisamy, 2018), including human glutathione-S- type (WT), and the variants S26G, S119P, and T214P. The radius of the secondary protein structure and the gradient color (from 0.5 to 4.5 Å) are transferase superfamily (Kitteringham et al., 2007). The structure of functions of the RMSF values calculated. As the most affected protein GSTO-1 is composed mainly of an N-terminal domain of the regions, arrowheads highlight the helix a10 (black) and the loop between thioredoxin type which consists of four central b-sheets helixes a6 and a7 (white), respectively. surrounded at each end by a-helices, and a C-terminal end which is mainly made of a-helices (Board et al., 2000). The latter helices, variations on sequences of interest in order to provide fold over the N-terminal domain generating a network of hydrogen semiquantitative information and primarily intended for its use bonds which define a continuous surface (b/a structure) (Board in the absence of biochemical characterization. The et al., 2000). This b/a structure could be destabilized by the SNP nonsynonymous SNPs were found only in the gsto-1 gene c.205A > G (S26 amino acid substitution) located in b2 sheet, which Frontiers in Genetics | www.frontiersin.org 13 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure could introduce a structural conformational freedom either in the likely to be determining in the stability or the catalytic b1 sheet as in the adjacent one (b2 sheet). In addition, on the N- capabilities of the enzyme. Since conformational dynamics are terminal end there is the glutathione binding site, also called G-site, an essential aspect of a protein function, these differences are in which the cysteine 32 (located at the N-terminal end of the a1 probably impacting the functioning of sGSTO-1, and helix) forms a disulfide bond in the presence of reduced glutathione consequently the fitness of the organism (Bhabha et al., 2013). (GSH) (Rossjohn et al., 1998). It has been described in other GST Enhanced rigidity is a trait commonly associated with more proteins that a H-bond between the GSH sulfur and the OH of the stable proteins (Radestock and Gohlke, 2011). Therefore, it is tyrosine or serine stabilizes the GSH thiolate anion (Kortemme and possible that these SNPs could alter the stability of GSTO-1. On Creighton, 1995). Mutations of these serine or tyrosine sites lead to a the other hand, although these variations are not located at the substantial or complete protein inactivation (Stenberg et al., 1991; binding site of the enzyme, it is possible that these SNPs could Board et al., 1995). Based on the predicted tridimensional have allosteric effects modulating the catalysis. This because S. salar GSTO-1, the S26 is located in the proximity of the G-site. conformational dynamics can play a role in multiple aspects of Taken together, the S26G substitution could affect the enzymatic activity like accessibility of the substrate to the active protein functionality, although more studies are needed to site, product release rate and the probability of finding establish the consequence at protein level of this candidate catalytically active populations within the conformational nonsynonymous SNP. ensemble (Henzler-Wildman et al., 2007; Kokkinidis et al., 2012). The predicted protein structure for sGSTO-1 S119P showed a Altogether, and based on the modeled sGSTO-1 structure, most modification on the a4 helix kink region, rendering it longer of the candidate nonsynonymous SNP identified showed changed compared with the wild-type. This probably because the proline in the secondary structure dynamics. Therefore, the strategy to substitution prevents the H-bonds network within the alpha evaluate nonsynonymous SNP based on homology modeling and helix. In a previous report has been described that proline all-atom MD simulations provide additional evidence that may help modifies the secondary protein structure, suggesting that the to rank them in order to subsequently validate those that are most replacement of Ser119Pro probably interferes with the a4 helix relevant according to their structural effect. Further studies should formation (Chiu et al., 2013). On the other hand, despite the high be focused in the in vitro characterization of these enzymes to SIFT+POLYPHEN score obtained for the SNP c.499T > C, no generate a complete biophysical picture of the effects of these SNPs appreciable changes in terms of the secondary structure were in these genes associated to immune function. observed in the predicted sGSTO-1 Y124H model. However, the physicochemical change may have impacts that are not easily predicted from a structural point of view. The homology model for the SNP c.769A > C (sGSTO-1 AUTHOR CONTRIBUTIONS T214P) showed a change in the tridimensional conformation at The conceptualization of the study was performed by EV-V, MI, the C-terminal end of the a7 helix of the protein. The a7 helix is and FER-L. The methodology was originally proposed by EV-V, directly involved in the H site of the protein, which has a KM, AMS, MI, and FER-L. The SNP search and identification hydrophobic motif and is adjacent to the G binding site. This was carried out by EV-V, SR-C, JY, and FER-L. The codon usage site is composed of both the N-terminal and the C-terminal ends analysis was carried out by EV-V, HV, and FER-L. The and has the putative function of being a binding site for GSH or homology modeling was conducted by EV-V, JAR-P, KM, and other target molecules. In the sGSTO-1, the threonine 214 FER-L. The molecular dynamics simulations were performed by (hydrophilic nature) forms part of a hydrophobic motif, JAR-P, PC, VC-F, and FER-L. EV-V, SR-C, JAR-P, KM, JY, HV, contributing to the formation of H bonds, thus allowing the PC, VC-F, LT, AMS, MI, and FER-L participated actively in the union of a greater variety of substrates (Board et al., 2000). data analysis and interpretation. EV-V, SR-C, and FER-L wrote Therefore, the T214 substitution could decrease the efficiency of the original draft. All the authors corrected, read, and approved the protein, since having an exclusively apolar character the final manuscript. decreases the range of substrates that can bind to this active site. MD simulations of those SNPs with a relevant impact based on the predicted GSTO-1 structure (S26G; S119P; T214P) were FUNDING also assessed. The results show that these variations affect the dynamics of common regions within the GSTO-1 structure, This study was supported by INNOVA-CORFO (No. 09MCSS-6691 turning the overall structure of the enzyme more rigid than the and 09MCSS-6698), FONDECYT (No. 1161015; 11150807; wild-type variant. From the three variants above-mentioned 11180705; 11181133), DICYT- USACH, VRIDEI-USACH assayed, only T214P is located directly in the main affected (USA1899 VRIDEI 021943IB-PAP), and Universidad Mayor helix, while the other two (S26G; S119P) are located closer to the startup funds (No. OI101205; SR-C). The authors also thank to the N-terminal region. Since the protein structure is a complex grants from CONICYT-BCH (Chile) Postdoctoral fellowship (No. network of covalent and noncovalent interactions, the result 74170091; EV-V), International postdoctoral stay 2019 (Universidad obtained by MD simulations implies that the individual local de Chile, UCH1566; EV-V) and VRIDEI-USACH (FR-L). The effect of each one of these substitutions is being transduced funders had no role in study design, data collection and analysis, throughout the protein structure and affecting regions that are decision to publish, or preparation of the manuscript. Frontiers in Genetics | www.frontiersin.org 14 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure SUPPLEMENTARY MATERIAL (hHSP70) structure. (B) Predicted tridimensional S. salar HSP70 (sHSP70) structure. (C) hHSP70 and sHSP70 overlay. (D) Ramachandran plot for the The Supplementary Material for this article can be found online predicted sHSP70 structure. The amino acid distribution into most favored regions (A,B,L), additional allowed regions (a,b,l,p), generously allowed regions (~a,~b,~l, at: https://www.frontiersin.org/articles/10.3389/fgene.2019. ~p) and disallowed regions (GLN83 and ALA540) is represented. 01406/full#supplementary-material FIGURE S4 | Homology modeling for predicted S. salar MHC class I based on FIGURE S1 | Homology modeling for predicted S. salar CCL19 based on tridimensional H. sapiens MHC class I structure. (A) Tridimensional H. sapiens tridimensional H. sapiens CCL19 structure. (A) Tridimensional H. sapiens CCL19 MHC class I (hMHC class I) structure. (B) Predicted tridimensional S. salar MHC (hCCL19) structure. (B) Predicted tridimensional S. salar CCL19 (sCCL19) structure. class I (sMHC class I) structure. (C) hMHC class I and sMHC class I overlay. (C) hCCL19 and sCCL19 overlay. (D) Ramachandran plot for the predicted sCCL19 (D) Ramachandran plot for the predicted sMHC class I structure. The amino structure. The amino acid distribution into most favored regions (A,B,L), additional acid distribution into most favored regions (A,B,L), additional allowed regions (a,b,l, allowed regions (a,b,l,p), generously allowed regions (~a,~b,~l,~p) is indicated. p), generously allowed regions (~a,~b,~l,~p) and disallowed regions (ASP39 and ILE40) is represented. FIGURE S2 | Homology modeling for predicted S. salar ITB2 based on tridimensional H. sapiens ITB2 structure. (A) Tridimensional H. sapiens ITB2 (hITB2) structure. (B) Predicted tridimensional S. salar ITB2 (sITB2) structure. (C) hITB2 and FIGURE S5 | Molecular dynamics simulation for the wild-type GSTO-1 (A) and sITB2 overlay. (D) Ramachandran plot for the predicted sITB2 structure. The amino the variants S26G (B), S119P (C), and T214P (D). The root mean square deviation acid distribution into most favored regions (A,B,L), additional allowed regions (a,b,l, (RMSD) of the Ca atoms on a time windows explored of 100 ns is shown for p), generously allowed regions (~a,~b,~l,~p) is represented. each replicate. FIGURE S3 | Homology modeling for predicted S. salar HSP70 based on FIGURE S6 | Root mean square fluctuation (RMSF) profiles for each mutant tridimensional H. sapiens HSP70 structure. (A) Tridimensional H. sapiens HSP70 analyzed (S26G; S119P; T214P) compared to the GSTO-1 wild-type. REFERENCES Bruno, D. W. (2004). Changes in prevalence of clinical infectious pancreatic necrosis among farmed Scottish Atlantic salmon, Salmo salar L. between 1990 Abdelrahman, H., ElHady, M., Alcivar-Warren, A., Allen, S., Al-Tobasei, R., Bao, and 2002. Aquaculture 235, 13–26. doi: 10.1016/j.aquaculture.2003.11.035 L., et al. (2017). Aquaculture genomics, genetics and breeding in the United Burmeister, C., Luersen, K., Heinick, A., Hussein, A., Domagalski, M., Walter, R. States: current status, challenges, and priorities for future research. BMC D., et al. (2008). Oxidative stress in Caenorhabditis elegans: protective effects of Genomics 18, 191. doi: 10.1186/s12864-017-3557-1 the Omega class glutathione transferase (GSTO-1). FASEB J. 22, 343–354. doi: Ahmad, T., Valentovic, M. A., and Rankin, G. O. (2018). Effects of cytochrome 10.1096/fj.06-7426com P450 single nucleotide polymorphisms on methadone metabolism and Cai, W., Tao, J., Zhang, X., Tian, X., Liu, T., Feng, X., et al. (2014). Contribution of pharmacodynamics. Biochem. Pharmacol. 153, 196–204. doi: 10.1016/ homeostatic chemokines CCL19 and CCL21 and their receptor CCR7 to j.bcp.2018.02.020 coronary artery disease. Arter. Thromb. Vasc. Biol. 34, 1933–1941. doi: Ali Mohamoud, H. S., Manwar Hussain, M. R., El-Harouni, A. A., Shaik, N. A., 10.1161/ATVBAHA.113.303081 Qasmi, Z. U., Merican, A. F., et al. (2014). First comprehensive in silico analysis Capasso, M., Ayala, F., Russo, R., Avvisati, R. A., Asci, R., and Iolascon, A. (2009). A of the functional and structural consequences of SNPs in human GalNAc-T1 predicted functional single-nucleotide polymorphism of bone morphogenetic gene. Comput. Math. Methods Med. 2014, 904052. doi: 10.1155/2014/904052 protein-4 gene affects mRNA expression and shows a significant association Bachmann, M. F., McKall-Faienza, K., Schmits, R., Bouchard, D., Beach, J., Speiser, with cutaneous melanoma in Southern Italian population. J. Cancer Res. Clin. D. E., et al. (1997). Distinct roles for LFA-1 and CD28 during activation of Oncol. 135, 1799–1807. doi: 10.1007/s00432-009-0628-y naive T cells: adhesion versus costimulation. Immunity 7, 549–557. doi: Carlson, C. S., Eberle, M. A., Rieder, M. J., Smith, J. D., Kruglyak, L., and 10.1016/S1074-7613(00)80376-3 Nickerson, D. A. (2003). Additional SNPs and linkage-disequilibrium Bhabha, G., Ekiert, D. C., Jennewein, M., Zmasek, C. M., Tuttle, L. M., Kroon, G., analyses are necessary for whole-genome association studies in humans. Nat. et al. (2013). Divergent evolution of protein conformational dynamics in Genet. 33, 518–521. doi: 10.1038/ng1128 dihydrofolate reductase. Nat. Struct. Mol. Biol. 20, 1243–1249. doi: 10.1038/ Cartegni, L., Wang, J., Zhu, Z., Zhang, M. Q., and Krainer, A. R. (2003). ESEfinder: nsmb.2676 A web resource to identify exonic splicing enhancers. Nucleic Acids Res. 31, Bhatti, P., Church, D. M., Rutter, J. L., Struewing, J. P., and Sigurdson, A. J. (2006). 3568–3571. doi: 10.1093/nar/gkg616 Candidate single nucleotide polymorphism selection using publicly available Cepeda, V., Cofre, C., González, R., MacKenzie, S., and Vidal, R. (2011). tools: a guide for epidemiologists. Am. J. Epidemiol. 164, 794–804. doi: 10.1093/ Identification of genes involved in immune response of Atlantic salmon aje/kwj269 (Salmo salar) to IPN virus infection, using expressed sequence tag (EST) Board, P. G., Coggan, M., Wilce, M. C., and Parker, M. W. (1995). Evidence for an analysis. Aquaculture 318, 54–60. doi: 10.1016/j.aquaculture.2011.04.045 essential serine residue in the active site of the Theta class glutathione Cepeda, V., Cofre, C., González, R., MacKenzie, S., and Vidal, R. (2012). transferases. Biochem. J. 311, 247–250. doi: 10.1042/bj3110247 Corrigendum to “Identification of genes involved in immune response of Board, P. G., Coggan, M., Chelvanayagam, G., Easteal, S., Jermiin, L. S., Schulte, G. Atlantic salmon (Salmo salar) to IPN virus infection, using expressed sequence K., et al. (2000). Identification, characterization, and crystal structure of the tags (EST) analysis” [Aquaculture 318 (2011) 54–60]. Aquaculture 338–341, Omega class glutathione transferases. J. Biol. Chem. 275, 24798–24806. doi: 314. doi: 10.1016/j.aquaculture.2012.01.002 10.1074/jbc.M001706200 Chapman, S. J., and Hill, A. V. S. (2012). Human genetic susceptibility to Bootland, L. M., Dobos, P., and Stevenson, R. M. W. (1991). The IPNV carrier infectious disease. Nat. Rev. Genet. 13, 175–188. doi: 10.1038/nrg3114 state and demonstration of vertical transmission in experimentally infected Cheng, Q., Wu, B., Kager, L., Panetta, J. C., Zheng, J., Pui, C. H., et al. (2004). A brook trout. Dis. Aquat. Organ. 10, 13–21. doi: 10.3354/dao010013 substrate specific functional polymorphism of human gamma-glutamyl Brest, P., Lapaquette, P., Souidi, M., Lebrigand, K., Cesaro, A., Vouret-Craviari, V., hydrolase alters catalytic activity and methotrexate polyglutamate et al. (2011). A synonymous variant in IRGM alters a binding site for miR-196 accumulation in acute lymphoblastic leukaemia cells. Pharmacogenetics 14, and causes deregulation of IRGM-dependent xenophagy in Crohn’s disease. 557–567. doi: 10.1097/01.fpc.0000114761.78957.7e Nat. Genet. 43, 242–245. doi: 10.1038/ng.762 Chiu, C. C., Singh, S., and De Pablo, J. J. (2013). Effect of proline mutations on the Brookes, A. J. (1999). The essence of SNPs. Gene 234, 177–186. doi: 10.1016/ monomer conformations of amylin. Biophys. J. 105, 1227–1235. doi: 10.1016/ S0378-1119(99)00219-X j.bpj.2013.07.029 Frontiers in Genetics | www.frontiersin.org 15 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure Correa, K., Lhorente, J. P., Bassini, L., López, M. E., Di Genova, A., Maass, A., et al. Hayes, B., Laerdahl, J. K., Lien, S., Moen, T., Berg, P., Hindar, K., et al. (2007b). An (2017). Genome wide association study for resistance to Caligus rogercresseyi extensive resource of single nucleotide polymorphism markers associated with in atlantic salmon (Salmo salar L.) using a 50k SNP genotyping array. Atlantic salmon (Salmo salar) expressed sequences. Aquaculture 265, 82–90. Aquaculture 472, 61–65. doi: 10.1016/j.aquaculture.2016.04.008 doi: 10.1016/j.aquaculture.2007.01.037 Cree, B. A. C., Rioux, J. D., McCauley, J. L., Pierre-Antoine, F. D. G., Goyette, P., He, Q.-Y., Ding, L.-J., Sun, H.-X., Yu, B., Dai, A.-Y., Zhang, N.-Y., et al. (2013). McElroy, J., et al. (2010). A major histocompatibility class I locus contributes to FSH modulates the expression of inhibin-alpha and the secretion of inhibins multiple sclerosis susceptibility independently from HLA-DRB1*15:01. PloS via orphan nuclear receptor NUR77 in ovarian granulosa cells. Mol. Reprod. One 5, e11296. doi: 10.1371/journal.pone.0011296 Dev. 80, 734–743. doi: 10.1002/mrd.22206 Den Dunnen, J. T., and Antonarakis, E. (2001). Nomenclature for the description Heinemeyer, T., Wingender, E., Reuter, I., Hermjakob, H., Kel, A. E., Kel, O. V., of human sequence variations. Hum. Genet. 109, 121–124. doi: 10.1007/ et al. (1998). Databases on transcriptional regulation: TRANSFAC, TRRD and s004390100505 COMPEL. Nucleic Acids Res. 26, 362–367. doi: 10.1093/nar/26.1.362 Dijk, E. L. van, Auger, H., Jaszczyszyn, Y., and Thermes, C. (2014). Ten years of Henzler-Wildman, K. A., Lei, M., Thai, V., Kerns, S. J., Karplus, M., and Kern, D. next-generation sequencing technology. Trends Genet. 30, 418–426. doi: (2007). A hierarchy of timescales in protein dynamics is linked to enzyme 10.1016/j.tig.2014.07.001 catalysis. Nature 450, 913–916. doi: 10.1038/nature06407 Dobashi, K., Aihara, M., Araki, T., Shimizu, Y., Utsugi, M., Iizuka, K., et al. (2001). Hill, A. V. S. (2001). The genomics and genetics of human infectious disease Regulation of LPS induced IL-12 production by IFN-gamma and IL-4 through susceptibility. Annu. Rev. Genomics Hum. Genet. 2, 373–400. doi: 10.1146/ intracellular glutathione status in human alveolar macrophages. Clin. Exp. annurev.genom.2.1.373 Immunol. 124, 290–296. doi: 10.1046/j.1365-2249.2001.01535.x Hirschhorn, J. N., and Daly, M. J. (2005). Genome-wide association studies for Elton, T. S., Sansom, S. E., and Martin, M. M. (2010). Cardiovascular disease, single common diseases and complex traits. Nat. Rev. Genet. 6, 95–108. doi: 10.1038/ nucleotide polymorphisms; and the renin angiotensin system: is there a MicroRNA nrg1521 connection? Int. J. Hypertens. 2010, 281692. doi: 10.4061/2010/281692 Holborn, M. K., Ang, K. P., Elliott, J. A. K., Powell, F., and Boulding, E. (2018). Genome Ewing, B., and Green, P. (1998). Base-Calling of automated sequencer traces using wide association analysis for bacterial kidney disease resistance in a commercial phred. ii. error probabilities. Genome Res. 8, 175–185. doi: 10.1101/gr.8.3.175 North American Atlantic salmon (Salmo salar) population using a 50 K SNP panel. Fairbrother, W. G., Yeh, R., Sharp, P. A., and Burge, C. B. (2002). Predictive Aquaculture 495, 465–471. doi: 10.1016/j.aquaculture.2018.06.014 identification of exonic splicing enhancers in human genes. Science 297 (5583), Houston, R. D., Davey, J. W., Bishop, S. C., Lowe, N. R., Mota-Velasco, J. C., 1007–1013. doi: 10.1126/science.1073774 Hamilton, A., et al. (2012). Characterisation of QTL-linked and genome-wide Fairbrother, W. G., Yeo, G. W., Yeh, R., Goldstein, P., Mawson, M., Sharp, P. A., restriction site-associated DNA (RAD) markers in farmed Atlantic salmon. et al. (2004). RESCUE-ESE identifies candidate exonic splicing enhancers in BMC Genomics 13, 244. doi: 10.1186/1471-2164-13-244 vertebrate exons. Nucleic Acids Res. 32, 187–190. doi: 10.1093/nar/gkh393 Irizarry, K., Kustanovich, V., Li, C., Brown, N., Nelson, S., Wong, W., et al. (2000). Fellay, J., Ge, D., Shianna, K. V., Colombo, S., Ledergerber, B., Cirulli, E. T., et al. Genome-wide analysis of single-nucleotide polymorphisms in human (2009). Common genetic variation and the control of HIV-1 in humans. PloS expressed sequences. Nat. Genet. 26, 233–236. doi: 10.1038/79981 Genet. 5, e1000791. doi: 10.1371/journal.pgen.1000791 Jang, J.-H., Wei, J.-D., Kim, M., Kim, J.-Y., Cho, A. E., and Kim, J.-H. (2017). Flanagan, S. E., Patch, A. M., and Ellard, S. (2010). Using SIFT and PolyPhen to Leukotriene B4 receptor 2 gene polymorphism (rs1950504, Asp196Gly) leads predict loss-of-function and gain-of-function mutations. Genet. Test Mol. to enhanced cell motility under low-dose ligand stimulation. Exp. Mol. Med. Biomarkers 14, 533–537. doi: 10.1089/gtmb.2010.0036 49, e402. doi: 10.1038/emm.2017.192 Fu, G. H., Bai, Z. Y., Xia, J. H., Liu, X. J., Liu, F., Wan, Z. Y., et al. (2014a). Jia, Z.-Y., Ge, H.-Z., Bai, Y.-Y., Li, C.-T., and Shi, L.-Y. (2015). Characterization and Characterization of the LECT2 gene and its associations with resistance to the SNP variation of the interleukin-1 b gene of bighead carp (Cyprinus pellegrini big belly disease in asian seabass. Fish Shellfish Immunol. 37, 131–138. doi: Tchang, 1933) and five strains of common carp [Cyprinus carpio Linnaeus (1758)] 10.1016/j.fsi.2014.01.019 in China. J. Appl. Ichthyol. 31, 1102–1106. doi: 10.1111/jai.12862 Fu, G. H., Wan, Z. Y., Xia, J. H., Liu, F., Liu, X. J., and Yue, G. H. (2014b). The Kitteringham, N. R., Palmer, L., Owen, A., Lian, L. Y., Jenkins, R., Dowdall, S., et al. MCP-8 gene and its possible association with resistance to Streptococcus (2007). Detection and biochemical characterisation of a novel polymorphism agalactiae in tilapia. Fish Shellfish Immunol. 40, 331–336. doi: 10.1016/ in the human GSTP1 gene. Biochim. Biophys. Acta 1770, 1240–1247. doi: j.fsi.2014.07.019 10.1016/j.bbagen.2007.05.001 Glover, K. A., Grimholt, U., Bakke, H. G., Nilsen, F., Storset, A., and Skaala, O. Kokkinidis, M., Glykos, N. M., and Fadouloglou, V. E. (2012). Protein flexibility (2007). Major histocompatibility complex (MHC) variation and susceptibility and enzymatic catalysis Adv. Protein Chem. Struct. Bio. 87, 181–218. doi: to the sea louse Lepeophtheirus salmonis in Atlantic salmon Salmo salar. Dis. 10.1016/B978-0-12-398312-1.00007-X Aquat Organ 76, 57–65. doi: 10.3354/dao076057 Kongchum, P., Sandel, E., Lutzky, S., Hallerman, E. M., Hulata, G., David, L., et al. Gowin, E., Świątek-Kościelna, B., Kałużna, E., Strauss, E., Wysocki, J., Nowak, J., (2011). Association between IL-10a single nucleotide polymorphisms and et al. (2018). How many single-nucleotide polymorphisms (SNPs) must be resistance to cyprinid herpesvirus-3 infection in common carp (Cyprinus tested in order to prove susceptibility to bacterial meningitis in children? carpio). Aquaculture 315, 417–421. doi: 10.1016/j.aquaculture.2011.02.035 analysis of 11 SNPs in five genes involved in the immune response and their Kortemme, T., and Creighton, T. E. (1995). Ionisation of cysteine residues at the effect on the susceptibility to bacterial meningitis in children. Innate Immun. termini of model alpha-helical peptides. Relevance to unusual thiol pKa values 24, 163–170. doi: 10.1177/1753425918762038 in proteins of the thioredoxin family. J. Mol. Biol. 253, 799–812. doi: 10.1006/ Grimholt, U., Drabløs, F., Jørgensen, S. M., Høyheim, B., and Stet, R. J. M. (2002). jmbi.1995.0592 The major histocompatibility class I locus in Atlantic salmon (Salmo salar L.): Kumar, G., and Kocour, M. (2017). Applications of next-generation sequencing in Polymorphism, linkage analysis and protein modelling. Immunogenetics 54, fisheries research: a review. Fish. Res. 186, 11–22. doi: 10.1016/ 570–581. doi: 10.1007/s00251-002-0499-8 j.fishres.2016.07.021 Grimholt, U., Larsen, S., Nordmo, R., Midtlyng, P., Kjoeglum, S., Storset, A., et al. Laskowski, R. A., Rullmannn, J. A., MacArthur, M. W., Kaptein, R., and Thornton, (2003). MHC polymorphism and disease resistance in Atlantic salmon (Salmo J. M. (1996). AQUA and PROCHECK-NMR: programs for checking the salar); facing pathogens with single expressed major histocompatibility class I quality of protein structures solved by NMR. J. Biomol NMR 8, 477–486. and class II loci. Immunogenetics 55, 210–219. doi: 10.1007/s00251-003-0567-8 doi: 10.1007/BF00228148 Guryev, V., Berezikov, E., Malik, R., Plasterk, R. H. A., and Cuppen, E. (2004). Liao, Z., Wan, Q., Shang, X., and Su, J. (2017). Large-scale SNP screenings identify Single nucleotide polymorphisms associated with rat expressed sequences. markers linked with GCRV resistant traits through transcriptomes of Genome Res. 14, 1438–1443. doi: 10.1101/gr.2154304 individuals and cell lines in Ctenopharyngodon idella. Sci. Rep. 7, 1184. doi: Hayes, B. J., Nilsen, K., Berg, P. R., Grindflek, E., and Lien, S. (2007a). SNP 10.1038/s41598-017-01338-7 detection exploiting multiple sources of redundancy in large EST collections Liu, Z. J., and Cordes, J. F. (2004). DNA marker technologies and their improves validation rates. Bioinformatics 23, 1692–1693. doi: 10.1093/ applications in aquaculture genetics. Aquaculture 238, 1–37. doi: 10.1016/ bioinformatics/btm154 j.aquaculture.2004.05.027 Frontiers in Genetics | www.frontiersin.org 16 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure Liu, Y. H., Li, C. G., and Zhou, S. F. (2009). Prediction of deleterious functional Olsson, M. H. M., Søndergaard, C. R., Rostkowski, M., and Jensen, J. H. (2011). effects of non-synonymous single nucleotide polymorphisms in human nuclear PROPKA3: Consistent Treatment of Internal and Surface Residues in receptor genes using a bioinformatics approach. Drug Metab. Lett. 3, 242–286. Empirical pKa predictions. J. Chem. Theory Comput. 7, 525–537. doi: doi: 10.2174/187231209790218145 10.1021/ct100578z Liu, J. M., Cui, Y. Z., Zhang, G. L., Zhou, X. Y., Pang, J. X., Wang, X. Z., et al. Panda, R., and Suresh, P. K. (2014). Computational identification and analysis of (2016). Association between dentin matrix protein 1 (Rs10019009) functional polymorphisms involved in the activation and detoxification genes polymorphism and ankylosing spondylitis in a Chinese Han population implicated in endometriosis. Gene 542, 89–97. doi: 10.1016/j.gene.2014.03.058 from Shandong province. Chin. Med. J. (Engl). 129, 657–664. doi: 10.4103/ Parmley, J. L., and Hurst, L. D. (2007). Exonic splicing regulatory elements skew 0366-6999.177972 synonymous codon usage near intron-exon boundaries in mammals. Mol. Biol. Maier, J. A., Martinez, C., Kasavajhala, K., Wickstrom, L., Hauser, K. E., and Evol. 24, 1600–1603. doi: 10.1093/molbev/msm104 Simmerling, C. (2015). ff14SB: improving the accuracy of protein side chain Phani, N. M., Acharya, S., Xavy, S., Bhaskaranand, N., Bhat, M. K., Jain, A., et al. and backbone parameters from ff99SB. J. Chem. Theory Comput. 11, 3696– (2014). Genetic association of KCNJ10 rs1130183 with seizure susceptibility 3713. doi: 10.1021/acs.jctc.5b00255 and computational analysis of deleterious non-synonymous SNPs of KCNJ10 Maisey, K., Toro-Ascuy, D., Montero, R., Reyes-López, F. E., and Imarai, M. gene. Gene 536, 247–253. doi: 10.1016/j.gene.2013.12.026 (2011). Identification of CD3ϵ, CD4, CD8b splice variants of Atlantic salmon. Radestock, S., and Gohlke, H. (2011). Protein rigidity and thermophilic Fish Shellfish Immunol. 31, 815–822. doi: 10.1016/j.fsi.2011.07.022 adaptation. Proteins Struct. Funct. Bioinforma. 79, 1089–1108. doi: 10.1002/ Majumdar, I., Nagpal, I., and Paul, J. (2017). Homology modeling and in silico prot.22946 prediction of Ulcerative colitis associated polymorphisms of NOD1. Mol. Cell. Rajasekaran, R., and Sethumadhavan, R. (2010). Exploring the structural and Probes 35, 8–19. doi: 10.1016/j.mcp.2017.05.009 functional effect of pRB by significant nsSNP in the coding region of RB1 gene Martin, M. P., Dean, M., Smith, M. W., Winkler, C., Cerrard, B., Michael, N. L., causing retinoblastoma. Sci. China Life Sci. 53, 234–240. doi: 10.1007/s11427- et al. (1998). Genetic acceleration of AlDS. Progression by a promoter variant 010-0039-y of ccr5. Science 282 (5395), 1907–1911. doi: 10.1126/science.282.5395.1907 Rebbeck, T. R., Spitz, M., and Wu, X. (2004). Assessing the function of genetic McAllister, P. E., Schill, W. B., Owens, W. J., and Hodge, D. L. (1993). Determining the variants in candidate gene association studies. Nat. Rev. Genet. 5, 589–597. doi: prevalence of infectious pancreatic necrosis virus in asymptomatic brook trout 10.1038/nrg1403 Salvelinus fontinalis : a study of clinical samples and processing methods. Dis. Reyes-Cerpa, S., Reyes-López, F., Toro-Ascuy, D., Montero, R., Maisey, K., Acuña- Aquat. Organims 15, 157–162. doi: 10.3354/dao015157 Castillo, C., et al. (2014). Induction of anti-inflammatory cytokine expression Melzer, A. M., and Palanisamy, N. (2018). Deleterious single nucleotide by IPNV in persistent infection. Fish Shellfish Immunol. 41, 172–182. doi: polymorphisms of protein kinase R identified by the computational 10.1016/j.fsi.2014.08.029 approach. Mol. Immunol. 101, 65–73. doi: 10.1016/j.molimm.2018.05.026 Reyes-López, F. E., Romeo, J. S., Vallejos-Vidal, E., Reyes-Cerpa, S., Sandino, Miller, K. M., Winton, J. R., Schulze, a D., Purcell, M. K., and Ming, T. J. (2004). A. M., Tort, L., et al. (2015). Differential immune gene expression profiles in Major histocompatibility complex loci are associated with susceptibility of susceptible and resistant full-sibling families of Atlantic salmon (Salmo salar) Atlantic salmon to infectious hematopoietic necrosis virus. Environ. Biol. challenged with infectious pancreatic necrosis virus (IPNV). Dev. Comp. Fishes 69, 307–316. doi: 10.1023/B:EBFI.0000022874.48341.0f Immunol. 53, 210–221. doi: 10.1016/j.dci.2015.06.017 Moen, T., Baranski, M., Sonesson, A. K., and Kjoglum, S. (2009). Confirmation Roe, D. R., and Cheatham, T. E. (2013). PTRAJ and CPPTRAJ: Software for and fine-mapping of a major QTL for resistance to infectious pancreatic processing and analysis of molecular dynamics trajectory data. J. Chem. Theory necrosis in Atlantic salmon (Salmo salar): population-level associations Comput. 9, 3084–3095. doi: 10.1021/ct400341p between markers and trait. BMC Genomics 10, 368. doi: 10.1023/B: Rossjohn, J., Polekhina, G., Feil, S. C., Allocati, N., Masulli, M., De Illio, C., et al. EBFI.0000022874.48341.0f (1998). A mixed disulfide bond in bacterial glutathione transferase: functional Moen, T., Torgersen, J., Santi, N., Davidson, W. S., Baranski, M., Ødegård, J., et al. and evolutionary implications. Structure 6, 721–734. doi: 10.1016/S0969-2126 (2015). Epithelial cadherin determines resistance to infectious pancreatic (98)00074-4 necrosis virus in Atlantic salmon. Genetics 200, 1313–1326. doi: 10.1534/ Sali, A., and Blundell, T. L. (1993). Comparative protein modelling by satisfaction genetics.115.175406 of spatial restraints. J. Mol. Biol. 234, 779–815. doi: 10.1006/jmbi.1993.1626 Mohamed, M. F., Frye, R. F., and Langaee, T. Y. (2008). Interpopulation variation Salkowitz, J. R., Bruse, S. E., Meyerson, H., Valdez, H., Mosier, D. E., Harding, frequency of human inosine 5’-monophosphate dehydrogenase type II C. V., et al. (2003). CCR5 promoter polymorphism determines macrophage (IMPDH2) genetic polymorphisms. Genet. Test 12, 513–516. doi: 10.1089/ CCR5 density and magnitude of HIV-1 propagation in vitro. Clin. Immunol. gte.2008.0049 108, 234–240. doi: 10.1016/S1521-6616(03)00147-5 Montero, J., Chaves-Pozo, E., Cuesta, A., and Tafalla, C. (2009). Chemokine Salomon-Ferrer, R., Case, D. A., and Walker, R. C. (2013). An overview of the transcription in rainbow trout (Oncorhynchus mykiss) is differently Amber biomolecular simulation package. Wiley Interdiscip. Rev. Comput. Mol. modulated in response to viral hemorrhagic septicaemia virus (VHSV) or Sci. 3, 198–210. doi: 10.1002/wcms.1121 infectious pancreatic necrosis virus (IPNV). Fish Shellfish Immunol. 27, 661– Sauna, Z. E., and Kimchi-Sarfaty, C. (2011). Understanding the contribution of 669. doi: 10.1016/j.fsi.2009.08.003 synonymous mutations to human disease. Nat. Rev. Genet. 12, 683–691. doi: Mooney, S. D., Krishnan, V. G., and Evani, U. S. (2010). Bioinformatic tools for 10.1038/nrg3051 identifying disease gene and SNP candidates. In: Barnes, M., Breen, G. (eds) Sharp, M. P., and Li, W.-H. (1987). The codon adaptation index - a measure of Genetics Variation. Methods Mol. Biol. 628, 307–319. Totowa, NJ: Humana directional synonymous codon usage bias, and its potential applications. Press. doi: 10.1007/978-1-60327-367-1_17 Nucleic Acids Res. 15, 1281–1295. doi: 10.1093/nar/15.3.1281 Mueller, K. L., Daniels, M. A., Felthauser, A., Kao, C., Jameson, S. C., and Shimizu, Skevaki, C., Pararas, M., Kostelidou, K., Tsakris, A., and Routsias, J. G. (2015). Y. (2004). Cutting edge: LFA-1 integrin-dependent T cell adhesion is regulated Single nucleotide polymorphisms of Toll-like receptors and susceptibility to by both ag specificity and sensitivity. J. Immunol. 173, 2222–2226. doi: 10.4049/ infectious diseases. Clin. Exp. Immunol. 180, 165–177. doi: 10.1111/cei.12578 jimmunol.173.4.2222 Srinivasan, S., Clements, J. A., and Batra, J. (2016). Single nucleotide Muinos-Gimeno, M., Guidi, M., Kagerbauer, B., Martin-Santos, R., Navines, R., polymorphisms in clinics: Fantasy or reality for cancer? Crit. Rev. Clin. Lab. Alonso, P., et al. (2009). Allele variants in functional MicroRNA target sites of Sci. 53, 29–39. doi: 10.3109/10408363.2015.1075469 the neurotrophin-3 receptor gene (NTRK3) as susceptibility factors for anxiety Stenberg, G., Board, P. G., and Mannervik, B. (1991). Mutation of an disorders. Hum. Mutat. 30, 1062–1071. doi: 10.1002/humu.21005 evolutionarily conserved tyrosine residue in the active site of a human class Ng, P. C., and Henikoff, S. (2003). SIFT: Predicting amino acid changes that affect Alpha glutathione transferase. FEBS Lett. 293, 153–155. doi: 10.1016/0014- protein function. Nucleic Acids Res. 31, 3812–3814. doi: 10.1093/nar/gkg509 5793(91)81174-7 Nielsen, M., Lundegaard, C., Lund, O., and Petersen, T. N. (2010). CPHmodels- Stergachis, A. B., Haugen, E., Shafer, A., Fu, W., Vernot, B., Reynolds, A., et al. 3.0–remote homology modeling using structure-guided sequence profiles. (2013). Exonic transcription factor binding directs codon choice and affects Nucleic Acids Res. 38, 576–581. doi: 10.1093/nar/gkq535 protein evolution. Science (80), 1367–1372. doi: 10.1126/science.1243490 Frontiers in Genetics | www.frontiersin.org 17 February 2020 | Volume 10 | Article 1406 Vallejos-Vidal et al. SNPs Mining’s Effect on 3D-Protein Structure Tang, J., Leunissen, J. A., Voorrips, R. E., van der Linden, C. G., and Vosman, B. population screenings are predicted to impact protein function. Genomics (2008). HaploSNPer: a web-based allele and SNP detection tool. BMC Genet. 9, 83, 970–979. doi: 10.1016/j.ygeno.2003.12.016 23. doi: 10.1186/1471-2156-9-23 Yu, Y., Chen, Y., Wang, F.-L., Sun, J., Li, H.-J., and Liu, J.-M. (2017). Cytokines Tarailo-Graovac, M., and Chen, N. (2009). Using RepeatMasker to identify Interleukin 4 (IL-4) and Interleukin 10 (IL-10) gene polymorphisms as repetitive elements in genomic sequences. Curr. Protoc. Bioinforma. 4.10.1– potential host susceptibility factors in virus-induced encephalitis. Med. Sci. 4.10.14. doi: 10.1002/0471250953.bi0410s25 Monit. 23, 4541–4548. doi: 10.12659/MSM.904364 Voisey, J., and Morris, C. P. (2008). SNP technologies for drug discovery: a current Yuan, H. Y., Chiou, J. J., Tseng, W. H., Liu, C. H., Liu, C. K., Lin, Y. J., et al. (2006). review. Curr. Drug Discovery Technol. 5, 230–235. doi: 10.2174/ FASTSNP: an always up-to-date and extendable service for SNP function analysis 157016308785739811 and prioritization. Nucleic Acids Res. 34, 635–641. doi: 10.1093/nar/gkl236 Wang, D. G., Fan, J., Siao, C., Berno, A., Young, P., Sapolsky, R., et al. (1998). Zhang, G., Chen, X., Chan, L., Zhang, M., Zhu, B., Wang, L., et al. (2011). An SNP Large-scale identification, mapping, and genotyping of single-nucleotide selection strategy identified IL-22 associating with susceptibility to tuberculosis polymorphisms in the human genome. Science 280 (5366), 1077–1082. doi: in Chinese. Sci. Rep. 1, 20. doi: 10.1038/srep00020 10.1126/science.280.5366.1077 Zhang, G., Zhou, B., Wang, W., Zhang, M., Zhao, Y., Wang, Z., et al. (2012). A Wang, Z., Rolish, M. E., Yeo, G., Tung, V., Mawson, M., and Burge, C. B. (2004). functional single-nucleotide polymorphism in the promoter of the gene Systematic identification and analysis of exonic splicing silencers. Cell 119, encoding interleukin 6 is associated with susceptibility to tuberculosis. 831–845. doi: 10.1016/j.cell.2004.11.010 J. Infect. Dis. 205, 1697–1704. doi: 10.1093/infdis/jis266 Wei, Q., Wang, L., Wang, Q., Kruger, W. D., and Dunbrack, R. L. (2010). Testing Zhao, R., Song, Z., Dong, R., Li, H., Shen, C., and Zheng, S. (2013). Polymorphism computational prediction of missense mutation phenotypes: functional of ITGB2 gene 3’-UTR+145C/A is associated with biliary atresia. Digestion 88, characterization of 204 mutations of human cystathionine beta synthase. 65–71. doi: 10.1159/000352025 Proteins 78, 2058–2074. doi: 10.1002/prot.22722 Zmorzynski, S., Swiderska-Kolacz, G., Koczkodaj, D., and Filip, A. A. (2015). Welter, D., MacArthur, J., Morales, J., Burdett, T., Hall, P., Junkins, H., et al. Significance of polymorphisms and expression of enzyme-encoding genes (2014). The NHGRI GWAS Catalog, a curated resource of SNP-trait related to glutathione in hematopoietic cancers and solid tumors. BioMed. associations. Nucleic Acids Res. 42, 1001–1006. doi: 10.1093/nar/gkt1229 Res. Int. 2015, 853573. doi: 10.1155/2015/853573 Wenne, R. (2018). Single nucleotide polymorphism markers with applications in aquaculture and assessment of its impact on natural populations. Aquat. Living Conflict of Interest: The authors declare that the research was conducted in the Resour. 31, 2. doi: 10.1051/alr/2017043 absence of any commercial or financial relationships that could be construed as a Wijmenga, C., and Zhernakova, A. (2018). The importance of cohort studies potential conflict of interest. in the post-GWAS era. Nat. Genet. 50, 322–328. doi: 10.1038/s41588- 018-0066-3 Copyright © 2020 Vallejos-Vidal, Reyes-Cerpa, Rivas-Pardo, Maisey, Yáñez, Valenzuela, Wynne, J. W., Cook, M. T., Nowak, B. F., and Elliott, N. G. (2007). Major Cea, Castro-Fernandez, Tort, Sandino, Imarai and Reyes-López. This is an open-access histocompatibility polymorphism associated with resistance towards amoebic article distributed under the terms of the Creative Commons Attribution License (CC BY). gill disease in Atlantic salmon (Salmo salar L.). Fish Shellfish Immunol. 22, 707– The use, distribution or reproduction in other forums is permitted, provided the original 717. doi: 10.1016/j.fsi.2006.08.019 author(s) and the copyright owner(s) are credited and that the original publication in this Xi, T., Jones, I. M., and Mohrenweiser, H. W. (2004). Many amino acid journal is cited, in accordance with accepted academic practice. No use, distribution or substitution variants identified in DNA repair genes during human reproduction is permitted which does not comply with these terms. Frontiers in Genetics | www.frontiersin.org 18 February 2020 | Volume 10 | Article 1406