Radial MR:两样本孟德尔随机化的radial plot



Radial MR: A package for implementing radial inverse variance weighted and MR-Egger methods

导入数据
library(RadialMR)
View(data_radial)
1. 分析
1.1 egger_radial

Example using format_radial data

ldl.dat <- data_radial[data_radial[,10]<5*10^-8,]

ldl.fdat <- format_radial(ldl.dat[,6], ldl.dat[,9],
                          ldl.dat[,15], ldl.dat[,21],
                          ldl.dat[,1])

egger_radial(ldl.fdat, 0.05, 1, TRUE)

# Radial MR-Egger

#               Estimate Std. Error   t value
# (Intercept) -0.4303552  0.3244151 -1.326557
# Wj           0.6129126  0.1109369  5.524877
#                 Pr(>|t|)
# (Intercept) 1.884296e-01
# Wj          3.988423e-07

# Residual standard error: 1.537 on 80 degrees of freedom

# F-statistic: 30.52 on 1 and 80 DF, p-value: 3.99e-07
# Q-Statistic for heterogeneity: 188.9285 on 80 DF , p-value: 8.445607e-11

# Outliers detected 

Example using TwoSampleMR format data
Example with one exposure-outcome pair

bmi_exp_dat <- TwoSampleMR::extract_instruments(outcomes = 'ieu-a-2')
chd_out_dat <- TwoSampleMR::extract_outcome_data(
  snps = bmi_exp_dat$SNP,
  outcomes = 'ieu-a-7')
tsmrdat <- TwoSampleMR::harmonise_data(exposure_dat = bmi_exp_dat,
                                       outcome_dat = chd_out_dat)
egger_radial(r_input = tsmrdat, alpha = 0.05,
             weights = 1, summary = TRUE)

# Radial MR-Egger

#               Estimate Std. Error   t value
# (Intercept) -0.4523432  0.3949287 -1.145379
# Wj           0.6071217  0.1525641  3.979452
#                 Pr(>|t|)
# (Intercept) 0.2555995047
# Wj          0.0001550786

# Residual standard error: 1.354 on 77 degrees of freedom

# F-statistic: 15.84 on 1 and 77 DF, p-value: 0.000155
# Q-Statistic for heterogeneity: 141.2444 on 77 DF , p-value: 1.129761e-05

#  Outliers detected 

Example using MendelianRandomization format data

dat <- data_radial[data_radial[,10] < 5e-8,]
mrdat <- MendelianRandomization::mr_input(bx = dat$ldlcbeta,
                                          bxse = dat$ldlcse,
                                          by = dat$chdbeta,
                                          byse = dat$chdse,
                                          snps = dat$rsid)
egger_radial(r_input = mrdat, alpha = 0.05,
             weights = 1, summary = TRUE)

# Radial MR-Egger

#               Estimate Std. Error
# (Intercept) -0.4303552  0.3244151
# Wj           0.6129126  0.1109369
#               t value     Pr(>|t|)
# (Intercept) -1.326557 1.884296e-01
# Wj           5.524877 3.988423e-07

# Residual standard error: 1.537 on 80 degrees of freedom

# F-statistic: 30.52 on 1 and 80 DF, p-value: 3.99e-07
# Q-Statistic for heterogeneity: 188.9285 on 80 DF , p-value: 8.445607e-11

#  Outliers detected 
1.2 format_radial
ldl.dat <- data_radial[data_radial[,10]<5*10^-8,]
ldl.fdat <- format_radial(ldl.dat[,6], ldl.dat[,9],
                          ldl.dat[,15], ldl.dat[,21], ldl.dat[,1])
head(ldl.fdat)
#          SNP beta.exposure beta.outcome se.exposure se.outcome
# 1 rs10903129        -0.033       -0.012 0.003692528 0.01366904
# 2  rs1998013        -0.380       -0.150 0.021953470 0.09647707
# 3  rs4587594        -0.049        0.017 0.003842235 0.01509245
# 4  rs6603981         0.034        0.012 0.004444080 0.01698989
# 5   rs646776         0.160        0.094 0.004375672 0.01724356
# 6  rs1010167        -0.025       -0.028 0.003969023 0.01897288
class(ldl.fdat)
# [1] "data.frame" "rmr_format"
1.3 ivw_radial
# Example using format_radial data
ldl.dat <- data_radial[data_radial[,10]<5*10^-8,]
ldl.fdat <- format_radial(ldl.dat[,6], ldl.dat[,9],
                          ldl.dat[,15], ldl.dat[,21], ldl.dat[,1])
ivw_radial(ldl.fdat, 0.05, 1, 0.0001, TRUE)

# Example using TwoSampleMR format data
## Not run: 
# Example with one exposure-outcome pair
bmi_exp_dat <- TwoSampleMR::extract_instruments(outcomes = 'ieu-a-2')
chd_out_dat <- TwoSampleMR::extract_outcome_data(
                               snps = bmi_exp_dat$SNP,
                               outcomes = 'ieu-a-7')
tsmrdat <- TwoSampleMR::harmonise_data(exposure_dat = bmi_exp_dat,
                                   outcome_dat = chd_out_dat)
ivw_radial(r_input = tsmrdat, alpha = 0.05,
           weights = 1, tol = 0.0001, summary = TRUE)

## End(Not run)

# Example using MendelianRandomization format data
dat <- data_radial[data_radial[,10] < 5e-8,]
mrdat <- MendelianRandomization::mr_input(bx = dat$ldlcbeta,
                                          bxse = dat$ldlcse,
                                          by = dat$chdbeta,
                                          byse = dat$chdse,
                                          snps = dat$rsid)
ivw_radial(r_input = mrdat, alpha = 0.05,
           weights = 1, tol = 0.0001, summary = TRUE)
2. mrinput_to_rmr_format

Convert an object of class MRInput from the MendelianRandomization package to the RadialMR rmr_format class

dat <- data_radial[data_radial[,10] < 5e-8,]
dat <- MendelianRandomization::mr_input(bx = dat$ldlcbeta,
                                        bxse = dat$ldlcse,
                                        by = dat$chdbeta,
                                        byse = dat$chdse,
                                        snps = dat$rsid)
dat <- mrinput_to_rmr_format(dat)
head(dat)
#          SNP beta.exposure beta.outcome se.exposure se.outcome
# 1 rs10903129        -0.033       -0.012 0.003692528 0.01366904
# 2  rs1998013        -0.380       -0.150 0.021953470 0.09647707
# 3  rs4587594        -0.049        0.017 0.003842235 0.01509245
# 4  rs6603981         0.034        0.012 0.004444080 0.01698989
# 5   rs646776         0.160        0.094 0.004375672 0.01724356
# 6  rs1010167        -0.025       -0.028 0.003969023 0.01897288
class(dat)
# [1] "data.frame" "rmr_format"
3. 绘图

plotly_radial {RadialMR}

ldl.dat <- data_radial[data_radial[,10]<5*10^-8,]
ldl.fdat <- format_radial(ldl.dat[,6], ldl.dat[,9],
                          ldl.dat[,15], ldl.dat[,21], ldl.dat[,1])
ivw.object <- ivw_radial(ldl.fdat, 0.05, 1, 0.0001, TRUE)

# Radial IVW

#               Estimate  Std.Error   t value     Pr(>|t|)
# Effect (1st) 0.4874900 0.05830409  8.361163 6.210273e-17
# Iterative    0.4873205 0.05827885  8.361874 6.172955e-17
# Exact (FE)   0.4958973 0.03804168 13.035630 7.673061e-39
# Exact (RE)   0.4910400 0.05451396  9.007602 7.682743e-14


# Residual standard error: 1.544 on 81 degrees of freedom

# F-statistic: 69.91 on 1 and 81 DF, p-value: 1.46e-12
# Q-Statistic for heterogeneity: 193.0843 on 81 DF , p-value: 3.827332e-11

#  Outliers detected 
# Number of iterations = 3

plotly_radial(ivw.object)
ldl.dat <- data_radial[data_radial[,10]<5*10^-8,]
ldl.fdat <- format_radial(ldl.dat[,6], ldl.dat[,9],
                          ldl.dat[,15], ldl.dat[,21], ldl.dat[,1])
ivw.object <- ivw_radial(ldl.fdat, 0.05, 1, 0.0001, TRUE)

# Radial IVW

#               Estimate  Std.Error   t value     Pr(>|t|)
# Effect (1st) 0.4874900 0.05830409  8.361163 6.210273e-17
# Iterative    0.4873205 0.05827885  8.361874 6.172955e-17
# Exact (FE)   0.4958973 0.03804168 13.035630 7.673061e-39
# Exact (RE)   0.4910400 0.05483473  8.954910 9.769963e-14


# Residual standard error: 1.544 on 81 degrees of freedom

# F-statistic: 69.91 on 1 and 81 DF, p-value: 1.46e-12
# Q-Statistic for heterogeneity: 193.0843 on 81 DF , p-value: 3.827332e-11

#  Outliers detected 
# Number of iterations = 3

plot_radial(ivw.object)
4. tsmr_to_rmr_format

Convert a data.frame containing a single exposure - outcome pair generated by TwoSampleMR::harmonise_data() to the RadialMR rmr_format class

## Not run: 
# Example with one exposure-outcome pair
bmi_exp_dat <- TwoSampleMR::extract_instruments(outcomes = 'ieu-a-2')
chd_out_dat <- TwoSampleMR::extract_outcome_data(
                               snps = bmi_exp_dat$SNP,
                               outcomes = 'ieu-a-7')
dat <- TwoSampleMR::harmonise_data(exposure_dat = bmi_exp_dat,
                                   outcome_dat = chd_out_dat)
dat <- tsmr_to_rmr_format(dat)
class(dat)
# [1] "data.frame" "rmr_format"
head(dat)
#          SNP beta.exposure beta.outcome se.exposure se.outcome
# 1  rs1000940        0.0184    -0.014538      0.0033  0.0098395
# 2 rs10132280       -0.0221    -0.012169      0.0033  0.0106640
# 3  rs1016287       -0.0228    -0.014087      0.0033  0.0105738
# 4 rs10182181        0.0309     0.018295      0.0029  0.0092790
# 5 rs10733682       -0.0188    -0.004541      0.0030  0.0095245
# 6 rs10840100        0.0206     0.014599      0.0030  0.0095234

## End(Not run)
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
禁止转载,如需转载请通过简信或评论联系作者。
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 212,332评论 6 493
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 90,508评论 3 385
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 157,812评论 0 348
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 56,607评论 1 284
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 65,728评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,919评论 1 290
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 39,071评论 3 410
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,802评论 0 268
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,256评论 1 303
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,576评论 2 327
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,712评论 1 341
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,389评论 4 332
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 40,032评论 3 316
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,798评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,026评论 1 266
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,473评论 2 360
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,606评论 2 350

推荐阅读更多精彩内容