library(nhanesR)
# 2.安装nhanesR
nhshelp::install_nhanesR('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')
config_path('E/nhanesR')
library(nhanesR)
config_path('E:/NHANES')
config_years()
config_items()
library(nhanesR)
# 2.安装nhanesR
nhshelp::install_nhanesR('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')
library(nhanesR)
alq<- nhs_tsv('alq','!~alqy')
config_path(E:/nhanesR)
# 1.安装nhshelp
devtools::install_github('yikeshu0611/nhshelp',force=TRUE)
# 2.安装nhanesR
nhshelp::install_nhanesR('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')
# 2.安装nhanesR
nhshelp::install_nhanesR('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')
library(nhanesR)
config_path('E:/nhanesR')
# 1.安装nhshelp
devtools::install_github('yikeshu0611/nhshelp',force=TRUE)
# 2.安装nhanesR
nhshelp::install_nhanesR('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')
library(nhanesR)
install.packages("Matrix")
##1.提取暴露[exp]数据------
exp <- extract_instruments(outcomes = 'ebi-a-GCST90016671',p1 = 5e-08,
clump = TRUE,r2 = 0.01,kb = 10000,
access_token = NULL)#提取暴露数据
library(TwoSampleMR)
library(MRPRESSO)
##1.提取暴露[exp]数据------
exp <- extract_instruments(outcomes = 'ebi-a-GCST90016671',p1 = 5e-08,
clump = TRUE,r2 = 0.01,kb = 10000,
access_token = NULL)#提取暴露数据
save(exp,file="exp_Visceral adipose tissue volume.rda")
##1.提取暴露[exp]数据------
exp <- extract_instruments(outcomes = 'ebi-a-GCST90016672',p1 = 5e-08,
clump = TRUE,r2 = 0.01,kb = 10000,
access_token = NULL)#提取暴露数据
colnames(exp)
save(exp,file="exp_Abdominal subcutaneous adipose tissue volume.rda")
getwd()
if (!require("BiocManager")) install.packages("BiocManager")
BiocManager::install("MungeSumstats")
BiocManager::install("BSgenome.Hsapiens.1000genomes.hs37d5")
install.packages("path/to/BSgenome.Hsapiens.1000genomes.hs37d5_0.99.1.tar.gz", repos = NULL, type = "source")
install.packages("~/Desktop/BSgenome.Hsapiens.1000genomes.hs37d5_0.99.1.tar.gz", repos = NULL, type = "source")
BiocManager::install("BSgenome.Hsapiens.1000genomes.hs37d5")
install.packages("C:/Users/12543/Downloads", repos = NULL, type = "source")
install.packages("D:/LenovoSoftstore/R-4.3.1/library/BSgenome.Hsapiens.1000genomes.hs37d5_0.99.1.tar.gz", repos = NULL, type = "source")
install.packages("BiocManager")
install.packages("BiocManager")
BiocManager::install("SNPlocs.Hsapiens.dbSNP144.GRCh37")
install.packages("C:/Users/12543/Desktop/3.1chr转RS/CHR转SNP本地安装包/SNPlocs.Hsapiens.dbSNP144.GRCh37.zip", repos = NULL, type = "win.binary")
install.packages("D:/LenovoSoftstore/R-4.3.1/library/SNPlocs.Hsapiens.dbSNP144.GRCh37.zip", repos = NULL, type = "win.binary")
install.packages("D:/LenovoSoftstore/R-4.3.1/library/SNPlocs.Hsapiens.dbSNP144.GRCh37.zip", repos = NULL, type = "win.binary")
install.packages("D:/LenovoSoftstore/R-4.3.1/library/SNPlocs.Hsapiens.dbSNP144.GRCh37.zip", repos = NULL, type = "win.binary")
library(BSgenome.Hsapiens.1000genomes.hs37d5)
library(SNPlocs.Hsapiens.dbSNP144.GRCh37)
library(MungeSumstats)
library(data.table)
eduAttainOkbayPth <- system.file("extdata", "data1.csv",
package = "MungeSumstats")
data <- fread(eduAttainOkbayPth)
eduAttainOkbayPth <- system.file("extdata", "UKBB.csv",
package = "MungeSumstats")
data <- fread(eduAttainOkbayPth)
library(dplyr)
new_data <- data %>%
rename( A1 = V4, A2 = V5 , P = ES , EAF =V8 , Beta = LP)
##1.提取暴露[exp]数据------
exp <- extract_instruments(outcomes = 'ieu-b-5106',p1 = 5e-20,
clump = TRUE,r2 = 0.01,kb = 10000,
access_token = NULL)#提取暴露数据
library(TwoSampleMR)
##1.提取暴露[exp]数据------
exp <- extract_instruments(outcomes = 'ieu-b-5106',p1 = 5e-20,
clump = TRUE,r2 = 0.01,kb = 10000,
access_token = NULL)#提取暴露数据
##2.提取结局[out]数据------
out<-extract_outcome_data(snps = exp$SNP,
outcomes = 'ebi-a-GCST90014022',
proxies=TRUE)#提取结局数据
##3.等位基因对齐 + 去除回文序列------
####action = 2，意味着剔除中等等位基因频率的回文序列
####action = 3，意味着剔除所有的回文序列
dat<-harmonise_data(exp,out,action = 2)
##如上数据缺失算不出F值后，使用b^2/se^2
dat$`F`<-((dat$beta.exposure)^2)/((dat$se.exposure)^2)
##5.MR分析------
#查看TwoSample内置的MR方法
mr_method_list()
results<-mr(dat)
results
##7.异质性检验------
#异质性分析啦-p>0.05表明没有异质性，p<0.05表明有异质性,有异质性，但是我使用了随机效应模型
heterogeneity <- mr_heterogeneity(dat)
heterogeneity
#if (!require("devtools")) { install.packages("devtools") } else {}
#devtools::install_github("rondolab/MR-PRESSO")
library(MRPRESSO)   #1水平多效性  2.挑出离群snp   3.计算跳出  /不跳出 离群snp 结果是否差异
mr_presso(BetaOutcome = 'beta.outcome',
BetaExposure = 'beta.exposure',
SdOutcome = 'se.outcome',
SdExposure = 'se.exposure',
data = dat, OUTLIERtest = TRUE,
DISTORTIONtest = TRUE, SignifThreshold = 0.05, NbDistribution = 1000, seed = NULL)
write.csv(dat,'dat.csv')
getwd()
library(vroom)
dat<-vroom('dat.csv')
##5.MR分析------
#查看TwoSample内置的MR方法
mr_method_list()
results<-mr(dat)
results
write.csv(results,'menarche_results.csv')
##7.异质性检验------
#异质性分析啦-p>0.05表明没有异质性，p<0.05表明有异质性,有异质性，但是我使用了随机效应模型
heterogeneity <- mr_heterogeneity(dat)
heterogeneity
write.csv(heterogeneity,'menarche_heterogeneity.csv')
##8.水平多效性【MR-Egger截距】---p<0.05有水平多效性，p>0.05没有水平多效性
pleio <- mr_pleiotropy_test(dat)
pleio
write.csv(pleio,'menarche_pleio.csv')
library(ggsci)
library(ggplot2)
##10.留一法+作图------
###逐个剔除每个SNP，计算剩余SNP的效应
single <- mr_leaveoneout(dat)
mr_leaveoneout_plot(single)
library(TwoSampleMR)
#save(dat,file="dat.rda")
##excel里操作
#①phenoscanner【二三假设】
#②如缺失暴露的样本量，需导出补齐samplesize.exposure
dat <- read.csv('dat.csv')
##5.MR分析------
#查看TwoSample内置的MR方法
mr_method_list()
results<-mr(dat)
results
##7.异质性检验------
#异质性分析啦-p>0.05表明没有异质性，p<0.05表明有异质性,有异质性，但是我使用了随机效应模型
heterogeneity <- mr_heterogeneity(dat)
heterogeneity
#if (!require("devtools")) { install.packages("devtools") } else {}
#devtools::install_github("rondolab/MR-PRESSO")
library(MRPRESSO)   #1水平多效性  2.挑出离群snp   3.计算跳出  /不跳出 离群snp 结果是否差异
mr_presso(BetaOutcome = 'beta.outcome',
BetaExposure = 'beta.exposure',
SdOutcome = 'se.outcome',
SdExposure = 'se.exposure',
data = dat, OUTLIERtest = TRUE,
DISTORTIONtest = TRUE, SignifThreshold = 0.05, NbDistribution = 1000, seed = NULL)
getwd()
View(dat)
setwd("D:/MR/endometriosis/E-selectin_finn-b-N14_ENDOMETRIOSIS")
library(TwoSampleMR)
dat <- read.csv('prot-a-2665_E_selectin_finn-b-N14_ENDOMETRIOSIS_dat.csv')
results<-mr(dat)
results
##10.留一法+作图------
###逐个剔除每个SNP，计算剩余SNP的效应
single <- mr_leaveoneout(dat)
mr_leaveoneout_plot(single)
##导出图
##11.散点图------
# 散点图
scatter_plot<-mr_scatter_plot(results,dat)[[1]]+
scale_color_lancet()+
scale_fill_lancet()+
theme_bw()
# 更改散点图的点的颜色
scatter_plot[["layers"]][[3]][["aes_params"]]$colour<-"black"
# 更改散点图的点的透明度
scatter_plot[["layers"]][[3]][["aes_params"]]$alpha<-0.5
library(ggsci)
library(ggplot2)
##导出图
##11.散点图------
# 散点图
scatter_plot<-mr_scatter_plot(results,dat)[[1]]+
scale_color_lancet()+
scale_fill_lancet()+
theme_bw()
# 更改散点图的点的颜色
scatter_plot[["layers"]][[3]][["aes_params"]]$colour<-"black"
# 更改散点图的点的透明度
scatter_plot[["layers"]][[3]][["aes_params"]]$alpha<-0.5
# 更改散点图的横坐标名字
scatter_plot[["labels"]][["x"]]<-"SNP effect on E_selectin"
# 更改散点图的纵坐标名字
scatter_plot[["labels"]][["y"]]<-"SNP effect on Endometriosis"
scatter_plot
##12.森林图---------
forest_plot<-mr_forest_plot(mr_singlesnp(dat))[[1]]+
scale_color_lancet()+
scale_fill_lancet()+
theme_bw()+
theme(legend.position = 'none')
forest_plot
##13.漏斗图------
funnel_plot<-mr_funnel_plot(mr_singlesnp(dat,all_method=c("mr_egger_regression","mr_weighted_median","mr_ivw","mr_simple_mode","mr_weighted_mode")))[[1]]+
theme_bw()+
scale_color_lancet()+
scale_fill_lancet()
# 更改漏斗图点的颜色
funnel_plot[["layers"]][[1]][["aes_params"]]$colour <- "black"
# 更改点的透明度
funnel_plot[["layers"]][[1]][["aes_params"]]$alpha<-0.5
funnel_plot
