The effect of QTL-rich region polymorphisms identified by targeted DNA-seq on pig production traits

The aim of the present study was to analyse the effect of PLCD4, PECR, FN1 and PNKD mutations on pig productive traits and tested the usefulness of targeted enrichment DNA sequencing method as tool for preselection of genetic markers. The potential genetic markers for pig productive traits were identified by using targeted enrichment DNA sequencing of chromosome 15 region that is QTL-rich. The selected mutations were genotyped by using HRM, Sanger sequencing and PCR-ACRS methods. The association study was performed by using GLM model in SAS program and included over 500 pigs of 5 populations maintained in Poland. The variation (C/T) of PLCD4 gene affected feed conversion, intramuscular fat and water exudation. The PNKD mutations were associated with texture parameters measured after cooking. In turn, the variation rs792423408 (C/T) in the FN1 gene influenced toughness measured in semimembranosus muscle and growth traits that was observed particularly in Duroc breed. Summarizing, the investigated gene variants delivered valuable information that could be used during developing SNP microarray for genomic estimated breeding value procedure in pigs. Moreover, the study showed that the TEDNA-seq method could be used to preselect the molecular markers associated with pig traits. However, the further association study that included large number animal populations is necessary.


Introduction
One of the approaches based on next-generation sequencing (NGS) is targeted enrichment of genomic DNA (TEDNA-seq) sequencing, which allows the identification of a directly selected genomic region by using hybridization probes during DNA library preparation. Compared with a whole-genome strategy, the cost and effort of TEDNAseq analysis can be lower because of selective recovery and subsequent sequencing. Targeted enrichment sequencing is mainly applied to whole-exome analysis, which enables the identification of many diseases or metabolic pathways [1]. However, the bioinformatics of TEDNA-seq is a big challenge due to the complex analysis that includes choosing the most suitable software, algorithm, and filtering parameters to avoid false positive gene variant calls [2].
On the other hand, based on the current release of the Pig QTL database, there are 17 955 QTLs (quantitative trait loci)/associations identified from 576 publications that represent 635 different traits associated with meat and carcass traits, health, exterior, production and reproductive traits (http://www.anima lgeno me.org/cgi-bin/ QTLdb /SS/index ). In pigs, chromosome 15 contains a region between the microsatellites SW1683 and SW906 (79.3-89.3 cM, 127-135 Mbp) which is highly rich in QTLs associated with meat quality, fat content and growth traits. In this region, the abundance of QTLs related to meat quality traits, in particular meat texture, is the highest in comparison to other pig chromosomes. That SSC15 1 3 region contains more than 60 genes, including ABCA12, PRKAG3, MARCH4, DES, CYP27A1, OBSL1, IGFBP2 and IGFBP5, which have been already investigated for their relationship with meat quality and growth traits [3]. The PRKAG3 encode a muscle-specific isoform of the regulatory γ subunit of adenosine monophosphate-activated protein kinase (AMPK) [4]. A missense mutation in the PRKAG3 gene was detected that affected the pH and meat colour in pigs [5]. In turn, the DES gene encodes a cytoskeletal protein called desmin, which is located in the outer periphery of the Z-disk in skeletal muscle fibres. This protein plays an important role in joining the neighbouring myofibrils and adhering them to the cell membrane [6]. Starkey et al. [7] postulate that there is a strong correlation between desmin degradation and increase in tenderness and water exudation in ovine meat products, but the desmin degradation was not researched in relation to pork. And the study of Piórkowska et al. [8] showed that the ABCA12 promoter insertion affected daily gain, feed conversion and meat brightness. Although, that numerous genes located in the SSC15 region of interest has been analysed for their effect on pig traits, this region is still highly curious. Because, it probably contains unknown genetic markers associated with pork texture parameters and fat content.
Therefore, the aim of the present study was the analysis of PECR, FN1, PNKD and PLCD4 mutation effect on pig productive traits with particular emphasis on pork texture parameters. These mutations are located in the region on chromosome 15 and they were selected using TEDNA-seq method. The preliminary study was conducted by using two pig breeds Polish Landrace and Puławska that are significantly different in fat content, meat quality and growth traits. In turn, the main association study was performed by using over 500 pigs belonging to 5 pig populations maintained in Poland.

Animals
All animals used in the investigation originated from different farms, were female and unrelated. They were delivered to Pig Test Station of National Research Institute of Animal Production located in Chorzelów, Mełno, Pawłowice, and Rossocha. All pigs were maintained under the same environmental and diet condition according to Pig Test Station procedure. They were fed ad libitum by initial weight at 30 kg and finished at 100 ± 2 kg. The growth traits were measured during Test and the body composition, fat content and meat quality parameters were assessed after slaughter according to Tyra and Żak et al. [9] and meat texture parameters as was described by Ropka-Molik et al. [10].

DNA isolation
The blood samples were collected during slaughter and stabilized by EDTA. DNA was isolated from whole blood by Sherlock AX (A&A Biotechnology) kit.

Targeted enrichment DNA sequencing
The TEDNA-seq analysis, performed for another study [11], included 16 pigs of the Polish Landrace (n = 8) and Puławska (n = 8) breeds, which differ in growth, meat quality and fat content traits. The enrichment technique applied in this approach was RNA hybrid capture (1× tiling). The hybridization probes for the region of interest with exclusion of repeat-masked elements [12], which reduce the recovery of undesirable products, were designed by the Agilent team. Each library was indexed with a unique adaptor that enabled identification of most mutations for specific animals and then made it possible to perform an association study on chosen pigs. Filtered sequences were aligned to the Sus scrofa genome (Sscrofa10.2 assembly).

Selection of potential biomarkers
After performing TEDNA-seq on 16 chosen pigs, the analysis provided information about all gene variants in the SSC15 region, our region of interest, which is QTL-rich. The preliminary association analysis was conducted using LRT test in R-project. The LRT test was used to compare likelihood of linear null model versus alternative mode with SNP fixed effect.
The difference in log-likelihoods was compared with a Chi-squared distribution. The selection criteria of gene variants for further analysis were: P value ≤ 0.05 of LRT test obtained for at least few investigated pig productive traits and the relationship of identified mutations with genes likely involved in shaping the pig phenotype. The functional analysis of genes having significant mutations was carried out by Panther, STRING and KEGG. All selected gene variants were validated by Sanger sequencing.
Genotyping of selected mutations with the potential to become biomarkers Five mutations with high potential were selected. The genotyping was performed for 535 pigs of 5 breeds: Puławska, Polish Landrace, Polish Large White, Pietrain and Duroc by using different molecular techniques. The missense variant rs792423408 (C/T) of the FN1 gene and the upstream gene variant rs343851532 (C/T) of the PECR gene were genotyped by the high-resolution melting (HRM) technique. The KAPA™ HRM FAST PCR Kit (Kapa Biosystems, USA) was used to perform the HRM method on a QuantStudio 7 Flex Real-Time PCR System (Thermo Scientific, USA). The mutation rs324680963 (C/T), located upstream of the PLCD4 gene, was analysed by the PCR-ACRS method. PCR was performed using the AmpliTaq® 360 Gold Master Mix Kit (Thermo Scientific, USA) on a Mastercycler® Nexus Gradient (Eppendorf), and then the PCR product was digested in 37 °C with a mixture of 4U restriction enzyme RsaI (New England Biolabs, USA). The missense rs329501722 (C/A) and the frameshift rs792243103 (-/C) in the first exon of the PNKD gene variants were sequenced by Sanger sequencing using the Big-Dye® Terminator v3.1 Sequencing Kit (Thermo Scientific, USA) on a 3500XL Genetic Analyzer (Applied Biosystems, USA).

Statistical analysis
The gene variant effects on 535 gilts were estimated by using the GLM procedure (SASv.8.02). The final model was where, Y ijk is the observation, µ is the overall mean of the trait, d i is the fixed effect of the k th genotypes of gene, b j is the fixed effect of breed, s k is the fixed effect of pig station, e ijk is the random error.
To estimate significance differences between means, the Least Mean Squares (LSM) test was used. All results are shown as LSM ± SE.

Results
The selection of interesting gene variants after targeted enrichment sequencing of DNA By using LRT test and functional analysis five mutations of PECR, FN1, PNKD and PLCD4 genes ( Table 1) showing possible association with pig production traits were selected. The whole gene variant list with P values after LRT test is available by link goo.gl/dXenef.

Validation of significant gene variants by using the Sanger sequencing method
The validation of identified significant gene variants was performed by Sanger sequencing method using a 3500XL Genetic Analyzer (Applied Biosystems, USA) to verify their position on porcine chromosome 15. Chromatograms that confirmed the presence of variations in chosen genes were shown in Fig. 1.

PLCD4 gene
The PCR-ACRS method was used to genotype the mutation (rs324680963) in the gene PLCD4. In was observed that alternate T allele occurred with 13% frequency and homozygotes TT constituted only 4%. The highest number of pig individuals with the TT genotype has been observed in Polish Landrace (7%). In turn, in the Duroc, TT pigs were not present (Table S1).
The variation (C/T) of the PLCD4 gene significantly influenced feed conversion (kg/kg), intramuscular fat (IMF) and water exudation. The individuals with the CC genotype presented the highest feed conversion values in comparison to pigs with the other genotypes. Moreover, TT Polish Landrace pigs showed higher IMF (P < 0.05) than the other investigated pigs, which revealed similar values of IMF. In  (Table 2).

FN1 gene
To genotype the variation rs792423408 (C/T) in the FN1 gene, the HRM method was applied. It was observed that the majority of investigated pigs belonged to CC homozygotes, CT heterozygous genotype was identified in just 10% pigs and TT individuals were absent. And in the Pietrain, only one CC genotype was detected. The frequencies of the genotypes are shown in (Table S2). The influence of the FN1 gene on the toughness of cooked semimembranosus muscle was observed with a significance of P < 0.05 (Table 3). In all investigated breeds, the CT heterozygotes were characterized by higher values of toughness (+ 16 N/mm/s) than the CC homozygotes. The overall analysis showed a tendency that the CC pigs characterized by higher daily gain and achieved faster slaughter mass than CT pigs. Significant results were observed for the Duroc pigs, where CC pigs had 83 g higher daily gain (P < 0.05) and achieved a faster finish day (P < 0.01, 26 days earlier) than  In regards to meat quality traits, the CT Puławska pigs had higher values of both meat lightness and yellowness than the CC homozygotes (P < 0.05). Moreover, it was observed that rs792423408 (C/T) FN1 mutation probably affected meat texture parameters. Significant results were detected in Polish Landrace where heterozygotes showed 15% higher values of toughness than the CC pigs (P < 0.01).

PECR gene
The upstream PECR rs343851532 (C/T) variant was genotyped by the same as FN1 mutation the HRM method. The analysis of five Polish pig breeds showed low frequency of C allele (13%) in these populations. In all breeds only two genotypes TT and CT were observed (Table S3).
The conducted association study showed that PECR rs343851532 (C/T) mutation affected meat content. The CC pigs characterized by higher carcass yield and weight of loin than the CT individuals. It was particularly observed in the Duroc pigs, where CC Duroc showed 300 g higher weight of loin than those with the CT genotype. In turn, in the Pietrain breed CC pigs presented 2% higher meat percentage than the heterozygotes (Table 4). In regards to meat texture parameters just tendencies were identified.
The alternate A allele frequency of rs329501722 PNKD mutation was low and accounted only 13%. The 5% of investigated pigs revealed AA genotype and 24% belonged to heterozygotes (Table S4).
The association analysis showed that the PNKD rs329501722 (C/A) mutation influenced feed conversion, water exudation, resilience, cohesiveness for semimembranosus (ham) muscle and weight of loin. All investigated AA pigs showed higher water exudation (P < 0.01), resilience and cohesiveness values measured in cooked semimembranosus, and higher weight of loin (approx. 200 g, P < 0.05). It was particularly observed in PLW and PUL pigs, where AA pigs characterized by 5 and 10 units higher values of water exudation (P < 0.01) than the CC individuals. In the Duroc and Pietrain breeds the AA pigs were absent. Nevertheless, in the Duroc was observed that AC pigs showed lower values of daily gains (close to 100 g, P < 0.05) than the CC pigs. In the Polish Landrace, opposite tendencies were observed. The AA PL pigs exhibited the lowest growth performance parameters such as the 70 g lower daily gain (P < 0.05) and higher feed/gain ratio (P < 0.01), but their meat had higher IMF content (P < 0.05) than pigs with the other genotypes ( Table 5).
The second identified in PNKD gene rs792243103 (-/C) mutation was more frequent. In investigated pig populations prevailed heterozygotes (46%) and the C deletion in both alleles occurred with 17% frequency (Table S5). The present study showed that this INDEL mutation affected the growth traits particularly in PUL pigs, which belong to native breeds and were not under selective pressure. In the PUL breed, the absence of C in position 202 of the PNKD transcript caused a slight decrease in daily feed intake. Those pigs were also characterized by higher toughness values measure in the semimembranosus muscle (P < 0.05) and resilience parameters. Moreover, this mutation affected meat yellowness and carcass yield (P < 0.05). The pigs with C insertions in both alleles presented higher values of meat yellowness (P < 0.05) and carcass yield and lower toughness values (− 12 N/mm/s, P < 0.05) than pigs with other genotypes. The CC PL individuals showed the lowest carcass yield (− 3%, P < 0.05), weight of loin (− 300 g, P < 0.05) and primer cuts (− 700 g, P < 0.05). On the other hand, in the CC Pietrain pigs significantly lower pH value measured in longissimus dorsi muscle in comparison to -/-pigs (Table 6).

Discussion
The present study examined the targeted enrichment DNA sequencing method as tool for preselection of potential genetic markers. The positive results of the analysis could lead to the discovery of new possible ways to use this novel, fast and cost-effective technique. By using this method and functional analysis, five potential genetic markers were selected and included into downstream association study to show their effects on porcine productive traits. In position − 123 bp of the pig PLCD4 gene, one mutation (rs324680963) was identified that can introduce changes to the PLCD4 promoter region or in the sequence that enhances or silences gene expression. The analysis using the PROMO3 freeware showed that this mutation can affect the binding of numerous transcription factors (TFs). The PLCD4 gene encodes phospholipase C delta 4, which activates Ca 2+ ions during the leptin signalling pathway and induces the activity of NPY neurons and possibly affects feed intake. Additionally, phospholipase C participates indirectly in the thyroid hormone signalling pathway, which controls the metabolism of carbohydrates, proteins and lipids, thus having a potential influence on the growth traits [13]. In the present study, the upstream gene variant in the PLCD4 gene showed a significant impact on the feed conversion, IMF, water exudation and carcass yield. The CC homozygotes presented lower feed conversion than the pigs with other genotypes. However, the influence on feed intake was not observed. These results suggest that the identified mutation has not altered PLCD4 expression or that the phospholipase C delta 4 does not play such an important role in leptin signalling pathways, which was reflected in the phenotype. Nevertheless, it was shown that there is some significance in relation to feed conversion that is connected to the metabolism rate, thus implicating that this gene could play a meaningful role in the regulation of metabolism via the thyroid signalling pathway. Such a hypothesis could be explained by the research carried out by Onteru et al. [14] that showed the relationship between lipid metabolism and residual feed intake in Yorkshire pigs. In turn, Piórkowska et al. [15] proposed the PLCD4 gene as a candidate gene for meat quality traits such as IMF and fat content because of the probable connection with fatty acid metabolism. Our research also showed a high impact of the PLCD4 gene on IMF and water exudation, which should be further developed in studies included other pig breeds.
Another interesting gene variant that was located − 23 bp upstream of the PECR gene seemed to be highly promising and was tested in the present study as a potential biomarker. Because the analysis of predicting the transcription factor binding sites by PROMO3 determined that the replacement of cytosine to thymine causes the loss of binding sites for 14 transcription factors and generates new binding sites for 5 novel TFs, this variant should significantly affect the PECR gene expression. The PECR encodes peroxisomal trans-2-enoyl-CoA reductase that is involved in fatty acid chain elongation and biosynthesis of unsaturated fatty acids. This reductase has high affinity to fatty acids with chain lengths up to 16 carbons [16]. A recent study by Huttlin et al. [17] postulated that PECR interacts with actin alpha 2 (ACTA2), which is involved in actin cytoskeleton and calcium signalling pathways, among others. It is possible that PECR, by interacting with ACTA2, can indirectly influence the growth traits, or more likely, it can have an impact on meat quality traits such as IMF and fat content because of its main function. This hypothesis was also confirmed by Sadkowski et al. [18], who researched novel genes that can be correlated with higher IMF. They found that the PECR gene is one of the novel genes that can be implicated in beef marbling. Another study supporting the participation of the PECR gene in fat deposition was carried out by Kärst et al. [19], where they investigated the effect of the PECR gene on shaping the intramuscular fat level and water holding capacity values in mice with high muscle mass. In turn, Piórkowska et al. [20] analysed the PECR145T>C mutation effect on pork quality in Polish pig breeds and found thicker backfat, lower loin intramuscular fat content, lower texture parameter values, and higher PECR expression measured in the longissimus dorsi muscle in pigs with the TT genotype. Therefore, they recommended the PECR gene as one of the candidate genes influencing fat mass. In the present study, − 23 bp upstream of the PECR gene mutation affected only meat content traits. The effect on fat or IMF content were not observed. It suggests that this mutation does not regulate PECR expression. However, according to the latest pig genome release (Sscrofa11.1 assembly), this mutation is located also in 5` upstream regions of XRCC5 and FAM169 genes, therefore it could be associated with controlling their expression. Nevertheless, the abundance of publications with significant results focused on the function of the PECR gene in fat deposition in different animals suggests that the PECR polymorphisms could negatively influence pig productive traits by increasing the fat content. Therefore PECR mutations should be further investigated. Another gene variation investigated in the present study was a missense variant in the FN1 gene (rs792423408). This type of mutation is highly interesting because it leads to an alteration in the amino acid sequence, which could influence the protein function. The SIFT tool, available in the Ensembl database, predicted that the identified missense mutation belongs to a deleterious gene variant located in the DNA sequence encoding important protein domains. The FN1 gene encodes fibronectin 1, which is involved in multiple pathways in the organism, the most interesting of which are the integrin pathway, extracellular matrix-receptor interaction, and regulation of actin cytoskeleton. The main function of fibronectin 1 is to bind cell surfaces and various compounds including collagen, fibrin, and actin as they participate in creating the extracellular matrix and basement membrane [21]. Fibronectins take part in cell adhesion and motility but also in opsonization, wound healing, and maintenance of cell shape [22,23]. In the present study, it was observed that the missense variant of the FN1 gene influenced the toughness values measured in cooked semimembranosus muscle (P < 0.05). In turn, Cassar-Malek et al. [24] investigated the effect of myostatin in double-muscled cattle, and its connection with FN1 indicated that the fibronectin 1 gene is down-regulated in double-muscled cattle. This result suggests that FN1 expression is switched off when massive muscle growth occurs. Moreover, Ponsuksili et al. [25], searching for the genes that can influence the water holding capacity in Duroc and Pietrain pigs, found that the FN1 gene has a negative impact on water holding capacity, which can indicate that this protein has an essential role in maintaining the tight structure of muscles. However, in the present study, similar FN1 gene effects were not observed.
The last gene variations that were investigated in the present study were missense and frameshift variants of the PNKD gene. However, according to the latest pig genome release (SusScrofa assembly 11.1) it was shown that both polymorphisms are located in intron region of PNKD gene. Thus, they are rather involved in the regulation of PNKD expression than they alter amino acid sequence. Nevertheless, the function of the PNKD gene is still unknown, but there are some assumptions that the protein encoded by this gene could be a possible hydrolase [23]. However, STRING v.10 predicted a connection between the PNKD and MAPK3, ENO3, MYL12B, ILK and UBC proteins. MAPK3 participates in the MAPK/ERK pathway, which is essential for muscle cell proliferation and regeneration [26]. ENO3 is a protein called enolase 3 that is involved in muscle development in animals [27]. ILK, or integrin-linked kinase, takes part in the mediation of protein-protein interactions, thus connecting the ECM with intracellular cytoskeleton and signalling proteins [28]. In turn, UBC is a protein called ubiquitin C that participates principally in protein degradation. As shown above, all proteins that potentially interact with PNKD are connected in different ways with muscle growth. In the present study, It was shown that rs329501722 variant affected water exudation, weight of loin, ham resilience and cohesiveness (P < 0.01), which was particularly apparent in Puławska pigs, a breed that was not under selective pressure. The pigs with the CC genotype showed lower water exudation and better texture parameters; thus, selection for this allele can affect the pork quality. In turn, the frameshift variant of PNKD seems to influence meat yellowness, ham resilience and meat content traits (P < 0.05). Moreover, in the Polish Landrace breed, the effect on weight of loin, primary cuts and carcass yield was observed. Unfortunately, no publication has been written about the PNKD gene and its effect on any pig production traits, which makes it difficult to compare the results. However, in this manuscript were shown the interesting dependency between PNKD mutations and pig productive traits, thus PNKD effects should be further studied.
To summarize, the targeted enrichment DNA sequencing method could be used to preselection of molecular markers, but downstream association studies including large number animal populations are necessary. On the other hand, the investigated gene variants provided valuable information that could be used for the development of SNP microarrays to improve genomic estimated breeding value procedures in pigs.

Compliance with ethical standards
Conflict of interest None of the authors has a financial or other relationship with other people or organizations that may inappropriately influence this work.
Ethical approval The research will be performed on biological material derived from pigs maintained and slaughtered in the Test Pig Station (National Research Institute of Animal Production). In the Test Station pigs are slaughtered, dissected and after carcass evaluation, meat is standard intended for consumption. Therefore, our research does not require the approval of Animal Experimentation committee.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creat iveco mmons .org/licen ses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.