Mapping splice QTLs reveals distinct transcriptional and post-transcriptional regulatory variation of gene expression and identifies putative alternative splicing variation mediating complex trait variation in pigs

BACKGROUND

Alternative splicing is an important step in gene expression, generating multiple isoforms for the same genes and greatly expanding the diversity of proteomes. Genetic variation in alternative splicing contributes to phenotypic diversity in natural populations. However, the genetic basis of variation in alternative splicing in livestock including pigs remains poorly understood.

RESULTS

In this study, using a Duroc x Pietrain F2 pig population, we performed genome-wide analysis of alternative splicing estimated from stranded RNA-Seq data in skeletal muscle. We characterized the genetic architecture of alternative splicing and compared its basic features with those of overall gene expression. We detected a large number of novel alternative splicing events that were not previously annotated. We found heritability of quantitative alternative splicing scores (percent spliced in or PSI) to be lower than that of overall gene expression. In addition, heritabilities showed little correlation between alternative splicing and overall gene expression. We mapped expression QTLs (eQTLs) and splice QTLs (sQTLs) and found them to be largely non-overlapping. Finally, we integrated sQTL mapping with phenotype QTL (pQTL mapping to identify potential mediator of pQTL effect by alternative splicing.

CONCLUSIONS

Our results suggest that regulatory variation exists at multiple levels and that their genetic controls are distinct, offering opportunities for genetic improvement.

Pig fetal skeletal muscle development is associated with genome-wide DNA hypomethylation and corresponding alterations in transcript and microRNA expression

Fetal myogenesis represents a critical period of porcine skeletal muscle development and requires coordinated expression of thousands of genes. Epigenetic mechanisms, including DNA methylation, drive transcriptional regulation during development; however, these processes are understudied in developing porcine tissues. We performed bisulfite sequencing to assess DNA methylation in pig longissimus dorsi muscle at 41- and 70-days gestation (dg), as well as RNA- and small RNA-sequencing to identify coordinated changes in methylation and expression between myogenic stages. We identified 45 739 differentially methylated regions (DMRs) between stages, and the majority (N = 34 232) were hypomethylated at 70 versus 41 dg. Integration of methylation and transcriptomic data revealed strong associations between differential gene methylation and expression. Differential miRNA methylation was significantly negatively correlated with abundance, and dynamic expression of assayed miRNAs persisted postnatally. Motif analysis revealed significant enrichment of myogenic regulatory factor motifs among hypomethylated regions, suggesting that DNA hypomethylation may function to increase accessibility of muscle-specific transcription factors. We show that developmental DMRs are enriched for GWAS SNPs for muscle- and meat-related traits, demonstrating the potential for epigenetic processes to influence phenotypic diversity. Our results enhance understanding of DNA methylation dynamics of porcine myogenesis and reveal putative cis-regulatory elements governed by epigenetic processes.

Transcriptome annotation of 17 porcine tissues using nanopore sequencing technology

The annotation of animal genomes plays an important role in elucidating molecular mechanisms behind the genetic control of economically important traits. Here, we employed long-read sequencing technology, Oxford Nanopore Technology, to annotate the pig transcriptome across 17 tissues from two Yorkshire littermate pigs. More than 9.8 million reads were obtained from a single flow cell, and 69 781 unique transcripts at 50 108 loci were identified. Of these transcripts, 16 255 were found to be novel isoforms, and 22 344 were found at loci that were novel and unannotated in the Ensembl (release 102) and NCBI (release 106) annotations. Novel transcripts were mostly expressed in cerebellum, followed by lung, liver, spleen, and hypothalamus. By comparing the unannotated transcripts to existing databases, there were 21 285 (95.3%) transcripts matched to the NT database (v5) and 13 676 (61.2%) matched to the NR database (v5). Moreover, there were 4324 (19.4%) transcripts matched to the SwissProt database (v5), corresponding to 11 356 proteins. Tissue-specific gene expression analyses showed that 9749 transcripts were highly tissue-specific, and cerebellum contained the most tissue-specific transcripts. As the same samples were used for the annotation of cis-regulatory elements in the pig genome, the transcriptome annotation generated by this study provides an additional and complementary annotation resource for the Functional Annotation of Animal Genomes effort to comprehensively annotate the pig genome.

1 Associations between Resilience to Early-Life Stress and Behavior in Finishing Gilts

Weaning is an acute early-life stressor and therefore, an excellent model for examining lasting impacts of stress resilience. Our objective was to identify gilts exhibiting resilience or vulnerability to weaning stress and characterize long-term impact on behavior. At weaning, blood samples were collected from all gilt piglets of 17 crossbred litters (n=112) at -1d, 0d, and +4d pre- and post-weaning. Serum cortisol concentrations were quantified using a commercial porcine-specific ELISA kit. For each litter, serum cortisol was used to identify the gilt most capable to return to baseline concentrations by +4d as stress resilient (SR, n=17) and the gilt least capable as stress vulnerable (SV, n=17). SR and SV gilts were followed and mixed into new social cohorts at approximately 8 wk-of-age. Skin wounds counts were recorded at -1d, +1d, and +4d surrounding the mix event. At approximately 12 wk-of-age, SR and SV gilts underwent a novel object test (NOT) to assess reactivity to a novel situation. Measures collected during the NOT included latency to cross 1m line from object, latency to cross 0.5m line from object, latency to touch object, total vocalizations, and percent of vocalizations that were high pitched. A paired t-test was conducted to analyze differences between SR and SV gilts. SR gilts had significantly more skin wounds at 1d post-mixing than SV gilts (P=0.007). No significant differences in vocalizations or latency to cross 1m line from object were found between the two groups. However, SR gilts had significantly shorter latency to cross the 0.5m line (P=0.010) and marginally shorter latency to touch the novel object (P=0.088). We found behavioral differences associated with resilience to weaning stress. Gilts resilient to weaning stress engaged in more agonistic behavior at mixing and were bolder when faced with a novel situation.

Assessment of DNA methylation in porcine immune cells reveals novel regulatory elements associated with cell-specific gene expression and immune capacity traits

Background

Genetics studies in the porcine immune system have enhanced selection practices for disease resistance phenotypes and increased the efficacy of porcine models in biomedical research; however limited functional annotation of the porcine immunome has hindered progress on both fronts. Among epigenetic mechanisms that regulate gene expression, DNA methylation is the most ubiquitous modification made to the DNA molecule and influences transcription factor binding as well as gene and phenotype expression. Human and mouse DNA methylation studies have improved mapping of regulatory elements in these species, but comparable studies in the pig have been limited in scope.

Results

We performed whole-genome bisulfite sequencing to assess DNA methylation patterns in nine pig immune cell populations: CD21⁺ and CD21⁻ B cells, four T cell fractions (CD4⁺, CD8⁺, CD8⁺CD4⁺, and SWC6γδ⁺), natural killer and myeloid cells, and neutrophils. We identified 54,391 cell differentially methylated regions (cDMRs), and clustering by cDMR methylation rate grouped samples by cell lineage. 32,737 cDMRs were classified as cell lowly methylated regions (cLMRs) in at least one cell type, and cLMRs were broadly enriched in genes and regions of intermediate CpG density. We observed strong correlations between differential methylation and expression across immune cell populations, with cell-specific low methylation disproportionately impacting genes exhibiting enriched gene expression in the same cell type. Motif analysis of cLMRs revealed cell type-specific enrichment of transcription factor binding motifs, indicating that cell-specific methylation patterns may influence accessibility by trans-acting factors. Lastly, cDMRs were enriched for immune capacity GWAS SNPs, and many such overlaps occurred within genes known to influence immune cell development and function (CD8B, NDRG1).

Conclusion

Our DNA methylation data improve functional annotation of the porcine genome through characterization of epigenomic regulatory patterns that contribute to immune cell identity and function, and increase the potential for identifying mechanistic links between genotype and phenotype.

Supplementary Information

The online version contains supplementary material available at 10.1186/s12864-022-08773-5.

Pig genome functional annotation enhances the biological interpretation of complex traits and human disease

The functional annotation of livestock genomes is crucial for understanding the molecular mechanisms that underpin complex traits of economic importance, adaptive evolution and comparative genomics. Here, we provide the most comprehensive catalogue to date of regulatory elements in the pig (Sus scrofa) by integrating 223 epigenomic and transcriptomic data sets, representing 14 biologically important tissues. We systematically describe the dynamic epigenetic landscape across tissues by functionally annotating 15 different chromatin states and defining their tissue-specific regulatory activities. We demonstrate that genomic variants associated with complex traits and adaptive evolution in pig are significantly enriched in active promoters and enhancers. Furthermore, we reveal distinct tissue-specific regulatory selection between Asian and European pig domestication processes. Compared with human and mouse epigenomes, we show that porcine regulatory elements are more conserved in DNA sequence, under both rapid and slow evolution, than those under neutral evolution across pig, mouse, and human. Finally, we provide biological insights on tissue-specific regulatory conservation, and by integrating 47 human genome-wide association studies, we demonstrate that, depending on the traits, mouse or pig might be more appropriate biomedical models for different complex traits and diseases.

Time budgets of group-housed pigs in relation to social aggression and production

Commercial producers house growing pigs by sex and weight to allow for efficient use of resources and provide pigs the welfare benefits of interacting with their conspecifics and more freedom of movement. However, the introduction of unfamiliar pigs can cause increased aggression for 24 to 48 h as pigs establish social relationships. To address this issue, a better understanding of pig behavior is needed. The objectives of this study were to quantify time budgets of pigs following introduction into a new social group and how these changed over time and to investigate how social aggression influences the overall time budgets and production parameters. A total of 257 grow-finish Yorkshire barrows across 20 pens were introduced into new social groups at 10 wk of age (~23 kg) and observed for aggression and time budgets of behavior at four periods: immediately after introduction and 3, 6, and 9 wk later. Pigs were observed for the duration of total aggression and initiated aggression (s) for 9 h after introduction and for 4 h at 3, 6, and 9 wk later. Time budgets were created by scan sampling inactive, movement, ingestion, social, and exploration behaviors every 2 min for 4 h in the afternoon and summarizing the proportion of time each behavior was performed by period. The least square means of each behavior were compared across time points. Pigs spent most of their time inactive. In general, the greatest change in pig behavior was observed between introduction and week 3 (P < 0.003), with gradual changes throughout the study period as pigs became more inactive (week 3 vs. week 6: P = 0.209; week 6 vs. week 9: P = 0.007) and spent less time on other behaviors. Pigs' nonaggressive behavior and production parameters were compared with aggression using generalized linear mixed models. The time pigs spent on nonaggressive behaviors was negatively related to aggression (P < 0.045) with few exceptions. Initiated aggression after introduction was negatively related to loin muscle area (P = 0.003). These results show how finishing pigs spend their time in commercial facilities and indicate that behavior continues to change for up to 9 wk after introduction into a new social group. Efforts to reduce chronic levels of aggression should focus on promoting nonaggressive behaviors, such as exploration and movement, after the initial fighting that occurs immediately after introduction has waned, and should be implemented for up to 9 wk after introduction into new social groups.

Estimation of direct and social effects of feeding duration in growing pigs using records from automatic feeding stations

Automatic feeding systems in pig production allow for the recording of individual feeding behavior traits, which might be influenced by the social interactions among individuals. This study fitted mixed models to estimate the direct and social effects on visit duration at the feeder of group-housed pigs. The dataset included 74,413 records of each visit duration time (min) event at the automatic feeder from 135 pigs housed in 14 pens. The sequence of visits at the feeder was employed as a proxy for the social interaction between individuals. To estimate animal effects, the direct effect was apportioned to the animal feeding (feeding pig), and the social effect was apportioned to the animal that entered the feeder immediately after the feeding pig left the feeding station (follower). The data were divided into two subsets: "non-immediate replacement" time (NIRT, N = 6,256), where the follower pig occupied the feeder at least 600 s after the feeding pig left the feeder, and "immediate replacement" time (IRT, N = 58,255), where the elapsed time between replacements was less than or equal to 60 s. The marginal posterior distribution of the parameters was obtained by Bayesian method. Using the IRT subset, the posterior mean of the proportion of variance explained by the direct effect (Prpσ^d2) was 18% for all models. The proportion of variance explained by the follower social effect (Prpσ^f2) was 2%, and the residual variance (σ^e2) decreased, suggesting an improved model fit by including the follower effect. Fitting the models with the NIRT subset, the estimate of Prpσ^d2 was 20% but the Prpσ^f2 was almost zero and σ^e2 was identical for all models. For the IRT subset, the predicted best linear unbiased predictor (BLUP) of direct (Direct BLUP) and social (Follower BLUP) random effects on visit duration at the feeder of an animal was calculated. Feeder visit duration time was not correlated with traits, such as weight gain or average feed intake (P > 0.05), whereas for the daily feeder occupation time, the estimated correlation was positive with the Direct BLUP (r^ = 0.51, P < 0.05) and negative with the Follower BLUP (r^= -0.26, P < 0.05). The results suggest that the visit duration of an animal at the single-space feeder was influenced by both direct and social effects when the replacement time between visits was less than 1 min. Finally, animals that spent a longer time per day at the feeder seemed to do so by shortening the meal length of the preceding individual at the feeder.

Integrated Genome-Wide Analysis of MicroRNA Expression Quantitative Trait Loci in Pig Longissimus Dorsi Muscle

Determining mechanisms regulating complex traits in pigs is essential to improve the production efficiency of this globally important protein source. MicroRNAs (miRNAs) are a class of non-coding RNAs known to post-transcriptionally regulate gene expression affecting numerous phenotypes, including those important to the pig industry. To facilitate a more comprehensive understanding of the regulatory mechanisms controlling growth, carcass composition, and meat quality phenotypes in pigs, we integrated miRNA and gene expression data from longissimus dorsi muscle samples with genotypic and phenotypic data from the same animals. We identified 23 miRNA expression Quantitative Trait Loci (miR-eQTL) at the genome-wide level and examined their potential effects on these important production phenotypes through miRNA target prediction, correlation, and colocalization analyses. One miR-eQTL miRNA, miR-874, has target genes that colocalize with phenotypic QTL for 12 production traits across the genome including backfat thickness, dressing percentage, muscle pH at 24 h post-mortem, and cook yield. The results of our study reveal genomic regions underlying variation in miRNA expression and identify miRNAs and genes for future validation of their regulatory effects on traits of economic importance to the global pig industry.

Weaning Induces Stress-Dependent DNA Methylation and Transcriptional Changes in Piglet PBMCs

Changes to the epigenome, including those to DNA methylation, have been proposed as mechanisms by which stress can induce long-term physiological changes in livestock species. Pig weaning is associated with dietary and social stress, both of which elicit an immune response and changes to the hypothalamic–pituitary–adrenal (HPA) axis. While differential methylation following stress has been assessed in model organisms, it remains poorly understood how the pig methylome is altered by stressors in production settings. We quantified changes in CpG methylation and transcript abundance in piglet peripheral blood mononuclear cells (PBMCs) following weaning and also assessed differential patterns in pigs exhibiting high and low stress response as measured by cortisol concentration and lesion scores. Blood was collected from nine gilt piglets 24 h before and after weaning, and whole-genome bisulfite sequencing (WGBS) and RNA-sequencing were performed on six and nine animals, respectively, at both time points. We identified 2,674 differentially methylated regions (DMRs) that were enriched within promoters of genes associated with lymphocyte stimulation and transcriptional regulation. Stress groups displayed unique differential methylation and expression patterns associated with activation and suppression of T cell immunity in low and high stress animals, respectively. Differential methylation was strongly associated with differential expression; specifically, upregulated genes were enriched among hypomethylated genes. We observed post-weaning hypermethylation of the glucocorticoid receptor (NR3C1) promoter and a significant decrease in NR3C1 expression (n = 9, p = 6.1 × 10^–3). Our results indicate that weaning-associated stress elicits genome-wide methylation changes associated with differential gene expression, reduced T cell activation, and an altered HPA axis response.

Genome-Wide Assessment of DNA Methylation in Chicken Cardiac Tissue Exposed to Different Incubation Temperatures and CO2 Levels

Temperature and CO₂ concentration during incubation have profound effects on broiler chick development, and numerous studies have identified significant effects on hatch heart weight (HW) as a result of differences in these parameters. Early life environment has also been shown to affect broiler performance later in life; it has thus been suggested that epigenetic mechanisms may mediate long-term physiological changes induced by environmental stimuli. DNA methylation is an epigenetic modification that can confer heritable changes in gene expression. Using reduced-representation bisulfite sequencing (RRBS), we assessed DNA methylation patterns in cardiac tissue of 84 broiler hatchlings incubated at two egg shell temperatures (EST; 37.8°C and 38.9°C) and three CO₂ concentrations (0.1%, 0.4%, and 0.8%) from day 8 of incubation onward. We assessed differential methylation between EST treatments and identified 2,175 differentially methylated (DM) CpGs (1,121 hypermethylated, 1,054 hypomethylated at 38.9° vs. 37.8°) in 269 gene promoters and 949 intragenic regions. DM genes (DMGs) were associated with heart developmental processes, including cardiomyocyte proliferation and differentiation. We identified enriched binding motifs among DM loci, including those for transcription factors associated with cell proliferation and heart development among hypomethylated CpGs that suggest increased binding ability at higher EST. We identified 9,823 DM CpGs between at least two CO₂ treatments, with the greatest difference observed between 0.8 and 0.1% CO₂ that disproportionately impacted genes involved in cardiac muscle development and response to low oxygen levels. Using HW measurements from the same chicks, we performed an epigenome-wide association study (EWAS) for HW, and identified 23 significantly associated CpGs, nine of which were also DM between ESTs. We found corresponding differences in transcript abundance between ESTs in three DMGs (ABLIM2, PITX2, and THRSP). Hypomethylation of an exonic CpG in PITX2 at 38.9°C was associated with increased expression, and suggests increased cell proliferation in broiler hatchlings incubated at higher temperatures. Overall, these results identified numerous epigenetic associations between chick incubation factors and heart development that may manifest in long-term differences in animal performance.

A comparative analysis of chromatin accessibility in cattle, pig, and mouse tissues

Background

Although considerable progress has been made towards annotating the noncoding portion of the human and mouse genomes, regulatory elements in other species, such as livestock, remain poorly characterized. This lack of functional annotation poses a substantial roadblock to agricultural research and diminishes the value of these species as model organisms. As active regulatory elements are typically characterized by chromatin accessibility, we implemented the Assay for Transposase Accessible Chromatin (ATAC-seq) to annotate and characterize regulatory elements in pigs and cattle, given a set of eight adult tissues.

Results

Overall, 306,304 and 273,594 active regulatory elements were identified in pig and cattle, respectively. 71,478 porcine and 47,454 bovine regulatory elements were highly tissue-specific and were correspondingly enriched for binding motifs of known tissue-specific transcription factors. However, in every tissue the most prevalent accessible motif corresponded to the insulator CTCF, suggesting pervasive involvement in 3-D chromatin organization. Taking advantage of a similar dataset in mouse, open chromatin in pig, cattle, and mice were compared, revealing that the conservation of regulatory elements, in terms of sequence identity and accessibility, was consistent with evolutionary distance; whereas pig and cattle shared about 20% of accessible sites, mice and ungulates only had about 10% of accessible sites in common. Furthermore, conservation of accessibility was more prevalent at promoters than at intergenic regions.

Conclusions

The lack of conserved accessibility at distal elements is consistent with rapid evolution of enhancers, and further emphasizes the need to annotate regulatory elements in individual species, rather than inferring elements based on homology. This atlas of chromatin accessibility in cattle and pig constitutes a substantial step towards annotating livestock genomes and dissecting the regulatory link between genome and phenome.

Deciphering Sex-Specific Genetic Architectures Using Local Bayesian Regressions

Many complex human traits exhibit differences between sexes. While numerous factors likely contribute to this phenomenon, growing evidence from genome-wide studies suggest a partial explanation: that males and females from the same population possess differing genetic architectures. Despite this, mapping gene-by-sex (G×S) interactions remains a challenge likely because the magnitude of such an interaction is typically and exceedingly small; traditional genome-wide association techniques may be underpowered to detect such events, due partly to the burden of multiple test correction. Here, we developed a local Bayesian regression (LBR) method to estimate sex-specific SNP marker effects after fully accounting for local linkage-disequilibrium (LD) patterns. This enabled us to infer sex-specific effects and G×S interactions either at the single SNP level, or by aggregating the effects of multiple SNPs to make inferences at the level of small LD-based regions. Using simulations in which there was imperfect LD between SNPs and causal variants, we showed that aggregating sex-specific marker effects with LBR provides improved power and resolution to detect G×S interactions over traditional single-SNP-based tests. When using LBR to analyze traits from the UK Biobank, we detected a relatively large G×S interaction impacting bone mineral density within ABO, and replicated many previously detected large-magnitude G×S interactions impacting waist-to-hip ratio. We also discovered many new G×S interactions impacting such traits as height and body mass index (BMI) within regions of the genome where both male- and female-specific effects explain a small proportion of phenotypic variance (R2 < 1 × 10-4), but are enriched in known expression quantitative trait loci.

Estimation of indirect social genetic effects for skin lesion count in group-housed pigs by quantifying behavioral interactions1

Mixing of pigs into new social groups commonly induces aggressive interactions that result in skin lesions on the body of the animals. The relationship between skin lesions and aggressive behavioral interactions in group-housed pigs can be analyzed within the framework of social genetic effects (SGE). This study incorporates the quantification of aggressive interactions between pairs of animals in the modeling of SGE for skin lesions in different regions of the body in growing pigs. The dataset included 792 pigs housed in 59 pens. Skin lesions in the anterior, central, and caudal regions of the body were counted 24 h after pig mixing. Animals were video-recorded for 9 h postmixing and trained observers recorded the type and duration of aggressive interactions between pairs of animals. The number of seconds that pairs of pigs spent engaged in reciprocal fights and unilateral attack behaviors were used to parametrize the intensity of social interactions (ISI). Three types of models were fitted: direct genetic additive model (DGE), traditional social genetic effect model (TSGE) assuming uniform interactions between dyads, and an intensity-based social genetic effect model (ISGE) that used ISI to parameterize SGE. All models included fixed effects of sex, replicate, lesion scorer, weight at mixing, premixing lesion count, and the total time that the animal spent engaged in aggressive interactions (reciprocal fights and unilateral attack behaviors) as a covariate; a random effect of pen; and a random direct genetic effect. The ISGE models recovered more direct genetic variance than DGE and TSGE, and the estimated heritabilities (h^D2) were highest for all traits (P < 0.01) for the ISGE with ISI parametrized with unilateral attack behavior. The TSGE produced estimates that did not differ significantly from DGE (P > 0.5). Incorporating the ISI into ISGE, even in a small dataset, allowed separate estimation of the genetic parameters for direct and SGE, as well as the genetic correlation between direct and SGE (r^ds), which was positive for all lesion traits. The estimates from ISGE suggest that if behavioral observations are available, selection incorporating SGE may reduce the consequences of aggressive behaviors after mixing pigs.

Genetic control of longissimus dorsi muscle gene expression variation and joint analysis with phenotypic quantitative trait loci in pigs

Background

Economically important growth and meat quality traits in pigs are controlled by cascading molecular events occurring during development and continuing throughout the conversion of muscle to meat. However, little is known about the genes and molecular mechanisms involved in this process. Evaluating transcriptomic profiles of skeletal muscle during the initial steps leading to the conversion of muscle to meat can identify key regulators of polygenic phenotypes. In addition, mapping transcript abundance through genome-wide association analysis using high-density marker genotypes allows identification of genomic regions that control gene expression, referred to as expression quantitative trait loci (eQTL). In this study, we perform eQTL analyses to identify potential candidate genes and molecular markers regulating growth and meat quality traits in pigs.

Results

Messenger RNA transcripts obtained with RNA-seq of longissimus dorsi muscle from 168 F2 animals from a Duroc x Pietrain pig resource population were used to estimate gene expression variation subject to genetic control by mapping eQTL. A total of 339 eQTL were mapped (FDR ≤ 0.01) with 191 exhibiting local-acting regulation. Joint analysis of eQTL with phenotypic QTL (pQTL) segregating in our population revealed 16 genes significantly associated with 21 pQTL for meat quality, carcass composition and growth traits. Ten of these pQTL were for meat quality phenotypes that co-localized with one eQTL on SSC2 (8.8-Mb region) and 11 eQTL on SSC15 (121-Mb region). Biological processes identified for co-localized eQTL genes include calcium signaling (FERM, MRLN, PKP2 and CHRNA9), energy metabolism (SUCLG2 and PFKFB3) and redox hemostasis (NQO1 and CEP128), and results support an important role for activation of the PI3K-Akt-mTOR signaling pathway during the initial conversion of muscle to meat.

Conclusion

Co-localization of eQTL with pQTL identified molecular markers significantly associated with both economically important phenotypes and gene transcript abundance. This study reveals candidate genes contributing to variation in pig production traits, and provides new knowledge regarding the genetic architecture of meat quality phenotypes.

Electronic supplementary material

The online version of this article (10.1186/s12864-018-5386-2) contains supplementary material, which is available to authorized users.

Evidence of evolutionary history and selective sweeps in the genome of Meishan pig reveals its genetic and phenotypic characterization

Background

Meishan is a pig breed indigenous to China and famous for its high fecundity. The traits of Meishan are strongly associated with its distinct evolutionary history and domestication. However, the genomic evidence linking the domestication of Meishan pigs with its unique features is still poorly understood. The goal of this study is to investigate the genomic signatures and evolutionary evidence related to the phenotypic traits of Meishan via large-scale sequencing.

Results

We found that the unique domestication of Meishan pigs occurred in the Taihu Basin area between the Majiabang and Liangzhu Cultures, during which 300 protein-coding genes have underwent positive selection. Notably, enrichment of the FoxO signaling pathway with significant enrichment signal and the harbored gene IGF1R were likely associated with the high fertility of Meishan pigs. Moreover, NFKB1 exhibited strong selective sweep signals and positively participated in hyaluronan biosynthesis as the key gene of NF-kB signaling, which may have resulted in the wrinkled skin and face of Meishan pigs. Particularly, three population-specific synonymous single-nucleotide variants occurred in PYROXD1, MC1R, and FAM83G genes; the T305C substitution in the MCIR gene explained the black coat of the Meishan pigs well. In addition, the shared haplotypes between Meishan and Duroc breeds confirmed the previous Asian-derived introgression and demonstrated the specific contribution of Meishan pigs.

Conclusions

These findings will help us explain the unique genetic and phenotypic characteristics of Meishan pigs and offer a plausible method for their utilization of Meishan pigs as valuable genetic resources in pig breeding and as an animal model for human wrinkled skin disease research.

Genome-wide association analyses of lesion counts in group-housed pigs

Aggression in group-housed pigs is a welfare concern and can negatively affect production. Skin lesions are reliable indicators of aggression and are moderately heritable, suggesting that selective breeding may reduce aggression. To further understand the genetic control of behavioral traits, such as the aggressive response to regrouping, associated single nucleotide polymorphisms (SNPs) can be identified within the genome, and the region in which these SNPs are located can be related to known genes. To investigate SNPs associated with aggression, 1093 purebred Yorkshire pigs were strategically remixed into new groups of familiar and unfamiliar animals at three life stages and lesion counts were recorded. Genomic best linear unbiased prediction (GBLUP) models were fitted for each trait. The genetic additive effect was obtained from a genetic relationship matrix constructed from the 50 924 SNPs. SNP effects and their variances were estimated from the GBLUP objects. SNPs that were associated with a significant portion of the trait variance were identified for lesions to the anterior (three SNPs, FDR <5%) and central (one SNP, FDR <5%) portions of the body in grow-finish pigs. These SNPs were located on chromosome 11, suggesting that chromosome 11 contains a region explaining variation in lesion counts that should be further explored to identify genes underlying biological control of aggression.

Estimation of genetic parameters for lesion scores and growth traits in group-housed pigs

Pigs housed in groups are remixed with unfamiliar individuals, which can trigger aggressive interactions, potentially compromising animal welfare. Skin lesions are a reliable indicator trait of aggression and are moderately heritable, suggesting that aggression may be reduced through selection. This study estimated genetic parameters of skin lesions of pigs at multiple life stages, explored genetic correlations of skin lesions between age groups and body location, and studied the relationship between skin lesions and production traits of commercial importance. A population of 1,079 Yorkshire pigs was strategically remixed into new groups of familiar and unfamiliar animals at 3 life stages (weaning, grow-finish, and mature gilts). Skin lesions (fresh, bright red cuts) were counted immediately prior to mixing and 24 h and 3 wk after mixing across 3 body regions: anterior, central, and caudal. Weights were recorded prior to each mixing event. Prior to slaughter, backfat thickness and loin muscle area were determined using ultrasound. Univariate analyses were performed to obtain heritability estimates of lesion scores. Bivariate analyses were performed with response variables being skin lesions, weight gain per life stage, backfat thickness, or loin muscle area, depending on the relationship of interest, to obtain correlations. Lesion score heritabilities ranged from 0.10 to 0.40 and were significant ( < 0.05). Heritability was highest for lesions on the anterior region of the body for 24 h and 3 wk after mixing. Lesions to the central and caudal areas showed the highest genetic correlation at each stage of production, whereas those to the anterior and caudal regions had the lowest correlation. The highest genetic correlation was found between the mature gilt and grow-finish stages, whereas the weaning and mature gilt stages had the lowest correlations. Genetic correlations between lesions and production traits were not significantly different from 0 for weight gain and backfat thickness, but loin muscle area was negatively correlated with lesions ( = 1.17 × 10, = 2.30 × 10, and = 6.08 × 10 for anterior, central, and caudal lesions, respectively). These results are promising for the industry because they suggest that pigs selected for reduced lesions will show increased loin muscle area without negative effects on growth. Alternatively, selection for these production traits would not increase lesions.

Genome-wide association study in an F2 Duroc x Pietrain resource population for economically important meat quality and carcass traits

Meat quality is essential for consumer acceptance, it ultimately impacts pork production profitability and it is subject to genetic control. The objective of this study was to map genomic regions associated with economically important meat quality and carcass traits. We performed a genome-wide association (GWA) analysis to map regions associated with 38 meat quality and carcass traits recorded for 948 F2 pigs from the Michigan State University Duroc × Pietrain resource population. The F0, F1, and 336 F2 pigs were genotyped with the Illumina Porcine SNP60 BeadChip, while the remaining F2 pigs were genotyped with the GeneSeek Genomic Profiler for Porcine Low Desnisty (LD) chip, and imputed with high accuracy ( = 0.97). Altogether the genomic dataset comprised 1,019 animals and 44,911 SNP. A Gaussian linear mixed model was fitted to estimate the breeding values and the variance components. A linear transformation was performed to estimate the marker effects and variances. Type I error rate was controlled at a False Discovery Rate of 5%. Seven putative QTL found in this study were previously reported in other studies. Two novel QTL associated with tenderness (TEN) were located on SSC3 [135.6:137.5Mb; False Discovery rate (FDR) < 0.03] and SSC5 (67.3:69.1Mb; FDR < 0.02). The QTL region identified on SSC15 includes Protein Kinase AMP-activated ɣ 3-subunit gene (), which has been associated with 24-h pH (pH24), drip loss (DL) and cook yield (CY). Also, novel candidate genes were identified for TEN in the region on SSC5 [A Kinase (PRKA) Anchor Protein 3 (], and for tenth rib backfat thickness (BF10) [Carnitine O-Acetyltransferase ()] on SSC1. The association of gene polymorphisms with pork quality traits has been reported for several pig populations. However, there are no SNP for this gene on the chip used, thus we genotyped the animals for 2 non-synonymous variants ( and ). We then performed a GWA conditioning on the genotype of both SNP and was associated with pH24, DL, protein content (PRO) and CY ( < 0.004) and T30N with Juiciness, TEN, shear force, pH24, PRO, and CY < 0.04). Finally, we performed a GWA conditioning on the genotype of the SNP peak detected in this study, and T30N remained associated only with PRO ( < 0.02). Therefore, in this study we identified 2 novel QTL regions, suggest 2 novel candidate genes, and conclude that other SNP in PRKAG3 or nearby gene(s) explain the observed associations on SSC15 in this population.

Profiling and characterization of a longissimus dorsi muscle microRNA dataset from an F2 Duroc × Pietrain pig resource population

To elucidate the effects of microRNA (miRNA) regulation in skeletal muscle of adult pigs, miRNA expression profiling was performed with RNA extracted from longissimus dorsi (LD) muscle samples from 174 F₂ pigs (~ 5.5 months of age) from a Duroc × Pietrain resource population. Total RNA was extracted from LD samples, and libraries were sequenced on an Illumina HiSeq 2500 platform in 1 × 50 bp format. After processing, 232,826,977 total reads were aligned to the Sus scrofa reference genome (v10.2.79), with 74.8% of total reads mapping successfully. The miRDeep2 software package was utilized to quantify annotated Sus scrofa mature miRNAs from miRBase (Release 21) and to predict candidate novel miRNA precursors. Among the retained 295 normalized mature miRNA expression profiles ssc­miR­1, ssc­miR­133a­3p, ssc­miR­378, ssc­miR­206, and ssc­miR­10b were the most abundant, all of which have previously been shown to be expressed in pig skeletal muscle. Additionally, 27 unique candidate novel miRNA precursors were identified exhibiting homologous sequence to annotated human miRNAs. The composition of classes of small RNA present in this dataset was also characterized; while the majority of unique expressed sequence tags were not annotated in any of the queried databases, the most abundantly expressed class of small RNA in this dataset was miRNAs. This data provides a resource to evaluate miRNA regulation of gene expression and effects on complex trait phenotypes in adult pig skeletal muscle. The raw sequencing data were deposited in the Sequence Read Archive, BioProject PRJNA363073.

Refining genomewide association for growth and fat deposition traits in an F pig population

The identification of genomic regions that affect additive genetic variation and contain genes involved in controlling growth and fat deposition has enormous impact in the farm animal industry (e.g., carcass merit and meat quality). Therefore, a genomewide association study was implemented in an F pig population using a 60,000 SNP marker panel for traits related to growth and fat deposition. Estimated genomic EBV were linearly transformed to calculate SNP effects and to identify genomic positions possibly associated with the genetic variability of each trait. Genomic segments were then defined considering the markers included in a region 1 Mb up- and downstream from the SNP with the smallest -value and a false discovery rate < 0.05 for each trait. The significance for each 2-Mb segment was tested using the Bonferroni correction. Significant SNP were detected on SSC2, SSC3, SSC5, and SSC6, but 2-Mb segment significant effects were observed on SSC3 for weight at birth (wt_birth) and on SSC6 for 10th-rib backfat and last-rib backfat measured by ultrasound at different ages. Furthermore, a 6-Mb segment on SSC6 was also considered because the 2-Mb segments for 10 different fat deposition traits were overlapped. Although the segment effects for each trait remain significant, the proportion of additive variance explained by this larger segment was slightly smaller in some traits. In general, the results confirm the presence of genetic variability for wt_birth on SSC3 (18.0-20.2 Mb) and for fat deposition traits on SSC6 (133.8-136.0 Mb). Within these regions, fibrosin () and myosin light chain, phosphorylatable, fast skeletal muscle () genes could be considered as candidates for the wt_birth signal on SSC3, and the SERPINE1 mRNAbinding protein 1 gene () may be a candidate for the fat deposition trait signals on SSC6.

Evidence for transcriptome-wide RNA editing among Sus scrofa PRE-1 SINE elements

Background

RNA editing by ADAR (adenosine deaminase acting on RNA) proteins is a form of transcriptional regulation that is widespread among humans and other primates. Based on high-throughput scans used to identify putative RNA editing sites, ADAR appears to catalyze a substantial number of adenosine to inosine transitions within repetitive regions of the primate transcriptome, thereby dramatically enhancing genetic variation beyond what is encoded in the genome.

Results

Here, we demonstrate the editing potential of the pig transcriptome by utilizing DNA and RNA sequence data from the same pig. We identified a total of 8550 mismatches between DNA and RNA sequences across three tissues, with 75% of these exhibiting an A-to-G (DNA to RNA) discrepancy, indicative of a canonical ADAR-catalyzed RNA editing event. When we consider only mismatches within repetitive regions of the genome, the A-to-G percentage increases to 94%, with the majority of these located within the swine specific SINE retrotransposon PRE-1. We also observe evidence of A-to-G editing within coding regions that were previously verified in primates.

Conclusions

Thus, our high-throughput evidence suggests that pervasive RNA editing by ADAR can exist outside of the primate lineage to dramatically enhance genetic variation in pigs.

Estimation of genome-wide and locus-specific breed composition in pigs

Advances in pig genomic technologies enable implementation of new methods to estimate breed composition, allowing innovative and efficient ways to evaluate and ensure breed and line background. Existing methods to test for homozygosity at key loci involve test mating the animal in question and observing phenotypic patterns among offspring, requiring extensive resources. In this study, whole-genome pig DNA microarray data from over 8,000 SNP was used to profile the composition of U.S. registered purebred pigs using a refined linear regression method that enhances the interpretation of coefficients. In a simulation analysis, a strong correlation between true and estimated breed composition was observed (R² = 0.94). Applying these methods to 930 Yorkshire animals registered with the National Swine Registry, 95% were estimated to have a “genome-wide” Yorkshire breed composition of at least 0.825 or 82.5%, with similar performance for evaluating datasets of registered Duroc (n = 88) Landrace (n = 129), and Hampshire (n = 17) breeds. We also developed new methods to evaluate locus-based breed probabilities. Such methods have been applied to multi-locus SNP genotypes flanking the KIT gene known to predominantly control coat color, thereby inferring the probability that an animal has haplotypes in the KIT region that are predominant in white breeds. These methods have been adopted by the National Swine Registry as a means to identify purebred Yorkshire animals.

Impact of improving dietary amino acid balance for lactating sows on efficiency of dietary amino acid utilization and transcript abundance of genes encoding lysine transporters in mammary tissue

Lactating multiparous Yorkshire sows ( = 64) were used in 2 experiments to test the hypothesis that reducing dietary CP intake and improving AA balance through crystalline AA (CAA) supplementation improves apparent dietary AA utilization efficiency for milk production and increases transcript abundance of genes encoding Lys transporter proteins in mammary tissue. In Exp. 1, 40 sows were assigned to 1 of 4 diets: 1) high CP (HCP; 16.0% CP, as-fed basis; analyzed concentration), 2) medium-high CP (MHCP; 15.7% CP), 3) medium-low CP (MLCP; 14.3% CP), and 4) low CP (LCP; 13.2% CP). The HCP diet was formulated using soybean meal and corn as the only Lys sources. The reduced-CP diets contained CAA to meet estimated requirements for essential AA that became progressively limiting with reduction in CP concentration, that is, Lys, Ile, Met + Cys, Thr, Trp, and Val. Dietary standardized ileal digestible (SID) Lys concentration was 80% of the estimated requirement. In Exp. 2, 24 sows were assigned to the HCP or LCP diets. In Exp. 1, blood samples were postprandially collected 15 h on d 3, 7, 14, and 18 of lactation and utilization efficiency of dietary AA for milk production was calculated during early (d 3 to 7) and peak (d 14 to 18) lactation. Efficiency values were estimated from daily SID AA intakes and milk AA yield, with corrections for maternal AA requirement for maintenance and AA contribution from body protein losses. In Exp. 2, mammary tissue was biopsied on d 4 and 14 of lactation to determine the mRNA abundance of genes encoding Lys transporter proteins. In peak lactation, Lys, Thr, Trp, and Val utilization efficiency increased with decreasing dietary CP (linear for Trp and Val, < 0.05; in sows fed the MHCP diet vs. sows fed the HCP diet for Lys and Thr, < 0.05). Total essential and nonessential 15-h postprandial serum AA concentrations increased with decreasing dietary CP (linear, = 0.09 and < 0.05, respectively), suggesting increased maternal body protein mobilization. Transcript abundance of several genes involved in Lys transport in mammary tissue did not differ between sows fed the LCP and HCP diets. Feeding lactating sows low-CP diets supplemented with CAA increases the efficiency of utilizing dietary Lys, Thr, Trp, and Val for milk protein production but is unrelated to abundance in mRNA of genes encoding Lys transport proteins in the mammary gland. Dietary Lys utilization for milk protein production in lactating sows appears to be optimized when crystalline Lys is included at a minimum of 0.10% in a diet containing 15.70% CP.

Estimates of the actual relationship between half-sibs in a pig population

Genomic relationships based on markers capture the actual instead of the expected (based on pedigree) proportion of genome shared identical by descent (IBD). Several methods exist to estimate genomic relationships. In this research, we compare four such methods that were tested looking at the empirical distribution of the estimated relationships across 6704 pairs of half-sibs from a cross-bred pig population. The first method based on multiple marker linkage analysis displayed a mean and standard deviation (SD) in close agreement with the expected ones and was robust to changes in the minor allele frequencies (MAF). A single marker method that accounts for linkage disequilibrium (LD) and inbreeding came second, showing more sensitivity to changes in the MAF. Another single marker method that considers neither inbreeding nor LD showed the smallest empirical SD and was the most sensible to changes in MAF. A higher mean and SD were displayed by VanRaden's method, which was not sensitive to changes in MAF. Therefore, the method based on multiple marker linkage analysis and the single marker method that considers LD and inbreeding performed closer to theoretical values and were consistent with the estimates reported in literature for human half-sibs.

Implementing meta-analysis from genome-wide association studies for pork quality traits

Pork quality plays an important role in the meat processing industry. Thus, different methodologies have been implemented to elucidate the genetic architecture of traits affecting meat quality. One of the most common and widely used approaches is to perform genome-wide association (GWA) studies. However, a limitation of many GWA in animal breeding is the limited power due to small sample sizes in animal populations. One alternative is to implement a meta-analysis of GWA (MA-GWA) combining results from independent association studies. The objective of this study was to identify significant genomic regions associated with meat quality traits by performing MA-GWA for 8 different traits in 3 independent pig populations. Results from MA-GWA were used to search for genes possibly associated with the set of evaluated traits. Data from 3 pig data sets (U.S. Meat Animal Research Center, commercial, and Michigan State University Pig Resource Population) were used. A MA was implemented by combining -scores derived for each SNP in every population and then weighting them using the inverse of estimated variance of SNP effects. A search for annotated genes retrieved genes previously reported as candidates for shear force (calpain-1 catalytic subunit [] and calpastatin []), as well as for ultimate pH, purge loss, and cook loss (protein kinase, AMP-activated, γ 3 noncatalytic subunit []). In addition, novel candidate genes were identified for intramuscular fat and cook loss (acyl-CoA synthetase family member 3 mitochondrial []) and for the objective measure of muscle redness, CIE a* (glycogen synthase 1, muscle [] and ferritin, light polypeptide []). Thus, implementation of MA-GWA allowed integration of results for economically relevant traits and identified novel genes to be tested as candidates for meat quality traits in pig populations.

Searching for causal networks involving latent variables in complex traits: Application to growth, carcass, and meat quality traits in pigs

Structural equation models (SEQM) can be used to model causal relationships between multiple variables in multivariate systems. Among the strengths of SEQM is its ability to consider causal links between latent variables. The use of latent variables allows modeling complex phenomena while reducing at the same time the dimensionality of the data. One relevant aspect in the quantitative genetics context is the possibility of correlated genetic effects influencing sets of variables under study. Under this scenario, if one aims at inferring causality among latent variables, genetic covariances act as confounders if ignored. Here we describe a methodology for assessing causal networks involving latent variables underlying complex phenotypic traits. The first step of the method consists of the construction of latent variables defined on the basis of prior knowledge and biological interest. These latent variables are jointly evaluated using confirmatory factor analysis. The estimated factor scores are then used as phenotypes for fitting a multivariate mixed model to obtain the covariance matrix of latent variables conditional on the genetic effects. Finally, causal relationships between the adjusted latent variables are evaluated using different SEQM with alternative causal specifications. We have applied this method to a data set with pigs for which several phenotypes were recorded over time. Five different latent variables were evaluated to explore causal links between growth, carcass, and meat quality traits. The measurement model, which included 5 latent variables capturing the information conveyed by 19 different phenotypic traits, showed an acceptable fit to data (e.g., χ/df = 1.3, root-mean-square error of approximation = 0.028, standardized root-mean-square residual = 0.041). Causal links between latent variables were explored after removing genetic confounders. Interestingly, we found that both growth (-0.160) and carcass traits (-0.500) have a significant negative causal effect on quality traits (-value ≤ 0.001). This result may have important implications for strategies for pig production improvement. More generally, the proposed method allows further learning regarding phenotypic causal structures underlying complex traits in farm species.

A comparison of methods to estimate genomic relationships using pedigree and markers in livestock populations

Accurate prediction of breeding values depends on capturing the variability in genome sharing of relatives with the same pedigree relationship. Here, we compare two approaches to set up genomic relationship matrices for precision of genomic relationships (GR) and accuracy of estimated breeding values (GEBV). Real and simulated data (pigs, 60k SNP) were analysed, and GR were estimated using two approaches: (i) identity by state, corrected with either the observed (G_VR-O ) or the base population (G_VR-B ) allele frequencies and (ii) identity by descent using linkage analysis (G_IBD-L ). Estimators were evaluated for precision and empirical bias with respect to true pedigree IBD GR. All three estimators had very low bias. G_IBD-L displayed the lowest sampling error and the highest correlation with true genome-shared values. G_VR-B approximated G_IBD-L 's correlation and had lower error than G_VR-O . Accuracy of GEBV for selection candidates was significantly higher when G_IBD-L was used and identical between G_VR-O and G_VR-B . In real data, G_IBD-L 's sampling standard deviation was the closest to the theoretical value for each pedigree relationship. Use of pedigree to calculate GR improved the precision of estimates and the accuracy of GEBV.

Meta-analysis of genome-wide association from genomic prediction models

Genome-wide association (GWA) studies based on GBLUP models are a common practice in animal breeding. However, effect sizes of GWA tests are small, requiring larger sample sizes to enhance power of detection of rare variants. Because of difficulties in increasing sample size in animal populations, one alternative is to implement a meta-analysis (MA), combining information and results from independent GWA studies. Although this methodology has been used widely in human genetics, implementation in animal breeding has been limited. Thus, we present methods to implement a MA of GWA, describing the proper approach to compute weights derived from multiple genomic evaluations based on animal-centric GBLUP models. Application to real datasets shows that MA increases power of detection of associations in comparison with population-level GWA, allowing for population structure and heterogeneity of variance components across populations to be accounted for. Another advantage of MA is that it does not require access to genotype data that is required for a joint analysis. Scripts related to the implementation of this approach, which consider the strength of association as well as the sign, are distributed and thus account for heterogeneity in association phase between QTL and SNPs. Thus, MA of GWA is an attractive alternative to summarizing results from multiple genomic studies, avoiding restrictions with genotype data sharing, definition of fixed effects and different scales of measurement of evaluated traits.

Exploring causal networks underlying fat deposition and muscularity in pigs through the integration of phenotypic, genotypic and transcriptomic data

BACKGROUND

Joint modeling and analysis of phenotypic, genotypic and transcriptomic data have the potential to uncover the genetic control of gene activity and phenotypic variation, as well as shed light on the manner and extent of connectedness among these variables. Current studies mainly report associations, i.e. undirected connections among variables without causal interpretation. Knowledge regarding causal relationships among genes and phenotypes can be used to predict the behavior of complex systems, as well as to optimize management practices and selection strategies. Here, we performed a multistep procedure for inferring causal networks underlying carcass fat deposition and muscularity in pigs using multi-omics data obtained from an F2 Duroc x Pietrain resource pig population.

RESULTS

We initially explored marginal associations between genotypes and phenotypic and expression traits through whole-genome scans, and then, in genomic regions with multiple significant hits, we assessed gene-phenotype network reconstruction using causal structural learning algorithms. One genomic region on SSC6 showed significant associations with three relevant phenotypes, off-midline10th-rib backfat thickness, loin muscle weight, and average intramuscular fat percentage, and also with the expression of seven genes, including ZNF24, SSX2IP, and AKR7A2. The inferred network indicated that the genotype affects the three phenotypes mainly through the expression of several genes. Among the phenotypes, fat deposition traits negatively affected loin muscle weight.

CONCLUSIONS

Our findings shed light on the antagonist relationship between carcass fat deposition and lean meat content in pigs. In addition, the procedure described in this study has the potential to unravel gene-phenotype networks underlying complex phenotypes.

Association of Insulin-like growth factor binding protein 2 genotypes with growth, carcass and meat quality traits in pigs

Background

This study was conducted to investigate the potential association of variation in the insulin-like growth factor binding protein 2 (IGFBP2) gene with growth, carcass and meat quality traits in pigs. IGFBP2 is a member of the insulin-like growth factor binding protein family that is involved in regulating growth, and it maps to a region of pig chromosome 15 containing significant quantitative trait loci that affect economically important trait phenotypes.

Results

An IGFBP2 polymorphism was identified in the Michigan State University (MSU) Duroc × Pietrain F₂ resource population (n = 408), and pigs were genotyped by MspI PCR-RFLP. Subsequently, a Duroc pig population from the National Swine Registry, USA, (n = 326) was genotyped using an Illumina Golden Gate assay. The IGFBP2 genotypic frequencies among the MSU resource population pigs were 3.43, 47.06 and 49.51 % for the AA, AB and BB genotypes, respectively. The genotypic frequencies for the Duroc pigs were 9.82, 47.85, and 42.33 % for the AA, AB and BB genotypes, respectively. Genotype effects (P < 0.05) were found in the MSU resource population for backfat thickness at 10^th rib and last rib as determined by ultrasound at 10, 13, 16 and 19 weeks of age, ADG from 10 to 22 weeks of age, and age to reach 105 kg. A genotype effect (P < 0.05) was also found for off test Longissimus muscle area in the Duroc population. Significant effects of IGFBP2 genotype (P < 0.05) were found for drip loss, 24 h postmortem pH, pH decline from 45 min to 24 h postmortem, subjective color score, CIE L* and b*, Warner-Bratzler shear force, and sensory panel scores for juiciness, tenderness, connective tissue and overall tenderness in MSU resource population pigs. Genotype effects (P < 0.05) were found for 45-min pH, CIE L* and color score in the Duroc population.

Conclusions

Results of this study revealed associations of the IGFBP2 genotypes with growth, carcass and meat quality traits in pigs. The results indicate IGFBP2 as a potential candidate gene for growth rate, backfat thickness, loin muscle area and some pork quality traits.

Whole blood microarray analysis of pigs showing extreme phenotypes after a porcine reproductive and respiratory syndrome virus infection

Background

The presence of variability in the response of pigs to Porcine Reproductive and Respiratory Syndrome virus (PRRSv) infection, and recent demonstration of significant genetic control of such responses, leads us to believe that selection towards more disease resistant pigs could be a valid strategy to reduce its economic impact on the swine industry. To find underlying molecular differences in PRRS susceptible versus more resistant pigs, 100 animals with extremely different growth rates and viremia levels after PRRSv infection were selected from a total of 600 infected pigs. A microarray experiment was conducted on whole blood RNA samples taken at 0, 4 and 7 days post infection (dpi) from these pigs. From these data, we examined associations of gene expression with weight gain and viral load phenotypes. The single nucleotide polymorphism (SNP) marker WUR10000125 (WUR) on the porcine 60 K SNP chip was shown to be associated with viral load and weight gain after PRRSv infection, and so the effect of the WUR10000125 (WUR) genotype on expression in whole blood was also examined.

Results

Limited information was obtained through linear modeling of blood gene differential expression (DE) that contrasted pigs with extreme phenotypes, for growth or viral load or between animals with different WUR genotype. However, using network-based approaches, molecular pathway differences between extreme phenotypic classes could be identified. Several gene clusters of interest were found when Weighted Gene Co-expression Network Analysis (WGCNA) was applied to 4dpi contrasted with 0dpi data. The expression pattern of one such cluster of genes correlated with weight gain and WUR genotype, contained numerous immune response genes such as cytokines, chemokines, interferon type I stimulated genes, apoptotic genes and genes regulating complement activation. In addition, Partial Correlation and Information Theory (PCIT) identified differentially hubbed (DH) genes between the phenotypically divergent groups. GO enrichment revealed that the target genes of these DH genes are enriched in adaptive immune pathways.

Conclusion

There are molecular differences in blood RNA patterns between pigs with extreme phenotypes or with a different WUR genotype in early responses to PRRSv infection, though they can be quite subtle and more difficult to discover with conventional DE expression analyses. Co-expression analyses such as WGCNA and PCIT can be used to reveal network differences between such extreme response groups.

Electronic supplementary material

The online version of this article (doi:10.1186/s12864-015-1741-8) contains supplementary material, which is available to authorized users.

A comparison of alternative random regression and reaction norm models for whole genome predictions

Whole genome prediction (WGP) based on high density SNP marker panels is known to improve the accuracy of breeding value (BV) prediction in livestock. However, these accuracies can be compromised when genotype by environment interaction (G×E) exists but is not accounted for. Reaction norm (RN) and random regression (RR) models have proven to be useful in accounting for G×E in pre-WGP evaluations by modeling BV as linear or higher order functions of environmental or temporal covariates. We extend these RR/RN models based on several alternative specifications for SNP-specific intercepts and linear slopes on environmental covariates. One specification is based on bivariate normality (BVN) of SNP-specific intercepts and slopes, whereas 2 others, IW-BayesA and based on inverted Wishart (IW) extensions IW-BayesB, are, respectively, bivariate Student t extensions of currently popular models without (BayesA) or with (BayesB) variable selection. We highlight alternative specifications based on the square root free Cholesky decomposition (CD) of SNP-specific variance-covariance (VCV) matrices in an attempt to better differentially model environmentally sensitive from environmentally robust QTL. Two CD specifications were considered with (CD-BayesB) or without (CD-BayesA) any variable selection on intercept and slope effects. We compared each of the 5 models based on an RN simulation study. Six scenarios were considered based on differences in overall genetic correlations between SNP-specific intercept and slope effects as well as on heritabilities and numbers of environmentally robust versus sensitive QTL. In most scenarios, IW-BayesA had the greatest accuracy, whereas CD-BayesB exhibited the greatest accuracy in low complexity architectures (i.e., low number of QTL). In an RR application of a Duroc × Pietrain resource population at Michigan State University, 5,271 SNP markers and 928 F2 animals with known pedigree were analyzed for backfat thickness at wk 10, 13, 16, 19, and 22. SNP-based RR methods had a 2.5% greater (P < 0.0001) cross-validation accuracy for predicting phenotypes than the SNP-based conventional BayesA/BayesB and/or pedigree based RR BLUP; however, none of the proposed RR models had performances that were different from each other.

A comparison of alternative random regression and reaction norm models for whole genome predictions

Whole genome prediction (WGP) based on high density SNP marker panels is known to improve the accuracy of breeding value (BV) prediction in livestock. However, these accuracies can be compromised when genotype by environment interaction (G×E) exists but is not accounted for. Reaction norm (RN) and random regression (RR) models have proven to be useful in accounting for G×E in pre-WGP evaluations by modeling BV as linear or higher order functions of environmental or temporal covariates. We extend these RR/RN models based on several alternative specifications for SNP-specific intercepts and linear slopes on environmental covariates. One specification is based on bivariate normality (BVN) of SNP-specific intercepts and slopes, whereas 2 others, IW-BayesA and based on inverted Wishart (IW) extensions IW-BayesB, are, respectively, bivariate Student t extensions of currently popular models without (BayesA) or with (BayesB) variable selection. We highlight alternative specifications based on the square root free Cholesky decomposition (CD) of SNP-specific variance-covariance (VCV) matrices in an attempt to better differentially model environmentally sensitive from environmentally robust QTL. Two CD specifications were considered with (CD-BayesB) or without (CD-BayesA) any variable selection on intercept and slope effects. We compared each of the 5 models based on an RN simulation study. Six scenarios were considered based on differences in overall genetic correlations between SNP-specific intercept and slope effects as well as on heritabilities and numbers of environmentally robust versus sensitive QTL. In most scenarios, IW-BayesA had the greatest accuracy, whereas CD-BayesB exhibited the greatest accuracy in low complexity architectures (i.e., low number of QTL). In an RR application of a Duroc × Pietrain resource population at Michigan State University, 5,271 SNP markers and 928 F2 animals with known pedigree were analyzed for backfat thickness at wk 10, 13, 16, 19, and 22. SNP-based RR methods had a 2.5% greater (P < 0.0001) cross-validation accuracy for predicting phenotypes than the SNP-based conventional BayesA/BayesB and/or pedigree based RR BLUP; however, none of the proposed RR models had performances that were different from each other.

Evaluation of housekeeping genes for normalizing real-time quantitative PCR assays in pig skeletal muscle at multiple developmental stages

Quantitative PCR (qPCR) requires a constantly expressed housekeeping gene as an internal control, the expression of which is similar in different biological samples. In the present study, we evaluated the applicability and compared the consistency of the gene expression of 16 reference genes, i.e., 10 common and 6 candidate genes, through qPCR assays in pig skeletal muscles at multiple developmental stages. The stability of these 16 potential reference genes was investigated using the geNorm and NormFinder methods. Our results indicated that DRAP1 and RNF7 were the most appropriate combination to normalize gene expression in the Yorkshire samples, the combination of DRAP1 and WSB2 were appropriate in the Tongcheng samples, H3F3A and DRAP1 in prenatal periods, DRAP1 and RNF7 in postnatal periods, and the combination of DRAP1 and WSB2 was most suitable for accurate normalization in whole samples. This result provides a strong basis for similar studies in pigs that explore the molecular mechanisms underlying skeletal muscle development.

Rapid screening for phenotype-genotype associations by linear transformations of genomic evaluations

Background

Currently, association studies are analysed using statistical mixed models, with marker effects estimated by a linear transformation of genomic breeding values. The variances of marker effects are needed when performing the tests of association. However, approaches used to estimate the parameters rely on a prior variance or on a constant estimate of the additive variance. Alternatively, we propose a standardized test of association using the variance of each marker effect, which generally differ among each other. Random breeding values from a mixed model including fixed effects and a genomic covariance matrix are linearly transformed to estimate the marker effects.

Results

The standardized test was neither conservative nor liberal with respect to type I error rate (false-positives), compared to a similar test using Predictor Error Variance, a method that was too conservative. Furthermore, genomic predictions are solved efficiently by the procedure, and the p-values are virtually identical to those calculated from tests for one marker effect at a time. Moreover, the standardized test reduces computing time and memory requirements.

The following steps are used to locate genome segments displaying strong association. The marker with the highest − log(p-value) in each chromosome is selected, and the segment is expanded one Mb upstream and one Mb downstream of the marker. A genomic matrix is calculated using the information from those markers only, which is used as the variance-covariance of the segment effects in a model that also includes fixed effects and random genomic breeding values. The likelihood ratio is then calculated to test for the effect in every chromosome against a reduced model with fixed effects and genomic breeding values. In a case study with pigs, a significant segment from chromosome 6 explained 11% of total genetic variance.

Conclusions

The standardized test of marker effects using their own variance helps in detecting specific genomic regions involved in the additive variance, and in reducing false positives. Moreover, genome scanning of candidate segments can be used in meta-analyses of genome-wide association studies, as it enables the detection of specific genome regions that affect an economically relevant trait when using multiple populations.

Electronic supplementary material

The online version of this article (doi:10.1186/1471-2105-15-246) contains supplementary material, which is available to authorized users.

Genotype imputation accuracy in a F2 pig population using high density and low density SNP panels

Background

F₂ resource populations have been used extensively to map QTL segregating between pig breeds. A limitation associated with the use of these resource populations for fine mapping of QTL is the reduced number of founding individuals and recombinations of founding haplotypes occurring in the population. These limitations, however, become advantageous when attempting to impute unobserved genotypes using within family segregation information. A trade-off would be to re-type F₂ populations using high density SNP panels for founding individuals and low density panels (tagSNP) in F₂ individuals followed by imputation. Subsequently a combined meta-analysis of several populations would provide adequate power and resolution for QTL mapping, and could be achieved at relatively low cost. Such a strategy allows the wealth of phenotypic information that has previously been obtained on experimental resource populations to be further mined for QTL identification. In this study we used experimental and simulated high density genotypes (HD-60K) from an F₂ cross to estimate imputation accuracy under several genotyping scenarios.

Results

Selection of tagSNP using physical distance or linkage disequilibrium information produced similar imputation accuracies. In particular, tagSNP sets averaging 1 SNP every 2.1 Mb (1,200 SNP genome-wide) yielded imputation accuracies (IA) close to 0.97. If instead of using custom panels, the commercially available 9K chip is used in the F₂, IA reaches 0.99. In order to attain such high imputation accuracy the F₀ and F₁ generations should be genotyped at high density. Alternatively, when only the F₀ is genotyped at HD, while F₁ and F₂ are genotyped with a 9K panel, IA drops to 0.90.

Conclusions

Combining 60K and 9K panels with imputation in F₂ populations is an appealing strategy to re-genotype existing populations at a fraction of the cost.

Methods of tagSNP selection and other variables affecting imputation accuracy in swine

BACKGROUND

Genotype imputation is a cost efficient alternative to use of high density genotypes for implementing genomic selection. The objective of this study was to investigate variables affecting imputation accuracy from low density tagSNP (average distance between tagSNP from 100kb to 1Mb) sets in swine, selected using LD information, physical location, or accuracy for genotype imputation. We compared results of imputation accuracy based on several sets of low density tagSNP of varying densities and selected using three different methods. In addition, we assessed the effect of varying size and composition of the reference panel of haplotypes used for imputation.

RESULTS

TagSNP density of at least 1 tagSNP per 340kb (~7000 tagSNP) selected using pairwise LD information was necessary to achieve average imputation accuracy higher than 0.95. A commercial low density (9K) tagSNP set for swine was developed concurrent to this study and an average accuracy of imputation of 0.951 based on these tagSNP was estimated. Construction of a haplotype reference panel was most efficient when these haplotypes were obtained from randomly sampled individuals. Increasing the size of the original reference haplotype panel (128 haplotypes sampled from 32 sire/dam/offspring trios phased in a previous study) led to an overall increase in imputation accuracy (IA = 0.97 with 512 haplotypes), but was especially useful in increasing imputation accuracy of SNP with MAF below 0.1 and for SNP located in the chromosomal extremes (within 5% of chromosome end).

CONCLUSION

The new commercially available 9K tagSNP set can be used to obtain imputed genotypes with high accuracy, even when imputation is based on a comparably small panel of reference haplotypes (128 haplotypes). Average imputation accuracy can be further increased by adding haplotypes to the reference panel. In addition, our results show that randomly sampling individuals to genotype for the construction of a reference haplotype panel is more cost efficient than specifically sampling older animals or trios with no observed loss in imputation accuracy. We expect that the use of imputed genotypes in swine breeding will yield highly accurate predictions of GEBV, based on the observed accuracy and reported results in dairy cattle, where genomic evaluation of some individuals is based on genotypes imputed with the same accuracy as our Yorkshire population.

Characterizing differential individual response to porcine reproductive and respiratory syndrome virus infection through statistical and functional analysis of gene expression

We evaluated differences in gene expression in pigs from the Porcine Reproductive and Respiratory Syndrome (PRRS) Host Genetics Consortium initiative showing a range of responses to PRRS virus infection. Pigs were allocated into four phenotypic groups according to their serum viral level and weight gain. RNA obtained from blood at 0, 4, 7, 11, 14, 28, and 42 days post-infection (DPI) was hybridized to the 70-mer 20K Pigoligoarray. We used a blocked reference design for the microarray experiment. This allowed us to account for individual biological variation in gene expression, and to assess baseline effects before infection (0 DPI). Additionally, this design has the flexibility of incorporating future data for differential expression analysis. We focused on evaluating transcripts showing significant interaction of weight gain and serum viral level. We identified 491 significant comparisons [false discovery rate (FDR) = 10%] across all DPI and phenotypic groups. We corroborated the overall trend in direction and level of expression (measured as fold change) at 4 DPI using qPCR (r = 0.91, p ≤ 0.0007). At 4 and 7 DPI, network and functional analyses were performed to assess if immune related gene sets were enriched for genes differentially expressed (DE) across four phenotypic groups. We identified cell death function as being significantly associated (FDR ≤ 5%) with several networks enriched for DE transcripts. We found the genes interferon-alpha 1(IFNA1), major histocompatibility complex, class II, DQ alpha 1 (SLA-DQA1), and major histocompatibility complex, class II, DR alpha (SLA-DRA) to be DE (p ≤ 0.05) between phenotypic groups. Finally, we performed a power analysis to estimate sample size and sampling time-points for future experiments. We concluded the best scenario for investigation of early response to PRRSV infection consists of sampling at 0, 4, and 7 DPI using about 30 pigs per phenotypic group.

Evaluation of BTA1 and BTA5 QTL Regions for Growth and Carcass Traits in American and Korean Cattle

Previously identified QTL regions on BTA1 and BTA5 were investigated to validate the QTL regions and to identify candidate genes for growth and carcass traits in commercial cattle populations from the USA and Korea. Initially, a total of 8 polymorphic microsatellite (MS) markers in the BTA1 and 5 QTL regions were used for Chi-square tests to compare the frequencies of individual alleles between high and low phenotypic groups for the US (Michigan Cattleman’s Association/Michigan State University; MCA/MSU) cattle. For a subsequent study, 24 candidate genes containing missense mutations and located within the QTL regions based on bovine genome sequence data were analyzed for genotyping in the two commercial cattle populations. Re-sequencing analyses confirmed 18 public missense SNPs and identified 9 new SNPs. Seventeen of these SNPs were used for genotyping of the MCA/MSU cattle (n = 98) and Korean native cattle (n = 323). On BTA1, UPK1B, HRG, and MAGEF1 polymorphisms residing between BM1312 and BMS4048 were significantly associated with growth and carcass traits in one or both of the MCA/MSU and Korean populations. On BTA5, ABCD2, IL22 and SNRPF polymorphisms residing between BL4 and BR2936 were associated with marbling and backfat traits in one or both of the MCA/MSU and Korean cattle populations. These results suggested that BTA 1 and 5 QTL regions may be segregating in both Korean Hanwoo and USA commercial cattle populations and DNA markers tested in this study may contribute to the identification of positional candidate genes for marker-assisted selection programs.

Versican, matrix Gla protein, and death-associated protein expression affect muscle satellite cell proliferation and differentiation

Our previous transcriptional profiling study using a turkey skeletal muscle-specific oligonucleotide microarray revealed over 3,000 genes that were differentially expressed at 3 critical stages of muscle development: 18 d embryonic, 1 d posthatch, and 16 wk of age. The genes versican, matrix Gla protein (MGP), and death-associated protein (DAP) were selected to study for their potential effects on muscle satellite cell proliferation and differentiation, as their functions in other tissues are suggestive of possible key roles in the regulation of myogenesis and they are differentially expressed throughout muscle development in the turkey. Using small interfering RNA to knockdown the expression of these genes during proliferation and differentiation, each of the genes was found to differentially affect proliferation and differentiation. Versican and MGP predominantly affected proliferation with line effects, but later stages of differentiation were affected by the knockdown of versican and MGP. The underexpression of DAP inhibited myotube formation, which is a necessary stage in the development of muscle fibers. Without myotube development, muscle fiber formation will be inhibited or abolished. This is the first report that these genes with no previously documented functions with regard to muscle development play a critical role in muscle cell proliferation and differentiation.

Evaluation of QTL for carcass merit and meat quality traits in a US commercial Duroc population

Putative quantitative trait loci (QTL) regions on 5 chromosomes (SSC3, 6, 12, 15, and 18) were selected from our previous genome scans of a Duroc×Pietrain F(2) resource population for further evaluation in a US commercial Duroc population (n=331). A total of 81 gene-specific single nucleotide polymorphism (SNP) markers were genotyped and 33 markers were segregating. The MDH1 SNP on SSC3 was associated with 45-min and ultimate pH (pHu), and pH decline. PRKAG3 on SSC15 was associated with pHu. The HSPG2 SNP on SSC6 was associated with marbling score and days to 113kg. Markers for NUP88 and FKBP10 on SSC12 were associated with 45-min pH and L*, respectively. The SSC15 marker SF3B1 was associated with L* and LMA, and the SSC18 marker ARF5 was associated with pHu and color score. These results in a commercial Duroc population showed a general consistency with our previous genome scan.

Estimation of linkage disequilibrium in four US pig breeds

BACKGROUND

The success of marker assisted selection depends on the amount of linkage disequilibrium (LD) across the genome. To implement marker assisted selection in the swine breeding industry, information about extent and degree of LD is essential. The objective of this study is to estimate LD in four US breeds of pigs (Duroc, Hampshire, Landrace, and Yorkshire) and subsequently calculate persistence of phase among them using a 60 k SNP panel. In addition, we report LD when using only a fraction of the available markers, to estimate persistence of LD over distance.

RESULTS

Average r2 between adjacent SNP across all chromosomes was 0.36 for Landrace, 0.39 for Yorkshire, 0.44 for Hampshire and 0.46 for Duroc. For markers 1 Mb apart, r2 ranged from 0.15 for Landrace to 0.20 for Hampshire. Reducing the marker panel to 10% of its original density, average r2 ranged between 0.20 for Landrace to 0.25 for Duroc. We also estimated persistence of phase as a measure of prediction reliability of markers in one breed by those in another and found that markers less than 10 kb apart could be predicted with a maximal accuracy of 0.92 for Landrace with Yorkshire.

CONCLUSIONS

Our estimates of LD, although in good agreement with previous reports, are more comprehensive and based on a larger panel of markers. Our estimates also confirmed earlier findings reporting higher LD in pigs than in American Holstein cattle, especially at increasing marker distances (> 1 Mb). High average LD (r2 > 0.4) between adjacent SNP found in this study is an important precursor for the implementation of marker assisted selection within a livestock species.Results of this study are relevant to the US purebred pig industry and critical for the design of programs of whole genome marker assisted evaluation and selection. In addition, results indicate that a more cost efficient implementation of marker assisted selection using low density panels with genotype imputation, would be feasible for these breeds.