Loading [MathJax]/jax/element/mml/optable/BasicLatin.js

6 votos

Combinación lineal de variables discretas Ti con P(Ti=1)=P(Ti=1)=1/2

Dejemos que T1,...,Tn sea iid con una distribución Rademacher; es decir , P(Ti=1)=P(Ti=1)=1/2 ; y dejar que w=(w1,...,wn)Rn sin más limitaciones en w .

¿Existe una forma de calcular P(T1ni=1Tiwi0) sin pasar por todos los 2n posibles resultados de la T 's?

0 votos

¿Se asumen algunas condiciones adicionales, como wi0 i y ni=1=1 ? ¿Qué pasa con la dependencia entre los Ti -- ¿se supone que son iid?

0 votos

No hay más condiciones. Sí, iid.

0 votos

No creo que este problema tenga una solución interesante si no se ponen al menos algunas condiciones a los posibles valores de los términos wi , i=1,,n . Poner Y=T1(ni=1wiTi)=w1+T1(ni=2wiTi) (ya que T211 ). Si w1 es (digamos) un número positivo grande s.t. w1>ni=2|wi| entonces P(Y0)=1 . Para obtener un resultado no trivial, creo que se necesita, como mínimo, algún tipo de límite en max .

6voto

jldugger Puntos 7490

Sí.

Empecemos por simplificar la cuestión. El evento

T_1\sum_{i=1}^n T_i w_i \geq 0

es la unión de los eventos disjuntos (determinados por T_1=\pm 1 )

w_1 + \sum_{i=2}^n T_i w_i \geq 0,\ -w_1 + \sum_{i=2}^n T_i w_i \geq 0.

Tras una sencilla manipulación algebraica, y utilizando el hecho de que cada T_i tiene la misma distribución que -T_i (es "simétrico"), estos eventos son ambos equivalentes a los de la forma

\sum_{i=2}^n T_i w_i \leq w

para w=\pm w_1 . Aprovechando de nuevo las simetrías del T_i (para que todos los w_i no negativo) y reordenándolos, se puede escribir

\sum_{i=2}^n T_i |w_i| \leq w

donde |w_2| \ge |w_3| \ge \cdots \ge |w_n| .

Necesitamos calcular las probabilidades de que se produzcan dos eventos de este tipo y deseamos hacerlo de forma eficiente. Asumiré que la probabilidad debe calcularse con gran precisión, de modo que las aproximaciones serían inaceptables. (Para muchas configuraciones del |w_i| De lo contrario, las aproximaciones normales a las distribuciones podrían funcionar bien. Se puede realizar un análisis más profundo utilizando la función característica \phi de \sum_{i=2}^n T_i |w_i| que por definición es \phi(t) = \prod_{i=2}^n\left(\frac{\exp(-j |w_i|t)}{2} + \frac{\exp(j |w_i|t)}{2}\right)=\prod_{i=2}^n\cos(w_it); j = \sqrt{-1} .)


A partir de ahora, reduzca n por 1 , inicie la indexación con i=1 y asumir todos los w_i son positivos. (Obviamente, cualquier valor cero puede ser ignorado). 1 \le k \lt n . Considere la distribución de X = \sum_{i=1}^n T_i w_i condicionada a la primera k valores de T_i :

\Pr(X \le w\,|\, T_1, \ldots, T_k) = {\Pr}_{(T_{k+1},\ldots,T_n)}\left(\sum_{i=k+1}^n T_i w_i \le w - \sum_{i=1}^k T_i w_i = w_0\right).

Desde |T_i| \le 1 la suma de la izquierda está acotada:

-u_{k+1} = -(w_{k+1} + \cdots + w_n) \le \sum_{i=k+1}^n T_i w_i \le w_{k+1} + \cdots + w_n = u_{k+1}.

Dejemos que u_{k+1} sea el lado derecho: es una suma acumulativa de los w_i acumulado a partir de los valores más bajos de la derecha. Obviamente, ahora si u_{k+1} \le w_0 entonces X \le w es seguro; y si -u_{k+1} \gt w_0 entonces X \le w tiene una probabilidad nula. Esto lleva a una simple rama y límite algoritmo, porque no necesitamos buscar más para evaluar la distribución. Mientras que una enumeración exhaustiva habría tenido que examinar 2^{n-k} casos posibles, hemos hecho una determinación de su contribución a la distribución en O(1) tiempo.

No hay muchas esperanzas de que esta mejora conduzca a una el peor de los casos algoritmo que es mejor que O(2^n) en el rendimiento (aunque creo que se puede reducir a O(2^{n/2}) que -aunque mucho mejor- sigue sin ser polinómico). Sin embargo, la mejora es lo suficientemente buena como para que valga la pena considerarla, especialmente para valores moderados de n donde la enumeración exhaustiva empieza a ser impracticable (en algún punto por encima de 20 y ciertamente por encima de 40 ). Por lo tanto, pasemos a evaluar el algoritmo. ¿Cuál es su eficacia?

El peor caso es cuando la mayoría de los w_i tienen tamaños comparables, ya que entonces la heurística branch-and-bound rara vez logra algo. Afortunadamente, esta es exactamente la situación en la que el Teorema Central del Límite puede proporcionar excelentes aproximaciones. Por lo tanto, quizás sea aún más interesante explorar situaciones en las que los tamaños de las w_i están muy repartidos.

Para ello, he creado cuatro conjuntos de datos con n=15 de cuatro distribuciones con formas radicalmente diferentes: Beta (1/5,1/5) (fuertemente bimodal), Exponencial (modo alto cerca de cero), Uniforme y Gamma (20) (casi Normal, con valores muy cercanos a 20 ). Cada conjunto de datos se normalizó a un máximo de 1 . Utilizando el método de rama y límite, \Pr(X \le w) se calculó para 19 valores de w que van desde 0 hasta 2\sqrt{n} . (Los valores negativos de w no es necesario mostrarlo porque las distribuciones de X son casi simétricos). La figura muestra los resultados, graficando las probabilidades (la función de distribución acumulativa empírica de X ) en la fila superior y el eficiencia (en escalas log-lineales) en la fila inferior. La "eficiencia" es la relación entre el número de configuraciones necesarias para la enumeración exhaustiva ( 2^n ) al número de configuraciones consideradas por el algoritmo branch-and-bound.

Figure

Para la distribución Beta, cuya fuerte bimodalidad limita esencialmente los cálculos a la mitad superior de los datos, las eficiencias son mayores. Como era de esperar, son menores para la distribución Gamma, cuyos valores son todos comparables. Sin embargo, incluso en este caso difícil, todas las eficiencias superan 3 .

Las eficiencias más pequeñas suelen ser para w=0 . Con el tiempo, la eficiencia aumentará a medida que |w| aumenta. En conjunto, estos resultados son una fuerte evidencia de que el enfoque descrito aquí no sólo es computacionalmente más eficiente que la enumeración exhaustiva, sino que tiende a ser mucho más eficiente.


La aplicación del algoritmo es directa y sencilla. R El código es el siguiente. Incluye una sección que prueba el algoritmo (comparando su salida con una enumeración exhaustiva) y otra sección para reproducir la figura.

#
# The algorithm.
#
f <- function(w0, w) {
  w <- sort(abs(w), decreasing=TRUE)
  w.sum <- c(rev(cumsum(rev(w)))[-1], 0)
  count <- 0                              # Counts calls to f()
  f <- function(u, u.sum, w0) {
    count <<- count + 1
    y <- sapply(w0 + c(-1,1)*u[1], function(w1) {
      if (w1 < -u.sum[1]) return(0)
      if (w1 >= u.sum[1]) return(1)
      return(f(u[-1], u.sum[-1], w1))
    })
    return(mean(y)) # (This could easily be changed to accommodate 
                    # other probabilities for the T[i].)
  }
  list(Value=f(w, w.sum, w0), Count=count)
}
#
# Test with complete enumeration.  The plot pairs should exactly coincide.
#
binary <- function(a, b, zero=-1, one=1) rep(c(rep(zero, 2^a), rep(one, 2^a)), 2^b)
n <- 9
b <- sapply(1:n, function(a) binary(n-a, a-1))
par(mfrow=c(2,2))
set.seed(17)
for (i in 1:4) {
  w <- rexp(n)
  x <- b %*% w
  y <- sapply(x, function(w0) f(w0, w)$Value) #$
  plot(ecdf(x))
  points(x, y, pch=16, cex=1/2, col="Red")
}
#
# Explore the efficiencies actually achieved.
#
n <- 15
qb <- function(q) qbeta(q, 1/5, 1/5); sigma <- sqrt(1/(1+2*1/5))/2
qg <- function(q) qgamma(q, 20)
dist <- list(Gamma=qg, Normal=qnorm, Uniform=qunif, Exponential=qexp, Beta=qb)
par(mfcol=c(2,4))
for (s in c("Beta", "Exponential", "Uniform", "Gamma")) {
  d <- dist[[s]]
  w <- d(((1:n)-1/2)/n)
  w <- w / max(abs(w))
  print(c(system.time({
    y <- sapply(x <- seq(0, 2*sqrt(n), length.out=19), 
                function(w0) {x <- f(w0, w); c(x$Value, x$Count)})
  })["elapsed"], count=sum(y[2, ])))

  plot(x, y[1, ], ylim=c(1/2, 1), type="b", ylab="Probability", main=s)
  plot(x, 2^n/y[2, ], ylim=c(1, 2^n), type="b", log="y", ylab="Efficiency")
}

1voto

Jeff Bauer Puntos 236

T_1\sum_{i=1}^n T_i w_i = w_1 + T_1\sum_{i=2}^n T_i w_i

desde T_1^2 = 1 .

Escribe la suma como S_{2w} . Así que estamos viendo la probabilidad

P\left(T_1S_{2w} \geq -w_1\right) = \\P\left(T_1S_{2w} \geq -w_1 \mid T_1 =1\right)\cdot P(T_1 =1) + P\left(T_1S_{2w} \geq -w_1 \mid T_1 =-1\right)\cdot P(T_1 =-1)

\Rightarrow P\left(T_1S_{2w} \geq -w_1\right)=\frac 12 P\left(S_{2w} \geq -w_1 \right)+ \frac 12P\left(-S_{2w} \geq -w_1 \right)

=\frac 12 \Big[1-P\left(S_{2w} \leq -w_1 \right)+P\left(S_{2w} \leq w_1 \right)\Big] \tag{1}

Ahora defina Z_i \equiv (T_i + 1)/2 . Entonces el Z_i son Bernoulli (1/2) variables aleatorias.

Así que

S_{2w}=\sum_{i=2}^n (2Z_i-1) w_i = 2\sum_{i=2}^n Z_iw_i - \sum_{i=2}^n w_i \tag{2}

Inserción en (1) y reordenando obtenemos

P\left(T_1S_{2w} \geq -w_1\right) = \frac 12 \left[1-P\left(\sum_{i=2}^n Z_iw_i \leq (1/2)\sum_{i=2}^n w_i-(w_1/2) \right)+P\left(\sum_{i=2}^n Z_iw_i \leq (1/2)\sum_{i=2}^n w_i+(w_1/2) \right)\right]

Así que necesitamos probabilidades relacionadas con S_{Zw}\equiv \sum_{i=2}^n Z_iw_i . El Función generadora de probabilidad , PGF de un Bernoulli (1/2) es

G_Z(z) = \frac 12(1+z)

y de la combinación lineal de Bernoullis independientes

G_{S_{Zw}}(z) = \prod_{i=2}^{n}\left(\frac 12(1+z^{w_i})\right) = \frac 1{2^{n-1}}\prod_{i=2}^{n}(1+z^{w_i})

Desde el PGF podemos recuperar las probabilidades de la función de masa de probabilidad y, por tanto, de la función de distribución acumulativa que es nuestro objetivo.

UN EJEMPLO

A efectos de exposición, supongamos que n=3 y (w_1 =2, w_2 =2, w_3 = 4) .

Entonces (1/2)\sum_{i=2}^n w_i-(w_1/2) = 2, \;\;\; (1/2)\sum_{i=2}^n w_i+(w_1/2) = 4

Así que queremos evaluar

P\left(T_1S_{2w} \geq -1\right) = \frac 12 \left[1-P\left(S_{Zw} \leq 2 \right)+P\left(S_{Zw} \leq 4 \right)\right]

=\frac 12 \left[1+P\left(S_{Zw} =3 \right)+P\left(S_{Zw} = 4 \right)\right]

En un ejemplo tan simple, no necesitamos realmente el PGF Dado los valores supuestos de los pesos en el Bernoullis, es fácil ver que P\left(S_{Zw} =3 \right) = 0 y que P\left(S_{Zw} = 4\right) = 1/4 . Así obtenemos la respuesta, para el ejemplo concreto

P\left(T_1S_{2w} \geq -1\right) = 5/8

Verificación a través del FGP
Para el ejemplo numérico concreto

G_{S_{Zw}}(z) = \frac 1{2^{n-1}}\prod_{i=2}^{n}(1+z^{w_i}) = \frac 1{4}(1+z^2)(1+z^4) = \frac 1{4}(1+z^2+z^4+z^6)

La relación que une el PGF con probabilidades es ( (k) denotando el orden de la derivada tomada

P(S_{Zw} = k) =\frac{ G_{S_{Zw}}^{(k)}(0)}{k!}

Tenemos

G_{S_{Zw}}^{(1)}(z) = \frac 1{4}(2z+4z^3+6z^5) G_{S_{Zw}}^{(2)}(z) = \frac 1{4}(2+12z^2+30z^4) G_{S_{Zw}}^{(3)}(z) = \frac 1{4}(24z+120z^3) G_{S_{Zw}}^{(4)}(z) = \frac 1{4}(24+360z^2)

Así que

P\left(S_{Zw} =3 \right) = \frac{ G_{S_{Zw}}^{(3)}(0)}{3!} = 0 P\left(S_{Zw} =4 \right) = \frac{ G_{S_{Zw}}^{(4)}(0)}{4!} = \frac 1{4}

Por lo tanto,

P\left(T_1S_{2w} \geq -1\right) = \frac 12 \left[1+0+\frac 14 \right] = 5/8

Ahora imagina lo que supondría si el ejemplo no fuera tan sencillo.

1 votos

Buena respuesta, pero Z_i = T_i + 1/2 no es una variable aleatoria Bernoulli(1/2). Se necesita Z_i = (T_i + 1) / 2 por eso. Aunque tu enfoque debería seguir funcionando. Sólo cambia el álgebra al final.

0 votos

@goodepic Ah, cierto, confundí las probabilidades con el soporte, y por error tuve en cuenta T_i \in \{-1/2,1/2\} . Gracias por detectar esto, comprobaré si el enfoque sobrevive de todos modos.

0 votos

@goodepic Creo que ya está bien.

i-Ciencias.com

I-Ciencias es una comunidad de estudiantes y amantes de la ciencia en la que puedes resolver tus problemas y dudas.
Puedes consultar las preguntas de otros usuarios, hacer tus propias preguntas o resolver las de los demás.

Powered by:

X