As sondagens americanas falharam - Como se avaliam sondagens?

Tenho as sondagens e os resultados e agora vou comparar as coisas.

É normal ouvir dizer que a estatística não é uma ciência exacta mas, mesmo que só dê probabilidades, os seus resultados podem ser avaliados. 

Claro que um aluno que fez uma cadeira de estatística não ficou com ferramentas conceptuais que o tornem capaz de fazer esta avaliação, limitação que não se coloca a mim :-)

Tenho que desenhar a experiência estatística.

Começo por identificar H0 e, por oposição, H1.

    H0: As sondagens foram retiradas da população (que são os resultados);

    H1: As sondagens não foram retiradas da população (estavam erradas).

Depois, repito a experiência.

O resultado deu 0,00% de probabilidade de se observar H0, isto é, de as amostras na sua globabilidade estarem correctas.

Código R utilizado.

dados <- data.frame(

  Estado = c("Arizona","Georgia","North Carolina","Pensilvania","Nevada","Wisconsin","Michigan"),

  Respostas = c(34311, 37616, 38732, 57342, 60630, 19176,43634),

  Harris = c(0.4666, 0.4851, 0.4780, 0.4856, 0.4748, 0.4897, 0.4848),

  Trump = c(0.5230, 0.5072, 0.5117, 0.5048, 0.5061, 0.4984, 0.4986),

  Outros = c(0.0104, 0.0077, 0.0103, 0.0096, 0.0192, 0.0119, 0.0166),

  dif_sond=c(0.0190,0.0127,0.0112,0.0010,-0.0001,-0.0021,-0.0047)

  )

# Função para realizar o sorteio e calcular a soma dos deputados de quem ganhar em cada estado

simula_votos <- function(dados) {

   #Avalia se poderiam sair sondagens iguais ou piores que as publicadas 

   avalia=TRUE

   # Percorre os 7 estados

   for (i in 1:nrow(dados)) {

      # Sorteia amostras dos resultados para Harris, Trump e Outros

      votos= sample(c(1,2,3), dados$Respostas[i], replace = TRUE, prob = c(dados$Harris[i],dados$Trump[i],dados$Outros[i]))

      votos = table(votos)

      # diferença na amostra sorteada

      def_TH = (votos[2]-votos[1])/dados$Respostas[i]

      avalia = (avalia) & (def_TH <=  dados$dif_sond[i])

      }

    if(avalia) 1 else 0

    }

niter=10000

avalia=rep(0,niter)

for (i in 1:niter){

    avalia[i]=simula_votos(dados)

    print(i)    

    }

mean(avalia)