Calculating accuracy and reliability of random effect estimates with ASReml-R

An important aim when fitting linear mixed models (LMM) is the exploration of the random effect estimates. In some analyses, such as genetic evaluations, the main objective of the analysis is to obtain these estimates. Random effects, also known as BLUPs (Best Linear Unbiased Predictions), are obtained after fitting a LMM, where their variance components are estimated by residual maximum likelihood (REML) and these components are then used to calculate BLUPs.

Consider the following LMM for a genetic experiment, where we have a randomized complete block design with blocks and a total of genotypes evaluated:

where is the phenotypic observation of the genotype on the block; is the overall mean; is the fixed effect of the block; is the random effect of the genotype, with that is associated with some information about pedigree, from which we can obtain a pedigree-based or genomic-based relationship matrix ; and is the random residual with

In the above model, also known as Animal Model, the random effects correspond to the breeding value (BV) from which selections of future progenitors or individuals will be done. As these effects are based on estimated variance components, then they are often known as EBVs (estimated breeding values). Statistically, the matrix  is critical as it is the place where we specify the correlations or dependencies that exist between random effects.

As many breeding and commercial decisions depend on the EBVs it is important to have an assessment of our confidence in these values, and this is done by calculating a statistical measure of precision. Probably the most common statistic used to assess these effects is the calculation of the correlation between true and predicted random effects, , also known as accuracy. This can be calculated using the following expression for a given genotype as:

This formula requires the PEV (predictor error variance) of , which corresponds to the square of the standard error (SE) of the random effect estimate, i.e. . Most statistical software, including SAS and ASReml-R, will report these SE values. In addition, we have , which is the population estimate of the variance associated with the EBVs, and this is our genetic additive variance, which is used to calculate narrow-sense heritability.

According to statistical theory, the values of PEV will range from zero to . In cases where we have full information about a given genotype, then the PEV will be close to zero as we have little uncertainty in the true additive value of a genotype. However, when there is very little information, then the PEV will approximate to . Note that these extreme values of PEV translate in values of 1.0 and 0.0, respectively. Therefore, values closer to 1.0 are an indication of very good quality of a given additive effect estimate.

Sometimes, you will see reliability instead of accuracy reported, where the former is the square of the later. Both statistics are commonly used when reporting the quality of a random effect estimate in genetic studies.

To illustrate the above calculations, we used genetic data originating from a trial on a population of cultivated two-row spring barley. The original trial measured height (cm) of plants grown in pots during 2011 for a total of 856 lines; however, the available data for this analysis consists of the adjusted mean values of only 478 lines, together with their molecular information. Further details can be found in Oakley et al. (2016).

In our example, we obtained the realized relationship matrix from the single nucleotide polymorphism (SNP) marker information available for these 478 lines, and an Animal Model was fitted considering this genomic matrix. This was done in the package R using the library ASReml-R, and the code used to fit this model was:

modelGBLUP <- asreml(fixed=Hadj11~1,
                     random=~vm(lines, ainv),
                     na.action=na.method(y='include'),
                     data=datagFULL)

Here, the response variable was and our single random effect is the factor that represents our additive effects, and this is associated with the inverse of the relationship matrix . In addition, in the above code it was requested that observations with missing data (i.e., ) in their response variable will be included in the analyses. This allowed us to obtain EBVs for all individuals, even those that were not phenotyped but that were included in the relationship matrix.

The above analysis, using the following code

vc <- summary(modelGBLUP)$varcomp

will output the following estimation of variance components.

These values correspond to a narrow-sense heritability of for the adjusted mean of height. Hence, the trait has a strong genetic control. Now, we can proceed to request the BLUP (or EBV) values together with their standard errors, this is done using the code:

BLUP <- summary(modelGBLUP, coef=TRUE)$coef.random

Below we can see a subset of this list with the top 10 and bottom 10 genotypes sorted by their SE values.

The EBV for this data are positive and negative, as they represent deviations from the overall height mean. For example, the genotype has an EBV of indicating that, if this parent is crossed with a parent of equal genetic worth, its progeny will have an expected deviation from the mean of cm. The full range of EBVs is between cm to cm, a wide interval considering that the mean of the population is cm, and this is a result of the relatively large heritability value found for this trait.

Also, there is an interesting range in SE values from to . The width of this range is related to the level of genetic relationships between individuals, that is specified through the matrix . Genotypes with many relatives will share or collect information from many phenotypic responses, and therefore they will have lower SE values. In contrast, individuals with few, or no relatives, will tend to have much higher SE values. However, it is also possible that individuals with no phenotypic data will have small SE values. This occurs when they have many direct relatives to help with the estimation of their EBVs.

Now, we can proceed to calculate the accuracy and reliability using the expressions shown before, with the code

Va <- vc[1,1]
BLUP$PEV <- BLUP$std.error^2
BLUP$accuracy <- sqrt(1 - BLUP$PEV/Va)
BLUP$reliability <- 1 - BLUP$PEV/Va

And again we can explore these calculated values for a subset of the top 10 and bottom 10 individuals.

For this trait, the range of correlation between true and predicted breeding values, , is from to , and the lowest accuracy and reliability values correspond to individuals with the largest SE values. Note that the genotype has a for accuracy, which originates from a negative square root, as this BLUP has an SE that is larger than .

In the figure below we present a boxplot of the accuracy values (left plot) and a scatterplot of the PEV against accuracy (right plot). It is clear from these plots that there is a group of three observations with very poor accuracy that have values lower than . However, the majority of the observations are found with accuracy values of or higher, where some of them reach values as high as .

It is interesting to speculate what the reasons could be for these three lower observations (genotypes , and ) together with the observation (genotype ), to be so different from the rest of the population. The unusual behaviour of these four observations is likely to reflect errors in their information, including phenotypic response (maybe a potential outlier), or incorrect genetic relationships (due to pedigree or genotyping errors).

In summary, we have calculated values of accuracy and reliability for the estimation of random effects from an Animal Model. For the data analyzed, these represented a wide range of values, reflecting different levels of precision. For example, a large value of accuracy , such as , will indicate that there is little expected change on the estimation of the EBV if additional information is collected from a given genotype or from its relatives (e.g., its offspring). Therefore, values provide us with a good measure of our confidence in the genetic evaluation performed.

In our analysis, and in many other genetic analyses, it is critical to consider the incorporation of the matrix of genetic relationships, as it will allow us to improve the estimation of the random effects because we have changed our assumption from independent random effects to a model that considers correlations or dependencies between these random effects.

A routine check of the accuracy and reliability (and indirectly of PEV values), can help to make future improvements on the quality of the information obtained from statistical analyses. For example, we can assess the consequences of incorporating additional relatives or replicates in future evaluations. Alternatively, these verifications can also help to assess the benefits of using more complex linear models, such as those that combine multiple traits, measurements or sites.

Finally, the above calculations can be extended to any other random effect on all linear mixed models fitted, and their use is recommended to assess the quality of the estimates beyond what a simple standard error will indicate.

About the author

Dr. Salvador Gezan is a statistician/quantitative geneticist with more than 20 years’ experience in breeding, statistical analysis and genetic improvement consulting. He currently works as a Statistical Consultant at VSN International, UK. Dr. Gezan started his career at Rothamsted Research as a biometrician, where he worked with Genstat and ASReml statistical software. Over the last 15 years he has taught ASReml workshops for companies and university researchers around the world. 

Dr. Gezan has worked on agronomy, aquaculture, forestry, entomology, medical, biological modelling, and with many commercial breeding programs, applying traditional and molecular statistical tools. His research has led to more than 100 peer reviewed publications, and he is one of the co-authors of the textbook Statistical Methods in Biology: Design and Analysis of Experiments and Regression.

Reference

Oakey, H; Cullis, B; Thompson, R; Comadran, J; Halpin, C; Waugh, R. 2016. Genomic selection in multi-environment crop trials. Genes | Genomes | Genetics 6:1313-1326

Notes: SAG Jan-2020