2020-05-04 eigenMT

eigenMT:解释变异间连锁不平衡的eQTL研究的有效多重检验调整

简介

特征最小二乘法是顺式eQTL研究中一种计算效率较高的多重检验校正方法。通常,现代顺式eQTL研究可以正确地使用排列对单个基因进行跨越数千个变异位点的多重测试。这种校正方法虽然准确,但对于大规模队列测序/大样本测序,需要数天到数周的高计算成本。与其相比,我们的方法减少了数量级的计算负担,同时保持了较高的精度。为了做到这一点,我们的方法试图通过直接考虑被测试变异位点之间的连锁不平衡LD关系来估计对给定基因执行的有效测试的数量。

下载安装

首先下载或克隆此存储库https://github.com/joed3/eigenMT。eigurnMT作为独立的python脚本运行。它需要安装Python(版本2.7或更高版本)才能安装。要下载并安装Python发行版,有几个方便的选项:

这些捆绑安装已经包括许多运行EigurnMT所需的模块,包括:

同样,这些模块应该预先打包在上面的某个捆绑的python安装中。如果您决定自己安装这些软件包,请查看列出的网站获取详细说明。

最后,如果您希望运行Matrix eQTL示例脚本,则需要安装Matrix eQTL包,这只需打开R会话并运行以下命令即可完成:

install.packages('MatrixEQTL')

关于MatrixEQTL,技能树的公众号已经做了详细的介绍:
用 Matrix eQTL 进行 eQTL 分析https://mp.weixin.qq.com/s/bAEdj3xArKAv2yG3YSzIrQ
我这里就不多说了,那么用python怎么来完成呢?

输入数据准备

可以用如下命令行:

python eigenMT.py --CHROM <chromosome ID>
    --QTL <Matrix eQTL SNP-gene tests file> \
    --GEN <Matrix eQTL genotype matrix> \
    --POS <Matrix eQTL genotype position file> \
    --var_thresh [variance explained threshold, default 0.99] \
    --OUT <output filename> \
    --window [window size, default 200]

eigenMT 的设计是按照以往发表的流程pipeline 出来的数据 Matrix eQTL (Shabalin, Bioinformatics 2012) 来进行 eQTL calling. 输入的参数具体解释如下:

Argument Description
CHROM 染色体ID。指示将在哪个染色体上执行分析。必须与Matrix eQTL SNP基因测试文件中使用的ID匹配。
QTL Matrix eQTL SNP-基因测试文件的文件名。有关示例,请参见Example/cis.qtls.txt。
GEN 矩阵eQTL格式的基因型矩阵。可以接受硬编码或基于剂量的基因型。丢失的基因型必须编码为NA,并将在校正过程中归入平均基因型。有关示例,请参见Example/genetyes.txt。
GENPOS 矩阵eQTL基因型位置文件。有关示例,请参见example/position tions.txt。
var_thresh 基因型相关矩阵中解释的变异量的阈值。默认值为99%解释的差异。增加这一阈值将增加对有效试验次数(M_eff)的估计,并降低近似的精度。
OUT 输出文件名。输出文件格式如下描述。
window 窗口大小参数。确定要将每个基因的基因型矩阵拆分成的不相交窗口的大小。默认值为200个SNP。我们建议使用至少50个SNP到200个SNP的窗口大小,以平衡精度和速度。
cis_dist bp 距离TSS 进行多重检验较正。默认1e6。
external 逻辑指示所提供的基因型矩阵是否不同于最初用于测试顺-eQTL的基因型矩阵。

虽然我们将软件设计为直接利用Matrix eQTL的输出,但eignMT可以处理更通用的QTL输入文件。这种更通用的输入格式只需要三列:1.SNP ID、2.基因ID和3.cis-eQTL p值。此文件必须包含标题行,其中 p值列由 p-value指示。

输出文件

输出文件为制表符分隔格式,包含以下列。

Column Description
1 variant ID
2 Gene ID
3 estimate of effect size BETA from Matrix eQTL
4 T-statistic from Matrix eQTL
5 p-value from Matrix eQTL
6 FDR estimate from Matrix eQTL
7 eigenMT corrected p-value
8 estimated number of independent tests for the gene

注意:输出的前5列将对应于典型的Matrix eQTL输出。此外,每个测试的基因将在输出文件中出现一次,其中包含其最重要的SNP和特征MT校正的p-value

示例

我们提供了一个用于测试的小示例数据集。我们以Matrix eQTL格式提供基因型和归一化基因表达矩阵以及相应的位置文件。该基因型矩阵是作为GEUVADIS研究的一部分的19号染色体的EUR373样本。归一化表达矩阵覆盖了19号染色体上的前100个基因,并在GEUVADIS数据集中进行了量化。在DESeq2进行方差稳定并去除30个同级因子之后,表达矩阵已被分位数归一化为标准正态分布。我们还提供了使用Matrix eQTL对这些样本进行顺-eQTL测试的样本。这个样本变体基因测试文件代表100个基因。要提取示例数据集,请运行

tar -zxvf example.tar.gz

我们提供了一个运行Matrix eQTL的示例R脚本。要了解有关使用Matrix eQTL的更多信息,请参阅CRAN上的软件包文档package documentation on CRANoriginal paper。要生成cis-eQTL测试文件,请使用以下命令:

Rscript example/MatrixEQTL.R --SNP example/genotypes.txt \
    --GE example/phenotypes.txt \
    --snploc example/gen.positions.txt \
    --geneloc example/phe.positions.txt \
    --OutputFilePath example/cis.eqtls.txt \
    --cisDist 1e6 \
    --ANOVA FALSE \
    --pValue 1 > example/matrix.eqtl.log

要对示例数据运行eogen MT,请使用以下命令:

python eigenMT.py --CHROM 19 \
    --QTL example/cis.eqtls.txt \
    --GEN example/genotypes.txt \
    --GENPOS example/gen.positions.txt \
    --PHEPOS example/phe.positions.txt \
    --OUT example/exampleOut.txt

注:本例使用窗口大小(200个SNP)和差异解释阈值(0.99)的默认设置。

要可视化更正结果,请使用以下命令:

Rscript example/compareToEmpirical.R

此R脚本将生成一个包含三个数字的PDF:
(1)特征MT校正的P值(y轴)与来自10000个排列(x轴)的经验P值的曲线图。
(2)绘图与(1)相同,但x轴和y轴的变换为-log10。
(3)与(2)中相同的曲线图被放大以排除经验P值=1。在这个最后的曲线图中,点应该是强相关的,并且应该略低于对角线(大部分)。

外部基因型数据

我们的方法提供了使用来自独立样本群体的基因型矩阵进行多重测试校正的能力,而不是最初用于顺式eQTL测试的基因型矩阵。重要的是要注意,这个外部基因型矩阵应该来自与正在研究的背景人群相同的背景人群。例如,如果在欧洲血统的样本中进行顺式eQTL测试,那么使用的任何外部基因型矩阵都应该来自类似的欧洲人群。我们已经证明,与使用基因型数据进行研究相比,使用具有更大样本量的研究中的基因型数据可以提高我们方法的准确性。来自大型研究的基因型数据将为每个基因周围的变异提供更好的连锁不平衡LD结构估计,从而改进我们对有效测试次数的估计。

人口分层和其他协变量分析

我们建议首先从表达式矩阵中去除群体分层(基因型PC、群体或祖先分配)和其他协变量(年龄、性别、同伴因素等)的影响。我们已经证明,通过首先消除这些影响,然后对逆秩归一化残差执行cis-eQTL调用,该软件保持了保守性和准确性。

该软件引用信息

Davis JR, Fresard L, Knowles DA, Pala M, Bustamante CD, Battle A, Montgomery SB (2015). An Efficient Multiple-Testing Adjustment for eQTL Studies that Accounts for Linkage Disequilibrium between Variants. The American Journal of Human Genetics, 98(1), 216–224. doi: (http://doi.org/10.1016/j.ajhg.2015.11.021)

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 204,293评论 6 478
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 85,604评论 2 381
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 150,958评论 0 337
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,729评论 1 277
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,719评论 5 366
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,630评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,000评论 3 397
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,665评论 0 258
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,909评论 1 299
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,646评论 2 321
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,726评论 1 330
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,400评论 4 321
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,986评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,959评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,197评论 1 260
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 44,996评论 2 349
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,481评论 2 342