多基因风险分数 PRS( Polygenic risk score)系列之五:使用PRS-CS计算PRS(beta-shrinkage方法)

本文内容

  1. PRS系列回顾
  2. PRS-CS简介
  3. 概念框架
  4. 使用方法
  5. 参考

系列回顾:

PRS-CS简介

本文将介绍一种使用GWAS数据与外部LD参考面板来计算SNP后验效应量的方法, PRS-CS。

该方法利用了高维贝叶斯回归的框架,与之前类似研究最大的不同点便是对于beta使用了连续收缩(continuous shrinkage, CS)先验分布,适应更多样的遗传结构,并大幅提高了计算效率,使局部LD特征的多变量建模成为可能。

概念框架

  1. 考虑一个贝叶斯高维回归框架
  2. 常用的beta的先验分布通常可以表示为多个正态分布的比例混合
  3. 为了构建更为灵活的模型以更好反映遗传结构,一些方法采用了两个或多个质点或分布的离散混合,相比单纯的正太先验分布,这类方法通常能够构建更广的效应量分布。例如LDpred。
  4. 但此类方法后验推断需要极大的计算量,可能不能准确的建模出局部的LD结构。
  5. 为了解决这些问题,PRS-CS的作者采用了一种与先前离散混合概念不同的先验分布,也就是连续收缩先验分布,也就是正态分布的全局-局部比例混合
  6. PRS-CS对于全局比例系数采取两种不同算法进行估计。一是搜索固定的系数,而是通过指定全局系数的先验分布来在给定数据中进行估计。

详细推导请参考原论文:T Ge, CY Chen, Y Ni, YCA Feng, JW Smoller. Polygenic Prediction via Bayesian Regression and Continuous Shrinkage Priors. Nature Communications, 10:1776, 2019.

使用方法

PRS-CS由Python语言写成,是python2与3同时兼容的。需要安装scipy与h5py这两个包。

PRS-CS 本体与外部LD参考面板的下载:

https://github.com/getian107/PRScs

git clone <https://github.com/getian107/PRScs.git>

下载完成后,检验是否安装成功:

./PRScs.py --help

下载对应人群的LD参考面板,作者提供了千人基因组以及UKB数据的参考面板:

https://github.com/getian107/PRScs#getting-started

输入准备

  1. 上面下载的LD参考面板 plink格式
  2. 目标群体的bim文件:该文件只是用来提供target数据集中的SNP列表
  3. GWAS sumstats

sumstats的格式如下

SNP          A1   A2   BETA      P
rs4970383    C    A    -0.0064   4.7780e-01
rs4475691    C    T    -0.0145   1.2450e-01
rs13302982   A    G    -0.0232   2.4290e-01

或是

SNP          A1   A2   OR        P
rs4970383    A    C    0.9825    0.5737                 
rs4475691    T    C    0.9436    0.0691
rs13302982   A    G    1.1337    0.0209
...

注意A1是effect allele, A2是non-effect allele

PRS模型计算

使用方法如下

python PRScs.py \
  --ref_dir=PATH_TO_REFERENCE \ # LD参考面板的文件夹路径
  --bim_prefix=VALIDATION_BIM_PREFIX \ #目标群体plink的bim格式文件的路径
  --sst_file=SUM_STATS_FILE \  # GWAS sumstats的路径
  --n_gwas=GWAS_SAMPLE_SIZE \ # GWAS的样本量大小
  --out_dir=OUTPUT_DIR  #输出文件夹

#以下为可选项
  #--a=PARAM_A  \ gamma-gamma prior中的参数a,默认为1
  #--b=PARAM_B  \ gamma-gamma prior中的参数b,默认为0.5
  #--phi=PARAM_PHI \ 全局比例系数,不指定时自动估计(PRS-CS-auto),也可以小规模网格搜索(phi=1e-6, 1e-4, 1e-2, 1)
  #--n_iter=MCMC_ITERATIONS \ MCMC迭代次数,默认1000
  #--n_burnin=MCMC_BURNIN \ MCMC中burnin的次数,默认500
  #--thin=MCMC_THINNING_FACTOR \ 马尔科夫链的thinning factor,默认为5
  #--chrom=CHROM \ 可以只单独计算一条染色体
  #--beta_std=BETA_STD \ 若为True,则输出标准化的后验SNP效应量 
  #--seed=SEED 随机数种子

PRS-CS使用scipy,会自动占用所有可用的cpu核心,使用服务器时可能会出现干扰,可以通过以下shell脚本指定使用核心数:

export MKL_NUM_THREADS=$N_THREADS
export NUMEXPR_NUM_THREADS=$N_THREADS
export OMP_NUM_THREADS=$N_THREADS

例如,只用一个核时,N_THREADS=1

输出结果

输出文件包含五列,rsID,碱基位置,A1,A2 与后验效应量估计值

然后使用该文件和plink可以对目标人群计算PRS: https://www.cog-genomics.org/plink/1.9/score

分染色体计算时,最后加和各个染色体的PRS即可

测试例子

使用eur的参考面板与test_data文件夹里的gwas数据,以及一个有22号染色体上1000个SNP的bim文件:

python PRScs.py \
--ref_dir=path_to_ref/ldblk_1kg_eur \
--bim_prefix=path_to_bim/test \
--sst_file=path_to_sumstats/sumstats.txt \
--n_gwas=200000 \
--chrom=22 \
--phi=1e-2 \
--out_dir=path_to_output/eur

参考

T Ge, CY Chen, Y Ni, YCA Feng, JW Smoller. Polygenic Prediction via Bayesian Regression and Continuous Shrinkage Priors. Nature Communications, 10:1776, 2019.

《多基因风险分数 PRS( Polygenic risk score)系列之五:使用PRS-CS计算PRS(beta-shrinkage方法)》有2个想法

发表评论

Fill in your details below or click an icon to log in:

WordPress.com 徽标

您正在使用您的 WordPress.com 账号评论。 注销 /  更改 )

Twitter picture

您正在使用您的 Twitter 账号评论。 注销 /  更改 )

Facebook photo

您正在使用您的 Facebook 账号评论。 注销 /  更改 )

Connecting to %s