聊聊Nhanes数据怎么做中介效应分析

最近不少粉丝私信问我NANESA数据怎么做中介效应分析，我也查了10余篇文献，全都是普通中介分析的做法，没有涉及加权。
先来看这篇8分的文章
Association of volatile organic compounds co-exposure with bone health indicators and potential mediators（共同接触挥发性有机化合物与骨骼健康指标和潜在介质的关系）

这个研究是探索挥发性有机物的个体代谢产物和组合代谢产物之间的关系尿液中的（mVOCs），代表VOCs的内部暴露水平，具有骨密度（BMD），骨质疏松症（OP）和骨折以及潜在的介质。
咱们来看下它的统计学部分

看到了把，使用的是mediation包，这是个常用的中介分析的包，我已经有多篇文章进行介绍了。

再来看一篇高分的，

咱们看下统计学部分，也是使用mediation包

这篇也是一样，Prospective evaluation of cardiovascular risk andmortality in patients with psoriasis: A population-base study（银屑病患者心血管风险和死亡率的前瞻性评估：一项基于人群的研究）

研究的是银屑病和冠心病发病的中介因素，咱们来看下统计学部分，这篇文章介绍缺失超过30%的变量被删除，以后咱们也可以这样做。采用链式插补的方法进行插补，咱们在文章《一篇来自JAMA（IF=51）的NHANES数据插补方法》已经介绍了链式插补，看来高分文章都喜欢链式插补。
咱们看下它的统计学部分

使用的是mma包，这个包我查过了和mediation包差不多。

最后来介绍一篇7分的文献，为什么放在最后介绍呢，因为它是个统计专业的杂志，应该可靠性很强，而且它还附带了代码

作者后面还很牛的表示，它的统计经过了大佬的肯定，让咱们不敢怀疑。

最后作者附带了代码，虽然作者扯到了什么前门公式，但是我看了一下并没有涉及权重，和咱们既往介绍的方法很像。代码如下

set.seed(12345)
#Set variables
d$exp <- d$pain #exposure
d$med <- d$opioid #mediator
d$out <- d$mortality #outcome
#Create temporary exposure and mediator
d$exptemp <- d$exp
d$medtemp <- d$med
#Run a logistic regression model (exposure -> mediator)
reg1 <- glm(med ~ as.factor(exptemp) + covariates, data=d, family=binomial(link="logit"))
#Run a logistic regression model (mediator -> outcome)
reg2 <- glm(out ~ as.factor(medtemp) + as.factor(exp) + covariates, data=d, family=binomial(link="logit"))
#Create a potential outcome model (exposure -> mediator)
#Create combined datasets: 1) everyone set to exp=0, and 2) everyone set to exp=1
levelsofexp <- unique(d$exp)
d1 <- d
d2 <- d
d1$expstar <- levelsofexp[1] #exp=0
d2$expstar <- levelsofexp[2] #exp=1
newmyd1 <- rbind(d1, d2)
N1 <- nrow(newmyd1)
newmyd1$exptemp <- newmyd1$expstar
newmyd1$med_po <- rbinom(N1, size = 1, prob=1/(1+exp(-predict(reg1,newdata=newmyd1))))
reg1_po <- glm(med_po ~ as.factor(exptemp), data=newmyd1, family=binomial(link="logit"))
#Create a potential outcome model (mediator -> outcome)
#Create combined datasets: 1) everyone set to med=0, and 2) everyone set to med=1
levelsofmed<-unique(d$med)
d3 <- d
d4 <- d
d3$medstar<- levelsofmed[1] #med=0
25
d4$medstar<- levelsofmed[2] #med=1
newmyd2<-rbind(d3, d4)
N2<-nrow(newmyd2)
newmyd2$medtemp <- newmyd2$medstar
newmyd2$out_po <- rbinom(N2, size = 1, prob=1/(1+exp(-predict(reg2, newdata=newmyd2))))
reg2_po <- glm(out_po~ as.factor(medtemp), data=newmyd2, family=binomial(link="logit"))
##Front-door adjustment
#Y=death, M=Opioid, X=Pain
#Path-specific front-door effect
newmyd1$medtemp <- newmyd1$med_po
newmyd1$po <- rbinom(N1, size = 1, prob=1/(1+exp(-predict(reg2, newdata=newmyd1))))
summary(glm(po ~ as.factor(exptemp), data=newmyd1, family=binomial(link="logit")))
#Bootstrap can be used to obtain 95% confidence intervals.

本期就这样结束了，祝大家新年快乐。