조건부간접효과중 가장 인기 있는 model 7을 해보장
x<-rnorm(100)
w<- rnorm(100)
me<-rnorm(100)+ 2*x*w + x
y<-rnorm(100, 0,1) +me +x
co1<-rnorm(100)
d<-data.frame(x,me,w,y,co1)#개수 구하기
boot7_1<-function(xxx,mmm,www,yyy,d){
n<-sample(1:nrow(d),replace = T)
nnk<-d[n,]
nnk<-as.data.frame(nnk)
k1<-lm(nnk[,mmm]~ nnk[,xxx], data=nnk)
s1<-summary(k1)
coem<-s1$coefficients
eff<-as.data.frame(coem)
eff<-eff[nrow(eff),1]
k2<-lm(nnk[,yyy] ~ nnk[,xxx]+ nnk[,mmm], data = nnk)
s2<-summary(k2)
coem2<-s2$coefficients
eff2<-as.data.frame(coem2)
eff5<-eff2[nrow(eff2),1]
eff*eff5 #indirect
k3<-lm(nnk[,mmm] ~ nnk[,xxx]+ nnk[,www] + nnk[,xxx]*nnk[,www], data = nnk)
s3<-summary(k3)
coem3<-s3$coefficients
eff2_1<-as.data.frame(coem3)
eff7<-eff2_1[nrow(eff2_1),1]
nnk2<-nnk
nnk2[,www]<-ifelse(nnk2[,www]> mean(nnk2[,www])+sd(nnk2[,www]), nnk2[,www], NA)
nnk2<-na.omit(nnk2)
k3_1<-lm(nnk2[,yyy] ~ nnk2[,xxx], data = nnk2)
s3_1<-summary(k3_1)
coem3_1<-s3_1$coefficients
effhx<-as.data.frame(coem3_1)
effhx<-effhx[nrow(effhx),1]
k3_2<-lm(nnk2[,yyy] ~ nnk2[,xxx] + nnk2[,mmm], data = nnk2)
s3_2<-summary(k3_2)
coem3_2<-s3_2$coefficients
effhm<-as.data.frame(coem3_2)
effhm<-effhm[nrow(effhm),1]
effh<-effhx*effhm
nnk3<-nnk
nnk3[,www]<-ifelse(nnk3[,www]< mean(nnk3[,www]) - sd(nnk3[,www]), nnk3[,www], NA)
nnk3<-na.omit(nnk3)
k4_1<-lm(nnk3[,yyy] ~ nnk3[,xxx], data = nnk3)
s4_1<-summary(k4_1)
coem4_1<-s4_1$coefficients
efflx<-as.data.frame(coem4_1)
efflx<-efflx[nrow(efflx),1]
k4_2<-lm(nnk3[,yyy] ~ nnk3[,xxx] + nnk3[,mmm], data = nnk3)
s4_2<-summary(k4_2)
coem4_2<-s4_2$coefficients
efflm<-as.data.frame(coem4_2)
efflm<-efflm[nrow(efflm),1]
effl<-efflx*efflm
meff<-eff*eff5
moff<-eff7
index_mo<-meff*moff
efff<-c(meff, moff, index_mo, effh, effl)
efff<-matrix(efff, ncol = 5)
efff
}
boot7_1(1,2,3,4,d)## [,1] [,2] [,3] [,4] [,5]
## [1,] 0.79773 2.05952 1.642941 4.663347 -2.620161#boot결과 써보기 솔직히 너무 정보가 많아서 쓰기 귀찮아진다. 다른거 해야겠다.
boot7<-function(xxx,mmm,www,yyy,d,bootnum){
###estimate a*m
boot7_1<-function(xxx,mmm,www,yyy,d){
n<-sample(1:nrow(d),replace = T)
nnk<-d[n,]
nnk<-as.data.frame(nnk)
k1<-lm(nnk[,mmm]~ nnk[,xxx], data=nnk)
s1<-summary(k1)
coem<-s1$coefficients
eff<-as.data.frame(coem)
eff<-eff[nrow(eff),1]
k2<-lm(nnk[,yyy] ~ nnk[,xxx]+ nnk[,mmm], data = nnk)
s2<-summary(k2)
coem2<-s2$coefficients
eff2<-as.data.frame(coem2)
eff5<-eff2[nrow(eff2),1]
eff*eff5 #indirect
k3<-lm(nnk[,mmm] ~ nnk[,xxx]+ nnk[,www] + nnk[,xxx]*nnk[,www], data = nnk)
s3<-summary(k3)
coem3<-s3$coefficients
eff2_1<-as.data.frame(coem3)
eff7<-eff2_1[nrow(eff2_1),1]
nnk2<-nnk
nnk2[,www]<-ifelse(nnk2[,www]> mean(nnk2[,www])+sd(nnk2[,www]), nnk2[,www], NA)
nnk2<-na.omit(nnk2)
k3_1<-lm(nnk2[,yyy] ~ nnk2[,xxx], data = nnk2)
s3_1<-summary(k3_1)
coem3_1<-s3_1$coefficients
effhx<-as.data.frame(coem3_1)
effhx<-effhx[nrow(effhx),1]
k3_2<-lm(nnk2[,yyy] ~ nnk2[,xxx] + nnk2[,mmm], data = nnk2)
s3_2<-summary(k3_2)
coem3_2<-s3_2$coefficients
effhm<-as.data.frame(coem3_2)
effhm<-effhm[nrow(effhm),1]
effh<-effhx*effhm
nnk3<-nnk
nnk3[,www]<-ifelse(nnk3[,www]< mean(nnk3[,www]) - sd(nnk3[,www]), nnk3[,www], NA)
nnk3<-na.omit(nnk3)
k4_1<-lm(nnk3[,yyy] ~ nnk3[,xxx], data = nnk3)
s4_1<-summary(k4_1)
coem4_1<-s4_1$coefficients
efflx<-as.data.frame(coem4_1)
efflx<-efflx[nrow(efflx),1]
k4_2<-lm(nnk3[,yyy] ~ nnk3[,xxx] + nnk3[,mmm], data = nnk3)
s4_2<-summary(k4_2)
coem4_2<-s4_2$coefficients
efflm<-as.data.frame(coem4_2)
efflm<-efflm[nrow(efflm),1]
effl<-efflx*efflm
meff<-eff*eff5
moff<-eff7
index_mo<-meff*moff
efff<-c(meff, moff, index_mo, effh, effl)
efff<-matrix(efff, ncol = 5)
efff
}
k<-1
l<-matrix(rep(NA,bootnum*5),ncol = 5)
l<-as.data.frame(l)
repeat{
l[k,]<-boot7_1(xxx,mmm,www,yyy,d)
k<-k+1
if(k>=bootnum+1) break
}
estimates<-list(l)
ci1<-quantile(l[,1],probs = c(.001,0.01,0.05,0.10,0.90,0.95,0.99,.999))
ci2<-quantile(l[,2],probs = c(.001,0.01,0.05,0.10,0.90,0.95,0.99,.999))
ci3<-quantile(l[,3],probs = c(.001,0.01,0.05,0.10,0.90,0.95,0.99,.999))
ci4<-quantile(l[,4],probs = c(.001,0.01,0.05,0.10,0.90,0.95,0.99,.999))
ci5<-quantile(l[,5],probs = c(.001,0.01,0.05,0.10,0.90,0.95,0.99,.999))
r1<-summary(lm(d[,mmm]~ d[,xxx] + d[,www] + d[,xxx]*d[,www], data = d))
r11<-r1$adj.r.squared
r2<-summary(lm(d[,yyy]~ d[,xxx] + d[,mmm] + d[,www] + d[,xxx]*d[,www], data = d))
r12<-r2$adj.r.squared
kmkmkmkm<-list(c(mean(l[,1]),sd(l[,1])),ci1, c(mean(l[,2]),sd(l[,2])),ci2, c(mean(l[,3]),sd(l[,3])),ci3, c(mean(l[,4]),sd(l[,4])),ci4, c(mean(l[,5]),sd(l[,5])),ci5, r11,r12)
names(kmkmkmkm)<-c("mediation_mean_BootSE", "mediation_CI","moderation_mean_BootSE", "moderation_CI","index_mean_BootSE", "index_CI","+1sd_indirect_mean_BootSE", "+1sd_indirect_CI","-1sd_indirect_mean_BootSE", "-1sd_indirect_CI",
"mediation r^2", "total r^2")
kmkmkmkm
}
boot7(1,2,3,4,d,5000)## $mediation_mean_BootSE
## [1] 1.2735563 0.3653454
##
## $mediation_CI
## 0.1% 1% 5% 10% 90% 95% 99% 99.9%
## 0.2509987 0.4695889 0.6928572 0.8094277 1.7451325 1.8786458 2.1727903 2.4209386
##
## $moderation_mean_BootSE
## [1] 1.8592766 0.1394208
##
## $moderation_CI
## 0.1% 1% 5% 10% 90% 95% 99% 99.9%
## 1.330371 1.502711 1.611950 1.675343 2.027867 2.070595 2.146376 2.239904
##
## $index_mean_BootSE
## [1] 2.3841682 0.7621835
##
## $index_CI
## 0.1% 1% 5% 10% 90% 95% 99% 99.9%
## 0.4181733 0.8150223 1.2302829 1.4286278 3.3949870 3.6937418 4.2901225 5.0475814
##
## $`+1sd_indirect_mean_BootSE`
## [1] 5.126119 1.344359
##
## $`+1sd_indirect_CI`
## 0.1% 1% 5% 10% 90% 95%
## 0.03462165 1.58517410 3.00820516 3.62764011 6.63770027 7.22775884
## 99% 99.9%
## 8.87115018 11.89170854
##
## $`-1sd_indirect_mean_BootSE`
## [1] -0.9530577 1.1246842
##
## $`-1sd_indirect_CI`
## 0.1% 1% 5% 10% 90% 95% 99% 99.9%
## -4.320002 -3.136619 -2.409837 -2.127015 0.788888 1.225049 1.762490 2.420024
##
## $`mediation r^2`
## [1] 0.8062191
##
## $`total r^2`
## [1] 0.9098085