Skip to content

Commit

Permalink
Correcciones menores
Browse files Browse the repository at this point in the history
  • Loading branch information
felipegonzalez committed Jan 17, 2024
1 parent 67b718a commit 64cbc8e
Show file tree
Hide file tree
Showing 3 changed files with 28 additions and 25 deletions.
31 changes: 12 additions & 19 deletions notas/02-flujo-basico-2.qmd
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@ En esta parte veremos cómo continuaríamos refinando el modelo tomando
en cuenta otros aspectos importantes del problema de estimar la seropositividad
de una población.

El primer aspecto importante es persistir **seguiendo el flujo de trabajo**. En
El primer aspecto importante es persistir **siguiendo el flujo de trabajo**. En
segundo lugar, este ejemplo ilustra que la construcción de modelos preferiblemente
se hace de manera iterativa: una vez que probamos nuestro código, entendemos cómo
funciona nuesto modelo, podemos continuar agregando componentes para hacerlo
Expand Down Expand Up @@ -48,7 +48,7 @@ digraph {
node [shape=plaintext]
N
Npos [label = <N<SUB>+</SUB>>]
Nobs
Nobs [label = <N<SUB>obs</SUB>>]
#Nneg [label = <N<SUB>-</SUB>>]
#sens
#esp
Expand Down Expand Up @@ -126,7 +126,8 @@ un positivo ahora es:
$\theta_{obs} = P(Positivo | \theta, sens, esp) = \theta \cdot sens + (1-\theta) \cdot (1-esp)$

Si llamamos a esta cantidad $\theta_{obs}$, de forma que dada una muestra
de 0's y 1's, tenemos
de 0's y 1's, tenemos que la verosimilitud de la muestra dada cada
conjetura $\theta$, y con $sens$ y $esp$ fijas, es:

$$p(D|\theta, sens, esp) = \theta_{obs}^{N_{+}}(1-\theta_{obs})^{N_{-}}$$
Suponiendo que la distribución apriori de $\theta$ es uniforme, tenemos entonces
Expand Down Expand Up @@ -215,6 +216,8 @@ Ahora usamos histogramas por ejemplo para mostrar cómo luce la posterior,
y comparamos con los valores de la simulación:

```{r}
#| label: fig-apriori-falla
#| fig-cap: "Verificación a priori"
ggplot(simulacion_rep_error, aes(x = theta)) +
geom_histogram(bins = 50) +
labs(x = "theta", y = "Prob posterior") +
Expand Down Expand Up @@ -245,6 +248,8 @@ Ahora usamos histogramas por ejemplo para mostrar cómo luce la posterior,
y comparamos con los valores de la simulación:

```{r}
#| label: fig-apriori-falla
#| fig-cap: "Verificación a priori fallida (modelo incorrecto)"
ggplot(simulacion_rep, aes(x = theta)) +
geom_histogram(bins = 50) +
labs(x = "theta", y = "Prob posterior") +
Expand Down Expand Up @@ -378,8 +383,8 @@ datos_lista <- list(N = 3300, n = 50,
kit_pos = 103, n_kit_pos = 122,
kit_neg = 399, n_kit_neg = 401)
ajuste <- mod_sc$sample(data = datos_lista, refresh = 1000)
sims <- ajuste$draws(c("p", "sens", "esp"), format = "df")
resumen <- ajuste$summary(c("p"))
sims <- ajuste$draws(c("theta", "sens", "esp"), format = "df")
resumen <- ajuste$summary(c("theta"))
```

```{r}
Expand All @@ -389,7 +394,7 @@ resumen |> select(variable, mean, q5, q95)
Y podemos graficar la posterior de la seroprevalencia:

```{r, fig.width = 5, fig.height=3.5}
ggplot(sims, aes(x = p)) +
ggplot(sims, aes(x = theta)) +
geom_histogram()
```

Expand All @@ -413,18 +418,6 @@ para explicar la prevalencia observada, pero si no pensamos con cuidado podríam
concluir que los falsos positivos no deberían ser problema por que la especificidad
para ser muy buena.













Y notamos que aún con una muestra relativamente grande,
el rango de $p$ es considerable: va desde valores cercanos a 0 hasta valores
el rango de $\theta$ es considerable: va desde valores cercanos a 0 hasta valores
alrededor de 0.025-0.03.
16 changes: 13 additions & 3 deletions notas/02-flujo-basico.qmd
Original file line number Diff line number Diff line change
Expand Up @@ -701,7 +701,7 @@ ggplot(datos_graf, aes(x=theta, y = densidad, group = n)) +

---

En este punto, podríamos pasar al siguiente paso, que es
En este punto, podríamos ir al siguiente paso, que es
escribir una función para calcular la posterior. En realidad ya
sabemos su función de densidad, pero cualquier resumen que hagamos de esta
distribución requerirá de integrales (¿por qué? piensa en cómo calcular la probabilidad
Expand Down Expand Up @@ -784,13 +784,23 @@ cuando $M\to \infty$. De este modo, podemos aproximar con la precisión que
requiramos la integral $I$.
:::

**Nota importante**:
Notamos, sin embargo, que sin más información del proceso de simulación,
**Nota 1**: Sin más información del proceso de simulación,
no es posible *demostrar* que una aproximación es "suficientemente" buena, no
importa que tan grande sea $M$. Más
adelante veremos una batería de diagnósticos para al menos excluir los casos comunes
en los que la aproximación es mala.

**Nota 2**: En nuestro caso, las integrales de interés usualmente son de
la forma
$$I = \int f(\theta)p(\theta|D) d\theta,$$
donde $D$ es la información de la muestra, $\theta$ es un vector de parámetros
del modelo, y $f(\theta)$ es una función de $\theta$ que nos interesa. Por ejemplo,
para la media posterior de $\theta_i$, usaríamos $f(\theta) = \theta_i$. Podemos
aproximar cualquier integral si tenemos simulaciones de la posterior:

$$\theta_i \sim p(\theta|D) \implies \frac{1}{M} \sum_{i=1}^{M} f(\theta_i) \to I.$$

---

Finalmente, checamos todo nuestra construcción de estimación como hicimos
arriba, la diferencia es que ahora usamos simulaciones para entender el comportamiento
Expand Down
6 changes: 3 additions & 3 deletions notas/src/sclara.stan
Original file line number Diff line number Diff line change
Expand Up @@ -8,15 +8,15 @@ data {
}

parameters {
real<lower=0, upper=1> p; //seroprevalencia
real<lower=0, upper=1> theta; //seroprevalencia
real<lower=0, upper=1> sens; //sensibilidad
real<lower=0, upper=1> esp; //especificidad
}

transformed parameters {
real<lower=0, upper=1> prob_pos;

prob_pos = p * sens + (1 - p) * (1 - esp);
prob_pos = theta * sens + (1 - theta) * (1 - esp);

}
model {
Expand All @@ -26,7 +26,7 @@ model {
kit_pos ~ binomial(n_kit_pos, sens);
kit_neg ~ binomial(n_kit_neg, esp);
// iniciales para cantidades no medidas
p ~ beta(1.0, 10.0);
theta ~ beta(1.0, 10.0);
sens ~ beta(2.0, 1.0);
esp ~ beta(2.0, 1.0);
}

0 comments on commit 64cbc8e

Please sign in to comment.