AFEchidna示例9--分析阈性状

阈性状 (threshold trait )：性状数值达到某一特定值时表现为正常，达不到则为不正常，如血压，血糖含量、生物的抗病力等，在数据方面以0、1表示，属于二元数据分布。AFEchidna可以轻松应付阈性状，通过family参数选择binomial函数来分析。

简单示例：

res<-echidna(es0.file="dfm2.es0",
             fixed=lt~1+Rep,
             family=esr_binomial(link='marginal'),
             random=~nrm(TreeID),
             residual=~units,
             predict=c('Rep')
            )

运行过程如下：

> res<-echidna(es0.file="dfm2.es0",
+                fixed=lt~1+Rep,
+                family=esr_binomial(link='marginal'),
+                random=~nrm(TreeID),
+                residual=~units,
+                trace=T,
+                predict=c('Rep')
+              )

Running Echidna for analysis:  lt 

 Tue Mar 16 16:52:15 2021 

    1 LogL= -673.01  0.9767           554 DF 
    2 LogL= -672.07  0.8883           554 DF 
    3 LogL= -672.02  0.8606           554 DF 
    4 LogL= -672.02  0.8591           554 DF 
    5 LogL= -672.02  0.8591           554 DF
 Tue Mar 16 16:52:15 2021 LogL Converged 

运行过程如下：

> Var(res)
         Term   Sigma      SE  Z.ratio
1    Residual 0.85911 0.10752 7.990234
2 nrm(TreeID) 0.57500 0.40652 1.414445
> pin(res,mulp=c(h2~V2/(V1+V2),Vp~V1+V2))
variance components are as following:
         Term   Sigma      SE vcS
1    Residual 0.85911 0.10752  V1
2 nrm(TreeID) 0.57500 0.40652  V2

pin formula: 
h2 ~ V2/(V1 + V2)
Vp ~ V1 + V2

   Term Estimate    SE
1    h2    0.401 0.196
11   Vp    1.434 0.322
> prd<-predict(res)
> names(prd)
[1] "heads" "pred"  "ased" 
> prd$pred
$pred1
  Prediction Stnd_Error Ecode Rep
1   -0.08508    0.20736     E   1
2   -0.13336    0.19901     E   2
3   -0.12722    0.18984     E   3
4    0.04244    0.18837     E   4
5   -0.29345    0.18827     E   5

当然还可以进行批量分析，以之前例子类似，不做演示了。