Mostrando que 100 mediciones para 5 sujetos proporcionan mucha menos información que 5 mediciones para 100 sujetos

21

En una conferencia escuché la siguiente declaración:

100 mediciones para 5 sujetos proporcionan mucha menos información que 5 mediciones para 100 sujetos.

Es bastante obvio que esto es cierto, pero me preguntaba cómo se podría probar matemáticamente ... Creo que se podría usar un modelo lineal mixto. Sin embargo, no sé mucho acerca de las matemáticas utilizadas para estimarlos (solo corro lmer4para LMM y bmrsGLMM :) ¿Podría mostrarme un ejemplo donde esto sea cierto? Prefiero una respuesta con algunas fórmulas, que solo un código en R. Siéntase libre de asumir una configuración simple, como por ejemplo un modelo mixto lineal con intercepciones aleatorias y pendientes distribuidas normalmente.

PD: una respuesta basada en matemáticas que no implique LMM también estaría bien. Pensé en LMM porque me parecían la herramienta natural para explicar por qué menos medidas de más sujetos son mejores que más medidas de pocos sujetos, pero es posible que me equivoque.

DeltaIV
fuente
3
+1. Supongo que la configuración más simple sería considerar una tarea de estimar la media de la población μ donde cada sujeto tiene su propia media anorte(μ,σuna2) y cada medición de este sujeto se distribuye como Xnorte(una,σ2) . Si tomamos norte medidas de cada uno de los metro sujetos, entonces ¿cuál es la forma óptima de establecer norte y metro dado un producto constante nm=N ?
ameba dice Restablecer Monica
"Óptimo" en el sentido de minimizar la varianza de la media muestral de los puntos de datos adquiridos. N
ameba dice Reinstate Monica
1
Sí. Pero para su pregunta no necesitamos preocuparnos por cómo estimar las variaciones; su pregunta (es decir, la cita en su pregunta) es que solo creo en estimar la media global y parece obvio que el mejor estimador viene dado por la gran media ˉ x de todos los puntos N = n m de la muestra. La pregunta entonces es: dado μ , σ 2 , σ 2 a , n y m , ¿cuál es la varianza de ˉ x ? Si sabemos eso, podremos minimizarlo con respecto a n dado el n mμx¯N=nmμσ2σa2nmx¯n restricción. nm=N
ameba dice Reinstate Monica
1
No sé cómo derivar nada de eso, pero estoy de acuerdo en que parece obvio: para estimar la varianza del error, sería mejor tener todas las mediciones de un solo sujeto; y para estimar la varianza del sujeto (¿probablemente?) sería mejor tener N sujetos diferentes con 1 medición cada uno. Sin embargo, no está tan claro sobre la media, pero mi intuición me dice que tener N sujetos con 1 medición cada uno también sería lo mejor. Me pregunto si eso es cierto ...NNN
ameba dice Reinstate Monica
2
Tal vez algo así: la varianza de las medias de muestra por sujeto debería ser , donde el primer término es la varianza del sujeto y el segundo es la varianza de la estimación de la media de cada sujeto. Entonces la varianza de la media de sobre-sujetos (es decir, gran media) será ( σ 2 a + σ 2 / n ) / m = σ 2 a / m + σ 2 / ( n m ) = σ 2 a / mσa2+σ2/n que se minimiza cuando m = N .
(σa2+σ2/n)/m=σa2/m+σ2/(nm)=σa2/m+σ2/N=σa2/m+const,
m=N
ameba dice Reinstate Monica

Respuestas:

25

La respuesta corta es que su conjetura es verdadera cuando y solo cuando hay una correlación positiva dentro de la clase en los datos . Hablando empíricamente, la mayoría de los conjuntos de datos agrupados la mayor parte del tiempo muestran una correlación positiva dentro de la clase, lo que significa que en la práctica su conjetura suele ser cierta. Pero si la correlación intraclase es 0, entonces los dos casos que mencionó son igualmente informativos. Y si la correlación intraclase es negativa , en realidad es menos informativo tomar menos medidas en más sujetos; en realidad preferiríamos (en lo que respecta a la reducción de la varianza de la estimación del parámetro) tomar todas nuestras mediciones en un solo tema.

Estadísticamente, hay dos perspectivas desde las cuales podemos pensar en esto: un efecto aleatorio (o mixto ) modelo , que usted menciona en su pregunta, o un modelo marginal , que termina siendo un poco más informativo aquí.

Modelo de efectos aleatorios (mixto)

Digamos que tenemos un conjunto de sujetos de los cuales hemos tomado m mediciones cada uno. Entonces, un modelo simple de efectos aleatorios de la medida j del sujeto i podría ser y i j = β + u i + e i j , donde β es la intersección fija, u i es el efecto aleatorio del sujeto (con varianza σ 2 u ), e i j es el término de error de nivel de observación (con varianza σ 2 enmji

yij=β+ui+eij,
βuiσu2eijσe2), y los dos últimos términos aleatorios son independientes.

En este modelo, representa la media de la población, y con un conjunto de datos equilibrado (es decir, un número igual de mediciones de cada sujeto), nuestra mejor estimación es simplemente la media de la muestra. Entonces, si tomamos "más información" para significar una varianza menor para esta estimación, entonces básicamente queremos saber cómo la varianza de la media muestral depende de n y m . Con un poco de álgebra podemos resolver esa var ( 1βnm Al examinar esta expresión, podemos ver quecada vezquehay alguna variación de sujeto(es decir,σ2u>0), al aumentar el número de sujetos (n), ambos términos serán más pequeños, al tiempo que aumenta el número de mediciones por sujeto (m) solo hará que el segundo término sea más pequeño. (Para una implicación práctica de esto para el diseño de proyectos de replicación de sitios múltiples, veaesta publicación de blog que escribí hace un tiempo).

var(1nmijyij)=var(1nmijβ+ui+eij)=1n2m2var(ijui+ijeij)=1n2m2(m2ivar(ui)+ijvar(eij))=1n2m2(nm2σu2+nmσe2)=σu2n+σe2nm.
σu2>0nm

mnnm

σu2n+constant,
nn=nmm=1

ρ=σu2σu2+σe2
(sketch of a derivation here). So we can write the variance equation above as
var(1nmijyij)=σu2n+σe2nm=(ρn+1ρnm)(σu2+σe2)
This doesn't really add any insight to what we already saw above, but it does make us wonder: since the intra-class correlation is a bona fide correlation coefficient, and correlation coefficients can be negative, what would happen (and what would it mean) if the intra-class correlation were negative?

In the context of the random-effects model, a negative intra-class correlation doesn't really make sense, because it implies that the subject variance σu2 is somehow negative (as we can see from the ρ equation above, and as explained here and here)... but variances can't be negative! But this doesn't mean that the concept of a negative intra-class correlation doesn't make sense; it just means that the random-effects model doesn't have any way to express this concept, which is a failure of the model, not of the concept. To express this concept adequately we need to consider the marginal model.

Marginal model

For this same dataset we could consider a so-called marginal model of yij,

yij=β+eij,
where basically we've pushed the random subject effect ui from before into the error term eij so that we have eij=ui+eij. In the random-effects model we considered the two random terms ui and eij to be i.i.d., but in the marginal model we instead consider eij to follow a block-diagonal covariance matrix C like
C=σ2[R000R000R],R=[1ρρρ1ρρρ1]
In words, this means that under the marginal model we simply consider ρ to be the expected correlation between two es from the same subject (we assume the correlation across subjects is 0). When ρ is positive, two observations drawn from the same subject tend to be more similar (closer together), on average, than two observations drawn randomly from the dataset while ignoring the clustering due to subjects. When ρ is negative, two observations drawn from the same subject tend to be less similar (further apart), on average, than two observations drawn completely at random. (More information about this interpretation in the question/answers here.)

So now when we look at the equation for the variance of the sample mean under the marginal model, we have

var(1nmijyij)=var(1nmijβ+eij)=1n2m2var(ijeij)=1n2m2(n(mσ2+(m2m)ρσ2))=σ2(1+(m1)ρ)nm=(ρn+1ρnm)σ2,
which is the same variance expression we derived above for the random-effects model, just with σe2+σu2=σ2, which is consistent with our note above that eij=ui+eij. The advantage of this (statistically equivalent) perspective is that here we can think about a negative intra-class correlation without needing to invoke any weird concepts like a negative subject variance. Negative intra-class correlations just fit naturally in this framework.

(BTW, just a quick aside to point out that the second-to-last line of the derivation above implies that we must have ρ1/(m1), or else the whole equation is negative, but variances can't be negative! So there is a lower bound on the intra-class correlation that depends on how many measurements we have per cluster. For m=2 (i.e., we measure each subject twice), the intra-class correlation can go all the way down to ρ=1; for m=3 it can only go down to ρ=1/2; and so on. Fun fact!)

So finally, once again considering the total number of observations nm to be a constant, we see that the second-to-last line of the derivation above just looks like

(1+(m1)ρ)×positive constant.
So when ρ>0, having m as small as possible (so that we take fewer measurements of more subjects--in the limit, 1 measurement of each subject) makes the variance of the estimate as small as possible. But when ρ<0, we actually want m to be as large as possible (so that, in the limit, we take all nm measurements from a single subject) in order to make the variance as small as possible. And when ρ=0, the variance of the estimate is just a constant, so our allocation of m and n doesn't matter.
Jake Westfall
fuente
3
+1. Great answer. I have to admit that the second part, about ρ<0, is quite unintuitive: even with a huge (or infinite) total number nm of observations the best we can do is to allocate all observations to one single subject, meaning that the standard error of the mean will be σu and it's not possible in principle to reduce it any further. This is just so weird! True β remains unknowable, whatever resources one puts into measuring it. Is this interpretation correct?
amoeba says Reinstate Monica
3
Ah, no. The above is not correct because as m increases to infinity, ρ cannot stay negative and has to approach zero (corresponding to zero subject variance). Hmm. This negative correlation is a funny thing: it's not really a parameter of the generative model because it's constrained by the sample size (whereas one would normally expect a generative model to be able to generate any number of observations, whatever the parameters are). I am not quite sure what is the proper way to think about it.
amoeba says Reinstate Monica
1
@DeltaIV What is "the covariance matrix of the random effects" in this case? In the mixed model written by Jake above, there is only one random effect and so there is no "covariance matrix" really, but just one number: σu2. What Σ are you referring to?
amoeba says Reinstate Monica
2
@DeltaIV Well, the general principle is en.wikipedia.org/wiki/Inverse-variance_weighting, and the variance of each subject's sample mean is given by σu2+σe2/mi (that's why Jake wrote above that the weights have to depend on the estimate of between-subject variance). The estimate of within-subject variance is given by the variance of the pooled within-subject deviations, the estimate of between-subject variance is the variance of subjects' means, and using all that one can compute the weights. (But I am not sure if this is 100% equivalent to what lmer will do.)
amoeba says Reinstate Monica
1
Jake, yes, it's exactly this hard-coding of m that was bothering me. If this is "sample size" then it cannot be a parameter of the underlying system. My current thinking is that negative ρ should actually indicate that there is another within-subject factor that is ignored/unknown to us. E.g. it could be pre & post of some intervention and the difference between them is so large that the measurements are negatively correlated. But this would mean that m is not really a sample size, but the number of levels of this unknown factor, and that can certainly be hard coded...
amoeba says Reinstate Monica