Análise estatística 2

Jacques Wainer

Erros de amostragem (ou diferenças devido a sorte apenas)

Dados retirados (amostrados) de uma mesma fonte podem ter médias diferentes. A diferença é apenas por causa da sorte ou do azar. Isso é chamado de “erro de amostragem”. O nome erro não é um bom nome, pense como “ruído de amostragem”

Testes estatisticos tentam avaliar qual a probabilidade que a differeça das 2 (ou mais) amostras são apenas devido ao erro de amostragem.

Esse é o p-valor: a probabilidade que dado tao diferentes como as 2 amostras tenham vindo da mesma fonte de dados e que as diferenças são apenas devido ao erro de amostragem.

P-valor como decisão: p<0.05

Tradicionalmente, na maioria das Ciências e em particular em Computação usa-se o valor de 0.05 de p-valor para afirmar que a diferença é estatisticamente significante

Na prática, se o p-valor do seu teste der 0.051 vc não tem um paper para publicar!!!

Há criticas sobre usar o p-valor como decisão, e veremos isso numa próxima aula, mas tradicionalmente essa é a prática.

(Em física de particulas, o p-valor é 0.000 000 3 ou 3x10^{-7} ou a regra do 5 sigma)

Resumo: teste estatísticos

Quando você precisa e não precisa de testes estatísticos?

Se não há custos para tomar a decisão - não use

Se há custos diferentes para a decisão, é preciso mais evidência - use

O que testes estatísticos não fazem

Existem muitos testes estatísticos!

Poder de um teste

Material novo

Exemplos

set.seed(1234)
fonte=rnorm(10000,100,30)
a=fonte[sample(1:10000,7)]
b=fonte[sample(1:10000,12)]
d = rnorm(12,150,30)


> t.test(a,b)

    Welch Two Sample t-test

data:  a and b
t = 0.78595, df = 8.0367, p-value = 0.4544
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -28.82529  58.66978
sample estimates:
mean of x mean of y 
110.62933  95.70708 

> wilcox.test(a,b)

    Wilcoxon rank sum test

data:  a and b
W = 51, p-value = 0.4824
alternative hypothesis: true location shift is not equal to 0

wilcox.test roda o rank sum (dados não pareados), que é um teste não paramétrico

veja que o teste T deu um p-valor menor, e confirma o que foi dito que teste T é mais poderoso. Mas veja:

> t.test(a,d)

    Welch Two Sample t-test

data:  a and d
t = -2.5849, df = 7.7994, p-value = 0.03305
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -92.287107  -5.055969
sample estimates:
mean of x mean of y 
 110.6293  159.3009 

> wilcox.test(a,d)

    Wilcoxon rank sum test

data:  a and d
W = 12, p-value = 0.009764
alternative hypothesis: true location shift is not equal to 0

Na verdade a,b, e d não satisfazem as pressuposições do teste T (n>30) e então não deveríamos ter usado e teste T em nenhum dos exemplos.

Dados (e testes) pareados e não pareados

Testes Paramétricos e não paramétricos

Tabela (parcial) de testes

Dois conjuntos de dados

Para comparar 2 conjuntos de amostras/dados

pareado não pareado
paramétrico teste T pareado teste T não pareado
não paramétrico Wilcoxon signed rank Wilcoxon rank sum

Mais de dois conjuntos de dados

pareado não pareado
paramétrico repeated ANOVA ANOVA
não paramétrico Friedman Kruskal-Wallis

Testes para mais de um conjunto de dados são um pouco complicados. Os testes acima se baseiam não hipótese nula que todos os grupos vieram da mesma fonte. E portanto um p-valor pequeno apenas indica que não é provável que todas as médias sejam iguais. Mas o teste não diz quais grupos são estatisticamente diferentes dos outros.

Para determinar isso é preciso fazer testes post hoc Este blog tem alguma explicação sobre post hoc tests mas não vai muito a fundo. Na aula de problemas com o p-value veremos porque isso testes post hoc não são tão simples.

Veremos mais sobre multipas comparações na aula sobre tecnicas particulares para Aprendizado de Maquina.

Outros tipos de dados

Há outros testes para outras situações (dados não numéricos, dados 0/1 pareados e não pareados, etc)

Há testes para outras perguntas:

Testes estatísticos em computação

Testes para um só conjunto de dados

?t.test

Note o parâmetro mu

Testes em R

Testes em python

https://machinelearningmastery.com/statistical-hypothesis-tests-in-python-cheat-sheet/

Links

A maioria das aulas de estatistica no Youtube que eu conheco nao seguem a abordagem de fonte de dados que eu usei acima. Eles usam mais frequentemente a ideia de comparar um conjunto de medidas com 1 valor apenas (que eu menciono apenas brevemente). Mas há vários videos sobre isso e voce certamente aprenderá algo útil neles.

Estudo dos fatores que influenciam o p-valor.

Vamos estudar o impacto do tamanho das amostras (N), da diferença das medias (media1 - media2), e do ruido das amostras (desv) no p-valor

  1. Com media1 = 10, media2 = 12, e desv=5 fixos, plote a media dos p-valores com N=20 a N=100 com passo 20

  2. Com N=30, media1 = 10, e desv=3, plote plote a media dos p-valores com media2 de 12 a 20 com passo 2

  3. Com media1 = 10, media2 = 12, e N=30 fixos, plote a media dos p-valores com desv de 1 a 6 com passo 1

Verifique que o p-valor cai quando aumentamos o número de dados nos conjuntos, dai quando aumentamos a diferença entre os conjuntos e sobe quando aumentamos o ruido dos dados.

ave_pval <- function(N, media1, media2, desv){
    rep = 40
    pvals = NULL
    for (i in 1:rep){
        mod = t.test(rnorm(N,media1,desv), rnorm(N,media2,desv))
        pvals = c(mod$p.value, pvals)
    }
    mean(pvals)
}
library(ggplot2)
exp1 <- function(){
    Ns = seq(20,100,20)
    nn = length(Ns)
    out = data.frame(N = Ns, p.val = NA)
    for (i in 1:nn)
        out$p.val[i] = ave_pval(Ns[i],10,12,5)
    ggplot(out,aes(x=N,y=p.val))+ylim(c(0,0.5))+
       geom_hline(yintercept=0.05, linetype="dotted") +
       geom_point()
   }
exp2 <- function(){
    mus = seq(12,20,2)
    nn = length(mus)
    out = data.frame(mu2 = mus, p.val = NA)
    for (i in 1:nn)
        out$p.val[i] = ave_pval(30,10,mus[i],5)
    ggplot(out,aes(x=mu2,y=p.val))+ylim(c(0,0.5))+
       geom_hline(yintercept=0.05, linetype="dotted") +
       geom_point()
    }
       
       
exp3 <- function(){
    sds = seq(1,6)
    nn = length(sds)
    out = data.frame(sd = sds, p.val = NA)
    for (i in 1:nn)
        out$p.val[i] = ave_pval(30,10,12, sds[i])
    ggplot(out,aes(x=sd,y=p.val))+ylim(c(0,0.5))+
       geom_hline(yintercept=0.05, linetype="dotted") +
       geom_point()
    }

Na verdade o p-valor depende diretamente de medida (mu1-mu2/desv) que é chamado de tamanho de efeito.

Mantendo o tamanho de efeito fixo em 1/3

exp4 <- function(){
    mus = seq(12,20,2)
    delta = mus-10
    nn = length(mus)
    desv = 3*delta
    out = data.frame(mu2 = mus, p.val = NA)
    for (i in 1:nn)
        out$p.val[i] = ave_pval(30,10,mus[i],desv[i])
    ggplot(out,aes(x=mu2,y=p.val))+ylim(c(0,0.5))+
       geom_hline(yintercept=0.05, linetype="dotted") +
       geom_point()
    }

Veremos mais sobre tamanho de efeito numa próxima aula.

Esse programa em Python, mas sem a parte grafica

import numpy as np
from scipy.stats import ttest_ind

def ave_pvalue(N, mu1, mu2, sd):
    REP=40
    array1 = np.random.normal(mu1, sd, size=(REP, N))
    array2 = np.random.normal(mu2, sd, size=(REP, N))
    p_values = []
    for i in range(REP):
      t_statistic, p_value = ttest_ind(array1[i], array2[i])
      p_values.append(p_value)
    return np.mean(p_values)
  
def exp1():
  pval = []
  Ns = list(range(20,101,5))
  for N in Ns
    pval.append(ave_pvalue(N,10,12,5))
  return Ns,pval

def exp2():
  pval = []
  mus = list(range(11,21,1))
  for mu in mus
    pval.append(ave_pvalue(20,10,mu,5))
  return mus,pval

def exp3():
  pval = []
  sds = list(range(11,21,1))
  for sd in sds
    pval.append(ave_pvalue(20,10,12,sd))
  return sds,pval