孟德尔随机化:根据孟德尔遗传规律,亲代的等位基因随机分配给子代,此过程相当于随机对照研究(RCT)的随机分组过程:不受混杂因素(社会地位、行为等)的影响;满足时间顺序合理性(遗传变异继承于父母,且保持不变)
基本思想:利用与暴露因素具有强相关的遗传变异作为工作变量,借此推断暴露因素与研究结局之间的因果效应。
现有全基因组关联研究(GWAS)数据大部分来源于欧洲,少量来源于亚洲、非洲、混合人种(暴露与结局最好来自同一个地方,消除人群异质性)
暴露与结局数据不存在或存在少量可接受范围内的样本重叠(最好不用同一库的)
Rs号(pos,chr可以转化为rs号)
effect_allele(A1;效应等位基因)
other_allele(A2)
beta(OR;效应值):表示SNP对表型起到的作用(结局是连续型变量用beta,大于0则增加风险;结局是二分类变量用OR,OR大于1则增加风险)
se
P
EAF(Freq;MAF;效应等位基因频率)
注:beta与OR可以相互转化:beta <- ln(OR)
P、beta、se知道其中两个可以计算另一个:se <- abs(beta/qnorm(p/2,lower.tail=F)
p<5e-08(视条件而定);F>10代表强相关; 注:p最低可放宽到5e-06(可能存在弱工具变量) ; P值越小,F越大
Clumping:LD r2<0.001(r2越小,SNP越少),kb=10000(视条件而定)#连锁不平衡r2小于0.001(越小SNP之间越独立),kb(指连锁不平衡的区域长度,在遗传学上我们认为在染色体上距离很近的遗传位点通常是“捆绑”在一起遗传给后代的)
R2 <- 2*(1-MAF)*MAF*(beta)2 #计算每个SNP的R2将其相加,带入下面公式
F=(R2/(1-R2))((n-k-1)/k)# n为sample size , k为最终的SNP个数
PhenoScanner (cam.ac.uk)#假设2和假设3(将假设1筛选出的SNP逐个输入该网站验证)
解决方案:1.增加样本量 2.增加表型解释度:相对于单个遗传变异,多个遗传变异能解释更大比例的表型变异
解决方案:使用生物学功能明确的遗传变异作为工具变量
解决方案:使用生物学功能明确的遗传变异作为工具变量
解决方案:纳入具有相同遗传背景的人群开展遗传关联研究
dat <- harmonise_data(exposure_dat = RI_exp_dat_clumped,outcome_dat = outcome_dat)
mRnd: Power calculations for Mendelian Randomization (cnsgenomics.com)
将数据下载到桌面,打开RStudio,导入暴露数据exposure.txt:
setwd("C:/Users/76325/Desktop")
exposure1 <- read.table("exposure.txt",header=T)
exposure2 <- (exposure1,p<5e-08)#p为表中的变量名
write.csv(exposure2,file="exposure_RI.csv")
RI_exp_dat <- read_exposure_data(filename=exposure_RI,sep=",",snp_col="SNP",bata_col="bata",
se_col="se",effect_allele_col=""effect_allele,other_allele_col="other_allele",eaf_col=“eaf”,
pval_col="p")
RI_exp_dat_clumped <- clump_data(RI_exp_dat,clump_kb= 10000,clump_r2=0.01,clump_p1=1,
clump_p2=1,pop="EUR")#EUR表示样本来自欧洲
exposure_outcome <- merge(RI_exp_dat_clumped,PCOS_outcome,by.x="SNP",by.y="MarkerName")
outcome_dat <- read_outcome_data(snps=RI_exp_dat_clumped,filename=outcome.csv,sep=",",snp_col="SNP",bata_col="bata",
se_col="se",effect_allele_col=""effect_allele,other_allele_col="other_allele",eaf_col=“eaf”,
pval_col="p")
dat <- harmonise_data(exposure_dat = RI_exp_dat_clumped,outcome_dat = outcome_dat)
> mr(dat)#默认用五种方法分析
mr_method_list()#查看总共有多少种方法
#mr(dat,method_list=c("mr_ivw","mr_raps"))#指定方法分析
Analysing 'ieu-a-2' on 'ieu-a-7'
id.exposure id.outcome outcome
1 ieu-a-2 ieu-a-7 Coronary heart disease || id:ieu-a-7
2 ieu-a-2 ieu-a-7 Coronary heart disease || id:ieu-a-7
3 ieu-a-2 ieu-a-7 Coronary heart disease || id:ieu-a-7
4 ieu-a-2 ieu-a-7 Coronary heart disease || id:ieu-a-7
5 ieu-a-2 ieu-a-7 Coronary heart disease || id:ieu-a-7
exposure method nsnp b se
1 Body mass index || id:ieu-a-2 MR Egger 79 0.5024935 0.14396056
2 Body mass index || id:ieu-a-2 Weighted median 79 0.3870065 0.07201409
3 Body mass index || id:ieu-a-2 Inverse variance weighted 9 79 0.4459091 0.05898302
4 Body mass index || id:ieu-a-2 Simple mode 79 0.3401554 0.15562464
5 Body mass index || id:ieu-a-2 Weighted mode 70.3888249 0.09681892
pval
1 8.012590e-04 #在MR-Egger的intercept分析时,P必须大于0.05,否则IVW结果不可靠
2 7.699254e-08
3 4.032020e-14 #IVW的P值<0.05,其他的大于0.05也可以当作阳性结果,但其他方法的beta值与IVW方向必须一致
4 3.183437e-02 #若方向不一致,将p值逐渐缩小 5e-9
5 1.352008e-04
generate_odds_ratios(mr_res=mr(dat))
mr_scatter_plot(mr_results = mr(dat,method_list = c("mr_ivw","mr_egger_regression",
"mr_weighted_median")),dat)
> mr_heterogeneity(dat)#Q值小于0.05说明存在异质性
id.exposure id.outcome outcome
1 ieu-a-2 ieu-a-7 Coronary heart disease || id:ieu-a-7
2 ieu-a-2 ieu-a-7 Coronary heart disease || id:ieu-a-7
exposure method Q Q_df Q_pval
1 Body mass index || id:ieu-a-2 MR Egger 143.3046 77 6.841585e-06
2 Body mass index || id:ieu-a-2 Inverse variance weighted 143.6508 78 8.728420e-06
run_mr_presso(dat,NbDistribution = 5000)#数字越大,精度越高
mr_funnel_plot(singlesnp_results =mr_singlesnp(dat))
mr_pleiotropy_test(dat)#若p<0.05说明存在多效性(有截距),则不做了
mr_leaveoneout(dat)
mr_leaveoneout_plot(leaveoneout_results = mr_leaveoneout(dat))
library(devtools)
devtools::install_github("MRCIEU/TwoSampleMR")
library(TwoSampleMR)
exposure <- extract_instruments(outcomes = "ieu-a-2",p1=5e-08)
outcome <- extract_outcome_data(snps = exposure$SNP,outcomes = "ieu-a-7")
dat <- harmonise_data(exposure_dat = exposure,outcome_dat = outcome)
mr(dat)
#异质性检测
mr_heterogeneity(dat)#Q值大于0.05说明不存在异质性
run_mr_presso(dat,NbDistribution = 5000)
#异质性可视化
mr_funnel_plot(singlesnp_results =mr_singlesnp(dat))
#多效性检测
mr_pleiotropy_test(dat)#P<0.05则别做了
#leaveoneout analysis
mr_leaveoneout_plot(leaveoneout_results = mr_leaveoneout(dat))
#可视化
mr_scatter_plot(mr_results = mr(dat,method_list = c("mr_ivw","mr_egger_regression","mr_weighted_median")),dat)
Briefly describe the study design and the corresponding fundamental assumption
Exposure & Outcome
Association threshold,LD cutoff,exclude SNP associated with outcome,Harmonization to exclude palindromic and incompatible SNPs, F>10 ,and so on
MR estimations like IVW,weighted median, MR-Egger,RAPS,and so on; Sensitivity analysis like Cochran Q test(异质性检测),Egger-intercept test(多效性检测), leave-one-out analysis
原创声明:本文系作者授权腾讯云开发者社区发表,未经许可,不得转载。
如有侵权,请联系 cloudcommunity@tencent.com 删除。
原创声明:本文系作者授权腾讯云开发者社区发表,未经许可,不得转载。
如有侵权,请联系 cloudcommunity@tencent.com 删除。