scDRS 单细胞疾病相关分数

基因组学与单细胞RNA测序的结合

目前复杂疾病研究方法热点之一就是多组学多方法的结合,近来多种新的结合方法中,scDRS是比较有代表性的将基因组学与单细胞RNA测序的结合的方法,类似的方法还有sc-linker等,这类方法也属于一个新的交叉领域单细胞遗传学(Single-cell Genetics)

传统的基因组学与单细胞测序的结合的方法,例如MAGMA或是LDSC-SEG本质上还是基于细胞特异表达的基因集而进行的富集检测,而近来以scDRS为代表的方法则利用scRNA-seq的表达矩阵进一步深入至单个细胞的层面,能够或得更高的分辨率,这对解析疾病异质性,疾病关联的细胞亚群等方面可以发挥巨大作用 (就是让精准医疗更精准)。

scDRS的全称是 single-cell disease relevance score, 单细胞疾病相关分数,正如其名称类似于多基因分风险分数的PRS,该方法结构上或多或少也类似PRS的构建和计算方法,不过这里的个体是一个一个的细胞。PRS用于评估个体疾病的风险,而scDRS则评估细胞是否高表达疾病相关的基因。注意这里的R,一个是risk风险 (有方向),一个是relevance相关(没有方向),概念上的差异。

scDRS的方法概况

  1. 首先从GWAS结果构建疾病基因集:使用MAGMA和目标疾病的GWAS sumstats进行基因水平的关联检验,得到每个基因与疾病关联的Z分数,选定前1000个基因作为假定的疾病基因 (这个步骤不是方法的重点,除了MAGMA也可以使用其他类似方法构建疾病基因集)
  2. 然后计算单个细胞的疾病分数:scDRS会对每个细胞进行计算,量化假定疾病基因的整体表达。为了最大化检验效能,会根据MAGMA所得Z分数进行加权,同时根据每个基因在单细胞测序中特异的技术性噪音进行进行加权。
  3. 最后scDRS对所有基因集和细胞,标准化其原始的疾病分数以及原始的对照分数。然后基于所有基因集和细胞的标准化后分数的经验分布,计算每个细胞的P值
  4. 利用得到的分数可以进行多种下游分析,分析包括(1)单个细胞层面的关联性检验(2)细胞类型关联性检验 以及(3)基因关联性检验等, 这些分析的P值均通过MC检验Monte Carlo 蒙特卡洛检验)获得。

用例(下游分析)

分析与表型相关的细胞类型

相比于传统方法LDSC或MAGMA,scDRS可以检验细胞类型的相关性,还可以检验同一细胞类型内的异质性。

作者列举了 22个表型与19种细胞的关联结果的热力图,Y轴为表型,X轴为细胞类型,方框表示显著相关,×表示细胞类型内的异质性,颜色深浅表示显著关联细胞的比例。

发现存在异质性的细胞亚群

作者研究了自身免疫性疾病中T细胞的异质性,11个T细胞群中(a),与IBD相关的T细胞构成了4个新的群(b)。

基因-分数的关联分析

scDRS检验基因是否与GWAS由来的基因集中的基因共表达。

相比于MAGMA,scDRS能更准确地识别出疾病相关的基因

典型的分析流程

只需要GWAS的Sumstats和单细胞RNA测序的数据,好消息是这两个都可以很容易从公开数据库中获得。

具体流程官方文档以及大牛博客已经写得很详细了,

分析代码可以参考 https://martinjzhang.github.io/scDRS/

以及 https://zhuanlan.zhihu.com/p/592128325

注意的点

  1. scDRS的假设是基因集里的基因与疾病有关,与疾病有关的基因会在疾病相关的细胞群体(可以是疾病细胞或健康细胞)里高表达,与基因的方向无关。 scDRS并不是假设基因集里的基因会在疾病细胞里高表达。注意这里概念上细节的差别。(https://github.com/martinjzhang/scDRS/issues/42)
  2. 基本数据要求:为了能够得到足够的检验效能,GWAS的heritability z-score最好大于5,或样本量大于10万。(https://martinjzhang.github.io/scDRS/faq.html#which-gwas-and-scrna-seq-data-to-use)
  3. Seurat格式的数据需要转换为scanpy使用的h5ad格式(表达矩阵不能有负值,Seurat的scaled.data里表达矩阵会有负值,转换时要注意)(https://github.com/martinjzhang/scDRS/issues/44);转换可以使用 SeuratDisk https://mojaveazure.github.io/seurat-disk/articles/convert-anndata.html
  4. 为了增加检验效能,单细胞RNA测序可以事先进行 imputationhttps://github.com/martinjzhang/scDRS/issues/32)
  5. -adj-prop 可以调整细胞类别的比例 (当某些细胞种类比例过高的时候) (https://github.com/martinjzhang/scDRS/issues/32)
  6. 显著细胞的数太少:正常现象,仍然可以进行group分析(检验效能通常更高) (https://martinjzhang.github.io/scDRS/faq.html#scdrs-detected-few-significant-cells-fdr-0-2)

参考

Zhang, M. J., Hou, K., Dey, K. K., Sakaue, S., Jagadeesh, K. A., Weinand, K., … & Price, A. L. (2022). Polygenic enrichment distinguishes disease associations of individual cells in single-cell RNA-seq data. Nature genetics54(10), 1572-1580.

MAGMA 基因及通路分析

目录

1 背景
2 算法简介
3 使用教程
 3.1 magma软件与参考数据下载
 3.2 注释
 3.3 基因关联分析
 3.4 通路分析
4 参考

背景

通过以有生物学意义的方式整合复杂疾病的信息的基因以及通路分析,是单变异GWAS的有效补充。

本文将介绍MAGMA,一个基于多变量回归模型的基因以及通路分析工具。

算法简介

Gene-based分析

magma的gene-based分析模型采用了多元线性主成分回归(multiple linear principal components regression)

α0g为截距,Xg* 为PC矩阵,W与βg为协变量变量的效应(可选),αg为遗传效应,εg为残差

计算时:

首先根据染色体位置提取某个基因SNP矩阵,计算PC

去除掉特征值过小的PC

对剩下的主成分进行回归

最后通过F检验获得p值 (H0: αg =0 )

Gene set /Pathway 基因集 通路分析

进行通路分析时,magma首先将上一步所得到的每个基因的p值通过probit方程转化为z值,

这里的zg整体上近似服从正态分布,zg反映了该基因关联的强度。

Self-contained gene-set analysis检验一个通路上的基因的整体上是否与表型相关联。使用通路里基因的z值,可以构建一个只有截距项的回归,然后检验β是否大于0:

Competitive gene-set analysis检验在一个通路里的基因,是否比这个通路里其他基因的相关性更高。

因为通路分析中的回归模型为标准线性回归,假设残差项独立正态分布,即

但实际上,由于LD,相邻的基因可能会互相关联,而打破上述假设。所以magma采用了广义最小二乘法,并假设

R为基因-基因的关联矩阵。

概括性数据分析

在没有原始基因型数据的情况下,magma也提供分析概括性数据SNP的p值的方法,即逐个分析某个基因中的SNP,然后将SNP的p值合并成该基因的检验统计量,但是该方法需要提供一个相近群体的LD参考面板。

使用教程

magma软件与参考数据下载

https://ctg.cncr.nl/software/magma

1.软件下载,下载对应系统的版本即可

2.基因注释参考,magma提供了多种版本,下载对因自己数据的版本

3.LD参考面板,来自1000 genome,下载相应人群的文件即可

magma三步骤:

  1. 注释:将自己的SNP根据染色体位置注释到基因上
  2. Gene-based 分析
  3. 基因集/通路分析

SNP注释 Annotation

第一步需要对gwas中所包含的SNP进行注释,在这里就是将SNP根据染色体上的位置对应到相应基因上。

示例代码如下:

magma \
	--annotate \
	--snp-loc [SNPLOC_FILE] \
	--gene-loc [GENELOC_FILE] \
	--out [ANNOT_PREFIX]

snp-loc 文件应当包含三列,SNPID,CHR以及POS,这里也可以直接使用plink的bim文件

rs540836310	1 861725
rs139858754	1 861728

gene-loc 可以直接使用magma提供的基因注释参考文件 (注意版本)

注释完成后会得到 .genes.annot文件,内容为第一例为geneID 之后为在这个基因内的SNP

816	rs540836310	rs139858754	rs188152259	rs13302982	rs76842830	rs13303101	rs74442310	rs6680268

基因关联分析 Gene-based analysis

第二部,进行基于基因的关联分析,

示例代码如下所示:

magma \
	--bfile [REFDATA] \
	--pval [PVAL_FILE] \
	N=[N] \
	--gene-annot [ANNOT_PREFIX].genes.annot \
	--out [GENE_PREFIX]

bfile为原始数据或是参考LD面板,如果数据量不大可以直接使用自己的plink的bed格式原始数据,在原始数据无法获得的时候可以使用magma提供的1000 genome参考数据,biobank级别的数据的情况下,可以随机抽取某个族裔无亲缘关系的一定人数(例如20000人)来构建自己的参考面板。

pval为SNP的p值文件,包含两列 SNP 以及 P, 这里N为样本量,

gene-annot为上一步注释后得到的文件。

完成后有两个文件输出;1 .genes.out 2 ..genes.raw

1 .genes.out 基因关联分析的结果

GENE       CHR      START       STOP  NSNPS  NPARAM       N        ZSTAT            P
148398       1     859993     879961    125      27  157848      -1.7291       0.9581

2 ..genes.raw 在第三步通路分析时会使用到这个文件

# VERSION = 109
# COVAR = NSAMP MAC
148398 1 859993 879961 125 27 157848 133.304 -1.72907

基因集 通路分析 gene-set analysis 、 pathway

第三步,通路分析

这一步是对显著的基因是否富集于某个基因集或通路进行检验,

如果是对经典的通路进行检验,可以在MSigDB(GSEA)下载通路文件 msigdb.v7.4.entrez.gmt

MSigDB:https://www.gsea-msigdb.org/gsea/downloads.jsp#msigdb

示例代码如下:

magma \
	--gene-results [GENE_PREFIX].genes.raw \
	--set-annot [SET_FILE] \
	--out [GS_PREFIX]

gene-results 为上一步基因关联分析所得的.genes.raw文件

set-annot 为基因集或通路的定义文件,可以直接使用MSigDB下载的文件(注意基因ID要与之前相对应)

默认输出三个.gsa.out .gsa.genes.out, .gsa.sets.genes.out

.gsa.out 为最主要的输出文件,包含了各个通路的检验结果

# TOTAL_GENES = 17198
# TEST_DIRECTION = one-sided, positive (set), two-sided (covar)
# CONDITIONED_INTERNAL = gene size, gene density, inverse mac, log(gene size), log(gene density), log(inverse mac)
VARIABLE                                TYPE  NGENES         BETA     BETA_STD           SE            P FULL_NAME
chr1p12                                  SET      18    -0.096902     -0.00381      0.30843      0.71093 chr1p12


.gsa.genes.out 为在此分析中使用到的基因关联结果

GENE       CHR      START       STOP  NSNPS  NPARAM       N        ZSTAT  ZFITTED_BASE  ZRESID_BASE
148398       1     859993     879961    125      42  12345      -1.62             0      -1.62

.gsa.sets.genes.out 为多重检验调整后, 显著的通路里的基因的信息

参考

MAGMA: Generalized Gene-Set Analysis of GWAS Data. de Leeuw CA, Mooij JM, Heskes T, Posthuma D (2015) MAGMA: Generalized Gene-Set Analysis of GWAS Data. PLOS Computational Biology 11(4): e1004219. https://doi.org/10.1371/journal.pcbi.1004219

https://ctg.cncr.nl/software/magma