奶牛乳腺炎关键SNPs位点rs20438858及2b-RAD基因分型和分析方法与流程

本发明涉及一种奶牛乳腺炎关键snps位点rs20438858及2b-rad基因分型和分析方法。
背景技术：
：限制性酶切位点关联dna测序(radseq)技术是利用限制性内切酶对基因组进行酶切，产生一定大小的dna片段，然后通过构建测序文库对酶切后产生的rad标记进行高通量测序。在过去的十年里，radseq被认为是最重要的科学突破之一，在全基因组中通过单一、简单且成本效益高的方法，一次能检测到成千上万个基因组内的单核苷酸多态性标记(singlenucleotidepolymorphism，snp)，从而推动基因组学的研究。与其它测序技术相比较，该技术具有通量高、准确性好、实验周期短、性价比高和不受有无参考基因组序列的限制等优点。目前已经成功应用于种群群体遗传结构和系统进化分析、动植物重要经济性状的数量性状位点(qtl)定位和辅助遗传育种、遗传图谱的构建及snp标记检测等研究领域。radseq技术流程包括：基因组dna的酶切(1种内切酶酶)，构建文库(适配体连接，片段大小的筛选，片段端部修饰，末端添加y型适配器，pcr扩增)，上机测序(主要是illuminagaii或hiseq测序平台)，生物信息学分析(常用分析软件：stacks，pyrad和uneak等)。其具体流程图如图1。现有技术的缺点：1、酶切片段的长短大小不一，需要筛选；2、酶切片段端部需要两次添加不同的接头；3、酶切片段需要添加特殊的a-尾部和“y”型接头；4、步骤比较繁琐，技术要求高并且耗时；5、每个样本测序费用较高。技术实现要素：为了克服上述缺陷，本发明提供一种核酸内切酶dna片段长短均一，免除后续筛选、不需要多次添加接头、步骤简单缩短测序时间；降低每个样本的测序成本的2b-rad基因分型和分析方法。本发明还提供一个奶牛乳腺炎关键snps位点，该关键snps位点rs20438858位于基因tnfrsf21内含子区，snps为g>a，涉及染色体ac_000180.1。筛选出前述的奶牛乳腺炎关键snps位点的2b-rad基因分型和分析方法，包括如下步骤：1)建库测序：酶切：≥200ng基因组dna采用iib型限制性内切酶进行酶切；加接头：酶切产物分别加入5组不同的接头，t4脱氧核苷酸连接酶连接；扩增；串联；混库；测序：质检合格的dna文库上机测序；2)生物信息学分析：(1)数据过滤：对cleanreads进行质控；(2)酶切序列提取：提取含有酶切识别位点的序列，用于后续分析；(3)数据比对：利用soap软件将酶切序列比对到构建好的参考序列上；(4)snp分型：根据比对结果，利用最大似然法(ml)进行分型；(5)分析：构建进化树、主成分分析、群体遗传结构分析或全基因组关联分析。利用soap软件将酶切序列比对到参考序列后利用最大似然法(ml)进行snp标记分型，分型工作完成后采用下述的1)-5)步骤对分型结果进一步过滤：1)剔除所有样品中低于80％个体可以分型的位点；2)剔除maf低于0.01的位点；3)剔除含有1种或4种碱基型的单核苷酸多态(snp)位点；4)剔除标签内多于1个snp的位点；5)剔除标签内低于2个基因型的位点。采用bayesa模型和logistic回归模型对奶牛临床乳腺炎表型性状进行全基因组关联分析(gwas)；在进行全基因组关联分析(gwas)之前，首先构建基于奶牛乳腺炎表型性状的线性回归模型方程，其中，yi表示第i个体的表型特征向量；m为总snps数；μ为总表型性状平均值的特征向量；αk是第k个snp的加性相关性效应向量；xik为第i个体的第k个snp的基因型；e是残差效应的矢量；k指snp位点的个数。bayesa模型假定snps效应符合先验正态分布，其“零均值”和“snps方差”(“零均值”和“snps方差”等同，仅文字描述不同)以σk2表示，其中，k＝1,2……，m，k指snp位点的个数；snps效应方差是相互独立的，每个方差的独立分布iid与逆的卡方先验正态分布相同：其中v是自由度的参数，s2是尺度参数，p表示每个方差的独立分布(iid)与逆的卡方先验正态分布，χ-2为“逆卡方”；每个snp效应的临界度的先验分布符合t-分布：其中n指“当概率为п时，snps为零效应，或符合正态分布且概率分布为(1-п),”，p(αk│v,s2)表示为每个snp效应的临界度的先验分布，αk表示第k个snp的加性相关性效应向量，αk的先验取决于每个snp的方差，而每个snp的方差都有一个逆的卡方；当概率为п时，snps为零效应，或符合正态分布且概率分布为(1-п),αk│п,其中，代表所有非零snps效应的共同方差，它按比例分配了符合卡方检验的先验分布：模型中未知的п值由其先验分布(在0和1之间被认为是均匀的)或п-一致(0，1)预测。va被指定为4，由加性方差计算：和其中，pk表示为第k个snps的等位基因频率；为给定标记的差异；通过snps对加性遗传方差进行解释或阐明；为卡方检验的先验分布；pk表示第k个snps的等位基因频率；k为总snps数。logistic回归分析模型：假设单核苷酸多态性对奶牛乳腺炎的临床表型性状有影响，建立逻辑(logistic)回归模型来预测奶牛临床乳腺炎发生的可能性，首先构建拟合的logistic回归方程，其中，其中pj是在条件xj下乳腺炎临床表现型的概率，(1-pj)是在条件xj下临床乳腺炎表型不发生的概率，j表示第j个snp位点，xij＝(x1j,x2j,x3j……xmj)为第i个个体在j位点的基因型(0,1和2)，βj是第j个snp的影响，m是样本数量，μ为总表型性状平均值的特征向量；在逻辑回归分析模型中，y＝(μ+σβixi)方程转化成另一种形式：其中y表示为第i个个体的乳腺炎表型，p代表临床乳腺炎表型概率；xi为第i个个体的基因型；βi是优势比or；p和可变量之间表达的方程通过方程变换：95％置信区间(ci)＝exp(βi±1.96se(βi))，p1表示的是病例组某个snp位点发生的概率，p0表示的是对照组对应位点发生的概率；se(βi)表示为：βi的标准误。本发明通过两种分析模型得到1个奶牛乳腺炎关键snps位点，如表1和2：表1bayesa分析模型结果表2逻辑回归分析模型结果相对于现有技术，本发明的有益效果为：相对于radseq，2b-rad测序技术具有以下几点优点：1、酶切片段长短均一，不需要后续筛选；2、酶切片段不需要添加“y”型接头；3、步骤简单；4、每个样本测序成本低；5、测序耗时短。本发明还构建两种全基因组关联分析模型(bayesa和logistics)；3、筛选到一个中国荷斯坦奶牛乳腺炎关键snps位点及对应基因(tnfrsf21)。附图说明图1为现有技术的radseq测序技术流程图；图2为本发明的2b-rad测序流程图；图3.pcr扩增片段直接测序序列与ncbi参考序列比对图，(a)和(b)为pcr扩增片段直接测序chromas图；(c)1为ncbi参考序列，a和b为直接测序序列；灰色方框为单核苷酸多态标记位点。具体实施方式下面结合具体实施例和附图对本发明作进一步说明。2b-rad是一种基于iib型限制性内切酶的、简化的rad基因分型方法，为研究种群基因组遗传学提供了一种强有力的技术和方法。本研究中我们以中国荷斯坦奶牛为研究对象，构建中国荷斯坦奶牛临床乳腺炎和正常健康对照组牛群，提取构建牛群奶牛的全基因组，利用bael核酸内切酶对所有奶牛样本全基因组dna进行酶切，获得标准的酶切片段，然后进行上机测序并分析，具体建库测序流程为(图2)：(1)酶切：≥200ng基因组dna采用iib型限制性内切酶进行酶切；(2)加接头：酶切产物分别加入5组不同的接头，t4脱氧核苷酸连接酶(t4dnaligase)连接；(3)扩增：聚合酶链式反应(pcr)扩增连接产物；(4)串联：根据5组接头信息，将五个标签按顺序串联；(5)混库(pooling)：连接产物添加条形码(barcode)序列，混库；(6)测序：质检合格的高质量文库上机测序。上述的建库测序流程参见serialsequencingofisolengthradtagsforcost-efficientgenome-wideprofilingofgeneticandepigeneticvariations，作者为shiwang等人，2016年10月6号在线公开。生物信息学分析：本发明以牛属(https://www.ncbi.nlm.nih.gov/genome/？term＝bos+taurus)基因组作为参考基因组，利用soap软件(version2.21)将测序数据比对到参考序列，利用最大似然法(ml)进行位点的分型。分析流程如下：(1)数据过滤：对cleanreads进行质控；(2)酶切序列(enzymereads)提取：提取含有酶切识别位点的序列(reads)，我们称之为enzymereads，用于后续分析；(3)数据比对：利用soap软件将enzymereads比对到构建好的参考序列上；(4)snp分型：根据比对结果，利用最大似然法(ml)进行分型；(5)分析内容：构建进化树、主成分分析、群体遗传结构分析、全基因组关联分析等。利用soap软件将enzymereads比对到参考序列后利用最大似然法(ml)进行snp标记分型。过程中使用的rad分型软件包(radtyping)，包含10余个软件组分，覆盖了从数据预处理至最终分型结果输出的全过程。为保证后续分析的准确性，分型工作完成后会通过以下指标对分型结果进一步过滤：1)剔除所有样品中低于80％个体可以分型的位点；2)剔除maf低于0.01的位点；3)剔除含有1种或4种碱基型的单核苷酸多态(snp)位点；4)剔除标签内多于1个snp的位点；5)剔除标签内低于2个基因型的位点；所有样品共得到snp标记10058个。统计学分析模型本研究采用bayesa模型和logistic回归模型对奶牛临床乳腺炎表型性状进行全基因组关联分析(gwas)。我们首先构建了基于奶牛乳腺炎表型性状的线性回归模型方程，其中，yi表示第i个体的表型特征向量；m为总snps数；μ为总表型性状平均值的特征向量；αk是第k个snp的加性相关性效应向量；xik为第i个体的第k个snp的基因型(0,1和2)；e是残差效应的矢量。bayesa模型假定snps效应符合先验正态分布，其“零均值”和“snps方差”以σk2表示，其中，k＝1,2……，m；snps效应方差是相互独立的，每个方差的独立分布(iid)与逆的卡方先验正态分布相同，其中v是自由度的参数；s2是尺度参数：每个snp效应的临界度的先验分布符合t-分布：αk的先验取决于每个snp的方差，而每个方差都有一个逆的卡方,。当概率为п时，snps为零效应，或符合正态分布且概率分布为(1-п),αk│п,其中，代表所有非零snps效应的共同方差，它按比例分配了符合卡方检验的先验分布：从先验分布预测模型中的未知п值(在0和1之间被认为是均匀的)或п-一致(0，1)预测。va被指定为4，由加性方差计算：和其中，pk表示为第k个snps的等位基因频率；为给定标记的差异；通过snps对加性遗传方差进行解释或阐明。逻辑回归分析模型，假设单核苷酸多态性对奶牛乳腺炎的临床表型性状有影响，我们建立了逻辑(logistic)回归模型来预测奶牛临床乳腺炎发生的可能性，并建立了一个拟合的logistic回归方程，其中，其中pj是在条件xj下乳腺炎临床表现型的概率，(1-pj)是临床乳腺炎表型不发生的概率；xij＝(x1j,x2j,x3j……xmj)为第i个个体在j位点的基因型(0,1和2)，例如，aa表示为0，tt表示为2，at表示为1；也可以是这样：cc表示为0，gg表示为2，cg表示为1；也可以aa表示为0，cc表示为2，ac表示为1…；βj是第j个snp的影响；m是样本数量，μ为总表型性状平均值的特征向量。在逻辑回归分析模型中，y＝(μ+σβixi)方程可以转化成另一种形式：其中y表示为第i个个体的乳腺炎表型，p代表临床乳腺炎表型概率；xi为第i个个体的基因型；βi是优势比(or)；p和可变量之间表达的方程可以通过方程变换：95％置信区间(ci)＝exp(βi±1.96se(βi))。本研究通过两种分析模型得到1个奶牛乳腺炎关键snps位点，如表1和2：表1bayesa分析模型结果表2逻辑回归分析模型结果注：*表示由卡方(<0.05)计算的p-值；**是逻辑回归模型的t-统计p值(<0.05)；chisq是卡方检验下的卡方值。stat是logistic回归模型下的t-统计系数。or：优势比。l95：95％置信区间的概率比95％的下限。u95：95％概率置信区间95％的上限。为验证snp标记与奶牛乳腺炎的相关性，采用病例对照研究的方法，对病例组和对照组的关键snp位点暴露率进行了比较分析。经统计学检验，如果两组间存在显着性差异，可以认为是与奶牛乳房炎性状相关snp位点。在比较中排除外界匹配因素的干扰，仅考虑了snps与乳腺炎的关联关系。我们采用匹配设计和案例控制不相等(case/control＝1/h)来确定验证样本的数量。or＝ad/bcn为验证群体中所需临床乳腺炎数量，n为验证群体奶牛总数量。p0为正常对照群体snp位点突变的暴露率，p1为临床乳腺炎群体中snp位点突变的暴露率，or为比值比(预期该snp位点的关联强度)，α为假设检验第i类错误的概率(期望达到的检验显著性水平)，β为假设检验第ii类错误的概率，(1-β)为期望达到的检验把握度，or95％ci为95％置信区间，χ2为关键snp位点卡方检验。a为临床乳腺炎群体中snp位点突变个体数量，b为正常对照群体中snp位点突变个体数量，c为临床乳腺炎群体中snp位点非突变个体数量，d为正常对照群体中snp位点非突变个体数量，见表3。rs20438858snp位点碱基临床乳腺炎正常对照合计a17(a)142(b)159g56(c)168(d)224合计73310383表3snp标记与奶牛乳腺炎的相关性验证自由度df＝1，or＝ad/bc＝0.359，or值＜1说明中国荷斯坦奶牛临床乳腺炎的危险度因rs20438858位点g>a而减少，即a与乳腺炎之间为“负”关联；卡方χ2＝12.34≥10.828，p＜0.001，结论为拒绝无效假设，即snp位点rs20438858差异有统计学显著性。本发明所述的实例是对本发明的说明而不能限制本发明，在与本发明相当的含义和范围内的任何改变和调整，都应认为是在本发明的范围内。当前第1页12