Cellular Automata (CA) have been used for modeling physical systems and for solving differential equations numerically [1]. In this abstract, we show that the MRI signal progression for repeating sequences can be modeled as a sum of CA, opening up the possibility that the theorems of CA can help solve MRI design problems. The sum naturally splits across excited magnetization and recovered magnetization providing insight into MRI signal progression through simulation.

Extended Phase Graphs (EPG) are extremely useful for understanding MRI signal progression [2]. They represent a spin system using a Fourier basis assuming gradients induce an integer number of cycles of phase across a voxel. A pictorial depiction of the EPG basis functions is shown in figure 1.

Given this assumption, the magnetization using the basis of figure 1 can be represented by the EPG coefficients defined by the analysis equations presented below.

$$\begin{align}F_n^+ &=\int_0^1 M_{xy}(z)\exp\left(-i\,2\pi\,n\,z\right) dz\\ F_n^- &= \int_0^1 M_{xy}^*(z)\exp\left(-i\,2\pi\,n\,z\right) dz \\ Z_n &= \int_0^1 M_z(z)\exp\left(-i\,2\pi\,n\,z\right) dz \end{align}$$

Note that $F_n^-=\left(F_{-n}^+\right)^*$, thus the synthesis equations can use only the coefficients for $n\geq 0$. With the EPG coefficients, one can perform synthesis with the equations below.

$$\begin{align}M_{xy}(z) &= F_0^+ + \sum_{n=1}^\infty F_n^+\,\exp\left(i\,2\pi\,n\,z\right) + \sum_{n=1}^\infty (F_n^-)^*\,\exp\left(-i\,2\pi\,n\,z\right) \\ M_z(z) &= \text{Real}\left(Z_0+2\sum_{n=1}^\infty Z_n\exp(i\,2\pi\,n\,z) \right)\end{align}$$

In the EPG domain, the effects of RF pulses, relaxations, gradients, and recovery are simple functions of the EPG coefficients as summarized in figure 2.

The rule of a discrete CA can be though of as a generalization of discrete convolution. A kernel of finite size is slide over the data. Whereas the output value for each location of the kernel with convolution is the result of a pointwise multiplication and a summation, the output of the CA’s rule is any function of the input values.

As shown in the equations below, longitudinal recovery only affects the $Z_0$ term (refer to figure 2); this special case makes it challenging to represent the magnetization over time as a CA since the kernel of a CA is unchanged for all state coefficients.

$$\begin{align} Z_n' &= \exp(-t/T1) \, Z_n \\ Z_0' &= \left(1-\exp(-t/T1)\right)M_0 + \exp(-t/T1)\, Z_0 \end{align}$$

Initially, we will let $Z_n’ = \exp(-t/T1)\,Z_n$ for all $n$; the recovery term for $Z_0$ will be taken into account later. The magnetization immediately after the $\alpha_0$ pulse is the initial condition of our CA: $Q_0=R_{\alpha_0}\left( Q_0^- \right)$, where $R_{\alpha_0}$ is the rotation matrix corresponding to flip angle $\alpha_0$. At time $t_1+t_2$, the magnetization is represented by $Q_1$:

$$Q_1(n) = \begin{bmatrix} e^{-t_2/T2} \, R_{\alpha,x}\left( e^{-t_1/T2} F_{n-2}^+, e^{-t_1/T2} F_{n}^-, e^{-t_1/T1} Z_{n-1} \right) \\ e^{-t_2/T2} \, R_{\alpha,y}\left( e^{-t_1/T2} F_n^+, e^{-t_1/T2} F_{n+2}^-, e^{-t_1/T1} Z_{n+1} \right) \\ e^{-t_2/T1} \, R_{\alpha,z}\left( e^{-t_1/T2} F_{n-1}^+, e^{-t_1/T2} F_{n+1}^-, e^{-t_1/T1} Z_n \right)\end{bmatrix}$$

This is a compact representation of the EPG propagation. Note that $Q_1(n)$ is a function of $\left\{ Q_0(n-2), \ldots, Q_0(n+2) \right\}$ for all $n$. Additionally, note that the magnetization $Q_m$ (without taking recovery into account) can be determined by applying this same function recursively. This is a Cellular Automata; we will denote this process as $C_{Q_0}(m,n)$ where $m$ is the time index and $n$ is the basis function index.

It is assumed that enough coefficients are stored such that the magnetization for $n>n_{\text{max}}$ is very small, and thus we impose a 0 boundary condition for large $n$. For small $n$, we exploit the symmetry in the $F$ coefficients and use the form of symmetric boundary conditions shown below.

$$\begin{array}{cc} F_{-1}^+ = (F_1^-)^* & F_{-1}^- = (F_1^+)^* \\ F_{-2}^+ + (F_2^-)^* & F_{-2}^- = (F_2^+)^* \end{array}$$

The magnetization $Q_m(n)$ without recovery is linear in the EPG coefficients. Thus, the total magnetization $Q_1$ is the sum of the CA process $C_{Q_0}(1,\cdot)$ and the magnetization due to recovery. At time $t_1+t_2$ the magnetization due to recovery has all EPG coefficients equal to 0 except the following.

$$\begin{align}{F_1^+}’ &= e^{-t_2/T2}\left(1-e^{-t_1/T1}\right)\,M_0\,\left(-ie^{i\phi}\right)\sin\alpha \\ Z_0' &= \left(1-e^{-t_1/T1}\right)\,M_0 + e^{-t_2/T1}\,\cos\alpha\,\left(1-e^{-t_1/T1}\right)\,M_0\end{align}$$

Let $\tilde{Q}$ have all zeros except for $F_1^+$ and $Z_0$, which are defined as written above. Then the magnetization at time $t_1+t_2$ is $Q_1(n)=\tilde{Q}+C_{Q_0}(1,n)$. The magnetization $\tilde{Q}$ is an initial condition for the same CA process that starts at time $t_1+t_2$ as shown in figure 4.

Thus, the magnetization $Q_m$ is the sum of CA processes as written below.

$$Q_m(n) = C_{Q_0}(m,n) + \sum_{i=1}^m C_{\tilde{Q}}(m-i,n)$$

Figure 5 shows the results of a CA simulation. In this simulation, we see the contribution from the initial excitation and from each recovery period separately.We have shown that the signal progression in MRI with repeated sequences can be represented by a sum of CA.