Problema locomotor con empresas de varios tamaños

Estoy trabajando a través de Think Bayes (gratis aquí: http://www.greenteapress.com/thinkbayes/ ) y estoy en el ejercicio 3.1. Aquí hay un resumen del problema:

"Un ferrocarril numera sus locomotoras en el orden 1..N. Un día ves una locomotora con el número 60. Estima cuántas locomotoras tiene el ferrocarril".

Esta solución se encuentra con la función de verosimilitud y la exponencial anterior de la siguiente manera:

class Train(Suite):
  def __init__(self, hypos, alpha=1.0):
    # Create an exponential prior
    Pmf.__init__(self)
    for hypo in hypos:
      self.Set(hypo, hypo**(-alpha))
    self.Normalize()
  def Likelihood(self, data, hypo):
    if hypo < data:
      return 0
    else:
      return (1.0/hypo)

Conceptualmente esto es decir, si vemos un número de tren mayor que una de nuestras hipótesis (1 ... 1000), entonces cada hipótesis más pequeña tiene cero posibilidades de ser correcta. El resto de las hipótesis tienen una probabilidad de 1 / number_of_trains de mostrarnos un tren con este número.

En el ejercicio en el que estoy trabajando, el autor agrega un poco más. Esto supone que solo hay una compañía. Sin embargo, en la vida real tendrías una mezcla de compañías grandes y pequeñas y compañías más grandes (ambas igualmente probables). Sin embargo, esto significaría que es más probable que vea un tren de una compañía más grande ya que tendrían más trenes.

Ahora la pregunta es ¿cómo reflejar esto en la función de probabilidad?

Esto no es Stack Overflow, así que realmente no estoy pidiendo ayuda para la codificación, sino que tal vez solo ayude sobre cómo podría pensar sobre este problema en términos de una función de probabilidad.

Primero esbozo un enfoque para dos compañías en detalle, la extensión a incluso más compañías debería ser intuitiva (al menos para la probabilidad, lo anterior podría ser más complicado).

Imagine que hay dos compañías A y B , donde A tiene locomotoras y B tiene locomotoras . Asumimos (siempre puede cambiar A y B para hacer que esto se mantenga). El número total para esa hipótesis de locomotoras es .$N_A$$N_B$$N_A \ge N_B$$N_{tot} = N_A + N_B$

Imagina que ves una locomotora con el número . Hay tres casos para la probabilidad:$n$

1. $N_A < n$ : Esto no puede suceder, por lo que la probabilidad es cero.
2. $N_B < n \le N_A$ : esta locomotora debe provenir de la compañía A , por lo que solo hay una locomotora con este número. Por lo tanto, la probabilidad es$1/N_{tot}$
3. $n \le N_B$ : esta locomotora puede ser de A o de B , por lo que hay dos locomotoras con este número. La probabilidad de ver uno de ellos es .$2/N_{tot}$

Como una verificación rápida de la cordura: La probabilidad de ver cualquier número es .

$$\sum_{i=1}^\infty L(i) = \sum_{i=1}^{N_B} \frac{2}{N_{tot}} + \sum_{i=N_B+1}^{N_A} \frac{1}{N_{tot}} \\ = \frac{2\cdot N_B}{N_{tot}} + \frac{N_A-N_B}{N_{tot}} = \frac{N_A+N_B}{N_{tot}} = 1$$

---

En general, habrá (número de empresas + 1) casos, uno para cada intervalo . Afortunadamente, podemos ver el problema desde un ángulo diferente y ver que lo que necesitamos para la probabilidad son en realidad solo dos números: , el número total de locomotoras; y , el número de locomotoras que tienen el número . ¿Qué tan probable es que veamos una de las locomotoras , fuera de las locomotoras ? Esto sucederá en de todos los casos, por lo que esta fracción es la probabilidad. En Python, puede calcular esto con dos generadores de suma (y ni siquiera tiene que ordenar las empresas por tamaño). Si$N_i < n \le N_{i+1}$$N_{tot}$$N_n$$n$$N_n$$N_{tot}$$\frac{N_n}{N_{tot}}$Nscontiene una lista (o tupla) de tamaños de compañía de acuerdo con su hipótesis, entonces esto le dará la posibilidad de ver una locomotora con número n:

total_number_of_locomotives = sum(N for N in Ns)
number_of_locomotives_with_that_number = sum(1 for N in Ns if n<=N)
likelihood = (number_of_locomotives_with_that_number / total_number_of_locomotives)

Tenga en cuenta que el caso trivial con una compañía también se maneja con este código (la primera suma será , la segunda suma será 0 o 1, dependiendo de si ).$N$$n\le N$

---

Para los anteriores, la ley de Zipf podría ser un buen punto de partida para una distribución realista de los tamaños de las empresas.

No voy a analizar el código, pero a continuación está la solución.

Dejar

• P (loc60) sea la probabilidad de que una locomotora aleatoria tenga el número 60
• P (N) sea la probabilidad previa de que existan exactamente N locomotoras
• P (loc60 | N) sea la probabilidad de que una locomotora aleatoria tenga el número 60, si el número total de locomotoras es N,
• P (N | loc60) sea la probabilidad de que existan exactamente N locomotoras, si una locomotora aleatoria tiene el número 60

Entonces

$$P(N|\text{loc60}) = {P(\text{loc60}|N) P(N) \over P(\text{loc60})} = {P(\text{loc60}|N) P(N) \over \sum_M P(\text{loc60}|M)}$$

Pero

$$P(\text{loc60}|N) = \cases {1/N & if $N\ge 60$ \\ 0 & otherwise }$$

De ahora en adelante, asumimos que .$N \ge 60$

$$P(N|\text{loc60}) = {P(N)/N \over \sum_{M=60}^\infty P(M)/M}$$

Ahora debemos seleccionar P (N), de lo contrario estamos atascados. Como no conocemos ni siquiera el orden de magnitud de P (N), es razonable suponer que está distribuido uniformemente entre 0 y algunos (es decir, la probabilidad de que es lo mismo que la probabilidad de que ). Guestimating es una tarea difícil, pero por mi conocimiento previo sobre ferrocarriles y locomotoras, puedo suponer que .$\log N$$\log N_\max$$10^2\le N<10^3$$10^3\le N<10^4$$N_\max$$N_\max \gg 60$

La distribución uniforme de significa que , donde c es una constante independiente de N.$\log N$

$$P(N) = c(\log (N+1)-\log N) \approx c/N$$

Sustituyendo esto por la fórmula anterior, tenemos:

$$P(N|\text{loc60}) \approx {c/N^2 \over \sum_{M=60}^{N_\max} c/M^2}$$

Pero

$$\sum_{M=60}^{N_\max} c/M^2 \approx \int_{60}^{N_\max} {c\over M^2}dM = {c \over 60} - {c \over N_\max} \approx {c\over60}$$

Ahora tenemos

$$P(N|\text{loc60}) \approx {60/N^2}$$

¿Cuál es el valor medio de N? Deja que sea , luego$N_\text{med}$

$$\int_{60}^{N_\text{med}} {60 \over N^2} dN = 1/2$$

$$60/N - {60 \over N_\text{med}} = 1/2$$

$$N_\text{med} = 120$$

Si lo que necesitamos es expectativa matemática en lugar de mediana, entonces

$$E(N) = \int_{60}^{N_\max} {60\over N^2} N dN = 60 \log {N_\max \over 60}$$

Por lo que sé sobre los ferrocarriles, debería estar entre y , por lo que E (N) está en algún lugar entre 170 y 600.$N_\max$$10^3$$10^6$