Authors: Frank Chen, Arvin Nguyen

Perfusion Angiography in Acute Ischemic Stroke

Papers Analysis:

Digital Subtraction Angiography (DSA) is the gold standard for rapid imaging.

The parameters, including cerebral blood flow (CBF) and cerebral blood volume (CBV), mean transit time (MTT), time-to-peak (TTP), and $T_{max}$, are computed using a bolus tracking method based on the deconvolution of the time-density curve on a pixel-by-pixel basis

What is Angiography? It is a technique that relies on X-ray imaging of iodinated radioopaque contrast agent previously injected into the blood stream, "The popularity of DSA can be attributed to its good spatiotemporal resolution which is not easily matched by other acquisition techniques such as magnetic resonance imaging (MRI) and computed tomography (CT)."

2.3

Section 2.3 discusses Video Densitometric Theory. It can be estimated through DSA as the intensity observed in the image is proportional to the contrast concentration. Given by the equation $$I(t)_{(x,y)} = k\mu C(t)_{(x,y)}\rho_{(x,y)}$$ where:

• $I(t)_{(x,y)}$ is the DSA image intensity value for a given pixel $(x,y)$ at time $t$
• $\mu$ is the mass attenuation coefficient of the contrast agent which is proportional to the X-ray energy
• $\rho_{(x,y)}$ is the thickness of the vessel
• $C(x,y)$ is the contrast concentration
• $k$ is a constant that accounts for the X-ray imaging system acquisition and amplification

We can compute the vessel thickness $\rho_{(x,y)}$ by doing the following:

• Apply a vessel detector based on vesselness filtering and thresholding.
• Centerlines are then obtained via skeletization.
• Finally, a perpendicular segment (computed along each point of the centerline) is used to measure the distance to the edges of the vessel and derive the thickness assuming cylindrical volume.
• The thickness is then applied on a cross-sectional basis to every point within the vessel using bicubic interpolation.

2.4

Section 2.4 discusses Perfusion Parameters from DSA using Bolus Tracking. First, to calculate Cerebral Blood Volume (CBV), we use the following equation: $$CBV = \frac{\int_{t=0}^{\infty}C_u(t)dt}{\int_{t=0}^{\infty}C_a(t)dt}$$ Where $u$ is any location in the image, $C_u$ is the amount of contrast agent that has passed through $u$ with respect to the total amount of contrast measured at the feeding arterial vessel $C_a$

The temporal relationship between $C_a$ and and $C_u$ is $$C_u(t) = C_a(t) \otimes h(t)$$ where $\otimes$ is the symbol for the convolution and $h$ is the distribution of the transit times, which is related to the fraction of injected contrast agent still present in the vasculature at any given time $t$. We define the residue function $R(t)$ as $$R(t) = 1 - \int_{\tau = 0}^{t} h(\tau)d\tau$$ We can obtain a more rigorous relationship between $C_u$ and $C_a$ as: $$C_u(t) = CBF(C_a \otimes R) (t)$$ where CBF is Cerebral Blood Flow. $C_a$ and $C_u$ can be estimated from DSA, but $R$ and $CBF$ will need more complex computations.

We can sample $C_u$ and $C_a$ at discrete time points $t_j \in [0, N-1]$: $$C_u(t_j) = \Delta CBF \sum^{N-1}_{i=0}C_a(t_i)R(t_j - t_i)$$ We can actually rewrite this as a matrix-vector notation: $$C_u = \Delta t CBFC_AR$$ where

• $C_u, R \in \mathbb{R}^N$
• $C_A$ is expanded to a Toeplitz matrix.

Now, we want to recover $R$. We can perform SVD on $C_A$: $C_A = UWV^T$. From here, we can find the solution $$R = V\widehat{W}^{-1}U^TC_u$$ We can extract the parameters ($CBF$, $CBV$, $MTT$, $TTP$, and $T_{max}$), see figure below.

2.5

A challenging aspect of the estimation of perfusion parameters is when there is overlap of the vessels (see figure below). We can solve this problem by representing the concentration over time with a Gamma-mixture distribution that is created by using an EM algorithm.

2.5.1 Gamma-Variate Fitting

The density function $\gamma_{\alpha, \beta}$ can be written as: $$ \gamma_{\alpha, \beta} (x) = \begin{cases} \frac{\beta^{\alpha}}{\Gamma(\alpha)}\text{exp}^{-(x-\mu)\beta}(x-\mu)^{\alpha - 1}, & \text{if}\ x - \mu \ge \Delta_{min} \\ 0, & \text{otherwise} \end{cases} $$ where $\alpha$, $\beta$, and $\mu$ are the shape, scale, and location parameters. We define the Gamma function $\Gamma(\alpha)$ as $$\Gamma(\alpha) = \int_0^{\infty}t^{\alpha-1} \text{exp}^{-t}dt$$ We like the Gamma distribution because of it can accomodate different concentration-time curves.

2.5.2 Mixture of Gamma-Variate Distribution

Our mixture equation is defined as the following: $$\mathscr{M}(x, \Theta) = \sum^K_{j=1}\tau_j \gamma_{\alpha_j, \beta_j} (x)$$ where $\gamma_{\alpha_j, \beta_j}(x)$ is the Gamma-variate distribution of the $j^{th}$ component.

2.5.3 Parameter Estimation

We use good-old MLE to estimate our parameter $\Theta$. Our log-likelihood is: $$\mathscr{L}(\Theta) = \sum^N_{i=0}\text{log}\mathscr{M}(x_i, \Theta)$$ and we use Expectation-Maximization (EM) algorithm to estimate our param.

Expectation Step: We calculate expected value $Q(\Theta, \Theta^m)$ of the log-likelihood given current parameteres $\Theta^m$ and $$Q(\Theta, \Theta^m) = \sum^N_{i=1}\sum^K_{j=1} z_{ij}^m \text{log} (\tau_j) + C$$ where $$z_{ij}^m = \frac{\tau_{j}^m \gamma_j(x_i\alpha_j^m, \beta_j^m)}{\mathscr{M}(x_i, \Theta^m)}$$

Maximization Step: We maximize $Q(\Theta, \Theta^m)$ w.r.t. $\Theta$ using numerical optimization (Newton's Method, Gradient Descent). $$\Theta^{m+1} = \mathop{\arg\,\max}\limits_{\Theta} Q(\Theta, \Theta^m)$$ We execute this iterative procedure until the convergence criterion $| \Theta^{m+1} - \Theta^m | < \tau_{EM}$ is satisfied.