Sensitivity / specificity with dependent observations

Question

Suppose we have a classification scenario, where each observation is comprised from a reader, a case and the status (0 or 1). Each reader get several cases (but each case is shown to only one reader). Our goal is to estimate the probability of the readers to predict the status.

That is:

$$ \mu = \beta_0 + \beta_1 X_{status} + \beta_2 Z_{reader}, $$

where $Z_{reader} \sim N(0, \sigma^2)$ is a random effect, and we assume $Y \sim Ber(logistic(\mu))$. After estimating the effect of $X_{status}$, we want to give confidence intervals for the specificity and sensitivity. Now, usually, we would employ a standard proportion interval, however the observations are dependent.

What is the standard method to overcome this issue? Is it ok to ignore it, since we are using the logistic model fit and there we account for the dependency?

Answer 1

No, you probably cannot ignore the random effect. The problem is that this is a rather non-linear model and the random effects (and the uncertainty around them) don't cancel out on average. E.g. think about an intercept of about 2.19 on the logit scale (with no uncertainty around that) and a random effect around that's estimated to be about 0, but with a standard error of 1.0. Then, the if the random effect were exactly 0, the corresponding probability is 0.9, but with the uncertainty around the random effect our estimate is 0.866 (R code: plogis(qlogis(0.9)) = 0.9 vs. mean( plogis(qlogis(0.9) + rnorm(1000000, mean=0, sd=1)) $\approx 0.866$).

There's multiple questions you could answer: 1) What is the sensitivity & specificity of reader $i$ (probably not really of interest to you?)? 2) What is the sensitivity & specificity of a reader randomly drawn from the population (you might think of this as "What would be the sensitivity and specificity of a future reader that is similar (=comes from the same population) as the ones in my study and works in a similar set-up?").

Perhaps you can think of further questions that may make sense, in which case you may need to think through what those need similar to the below.

For the first question, if you have the estimates of $\beta_0$, $\beta_1$, $\beta_2$, and $Z_i$, and the covariance matrix of them (we can probably reasonably assume a normal approximation to the sampling distribution of the estimates), you can sample from that normal distribution, calculate the sensitivity and specificity for each of those samples, and thus, get estimate, SE and CIs for those metrics. Maybe you can even do this analytically. Your main challenge is probably that your software may not even give you a SE for $\hat{Z}_i$, nor the covariance with other model parameters (I'm not 100% sure, but it might be approx. okay to assume 0 correlation, in which case you'd just need the SE).

For the second question you take a similar approach, except that you sample a new $Z_* \sim N(0, \sigma^2)$. Here you'd want the estimates of $\beta_0$, $\beta_1$, $\beta_2$ and $\sigma$, and their covariance matrix. Based on samples from those you then sample $Z_*$ and sensitivities & specificities. Again, your main challenge is getting informtion on the uncertainty around $\hat{\sigma}$ and how your estimate of $\sigma$ is correlated with the other parameter estimates.

The solution to this problem that I've used (in a slightly more complicated case) is to use a Bayesian model with weakly informative (or default "objective") priors. That way, you can quite easily get MCMC samples from the joint posterior distribution of the model parameters and do any simulations you like based on those. E.g. in R you can probably use the brms (or rstanarm) packages to do something like this (arguably you don't really need $\beta_2$, so I omit it below):

library(brms)
brm(formula = Y ~ (1|reader) + status,  
                   data=mydata, 
                   family = bernoulli(link = "logit"),
                   seed = 123)

Then you can calculate the quantities you like afterwards. E.g. the tidybayes package is nice for working with posterior samples, but it's probably easier to just ask the predict method for the brmfit object into giving you predictions for a previously unseen value for reader and for all possible categories for status.