Incorporating Measurement Error for Death Event In Survival Model

Modeling
1

Hello,

I am newer to bayesian modeling and even newer to Stan specifically. I am working on a power survival model and would like to incorporate measurement error - specifically measurement error related to the observation of the death event, not measurement error of any covariates. The specific use-case is when the death event is observed by another model (statistical or machine learning) leading to error in the recorded death event. I have started from the Weibull model on the Stan user guide:

data {
    int<lower=0> N;
    array[N] <lower=0> observed_times;
}

parameters {
    real<lower=0> alpha;
    real<lower=0> sigma;
}

model {
    alpha ~ lognormal(0,1);
    sigma ~ lognormal(0,1);
    obs_times ~ weibull(alpha, sigma);
}

I am simulating data from an inverse survival function to try to fit using this code;

functions {
    real inv_survival(real u, real alpha, real sigma) {
        return sigma * pow(-log(u), 1 / alpha)
    }
}

data {
    int<lower=1> N;
    real<lower=0> alpha;
    real<lower=0> beta;
}

generated quantities {
    array[N] real<lower=0> observed_times;
    for (n in 1:N) {
        real_u = uniform_rng(0,1);
        observed_times[n] = inv_survival(u, alpha, sigma);
    }
}

All of which works as expected. I can properly fit my simulated data and extract the parameters used. What is lacking is my ability to model any error in my data. After generating some “ground truth” survival observation times, I add some random noise to the observed times to simulate the effect of uncertainty in the event-observing model (while still restricting observed_time < 0). How can I add to this model to incorporate this error specifically? My mental model suggests that when you know generally how much measurement error is added, you should be able to recover the same parameters with wider HDI bounds. Instead, I get a totally different parameter back (when fitting alpha = 2.5, sigma = 3, adding noise returned alpha = 2.9, sigma = … 9.5). Any ideas how to work measurement error in the event observation into my model?

2

Hi, @outStanding, and welcome to the Stan Forums!

Measurement error models all tend to look similar. There’s a chapter in the User’s Guide: Measurement Error and Meta-Analysis

The basic idea is that you introduce a parameter for the real time and then treat the observed time as being observed with error. Given that everything’s positive here, let’s say you have lognormal measurement error with scale tau. You’d code that up in a case like you have as:

data {
    int<lower=0> N;
    array[N] <lower=0> obs_times;  // renamed for consistency
}

parameters {
    real<lower=0> alpha;
    real<lower=0> sigma;
    vector<lower=0>[N] actual_times;
    real<lower=0> tau;
}

model {
    alpha ~ lognormal(0,1);
    sigma ~ lognormal(0,1);
    tau ~ exponential(1);  // or whatever makes sense for scale of error prior
    obs_times ~ lognormal(actual_times, tau);
    actual_times ~ weibull(alpha, sigma);
}

If you know more about the measurement error process, you can model that. For instance, there might be covariates that affect the error, or the error might not be unbiased so you want to add a bias term or it might not be multiplicative like the lognormal model implies, etc.