Identification of the putative regulatory regions and candidate genes associated with backfat thickness and intramuscular fat content traits in Xiang pigs

Abstract

Backfat thickness (BFT) and intramuscular fat (IMF) content reflect fat deposition and meat quality in pigs. Numerous genetic variations contribute to differences in BFT and IMF for different pig breeds. However, limited information is available on the genetic variants associated with IMF deposition in Xiang pigs. To identify the putative regulatory regions and candidate genes associated with BFT and IMF traits in Xiang pigs, we integrated the whole genomic DNA genotyping and RNA-sequencing data in Longissimus thoracis muscle from F2 individuals of Large White × Xiang pig cross and performed the expression quantitative trait loci (eQTL) analyses. We identified a total of 64,419 significant cis-variant-eGene, including 33,288 eSNPs, 7,910 eInDels, and 5,940 eSVs. Among them, 1,156 cis-eQTL loci overlapped with QTL for IMF and BFT, affecting the expressions of 541 genes, including 54 transcription factors. Enrichment analysis revealed that key candidate genes were significantly associated with lipid metabolism pathways, such as mTOR and PPAR signaling. Overlap analysis demonstrated that several eQTLs exhibited dual regulatory effects on gene expression. Notable eQTL hotspots included eSNP 12:52760182_T, eInDel 2:6807484-6,807,485, and three eSV sites, with the P2RY6 gene significantly associated with multiple variants. Our results indicated that many regions across the pig genome that affected gene expression level (eQTL) and overlapped with QTL regions associated with the fat store and lipid metabolism for IMF and BFT. These findings enhance our understanding of the genetic mechanisms underlying fat deposition in Xiang pigs and provide a theoretical framework for future functional studies and genomic selection.

Keywords: EQTL; Gene expression; Intramuscular fat deposition; Longissimus thoracis; Xiang pigs.

Fig. 1

Features of cis-eVariants associations. A: Venn diagram showing the overlapped genes and specific genes between three type eQTLs. B: Disaggregated statistical pie chart of gene consequence by 47138 cis-eVariants. C: Box plots of distances between eQTLs sites of eSNPs, eInDel, and eSVs using logarithmic coordinate (Log10). D: Box plots denoted the distances between eQTLs (eSNPs, eInDels, and eSVs) and the TSS. E: Number proportion of the associated eGenes with three types of eQTLs such as eSNPs, eInDels, and eSVs. The legend “eGene_counts=1” in gray color showed that the eQTL was associated with one eGene, while “eGene_counts > 1” in orange color indicated that the eQTL was associated with more than one eGene

Fig. 2

Associations between eGenes and beta values at specific loci. The chord diagram represented the eGenes and their corresponding beta values at loci eSNP at 12:52760182_T, eInDel at 2:6807484-6807485, eSVs at 2:6040837-6041134, 2:6682987-6683057 and 2:6744727-6745016

Fig. 3

Associations of eGene and cis-eQTLs for BFT and IMF A: The chromosomal distribution of the detected eQTLs associated with IMF and BFT in F2 individuals from the cross between Large White and Xiang pigs. B: The Venn diagram illustrated the overlapping eQTLs and eGenes among eSNPs, eInDels, eSVs, as well as positive and negative effects, along with the specific eQTLs and eGenes unique to each category. C: The mRNA expression of the CENPQ, LRP2, PAFAH1B2, and PIK3C3 genes were regulated by three types of eVariants (eSNPs, eInDels, and eSVs) at distinct genotypes (0, 1, and 2)

Fig. 4

Functional annotations of candidate eQTLs and eGenes for BFT and IMF. A: A Sankey diagram depicting the eGenes associated with the eSNPs at loci 7:30120931_T, 6:7364793_A, 4:564432_G, 4:41466_C, 4:34519_A, 4:235130_C, 2:899_A, 2:6191_C, and 14:40473715_C, under different beta values. B: The protein–protein interaction (PPI) among the 493 eGenes which mapped on the QTL regions for the BFT and IMF content traits. C: GO categories of the 541 genes associated with BFT and IMF. D: Top 30 of KEGG enrichment of the 541 genes associated with BFT and IMF

Fig. 5

The PPI network among the 57 candidate genes associated with the BFT and IMF content traits

Fig. 6

Validation of candidate cis-eQTLs and eGenes. A: The eSNP rs319855910 (C>G) on chromosome 11 was validated in 19 F2 individuals. This locus is associated with the eGene MIPEP. We employed an RFLP-PCR assay, in which the allele G introduced an EcoR II site. The enzyme cleaved PCR-products containing rs319855910 site into two fragments (159 bp and 77 bp). After detection by agarose gel electrophoresis, three distinct genotype patterns were observed: (1) a single 236 bp band (CC homozygous), (2) three bands (236 bp, 159 bp, and 77 bp) (CG heterozygous), and (3) two bands (159 bp and 77 bp) (GG homozygous). B: The 100 bp deletion type of eSV on chromosome 6 (6:70310493-70310594) was validated bty PCR method in 19 F2 individuals. This locus is associated with the eGene NMNAT1. PCR analysis revealed three distinct genotypes: (1) a single 284 bp band in wild-type individuals (WW), (2) both 284 bp and 184 bp bands in heterozygous deletions (WD), and (3) a single 184 bp band in homozygous deletions (DD). C: The PCR validation in 19 F2 individuals confirmed the 164 bp deletion of eSV on chromosome 1 (1:215901356-215901521), which is associated with the eGene IL33. The analysis revealed two distinct genotypes: (1) a single 434 bp band for the wild-type (WW) and (2) both 434 bp and 298 bp bands for the heterozygous deletion (WD). Multi-sample genotyping (Supplementary Fig. S1) and Sanger sequencing (Supplementary Fig. S2) were performed to validate these genotypes. D-F: The box plots show the allelic differential expression in the Longissimus thoracis muscle. The horizontal lines within the boxes denote the median normalized expression, the box boundaries represent the 1st and 3rd quartiles, and the whiskers extend to 1.5 times the interquartile range. G-I Histograms presented the relative levels of eGene expression under different genotypes measured using RT-qPCR. Statistical significance was determined by one-way ANOVA with Tukey’s multiple comparisons test *: P < 0.05, all data were indicated as mean ± SEM