Patient-specific Conditional Joint Models of Shape, Image Features and Clinical Indicators

by   Bernhard Egger, et al.

We propose and demonstrate a joint model of anatomical shapes, image features and clinical indicators for statistical shape modeling and medical image analysis. The key idea is to employ a copula model to separate the joint dependency structure from the marginal distributions of variables of interest. This separation provides flexibility on the assumptions made during the modeling process. The proposed method can handle binary, discrete, ordinal and continuous variables. We demonstrate a simple and efficient way to include binary, discrete and ordinal variables into the modeling. We build Bayesian conditional models based on observed partial clinical indicators, features or shape based on Gaussian processes capturing the dependency structure. We apply the proposed method on a stroke dataset to jointly model the shape of the lateral ventricles, the spatial distribution of the white matter hyperintensity associated with periventricular white matter disease, and clinical indicators. The proposed method yields interpretable joint models for data exploration and patient-specific statistical shape models for medical image analysis.



There are no comments yet.


page 1

page 2

page 3

page 4


Dynamic multi-object Gaussian process models: A framework for data-driven functional modelling of human joints

Statistical shape models (SSMs) are state-of-the-art medical image analy...

A Wide and Deep Neural Network for Survival Analysis from Anatomical Shape and Tabular Clinical Data

We introduce a wide and deep neural network for prediction of progressio...

MedMNIST Classification Decathlon: A Lightweight AutoML Benchmark for Medical Image Analysis

We present MedMNIST, a collection of 10 pre-processed medical open datas...

Dynamic multi feature-class Gaussian process models

In model-based medical image analysis, three features of interest are th...

Monte Carlo dropout increases model repeatability

The integration of artificial intelligence into clinical workflows requi...

A joint model for multiple dynamic processes and clinical endpoints: application to Alzheimer's disease

Alzheimer's disease, the most frequent dementia in the elderly, is chara...

On the Evaluation and Validation of Off-the-shelf Statistical Shape Modeling Tools: A Clinical Application

Statistical shape modeling (SSM) has proven useful in many areas of biol...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

In medical imaging almost every image is at least weakly labeled with clinical indicators. A minimal set of such labels usually includes age, sex, demographics, and diagnostic or prognostic clinical values such as blood pressure, history of smoking or outcome scores. The labels are of interest for clinical decisions or clinical research, however, they are rarely included in the analysis of shape or image features. We propose to construct joint statistical shape and attribute models to better understand and visualize the interplay between shape deformations and the variations of clinical indicators. The proposed model promises to reveal associations between anatomical shapes, image features, and clinical phenotypes, to help understand disease mechanisms and risk factors, and to improve predictions and patient stratification in clinical trials.

We demonstrate our approach in the context of a clinical cohort of ischemic stroke patients, where we investigate the relationships between the patterns of expansion of the lateral ventricles, the spatial patterns of white matter hyperintensity (WMH) indicative of the white matter disease, and patient data that include stroke outcome scores. Our model is built from segmentations of ventricles and white matter hyperintensities in T2-FLAIR magnetic resonance images.

We integrate copula [13] and Gaussian processes [11] to jointly analyze shape and other variables. Copula methods provide a decomposition of every continuous and multivariate distribution into its marginal distributions and a copula that captures the dependency structure [13]

. The separate modeling of the marginal distributions and of the copula offers better flexibility than the classical statistical shape models. The marginal distributions can be modeled as Gaussian, or other parametric and non-parametric distributions. The choice of the model for each marginal distribution is separate from all others, e.g., a Gaussian distribution for the shape components and the empirical distributions for the clinical indicators. In this work, we use a Gaussian copula, which means we transform the data into a space where all marginal distributions are Gaussian. In this latent space, we then apply standard shape modeling techniques, such as the principal component analysis (PCA) or the Gaussian processes, to model the dependency structure between 3D coordinates of the shape and other indicators. The transformation to the latent space is invertible and our shape model remains generative. We can, therefore, construct the conditional distribution given clinical indicators using Gaussian process regression. Here we focus on the setting where some clinical indicators are given and we capture the conditional distribution of shapes, i.e., a patient-specific statistical shape model.

The main contributions of this paper are: (i) a novel approach for generating conditional shape models based on clinical indicators or image features; (ii) copula-based shape models that include binary, discrete, ordinal and continuous variables; (iii) novel associations between ventricle growth, WMH patterns and clinical indicators in the ischemic stroke patient cohort.

Related Work Although many clinical indicators are collected and available with medical images, the idea of combining statistical shape models with clinical indicators is rarely explored. The work of Blanc et al. [2] includes patient attributes in femur modeling. The method assumes the data is Gaussian distributed and cannot handle binary, discrete, or ordinal variables. Pereanez et al. [10] proposed constraining shape models based on clinical indicators. They demonstrated that such constrained models are beneficial for model-based image segmentation. This approach is limited to a model prediction given a complete set of clinical indicators and handles ordinal variables as a set of binary attributes. Our method combines these ideas to employ the copula and uses the conditional distribution to capture the remaining variability given partial data. Han and Liu [7]

proposed a variant of PCA called Copula Component Analysis that is inherently scale-invariant, more robust to outliers and handles arbitrary discrete ordinal marginal distributions. Hoff


presented the idea of an extended rank likelihood which can also account for non-continuous variables. In computer vision, Egger

et al. [4, 5] used copulas for shape and appearance modeling to account for the non-Gaussian appearance of human faces, however, they did not handle discrete or ordinal variables.

Conditional shape models without clinical indicators are well studied. A Bayesian approach was proposed by Albrecht et al. [1] to model the remaining variability of the femur given partial information. They showed that (probabilistic) PCA and Gaussian process regression with a statistical kernel are equivalent for such modeling. We build upon this model and include clinical indicators to construct patient-specific conditional shape models based on image features and/or clinical indicators. In our experiments, we visualize the conditional models and demonstrate superior performance in a reconstruction task when including all data in the statistical model.

2 Methods

We propose a joint model for anatomical shape, image features and clinical indicators based on the copula [13] and combine it with the idea of conditional shape models based on Gaussian processes [1]. To enable the use of not only continuous, but also binary, discrete and ordinal variables, we use a sampling strategy similar to [8].

We define an instance of our model to be a combination of 3D shape coordinates and image-based feature associated with surface point . All surface points are in correspondence across training data, achieved via shape registration to an atlas in our experiments. In addition each instance includes binary, discrete, ordinal or continuous variables . Let

be the vector representation of the instance,

with . The training set of shapes with image features and clinical indicators forms the data matrix

. Let random variable

represent component of the instance and represents its marginal distribution.

We use a Gaussian copula to capture the dependency pattern in combination with marginal distributions [6, 14]. According to Sklar’s theorem [13] the copula

provides the dependency structure in the joint cumulative distribution function (CDF)

such that


Where is the inverse standard normal CDF and is the joint CDF of a zero-mean multivariate Gaussian distribution with covariance matrix . The multivariate latent representation

distributed according to the standard normal distribution.

Marginal Distributions

To learn the model above from training data we first estimate the marginal distributions with parametric or nonparametric methods. For properties that cannot be naturally captured by a Gaussian distribution (e.g., image features, sex or outcome scores), we employ the empirical marginal distribution for all the components as proposed in


. If only limited data is available, various assumptions on the marginal distribution must be made, separately for each component. In statistical shape modeling, a Gaussian marginal distribution is usually assumed. To transform our data into the latent space, we first transform the empirical marginal distribution into a uniform distribution by replacing every entry of the data by the CDF of the empirical distribution. Second, we use the inverse CDF of a standard normal distribution to map to the latent space


While the original copula framework is limited to continuous distributions, multiple solutions have been proposed to overcome this limitation. The first step, the transformation to a uniform distribution, relies on unique ranks (values of CDF). For binary, discrete and ordinal variables this sorting is not unique. In our experiments, we observed that the latent representation is stable when permuting elements with the same values for ranking, in the setting where we have many more variables than samples (). To construct a consistent dependency structure we generate multiple (50) random rankings for non-unique values and average the resulting estimates of the model parameters. Once the data is transformed to its latent representation we estimate the covariance matrix and represent it by its eigenvectors (i.e. principal components).

Conditional modeling The Gaussian copula represents the dependency structure of our data based on the latent representation . We use Gaussian processes to capture correlations in the latent space [11]. We define function on a finite domain such that refers to the th element of the th training example in our latent space, , . The Gaussian process is fully defined by the mean and covariance functions estimated from the training data:


Given a subset of the vector as an observation, we construct component-wise its latent representation . The prediction of the missing parts reduces to a standard Gaussian process regression problem. We assume observations are corrupted by uncorrelated Gaussian noise and obtain the prediction for the distribution of the complete vector from the observed vector . The conditional distribution is again a Gaussian process with mean:


and covariance



is the identity matrix.

The parameter has an impact on the accuracy of the prediction and varies between different observations. The different modalities in our model have a different amount of uncertainty. For shape observations, is an estimate of the segmentation accuracy and the registration quality, for clinical indicators, it is a measurement of how accurate they are measured. We estimated the optimal parameters for observation uncertainty for each set of observations via cross-validation.

3 Experiments

Figure 1: Joint model of ventricle shape, WMH burden, and clinical indicators. Mean and variation along the first principal mode is shown. a) the original model, b), c), d) the model conditioned on , and respectively. The color indicates the spatially varying amount of WMH burden. The variables per instance are ordered as age, sex, hypertension, hyperlipidemia, atrial fibrillation, smoking, NIHSS, and mRS.

We produce a statistical model based on a dataset of T2-FLAIR magnetic resonance images from 793 acute ischemic stroke patients, without gross pathologies, such as mass effects, from the Genes Affecting Stroke Risk and Outcomes Study (GASROS), for which patients were enrolled at Massachusetts General Hospital between 2003 and 2011 [15]. Manual WMH segmentation was performed using MRIcro [12], and ventricles were automatically segmented (Dice: 0.89), using a 3D U-Net [3]. Subsequently, we isolated periventricular from subcortical WMH, based on connectivity of individual WMH lesions to the ventricles. After segmentation, we establish correspondence between every ventricle segmentation and the atlas ventricles by fitting a Gaussian process deformation model estimated via the iterative closest point algorithm [9]. The ventricles are represented as a surface triangle mesh with 1504 vertices. The WMH is modeled as a feature based on the vertices of the mesh. Every voxel labeled as periventricular WMH is assigned to the closest vertex point. This results in a WMH representation as a surface feature. In addition to the shape of the ventricles and the WMH burden feature, we include clinical indicators that represent age, sex, hypertension (binary), hyperlipidemia (binary), atrial fibrillation (binary), smoking (ordinal, 1-6) and NIHSS (stroke score, ordinal 0-42), mRS (outcome score, ordinal 0-6) and the overall WMH volume. With 1504 3D ventricle vertices, 1504 WMH burden features and 9 clinical indicators, the total dimensionality of a single data point is .

Fig. 1 presents the joint model of the ventricle shape, WMH burden features, and clinical indicators and also three conditional models for different clinical indicators, like age and sex. The conditional models have similar mean and variability as the unconditioned model when conditioned on features with low correlations to the remaining components Fig. 1b and show different means and less variability when conditioned on stronger correlated Fig. 1c or combined features Fig. 1d. The supplementary material contains a video111 of our graphical user interface for data exploration. It enables the user to fix clinical indicators and visualize the remaining variability in the conditional model.

mean ventricles WMH indicators ind + vol combined
Age 15.77 9.49 11.27 6.62 12.06 7.58 13.21 8.19 11.57 7.04 9.85 5.76
mRS 2.01 1.50 1.88 1.37 1.96 1.45 1.65 1.20 1.65 1.20 1.65 1.19
Sex 60.6% 62.7% 61.7% 65.3% 64.8% 67.0%
Table 1:

Clinical indicator prediction. Mean error and standard deviation are reported for the prediction of age, stroke outcome (mRS) and sex. The columns of the matrix correspond to the prediction based on the mean of the marginal distribution and conditioned on the ventricle shape, WMH burden, all other clinical indicators, indicators combined with WMH volume (ind + vol), and all available components except for the one being predicted. For sex, the percent of correct predictions is reported.

mean ventricles WMH indicators ind + vol combined
ventricles 1582.2 34.5 - 132.9 35.0 146.9 45.1 144.8 43.0 132.3 34.8
WMH 3980.1 2830 2927 2270 - 3073 2419 1291 897 1254 866
Table 2: Ventricle Shape and WMH burden prediction. The setting is the same as in Table 1 but we predict the shape of the ventricles and the WMH burden from the remaining values. Mean and standard deviation of distances between mesh vertices in mm is reported for the ventricle shape and WMH feature distance in voxels.

We evaluate the joint model in a reconstruction task of estimating the full vector from partial data in three different tasks: attribute prediction (Table 1), shape prediction, and WMH burden prediction (Table 2). For all tasks, we use the mean of the conditional distribution in the latent space to generate predictions. The values are obtained by randomly splitting the dataset into a training set of 600 patients and a validation set of the remaining 193 patients. The attribute prediction task (Table 1) aims to estimate clinical indicators, e.g., stroke outcome, given the shape and WMH burden. We observe that the ventricle shape and WMH burden are strong indicators for age whilst they contribute little to the stroke outcome score. For age and sex, using the combination of all modalities provides the best prediction. The shape and WMH prediction tasks (Table 2) reveal how the imaging data from a patient relates to the clinical indicators. The shape of the ventricles can be predicted from the WMH burden. The prediction in the other direction is more challenging, due to strong variability in WMH burden in the population. Clinical indicators including WMH volume combined with the ventricle shape produce the best prediction.

Figure 2: Sample distributions of age, mRS (outcome score), hypertension and the log WMH volume for unconditional and conditional model conditioned on , and respectively.

In the third experiment, we demonstrate the difference between the unconditional and conditional models. We produce the sample distributions by drawing 1,000 random samples from the corresponding model. Fig. 2 reports those sample distributions before and after conditioning on the same three events as in Fig. 1. We observe that the conditioning on male has a weak effect whilst the conditioning on age or a combination of age and sex has a strong effect on the resulting marginal distributions.

We implemented the copula model based on the code provided in [4] in the Scalismo222Scalismo - A Scalable Image Analysis and Shape Modeling Software Framework
framework for statistical shape modeling based on Gaussian processes. The computation of the conditional shape models based on the clinical indicators takes a few seconds for our dataset on a single core. Reducing the model via PCA will speed it up to real-time performance, enabling interactions with the model.

Limitations The proposed model has two main limitations. First, we assume dense correspondence on the surface of the ventricles as well as correspondence across image features. This is especially limiting if some shape features have weak or no correspondence across the population. Second, the Gaussian model for the dependency structure we are using is limited to second-order dependency and does not reflect higher-order dependency. These two limitations and assumptions are very common in the statistical shape modeling, especially when working with small datasets.

4 Conclusion

We present a joint generative model of shape, image features, and clinical indicators. Conditional models for partially observed data can be straightforwardly derived. The copula can be implemented as a simple pre- and postprocessing step in existing pipelines relying on (probabilistic) PCA or Gaussian processes. In addition, conditional models can be constructed efficiently for interactive data exploration. Our work provides a highly effective approach to exploring high dimensional data, especially spatial correlations, for example between ventricle expansion and WMH burden patterns. Statistical shape models are commonly used as statistical priors for image analysis tasks, and the proposed conditional models can be interpreted as patient-specific statistical shape models for patient-specific medical image analysis. Our model includes binary, discrete, ordinal and continuous variables and hence can handle a vast variety of available clinical indicators.

Aknowledgments This research was funded by SNSF P2BSP2_178643, NIH NIBIB NAC P41EB015902, NIH NINDS R01NS086905, Horizon2020 753896, De Drie Lichten 24/18, ZonMw 104003005, Wistron Corporation, AWS, SIP and NVIDIA.


  • [1] Albrecht, T., Lüthi, M., Gerig, T., Vetter, T.: Posterior shape models. Medical image analysis 17(8), 959–973 (2013)
  • [2] Blanc, R., Seiler, C., Székely, G., Nolte, L., Reyes, M.: Statistical model based shape prediction from a combination of direct observations and various surrogates: application to orthopaedic research. Medical image analysis 16(6), 1156–1166 (2012)
  • [3] Dubost, F., de Bruijne, M., Nardin, M., Dalca, A.V., Donahue, K.L., et al.: Automated Image Registration Quality Assessment Utilizing Deep-learning based Ventricle Extraction in Clinical Data. arXiv e-prints arXiv:1907.00695 (Jul 2019)
  • [4] Egger, B., Kaufmann, D., Schönborn, S., Roth, V., Vetter, T.: Copula Eigenfaces with Attributes: Semiparametric Principal Component Analysis for a Combined Color, Shape and Attribute Model, pp. 95–112. Springer (2017)
  • [5] Egger, B., Kaufmann, D., Schönborn, S., Roth, V., Vetter, T.: Copula eigenfaces - semiparametric principal component analysis for facial appearance modeling. In: Proceedings of the 11th Joint Conference on Computer Vision, Imaging and Computer Graphics Theory and Applications: Volume 1: GRAPP. pp. 50–58. SCITEPRESS-Science and Technology Publications, Lda (2016)
  • [6] Genest, C., Ghoudi, K., Rivest, L.P.: A semiparametric estimation procedure of dependence parameters in multivariate families of distributions. Biometrika 82(3), 543–552 (1995)
  • [7] Han, F., Liu, H.: Semiparametric principal component analysis. In: Advances in Neural Information Processing Systems. pp. 171–179 (2012)
  • [8] Hoff, P.D.: Extending the rank likelihood for semiparametric copula estimation. The Annals of Applied Statistics pp. 265–283 (2007)
  • [9] Lüthi, M., Gerig, T., Jud, C., Vetter, T.: Gaussian process morphable models. IEEE PAMI 40(8), 1860–1873 (2018)
  • [10] Pereañez, M., Lekadir, K., Albà, X., Medrano-Gracia, P., Young, A.A., Frangi, A.: Patient metadata-constrained shape models for cardiac image segmentation. In: STATCOM, MICCAI. pp. 98–107. Springer (2015)
  • [11]

    Rasmussen, C.E.: Gaussian processes in machine learning. In: Summer School on Machine Learning. pp. 63–71. Springer (2003)

  • [12] Rost, N.S., Rahman, R.M., Biffi, A., et al.: White matter hyperintensity volume is increased in small vessel stroke subtypes. Neurology 75(19), 1670–1677 (2010)
  • [13] Sklar, M.: Fonctions de répartition à n dimensions et leurs marges. Université Paris 8 (1959)
  • [14] Tsukahara, H.: Semiparametric estimation in copula models. Canadian Journal of Statistics 33(3), 357–375 (2005)
  • [15] Zhang, C.R., Cloonan, L., Fitzpatrick, K.M., Kanakis, A.S., Ayres, A.M., Furie, K.L., Rosand, J., Rost, N.S.: Determinants of white matter hyperintensity burden differ at the extremes of ages of ischemic stroke onset. Journal of Stroke and Cerebrovascular Diseases 24(3), 649–654 (2015)