Title: COUNTDOWN REGRESSION: SHARP AND CALIBRATED SURVIVAL PREDICTIONS

Abstract: Personalized probabilistic forecasts of time to event (such as mortality) can be crucial in decision making, especially in the clinical setting. Inspired by ideas from the meteorology literature, we approach this problem through the paradigm of maximizing sharpness of prediction distributions, subject to calibration. In regression problems, it has been shown that optimizing the continuous ranked probability score (CRPS) instead of maximum likelihood leads to sharper prediction distributions while maintaining calibration. We introduce the Survival-CRPS, a generalization of the CRPS to the time to event setting, and present right-censored and interval-censored variants. To holistically evaluate the quality of predicted distributions over time to event, we present the scale agnostic Survival-AUPRC evaluation metric, an analog to area under the precision-recall curve. We apply these ideas by building a recurrent neural network for mortality prediction, using an Electronic Health Record dataset covering millions of patients. We demonstrate significant benefits in models trained by the Survival-CRPS objective instead of maximum likelihood.

1. INTRODUCTION:

Having patient-specific predictions of time to an event such as mortality or bone fracture allows caregivers to make better informed decisions around patient care. Historically, prognosis scores have served as simple tools to stratify patient risk within a predefined time window (Lau et al., 2006; Cardona-Morrell & Hillman, 2015) . However, such models tend to be too simplistic to be widely useful. They are often estimated from a large population of patients, and do not take into account patient-specific information to make individualized predictions (Yu et al., 2011) . Meanwhile, the adoption of Electronic Health Record (EHR) systems over the past few decades has resulted in the collection of observational data on millions of patients spanning multiple years. This data enables development of patient-specific prediction models using machine learning. Such models are applicable to the larger patient population without being specific to a disease type or demographic, and this makes it possible to develop novel workflows in care delivery. For example, a high predicted probability of 3-12 month mortality could proactively notify palliative care teams of otherwise overlooked patients with end-of-life needs (Avati et al., 2017) .

One way to obtain patient-specific survival predictions is to treat the problem as probabilistic classification; that is, training a binary classifier to predict outcomes of event by a particular time of interest (Avati et al., 2017; Rajkomar et al., 2018) . However, such an approach has drawbacks. First, the model is specific to the time of interest it was trained upon -it is not straightforward how to take a model that was trained to predict probabilities of 1-year mortality and obtain predictions of 6-month mortality from it. Second, it is not usually possible to use data on all patients -for example, if a patient has only 3 months of history in the EHR system, it is neither possible to include that patient as a positive case nor a negative case in the 1-year mortality prediction task. Third, the process of constructing the data set implicitly conditions on the future outcome to select prediction times -evaluation is performed only at times looking backward from the event of interest. It has been shown that evaluation metrics can be overly optimistic relative to real world performance as a result (Sherman et al., 2017) .

An alternative approach to the problem is survival prediction; that is, predicting time to event by estimating a distribution over future time. In this setting, traditional survival analysis methods such Figure 1 : Example of a patient's predicted distributions for age of death under different models. Our proposed techniques improve sharpness of predicted distributions, subject to calibration. Repeated interactions (indicated by darker color) between the patient and the EHR yield more confident predictions of time of death.

as the Cox proportional hazards model (R., 1972) or accelerated failure time models (J.) are capable of handling data with censored observations (cases in which the event was not observed, but we know that the event did not occur up to a certain time). This addresses concerns raised by the classification approach, but there are a few nuances. First, traditional models typically make strong assumptions, such as proportional hazards or linearity. Second, challenges of low prevalence often arise when these methods are applied to large-scale observational datasets with heavy censoring, which is the case in real EHR data. Third, these survival analysis methods are typically evaluated as point estimates of risk, such as 10-year probabilities of events, rather than holistic measures of quality of the predicted distributions (Goff et al., 2014; Ranganath et al., 2016; Lee et al., 2018) . Common metrics of evaluation include the C-statistic (Uno et al., 2007) , log-1 loss (Yu et al., 2011) , and mean-squared-error (Katzman et al., 2018) . While useful for the purposes of relative risk stratification, model comparisons made using point estimates leaves the quality of uncertainty in predicted distributions left unmeasured. If a point prediction is way off, it is penalized by the same amount whether the model was confident or not (that is, whether the predicted distribution had low or high variance).

In contrast, forecasts in the field of meteorology are typically made as full prediction distributions over all weather conditions given past and current observations (Gneiting et al., 2008) . Evaluation of predictive performance is assessed by the paradigm of maximizing the sharpness of the predictive distribution, subject to calibration (Gneiting & Katzfuss, 2014) . The intuition behind this paradigm is that probabilities have to be calibrated in order to be correct. However, that does not necessarily make them useful (one could always predict the marginal probability of an outcome without looking at the data, and still be well calibrated). The usefulness of a prediction distribution lies in its sharpness, or how well its mass concentrates. In summary, uncalibrated predictions (sharp or not) are useless, calibrated but non-sharp predictions are correct but less useful, and calibrated and sharp distributions are most useful.

To improve the sharpness of prediction distributions in the survival setting, we propose the use of proper scoring rules beyond maximum likelihood as the training objective. Proper scoring rules are known to measure calibration, and any model trained with a proper scoring rule will tend to maintain calibration (Gneiting & Katzfuss, 2014) . For our purposes, we focus on the continuous ranked probablity score (CRPS) which has been used as an objective in the regression setting (Gneiting et al., 2008; Mohammadi et al., 2016; 2015) . We generalize the CRPS for the survival setting, called Survival-CRPS, with right-censored and interval-censored extensions. Our work is one among many recent works Chapfuwa et al. (2018) in using non-MLE training objectives for survival prediction.

Summary of contributions.:

(1) We introduce the proper scoring rule Survival-CRPS, a generalization of CRPS, as an objective in survival prediction. We present its right-censored and interval-censored variants.

(2) We propose a new metric, Survival-AUPRC, inspired by the paradigm of maximizing sharpness subject to calibration, to holistically measure the quality of a prediction distribution with respect to a possibly censored outcome. (3) We give practical recommendations for the mortality prediction task, by recommending use of the log-normal parameterization and interval censoring when training. (4) We employ the above techniques and demonstrate their efficacy by training a deep recurrent neural network model for accurate survival prediction of patient mortality using EHR data.

Parametric survival prediction methods model the time to an event of interest with a family of probability distributions, indexed by the distribution parameters. The survival function, denoted S(t) : [0, ∞) → [0, 1], is a monotonically decreasing function over the positive reals with S(0) = 1 and lim t→∞ S(t) = 0. The survival function represents the probability of an individual not having the event of interest up to a given time. Every survival function has a corresponding cumulative density function (CDF), denoted F (t) = 1 -S(t), and probability density function (PDF), denoted f (t) = d dt F (t). The choice of the family of probability distributions implies assumptions made about the nature of the data generating process.

We denote the medical record of a patient i as {(x

(i) t , a (i) t )} T (i)

t=1 , d (i) , c (i) , where t ∈ {1 . . . T (i) } denotes the interaction number of this patient with the health record, x (i) t ∈ R D is the set of features corresponding to the t-th interaction, a (i) t ∈ R + is age at time t, d (i) ∈ R + is the age of death or age of last known (alive) encounter, and c (i) ∈ {0, 1} is a censoring indicator where c (i) = 0 means the age of death is d (i) , and c (i) = 1 means the age of death is at least d (i) . For each x (i) t we define the quantity y

(i) t = d (i) -a (i)

t which represents the corresponding time to event or time to censoring. Traditional methods in survival analysis are designed to handle right-censored outcomes, but we observe that in many common scenarios outcomes are actually interval-censored. In the context of mortality prediction, for example, we know that humans almost never live past 120 years of age. Therefore, we assume that the true age of death lies between d (i) and A = 120 years, implying that the true time to death lies between 0 and

T (i) t = A -a (i)

t . We omit patient superscripts i and interaction subscripts t for succinctness where possible. We note that although our notation focuses on the problem of mortality prediction, our techniques generalize to any time to event task of interest.

2.1. SURVIVAL-CRPS: PROPER SCORING RULE OBJECTIVES:

A scoring rule is a measure of the quality of a probabilistic forecast. A forecast over a continuous outcome is a probability density function over all possible outcomes, f with corresponding cumulative density function F . In reality, we observe some actual outcome, y. A scoring rule S takes a predicted distribution and an actual outcome, and returns a loss S( F , y). It is considered a proper scoring rule if for all possible distributions G,

E y∼ F [S( F , y)] ≤ E y∼ F [S(G, y)],

and strictly proper when equality holds if and only if F = G (Gneiting et al., 2008) . A proper scoring rule is one in which the expected score is minimized by the distribution with respect to which the expectation is taken. Intuitively, it encourages a model for being honest by predicting what it actually believes (Savage, 1971) . When a proper scoring rule is employed as a loss function, it naturally forces the model to output calibrated probabilities (Gneiting & Katzfuss, 2014) .

There are many commonly used proper scoring rules. Perhaps the most widely used is the logarithmic scoring rule, equivalent to the maximum likelihood objective:

S MLE ( F , y) = -log f (y).

In the presence of possibly censored data, we maximize the density for observed outcomes, and tail or interval mass for censored outcomes, and this is a proper scoring rule (Dawid & Musio, 2014) .

S MLE-RIGHT ( F , (y, c)) = -log (1 -c) f (y) + c Ŝ(y) S MLE-INTVL ( F , (y, c, T )) = -log (1 -c) f (y) + c( F (T ) -F (y))

However, the logarithmic scoring rule is asymmetric, and harshly penalizes predictions that are wrong yet confident. This results in the training process becoming sensitive to outliers, and in general conservative in prediction-making (that is, hesitant to make sharp predictions) (Gneiting & Raftery, 2007) . Another proper scoring rule for forecasts over continuous outcomes is the CRPS (Gneiting et al., 2007) , defined as

S CRPS ( F , y) = ∞ -∞ F (z) -1{z ≥ y} 2 dz = y -∞ F (z) 2 dz + ∞ y (1 -F (z)) 2 dz.

The CRPS has been used in regression as an objective function that yields sharper predicted distributions compared to maximum likelihood, while maintaining calibration (Gneiting et al., 2008) . Intuition for the CRPS is better understood by analyzing the latter expression and noting that the two integral terms correspond to the two shaded regions in Figure 2a . The CRPS score is completely reduced to zero when the predicted distribution places all the mass on the point of true outcome, or equivalently, when the shaded region completely vanishes.

In the context of time to event predictions we propose the Survival-CRPS which accounts for the possibility of right-censored or interval-censored data:

S CRPS-RIGHT ( F , (y, c)) = y 0 F (z) 2 dz + (1 -c) ∞ y (1 -F (z)) 2 dz, S CRPS-INTVL ( F , (y, c, T )) = y 0 F (z) 2 dz + (1 -c) T y (1 -F (z)) 2 dz + ∞ T (1 -F (z)) 2 dz.

Note that when c = 0, both of the above expressions are equivalent to the original CRPS. Again, the intuition behind the Survival-CRPS is better understood by mapping each of the integral terms to the corresponding shaded region in Figure 2b and Figure 2c . The Survival-CRPS behaves like the original CRPS when the time of event is uncensored. For censored outcomes, it penalizes the predicted mass that occurs before the time of censoring and, if interval censored, also the mass after time by which the event must have occurred.

Both variants of the Survival-CRPS are proper scoring rules. They are special cases of the threshold weighted CRPS (Gneiting & Ranjan, 2011) , where the weighting function is an indicator over the uncensored regions.

2.2. EVALUATION BY SHARPNESS SUBJECT TO CALIBRATION:

Calibration assesses how well forecasted event probabilities match up to observed event probabilities. It is crucial in development of useful predictive models, especially for clinical decision-making. In binary prediction tasks without censoring, the Hosmer-Lemeshow test statistic (Hosmer et al., 2011) is commonly used to assess goodness-of-fit by comparing observed versus predicted event probabilities at quantiles of predicted probabilities. Extensions to account for censoring have been proposed (Grønnesby & Borgan, 1996; D'Agostino & Nam, 2003; Demler et al., 2015) , but these methods apply only to predictions of dichotomous outcomes within a particular time frame (for example, 1-year risks of mortality).

There is no widely accepted method for evaluating the calibration of a set of entire prediction distributions, over multiple time frames, in the survival setting. D-calibration has been recently proposed as a method for holistic evaluation (Andres et al., 2018) , but relies on handling censored observations by assuming the true times to death are uniformly distributed past the times of censoring in the predicted distributions. When censored observations far outnumber the uncensored observations, this can lead to overly optimistic assessments of calibration. Another option is to evaluate observed event times on the cumulative density scale of predicted distributions, using a Kaplan-Meier estimate to account for censoring (Harrell, 2006) . Again, this method has limitations in the heavily censored setting, as the quantiles in the tail of predicted cumulative densities have few uncensored observations, and will rarely yield well calibrated values.

We instead employ the following method to measure calibration. We compare predicted cumulative densities against observed event frequencies, evaluated at quantiles of predicted cumulative density. Right-censored observations are removed from consideration in quantiles that correspond to times after their points of censoring. Interval-censored observations are similarly removed from consideration in quantiles that correspond to times after censoring, but are additionally re-introduced in quantiles that correspond to times past the time by which the event must have occurred (in the mortality prediction task, this corresponds to 120 years of age).

Subject to calibration, we strive for prediction distributions that are sharp (i.e, concentrated). There are several metrics that could be used for measuring sharpness, such as variance or entropy. In the context of time to event predictions, holding two distributions with vastly different means to the same standard of variance or entropy would be unfair (for example, we would want lower variance for a prediction distribution with a mean of a day, compared to a mean of a year). Instead, we use the coefficient of variation (CoV) as a reasonable measure of sharpness. The CoV is defined as the ratio of one standard deviation to the mean, CoV(

F ) = √ Var[ F ] E[ F ] .

2.3. SURVIVAL-AUPRC: HOLISTIC EVALUATION OF TIME TO EVENT PREDICTIONS:

Since sharpness is only a function of the predicted distributions, a measure of sharpness is only meaningful if the model is sufficiently calibrated. We now propose a metric that measures how concentrated the mass of the prediction distribution is around the true outcome, robust to miscalibration. The idea is similar to the area under a precision-recall curve, except here it is with respect to only one predicted distribution and one outcome. We first consider the uncensored case. As an analog to precision, we consider intervals relative to the true time of event, defined by ratios. For example, a region of precision 0.9 around an event that occurs at time y is the interval [0.9y, y/0.9]. Corresponding to this region of precision, the analogy to recall is the mass assigned by the predicted distribution over this interval, F (y/0.9) -F (0.9y). By exploring the full range of precision from 0 to 1, we obtain the Survival Precision Recall Curve. The area under this curve measures how quickly predicted mass concentrates around the true outcome as we expand the precision window.

Survival-AUPRC:

UNCENS ( F , y) = 1 0 ( F (y/t) -F (yt))dt

The highest possible score is 1, when the predicted distribution is a Dirac δ function centered over the time of outcome. The lowest possible score is 0, when the predicted distribution is infinitely dispersed. The mean of all Survival-AUPRC scores across examples provides an overall measure of the quality of the predictions.

The aforementioned metric only applies when the event outcome is uncensored. In the case of censored observations, we use the same analogy but with the right end of precision intervals defined with respect to the time by which the event must have occurred in the interval-censored case, or infinity in the right-censored case.

RIGHT ( F , y) = 1 0 (1 -F (yt))dt Survival-AUPRC INTVL ( F , y, T ) = 1 0 ( F (T /t) -F (yt))dt

2.4. RECURRENT NEURAL NETWORK MODEL:

We apply our techniques to the mortality prediction task by building a multilayer recurrent neural network (RNN) with parameters θ, denoted RNN θ , that takes as input a sequence of features (in our case, information about a patient recorded in the EHR, for each interaction they had with the hospital) to predict parameters of a parametric probability distribution F over time to death at each timestep (Figure 3 ). The network depends only on data from the current and previous timesteps, and not the future. The approach here is similar to the recently proposed Weibull time to event RNN (Martinsson, 2016) , though we generalize to any choice of noise distribution. The distributions that are output in each timestep are used to construct an overall loss,

L RIGHT = N i=1 T (i) t=1 S RIGHT FRNN θ x (i) 1:t , y (i) t , c (i) L INTVL = N i=1 T (i) t=1 S INTVL FRNN θ x (i) 1:t , y (i) t , c (i) , T (i) t ,

where N is the total number of patients in the training set, T (i) is the sequence length for patient i, and FRNN θ denotes the distribution parameterized by the output of the RNN. It is the sequential and monotonically decreasing predicted times to event that inspires the name Countdown Regression.

2.5. CHOICE OF LOG-NORMAL NOISE DISTRIBUTION:

Common parametric distributions over time to event used in traditional survival analysis models include the Weibull, log-normal, log-logistic, and gamma (in order to be sufficiently expressive in model space, we seek distributions with at least two parameters). We choose the log-normal distribution because other distributions either involve the Gamma function in their density, or involve the pattern (y/p 1 ) p2 , where p 1 and p 2 are parameters output from the neural network. We found these patterns to be highly sensitive to the inputs and to suffer from numerical instability issues.

For the log-normal distribution, a closed form expression for the CRPS is well known (Baran & Lerch, 2015) . However, a closed form expression for the Survival-CRPS does not exist. We perform a change of variable to express the integral terms as finite integrals, and numerically approximate with the trapezoid rule. When training, we then back-propagate through the trapezoidal approximation. Details are given in Appendix B and C. We note that the approximation formulas are themselves proper scoring rules, as they are just weighted sums of Brier scores. Closed form expressions for the log-normal Survival-AUPRC are also given in Appendix D, E, and F.

3. EXPERIMENTS:

We run experiments for the mortality prediction task to evaluate four different training objectives: maximum likelihood S MLE-RIGHT and S MLE-INTVL , and our scoring based loss S CRPS-RIGHT and S CRPS-INTVL . For interval censoring we assume a maximum lifespan of A = 120 years. 

RNN:

Step 1

Step 2

Step 3 (c) Alive (interval-cens) patients.



Step 1

Step 2

Step 3 0.233 ± 2e-4 0.319 ± 3e-4 0.343 ± 4e-4 0.366 ± 4e-4 Alive: mean Surv-AUPRC (intvl-cen) 0.407 ± 6e-5 0.963 ± 2e-5 0.977 ± 3e-5 0.976 ± 3e-5 Figure 4 : Calibration plots for each of the models. We compare predicted cumulative densities against observed event frequencies, evaluated at quantiles of predicted cumulative density. Right-censored observations are removed from consideration in quantiles past times of censoring, interval-censored observations are additionally re-introduced in quantiles corresponding to times past 120 years.

The neural network architecture is kept identical for all four experiments and implemented in PyTorch (Paszke et al., 2017) . The input at each timestep consists of both real valued (for example, age of patient) and discrete valued (for example, ICD codes) data. Discrete data is embedded into a trainable real-valued 126-dimensional vector space, and vectors corresponding to the codes recorded at a given timestep are combined into a weighted mean by a soft self-attention mechanism. All real valued inputs are appended to the averaged embedding vector. We also provide the real valued features to every layer by appending them to the output of previous layer. The input vector feeds into a fully connected layer, followed by multiple recurrent layers. We use the Swish activation function (Ramachandran et al., 2017) and layer normalization (Ba et al., 2016) at every layer. Recurrent layers are defined using GRU units (Chung et al., 2014) with layer normalization inside. After the set of recurrent layers, the network has multiple branches, one per parameter of the survival distribution (for the lognormal, µ and σ 2 ). The final layer in each branch has scalar output, optionally enforced positive with the softplus function, Softplus(z) = log(1 + exp(z)). We use Bernoulli dropout (Srivastava et al., 2014) at all fully connected layers, and Variational RNN dropout (Gal & Ghahramani, 2015) in the recurrent layers, with a dropout probability of 0.5. Optimization is performed using the Adam optimizer (Kingma & Ba, 2014) , with a fixed learning rate of 1e-3.

3.1. DATA:

We use electronic health records, with IRB approval, from the STARR Data Warehouse (previously known as STRIDE) for training and evaluation (Lowe et al., 2009) . The Warehouse contains deidentified data for over 3 million patients (about 2.6% having a recorded date of death), spanning approximately 27 years. Each timestep in the sequence for a patient corresponds to all the data in the EHR for a given day. Only days having any data have a corresponding timestep in the sequence for each patient. We use diagnostic codes, medication order codes, lab test order codes, encounter type codes, and demographics (age and gender). Each code has a randomly initialized embedding vector as a trainable parameter. The set of 3 million patients, correspond to 51 million overall timesteps, and was randomly split in the ratio 8:1:1 into train, validation and test splits.

3.2. RESULTS:

We first verify that all models are reasonably well-calibrated (Figure 4 ). Both the coefficient of variation and the Survival-AUPRC metrics suggest that the Survival-CRPS with interval censoring yields the sharpest prediction distributions (Table 1 ). Inspecting the mass past 120 years of age shows that a naively trained prediction model with maximum likelihood can assign more than 75% of the mass to unreasonable regions, which is highly undesirable for the purpose of prediction. We note that this behavior is largely due to low prevalence of uncensored examples, which is typical in real world EHR data sets. As a result, the loss for the censored examples, which can be minimized by pushing mass as far away to the right as possible, dominates the small number of uncensored examples.

By predicting an entire distribution over time to death, the same model can be used to make classification predictions at various time points, highlighting the flexibility of our approach. When evaluated at 6 month, 1 year, and 5 year probabilistic predictions of mortality, our model remains well-calibrated with high discriminative ability (Appendix G, Figure 6 ).

In the interest of reproducibility, we run similar experiments on the publically available MIMIC-III dataset (Johnson et al., 2016) (Appendix H), and have published our corresponding source code1 .

4. RELATED WORK:

Recent works have demonstrated potential to significantly improve patient care by making predictions with deep learning models on EHR data (Avati et al., 2017; Rajkomar et al., 2018) , but these have been limited in treating the task as binary classification over a fixed time frame. Predicting survival curves instead of dichotomous outcomes has been explored (Yu et al., 2011; Lee et al., 2018) , but only over finite length horizons. Deep survival analysis (Ranganath et al., 2016) has been proposed, but is limited to a fixed shape Weibull (bypassing the concerns we raised about stability, but limited in expressivity). The work by Yang et al. (2017) is similar to ours in terms of using an RNN and log-normal noise distribution, but limited to MLE training. DeepSurv (Katzman et al., 2018) uses a Cox proportional hazards model, which similarly makes a set of inflexible assumptions. The WTTE-RNN (Martinsson, 2016) model has a similar network architecture to ours, but is also limited to a Weibull distribution. All aforementioned models have only been optimized for maximum likelihood, instead of more robust proper scoring rules. The CRPS scoring rule has been used with Neural Networks in (Rasp & Lerch, 2018) . Work in (Alaa & van der Schaar, 2017 ) also predicts full survival curves specific to a patient, but the use of GPs makes it difficult to scale to millions of patients. The work in Miscouridou et al. (2018) 

5. CONCLUSION:

Better survival prediction models can be built by exploring objectives beyond maximum likelihood and evaluation metrics that assess the holistic quality of predicted distributions, instead of point estimates. We introduce the Survival-CRPS objective, motivated by the fact that the CRPS scoring rule is known to yield sharp prediction distributions while maintaining calibration. There are perhaps others scoring rules that work better, leaving avenues for future work. To evaluate, we introduce the Survival-AUPRC metric, which captures the degree to which a prediction distribution concentrates around the observed time of event. We demonstrate success in large-scale survival prediction by using a deep recurrent model employing a log-normal parameterization. By predicting an entire distribution for time-to-event, we circumvent issues associated with binary classification. Meanwhile, our model still yields accurate predictions when evaluated as dichotomous outcomes at particular times. The impact of having meaningfully accurate survival models is tremendous, especially in healthcare. We hope our work will be useful to those looking to build and deploy such models. 

:

http://github.com/anonymoususer/anonymous.git

Figures and Tables:

Figure2: Graphical intuition for the Survival-CRPS scoring rule. For uncensored observations, we minimize mass before and after the observed time of event. For right-censored observations, we minimize mass before observed time of censoring. For interval-censored observations, we minimize mass before observed time of censoring, and mass after the time by which event must have occurred.

Figure3: RNN model overview. For each interaction, we attend over recorded ICD codes at that timestep and predict parameters µ, σ 2 of a log-normal distribution, minimizing a proper scoring rule.

predict survival curves (both non-parametric, and flexible flow based parametric curves) while also handling missing covariates. Another recent work Chapfuwa et al. (2018) uses adversarial training for survival prediction. It has been shown that modern neural networks can be uncalibrated, and the work by Guo et al. (2017) suggests ways to improve calibration (though the work focuses on classification).

Figure 5: Median predicted time to death (with 95% intervals) for individual patients from the interval-censored Survival-CRPS model. Our model gives more confident predictions upon repeated interactions between patients and the EHR. True times to death generally lie within predicted intervals.

N (µ,σ 2 ) ( Ũ -log t) -F N (µ,σ 2 ) ( L + log t) dt(substituting s = log t) = 0 -∞ F N (µ,σ 2 ) ( Ũ -s) -F N (µ,σ 2 ) ( L + s) e s ds = F N (µ,σ 2 ) ( Ũ -s) -F N (µ,σ 2 ) ( L + s) e s (µ,σ 2 ) ( Ũ -s) -f N (µ,σ 2 ) ( L + s) e s ds = F N (µ,σ 2 ) ( Ũ ) -F N (µ,σ 2 ) µ,σ 2 ) ( Ũ -s) + f N (µ,σ 2 ) ( L + s) e s ds = F N (µ,σ 2 ) ( Ũ ) -F N (µ,σ 2 ) µ,σ 2 ) ( Ũ -s)e s ds + 0 -∞ f N (µ,σ 2 ) ( L + s)e s ds = F N (µ,σ 2 ) ( Ũ ) -F N (µ,σ 2 ) substituting u = Ũ -s -µ σ = F N (µ,σ 2 ) ( Ũ ) -F N (µ,σ 2 ) ( L) + Ũ -µ σ ∞ 1 σ φ (u) e Ũ -σu-µ (-σ)du + (µ,σ 2 ) ( Ũ ) -F N (µ,σ 2 ) e vσ-L+µ σdv = F N (µ,σ 2 ) ( Ũ ) -F N (µ,σ 2 ) ( L) -e Ũ -µ Ũ -µ σ ∞ φ (u) e -σu du + e -L+µ x -c) = F N (µ,σ 2 ) ( Ũ ) -F N (µ,σ 2 ) F N (µ,σ 2 ) ( Ũ ) -F N (µ,σ 2 ) F N (µ,σ 2 ) ( Ũ ) -F N (µ,σ 2 ) ,σ 2 (log U ) -Φ µ,σ 2 (log L) -AUPRC FOR LOG-NORMAL (RIGHT-CENSORED)For a general continuous prediction distribution F with an interval outcome [L, ∞), we define Survival-AUPRC asSurvival-AUPRC(F, L) = 1 -F (Lt)] dt.Specifically for the case of log-normal, where Φ is the CDF of N (0, 1), andL = log L (following Appendix-D), Survival-AUPRC(F LN(µ,σ 2 ) , L) = 1 -F LN(µ,σ 2 ) (Lt) dt = 1 -Φ µ,σ 2 ( -AUPRC FOR LOG-NORMAL (UNCENSORED)For a general continuous prediction distribution F with a point outcome y, we define Survivaly/t) -F (yt)] dt.Specifically for the case of log-normal, where Φ is the CDF of N (0, 1), and ỹ = log y (following Appendix-D),Survival-AUPRC(F LN(µ,σ 2 ) , y) = 1 LN(µ,σ 2 ) (y/t) -F LN(µ,σ 2 ) (yt) dt

Metrics measuring sharpness and calibration for models trained on the right-censored and interval-censored variants of the maximum likelihood and Survival-CRPS objectives.

Table content: <table><tr><td>Metric</td><td colspan="3">MLE-RIGHT MLE-INTVL CRPS-RIGHT CRPS-INTVL</td></tr><tr><td>Calibration slope</td><td>1.125 ± 3e-4 1.139 ± 3e-4</td><td>1.003 ± 3e-4</td><td>0.959 ± 5e-4</td></tr><tr><td>Mean coefficient of variation</td><td>18.42 ± 5e-3 0.911 ± 4e-4</td><td>0.332 ± 1e-4</td><td>0.301 ± 1e-4</td></tr><tr><td>Mean prob of survival to age 120 yrs</td><td>0.754 ± 2e-5 0.045 ± 3e-5</td><td>0.015 ± 3e-5</td><td>0.005 ± 1e-6</td></tr><tr><td>Dead: mean Surv-AUPRC (uncen)</td><td/><td/><td/></tr></table>

Peer review (iclr2019):
 My main concern is that the authors fail to compare their appproach to any of the modelling approaches discussed in the related works section. In particular, as mentioned by the authors the WTTTE-RNN has a similar architecture and thus would have been a crucial baseline for comparisons. Furthermore, I would have liked to see an evaluation on more datasets, especially since the data in Appendix H indicate that the proposed approach is only marginally better than MLE-based model fitting. Finally, in addition to the metrics presented, conventional metrics such as the C-statistic would have been interesting. I further miss a discussion of alternative approaches to achieve well calibrated scores, especially posthoc calibration using the validation set as discussed in Guo et al, ICML 2017. Related work is incomplete, for example the use of tensor-trains in RNNs to model EHR data (Yang et al) - would the proposed approach not benefit for the use of such tensorization to better model the high-dimensional, sparse EHR data? references: Guao et al, On Calibration of Modern Neural Networks, ICML 2017 Yang et al, Modeling progression free survival in breast cancer with tensorized recurrent neural networks and accelerated failure time models, Machine Learning for Healthcare Conference 2017