有何AI與醫(yī)學(xué):孟德爾隨機(jī)化準(zhǔn)備從大家問題多地方以及簡單介紹寫起來,然后補(bǔ)全整個(gè)分析流程!歡迎轉(zhuǎn)發(fā)關(guān)注哦!
之前的文章介紹了非線性孟德隨機(jī)化的基本原理姥饰,以及多基因風(fēng)險(xiǎn)得分的計(jì)算(GRS),今天來介紹使用nlmr-R包完成非線孟德爾隨機(jī)化分析孝治。
nlmr包實(shí)踐
1#library(devtools)
2#install_github("jrs95/nlmr")
3library(nlmr)
4###?IV?(g),?exposure?(x)?&?outcome?(y)
5epsx?<-?rexp(10000)#生成10000個(gè)均值為1的指數(shù)分布隨機(jī)數(shù)
6u?<-?runif(10000,?0,?1)#最小0最大1
7g?<-?rbinom(10000,?2,?0.3)#0?1?2對應(yīng)三個(gè)基因型
8epsy?<-?rnorm(10000)#均值為0方差為1
9ag?<-?0.25
10x?<-?1?+?ag?*?g?+?u?+?epsx#這里可以找最行的GRS算
11y?<-?0.15?*?x^2?+?0.8?*?u?+?epsy#這里直接是我們對應(yīng)的結(jié)果,疾病就二分類就行列粪,family對應(yīng)binomial
12
13###?Covariates?(c)
14c1?=?rnorm(10000)#模擬連續(xù)協(xié)變量,例如年齡
15c2?=?rnorm(10000)#同上
16c3?=?rbinom(10000,?2,?0.33)#模擬離散協(xié)變量谈飒,例如吸煙
17c?=?data.frame(c1?=?c1,?c2?=?c2,?c3?=?as.factor(c3))
18
19###?Analyses
20#通過使用回歸將分?jǐn)?shù)多項(xiàng)式模型擬合到局部平均因果效應(yīng)來進(jìn)行工具變量分析岂座,光滑曲線。
21#y結(jié)局x暴露g基因型c協(xié)變量舉證杭措,family方法费什,q分段數(shù),d擬合次數(shù)次1或者2默認(rèn)1不重要手素,ci置信區(qū)間
22fp?=?fracpoly_mr(y,?x,?g,?c,?family?="gaussian",?q?=?10,?d?="both",?ci?="model_se",?fig?=?TRUE)
23summary(fp)
24#分段節(jié)點(diǎn)式
25plm?=?piecewise_mr(y,?x,?g,?c,?family?="gaussian",?q?=?10,?nboot?=?50,?fig?=?TRUE)
26summary(plm)
markdown結(jié)果
1#library(devtools)
2#install_github("jrs95/nlmr")
3library(nlmr)
4###?IV?(g),?exposure?(x)?&?outcome?(y)
5epsx?<-?rexp(10000)
6u?<-?runif(10000,?0,?1)
7g?<-?rbinom(10000,?2,?0.3)
8epsy?<-?rnorm(10000)
9ag?<-?0.25
10x?<-?1?+?ag?*?g?+?u?+?epsx
11y?<-?0.15?*?x^2?+?0.8?*?u?+?epsy
12
13#
##?Covariates?(c)
14c1?=?rnorm(10000)
15c2?=?rnorm(10000)
16c3?=?rbinom(10000,?2,?0.33)
17c?=?data.frame(c1?=?c1,?c2?=?c2,?c3?=?as.factor(c3))
18
19#
##?Analyses
20fp?=?fracpoly_mr(y,?x,?g,?c,?family?=?"gaussian",?q?=?10,?d?=?"both",?ci?=?"model_se",?fig?=?TRUE)
21summary(fp)
22##?Call:?fracpoly_mr
23##?
24##?Number?of?individuals:?10000;?Quantiles:?10;?95%CI:?Model?based?SEs
25##?
26##?Powers:?2
27##?
28##?Coefficients:
29##???Estimate?Std.?Error?95%CI?Lower?95%CI?Upper???p.value????
30##?2?0.144987???0.012498????0.120491??????0.1695?<?2.2e-16?***
31##?---
32##?Signif.?codes:??0?'***'?0.001?'**'?0.01?'*'?0.05?'.'?0.1?'?'?1
33##?
34##?Non-linearity?tests
35##?Fractional?polynomial?degree?p-value:?0.924
36##?Fractional?polynomial?non-linearity?p-value:?0.000535
37##?Quadratic?p-value:?0.000532
38##?Cochran?Q?p-value:?0.103
39##?
40##?Heterogeneity?tests
41##?Cochran?Q?p-value:?0.937
42##?Trend?p-value:?0.215
43
44plm?=?piecewise_mr(y,?x,?g,?c,?family?=?"gaussian",?q?=?10,?nboot?=?50,?fig?=?TRUE)
45summary(plm)
46##?Call:?piecewise_mr
47##?
48##?Number?of?individuals:?10000;?Quantiles:?10;?Number?of?bootstrap?replications:?50
49##?
50##?LACE:
51##?????Estimate?Std.?Error?95%CI?Lower?95%CI?Upper???p.value????
52##?1???0.370222???0.194433???-0.010866??????0.7513?0.0568951?.??
53##?2???0.324060???0.194220???-0.056611??????0.7047?0.0952126?.??
54##?3???0.728189???0.196023????0.343984??????1.1124?0.0002034?***
55##?4???0.674931???0.191170????0.300237??????1.0496?0.0004147?***
56##?5???0.767218???0.207506????0.360507??????1.1739?0.0002179?***
57##?6???0.766719???0.191298????0.391775??????1.1417?6.124e-05?***
58##?7???0.673263???0.190727????0.299438??????1.0471?0.0004156?***
59##?8???0.888327???0.199124????0.498044??????1.2786?8.151e-06?***
60##?9???1.096948???0.196101????0.712590??????1.4813?2.222e-08?***
61##?10??1.389435???0.416601????0.572897??????2.2060?0.0008525?***
62##?---
63##?Signif.?codes:??0?'***'?0.001?'**'?0.01?'*'?0.05?'.'?0.1?'?'?1
64##?
65##?Non-linearity?tests
66##?Quadratic?p-value:?0.000532
67##?Cochran?Q?p-value:?0.103
68##?
69##?Heterogeneity?tests
70##?Cochran?Q?p-value:?0.937
71##?Trend?p-value:?0.215