首先是Windows
一个做孟德尔随机化的过程遇到的报错:
bmi_exp_dat <- clump_data(bmi_exp_dat,clump_r2=0.01,pop = "EUR")
Please look at vignettes for options on running this locally if you need to run many instances of this command.
Clumping C5nTuK, 5340156 variants, using EUR population reference
Error in api_query("ld/clump", query = list(rsid = dat[["rsid"]], pval = dat[["pval"]], :
The query to MR-Base exceeded 300 seconds and timed out. Please simplify the query
以下是作者给出的解决办法:
- Extracting without clumping
- Performing clumping on each chromosome separately
Alternatively, you could try downloading the VCF from IEU OpenGWAS project and using the gwasvcf to extract based on p-value, and use the ieugwasr package to do clumping locally. We're trying to get to a point where it's easy to do heavier computation on these data locally.
LD reference files listed in the gwasvcf page
看到作者建议可以用ieugwasr这个包做本地clump。
下面尝试做本地clump
首先安装ieugwasr包
devtools::install_github("mrcieu/ieugwasr")
然后我们可以看到
|
Wrapper for clump function using local plink binary and ld reference dataset |
ld_clump_local()这个函数是用来做本地clump的,用法和参数如下:
ld_clump_local(dat, clump_kb, clump_r2, clump_p, bfile, plink_bin)
Arguments
dat | Dataframe. Must have a variant name column ("variant") and pval column called "pval". If id is present then clumping will be done per unique id. |
---|---|
clump_kb | Clumping kb window. Default is very strict, 10000 |
clump_r2 | Clumping r2 threshold. Default is very strict, 0.001 |
clump_p | Clumping sig level for index variants. Default = 1 (i.e. no threshold) |
bfile | If this is provided then will use the API. Default = NULL |
plink_bin | Specify path to plink binary. Default = NULL. See https://github.com/explodecomputer/plinkbinr for convenient access to plink binaries |
需要plink包,安装一下plink包,获得可执行二进制plink地址
devtools::install_github("explodecomputer/plinkbinr")
library(plinkbinr)
get_plink_exe()
#[1] "D:/R-4.1.1/library/plinkbinr/bin/plink_Windows.exe"
然后我们下载bfile,我的GWAS是欧洲人种所以用EUR*
wget http://fileserve.mrcieu.ac.uk/ld/1kg.v3.tgz
#ld_clump_local()报错,作者是把两个函数合并了,在ld_clump函数里包含bfile参数的话,会自动调用ld_clump_local
所以改用ld_clump
ld_clump( dplyr::tibble(rsid=dat$rsid, pval=dat$pval, id=dat$trait_id), plink_bin = genetics.binaRies::get_plink_binary(), bfile = "/path/to/reference/EUR" )
"
b <- ld_clump(
dplyr::tibble(rsid=a$rsid, pval=a$p, id=a$id),
#get_plink_exe()
plink_bin = "D:/R-4.1.1/library/plinkbinr/bin/plink_Windows.exe",
#欧洲人群参考基因组位置
bfile = "D:/EUR_ref/EUR"
)
注意a的列名,必须要有:(with the following columns:)
rsid
pval
trait_id
然后就完成了,可以看到跟在线方法去除的LD SNP个数是一样的
PLINK v1.90b6.10 64-bit (17 Jun 2019) www.cog-genomics.org/plink/1.9/
(C) 2005-2019 Shaun Purcell, Christopher Chang GNU General Public License v3
Logging to C:\Users\Lenovo\AppData\Local\Temp\RtmpS8aq1j\file2e646a65362c.log.
Options in effect:
--bfile D:/EUR_ref/EUR
--clump C:\Users\Lenovo\AppData\Local\Temp\RtmpS8aq1j\file2e646a65362c
--clump-kb 10000
--clump-p1 0.99
--clump-r2 0.001
--out C:\Users\Lenovo\AppData\Local\Temp\RtmpS8aq1j\file2e646a65362c
32549 MB RAM detected; reserving 16274 MB for main workspace.
8550156 variants loaded from .bim file.
503 people (0 males, 0 females, 503 ambiguous) loaded from .fam.
Ambiguous sex IDs written to
C:\Users\Lenovo\AppData\Local\Temp\RtmpS8aq1j\file2e646a65362c.nosex .
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 503 founders and 0 nonfounders present.
Calculating allele frequencies... done.
8550156 variants and 503 people pass filters and QC.
Note: No phenotypes present.
--clump: 78 clumps formed from 2011 top variants.
Results written to
C:\Users\Lenovo\AppData\Local\Temp\RtmpS8aq1j\file2e646a65362c.clumped .
Warning: 'rs9930333' is missing from the main dataset, and is a top variant.
Warning: 'rs8083289' is missing from the main dataset, and is a top variant.
Warning: 'rs2683992' is missing from the main dataset, and is a top variant.
27 more top variant IDs missing; see log file.
Removing 1963 of 2041 variants due to LD with other variants or absence from LD reference panel
Linux同理
expo_dat <- expo_dat[which(expo_dat$pval.exposure<1e-5),]
b <- ld_clump(
dplyr::tibble(rsid=expo_dat$SNP, pval=expo_dat$pval.exposure, id=expo_dat$id.exposure),
#plink位置
plink_bin = "/GM_GWAS_LD_clumped_snps/plink",
bfile = "/GM_GWAS_LD_clumped_snps/EUR_ref/EUR",
clump_kb = 1000,clump_r2 = 0.1
)
expo_dat <- expo_dat[which(expo_dat$SNP %in% b$rsid),]
Clumping A9o7Sb, 50 variants, using EUR population reference
PLINK v1.90b6.21 64-bit (19 Oct 2020) www.cog-genomics.org/plink/1.9/
(C) 2005-2020 Shaun Purcell, Christopher Chang GNU General Public License v3
Logging to /tmp/RtmpdrrA6m/file4bc2f27623281.log.
Options in effect:
--bfile /share2/pub/zhenggw/zhenggw/GM_GWAS_LD_clumped_snps/EUR_ref/EUR
--clump /tmp/RtmpdrrA6m/file4bc2f27623281
--clump-kb 1000
--clump-p1 0.99
--clump-r2 0.1
--out /tmp/RtmpdrrA6m/file4bc2f27623281
515328 MB RAM detected; reserving 257664 MB for main workspace.
Allocated 61144 MB successfully, after larger attempt(s) failed.
8550156 variants loaded from .bim file.
503 people (0 males, 0 females, 503 ambiguous) loaded from .fam.
Ambiguous sex IDs written to /tmp/RtmpdrrA6m/file4bc2f27623281.nosex .
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 503 founders and 0 nonfounders present.
Calculating allele frequencies... done.
8550156 variants and 503 people pass filters and QC.
Note: No phenotypes present.
Warning: 'rs55654746' is missing from the main dataset, and is a top variant.
--clump: 8 clumps formed from 49 top variants.
Results written to /tmp/RtmpdrrA6m/file4bc2f27623281.clumped .
Removing 42 of 50 variants due to LD with other variants or absence from LD reference panel
更新了 更新了
新增加In file(file, "rt") : cannot open file '/tmp/Rtmp8ur92k/xxx': No such file or directory 解决办法
好多人出现了这个问题,我今天做LD matrix也遇到了,来解决一下,报错信息是找不到文件,我考虑是tmp文件在程序运行过程中被自动删除了。尝试手动添加tmp文件名
ld_matrix_local <- function(variants, bfile, plink_bin, with_alleles=TRUE)
{
# Make textfile
shell <- ifelse(Sys.info()['sysname'] == "Windows", "cmd", "sh")
fn <- tempfile()
ld_matrix_local <- function(variants, bfile, plink_bin, with_alleles=TRUE)
{
# Make textfile
shell <- ifelse(Sys.info()['sysname'] == "Windows", "cmd", "sh")
fn <- "~/plinkLD/tmpfile/tmp"
发现还是同样error
但是提示里面发现Using up to 39 threads (change this with --threads),考虑可能是线程冲突,导致文件没有产生,所以在源代码加入线程变量
fun2 <- paste0(
shQuote(plink_bin, type=shell),
" --bfile ", shQuote(bfile, type=shell),
" --extract ", shQuote(fn, type=shell),
" --r square ",
" --keep-allele-order ",
" --threads 1 ",
" --out ", shQuote(fn, type=shell)
)
515328 MB RAM detected; reserving 257664 MB for main workspace.
Allocated 61144 MB successfully, after larger attempt(s) failed.
8550156 variants loaded from .bim file.
503 people (0 males, 0 females, 503 ambiguous) loaded from .fam.
Ambiguous sex IDs written to /tmp/Rtmp8ur92k/file557ea8a4015e.nosex .
--extract: 3839 variants remaining.
Using 1 thread.
Before main variant filters, 503 founders and 0 nonfounders present.
Calculating allele frequencies... done.
3839 variants and 503 people pass filters and QC.
Note: No phenotypes present.
--r square to /tmp/Rtmp8ur92k/file557ea8a4015e.ld ... done.
顺利解决!!!文章来源:https://www.toymoban.com/news/detail-411528.html
所以原函数加一行代码就好 --threads 1 ,这个线程1好像速度也蛮快的,不懂有什么影响文章来源地址https://www.toymoban.com/news/detail-411528.html
到了这里,关于TwoSampleMR:local clump(MR-Base exceeded 300 seconds) 包括Windows和Linux R解决办法的文章就介绍完了。如果您还想了解更多内容,请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章,希望大家以后多多支持TOY模板网!