¿Cómo modelar la suma de variables aleatorias de Bernoulli para datos dependientes?

Tengo casi las mismas preguntas como esta: ¿Cómo puedo modelar eficientemente la suma de variables aleatorias de Bernoulli?

Pero la configuración es bastante diferente:

1. $S=\sum_{i=1,N}{X_i}$ , , ~ 20, ~ 0.1$P(X_{i}=1)=p_i$$N$$p_i$

2. Tenemos los datos para los resultados de las variables aleatorias de Bernoulli: ,$X_{i,j}$$S_j=\sum_{i=1,N}{X_{i,j}}$

3. Si estimamos el con la estimación de máxima verosimilitud (y obtenemos ), resulta que es mucho más grande que esperado por los otros criterios:$p_i$$\hat p^{MLE}_i$$\hat P\{S=3\} (\hat p^{MLE}_i)$$\hat P\{S=3\} (\hat p^{MLE}_i) - \hat P^{expected} \{S=3\}\approx 0.05$

4. Por lo tanto, y no pueden tratarse como independientes (tienen una pequeña dependencia).$X_{i}$$X_{j}$ $(j>k)$

5. Hay algunas restricciones como estas: y (conocido), lo que debería ayudar con la estimación de .$p_{i+1} \ge p_i$$\sum_{s \le 2}\hat P\{S=s\}=A$$P\{S\}$

¿Cómo podríamos intentar modelar la suma de las variables aleatorias de Bernoulli en este caso?

¿Qué literatura podría ser útil para resolver la tarea?

ACTUALIZADO

Hay algunas ideas más:

(1) Es posible suponer que la dependencia desconocida entre comienza después de 1 o más éxitos en serie. Entonces, cuando , y .${X_i}$$\sum_{i=1,K}{X_i} > 0$$p_{K+1} \to p'_{K+1}$$p'_{K+1} < p_{K+1}$

(2) Para usar MLE necesitamos el modelo menos cuestionable. Aquí hay una variante:

$P\{X_1,...,X_k\}= (1-p_1) ... (1-p_k)$ si para cualquier k if y , y para cualquier k.$\sum_{i=1,k}{X_i} = 0$$P\{X_1,...,X_k,X_{k+1},...,X_N\}= (1-p_1) ... p_k P'\{X_{k+1},...,X_N\}$$\sum_{i=1,k-1}{X_i} = 0$$X_k = 1$$P'\{X_{k+1}=1,X_{k+2}=1,...,X_N=1\} \le p_{k+1} p_{k+2} ... p_N$

(3) Dado que solo nos interesa , podemos establecer (la probabilidad de éxitos para N- (k + 1) +1 sumandos desde la cola). Y use la parametrización$P\{S\}$$P'\{X_{k+1},...,X_N\} \approx P''\{\sum_{i=1,k}{X_i}=s' ; N-(k+1)+1=l\}$$\sum_{i=k+1,N}{X_i}$$P''\{\sum_{i=k,N}{X_i}=s' ; N-k+1=l\}= p_{s',l}$

(4) Use MLE para el modelo basado en los parámetros y con para (y cualquier ) y algunas otras restricciones nativas .$p_1,...,p_N$$p_{0,1}, p_{1,1}; p_{0,2}, p_{1,2}, p_{2,2};...$$p_{s',l}=0$$s' \ge 6$$l$

¿Está todo bien con este plan?

ACTUALIZADO 2

Algunos ejemplos de distribución empírica (rojo) en comparación con la distribución de Poisson (azul) (las medias de poisson son 2.22 y 2.45, los tamaños de muestra son 332 y 259):$P\{S\}$

Para muestras (A1, A2) con las medias de Poisson 2.28 y 2.51 (los tamaños de muestra son 303 y 249):

Para samlpe A1 + A2 unido (el tamaño de la muestra es 552):

Parece que alguna corrección a Poisson debería ser el mejor modelo :).

Un enfoque sería modelar las con un modelo lineal generalizado (GLM). Aquí, formularía , la probabilidad de éxito en la prueba como prueba (función logística lineal) del historial de observación reciente. Por lo tanto, básicamente está ajustando un GLM autorregresivo donde el ruido es Bernoulli y la función de enlace es logit. La configuración es:$X$$p_i$$i$

$p_i = f(b + a_1 X_{i-1} + a_2 X_{i-2} + \ldots a_k X_{i-k})$ , donde

$f(x) = \frac{1}{1+\exp(x)}$ , y

$X_i \sim Bernoulli(p_i)$

Los parámetros del modelo son , que pueden estimarse por regresión logística. (Todo lo que tiene que hacer es configurar su matriz de diseño utilizando la parte relevante del historial de observación en cada ensayo, y pasarla a una función de estimación de regresión logística; la probabilidad de registro es cóncava, por lo que hay un máximo global único para los parámetros). Si los resultados son realmente independientes, entonces los se establecerán en cero; positivo significa que los posteriores aumentan cada vez que se observa un éxito.$\{b, a_1, \ldots a_k\}$$a_i$$a_i$$p_i$

El modelo no proporciona una expresión simple para la probabilidad sobre la suma de las 's, pero esto es fácil de calcular por simulación (filtrado de partículas o MCMC) ya que el modelo tiene una estructura Markovian simple.$X_i$

Este tipo de modelo se ha utilizado con gran éxito para modelar las dependencias temporales entre los "picos" de las neuronas en el cerebro, y existe una extensa literatura sobre modelos de procesos puntuales autorregresivos. Ver, por ejemplo, Truccolo et al 2005 (aunque este documento usa una probabilidad de Poisson en lugar de Bernoulli, pero el mapeo de uno a otro es sencillo).

Si la dependencia se debe a la aglomeración, un modelo compuesto de Poisson podría ser la solución como modelo de . Una referencia algo aleatoria es esta de Barbour y Chryssaphinou.$S_j$

En una dirección completamente diferente, ya que indica que es 20, y por lo tanto relativamente pequeño, podría ser construir un modelo gráfico de los , pero no sé si su configuración y sus datos lo hacen posible. Como comenta @chl, será útil si describe cuáles son las .$N$$X_{ij}$$X_{i,j}$

Si las representan mediciones secuenciales, por ejemplo, a lo largo del tiempo, y la dependencia está relacionada con esto, una tercera posibilidad, y hasta cierto punto un compromiso entre las dos sugerencias anteriores, es usar un modelo oculto de Markov de las 's.$X_{i,j}$$X_{i,j}$