SNP- and pedigree-based estimation of genetic parameters for back fat thickness traits in an F2 intercross between Landrace and Jeju native black pigs

Research Article
박 희복  Hee-Bok Park1이 재봉  Jae-Bong Lee2신 문철  Moon-Cheol Shin 1김 남영  Nam-Young Kim1손 준규  Jun-Kyu Son우 제훈  Jae-Hoon Woo1유 지현  Ji-Hyun Yoo1신 상민  Sang-Min Shin1박 남건  Nam-Geon Park1조 인철  In-Cheol Cho3양 병철  Byoung-Chul Yang1*

Abstract

Backfat thickness, an important quantitative economic trait, directly affects productivity in the swine industry. In this study, we estimated genetic parameters for four backfat thickness traits using information from both genome-wide single nucleotide polymorphism (SNP) chip and pedigree data. Four backfat thickness phenotypes were measured in 1,105 F2 progeny from an intercross between Landrace and Jeju native black pigs. All animals in the F2 intercross were subjected to genotypic analysis using PorcineSNP60K BeadChip platform and 39,992 SNP markers on autosomes filtered by quality control criteria were used to construct genomic relationship matrix for genetic parameter estimation. Restricted maximum likelihood (REML) and genomic restricted maximum likelihood (GREML) estimates of genetic parameters were obtained using both genomic- and pedigree- relationship matrix in linear mixed models. The heritability estimates showed a high level of concordance between the REML (0.48-0.57) and GREML (0.48-0.60) estimates. However, the standard errors of the heritability estimated using GREML (0.04-0.05) were smaller than those which were estimated using REML (0.10-0.11). The backfat thickness traits showed high positive phenotypic correlations with each other, ranging from 0.62 to 0.83. These traits showed higher positive genetic correlation coefficients (REML: 0.80-0.99, GREML: 0.93-0.99) compared to the phenotypic correlation coefficients. The standard errors of the genetic correlation coefficient estimated using GREML (0.003-0.02) were smaller than those which were estimated using REML (0.01-0.09). The results of this study showed that the backfat thickness traits had high heritability and genetic correlation. The accuracy of genetic parameters estimated by GREML was higher than those which were estimated using REML. However, additional studies may be necessary to verify the results obtained from this study.

Keyword



Introduction

양돈산업에서 등지방 두께는 주요 양적경제형질(quantitative economic trait)로서 정육생산량에 결정적인 영향을 미친다. 등지방 두께가 줄어 들수록 정육생산량이 늘어남으로써 수익률 증대로 이어 질 수 있기 때문이다(Cliplef and Mckay, 1993). 따라서, 양돈산업에서 등지방 두께형질의 유전적 능력 향상은 경제적 효율성과 생산성에 직결되어있다고 할 수 있다.

돈군에서 등지방 두께형질의 유전적 능력 향상을 위해서는 돈군의 유전적 특성을 나타내는 유전력 (heritability, h2) 이나 유전상관계수(genetic correlation coefficient, rg)와 같은 유전모수(genetic paramenter)의 정확한 추정이 필수적이다. 이제까지 경제동물의 유전적 개량을 위한 유전모수의 추정은 가계(pedigree) 정보를 이용하여, 그 가계 구성원 간의 상가적(additive)인 수리적 관계를 계산하여 구축한 혈연관계행렬(numeric relationship matrix, NRM)이 사용되어 왔다(Henderson, 1984). 최근 들어서 유전체학 기법의 급속한 발달로 단일염기서열(single nucleotide polymorphism, SNP) 마커를 고밀도로 집적한 chip 기반의 대용량 유전자형 결정법이 실용화 되었다. 이에 따라 집단(population) 구성원 간 유전적 관계를 SNP마커 정보를 이용하여 SNP chip 기반 유전체관계행렬(genomic relationship matrix, GRM)을 구축할 수 있는 방법들이 개발되었고, GRM은 유전모수 추정을 위한 혼합선형모형 분석에 있어서 NRM을 대치 할 수 있는 대안으로 대두되었다(VanRaden, 2008; Yang 등 2011; Goddard 등 2011).

국민소득증가와 이에 따른 생활 수준의 향상은 돈육소비 패턴의 변화를 가져와 양적인 소비에서 맛과 같은 질적인 소비 위주로 바뀌고 있으며, 이와 함께 양돈산업 및 일반소비자들에게 제주재래흑돼지(Jeju native black pig, JBP)에 대한 관심이 증대되고 있다. JBP는 제주도에서 길러져 오던 재래돼지품종으로써, 요크셔나 랜드레이스와 같은 외래돼지품종에 비하여 적색의 육색, 백색의 단단한 지방과 그리고 우수한 상강도(marbling degree)를 보이지만(Kim 등, 2009), 체구가 작고, 성장이 느리며, 산자수가 적은 단점이 있다. JBP는 경제성이 우수한 외래돼지품종의 제주도 도입으로 멸종의 위기를 맞았지만, 1980년대말 시작된 JBP 복원 및 보존사업으로 순종의 보호와 육성이 진행되어 왔으며, 2015년도에 천연기념물 제550호로 지정되었다. 하지만, 성장형질 및 번식형질에 관련된 단점들로 인해 순종 JBP의 산업화는 아직 미흡한 상황이다. 2006년 JBP를 소재로 이용한 본격적인 산업화를 위해, 단점인 성장형질을 개량하고, 장점인 육질형질을 보존하고자 하는 연구들이 시작되었다. JBP 품종과 랜드레이스 품종 간 F2 교배축군이 조성되어 육질형질에 영향을 미치는 양적 형질 좌위(quantitative trait loci, QTLs)를 전장유전체검색 방법을 적용하여 검출해 낼 수가 있었다(Cho 등, 2015; Yoo 등, 2012). 본 연구에서는 우리나라 재래돼지품종인 JBP와 외래돼지품종인 랜드레이스 품종 간 교배축군(F2세대)에서 얻어진 등지방 두께형질에 대한 유전모수를 SNP chip 및 가계 기반으로 추정하고 이를 비교 분석하였다.

Materials and methods

공시 동물

본 연구에 이용된 F2 집단은 농촌진흥청 국립축산과학원 난지축산연구소에서 JBP와 랜드레이스의 상호교차교배(reciprocal intercross)에 의해서 생산되고 육성되었다(Kang 등, 2016). P세대 돼지들은 랜드레이스가 17두, JBP가 19두가 F1세대 생산에 이용되었다. 상호교차교배를 통하여 생산된 F1세대 91두를 다시 상호교차교배를 통하여 약 1200두의 F2세대 돼지를 생산하였다. 생산된 F2세대 자손 중 총 1,105두(암 537, 수 568)를 유전자형 결정과 표현형 수집에 이용하였다.

분석 자료

모든 F2 세대 돼지들은 평균일령 199일에 동일한 도축장에서 도축되었으며 도축 전에 돼지들은 24시간 공복 시켰으나, 물은 자유급여 하였다. 등지방두께는 좌측 반도체의 세 부위에서 직접 측정하였다. 첫번째 등지방두께는 네번째와 다섯번째 갈비뼈 사이에서 측정하였고(bf_tho4_5, mm), 두번째 등지방 두께는 11번째와 12번째 갈비뼈 사이에서 측정하였으며(bf_tho11_12, mm), 마지막 갈비뼈와 첫번째 요추뼈(first lumbar vertebrae) 사이에서도 측정하였다(bf_tho_lum, mm). 이에 더하여 축산물품질평가원에서 측정한 등지방두께(BF, mm) [(11번째와 12번째 갈비뼈 사이 등지방 두께)+(마지막 갈비뼈와 첫번째 요추뼈 사이 등지방 두께)]/2 도 확보하여 본 연구에 이용하였다. SNP 마커의 유전자형 결정은 PorcineSNP60K BeadChip (Illumina Inc., USA)을 사용하여 수행하였다. 결정된 SNP 마커의 유전자형에 대하여 minor allele frequency(MAF)가 0.05 이하, 유전자형 결측 비율이 5% 이상, 하디-바인베르그평형(Hardy-Weinberg equilibrium)을 극도로 벋어 나는 경우(P<0.000001)를 제외한 돼지 18개 상염색체상의 SNP 마커 39,992개의 분석결과를 이용하였다. 이와 같은 작업은 PLINK 프로그램(Purcell 등, 2007)을 사용하였다.

통계 분석

측정된 등지방 두께 단일 형질들에 대한 유전력을 추정하기 선형혼합모형은 다음과 같다:

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/1.jpg (모형 1)

여기서,y는 등지방 두께 관측치들에 대한 벡터이고, b는 전체평균, 성별, 배치 그리고 도체중의 고정효과 벡터이며, X는 관측치에 각 고정 효과를 연관시키는 빈도행렬(incidence matrix)이며, u는 개체의 상가적 임의 유전효과 벡터이고, Z는 u에 연관이 된 빈도행렬이고, e는 임의잔차효과 벡터이다. e에 대한 평균과 분산은 e ~N(0, Ise2) 로 설정하였으며, u에 대한 평균과 분산은 u ~N(0, Asu2) 또는 u ~N(0, Gsu2)과 같이 설정하였다. 여기서 I는 단위행렬(identity matrix)이고, A는 혈연관계행렬(NRM)이며(Henderson, 1984), G는 유전체관계행렬(GRM)을 나타낸다. 앞서 얻은 39,992개의 SNP 마커를 이용해서 G를 구하기 위해 SNP 마커의 대립유전자의 빈도와 유전자형을 이용하는 아래와 같은 행렬식이 이용되었다(Yang 등, 2011):

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/2.jpg

여기서, xijxiki번째 존재하는 SNP마커에서 j번째와 k번째 개체의 소수대립유전자 개수(0, 1, 2)이며, 는 소수대립유전자의 빈도이다. GRM 구축에는 Yang 등(2011)이 개발한 GCTA(Genome-wide Complex Trait Analysis) 프로그램을 이용하였다. 등지방 두께 형질들간의 유전상관계수를 추정하기 위한 이변량선형혼합모형은 다음과 같다:

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/3.jpg

여기서, y1와 y2 는 분석에 고려된 두 표현형에 대한 벡터이고, b1와 b2는 전체평균, 성별, 배치 그리고 도체중의 고정효과 벡터이며, X1과 X2는 표현형관측치에 각 고정 효과를 연관시키는 빈도행렬(incidence matrix)이다. u1와 u2는 두 형질에 대한 개체의 상가적 임의 유전효과 벡터이고, Z1와 Z2는 u1와 u2에 연관된 빈도행렬이고, e1와 e2는 임의잔차효과 벡터이다. 이변량선형혼합모형의 기대치는 다음과 같다.

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/4.jpg

NRM(A)을 사용한 이변량선형혼합모형의 분산은 다음과 같으며,

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/5.jpg

GRM(G)을 사용한 이변량선형혼합모형의 분산은 다음과 같다.

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/6.jpg

여기서, σu112, σu122, σu212, σu222 는 고려된 두 형질들의 상가적 유전(공)분산이며, σe112, σe122, σe212, σe222 두 형질의 잔차(공)분산이다. I는 단위행렬(identity matrix)이다.

일변량에 대한 유전분산성분과 그리고 환경분산성분의 추정에는 NRM을 이용한 경우 restricted maximum likelihood (REML) 방법에 기반한 ASReml 프로그램(VSN international, UK)을 이용하였고, GRM을 이용한 경우 genomic restricted maximum likelihood (GREML)방법에 기반한 GCTA프로그램(Yang 등 2011)을 이용하였다. 모형 1을 이용해서 얻은 분산 성분 추정을 위한 결과를 이용해서는 아래의 식을 이용하여 유전력(h2) 및 SNP-유전력(SNP-h2)을 계산하였다:

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/7.jpg

이변량에 대한 유전(공)분산성분과 그리고 환경(공)분산성분의 추정에는 NRM을 이용한 경우 REML 방법에 기반한 ASReml 프로그램(VSN international, UK)을 이용하였고, GRM을 이용한 경우 GREML방법에 기반한 MTG2 프로그램(Lee and van der Werf, 2016)을 이용하였다.

표현형 상관계수 및 유전상관계수는 아래의 식을 이용하여 계산하였다:

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/8.jpg

각 모수들의 표준오차는 테일러전개식의 1차 근사(first order Talyor series approximation)를 이용한 델타방법(delta method)을 이용하여 구하였다.

Results and discussion

등지방 두께 측정형질의 기초 통계치들은 Yoo 등의 보고(2012)에 이미 제시되어 있으며, bf_tho4_5, bf_tho11_12, bf_tho_lum 및 BF의 평균은 각각 34.0 mm, 28.0 mm, 26.2 mm, 그리고 22.9 mm 이었다. Porcine SNP 60K BeadChip을 이용하여 유전자형이 분석된 총 62,163개의 SNP 표지들 중에서 quality control 작업을 거친 전체 39,992개 SNP의 상염색체 상의 위치에 분포에 대한 정보는 Park등의 보고(2016)에 제시되어 있다.

Figure 1 에서는 F2 교배축군의 돼지들의 SNP 정보를 이용하여 구성된 GRM을 heatmap으로 나타내었고 GRM의 대각선 원소(diagonal element) 및 비대각선 원소(off-diagonal element)들의 분포를 나타내는 히스토그램이 제시되어 있다. GRM heatmap을 살펴보면 F2 교배축군의 15 전형매(fullsib) 가계를 나타내는 블록들을 관찰 할 수 있다. 이 블록들을 구성하는 개체들간의 genomic relationship은 대략 0.25임을 알 수 있다 (Figure 1a). 대각선원소 분포의 평균은 기대치와 같이 1이었으나 표준편차가 0.085로 계산되어 평균을 중심으로 한 변이가 존재함을 알 수 있다 (Figure 1b). 비대각선원소 분포의 평균은 0에 근사하였으나 역시 상당한 변이가 존재하였다 (Figure 1c).

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/Figure_JABG_02_02_04_F1_0323.jpg

Figure 1 (a) Heatmap of the genomic relationship matrix G used in this study. (b) Histogram of offdiagonal elements of the genomic relationship matrix. (c) Histogram of diagonal elements of the genomic relationship matrix.

Table 1는 공시축군에서 bf_tho4_5, bf_tho11_12, bf_tho_lum 및 BF에 대한 유전 분산성분과 환경 분산성분, 그리고 유전력을 모형 1을 기반으로 REML과 GREML 방법을 이용하여 추정한 결과를 요약하였다. 모형 1을 이용해서 REML 방법으로 구한 유전력 추정치는 0.48에서 0.57이었며, GREML 방법으로 구한 유전력 추정치는 0.48에서 0.60으로 두 방법사이의 추정치가 거의 유사하였다. 하지만, 이들 추정치의 표준오차를 살펴보면 REML의 경우는 0.10에서 0.11이으나, GREML의 경우는 0.04에서 0.05으로서 REML에 의해 추정된 표준오차에 비해 거의 반 수준을 보였다(Table 1).

Table 1. Estimates of heritability (h2) using Model 1 http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/Table_JABG_02_02_04_T1.jpg

Figure 2는 공시축군에서 네 가지 등지방 두께 형질에 대해서 이 들의 표현형상관계수와 유전상관계수를 추정한 결과이다. 표현형 상관계수의 추정 범위는 0.62에서 0.83으로 고도의 양의 상관관계를 나타냈다. 모형 2에서 REML 방법으로 구한 유전상관계수의 추정치는 0.80에서 0.99로 추정되어 표현형 상관계수 보다 높은 양의 상관관계를 나타냈으나, GREML 방법으로 구한 유전상관계수추정치는 0.93에서 0.99으로 추정되어 더 높은 양의 상관관계를 나타냈다. 이들 추정치의 표준오차를 살펴보면 REML의 경우는 0.01에서 0.09이었으나, GREML의 경우는 0.003에서 0.02로서 REML에 의해 추정된 표준오차에 비해 낮은 수준을 보였다. 이와 같이 높은 수준에서 추정된 양의 유전상관계수는 비록 등지방 두께를 측정한 부위는 다르지만 유사한 유전요인들이 등지방 두께 변이에 영향을 주는 것으로 사료된다.

http://dam.zipot.com:8080/sites/jabg/images/N0270020204_image/Figure_JABG_02_02_04_F2.jpg

Figure 2 (a) Phenotypic correlation coefficients of the backfat thickness traits. (b) Phenotypic correlation coefficients (REML estimates) of the backfat thickness traits. (c) Phenotypic correlation coefficients (GREML estimates) of the backfat thickness traits. Numbers in parentheses represent standard error of each estimated coefficient.

동일한 통계 모형과 같은 공시동물에서 비롯된 표현형 자료가 사용되었음에도 불구하고, 이상과 같이 REML 방법과 GREML 방법의 유전모수 추정에 있어서 정확도에서 차이를 보이는 이유를 살펴본다면, 기존의 NRM을 이용하는 REML방법으로 분산성분을 추정하게 되면, Mendelian sampling variance이 이용하지 못하게 된다. 예를 들어 전형매(full-sib)의 경우 동일한 정보를 분산성분추정에 제공 할 수 밖에 없다. 반면에, GRM을 이용하는 GREML방법을 이용하게 된다면 Mendelian sampling variance를 이용 할 수 있게 되어 전형매(full-sib)의 경우라도 상이한 정보를 분산성분추정에 제공 할 수 있어 좀 더 정확한 분산성분추정치를 기대 할 수 있다(Hill and Weir, 2011). 하지만, 본 연구에서 사용된 F2 교배축군의 경우 유전체 구조가 이질적(heterogeneous)이므로 추정 된 유전모수의 해석에 조심을 기울일 필요가 있다.

본 연구에서는 제주재래흑돼지와 랜드레이스의 교배축군으로부터 등지방 두께형질의 유전모수를 추정을 위하여 기존의 가계정보를 이용한 방법과 SNP chip 분석에서 산출된 유전자형 정보를 활용한 방법을 비교하였고, 그 차이를 “Mendelian sampling variance” 를 이용 하여 고찰하였다. 본 연구의 결과는 제주재래돼지를 기반으로 계통집단 구축에 있어 등지방 두께형질 형질 개량 연구에 기초자료를 활용될 것으로 기대된다.

Acknowledgements

본 연구는 농촌진흥청 경상과제 “재래돼지 집단간 유전특성 및 번식형질 유전자 특성분석 연구(과제번호 PJ01198501)”의 지원과 2017년도 농촌진흥청 국립축산과학원 전문연구원 연수과정 지원사업에 의해 이루어진 것임.

References

1  1. Cliplef RL and Mckay RM.1993. Carcass quality characteristics of swine selected for reduced backfat thickness and increased growth rate. Can. J. Anim. Sci. 73: 483-494.  

2  2. Cho IC, Yoo CK, Lee JB, Jung EJ, Han SH, Lee SS, Ko MS, Lim HT and Park HB. 2015. Genome-wide QTL analysis of meat quality-related traits in a large F2 intercross between Landrace and Korean native pigs. Genet. Sel. Evol. 47:7. 

3  3. Goddard ME, Hayes BJ, Meuwissen TH. 2011. Using the genomic relationship matrix to predict the accuracy of genomic selection. J. Anim. Breed. Genet. 128:409-421 

4  4. Henderson C. 1984. Applications of linear models in animal breeding. University of Guelph, Guelph. Canada. 

5  5. Hill WG and Weir BS. 2011. Variation in actual relationship as a consequence of Mendelian sampling and linkage. Genet. Res. 93: 47-64 

6  6. Kang YJ, Cho SR, Jeong DK, Lee JB, Park HB, Cho IC, Han SH. 2016. Effect of Mating Types on the Growth Traits of F2 Population between Landrace and the Jeju Native Black Pigs. J. Emb. Trans. 31: 65-72. 

7  7. Kim DH, Seong PN, Cho SH, Kim JH, Lee JM, Jo C and Lim DG. 2009. Fatty acid composition and meat quality traits of organically reared Korean native black pigs. Livest. Sci. 120:96-102. 

8  8. Lee SH and van der Werf JHJ. 2016. MTG2: an efficient algorithm for multivariate linear mixed model analysis based on genomic information. Bioinformatics 32: 1420-1422. 

9  9. Park HB, Han SH, Lee JB, Kim SG, Kang YJ, Shin HS, Shin SM, Kim JH, Son JK, Baek KS, Cho SR, and Cho IC. 2016. SNP-based and pedigree-based estimation of heritability and maternal effect for body weight traits in an F2 intercross between Landrace and Jeju native black pigs. J Emb. Trans. 31: 243-247. 

10  10. Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira MA, Bender D, Maller J, Sklar P, de Bakker PI, Daly MJ and Sham PC .2007. PLINK: a tool set for whole-genome association and population-based linkage analyses. Am. J. Hum. Genet. 81:559-575. 

11  11. VanRaden PM. 2008. Efficient methods to compute genomic predictions. J. Dairy Sci. 91:4414-4423. 

12  12. Yang J, Lee SH, Goddard ME and Visscher PM. 2011. GCTA: a tool for genome-wide complex trait analysis. Am. J. Hum. Genet. 88:76-82. 

13  13. Yoo CK, Lim HT, Han SH, SS Lee, Ko MS, Kang T, Lee JH, Park HB, Cho IC. 2012. QTL analysis of backfat thickness and carcass pH in an F2 intercross between Landrace and Korean native pigs. Mol. Biol. Rep. 39: 8327-8333.