赞
踩
- setwd("D:/R/MR")
- library(TwoSampleMR)
- #-----导入菌群的数据
- all_gut <- read.table('MBG.allHits.p1e4.txt',header = T)
- all_gut <- subset(all_gut,P.weightedSumZ<1e-05) #过滤一
- write.csv(all_gut,"exposure_all_gut.csv")
- #-----读取exposure
- exposure_data <- read_exposure_data(filename = "exposure_all_gut.csv", sep = ",", snp_col = "rsID",
- beta_col = "beta", se_col = "SE", phenotype_col = "bac",
- effect_allele_col = "eff.allele", other_allele_col = "ref.allele",
- chr_col = "chr", pos_col = "bp", clump = FALSE)
- #clump这一步需要联网
- exposure_data <- clump_data(exposure_data, clump_r2 = 0.001, pop = "EUR", clump_kb = 10000) #过滤二
- write.csv(exposure_data,"exposure_all_gut_clumped.csv")
- #-----读取outcome_data
- outcome_data <- read_outcome_data(filename = "outcome_eur_rsid_reformatted_new.csv",
- snps = exposure_data$SNP, snp_col = "SNP",sep = ",",beta_col = "BETA",
- se_col = "SE",effect_allele_col = "A1",other_allele_col = "A2",
- pval_col = "p",eaf_col = "FRQ",chr_col = "CHR",pos_col = "BP")
- #-----预处理数据
- dat <- harmonise_data(exposure_data,outcome_data)
- write.csv(dat,"dat_harmonised_gut_NC.csv")
- #-----自选方法进行MR分析
- res <- mr(dat,method_list = c("mr_ivw","mr_two_sample_ml","mr_egger_regression","mr_weighted_median","mr_weighted_mode"))
- head(res)
Copyright © 2003-2013 www.wpsshop.cn 版权所有,并保留所有权利。