Last updated: 2018-10-19
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.
✔ Repository version: 1476371
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/Test.Rmd/
| File | Version | Author | Date | Message |
|---|---|---|---|---|
| Rmd | 1476371 | zouyuxin | 2018-10-19 | wflow_publish(“analysis/MouthwashMash.Rmd”) |
| html | 34d5339 | zouyuxin | 2018-10-19 | Build site. |
| Rmd | 2b25b41 | zouyuxin | 2018-10-19 | wflow_publish(“analysis/MouthwashMash.Rmd”) |
| html | 7fbf98a | zouyuxin | 2018-10-19 | Build site. |
| Rmd | e8f7e50 | zouyuxin | 2018-10-19 | wflow_publish(“analysis/MouthwashMash.Rmd”) |
| html | f58230a | zouyuxin | 2018-10-19 | Build site. |
| Rmd | ff8c336 | zouyuxin | 2018-10-19 | wflow_publish(“analysis/MouthwashMash.Rmd”) |
| html | 9156912 | zouyuxin | 2018-10-19 | Build site. |
| Rmd | 2c4622c | zouyuxin | 2018-10-19 | wflow_publish(“analysis/MouthwashMash.Rmd”) |
Suppose we have individual level data for each sample n, gene j, at tissue r \[ \mathbf{Y}_{N\times J}^{(r)} = \mathbf{x}_{N\times 1}^{(r)}\boldsymbol{\beta}_{1\times J}^{(r)} + \mathbf{Z}_{N\times q}^{(r)}\boldsymbol{\alpha}_{q\times J}^{(r)} + \mathbf{E}_{N\times J}^{(r)} \]
where \(y_{ij}^{(r)}\) is the normalized expression level of gene j in sample i, tissue r; \(\mathbf{x}\) is a matrix containing observed covariates, with \(\boldsymbol{\beta}\) a matrix of corresponding effects; \(\mathbf{Z}\) is a matrix of unobserved factors causing unwanted variation, with \(\boldsymbol{\alpha}\) a matrix of corresponding effects; and \(\mathbf{E}\) has independent (Gaussian) errors with means 0 and column-specific variances var(\(e_{ij}\)) = \(\sigma_{j}^{2}\). Only \(\mathbf{Y}\) and \(\mathbf{x}\) are known; other quantities are to be estimated.
For each gene j (\(j = 1, \cdots, J\)) at each tissue r (\(r = 1, \cdots, R\)), obtain an initial estimate \(\hat{\beta}_{j}^{(r)}\) for \(\beta_{j}^{(r)}\) ignoring unwanted variation by using ordinary least squares (OLS) regression of the jth column of \(Y^{(r)}\) on \(x^{(r)}\). \[ \hat{\beta}_{j}^{(r)} = \frac{\mathbf{x}^{(r)T}\mathbf{y}_{\cdot j}^{(r)} }{\mathbf{x}^{(r)T}\mathbf{x}^{(r)}} \] We have \(\hat{\boldsymbol{\beta}}^{(r)} \in \mathbb{R}^{J}\), \(\hat{B} \in \mathbb{R}^{J\times R}\).
For each tissue r, perform FA on the residuals \(\tilde{\mathbf{Y}}^{(r)} = \mathbf{Y}^{(r)} - \mathbf{x}^{(r)}\hat{\boldsymbol{\beta}}^{(r)T}\). The FA model is \[ \tilde{\mathbf{Y}}^{(r)} = \tilde{\mathbf{Z}}^{(r)}\tilde{\boldsymbol{\alpha}}^{(r)} + \tilde{\mathbf{E}}^{(r)} \] We obtain \(\hat{\boldsymbol{\alpha}}^{(r)} \in \mathbb{R}^{q\times p}\) and \(\hat{\sigma}_{j(r)}^{2}\).
A simplified model: \[ \hat{\beta}_{j}^{(r)} \sim N(\beta_{j}^{(r)} + \hat{\boldsymbol{\alpha}^{(r)}}_{j}^{T}\mathbf{z}^{(r)}, \hat{s}_{j(r)}^{2}) \] \(\mathbf{z}^{(r)} \in \mathbb{R}^{q}\), \(\mathbf{S}^{(r)} = diag(\hat{s}_{1(r)}^{2}, \cdots, \hat{s}_{J(r)}^{2})\) \[ \hat{s}_{j(r)}^{2} = \frac{\hat{\sigma}_{j(r)}^{2}}{\mathbf{x}^{T}\mathbf{x}} \] We can formulate it in the mash framework, \(\mathbf{V}\) is the correlation matrix which is unknown \[ \hat{\boldsymbol{\beta}}_{j} \sim N_{R}(\boldsymbol{\beta}_{j} + \hat{\boldsymbol{\alpha}}_{j}^{T} \mathbf{z}, \mathbf{S}_{j}\boldsymbol{\xi}^{1/2}\mathbf{V}\boldsymbol{\xi}^{1/2}\mathbf{S}_{j}) = N_{R}(\boldsymbol{\beta}_{j} + \hat{\boldsymbol{\alpha}}_{j}^{T} \mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}) \] where \[ \hat{\boldsymbol{\alpha}}_{j} = \left(\begin{matrix} \hat{\mathbf{\alpha}}_{j}^{(1)} & & \\ & \ddots & \\ & & \hat{\mathbf{\alpha}}_{j}^{(R)}\end{matrix}\right)_{Rq\times R} \quad \mathbf{z} = \left(\begin{matrix} \mathbf{z}^{(1)} \\ \vdots \\ \mathbf{z}^{(R)} \end{matrix}\right)_{Rq\times 1} \quad \boldsymbol{\xi} = diag(\xi^{(1)}, \cdots, \xi^{(R)}) \quad \mathbf{S}_{j} = diag(\hat{s}_{j}^{(1)}, \cdots, \hat{s}_{j}^{(R)}) \] \(\hat{s}_{j}^{(r)} = \frac{\hat{\sigma}_{j}^{(r)}}{\mathbf{x}^{(r)T}\mathbf{x}^{(r)}}\), \(\mathbf{C} = \boldsymbol{\xi}^{1/2}\mathbf{V}\boldsymbol{\xi}^{1/2}\). \ Estimate the unimodel effects distribution g and the unwanted variation effects \(\mathbf{z}\), multiplicative parameter \(\boldsymbol{\xi}\) by maximum (marginal) likelihood \[ \begin{align*} (\hat{g}, \hat{\mathbf{z}}, \hat{\boldsymbol{\xi}}) &= \arg \max_{(g,\mathbf{z}, \boldsymbol{\xi})\in \mathcal{U}\times \mathbb{R}^{Rq} \times \mathbb{R}^{R+}} p(\hat{\mathbf{B}}|g,\mathbf{z},\hat{\boldsymbol{\alpha}}, \hat{\mathbf{S}}, \boldsymbol{\xi})\\ (\hat{g}, \hat{\mathbf{z}}, \hat{\mathbf{C}}) &= \arg \max_{(g,\mathbf{z}, \mathbf{C})\in \mathcal{U}\times \mathbb{R}^{Rq} \times \mathbb{R}_{+}^{R\times R}} \prod_{j=1}^{J}\int_{\boldsymbol{\beta}_{j}} N_{R}(\hat{\boldsymbol{\beta}}_{j}|\boldsymbol{\beta}_{j}+\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}) g(d\boldsymbol{\beta}_{j}) \end{align*} \] Suppose \[ g(\boldsymbol{\beta}_{j}) \sim \sum_{l,k} \pi_{l,k} N_{R}(\boldsymbol{\beta}_{j}; \mathbf{0}, \omega_{l}U_{k}) \] \[ \begin{align*} (\hat{g}, \hat{\mathbf{z}}, \hat{\mathbf{C}}) &= \arg \max_{(g,\mathbf{z}, \mathbf{C})\in \mathcal{U}\times \mathbb{R}^{Rq} \times \mathbb{R}_{+}^{R\times R}} \prod_{j=1}^{J} \sum_{l,k} \pi_{l,k} N_{R}(\hat{\boldsymbol{\beta}}_{j}|\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j} + \omega_{l}U_{k}) \\ &= \arg \max_{(g,\mathbf{z}, \mathbf{C})\in \mathcal{U}\times \mathbb{R}^{Rq} \times \mathbb{R}_{+}^{R\times R}} \prod_{j=1}^{J} \sum_{p=0}^{P} \pi_{p} N_{R}(\hat{\boldsymbol{\beta}}_{j}|\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j} + \Sigma_{p}) \end{align*} \]
Compute posterior distributions \(p(\boldsymbol{\beta}_{j}|\hat{g}, \hat{\mathbf{z}}, \hat{\mathbf{B}}, \hat{\mathbf{S}}, \hat{\mathbf{C}})\), and return posterior summaries.
\[ p(\hat{\mathbf{B}}, \mathbf{W}|\mathbf{z}, \boldsymbol{\pi}, \mathbf{C}) h(\boldsymbol{\pi}|\boldsymbol{\lambda}) = \left(\prod_{p=0}^{P} \pi_{p}^{\lambda_{p}-1} \right) \prod_{j=1}^{J}\prod_{p=0}^{P} \left[\pi_{p} N_{R}(\hat{\boldsymbol{\beta}}_{j}|\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j} + \Sigma_{p})\right]^{\mathbb{I}(w_{jp}=1)} \] The complete log likelihood is \[ l(\mathbf{z}, \boldsymbol{\pi}, \mathbf{C}; \hat{\mathbf{B}}, \mathbf{W}) = \sum_{j=1}^{J} \sum_{p=0}^{P} \mathbb{I}(w_{jp}=1) \left[\log\pi_{p} - \frac{R}{2}\log 2\pi - \frac{1}{2}\log|\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j} + \Sigma_{p}| - \frac{1}{2}(\hat{\boldsymbol{\beta}}_{j}-\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T}(\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j} + \Sigma_{p})^{-1}(\hat{\boldsymbol{\beta}}_{j}-\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}) \right] + \sum_{p=0}^{P} (\lambda_{p}-1)\log \pi_{p} \]
Let \(\mathbf{z}^{(old)}\), \(\boldsymbol{\pi}^{(old)}\), \(\mathbf{C}^{(old)}\) be the current value of the parameters, \[ p(w_{jp} = 1|\hat{\boldsymbol{\beta}}_{j}, \mathbf{z}^{(old)}, \mathbf{C}^{(old)}, \boldsymbol{\xi}^{(old)}) = \frac{\pi_{p}^{(old)} N_{R}(\hat{\boldsymbol{\beta}}_{j}|\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}^{(old)}, \mathbf{S}_{j}\mathbf{C}^{(old)}\mathbf{S}_{j} + \Sigma_{p})}{\sum_{p'=0}^{P}\pi_{p'}^{(old)} N_{R}(\hat{\boldsymbol{\beta}}_{j}|\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}^{(old)}, \mathbf{S}_{j}\mathbf{C}^{(old)}\mathbf{S}_{j} + \Sigma_{p'}) } = \gamma_{jp} \]
Replace \(\mathbb{I}(w_{jp}=1)\) with \(\gamma_{jp}\)
\[ \mathbb{E} = \sum_{j=1}^{J} \sum_{p=0}^{P} \gamma_{jp} \left[\log\pi_{p} - \frac{R}{2}\log 2\pi - \frac{1}{2}\log|\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j} + \Sigma_{p}| - \frac{1}{2}(\hat{\boldsymbol{\beta}}_{j}-\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T}(\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j} + \Sigma_{p})^{-1}(\hat{\boldsymbol{\beta}}_{j}-\hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}) \right] + \sum_{p=0}^{P} (\lambda_{p}-1)\log \pi_{p} \]
\(\boldsymbol{\pi}\): \[ \pi_{m} \leftarrow \frac{(\sum_{j=1}^{J} \gamma_{jp}) + \lambda_{p} - 1}{J - (P+1) + \sum_{p=0}^{P}\lambda_{p}} \]
Then we perform a few iterative updates on \(\mathbf{C}\) and \(\mathbf{z}\).
To update \(\mathbf{z}\) given \(\mathbf{C}\), \[ \hat{\mathbf{z}} = \left[\sum_{j=1}^{J} \hat{\boldsymbol{\alpha}_{j}}\left(\sum_{p=0}^{P} \gamma_{jp}(\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j} + \Sigma_{p})^{-1}\right) \hat{\boldsymbol{\alpha}_{j}}^{T} \right]^{-1} \left[\sum_{j=1}^{J} \hat{\boldsymbol{\alpha}_{j}}\left(\sum_{p=0}^{P} \gamma_{jp}(\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j} + \Sigma_{p})^{-1}\right) \hat{\boldsymbol{\beta}_{j}}\right] \]
To update \(\mathbf{C}\) given \(\mathbf{z}\), \[ \hat{\mathbf{C}} = \frac{1}{J} \sum_{j=1}^{J} \mathbf{S}_{j}^{-1} \mathbb{E}\left[(\hat{\boldsymbol{\beta}}_{j} -\boldsymbol{\beta}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z} )(\hat{\boldsymbol{\beta}}_{j} - \boldsymbol{\beta}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T}|\hat{\boldsymbol{\beta}}_{j}\right]\mathbf{S}_{j}^{-1} \] Details: \[ p(\hat{\mathbf{B}}, \mathbf{B}|\mathbf{z}, \boldsymbol{\pi}, \mathbf{C}) h(\boldsymbol{\pi}|\boldsymbol{\lambda}) = \left(\prod_{p=0}^{P} \pi_{p}^{\lambda_{p}-1} \right) \prod_{j=1}^{J} \left[ N_{R}\left(\hat{\boldsymbol{\beta}}_{j}|\boldsymbol{\beta}_{j} + \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}\right)\sum_{p=0}^{P} \pi_{p} N_{R}(\boldsymbol{\beta}_{j}|\mathbf{0}, \Sigma_{p})\right] \] \[ \mathbb{E}_{\mathbf{B}|\hat{\mathbf{B}}} \log p(\hat{\mathbf{B}}, \mathbf{B}|\mathbf{z}, \boldsymbol{\pi}, \boldsymbol{\xi}) h(\boldsymbol{\pi}|\boldsymbol{\lambda}) = \sum_{p=0}^{P} (\lambda_{p}-1)\log \pi_{p} + \sum_{j=1}^{J} \mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}} \log N_{R}\left(\hat{\boldsymbol{\beta}}_{j}|\boldsymbol{\beta}_{j} + \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}\right) + \sum_{j=1}^{J} \mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}} \sum_{p=0}^{P} \pi_{p} N_{R}(\boldsymbol{\beta}_{j}|\mathbf{0}, \Sigma_{p}) \]
We find \(\mathbf{C}\) that maximize \(\sum_{j=1}^{J} \mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}} \log N_{R}\left(\hat{\boldsymbol{\beta}}_{j}|\boldsymbol{\beta}_{j} + \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}\right)\)
\[ \begin{align*} \log N_{R}\left(\hat{\boldsymbol{\beta}}_{j}|\boldsymbol{\beta}_{j} + \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}\right) &= -\frac{R}{2}\log 2\pi -\frac{1}{2}\log|\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}| - \frac{1}{2} (\hat{\boldsymbol{\beta}}_{j}-\boldsymbol{\beta}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1} (\hat{\boldsymbol{\beta}}_{j}-\boldsymbol{\beta}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}) \\ &= -\frac{R}{2}\log 2\pi -\frac{1}{2}\log|\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}| - \frac{1}{2} (\hat{\boldsymbol{\beta}}_{j}- \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1} (\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}) \\ & - \frac{1}{2} \boldsymbol{\beta}_{j}^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1} \boldsymbol{\beta}_{j} + \frac{1}{2}\boldsymbol{\beta}_{j}^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1}(\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}) +\frac{1}{2}(\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1}\boldsymbol{\beta}_{j} \\ \mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}} \log N_{R}\left(\hat{\boldsymbol{\beta}}_{j}|\boldsymbol{\beta}_{j} + \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}\right) &= -\frac{R}{2}\log 2\pi -\frac{1}{2}\log|\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}| - \frac{1}{2} (\hat{\boldsymbol{\beta}}_{j}- \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1} (\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}) \\ & + \frac{1}{2}\boldsymbol{\mu}_{j}^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1}(\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}) +\frac{1}{2}(\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1}\boldsymbol{\mu}_{j} - \frac{1}{2} tr\left((\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1} \mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}}\left[\boldsymbol{\beta}_{j}\boldsymbol{\beta}_{j}^{T}\right] \right) \end{align*} \]
\[ \begin{align*} f(\mathbf{C}) &= \sum_{j=1}^{J} \mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}} \log N_{R}\left(\hat{\boldsymbol{\beta}}_{j}|\boldsymbol{\beta}_{j} + \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}, \mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j}\right) \\ &= \sum_{j=1}^{J} -\frac{R}{2}\log 2\pi -\log|\mathbf{S}_{j}| - \frac{1}{2} \log|C| - \frac{1}{2} (\hat{\boldsymbol{\beta}}_{j}- \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1} (\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}) \\ & + \frac{1}{2}\boldsymbol{\mu}_{j}^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1}(\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z}) +\frac{1}{2}(\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} (\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1}\boldsymbol{\mu}_{j} - \frac{1}{2} tr\left((\mathbf{S}_{j}\mathbf{C}\mathbf{S}_{j})^{-1} \mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}}\left[\boldsymbol{\beta}_{j}\boldsymbol{\beta}_{j}^{T}\right] \right) \end{align*} \]
\[ \begin{align*} f(\mathbf{C})' &= \sum_{j=1}^{J} -\frac{1}{2}\mathbf{C}^{-1} + \frac{1}{2}\mathbf{C}^{-1}\mathbf{S}_{j}^{-1}(\hat{\boldsymbol{\beta}}_{j}- \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})(\hat{\boldsymbol{\beta}}_{j}- \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T}\mathbf{S}_{j}^{-1}\mathbf{C}^{-1} \\ & - \frac{1}{2} \mathbf{C}^{-1}\mathbf{S}_{j}^{-1} \boldsymbol{\mu}_{j} (\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} \mathbf{S}_{j}^{-1}\mathbf{C}^{-1} - \frac{1}{2} \mathbf{C}^{-1}\mathbf{S}_{j}^{-1} (\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})\boldsymbol{\mu}_{j}^{T} \mathbf{S}_{j}^{-1}\mathbf{C}^{-1} + \frac{1}{2} \mathbf{C}^{-1} \mathbf{S}_{j}^{-1}\mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}}\left[\boldsymbol{\beta}_{j}\boldsymbol{\beta}_{j}^{T}\right]\mathbf{S}_{j}^{-1} \mathbf{C}^{-1} = 0 \\ 0 &= \sum_{j=1}^{J} -\mathbf{C} + \mathbf{S}_{j}^{-1}(\hat{\boldsymbol{\beta}}_{j}- \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})(\hat{\boldsymbol{\beta}}_{j}- \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T}\mathbf{S}_{j}^{-1} \\ & - \mathbf{S}_{j}^{-1} \boldsymbol{\mu}_{j} (\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} \mathbf{S}_{j}^{-1} - \mathbf{S}_{j}^{-1} (\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})\boldsymbol{\mu}_{j}^{T} \mathbf{S}_{j}^{-1} + \mathbf{S}_{j}^{-1}\mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}}\left[\boldsymbol{\beta}_{j}\boldsymbol{\beta}_{j}^{T}\right]\mathbf{S}_{j}^{-1} \\ \hat{\mathbf{C}} &= \frac{1}{J} \sum_{j=1}^{J} \mathbf{S}_{j}^{-1} \left[ (\hat{\boldsymbol{\beta}}_{j}- \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})(\hat{\boldsymbol{\beta}}_{j}- \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} - \boldsymbol{\mu}_{j} (\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T} - (\hat{\boldsymbol{\beta}}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})\boldsymbol{\mu}_{j}^{T} + \mathbb{E}_{\boldsymbol{\beta}_{j}|\hat{\boldsymbol{\beta}}_{j}}\left[\boldsymbol{\beta}_{j}\boldsymbol{\beta}_{j}^{T}\right] \right] \mathbf{S}_{j}^{-1} \\ &= \frac{1}{J} \sum_{j=1}^{J} \mathbf{S}_{j}^{-1} \mathbb{E}\left[(\hat{\boldsymbol{\beta}}_{j} -\boldsymbol{\beta}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z} )(\hat{\boldsymbol{\beta}}_{j} - \boldsymbol{\beta}_{j} - \hat{\boldsymbol{\alpha}_{j}}^{T}\mathbf{z})^{T}|\hat{\boldsymbol{\beta}}_{j}\right]\mathbf{S}_{j}^{-1} \end{align*} \]
This reproducible R Markdown analysis was created with workflowr 1.1.1