Last updated: 2018-09-04
workflowr checks: (Click a bullet for more information) ✔ R Markdown file: up-to-date
Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.
✔ Environment: empty
Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.
✔ Seed:
set.seed(20180529)
The command set.seed(20180529)
was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.
✔ Session information: recorded
Great job! Recording the operating system, R version, and package versions is critical for reproducibility.
✔ Repository version: ce64c26
wflow_publish
or wflow_git_commit
). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:
Ignored files:
Ignored: .Rhistory
Ignored: .Rproj.user/
Ignored: analysis/.Rhistory
Ignored: docs/.DS_Store
Ignored: docs/figure/
Untracked files:
Untracked: analysis/SERmodel.Rmd
Unstaged changes:
Modified: analysis/_site.yml
Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.
File | Version | Author | Date | Message |
---|---|---|---|---|
Rmd | ce64c26 | zouyuxin | 2018-09-04 | wflow_publish(“analysis/BayesianLinearRegression.Rmd”) |
html | 1c9d4bc | zouyuxin | 2018-09-04 | Build site. |
Rmd | 29b2082 | zouyuxin | 2018-09-04 | wflow_publish(“analysis/BayesianLinearRegression.Rmd”) |
html | fe9dbcb | zouyuxin | 2018-09-04 | Build site. |
Rmd | 7b41f36 | zouyuxin | 2018-09-04 | wflow_publish(“analysis/BayesianLinearRegression.Rmd”) |
html | cc28cde | zouyuxin | 2018-09-04 | Build site. |
Rmd | 44601f6 | zouyuxin | 2018-09-04 | wflow_publish(“analysis/BayesianLinearRegression.Rmd”) |
html | e8fbbc2 | zouyuxin | 2018-08-20 | Build site. |
Rmd | 225c5c2 | zouyuxin | 2018-08-20 | wflow_publish(“analysis/BayesianLinearRegression.Rmd”) |
html | b8471c6 | zouyuxin | 2018-05-29 | Build site. |
Rmd | 2b20f8d | zouyuxin | 2018-05-29 | wflow_publish(“analysis/BayesianLinearRegression.Rmd”) |
html | 44284bd | zouyuxin | 2018-05-29 | Build site. |
Rmd | 751426d | zouyuxin | 2018-05-29 | wflow_publish(c(“analysis/BayesianLinearRegression.Rmd”, |
We have linear model, \(\mathbf{y} = (y_{1}, \cdots, y_{n})^{T}\) \[ \mathbf{y} = \mathbf{x}\beta + \mathbf{\epsilon} \quad \mathbf{\epsilon} \sim N_{n}(0, \frac{1}{\tau}I) \] Here \(\mathbf{y}\) is the n-vector of response data (centered to have mean 0), \(\mathbf{x}\) is an n-vector of covariates (similarly centered), and \(\epsilon\) is an n-vector of independent error terms with precision (inverse variance) \(\tau\).
Suppose we have prior on \(\beta\) \[ \beta\sim N(0,1/\tau_{0}) \] The likelihood is \[ p(\mathbf{y}|\mathbf{x}, \beta, \sigma^2) = \left(\frac{1}{\sqrt{2\pi/\tau}}\right)^{n} \exp\left( -\frac{\tau}{2} \|\mathbf{y} - \mathbf{x}\beta\|^2 \right) \] The posterior for \(\beta\) is \[ \begin{align*} p(\beta|\mathbf{y}, \mathbf{x}, \tau, \tau_{0}) &\propto p(\mathbf{y}|\mathbf{x}, \beta, \tau)p(\beta|\tau_{0}) \\ &\propto \exp\left( -\frac{\tau}{2} \|\mathbf{y} - \mathbf{x}\beta\|^2 - \frac{1}{2}\beta^2 \tau_{0} \right) \\ &\propto \exp\left( -\frac{\tau}{2} \left[ -2\mathbf{y}^{T}\mathbf{x}\beta + \beta^2 \mathbf{x}^{T}\mathbf{x}\right] - \frac{1}{2}\beta^2 \tau_{0} \right) \\ &= \exp\left( -\frac{\tau}{2}\left[ -2\mathbf{y}^{T}\mathbf{x}\beta + (\mathbf{x}^{T}\mathbf{x} + \frac{\tau_{0}}{\tau})\beta^2 \right] \right) \\ &= \exp\left( -\frac{1}{2}\left[(\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0})\beta^2 -2\tau\mathbf{y}^{T}\mathbf{x}\left(\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0}\right)^{-1} \left( \tau\mathbf{x}^{T}\mathbf{x} + \tau_{0} \right) \beta \right] \right) \\ &\propto \exp\left( -\frac{1}{2} \left( \tau\mathbf{x}^{T}\mathbf{x} + \tau_{0} \right) (\beta - \tau\mathbf{x}^{T}\mathbf{y}\left( \tau\mathbf{x}^{T}\mathbf{x} + \tau_{0} \right)^{-1})^2 \right) \end{align*} \]
\[ \color{red}{\beta|\mathbf{y},\mathbf{x}, \sigma^2, \tau^2 \sim N(\tau\mathbf{x}^{T}\mathbf{y}\left( \tau\mathbf{x}^{T}\mathbf{x} + \tau_{0} \right)^{-1}, \left( \tau\mathbf{x}^{T}\mathbf{x} + \tau_{0} \right)^{-1})} \]
Suppose we have prior on \(\beta\), \[ \beta|\gamma = 1 \sim N(0,1/\tau_{0}) \] \(P(\beta = 0|\gamma = 0) = 1\), \(P(\gamma=1) = \pi\)
\(\Rightarrow\) \[ \beta \sim \pi N(0, 1/\tau_{0}) + (1-\pi) \delta_{0} \]
Or, we can formulate the problem as \[ \mathbf{y} = \mathbf{x}\tilde{\beta} + \mathbf{\epsilon} \quad \mathbf{\epsilon} \sim N_{n}(0, \frac{1}{\tau}I) \\ \tilde{\beta} = \gamma \beta \quad \gamma \sim Bern(\pi) \quad \beta \sim N(0,1/\tau_{0}) \] The likelihood is \[ p(\mathbf{y}|\mathbf{x}, \beta, \tau) = \left(\frac{1}{\sqrt{2\pi/\tau}}\right)^{n} \exp\left( -\frac{\tau}{2} \|\mathbf{y} - \mathbf{x}\beta\|^2 \right) \]
The posterior for \(\beta\) is \[ \begin{align*} p(\beta|\mathbf{y}, \mathbf{x}, \tau, \tau_{0}, \gamma = 1) &\propto p(\mathbf{y}|\mathbf{x}, \beta, \tau)p(\beta|\tau_{0}, \gamma=1) \\ &\propto \exp\left( -\frac{\tau}{2} \|\mathbf{y} - \mathbf{x}\beta\|^2 - \frac{1}{2}\beta^2 \tau_{0} \right) \end{align*} \]
\[ \beta|\mathbf{y},\mathbf{x}, \tau,\tau_{0},\gamma=1 \sim N(\tau\mathbf{x}^{T}\mathbf{y}\left( \tau\mathbf{x}^{T}\mathbf{x} + \tau_{0} \right)^{-1}, \left( \tau\mathbf{x}^{T}\mathbf{x} + \tau_{0} \right)^{-1}) \]
The posterior probability of \(\gamma\) \[ \begin{align*} p = P(\gamma=1|\mathbf{y}, \mathbf{x}, \tau, \tau_{0}) &\propto \pi \int \left(\frac{1}{\sqrt{2\pi/\tau}}\right)^{n} \exp\left( -\frac{\tau}{2} \|\mathbf{y} - \mathbf{x}\beta\|^2 \right) N(\beta; 0, 1/\tau_{0}) d\beta \\ & \propto \pi N(0, (\tau\mathbf{x}^{T}\mathbf{x})^{-1} + \tau^{-1}) \\ P(\gamma=0|\mathbf{y}, \mathbf{x}, \tau, \tau_{0}) & \propto (1-\pi) N(0, (\tau\mathbf{x}^{T}\mathbf{x})^{-1}) \end{align*} \]
Note: \(\mathbf{y} \sim N(\mathbf{x}\beta, \frac{1}{\tau}I) \rightarrow \frac{\mathbf{x}^{T}\mathbf{y}}{\mathbf{x}^{T}\mathbf{x}} \sim N(\beta, (\tau\mathbf{x}^{T}\mathbf{x})^{-1})\)
\[ \begin{align*} \frac{p}{1-p} = \frac{P(\gamma=1|\mathbf{y}, \mathbf{x}, \tau, \tau_{0})}{P(\gamma=0|\mathbf{y}, \mathbf{x}, \tau, \tau_{0})} = \frac{\pi}{1-\pi} \frac{N(0, (\tau\mathbf{x}^{T}\mathbf{x})^{-1} + \tau_{0}^{-1})}{N(0, (\tau\mathbf{x}^{T}\mathbf{x})^{-1})} = \frac{\pi}{1-\pi} BF(\mathbf{y}, \mathbf{x}) \end{align*} \]
\[ p = \frac{\pi BF(\mathbf{y}, \mathbf{x})}{1-\pi+\pi BF(\mathbf{y}, \mathbf{x})} \] \[ \color{red}{\beta|\mathbf{y},\mathbf{x}, \tau, \tau_{0} \sim p N(\frac{\tau\mathbf{x}^{T}\mathbf{y}}{\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0}}, \left( \tau\mathbf{x}^{T}\mathbf{x} + \tau_{0} \right)^{-1}) + (1-p) \delta_{0}} = pN(\mu, 1/\tau_{1}) + (1-p)\delta_{0} \]
\[ \begin{align*} BF(\mathbf{y}, \mathbf{x}) &= \frac{N(0, (\tau\mathbf{x}^{T}\mathbf{x})^{-1} + \tau_{0}^{-1})}{N(0, (\tau\mathbf{x}^{T}\mathbf{x})^{-1})} \\ &= \sqrt{\frac{1/\tau\mathbf{x}^{T}\mathbf{x}}{1/\tau\mathbf{x}^{T}\mathbf{x} + 1/\tau_{0}}} \exp\left( \frac{1}{2}(\frac{\mathbf{x}^{T}\mathbf{y}}{\mathbf{x}^{T}\mathbf{x}})^2 \left( \tau \mathbf{x}^{T}\mathbf{x} - \frac{\tau_{0}\tau\mathbf{x}^{T}\mathbf{x}}{\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0}} \right) \right) \\ &= \sqrt{\frac{\tau_{0}}{\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0}}} \exp\left( \frac{1}{2}(\frac{\mathbf{x}^{T}\mathbf{y}}{\mathbf{x}^{T}\mathbf{x}})^2 \left( \tau \mathbf{x}^{T}\mathbf{x} - \frac{\tau_{0}\tau\mathbf{x}^{T}\mathbf{x}}{\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0}} \right) \right) \\ &= \sqrt{\frac{\tau_{0}}{\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0}}} \exp\left( \frac{1}{2}(\frac{\mathbf{x}^{T}\mathbf{y}}{\mathbf{x}^{T}\mathbf{x}})^2 \frac{\tau^2 (\mathbf{x}^{T}\mathbf{x})^2}{\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0}} \right) \\ &= \sqrt{\frac{\tau_{0}}{\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0}}} \exp\left( \frac{1}{2} \frac{\tau^2 (\mathbf{x}^{T}\mathbf{y})^2 }{\tau\mathbf{x}^{T}\mathbf{x} + \tau_{0}} \right) \\ \end{align*} \]
\[ \color{red}{BF(\mathbf{y}, \mathbf{x}) = \sqrt{\frac{\tau_{0}}{\tau_{1}}} \exp\left( \frac{1}{2}\mu^2\tau_{1} \right)} \]
sessionInfo()
R version 3.5.1 (2018-07-02)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS High Sierra 10.13.6
Matrix products: default
BLAS: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRlapack.dylib
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
attached base packages:
[1] stats graphics grDevices utils datasets methods base
loaded via a namespace (and not attached):
[1] workflowr_1.1.1 Rcpp_0.12.18 digest_0.6.15
[4] rprojroot_1.3-2 R.methodsS3_1.7.1 backports_1.1.2
[7] git2r_0.23.0 magrittr_1.5 evaluate_0.11
[10] stringi_1.2.4 whisker_0.3-2 R.oo_1.22.0
[13] R.utils_2.6.0 rmarkdown_1.10 tools_3.5.1
[16] stringr_1.3.1 yaml_2.2.0 compiler_3.5.1
[19] htmltools_0.3.6 knitr_1.20
This reproducible R Markdown analysis was created with workflowr 1.1.1