Estadística

Primero, una simulación:

n <- 100
delta <- 0.2
n.iter <- 10000

p_valores <- function(n, delta){
  tmp <- replicate(n.iter, {
    x <- rnorm(n)
    y <- rnorm(n, mean = delta)
    t.test(x, y)$p.value
  })

  res <- tmp[tmp < 0.05]

  hist(res, freq = FALSE, xlab = "p value", ylab = "", col = "gray", main = "histograma de p-valores publicables")

  res
}

null_effect_p_values <- p_valores(n, 0)
some_effect_p_values <- p_valores(n, delta)

Lo que simula son n.iter experimentos en los que se comparan n valores N(0,1) con otros n valores N(delta, 1) y se extrae el correspondiente p-valor. Luego se grafican los publicables (<0.05).

Alguno de mis lectores, supongo, estará metido en ese mundo de ir escribiendo cosas y cosechado méritos, impactos y anecosas para salir del precariado y pillar moscosos. Que dejen de leer. Es una orden.

A aquellos que tengan tiempo y talento los invito a escribir el artículo titulado Temperaturas umbrales de disparo de la mortalidad atribuible al frío y al calor en España en el periodo 2007-2017.

Se trata, esencialmente, de aggiornar metodológica, gráfica y sintácticamente esta cosa viejuna y manifiestamente mejorable en todas las dimensiones concebibles.

Cuando tienes una serie temporal al uso (sin entrar a definir qué es eso), uno puede aplicar descomposiciones tmabién al uso, como stl, para extraer tendencia y estacionalidad, de la forma

como en esta entrada previa.

Lluvia.

La serie de la lluvia es otra cosa. Uno ve si llueve o no llueve (típicamente no). Lo que uno no ve es la probabilidad, que varía a lo largo del año, de que llueva. Pasa lo mismo con monedas (uno ve caras o cruces, no la probabilidad de cara), clientes que compran (uno ve si compra o no, no la probabilidad de compra), etc. Pero la probabilidad existe y en estimarla consiste frecuentemente el trabajo de algunos.

Hoy toca esto:

Esto es lo que provoca la contaminación: los picos de contaminación coinciden con un aumento radical en los ingresos de los hospitales https://t.co/GpEBg6hqko pic.twitter.com/tvwS1r3Ldi

— Ignacio Escolar (@iescolar) November 23, 2017

Se trata de una invitación para leer el artículo Los picos de contaminación coinciden con un aumento radical en los ingresos hospitalarios, un cúmulo de desafueros epilogados por el ya habitual

Los resultados de esta investigación tienen puntos en común con la metodología científica aunque en ningún momento tendrán la misma validez ni tampoco es su intención que la tenga.

Se ve que hay arqueólogos bayesianos. Un problema con el que se encuentran es que tropiezan con cacharros antiguos y quieren estimar su antigüedad.

Así que prueban distintos métodos (¿químicos?), cada uno de los cuales con su precisión, y acaban recopilando una serie de estimaciones y errores. Obviamente, tienen que combinarlas de alguna manera.

El modelo más simple es

$$ M_i \sim N(\mu, \sigma_i)$$

donde $latex \mu$ es la antigüedad (desconocida) del artefacto y los $latex \sigma_i$ son las varianzas distintas de los distintos métodos de medida, que arrojan las estimaciones $latex M_i$.

Hay disciplinas que parecen puras colecciones de anécdotas, recetarios ad hoc y listas de contraejemplos. Tal se ha predicado, por ejemplo, de la economía conductual.

Pero, ¿pueden reconstruirse a partir de primeros principios? Si se ha ensayado con la economía conductual, ¿por qué no intentarlo con nuestra modestísima estadística descriptiva?

Un caso particular: cuando de una variable aleatoria calculo y escribo o represento su media y su desviación estándar, de alguna manera estoy modelizándola como una distribución normal. Esta modelización puede ser explícita, aunque casi siempre es implícita. Si la variable aleatoria tiene una distribución muy alejada de la normal, habrá quien proteste: que si la media es engañosa, que si… Pero, ¿por qué habría de ser engañosa en este caso y no en otro? Precisamente por la (incorrecta) modelización implícita: estaría usando lo de la normal donde no aplica.

El objeto único de la estadística es informar decisiones.

V.g, si conceder un préstamo, proceder a la quimio, construir una línea de AVE entre Calatayud y Soria o permitir aparcar mañana en el centro de Madrid a los de Móstoles.

Pero quienes toman decisiones y quienes analizan datos suelen ser personas distintas. Típicamente, ni se conocen. Lo cual es tanto pésimo como tema para otra entrada distinta de esta.

Lo fundamental es que estas personas se comunican a través de, metafóricamente, APIs. Unas de las más usadas son los p-valores. Que son tan pésismos como tema para otra entrada distinta de esta.

He usado el vídeo

en un curso de estadística básica para ilustrar a través de experimentos se construyen histogramas y estos convergen a y, en última instancia, justifican el uso de distribuciones de probabilidad.

Es decir,

experimentos -> histogramas -> funciones de distribución.

Y de ahí, el resto.

Envalentonado por el comentario de Iñaki Úcar a mi entrada del otro día, que me remitía a este artículo, decidí rizar el rizo y crear intervalos de confianza no ya discontinuos sino con otra propiedad topológica imposible: homeomorfos con un toro.

Y aquí está:

El modelo, el código y demás,

library(rstan)
library(ggplot2)

n <- 100

a1 <- 1
a2 <- 1
sigma <- 0.4

datos <- data.frame(x1 = rnorm(n, 2, 0.1),
                    x2 = rnorm(n, 2, 0.1))

datos$y <- a1^datos$x1 + a2^datos$x2 + rnorm(n, 0, sigma)

codigo <- "
data {
  int<lower=1> N;
  real y[N];
  real x1[N];
  real x2[N];
}

parameters {
  real<lower=-3, upper="3"> a1;
  real<lower=-3, upper="3"> a2;
  real<lower=0, upper="3"> sigma;
}

model {
  for (n in 1:N)
    y[n] ~ normal(fabs(a1)^x1[n] +
      fabs(a2)^x2[n], sigma);
}"

fit <- stan(model_code = codigo,
    data = list(N = length(datos$y), y = datos$y,
                x1 = datos$x1, x2 = datos$x2),
    iter=40000, warmup=2000,
    chains=1, thin=10)

res <- as.data.frame(fit)

ggplot(res, aes(x = a1, y = a2)) + geom_point(alpha = 0.1)

De nuevo, no son intervalos propiamente dichos, lo convengo. Pero son configuraciones más fieles al espíritu de lo que un intervalo de confianza es y representa que su(s) letra(s) I N T E R V A L O.

Me sorprende haber leído tantos de los mejores 100 libros del siglo XX según Le Monde. Sobre todo porque no leo ficción casi en lo que va de siglo y porque, carajo, los libros estupendos que he leído de tapa, como el Análisis Real de Folland o la Introducción a la Teoría de la Probabilidad de Feller parece que no cualifican para esa listeja de textos sin una mala integral preparada por gentecilla de letras.