Table of Contents

Example 1 for WOMBAT

Download Example1

One of the simplest analysis in quantitative genetics is the case where we have one trait of interest which is recorded only once, and where we have no additional random effects. Assuming a general pedigree (i.e. information from various types of relatives), we then fit a simple animal model, with additive genetic effects due to the individual to estimate the additive genetic and residual variance components and calculate the heritability of the trait from them.

Data and pedigree files

This is illustrated in Example 1 using some simulated data. The first few lines of the file tstmo1.d are shown in this box:

   10025       1       2   1  220.0000
   10026       1       2   1  212.0000
   10027       1       2   1  221.0000
   10028       1       2   1  207.0000
   10029       1       2   1  218.0000
   10030       1       2   1  201.0000
   10031       1       2   1  214.0000
   10032       1       2   1  229.0000

Column 1 contains the animal identity and columns 2 and 3 give the codes for its sire and dam, respectively. Note that all are numeric codes and that animal codes must be higher than both sire and dam codes. In this case, we have a single fixed effect – generation – with 2 levels which is listed in column 4. Column 5 finally contains the trait value. As the data file has only one record per individual and as it contains the pedigree information for all individuals (i.e. there is no additional pedigree information on parents without records), we can use the same file as pedigree file.

Parameter file

The parameter file to analyse these data then looks like this (note that line numbers are NOT part of the parameter files and added here for ease of referencing only!):

  1. COM Example 1 from DFREML : Simple univariate analysis
  2. PED ../tstmo1.d
  3. DAT ../tstmo1.d
  4. animal
  5. sire
  6. dam
  7. fixeff 9
  8. record
  9. end
  10. ANAL UNI
  11. MODEL
  12. RAN animal NRM
  13. FIX fixeff
  14. TR record
  15. END
  16. VAR animal 1
  17. 60
  18. VAR error 1
  19. 40

Lines within a `block' are indented here – this is for clarity (and recommended), but not actually required by the program.

In this example, we call the parameter file wombat.par.

Running WOMBAT

Running WOMBAT then requires typing in the program name, possibly with the appropriate path, together with any run time options desired on the `command line' Here we use

wombat -v

N.B. This requires a set-up where the operating system knows where to look for the executable file! Otherwise, you may need to specify the complete path, e.g. /home/kmeyer/bin/wombat. The run time option -v requests `verbose' output to the screen. Since our parameter file has the default name wombat.par, we don't need to give its name.

Screen output

This generates the screen output as shown in the following:

 Parameter file opened : "wombat.par"
**COM Example 1 from DFREML : Simple univariate analysis**
**PED ../tstmo1.d**
  Pedigree file : "../tstmo1.d"
**DAT ../tstmo1.d**
  Data file : "../tstmo1.d"
**animal**
  Data - column :   1  "animal"    -1
**sire**
  Data - column :   2  "sire"    -1
**dam**
  Data - column :   3  "dam"    -1
**fixeff 9**
  Data - column :   4  "fixeff"     9
**record**
  Data - column :   5  "record"    -1
**end**
**ANAL UNI**
  Analysis type = 1  "UNI"
  No. of traits =  1
**MODEL**
**RAN animal NRM**
**FIX fixeff**
**TR record**
  Trait no.  1  "record"  Missing value  -123456789
**END**
**VAR animal 1**
**60**
**VAR error 1**
**40**

 End of parameter file reached :  19 lines read
 Data file appers to be sorted correctly
 No. of records                             =       282
 No. of animal IDs in data                  =       282
 No. of animal IDs in data                  =       282
 Total no. of IDs in reduced ped file       =       306
 Pedigree file                  =   ../tstmo1.d                   
 Total no. of animal IDs found  = 306
 no. of NRM matrices  1 1
 Pruned pedigree for random effect no. 1
 Original no. of animals = 306
 New no. of animals      = 306

Inbreeding coefficients for random effect   1 computed
No. of inbred animals          =      38
Average inbreeding coefficient =  2.2059 (in %)
... amongst inbred animals     = 17.7632 (in %)
 NRM Inverse no. : 1   "animal"
 No. of elements = 906   log determinant = -195.46750491790496
 No. of "subjects" in data 282
  maxnr = 2   ncombi = 1
 Parameters to be estimated & their starting values
    1  COVS Z 1 1           40.000000    
    2  COVS A 1 1           60.000000    
 END data
 Coefficient matrix fixed effects - rank =     2  no. of rows =     2
 END "lsqfixed"
 No. of equations in full rank submatrix = 308
 Elements in MMM - after data                 282
 Elements in MMM - after RE    1              882
 End of "GENMMD"      308  nofsub =           950
 Elements in MMM - after RHS                 1166
 End of "SMBFCT": maxsub                      318
                  maxlnz                     1274
                  nnops                      4599.  999999999999999.
 END "eqnorders" 
***log L    1          -760.345   Time       0
 search 1
***AI iterate     0  log L =  -760.344597      D =   0.0000       0.0000      Time      0
 AI matrix : min. & max. eigenvalue 0.580946      6.50567      tau =   0.0000    
***log L    2          -759.598   Time       0
 Tau =  0.0000     Alpha =  1.0000      log L =        -759.598
***AI iterate     1  log L =  -759.597721      D =  0.74688      0.18539      Time      0
 Estimates for "tau" =  0.0000       "alpha" =  1.0000       "norm" =  1.6662    
   1   CHOL Z 1 1       Old =    6.324555     New =    7.396408                   
   2   CHOL A 1 1       Old =    7.745967     New =    6.298573                   
 AI matrix : min. & max. eigenvalue 0.689232      6.18917      tau =   0.0000    
***log L    3          -759.507   Time       0
 Tau =  0.0000     Alpha =  1.0000      log L =        -759.507
***AI iterate     2  log L =  -759.506748      D =  0.90973E-01  0.52339E-01  Time      0
 Estimates for "tau" =  0.0000       "alpha" =  1.0000       "norm" = 0.64462    
   1   CHOL Z 1 1       Old =    7.396408     New =    7.097490                   
   2   CHOL A 1 1       Old =    6.298573     New =    6.713435                   
 AI matrix : min. & max. eigenvalue 0.637744      6.31627      tau =   0.0000    
***log L    4          -759.504   Time       0
 Tau =  0.0000     Alpha =  1.0000      log L =        -759.504
***AI iterate     3  log L =  -759.503659      D =  0.30884E-02  0.11687E-01  Time      0
 Estimates for "tau" =  0.0000       "alpha" =  1.0000       "norm" = 0.72789E-01
   1   CHOL Z 1 1       Old =    7.097490     New =    7.148321                   
   2   CHOL A 1 1       Old =    6.713435     New =    6.611619                   
 AI matrix : min. & max. eigenvalue 0.652444      6.34003      tau =   0.0000    
***log L    5          -759.503   Time       0
 Tau =  0.0000     Alpha =  1.0000      log L =        -759.503
***AI iterate     4  log L =  -759.503454      D =  0.20554E-03  0.29391E-02  Time      0
 Estimates for "tau" =  0.0000       "alpha" =  1.0000       "norm" = 0.19922E-01
   1   CHOL Z 1 1       Old =    7.148321     New =    7.134470                   
   2   CHOL A 1 1       Old =    6.611619     New =    6.636684                   
*** AIREML has converged
    Last change in log L 2.05541655532215373E-4 5.0000000000000001E-4
   1  animal                   306  zero soln.s       0
 SS residuals 9470.5877529900336
 "WOMBAT" has finished !

Results files

Summary of pedigree information

The first output file generated is a summary of the pedigree information found. This file is named SumPedigree.out.

======= Version 15-09-2009 ======================================= **KM** ====
 
        Program WOMBAT : Summary of Pedigree Information
==============================================================================
 
 Example 1 from DFREML : Simple univariate analysis                         
 
 Analysis type          :   "UNI"
 Data file              :   "../tstmo1.d"
 Pedigree file          :   "../tstmo1.d"
 Parameter file         :   "wombat.par"
 No. of animal IDs in data file =           =       282
 No. of animal IDs in total     =           =       306
 
*****Pedigree Structure for random effect :    1  ****************************
 Original no. of animals                    =       306
 No. of animals after pruning               =       306
 ... proportion (%) remaining               =     100.0
 
 No. of levels w/out records                =        24
 No. of levels with records                 =       282  100.0%
  ...  1 record(s)                          =       282  100.0%
 
 No. of animals w/out offspring             =       258   84.3%
 No. of animals with offspring              =        48   15.7%
  ... and records                           =        24    7.8%
 
 No. of animals with unknown sire           =        24
 No. of animals with unknown dam            =        24
 No. of animals with both parents unknown   =        24
 No. of animals with records                =
  ... and unknown sire                      =         0
  ... and unknown dam                       =         0
  ... and both parents unknown              =         0
 No. of sires                               =        12
  ... with progeny in the data              =        12
  ... with records & progeny in data        =         6
 No. of dams                                =        36
  ... with progeny in the data              =        36
  ... with records & progeny in data        =        18
 
 No. of animals with known/unpruned grand-parents
  ... with paternal grandsire               =       144
  ... with paternal granddam                =       144
  ... with maternal grandsire               =       144
  ... with maternal granddam                =       144
 
Inbreeding coefficients for random effect   1 computed
No. of inbred animals          =      38
Average inbreeding coefficient =  2.2059 (in %)
... amongst inbred animals     = 17.7632 (in %)
  random effect no. = 1 NRM
  no. of elements in NRM/GIN inverse  906
  log determinant   = -195.46750491790496
======== end of file ============================21-09-2009==========15:40====

Summary for data file and model specified

The second summary file generated at the beginning of a run is SumModel.out. It gives a summary of the number of records, means and standard deviations per trait and the number of levels for each effect fitted that WOMBAT has found.

======= Version 15-09-2009 ======================================= **KM** ====
 
        Program WOMBAT : Summary of information from Set-up step
==============================================================================
 
 Example 1 from DFREML : Simple univariate analysis                         
 
 Analysis type          :   "UNI"
 Data file              :   "../tstmo1.d"
 Pedigree file          :   "../tstmo1.d"
 Parameter file         :   "wombat.par"
 
 
 No. of traits          =   1
                      nrec       mean         sdev         min.         max.
   1  "record"         282    229.787      13.4938      194.000      269.000    
 
 Fixed effects
 1 "record"           nlev
   1  "fixeff"           2
 
 Random effects       nlev
   1  "animal"         306    NRM
======== end of file ============================21-09-2009==========15:40====

Summary of variance component estimates

The main results file is named SumEstimates.out.

======= Version 15-09-2009 ======================================= **KM** ====
 
        Program WOMBAT : Estimates of covariance components
==============================================================================
 
 Example 1 from DFREML : Simple univariate analysis                         
 
 Analysis type          :   "UNI"
 Data file              :   "../tstmo1.d"
 Pedigree file          :   "../tstmo1.d"
 Parameter file         :   "wombat.par"
 
 
 No. of traits          =      1                record
 No. of records         =    282           282
 No. of parameters      =      2
 Maximum log L          =          -759.503
 -1/2 AIC & AICC        =          -761.503          -761.525
 -1/2 BIC               =          -765.138    "Penalty factor" =   2.817
 
 Parameter estimates with approx. sampling erors
   1   CHOL Z 1 1              7.13447          0.662468              
   2   CHOL A 1 1              6.63668           1.12245              
 
 Convergence criteria for last 3 iterates
 Change in log likelihood    =    0.090973    0.003088    0.000206
 Change in parameter vector  =    0.052339    0.011687    0.002939
 Norm of gradient vector     =      0.6446      0.0728      0.0199
 Newton decrement            =     -0.2226     -0.0083     -0.0005
***** Estimates of residual covariances   ************************************
      Order of fit         =       1
 Covariance matrix
   1     50.901    
 Matrix of correlations and variance ratios
   1      0.5361
 Covariances & correlations with approximate sampling errors
   1  COVS Z 1 1        50.9007        9.45272       vrat     0.536   0.123
 
***** Estimates for RE  1   "animal"   ***************************************
      No. of levels        =     306
      Covariance structure =    NRM
      Order of fit         =       1
 Covariance matrix
   1     44.046    
 Matrix of correlations and variance ratios
   1      0.4639
 Covariances & correlations with approximate sampling errors
   2  COVS A 1 1        44.0456        14.8987       vrat     0.464   0.123
 
***** Estimates of phenotypic covariances  ***********************************
 Covariance matrix
   1     94.946    
 Covariances & correlations with approximate sampling errors
   3  COVS T 1 1        94.9462        10.0250    
======== end of file ============================21-09-2009==========15:40====

Other output to look at