Repeated measurements on the right hand side of the regression equation

Question

Suppose I observe an outcome $Y$ (e.g., survival binary) and I have variables $X_1, X_2, \dots, X_m$ that all refer to the same entity (e.g., health condition assessed by rater $1, 2, \dots, m$). So for each $X_i$ I have the same $Y$, but $X_1, X_2, \dots, X_m$ may be different. I thought about the following model: $$y_{ij} = \beta_0 + \beta_1 x_{ij} + u_{ij}$$ for $i = 1,2,\dots,n$ and $j = 1,2,\dots,m$. Note that $y_{ij} =y_i$ for all $j = 1,2,\dots,m$. That is, the dependent variable $y_i$ for individual $i$ is always the same. The variance-covariance matrix of the OLSE $\hat{\boldsymbol\beta} = (\hat\beta_0,\hat\beta_1)$ is given by $$\operatorname{Var}[\hat{\boldsymbol\beta}|\boldsymbol X] = (\boldsymbol X^\top\boldsymbol X)^{-1}(\boldsymbol X^\top\boldsymbol\Sigma\boldsymbol X)(\boldsymbol X^\top\boldsymbol X)^{-1},$$ where $$\boldsymbol X = \begin{pmatrix} 1 & x_{11} \\\ 1 & x_{12} \\ \vdots & \vdots \\\ 1 & x_{1m} \\\ 1 & x_{21} \\\ 1 & x_{22} \\\ \vdots & \vdots \\ 1 & x_{nm}\end{pmatrix}\quad\text{and}\quad\boldsymbol\Sigma = \operatorname{Var}\left[\begin{pmatrix} u_{11} \\ u_{12} \\ \vdots \\ u_{1m} \\ u_{21} \\ u_{22} \\ \vdots \\ u_{nm} \end{pmatrix}\middle|\boldsymbol X\right].$$ Now my questions concerns $\boldsymbol\Sigma$:

1. Is $\boldsymbol\Sigma$ a diagonal matrix, i.e., $u_{ij}$ and $u_{kl}$ are uncorrelated for $i\neq k$ and $j\neq l$?
2. Is $\boldsymbol\Sigma$ a block-digonal matrix, i.e., $u_{ij}$ and $u_{kl}$ are correlated for $i\neq k$ but there might be non-zero correlations for $i = k$.

I find 1. not very realistic because all raters rate the same individual, so there should be some dependency (driven by the, say, actual health condition of the individual. Of course, this implicitly assumes the raters share a common perspective, but I think this is a reasonable assumption). On the other hand, I find that 2. can't be true as $\operatorname{Var}[y_{ij}] = \operatorname{Var}[u_{ij}]$ in the linear model. Hence we would consider the covariance of $y_{ij}$ and $y_{il}$. But since for any $i = 1,2,\dots,n$, the outcome $y_{ij}$ is the same for all $j = 1,2,\dots,m$, we would have the covariance of a random variable with itself.

What would you think is a better model? I know that this will likely depend on the specific use-case at hand, but I would like to get a general understanding. The reason I am asking this question is whether there is a need to consider a method that corrects the standard errors (e.g., mixed-model, robust standard errors, clustered standard errors, ...).

Answer 1

What about the following joint-modeling perspective:

1. There is a latent quantity (you might think of it as the true underlying health state of the patient). Let's call that $x_i$ for subject $i$.
2. The outcome (or survival) model get specified as relating the outcome with the not directly observable true latent state.
3. There is a second component to the model that describes the observations $X_{ij}$ (by different raters indexed by $j$). E.g. this might just be that there's both a random subject effect (that allows severity vary between patients) and random rater effects (that allows different raters to have different systematic biases).

This would be a common approach to this kind of problem. E.g. Stan or if you don't want to write all of your own code the brms R package could implement such a model in a straightforward manner.

Update:

Model 1: To write out model equations in the case of random effects on both subject and rater: $$Y_i | (X_i = x_i) = \beta_0 + \beta_1 \times x_i + \epsilon_i$$ where $\epsilon_i \sim N(0, \sigma^2)$, $$X_{ij} | (X_i = x_i) = x_i + \nu_j,$$ where $\nu_j \sim N(0, \rho^2)$ (arguably you could also have a second error term $\psi_{ij}$ that is independent between observations $X_{ij}$) and $$X_i = \mu + \eta_i,$$ where $\eta_i \sim N(0, \tau^2)$.

Model 2: Alternatively, you could consider $x_i$ for each subject a fixed effect where we don't assume a distribution across subjects, in which case this all becomes: $$Y_i = \beta_0 + \beta_1 \times x_i + \epsilon_i$$ where $\epsilon_i \sim N(0, \sigma^2)$ and $$X_{ij} = x_i + \nu_j,$$ where $\nu_j \sim N(0, \rho^2)$ (again, you could introduce a second error term $\psi_{ij}$ independent between observations $X_{ij}$). Model 2 makes fewer assumptions than model 1 (about how the $X_i$ vary across subjects), but needs to estimate a separate parameter $X_i$ for each subject with not sharing of information across subjects. So in some situations parameters might be non-identifiable when in the first formulation they would have been.

Answer 2

The cluster sandwich covariance estimator works quite well in situations like this. As a test, you can just duplicate each row of data in a standard study and you’ll find the cluster covariance estimate is virtually the same as the standard covariance matrix from the original non-duplicated dataset. So arranging the data in a tall and thin format and making each patient a cluster should work fine.

One thing is lost in this layout though. It’s possible that what predicts survival is how difficult it is to rate a condition. The original short and wide layout allows you to create summary variables such as the variance across raters or the proportion of raters who deemed the condition as positive. Perhaps it’s possible to compute summary measures and smear them over the tall and thin records, but watch out for collinearities with individual ratings. You might think of doing a cluster sandwich-based Wald chunk test for the added value of individual ratings after adjusting for inter-rater summaries.

Slightly related: https://hbiostat.org/bbr/obsvar