Bayesian Estimating of variance components on visit lenght at the feeder in pigs
Results bayesian estimating of variance components with Stan program, on 66.788 records of visit length time at the feeder from 6 trials, including locations (2 per trial), median weight as fixed effects and Eartag (118) and followers (118) as random effects in the models. iter = 12000, chain= 3, burn-in = 2000, Thin=1.
1. Mixed model: Location + Median wt + Eartag
1.1.Convergence Diagnostic
Check autorrelation, effective sample size, traceplot
ggs_autocorrelation(ggs(outp1)%>%filter(Parameter%in%lt1))
ggs_autocorrelation(ggs(outp1)%>%filter(Parameter%in%lt2))
ggs_autocorrelation(ggs(outp1)%>%filter(Parameter%in%c("Median Weight",
"var_eartag","var_error","prp_var_eartag","prp_var_error")))
autocorr.diag(outp1)## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4
## Lag 0 1.0000000000 1.0000000000 1.0000000000 1.000000000
## Lag 1 -0.1471087846 -0.0732679942 -0.0903091327 -0.168101632
## Lag 5 -0.0013865537 0.0064696105 -0.0052886663 -0.003085494
## Lag 10 -0.0004056016 0.0002101426 -0.0064972652 0.006270659
## Lag 50 0.0095696977 -0.0115891740 -0.0006608415 -0.003299075
## Loc1_t_5 Loc1_t_6 Loc2_t_1 Loc2_t_2 Loc2_t_3
## Lag 0 1.0000000000 1.000000000 1.0000000000 1.000000000 1.000000000
## Lag 1 -0.1576735352 -0.167466817 -0.1249448175 -0.091310126 -0.110954905
## Lag 5 0.0004148363 -0.002615544 -0.0007691954 0.002380314 0.002054322
## Lag 10 -0.0041770850 -0.005311132 -0.0008523821 0.005635538 0.003587318
## Lag 50 -0.0068091912 -0.009599190 0.0033201149 0.005839662 -0.001348789
## Loc2_t_4 Loc2_t_5 Loc2_t_6 Median Weight var_eartag
## Lag 0 1.000000000 1.000000000 1.000000000 1.0000000000 1.0000000000
## Lag 1 -0.162096122 -0.197137328 -0.139462238 0.0005343875 0.0628949588
## Lag 5 -0.003185497 0.010059382 0.005479880 -0.0065520315 -0.0003994407
## Lag 10 -0.008138205 -0.005912223 -0.006605805 -0.0004538409 -0.0013517849
## Lag 50 -0.001796362 -0.012676014 -0.003613669 -0.0006002109 -0.0046814530
## var_error prp_var_eartag prp_var_error lp__
## Lag 0 1.000000000 1.0000000000 1.0000000000 1.000000000
## Lag 1 -0.010628532 0.0580406038 0.0580406038 0.489826650
## Lag 5 -0.005053202 -0.0007058289 -0.0007058289 0.027630730
## Lag 10 -0.003881897 -0.0016326874 -0.0016326874 -0.010112585
## Lag 50 0.002083661 -0.0033205231 -0.0033205231 -0.007768452effectiveSize(outp1)## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4 Loc1_t_5
## 43936.58 36760.72 36655.08 42869.75 39172.92
## Loc1_t_6 Loc2_t_1 Loc2_t_2 Loc2_t_3 Loc2_t_4
## 42521.17 38855.17 36748.43 36963.68 43947.09
## Loc2_t_5 Loc2_t_6 Median Weight var_eartag var_error
## 45641.55 40864.56 30584.68 24626.22 30460.43
## prp_var_eartag prp_var_error lp__
## 25174.92 25174.92 10127.43geweke.diag(outp1)## [[1]]
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4 Loc1_t_5
## -0.27922 -1.72199 -0.96672 -1.74382 0.40273
## Loc1_t_6 Loc2_t_1 Loc2_t_2 Loc2_t_3 Loc2_t_4
## 0.41265 0.32945 1.43042 -0.58987 0.13493
## Loc2_t_5 Loc2_t_6 Median Weight var_eartag var_error
## -0.31106 -0.62678 -0.32678 0.32195 -1.28358
## prp_var_eartag prp_var_error lp__
## 0.35320 -0.35320 0.08677
##
##
## [[2]]
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4 Loc1_t_5
## 1.5796 0.5326 -0.5708 0.3733 0.7997
## Loc1_t_6 Loc2_t_1 Loc2_t_2 Loc2_t_3 Loc2_t_4
## 0.7577 -0.0771 0.1930 1.4573 1.4353
## Loc2_t_5 Loc2_t_6 Median Weight var_eartag var_error
## -1.1317 -0.8819 -0.7057 0.5718 -1.3023
## prp_var_eartag prp_var_error lp__
## 0.5971 -0.5971 -0.5174
##
##
## [[3]]
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4 Loc1_t_5
## -1.29541 0.40753 0.08778 0.09720 -0.32688
## Loc1_t_6 Loc2_t_1 Loc2_t_2 Loc2_t_3 Loc2_t_4
## 1.05565 -0.28786 -0.50399 -0.55452 -1.78947
## Loc2_t_5 Loc2_t_6 Median Weight var_eartag var_error
## 0.72070 -0.81843 0.42240 -1.06572 1.52031
## prp_var_eartag prp_var_error lp__
## -1.18914 1.18914 -0.30042#gelman.diag(outp1)
traplot(outp1,col =c("red1","blue4","purple3"))
denplot(outp1)
1.2. Summary Posterior Distribution
summary(outp1)##
## Iterations = 2001:12000
## Thinning interval = 1
## Number of chains = 3
## Sample size per chain = 10000
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## Loc1_t_1 1.354e+01 0.912104 5.266e-03 4.614e-03
## Loc1_t_2 7.653e+00 0.930848 5.374e-03 5.233e-03
## Loc1_t_3 1.049e+01 0.929341 5.366e-03 5.130e-03
## Loc1_t_4 8.459e+00 1.071551 6.187e-03 5.403e-03
## Loc1_t_5 9.494e+00 1.082744 6.251e-03 5.820e-03
## Loc1_t_6 5.357e+00 1.004754 5.801e-03 4.997e-03
## Loc2_t_1 1.112e+01 0.890776 5.143e-03 4.754e-03
## Loc2_t_2 8.372e+00 0.923649 5.333e-03 5.149e-03
## Loc2_t_3 1.531e+01 0.922397 5.325e-03 5.045e-03
## Loc2_t_4 7.582e+00 1.078935 6.229e-03 5.494e-03
## Loc2_t_5 1.005e+01 1.085097 6.265e-03 5.314e-03
## Loc2_t_6 7.275e+00 1.005811 5.807e-03 5.219e-03
## Median Weight 6.087e-03 0.001561 9.011e-06 8.928e-06
## var_eartag 9.203e+00 1.307650 7.550e-03 8.340e-03
## var_error 4.376e+01 0.239816 1.385e-03 1.375e-03
## prp_var_eartag 1.733e-01 0.020185 1.165e-04 1.274e-04
## prp_var_error 8.267e-01 0.020185 1.165e-04 1.274e-04
## lp__ -1.598e+05 8.215213 4.743e-02 8.167e-02
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## Loc1_t_1 1.174e+01 1.293e+01 1.354e+01 1.415e+01 1.533e+01
## Loc1_t_2 5.815e+00 7.036e+00 7.652e+00 8.275e+00 9.475e+00
## Loc1_t_3 8.657e+00 9.883e+00 1.049e+01 1.111e+01 1.232e+01
## Loc1_t_4 6.347e+00 7.740e+00 8.459e+00 9.187e+00 1.057e+01
## Loc1_t_5 7.339e+00 8.776e+00 9.487e+00 1.022e+01 1.162e+01
## Loc1_t_6 3.377e+00 4.683e+00 5.356e+00 6.036e+00 7.315e+00
## Loc2_t_1 9.356e+00 1.051e+01 1.111e+01 1.172e+01 1.286e+01
## Loc2_t_2 6.544e+00 7.757e+00 8.366e+00 8.989e+00 1.020e+01
## Loc2_t_3 1.350e+01 1.470e+01 1.531e+01 1.593e+01 1.714e+01
## Loc2_t_4 5.472e+00 6.857e+00 7.583e+00 8.299e+00 9.710e+00
## Loc2_t_5 7.921e+00 9.322e+00 1.005e+01 1.078e+01 1.219e+01
## Loc2_t_6 5.314e+00 6.601e+00 7.268e+00 7.956e+00 9.233e+00
## Median Weight 3.014e-03 5.041e-03 6.096e-03 7.148e-03 9.143e-03
## var_eartag 6.984e+00 8.282e+00 9.082e+00 9.993e+00 1.211e+01
## var_error 4.330e+01 4.360e+01 4.376e+01 4.392e+01 4.424e+01
## prp_var_eartag 1.377e-01 1.591e-01 1.719e-01 1.859e-01 2.168e-01
## prp_var_error 7.832e-01 8.141e-01 8.281e-01 8.409e-01 8.623e-01
## lp__ -1.598e+05 -1.598e+05 -1.598e+05 -1.598e+05 -1.598e+05print(Et.model)## Inference for Stan model: Eartagtrp_model.
## 3 chains, each with iter=12000; warmup=2000; thin=1;
## post-warmup draws per chain=10000, total post-warmup draws=30000.
##
## mean se_mean sd 2.5% 25% 50%
## beta[1] 13.54 0.00 0.91 11.74 12.93 13.54
## beta[2] 7.65 0.01 0.93 5.81 7.04 7.65
## beta[3] 10.49 0.01 0.93 8.66 9.88 10.49
## beta[4] 8.46 0.01 1.07 6.35 7.74 8.46
## beta[5] 9.49 0.01 1.08 7.34 8.78 9.49
## beta[6] 5.36 0.00 1.00 3.38 4.68 5.36
## beta[7] 11.12 0.00 0.89 9.36 10.51 11.11
## beta[8] 8.37 0.01 0.92 6.54 7.76 8.37
## beta[9] 15.31 0.01 0.92 13.50 14.70 15.31
## beta[10] 7.58 0.01 1.08 5.47 6.86 7.58
## beta[11] 10.05 0.01 1.09 7.92 9.32 10.05
## beta[12] 7.28 0.01 1.01 5.31 6.60 7.27
## beta[13] 0.01 0.00 0.00 0.00 0.01 0.01
## var_eartag 9.20 0.01 1.31 6.98 8.28 9.08
## var_error 43.76 0.00 0.24 43.30 43.60 43.76
## prp_var_eartag 0.17 0.00 0.02 0.14 0.16 0.17
## prp_var_error 0.83 0.00 0.02 0.78 0.81 0.83
## lp__ -159767.12 0.08 8.22 -159783.94 -159772.56 -159766.73
## 75% 97.5% n_eff Rhat
## beta[1] 14.15 15.33 37302 1
## beta[2] 8.27 9.47 30940 1
## beta[3] 11.11 12.32 33375 1
## beta[4] 9.19 10.57 35610 1
## beta[5] 10.22 11.62 34484 1
## beta[6] 6.04 7.32 40489 1
## beta[7] 11.72 12.86 34982 1
## beta[8] 8.99 10.20 30562 1
## beta[9] 15.93 17.14 32164 1
## beta[10] 8.30 9.71 38432 1
## beta[11] 10.78 12.19 40682 1
## beta[12] 7.96 9.23 37120 1
## beta[13] 0.01 0.01 29969 1
## var_eartag 9.99 12.11 24537 1
## var_error 43.92 44.24 29738 1
## prp_var_eartag 0.19 0.22 24881 1
## prp_var_error 0.84 0.86 24881 1
## lp__ -159761.38 -159752.03 10326 1
##
## Samples were drawn using NUTS(diag_e) at Sun Sep 29 09:04:11 2019.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at
## convergence, Rhat=1).Posterior Correlation of model parameters
parcorplot(outp1,col = terrain.colors(15,0.5,T), cex.axis=0.75)
2. Mixed model: Location + Median wt + Eartag + Follower
## [1] "mcmc.list"2.1.Convergence Diagnostic
Check autorrelation, effective sample size, traceplot
ggs_autocorrelation(ggs(outp2)%>%filter(Parameter%in%lt1))
ggs_autocorrelation(ggs(outp2)%>%filter(Parameter%in%lt2))
ggs_autocorrelation(ggs(outp2)%>%filter(Parameter%in%c("Median Weight",
"var_eartag", "var_follower", "var_error",
"prp_var_eartag","prp_var_follower","prp_var_error")))
autocorr.diag(outp2)## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4 Loc1_t_5
## Lag 0 1.000000000 1.000000000 1.000000000 1.00000000 1.00000000
## Lag 1 0.594002608 0.618300364 0.639921465 0.58903279 0.59315420
## Lag 5 0.143892653 0.158340994 0.163086715 0.11255210 0.12902442
## Lag 10 0.032271950 0.041433346 0.036490151 0.01985852 0.01850784
## Lag 50 0.004215207 -0.002853048 0.006197219 -0.01276794 -0.01263318
## Loc1_t_6 Loc2_t_1 Loc2_t_2 Loc2_t_3 Loc2_t_4
## Lag 0 1.000000000 1.000000000 1.000000000 1.000000000 1.000000000
## Lag 1 0.639825327 0.606355419 0.623261740 0.617011771 0.598034082
## Lag 5 0.161308123 0.135673113 0.158599043 0.154108266 0.120278300
## Lag 10 0.035861957 0.025648976 0.051095396 0.036564891 0.024250692
## Lag 50 -0.008846938 -0.005996926 -0.003731802 0.001050777 -0.002535053
## Loc2_t_5 Loc2_t_6 Median Weight var_eartag var_follower
## Lag 0 1.000000000 1.000000000 1.000000000 1.000000000 1.000000000
## Lag 1 0.571115813 0.622484311 -0.325457757 0.007039815 0.032155487
## Lag 5 0.114004812 0.154890123 0.003143067 0.026226007 0.028048595
## Lag 10 0.023716365 0.033294392 0.005613138 -0.001739518 -0.006046875
## Lag 50 0.002476578 0.002289386 -0.001579031 0.002985799 -0.005553047
## var_error prp_var_eartag prp_var_follower prp_var_error
## Lag 0 1.0000000000 1.0000000000 1.000000000 1.0000000000
## Lag 1 -0.2722801360 0.0001802698 0.032527274 -0.0001569912
## Lag 5 -0.0020284288 0.0257939502 0.026363432 0.0270425119
## Lag 10 0.0004952582 -0.0013706333 -0.006734379 -0.0007329823
## Lag 50 -0.0011531812 0.0028527428 -0.005571189 0.0029827281
## lp__
## Lag 0 1.000000000
## Lag 1 0.505119433
## Lag 5 0.048943214
## Lag 10 0.011836083
## Lag 50 0.009907233effectiveSize(outp2)## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4
## 6401.532 5848.650 5967.566 7084.777
## Loc1_t_5 Loc1_t_6 Loc2_t_1 Loc2_t_2
## 6699.894 5764.725 6355.171 5953.901
## Loc2_t_3 Loc2_t_4 Loc2_t_5 Loc2_t_6
## 6107.235 6908.411 7292.276 6155.945
## Median Weight var_eartag var_follower var_error
## 58949.799 20831.665 19326.915 52492.818
## prp_var_eartag prp_var_follower prp_var_error lp__
## 21151.015 19473.568 21106.489 9864.223geweke.diag(outp2)## [[1]]
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4
## 0.5777 -0.4235 1.3525 -1.5509
## Loc1_t_5 Loc1_t_6 Loc2_t_1 Loc2_t_2
## 0.5631 -1.3873 0.3729 -0.7813
## Loc2_t_3 Loc2_t_4 Loc2_t_5 Loc2_t_6
## 1.1437 0.2879 -2.0518 0.9542
## Median Weight var_eartag var_follower var_error
## 0.7693 -0.5184 0.5412 2.1732
## prp_var_eartag prp_var_follower prp_var_error lp__
## -0.6573 0.6127 0.5788 0.2980
##
##
## [[2]]
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4
## 1.09723 1.76390 -1.35289 0.67903
## Loc1_t_5 Loc1_t_6 Loc2_t_1 Loc2_t_2
## 0.33461 1.24954 0.30048 0.84359
## Loc2_t_3 Loc2_t_4 Loc2_t_5 Loc2_t_6
## -0.08648 -1.02497 1.04759 -0.88728
## Median Weight var_eartag var_follower var_error
## -1.60636 1.72855 1.02087 -0.43271
## prp_var_eartag prp_var_follower prp_var_error lp__
## 1.73276 0.76178 -1.94014 -0.34260
##
##
## [[3]]
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4
## 0.57982 1.09614 0.19825 -0.23838
## Loc1_t_5 Loc1_t_6 Loc2_t_1 Loc2_t_2
## -0.66364 0.05591 1.28029 -0.21406
## Loc2_t_3 Loc2_t_4 Loc2_t_5 Loc2_t_6
## -0.13724 0.94658 0.71199 -0.35630
## Median Weight var_eartag var_follower var_error
## -0.47151 -0.98352 0.22910 0.45921
## prp_var_eartag prp_var_follower prp_var_error lp__
## -1.03208 0.41784 0.98461 1.32510gelman.diag(outp2)## Potential scale reduction factors:
##
## Point est. Upper C.I.
## Loc1_t_1 1 1.00
## Loc1_t_2 1 1.01
## Loc1_t_3 1 1.00
## Loc1_t_4 1 1.00
## Loc1_t_5 1 1.00
## Loc1_t_6 1 1.00
## Loc2_t_1 1 1.00
## Loc2_t_2 1 1.00
## Loc2_t_3 1 1.00
## Loc2_t_4 1 1.00
## Loc2_t_5 1 1.01
## Loc2_t_6 1 1.00
## Median Weight 1 1.00
## var_eartag 1 1.00
## var_follower 1 1.00
## var_error 1 1.00
## prp_var_eartag 1 1.00
## prp_var_follower 1 1.00
## prp_var_error 1 1.00
## lp__ 1 1.00
##
## Multivariate psrf
##
## 1.01traplot(outp2, col = c("red1","purple1","blue4"))
denplot(outp2,col = c("red1","purple1","blue4"))
## 2.2. Summary Posterior Distribution
summary(outp2)##
## Iterations = 2001:12000
## Thinning interval = 1
## Number of chains = 3
## Sample size per chain = 10000
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## Loc1_t_1 1.358e+01 0.934862 5.397e-03 1.188e-02
## Loc1_t_2 7.621e+00 0.959333 5.539e-03 1.287e-02
## Loc1_t_3 1.040e+01 0.961302 5.550e-03 1.290e-02
## Loc1_t_4 8.473e+00 1.128124 6.513e-03 1.393e-02
## Loc1_t_5 9.546e+00 1.125146 6.496e-03 1.414e-02
## Loc1_t_6 5.407e+00 1.080693 6.239e-03 1.462e-02
## Loc2_t_1 1.142e+01 0.922965 5.329e-03 1.175e-02
## Loc2_t_2 8.395e+00 0.964322 5.568e-03 1.307e-02
## Loc2_t_3 1.539e+01 0.953857 5.507e-03 1.254e-02
## Loc2_t_4 7.503e+00 1.120811 6.471e-03 1.385e-02
## Loc2_t_5 9.987e+00 1.125537 6.498e-03 1.347e-02
## Loc2_t_6 7.529e+00 1.068002 6.166e-03 1.424e-02
## Median Weight 6.026e-03 0.001554 8.971e-06 6.402e-06
## var_eartag 9.093e+00 1.295966 7.482e-03 9.069e-03
## var_follower 8.662e-01 0.136249 7.866e-04 9.802e-04
## var_error 4.300e+01 0.235474 1.360e-03 1.028e-03
## prp_var_eartag 1.712e-01 0.020064 1.158e-04 1.393e-04
## prp_var_follower 1.636e-02 0.002564 1.480e-05 1.837e-05
## prp_var_error 8.124e-01 0.019756 1.141e-04 1.374e-04
## lp__ -1.592e+05 11.483490 6.630e-02 1.157e-01
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## Loc1_t_1 1.172e+01 1.296e+01 1.358e+01 1.421e+01 1.538e+01
## Loc1_t_2 5.743e+00 6.985e+00 7.621e+00 8.251e+00 9.524e+00
## Loc1_t_3 8.526e+00 9.755e+00 1.040e+01 1.105e+01 1.230e+01
## Loc1_t_4 6.260e+00 7.720e+00 8.476e+00 9.215e+00 1.069e+01
## Loc1_t_5 7.329e+00 8.795e+00 9.543e+00 1.030e+01 1.176e+01
## Loc1_t_6 3.261e+00 4.682e+00 5.414e+00 6.132e+00 7.526e+00
## Loc2_t_1 9.607e+00 1.080e+01 1.143e+01 1.204e+01 1.321e+01
## Loc2_t_2 6.479e+00 7.745e+00 8.403e+00 9.046e+00 1.026e+01
## Loc2_t_3 1.352e+01 1.476e+01 1.539e+01 1.604e+01 1.727e+01
## Loc2_t_4 5.285e+00 6.760e+00 7.490e+00 8.262e+00 9.719e+00
## Loc2_t_5 7.765e+00 9.231e+00 9.998e+00 1.074e+01 1.218e+01
## Loc2_t_6 5.452e+00 6.816e+00 7.525e+00 8.241e+00 9.631e+00
## Median Weight 2.972e-03 4.982e-03 6.019e-03 7.073e-03 9.077e-03
## var_eartag 6.890e+00 8.176e+00 8.979e+00 9.884e+00 1.195e+01
## var_follower 6.337e-01 7.704e-01 8.536e-01 9.485e-01 1.167e+00
## var_error 4.254e+01 4.284e+01 4.300e+01 4.316e+01 4.347e+01
## prp_var_eartag 1.357e-01 1.571e-01 1.699e-01 1.839e-01 2.141e-01
## prp_var_follower 1.198e-02 1.456e-02 1.613e-02 1.793e-02 2.199e-02
## prp_var_error 7.702e-01 8.000e-01 8.137e-01 8.263e-01 8.476e-01
## lp__ -1.593e+05 -1.592e+05 -1.592e+05 -1.592e+05 -1.592e+05print(Foll.model)## Inference for Stan model: EFtrp_model.
## 3 chains, each with iter=12000; warmup=2000; thin=1;
## post-warmup draws per chain=10000, total post-warmup draws=30000.
##
## mean se_mean sd 2.5% 25% 50%
## beta[1] 13.58 0.01 0.93 11.72 12.96 13.58
## beta[2] 7.62 0.01 0.96 5.74 6.98 7.62
## beta[3] 10.40 0.01 0.96 8.53 9.75 10.40
## beta[4] 8.47 0.01 1.13 6.26 7.72 8.48
## beta[5] 9.55 0.01 1.13 7.33 8.80 9.54
## beta[6] 5.41 0.01 1.08 3.26 4.68 5.41
## beta[7] 11.42 0.01 0.92 9.61 10.80 11.43
## beta[8] 8.40 0.01 0.96 6.48 7.75 8.40
## beta[9] 15.39 0.01 0.95 13.52 14.76 15.39
## beta[10] 7.50 0.01 1.12 5.29 6.76 7.49
## beta[11] 9.99 0.01 1.13 7.76 9.23 10.00
## beta[12] 7.53 0.01 1.07 5.45 6.82 7.53
## beta[13] 0.01 0.00 0.00 0.00 0.00 0.01
## var_eartag 9.09 0.01 1.30 6.89 8.18 8.98
## var_follower 0.87 0.00 0.14 0.63 0.77 0.85
## var_error 43.00 0.00 0.24 42.54 42.84 43.00
## prp_var_eartag 0.17 0.00 0.02 0.14 0.16 0.17
## prp_var_follower 0.02 0.00 0.00 0.01 0.01 0.02
## prp_var_error 0.81 0.00 0.02 0.77 0.80 0.81
## lp__ -159229.37 0.12 11.48 -159252.76 -159236.88 -159229.09
## 75% 97.5% n_eff Rhat
## beta[1] 14.21 15.38 5980 1
## beta[2] 8.25 9.52 5393 1
## beta[3] 11.05 12.30 5178 1
## beta[4] 9.22 10.69 6609 1
## beta[5] 10.30 11.76 6173 1
## beta[6] 6.13 7.53 5515 1
## beta[7] 12.04 13.21 5906 1
## beta[8] 9.05 10.26 5426 1
## beta[9] 16.04 17.27 5594 1
## beta[10] 8.26 9.72 6520 1
## beta[11] 10.74 12.18 6847 1
## beta[12] 8.24 9.63 5630 1
## beta[13] 0.01 0.01 57242 1
## var_eartag 9.88 11.95 20155 1
## var_follower 0.95 1.17 19172 1
## var_error 43.16 43.47 53180 1
## prp_var_eartag 0.18 0.21 20446 1
## prp_var_follower 0.02 0.02 19266 1
## prp_var_error 0.83 0.85 20360 1
## lp__ -159221.53 -159207.60 9068 1
##
## Samples were drawn using NUTS(diag_e) at Sun Sep 29 03:45:39 2019.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at
## convergence, Rhat=1).Posterior Correlation of model parameters
parcorplot(outp2,col = terrain.colors(15,0.5,T), cex.axis=0.75)
3. Mixed model: Location + Median wt + Eartag + Follower + covariance between eartag and follower random effects
3.1.Convergence Diagnostic
Check autorrelation, effective sample size, traceplot
ggs_autocorrelation(ggs(outp3)%>%filter(Parameter%in%lt1))
ggs_autocorrelation(ggs(outp3)%>%filter(Parameter%in%lt2))
ggs_autocorrelation(ggs(outp3)%>%filter(Parameter%in%c("Median Weight",
"rho", "var_eartag","var_follower", "var_error",
"prp_var_eartag","prp_var_follower", "prp_var_error")))
autocorr.diag(outp3)## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4 Loc1_t_5
## Lag 0 1.00000000 1.000000000 1.00000000 1.00000000 1.00000000
## Lag 1 0.61860344 0.670766419 0.66216005 0.62488624 0.60780882
## Lag 5 0.17781019 0.225939350 0.19704306 0.15431546 0.15758009
## Lag 10 0.04138631 0.069112818 0.03683993 0.03414666 0.03494142
## Lag 50 -0.01135936 0.007389174 0.02682615 0.01335822 0.01071376
## Loc1_t_6 Loc2_t_1 Loc2_t_2 Loc2_t_3 Loc2_t_4
## Lag 0 1.000000000 1.000000000 1.0000000000 1.00000000 1.00000000
## Lag 1 0.656078546 0.617990719 0.6553304565 0.64454478 0.61432343
## Lag 5 0.186447136 0.176769406 0.2138857890 0.18831438 0.15624435
## Lag 10 0.062029492 0.059820956 0.0692007090 0.04507550 0.03089223
## Lag 50 0.001010602 -0.007224273 0.0003442391 -0.00942489 0.01212918
## Loc2_t_5 Loc2_t_6 Median Weight rho var_eartag
## Lag 0 1.000000000 1.000000000 1.000000000 1.000000000 1.000000000
## Lag 1 0.606332923 0.642690825 -0.337047152 0.035353019 0.011332317
## Lag 5 0.132409011 0.187686625 -0.009086360 0.015448157 0.026623969
## Lag 10 0.026056938 0.054659868 0.001830421 0.005032227 -0.002153371
## Lag 50 0.005396061 -0.003549437 0.013549159 0.007607342 -0.006895332
## var_follower var_error prp_var_eartag prp_var_follower
## Lag 0 1.000000000 1.000000000 1.000000000 1.000000000
## Lag 1 0.044131219 -0.287035420 0.001622168 0.038775791
## Lag 5 0.016251067 -0.022551156 0.026483227 0.017650602
## Lag 10 0.002443573 0.001646444 -0.001656122 0.003473616
## Lag 50 0.001775547 -0.006015079 -0.007437080 0.001512490
## prp_var_error lp__
## Lag 0 1.000000000 1.0000000000
## Lag 1 0.005168091 0.5009938120
## Lag 5 0.025515323 0.0427589382
## Lag 10 -0.002170438 -0.0010363361
## Lag 50 -0.007339322 0.0009184027effectiveSize(outp3)## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4
## 6193.227 5194.153 5729.208 6276.004
## Loc1_t_5 Loc1_t_6 Loc2_t_1 Loc2_t_2
## 6489.591 5751.292 6221.700 5487.164
## Loc2_t_3 Loc2_t_4 Loc2_t_5 Loc2_t_6
## 5920.647 6732.954 6914.040 6179.950
## Median Weight rho var_eartag var_follower
## 60556.214 21775.197 20090.471 20307.905
## var_error prp_var_eartag prp_var_follower prp_var_error
## 54037.431 20435.764 20420.498 20217.172
## lp__
## 9858.165traplot(outp3,col =c("red1","blue4","purple3"))
geweke.diag(outp3)## [[1]]
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4
## 0.7466 -0.3187 -4.0814 1.9540
## Loc1_t_5 Loc1_t_6 Loc2_t_1 Loc2_t_2
## -1.1884 -1.4185 1.9389 0.2663
## Loc2_t_3 Loc2_t_4 Loc2_t_5 Loc2_t_6
## -2.5897 -0.4564 1.1675 0.3609
## Median Weight rho var_eartag var_follower
## -1.7180 -2.1725 0.3222 0.9613
## var_error prp_var_eartag prp_var_follower prp_var_error
## -2.1611 0.4195 0.9206 -0.5420
## lp__
## -0.6094
##
##
## [[2]]
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4
## 1.30399 1.74579 -0.44428 -1.33039
## Loc1_t_5 Loc1_t_6 Loc2_t_1 Loc2_t_2
## 0.32080 0.13947 -0.04603 -1.21169
## Loc2_t_3 Loc2_t_4 Loc2_t_5 Loc2_t_6
## -0.99788 0.91100 -0.13264 -0.39034
## Median Weight rho var_eartag var_follower
## 0.33713 -1.93620 -1.24825 -0.24059
## var_error prp_var_eartag prp_var_follower prp_var_error
## 0.25971 -1.26234 -0.03771 1.28407
## lp__
## -0.99320
##
##
## [[3]]
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## Loc1_t_1 Loc1_t_2 Loc1_t_3 Loc1_t_4
## 0.22681 -1.49278 0.20749 -0.23288
## Loc1_t_5 Loc1_t_6 Loc2_t_1 Loc2_t_2
## 0.15888 0.48129 -0.63671 1.09432
## Loc2_t_3 Loc2_t_4 Loc2_t_5 Loc2_t_6
## -0.53911 1.43280 -0.30301 0.02446
## Median Weight rho var_eartag var_follower
## 0.85535 0.12719 0.37711 -0.43710
## var_error prp_var_eartag prp_var_follower prp_var_error
## -0.00472 0.45078 -0.45093 -0.39220
## lp__
## 0.33149#gelman.diag(outp3)
denplot(outp3,col = c("red1","blue4","purple3"))
3.2. Summary Posterior Distribution
summary(outp3)##
## Iterations = 2001:12000
## Thinning interval = 1
## Number of chains = 3
## Sample size per chain = 10000
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## Loc1_t_1 1.354e+01 0.960095 5.543e-03 1.344e-02
## Loc1_t_2 7.635e+00 0.981886 5.669e-03 1.488e-02
## Loc1_t_3 1.041e+01 0.978042 5.647e-03 1.392e-02
## Loc1_t_4 8.481e+00 1.143923 6.604e-03 1.531e-02
## Loc1_t_5 9.562e+00 1.158813 6.690e-03 1.530e-02
## Loc1_t_6 5.414e+00 1.081300 6.243e-03 1.533e-02
## Loc2_t_1 1.146e+01 0.949496 5.482e-03 1.300e-02
## Loc2_t_2 8.376e+00 0.986848 5.698e-03 1.459e-02
## Loc2_t_3 1.537e+01 0.984655 5.685e-03 1.367e-02
## Loc2_t_4 7.525e+00 1.133340 6.543e-03 1.493e-02
## Loc2_t_5 9.992e+00 1.145833 6.615e-03 1.467e-02
## Loc2_t_6 7.549e+00 1.094912 6.321e-03 1.548e-02
## Median Weight 6.045e-03 0.001574 9.086e-06 6.402e-06
## rho 7.103e-02 0.105076 6.067e-04 7.193e-04
## var_eartag 9.169e+00 1.321209 7.628e-03 9.399e-03
## var_follower 8.753e-01 0.138160 7.977e-04 9.751e-04
## var_error 4.300e+01 0.236620 1.366e-03 1.019e-03
## prp_var_eartag 1.723e-01 0.020352 1.175e-04 1.436e-04
## prp_var_follower 1.650e-02 0.002588 1.494e-05 1.825e-05
## prp_var_error 8.111e-01 0.020090 1.160e-04 1.424e-04
## lp__ -1.592e+05 11.427947 6.598e-02 1.151e-01
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## Loc1_t_1 1.166e+01 1.291e+01 1.355e+01 1.419e+01 1.541e+01
## Loc1_t_2 5.724e+00 6.970e+00 7.633e+00 8.303e+00 9.540e+00
## Loc1_t_3 8.448e+00 9.758e+00 1.042e+01 1.106e+01 1.231e+01
## Loc1_t_4 6.230e+00 7.714e+00 8.503e+00 9.258e+00 1.069e+01
## Loc1_t_5 7.288e+00 8.780e+00 9.558e+00 1.033e+01 1.187e+01
## Loc1_t_6 3.284e+00 4.698e+00 5.404e+00 6.141e+00 7.550e+00
## Loc2_t_1 9.598e+00 1.082e+01 1.145e+01 1.211e+01 1.332e+01
## Loc2_t_2 6.464e+00 7.705e+00 8.380e+00 9.035e+00 1.032e+01
## Loc2_t_3 1.346e+01 1.471e+01 1.536e+01 1.602e+01 1.732e+01
## Loc2_t_4 5.275e+00 6.763e+00 7.540e+00 8.297e+00 9.712e+00
## Loc2_t_5 7.760e+00 9.222e+00 9.992e+00 1.077e+01 1.225e+01
## Loc2_t_6 5.384e+00 6.815e+00 7.552e+00 8.282e+00 9.692e+00
## Median Weight 2.948e-03 4.993e-03 6.047e-03 7.100e-03 9.123e-03
## rho -1.380e-01 9.629e-04 7.127e-02 1.430e-01 2.721e-01
## var_eartag 6.924e+00 8.240e+00 9.049e+00 9.960e+00 1.212e+01
## var_follower 6.442e-01 7.783e-01 8.620e-01 9.582e-01 1.184e+00
## var_error 4.254e+01 4.284e+01 4.300e+01 4.316e+01 4.347e+01
## prp_var_eartag 1.363e-01 1.581e-01 1.710e-01 1.850e-01 2.167e-01
## prp_var_follower 1.216e-02 1.468e-02 1.627e-02 1.806e-02 2.228e-02
## prp_var_error 7.678e-01 7.987e-01 8.125e-01 8.252e-01 8.468e-01
## lp__ -1.593e+05 -1.592e+05 -1.592e+05 -1.592e+05 -1.592e+05print(EFcov.model)## Inference for Stan model: covEartFolltrp_model.
## 3 chains, each with iter=12000; warmup=2000; thin=1;
## post-warmup draws per chain=10000, total post-warmup draws=30000.
##
## mean se_mean sd 2.5% 25% 50%
## beta[1] 13.54 0.01 0.96 11.66 12.91 13.55
## beta[2] 7.63 0.01 0.98 5.72 6.97 7.63
## beta[3] 10.41 0.01 0.98 8.45 9.76 10.42
## beta[4] 8.48 0.02 1.14 6.23 7.71 8.50
## beta[5] 9.56 0.02 1.16 7.29 8.78 9.56
## beta[6] 5.41 0.02 1.08 3.28 4.70 5.40
## beta[7] 11.46 0.01 0.95 9.60 10.82 11.45
## beta[8] 8.38 0.01 0.99 6.46 7.70 8.38
## beta[9] 15.37 0.01 0.98 13.46 14.71 15.36
## beta[10] 7.53 0.02 1.13 5.28 6.76 7.54
## beta[11] 9.99 0.01 1.15 7.76 9.22 9.99
## beta[12] 7.55 0.02 1.09 5.38 6.82 7.55
## beta[13] 0.01 0.00 0.00 0.00 0.00 0.01
## rho 0.07 0.00 0.11 -0.14 0.00 0.07
## var_eartag 9.17 0.01 1.32 6.92 8.24 9.05
## var_follower 0.88 0.00 0.14 0.64 0.78 0.86
## var_error 43.00 0.00 0.24 42.54 42.84 43.00
## prp_var_eartag 0.17 0.00 0.02 0.14 0.16 0.17
## prp_var_follower 0.02 0.00 0.00 0.01 0.01 0.02
## prp_var_error 0.81 0.00 0.02 0.77 0.80 0.81
## lp__ -159230.53 0.12 11.43 -159253.75 -159238.00 -159230.30
## 75% 97.5% n_eff Rhat
## beta[1] 14.19 15.41 5311 1
## beta[2] 8.30 9.54 4399 1
## beta[3] 11.06 12.31 4399 1
## beta[4] 9.26 10.69 5718 1
## beta[5] 10.33 11.87 5598 1
## beta[6] 6.14 7.55 4714 1
## beta[7] 12.11 13.32 5110 1
## beta[8] 9.03 10.32 4482 1
## beta[9] 16.02 17.32 5181 1
## beta[10] 8.30 9.71 5707 1
## beta[11] 10.77 12.25 5970 1
## beta[12] 8.28 9.69 4553 1
## beta[13] 0.01 0.01 61991 1
## rho 0.14 0.27 20727 1
## var_eartag 9.96 12.12 19063 1
## var_follower 0.96 1.18 19549 1
## var_error 43.16 43.47 55476 1
## prp_var_eartag 0.19 0.22 19352 1
## prp_var_follower 0.02 0.02 19752 1
## prp_var_error 0.83 0.85 19159 1
## lp__ -159222.57 -159209.16 9539 1
##
## Samples were drawn using NUTS(diag_e) at Sun Sep 29 04:04:44 2019.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at
## convergence, Rhat=1).Posterior Correlation of model parameters
parcorplot(outp3,col = terrain.colors(15,0.5,T), cex.axis=0.75 )