孟德尔随机化分析流程

这个只是留给自己看

a <- read.table('20002_1464.gwas.imputed_v3.both_sexes.tsv', header = T)
b <- subset(a, pval<10e-08)
library('TwoSampleMR')
write.csv(b, file = 'exposure.csv')


bmi <- system.file("/home/data/t050312/aa/exposure.csv", package = "TwoSampleMR")

bmi_exp_dat_clumped<-read_exposure_data(filename = '/home/data/t050312/aa/exposure.csv',sep = ",",snp_col = "SNP",beta_col = "beta",se_col = "se",effect_allele_col = "effect_allele",other_allele_col = "other_allele",eaf_col = "eaf",clump = TRUE)

write.csv(bmi_exp_dat_clumped, file = "/home/data/t050312/aa/outcome.csv")

outcome_dat<-read_outcome_data(snps = bmi_exp_dat_clumped$SNP,filename = "/home/data/t050312/aa/MR/outplose/abnormalities_of_heartbeat/outcome2.txt",sep = "\t",snp_col = "SNP",beta_col = "beta",se_col = "se",effect_allele_col = "effect_allele",other_allele_col = "other_allele",pval_col = "pval")
dat<-harmonise_data(exposure_dat = bmi_exp_dat_clumped,outcome_dat = outcome_dat)
write.csv(dat, file = "/home/data/t050312/aa/harmonised_data.csv")


mr(dat)
generate_odds_ratios(mr_res = mr(dat))
mr(dat,method_list = c("mr_ivw","mr_egger_regression","mr_weighted_median"))
mr_scatter_plot(mr_results = mr(dat,method_list = c("mr_ivw","mr_egger_regression","mr_weighted_median")),dat)
mr_heterogeneity(dat)

mr_funnel_plot(singlesnp_results = mr_singlesnp(dat))
mr_pleiotropy_test(dat)    #多效性检测
mr_leaveoneout_plot(leaveoneout_results = mr_leaveoneout(dat))
















第二次
a <- read.table('20002_1464.gwas.imputed_v3.both_sexes.tsv', header = T)
b <- subset(a, pval<5e-08)
library('TwoSampleMR')


exposure <- extract_instruments(outcomes = 'ukb-e-20128_AFR')


exposure <- system.file("/home/zhang123/LBD/Sex_hormone_binding_globulin_levels_adjusted_for_BMI.csv", package = "TwoSampleMR")

outcome_dat<-read_outcome_data(snps = exposure$SNP,filename = "/home/zhang123/LBD/GCST90001390_buildGRCh38.tsv",sep = "\t",snp_col = "variant_id",beta_col = "beta",se_col = "standard_error",effect_allele_col = "effect_allele",other_allele_col = "other_allele",pval_col = "p_value", eaf_col = "effect_allele_frequency")
dat<-harmonise_data(exposure_dat = exposure,outcome_dat = outcome_dat)
write.csv(dat, file = "/home/zhang123/LBD/harmonise/Sex_hormone_binding_globulin_levels_adjusted_for_BMI.csv")



generate_odds_ratios(mr_res = mr(dat))

mr_heterogeneity(dat)
mr_pleiotropy_test(dat)    #多效性检测

pdf("/home/zhang123/LBD/plot/Feeling_tense.pdf", onefile = TRUE)
mr_scatter_plot(mr_results = mr(dat,method_list = c("mr_ivw","mr_egger_regression","mr_weighted_median")),dat)

mr_funnel_plot(singlesnp_results = mr_singlesnp(dat))

mr_leaveoneout_plot(leaveoneout_results = mr_leaveoneout(dat))

dev.off()



2023.06.22记录：

a <- read.table('20002_1464.gwas.imputed_v3.both_sexes.tsv', header = T)
b <- subset(a, pval<5e-08)
library('TwoSampleMR')


exposure <- extract_instruments(outcomes = 'ukb-b-2988')



outcome_dat<-read_outcome_data(snps = exposure$SNP,filename = "/home/zhang123/LBD/GCST90001390_buildGRCh38.tsv",sep = "\t",snp_col = "variant_id",beta_col = "beta",se_col = "standard_error",effect_allele_col = "effect_allele",other_allele_col = "other_allele",pval_col = "p_value", eaf_col = "effect_allele_frequency")
dat<-harmonise_data(exposure_dat = exposure,outcome_dat = outcome_dat)
write.csv(dat, file = "/home/zhang123/LBD/harmonise/Intelligence.csv")



generate_odds_ratios(mr_res = mr(dat))

mr_heterogeneity(dat)
mr_pleiotropy_test(dat)    #多效性检测

pdf("/home/zhang123/LBD/plot/Intelligence.pdf", onefile = TRUE)
mr_scatter_plot(mr_results = mr(dat,method_list = c("mr_ivw","mr_egger_regression","mr_weighted_median")),dat)

mr_funnel_plot(singlesnp_results = mr_singlesnp(dat))

mr_leaveoneout_plot(leaveoneout_results = mr_leaveoneout(dat))

dev.off()