||
本人在ASReml-R程序包的使用中,编写了一个程序包AAfun,增加一些额外的函数,作为ASReml-R的辅助。主要功能如下:
功能一、asreml.pin,用于计算遗传力、遗传相关等参数,并得到遗传相关的显著性水平。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 | library(AAfun) ### function 1 asreml.pin #### ## pin(): 1. count standard error(se) for heritability(h2) and corr value; ## 2. output significent level for corr value ## work for asreml package data(PrSpa) df<-PrSpa library(asreml) # exmaple 1 for sigle trait model fm<-asreml(h5~1+Rep, random=~Fam, subset=Spacing=='3',data=df) summary(fm)$varcomp[,1:3] pin(fm, h2 ~4*V1/(V1+V2)) # asreml result, h2 formlua # exmaple 2 for us model for bi-trait fm2<-asreml(cbind(dj,h5)~ trait+trait:Rep, random=~ us(trait):Fam, rcov=~units:us(trait), subset=Spacing=='3',data=df,maxit=40) summary(fm2)$varcomp[,1:3] pin(fm2, h2_A ~ 4 * V1/(V1+V5)) # heritability for trait A pin(fm2, h2_B ~ 4 * V3/(V3+V7)) # heritability for trait B # genetic corr: asreml result, corr formlua, 1 for signif pin(fm2, gCORR ~ V2/sqrt(V1*V3),1) pin(fm2, pCORR ~ (V2+V6)/sqrt((V1+V5)*(V3+V7)),1) # phenotype corr # exmaple 3 for corr model for multi-trait fm3<-asreml(cbind(dj,h3,h5)~ trait+trait:Rep, random=~ corgh(trait):Fam, rcov=~units:us(trait), subset=Spacing=='3',data=df,maxit=40) summary(fm3)$varcomp[,1:3] pin(fm3,N=3) # asreml result, corr coef number pin(fm3) # corr coef default N is 1, for the first corr |
运行结果如下:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 | > pin(fm, h2 ~4*V1/(V1+V2)) # asreml result, h2 formlua Estimate SE h2 0.306 0.124 > pin(fm2, h2_A ~ 4 * V1/(V1+V5)) # heritability for trait A Estimate SE h2_A 0.455 0.145 > pin(fm2, gCORR ~ V2/sqrt(V1*V3),1) Estimate SE sig.level gCORR 0.442 0.257 * --------------- Sig.level: 0'***' 0.001 '**' 0.01 '*' 0.05 'Not signif' 1 > pin(fm3,N=3) Estimate SE sig.level trait:Fam!trait.h3:!trait.dj.cor 0.751 0.233 *** trait:Fam!trait.h5:!trait.dj.cor 0.448 0.257 * trait:Fam!trait.h5:!trait.h3.cor 0.798 0.123 *** --------------- Sig.level: 0'***' 0.001 '**' 0.01 '*' 0.05 'Not signif' 1 |
功能二、asreml.met,用于辅助ASReml在多地点试验(MET)的因子分析。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 | ### function 2 asreml.met #### ## met.plot(): plots MET data ## met.corr(): calculate var/cov/corr from asreml MET factor analytic results ## met.biplot(): biplots MET factor analytic results from asreml ## work for asreml package data(MET) ## plot MET data -- example 1 # variable order: genotype,yield,site,row,col MET2<-MET[,c(1,9,2,4:5)] met.plot(MET2) ## plot MET data -- example 2 MET3<-MET[,c(1,9,2,4:7)] # add variable order on MET2: Rep, Block met.plot(MET3,"My met trials") MET$yield <- 0.01* MET$yield met.asr<-asreml(yield~Loc, random=~ Genotype:fa(Loc,2), rcov=~ at(Loc):ar1(Col):ar1(Row), data=MET, maxiter=50) ## count var/cov/corr matrix, etc. met.corr(met.asr,MET$Loc,2) met.corr(met.asr,MET$Loc,kn=2) # asreml met-result, Site, Site cluster number ## biplot asreml-met results met.biplot(met.asr,6,36,2) met.biplot(met2.asr,6,36,3) met.biplot(met.asr,siteN=6,VarietyN=36,faN=2,dSco.u=1.8,dLam.u=1.0) # dLam.u -- least distance from center # dSco.u -- least Variety breeding value # if can not draw fig 3, try multiplying or being devided by 10 for aim trait data. |
运行结果如下:
1. 绘制MET数据的图形
2. 输出MET分析的相关矩阵
1 2 3 4 5 6 7 8 9 10 11 12 13 14 | > met.corr(met.asr,MET$Loc,2) Site cluster results: S1 S2 S3 S4 S5 S6 1 1 1 2 1 2 CovVarCorr matrix S1 S2 S3 S4 S5 S6 S1 5.38 0.816 0.875 0.520 0.979 0.130 S2 3.90 4.246 0.688 0.472 0.811 0.161 S3 4.97 3.469 5.988 0.042 0.759 -0.366 S4 2.84 2.288 0.243 5.535 0.682 0.914 S5 4.31 3.169 3.525 3.044 3.598 0.328 S6 0.32 0.353 -0.954 2.291 0.662 1.135 |
3. 绘制MET分析结果的biplot图
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 | > met.biplot(met.asr,6,36,2) Scores.u is: 4 Fa1 Fa2 Scores V29 -1.164 -5.99 6.11 V30 0.470 -4.49 4.51 V31 -0.739 -5.02 5.08 V32 0.541 -4.03 4.06 V33 -0.218 -4.96 4.97 distFC.u is: 1.2 FA1 FA2 distFC S3 -2.03 1.37 1.68 S4 -1.39 -1.90 1.65 S6 -0.23 -1.04 1.57 |
未完待续……
程序包AAfun的下载地址:
1. 网盘下载 http://yzhlin-asreml.ys168.com/
2. 本地下载:
AAfun_1.0.zip
3. AAfun的例子脚本
AAfun.R
Archiver|手机版|科学网 ( 京ICP备07017567号-12 )
GMT+8, 2024-9-27 07:43
Powered by ScienceNet.cn
Copyright © 2007- 中国科学报社