METHOD FOR DETECTING ADVERSE CARDIAC EVENTS

Abstract

A method (1) is described for training a machine learning model (2) to receive as input a time-resolved three-dimensional model (4) of a heart or a portion of a heart, and to output (3) a predicted time-to-event or a measure of risk for an adverse cardiac event. The method includes receiving a training set (5). The training set (5) includes a number of time-resolved three-dimensional models (4.sub.1, . . . , 4.sub.N) of a heart or a portion of a heart. The training set (5) also includes, for each time-resolved three-dimensional model (4.sub.1, . . . , 4.sub.N), corresponding outcome data (7.sub.1, . . . , 7.sub.N) associated with the time-resolved three-dimensional model (4.sub.1, . . . , 4.sub.N). The method (1) of training a machine learning model (2) also includes, using the training set (5) as input, training the machine learning model (2) to recognise latent representations (12) of cardiac motion which are predictive of an adverse cardiac event. The method (1) of training a machine learning model (2) also includes storing the trained machine learning model (2).

Claims

1. A method of training a machine learning model to: receive as input a time-resolved three-dimensional model of a heart or a portion of a heart; and output a predicted time-to-event or a measure of risk for an adverse cardiac event; the method comprising: receiving a training set which comprises: a plurality of time-resolved three-dimensional models of a heart or a portion of a heart, for each time-resolved three-dimensional model, corresponding outcome data associated with the time-resolved three-dimensional model; using the training set as input, training the machine learning model to recognise latent representations of cardiac motion which are predictive of an adverse cardiac event; storing the trained machine learning model.

2. A method according to any claim 1, wherein each time-resolved three-dimensional model comprises a plurality of vertices, each vertex comprising a coordinate for each of a plurality of time points; wherein each time-resolved three-dimensional model is input to the machine learning model as an input vector which comprises, for each vertex, the relative displacement of the vertex at each time point after an initial time point.

3. A method according to claim 1, wherein the machine learning model comprises an encoding layer configured to encode latent representations of cardiac motion.

4. A method according to claim 3, wherein the machine learning model is configured so that the output predicted time-to-event or measure of risk for an adverse cardiac event is determined using a prediction branch which receives as input the latent representation of cardiac motion encoded by the encoding layer.

5. A method according to claim 1, wherein the machine learning model comprises a de-noising autoencoder.

6. A method according to claim 3, wherein the machine learning model is trained according to a hybrid loss function which comprises a weighted sum of: a first contribution determined based on the input time-resolved three-dimensional models and corresponding reconstructed models of cardiac motion, each reconstructed model determined based on the latent representations of cardiac motion encoded by the encoding layer; and a second contribution determined based on the outcome data and the corresponding outputs of predicted time-to-event or measure of risk for an adverse cardiac event.

7. A method according to claim 1 wherein training the machine learning model comprises optimising one or more hyperparameters selected from the group consisting of: a predetermined fraction of inputs to the machine learning model which are set to zero at random; a number of nodes included in a hidden layer of the machine learning model; the dimensionality of an encoding layer which encodes a latent representation of cardiac motion; weights of the first and second contributions to the hybrid loss function; a learning rate for training the machine learning model; and an l.sub.1 regularization penalty used for training the machine learning model.

8. A method according to claim 7, wherein optimising one or more hyperparameters comprises particle swarm optimisation.

9. A method according to claim 1, wherein the machine learning model is trained to output a predicted time-to-event or a measure of risk for an adverse cardiac event associated with heart dysfunction.

10. (canceled)

11. A non-transient computer-readable storage medium storing a machine learning model trained to receive as input a time-resolved three-dimensional model of a heart or a portion of a heart, and to output a predicted time-to-event or a measure of risk for an adverse cardiac event.

12. A method comprising: receiving a time-resolved three-dimensional model of a heart or a portion of a heart; providing the time-resolved three-dimensional model to a trained machine learning model, the trained machine learning model configured to recognise latent representations of cardiac motion which are predictive of an adverse cardiac event; obtaining, as output of the trained machine learning model, a predicted time-to-event or a measure of risk for an adverse cardiac event.

13. A method according to claim 12, wherein the time-resolved three-dimensional model comprises a plurality of vertices, each vertex comprising a coordinate for each of a plurality of time points; wherein the time-resolved three-dimensional model is input to the trained machine learning model as an input vector which comprises, for each vertex, the relative displacement of the vertex at each time point after an initial time point.

14. A method according to claim 12, wherein the trained machine learning model comprises an encoding layer configured to encode a latent representation of cardiac motion.

15. A method according to claim 14, wherein the trained machine learning model is configured so that the output predicted time-to-event or measure of risk for an adverse cardiac event is determined using a prediction branch which receives as input the latent representation of cardiac motion encoded by the encoding layer.

16. A method according to claim 12, wherein the machine learning model further outputs a reconstructed model of cardiac motion.

17. A method according to claim 12, wherein the trained machine learning model comprises a de-noising autoencoder.

18. A method according to claim 12, wherein the trained machine learning model is configured to output a predicted time-to-event or a measure of risk for an adverse cardiac event associated with heart dysfunction.

19. (canceled)

20. A method according to claim 12, further comprising: obtaining a plurality of images of a heart or a portion of a heart, each image corresponding to a different time or a different point within a cycle of the heart; generating the time-resolved three-dimensional model of the heart or the portion of the heart by processing the plurality of images using a second machine learning model.

21. A method according to claim 12, wherein the trained machine learning model is a machine learning model trained using steps comprising: receiving a training set which comprises: a plurality of time-resolved three-dimensional models of a heart or a portion of a heart, for each time-resolved three-dimensional model, corresponding outcome data associated with the time-resolved three-dimensional model; using the training set as input, training the machine learning model to recognise latent representations of cardiac motion which are predictive of an adverse cardiac event.

22. A method according to claim 5, wherein the machine learning model is trained according to a hybrid loss function which comprises a weighted sum of: a first contribution determined based on the input time-resolved three-dimensional models and corresponding reconstructed models of cardiac motion, each reconstructed model determined based on the latent representations of cardiac motion encoded by the encoding layer; and a second contribution determined based on the outcome data and the corresponding outputs of predicted time-to-event or measure of risk for an adverse cardiac event.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0060] Certain embodiments of the present invention will now be described, by way of example, with reference to the accompanying drawings in which:

[0061] FIG. 1 illustrates a method of training a machine learning model;

[0062] FIG. 2 illustrates a method of using a machine learning model;

[0063] FIG. 3A shows examples of automatically segmented cardiac images;

[0064] FIG. 3B shows examples of time resolved three-dimensional models;

[0065] FIG. 4A shows Kaplan-Meier plots of survival probabilities for subjects in a clinical study, obtained using a conventional parameter model;

[0066] FIG. 4B shows Kaplan-Meier plots of survival probabilities for subjects in a clinical study, obtained using an exemplary machine learning model (herein termed the 4Dsurvival network);

[0067] FIG. 5A shows a 2-dimensional projection of latent representations 12 of cardiac motion derived and used by the 4Dsurvival network;

[0068] FIG. 5B shows saliency maps derived for the 4D survival network;

[0069] FIG. 6 is a flow diagram of the clinical study;

[0070] FIG. 7 illustrates the architecture of an exemplary second machine learning model for processing image data;

[0071] FIG. 8 illustrates the architecture of the 4Dsurvival network

[0072] FIG. 9 illustrates automated segmentation of the left and right ventricles in a patient with left ventricular failure; and

[0073] FIG. 10 shows a three-dimensional model of the left and right ventricles of a patient with left ventricular failure.

DETAILED DESCRIPTION OF CERTAIN EMBODIMENTS

[0074] In the following, like parts are denoted by like reference numbers. The interpretation of dynamic biological systems requires accurate and precise motion tracking, as well as efficient representations of high-dimensional motion trajectories in order to enable use for prediction and/or risk classification tasks. Such motion information may be important in biological systems which exhibit complex spatio-temporal behaviour in response to stimuli or as a consequence of disease processes. In the present specification, methods are described which provide a generalisable approach for modelling time-to-event outcomes and/or event risk classification from time-resolved three-dimensional model data.

[0075] The present specification is concerned with the task of predicting, for a particular subject (also referred to as a patient), a time-to-event for an adverse cardiac event, and/or a measure of risk for an adverse cardiac event. The general methods described in this specification have also been assessed in a clinical study described herein.

[0076] The motion dynamics of the beating heart are a complex rhythmic pattern of non-linear trajectories regulated by molecular, electrical and biophysical processes. Heart failure is a disturbance of this coordinated activity characterised by adaptations in cardiac geometry and motion that often leads to impaired organ perfusion.

[0077] A major challenge in medical image analysis has been to automatically derive quantitative and clinically-relevant information in patients with disease phenotypes such as, for example, heart failure. The present specification describes methods to solve such problems by training a machine learning model to learn latent representations of cardiac motion which are both robust against noise and also relevant for survival prediction and/or risk estimation.

[0078] Method of Training a Machine Learning Model

[0079] Referring to FIG. 1, a block diagram of a method 1 of training a machine learning model 2 is shown.

[0080] The method is used to train the machine learning model 2 to calculate output data 3 in the form of a predicted time-to-event of an adverse cardiac event, and/or a measure of risk for an adverse cardiac event. The machine learning model 2 receives as input a time-resolved three-dimensional model 4 of a heart, or a portion of a heart. An adverse cardiac event may include death from heart disease, heart failure and so forth. An adverse cardiac event may include death from any cause. The adverse cardiac event may be associated with cardiovascular disease and/or heart dysfunction. Cardiovascular disease and/or heart dysfunction may affect one or more of the left ventricle, right ventricle, left atrium, right atrium and/or myocardium. One example of cardiovascular disease is pulmonary hypertension, such as pulmonary hypertension characterised by right and/or left ventricular dysfunction. Another example of cardiovascular disease is left ventricular failure, sometimes also referred to as dilated cardiomyopathy.

[0081] The method of training utilises a training set 5. The training set 5 may be either pre-prepared or generated at the point of training, and includes training data 6.sub.1, . . . , 6.sub.n, . . . , 6.sub.N corresponding to a number, N, of distinct subjects (also referred to as patients). Each subject for whom data 6.sub.n is included in the training set 5 has had a scan performed from which a time resolved three-dimensional model 4.sub.n has been generated. Each time resolved three-dimensional model 4.sub.n may include a representation of the whole or any part of the subject's heart, such as, for example, the right ventricle, left ventricle, right atrium, left atrium, myocardium, and so forth. Each time resolved three-dimensional model 4.sub.n may be generated from a sequence of images obtained at different time points, or different points within a cycle of the heart of the n.sup.th of N subjects. Each time resolved three-dimensional model 4.sub.n may be generated from a sequence of gated images of the subject's heart. A gated image may be built up across a number of heartbeat cycles of the subject's heart, by capturing data from the same relative time within numerous successive heartbeat cycles. For example, gated imaging may be synchronised to electro-cardiogram measurements. Each time-resolved three-dimensional model 4.sub.n may span at least one heartbeat cycle of the corresponding subject.

[0082] The time resolved three-dimensional models 4.sub.1, . . . , 4.sub.n, . . . , 4.sub.N included in the training set 5 may include or be derived from magnetic resonance (MR) imaging data. MR imaging data is typically acquired by means of gated imaging. Additionally or alternatively, some or all of the time resolved three-dimensional models 4.sub.n, . . . , 4.sub.n, . . . , 4.sub.N included in the training set 5 may include or be derived from ultrasound data. Although ultrasound data may typically have relatively lower resolution compared to MR imaging data, ultrasound data is easier and quicker to obtain, and the required equipment is significantly less expensive and more portable than an MR imaging scanner. In general, the time resolved three-dimensional models 4.sub.1, . . . , 4.sub.n, . . . , 4.sub.N included in the training set 5 may be derived from a single type of image data 23 (FIG. 2) or from a variety of types of image data 23 (FIG. 2). The machine learning methods 1, 22 of the present specification are based on latent representations 12.sub.n of cardiac motion which are robust against noise, and consequently the machine learning methods 1, 22 merely require that it is possible to acquire the necessary data to produce the time resolved three-dimensional models 4.sub.n, . . . , 4.sub.n, . . . , 4.sub.N used as input. The training data 6.sub.n for the n.sup.th of N subjects also includes corresponding outcome data 7.sub.n for that subject. Outcome data 7.sub.n may indicate the timing and nature of any adverse cardiac events associated with the subject, and hence also associated with the corresponding time-resolved three-dimensional model 4.sub.n. Outcome data 7.sub.n is obtained from long term follow-up of subjects following the scan from which the data for the time-resolved three-dimensional model 4.sub.n is obtained. The follow-up period may be as short as a few months, or may be up to several decades, depending on the subject.

[0083] According to the method 1, the machine learning model 2 is trained to recognise latent representations 12.sub.1, . . . , 12.sub.n, . . . , 12.sub.N of cardiac motion which are predictive of either the time to an adverse cardiac event and/or the risks of an adverse cardiac event. Once trained, the machine learning model 2 may be used to encode a latent representation 12 for a new subject, and use the latent representation 12 to calculate output data 3 in the form of a predicted time-to-event of an adverse cardiac event, and/or a measure of risk for an adverse cardiac event.

[0084] Once the machine learning model 2 has been trained, for example once the predictive accuracy of the machine learning model 2 when applied to a validation set (not shown) shows no further improvement, the trained machine learning model 2 is stored. For example, when the trained machine learning model 2 (FIG. 2) takes the form of a neural network, the trained machine learning model 2 may be stored by recording the weights of each interconnection between a pair of nodes. In some examples, the numbers of nodes and the connectivity of each node may be varied. In such examples, storing the trained machine learning model 2 may also include storing the number and connectivity of nodes forming one or more layers of the trained machine learning model 2. The validation set (not shown) is structurally identical to the training set 5, except that the time resolved three-dimensional models 4 and outcome data 7 included in the validation set (not shown) correspond to subjects who are not included in the training set 5. The sampling of subjects to form the training set 5 and the validation set (not shown) should be performed at random from the pool of available subjects.

[0085] In some examples, a validation set need not be used. This may be the case when the pool of potential subject is small. When a validation set is not used or not available, the predictive accuracy of the machine learning model 2 may be confirmed using a bootstrap internal validation procedure described hereinafter in relation to a clinical study.

[0086] Structure of the Machine Learning Model

[0087] The machine learning model 2 includes an input layer 9 and an output layer 10. The input layer 9 receives a time-resolved three-dimensional model 4.sub.n. Each time-resolved three-dimensional model 4.sub.n takes the form of a plurality of vertices N.sub.v. The v.sup.th of N.sub.v vertices takes the form of a three-dimensional coordinate, for example, (x.sub.v, y.sub.v, z.sub.v) in Cartesian coordinates. The vertices are mapped to features of the subject's heart to ensure that the same vertex corresponds to the same portion of the subject's heart at each time of the time-resolved three-dimensional model 4.sub.n.The time-resolved three-dimensional models may all have an equal number of vertices (x.sub.v, y.sub.v, z.sub.v). The time-resolved three-dimensional models may also include connectivity data defining which vertices are connected to which other vertices to define faces used for rendering the time-resolved three-dimensional model 4.sub.n.Although some examples of the machine learning model 2 may additionally make use of such connectivity data, this is not required.

[0088] The N.sub.v vertices of the time-resolved three-dimensional models 4.sub.1, . . . , 4.sub.n, . . . , 4.sub.N may be co-registered. In other words, there may be a spatial correspondence between the position of the N.sub.v vertices in each of the time-resolved three-dimensional model 4.sub.1, . . . , 4.sub.n, . . . , 4.sub.N. The mapping of vertices to features of subject's hearts may be used to provide such co-registration of vertex locations across different subjects.

[0089] The vertex positions (x.sub.v, y.sub.v, z.sub.v) are functions of time, i.e. x.sub.v(t.sub.o+(k−1)δt), y.sub.v(t.sub.o+(k−1)δt), z.sub.v(t.sub.o+(k−1)δt), in which t.sub.o is an initial time within the heartbeat cycle, for example t.sub.o=0, and δt is the interval between sampling times for the image sequence used to generate the time-resolved three-dimensional model 4.sub.n. A more concise notation for the vertex coordinates is used hereinafter, wherein x.sub.vk=x.sub.v(t.sub.o+(k−1)δt), y.sub.vk=y.sub.v(t.sub.o+(k−1)δt) and z.sub.vk=z.sub.v(t.sub.o+(k−1)δt). Although explained with reference to Cartesian coordinates for convenience, any suitable three-dimensional coordinate system may be used. The total number of sampling times (or gated times) may be denoted N.sub.t so that 1≤k≤N.sub.t.

[0090] Each time-resolved three-dimensional model 4.sub.n may be input to the machine learning model 2 as an input vector x which includes, for each vertex (x.sub.vk, y.sub.vk, z.sub.vk), the relative displacement of the vertex (x.sub.vk, y.sub.vk, z.sub.vk) at each time point after an initial time point. For each vertex of a given time-resolved three-dimensional model 4.sub.n, the relative displacements for the input vector x may be calculated with respect to an initial coordinate (x.sub.v1, y.sub.v1, z.sub.v1) of the vertex (x.sub.vk, y.sub.vk, z.sub.vk). For example, the input vector x may be formulated as:

x=(x.sub.vk−x.sub.v1,y.sub.vk−y.sub.v1, z.sub.vk−z.sub.v1) for all 1≤v≤N.sub.v,2≤k≤N.sub.t   (1)

[0091] Each time-resolved three-dimensional model 4.sub.n is separately converted to a corresponding input vector x.sub.n, and the time-resolved three-dimensional models 4.sub.1, . . . , 4.sub.n, . . . , 4.sub.N are processed one at a time or in batches, i.e. sequentially and not in parallel. The input layer 9 includes a number of nodes equal to the length (number of entries) of the input vectors x.sub.n, and each input vector x.sub.n in a given training set 5 is of equal length.

[0092] The machine learning model may include an encoding layer 11 which encodes a latent representation 12 of cardiac motion. In other words, the machine learning model 2 takes an input vector x.sub.n corresponding to the n.sup.th of N subjects and converts it into the latent representation 12.sub.n, which may be encoded in the values of the encoding layer 11. Each latent representation 12.sub.n is a dimensionally reduced representation of the same information as the input vector x.sub.n. Thus, the number of nodes, or dimensionality d.sub.h, of the encoding layer 11 is less than, preferably significantly less than, the number of nodes, or dimensionality d.sub.in, of the input layer 9 (equal to the length of x.sub.n). In some examples, the dimensionality d.sub.h of the encoding layer 11 may be a hyperparameter of the machine learning model 2, which may be optimised during the method 1 of training the machine learning model 2. The conversion of the input vector x.sub.n into the latent representation 12 may be performed by one or more encoding hidden layers 13 of the machine learning model 2, connected in order of decreasing dimensionality d (number of nodes) between the input layer 9 and the encoding layer 11.

[0093] The machine learning model 2 may be configured so that an output 3.sub.n in the form of a predicted time-to-event of an adverse cardiac event, or a measure of risk for an adverse cardiac event, is determined using a prediction branch 14 which receives as input the latent representation 12 of cardiac motion encoded by the encoding layer 11. The prediction branch 14 may be based on a Cox proportional hazards model, or any other suitable predictive model for adverse cardiac events. The output 3.sub.n in the form of a predicted time-to-event of an adverse cardiac event, or a measure of risk for an adverse cardiac event, is provided at one or more nodes of the output layer 10.

[0094] Additionally, the output layer 10 also provides a reconstructed model 15n of the cardiac motion, which is generated based on the latent representation 12.sub.n, for example as encoded by an encoding layer ii. The reconstructed model 15.sub.n may be determined from the latent representation 12.sub.n by one or more decoding hidden layers 16. The decoding hidden layers 16 may be symmetric with the encoding hidden layers 13, in terms of dimensionality d and connectivity.

[0095] In one example, the machine learning model 2 may include hidden layers 13, 16 and an encoding layer 11 which form a de-noising autoencoder. Such a de-noising auto-encoder may be symmetric about the central, encoding layer 11. When the machine learning model 2 includes a de-noising autoencoder, the input layer 9 and/or one or more encoding hidden layers 13 may implement a mask configured to apply stochastic noise to the inputs. For example, the input layer 9 and/or one or more encoding hidden layers 13 may be configured to set a predetermined fraction, f, of entries (i.e. inputs to the machine learning model 2) of each input vector x.sub.n to zero, the specific entries being selected at random. Herein, the term random encompasses pseudo-random numbers and processes. The predetermined fraction f may be a hyperparameter of the machine learning model 2 which may be optimised during the method 1 of training the machine learning model 2.

[0096] Alternatively, the input layer 9 and/or one or more encoding hidden layers 13 may be configured to add a random amount of noise to a predetermined fraction, f, of entries (i.e. inputs to the machine learning model) of each input vector x.sub.n, and so forth.

[0097] Updating the Machine Learning Model

[0098] Each time-resolved three-dimensional model 4.sub.n in the training set 5 is processed in sequence, and the corresponding output data 3n and reconstructed model 15n are used as input to a loss function 16 for training the machine learning model 2. The loss function provides error(s) 17 (also referred to as discrepancies or losses) to a weight adjustment process 18.

[0099] For example, the error 17 may take the form of a hybrid loss function which is a weighted sum of [0100] 1. A first contribution in the form of a reconstruction loss 19, determined based on the input time-resolved three-dimensional model 4.sub.n and the corresponding reconstructed model 15.sub.n of cardiac motion; and [0101] 2. A second contribution in the form of a prediction loss 20, determined based on the outcome data 7.sub.n obtained by clinical follow-up of the n.sup.th subject and the corresponding output data 3.sub.n.

[0102] The reconstruction loss 19 may be determined based on differences between the input time-resolved three-dimensional model 4.sub.n and the corresponding reconstructed model 15.sub.n of cardiac motion. In some examples, the prediction loss 20 may be determined based on differences between the outcome data and the corresponding outputs of predicted time-to-event or measure of risk for an adverse cardiac event.

[0103] Training the machine learning model 2 based on a loss function 16 having contributions from a reconstruction loss 19 and also a prediction loss 20 may help to ensure that the machine learning model 2 is trained to recognise latent representations 12 which are indicative of the most important geometric/dynamic aspects of a time resolved three-dimensional model 4. Use of a hybrid loss function may help to enforce that said geometric/dynamic aspects are relevant to the prediction task of estimating output data 3 in the form of a predicted time-to-event of an adverse cardiac event, and/or a measure of risk for an adverse cardiac event.

[0104] The relative weightings of the reconstruction loss 19 and the prediction loss 20 may each be hyperparameters of the machine learning model 2 which may be optimised during the method 1 of training the machine learning model 2.

[0105] In one example, the loss function 16 used to train the machine learning model 2 may take the form of a hybrid loss function, L.sub.hybrid according to:

[00002]$\begin{array}{cc}{L}_{\mathrm{hybrid}}=\alpha L_r+\gamma L_s{L}_r=\frac{1}{N}\underset{n=1}{\overset{N}{.Math.}}{.Math.x_n-\psi \left(\phi \left(x_n\right)\right).Math.}^2{L}_s=-\underset{n=1}{\overset{N}{.Math.}}{\delta }_n\left[W^{\prime }\phi \left(x_n\right)-\log \underset{j\in R\left(t_n\right)}{.Math.}\exp \left(W^{\prime }\phi \left(x_j\right)\right)\right]& \left(2\right)\end{array}$

[0106] In which: [0107] a is a weighting coefficient of the reconstruction loss, L.sub.r, [0108] y is a weighting coefficient of the prediction loss, L.sub.s, [0109] N is sample size, in terms of the number of subjects, [0110] x.sub.n is the n.sup.th of N input vectors to the machine learning model 2, [0111] δ.sub.n is an indicator of the status of the n.sup.th of N subjects (0=Alive, 1=Dead), [0112] W′ denotes a (1×d.sub.h) vector of weights, which when multiplied by the d.sub.h-dimensional latent code 12, φ(x) yields a single scalar W′φ(x.sub.i) representing the survival prediction for the n.sup.th of N subjects, [0113] ψ(φ(x.sub.n)) is the reconstructed model 15.sub.n for the n.sup.th of N subjects, expressed in an equivalent way to the input vector n.sup.th (and having dimensionality equal to input vector x.sub.n), [0114] R(t.sub.n) represents the risk set for the n.sup.th of N subjects, i.e. subjects still alive (and thus at risk) at the time the n.sup.th of N subjects died or became censored ({j:>t.sub.j>t.sub.n}), herein censored refers to the subjects outcome being only partially known because, for example, the patient underwent surgery, and [0115] n and j are summation indices.

[0116] The weight adjustment process 18 calculates updated weights/adjustments 21 for each node of the machine learning model 2 and/or connections between the nodes, and updates the machine learning model 2. For example, the updating may utilise back-propagation of errors. The updating of the machine learning model 2 is typically performed using a learning rate to avoid over-fitting to the most recently processed time resolved three-dimensional model 4.sub.n.In accordance with common practices, training of the machine learning model 2 may take place across two or more epochs. In some example, the size of the training set 5 may be expanded using suitable data augmentation strategies.

[0117] The method 1 of training the machine learning model 2 may include optimising one or more hyperparameters selected from the group of: [0118] a predetermined fraction f of entries in the input vector x.sub.n which are randomly set to zero, or otherwise modified at random; [0119] a dimensionality d (number of nodes) of one or more hidden layers 13, 16 of the machine learning model 2; [0120] a dimensionality d.sub.h of the encoding layer 11 which encodes the latent representation 12 of cardiac motion; [0121] weights a, y of the reconstruction loss 19 and/or the prediction loss 20; [0122] a learning rate for training the machine learning model 2; and [0123] an l.sub.1 regularization penalty used for training the machine learning model 2.

[0124] Depending upon the structure of machine learning model 2, not all of these hyperparameters will be used in every example of the machine learning model 2. Some examples of the machine learning model 2 may not use any hyperparameters, or may use different hyperparameters to those listed herein. Optimising one or more hyperparameters of the machine learning model 2 may be performed using any suitable technique such as, for example, particle swarm optimisation.

[0125] Each of the time resolved three-dimensional models 4.sub.1, . . . , 4.sub.n, . . . , 4.sub.N may be generated from original image data 23 (FIG. 2) using a second machine learning model 24 (FIGS. 2, 7). The second trained machine learning model 24 (FIGS. 2, 7) may be a convolutional neural network trained to identify one or more anatomical boundaries and/or features of a subject's heart. The second machine learning model 24 (FIGS. 2, 7) may generate segmentations of image date 23 (FIG. 2) in the form of a plurality of images corresponding to one or more anatomical boundaries and/or features of the subject's heart. The second machine learning model 24 (FIGS. 2, 7) may employ image registration to track and correlate one or more anatomical features within the plurality of images. An example of second machine learning model 24 (FIGS. 2, 7) is explained hereinafter.

[0126] Once the method 1 is complete, the trained machine learning model 2, or at least the portions of the trained machine learning model 2 necessary for obtaining output data 3 from an input time resolved three-dimensional model 4, may be stored on a non-transient computer-readable storage medium (not shown). For example, when a reconstructed model 15 is not needed in use, it may be considered to store only the input layer 9, the encoding hidden layers 13, the encoding layer ii, the prediction branch 14 and the part of the output layer 10 providing output data 3. However, in practice, the entire machine learning model 2 would typically be stored for convenience and also to allow inspection of the reconstructed models 15 to enable checking that output data 3 has been derived from a sensible latent representation 12. For example, if the reconstructed model 15 does not look like a heart, then the corresponding output data 3 may be regarded as questionable.

[0127] Method of Estimating a Predicted Time-to-Event of an Adverse Cardiac Event, and/or a Measure of Risk for an Adverse Cardiac Event

[0128] Referring also to FIG. 2, a block diagram of a method 22 of using a machine learning model 2 trained according to the method 1 is shown.

[0129] The method 22 includes receiving a time-resolved three-dimensional model 4 of a heart or a portion of a heart, and providing the time-resolved three-dimensional model 4 to the trained machine learning model 2. As explained hereinbefore, the trained machine learning model 2 is configured to recognise latent representations 12 of cardiac motion which are predictive of an adverse cardiac event and/or indicative of a measure of risk for an adverse cardiac event. The method 22 also includes obtaining output data 3 from the trained machine learning model 2 in the form of a predicted time-to-event of an adverse cardiac event, and/or a measure of risk for an adverse cardiac event. The time resolved three-dimensional model 4, the trained machine learning model 2, and the output data 3 are all the same as described in relation to the method 1 of training a machine learning model 2. The trained machine learning model 2 is the product of the method 1 of training a machine learning model 2.

[0130] Although not essential, the method 22 may also include obtaining a reconstruction 15 of the input time-resolved three-dimensional model 4. Obtaining the reconstruction 15 may be useful for visualisation purposes, for example to allow inspection of the reconstructed models 15 to check that output data 3 has been derived from a sensible latent representation 12. For example, if the reconstructed model 15 does not look like a heart, then the corresponding output data 3 may be regarded as questionable. Optionally, the method 22 may also include obtaining or receiving image data 23 of a subject's heart, or a portion thereof. The image data 23 may take the form of a sequence of images corresponding to different time points throughout one or more complete cardiac cycles. In general, the image data 23 will include a number of images for each time point, for example, a stack of images for each time point, each image corresponding to a slice through a cross-section of the subject's heart which is offset from each other image. The image data 23 may be obtained using any suitable technique such as, for example, magnetic resonance imaging, ultrasound, and so forth.

[0131] The method may also include processing the image data 23 to generate segmented images, then using the segmented images to generate a corresponding time-resolved three-dimensional model 4 of the subject's heart or a portion thereof, using a second machine learning model 24. The second trained machine learning model 24 may be a convolutional neural network trained to identify one or more anatomical boundaries and/or features of a subject's heart. The second machine learning model 24 may generate segmentations of a plurality of images corresponding to one or more anatomical boundaries and/or features of the subject's heart. The second machine learning model 24 may employ image registration to track and correlate one or more anatomical features within the plurality of images. An example of second machine learning model 24 is detailed hereinafter.

[0132] Although it has been described to optionally process image data 23 using the second machine learning model 24 in order to generate a time resolved three-dimensional model 4, this is not essential. The trained machine learning model 8 may generate the output data 3 by processing any suitable time resolved three-dimensional model 4, however it is originally obtained.

[0133] Experimental Study

[0134] The methods 1, 22 of the present specification have been investigated in a clinical study, the results and methods of which shall be described and discussed hereinafter in order to provide relevant context. The clinical study relates to one exemplary implementation of the general methods 1, 22 of the present specification. Although details of the exemplary machine learning model 2 used in the clinical study, termed 4Dsurvival network, provide context and verification of the methods 1, 22, the methods 1, 22 and the appended claims should not be construed as being limited by or to any specific details of the clinical study or the 4Dsurvival network described hereinafter.

[0135] The clinical study used image data 23 corresponding to the hearts of 302 subjects (patients), acquired using cardiac magnetic resonance (MR) imaging, to create time-resolved three-dimensional models 4.sub.1, . . . , 4.sub.n, . . . , 4.sub.N, which were generated using an exemplary second machine learning model 24 in the form of a fully convolutional network trained on anatomical shape priors. The time-resolved three-dimensional models 4.sub.1, . . . , 4.sub.n, . . . , 4N so generated formed the input to an exemplary machine learning model 2 in the form of a supervised denoising autoencoder, herein referred to as the 4Dsurvival network, which took the form of a hybrid network including an autoencoder configured to learn a task-specific latent code representations 12 trained on observed outcome data 7.sub.1, . . . , 7.sub.n, . . . , 7.sub.N. In this way, the trained machine learning model 2, i.e. the trained 4D survival network, was able to generate latent representations 12 optimised for survival prediction.

[0136] In order to handle right-censored survival outcomes, the 4D survival network 2 used for the clinical study was trained using a loss function 16 based on a Cox partial likelihood loss function. The clinical study included 302 subject (patients), and the predictive accuracy (quantified by the C-index, see Equation (8)) was significantly higher (p<0.0001) for the 4D survival network 2, with C=0.59 (95% confidence interval, CI: 0.68-0.78), than for a comparison human benchmark of C=0.59 (95% C1: 0.53-0.65). The clinical study provides evidence of how the methods 1, 22 of the present specification may be used to efficiently and accurately predict human survival by estimating a time-to-event for an adverse cardiac event and/or a measure of risk for an adverse cardiac event

[0137] For the clinical study, the 302 subjects (patients) studied had been diagnosed with pulmonary hypertension (PH), characterised by right ventricular (RV) dysfunction. This group was chosen as this is a disease with high mortality where the choice of treatment depends on individual risk stratification.

[0138] The training set 5 used for the clinical study was derived from cardiac magnetic resonance (CMR), which acquires imaging of the heart in any anatomical plane for dynamic assessment of function. A separate validation set was not used. Instead, a bootstrap internal validation procedure described hereinafter was used. While conventional, explicit measurements of performance obtained from myocardial motion tracking may be used to detect early contractile dysfunction and may act as discriminators of different pathologies, one outcome of the clinical study has been to demonstrate that learned features of complex three-dimensional cardiac motion, as learned by a trained machine learning model 2 in the form of the 4Dsurvival network 2, may provide enhanced prognostic accuracy.

[0139] A major challenge for medical image analysis has been to automatically derive quantitative and clinically-relevant information in patients with disease phenotypes. The methods 1, 22 of the present specification provide one solution to such challenges.

[0140] An example of a second machine learning model 24 was used, in the form a fully convolutional network (FCN), to learn a cardiac segmentation task from manually-labelled priors. The outputs of the exemplary second machine learning model 24 were time resolved three-dimensional models 4, in the form of smooth 3D renderings of frame-wise cardiac motion. The generated time resolved three-dimensional models 4 were used as part of a training set 5 for training the 4Dsurvival network 2, which took the form of a denoising autoencoder prediction network. The 4Dsurvival network was trained to learn latent representations 12 of cardiac motion which are robust against noise, and also relevant for estimating output data 3 in the form of a predicted time-to-event of an adverse cardiac event in the form of subject death. The performance of the trained 4Dsurvival network (which is only one example of a trained machine learning model 2 according to the present specification) was also compared against a benchmark in the form of conventional human-derived volumetric indices used for survival prediction.

[0141] The 4Dsurvival network 2 included an autoencoder. Autoencoding is a dimensionality reduction technique in which an encoder (e.g. encoding hidden layers 13) takes an input (e.g. vector x representing a time resolved three-dimensional model 4) and maps it to a latent representation 12 (lower-dimensional space) which is in turn mapped back to the space of the original input (e.g. reconstructed model 15). The latter step represents an attempt to ‘reconstruct’ the input time resolved three-dimensional model 4 from the compressed (latent) representation 12, and this is done in such a way as to minimise the reconstruction loss 19, i.e. the degree of discrepancy between the input time resolved three-dimensional model 4 and the corresponding reconstructed model 15 (alternatively, between input vector x and a corresponding reconstructed output vector, denoted ψ(φ(x.sub.n)) and further described hereinafter).

[0142] The 4Dsurvival network 2 was based on a denoising autoencoder (DAE), which is a type of autoencoder which aims to extract more robust latent representations 12 by corrupting the input, for example vector x representing a time resolved three-dimensional model 4 with stochastic noise. The denoising autoencoder used in the 4Dsurvival network 2 was augmented with a prediction branch 14, in order to allow training the 4Dsurvival network 2 to learn latent representations 12 which are both reconstructive and discriminative. A loss function 16 was used in the form of a hybrid loss function having a contribution from a reconstruction loss 19 and a contribution from a prediction loss 20. The prediction loss 20 for training the exemplary machine learning model 2 was inspired by the Cox proportional hazards model. A hybrid loss function 16, L.sub.hybrid was used in order to permit optimisation of the trade-off between accuracy of the output data 3 and accuracy of the reconstructed model 15, and the balance between these aspects was calibrated during training by adjusting the relative weightings α, γ of the contributions 19, 20 to the overall loss function 16. As described hereinafter, the output data 3 from the 4Dsurvival network 2, based on latent representations 12 of cardiac motion, may be observed to predict survival more accurately than a composite measure of conventional manually-derived parameters measured on the same image data 23. To safeguard against overfitting on the training set 5, dropout and L.sub.1 regularization were used in order to yield a robust prediction model.

[0143] Baseline Characteristics

[0144] Data from all 302 subjects with incident PH were included for analysis. Objective is diagnosis was made according to haemodynamic criteria. Subjects were investigated between 2004 and 2017, and were followed-up until Nov. 27, 2017 (median 371 days). All-cause mortality was 28% (85 of 302). Table 1 summarizes characteristics of the study sample at the date of diagnosis. No subjects' data were excluded.

[0145] MR Image Processing

[0146] Automatic segmentation of the ventricles from image data 23 in the form of gated CMR images was performed for each slice position at each of 20 temporal phases producing a total of 69,820 label maps for the cohort.

[0147] Referring also to FIG. 3A, an example is shown of an automatic cardiac image segmentation of each short-axis cine image from apex (slice 1) to base (slice 9) across 20 time points.

[0148] Data were aligned to a common reference space to build a population model of cardiac motion. In each image, the right ventricular wall 25, the left ventricular wall 26, the right ventricular blood pool 27 and the left ventricular blood pool 28 may be observed to have been clearly segmented.

[0149] Image registration was used to track the motion of corresponding anatomic points. Segmented image data 23 for each subject was aligned producing a dense time resolved three-dimensional model 4 of cardiac motion, which was then used as an input for training or validating the 4Dsurvival network.

[0150] Referring also to FIG. 3B, examples of time resolved three-dimensional models 4 are shown for the freewall 29 and septum 30 of the subject's hearts, averaged across the study population. The time resolved three-dimensional models 29, 30 shown in FIG. 3B were generated by averaging vertex-wise, time-resolved displacement values (along x, y and z coordinates) across all subjects.

[0151] Trajectories of right ventricular contraction and relaxation averaged across the study population are also plotted in FIG. 3B as looped pathlines for a sub-sample of 100 points (vertices) on the heart, using a magnification factor of 4 times. The greyscale shading represents relative myocardial velocity at each phase of the cardiac cycle. The surface-shaded models 29, 30 are shown at the end-systole point of a heartbeat cycle. Such dense myocardial motion fields for each subject, for example represented in the form of an input vector x, were used as the inputs to the 4Dsurvival network.

[0152] Predictive Performance

[0153] Bootstrapped internal validation was applied to the 4Dsurvival network, and also to the benchmark conventional parameter models.

[0154] Referring also to Table 1, Patient characteristics are tabulated at baseline (date of MRI scan). The acronyms in Table 1 have the following correspondences: WHO, World Health Organization; BP, Blood pressure; LV, left ventricle; RV, right ventricle.

[0155] Referring also to FIG. 4A, Kaplan-Meier plots are shown for a conventional parameter model using a composite of manually-derived volumetric measures.

[0156] Referring also to FIG.ure 4B, Kaplan-Meier plots are shown for the 4Dsurvival network, using the time resolved three-dimensional models 4 of cardiac motion as input.

[0157] For both models, subjects were divided into a low risk group 32 and a high-risk group 31 by median risk score. Survival function estimates for each group 31, 32 (with 95% confidence intervals as error bars) are shown. For the data shown in FIGS. 4A and 4B, the Logrank test was performed to compare survival curves between the risk groups 31, 32. For the conventional parameter model: x.sup.2=5.7, p=0.0173; for the 4Dsurvival network: x.sup.2=20.7, p<0.0001).

[0158] The apparent predictive accuracy for the 4Dsurvival network was C=0.85 and the optimism-corrected value was C=0.73 (95% CI: 0.68-0.78). For the benchmark conventional parameter model, the apparent predictive accuracy was C=0.61 with the corresponding optimism-adjusted value being C=0.59 (95% CI: 0.53-0.65). The accuracy for the 4Dsurvival network was significantly higher than that of the conventional parameter model (p<0.0001). After bootstrap validation, a final model was created using the training and optimization procedure outlined hereinafter, with the Kaplan-Meier plots shown in FIGS. 4A and 4B showing the survival probability estimates over time, stratified by risk groups 31, 32 defined by each model's predictions. Further details of the methods used to validate the 4Dsurvival model are described hereinafter.

[0159] Referring also to FIG. 5A, a 2-dimensional projection is shown of latent representations 12 of cardiac motion derived and used by the 4Dsurvival network. Visualisations of right ventricular motion are also shown for two subjects with contrasting risks.

[0160] To assess the ability of the 4Dsurvival network (i.e. one example of a machine learning model 2) to learn discriminative features from the data, the encoded latent representations 12 were examined by projection to 2D space using Laplacian Eigenmaps, as shown in FIG. 5A. In FIG. 5A, each subject is represented by a point, the greyscale shade of which is based on the subject's survival time, i.e. time elapsed from baseline (date of MR imaging scan) to death (for uncensored patients), or to the most recent follow-up date (for censored patients surviving beyond .sub.7 years).

[0161] Survival time was truncated at 7 years for ease of visualization. As may be observed from FIG. 5A, the 4Dsurvival network's latent representations 12 of cardiac motion show distinct patterns of clustering according to survival time. FIG. 5A also shows visualizations of right ventricular motion for a pair of exemplar subjects at opposite ends of the risk spectrum.

[0162] The extent to which motion in various regions of the right ventricle contributed to overall survival prediction was also assessed.

[0163] Referring also to FIG. 5B, saliency maps are shown for freewall 33 and septum 34, each showing regional contributions to the survival prediction (output data 3) by right ventricular motion. The greyscale shading corresponds to absolute regression coefficients which are expressed on a log-scale. For each saliency map 33, 34, a region of relatively high saliency 35, a region of relatively low saliency 36, and a region of intermediate saliency 37 are indicated in FIG. 5 for reference.

[0164] Fitting univariate linear models to each vertex in the mesh making up a time resolved three-dimensional model 4, the association between the magnitude of cardiac motion and the 4Dsurvival network's predicted risk score was computed, yielding the saliency maps 33, 34 shown in FIG. 5B. It may be observed from the saliency maps 33, 34 that contributions from spatially distant but functionally synergistic regions of the right ventricle may influence survival of subjects suffering from pulmonary hypertension.

[0165] Methods of the Clinical Study

[0166] Referring also to FIG. 6, a flowchart of the clinical study is shown.

[0167] The clinical study was a single-centre observational study. The analysed data were collected from subjects referred to the National Pulmonary Hypertension Service at the Imperial College Healthcare NHS Trust between May 2004 and October 2017. The study was approved by the Heath Research Authority and all subjects gave written informed consent. Criteria for inclusion were a documented diagnosis of Group 4 pulmonary hypertension investigated by right heart catheterization (RHC) and non-invasive imaging. All subjects were treated in accordance with current guidelines including medical and surgical therapy as clinically indicated.

[0168] In total 302 subject had cardiac magnetic resonance imaging, and the corresponding image data 23 was used for both manual volumetric analysis to generate manual segmentations 38, and also for automated image segmentation encompassing the right ventricle 39 and the left ventricle 40, across N.sub.t=20 time points (k=1, . . . , 20). Internal validity of the predictive performance of a conventional parameter model and a deep learning motion model was assessed using a bootstrapped internal validation procedure described hereinafter.

[0169] MR Image Acquisition, Processing and Computational Image Analysis

[0170] Cardiac magnetic resonance imaging was performed on a 1.5T Achieva (Philips, Best, Netherlands), using a standard clinical protocol based on international guidelines. The specific images analysed in the clinical study were retrospectively-gated cine sequences, in the short axis plane of the subject's heart, with a reconstructed spatial resolution of 1.3×1.3×10.0 mm and a typical temporal resolution of 29 ms.

[0171] Manual volumetric analysis of the images was independently performed by accredited physicians, according to international guidelines with access to all available images for each subject and no analysis time constraint. The derived parameters included the strongest and most well-established CMR findings for prognostication reported in a disease-specific meta-analysis.

[0172] Referring also to FIG. 7, the architecture of an exemplary second machine learning model 24 used for segmenting image data 23 is illustrated.

[0173] Briefly, the exemplary second machine learning model 24 took the form of a fully convolutional neural network (CNN), which takes each stack of cine images as an input, applies a branch of convolutions, learns image features from fine to coarse levels, concatenates multi-scale features and finally predicts the segmentation and landmark location probability maps simultaneously. These maps, together with the ground truth landmark locations and label maps, are then used in a loss function which is minimised via back-propagation stochastic gradient descent. Further details of the exemplary second machine learning model 24 used for the clinical study are described hereinafter.

[0174] The exemplary second machine learning model 24 was developed as a CNN combined with image registration for shape-based biventricular segmentation of the CMR images forming the image data 23 for each subject. The pipeline method has three main components: segmentation, landmark localisation and shape registration. Firstly, a 2.5D multi-task fully convolutional network (FCN) is trained to effectively and simultaneously learn segmentation maps and landmark locations from manually labelled volumetric CMR images. Secondly, multiple high-resolution three-dimensional atlas shapes are propagated onto the network segmentation to form a smooth segmentation model. This step effectively induces a hard anatomical shape constraint and is fully automatic due to the use of predicted landmarks from the exemplary second machine learning model 24.

[0175] The problem of predicting segmentations and landmark locations was treated as a multi-task classification problem. First, the learning problem may be formulated as follows: denoting the input training dataset by S={(U.sub.n, R.sub.n, L.sub.n), n=1, . . . , N }, where N is the sample size of the training data, U.sub.n={u.sup.n.sub.m, m=1, . . . , |U.sub.n|} is the raw input CMR volume for the n.sup.th of N subjects, R.sub.n={r.sup.n.sub.m, m=1, . . . , |R.sub.n|}, r.sup.n.sub.m∈{1, . . . , N.sub.r} are the ground truth region labels for volume U.sub.n(N.sub.r=5 representing 4 regions and background), and L.sub.n={l.sup.n.sub.m, m=1, . . . , |L.sub.n|}, l.sup.n.sub.m∈{1, . . . , N.sub.l} are the labels representing ground truth landmark locations for U.sub.n (N.sub.l=7 representing 6 landmark locations and background). Note that |U.sub.n|=|R.sub.n|=|L.sub.n| stands for the total number of voxels in a CMR volume. Let W denote the set of all network layer parameters. In a supervised setting, the following objective function is minimised via standard (backpropagation) stochastic gradient descent (SGD):

L(W)=L.sub.S(W)+aL.sub.d(W)+bL.sub.L(W)+c∥W∥.sub.F.sup.2   (3)

[0176] in which a, b and c are weight coefficients balancing the four terms. L.sub.S(W) and L.sub.D(W) are the region-associated losses that enable the network to predict segmentation maps. L.sub.L(W) is the landmark-associated loss for predicting landmark locations. ∥W∥.sub.F.sup.2, known as the weight decay term, represents the Frobenius norm on the weights W. This term is used to prevent the network from overfitting. The training problem is therefore to estimate the parameters W associated with all the convolutional layers. By minimising Equation (3), the exemplary second machine learning model 24 is able to simultaneously predict segmentation maps and landmark locations. The definitions of the loss functions L.sub.S(W), L.sub.D(W) and L.sub.L(W), used for predicting landmarks and segmentation labels, have been described previously, see Duan, J. et al. “Automatic 3D bi-ventricular segmentation of cardiac images by a shape-constrained multi-task deep learning approach.” ArXiv 1808.08578 (2018).

[0177] The FCN segmentations are used to perform a non-rigid registration using cardiac atlases built from >1000 high resolution images, allowing shape constraints to be inferred. This approach produces accurate, high-resolution and anatomically smooth segmentation results from input images with low through-slice resolution thus preserving clinically-important global anatomical features. Motion tracking was performed for each subject using a four-dimensional spatio-temporal B-spline image registration method with a sparseness regularisation term. The motion field estimate is represented by a displacement vector at each voxel and at each time frame k=1, . . . , 20. Temporal normalisation was performed before motion estimation to ensure consistency across the cardiac cycle.

[0178] Spatial normalisation of each subject's data was achieved by registering the motion fields to a template space. A template image was built by registering the high-resolution atlases at the end-diastolic frame and then computing an average intensity image. In addition, the corresponding ground-truth segmentations for these high-resolution images were averaged to form a segmentation of the template image. A template surface mesh was then reconstructed from its segmentation using a three-dimensional surface reconstruction algorithm. The motion field estimate lies within the reference space of each subject, and so to enable inter-subject comparison all the segmentations were aligned to this template space by non-rigid B-spline image registration. The template mesh was then warped using the resulting non-rigid deformation and mapped back to the template space. Twenty surface meshes, one for each temporal frame, were subsequently generated by applying the estimated motion fields to the warped template mesh accordingly. Consequently, the surface mesh of each subject at each frame contained the same number of vertices (18, 028), which maintained their anatomical correspondence across temporal frames, and across subjects (FIG. 7).

[0179] Characterization of Right Ventricular Motion

[0180] The time-resolved three-dimensional models 4 generated as described in the previous section were used to produce a relevant representation of cardiac motion—in this example of right-side heart failure limited to the RV. For this purpose, a sparser version of the meshes was utilized (down-sampled by a factor of ˜90) with 202 vertices. Anatomical correspondence was preserved in this process by utilizing the same vertices across all meshes.

[0181] This approach was used to produce a simple numerical representation of the trajectory of each vertex, i.e. the path each vertex traces through space during a cardiac cycle (FIG. 3B). The vertex positions (x.sub.v, y.sub.v, z.sub.v) are functions of time, i.e. x.sub.v(t.sub.0+(k−1)δt), y.sub.v(t.sub.0+(k−1)δt), z.sub.v(t.sub.0+(k−1)δt), in which t.sub.0 is an initial time within the heartbeat cycle, for example t.sub.0=0, and δt is the interval between sampling times for the image sequence used to generate the time-resolved three-dimensional model 4.sub.n. A more concise notation for the vertex coordinates is used hereinafter wherein x.sub.vk=x.sub.v(t.sub.0+(k−1)δt), y.sub.vk=y.sub.v(t.sub.0+(k−1)δt) and z.sub.vk=z.sub.v(t.sub.0+(k−1)δt). The total number of sampling times may be denoted N.sub.t so that 1≤k≤N.sub.t. For the clinical study, N.sub.v=202 and N.sub.t=20. The input vectors x are formulated according to Equation (1):

x=(x.sub.vk−x.sub.v1,y.sub.vk−y.sub.v1, z.sub.vk−z.sub.v1) for all 1≤v≤N.sub.v,2≤k≤N.sub.t   (1)

[0182] For the data in the clinical study, input vector x has length 11,514 (3×19×202), and was used as input to the 4Dsurvival network.

[0183] 4Dsurvival Network Design and Training

[0184] Referring also to FIG. 8, the architecture of the 4Dsurvival network is shown (i.e. one example of a machine learning model 2).

[0185] The 4Dsurvival network includes a denoising autoencoder that takes time-resolved three-dimensional models 4 of cardiac motion meshes as its input. The time-resolved is three-dimensional models 4 include representations of the right ventricle 39 and the left ventricle 40. For the sake of simplicity two hidden layers 13, 16, one immediately preceding and the other immediately following the central encoding layer 11, are not shown in FIG. 8. The autoencoder learns a task-specific latent code representation trained on observed outcome data 7, yielding a latent representation 12 optimised for survival prediction that is robust to noise. The actual number of latent factors is treated as an optimisable parameter.

[0186] The 4Dsurvival network provides an architecture capable of learning a low-dimensional latent representation 12 of right ventricular motion that robustly captures prognostic features indicative of poor survival. The hybrid design of the 4Dsurvival network combines a denoising autoencoder with an example of a prediction branch 14 which is based on a Cox proportional hazards model (described hereinafter). Again denoting the input vector by x∈.sup.d.sup.p, where d.sub.p=11,514, is the input dimensionality.

[0187] The 4Dsurvival network is based on a denoising autoencoder (DAE), an autoencoder variant which learns features robust to noise. The input vector x feeds directly into the encoder 41, the first layer of which is a stochastic masking filter that produces a corrupted version of x. The masking is implemented using random dropout, i.e. a predetermined fraction f of the elements of input vector x were set to zero (the value of f is treated as an optimizable parameter of the 4Dsurvival network). The corrupted input from the masking filter is then fed into a hidden layer 13, the output of which is in turn fed into a central, encoding layer ii. This central, encoding layer 11 represents the latent code, i.e. the encoded/compressed latent representation 12 of the input vector x.

[0188] This central encoding layer 11 is sometimes also referred to as the ‘code’, or ‘bottleneck’ layer. Therefore the encoder 41 may be considered as a function φ(.) mapping the input vector x∈.sup.d.sup.p to a latent code φ(x)∈.sup.d.sup.h, where d.sub.hd.sub.p (for notational convenience we consider the corruption, or dropout, step as part of the encoder 41). This produces a compressed latent representation 12 having a dimensionality which is lower than that of the input vector x (an undercomplete representation). Note that the number of units in the encoder's hidden layer 13, and the dimensionality d.sub.h of the latent code are not predetermined but, rather, treated as optimisable parameters of the 4Dsurvival network.

[0189] The latent representation 12, φ(x) is then fed into the second component of the denoising autoencoder, a multilayer decoder network 42 that upsamples the code back to the original input dimension d.sub.p. Like the encoder 41, the decoder 42 has one intermediate hidden layer 16 that feeds into the final, output layer 10, which in turn outputs a decoded representation (with dimension d.sub.p matching that of the input). In the 4Dsurvival network, this decoded representation corresponds to the reconstructed model 15.

[0190] The size of the decoder's 42 intermediate hidden layer 16 is constrained to match that of the encoder 41 networks hidden layer 13, to give the autoencoder a symmetric architecture. Dissimilarity between the original (uncorrupted) input vector x and the decoder's reconstructed model 15 (denoted here by ψ(φ(x))) is penalized by minimizing a loss function of general form L(x, ψ(φ(x))). Herein, a simple mean squared error form is chosen for L:

[00003]$\begin{array}{cc}{L}_r=\frac{1}{N}\underset{n=1}{\overset{N}{.Math.}}{.Math.x_n-\psi \left(\phi \left(x_n\right)\right).Math.}^2& \left(4\right)\end{array}$

[0191] in which N again represents the sample size in terms of the number of subjects. Minimizing this loss reconstruction loss 19, L.sub.r forces the autoencoder 41, 42 to reconstruct the input x from a corrupted/incomplete version, thereby facilitating the generation of a latent representation 12 with robust features.

[0192] As explained hereinbefore, in order to ensure that learned latent representations 12 are actually relevant for estimating output data 3, in this instance in the form of a survival prediction, the autoencoder 41, 42 of the 4Dsurvival network was augmented by adding a prediction branch 14. The latent representation 12 learned by the encoder 41, φ(x) is therefore linked to a linear predictor of survival (see Equation (5)), in addition to the decoder 42. This encourages the latent representation 12, φ(x) to contain features which are simultaneously robust to noisy input and salient for survival prediction. The prediction branch 14 of the 4Dsurvival network is trained with observed outcome data 7, in this instance survival/follow-up time. For each subject, this is time elapsed from MRI acquisition until death (all-cause mortality), or if the subject is still alive, the last date of follow-up. Also, patients receiving surgical interventions were censored at the date of surgery. This type of outcome is called a right-censored time-to-event outcome, and is typically handled using survival analysis techniques, the most popular of which is Cox's proportional hazards regression model:

[00004]$\begin{array}{cc}\log \frac{h_n\left(t\right)}{h_o\left(t\right)}={\beta }_1{z}_{n1}+{\beta }_2{z}_{n2}+\mathrm{.Math.}+{\beta }_p{z}_{np}& \left(5\right)\end{array}$

[0193] in which h.sub.n(t) represents the hazard function for subject n, i.e the ‘chance’ (normalized probability) of subject n dying at time t. The term h.sub.0(t) is a baseline hazard level to which all subject-specific hazards h.sub.n(t) (n=1, . . . , N) are compared. The key assumption of the Cox survival model is that the hazard ratio h.sub.n(t)/h.sub.0(t) is constant with respect to time (which is termed the proportional hazards assumption). The natural logarithm of this ratio is modelled as a weighted sum of a number of predictor variables (denoted here by z.sub.n1, . . . , z.sub.np), where the weights/coefficients are unknown parameters denoted by β.sub.1, . . . , ⊕.sub.p. These parameters are estimated via maximization of the Cox proportional hazards partial likelihood function:

[00005]$\begin{array}{cc}\log L\left(\beta \right)=\underset{n=1}{\overset{N}{.Math.}}{\delta }_n\left\{{\beta }^{\prime }{z}_n-\log \underset{jϵR\left(t_n\right)}{.Math.}\exp \left({\beta }^{\prime }{z}_j\right)\right\}& \left(6\right)\end{array}$

[0194] in which, z.sub.n is the vector of predictor/explanatory variables for subject n, δ.sub.n is an indicator of subject n's status (0=Alive, 1=Dead) and R(t.sub.n) represents subject n's risk set, i.e. subjects still alive (and thus at risk) at the time subject n died or became censored ({j:t.sub.j>t.sub.n}). This loss function was adapted to provide the prediction loss 20 for the 4Dsurvival network architecture as follows:

[00006]$\begin{array}{cc}{L}_s=-\underset{n=1}{\overset{N}{.Math.}}{\delta }_n\left[W^{\prime }\phi \left(x_n\right)-\log \underset{j\in R\left(t_n\right)}{.Math.}\exp \left(W^{\prime }\phi \left(x_j\right)\right)\right]& \left(7\right)\end{array}$

[0195] The term W′ denotes a (1×d.sub.h) vector of weights, which when multiplied by the d.sub.h-dimensional latent code φ(x) yields a single scalar (W′φ(x.sub.i)) representing the survival prediction (specifically, natural logarithm of the hazard ratio) for subject n. Note that this makes the prediction branch 14 of the 4Dsurvival network essentially a simple linear Cox proportional hazards model, and the predicted output data 3 may be seen as an estimate of the log hazard ratio (see Equation (5)).

[0196] For the 4Dsurvival network, the prediction loss 20 (Equation (7)) is combined with the reconstruction loss 19 (Equation (5)) to form the hybrid loss function 16 of Equation (2), reproduced for convenience:

[00007]$\begin{array}{cc}{L}_{\mathrm{hybrid}}=\alpha L_r+\gamma L_s{L}_r=\frac{1}{N}\underset{n=1}{\overset{N}{.Math.}}{.Math.x_n-\psi \left(\phi \left(x_n\right)\right).Math.}^2{L}_s=-\underset{n=1}{\overset{N}{.Math.}}{\delta }_n\left[W^{\prime }\phi \left(x_n\right)-\log \underset{j\in R\left(t_n\right)}{.Math.}\exp \left(W^{\prime }\phi \left(x_j\right)\right)\right]& \left(2\right)\end{array}$

[0197] in which the weighting coefficients α and γ are used to calibrate the contributions of each term 19, 20 to the overall loss function 16, i.e. to control the tradeoff between accuracy of the output data 3 in the form of a survival prediction versus accuracy of the reconstructed model 15. During training of the 4Dsurvival network, the weights α and γ are treated as optimisable network hyperparameters. For the clinical study, γ was chosen to equal (1−α) for convenience.

[0198] The loss function 16 was minimized via backpropagation. To avoid overfitting and to encourage sparsity in the encoded representation, we applied L.sub.1 regularization. The rectified linear unit (ReLU) activation function was used for all layers, except the prediction output layer (linear activation was used for this layer). Using the adaptive moment estimation (Adam) algorithm, the 4Dsurvival network was trained for 100 epochs with a batch size of 16 subjects. The learning rate was also treated as a hyperparameter (see Table 2). During training of the 4Dsurvival network, the random dropout (input corruption) was repeated at every backpropagation pass. The entire training process, including hyperparameter optimisation and bootstrap-based internal validation (described hereinafter) took a total of 76 hours.

[0199] Hyperparameter Tuning

[0200] To determine optimal hyperparameter values, particle swarm optimization (PSO) was used. Particle swarm optimization is a gradient-free meta-heuristic approach for finding optima of a given objective function. Inspired by the social foraging behavior of birds, particle swarm optimization is based on the principle of swarm intelligence, which refers to problem-solving ability that arises from the interactions of simple information-processing units. In the context of hyperparameter tuning, it can be used to maximize the prediction accuracy of a model with respect to a set of potential hyperparameters. Particle swarm optimization was utilised to choose the optimal set of hyperparameters from among predefined ranges of values, summarized in Table 2. The particle swarm optimization algorithm was run for 50 iterations, at each step evaluating candidate hyperparameter configurations using 6-fold cross-validation. The hyperparameters at the final iteration were chosen as the optimal set.

[0201] Model Validation and Comparison

[0202] Discrimination was evaluated using Harrell's concordance index, an extension of area under the receiver operating characteristic curve (AUC) to censored time-to-event data:

[00008]$\begin{array}{cc}C=\frac{{.Math.}_{n1,\mathrm{n2}}{\delta }_{n1}I\left({\eta }_{n1}>{\eta }_{n2}\right)I\left(t_{n1}<t_{n2}\right)}{{.Math.}_{n1,\mathrm{n2}}{\delta }_{n1}I\left(t_{n1}<t_{n2}\right)}& \left(8\right)\end{array}$

[0203] in which the indices n1 and n2 refer to pairs of subjects in the sample and I(.Math.) denotes an indicator function that evaluates to 1 if its argument is true (and o otherwise). Symbols η.sub.n1 and η.sub.n2 denote the predicted risks for subjects n1 and n2. The numerator tallies the number of subject pairs (n1, n2) where the pair member with greater predicted risk has shorter survival, representing agreement (concordance) between the model's risk predictions and ground-truth survival outcomes. Multiplication by δ.sub.n1 restricts the sum to subject pairs where it is possible to determine who died first (i.e. informative pairs). The C index therefore represents the fraction of informative pairs exhibiting concordance between predictions and outcomes. In this sense, the index has a similar interpretation to the AUC (and consequently, the same range).

[0204] Internal Validation In

[0205] order to get a sense of how well the 4Dsurvival network would generalize to an external validation cohort, its predictive accuracy was assessed within the training sample using a bootstrap-based procedure recommended in the guidelines for Transparent Reporting of a multivariable model for Individual Prognosis Or Diagnosis (TRIPOD)—see Moons, K. et al. Transparent reporting of a multivariable prediction model for Individual Prognosis Or Diagnosis (TRIPOD): Explanation and elaboration. Ann Intern Med 162, W1-W73 (2015).

[0206] This procedure attempts to derive realistic, ‘optimism-adjusted’ estimates of the model's generalization accuracy using the training sample.

[0207] (Step 1) A prediction model was developed on the full training sample (size N), utilizing the hyperparameter search procedure discussed above to determine the best set of hyperparameters. Using the optimal hyperparameters, a final model was trained on the full sample. Then the Harrell's concordance index (C) of this model was computed on the full sample, yielding the apparent accuracy, i.e. the inflated accuracy obtained when a model is tested on the same sample on which it was trained/optimized.

[0208] (Step 2) A bootstrap sample was generated by carrying out N random selections (with replacement) from the full sample. On this bootstrap sample, a model was developed (applying exactly the same training and hyperparameter search procedure used in Step 1) and computed C for the bootstrap sample (henceforth referred to as bootstrap performance). Then the performance of this bootstrap-derived model on the original data (the full training sample) was also computed (henceforth referred to as test performance)

[0209] (Step 3) For each bootstrap sample, the optimism was computed as the difference between the bootstrap performance and the test performance.

[0210] (Step 4) Steps 2 to 3 were repeated B times (where B=100).

[0211] (Step 5) The optimism estimates derived from Steps 2 to 4 were averaged across the B=100 bootstrap samples and the resulting quantity was subtracted from the apparent predictive accuracy from Step 1.

[0212] This procedure yields an optimism-corrected estimate of the model's concordance index:

[00009]$\begin{array}{cc}{C}_{cor\mathrm{rected}}=C_{full}^{full}-\frac{1}{B}\underset{b=1}{\overset{B}{.Math.}}\left(C_b^b-C_b^{full}\right)& \left(9\right)\end{array}$ [0213] Above, symbol C.sub.s1.sup.s1 refers to the concordance index of a model trained on sample s.sub.1 and tested on sample s.sub.2. The first term refers to the apparent predictive accuracy, i.e. the (inflated) concordance index obtained when a model trained on the full sample is then tested on the same sample. The second term is the average optimism (difference between bootstrap performance and test performance) over the B=100 bootstrap samples. It has been demonstrated that this sample-based average is a nearly unbiased estimate of the expected value of the optimism that would be observed in external validation. Subtraction of this optimism estimate from the apparent predictive accuracy gives the optimism-corrected predictive accuracy.

[0214] Conventional Parameter Model

[0215] As a benchmark comparison to the 4Dsurvival motion model, a Cox proportional hazards model was trained using conventional right ventricular (RV) volumetric indices including right ventricular end-diastolic volume (RVEDV), right ventricular end-systolic volume (RVESV) and the difference between these measures expressed as a percentage of RVEDV, right ventricular ejection fraction (RVEF) as survival predictors. To account for collinearity among these predictor variables, an L.sub.2-norm regularization term was added to the Cox partial likelihood function:

[00010]$\begin{array}{cc}\log L\left(\beta \right)=\underset{n=1}{\overset{N}{.Math.}}{\delta }_n\left\{{\beta }^{\prime }{z}_n-\log \underset{jϵR\left(t_n\right)}{.Math.}\exp \left({\beta }^{\prime }{z}_j\right)\right\}+\frac{1}{2}\lambda {.Math.\beta .Math.}^2& \left(10\right)\end{array}$

[0216] in which λ is a parameter that controls the strength of the penalty term. The optimal value of λ was selected via cross-validation.

[0217] Interpretation of the 4Dsurvival Model

[0218] To facilitate interpretation of the 4Dsurvival network, Laplacian Eigenmaps were used to project the learned latent representations 12 into two dimensions (FIG. 5A), allowing latent space visualization. Neural networks derive predictions through multiple layers of nonlinear transformations on the input data. This complex architecture does not lend itself to straightforward assessment of the relative importance of individual input features. In order to analyse this, a simple regression-based inferential mechanism was used to evaluate the contribution of motion in various regions of the RV to the model's predicted risk (FIG. 5B). For each of the 202 vertices in the time resolved three-dimensional models 4 used in the clinical study, a single summary measure of motion was computed by averaging the displacement magnitudes across 20 frames. This yielded one mean displacement value per vertex. This process was repeated across all subjects. Then the predicted risk scores were regressed onto these vertex-wise mean displacement magnitude measures using a mass univariate approach, i.e. for each vertex v (v=1, . . . , 202), a linear regression model was fitted where the dependent variable was predicted risk score, and the independent variable was average displacement magnitude of vertex v.

[0219] Each of these 202 univariate regression models was fitted on all subjects and yielded one regression coefficient representing the effect of motion at a vertex on predicted risk. The absolute values of these coefficients, across all vertices, were then mapped onto a template RV mesh to provide a visualization (FIG. 5B) of the differential contribution of various anatomical regions to predicted risk.

[0220] From the results of the clinical study, it may be observed that the generalised methods of the present specification permit learning of meaningful latent representations 12 of cardiac motion, which encode information useful for estimating output data 3 in the form of a predicted time to event for an adverse cardiac event and/or an estimate of risk for an adverse cardiac event.

[0221] Modifications

[0222] It will be appreciated that many modifications may be made to the embodiments hereinbefore described. Such modifications may involve equivalent and other features which are already known in the design, training and application of machine-learning methods for image processing, and which may be used instead of or in addition to features already described herein. Features of one embodiment may be replaced or supplemented by features of another embodiment.

[0223] Although the clinical study presented hereinbefore related to a particular type of heart failure, the methods of the present specification are equally applicable to similar analysis of any other heart condition and/or irregularity. This is expected to be the case because any heart condition will, inherently, have an effect on cardiac motion, and the methods of the present specification have been demonstrated, through the clinical study, to be capable of learning robust and meaningful latent representations of cardiac motion.

[0224] The same methods described hereinbefore may be applied to groups of patients experiencing different type of cardiac dysfunction. For example, the methods of the present specification may be applied a training set 5 corresponding to patients with left ventricular failure (also known as dilated cardiomyopathy).

[0225] Referring also to FIG. 9, automated segmentation of the left and right ventricles in a patient with left ventricular failure is shown. Referring again to FIG. 3A, further examples of segmenting the left ventricular wall 26 and left ventricular blood pool 28 may be seen (though the data of FIG. 3A relates to patients with pulmonary hypertension rather than left ventricular failure as shown in FIG. 9). The segmented images may be used to create a time-resolved three-dimensional model 4.

[0226] Referring also to FIG. 10, a three-dimensional model of the left and right ventricles describing cardiac motion trajectory is shown for a patient with left ventricular failure.

[0227] Such a time-resolved three-dimensional model may be used as input for training a machine learning model, for example the 4Dsurvival network described hereinbefore. The input to the machine learning model 2 may take the form of the time-resolved three-dimensional model 4, or time-resolved trajectories of three-dimensional contraction and relaxation extracted therefrom. The loss function used to the train the machine learning model 2, for example including a reconstruction loss 19 and a prediction loss 20, may be the same as described hereinbefore. Once a trained machine learning model 2 has been obtained, this may be used as described hereinbefore to obtain predictions of outcomes for patients with left ventricular failure (or any other type of cardiac dysfunction).

[0228] Although claims have been formulated in this application to particular combinations of features, it should be understood that the scope of the disclosure of the present invention also includes any novel features or any novel combination of features disclosed herein either explicitly or implicitly or any generalization thereof, whether or not it relates to the same invention as presently claimed in any claim and whether or not it mitigates any or all of the same technical problems as does the present invention. The applicant hereby gives notice that new claims may be formulated to such features and/or combinations of such features during the prosecution of the present application or of any further application derived therefrom.