AAfun4包的新函数

主要针对原有的asreml.batch()，增加了两个函数：repin和update，用于计算新的遗传参数，或者重新运行新的批量分析。最新版的asreml.batch()，基本上能涵盖常见所有的可用模型，比如家系模型，个体模型，空间模型，基因组选择模型，及其多变量模型。

repin用于计算批量分析结果的新参数，无需重新调用Asrem程序。

update用于重新运行批量分析，只需输入相应的新模型项，比如固定效应、随机效应和误差效应等，无需重复输入原有的其它参数。同时，运行的新结果，可直接传递给repin计算遗传参数。

1 repin

用法：

repin(object, mulp, asrV = 4)
     -- object：批量分析运行的结果；
     --mulp：新的参数计算公式；
     --asrV：asreml版本，默认为最新版4.1。

简单示范如下：

// example
pkgs <- c('asreml','AAfun4')
lapply(pkgs,library,character.only=TRUE)
data(PrSpa)
df1<-subset(PrSpa,Spacing==3)

mulp1<-c(h2 ~ 4 * V1/(V1+V2))
run1<-asreml.batch(data=df1,factorN=1:5,traitN=c(9:13),
             FMod=y~1+Rep+Plot,
             RMod=~Fam,
             asrV=4,all.result=T,
             #family=asr_gaussian(),
             mulp=mulp1)
repin(run1,mulp=c(h2 ~ 4 * V1/(V1+V2),
                             H2 ~  V1/(V1+V2/4)))

结果如下：

// result
> run1<-asreml.batch(data=df1,factorN=1:5,traitN=c(9:13),
+              FMod=y~1+Rep+Plot,
+              RMod=~Fam,
+              asrV=4,all.result=T,
+              #family=asr_gaussian(),
+              mulp=mulp1)
 
Program starts running batch analysis ------ 
ASReml-R batch analysis results:
Fixed Factors -- Rep, Plot 
Randomed Factors -- Fam, units 
Index formula -- 
	h2 ~ 4 * V1/(V1 + V2)
Variance order: V1-Fam, V2-units
  Trait       V1        V2    V1.se    V2.se     h2  h2.se Maxit Converge
1    h1  11.9614   43.0035   3.1371   2.7240 0.8705 0.1872     6     TRUE
2    h2  53.7874  623.9830  22.4836  39.4926 0.3174 0.1266     4     TRUE
3    h3 132.2286 1591.7955  56.4882 100.7596 0.3068 0.1253     4     TRUE
4    h4 239.4825 3618.3452 115.6495 229.4253 0.2483 0.1161     5     TRUE
5    h5 445.0040 5345.0410 188.3244 338.9637 0.3074 0.1243     4     TRUE

Total running times: 0.38 sec elapsed
	
> repin(run1,mulp=c(h2 ~ 4 * V1/(V1+V2),
+                              H2 ~  V1/(V1+V2/4)))

Re-count new parameters for run1 :
       h2     H2  h2.se  H2.se
h1 0.8705 0.5266 0.1872 0.0685
h2 0.3174 0.2564 0.1266 0.0826
h3 0.3068 0.2494 0.1253 0.0828
h4 0.2483 0.2093 0.1161 0.0825
h5 0.3074 0.2498 0.1243 0.0821

2 update

用法：

update(object, FMod = NULL, RMod = NULL, EMod = NULL,
                          maxit = 30, asrV = 4)
      -- object：批量分析运行的结果；
      --FMod：新的固定效应；
      --RMod：新的随机效应；
      --EMod：新的误差效应；
      --asrV：asreml版本，默认为最新版4.1。

同时新运行的结果，可直接通过repin计算遗传参数。简单示范如下：

// An highlighted block
> run1a<-update(run1,FMod=y~1+Rep)

Re-run new batch for  run1 :

Variance order: V1-Fam, V2-units

         Fam      units     Fam.se   units.se Maxit Converge
h1  12.04852   43.06196   3.150461   2.719396     6     TRUE
h2  54.21185  621.04578  22.486323  39.191850     4     TRUE
h3 132.58582 1588.97922  56.464010 100.290167     4     TRUE
h4 241.38466 3599.02926 115.528283 227.535259     5     TRUE
h5 441.97490 5333.22697 187.241317 337.202806     4     TRUE

FMod :	y ~ 1 + Rep
RMod :	~Fam
EMod :	~units
<environment: 0x000001d05a0fcb78>
Warning message:
In if (grepl("!cor", Nvar2)) { :
  the condition has length > 1 and only the first element will be use
> repin(run1a,mulp=c(h2 ~ 4 * V1/(V1+V2),
+                          H2 ~  V1/(V1+V2/4)))

Re-count new parameters for run1a :
       h2     H2  h2.se  H2.se
h1 0.8745 0.5281 0.1872 0.0683
h2 0.3211 0.2588 0.1269 0.0824
h3 0.3081 0.2502 0.1254 0.0827
h4 0.2514 0.2115 0.1164 0.0824
h5 0.3061 0.2490 0.1240 0.0820

3 coef2

coef2函数主要用于获取随机效应、固定效应等的标准误。解决原函数coef不太容易获得标准误。
用法：

coef2(object, pattern = NULL, asrV = 4, ...)
-- object：asreml模型运行结果；
--pattern：固定效应或随机效应名，或缺失值mv
--asrV：asreml版本

简单示例如下：

// example
data(oats)
oats.asr <- asreml(yield ~ Variety*Nitrogen, 
                   random = ~ Blocks/Wplots, data=oats)
coef(oats.asr, pattern="Blocks:Wplots")
coef2(oats.asr, pattern="Blocks:Wplots")

运行结果：

> coef(oats.asr, pattern="Blocks:Wplots")
                     effect
Blocks_1:Wplots_1 -3.851187
Blocks_1:Wplots_2 14.080632
Blocks_1:Wplots_3  2.351459
Blocks_2:Wplots_1 -7.116152
Blocks_2:Wplots_2 -1.001696
Blocks_2:Wplots_3  5.788877
Blocks_3:Wplots_1  4.299037
Blocks_3:Wplots_2  6.209805
Blocks_3:Wplots_3 -9.193921
Blocks_4:Wplots_1 -9.849414
Blocks_4:Wplots_2  1.115452
Blocks_4:Wplots_3  3.496562
Blocks_5:Wplots_1 -7.916764
Blocks_5:Wplots_2 -6.064789
Blocks_5:Wplots_3 10.749965
Blocks_6:Wplots_1 -1.316788
Blocks_6:Wplots_2 -5.638062
Blocks_6:Wplots_3  3.856983
> coef2(oats.asr, pattern="Blocks:Wplots")
$`Blocks:Wplots`
                   solution std error    z.ratio
Blocks_1:Wplots_1 -3.851187    7.7817 -0.4949030
Blocks_1:Wplots_2 14.080632    7.7817  1.8094544
Blocks_1:Wplots_3  2.351459    7.7817  0.3021780
Blocks_2:Wplots_1 -7.116152    7.7817 -0.9144727
Blocks_2:Wplots_2 -1.001696    7.7817 -0.1287246
Blocks_2:Wplots_3  5.788877    7.7817  0.7439091
Blocks_3:Wplots_1  4.299037    7.7817  0.5524548
Blocks_3:Wplots_2  6.209805    7.7817  0.7980011
Blocks_3:Wplots_3 -9.193921    7.7817 -1.1814797
Blocks_4:Wplots_1 -9.849414    7.7817 -1.2657149
Blocks_4:Wplots_2  1.115452    7.7817  0.1433429
Blocks_4:Wplots_3  3.496562    7.7817  0.4493314
Blocks_5:Wplots_1 -7.916764    7.7817 -1.0173566
Blocks_5:Wplots_2 -6.064789    7.7817 -0.7793655
Blocks_5:Wplots_3 10.749965    7.7817  1.3814417
Blocks_6:Wplots_1 -1.316788    7.7817 -0.1692159
Blocks_6:Wplots_2 -5.638062    7.7817 -0.7245283
Blocks_6:Wplots_3  3.856983    7.7817  0.4956479

注：新版本AAfun4不再免费。