A Study of Coronary Bifurcation Shape in a Normal Population.

During percutaneous coronary intervention, stents are placed in narrowings of the arteries to restore normal blood flow. Despite improvements in stent design, deployment techniques and drug-eluting coatings, restenosis and stent thrombosis remain a significant problem. Population stent design based on statistical shape analysis may improve clinical outcomes. Computed tomographic (CT) coronary angiography scans from 211 patients with a zero calcium score, no stenoses and no intermediate artery, were used to create statistical shape models of 446 major coronary artery bifurcations (left main, first diagonal and obtuse marginal and right coronary crux). Coherent point drift was used for registration. Principal component analysis shape scores were tested against clinical risk factors, quantifying the importance of recognised shape features in intervention including size, angles and curvature. Significant differences were found in (1) vessel size and bifurcation angle between the left main and other bifurcations; (2) inlet and curvature angle between the right coronary crux and other bifurcations; and (3) size and bifurcation angle by sex. Hypertension, smoking history and diabetes did not appear to have an association with shape. Physiological diameter laws were compared, with the Huo-Kassab model having the best fit. Bifurcation coronary anatomy can be partitioned into clinically meaningful modes of variation showing significant shape differences. A computational atlas of normal coronary bifurcation shape, where disease is common, may aid in the design of new stents and deployment techniques, by providing data for bench-top testing and computational modelling of blood flow and vessel wall mechanics.

Introduction

Cardiovascular disease is the leading global cause of death, accounting for an estimated 17.5 million deaths in 2012, and coronary artery disease represented 42 % of these deaths [1]. Percutaneous coronary intervention (PCI) is a well-established non-surgical medical procedure to deliver a stent to hold open occluded arteries.

However, in-stent restenosis of bare metal stents (BMS) and late thrombosis of drug-eluting stents (DES) are important limitations of coronary stenting today with poor stent expansion being a predictor of restenosis [2].

Both stent over-sizing and under-sizing can have detrimental effects. Larger stents induce more trauma to vessels and can lead to coronary dissection and/or rupture [3]. Underexpanded or undersized stents increase the risk of both restenosis and the likelihood of stent thrombosis [2]. It is therefore important for the cardiovascular community to understand and quantify the anatomy of the coronary vessels, in particular, that of bifurcations where atheromatous narrowing is common. In the USA, there were an estimated 954,000 PCI procedures [4] in 2010. In a large European study with nearly 100,000 stent deployments, restenosis rates were estimated at 7.4 % for BMS and up to 5.8 % for DES [5]. With the advent of computational atlasing and associated tools, we can deliver an innovative and accurate database to assist in the informed improvement of stent technology. The aims of this study are therefore to (1) characterise the coronary anatomy for future stent design, and to (2) establish a large normal database of cases to which future models with disease can be compared.

Most studies in the literature rely on anatomical measures derived from patients with established disease, use sub-optimal two-dimensional projection imaging techniques, and suffer from high inter-observer variability (e.g. [6]). Previously, we published selected basic measures such as angles and diameters [7], and we investigated the use of centrelines and cross-sectional shape data to create an average coronary artery model using tree-spaces [8]. In this work, we present a robust methodology to register any computer-aided design (CAD) to a correspondenceless population of shapes, and apply it to build a large expert-reviewed fully 3D atlas of vessel bifurcations. This atlas allows us to test for shape differences in a holistic fashion by using clinical risk factors.

The contributions of this paper are thus as follows: (1) a robust and semi-automatic correspondenceless or landmark-free pipeline for coronary bifurcations, (2) a high-resolution statistical shape model of the major coronary bifurcations built from 446 instances, (3) an association analysis of clinical risk factors with shape and (4) a physiological analysis of diameter relationships in bifurcations.

Methods

Sample Data

Computed tomographic coronary angiography (CTCA) is routinely used in the assessment of coronary heart disease due to its superior spatial resolution (0.4 mm for 64 multi-slice scanners) compared with other imaging modalities [9]. A coronary calcium score is typically included in CTCA and has been shown to be highly specific for predicting prognosis, with zero-calcium patients (mean age 58 years, sample size 35,765) having an estimated 10-year event rate of only 0.3 % [10].

A population of cases with zero calcium, no stenoses (as confirmed by experienced cardiologists) and no intermediate artery (no trifurcations) was retrospectively selected from our CTCA database. All patients gave written consent under locality approval by the University of Auckland (reg. no. 8605). The population demographics were 67 % female, age years, weight kg and height cm.

Segmentation

The left coronary artery originates from the aortic sinus and splits into two main vessels, the left circumflex (LCX) and the left anterior descending (LAD). This bifurcation is referred to as the left main (LMB). In our database, containing over 300 cases, approximately 30 % of subjects have a third branch (trifurcation) termed intermediate, but these were excluded from this study. Typically, the LAD presents one or more child vessels named diagonals. Analogously, the LCX gives rise to obtuse marginals. We only considered the first diagonal (D1) and first obtuse marginal (OM1) arteries. A standard operational procedure ensured that the largest branch was labelled as the first obtuse marginal or diagonal. Subsequent branches were identified and labelled but were not included in this bifurcation study. In the right coronary artery, the most important bifurcation from an interventional point of view typically occurs at the crux of the heart, splitting into the posterior descending and the postero-lateral branches; we refer to this bifurcation as simply the (right-coronary) crux.

CTCA scans were semi-automatically segmented as previously described [8]. Briefly, a virtual catheter was placed in the left coronary ostium in the CT volumes (Fig. 1 left). Centrelines were then automatically created and further refined manually. Luminal meshes were generated from these centrelines and the CTCA scans using customised level-set software [11]. All the cases were visually checked for quality and manually labelled. An example of the generated mesh can be seen in Fig. 1 right, typically comprising 150,000 vertices and 50,000 triangles.

Fig. 1

Left: Screen-shot of the segmentation tool MIA Lite used to extract the luminal mesh from the centrelines and the CT volume [11]. Right: Example of the resulting segmentation with the left-main bifurcation highlighted in red. Units in mm

The European Bifurcation Club published in 2014 a naming convention for the planar angles found in any coronary bifurcation [12]. These are defined according to the hierarchy of each vessel (Fig. 2). Blood flows from the primary main vessel (PMV) through to the distal main vessel (DMV), giving rise to a smaller branch termed the sub-branch (SB). In the case of the LAD, the SB corresponds to the diagonal, and in the LCX, to the obtuse marginal. The angle between the DMV and the SB is named bifurcation angle or simply B (also ). The angle between the PMV and SB is angle A, and between the PMV and DMV, angle C. Other out-of-plane angles are also present, as the heart is curved, but do not have a consensus name yet.

Fig. 2

Angle definitions as per the European Bifurcation Club recommendation [12]. PMV primary main vessel, DMV distal main vessel, SB sub-branch

Registration

An idealised computer-aided design (CAD) model was generated using the average diameters at the bifurcation entrance and two outlets, and the average bifurcation angle from the dataset, analogous to our previous work [8]. This template comprised 3,057 vertices and 6,110 faces (Fig. 3) and defined our atlas coordinates.

Fig. 3

CAD template shown at the bottom in atlas coordinates. The template (shown in blue) was registered rigidly and then non-rigidly sequentially to each bifurcation (red) in DICOM coordinates using coherent point drift

Bifurcation points were computed automatically from the centrelines. Luminal meshes were cut using a sphere with a radius of 10 mm centred at the bifurcation point (Fig. 1 right). Our dataset shows that 10 mm typically includes the left-main trunk and excludes the next bifurcation. These were considered the ‘bifurcation models’. If another bifurcation was present in the sphere, the bifurcation was discarded from the sample (% cases).

Coherent point drift (CPD) is a probabilistic registration algorithm capable of registering two point sets with different cardinality while preserving their topological structure by fitting Gaussian mixture models (GMM) to their centroids and maximising their likelihood [13].

If we consider as the dataset of points from the template and (with cardinality ) as the i-th bifurcation model, CPD models the points in as the GMM centroids, the probability density function of which is: where is the number of points (defined at the vertices) in the template, (uniform probabilities) and a three-dimensional Gaussian distribution with isotropic covariances . To account for noise and outliers, the model also incorporates an additional weighted noise term giving rise to the mixture model: Parameters are then estimated by the expectation-maximisation algorithm. By imposing different regularisation constraints, two versions of the algorithm are publicly available, rigid and non-rigid. The non-rigid algorithm regularises the displacement field between point sets, from a motion coherence theory perspective [13].

First, the template was rigidly aligned (noise weight , maximum iterations 100, likelihood tolerance ) to each bifurcation mesh giving rise to an affine patient-to-template transformation (including isotropic scaling) mapping for the i-th patient. Second, the result was subsequently registered non-rigidly (additional regularisation parameters were , ) to establish fine local point correspondences, yielding the mapping . Finally, all template points were mapped bijectively to patient data (Fig. 3). Registration parameters were empirically determined to ensure convergence and visual quality. Registration error was measured by Euclidean distances between corresponding points.

Statistical Shape Analysis

The resulting registered models (Fig. 3) were aligned to atlas coordinates using the inverse rigid transformations ; hence, the final functional composition was .

All models were therefore represented in atlas space by 3057 points preserving their warped topology. These aligned and scaled shapes were statistically analysed by principal component analysis (PCA) [14], a matrix factorisation which linearly decomposes the covariance matrix to achieve maximum variance explained by its eigenvectors, whilst maintaining orthogonality and statistical independence.

By projection, a score can be calculated for each eigenvector, representing how much of that shape feature is present in a bifurcation model. Thus, each model can be quantified in each shape pattern or direction. These scores can be used as a compact alternative representation of the models’ points. Shape indices thus comprised the scaling factors arising from the rigid registration (s) and the PCA scores corresponding to the first three modes of variation (PC1, PC2 and PC3). These were tested for associations with bifurcation type (LMB/D1/OM1/Crux) and clinical variables: sex (M/F), hypertension (Y/N), smoking and diabetes status (Y/N) using one-way analysis of variance (ANOVA) with Scheffé’s post-hoc analysis with [15]. A p value of was considered significant.

In interventional cardiology, the relationship between the inlet () and outlet diameters (, ) in the bifurcation area is of importance for accurate percutaneous treatment [16]. The four diameter laws currently described in the literature were tested with the atlas.

Results

The segmentation of a typical case took 15–20 min for the centrelines (semi-manual) and 5 min for the mesh computation (automatic). Registration time was approximately 30 sec for the rigid part and 2 min for the non-rigid on a Dell Precision T5600 system with 16 GB of memory and two Intel Xeon E5-2630 processors.

The median registration error was 0.27 mm with a median standard deviation of 0.17 mm and a median maximum error of 1.36 mm. Maximum errors occurred in areas of high deformation such as unusually large or curved ostium shapes. In total, there were bifurcation models as follows: 166 LMB, 106 D1, 67 OM1 and 107 cruxes.

Having removed isotropic size (scaling factor), the average model and first three modes of variation are shown in Fig. 4, and animations can be seen in the supplementary data. The total variation explained by these three modes was 71 %.

Fig. 4

Top model represents the average bifurcation (N = 446). Below are modes of variation for [−2σ,+ 2σ] shown at the same scale. Total variation (power) is shown in (%)

Scaling factors and the first three PCA modes of variation are tabulated by bifurcation type and sex in Table 1. No significant associations were found between these shape indices and hypertension, smoking or diabetes status.

Table 1

Descriptive statistics (μ±σ) of shape indices by bifurcation type and sex. Significant ANOVA tests are marked with * (p < 0.01). Super-indices show homogeneous subsets using Scheffé’s post-hoc analysis with α = 0.01. LMB left main bifurcation, D1 first diagonal, OM1 first obtuse marginal, Crux right coronary crux

By type	LMB	D1	OM1	Crux
Scale factor s*	0.998±0.098^1	0.935±0.045^2	0.938±0.070^2	0.934±0.061^2
PC1 ‘Angle B’*	44.9±61.0^1	−30.3±40.2^2	−23.4±43.8^2	−25.4±51.9^2
PC2 ‘Angle A/C’*	12.4±32.0^1	3.2±36.0^1	2.5±41.4^1	−23.9±37.9^2
PC3 ‘Curvature’*	−8.9±26.2^1	−3.7±31.3^1	0.4±37.3^1	17.0±25.6^2
By sex	F	M
Scale factor s*	0.940±0.072	0.999±0.087
PC1 ‘Angle B’*	−5.7±60.3	11.7±65.6
PC2 ‘Angle A/C’	−0.1±41.0	0.5±32.9
PC3 ‘Curvature’	−2.6±33.7	5.4±23.0

A scatter plot of the first three principal components by bifurcation type is shown in Fig. 5.

Fig. 5

Score distributions. Points represent bifurcation scores according to the first three modes of variation, colour-coded by bifurcation type. Contour plots are kernel estimations of the joint probability density functions. PC1 ∼ ‘Angle B’; PC2 ∼ ‘Angle A/C’; PC3 ∼ ‘Curvature’; LMB left main bifurcation, D1 first diagonal, OM1 first obtuse marginal, RCA right coronary artery crux

Bifurcation Diameter Laws

The diameter data from our atlas were in close agreement with Finet’s ratio where , a common, yet generally inaccurate, practical reference [17]. This ratio was in our data sample. In comparison with other diameter ratio models based on flow [16], a least-square data fit to an exponential model yielded a minimum for , closest to the Huo-Kassab model (Fig. 6).

Fig. 6

Relationships between bifurcation diameters [16]. The vertical axis shows the increase in diameter between the proximal and two daughter vessels. The horizontal axis shows the ratio of the smaller over the larger daughter vessel. Light blue dots represent the bifurcations from the atlas. These were horizontally binned into 50 segments spanning their range where their median was taken (blue squares). Four geometric laws are compared: Finet’s, Huo-Kassab’s (HK), Murray’s and area-preservation (AP). The mathematical boundary () of the exponential models is also shown in light blue

Discussion

We have presented a registration methodology based on CPD which is applicable to any collection of shapes and have applied it to create a probabilistic model of the major coronary bifurcations from over 211 patient models, comprising 446 bifurcations in total.

Traditional clinical approaches to the quantification of vascular bifurcation shape focus almost entirely on simple geometric parameters such as diameter and angle, with their averages determined independently within various populations [6, 7].

Clinical questions increasingly require a more sophisticated description of shape including (i) the relationship between features such as angle, curvature and diameter, (ii) the definition of shape statistics such as average + / two standard deviations, and (iii) elucidation of the most important modes of variation.

Stent design benefits from this information by informing the sizes (and future shapes) required to span the population. Individual stent selection can be informed by placing the specific patient geometry within the statistical space (Fig. 5) where available stent designs have been mapped, and the deployment technique can be altered depending on the detailed geometry.

The dominant features of the first three PCA modes of variation can be clinically interpreted (Fig. 4). In mode 1 (43 % of the variation found in the sample), the dominant feature is the bifurcation angle B, a well-known parameter previously linked to clinical stenting outcomes [6]. The second mode primarily represents a ratio between angles A and C while maintaining a constant angle B, and mode 3 is possibly related to the curvature of the heart, with smaller hearts expected to have greater curvature.

When observing these modes numerically through their scores (Table 1) and by statistical analysis, PC1 or ‘Angle B’ was found to be significantly different between the LMB and the other bifurcations. Conversely, results for PC2 or ‘Angle A/C’ and the curvature mode (or PC3) separate the right-coronary crux from the rest of bifurcations. These observations can be visualised in Fig. 5 where the different distributions are encoded by bifurcation type. In the distribution of PC1 (bifurcation angle), the LMB group presents wider variation; for PC2, the distributions appear similar to PC3 (curvature) where the RCA group presents a more separate and compact cluster. The former can be explained by the wide variation in bifurcation angle present in LMB, which has to cover a wider myocardial perfusion area. The latter can be explained by its position in the heart. Imagining the heart to be a sphere, most right-coronary cruxes were found to score toward the positive side () of mode 3 (Fig. 4), following the curvature of the inferior side of the heart. Conversely, the rest of the bifurcations are typically found in the northern hemisphere of this sphere, showing negative scores (→−2σ) for mode 3.

These observations suggest that in the future, stents could be selected by mapping their regions of applicability onto PCA mode figures, and then adding the individual patient geometry to be treated to make the best selection from the available alternatives. These statistics also present an opportunity for bench-top testing of stents or hemodynamic simulations in population-representative phantom models [18].

One potential limitation of this work is that the anatomy has been modelled from a population of zero-calcium score with no stenoses which will have a very low incidence of coronary heart disease [19]. This does however establish a statistical map of normal anatomy for comparison with pathological groups, and may be used to address the questions of whether certain shapes (or modes) are associated with a greater incidence of disease. Determining which, if any, geometric features change in the presence of disease remains future work.

The concatenation of the rigid and non-rigid CPD algorithms provided highly accurate correspondences across point sets with disparate cardinality and permitted the alignment of the relatively featureless models for shape analysis. Regarding the scaling factors, post-hoc analysis of variance shows that the LMB required significantly larger values of isotropic scaling to fit the template than the rest (Table 1). This is expected as the diameters decrease along the vessels and the LMB is the first bifurcation by size order.

Although the effect of the template bias has been purposely reduced by using a down-sampled average CAD model, further work should address the minimisation of this bias by means of template iteration, e.g. [20]. Other non-linear statistical models could also be used. A future registration algorithm could incorporate other features such as centrelines, bifurcation point and tangential directions to optimise computation time and accuracy.

A closer analysis of the bifurcation diameter relationships revealed a physiological domain map where most atlas bifurcations clustered below Finet’s 0.678 boundary and approximately follow Huo-Kassab’s exponential model. A justification for the higher error found in Finet’s model has been previously theorised [16, 17]. Briefly, Finet’s law represents a linear relationship among the diameters in a bifurcation which can only approximate, by construction, the complex exponential relationships arising from vascular flow. Finet’s law is therefore only applicable in a small area of the domain shown in Fig. 6.

Conclusions

In conclusion, an atlas of normal coronary bifurcation anatomy delivers clinically important insights, defining range and distribution of shape. It also confirms, in a quantifiable and computational fashion, a long-standing question in interventional cardiology, i.e. there are statistically significant shape differences among the major bifurcations. This supports the case for dedicated stent designs and enables more accurate flow simulations.

Below is the link to the electronic supplementary material.

(GIF 546 KB)

(GIF 544 KB)

(GIF 589 KB)

Model	Geometric relation	RMS error (binned data)
Murray	DPMV3=DDMV3+DSB3	0.1693
Area-preservation (AP)	DPMV2=DDMV2+DSB2	0.1685
Huo-Kassab (HK)	DPMV7/3=DDMV7/3+DSB7/3	0.1650
Finet	D PMV=0.678(D DMV + D SB)	0.4564