RSS

workflowr

Ignored files:
    Ignored:    .DS_Store
    Ignored:    .Rhistory
    Ignored:    .Rproj.user/
    Ignored:    analysis/.Rhistory
    Ignored:    docs/.DS_Store
    Ignored:    docs/figure/Test.Rmd/

Untracked files:
    Untracked:  analysis/bibliography.bib

Zhu and Stephens (2017) introduce a ‘Regression with Summary Statistics’ (RSS) likelihood, which relates the multiple regression coefficients to univariate regression results that are often easily available. They focus on the simplest case of a quantitative trait (e.g., height) measured on random samples from a population. Here is the summary of the method.

Consider the multiple linear regression model: \[ \mathbf{y} = \mathbf{X} \mathbf{b} + \mathbf{\epsilon}, \] where \(\mathbf{y}\)is an \(n \times 1\) (centered) vector, \(\mathbf{X}\) is an \(n \times p\) (column-centered) matrix, \(\mathbf{b}\) is the \(p \times 1\) vector of multiple regression coefficients, and \(\mathbf{\epsilon}\) is the error term. Assuming that individual-level data are not available, but instead the summary statistics from p simple linear regression are provided: \[ \hat{b}_j = (X_j^\intercal X_j)^{-1} X_j^\intercal \mathbf{y} \quad \hat{s}_j^2 = (nX_j^\intercal X_j)^{-1} (\mathbf{y} - X_j\hat{b}_j)^\intercal (\mathbf{y} - X_j\hat{b}_j) \] The RSS likelihood is \[ L(\mathbf{b};\mathbf{\hat{b}}, \hat{D}, \hat{R}) = N(\mathbf{\hat{b}}; \hat{D} \hat{R} \hat{D}^{-1} \mathbf{b}, \hat{D} \hat{R} \hat{D}), \] where \(\hat{d}_j^2 = \hat{s}_j^2 + \frac{\hat{b}_j^2}{n}\), \(\hat{D} = diag(\hat{d}_j)\), and \(\hat{R}\) is the estimated LD among SNPs in the population from which the genotypes were sampled. Note that \(\hat{d}_j^2 = \mathbf{y}^\intercal \mathbf{y}/(nX_j^\intercal X_j)\).

The prior assumes that \(\mathbf{b}\) is independent of \(R\) apriori, with the prior on \(b_j\) being a mixture of two normal distributions \[ b_j \sim \pi N(0, \sigma_B^2 + \sigma_P^2) + (1-\pi) N(0, \sigma_P^2), \] \[ \log \pi \sim \mathcal{U}(\log(1/p), 0) \quad h \sim \mathcal{U}(0,1) \quad \rho \sim \mathcal{U}(0,1), \] \[ \sigma_B^2(S) = h\rho \left( \pi \sum_{j=1}^{p} n^{-1} \hat{d}_j^{-2}\right) \quad \sigma_P^{2}(S) = h(1-\rho) \left(\sum_{j=1}^{p} n^{-1} \hat{d}_j^{-2}\right), \] Since \(nd_j^2 = \mathbf{y}^\intercal \mathbf{y}/X_j^\intercal X_j\), the prior on hyperparameters ensure that the effect sizes of both components do not depend on n, and have the same measurement unit as \(\mathbf{y}\). \(\rho = \frac{\pi\sigma_B^2}{\pi\sigma_B^2 + \sigma_P^2}\) represents the expected proportion of total genetic variation explained by the sparse component. \(h\) represents the proportion of total variation in \(\mathbf{y}\) explained by X.