Single trait model aims to estimate variance components and its heritability of all random effects for a trait. Always remember that there is no need to adjust the order of individuals or markers to be consistent among different files, HIBLUP can do it automatically. Example command of running SSGLBUP for single trait model:
./hiblup --single-trait
--pheno demo.phe
--pheno-pos 8
--dcovar 2,3 #fixed effect
--qcovar 4,5 #covariates
--rand 6,7 #non-genetic (environmental) random effect
--pedigree demo.ped #genetic random effect
--bfile demo #genetic random effect
--add --dom
--threads 32
--out demo
The detailed format of input files please refer to phenotype, genotype. The way of how to set covariates (--qcovar), fixed (--dcovar) and random (--rand) effect for multiple traits, please refer to here.
Note that all the parameters (e.g., weighted GRM, mixed GRM, or HRM, …) supported for XRM construction can also be applied in the above command.
There are 5 types of algorithms for variance components estimation:
“AI”, “EM”, “HE”, “EMAI”, “HI”, where “AI” is AIREML,”EM” is EMREML, and “HE” represents HE regression, “HI” is “HE+AI”, which uses the outcome of HE regression as the start values of AIREML, users can choose one of them by flag --vc-method, and change the maximum iteration number for AIREML and EMREML by flag --ai-maxit and --em-maxit respectively. For “HE” regression, if covariates for fixed effects are provided, HIBLUP will first regress the phenotype on these, then perform HE Regression using the residuals, the standard error(SE) of HE is computed by Jackknife resampling strategy at 100 repeats. By default, HIBLUP doesn’t calculate the random effects for HE regression, user can add a flag --he-pred to output it, and only in this situation, the flag --pcg can be used for fast computing on very large dataset, which we call it as “HE+PCG” strategy.
However, above command maybe time-consuming when fitting trait by trait, as it repeatedly constructs XRM for each trait. We can construct XRM only once at the beginning, then assign it to --xrm when fitting multiple traits, which can give the same results, but would be much more efficient, the command can be referred as following:
# Step 1: construct XRMs
./hiblup --make-xrm
--pedigree demo.ped
--bfile demo
--add --dom
--out demo
# Step 2: link XRMs to fit model
./hiblup --single-trait
--pheno demo.phe
--pheno-pos 8
--dcovar 2,3
--qcovar 4,5
--rand 6,7
--xrm demo.HA,demo.HD
--add --dom
--threads 32
--out demo
Several files will be generated in the work directory:demo.vars: the estimated variance and heritability of all random effects:
Var Var_SE h2 h2_SE loc 12.2127 8.9109 0.1035 0.0682 dam 10.6268 4.8605 0.0901 0.041 HA 59.2007 11.7641 0.5019 0.0807 HD 28.9946 8.9294 0.2458 0.083 e 6.9232 1.6957 0.0587 0.0161
demo.beta: the estimated coefficients of all covariates and fixed effects:
Levels Estimation SE mu 32.5285 3.4758 sex_F 6.2933 1.7593 sex_M 26.2353 1.7634 season_Autumn 18.0905 0.9753 season_Spring -1.9168 0.9904 season_Summer 7.8471 1.0066 season_Winter 8.5078 0.9944 day 0.1547 0.0574 bornweight 1.5703 0.4614
demo.anova: analysis of variance table of all covariates and fixed effects:
Factors Df SumSq MeanSq F Pr(>F) sex 1 49710.3845 49710.3845 7180.311 <2e-16 season 3 24515.0941 8171.698 1180.344 <2e-16 day 1 925.3352 925.3352 133.658 <2e-16 bornweight 1 437.3432 437.3432 63.171 1.30E-14 e 493 3413.1141 6.9232
demo.rand: the estimated environmental random effects, genetic random effects (the column named “PA”, “GA”, or “HA” is the list of breeding values for phenotypic and non-phenotypic individuals) and residuals.
Users can specify the start values for all unknown variance to --vc-priors for fast convergence, the way of how to input the variance components can be found at here.