R语言在遗传统计学中的应用

遗传与疾病•人类的某些性状及部分疾病与人体的遗传因素密切相关•阐明遗传因素与人体疾病或健康状态的关系有非常重要的意义•遗传统计学在这其中起着至关重要的作用遗传与疾病疾病的易感基因研究研究的特点•收集的数据即包含一般的表型数据也包含基因型数据•数据分析时需要用到不同的遗传模型•需要一些遗传统计特有的分析方法：LD的计算，家系图的绘制等R在遗传统计中的应用•数据整理•获取位点的基本信息•Hardy-Weinberg平衡检验•连锁不平衡的计算•关联研究常用分析方法•家系图的绘制•……数据整理•R中的genetics包专门为基因型数据提供一个新的类—genotype•genotype函数是genetics包里最基本的函数，可以将以下四种形式的初始基因型数据转换成便于分析的带有genotype类的数据1.以一个字符分隔的向量g1-genotype(c('C-C','C-T','C-C','T-T','C-C',''),sep='-')2.可以按某一位置分隔的向量g2-genotype(c('DD','DI','DD','II',''),sep=1)3.两个分开的向量allele1-c('D','D','D','I','')allele2-c('D','I','D','I','')g3-genotype(allele1,allele2)4.数据框或矩阵中的两列data-data.frame(allele1=c('D','D','D','I',''),allele2=c('D','I','D','I',''))g4-genotype(data$allele1,data$allele2)或data1-cbind(allele1=c('D','D','D','I',''),allele2=c('D','I','D','I',''))g5-genotype(data1)获取位点的基本信息多态位点的基本信息包括：位点分型成功率(callrate)、等位基因频率、基因型频率、杂合度和多态信息含量(PIC)一个简单的例子：#载入popn数据data(popn,package=DGCgenetics)#获取A位点的基本信息summary(popn$A)Numberofsamplestyped:1489(96.9%)AlleleFrequency:(2alleles)CountProportion117860.6211920.4NA94NAGenotypeFrequency:CountProportion1/27040.472/22440.161/15410.36NA47NAHeterozygosity(Hu)=0.4802686Poly.Inf.Content=0.3648558Hardy-Weinberg定律•Hardy-Weinberg定律是由英国数学家哈迪(D.H.Hardy)和德国医生温伯格(W.Weinberg)于1908年分别独立发现的，也称遗传平衡定律~(geneticequilibriumlaw)•该定律可以简单描述为，遗传平衡群体的等位基因频率与基因型频率在世代间维持恒定•该定律的适用条件是：随机婚配，群体足够大，没有突变、选择、迁移和遗传漂变Hardy-Weinberg平衡检验•关联研究中Hardy-Weinberg平衡检验常被用来评价基因分型的质量。我们通常对病例和对照组分别进行Hardy-Weinberg平衡检验•如果某一位点在对照组中不符合Hardy-Weinberg平衡，我们通常会怀疑该位点的基因型鉴定的质量•如果该位点在对照组平衡而在病例组出现不平衡，则该位点可能和疾病有关Hardy-Weinberg平衡检验•genetics包里面提供两种不同的检验方法•一种是Pearson‘schi-squaretest，可以用HWE.chisq函数进行该检验，另一种是Fisherexacttest，对应于HWE.exact函数•HWE.chisq常用于MAF较高、样本量较大的场合；MAF较低的位点建议使用HWE.exact函数LD的计算•连锁不平衡则是指人群中两个位点处在同一个单体型的频率比期望值高•评价连锁不平衡程度的指标包括D'、r2等•genetics包提供计算LD各种指标的函数，并能以文字和图形两种形式显示位点间的连锁不平衡程度LD的计算#用LD函数计算位点间的LDldresult-LD(popn)#用文字显示D'值summary(ldresult,which=D')#用图形显示结果LDtable(ldresult,which=D')PairwiseLD-----------BCDAD'0.9790.9760.976BD'0.9980.991CD'0.997关联研究常用分析方法•卡方检验•Logistic回归•线性回归•……卡方检验data(popn,package=DGCgenetics)#首先载入popn数据(geno-table(popn$A,popn$affected))ControlCase1/12612801/24062982/216183chisq.test(geno)Pearson'sChi-squaredtestdata:genoX-squared=23.7385,df=2,p-value=7.003e-06(alle-allele.table(popn$A,popn$affected))popn$affectedpopn$AControlCase19288582728464chisq.test(alle)Pearson'sChi-squaredtestwithYates'continuitycorrectiondata:alleX-squared=23.6881,df=1,p-value=1.133e-06Logistic回归1.共显性模型summary(glm(affected~A+sex,family=binomial,data=popn))Call:glm(formula=affected~A+sex,family=binomial,data=popn)DevianceResiduals:Min1QMedian3QMax-1.4081-1.2428-0.65151.11341.8190Coefficients:EstimateStd.ErrorzvaluePr(|z|)(Intercept)-0.65890.1354-4.8681.13e-06***A1/2-0.37520.1234-3.0410.00236**A2/2-0.78320.1695-4.6203.84e-06***sexFemale1.18660.13358.8902e-16***Logistic回归2.加性模型summary(glm(affected~allele.count(A,'2')+sex,family=binomial,data=popn))Call:glm(formula=affected~allele.count(A,2)+sex,family=binomial,data=popn)DevianceResiduals:Min1QMedian3QMax-1.410-1.239-0.6551.1171.814Coefficients:EstimateStd.ErrorzvaluePr(|z|)(Intercept)-0.653960.13054-5.0105.45e-07***allele.count(A,2)-0.388170.08107-4.7881.69e-06***sexFemale1.186760.133498.8902e-16***Logistic回归3.显性或隐性模型summary(glm(affected~carrier(A,'2')+sex,family=binomial,data=popn))Call:glm(formula=affected~carrier(A,2)+sex,family=binomial,data=popn)DevianceResiduals:Min1QMedian3QMax-1.4078-1.1979-0.74651.15711.6817Coefficients:EstimateStd.ErrorzvaluePr(|z|)(Intercept)-0.65660.1352-4.8571.19e-06***carrier(A,2)TRUE-0.47880.1164-4.1153.87e-05***sexFemale1.18350.13328.8842e-16***家系图的绘制library(kinship)#载入kinship包p1-scan(nlines=6,what=list(0,0,0,0,0,0))1100101102111101001211020021110310110212110410110211110510110211ped-as.data.frame(p1)names(ped)-c(famid,id,dadid,monid,sex,aff)par(xpd=TRUE)ped1=pedigree(ped$id,ped$dadid,ped$monid,ped$sex,ped$aff)plot(ped1,symbolsize=1.2)谢谢各位