Assumptions and Properties of Limiting Pathway Models for Analysis of Epistasis in Complex Traits

For most complex traits, results from genome-wide association studies show that the proportion of the phenotypic variance attributable to the additive effects of individual SNPs, that is, the heritability explained by the SNPs, is substantially less than the estimate of heritability obtained by standard methods using correlations between relatives. This difference has been called the “missing heritability”. One explanation is that heritability estimates from family (including twin) studies are biased upwards. Zuk et al. revisited overestimation of narrow sense heritability from twin studies as a result of confounding with non-additive genetic variance. They propose a limiting pathway (LP) model that generates significant epistatic variation and its simple parametrization provides a convenient way to explore implications of epistasis. They conclude that over-estimation of narrow sense heritability from family data (‘phantom heritability’) may explain an important proportion of missing heritability. We show that for highly heritable quantitative traits large phantom heritability estimates from twin studies are possible only if a large contribution of common environment is assumed. The LP model is underpinned by strong assumptions that are unlikely to hold, including that all contributing pathways have the same mean and variance and are uncorrelated. Here, we relax the assumptions that underlie the LP model to be more biologically plausible. Together with theoretical, empirical, and pragmatic arguments we conclude that in outbred populations the contribution of additive genetic variance is likely to be much more important than the contribution of non-additive variance.     

Assumption-free estimation of heritability from genome-wide identity-by-descent sharing between full siblings

The study of continuously varying, quantitative traits is important in evolutionary biology, agriculture, and medicine. Variation in such traits is attributable to many, possibly interacting, genes whose expression may be sensitive to the environment, which makes their dissection into underlying causative factors difficult. An important population parameter for quantitative traits is heritability, the proportion of total variance that is due to genetic factors. Response to artificial and natural selection and the degree of resemblance between relatives are all a function of this parameter. Following the classic paper by R. A. Fisher in 1918, the estimation of additive and dominance genetic variance and heritability in populations is based upon the expected proportion of genes shared between different types of relatives, and explicit, often controversial and untestable models of genetic and non-genetic causes of family resemblance. With genome-wide coverage of genetic markers it is now possible to estimate such parameters solely within families using the actual degree of identity-by-descent sharing between relatives. Using genome scans on 4,401 quasi-independent sib pairs of which 3,375 pairs had phenotypes, we estimated the heritability of height from empirical genome-wide identity-by-descent sharing, which varied from 0.374 to 0.617 (mean 0.498, standard deviation 0.036). The variance in identity-by-descent sharing per chromosome and per genome was consistent with theory. The maximum likelihood estimate of the heritability for height was 0.80 with no evidence for non-genetic causes of sib resemblance, consistent with results from independent twin and family studies but using an entirely separate source of information. Our application shows that it is feasible to estimate genetic variance solely from within-family segregation and provides an independent validation of previously untestable assumptions. Given sufficient data, our new paradigm will allow the estimation of genetic variation for disease susceptibility and quantitative traits that is free from confounding with non-genetic factors and will allow partitioning of genetic variation into additive and non-additive components.     

Estimating Modifying Effect of Age on Genetic and Environmental Variance Components in Twin Models

Twin studies have been adopted for decades to disentangle the relative genetic and environmental contributions for a wide range of traits. However, heritability estimation based on the classical twin models does not take into account dynamic behavior of the variance components over age. Varying variance of the genetic component over age can imply the existence of gene–environment (G × E) interactions that general genome-wide association studies (GWAS) fail to capture, which may lead to the inconsistency of heritability estimates between twin design and GWAS. Existing parametric G × E interaction models for twin studies are limited by assuming a linear or quadratic form of the variance curves with respect to a moderator that can, however, be overly restricted in reality. Here we propose spline-based approaches to explore the variance curves of the genetic and environmental components. We choose the additive genetic, common, and unique environmental variance components (ACE) model as the starting point. We treat the component variances as variance functions with respect to age modeled by B-splines or P-splines. We develop an empirical Bayes method to estimate the variance curves together with their confidence bands and provide an R package for public use. Our simulations demonstrate that the proposed methods accurately capture dynamic behavior of the component variances in terms of mean square errors with a data set of >10,000 twin pairs. Using the proposed methods as an alternative and major extension to the classical twin models, our analyses with a large-scale Finnish twin data set (19,510 MZ twins and 27,312 DZ same-sex twins) discover that the variances of the A, C, and E components for body mass index (BMI) change substantially across life span in different patterns and the heritability of BMI drops to ∼50% after middle age. The results further indicate that the decline of heritability is due to increasing unique environmental variance, which provides more insights into age-specific heritability of BMI and evidence of G × E interactions. These findings highlight the fundamental importance and implication of the proposed models in facilitating twin studies to investigate the heritability specific to age and other modifying factors.     

A rapid generalized least squares model for a genome-wide quantitative trait association analysis in families

Genome-wide association studies (GWAS) using family data involve association analyses between hundreds of thousands of markers and a trait for a large number of related individuals. The correlations among relatives bring statistical and computational challenges when performing these large-scale association analyses. Recently, several rapid methods accounting for both within- and between-family variation have been proposed. However, these techniques mostly model the phenotypic similarities in terms of genetic relatedness. The familial resemblances in many family-based studies such as twin studies are not only due to the genetic relatedness, but also derive from shared environmental effects and assortative mating. In this paper, we propose 2 generalized least squares (GLS) models for rapid association analysis of family-based GWAS, which accommodate both genetic and environmental contributions to familial resemblance. In our first model, we estimated the joint genetic and environmental variations. In our second model, we estimated the genetic and environmental components separately. Through simulation studies, we demonstrated that our proposed approaches are more powerful and computationally efficient than a number of existing methods are. We show that estimating the residual variance-covariance matrix in the GLS models without SNP effects does not lead to an appreciable bias in the p values as long as the SNP effect is small (i.e. accounting for no more than 1% of trait variance).     

On the twin risk in autism

Autism is considered by many to be the most strongly genetically influenced multifactorial childhood psychiatric disorder. In the absence of any known gene or genes, the main support for this is derived from family and twin studies. Two recent studies (Greenberg et al. 2001; Betancur et al. 2002) suggested that the twinning process itself is an important risk factor in the development of autism. If true, this would have major consequences for the interpretation of twin studies. Both studies compared the number of affected twin pairs among affected sib pairs to expected values in two separate samples of multiplex families and reported a substantial and significant excess of twin pairs. Using data from our epidemiological study in Western Australia, we investigated the possibility of an increased rate of autism in twins. All children born between 1980 and 1995 with autism, Asperger syndrome, or pervasive developmental disorder not otherwise specified (PDD-NOS) were ascertained. Of the 465 children with a diagnosis, 14 were twin births (rate 30.0/1,000) compared to 9,640 children of multiple births out of a total of 386,637 births in Western Australia between 1980 and 1995 (twin rate weighted to number of children with autism or PDD per year 26.3/1,000). These data clearly do not support twinning as a substantial risk factor in the etiology of autism. We demonstrate that the high proportion of twins found in affected-sib-pair studies can be adequately explained by the high ratio of concordance rates in monozygotic (MZ) twins versus siblings and the distribution of family size in the population studied. Our results are in agreement with those of two similar studies by Croen et al. (2002) in California and Hultman et al. (2002) in Sweden.     

Analytic approaches to twin data using structural equation models

The classical twin study is the most popular design in behavioural genetics. It has strong roots in biometrical genetic theory, which allows predictions to be made about the correlations between observed traits of identical and fraternal twins in terms of underlying genetic and environmental components. One can infer the relative importance of these 'latent' factors (model parameters) by structural equation modelling (SEM) of observed covariances of both twin types. SEM programs estimate model parameters by minimising a goodness-of-fit function between observed and predicted covariance matrices, usually by the maximum-likelihood criterion. Likelihood ratio statistics also allow the comparison of fit of different competing models. The program Mx, specifically developed to model genetically sensitive data, is now widely used in twin analyses. The flexibility of Mx allows the modelling of multivariate data to examine the genetic and environmental relations between two or more phenotypes and the modelling to categorical traits under liability-threshold models.     

Bias of heritability estimates from twin data in presence of errors of zygosity determination

Simple heritability estimators of continuous as well as discrete traits from twin data are known to overestimate the degree of genetic determination of the measured traits for several reasons. Errors of zygosity determination will, however, underestimate the true heritability. The bias due to wrong assignment of dizygous twin pairs into monozygous type is evaluated here, and the results indicate that this negative bias has a compensatory effect on the estimate of the degree of genetic determination when other factors of similarity between twin pairs are taken into account. It is shown that when an estimate of zygosity error is available, the bias due to this factor can be evaluated quantitatively, and hence the adjustment for zygosity error can be incorporated in the estimation of the degree of genetic determination of a trait. Although this theory is explicitly developed here for twin studies, the general principle also applies for other types of errors of determining the degree of biological relationships for estimation of heritability, in which case this type of error may be more important than the simple zygosity error.     

Analyzing Clustered Data: Why and How to Account for Multiple Observations Nested within a Study Participant?

A conventional study design among medical and biological experimentalists involves collecting multiple measurements from a study subject. For example, experiments utilizing mouse models in neuroscience often involve collecting multiple neuron measurements per mouse to increase the number of observations without requiring a large number of mice. This leads to a form of statistical dependence referred to as clustering. Inappropriate analyses of clustered data have resulted in several recent critiques of neuroscience research that suggest the bar for statistical analyses within the field is set too low. We compare naïve analytical approaches to marginal, fixed-effect, and mixed-effect models and provide guidelines for when each of these models is most appropriate based on study design. We demonstrate the influence of clustering on a between-mouse treatment effect, a within-mouse treatment effect, and an interaction effect between the two. Our analyses demonstrate that these statistical approaches can give substantially different results, primarily when the analyses include a between-mouse treatment effect. In a novel analysis from a neuroscience perspective, we also refine the mixed-effect approach through the inclusion of an aggregate mouse-level counterpart to a within-mouse (neuron level) treatment as an additional predictor by adapting an advanced modeling technique that has been used in social science research and show that this yields more informative results. Based on these findings, we emphasize the importance of appropriate analyses of clustered data, and we aim for this work to serve as a resource for when one is deciding which approach will work best for a given study.     

Segregation analysis of phenotypic components of learning disabilities. I. Nonword memory and digit span

Dyslexia is a common and complex disorder with evidence for a genetic component. Multiple loci (i.e., quantitative-trait loci [QTLs]) are likely to be involved, but the number is unknown. Diagnosis is complicated by the lack of a standard protocol, and many diagnostic measures have been proposed as understanding of the component processes has evolved. One or more genes may, in turn, influence these measures. To date, little work has been done to evaluate the mode of inheritance of individual component-as opposed to composite-phenotypes, beyond family or twin correlation studies that initially demonstrate evidence for a genetic basis of such components. Here we use two approaches to segregation analysis in 102 nuclear families to estimate genetic models for component phenotypes associated with dyslexia: digit span and a nonword-repetition task. Both measures are related to phonological skills, one of the key component processes in dyslexia. We use oligogenic-trait segregation analysis to estimate the number of QTLs contributing to each phenotype, and we use complex segregation analysis to identify the most parsimonious inheritance models. We provide evidence in support of both a major-gene mode of inheritance for the nonword-repetition task, with approximately 2.4 contributing QTLs, and for a genetic basis of digit span, with approximately 1.9 contributing QTLs. Results obtained by reciprocal adjustment of measures suggest that genes contributing to digit span may contribute to the nonword-repetition score but that there are additional QTLs involved in nonword repetition. Our study adds to existing studies of the genetic basis of composite phenotypes related to dyslexia, by providing evidence for major-gene modes of inheritance of these single-measure component phenotypes.     

Genetics of Callous-Unemotional Behavior in Children

Callous-unemotional behavior (CU) is currently under consideration as a subtyping index for conduct disorder diagnosis. Twin studies routinely estimate the heritability of CU as greater than 50%. It is now possible to estimate genetic influence using DNA alone from samples of unrelated individuals, not relying on the assumptions of the twin method. Here we use this new DNA method (implemented in a software package called Genome-wide Complex Trait Analysis, GCTA) for the first time to estimate genetic influence on CU. We also report the first genome-wide association (GWA) study of CU as a quantitative trait. We compare these DNA results to those from twin analyses using the same measure and the same community sample of 2,930 children rated by their teachers at ages 7, 9 and 12. GCTA estimates of heritability were near zero, even though twin analysis of CU in this sample confirmed the high heritability of CU reported in the literature, and even though GCTA estimates of heritability were substantial for cognitive and anthropological traits in this sample. No significant associations were found in GWA analysis, which, like GCTA, only detects additive effects of common DNA variants. The phrase ‘missing heritability’ was coined to refer to the gap between variance associated with DNA variants identified in GWA studies versus twin study heritability. However, GCTA heritability, not twin study heritability, is the ceiling for GWA studies because both GCTA and GWA are limited to the overall additive effects of common DNA variants, whereas twin studies are not. This GCTA ceiling is very low for CU in our study, despite its high twin study heritability estimate. The gap between GCTA and twin study heritabilities will make it challenging to identify genes responsible for the heritability of CU.     

Conditional Random Fields for Fast, Large-Scale Genome-Wide Association Studies

Understanding the role of genetic variation in human diseases remains an important problem to be solved in genomics. An important component of such variation consist of variations at single sites in DNA, or single nucleotide polymorphisms (SNPs). Typically, the problem of associating particular SNPs to phenotypes has been confounded by hidden factors such as the presence of population structure, family structure or cryptic relatedness in the sample of individuals being analyzed. Such confounding factors lead to a large number of spurious associations and missed associations. Various statistical methods have been proposed to account for such confounding factors such as linear mixed-effect models (LMMs) or methods that adjust data based on a principal components analysis (PCA), but these methods either suffer from low power or cease to be tractable for larger numbers of individuals in the sample. Here we present a statistical model for conducting genome-wide association studies (GWAS) that accounts for such confounding factors. Our method scales in runtime quadratic in the number of individuals being studied with only a modest loss in statistical power as compared to LMM-based and PCA-based methods when testing on synthetic data that was generated from a generalized LMM. Applying our method to both real and synthetic human genotype/phenotype data, we demonstrate the ability of our model to correct for confounding factors while requiring significantly less runtime relative to LMMs. We have implemented methods for fitting these models, which are available at http://www.microsoft.com/science.     

A model of HIV-1 infection with two time delays: mathematical analysis and comparison with patient data

Mathematical models have made considerable contributions to our understanding of HIV dynamics. Introducing time delays to HIV models usually brings challenges to both mathematical analysis of the models and comparison of model predictions with patient data. In this paper, we incorporate two delays, one the time needed for infected cells to produce virions after viral entry and the other the time needed for the adaptive immune response to emerge to control viral replication, into an HIV-1 model. We begin model analysis with proving the positivity and boundedness of the solutions, local stability of the infection-free and infected steady states, and uniform persistence of the system. By developing a few Lyapunov functionals, we obtain conditions ensuring global stability of the steady states. We also fit the model including two delays to viral load data from 10 patients during primary HIV-1 infection and estimate parameter values. Although the delay model provides better fits to patient data (achieving a smaller error between data and modeling prediction) than the one without delays, we could not determine which one is better from the statistical standpoint. This highlights the need of more data sets for model verification and selection when we incorporate time delays into mathematical models to study virus dynamics.     

Estimating two-stage models for genetic influences on alcohol, tobacco or drug use initiation and dependence vulnerability in twin and family data.

Genetic research on risk of alcohol, tobacco or drug dependence must make allowance for the partial overlap of risk-factors for initiation of use, and risk-factors for dependence or other outcomes in users. Except in the extreme cases where genetic and environmental risk-factors for initiation and dependence overlap completely or are uncorrelated, there is no consensus about how best to estimate the magnitude of genetic or environmental correlations between Initiation and Dependence in twin and family data. We explore by computer simulation the biases to estimates of genetic and environmental parameters caused by model misspecification when Initiation can only be defined as a binary variable. For plausible simulated parameter values, the two-stage genetic models that we consider yield estimates of genetic and environmental variances for Dependence that, although biased, are not very discrepant from the true values. However, estimates of genetic (or environmental) correlations between Initiation and Dependence may be seriously biased, and may differ markedly under different two-stage models. Such estimates may have little credibility unless external data favor selection of one particular model. These problems can be avoided if Initiation can be assessed as a multiple-category variable (e.g. never versus early-onset versus later onset user), with at least two categories measurable in users at risk for dependence. Under these conditions, under certain distributional assumptions, recovery of simulated genetic and environmental correlations becomes possible. Illustrative application of the model to Australian twin data on smoking confirmed substantial heritability of smoking persistence (42%) with minimal overlap with genetic influences on initiation.     

Modelling multiple ovulation,fertilization, and embryo loss in human fertility studies

Models of human fertility that incorporate information on timing of intercourse have assumed that a single ovum is released each menstrual cycle. These models are misspecified if two or more viable ova are sometimes released in a single cycle, which is known to occur in dizygotic twin pregnancies. In this paper, we propose a model for multiple ovulation in humans. We assume that the unobservable number of viable ova in each cycle follows a multinomial distribution. Successful fertilization of each ovum depends on the ability of the cycle to support a pregnancy and on the aggregate of a set of unobservable Bernoulli trials representing the fertilizing effects of intercourse on various days. Our model accommodates general covariate effects, allows for heterogeneity among couples, and accounts for a sterile subpopulation of couples. Information on early detection of pregnancy can be incorporated to estimate the probability of embryo loss. We outline a Markov chain Monte Carlo algorithm for estimation of the posterior distributions of the parameters. The methods are applied to data from a North Carolina pregnancy study, and applications to studies of assisted reproduction are described.     

Biometrical modeling of twin and family data using standard mixed model software

Summary .   Biometrical genetic modeling of twin or other family data can be used to decompose the variance of an observed response or 'phenotype' into genetic and environmental components. Convenient parameterizations requiring few random effects are proposed, which allow such models to be estimated using widely available software for linear mixed models (continuous phenotypes) or generalized linear mixed models (categorical phenotypes). We illustrate the proposed approach by modeling family data on the continuous phenotype birth weight and twin data on the dichotomous phenotype depression. The example data sets and commands for Stata and R/S-PLUS are available at the Biometrics website.     

Stability and stabilization for models of chemostats with multiple limiting substrates

We study chemostat models in which multiple species compete for two or more limiting nutrients. First, we consider the case where the nutrient flow and species removal rates and input nutrient concentrations are all given as positive constants. In that case, we use Brouwer degree theory to give conditions guaranteeing that the models admit globally asymptotically stable componentwise positive equilibrium points, from all componentwise positive initial states. Then we use the results to develop stabilization theory for a class of controlled chemostats with two or more limiting nutrients. For cases where the dilution rate and input nutrient concentrations can be selected as controls, we prove that many different componentwise positive equilibria can be made globally asymptotically stable. This extends the existing control results for chemostats with one limiting nutrient. We demonstrate our methods in simulations.     

Stability and stabilization for models of chemostats with multiple limiting substrates

We study chemostat models in which multiple species compete for two or more limiting nutrients. First, we consider the case where the nutrient flow and species removal rates and input nutrient concentrations are all given as positive constants. In that case, we use Brouwer degree theory to give conditions guaranteeing that the models admit globally asymptotically stable componentwise positive equilibrium points, from all componentwise positive initial states. Then we use the results to develop stabilization theory for a class of controlled chemostats with two or more limiting nutrients. For cases where the dilution rate and input nutrient concentrations can be selected as controls, we prove that many different componentwise positive equilibria can be made globally asymptotically stable. This extends the existing control results for chemostats with one limiting nutrient. We demonstrate our methods in simulations.     

Population HIV-1 dynamics in vivo: applicable models and inferential tools for virological data from AIDS clinical trials

In this paper, we introduce a novel application of hierarchical nonlinear mixed-effect models to HIV dynamics. We show that a simple model with a sum of exponentials can give a good fit to the observed clinical data of HIV-1 dynamics (HIV-1 RNA copies) after initiation of potent antiviral treatments and can also be justified by a biological compartment model for the interaction between HIV and its host cells. This kind of model enjoys both biological interpretability and mathematical simplicity after reparameterization and simplification. A model simplification procedure is proposed and illustrated through examples. We interpret and justify various simplified models based on clinical data taken during different phases of viral dynamics during antiviral treatments. We suggest the hierarchical nonlinear mixed-effect model approach for parameter estimation and other statistical inferences. In the context of an AIDS clinical trial involving patients treated with a combination of potent antiviral agents, we show how the models may be used to draw biologically relevant interpretations from repeated HIV-1 RNA measurements and demonstrate the potential use of the models in clinical decision-making.