RNA-Seq Data for Reliable SNP Detection and Genotype Calling: Interest for Coding Variant Characterization and Cis-Regulation Analysis by Allele-Specific Expression in Livestock Species

Transcriptome-wide association mapping provides insights into the genetic basis and candidate genes governing flowering, maturity and seed weight in rice bean (Vigna umbellata)

BACKGROUND

Rice bean (Vigna umbellata), an underrated legume, adapts to diverse climatic conditions with the potential to support food and nutritional security worldwide. It is used as a vegetable, minor food crop and a fodder crop, being a rich source of proteins, minerals, and essential fatty acids. However, little effort has been made to decipher the genetic and molecular basis of various useful traits in this crop. Therefore, we considered three economically important traits i.e., flowering, maturity and seed weight of rice bean and identified the associated candidate genes employing an associative transcriptomics approach on 100 diverse genotypes out of 1800 evaluated rice bean accessions from the Indian National Genebank.

RESULTS

The transcriptomics-based genotyping of one-hundred diverse rice bean cultivars followed by pre-processing of genotypic data resulted in 49,271 filtered markers. The STRUCTURE, PCA and Neighbor-Joining clustering of 100 genotypes revealed three putative sub-populations. The marker-trait association analysis involving various genome-wide association study (GWAS) models revealed significant association of 82 markers on 48 transcripts for flowering, 26 markers on 22 transcripts for maturity and 22 markers on 21 transcripts for seed weight. The transcript annotation provided information on the putative candidate genes for the considered traits. The candidate genes identified for flowering include HSC80, P-II PsbX, phospholipid-transporting-ATPase-9, pectin-acetylesterase-8 and E3-ubiquitin-protein-ligase-RHG1A. Further, the WRKY1 and DEAD-box-RH27 were found to be associated with seed weight. Furthermore, the associations of PIF3 and pentatricopeptide-repeat-containing-gene with maturity and seed weight, and aldo-keto-reductase with flowering and maturity were revealed.

CONCLUSION

This study offers insights into the genetic basis of key agronomic traits in rice bean, including flowering, maturity, and seed weight. The identified markers and associated candidate genes provide valuable resources for future exploration and targeted breeding, aiming to enhance the agronomic performance of rice bean cultivars. Notably, this research represents the first transcriptome-wide association study in pulse crop, uncovering the candidate genes for agronomically useful traits.

Transcriptome analysis of primary adult B-cell lineage acute lymphoblastic leukemia identifies pathogenic variants and gene fusions, and predicts subtypes for in depth molecular diagnosis

BACKGROUND

B-cell acute lymphoblastic leukemia (B-ALL) is classified into subgroups based on known driver oncogenes and molecular lesions, including translocations and recurrent mutations. However, the current diagnostic tests do not identify subtypes or oncogenic lesions for all B-ALL samples, creating a heterogeneous B-ALL group of unknown subtypes.

METHODS

We sorted primary adult B-ALL cells and performed transcriptome analysis by bulk RNA sequencing (RNA-seq).

RESULTS

Transcriptomic analysis of an adult B-ALL cohort allowed the classification of four patient samples with subtypes that were not previously revealed by standard gene panels. The leukemia of two patients were of the DUX4 subtype and two were CRLF2+ Ph-like B-ALL. Furthermore, single nucleotide variant analysis detected the oncogenic NRAS-G12D, KRAS-G12D, and KRAS-G13D mutations in three of the patient samples, presenting targetable mutations. Additional oncogenic variants and gene fusions were uncovered, as well as multiple variants in the PDE4DIP gene across five of the patient samples.

CONCLUSION

We demonstrate that RNA-seq is an effective tool for precision medicine in B-ALL by providing comprehensive molecular profiling of leukemia cells, identifying subtype and oncogenic lesions, and stratifying patients for appropriate therapy.

Genotype prediction of 336,463 samples from public expression data

Tens of thousands of RNA-sequencing experiments comprising hundreds of thousands of individual samples have now been performed. These data represent a broad range of experimental conditions, sequencing technologies, and hypotheses under study. The Recount project has aggregated and uniformly processed hundreds of thousands of publicly available RNA-seq samples. Most of these samples only include RNA expression measurements; genotype data for these same samples would enable a wide range of analyses including variant prioritization, eQTL analysis, and studies of allele specific expression. Here, we developed a statistical model based on the existing reference and alternative read counts from the RNA-seq experiments available through Recount3 to predict genotypes at autosomal biallelic loci in coding regions. We demonstrate the accuracy of our model using large-scale studies that measured both gene expression and genotype genome-wide. We show that our predictive model is highly accurate with 99.5% overall accuracy, 99.6% major allele accuracy, and 90.4% minor allele accuracy. Our model is robust to tissue and study effects, provided the coverage is high enough. We applied this model to genotype all the samples in Recount 3 and provide the largest ready-to-use expression repository containing genotype information. We illustrate that the predicted genotype from RNA-seq data is sufficient to unravel the underlying population structure of samples in Recount3 using Principal Component Analysis.

Development of Gene-Based InDel Markers on Putative Drought Stress-Responsive Genes and Genetic Diversity of Durian (Durio zibethinus)

Insertion-deletion (InDel) markers are co-dominant, relatively abundant and practical for agarose gel genotyping. InDel polymorphism usually affects gene functions. Nucleotide sequences of durian (Durio zibethinus) are available, but InDel makers have not been well established. This study aimed to develop drought-related gene-based InDel markers for durian through bioinformatic analysis of RNA-Seq datasets. The polymorphism of the markers was verified in 24 durian genotypes local to Thailand. Bioinformatic analysis indicated 496 InDel loci having lengths more than 9 bp. To evaluate these InDel markers, 15 InDel loci were selected. Nine markers were successfully amplified a clear polymorphic band pattern on 2% agarose gel. The polymorphic information content (PIC) of these nine markers ranged from 0.1103 to 0.5808. The genetic distance between the 24 genotypes ranged from 0.222 to 0.889. The phylogeny based on the nine InDel loci distinguished the 24 genotypes and divided samples into four groups. This set of gene-based InDel markers on putative drought-responsive genes will be useful for genetic studies.

Identification and high-throughput genotyping of single nucleotide polymorphism markers in a non-model conifer (Abies nordmanniana (Steven) Spach)

Single nucleotide polymorphism (SNP) markers are powerful tools for investigating population structures, linkage analysis, and genome-wide association studies, as well as for breeding and population management. The availability of SNP markers has been limited to the most commercially important timber species, primarily due to the cost of genome sequencing required for SNP discovery. In this study, a combination of reference-based and reference-free approaches were used to identify SNPs in Nordmann fir (Abies nordmanniana), a species previously lacking genomic sequence information. Using a combination of a genome assembly of the closely related Silver fir (Abies alba) species and a de novo assembly of low-copy regions of the Nordmann fir genome, we identified a high density of reliable SNPs. Reference-based approaches identified two million SNPs in common between the Silver fir genome and low-copy regions of Nordmann fir. A combination of one reference-free and two reference-based approaches identified 250 shared SNPs. A subset of 200 SNPs were used to genotype 342 individuals and thereby tested and validated in the context of identity analysis and/or clone identification. The tested SNPs successfully identified all ramets per clone and five mislabeled individuals via identity and genomic relatedness analysis. The identified SNPs will be used in ad hoc breeding of Nordmann fir in Denmark.

Affordable, accurate and unbiased RNA sequencing by manual library miniaturization: A case study in barley

We present an easy-to-reproduce manual miniaturized full-length RNA sequencing (RNAseq) library preparation workflow that does not require the upfront investment in expensive lab equipment or long setup times. With minimal adjustments to an established commercial protocol, we were able to manually miniaturize the RNAseq library preparation by a factor of up to 1:8. This led to cost savings for miniaturized library preparation of up to 86.1% compared to the gold standard. The resulting data were the basis of a rigorous quality control analysis that inspected: sequencing quality metrics, gene body coverage, raw read duplications, alignment statistics, read pair duplications, detected transcripts and sequence variants. We also included a deep dive data analysis identifying rRNA contamination and suggested ways to circumvent these. In the end, we could not find any indication of biases or inaccuracies caused by the RNAseq library miniaturization. The variance in detected transcripts was minimal and not influenced by the miniaturization level. Our results suggest that the workflow is highly reproducible and the sequence data suitable for downstream analyses such as differential gene expression analysis or variant calling.

Genetic and molecular dissection of ginseng (Panax ginseng Mey.) germplasm using high-density genic SNP markers, secondary metabolites, and gene expressions

Genetic and molecular knowledge of a species is crucial to its gene discovery and enhanced breeding. Here, we report the genetic and molecular dissection of ginseng, an important herb for healthy food and medicine. A mini-core collection consisting of 344 cultivars and landraces was developed for ginseng that represents the genetic variation of ginseng existing in its origin and diversity center. We sequenced the transcriptomes of all 344 cultivars and landraces; identified over 1.5 million genic SNPs, thereby revealing the genic diversity of ginseng; and analyzed them with 26,600 high-quality genic SNPs or a selection of them. Ginseng had a wide molecular diversity and was clustered into three subpopulations. Analysis of 16 ginsenosides, the major bioactive components for healthy food and medicine, showed that ginseng had a wide variation in the contents of all 16 ginsenosides and an extensive correlation of their contents, suggesting that they are synthesized through a single or multiple correlated pathways. Furthermore, we pair-wisely examined the relationships between the cultivars and landraces, revealing their relationships in gene expression, gene variation, and ginsenoside biosynthesis. These results provide new knowledge and new genetic and genic resources for advanced research and breeding of ginseng and related species.

BBmix: a Bayesian beta-binomial mixture model for accurate genotyping from RNA-sequencing

MOTIVATION

While many pipelines have been developed for calling genotypes using RNA-sequencing (RNA-Seq) data, they all have adapted DNA genotype callers that do not model biases specific to RNA-Seq such as allele-specific expression (ASE).

RESULTS

Here, we present Bayesian beta-binomial mixture model (BBmix), a Bayesian beta-binomial mixture model that first learns the expected distribution of read counts for each genotype, and then deploys those learned parameters to call genotypes probabilistically. We benchmarked our model on a wide variety of datasets and showed that our method generally performed better than competitors, mainly due to an increase of up to 1.4% in the accuracy of heterozygous calls, which may have a big impact in reducing false positive rate in applications sensitive to genotyping error such as ASE. Moreover, BBmix can be easily incorporated into standard pipelines for calling genotypes. We further show that parameters are generally transferable within datasets, such that a single learning run of less than 1 h is sufficient to call genotypes in a large number of samples.

AVAILABILITY AND IMPLEMENTATION

We implemented BBmix as an R package that is available for free under a GPL-2 licence at https://gitlab.com/evigorito/bbmix and https://cran.r-project.org/package=bbmix with accompanying pipeline at https://gitlab.com/evigorito/bbmix_pipeline.

Identification of SNPs Related to Salmonella Resistance in Chickens Using RNA-Seq and Integrated Bioinformatics Approach

Potential single nucleotide polymorphisms (SNPs) were detected between two chicken breeds (Kashmir favorella and broiler) using deep RNA sequencing. This was carried out to comprehend the coding area alterations, which cause variances in the immunological response to Salmonella infection. In the present study, we identified high impact SNPs from both chicken breeds in order to delineate different pathways that mediate disease resistant/susceptibility traits. Samples (liver and spleen) were collected from Salmonella resistant (K. favorella) and susceptible (broiler) chicken breeds. Salmonella resistance and susceptibility were checked by different pathological parameters post infection. To explore possible polymorphisms in genes linked with disease resistance, SNP identification analysis was performed utilizing RNA seq data from nine K. favorella and ten broiler chickens. A total of 1778 (1070 SNPs and 708 INDELs) and 1459 (859 SNPs and 600 INDELs) were found to be specific to K. favorella and broiler, respectively. Based on our results, we conclude that in broiler chickens the enriched pathways mostly included metabolic pathways like fatty acid metabolism, carbon metabolism and amino acid metabolism (Arginine and proline metabolism), while as in K. favorella genes with high impact SNPs were enriched in most of the immune-related pathways like MAPK signaling pathway, Wnt signaling pathway, NOD-like receptor signaling pathway, etc., which could be a possible resistance mechanism against salmonella infection. In K. favorella, protein-protein interaction analysis also shows some important hub nodes, which are important in providing defense against different infectious diseases. Phylogenomic analysis revealed that indigenous poultry breeds (resistant) are clearly separated from commercial breeds (susceptible). These findings will offer fresh perspectives on the genetic diversity in chicken breeds and will aid in the genomic selection of poultry birds.

A deep-learning-based RNA-seq germline variant caller

Summary

RNA sequencing (RNA-seq) can be applied to diverse tasks including quantifying gene expression, discovering quantitative trait loci and identifying gene fusion events. Although RNA-seq can detect germline variants, the complexities of variable transcript abundance, target capture and amplification introduce challenging sources of error. Here, we extend DeepVariant, a deep-learning-based variant caller, to learn and account for the unique challenges presented by RNA-seq data. Our DeepVariant RNA-seq model produces highly accurate variant calls from RNA-sequencing data, and outperforms existing approaches such as Platypus and GATK. We examine factors that influence accuracy, how our model addresses RNA editing events and how additional thresholding can be used to facilitate our models' use in a production pipeline.

Supplementary information

Supplementary data are available at Bioinformatics Advances online.

RNA-Seq-based discovery of genetic variants and allele-specific expression of two layer lines and broiler chicken

Recent advances in the selective breeding of broilers and layers have made poultry production one of the fastest-growing industries. In this study, a transcriptome variant calling approach from RNA-seq data was used to determine population diversity between broilers and layers. In total, 200 individuals were analyzed from three different chicken populations (Lohmann Brown (LB), n = 90), Lohmann Selected Leghorn (LSL, n = 89), and Broiler (BR, n = 21). The raw RNA-sequencing reads were pre-processed, quality control checked, mapped to the reference genome, and made compatible with Genome Analysis ToolKit for variant detection. Subsequently, pairwise fixation index (F ST) analysis was performed between broilers and layers. Numerous candidate genes were identified, that were associated with growth, development, metabolism, immunity, and other economically significant traits. Finally, allele-specific expression (ASE) analysis was performed in the gut mucosa of LB and LSL strains at 10, 16, 24, 30, and 60 weeks of age. At different ages, the two-layer strains showed significantly different allele-specific expressions in the gut mucosa, and changes in allelic imbalance were observed across the entire lifespan. Most ASE genes are involved in energy metabolism, including sirtuin signaling pathways, oxidative phosphorylation, and mitochondrial dysfunction. A high number of ASE genes were found during the peak of laying, which were particularly enriched in cholesterol biosynthesis. These findings indicate that genetic architecture as well as biological processes driving particular demands relate to metabolic and nutritional requirements during the laying period shape allelic heterogeneity. These processes are considerably affected by breeding and management, whereby elucidating allele-specific gene regulation is an essential step towards deciphering the genotype to phenotype map or functional diversity between the chicken populations. Additionally, we observed that several genes showing significant allelic imbalance also colocalized with the top 1% of genes identified by the FST approach, suggesting a fixation of genes in cis-regulatory elements.

Robust identification of regulatory variants (eQTLs) using a differential expression framework developed for RNA-sequencing

BACKGROUND

A gap currently exists between genetic variants and the underlying cell and tissue biology of a trait, and expression quantitative trait loci (eQTL) studies provide important information to help close that gap. However, two concerns that arise with eQTL analyses using RNA-sequencing data are normalization of data across samples and the data not following a normal distribution. Multiple pipelines have been suggested to address this. For instance, the most recent analysis of the human and farm Genotype-Tissue Expression (GTEx) project proposes using trimmed means of M-values (TMM) to normalize the data followed by an inverse normal transformation.

RESULTS

In this study, we reasoned that eQTL analysis could be carried out using the same framework used for differential gene expression (DGE), which uses a negative binomial model, a statistical test feasible for count data. Using the GTEx framework, we identified 35 significant eQTLs (P < 5 × 10^-8) following the ANOVA model and 39 significant eQTLs (P < 5 × 10^-8) following the additive model. Using a differential gene expression framework, we identified 930 and six significant eQTLs (P < 5 × 10^-8) following an analytical framework equivalent to the ANOVA and additive model, respectively. When we compared the two approaches, there was no overlap of significant eQTLs between the two frameworks. Because we defined specific contrasts, we identified trans eQTLs that more closely resembled what we expect from genetic variants showing complete dominance between alleles. Yet, these were not identified by the GTEx framework.

CONCLUSIONS

Our results show that transforming RNA-sequencing data to fit a normal distribution prior to eQTL analysis is not required when the DGE framework is employed. Our proposed approach detected biologically relevant variants that otherwise would not have been identified due to data transformation to fit a normal distribution.

RNA-Seq-Pop: Exploiting the sequence in RNA sequencing-A Snakemake workflow reveals patterns of insecticide resistance in the malaria vector Anopheles gambiae

We provide a reproducible and scalable Snakemake workflow, called RNA-Seq-Pop, which provides end-to-end analysis of RNA sequencing data sets. The workflow allows the user to perform quality control, perform differential expression analyses and call genomic variants. Additional options include the calculation of allele frequencies of variants of interest, summaries of genetic variation and population structure, and genome-wide selection scans, together with clear visualizations. RNA-Seq-Pop is applicable to any organism, and we demonstrate the utility of the workflow by investigating pyrethroid resistance in selected strains of the major malaria mosquito, Anopheles gambiae. The workflow provides additional modules specifically for An. gambiae, including estimating recent ancestry and determining the karyotype of common chromosomal inversions. The Busia laboratory colony used for selections was collected in Busia, Uganda, in November 2018. We performed a comparative analysis of three groups: a parental G24 Busia strain; its deltamethrin-selected G28 offspring; and the susceptible reference strain Kisumu. Measures of genetic diversity reveal patterns consistent with that of laboratory colonization and selection, with the parental Busia strain exhibiting the highest nucleotide diversity, followed by the selected Busia offspring, and finally, Kisumu. Differential expression and variant analyses reveal that the selected Busia colony exhibits a number of distinct mechanisms of pyrethroid resistance, including the Vgsc-995S target-site mutation, upregulation of SAP genes, P450s and a cluster of carboxylesterases. During deltamethrin selections, the 2La chromosomal inversion rose in frequency (from 33% to 86%), supporting a previous link with pyrethroid resistance. RNA-Seq-Pop is hosted at: github.com/sanjaynagi/rna-seq-pop. We anticipate that the workflow will provide a useful tool to facilitate reproducible, transcriptomic studies in An. gambiae and other taxa.

Different impacts of TP53 mutations on cell cycle-related gene expression among cancer types

Functional properties caused by TP53 mutations are involved in cancer development and progression. Although most of the mutations lose normal p53 functions, some of them, gain-of-function (GOF) mutations, exhibiting novel oncogenic functions. No reports have analyzed the impact of TP53 mutations on the gene expression profile of the p53 signaling pathway across cancer types. This study is a cross-cancer type analysis of the effects of TP53 mutations on gene expression. A hierarchical cluster analysis of the expression profile of the p53 signaling pathway classified 21 cancer types into two clusters (A1 and A2). Changes in the expression of cell cycle-related genes and MKI67 by TP53 mutations were greater in cluster A1 than in cluster A2. There was no distinct difference in the effects between GOF and non-GOF mutations on the gene expression profile of the p53 signaling pathway.

Evolution of weedy giant ragweed (Ambrosia trifida): Multiple origins and gene expression variability facilitates weediness

Agricultural weeds may originate from wild populations, but the origination patterns and genetics underlying this transition remain largely unknown. Analysis of weedy-wild paired populations from independent locations may provide evidence to identify key genetic variation contributing to this adaptive shift. We performed genetic variation and expression analyses on transcriptome data from 67 giant ragweed samples collected from different locations in Ohio, Iowa, and Minnesota and found geographically separated weedy populations likely originated independently from their adjacent wild populations, but subsequent spreading of weedy populations also occurred locally. By using eight closely related weedy-wild paired populations, we identified thousands of unique transcripts in weedy populations that reflect shared or specific functions corresponding, respectively, to both convergently evolved and population-specific weediness processes. In addition, differential expression of specific groups of genes was detected between weedy and wild giant ragweed populations using gene expression diversity and gene co-expression network analyses. Our study suggests an integrated route of weedy giant ragweed origination, consisting of independent origination combined with the subsequent spreading of certain weedy populations, and provides several lines of evidence to support the hypothesis that gene expression variability plays a key role in the evolution of weedy species.

Genomic and transcriptomic analyses reveal genetic adaptation to cold conditions in the chickens

Under the pressure of natural and artificial selection, domestic animals, including chickens, have evolved unique mechanisms of genetic adaptations such as high-altitude adaptation, hot and arid climate adaptation, and desert adaptation. Here, we investigated the genetic basis of cold tolerance in chicken by integrating whole-genome and transcriptome sequencing technologies. Genome-wide comparative analyses of 118 chickens living in different latitudes showed 46 genes and several pathways that may be involved in cold adaptation. The results of the functional enrichment analysis of differentially expressed genes proved the important role of metabolic pathways and immune-related pathways in cold tolerance in chickens. The subsequent integration of whole genome and transcriptome sequencing technology further identified six genes - dnah5 (dynein axonemal heavy chain 5), ptgs2 (prostaglandin-endoperoxide synthase 2), inhba (inhibin beta A subunit), irx2 (iroquois homeobox 2), ensgalg00000054917, and ensgalg00000046652 - requiring more detailed studies. In addition, we also discovered different allele frequency distributions of five SNPs (single nucleotide polymorphisms) within ptgs2 and nine SNPs within dnah5 in chickens in different latitudes, suggesting strong selective pressure of these two genes in chickens. We provide a novel insight into the genetic adaptation in chickens to cold environments, and provide a reference for evaluating and developing adaptive chicken breeds in cold environments.

GDF15 and ACE2 stratify COVID-19 patients according to severity while ACE2 mutations increase infection susceptibility

Coronavirus disease 19 (COVID-19) is a persistent global pandemic with a very heterogeneous disease presentation ranging from a mild disease to dismal prognosis. Early detection of sensitivity and severity of COVID-19 is essential for the development of new treatments. In the present study, we measured the levels of circulating growth differentiation factor 15 (GDF15) and angiotensin-converting enzyme 2 (ACE2) in plasma of severity-stratified COVID-19 patients and uninfected control patients and characterized the in vitro effects and cohort frequency of ACE2 SNPs. Our results show that while circulating GDF15 and ACE2 stratify COVID-19 patients according to disease severity, ACE2 missense SNPs constitute a risk factor linked to infection susceptibility.

H3K27ac chromatin acetylation and gene expression analysis reveal sex- and situs-related differences in developing chicken gonads

BACKGROUND

Birds exhibit a unique asymmetry in terms of gonad development. The female left gonad generates a functional ovary, whereas the right gonad regresses. In males, both left and right gonads would develop into testes. How is this left/right asymmetry established only in females but not in males remains unknown. The epigenetic regulation of gonadal developmental genes may contribute to this sex disparity. The modification of histone tails such as H3K27ac is tightly coupled to chromatin activation and gene expression. To explore whether H3K27ac marked chromatin activation is involved in the asymmetric development of avian gonads, we probed genome-wide H3K27ac occupancy in left and right gonads from both sexes and related chromatin activity profile to the expression of gonadal genes. Furthermore, we validated the effect of chromatin activity on asymmetric gonadal development by manipulating the chromatin histone acetylation levels.

METHODS

The undifferentiated gonads from both sides of each sex were collected and subjected to RNA-Seq and H3K27ac ChIP-Seq experiments. Integrated analysis of gene expression and active chromatin regions were performed to identify the sex- and situs-specific regulation and expression of gonadal genes. The histone deacetylase inhibitor trichostatin A (TSA) was applied to the undifferentiated female right gonads to assess the effect of chromatin activation on gonadal gene expression and cell proliferation.

RESULTS

Even before sex differentiation, the gonads already show divergent gene expression between different sexes and between left/right sides in females. The sex-specific H3K27ac chromatin distributions coincide with the higher expression of male/female specification genes in each sex. Unexpectedly, the H3K27ac marked chromatin activation show a dramatic difference between left and right gonads in both sexes, although the left/right asymmetric gonadal development was observed only in females but not in males. In females, the side-specific H3K27ac occupancy instructs the differential expression of developmental genes between the pair of gonads and contributes to the development of left but not right gonad. However, in males, the left/right discrepancy of H3K27ac chromatin distribution does not drive the side-biased gene expression or gonad development. The TSA-induced retention of chromatin acetylation causes up-regulation of ovarian developmental genes and increases cell proliferation in the female right gonad.

CONCLUSIONS

We revealed that left/right asymmetry in H3K27ac marked chromatin activation exists in both sexes, but this discrepancy gives rise to asymmetric gonadal development only in females. Other mechanisms overriding the chromatin activation would control the symmetric development of male gonads in chicken.

Watch Out for a Second SNP: Focus on Multi-Nucleotide Variants in Coding Regions and Rescued Stop-Gained

Most single-nucleotide polymorphisms (SNPs) are located in non-coding regions, but the fraction usually studied is harbored in protein-coding regions because potential impacts on proteins are relatively easy to predict by popular tools such as the Variant Effect Predictor. These tools annotate variants independently without considering the potential effect of grouped or haplotypic variations, often called "multi-nucleotide variants" (MNVs). Here, we used a large RNA-seq dataset to survey MNVs, comprising 382 chicken samples originating from 11 populations analyzed in the companion paper in which 9.5M SNPs- including 3.3M SNPs with reliable genotypes-were detected. We focused our study on in-codon MNVs and evaluate their potential mis-annotation. Using GATK HaplotypeCaller read-based phasing results, we identified 2,965 MNVs observed in at least five individuals located in 1,792 genes. We found 41.1% of them showing a novel impact when compared to the effect of their constituent SNPs analyzed separately. The biggest impact variation flux concerns the originally annotated stop-gained consequences, for which around 95% were rescued; this flux is followed by the missense consequences for which 37% were reannotated with a different amino acid. We then present in more depth the rescued stop-gained MNVs and give an illustration in the SLC27A4 gene. As previously shown in human datasets, our results in chicken demonstrate the value of haplotype-aware variant annotation, and the interest to consider MNVs in the coding region, particularly when searching for severe functional consequence such as stop-gained variants.