Image-based view-angle independent cardiorespiratory motion gating and coronary sinus catheter tracking for x-ray-guided cardiac electrophysiology procedures

The CS catheter used here is composed of 10 electrodes, distributed in pairs along the catheter. The CS tracking technique (Ma et al 2010) uses a fast multi-scale blob detection method (Lindeberg 1993) to detect all electrode-like objects in an x-ray image and a cost function to discriminate the CS catheter from other catheters. After application of the tracking algorithm for all frames in the training view, ${{\theta}_{1}}$, the x and y positions of each of the 10 electrodes are concatenated into a single column vector for each time point. Hence, the data generated by the tracking process consists of:

$$\begin{eqnarray}&&{{\mathbf{s}}_{i}}={{\left({{z}_{i,1,x}},{{z}_{i,1,y}},\ldots ,{{z}_{i,10,x}},{{z}_{i,10,y}}\right)}^{T}},1\leqslant i\leqslant N\end{eqnarray} \tag{ 1 }$$

where z_i, e, x and z_i, e, y represent the x and y coordinates of the eth electrode in the ith frame. e  =  1 corresponds to the distal electrode, and N is the number of frames.

PCA transforms a multivariate dataset of possibly correlated variables into a new dataset of a smaller number of uncorrelated variables called principal components (PCs), without any loss of information (Jolliffe 2002). We first compute the mean vector, $\mathbf{\bar{s}}$, and the covariance matrix, $\mathbf{S}$. The eigenvectors ${{\mathbf{v}}_{m}}$, $1\leqslant m\leqslant M$ of $\mathbf{S}$ represent the PCs and the corresponding eigenvalues d_m, $1\leqslant m\leqslant M$ represent the variance of the data along the direction of the eigenvectors. For our application M  =  20 since the ${{\mathbf{s}}_{i}}$ are of length 20. The result of the PCA gives possible CS catheter electrode positions in ${{\theta}_{1}}$. It was found, by correlation with our gold standard cardiorespiratory motion, that the first and second modes of variation represented the cardiac and respiratory motions, respectively (Panayiotou et al 2013).

Epipolar reconstruction allows loss of depth information to be recovered from a 2D projection image if another 2D image is taken of the same target from a different view, assuming we know the complete camera parameters of the x-ray fluoroscopy projective modality (Ganapathy 1984, Strat 1987, Rougée et al 1993, Hartley and Zisserman 2004). The camera parameters include $C\left({{c}_{s}},{{l}_{s}},f,{{\Delta }_{u}},{{\Delta }_{v}},{{\delta}_{x}},{{\delta}_{y}},{{\delta}_{z}},\theta ,\phi ,\psi \right)$. $\left({{c}_{s}},{{l}_{s}}\right)$ are the coordinates in the image where the beam is normal to the projection plane, f is the distance between the source and the projection plane, ${{\Delta }_{u}}$ and ${{\Delta }_{v}}$ are the pixel spacing in the (u, v)-axes, respectively on the projection plane, $\delta =\left({{\delta}_{x}},{{\delta}_{y}},{{\delta}_{z}}\right)$ and ($\theta ,\phi ,\psi$) are the parameters defining the position and orientation of the camera in 3D space. The camera parameters are obtained from the DICOM header of the x-ray images. Using epipolar geometry, each of the 2D electrode positions obtained when forming our statistical model in ${{\theta}_{1}}$ is backward projected to form a 3D line, ${{\mathbf{l}}_{b}}$, which is then forward projected to generate a 2D epipolar line, ${{\mathbf{l}}_{e}}$, in the current view, ${{\theta}_{2}}$, that contains the corresponding electrode position. Hence, we generate 10 epipolar lines in ${{\theta}_{2}}$ for each frame in the current view.

We now determine the 2D position of the CS catheter in the new image in ${{\theta}_{2}}$ and simultaneously determine its state of cardiorespiratory motion. Although the real-time catheter tracking algorithm (Ma et al 2010) is used in the training sequence, it fails in rotational datasets at view-angles that distort the solid circular shape of the CS catheter electrodes. Additionally, the multi-scale blob detection algorithm misdetects blobs at very low doses. Consequently, we cannot solely use this method for 2D catheter tracking in rotational sequences at low doses. Instead, we determine the cardiorespiratory state in the PCA model of ${{\theta}_{1}}$ that produces epipolar lines that most closely match the catheter electrodes in the new image in ${{\theta}_{2}}$. First, the fast multi-scale blob detector (Lindeberg 1993) is used to detect all possible electrode-like objects in the new image. The positions of all detected blobs are concatenated into a blob list, $\mathbf{BL}$.

Using the PCA model, an instance of the CS catheter, $\mathbf{\hat{s}}$, can be generated in ${{\theta}_{1}}$ according to

$$\begin{eqnarray}&&\mathbf{\hat{s}}=\mathbf{\bar{s}}+\underset{m}{\mathop \sum}\,{{w}_{m}}{{\mathbf{v}}_{\mathbf{m}}}\end{eqnarray} \tag{ 2 }$$

where w_m are the weights for each eigenvector. If the first two eigenvectors are used then we can represent the instance by ${{\mathbf{\hat{s}}}_{{{w}_{1}}{{w}_{2}}}}$. The weights ${{\hat{w}}_{1}}$ and ${{\hat{w}}_{2}}$ that represent the cardiorespiratory phase of the image for ${{\theta}_{2}}$ are then estimated according to

$$\begin{eqnarray}&&\widehat{{{w}_{1}}},\widehat{{{w}_{2}}}=\underset{{{w}_{1}}{{w}_{2}}}{\mathop{\text{argmin}}}\,\left[D\left({{\mathbf{\hat{s}}}_{{{w}_{1}}{{w}_{2}}}},\mathbf{BL}\right)+A\left({{\mathbf{\hat{s}}}_{{{w}_{1}}{{w}_{2}}}}\right)\right]\end{eqnarray} \tag{ 3 }$$

where $D\left({{\mathbf{\hat{s}}}_{{{w}_{1}}{{w}_{2}}}},\mathbf{BL}\right)$ is a function that computes the sum of the minimum Euclidean distances between the epipolar lines generated by ${{\mathbf{\hat{s}}}_{{{w}_{1}}{{w}_{2}}}}$ and their nearest corresponding blobs from $\mathbf{BL}$. Its purpose is to position the epipolar lines so that they pass through the electrodes in the new view. $A\left({{\mathbf{\hat{s}}}_{{{w}_{1}}{{w}_{2}}}}\right)$ is a function that computes the sum of angles between line segments joining successive blobs, in order for the optimisation to favour a smoothly curving line of blobs representing the electrodes. Motivated by the $2\sigma$ rule, which states that for a normal distribution, about 95% of the values lie within 2 standard deviations of the mean, our search space for the weights varies from $-2\sqrt{{{d}_{m}}}$ to $+2\sqrt{{{d}_{m}}}$, where d_m is the mth eigenvalue, and $1\leqslant m\leqslant 2$. We minimized the function using a coarse-scale exhaustive search followed by an iterative optimisation using Matlab's lsqnonlin function with the trust-region-reflective algorithm (MATLAB Release 2013a, The MathWorks, Inc., Natick, Massachusetts, USA). For the application of our model on rotational sequences, where the CS catheter position deviation was larger, our search space for the weights varies from $-3\sqrt{{{d}_{m}}}$ to $+3\sqrt{{{d}_{m}}}$, to solve the optimisation. Additionally, for the phantom rotational datasets, the exact view angles for each frame were not available in the DICOM headers. Instead, initial estimates for the angles were computed by assuming a constant rotation speed. These were refined by including a search over the angles ($\pm {{5}^{{}^\circ}}$ from the initial estimate in steps of 1°) in the exhaustive search for the PCA weights.

The resulting weights indicate the cardiorespiratory phases and the blobs from $\mathbf{BL}$ nearest to the epipolar lines give the catheter location in the new image.

${{\hat{w}}_{1}}$ is used to detect systolic frames, ${{\Omega }_{\text{sys}}}$, of the image sequence. These are represented by the peaks of the variation of ${{\hat{w}}_{1}}$ with frame number, i.

$$\begin{eqnarray}&&{{\Omega }_{\text{sys}}}=\left\{i\mid {{w}_{1,i-1}} < {{w}_{1,i}}>{{w}_{1,i+1}}\right\}\end{eqnarray} \tag{ 4 }$$

The peaks of the variation of ${{\hat{w}}_{2}}$ over time represent end-inspiration (EI) respiratory frames, ${{\Omega }_{EI}}$, while the troughs represent end-expiration (EX) respiratory frames, ${{\Omega }_{EX}}$.

$$\begin{eqnarray}&&{{\Omega }_{EI}}=\left\{i\mid {{w}_{2,i-1}} < {{w}_{2,i}} > {{w}_{2,i+1}}\right\}\end{eqnarray} \tag{ 5 }$$

$$\begin{eqnarray}&&{{\Omega }_{EX}}=\left\{i\mid {{w}_{2,i-1}}>{{w}_{2,i}} < {{w}_{2,i+1}}\right\}\end{eqnarray} \tag{ 6 }$$

The accuracy of our gating and catheter tracking were validated at both normal and low dose, on sequential biplane sequence pairs, from clinical and phantom data. The model was formed in one view at normal dose, and used to gate and track the catheter in a second view at either normal or low dose. A comparison to our existing single-view method (Panayiotou et al 2013) was performed to show the effect of the extra step of changing view angle. Comparison to all possible single view methods was not done, as this has already been done in our previous work (Panayiotou et al 2013) and was not our objective in this paper.

Systole was chosen as opposed to diastole for validation of cardiac gating since the manual ground truth is more reliable for end-systole where rapid motion can be used as the visual cue. Although the majority of the experiments looked at only systole, one set of experiments was performed for diastolic accuracy in normal dose clinical images to ensure that the method is comparably accurate on diastolic frames.

Secondly, the accuracy of the technique was validated at normal dose on rotational angiography sequences, from clinical and phantom data. The model was formed on a training view at a constant angle at normal dose, then applied to a rotational sequence with continuously changing angle between RAO90° and LAO90°.

Finally, it is important to investigate the dependence of the phase determination accuracy on the view angle of the image sequence. Consequently, experiments were performed on phantom datasets with different view angles used for the formation of the statistical model in the training view. These were then applied to a rotational phantom dataset in the current view. This allowed comparison of different training angles and also investigation of how the accuracy varies with the current angle.

To validate our technique, gold standard cardiorespiratory signals were generated along with a gold standard for the electrode positions in each x-ray image. The latter was done by the use of the CS tracking technique (Ma et al 2010). Manual tracking was performed in cases where the technique failed to accurately detect electrodes. This was done by manually localizing the centre of any misdetected electrode of the CS catheter. The 2D detection error of the CS tracking technique is $0.39\pm 0.22$ mm over all electrodes (Ma et al 2010). The manual electrode detection gold standard was validated by two independent observers on 330 frames of phantom data (3300 electrode detections) and 122 frames of clinical data (1220 electrodes). The average Euclidean distance between the results of the two observers over all electrodes and frames of these randomly chosen phantom and rotational clinical sequences were 0.47 mm and 0.35 mm, respectively.

For the gating gold standard in clinical data, manual gating of the cardiac cycle at systole was performed by an experienced observer, by visually detecting the end of contraction of the left ventricle from the fluoroscopic left heart border shadow. Manual gating of the cardiac cycle at end-systole was further validated by three additional observers who were trained to identify the end-systolic frames throughout the x-ray sequences. All four chose identical systolic frames, indicating perfect precision of the gold standard. Manual gating of the cardiac cycle at end-diastole was also performed by three experienced observers, by visually detecting the end of relaxation of the left ventricle from the fluoroscopic left heart border shadow. No variation between the multiple observers was found. For the respiratory gating, a gold standard was obtained for the clinical data using either diaphragm or heart border tracking as described in Ma et al (2012). The respiratory signal obtained was cardiac gated using the manually identified systolic frames.

For the phantom sequences, two different clearly-visible regions of interest moving independently with cardiac and respiratory motion were automatically tracked using the diaphragm/heart border tracking method (Ma et al 2012). The gating accuracy of the automatic diaphragm/heart border tracking technique was found to be identical to manual gating based on landmarks attached to the phantom.

For both cardiac and respiratory gating, the absolute frame difference was computed between our technique and the gold standard techniques. Specifically, systolic, end-inspiration, end-expiration, and for one experiment end-diastolic, frames were recorded from our automatic and gold standard methods and their corresponding absolute frame differences were computed. Even though we do detect the whole cycle, for evaluation, we use only these peaks and troughs, where the gold standard is well-defined. Motion gating accuracy objectives were set based on potential clinical applications that the proposed method could tackle. For motion gating of 3D rotational x-ray angiography sequences and 3D catheter reconstruction, cardiac motion is the limiting factor. These applications will use the diastolic cardiac phase, where the heart is more relaxed, and it is expected that its shape will be more repeatable over several cycles. Additionally, when registering a 3D model, the model is almost always from a diastolic image. Since the heart will be relatively stationary in the diastolic phase for a period of about 0.3 s, the motion gating accuracy objective is set to 0.1 s. Faultless gating results are signified when the absolute frame difference is within this motion gating accuracy objective.

2D catheter tracking accuracy was evaluated simply as the distance between the measured electrode positions and gold standard positions. For cardiac catheterization procedures such as RF ablation, CS catheter tracking accuracy within 2 mm is an acceptable tolerance (Esteghamatian et al 2008).

Clinical images were acquired from 11 patients undergoing radiofrequency ablation procedures. Phantom images were acquired using a bespoke beating and breathing left ventricular phantom (Manzke et al 2010) with an inserted CS catheter, designed to emulate the clinical workflow of a typical catheterisation, as illustrated in figure 3(a). For this phantom, polyvinyl alcohol (PVA, Lenticats, GeniaLab, Braunschweig, Germany) (Mano et al 1986, Chu and Rutt 1997) was used as a tissue-mimicking material, to mimic the heart during cardiac interventions. A cylindrical PVA phantom with an inner radius of 6.6 cm, an outer radius of 8.2 cm and a length of about 9 cm, providing an initial wall thickness of about 1.6 cm, was generated using one freeze / thaw cycle at $-35{{~}^{{}^\circ}}\text{C}/\,+\,20{{~}^{{}^\circ}}\text{C}$. The PVA phantom, illustrated in figure 3(b), was placed in an MR-compatible air-pressured actuator (Timinger et al 2004) which compresses and rotates the PVA according to a preset heart rate. Along with the cardiac motion, the phantom actuator can translate the phantom up and down to simulate breathing motion. For the experiments demonstrated in this chapter, the cardiac rate varied between 70 and 100 cycles per minute (cpm), while the respiratory rate varied between 8 and 12 cpm. X-ray imaging was performed at between 3 and 30 frames per second with an image size of either 512² or 1024² pixels in resolution, with a pixel size, ${{R}_{\text{x}-\text{ray}}}$, of 0.25 mm pixel⁻¹. Included in this ratio is the typical magnification factor of the x-ray system.

Zoom In Zoom Out Reset image size

The effect of low x-ray doses was investigated by adding Poisson noise to both the clinical and phantom test images. The magnitude of the noise varies across the image, as it depends on the local image intensity. The image noise is therefore specified in terms of the median value of the Poisson mean (Panayiotou et al 2013). We used SNR values of $\sqrt{50}$, $\sqrt{10}$, $\sqrt{8}$, $\sqrt{6}$ and $\sqrt{5}$, where the SNR value of $\sqrt{5}$ represents a dose reduction of more than 10 times. These values were chosen by initially applying our CS catheter tracking algorithm on low noise levels that are almost indistinguishable from the noise-free image, and then gradually increasing the noise level until the algorithm started to fail.

The multiplane normal and low-dose experiments used 3 sequential biplane clinical sequences (244 frames, 6 runs in total, running from a minimum of 8.0 to a maximum of 19.3 s, covering at least 2 respiratory cycles) and 2 sequential multiplane (one biplane and one triplane) phantom sequences (1741 frames, 5 runs in total, running from a minimum of 9.8 to a maximum of 22.7 s). The PCA model was formed for each sequence individually and the method run on each possible combination of sequence pairs, giving a total of 6 clinical and 8 phantom experiments.

For the rotation experiments, 3 clinical rotational x-ray angiography sequences (341 frames) and 2 phantom rotational x-ray sequences (540 frames) were used. Each sequence covered the range RAO 90° to LAO 90°. For the patient sequences the statistical model was formed on a monoplane x-ray sequence acquired for each patient in the AP view. For the application of our technique on the phantom rotational sequences, the model was formed on monoplane x-ray sequences in AP and RAO30° views in turn.

Finally, the angular dependency experiments used the same phantom rotational sequence. The training views used were PA, LAO 30° and RAO 30°.

Figure 4 illustrates the results of the CS catheter tracking technique on the 1st image of an uncorrupted example x-ray sequence and the simulated low-dose images, with SNR values of $\sqrt{50}$, $\sqrt{10}$, $\sqrt{8}$, $\sqrt{6}$ and $\sqrt{5}$, respectively.

Zoom In Zoom Out Reset image size

Figure 5(a) and 6(a) illustrate the CS catheter tracking median errors per electrode for phantom and patient x-ray images, respectively, at normal dose and for the 5 SNR values. Errors were calculated for each electrode and are shown in different grey scale colours, starting from the proximal (black colour) and moving to the distal (white colour) electrode. For the technique to be acceptable in clinical practice, failure cases were considered to be the ones where median errors per electrode were above 2 mm (Esteghamatian et al 2008). Figure 5(b) and 6(b) show the range of the values in terms of the first quartile (25th percentile) median and third quartile (75th percentile) over all electrodes, for each noise level, for phantom and patient datasets, respectively. Success rates were also calculated and are shown as percentages (%) in figures 7(a) and (b), for the patient and phantom datasets.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

For cardiac gating validation, a plot of the first weight for the first 30 frames of an example clinical x-ray sequence is illustrated in figure 8(a). The result of respiratory gating for the same sequence is shown in figure 8(b).

Zoom In Zoom Out Reset image size

Figure 9 illustrates the frame difference errors in frequency distribution bar charts, i.e. the number of peaks or troughs with 0, 1 or 2 frames separation from the gold standard for the patient biplane sequences, for all gating tasks and noise levels. For both systolic and diastolic gating, the absolute frame difference was computed between the view-angle independent method and the gold standard methods on the normal dose images. The results showed no difference between end-systolic and end-diastolic gating so they are excluded from the figure. This indicates that the method is just as accurate on diastolic frames, and consequently the systolic accuracy results are relevant for applications that require diastolic gating. For a more comprehensive comparison, percentage (%) success rates for end-systolic, end-expiration (EX) and end-inspiration (EI) gating were also calculated and are shown in figure 10 for both the used phantom and clinical images. Using the model on unmodified clinical images for near-real time cardiorespiratory motion gating, end-systole, EX and EI gating success rates of 100.0%, 85.7% and 92.3%, respectively, were established, based on the accuracy objective of 0.1 s. Additionally, gating success rates of 80.3%, 71.4% and 69.2% were established, even at the low SNR value of $\sqrt{5}$. Using the model on unmodified and modified phantom sequences, for near-real time cardiorespiratory motion gating, success rates of 100.0% were established for all gating tasks. To comparatively validate the technique to the Tracked-PCA method (Panayiotou et al 2013), where the statistical model was formed from normal dose images during a calibration phase and applied to low dose images acquired from the same view-angle, % success rates from both techniques for the same sequences are illustrated on the figure. No false positives or negatives (i.e. extra/fewer detected peaks/troughs) over the processed sequences were found. Outcomes show that the technique is robust and accurate in both CS catheter tracking and motion extraction even at the lowest SNR values of $\sqrt{5}$.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Regarding our algorithm's performance on the different experiments, depending on the step size of the optimisation, the execution time was between 0.41 and 1.73 s per frame running in Matlab on Windows 7 with a 3.4 GHz Intel Core i7 CPU and 8 GB of RAM.

Figures 11 and 12 illustrate the results of the CS catheter tracking technique on normal dose x-ray images at different view angles from a patient and a phantom example x-ray sequence, respectively. The results show that, this novel view-angle independent technique can successfully track the CS catheter in x-ray sequences where the angulation of the scanner is changed between frames. The multi-scale blob detection (Lindeberg 1993) and real-time tracking technique (Ma et al 2010) do not reliably detect the CS catheter electrodes in all views, even at normal dose, because of the different appearance of the electrodes at some angles.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 13 illustrates (a) the CS catheter tracking technique median errors per frame per electrode and (b) the % success rates for the patient and phantom x-ray rotational sequences.

Zoom In Zoom Out Reset image size

For both cardiac and respiratory gating, the absolute frame difference was computed for the phantom procedures. For patient procedures the absolute frame difference was computed only for cardiac gating as patients were in breath-hold. The results can be seen in the frequency distribution bar charts in figure 14 for both patient and phantom rotational sequences. No false positives or negatives (i.e. extra/fewer detected peaks/troughs) over all processed sequences were found.

Zoom In Zoom Out Reset image size

The results demonstrate that our proposed technique is robust and accurate in tracking the CS catheter and is doing a faultless job in extracting cardiorespiratory motion from any view-angle, over the range of ${{180}^{{}^\circ}}$, in normal dose scenarios for rotational acquisitions.

For both cardiac and respiratory gating, the gating success rate and the median errors per frame over all CS catheter electrodes were computed for the phantom experiments where the model was formed on different training views, PA, LAO 30° and RAO 30°, and tested on the same rotational sequence, ranging from RAO 90° to LAO 90°. This examines whether the difference in the viewing angles between the training and the current image sequences has an effect on the algorithm's success. The results for gating success rates were 100% in all cases. Therefore, the median errors of CS catheter tracking at each angle were compared and are shown in figure 15. CS catheter tracking is related to motion gating in the sense that a correctly annotated CS catheter would correspond to the best motion gating results we could obtain using our algorithm. CS tracking errors therefore indicate a potential for gating errors. The results show that there is no obvious dependence on the current or training angle. A paired two-sided Wilcoxon signed-rank test was conducted under the null hypothesis that there was no difference in building the model on different training views. The test was done on the median errors over all CS catheter electrodes with respect to viewing angles and was conducted for each of the viewing angle pairs. For all three pairs the p-values were more than 0.05, indicating that there is no statistically significant difference when building the model on different viewing-angles.

Zoom In Zoom Out Reset image size