I've set up a Bayesian regression model in WinBUGS to determine values for the unknown parameters (b1, b2, b3, b4) and intercept value (b0) in a linear regression model. The code is as follows:
model {
for (i in 1:(J-1)) {
FC ~ dnorm(mu, tau)
mu <- b0 + b1*(Factor_b1-mean_Factor_b1) + b2*(Factor_b2-mean_Factor_b2) + b3*(Factor_b3-mean_Factor_b3) + b4*(Factor_b4-mean_Factor_b4)
}
b0 ~ dflat()
b1 ~ dflat()
b2 ~ dflat()
b3 ~ dflat()
b4 ~ dflat()
tau <- 1/sigma2
log(sigma2) <- 2*log.sigma
log.sigma ~ dflat()
}
Inits:
list(b0 =0,b1 = 0, b2 =0, b3 = 0, b4 =0, log.sigma=0)
Data1
list(J = 20, FC = c(1.87315166256848, 1.87315166256848,
1.87315166256848, 1.8708501655802, 1.8708501655802, 1.8708501655802,
1.93248104062608, 1.93248104062608, 1.93248104062608, 1.93248104062608,
1.80846914258265, 1.80846914258265, 1.80846914258265, 2.10555453929548,
2.10555453929548, 2.10555453929548, 2.10555453929548, 2.10555453929548,
2.12908503670568, 2.12908503670568), Factor_b1 = c(7.0057890192535,
7.0057890192535, 7.0057890192535, 7.05012252026906, 7.05012252026906,
7.05012252026906, 7.13329595489607, 7.13329595489607, 7.13329595489607,
7.13329595489607, 7.13329595489607, 7.13329595489607, 7.13329595489607,
7.11720550316434, 7.11720550316434, 7.11720550316434, 7.11720550316434,
7.11720550316434, 7.14124512235049, 7.14124512235049), mean_Factor_b1 = 7.09846620316814,
Factor_b2 = c(7.2211050981825, 7.2211050981825, 7.2211050981825,
7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825,
7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825,
7.2211050981825, 7.2211050981825, 7.37650812632622, 7.37650812632622,
7.37650812632622, 7.37650812632622, 7.37650812632622, 7.46565531013406,
7.46565531013406), mean_Factor_b2 = 7.28441087641358, Factor_b3 = c(2.37954613413017,
2.37954613413017, 2.37954613413017, 2.28238238567653, 2.28238238567653,
2.28238238567653, 2.28238238567653, 2.28238238567653, 2.28238238567653,
2.28238238567653, 2.28238238567653, 2.28238238567653, 2.28238238567653,
2.33214389523559, 2.33214389523559, 2.33214389523559, 2.33214389523559,
2.33214389523559, 2.33214389523559, 2.33214389523559), mean_Factor_b3 = 2.31437347629025,
Factor_b4 = c(2.06686275947298, 2.06686275947298, 2.06686275947298,
2.09186406167839, 2.09186406167839, 2.09186406167839, 2.10413415427021,
2.10413415427021, 2.10413415427021, 2.10413415427021, 2.06686275947298,
2.06686275947298, 2.06686275947298, 2.32238772029023, 2.32238772029023,
2.32238772029023, 2.32238772029023, 2.32238772029023, 2.2082744135228,
2.2082744135228), mean_Factor_b4 = 2.15608963937253)
This WinBUGS code returns the following outputs for the unknown intercept (b0) and parameter (b1, b2, b3, b4) values:
node mean sd MC error 2.5% median 97.5% start sample
b0 1.957 0.009337 3.764E-5 1.939 1.957 1.976 1001 56000
b1 0.1068 0.3296 0.001438 -0.5529 0.1072 0.7615 1001 56000
b2 0.5977 0.2758 0.001068 0.05286 0.5967 1.147 1001 56000
b3 0.1892 0.4394 0.001825 -0.6871 0.1899 1.061 1001 56000
b4 0.5757 0.1886 7.423E-4 0.1986 0.5765 0.9472 1001 56000
MY PROBLEM: When I compare these Bayesian estimates with results from a linear MLE regression in R, I seem to be getting a different result for the intercept value (b0). The code for the R linear regression is as follows:
FC = c(1.87315166256848, 1.87315166256848, 1.87315166256848, 1.8708501655802, 1.8708501655802, 1.8708501655802, 1.93248104062608, 1.93248104062608, 1.93248104062608, 1.93248104062608, 1.80846914258265, 1.80846914258265, 1.80846914258265, 2.10555453929548, 2.10555453929548, 2.10555453929548, 2.10555453929548, 2.10555453929548, 2.12908503670568, 2.12908503670568)
b1 = c(7.0057890192535, 7.0057890192535, 7.0057890192535, 7.05012252026906, 7.05012252026906, 7.13329595489607, 7.13329595489607, 7.13329595489607, 7.13329595489607, 7.13329595489607, 7.13329595489607, 7.13329595489607, 7.11720550316434, 7.11720550316434, 7.11720550316434, 7.11720550316434, 7.11720550316434, 7.14124512235049, 7.14124512235049)
b2 = c(7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.2211050981825, 7.37650812632622, 7.37650812632622, 7.37650812632622, 7.37650812632622, 7.37650812632622, 7.46565531013406, 7.46565531013406)
b3 = c(2.37954613413017, 2.37954613413017, 2.37954613413017, 2.28238238567653, 2.28238238567653, 2.28238238567653, 2.28238238567653, 2.28238238567653, 2.28238238567653, 2.28238238567653, 2.28238238567653, 2.28238238567653, 2.28238238567653, 2.33214389523559, 2.33214389523559, 2.33214389523559, 2.33214389523559, 2.33214389523559, 2.33214389523559, 2.33214389523559)
b4 = c(2.06686275947298, 2.06686275947298, 2.06686275947298, 2.09186406167839, 2.09186406167839, 2.09186406167839, 2.10413415427021, 2.10413415427021, 2.10413415427021, 2.10413415427021, 2.06686275947298, 2.06686275947298, 2.06686275947298, 2.32238772029023, 2.32238772029023, 2.32238772029023, 2.32238772029023, 2.32238772029023, 2.2082744135228, 2.2082744135228)
# ======================= Linear Model =======================
lmfit_Linear_Model_Test =lm(FC ~ (b1 + b2 + b3 + b4))
print (summary(lmfit_Linear_Model_Test))
And the results from this MLE regressions are as follows:
Call:
lm(formula = FC ~ (b1 + b2 + b3 + b4))
Residuals:
Min 1Q Median 3Q Max
-0.05837 -0.00823 -0.00044 0.01307 0.04593
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -5.247 2.305 -2.28 0.0379 *
b1 0.105 0.298 0.35 0.7280
b2 0.674 0.195 3.45 0.0035 **
b3 0.177 0.394 0.45 0.6599
b4 0.529 0.141 3.75 0.0019 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.0359 on 15 degrees of freedom
Multiple R-squared: 0.931, Adjusted R-squared: 0.913
F-statistic: 50.8 on 4 and 15 DF, p-value: 0.0000000153
SUMMARY: Why is the intercept value (b0) coming to -5.247 with the MLE model and 1.957 with the Bayesian model? Should they not be the same?
r regression bayesian winbugs
扫码加好友，拉您进群
全部版块
我的主页
收藏
