## 2 Materials and Methods
### The Hierarchical Bayesian MAP  
The hierarchical Bayesian MAP is a statistical model written in multiple levels that estimates the parameters of the posterior distribution using the MAP method. The estimation of the sparse representation \(\boldsymbol{\alpha}\) can be decomposed into the estimation of two independent parameters \(\mathbf{x}\) and \(\mathbf{V}\) as shown in (7). \[\mathbf{\hat{x}}\ =\ \operatorname*{argmax}_{\mathbf{x}}\mathbf{P}^{\mathrm{T}} \mathbf{M}^{-1}(\mathbf{x})\mathbf{P}+\log\det\mathbf{M}(\mathbf{x})+\mathbf{ x}^{\mathrm{T}}\sum\
olimits_{x}^{-1}\mathbf{x} \tag{9}\]  
where, \(\det(\mathbf{M}(\mathbf{x}))\) is the determinant of \(\mathbf{M}(\mathbf{x})\) and  
\[\mathbf{M}(\mathbf{x})\ =\ \mathbf{A}\mathbf{H}\psi\mathbf{G}(\mathbf{x})\sum \limits_{v}\mathbf{G}(\mathbf{x})\psi^{\mathrm{T}}\mathbf{H}^{\mathrm{T}} \mathbf{A}^{\mathrm{T}}+\sum\
olimits_{m^{\prime}}\quad\mathbf{M}(\mathbf{x}) \in\mathbb{R}^{m\times m} \tag{10}\]  
\[\mathbf{G}(\mathbf{x})\ =\ diag(g(\mathbf{x}))\in\mathbb{R}^{k\times k} \tag{11}\]  
where \(diag(\cdot)\) is a square diagonal. Assuming that \(\mathbf{x}\) is known, we can obtain the following solution by using the Bayesian formula  
\[\mathbf{\hat{V}}\ =\ \operatorname*{argmax}_{\mathbf{V}}\{\log p(\mathbf{P}| \mathbf{V})+\log p(\mathbf{V})\} \tag{12}\]  
Substituting (7) and (9) into (12), the following solution can be obtained  
\[\mathbf{\hat{V}}\ =\ \operatorname*{argmax}_{\mathbf{\hat{V}}}(\mathbf{P}- \mathbf{A}\mathbf{H}\psi\mathbf{G}(\mathbf{\hat{x}})\mathbf{V})^{\mathrm{T}} \sum\
olimits_{m}^{-1}(\mathbf{P}-\mathbf{A}\mathbf{H}\psi\mathbf{G}(\mathbf{ \hat{x}})\mathbf{V})+\mathbf{V}^{\mathrm{T}}\sum\
olimits_{v}^{-1}\mathbf{V} \tag{13}\]  
(12) is quadratic on the unknown variable \(\mathbf{V}\) and its solution is given,  
\[\mathbf{\hat{V}}\ =\ \left(\left(\mathbf{G}^{\mathrm{T}}(\mathbf{\hat{x}})\psi^{ \mathrm{T}}\mathbf{H}^{\mathrm{T}}\mathbf{A}^{\mathrm{T}}\right)\sum\
olimits_{ m}^{-1}(\mathbf{A}\mathbf{H}\psi\mathbf{G}(\mathbf{\hat{x}}))+\sum\
olimits_{v}^{-1} \right)^{-1}(\mathbf{A}\mathbf{H}\psi\mathbf{G}(\mathbf{\hat{x}}))^{\mathrm{T}} \mathbf{P} \tag{14}\]  
According to (7) and (14), the estimation of the sharp CT image is given by \(\mathbf{\hat{f}}=\psi\mathbf{\hat{\alpha}}\ =\ \psi\left(\mathbf{\hat{z}}\odot\mathbf{\hat{V}}\right)\). In (9) and (14), the blurring matrix \(\mathbf{H}\) is usually unknown, so the degradation matrix \(\mathbf{H}\) must be estimated. Generally, \(\mathbf{H}\) is a circular matrix, and it can be obtained by a cyclic shift of any row \(\mathbf{h}_{j},(j\ =\ 1,2,\cdots,m)\). Therefore, we can obtain the blurring matrix \(\mathbf{H}\) by estimating the point spread function. Here, the Fourier transform was implemented before we estimate \(\mathbf{h}\), \(\mathbf{\overline{h}}\) and \(\mathbf{\overline{h}}_{0}\) are the Fourier transform of \(\mathbf{h}\) and \(\mathbf{h}_{0}\), respectively. \(\mathbf{h}_{0}\) is an initial estimation of \(\mathbf{h}\). A multivariate Gaussian distribution is chosen for the prior of the point spread function \(\mathbf{\overline{h}}\)[23]. \[p\left(\mathbf{\overline{h}}\right)\ =\ \frac{1}{\left(2\pi\right)^{\frac{\beta}{2}}| \Sigma_{h}|^{\frac{1}{2}}}\exp\left(-\frac{1}{2}\left(\mathbf{\overline{h}}- \mathbf{\overline{h}}_{0}\right)^{\mathrm{T}}\Sigma_{\overline{h}}^{-1}\left( \mathbf{\overline{h}}-\mathbf{\overline{h}}_{0}\right)\right) \tag{15}\]where \(\Sigma_{\widetilde{h}}\) is the covariance matrix of the variable \(\widetilde{\mathbf{h}}\).