COMPUTER IMPLEMENTED METHOD AND APPARATUS FOR ANALYZING GENETIC DATA - Patent application

Abstract

A computer-implemented method of analyzing genetic data is provided, the method including receiving a polygenic risk score for a phenotype or phenotype combination of interest for a subject individual, receiving individual genetic data for the subject individual, the individual genetic data having information regarding the ancestry of the subject individual, determining an individual position in ancestry space using the individual genetic data, and calculating a genetic contribution to the subject individual's risk for the phenotype or phenotype combination of interest using the polygenic risk score and the individual position. Corresponding apparatus is also provided.

Description

[Technical field]

[0001] The present invention relates to the analysis of genetic and phenotypic data about an organism to obtain information about the organism, and in particular to enable improved estimation of an organism's risk of having a phenotype or combination of phenotypes of interest based on the organism's ancestry. [Background technology]

[0002] A polygenic risk score (PRS) is a quantitative summary of the contribution of DNA to the phenotypes that an organism may express from its genetic inheritance. A PRS may include in its calculation all DNA variants that are related (directly or indirectly) to the phenotype of interest, or it may use components if they are more relevant to a particular aspect of the organism's biology (including cells, tissues, other biological units, mechanisms, or processes). A PRS can be used directly, or as part of multiple measurements or records of an organism, to infer aspects of the organism's past, present, and future biology.

[0003] PRSs have gained attention as tools for disease prevention, stratification, and diagnosis. In the context of improving human health and medical care, PRSs have a wide range of practical applications, including but not limited to predicting the risk of developing a disease or phenotype, predicting the age of onset of a phenotype, predicting disease severity, predicting disease subtype, predicting response to treatment, selecting appropriate screening strategies for individuals, selecting appropriate drug interventions, and setting prior probabilities for other predictive algorithms.

[0004] PRS can have direct application as an input source in the use of artificial intelligence and machine learning methods to make predictions or classifications from other high-dimensional input data (e.g., imaging). They can be used to help train these algorithms to identify predictive measurements, for example, based on non-genetic data. In addition to being useful for making predictive statements about individuals, they can also be used to identify cohorts of individuals, including but not limited to the above applications, by calculating the PRS of a large number of individuals and then grouping the individuals based on their PRS.

[0005] A PRS can also aid in the selection of individuals for clinical trials, for example optimizing trial design by recruiting individuals who are more likely to develop the relevant disease or phenotype, thereby improving the assessment of the efficacy of new treatments. A PRS has information about the individual for whom it is calculated, but also about its relatives (with whom the individual shares some of the DNA that he or she inherited genetically). Information about the impact of an individual's DNA on its phenotype can be derived from an appropriate assessment of the potential impact of carrying a particular combination of DNA variants.

[0006] In what follows, we focus on the analysis of the wealth of information in recent years that has come from genetic association studies (GAS). These studies systematically assess the potential contribution of DNA variants to the genetic basis of a phenotype. Since the mid-2000s, GAS (genome-wide association studies in general, or association studies targeting single variants or variants in a region of the genome, i.e., GWAS restricted to specific regions of the genome) have been conducted on thousands of (mostly human) phenotypes in millions of individuals, generating billions of potential links between genotype and phenotype. The resulting raw data are often subsequently simplified to obtain summary statistics. GAS summary statistics consist, for each genetic variant (either imputed or observed), of the estimated effect size of the genetic variant on the GAS phenotype and the standard error of the estimated effect size. In other cases, individual-level data are directly available, consisting of the complete genetic profile of the individuals in a study and information on their phenotype. However, individual-level data is typically not widely available due to requirements for privacy of individual data.

[0007] A PRS aggregates a large number of appropriately weighted genetic variants to provide an individual-specific relative risk of a disease such as coronary artery disease or breast cancer. However, mapping from weighted genetic variants to relative risk is not straightforward. One reason is that the analytical strategies used to generate a PRS typically yield uncalibrated scores, i.e., the PRS is not easily interpretable even in the population from which it was derived. Furthermore, the effect size associated with each genetic variant varies as a function of ancestry, making the interpretation of the PRS variable across populations. Thus, the aggregated "effect size per unit of PRS" is not constant across human populations.

[0008] Given these limitations, one practical option is to investigate PRS in a collection of cohorts representing diverse human populations. The effect size of the PRS in each of those cohorts can then be estimated. However, this is not always possible because the relevant data (combining genetic and outcome data) may not exist or because not all individuals fit neatly into well-defined ancestry groups. Individuals of mixed ancestry, i.e., those from small or poorly studied populations, often do not fit into the few commonly used ancestry groups.

[0009] The challenge is extrapolation, i.e., the need to extrapolate PRS effect sizes based on limited case-control or prospective data sets from a particular population to individuals who either come from different populations or do not exactly fit the characteristics of the training set. Existing methods for calculating PRS and relative risk do not adequately take these factors into account, leading to scores that are inaccurate for a large proportion of individuals. Summary of the Invention

[0010] To address these and other limitations, according to a first aspect of the present invention there is provided a computer-implemented method of analysing genetic data comprising: receiving a polygenic risk score for a phenotype or phenotype combination of interest for a subject individual; receiving individual genetic data for the subject individual, the genetic data having information about the ancestry of the subject individual; determining an individual position in ancestry space using the individual genetic data; and calculating a genetic contribution to the subject individual's risk for the phenotype or phenotype combination of interest using the polygenic risk score and the individual position.

[0011] Using an individual's position in ancestry space to determine genetic contribution to risk allows the method to account for individuals who do not fit neatly into one of a small number of predefined ancestries, which can improve risk estimation and the ability to make appropriate interventions based on genetic risk.

[0012] In some embodiments, individual locations are represented by a combination of orderable variables, which allow locations to be consistently ranked or positioned relative to one another to infer how similar they are to a particular predefined ancestor.

[0013] In some embodiments, individual locations are represented by a combination of continuous or pseudo-continuous variables. The use of continuous or pseudo-continuous variables further improves the ability to rank and compare locations.

[0014] In some embodiments, an individual location includes an assignment to one or a weighted combination of multiple ancestors. Using a weighted combination of ancestors allows for improved predictions for individuals of mixed ancestry.

[0015] In some embodiments, the genetic contribution has a continuous or quasi-continuous dependence on individual location. Allowing the genetic contribution to vary at least quasi-continuously with location improves the resolution at which appropriate risk can be assigned to individuals of mixed ancestry.

[0016] In some embodiments, the genetic contribution comprises a sum of partial contributions corresponding to each axis of ancestry space, where each partial contribution is calculated using the coordinates of the individual's position along the respective axis of ancestry space. This allows for maximum flexibility in how genetic contributions vary with position by taking into account position relative to each axis of ancestry space.

[0017] In some embodiments, the ancestry space is non-isotropic, such that the dependence of each partial contribution on the respective coordinates of the individual position varies between partial contributions, further enhancing maximum flexibility in how genetic contributions vary with position in ancestry space.

[0018] In some embodiments, at least two of the partial contributions have dependencies on the respective coordinates of the individual position that are related by a shared prior distribution. Relating dependencies on different axes introduces constraints in determining the dependencies that can reduce the tendency for overfitting.

[0019] In some embodiments, the shared prior distribution is specified such that the dependencies of at least two partial contributions are sampled from the same distribution, which can be advantageous when it is known that the variations in different axes have similar functional forms, but their intrinsic values ​​may vary.

[0020] In some embodiments, the shared prior distribution is determined using training data from multiple training individuals and one or more predefined hyper-parameters. The use of hyper-parameters allows for varying the strength of constraints on fitting depending on available information.

[0021] In some embodiments, each partial contribution comprises the product of the polygenic risk score and the coordinates of the individual's location along the respective axis of ancestry space, which is a simple and effective way of combining the PRS with coordinates in ancestry space.

[0022] In some embodiments, calculating the genetic contribution comprises calculating a distance in ancestry space between the individual location and a reference location in ancestry space and using the distance to calculate the genetic contribution. By reducing variation to rely on a single distance, the risk of extrapolating outside of known regions of ancestry space is reduced.

[0023] In some embodiments, the reference position is the position in the ancestry space of the ancestors used to train the coefficients used to calculate the polygenic risk score. This is generally an appropriate reference point, as it represents the point at which the PRS is likely to have the largest effect size, decreasing in any direction.

[0024] In some embodiments, the ancestry space is defined using reference genetic data from multiple reference individuals with multiple different ancestries, and calculating the genetic contribution comprises scaling each axis of the ancestry space using the variance explained by the respective axis in the reference genetic data before calculating the distance. This scaling means that the most significant axes in terms of the variance they explain will contribute the most to the distance, thereby weighting them more heavily.

[0025] In some embodiments, the distance is the Euclidean distance in ancestry space, which is a simple and easily computed measure of distance.

[0026] In some embodiments, the genetic contribution comprises the product of the polygenic risk score and the distance, which is a simple and effective way of combining the PRS with the distance in ancestry space.

[0027] In some embodiments, calculating the genetic contribution comprises using a linear dependence on the individual location, which provides a simple way to model the dependence with consistent behavior.

[0028] In some embodiments, calculating the genetic contribution includes using a non-linear dependency on individual location, which allows for more complex dependencies that may be appropriate in some circumstances.

[0029] In some embodiments, the nonlinear dependence comprises a regularized function. Using a regularized function ensures that the dependence has a reasonable and smooth shape and helps to avoid overfitting, especially when the data is sparse. For example, the nonlinear dependence comprises a penalized B-spline.

[0030] In some embodiments, the nonlinear dependence is determined using a Gaussian process as a prior distribution for calculating the genetic contribution using Bayesian inference. Gaussian processes are suitable for determining nonlinear dependence, especially when data are obtained from stochastic processes with unknown functional dependence, such as in genetics.

[0031] In some embodiments, the Gaussian process has a mean vector of zero. This simplifies the analysis of the Gaussian process and does not affect the generality of the method since the effect of the mean can be added later.

[0032] In some embodiments, the Gaussian process has a mean vector corresponding to a priori estimates of the genetic contributions to the risk of the subject individual, which allows the method to control the odds ratios in regions of ancestry space away from the reference genetic data to take into account knowledge of the risks of various populations.

[0033] In some embodiments, the Gaussian process kernel function is a stationary function. This choice helps ensure that the variation in genetic contributions will be similar in different parts of ancestry space, as is typically expected, even if their absolute values ​​differ.

[0034] In some embodiments, the kernel function of the Gaussian process decays to zero as the similarity between samples decreases. This helps ensure that the genetic contribution is not distorted in areas of ancestry space where reference genetic data is relatively sparse. For example, in some embodiments, the kernel function is a radial basis function or a rational quadratic covariance function.

[0035] In some embodiments, the posterior distribution for Bayesian inference is determined using a Gaussian process and training data from multiple training individuals with different ancestries, which ties the posterior distribution to the particular dataset being used.

[0036] In some embodiments, determining the posterior distribution includes approximating the posterior distribution as a normal distribution, which helps maintain tractability since some implementations may result in non-normal distributions.

[0037] In some embodiments, the Gaussian process kernel function depends on one or more hyperparameters. Using the hyperparameters to control the kernel function prevents overfitting and enforces a meaningful notion of distance in ancestral space. For example, in some embodiments, the hyperparameters include hyperparameters associated with each of the polygenic risk score, the individual location, and the interaction between the polygenic risk score and location.

[0038] In some embodiments, the genetic contribution includes an ancestry-independent component that is calculated using an individual location-independent polygenic risk score, which allows the PRS to contribute to risk to some degree regardless of the variation due to individual location in ancestry space.

[0039] In some embodiments, the genetic contribution includes an ancestry-dependent component that is calculated based on individual position independent of the polygenic risk score, allowing the risk to take into account increased risk due to ancestry without regard to the genetic variation of other individuals.

[0040] In some embodiments, ancestry space is defined using the reference genetic data from multiple reference individuals with multiple different ancestors.In some embodiments, each reference individual is assigned to one of multiple ancestors.Using the reference genetic data from individuals with different ancestors allows the method to take into account a wide variety of ancestors when determining ancestry space and individual position.

[0041] In some embodiments, the coordinate system of the ancestry space is determined by applying dimensionality reduction to the reference genetic data, which is an efficient technique for characterizing high-dimensional data, such as genetic data, in a more compact and efficient manner. For example, in some embodiments, the dimensionality reduction comprises principal component analysis or independent component analysis.

[0042] In some embodiments, dimensionality reduction involves discretizing the ancestry space into a finite set of ancestors, and the individual location includes a continuous or quasi-continuous membership ratio for each ancestor in the finite set of ancestors, which provides a simple way of expressing individual locations in terms of defined ancestral groups.

[0043] In some embodiments, the coordinate system of the ancestry space is selected to maximize the variance of the reference genetic data explained by the ancestry space, which ensures that the genetic contributions account for as much variation due to ancestry as possible.

[0044] In some embodiments, the ancestry space has a lower dimensionality than the individual genetic data, and determining the individual position comprises projecting the individual genetic data into the ancestry space, which allows the genetic contribution for the new individual to be calculated.

[0045] In some embodiments, the dependence of genetic contributions on individual location and polygenic risk scores is determined using training data from multiple training individuals of multiple different ancestries, where the training data includes, for each training individual, genetic data and whether the training individual has a phenotype or phenotype combination of interest. The training data may or may not be the same as the reference genetic data, depending on the specific information available from various studies.

[0046] In some embodiments, the training data further comprises data for each of the training individuals having information about one or more non-genetic covariates, and the genetic contribution is jointly estimated in the presence of the non-genetic covariates, and the method further comprises receiving individual covariate data for the subject individuals, the individual covariate data having information about additional non-genetic covariates for the subject individuals. This allows the genetic contribution to take into account other factors that may indirectly affect the genetic contribution and that may themselves be correlated with ancestry. For example, in some embodiments, the non-genetic covariates comprise one or more of weight, height, behavioral characteristics, medical traits, and other biomarkers such as blood or urine-based measurements.

[0047] In some embodiments, the method further comprises outputting the genetic contribution to risk, which allows for downstream utilization of the genetic contribution.

[0048] In some embodiments, the risk is a relative risk relative to an individual with an average estimated genetic contribution, and the method further comprises calculating the subject individual's relative risk for the phenotype or phenotype combination of interest using the genetic and non-genetic contributions to the relative risk, and outputting the relative risk. Relative risk is a measure that can be more directly correlated to the likelihood that an individual will develop a phenotype.

[0049] In some embodiments, calculating the relative risk includes determining a value of the relative risk for the individual from a distribution of relative risks for the individual using a loss function. The loss function determines an appropriate selection of a single value from the distribution as the expected value. For example, in some embodiments, the loss function is a mean squared error function or an asymmetric exponential loss function.

[0050] In some embodiments, the method further comprises calculating the subject individual's absolute risk for the phenotype or phenotype combination of interest using the genetic contribution calculated by the method of any preceding claim, and outputting the absolute risk. Absolute risk is another useful measure for expressing the likelihood that an individual will develop a particular phenotype.

[0051] According to a second aspect of the present invention there is provided an apparatus for analysing genetic data comprising a processor configured to: receive a polygenic risk score for a phenotype or phenotype combination of interest for a subject individual; receive individual genetic data for the subject individual, the genetic data having information about the ancestry of the subject individual; determine an individual position in ancestry space using the individual genetic data; and calculate a genetic contribution to the subject individual's risk for the phenotype or phenotype combination of interest using the polygenic risk score and the individual position. The processor may be further configured to perform operations similar to those described with respect to the computer-implemented method above.

[0052] The invention may also be embodied as a computer program comprising instructions to cause a computer to perform the above method, or as a computer readable medium comprising instructions which, when executed by a computer, cause the computer to perform the above method.

[0053] Embodiments of the present invention will now be further described, by way of example only, with reference to the accompanying drawings, in which: [Brief description of the drawings]

[0054] [Figure 1] FIG. 1 shows the estimated odds ratio per unit of PRS according to the prior art method. [Diagram 2] 1 is a flow chart of a method according to one embodiment of the present invention. [Diagram 3] FIG. 13. Training data lacking South Asian individuals projected into ancestry space. [Figure 4] FIG. 4 shows effect size estimates for various ancestries under a linear, absolute position model of genetic contribution using the training data from FIG. 3 . [Diagram 5] FIG. 5 shows maximum likelihood estimates for various terms in the genetic contribution under the same model and conditions as in FIG. 4. [Figure 6] FIG. 5 shows effect size estimates similar to those in FIG. 4, but with genetic contributions using fewer dimensions in ancestry space. [Figure 7] FIG. 7 shows the estimated odds ratios per unit of PRS for various ancestries under a linear relative position model of genetic contribution for the same training data as in FIGS. 3-6. [Figure 8] FIG. 8 shows effect size estimates for various ancestries under the same model as in FIG. 7. [Figure 9] FIG. 9 shows effect size estimates for various ancestries under a hierarchical linear absolute position model of genetic contributions for the same training data as in FIGS. 3-8. [Figure 10] FIG. 10 shows effect size estimates for various ancestries under a nonlinear, absolute position model of genetic contributions for the same training data as in FIGS. 3-9. [Figure 11] FIG. 8 shows effect size estimates for various ancestries under a Gaussian process relative position model of genetic contribution for the same training data as in FIGS. 2-7. DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS

[0055] prior art As mentioned above, mapping weighted genetic variants in the form of PRS to individual relative risks is not straightforward: the challenge is extrapolation, i.e., the need to extrapolate PRS effect sizes based on limited case-control or prospective data sets from a particular population to individuals who either come from different populations or do not exactly fit the characteristics of the training set.

[0056] The simplest approach is to assume that the effect size of the standardized polygenic risk score (PRS) is constant across all populations, or equivalently, that the phenotypic variance explained by the standardized PRS is independent of ancestry. This is obviously incorrect, but may provide a suitable starting point when limited data are available. This results in the following model: Y i ~Bernoulli(π i ) (1) logit(π i )=β0+β PRS X PRS,i Here, Y i is a random variable indicating whether individual i has a disease (1) or not (0), and π i is the probability that individual i falls into this category, and X PRS,i is the PRS of individual i, and β PRS is the effect size of the PRS and β0 is a constant coefficient. Other covariates could be included but are omitted here for clarity. Some existing methods include principal components (PCs, explained in more detail below) or proportion of European ancestry as covariates in this model, but without any interaction with the PRS. Thus, this still yields a single PRS effect size (Amariuta, T., Ishigaki, K., Sugishita, H. et al. 2020 Methods, Eq. 2, Fritsche et al. 2021 Methods, Eq. 1 and Bitarello and Mathieson 2021).

[0057] A second, and more appropriate, approach that is comparable to the state of the art, consists of assigning each individual to one of a set of predefined ancestors and fitting the model of equation (1) on each ancestor. This results in the following model:

[0058]

number

[0059]

number

[0060] This approach is visualized in Figure 1, which shows samples in the test dataset projected onto the first two principal components (PCs) defined by the 1,000 Genomes dataset and colored by the estimated odds ratio (Odds Ratio: OR) associated with one standard deviation of the PRS. For the data shown in Figure 1, the PRS is standardized to have a standard distribution of 1, so the odds ratio per one standard deviation of the PRS simply corresponds to e β PRS k in model (2).

[0061] Figure 1 shows a projection of an individual into ancestry space using PCs defined by the 1000 Genomes dataset. This type of projection into ancestry space is described further below in connection with the present invention, but is not typically done in existing approaches. The projection shown in Figure 1 serves to demonstrate existing approaches in a form that can be easily compared with embodiments of the present invention described below.

[0062] As shown, each individual is rigidly assigned (or hard-called) to one of five predefined ancestry groups: European (EUR), South Asian (SAS), Native American (AMR_NAT), East Asian (EAS), and African (AFR_SS). In other words, each individual is simply assigned to the ancestry group to which it is most similar by some measure. Odds ratios are estimated separately for each hard-called group, which results in discontinuities along the continuous clines of ancestry between Europe and East Asian and Europe and Africa. Similar discontinuities exist along the clines to South Asian and Native American along the axis defined by PC3 and PC4, but are not shown here. Red dots represent cluster centers of samples within each of the five superpopulations in 1,000 Genomes, defined by the mean value of samples assigned to that group after removal of admixture samples and subsequent cluster assignment.

[0063] One previous study demonstrated a decay in predictive performance with distance in ancestry space, as predicted on theoretical grounds (Prive et al. 2021). However, the use of distance in ancestry space has not been associated with the significance of PRS effect sizes or used to correct risk scores based on ancestry. Although some studies have attempted to quantify the degree to which PRS accuracy declines with genetic distance (Prive et al. 2021), no methodology has been developed to account for this phenomenon.

[0064] Inaccurately estimating relative effect sizes results in uncalibrated datasets that either over- or under-estimate the role of genetics. This can be harmful because preventive or diagnostic measures may be applied to individuals who do not qualify. Estimating the correct effect size of the standardized PRS requires appropriate weighting of the genetics of the individuals, and the interpretation of the effect size on this scale is the change in log odds with a one standard deviation change in the PRS. The present invention addresses these issues. To move away from grouping all individuals into separate homogenous groups, the present invention instead allows the PRS effect size to vary with the individual's position in ancestry space and builds a model that uses this continuous definition of ancestry.

[0065] Introduction of the invention FIG. 2 illustrates one embodiment of a computer-implemented method for analyzing genetic data. The method includes receiving S10 a polygenic risk score (PRS) 10 for a phenotype or phenotype combination of interest for a subject individual. As described above, a PRS is a quantitative summary of the effect of an individual's genetic variants on the risk for a particular phenotype or phenotype combination. The PRS 10 may be calculated using any suitable method. Calculation of the PRS 10 may be performed immediately prior to the method, for example by the same computer system. Alternatively, the PRS 10 may be calculated elsewhere at another time and transmitted to a system performing the method. The PRS 10 may relate to any phenotype. For example, the phenotype of interest may be a disease, such as heart disease, cancer, diabetes, or any other disease of interest.

[0066] The method further comprises a step S20 of receiving individual genetic data 20 of the individual of interest. The individual genetic data 20 comprises information regarding the ancestry of the individual of interest. For example, the individual genetic data 20 may comprise data regarding a number of genetic variants known to be indicative of the ancestry of the individual. However, it is not necessary that the individual genetic data 20 comprises information regarding genetic variants that bear information regarding a phenotype of interest or a phenotype combination of interest. Information regarding the genetic disposition of the individual that is directly related to the phenotype of interest is already encoded in the PRS 10.

[0067] Ancestral Space The method further includes determining an individual position in the ancestry space using the individual genetic data S30. The ancestry space may be defined using reference genetic data from a plurality of reference individuals having a plurality of different ancestries. For example, the reference genetic data may be derived from GWAS or a publicly available database such as the 1000 Genomes database mentioned in connection with FIG. 1. The reference genetic data has information about the ancestry of the reference individuals. For example, the reference genetic data may include data for each reference individual regarding at least a plurality of genetic variants known to be indicative of the individual's ancestry.

[0068] The plurality of genetic variants for which the reference genetic data has information may be the same as or substantially overlap with the plurality of genetic variants for which the individual genetic data has information (e.g., at least 50% of the variants are the same between the two populations). Defining the ancestral space using the same genetic variants as used for the individual can improve the accuracy of the placement of the individual in the space.

[0069] Each reference individual may be assigned to one of multiple ancestries. For example, the reference individuals may be assigned to one of the predefined ancestry groups mentioned above, namely, European (EUR), South Asian (SAS), Native American (AMR_NAT), East Asian (EAS), and African (AFR_SS). Depending on the reference genetic data used, more or fewer predefined ancestry groups may be included. Labeling the reference individuals in this manner helps define the regions in ancestry space associated with each ancestry group. Alternatively, only a subset of the reference individuals may be assigned to one of multiple ancestries. This may be preferred when the reference genetic data includes mixed individuals that do not match well to any of the predefined ancestry groups.

[0070] It is preferable to include data from as many reference individuals as possible, and from individuals with a wide range of ancestry, in the reference genetic data, as this helps to better define ancestry space and allows interpolation for the subject individual to be made more reliably over a larger region of ancestry space.

[0071] In some embodiments, the coordinate system of ancestry space is determined by applying dimensionality reduction to reference genetic data. In human genetics, it is known to use dimensionality reduction statistical methodologies to map very high dimensional genetic space (with millions of variants) into a subspace with a smaller number of dimensions. This subspace is ancestry space. The present invention is primarily concerned with using the position in ancestry space to correct the effect size of PRS. However, the choice of ancestry space and how it is derived has downstream effects on the model used to estimate PRS effect size.

[0072] In some embodiments, the dimensionality reduction may include discretizing the ancestry space into a finite set of ancestors, and the individual location includes a continuous or pseudo-continuous membership ratio for each ancestor of the finite set of ancestors. For example, the individual location may include an assignment to one or a weighted combination of multiple ancestors. By allowing the individual location to represent a mixture of different distinct ancestors for an individual, a more accurate estimate of the ancestry of the individual can be obtained compared to the prior art "hard call" approach described in FIG. 1. This in turn allows for a more accurate estimate of the effect size of the PRS for the individual of interest.

[0073] The ancestral space is typically approximated through dimensionality reduction of the reference genetic data to k dimensions. Typically, linear dimensionality reduction techniques are used. For example, dimensionality reduction may include principal component analysis (PCA), independent component analysis, non-negative matrix decomposition, or factor analysis. Non-linear dimensionality reduction techniques may also be used, as long as it is possible to project new samples (i.e., individual genetic data from the subject individual) into the dimensionally reduced subspace that is the ancestral space. For example, the ancestral space may be defined by performing a principal component analysis (PCA) and taking the first k principal components (PCs). Typically, the ancestral space may have 2, 3, or 4 dimensions. However, higher dimensional ancestral spaces may be used, for example, having 5, 6, 7, 8, 9, 10, or more than 10 dimensions.

[0074] For the examples in the remainder of this application, the linear subspace occupied by the first four PCs is used.

[0075] In some embodiments, the coordinate system of the ancestry space is selected to maximize the variance of the reference genetic data explained by the ancestry space. For example, when the dimensionality reduction includes PCA, this may be achieved by selecting the PC that explains the largest variance in the reference genetic data.

[0076] This choice means that the axes of the ancestry space best correspond to variations that reflect ancestry differences. For example, if the reference genetic data includes data on multiple genetic variants for each reference individual, the variance of the reference genetic data can be the variance among genetic variants known to be indicative of an individual's ancestry, or among variants that are seen to vary among ancestry groups in the dataset.

[0077] In some embodiments, the ancestral space has a lower dimensionality than the individual genetic data, and determining the individual location includes projecting the individual genetic data into the ancestral space. Once the ancestral space is determined, the individual location can be determined by projecting the individual genetic data into the ancestral space. The individual location may be represented in the ancestral space in various ways. In some embodiments, the individual location may be represented by a combination of orderable variables, pseudo-continuous variables, or continuous variables. For example, the individual location may include orderable variables, pseudo-continuous variables, or continuous variables corresponding to each axis or dimension of the ancestral space. The type of variables may be the same for all dimensions of the ancestral space, or may differ between the dimensions of the ancestral space.

[0078] It is often useful to understand the PRS10 as a relative risk conditional on location in a genetically defined ancestry space. A PRS is called "centered" if the expectation (mean) of the PRS conditional on a location in ancestry space is zero. A (centered) PRS is called "standardized" if the variance, or population label, of the PRS conditional on a location in ancestry space is one. Although not required, it is preferred that a centered and standardized PRS is used so that the PRS reflects only an individual's relative risk and does not capture differences between ancestry groups.

[0079] A centered and standardized PRS can be considered, due to the central limit theorem, as a normal distribution with mean 0 and variance 1 at any point in the PC map (i.e., ancestry space). Different but related techniques can be used to obtain this standardization, but they are not the subject of this document.

[0080] Calculation of genetic contributions The method further includes calculating S40 a genetic contribution to the subject individual's risk for a phenotype or phenotype combination of interest using the polygenic risk score 10 and the individual location.

[0081] Once an individual's location in ancestry space has been determined, the relationship between the individual's location and the effect size per unit PRS can be used to determine the genetic contribution to risk for that individual. The relationship between individual location, PRS, and genetic contribution to risk can be modeled in a variety of ways. Two broad approaches are considered in this method. As mentioned above, the examples in this application use a linear subspace occupied by the first four PCs defined from the 1000 Genomes Project data to define the ancestry space.

[0082] The first broad approach is the multi-dimensional "absolute position" approach. In such an embodiment, genetic contribution comprises the sum of the partial contributions corresponding to each axis of ancestry space, and each partial contribution is calculated using the coordinates of the individual position along each axis of ancestry space. In other words, the coordinates of the individual position on each axis or dimension of ancestry space all contribute independently to the dependency of genetic contribution on individual position. To obtain the dependency of genetic contribution on individual position and PRS, a model is fitted that includes PRS, each axis of ancestry space (defined by truncated dimensionality reduction of reference genetic data) and their interactions to obtain ancestry-specific estimates of the effect size contribution of PRS to risk.

[0083] In such an embodiment of the "absolute position" approach, the ancestry space may be non-isotropic, such that the dependence of each partial contribution on the respective coordinates of the individual position differs between partial contributions. This allows maximum flexibility in fitting the dependence of genetic contributions to account for different variations in ancestry of individuals of interest, represented by changes in individual position along different axes in ancestry space.

[0084] A second broad approach is a unidimensional "relative position" approach. In such an embodiment, calculating the genetic contribution involves calculating a distance in ancestry space between the individual position and a reference position in ancestry space, and using the distance to calculate the genetic contribution. In other words, the genetic contribution depends only on the absolute distance in ancestry space between the individual position and the reference point. For example, the distance may be a Euclidean distance in ancestry space. The genetic contribution may include the product of the polygenic risk score and the distance.

[0085] In this set of "relative positions" of the model, the absolute position of each sample in ancestry space is not directly incorporated into the dependency of genetic contributions, but rather the relative distance between the individual position of the subject individual and some reference point is used as some function of the axes of the multidimensional ancestry space. In some embodiments, the reference positions are the positions in ancestry space of the ancestors used to train the coefficients used to calculate the polygenic risk score. In other words, the reference positions correspond to the ancestors of the training data used to generate the PRS, i.e., to train (e.g., by a machine learning algorithm) the set of coefficients used to obtain the PRS of the subject individual based on its individual genetic data.

[0086] This "relative position" approach is based on the observation that the variance explained by a PRS decays monotonically with ancestral distance from the population used to fit the PRS. This relative distance and its interaction with the PRS are then included as covariates when fitting a model of the dependence of genetic contributions. This allows us to estimate PRS effect sizes that vary continuously with ancestral distance from the population used to train the PRS.

[0087] In some embodiments, the ancestry space is defined using reference genetic data from multiple reference individuals with multiple different ancestries, and calculating the genetic contribution comprises scaling each axis of the ancestry space using the variance explained by the respective axis in the reference genetic data before calculating the distances. When using a relative position approach, scaling the axes of the ancestry space by the variance they explain ensures that more significant coordinates are weighted more heavily in the relative position.

[0088] In either the first (multidimensional) or second (unidimensional) approach, the genetic contribution may have a continuous or pseudo-continuous dependence on individual location. Furthermore, both the multidimensional and unidimensional approaches consider both linear and nonlinear dependence of the genetic contribution on PRS and individual location. Specifically, a logistic regression-based estimation is considered, which implicitly assumes a linear dependence between individual location and the logarithmic change of odds ratio, and a nonlinear-based estimation of ancestry-specific PRS using Gaussian processes and general additive models (GAMs).

[0089] Linear Models In some embodiments, calculating the genetic contribution involves using a linear dependence on individual location, which, as described above, can capture individual location either as a multi-dimensional absolute location or as a single-dimensional relative location.

[0090] Absolute Position In some embodiments, the genetic contribution comprises a sum of partial contributions corresponding to each axis of ancestry space, where each partial contribution is calculated using the coordinates of the individual's position along the respective axis of ancestry space. A simple way to incorporate an estimate of the ancestry-dependent effect size PRS is to include the individual's PRS, individual position (as a principal component in this example), and their interaction as covariates in the following logistic model:

[0091]

number

[0092]

number

[0093]

number

[0094]

number

[0095] In this model, we assume some total effect size in the absence of ancestry information, which is then modified by the individual's location in ancestry space. If the PRS is standardized to have mean 0 and variance 1 regardless of ancestry location, then the ancestry-specific PRS effect size β for individual i is PRS,i can be estimated, which is the effect size equivalent to one standard deviation of the PRS:

[0096]

number

[0097]

number

[0098]

number

[0099] Model (3) is simple and extremely fast, and offers a continuous or quasi-continuous extension to prior art methods that make discrete assignments to ancestry categories. Because the underlying model of potential risk is linear, the interpretation of one unit of the standardized PRS does not vary with position on the PRS scale.

[0100] In this and other embodiments, the dependence of genetic contributions on individual location and polygenic risk scores is determined using training data from a plurality of training individuals of different ancestries, the training data including, for each training individual, genetic data and whether the training individual has a phenotype or phenotype combination of interest.

[0101] It is worth noting that the training data used to train the model of the dependency of genetic contribution on PRS and individual location can be different from the reference genetic data used to define ancestry space.For example, the reference genetic data can only include information about genetic variants that are relevant to determining genetic ancestry, and do not include enough data related to the phenotype of interest to also be used to train the model.In some embodiments, the reference genetic data and the training data can be the same, but this is not always the case.

[0102] Care must be taken to interpolate rather than extrapolate, as extrapolating to regions of ancestry space not represented in the training data can lead to spurious results. If a superpopulation (e.g., substantially corresponding to one of the predefined ancestral groups) is not present in the training data but is present in the reference genetic data (which are assumed to capture global genetic diversity), then a single PC will distinguish that superpopulation from the others. Moreover, the remaining individuals from the other superpopulations are likely not to differ significantly along the axes of ancestry space defined by that PC. This results in instability of effect size estimates for that single PC, making it prone to extreme estimates.

[0103] For example, Figure 3 shows an example of training data projected onto the PC3 / PC4 plane of ancestry space (defined by the first four PCs from the 1000 Genomes), where the training data lacks individuals of SAS ancestry. The solid red dots and labels in Figure 3 indicate the centroids of each superpopulation from the 1000 Genomes data. As can be seen, the training completely lacks samples around the SAS centroids.

[0104] Figure 4 shows effect size estimates (estimated effect size per unit PRS) for individuals at the mean position of each of the five major super-populations in 1000 Genomes under model (3), which was trained using the training data in Figure 3. Effect size estimates are reported for each super-population using the first four PCs in a logistic regression (blue) and compared to the discrete estimates within that population in the test (red). Discrete estimates are estimates calculated using only the test data from that super-population, rather than using the full test data from all super-populations with a correction appropriate for that position in ancestry space. The blue circle ("1KG") refers to the odds ratio at the median PC position of the super-population in 1000 Genomes, and the blue triangle ("Empirical") refers to the odds ratio at the mean PC position of the test set individuals in a given super-population. The blue square point estimates ("Group Mean") represent the average effect size estimate across the training set individuals in each super-population label. Due to the missing SAS individual in the training data, the estimated effect size for individuals in the SAS position is below the minimum bound of 1.0 that could be reasonably expected.

[0105] Figure 5 shows the interaction coefficients (Model (3)

[0106]

number

[0107] To prevent this kind of unstable estimation, several modifications can be considered. One approach is to determine the PCs that separate superpopulations not present in the training data from other superpopulations in ancestry space, and then incorporate these PCs into model (3) as covariates (and their corresponding interaction coefficients with the PRS).

[0108]

number

[0109] This correction technique is illustrated in Figure 6, continuing to use the same data used in the previous examples in Figures 3-5. Figure 6 shows effect size estimates for individuals at the mean location of each superpopulation under model (3) in an ancestry space defined by only three PCs. Effect size estimates are reported for each of the five major superpopulations in 1000 Genomes using three PCs in a logistic regression according to model (3) after removing the PCs that distinguish South Asia from the remaining superpopulations in the model using the first four PCs (blue) and compared to the discrete estimates within that population in the test (red). The blue circle ("1KG") refers to the odds ratio at the median PC location of the superpopulation in 1000 Genomes, and the blue triangle ("Empirical") refers to the odds ratio at the mean PC location of the test set individuals in a given superpopulation. The blue square point estimates ("Group Mean") represent the average effect size estimate across individuals in each superpopulation label. As can be seen, the odds ratios for the SAS are higher than the discrete estimates in the superpopulation in the study.

[0110] An alternative approach is to introduce a regularization term into the model of the dependence of genetic contributions to prevent overfitting, for example using LASSO, elastic net, or ridge regression.

[0111] Relative Position In some embodiments, calculating the genetic contribution comprises calculating a distance in ancestry space between the individual position and a reference position in ancestry space, and using the distance to calculate the genetic contribution. In this case, the effect size can be modeled as continuously varying with ancestry using the following model: Y i ~Bernoulli(π i ) (4) logit(π i )=β0+β PRS X PRS,i +β EUR X EUR,i +β PRS×EUR X PRS,i X EUR,i Where X EUR,i is the genetic distance from individual i to the median position of the EUR superpopulation in ancestry space. In the following, this distance is defined as the Euclidean distance after scaling the first four PCs by the variance they explain in the reference genetic data. In other words, the ancestry space is defined using reference genetic data from multiple reference individuals with multiple different ancestries, and calculating the genetic contribution involves scaling each axis of the ancestry space using the variance explained by the respective axis in the reference genetic data before calculating the distance. However, the method is not limited to these choices, and more or fewer PCs can be included, and non-Euclidean distance measures can be used.

[0112] FIG. 7 shows the same exemplary training data used in FIGS. 3-6. Here, the data is projected onto a plane in ancestry space defined by PC1 and PC2. Each point is a training individual in the training data, colored by its estimated PRS effect size ("odds ratio") according to model (4). The linear model of the present method has no discontinuities, unlike the prior art approach shown in FIG.

[0113] Model (4) results in close agreement with the SAS effect size estimates for the test set, as shown in Figure 8. In Figure 8, effect size estimates are reported for each of the five major superpopulations in 1,000 Genomes using three PCs in a logistic regression according to model (4) and compared to the discrete estimates within that population in the test (red). The blue circle ("1KG") refers to the odds ratio at the median PC location for the superpopulation in 1000 Genomes, and the blue triangle ("Empirical") refers to the odds ratio at the mean PC location for individuals in the test set in a given superpopulation. The blue square point estimates ("Group Mean") represent the average effect size estimate across individuals in each superpopulation label.

[0114] However, study estimates may differ due to large uncertainties resulting from cohort / environment effects or small sample sizes (as in the case of AMR_NAT). Note that effect size estimates for EAS are slightly lower compared to the previous approach, and AFR_SS are slightly higher, which may reflect the lower flexibility of this model.

[0115] Model (4) (and other similar "relative position" models) assume that PRS effect sizes change in the same way as individual positions change, regardless of the axis in ancestry space along which the individual positions move. This is a reasonable assumption, since it has been previously observed that discrimination / precision of predictions falls off approximately linearly with distance in ancestry space. In addition, similar effect sizes are observed when binning individuals by distance.

[0116] Model (4) has the same advantages as the multidimensional linear model. However, it has the additional benefit that new individuals are much less likely to have locations in the extrapolation range (i.e., outside the area covered by the training data) because the ancestral space is effectively reduced to a single dimension for the purposes of modeling genetic contributions. Rather, the model interpolates between the training individuals in the individual training data, resulting in more stable estimates of the model parameters (i.e., the interaction coefficients). This effectively regularizes the interaction coefficients.

[0117] A limitation of model (4) is that it makes additional assumptions compared to the multidimensional model (exemplified by model (3)). This means that this model is less flexible and has fewer degrees of freedom. It is possible that the assumption that PRS effect sizes vary the same way in all directions in ancestry space may be inaccurate. In this case, the relative position model will underfit the training data, resulting in biased estimates. In addition, the relative position model also requires the specification of a reference point. This model has similar extrapolation problems when African samples are not present in the training data used to fit the model, but is likely to be more robust than the multidimensional model to the absence of such training data.

[0118] Hierarchical Model Models (3) and (4) can be seen as two extremes of a spectrum. At one end of the spectrum (model (3)), each PC interaction coefficient is estimated without sharing information between the axes of variation in ancestry space. In this case, the partial contributions corresponding to each axis of ancestry space are independent of each other. At the other end of the spectrum (model (4)), the axes of ancestry space are collapsed into a single distance. In effect, this distance shares information completely across PCs.

[0119] Between these two extremes are "hierarchical" models, which share information across different interaction terms and PCs while retaining the flexibility to allow different coefficients for each. This effectively provides some degrees of freedom between 1 (corresponding to a "relative distance" model) and the number of PCs (corresponding to an "absolute distance" model). In such an embodiment, at least two of the partial contributions have a dependency on the respective coordinates of the individual position that is related by a shared prior distribution. These models use a shared prior distribution across the interaction coefficients. This allows for the sharing of information and also represents a prior belief that the PRS effect size should change in the same way along the axis in ancestral space corresponding to each PC, while allowing that the change need not be exactly the same along each axis. In some embodiments, the model can learn from the training data how similar the changes along each axis are to each other.

[0120] Relating the dependence of partial contributions through a shared prior has some similarity to standard hierarchical linear models (also known as multilevel models) that share information between groups of "individuals". However, in this case, as the model moves into a continuous ancestry space, the information is shared across the interaction coefficients themselves (which represent the dependence of partial contributions on each coordinate).

[0121] The present hierarchical model also differs from standard hierarchical models in that, since each interaction coefficient has the same number of data points, the model relies on the different uncertainties of the interaction coefficients to inform how much each interaction coefficient should be regularized.

[0122] For PCs that distinguish superpopulations that are present in the reference genetic data used to define the ancestral space but missing in the training data, estimates are more uncertain due to the smaller range over which the gradients are trained. As a result, these estimates are more strongly regularized.

[0123] This model can be expressed as:

[0124]

number

[0125]

number

[0126] For model (3), each partial contribution involves the product of a polygenic risk score and the coordinates of the individual's position along the respective axis in ancestry space, but here the partial contributions are related to each other by a shared prior distribution.

[0127] For other models, the data and examples given in the figures use only the first four principal components, but in general any number can be included. Increasing this number to the total number of signal PCs in the reference gene data increases the dispersion parameter σ PRS×PC This can assist in learning.

[0128] Since the signs of the principal components are arbitrary, they need to be matched to the training data, otherwise the learned σ PRS×PC would be artificially inflated, for example because one coefficient is -0.1 and the other is +0.1. This can result in a loss of information sharing and failure to regularize sufficiently. Harmonization may be done by changing the sign of the principal components (e.g., by multiplying by minus 1) so that the signs of all learned interaction coefficients are the same. Since the signs of the principal components are arbitrary, it does not matter which sign (positive or negative) is chosen, as long as everything is consistent.

[0129]

number

[0130]

number

[0131] The shared prior distribution is determined using training data from multiple training individuals and one or more predefined hyperparameters. The strength of the prior belief that the interaction coefficients should be similar can be varied by varying the hyperparameter. In the particular implementation exemplified above, the hyperparameter d can be varied, with smaller values ​​of d corresponding to a stronger belief that they are similar.

[0132] Interaction coefficient

[0133]

number

[0134]

number

[0135]

number

[0136] Then, one can derive a point estimate for each individual given some loss function, typically the mean squared error, although the loss function could be another, e.g. an asymmetric exponential loss, if overestimation is considered worse than underestimation, or vice versa.

[0137] This model uses the model parameters, especially the grand σ PRS×PC , but this model can also be constructed in a "centered" version, which is equivalent and is shown as the suggested model.

[0138]

number

[0139]

number

[0140] The parameters in this Bayesian model (5) are fitted using Hamiltonian Monte Carlo to obtain the example results shown in the figure. However, other approximations can be used, such as variational inference or the integrated nested Laplace approximation. When fitted using Hamiltonian Monte Carlo, the model is slower than other linear methods, but faster than the variational Gaussian process, which is described further below.

[0141] Figure 9 shows effect size estimates for individuals at the mean position of each superpopulation under model (5). Effect size estimates are reported for each of the five major superpopulations in the 1,000 Genomes using model (5) (blue) and compared to the discrete estimates within that population in the study (red). The blue circle ("1KG") refers to the odds ratio at the median PC position of the superpopulation in the 1000 Genomes, and the blue triangle ("Empirical") refers to the odds ratio at the mean PC position of the test set individuals in the given superpopulation. The blue square point estimates ("Group Mean") represent the average effect size estimate across individuals in each superpopulation label.

[0142] Model (5), like the unidimensional models above, regularizes well enough to produce SAS effect sizes that are very close to the discrete estimates in the test set. But model (5) is also flexible enough to have EAS / AFR_SS effect size estimates that are similar to the discrete estimates in the test set, like other multidimensional models. Model (5) strikes a balance between unidimensional (relative position) and multidimensional (absolute position) models by learning how similar the interaction coefficients are to the training data, rather than hard-coding identical effects or allowing effects to be "free". If the training data suggests that the dependence of partial contributions along different axes of ancestry space is indeed different, then a large standard deviation will be learned and the dependence will not be regularized as strongly.

[0143] Nonlinear Models All models described so far have been linear models, where the dependence of the genetic contribution (or its partial contribution) on the individual position in ancestry space has been linear. This means that a change of a certain magnitude and direction in the individual position has the same effect on the genetic contribution regardless of where it occurs in ancestry space. However, non-linear models can also be considered. In such an embodiment, calculating the genetic contribution involves using a non-linear dependence on the individual position.

[0144] Generalized Additive Model (GAM) The first example of a nonlinear model is the generalized additive model (GAM). In the case of a GAM, instead of having a simple linear dependence between each coefficient and a covariate (PRS or individual location), this dependence can be modeled as an arbitrary function of each covariate: Y i ~Bernoulli(p i |X) (6) logit(p i |X)=β0+f1(X1)+f2(X2)+···+f M (X N )

[0145] The nonlinear dependence may include a regularized function. This ensures that the dependence is smooth and does not include sharp discontinuities or other features that are unlikely to be biologically plausible. A common functional form to choose is a penalized B-spline (also called a P-spline), which allows for modeling of nonlinear forms. In such an embodiment, the nonlinear dependence includes a penalized B-spline. In a penalized B-spline, a penalty ("lambda") is imposed on the curvature of the function, making the function perfectly linear as the value gets larger. This lambda value may be learned using cross-validation over a grid of values. One way to choose a preferred lambda value is to choose the value that minimizes the unbiased risk estimator.

[0146] As an example of this type of model, Figure 10 shows effect size estimates for individuals at the mean position of each superpopulation under model (6). A unidimensional relative position model is used with a single distance from EUR and its interaction with PRS. Effect size estimates are reported for each of the five major superpopulations in the 1,000 Genomes (blue) using model (6) and compared to the discrete estimates within that population in the test (red). The blue circle ("1KG") refers to the odds ratio at the median PC position of the superpopulation in the 1,000 Genomes, and the blue triangle ("Empirical") refers to the odds ratio at the mean PC position of the test set individuals in a given superpopulation. The blue square points ("Group Mean") represent the average effect size estimate across the test set individuals in each superpopulation label.

[0147] The results using model (6) are very similar to those of the unidimensional model shown in Figure 8. Model (6) is faster than the variational Gaussian process strategy described below, while still allowing for flexible nonlinear dependencies. However, the grid search over different lambda values ​​to fit the penalized B-spline example in Figure 10 is slow, which limits the grid spacing that can be used to test different lambda values.

[0148] Gaussian Process An alternative method for fitting the non-linear dependence is to use a Gaussian process, in such an embodiment, the non-linear dependence is determined using a Gaussian process as a prior distribution for calculating the genetic contribution using Bayesian inference.

[0149] For Gaussian process regression, the following setup is assumed: the training data is a set of noise-added vectors at positions [x1,x2,...,x n ] from some unknown function f:X→R in the set of observations [y1,y2,...,y n ] collection: y i =f(x i )+ε i where i=1,2,...,n Here, ε iis a noise term.

[0150] From this training data, we estimate the unknown function f for the incoming new data: y * =f(x * )+ε * is determined. To do this, it is assumed that the data is sampled from a stochastic process known as a Gaussian Process (GP). Samples taken from a GP have functions μ(.) and k(.,.) such that expectation function E[f(x)]=μ(x) and covariance function Cov[f(x),f(x')]=k(x,x'). A simple notation for this sampling procedure is GP(μ(.),k(·,·)).

[0151] Each finite subset of the GP (specifically, the subset of the data used for training) is sampled from a multivariate normal distribution. How to interpolate between the training data points is determined by the choice of functions μ(.) and k(.,.). In this situation, there is an additional complication that Y is not continuous but discrete case / control labels. We assume that Y is Bernoulli distributed and that the underlying latent distribution governing the risk of being a case is sampled from the GP. Y i ~Bernoulli(π i ) (7) logit(π i )=GP(μ(.),k(·,·))

[0152] Also, the standard method for evaluating the predicted mean and variance is computationally expensive, O(n 3 ), where n is the number of data points. This is not an issue when n is in the thousands, but makes standard techniques intractable when sample sizes approach hundreds of thousands, as can be the case for genetic data.

[0153] To estimate the effect size per PRS for an individual, we use the posterior distribution. The posterior distribution for Bayesian inference is determined using a Gaussian process and training data from multiple training individuals with different ancestries.

[0154] One of the problems with the introduction of a link function (Bernoulli function) that relates the latent space to the case-control state Y is that the posterior distribution is no longer normal. To maintain tractability, in some embodiments, the non-Gaussian posterior distribution is approximated using a normal distribution. Thus, determining the posterior distribution involves approximating the posterior distribution as a normal distribution. There are a series of methods that can be used to make this approximation, for example, Laplace approximation, expectation propagation, and Kullback-Leibler (KL) divergence minimization. In the present embodiment, the latter technique is used by approximating the exact posterior distribution p by a normal distribution q and minimizing KL(q||p). This minimization problem can be solved using Newton's method, which also takes O(n 3 ) The full details of the variational approximation are given in (Nickisch and Rasmussen 2008).

[0155] When using a Gaussian process, it is common to assume a zero mean function μ(.)=0, such that the Gaussian process is given by GP(0,k(·,·). This choice simplifies fitting, and because the shape of a Gaussian process is determined entirely by its covariance function, the mean can always be added back in later. In such embodiments, the Gaussian process has a zero mean function. However, in some cases it may be desirable to include a non-zero mean function, either before fitting or added back in after fitting. For example, a non-zero mean function can be used to control the per-unit OR (effect size) of a PRS in a region of ancestry space away from the training data, essentially as a prior distribution for the PRS contribution in the absence of genetic data at that ancestral position. For example, the effect sizes seen for individuals of African ancestry can be set as the smallest effect sizes across ancestry space, based on the fact that they represent a "worst case scenario" for ancestry decay. Thus, in some embodiments, the Gaussian process has a mean vector that corresponds to a prior estimate of the genetic contribution to risk for the individual of interest. This approach may make fitting more difficult due to the non-linear nature of GPs.

[0156] Most important to the behavior of a Gaussian process is the choice of the kernel function k(x,x'), which describes the functional form for the expected covariance between any pair of locations x and x'.

[0157] Optionally, the kernel function of the Gaussian process is a stationary function. A stationary function depends on the distance between two points, not on their absolute positions. This means that the behavior of the Gaussian process will be similar across the entire ancestral space. The kernel function may additionally be isotropic, in which case it depends only on the magnitude of the distance between two points. Optionally, the kernel function of the Gaussian process decays to zero as the similarity between samples decreases. Selecting a decaying kernel function ensures that the value of the function decays back to the mean value in regions of the ancestral space that are far from any training data. This can help control and reduce unexpected behavior of the function when the target individual has an individual position in the ancestral space that is not in the vicinity of the training data. For example, the kernel function may be a radial basis function or a rational quadratic covariance function. In the following exemplary embodiment, a radial basis function (RBF) is used as the kernel.

[0158] Optionally, the Gaussian process kernel function depends on one or more hyperparameters. For example, the hyperparameters may include hyperparameters related to each of the polygenic risk scores, the individual location, and the interaction between the polygenic risk scores and the location. In the following example, one variance hyperparameter θ, and three length scale hyperparameters l PRS , l PC , l PRS×PC are used for PRS, PC, and PRS × PC interactions, respectively.

[0159] The hyperparameters may include hyperparameters related to the coordinates of the individual location along each axis of the ancestry space, i.e., corresponding to each PC. The hyperparameters may include hyperparameters related to the interaction between the PRS and each coordinate of the individual location. However, in this example, a single length scale hyperparameter is used for the individual location (PC) and a single length scale hyperparameter is used for the interaction between the PRS and the individual location (PRS×PC). This is achieved by scaling the relative length scale of the PC with the variance that the PC explains in the reference genetic data. Thus, the ancestry space is defined using reference genetic data from multiple reference individuals with multiple different ancestries (as described above), and calculating the genetic contribution includes scaling each axis of the ancestry space using the variance explained by the respective axis in the reference genetic data before calculating the distance. This scaling and use of a reduced number of hyperparameters improves speed, prevents overfitting, and enforces a robust notion of "distance" in the ancestry space.

[0160] As with the other models described above, Gaussian processes can be used with either a multi-dimensional approach to absolute position or a uni-dimensional approach to relative position.

[0161] Absolute Position In the absolute position model, the genetic contribution comprises the sum of the partial contributions corresponding to each axis of ancestry space, where each partial contribution is calculated using the coordinates of the individual's position along the respective axis of ancestry space. To apply the Gaussian process framework to the prediction problem in this case, the samples corresponding to each of the training individuals are defined by the vector: x i =[x i,PRS ,x i,PC ,x i,PC×PRS ] where x i,PRS , x i,PC , and x i,PRS×PCare respectively the standardized PRS, the position in ancestral space, and the corresponding interaction for individual i. Then, the kernel function using radial basis functions is given by:

[0162]

number

[0163] Appealing to Gaussian process categorization allows the nonlinear effects of PCs to be naturally incorporated into the prediction of individual case status. Furthermore, depending on the functional form of the kernel, Gaussian processes can enforce that the effect size of the PRS goes to zero away from the training data. This is a desirable property; that is, in the absence of information (prescribed by the kernel), the estimate is expected to be the prevalence of the disease in the training data.

[0164] The benefits of using Gaussian process categorization are also its drawbacks. In this framework, changes in the PRS affect the case probability nonlinearly. Thus, the impact of the PRS on risk varies as a function of its value (unless a linear kernel is enforced for the terms that contain the PRS). Fitting a GP can also be slow. Speed ​​can be improved by using sparse variational Gaussian processes or by making an approximation via dimensionality reduction. The former approach can be difficult to optimize robustly. However, recent papers have proposed techniques to improve the robustness of Gaussian process regression, but not yet for categorization. The radial basis function kernel is a good choice because the trained length scale hyperparameter results in estimates that quickly revert to baseline levels in the absence of data.

[0165] Relative Position In the relative position approach, calculating the genetic contribution involves calculating the distance in ancestral space between the individual position and a reference position in ancestral space, and using that distance to calculate the genetic contribution. In this unidimensional implementation of the GP framework, a different set of predictors is used. The predictors in model (4), namely: Y i ~Bernoulli(π i ) (4) logit(π i )=β0+β PRS X PRS,i +β EUR X EUR,i +β PRS×EUR X PRS,i X EUR,i Instead, a Gaussian process with zero mean function and the following kernel function is used according to model (7):

[0166]

number

[0167] As with the multidimensional embodiment above, four hyperparameters are used in the kernel function. In this embodiment, these are θ, l PRS , l EUR and l PRS×EUR These hyperparameters govern the global scaling and length scale by which information is shared between pairs of points about the PRS, the distance from a reference point (in this example, the Euclidean distance from the European cluster center, determined as the median European location in the reference genetic data), and the interaction between the PRS and distance, respectively.

[0168] Figure 11 shows effect size estimates for individuals at the mean position of each superpopulation. Effect size estimates are derived from the GP framework using model (7), a zero mean function, and the kernel of (8). Effect size estimates are reported for each of the five major superpopulations in 1,000 Genomes using a single distance from the European cluster center (as defined above) and its interaction with the PRS in the GP framework. The fitted estimate (blue) at the centroid of the superpopulation is compared to the discrete estimate (red) within that population in the test. The blue circle ("1KG") refers to the odds ratio at the median PC position of the superpopulation in 1000 Genomes, and the blue triangle ("Empirical") refers to the odds ratio at the mean PC position of the test set individuals in a given superpopulation. The blue square point estimate ("Group Mean") represents the average effect size estimate across individuals in each superpopulation label.

[0169] This GP implementation of the unidimensional relative position model allows for nonlinear changes in PRS effect size as a function of distance from a reference point. That is, it captures the empirically observed relationship between heritability explained by PRS and genomic distance from a reference point, but does not constrain this relationship to be linear. Like other unidimensional (relative position) ancestry space models considered herein, this model assumes that the change in heritability explained by PRS is independent of the direction of movement from the reference point in ancestry space. The drawback of this is that this assumption may have a greater influence than necessary on effect size estimates for areas of ancestry space where training data are available, and may provide overly confident effect size estimates in areas where little training data is available.

[0170] Contributions that are not dependent on both PRS and position All the exemplary models presented herein can be expanded. For example, additional PRS trained with different or similar ancestry can be incorporated by adding new standardized PRS and interaction terms to the model. The aggregate genetic contribution of these PRS to the risk of the subject individual can then be evaluated.

[0171] The component of the genetic contribution (or its partial contribution) that depends on both the PRS and individual location is generally referred to above as the interaction term. In addition, as shown in the various models above, the genetic contribution can include a component that is not dependent on both the PRS and individual location. The genetic contribution may include an ancestry-independent component that is calculated using a polygenic risk score that is independent of individual location. For example, the ancestry-independent component is β in model (4) above. PRS X PRS,i Similarly, the genetic contribution may include an ancestry-dependent component that is calculated based on the individual position independent of the polygenic risk score. For example, the ancestry-dependent component is expressed as β EUR X EUR,i It is written as follows.

[0172] Other features As described above, the dependence of genetic contributions on individual location and polygenic risk scores is determined using training data from a plurality of training individuals having a plurality of different ancestries. The training data includes, for each of the training individuals, genetic data and whether the training individual has a phenotype or phenotype combination of interest. The training data may be different from the reference genetic data used to define the ancestry space.

[0173] Integrated Risk Tools In some embodiments, the training data further comprises data for each of the training individuals having information about one or more non-genetic covariates, and the genetic contribution is jointly estimated in the presence of the non-genetic covariates. For example, the non-genetic covariates may comprise one or more of weight, height, behavioral characteristics, medical traits, and other biomarkers, such as blood or urine-based measurements.

[0174] This allows the model to take into account the correlation between genetic covariates and non-genetic covariates when determining the dependency of genetic contribution, which can further improve accuracy. In such an embodiment, the method further comprises receiving individual covariate data of the subject individual, the individual covariate data having information about additional non-genetic covariates for the subject individual. Calculating the genetic contribution may further comprise using this individual covariate data. Alternatively or additionally, the individual covariate data may be used when calculating the non-genetic contribution to the risk of the subject individual.

[0175] output Once calculated, the genetic contribution to risk may be output. In some embodiments, the method further comprises outputting the genetic contribution to risk. Alternatively or additionally, the genetic contribution may be used as part of further calculations used to calculate other useful indices for assessing an individual's risk of expressing a particular phenotype or phenotype combination.

[0176] In some embodiments, the risk is a relative risk to an individual with an average estimated genetic contribution. The average estimated genetic contribution may be based on average genetic data of similar individuals, e.g., individuals with similar ancestry. The average estimated genetic contribution may be based on the average prevalence of a phenotype or phenotype combination in a population or similar individuals. In such an embodiment, the method further comprises calculating S50 the subject individual's relative risk for a phenotype or phenotype combination of interest using the genetic and non-genetic contributions to the relative risk calculated using any of the methods described in any preceding claim above, and outputting S60 the relative risk. As mentioned above, the non-genetic contribution may be calculated based on non-genetic covariates.

[0177] Calculating the relative risk may include determining a value of the relative risk for a subject from a distribution of the subject's relative risks using a loss function, for example, the loss function may be a mean squared error function or an asymmetric exponential loss function.

[0178] In addition to relative risk, the method may be used to calculate absolute risk. In some embodiments, the method further comprises calculating the subject individual's absolute risk for a phenotype or phenotype combination of interest using the genetic contribution calculated using the method of any preceding claim, and outputting the absolute risk. [Explanation of symbols]

[0179] 10 Polygenic Risk Score (PRS) 20 Individual genetic data 30 Non-genetic contributions 40 Relative Risk

Claims

1. 1. A computer-implemented method for analyzing genetic data, comprising: receiving a polygenic risk score for a phenotype or combination of phenotypes of interest for the subject individual; receiving individual genetic data for the individual of interest, the individual genetic data comprising information regarding the ancestry of the individual of interest; determining an individual position in ancestry space using the individual genetic data, wherein the individual position is represented by a combination of continuous or pseudo-continuous variables or comprises an assignment to one or a weighted combination of multiple ancestors; and calculating a genetic contribution to risk for the subject individual for the phenotype or phenotype combination of interest using the polygenic risk score and the individual location, wherein the genetic contribution has a continuous or quasi-continuous dependence on the individual location; The ancestry space is defined using reference genetic data from a plurality of reference individuals with a plurality of different ancestries. The computer-implemented method.

2. 2. The method of claim 1, wherein the genetic contribution comprises a sum of partial contributions corresponding to each axis of the ancestry space, each partial contribution being calculated using coordinates of the individual's position along a respective axis of the ancestry space.

3. The method of claim 2 , wherein the ancestral space is anisotropic, such that the dependence of each partial contribution on the respective coordinates of the individual positions differs between the partial contributions.

4. The method of claim 3 , wherein at least two of the partial contributions have dependencies on respective coordinates of the individual positions that are related by a shared prior distribution.

5. The method of claim 4 , wherein the shared prior distribution is specified such that the dependencies of the at least two partial contributions are sampled from the same distribution.

6. The method of claim 5 , wherein the shared prior distribution is determined using training data from multiple training individuals and one or more predetermined hyperparameters.

7. 3. The method of claim 2, wherein each partial contribution comprises the product of the polygenic risk score and the coordinate of the individual's position along a respective axis of the ancestry space.

8. 2. The method of claim 1 , wherein calculating the genetic contribution comprises calculating a distance in the ancestral space between the individual location and a reference location in the ancestral space, and calculating the genetic contribution using the distance.

9. 9. The method of claim 8, wherein the reference position is one or more of the following: a) the positions in the ancestry space of ancestors used to train the coefficients used to calculate the polygenic risk score; b) the ancestry space is defined using reference genetic data from multiple reference individuals with multiple different ancestries; calculating the genetic contribution includes scaling each axis of the ancestry space using the variance explained by the respective axis in the reference genetic data before calculating the distance. c) the distance is a Euclidean distance in the ancestral space; and d) the genetic contribution comprises the product of the polygenic risk score and the distance.

10. The method of claim 1 , wherein calculating the genetic contribution comprises using a linear dependence on the individual location.

11. The method of claim 1 , wherein calculating the genetic contribution comprises using a non-linear dependence on the individual location.

12. The method of claim 11 , wherein the nonlinear dependence comprises one or both of a regularized function and a penalized B-spline.

13. 12. The method of claim 11, wherein the nonlinear dependence is determined using a Gaussian process as a prior distribution for calculating the genetic contribution using Bayesian inference.

14. The method of claim 13, wherein the method is one or more of the following a) to f): a) the Gaussian process is i) has a zero mean value function, or ii) have a mean vector corresponding to a prior estimate of the genetic contribution to the risk for the subject individual; b) the kernel function of the Gaussian process is a stationary function; c) the Gaussian process kernel function decays to zero as the similarity between samples decreases; d) the kernel function is a radial basis function or a rational quadratic covariance function; e) a posterior distribution for the Bayesian inference is determined using the Gaussian process and training data from a plurality of training individuals having a plurality of different ancestries, and optionally, determining the posterior distribution includes approximating the posterior distribution as a normal distribution; and f) the Gaussian process kernel function depends on one or more hyperparameters, optionally including hyperparameters related to each of the polygenic risk score, the individual location, and an interaction between the polygenic risk score and the location.

15. The method according to any one of claims 1 to 9, wherein the genetic contribution is one or more of the following a) to d): a) comprising an ancestry-independent component calculated using the polygenic risk score that is independent of the individual location; b) the genetic contribution comprises an ancestry-dependent component that is calculated based on the individual location independent of the polygenic risk score; c) each reference individual is assigned to one of multiple ancestors; and d) the ancestry space has a lower dimensionality than the individual genetic data, and determining the individual location comprises projecting the individual genetic data into the ancestry space.

16. 10. The method of claim 1, wherein the coordinate system of the ancestry space is determined by applying dimensionality reduction to the reference genetic data.

17. The method of claim 16, wherein the method is one or both of the following a) and b): a) the dimensionality reduction is i) includes principal component analysis, independent component analysis, non-negative matrix decomposition, or factor analysis; or ii) discretizing the ancestry space into a finite set of ancestors, wherein the individual positions comprise continuous or quasi-continuous membership ratios for each ancestor of the finite set of ancestors; and b) The coordinate system of the ancestry space is selected to maximize the variance of the reference genetic data explained by the ancestry space.

18. 10. The method of claim 1, wherein the dependence of the genetic contribution on the individual location and the polygenic risk score is determined using training data from a plurality of training individuals with a plurality of different ancestries, the training data including, for each of the training individuals, genetic data and whether the training individual has the phenotype or phenotype combination of interest.

19. the training data further comprises, for each of the training individuals, data having information about one or more non-genetic covariates, and the genetic contribution is jointly estimated in the presence of the non-genetic covariates; 20. The method of claim 18, wherein the method further comprises receiving individual covariate data for the subject individual, the individual covariate data comprising information regarding additional non-genetic covariates for the subject individual.

20. 20. The method of claim 19, wherein the non-genetic covariates include one or more of weight, height, behavioral characteristics, medical traits, and other biomarkers such as blood or urine based measurements.

21. wherein the risk is a relative risk to an individual with an average estimated genetic contribution, and the method comprises: calculating the relative risk for the subject individual for the phenotype or phenotype combination of interest using the genetic and non-genetic contributions to the relative risk; outputting the relative risk; 10. The method of claim 1, further comprising:

22. The method of claim 21, wherein the method is one or both of the following a) and b): a) calculating the relative risk includes determining the value of the relative risk for the subject individual from a distribution of the relative risks for the subject individual using a loss function; and b) calculating the relative risk comprises determining a hazard ratio for the subject individual, wherein the hazard ratio is normalized using the polygenic risk score and the genetic contribution.

23. 23. The method of claim 22, wherein the loss function is a mean squared error function or an asymmetric exponential loss function.

24. The method of any one of claims 1 to 9, further comprising one or both of the following a) and b): a) outputting the genetic contribution to said risk; and b) using said genetic contribution to calculate the absolute risk of said subject individual for said phenotype or phenotype combination of interest, and outputting said absolute risk;

25. 25. The method of claim 24, wherein calculating the absolute risk comprises determining a hazard ratio for the subject individual, wherein the hazard ratio is normalized using the polygenic risk score and the genetic contribution.

26. A computer program or computer readable medium comprising instructions which, when the program is executed by a computer, cause the computer to perform a method according to any of claims 1 to 9.

27. 1. An apparatus for analyzing genetic data, comprising: a processor, said processor comprising: receiving a polygenic risk score for a phenotype or combination of phenotypes of interest for the subject individual; receiving individual genetic data for the individual of interest, the genetic data comprising information regarding the ancestry of the individual of interest; determining an individual's position within ancestry space using the individual's genetic data; using said polygenic risk score and said individual location to calculate a genetic contribution to said subject individual's risk for said phenotype or phenotype combination of interest; and wherein the individual location is represented by a combination of continuous or pseudo-continuous variables, or comprises an assignment to one or a weighted combination of multiple ancestors; the ancestry space is defined using reference genetic data from a plurality of reference individuals having a plurality of different ancestries; the genetic contribution has a continuous or quasi-continuous dependence on the individual location; The device.

Images